Bruno
Bruno
Bruno
Bruno Sportisse
Field : Geophysics
Air Pollution
Modeling and Simulation
Jury composed of
Acknowledgements
I would like to thank François-Xavier Le Dimet for his constant help: his report about my
work was only a small part of the large amount of sollicitations to him.
Thanks to Matthias Beekmann for having written a report about my work. I hope that an
improved knowledge of my (our) activities will be an opportunity for joint works in the future.
I have been honored by the report of Greg Carmichael and by his presence at my defense in
Paris.
Thanks to Jan Verwer for having attended to my defense, even if my research works have
switched from simulation to modeling.
I would like to thank my colleagues for many joint works, especially Marc Bocquet, Jaouad
Boutahar, Édouard Debry, Rafik Djouad, Olivier Isnard, Vivien Mallet, Luc Musson-Genon,
Denis Quélo, Yelva Roustan, Karine Sartelet and Marilyne Tombette.
4
Contents
1 Introduction 7
3 Research works 19
3.1 Parameterizations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
3.1.1 Turbulent diffusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
3.1.2 Scavenging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
3.1.3 Segregation effects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
3.1.4 Plume In Grid models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.1.5 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.2 Numerical analysis for reactive dispersion . . . . . . . . . . . . . . . . . . . . . . 25
3.2.1 Splitting methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.2.2 Time integration of chemical kinetics . . . . . . . . . . . . . . . . . . . . . 27
3.2.3 Advection and mass consistency . . . . . . . . . . . . . . . . . . . . . . . 30
3.2.4 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.3 Aerosol modeling and simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
3.3.1 The General Dynamics Equation for aerosols . . . . . . . . . . . . . . . . 32
3.3.2 Stochastic methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
3.3.3 Modal methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
3.3.4 Size-resolved methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
3.3.5 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
3.4 Inverse modeling and data assimilation . . . . . . . . . . . . . . . . . . . . . . . . 40
3.4.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
3.4.2 Sequential versus variational methods . . . . . . . . . . . . . . . . . . . . 42
3.4.3 Coming issues . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
3.4.4 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
3.5 Sensitivity analysis, ucertainty propagation and ensemble forecast . . . . . . . . . 48
5
6 CONTENTS
3.5.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
3.5.2 Some possible strategies for assessing the impact of uncertainties . . . . . 49
3.5.3 Multi-modeling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
3.5.4 Ensemble forecast . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
3.5.5 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
3.6 Proxy and reduced models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
3.6.1 Motivations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
3.6.2 Slow/fast models and dynamic reduction . . . . . . . . . . . . . . . . . . 52
3.6.3 Building look-up tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
3.6.4 POD/EOF methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
3.6.5 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
3.7 Chemistry-Transport Models and applications . . . . . . . . . . . . . . . . . . . . 55
3.7.1 Development of CTM . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
3.7.2 Applications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
3.7.3 Future works . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58
Chapter 1
Introduction
• in an operational context (for instance for risk assessment of chemical, biological or radio-
logical releases), the model has to be robust and will be used with a small amount of input
data;
• for environmental forecast (urban air quality and photochemistry), the model will be fitted
in order to provide good skills for target species;
• for impact studies (or for process studies), a comprehensive model describing the whole
processes will be used in order to be runned in many conditions (far from the current
condition).
Air pollution modeling is based on the so-called Chemistry Transport Models (CTM here-
after). The core of a CTM (see Seinfeld and Pandis [1998]) is the simulation of the time evolu-
tion of chemical concentrations. The forcing conditions (meteorological conditions, emissions)
are fixed. Let ci stand for the concentration of a chemical gas-phase species i (in mol l −1 , ppb,
µg m−3 ...). The evolution is then governed by a reaction-diffusion-advection equation (equation
7
8 Chapitre 1 - Introduction
gaseous chemistry
nucleation condensation
H2SO4 evaporation
rainout
HNO3
coagulation
SO2 NOx
turbulent VOC washout
diffusion
VOC VOC
NOx primary dust
particles
NH3
dry deposition
TRAFFIC AGRICULTURE INDUSTRY BIOGENIC
of reactive dispersion):
∂ci ci
+ div (V ci ) = div ρK∇ + χi (c, T, I) − Λi (x, t)ci + Si (x, t) (1.1)
∂t | {z } ρ | {z } | {z } | {z }
advection | {z } gas phase chemistry scavenging volumic emission
turbulent diffusion
~ is the wind velocity, K the eddy diffusion coefficient, T the temperature, I the solar actinic
V
flux and ρ the air density. These variables are supposed to be known. χi stands for the chemical
sources and depends on the solar flux (through photolysis) and on temperature (through thermal
reactions). Λ is a parameterization of scavenging processes (by rain for instance). S i describes
the point sources.
Some boundary conditions are specified, for instance at ground for emission and deposition:
∂ci
−Kz = Ei (~x, t) − videp (~x, t)ci (1.2)
∂z | {z } | {z }
surfacic emission dry deposition
Kz is the vertical eddy coefficient, Ei is the flux of surfacic emission and videp is the dry deposi-
tion velocity.
A key component of these models is the computation of the chemical sources (χ i ). Gas-
phase chemical kinetics is well understood for some chemical species. Current comprehensive
gas-phase mechanisms describe up to hundreds of chemical species. In practice, only surrogate
lumped species (the so-called model species) are taken into account in CTM.
9
The condensed matter (in liquid or solid forms, Table 1.2 for the typical scales) has a key
role for many processes. One usually distinguishes the condensed matter related to water (rain
drop, cloud drop) and the aerosols (particles in suspension composed of a mixture of chemical
species). The gas-phase models have therefore to be extended to multiphase processes by de-
scribing the coupling to cloud droplets (diphasic models: gas/cloud) and the aerosol dynamics.
Source terms have then to be added to Equation (1.1).
A CTM can then be used for simulating the chemical composition of the atmosphere (forward
mode). It can be used for forecast by using data assimilation methods (coupling to observational
data). Another use is the identification of some poorly known input data (inverse modeling of
emissions for instance). The assessment of the uncertainties is another key point that can be
performed through multi-model methods or ensemble runs.
Outlines
My research activity has been led in this framework. I have focused on the following topics:
5. Treatment of uncertainties:
Which are the uncertainties of CTM ? How to extend mono-model runs to multi-model
runs ?
6. Model reduction:
How to deal with the high dimension ? How to build reduced models, especially for
integrated modeling ?
My research works are detailed in the following articles. I have omitted technical reports (related
to my application works especially with industrialists) except those with methodological points
or interesting case studies.
2.1 Parameterizations
Articles
Sportisse, B. (2001b). Box models versus Eulerian models in air pollution modeling. Atmos.
Env., 35:173–178
Sportisse, B. and Djouad, R. (2003). Mathematical investigation of mass transfer for atmo-
spheric pollutants into a fixed droplet with aqueous chemistry. J. Geophys. Res., 108(D2):4073
Sportisse, B. (2007d). A review of parameterizations for modeling dry deposition and scavenging
of radionuclides. Atmos. Env., (41):2683–2698
Sportisse, B., Bencteux, G., and Plion, P. (2000). Method of Lines versus Operator Splitting
methods for reaction-diffusion systems in air pollution modelling. Env. Mod. Soft., 15(6-7):673–
679
Djouad, R., Sportisse, B., and Audiffren, N. (2002). Numerical simulation of aqueous-phase
atmospheric models : use of a non-autonomous Rosenbrock method. Atmos. Env., 36:873–879
11
12 Chapitre 2 - Bibliography of my research works
Djouad, R. and Sportisse, B. (2003). Solving reduced models in air pollution modelling. Appl.
Numer. Math., 44(1):49:61
Sportisse, B., Quélo, D., and Mallet, V. (2007). Impact of mass consistency errors for atmo-
spheric dispersion. Atmos. Env., 41:6132–6142
Mallet, V., Pourchet, A., Quélo, D., and Sportisse, B. (2007b). Investigation of some numerical
issues in a Chemistry-Transport Model: gas-phase simulations. J. Geophys. Res., 112(D15301).
doi:10.1029/2006JD008373
Verwer, J. and Sportisse, B. (1998). A note on operator splitting analysis in the stiff linear case.
Technical Report MAS-R9830, CWI
Sportisse, B. and Mallet, V. (2005). Calcul Scientifique pour l’Environnement. Cours ENSTA.
86 pages
Sartelet, K. N., Hayami, H., Albriet, B., and Sportisse, B. (2005). Development and preliminary
validation of a Modal Aerosol Model for tropospheric chemistry: MAM. Aerosol Sci. and Tech-
nol., 40(2):118–127
Debry, E. and Sportisse, B. (2007a). Numerical simulation of the General Dynamics Equation
(GDE) for aerosols with two collocation methods. Appl. Numer. Math., 57(8):885–898
Debry, E. and Sportisse, B. (2007b). Solving aerosol coagulation with size-binning methods.
Appl. Numer. Math., 47:1008–1020
Debry, E., Fahey, K., Sartelet, K., Sportisse, B., and Tombette, M. (2007). A new SIze REsolved
Aerosol Model: SIREAM. Atmos. Chem. Phys., 7(6):1537–1547
Tombette, M. and Sportisse, B. (2007). Aerosol modeling at regional scale: Model-to-data com-
parison and sensitivity analysis over Greater Paris. Atmos. Env. In press
Sartelet, K., Debry, E., Fahey, K., Tombette, M., Roustan, Y., and Sportisse, B. (2007a). Simula-
tion of aerosols and gas phase species over Europe with the Polyphemus system. Part I: model-
to-data comparison for year 2001. Atmos. Env., 41:6116–6131. doi:10.1016/j.atmosenv.2007.04.024
Data assimilation and inverse modeling 13
Sartelet, K., Hayami, H., and Sportisse, B. (2007c). MICS-Asia Phase II: sensitivity to the
aerosol module. Atmos. Env. doi:10.1016/j.atmosenv.2007.03.005
Fahey, K., Debry, E., Foudhil, H., and Sportisse, B. (2005). Formulation, development and pre-
liminary validation of the SIze REsolved Aerosol Model, SIREAM. In Proceedings GLOREAM
2004, pages 7–17
Sportisse, B., Debry, E., Fahey, K., Roustan, Y., Sartelet, K., and Tombette, M. (2006b). Rap-
port du projet Primequal-Predit PAM (Pollution Atmosphérique Multiphasique). Modélisation.
Technical Report 2006-1, CEREA. 171 pages
Sportisse, B., Sartelet, K., Debry, E., Fahey, K., Roustan, Y., and Tombette, M. (2006c). The
SIze REsolved Aerosol Model (SIREAM) and the Modal Aerosol Model (MAM). Technical
documentation. Technical Report 2006-8, CEREA. 115 pages
Sportisse, B. and Quélo, D. (2003). Data assimilation and inverse modeling of atmospheric
chemistry. Proc. of Indian National Science Academy. Part A Physical Sciences, 69
Quélo, D., Sportisse, B., and Isnard, O. (2005b). Data assimilation for short-range dispersion of
radionuclides: a case study for second-order sensitivity. J. Environ. Radioactivity, 84:393–408
Quélo, D., Mallet, V., and Sportisse, B. (2005a). Inverse modeling of NO x emissions at regional
scale over Northern France. Preliminary investigation of the second-order sensitivity. J. Geo-
phys. Res., 110(D24310)
Krysta, M., Bocquet, M., Sportisse, B., and Isnard, O. (2006). Data assimilation for short-range
dispersion of radionuclides: an application to wind tunnel data. Atmos. Env., 40(38):7267–7279
Bocquet, M. and Sportisse, B. (2007). Modélisation inverse pour la qualité de l’air : éléments
de méthodologie et exemples. Pollution Atmosphérique. Accepted for publication
Lacour, S. and Sportisse, B. (2007). Estimation of indoor deposition velocity for ozone with a
simplified reactive box model. Atmos. Env. In revision
Sportisse, B. and Mallet, V. (2007). Data assimilation and inverse modelling of slow/fast atmo-
spheric chemistry. Journal of Computational Geosciences. Submitted
14 Chapitre 2 - Bibliography of my research works
Quélo, D., Sportisse, B., Berroir, J., and Charpentier, I. (2002). Some remarks concerning in-
verse modeling and data assimilation for slow-fast atmospheric chemical kinetics. In Sportisse,
B., editor, APMS 2001, Geosciences, pages 499–513. Springer
Krysta, M., Bocquet, M., Isnard, O., Issartel, J.-P., and Sportisse, B. (2004). Data assimilation
of radionuclides at short and regional scales. In Farago, I., Georgiev, K., and Havasi, A., editors,
Proceedings of the NATO Advanced Research Workshop on Advances in Air Pollution Modeling
for Environmental security, pages 275–285. Kluiwer
Isnard, O., Krysta, M., Bocquet, M., Dubiau, P., and Sportisse, B. (2005). Data assimilation of
radionuclides atmospheric dispersion at small scale: a tool to assess the consequences of radio-
logical emergencies. In Proceedings of the IAEA Conference. Rio Conference
Articles
Djouad, R., Audiffren, N., and Sportisse, B. (2003a). Sensitivity analysis using automatic dif-
ferentiation applied to a multiphase chemical mechanism. Atmos. Env., 37(22):3029:3038
Mallet, V. and Sportisse, B. (2005). A comprehensive study of ozone sensitivity with respect to
emissions over Europe with a chemistry-transport model. J. Geophys. Res., 110(D22)
Tombette, M., Fahey, K., Sartelet, K., and Sportisse, B. (2005). Aerosol modelling at regional
scale: a sensitivity study with the Polyphemus platform. In Proceedings of GLOREAM 2005.
Also as Report CEREA 2005-50
Model reduction 15
Mallet, V. and Sportisse, B. (2006b). Ensemble-based air quality forecasts: a multi-model ap-
proach applied to ozone. J. Geophys. Res., 111(D18):18302
Mallet, V. and Sportisse, B. (2006a). Air quality modeling: from deterministic to stochastic
modeling. Computers and Mathematics with Application. In press
Djouad, R. and Sportisse, B. (2002a). Partitioning techniques for reduction in chemical kinetics.
APLA: an Automatic Partitioning and Lumping Algorithm. Appl. Numer. Math., 43(4):383:398
Djouad, R., Sportisse, B., and Audiffren, N. (2003b). Reduction of multiphase atmospheric
chemistry. Journal of Atmospheric Chemistry, 46:131–157
Sportisse, B. and Djouad, R. (2007). Use of Proper Orthogonal Decompositions for the reduction
of atmospheric chemistry. J. Geophys. Res., 112(D06303)
Larrouturou, B. and Sportisse, B. (1997). Some mathematical and numerical aspects of reduc-
tion in chemical kinetics. In Bristeau, M. and al, editors, Computational Science for the 21st
16 Chapitre 2 - Bibliography of my research works
Sportisse, B., Jaubertie, A., and Plion, P. (1997). Reducing mechanisms in chemical kinetics for
the simulation of reactive transport; an application to air pollution modelling. In Proceedings
of the St Venant Symposium
Jaubertie, A., Plion, P., Sportisse, B., Aumont, B., and Toupance, G. (1998). Reduction of
kinetic schemes for atmospheric chemistry. In Proceedings of the Conference APMS98. INRIA,
France
Djouad, R. and Sportisse, B. (2002b). Some reduction techniques for simplifying air pollution
models. In Sportisse, B., editor, APMS 2001, Geosciences, pages 235–245. Springer
Boutahar, J. and Sportisse, B. (2002). Reduction methods for atmospheric chemistry. In Global
and Regional Atmospheric Modeling Meeting, pages 133–139. Un. of Aveiro, Portugal
Boutahar, J. and Sportisse, B. (2005). Reduction methods and uncertainty propagation: appli-
cation to a Chemistry-Transport-Model. In Proceedings of the TAM-TAM conference
Boutahar, J., Lacour, S., Mallet, V., Quélo, D., Roustan, Y., and Sportisse, B. (2004). De-
velopment and validation of a fully modular platform for numerical modelling of air pollution:
Polair3D. Int. J. Env. and Pollution, 22(1-2)
Roustan, Y., Bocquet, M., Musson-Genon, L., and Sportisse, B. (2006). Modélisation du mer-
cure, du plomb et du cadmium à l’échelle du continent européen. Pollution Atmosphérique,
(191):317–327
Mallet, V., Korssakissok, I., Quélo, D., and Sportisse, B. (2007a). Polyphemus : un système
modulaire de modélisation pour la dispersion atmosphérique et l’évaluation des risques. Pollu-
tion Atmosphérique
Mallet, V., Quélo, D., Sportisse, B., Ahmed de Biasi, M., Debry, É., Korsakissok, I., Wu, L.,
Roustan, Y., Sartelet, K., Tombette, M., and Foudhil, H. (2007c). Technical Note: The air qual-
ity modeling system Polyphemus. Atmos. Chem. Phys. Discuss., 7(3):6,459–6,486
Quélo, D., Krysta, M., Bocquet, M., Isnard, O., Minier, Y., and Sportisse, B. (2007). Validation
of the polyphemus system: the ETEX, Chernobyl and Algeciras cases. Atmos. Env., 41:5300–
Books 17
5315
Sportisse, B. (2007c). A review of current issues in air pollution modeling and simulation. Com-
putational Geosciences, 11(2):159–181
Sartelet, K., Hayami, H., and Sportisse, B. (2007b). Dominant aerosol processes during high-
pollution episodes over Greater Tokyo. J. Geophys. Res. Accepted for publication
Quélo, D., Lagache, R., and Sportisse, B. (2004). Etude de l’impact qualité de l’air des scénarios
PDUs sur Lille à l’aide du modèle de Chimie-Transport Polair3D. Technical Report 2004-28,
CEREA. Rapport PREDIT
Taghavi, M., Musson-Genon, L., and Sportisse, B. (2004). Impact des rejets de la centrale
thermique de Martigues sur la qualité de l’air. Modélisation avec Polair3D. Technical Report
2004-27, CEREA. Rapport pour la Mission Thermique EDF
Roustan, Y., Bocquet, M., Musson-Genon, L., and Sportisse, B. (2005). Modeling of mercury at
the European scale with the chemistry-transport model Polair3D. In Proceedings GLOREAM
2004, pages 121–130. Also as Technical Report CEREA 2005-04
Tombette, M. and Sportisse, B. (2006). Rapport de contrat préliminaire: évaluation des impacts
du CPT Porcheville. Technical Report 2006-55, CEREA. Rapport pour la Mission Thermique
EDF
Lagache, R., Declercq, C., Quélo, D., Sportisse, B., Palmier, P., Quetelard, B., and Haziak, F.
(2006). Evaluation du PDU de Lille-Métropole sur le trafic, les concentrations de polluants at-
mosphériques et la mortalité. In Actes de la 15ième conférence sur les transports et la pollution
de l’air
Roustan, Y., Lasry, F., and Sportisse, B. (2007). Modélisation de l’impact des émissions des
centrales EDF et SNET en France sur le transport atmosphérique des PM 10 et PM2.5 . Technical
Report 2007-9, CEREA
2.8 Books
As Editor
Sportisse, B., editor (1998). Air Pollution Modeling and Simulation APMS98, INRIA/ENPC
APMS Conference. INRIA
Sportisse, B., editor (2001a). Air Pollution Modeling and Simulation APMS01. Springer Verlag
18 Chapitre 2 - Bibliography of my research works
Le Dimet, F. and Sportisse, B., editors (2002). Data Assimilation for Geophysical flows, INRIA-
CEA-EDF School on Applied Non-linear Problems. INRIA
As Author
Sportisse, B. (2007b). Pollutions atmosphériques. Enjeux, processus, modélisation. Springer.
350 pages
2.9 Miscellaneous
Dumas, L., Coron, F., Dupuis, D., Pallegoix, J., and Sportisse, B. (1994). Applications of DSMC
methods to aerodynamic problems in an engineering context. In Méthodes de Monte-Carlo et
applications aux gaz raréfiés. Ecole CEA-INRIA-EDF, 1994
Coron, F., Pallegoix, J., and Sportisse, B. (1994). DSMC computations of complex test cases and
industrial configuration. In 19th International Symposium on Rarefied Gas Dynamics, Oxford
Sportisse, B. (2003a). Ozone des villes, ozone des champs: qualité de l’air et pollution pho-
tochimique. PCM Le Pont
Research works
The following presentation is mainly based on the review article Sportisse [2007c]. My works
are referenced through underlined references. The general framework (atmospheric chemistry
and physics for air pollution) is also detailed in a book written on the basis of a 10 years course
at Ecole Nationale des Ponts et Chaussées (Sportisse [2007b]).
Most of the works have been done with PhD students or post-doctoral fellows I have su-
pervised. For PhD works, I can cite: Rafik Djouad for model reduction and aqueous-phase
chemistry, Jaouad Boutahar for model reduction and 3D modeling, Edouard Debry for aerosol
modeling, Denis Quélo for inverse modeling and dispersion of radionuclides (among many other
applied studies), Yelva Roustan (I have initiated his PhD work with Luc Musson-Genon and
the supervision has then be made jointly with Marc Bocquet) for mercury and heavy metals,
Monika Krysta (I have initiated her PhD work and the supervision has then be made jointly with
Marc Bocquet) for inverse modeling of radionuclides, and Marilyne Tombette (I have initiated
her PhD work and the supervision has then be made jointly with Patrick Chazette) for aerosol
modeling and radiative effects. For postdoctoral works, I can cite the works of Karine Sartelet
(for air quality modeling and modal aerosol modeling) and Kathleen Fahey (for 3D modeling
of aerosols). The main part of the works devoted to radionuclides have been made with my
colleague from the French Institute for Nuclear Safety and Radiological Protection, Olivier Is-
nard. Even if we have not published papers together, many works have been initiated through
discussions with Luc Musson-Genon, especially when applied to topics relevant for Electricité de
France. Many works have also been initiated and discussed with my colleague Stephanie Lacour
for applications relevant for the French Ministry of Transport.
To conclude, many works have been performed with Vivien Mallet, first with his PhD work
devoted to uncertainties and ensemble forecast, and then for many topics related to Chemistry
Transport Models within the Polyphemus project.
3.1 Parameterizations
Subgrid parameterization is a key component of geophysical models. Many unresolved processes
have indeed characteristic scales much finer than the resolved numerical scale (typically the size
of a cell in a Finite Differences framework).
hΨi will stand for the average value of a physical field Ψ (the average has to be precised
0
rigorously). In practice, this corresponds to the numerical value. The local fluctuation Ψ is
0
defined through Ψ = Ψ − hΨi.
For
D non-linear
E processes, given for instance by quadratic terms such as Ψ 1 Ψ2 , a correlation
0 0
term Ψ1 Ψ2 has then to be parameterized as a function of the resolved values hΨ1 i and hΨ2 i.
19
20 Chapitre 3 - Research works
The issue is similar to the one arising for turbulence modeling (closure scheme).
The applications for air pollution modeling are related to the vertical turbulent fluxes (pa-
rameterization of Kz ), the microphysical mass transfer, the coupling between chemistry and
turbulence and Plume-In-Grid modeling.
My works have focused on mass transfer (scavenging by rain drops and cloud drops). Any-
how, I have made the choice to briefly review the whole topics.
3.1.2 Scavenging
The mass transfer between gas phase, aqueous phase and condensed matter (cloud droplets, rain
drops and aerosols) requires appropriate parameterizations. The typical microphysical scales are
indeed the following ones: an aerosol diameter range from a few nanometers to a few micrometers,
a cloud droplet diameter from a few tens of micrometers to one hundred micrometers, etc. Such
scales cannot be solved in a numerical model.
The atmospheric residence time of radionuclides is stronly related to wet scavenging and dry
deposition. The processes to be taken into account may vary for the radionuclides (iodine and
cesium). For instance, cesium is mainly bound to aerosols, which emphazises the role of wet
scavenging for its atmospheric residence. A review of the available parameterizations is presented
in Sportisse [2007d]. A key point is the representation of the rain distribution. The scavenging
coefficients are non-linear functions of the aerosol size and of the rain characteristics (see Figure
3.1 for the collision efficiency, a key component for scavenging).
A size resolved aerosol model is usually not available for radionuclides. It is therefore
mandatory to deal with parameterized models (see Sportisse [2007d]). A reasonable model
is Λ = A × p0 B (p0 is the rain intensity in mm/h, A ' 10−5 typically, B ' 0.8).
Parameterizations 21
0.5
Dr=0.1 mm
Dr=0.5 mm
0.0 Dr=1 mm
Dr=2 mm
-0.5
-1.0
log10 E(Dr,dp)
-1.5
-2.0
-2.5
-3.0
-3.5
-4.0
-3.0 -2.5 -2.0 -1.5 -1.0 -0.5 0.0 0.5 1.0
log10 dp (in m)
Figure 3.1: Scavenging of an aerosol. Collision efficiency as a function of the aerosol size for
different representative raindrop diameters (Dr ). dp is the wet aerosol diameter.
22 Chapitre 3 - Research works
S(IV)
S(VI)
practice, one can define a lumped model by replacingR the microphysical concentrations c a and
cg by averaged concentrations, for instance c̄a = Ω1a Ωa ca (ω, t)dω. The time evolution of the
averaged concentrations is then given by a system of Ordinary Differential Equations (ODE):
dc̄g = χg (c̄g ) − L I(c̄a , c̄g )
dt (3.4)
dc̄a
= χa (c̄a ) + I(c̄a , c̄g )
dt
L is the liquid water content (ratio of the volume Ωa to the volume Ωg ). The mass transfer
c̄a
flux I is parameterized as I(c̄a , c̄g ) = kmt (c̄g − ), with kmt a mass transfer coefficient. A
HRT
classical formula is: 2 −1
a 4a
kmt = + (3.5)
3Dg 3v̄α
Such derivations are rigorously presented in Sportisse and Djouad [2003] with the use of
asymptotic expansions with respect to “small” parameters (defined as Peclet-Dahmköhler num-
bers). The mathematical framework can be found in Fife [1979] for reaction-diffusion systems
and Vasileva et al. [1994] for boundary layer methods.
It is for instance proven that the double film parameterization has to be preferred under
some conditions: 2 −1
a a2 4a
kmt = + + (3.6)
3Dg 3HRT Da 3vα
There are still many other challenging issues: how to extend this approach to polydisperse
models ? How to take into account surface reactions (not described in Eq. 3.3) ?
Below-cloud scavenging
A second example is given by below-cloud scavenging. In this case, the scavenging is related to
falling raindrops (of falling velocity U and of radius a). With notations similar to below, the
microphysical model is:
∂cg L ca
=− cg −
∂t τtransfer HRT
∂ca ∂ca 1 ca (3.7)
+ U = c −
g
∂t ∂z τtransfer HRT
τtransfer is a timescale associated to mass transfer (not detailed here), given as the inverse value
of a mass transfer coefficient.
This model can be replaced by the following parameterization after having applied an asymp-
totic expansion (the “small” parameter is in this case related either to a large solubility of to a
small liquid water content):
dcg
= − Λ g cg (3.8)
dt
with Λg (p0 , a, HRT, τtransf er ) a scavenging coefficient. One refers to Sportisse and Du Bois [2002]
for more details.
0.1
-0
-0.1
-0.2
-0.3
Is
-0.4
-0.5
-0.6
-0.7
0 10 20 30 40 50
Time (UTC)
Figure 3.3: Evolution of surface segregation intensity (15-16 September 1989), from observations
in the Netherlands (from Vila-Guerau de Arellano et al. [1993]).
assumes that this term is equal to χi (hci), which may be a crude approximation for the fastest
reactions. The underlying assumption is indeed that the averaging is quicker than any chemical
reactions. We refer for instance to the review Vila-Guerau de Arellano et al. [2004].
A typical example is given by the so-called segregation effect between an updraft of emitted
NO and a downdraft of O3 (entrainment) in the atmospheric boundary layer. The chemical
reaction of O3 with NO (titration) is:
O3 + NO → NO2 + O2
For the sake of clarity we consider two processes that do not depend on time. The evolution
equation is therefore:
dc
= A1 (c) + A2 (c) , c(0) = c0 (3.10)
dt
A1 (for instance diffusion) and A2 (for instance chemistry) are the source terms associated to
the processes.
Let ∆t be the splitting timestep (for a CTM, it typically ranges from 100 to 1000 seconds
at continental scales). We formally write LA1 +A2 (c0 , t) the solution of (3.10) at time t.
The simplest method is the first-order splitting: one first solve the process A 1 over [0, ∆t]
and then (after having modified the initial condition) the process A2 over [0, ∆t]. The solution
is then LA2 (LA1 (c0 , ∆t), ∆t). Notice that we can also reverse the sequence between A1 and A2 .
The so-called Strang splitting (Strang [1968]) eliminates the lack of symmetry of the previous
method. There are three steps in the sequence: A1 is integrated over [0, ∆t/2], then A2 is inte-
grated over [0, ∆t] and finally A1 over [0, ∆t/2]. The solution is therefore LA1 (LA2 (LA1 (c0 , ∆t/2), ∆t), ∆t/2).
Operator splitting methods induce a loss of accuracy. This is usually estimated with the
computation of the local splitting error in the linear case for Ai (c) = Ai c (Ai is a matrix).
In this case, LA (c0 , ∆t) = eA∆t c0 . The splitting error is related to the commutation default:
eA1 +A2 = eA2 eA1 if and only if A1 A2 = A2 A1 .
Asymptotic analysis with respect to ∆t → 0 prove that the Strang splitting method is of
second-order (and that the first-order method is well named!).
Some results are available for the reaction/diffusion/advection PDE with the extension to
the nonlinear case by using the Lie formalism (Lanser and Verwer [1998] or Sportisse [1999] for
direct calculations). The key results are the following ones (boundary conditions are not taken
into account):
• advection and chemistry commute if chemistry does not depend on space and if the wind
field is divergence free;
• advection and diffusion commute if the wind field and the diffusion matrix do not depend
on space;
• diffusion and chemical kinetics commute if chemical kinetics is linear and does not depend
on space.
It is clear that the assumption concerning the independence with respect to time and space is
not rigorously met for the eddy diffusivity matrix and the temperature field. This is however a
rather good local assumption. On the contrary, the assumption of linearity for chemistry is clearly
wrong (chemical reactions are supposed to describe... reactions between reactants!). However,
these results give an indication for the main source of errors related to operator splitting methods
in air pollution models: chemistry and vertical diffusion have to be handled with a special care
(Graf and Moussiopoulos [1991]).
A difficult point (not addressed here) is also related to the way the boundary conditions have
to be taken into account (see for instance Sommeijer et al. [1981] or Sportisse et al. [2000] for
an illustration).
For atmospheric chemistry, the previous results (especially the order analysis) cannot be
applied as such. A major point is indeed that one process (for instance chemistry) contain fast
dynamics, that quickly reach an equilibrium (see below for more details and the relation to
stiffness). In practice, this means that there exist processes whose timescales τ are much lower
Numerical analysis for reactive dispersion 27
than the splitting timestep ∆t. The asymptotic analysis ∆t → 0 is therefore no more justified.
We refer for instance to Sportisse [2000]; Verwer and Sportisse [1998]; Sportisse et al. [2000] for
a deeper understanding.
One key result is the impact of the splitting sequence: it is advocated to end the sequence
with the stiff processes (chemistry for instance in the case of gas-phase chemistry, mass transfer
for aerosol dynamics) in order to ensure that the system is still along the equilibrium manifold
at each step (see Figure 3.4). A second result is related to the so-called order reduction of the
numerical algorithms. Similar results have been obtained in other fields, for instance in combus-
tion (Yang and Pope [1998]). See also the discussion for oceanography in Higdon and Bennett
[1996] or for climate models in Williamson [2001].
Methods that avoid splitting are also used. Among them, one can cite implicit-explicit
methods that are based on a partitioning between the stiff processes (with an implicit treatment)
and the nonstiff processes (with an explicit treatment). There are many variations of this method
(for instance Verwer et al. [1996a]; Sun [1996]; Wolke and Knoth [1998]; Nguyen and Dabdub
[2003]).
For instance, let us suppose that A1 and A2 are the nonstiff process and the stiff process,
respectively. A possible implementation (usually referred as Source Splitting) is to first solve
LA1 (c0 , ∆t) − c0
A1 and then to solve A2 after having added a correction term Ã1 = without
∆t
modifying the initial condition. The result is then: LA2 +Ã1 (c0 , ∆t). The key advantage is that
the integration is always performed in the vicinity of the equilibrium manifold of the stiff process
(Sportisse [2000]).
Another approaches are the so-called Approximate Matrix Factorizations (also referred as
“internal splitting”, see Berkvens et al. [2002]). The idea is to perform splitting at the algebraic
level. Suppose that an implicit algorithm (let say Backward Euler) is used for solving equation
(3.10) in the linear case.
We have therefore (I − (A1 + A2 )∆t) cn+1 = cn . The method is based on approximations
such as (I − (A1 + A2 )∆t) ' (I − A2 ∆t)(I − A1 ∆t) + O(∆t2 ).
y
y
1 2 y = h(x) 1-2
M OL NT S
x x
y
y 2 2
1
1
(T − χ) (χ − T )
x
x
y
y 3 3
2
2
1
1
(χ − T − χ) (T − χ − T )
x
x
Figure 3.4: Impact of the splitting sequence. M OL stands for the Method of Lines, N T S for
No Time Splitting. T stands for transport and χ for chemical kinetics. x and y are the solw
and fast components of concentrations, respectively. The equilibrium manifold defined by fast
processes is defined by the constraint y = h(x).
Numerical analysis for reactive dispersion 29
10
−2
−4
−6
−8
−10
−12
0 10 20 30 40 50 60
Figure 3.5: Typical distribution of the eigenvalues for the Jacobian matrix associated to the
chemical mechanism RACM (logarithmic scale).
Among many candidates, the second-order Rosenbrock solver (Verwer et al. [1999] for gas-
phase chemistry, Djouad et al. [2002] for aqueous-phase chemistry) is a promising algorithm as it
realizes a good trade-off between a reasonable accuracy and a rather low computational burden.
The Rosenbrock method is defined in the following way from iteration n to iteration n + 1:
k1 + k 2
cn+1 = cn +
2 (3.12)
(I − γ∆t J)k1 = ∆tχ(cn )
n
(I − γ∆t J)k2 = ∆tχ(c + ∆t k1 ) − 2γ∆t J k1
√
J = ∂χ∂c is the Jacobian matrix for chemical kinetics. γ = 1 + 1/ 2 ensures the L-stability
of the method (see Verwer et al. [1998] and also Mallet et al. [2007b] for the non-autonomous
case). Preprocessing the computation of J with appropriate techniques that take into account
the sparse nature (a chemical species is actually involved in a small set of reactions) is a powerful
technique (Sandu et al. [2003]; Daescu et al. [2003] and Mallet and Sportisse [2004]).
Another candidates are provided by explicit and quasi-explicit methods. The so-called QSSA
(Quasi Steady State Assumption) methods (Young and Boris [1977] for the original paper and
then for instance Jay et al. [1997]; Verwer and Van Loon [1994]; Mott et al. [2000]) are still
widely used even if they do not perform well in fair benchmarks (with error/CPU diagrams).
Such algorithms are strongly related to the so-called production-loss form of (3.11):
with Pi and Li two nonnegative terms. This form is easy to get by partitioning the chemical
reactions into two sets: reactions in which the species i is a reactant (loss term), reactions
in which it is a product (production term). The fast species and the slow species are defined
by timescales Li 1 and Li 1, respectively. There are many variations on the basis of
this partitioning of species. I refer for instance to my work Djouad and Sportisse [2002a] for a
30 Chapitre 3 - Research works
precise definition of the partitioning. In the simplest form of the QSSA method, the fast species
are computed from iteration n to iteration n + 1 with the equilibrium cn+1
i =L Pi
i
while the slow
Pi n Pi
species are computed with ci = Li + (ci − Li ) exp(−Li ∆t).
Beyond QSSA methods, quasi-explicit methods are also very popular. The Two-Step method
(Verwer [1994]; Verwer and Simpson [1995]) is an interesting algorithm that can be seen as an
extension of the EBI method (Euler Backward Iterative, Preussner and Brand [1981]; Hertel
et al. [1993]). Using the Backward Euler method leads to the following nonlinear algebraic
system:
cn+1 − cni
i
= χi (cn+1 ) (3.14)
∆t
This can be solved by an iterative algorithm, for instance the Newton algorithm, which re-
quires to compute the Jacobian matrix J. Another approach is to approximate χ i (cn+1 ) '
Pi (cn ) − Li (cn )cn+1
i . The algebraic system can then be solved with:
Variations of this idea lead to the EBI (first-order) and Two-Step (second-order) algorithms.
Notice that P (.) and L(.) can be evaluated with the components of cn+1
i that have been already
computed (which is the starting point of many tailored techniques).
Let us mention for the sake of completeness the combined explicit/implicit methods (for
instance Gong and Cho [1993]), that are based on the partitioning of chemical species in a set
of slow species (integrated with an explicit solver) and a set of fast species (integrated with an im-
plicit solver). We will go back to this approach in Section 3.6 (see also Djouad and Sportisse [2003]).
To conclude, a specific issue for the time integration of chemical kinetics is related to the
positivity of concentrations. A usual tricky approach is the clipping (set the negative numerical
concentrations to 0 or to a small predefined parameter ). We refer to Sandu [2001b] for a rigorous
method. The impact for real cases is also evaluated in Mallet et al. [2007b].
The conservation of uniform mixing ratio in the case of passive tracers is then easy to
prove: if q(0) = 1, one gets for any time t, q(t) = 1. This property is usually referred as mass
consistency. The key point is the consistency between the wind field used for ρ and the wind field
used for c. The discrepancy may result from different numerical algorithms and discretization
strategies between the computation of the wind velocity (inside the meteorological model) and
the advection step (inside the CTM).
We briefly present three possible methods in order to circumvent these difficulties.
A first approach is the renormalization technique. Let cn+1 = S(cn , V n ) be the advection
scheme used by the CTM (supposed to be an explicit solver, which is usually the case in practice).
Let us define ρ̃n+1 = S(ρn , V n ) the density that is computed with this scheme (this is not the
“true” density ρn+1 if the meso-scale model does not use the same scheme). The renormalization
approach is then based on the correction:
ρn+1
cn+1 = S(cn , V n ) (3.16)
ρ̃n+1
It is of course easy to check that mass consistency is met. One refers for instance to Byun
and Lee [2002]; Lee et al. [2004]. A major drawback of this approach is that it may artificially
introduce mass.
A second approach is based on the use of mass mixing ratio in the computation of fluxes after
∂c c
having noticed that + div(ρ V ) = 0. Convergence results for this algorithm are however not
∂t ρ
proved. One refer to Petersen et al. [1998] for the special case of the upwind scheme (without
flux limiting).
A third approach is to modify the vertical component of wind velocity, w, in order to satisfy
mass consistency. Suppose that horizontal advection is defined by cH = SH (cn , un , v n ) (even-
tually with directional splitting). One defines ρH = SH (ρn , nn , v n ). Let the vertical advection
scheme be given by cn+1 = SV (cH , wn ). We therefore try to compute w n such that the con-
straint ρn+1 = SV (ρH , wn ) is met. Of course, this implies that this constraint may be solved in
wn , which implies that a nonlinear scheme (for instance with flux limiting) cannot be used for
vertical advection. One example of a scheme easy to “invert” is the upwind scheme (Hu et al.
[2006]; Odman and Russel [2000]).
These approaches are clearly defined and benchmarked in Sportisse et al. [2007] for two case
studies (air quality modeling over Europe and Chernobyl accident). The third approach is an
outlyer, especially due to the strong modification of the vertical wind profile (Figure 3.6).
3000
dst3
upwind
2500
2000
height (m)
1500
1000
500
0
-2.0 -1.5 -1.0 -0.5 0.0 0.5
x1e-3
vertical wind (m/s)
Figure 3.6: Profile of the vertical wind velocity: reference case (DST3) and vertical adjustment
method (from Sportisse et al. [2007]).
One cas also wonder if the next generation of CTM will not be based on unstructured
methods, similar to those used in CFD (Computational Fluid Dynamics).
• one has to deal with high-dimensional models (the number of chemical species inside the
aerosol mixture has to be multiplied by the number of degrees of freedom chosen to describe
the size),
Aerosol modeling and simulation 33
500
deliquescence branch
metastable branch
400
200
100
0
0 20 40 60 80 100
relative humidity
Figure 3.7: Hysteresis for a salt of ammonium nitrate (computed with Isorropia: T =300 K,
[HNO3 ]=[NH3 ]=2 µg).
• the processes that affect the size and chemical composition of aerosols are nonlinear and
discontinuous (thresholds for phase transition: see for instance Figure 3.7 for the hysteresis
of deliquescence/cristallisation).
The evolution of size distribution and of chemical composition for aerosols is governed by the
so-called General Dynamics Equation for aerosols (GDE in the sequel). It describes the impact
of some processes (Figure 3.8) such as nucleation, coagulation, and mass transfer (condensation
and evaporation).
Nucleation The formation of the smallest aerosols is induced by the agregation of gaseous
molecules through thermodynamically stable clusters. One usually takes into account binary
nucleation (H2 O-H2 SO4 ) or ternary nucleation (H2 O-H2 SO4 -NH3 ).
Coagulation Coagulation is driven for instance by the agitation of the fluid. In practice,
one only describes Brownian coagulation due to thermal agitation. Other effects (the so-called
phoretic effects due to gradients in fields (such as temperature, electric fields, etc) are not taken
into account. Coagulation may be neglected for the evolution of particles above a few µm.
Condensation/evaporation Many gaseous species with a low saturation vapour pressure may
condense onto or evaporate from existing particles. This mass transfer is governed by the
gradient between the concentration at the aerosol surface (supposed to be at local equilibrium
with the internal mixture) and the concentration far from the aerosol.
Condensation is a fast process for submicronic aerosols.
The composition of the continental inorganic aerosol is mainly governed by the interaction
between sulfate, nitrate and ammonium. The corresponding gaseous species are:
34 Chapitre 3 - Research works
heterogeneous condensation
reactions coalescence
evaporation
wet
scavenging
dry deposition
D < 0.01 0.01 < D < 0.1 0.1 < D < 1 1 < D < 10
gas molecules aerosols aerosols aerosols aerosols cloud droplets
• sulfur dioxide (SO2 ) to be oxided to sulfuric acid (H2 SO4 ) whose saturation vapour pressure
can be neglected (H2 SO4 is only under particulate form);
Nitrate and sulfate are in competition inside aerosols for the formation of ammonium nitrate
(NH4 NO3 =NH+ – + 2–
4 +NO3 ) and ammonium sulfate ((NH4 )2 SO4 =2NH4 +SO4 ). A key point is
that the second salt is preferred.
Mass transfer to cloud droplets Part of the aerosol distribution acts as Cloud Condensation
Nuclei. The chemical composition of the activated aerosols is then governed by aqueous-phase
chemistry.
Heterogeneous reactions Surface reactions at the surface of condensed matter (aerosols and
cloud droplets) may have a significant impact for some species like ozone. These processes are
usually described by first-order kinetics.
Loss processes Scavenging and dry deposition have also to be taken into accound. For the
coarsest aersols, gravitational settling (sedimentation) has also to be described.
The aerosols are complicated mixtures of organic species, inorganic species, soot (black car-
bon) and liquid water. A typical composition of the continental aerosol over Europe can be
found in Figure 3.9. Notice that part of the composition is not known and that there are many
Aerosol modeling and simulation 35
Nitrate
Organics
15.0%
20.0%
Ammonium
7.0%
Black Carbon
5.0%
13.0% 9.0%
Sulfate
Mineral Dust
4.0%
27.0%
Sea salt
Other (27 %)
Figure 3.9: Composition of the “European” aerosol (PM10 ) in urban conditions, from Putaud
et al. [2003]). ’Other’ refers to the unknown part.
Let n(m, t) be the number distribution (for instance with respect to aerosol mass) and q i (m, t)
the mass distribution for species i (we use the same notation as for the mass mixing ratio in
the previous section). m0 is the smallest aerosol dry mass (in practice given by the nucleation
threshold). The GDE is composed of the following equations:
• For the number distribution:
Z
∂n 1 m−m0
(m, t) = K(u, m − u)n(u, t)n(m − u, t) du
∂t 2 m0
| {z }
coagulation gain
Z ∞ (3.17)
∂(I0 n)
− n(m, t) K(m, u)n(u, t) du − + δ(m0 , m)J0 (t)
| m0
{z } | ∂m{z } | {z
nucleation
}
mass transfer
coagulation loss
P c
An aerosol of mass m is composed by nc components of mass mi (m), such that i=ni=1 mi (m) =
m. Notice that qi (m) = mi (m)n(m). An aerosol with a given mass is assumed to have a unique
composition (internal mixing).
36 Chapitre 3 - Research works
K(., .) is the coagulation kernel (due to Brownian motion, not detailed here). I i is the
condensation/evaporation rate for volatile species i:
where ai (m) stands for a parameterization of a mass transfer coefficient (a function of accommo-
dation, of molecular diffusion and of particle size), ci is the gas-phase concentration of species i
and ceq
i is the surface concentration of species i (supposed to be at equilibriumP with the aerosol
mixture). ηi is the Kelvin coefficient that describes the curvature effect. I0 = i Ii is the total
rate.
J0 is the nucleation rate, given by parameterizations as a function of thermodynamic con-
ditions (temperature and relative humidity) and of some gaseous concentrations (water, ammo-
nium, sulfate).
The internal chemistry is usually neglected.
This set of Partial Integro-Differential equations is strongly coupled, especially through the
condensation/evaporation term. One key component is the availability of a thermodynamic
model for computing ceq i (it is the most expensive part of the model, as far as CPU is concerned,
Zhang et al. [1998b]).
Notice that the GDE is usually not solved simultaneously with the other processes. We re-
fer to Müller [2001] for the investigation of splitting errors related to microphysical-multiphase
chemical systems.
One refers to Sportisse et al. [2006a,b] for a description of the parameterizations used in the
aerosol models that have been developed under my leadership: the Model Aerosol Model MAM
(Sartelet et al. [2005]) and the SIze REsolved Aerosol Model SIREAM (Debry et al. [2007]). I
will focus on numerical points thereafter.
6000
Mode 1
Mode 2
Mode 3
Total
Distribution en nombre (n) 5000
4000
3000
2000
1000
0
2 1 0 1
10 10 10 10
DiameÁtre ( m)
Figure 3.10: Modal aerosol distribution (measured at Leipzig, Germany, from Putaud et al.
[2003]).
The third moments are computed for each chemical species in order to follow the chemical
composition in each mode.
The resulting system of (stiff) ODEs can be solved with appropriate time integration solvers
(see previous section).
A specific issue for the 3D implementation of modal methods is related to the management
of mode merging and splitting. Merging between two modes is required to force modes to be
of distinct size ranges throughout the simulations (Binkowski and Roselle [2003]; Whitby et al.
[2002]). When the combined effect of nucleation and condensation and the effect of coagulation
are of same order of magnitude but act in opposite direction, a mode may split into two modes.
This is particularly likely to happen for the nucleation mode i during high nucleation episodes.
A splitting scheme has then to be applied (Binkowski and Roselle [2003]).
Variational methods
Let us write the GDE under the following form:
dn dqi
= f (n, x, t) , = fi (n, qi , x, t) (3.23)
dt dt
x = ln m is the logarithm of the dry mass (the variable that is classically used for numerical
algorithms). If V is a functional space (of functions of x) with an appropriate scalar product
(., .), one defines the weak formulation of the GDE (for instance for the number distribution) as:
dn
∀v ∈ V, ( , v) = (f (n, x, t), v) (3.24)
dt
A numerical method is then based on a projection onto a subspace Vk of V , of finite dimension
k:
j=k
X
Πn(x, t) = nj (t)Lj (x) (3.25)
j=1
Apart from the choice of the basis functions Lj , the test functions v have also to be chosen.
In the case of a collocation method, the test functions are Dirac functions at collocation points
(xi )i=1,..., k : vi (x) = δ(x − xi ). For a basis of Lagrange polynomials, one easily gets for the
number distribution (Sandu and Borden [2003] and Debry and Sportisse [2007a]):
d~n
= A(~n)~n − diag(B~n)~n − |{z}
C~n (3.26)
dt | {z }
coagulation c/e
where (~n)j = nj . The matrices A(n) (the line i of A(n) is under the form ~nT Ai~n with Ai a
matrix) and B are related to coagulation. The matrix C is related to condensation/evaporation.
The resulting system of ODEs has then to be integrated with appropriate solvers. We refer
to Debry and Sportisse [2007a] for the benchmark of two collocation methods.
The advantage of such methods is that one can solve the GDE without using splitting
methods. However, such methods may suffer of numerical instabilities and are computationally
demanding.
Sectional methods
Another approach is based on Finite Differences. The aerosol distribution is discretized in a
fixed number of bins (let say nb ) or sections. This justifies the name of sectional method. The
bins correspond to the cells in the classical terminology of numerical analysis.
The objective is then to compute Qji and N j , the integrated values over the bin j for the
mass distribution of species i and the number distribution, respectively.
Splitting of coagulation and condensation/evaporation is often used. For the sake of clarity,
we use splitting in our presentation.
Coagulation Here, we only present the equations for the number distribution. After integra-
tion of the coagulation equation over bins, one typically gets:
j j nb
dN j 1X X j j1 j2 j
X
(t) = X j1 j2 N N − N Xjj1 N j1 (3.27)
dt 2
j1 =1 j2 =1 j1 =1
The coefficients Xjj1 j2 are the key components of the method and describe the fraction of the
aerosol generated by the coagulation of bins j1 and j2 located in bin j. Xj1 j2 is defined by
Pj=nb j
summing over all possible resulting bins: Xj1 j2 = j=1 X j1 j2 .
A popular method is based on Jacobson et al. [1994]; Stratton et al. [1994]. We refer to
Debry and Sportisse [2007b] for a review and a rigorous derivation of such equations. The key
point is to choose the closure scheme inside each bin (the a priori local form of the aerosol
distribution).
The resulting equation is then integrated with an ODE solver. A usual approach is based
on a semi-implicit approach, for instance:
∆t Pj−1 Pj j
j n+1
(N j )n + 2 j1 =1 j2 =1 Xj1 j2 (N )
j1 n+1 (N j2 )n
(N ) = P nb j
(3.28)
1+ ∆t j1 =1 (Xjj1 − Xjj 1
)(N j1 )n
viewed as an advection process). Classical methods are used in this aim (for instance, semi-
Lagrangian solvers Dhaniyala and Wexler [1995]; Dabdub and Nguyen [2002]).
An alternative approach is to use Lagrangian methods by following the evolution of bins
according to condensation/evaporation. This method is also referred as a moving sectional
algorithm (Kim and Seinfeld [1990]; Jacobson and Turco [1995]). A motivation is that one
can only use a small number of bins in 3D applications. A crucial point is then to map the
resulting moving distribution onto a fixed grid after having solved the GDE because the aerosol
distribution is given on a fixed grid for 3D processes, such as advection or diffusion. One refers
to Debry et al. [2007] for a method that ensures mass and number conservation.
From (3.19), one gets a system of ODEs that gives the time evolution of Lagrangian variables
(Debry and Sportisse [2006]). For instance for mass distributions:
nb
dQji X 0
= N j ai (mj ) Ki − Qji − ηij ceq j j
i (Q1 , . . . , Qnc ) (3.29)
dt 0
j =1
The mass conservation law for the total mass (gas and aerosol) is written under the form c i +
Pn b 0
j
0
j =1
Q i = Ki .
This system is highly stiff due to the wide range of timescales for condensation/evaporation.
We refer to the previous section for the issue of time integration. In the SIREAM model (SIze
REsolved Aerosol Model, Debry et al. [2007]), a Rosenbrock method is used. Some specific flux
limitations can also be applied in order to reduce the stiffness (especially this due to the fluxes
of H+ , Pilinis et al. [2000]).
We refer to Zhang et al. [1999] for a comparative review of algorithms.
• the modeling of Secondary Organic Aerosol, a strong weakness of current aerosol models;
• the modeling of ultrafine particles at local scales (with the competition between nucleation
and condensation, Figure 3.12);
• the modeling of non-spherical soot aerosols (see for instance Naumann [2003]);
• the extension from internal mixing to external mixing, a key point for the assessment of
radiative effects (Tombette et al. [2007], Figure 3.11).
The numerical simulation of the GDE is also still a difficult problem, especially with the increase
in complexity to arise in near future. See for instance Mallet [2007] for some strange behaviours
of the GDE.
0.8
55
La_Hague 0.7
Cabauw Belsk
Oostende IFT-Leipzig
Lille
50 Mainz
Karlsruhe
Palaiseau Paris 0.6
Laegeren
Ispra
Venise
45 0.5
Modena
Carpentras
Villefranche
Toulon
Palencia
Barcelona 0.4
40
Cabo_da_Roca
Evora
Messina
0.3
El_Arenosillo Granada
Blida
Lampedusa
35
-10 -5 0 5 10 15 20
Figure 3.11: Average value of the AOT (Aerosol Optical Thickness) at 440 nm, simulated with
the Polyphemus system for year 2005, Tombette et al. [2007].
A field of growing interest is the use of data assimilation techniques. The objective is to couple
model outputs and observational data provided by monitoring networks in order to improve the
forecasts (data assimilation) and/or to reduce the uncertainties of some input parameters (inverse
modeling, for instance of emissions). We refer for instance to Tarantola [1987]; Bennett [2002];
Kalnay [2003] for a general presentation in geophysics.
Let use write the evolution system (1.1) under the form:
dc
= f (c, Ψ) , c(0) = c0 (3.30)
dt
Here, Ψ stands for some uncertain forcing parameters (for instance emissions or boundary con-
ditions). The uncertainties are related to f (through the parameterizations), to the initial
conditions c0 and to the forcing parameters Ψ. Notice that in meteorology, the key point is
the sensitivity with respect to initial conditions, due to the chaotic behavior of atmospheric
dynamics (for instance Kalnay [2003]). In atmospheric chemistry, the sensitivity with respect to
the initial conditions is much lower due to the dissipative nature of chemical kinetics. A typical
spin-up is of magnitude one week at continental scale. Notice that in atmospheric chemistry
the eigenvalues of J = ∂f ∂c have strictly negative real parts. This is not a wave problem as in
meteorology.
Some observational data may be provided by ground stations (for instance the local mon-
itoring networks) or satellital data related to earth observation. Let us assume that some
observations obs1 , . . . , obsno are available at discrete times t1 , . . . , tno . At time tn , the model
output is c(tn ) = F n (Ψ) (we use this notation in order to underline the input/output relation).
The observation is not necessarily the model output and we write H n the observation operator
at time tn : if there were no errors, one should get obsn = H n F n (Ψ).
Data assimilation techniques are based on the reduction of discrepancies between the model
outputs and the observational data, measured by the cost function of the model parameters Ψ
(when required, it could also depend on the initial conditions):
n=n
Xo
Jo (Ψ) = (obsn − H n F n (Ψ))T Rn−1 (obsn − H n F n (Ψ)) (3.31)
| {z }
n=1
Jn
42 Chapitre 3 - Research works
3
8 8
10 10
6 6
10 10
4 4
10 10
Reference (250 ppm) Reference (250 ppm)
2 2
10 Sulfur / 5 (50 ppm) 10 Sulfur / 5 (50 ppm)
Sulfur / 10 (25 ppm) Sulfur / 10 (25 ppm)
0 0
10
3
2
1 0 1 10
3
2
1 0 1
10 10 10 10 10 10 10 10 10 10
D i am e t r e ( m i c r om e t r e s ) Diametre (micrometres)
Figure 3.12: Simulation of the impact of the reduction in sulfur composition (gasoline) in the
near exhaust of a car. Computed with the CFD model Mercure Saturne and the MAM aerosol
model. Left: winter; right: summer. In Albriet [2007].
Rn is the covariance matrix for observation error (to be specified). A second part of the cost
function can be defined by taking into account the so-called background term (an a priori
knowledge of Ψ):
Jb (Ψ) = (Ψ − Ψb )T B −1 (Ψ − Ψb ) (3.32)
with Ψb the background value and B the background error covariance matrix. Adding J b to Jo
is a penalization technique (Tikonov regularization) that may be necessary for ill-posed (high-
dimensional) problems.
The problem is then to minimize J = Jo + Jb with respect to Ψ ∈ Vadm . This provides
an optimized value Ψ? = argminJ, that can be used to compute, when required, an optimal
forecast c? = F (Ψ? ). Vadm is the space of admissible parameters (for instance, some parameters
have nonnegative values).
There are many applications of this approach in inverse modeling of emissions: for CH 4
(for instance Bergamaschi et al. [2005]), for carbon (Chang et al. [1997] at regional scale and
Bergamaschi et al. [2000] at global scale), for NOx (Quélo et al. [2005a], Figure 3.13) and VOC
(Chang et al. [1996]). Another related problem is the so-called localization problem for point
sources of accidental releases (for instance, radionuclides: see the works of my colleague Bocquet
[2005a,b]).
A brief review of inverse modeling in air quality can also be found in Bocquet and Sportisse [2007].
Two families of methods are widely used: the variational algorithms, based on a minimization
procedure (3D-Var when time evolution is not taken into account, 4D-Var when it is), and the
sequential methods, based on the Kalman filter theory.
My works have not been devoted to the development of new methods but to the application
of the existing approaches to CTM, especially the 4D-Var approach (Sportisse and Quélo [2003]
and Sportisse and Quélo [2004]).
Inverse modeling and data assimilation 43
Figure 3.13: Optimized hourly emission factors for NOx (simulation over Lille in Northern
France, Quélo et al. [2005a]).
Variational methods
The variational methods are based on gradient-like minimization methods. They require to
compute ∇Ψ J. The Finite Differences techniques are not useable due to the computational
burden and one often advocates the use of the so-called adjoint models. The CPU cost required
for the adjoint model ranges typically from 3 to 7 times the CPU cost required for the model.
Let us detail the algorithm for a CTM. Here, we assume that the CTM (actually the com-
putation of the cost function) is given by the sequence:
1. initialization of the cost function: J = 0 (or a background term);
The forcing fields φn include the meteorological data, the results of parameterizations, etc.
A direct calculation of ∇J could be performed but the storage requirement would be too
high because (φ0 , . .P
. , φno −1 ) and the trajectory (c0 , . . . , cno ) should be stored. The idea is to
notice that ∇Ψ J = n ∇Ψ Jn , and:
X ∂cm T ∂cn ∂Jn X ∂cn ∂Jn
∇Ψ J = ( )|tm−1 . . . ( n−1 )T|tn−1 ( n )T|tn + ( )T|tn−1 ( n )T|tn (3.33)
∂Ψ ∂c ∂c ∂Ψ ∂c
1≤m<n≤no 1≤n≤no
n n
The adjoint model gives, for any vector v, ( ∂c T ∂c T
∂Ψ )|tn−1 v and ( ∂cn−1 )|tn−1 v. If the differentiation is
performed with respect to the initial conditions (Ψ = c0 ), the formula is:
P 1 n
∇c0 J = 1≤n≤no ( ∂c )T . . . ( ∂c∂cn−1 )T|tn−1 ( ∂J
∂c0 |t0
n T
∂cn )|tn
h
∂J −2 T no −1
T [( ∂Jno −1 )T + ( ∂cno )T ( ∂Jno )T ]
i (3.34)
= . . . × ( ∂cnnoo−2 ) + ( ∂c∂c on −2 ) n
∂c o −1 n
∂c o −1 n
∂c o
44 Chapitre 3 - Research works
8e-06
Computation
7e-06
(d) Measurements
6e-06
4e-06
3e-06
2e-06
1e-06
0
-500 -400 -300 -200 -100 0 100 200 300 400
Y [m]
Figure 3.14: Model-to-data comparison for a wind tunnel experiment (Bugey nuclear power
plant). The model outputs are simulated with data assimilation; the controlled inputs are either
emissions and/or “physical” parameters. Observations given by ’x’ and model outputs by ’+’.
In Krysta et al. [2006].
1. initialization: (cno )? = ( ∂J no T
∂cno ) ,
A specific issue is related to the availability of the adjoint model. There are two solutions:
the adjoint model may be computed on the basis of the discretization of the continuous adjoint
model or it may be given by the adjoint code of the discretized model. We refer for instance to
Sirkes and Tziperman [1997].
In many cases, automatic differentiation with appropriate tools is performed (Giering and
Kaminski [1998]; Griewank and Walther [1999] and Mallet and Sportisse [2004]). Notice that
discontinuities exist in the CTM (especially related to aerosol thermodynamics with thresholds)
and may require a special care.
We refer for instance to Elbern et al. [1997]; Elbern and Schmidt [1999]; Elbern et al. [2000];
Elbern and Schmidt [2001] and Quélo et al. [2005a] for the application of variational methods to
air quality modeling at regional scale. Another example is given by the inverse modeling of ra-
dionuclides with a Puff model (in the context of a wind tunnel experiment, Krysta et al. [2006],
Figure 3.14).
Inverse modeling and data assimilation 45
In the linear case, when the chemical production term is given by a linear term (χ(c) = χ c,
with χ a matrix), a simplification can be used. For the sake of clarity, we consider the continuous
RT
case. The cost function is then given over an assimilation window [0, T ] as J(c) = 0 j(c(t)).
For instance, we suppose that the control parameter is the source term S. Let us define the
adjoint variable c? governed by (we omit density for the sake of clarity):
∂c?i
− − V div(c?i ) = div(K∇c?i ) + (χT c? )i − Λi c? i + ∇c j , c?i (T ) = 0 (3.35)
∂t
When the wind field is divergence free, one gets:
∂c?i
− − div(V c?i ) = div(K∇c?i ) + (χT c? )i − Λi c? i + ∇c j , c?i (T ) = 0 (3.36)
∂t
The key result is that ∇S(t) J = c? (t). The initial model can be used in order to compute the
sensitivity by modifying the input parameters and performing a backward integration (compare
(1.1) and (3.36)). The so-called retroplume methods (for instance Hourdin and Issartel [2000]
and Krysta et al. [2004]) are based on this approach.
Sequential methods
The sequential methods are based on the Kalman filter theory (for instance Tarantola [1987]).
We distinguish the forecast and the analysis obtained by taking into account the observational
data. Both are supposed to be stochastic variables that follow normal laws, N(c f , Pf ) and
N(ca , Pa ), respectively. P is a covariance error matrix. The model is supposed to be given by
cn+1 = F (cn ) from tn to tn+1 , with an unbiased model error following N(0, Q).
The most basic algorithm (BLUE for Best Linear Unbiased Estimator and the extension to
the nonlinear case with the Extended Kalman Filter) is given by the following sequence:
3. Analysis at tn+1 :
cn+1
a = cn+1
f + K n+1
obs n+1
− H n+1 n+1
c f (3.39)
Mn,n+1 is the linear tangent model from time tn to tn+1 . normal PDF.
This method is easy to implement but a key issue is the large dimension for reactive problems
and the computational burden associated to the step (3.38). A powerful technique is given by
the use of reduced filters.
46 Chapitre 3 - Research works
Notice that the covariance error matrices can be written as P f = Sf SfT , Pa = Sa SaT , Q =
Sq SqT and R = Sr SrT . It is straightforward to check that Sf = (M Sa , Sq ), K = Sf A(AT A +
R)−1 and Sa = Sf (I − A(R + AT A)−1 AT )1/2 with AT = HSf . The advantage of using the
square matrices S rather than the covariance matrices P is that these matrices are symmetric
and strictly positive. Moreover, with S = (s1 . . . sns ) (ns depends on the matrix, si is a column),
one gets:
Xns
SS T = si si T (3.42)
i=1
Once the sTi si values have been ranked in a decreasing order, it is easy to define truncated
matrices related to S.
For the RRSQRT filter (Verlaan [1998]), the modes are defined from the leading eigenvalues
of SaT Sa (that has the same eigenvalues as Pa ). For the SEEK (Singular Evolutive Extended
Kalman, Pham et al. [1998]) and the SEIK (Singular Evolutive Interpolated Kalman, Pham
[1996]) filters, the truncation is based on the Empirical Orthogonal Functions (EOF) of a com-
puted trajectory.
An alternative approach is to use an ensemble filter by computing the covariance matrix on
the basis of a Monte Carlo method:
j=N
1 X (j)
Pf ' (c − c̄)(c(j) − c̄)T (3.43)
N −1
j=1
with an ensemble of N possible states c(j) of average c̄.
We refer for instance to Segers [2002] for a benchmark of Kalman filters applied to air
pollution modeling and also to Constantinescu et al. [2006] for the application of Ensemble
Kalman filters.
The Polyphemus system (Mallet et al. [2007c]) includes 4D-Var, optimal interpolation, the
RRSQRT and EnKF methods, applied to the generic class of models available in the platform.
These methods are currently benchmarked for the simulation of air quality over Europe (Figure
3.15).
Slow/fast behaviour
The impact of the slow/fast behaviour of chemical kinetics (Section 3.6) is investigated in
Sportisse and Mallet [2007]. A theoretical framework is given in order to study the dynami-
cal issues for data assimilation. It is proven that the data assimilation problem is ill-posed for
the fast components (fast species). One also refers to Section 3.6.
Inverse modeling and data assimilation 47
Montandon
140
reference
observation
Optimal Interpolation
EnKF
120
4DVar
RRSQRT
100
Concentration
80
60
40
20
0 10 20 30 40 50
time
Figure 3.15: Comparaison of different data assimilation methods for the ozone concentration:
without assimilation or with Optimal Interpolation, RRSQRT, EnKF or 4D-Var methods. Sta-
tion of Montandon (continental simulation over Europe with Polyphemus in summer 2001).
Credit: Lin Wu, CEREA/CLIME Project, INRIA.
Second-order sensitivity
The first issue is related to the so-called second-order sensitivity. Once the optimized parameter
Ψ? (for instance emissions) has been computed, one can investigate its robustness with respect
to other uncertainties. For instance, if Φ stands for uncertain parameters that are not optimized
(due the large dimension for instance), Ψ? is defined as a function of Φ:
∂Ψ ∂∇Ψ J
= − (Hess J)−1 (3.45)
∂Φ ∂Φ
with the use of the implicit function theorem. We refer for instance to Quélo et al. [2005b] with
the application to a simple Gaussian model for the atmospheric dispersion of radionuclides.
Computing the Hessian matrix of J may be a difficult task for high-dimensional systems. A
brute-force approach is then to investigate the sensitivity of Ψ? with respect to modified model
configurations or model parameters. An example is given by the inverse modeling of NO x at re-
gional scale and the sensitivity analysis with respect to parameterizations in Quélo et al. [2005a].
Notice that this point is crucial in order to use the optimized parameters for themselves and
not only for forecast.
Network design
An interesting question is related to what is usually called network design: how to optimize a
mobile or a fixed monitoring network (related to the observational operator H) ? This issue is a
classical issue, for instance in meteorology. Mathematically speaking, the quality of a monitoring
48 Chapitre 3 - Research works
network can be estimated by the condition number of the Hessian matrix of J (that depends on
the network through H).
We refer for instance to Sportisse and Quélo [2002]. A preliminary application to air pollu-
tion modeling can also be found in Sandu [2006] (Targeting) and in Krysta et al. [2006] (Gen-
eralized Cross Validation).
• many input data are poorly known (see for instance Table 3.1), especially emissions. More-
over, some forcing fields, such as meteorological fields, are computed with numerical models
that may be uncertain.
• the numerical algorithms and the discretization also induce uncertainties (remember that
the numerical resolution may be coarse, due to the computational burden, especially for
aerosols).
It is therefore not relevant to view outputs of CTMs as deterministic values (Mallet and Sportisse [2006a]).
Even if the models are “validated” (it is more rigorous to say that model-to-data comparisons
are performed, when possible), one must keep in mind that there are a large amount of degrees
of freedom (especially in parameterizations) and only a small number of model outputs can
be measured. For instance, most of the existing CTMs have been extensively tuned to meet
acceptable model-to-data error statistics for ozone peaks at ground. It does not ensure that the
Sensitivity analysis, ucertainty propagation and ensemble forecast 49
Table 3.1: Uncertainties for some input parameters of CTM (parameters of log-normal proba-
bility density functions, following Hanna et al. [2001].)
results are satisfactory for 3D fields and many other trace species. There is therefore a danger
of using “overtuned” models, especially for impact studies or long-term scenario studies.
where Ψ stands for the uncertain input parameter, Pi are orthogonal polynomial functions
(“chaos” expansion) and N is the order of the approximation.
The coefficients (αi ) are computed at collocation points based on the PDF for Ψ (we refer
to Tatang et al. [1997]). The essence of the algorithm is to simulate the model response surface
50 Chapitre 3 - Research works
Figure 3.16: An ensemble for hourly profiles of ozone (averaged over Europe, summer 2001),
generated by the Polyphemus system. The points indicate the result of the ensemble forecast
for an ensemble algorithm. In Mallet and Sportisse [2006b].
by using orthogonal polynomials of the uncertain input parameters. Using (3.46) is much more
efficient than using the exact model.
The application of this method to simulations over Europe for uncertainties in NO x , NH3
and SO2 emissions can be found in Boutahar [2004]; Boutahar and Sportisse [2002].
3.5.3 Multi-modeling
As pointed out before, a key point is to avoid “overtuning” of models due to the large amount
of output data and the small number of observed data. For instance, a typical version of a
sectional aerosol model may include up to hundreds of species per grid cell (let say 10 bins times
10-20 chemical species per grid cell). In comparison, observational data are routinely available
for aggregated values, such as PM2 .5 or PM1 0 (PMx is the mass of aerosols whose aerodynamic
diameter is less or equal to x micrometers) but not for size distribution or chemical composition.
A promising approach is not to rely on one single model but to use a set (an ensemble) of
models or of model configurations to deliver the forecast (Figures 3.16 et 3.17). We refer for
instance to Mallet and Sportisse [2006d] for the study of model spreads due to uncertainties in
model inputs, numerical algorithms and parameterizations.
The extension to multiphase models can be found in Tombette et al. [2005] (at regional scale,
over Paris) and in Sartelet et al. [2007c] (at continental scale, in the framework of the MICS
Asia exercice).
Figure 3.17: Ensemble relative variance (in %) for ozone, summer 2001 (4 months). In
Mallet and Sportisse [2006b].
at position x). The simplest approach is to use the ensemble mean (EM in the sequel):
1 X
EMt,x = Mm,t,x (3.47)
|E|
m∈E
We refer for instance to Galmarini et al. [2004] with the application to the forecast of radionu-
clides. Another similar approach is the use of the median model.
A promising approach is to take into account observations. Let obst,x be the observation
corresponding to the model output Mt,x . An improved forecast may be given as a linear com-
bination of model outputs as: X
ELSt,x = αm Mm,t,x (3.48)
m
P P 2
where α minimize the cost function t,x [obst,x − m αm Mm,t,x ] which measures a model-
to-data discrepancy. The minimization can be performed in the Least Square sense (ELS
stands for Ensemble Least Square). The applications to air quality modeling can be found
in Mallet and Sportisse [2006b] and Pagowski et al. [2005].
A challenging point is to forecast the weights α by only using past observational data.
Some techniques based on machine learning and/or expert selection (some popular tools in
statistics) can be used. One of the simplest methods uses a so-called loss function L t (αt ) =
!2
X
αm,t Mm,t − obst . The weights can be updated with a gradient algorithm (Cesa-Bianchi
m
et al. [1996]) αt = αt−1 − η∇Lt−1 (αt−1 ). We refer for instance to Mallet and Sportisse [2006b].
10000
NH3
NO3
1000 Cl
10
0
0.001 0.01 0.1 1 10
aerosol dry diameter (µm)
Figure 3.18: Timescales for the condensation of inorganic species (ammonium, nitrate, chloride)
as a function of the aerosol size (Debry and Sportisse [2006]).
Using ensemble forecast for the input meteorological fields could be also relevant.
There are many open issues at the theoretical level: how to define spatialized ensemble
forecast ? How to couple to data assimilation methods ? etc.
ref
no solver
4 bulk fixed
bulk moving
iter fixed
0
0 0.5 1 1.5
cpu time (sec)
Figure 3.19: CPU/RMSE (relative error) diagram for different reducing methods applied to
condensation/evaporation (Debry and Sportisse [2006]).
dcs
0 = f (cf , cs ) , = g(cf , cs ) (3.50)
dt
This is the theoretical framework for the QSSA methods (Section 3.2). The key point is of
course to partition the system under the slow/fast form. Simple criteria can be derived from the
production-loss form of chemical kinetics (Djouad and Sportisse [2002a]). The use of reduced
models (that are no more stiff) is investigated in Sportisse and Djouad [2000] and Lowe and
Tomlin [2000]; Heard et al. [1998]; Kalachev and Field [2001]; Neophytou et al. [2004]; Lovas
et al. [2006] for gas-phase chemistry and Djouad et al. [2003b] for multiphase chemistry. We
also refer to Djouad and Sportisse [2003] for the numerical solution of the resulting model.
For aerosol modeling, the fast dynamics is related to condensation/evaporation. It appears
that the finest aerosols quickly reach equilibrium for mass transfer. This leads to the definition
of the so-called hybrid methods (Capaldo et al. [2000]): the smallest aerosols whose diameter
is below a cut-off diameter are supposed to be at equilibrium (they define the fast part) while
the coarsest aerosols are solved with kinetic mass transfer (they define the slow part). We refer
to Debry and Sportisse [2006] for a benchmark of different methods and the application to the
GDE (Figure 3.19).
54 Chapitre 3 - Research works
i=n
X X
F (Ψ) = F0 + Fi (Ψi ) + Fij (Ψi , Ψj )
i=1 1≤i<j≤n
X (3.51)
+ Fijk (Ψi , Ψj , Ψk ) + . . .
1≤i<j<k≤n
+F1,2,...,n (Ψ1 , Ψ2 , ..., Ψn )
Fi describe the independent action of Ψi , Fij the cooperative action of Ψi and Ψj , etc. A reduced
model is then to use a truncated development (of second order or of third order in practice).
The functions Fi are computed such that the truncation error is minimal for a given norm, over
the space of functions that satisfy the Ansatz. This procedure can be repeated for each order.
A crucial point is the choice of the error metrics. We refer to Wang et al. [1999] and
Boutahar [2004]; Boutahar and Sportisse [2002] for the application to atmospheric chemistry.
i=N
Xr
c= (c, Φi )Φi (3.52)
i=1
where the POD vectors Φi are orthonormal with respect to a given scalar product (., .). The
POD basis is defined such that for each 1 ≤ k ≤ Nr , (Φ1 . . . , Φk ) minimizes the residual Jk for
0
any subspace of dimension k, span{Φi }ki=1 :
Z T i=k
X 0
0 0 0
Jk (Φ1 , . . . , Φk ) = kc(t) − Φi , c(t) Φi k2 dt (3.53)
0 i=1
along a trajectory c(t). The norm k.k is related to the scalar product (., .). A classical result is
that the POD basis can be computed on the basis of the eigenvalues of the correlation matrix
[(c(ti ), c(tj ))]ij . The choice of the scalar product is once more a challenging point.
In Sportisse and Djouad [2007], the application to comprehensive chemical mechanism for
tropospheric chemistry indicates that there are only 2 or 3 local degrees of freedom for a species
like ozone. This result is promising but the extension to 3D modeling is still an open question.
Chemistry-Transport Models and applications 55
3.7.2 Applications
I have been implied in many applications of the developed models.
I refer at regional scale to air pollution modeling over Paris (Sartelet et al. [2002] and
Tombette and Sportisse [2007]), Lille (Quélo et al. [2004] and Quélo et al. [2005a]), Marseille
(Taghavi et al. [2004]) and Tokyo (Sartelet et al. [2007b]). I have already cited many works de-
voted to the continental scale, to which I can add Sartelet et al. [2007c] for Asia and Sartelet et al. [2007a]
for Europe (Figures 3.2 and 3.3). For both cases it was related to a benchmark exercice.
I have also been strongly implied in the model development for the forecast of the dispersion
of radionuclides. This is the subject of a joint project with the Emergency Center of the French
56 Chapitre 3 - Research works
Table 3.2: Model-to-data comparison for Polyphemus over Europe in 2001 for gaseous species:
number of stations (from the networks EMEP, Airbase and BDQA), observed mean (µg m −3 ),
simulated mean (µg m−3 ), RMSE (µg m−3 ), correlation (%), MFB (%), MFE (%). In
Sartelet et al. [2007a].
Table 3.3: Model-to-data comparison for Polyphemus over Europe in 2001 for aerosols: num-
ber of stations (from the networks EMEP, Airbase and BDQA), observed mean (µg m −3 ),
simulated mean (µg m−3 ), RMSE (µg m−3 ), correlation (%), MFB (%), MFE (%). In
Sartelet et al. [2007a].
Chemistry-Transport Models and applications 57
1986-04-26T13:03:00
Chernobyl
Figure 3.20: Simulated cesium 137 for 26 April 1986 at 3:00 (TU) in the framework of the
Chernobyl simulation (Quélo et al. [2007]).
58 Chapitre 3 - Research works
102
Polyphemus [Bq.m^{-3}]
101
0
10
1
10
2
10
3
10
1
10 3
10 2
10 100 101 102
Measurements [Bq.m^{-3}]
Figure 3.21: Model-to-data comparison (scatter plot) for cesium 137 for the first two weeks of
the Chernobyl release (Quélo et al. [2007]).
Institute of Radiological Protection and Nuclear Safety (IRSN, France). I refer for instance to
Quélo et al. [2007] (Figures 3.20 and 3.21) for the evaluation of the Polyphemus system that
constitutes the core of the coming forecast model at IRSN.
To date, the developed models are used by EDF R&D for impact studies of power plants
(for instance Tombette and Sportisse [2006] or Roustan et al. [2007], Figure 3.22), by IRSN for
its emergency center and by other partners of CEREA (EDF Polska, Meteo Chile, Technical
Center of the French Ministry for Defence/DGA, some academic teams, etc).
55 1.25
1.00
50
0.75
lat
45
0.50
0.25
40
0.00
35
-10 -5 0 5 10 15 20
-0.25
lon
Figure 3.22: Impact for the year 2005 of the emissions from the French thermal power plants
(EDF and SNET). The target is PM2.5 . The impact is a relative impact computed in % by
comparing a reference simulation and a simulation without the emissions (Roustan et al. [2007]).
60 Chapitre 3 - Research works
Conclusion
This report summarizes my research works related to air quality modeling and simulation for
the period 2000-2006.
Many points are of course still open. Some have been briefly presented and are follow-up of
my previous works. I focus here on five challenging topics for Chemistry-Transport Models.
In near future, a first trend will be the increase of the spatial resolution. At continental scale,
the current model resolution range from 10 to 100 kilometers and will go down to 1-10 kilometers.
This point is of course related to the improvement of the meso-scale meteorological models and
of the emission inventory. This will induce many new open questions at the numerical level, at
the parameterization level and at the data level. One can even wonder if a CTM will not use
meshes similar to those used for CFD (with millions of cells and unstructured adaptative grids).
The second point is related to the improvement of the “chemical” resolution. The complexity
of physical models, especially for aerosols, will increase: a better understanding of Secondary
Organic Aerosol (SOA) implies a better representation of organic species (more variables, Griffin
et al. [2002]), a coupling between organic and inorganic thermodynamics (more computational
burden, Pun et al. [2002]) and new processes (surface heterogeneous chemistry, polymeriza-
tion, etc). Apart from modeling issues, this will require the search for appropriate numerical
algorithms.
A third point is related to the couplings and the feedbacks. This concerns for instance
the coupling between radiative transfer models, CTM and meso-scale models for the radiative
behavior of the atmosphere. On-line coupling will be the standart. Another example is provided
by detailed multi-media models in the estimation of environmental impacts through multi-media
models.
One must keep in mind that the final use of such models always needs a drastic reduction
of the CPU costs. Finding the appropriate trade-off between highly detailed physical models
(the “arms race”) and tailored models that can be used in integrated modeling for impact
assessment will still remain a key issue. The current situation is that one often uses low-level
models and sometimes detailed models for such applications. Models of intermediate complexity
or algorithms for using detailed models are therefore of great interest. The fourth point is
therefore the search for reduced models to be used for integrated modeling and to be derived
with a rigourous basis.
The last point is related to uncertainties. All-in-one models, tuned to a small set of target
species, should be replaced by platforms able to generate a large set of models and model con-
figurations, not only restricted to a small set of models that are based on the same assumptions.
A key issue will be the ability of computing Probability Density Functions for most outputs.
61
62 Chapitre 3 - Research works
Bibliography
Albriet, B. (2007). Modélisation numérique des aérosols à l’échelle locale. PhD thesis, ENPC.
Babovsky, H. (1999). On a Monte Carlo scheme for Smoluchowski’s coagulation equation. Monte
Carlo Methods and Applications, 5(1):1–18.
Bartnicki, J. (2000). Nonlinear effects in the source-receptor matrices computed with the EMEP
Eulerian acid deposition model. Technical Report 4/2000, EMEP/MSC-W.
Bennett, A. (2002). Inverse Modeling of the ocean and atmosphere. Cambridge University Press.
Bergamaschi, P., Hein, R., Heimann, M., and Crutzen, P. (2000). Inverse modeling of the global
CO cycle. 1. Inversion of CO mixing ratios. J. Geophys. Res., 105(D2):1909–1927.
Bergamaschi, P., Krol, M., Dentener, F., Vermeulen, A., Meinhardt, F., Graul, R., Ramonet,
M., Peteres, W., and Dlugokencky, E. (2005). Inverse modelling of national and european
CH4 emissions using the atmospheric zoom model TM5. Atmos. Chem. Phys., 5:2431–2460.
Berkooz, G., Holmes, P., and Lumley, J. (1993). The proper orthogonal decomposition in the
analysis of turbulent flows. Annu. Rev. Fluid Mech., pages 539–575.
Berkvens, P., Botchev, M., Krol, M., Peters, W., and Verwer, J. (2002). Solving vertical trans-
port and chemistry in air pollution models. In Chock, D. and Carmichael, G., editors, Atmo-
spheric Modeling, IMA volumes in Mathematics and its interactions. Springer Verlag.
Bocquet, M. and Sportisse, B. (2007). Modélisation inverse pour la qualité de l’air : éléments
de méthodologie et exemples. Pollution Atmosphérique. Accepted for publication.
Boutahar, J., Lacour, S., Mallet, V., Quélo, D., Roustan, Y., and Sportisse, B. (2004). Devel-
opment and validation of a fully modular platform for numerical modelling of air pollution:
Polair3D. Int. J. Env. and Pollution, 22(1-2).
63
64 BIBLIOGRAPHY
Boutahar, J. and Sportisse, B. (2002). Reduction methods for atmospheric chemistry. In Global
and Regional Atmospheric Modeling Meeting, pages 133–139. Un. of Aveiro, Portugal.
Boutahar, J. and Sportisse, B. (2005). Reduction methods and uncertainty propagation: appli-
cation to a Chemistry-Transport-Model. In Proceedings of the TAM-TAM conference.
Byrne, G. and Hindmarsh, A. (1987). Stiff ODE solvers: a review of current and coming
attractions. J. Comp. Phys., 70:1–62.
Byun, D. and Lee, S. (2002). Numerical solution of trace species advection under non-
uniform density distribution: experiment with two-dimensional linear flows. In Chock, D.
and Carmichael, G., editors, Atmospheric Modeling, IMA volumes in Mathematics and its
interactions, pages 109–151. Springer.
Capaldo, K., Pilinis, C., and Pandis, S. (2000). A computationally efficient hybrid approach for
dynamic gas/aerosol transfer in air quality models. Atmos. Env., 34:3617–3627.
Carmichael, G. and al (1997). Sensitivity analysis for atmospheric chemistry models via auto-
matic differentiation. Atmos. Env., 31(3):475–489.
Cesa-Bianchi, N., Long, P. M., and Warmuth, M. K. (1996). Worst-case quadratic loss bounds
for prediction using linear functions and gradient descent. IEEE Trans. Neural Net., 7(3):604–
619.
Chang, M., Hartley, D., Cardelino, C., and Chang, W.-L. (1996). Inverse modeling of biogenic
emissions. Geophys.Res.Lett., 23:3007.
Chang, M., Hartley, D., Cardelino, C., Haas-Laursen, D., and Chang, W.-L. (1997). On using
inverse methods for resolving emissions with large spatial inhomogeneities. J. Geophys. Res.,
102(D13):16023–16036.
Constantinescu, E., Sandu, S., Chai, T., and Carmichael, G. (2006). Investigation of ensemble
based chemical data assimilation in an idealized setting. Atmos. Env. Accepted for publication.
Coron, F., Pallegoix, J., and Sportisse, B. (1994). DSMC computations of complex test cases
and industrial configuration. In 19th International Symposium on Rarefied Gas Dynamics,
Oxford.
Dabdub, D. and Nguyen, K. (2002). Semi-lagrangian flux scheme for the solution of the aerosol
condensation/ evaporation equation. Aerosol Sci. and Technol., 36:407–418.
Daescu, D., Sandu, A., and Carmichael, G. (2003). Direct and adjoint sensitivity analysis of
chemical kinetics with KPP: II Numerical validation and applications. Atmos. Env.
Debry, E., Fahey, K., Sartelet, K., Sportisse, B., and Tombette, M. (2007). A new SIze REsolved
Aerosol Model: SIREAM. Atmos. Chem. Phys., 7(6):1537–1547.
Debry, E., Jourdain, B., and Sportisse, B. (2001). Modeling of aerosol dynamics: a stochastic
algorithm. In Sportisse, B., editor, APMS 2001, Geosciences, pages 308–319. Springer.
BIBLIOGRAPHY 65
Debry, E. and Sportisse, B. (2007a). Numerical simulation of the General Dynamics Equation
(GDE) for aerosols with two collocation methods. Appl. Numer. Math., 57(8):885–898.
Debry, E. and Sportisse, B. (2007b). Solving aerosol coagulation with size-binning methods.
Appl. Numer. Math., 47:1008–1020.
Debry, E., Sportisse, B., and Jourdain, B. (2003). A stochastic approach for the numerical
simulation of the General Dynamics Equation for aerosols. J. Comp. Phys., 184:649:689.
Dhaniyala, S. and Wexler, A. (1995). Numerical schemes to model condensation and evaporation
of aerosols. Atmos. Env., 30(6):919–927.
Djouad, R., Audiffren, N., and Sportisse, B. (2003a). Sensitivity analysis using automatic
differentiation applied to a multiphase chemical mechanism. Atmos. Env., 37(22):3029:3038.
Djouad, R. and Sportisse, B. (2002a). Partitioning techniques for reduction in chemical ki-
netics. APLA: an Automatic Partitioning and Lumping Algorithm. Appl. Numer. Math.,
43(4):383:398.
Djouad, R. and Sportisse, B. (2002b). Some reduction techniques for simplifying air pollution
models. In Sportisse, B., editor, APMS 2001, Geosciences, pages 235–245. Springer.
Djouad, R. and Sportisse, B. (2003). Solving reduced models in air pollution modelling. Appl.
Numer. Math., 44(1):49:61.
Djouad, R., Sportisse, B., and Audiffren, N. (2002). Numerical simulation of aqueous-phase
atmospheric models : use of a non-autonomous Rosenbrock method. Atmos. Env., 36:873–
879.
Djouad, R., Sportisse, B., and Audiffren, N. (2003b). Reduction of multiphase atmospheric
chemistry. Journal of Atmospheric Chemistry, 46:131–157.
Domilovskii, E., Lushnikov, A., and Piskunov, V. (1978). Monte Carlo simulation of coagulation
processes. Dockl. Akad. Nauk SSSR, Ser. Phys. Chem., 240(1).
Dumas, L., Coron, F., Dupuis, D., Pallegoix, J., and Sportisse, B. (1994). Applications of DSMC
methods to aerodynamic problems in an engineering context. In Méthodes de Monte-Carlo et
applications aux gaz raréfiés. Ecole CEA-INRIA-EDF, 1994.
Eibeck, A. and Wagner, W. (2000). Stochastic particle approximations for Smoluchowski’s coag-
ulation equation. Technical report, Weierstrass-Institut for Applied Analysis and Stochastics.
Preprint No. 585.
Elbern, H. and Schmidt (1999). A four dimensional variational chemistry data assimilation
scheme for Eulerian chemistry transport modeling. J. Geophys. Res., 104:18 583–18598.
Elbern, H. and Schmidt, H. (2001). Ozone episode analysis by four dimensional variational
chemistry data assimilation. J. Geophys. Res., 106:3569–3590.
Elbern, H., Schmidt, H., and Ebel, A. (1997). Variational data assimilation for tropospheric
chemistry modeling. J. Geophys. Res., 102:15967:15985.
66 BIBLIOGRAPHY
Elbern, H., Schmidt, H., Talagrand, O., and Ebel, E. (2000). 4D variational data assimilation
with an adjoint air quality model for emission analysis. Env. Mod. Soft., 15:539–548.
Fahey, K., Debry, E., Foudhil, H., and Sportisse, B. (2005). Formulation, development and pre-
liminary validation of the SIze REsolved Aerosol Model, SIREAM. In Proceedings GLOREAM
2004, pages 7–17.
Fife, P. (1979). Mathematical aspects of reacting and diffusing systems, volume 28. Springer
Verlag.
Fiore, A. and Jacob, D. (2003). Application of empirical orthogonal functions to evaluate ozone
simulations with regional and global models. J. Geophys. Res., 108(D14).
Frohn, L., Christensen, J., and Brandt, J. (2002). Development of a high-resolution nested air
pollution model. J. Comp. Phys., 179:68–94.
Galmarini, S., Bianconi, R., Klug, W., Mikkelsen, T., Addis, R., Andronopoulos, S., Astrup,
P., Baklanov, A., Bartniki, J., Bartzis, J. C., et al. (2004). Ensemble dispersion forecasting –
part I: concept, approach and indicators. Atmos. Env., 38(28):4,607–4,617.
Gelbard, F. and Seinfeld, J. (1977). Numerical solution of the dynamic equation for particulate
systems. J. Comp. Phys., 28:357–375.
Gelbard, F. and Seinfeld, J. (1979). The General Dynamics Equation for aerosols. J. Colloid
Interface Sci., 68(2):363–382.
Giering, R. and Kaminski, T. (1998). Recipes for adjoint code construction. ACM Trans. Math.
Software, 24(4):437–474.
Gong and Cho (1993). A numerical scheme for the integration of the gas-phase chemical rate
equations in 3D atmospheric models. Atmos. Env., 27A:2591–2611.
Graf, J. and Moussiopoulos, N. (1991). Intercomparison of two models for the dispersion of
chemically reacting pollutants. Beitr.Phys.Atmosph., 64(1):13–25.
Griewank, A. and Walther, A. (1999). Revolve: An implementation of checkpoint for the reverse
or adjoint mode of computational differentiation. ACM Trans. Math. Software.
Griffin, R., Dabdub, D., and Seinfeld, J. (2002). Secondary organic aerosol 1. Atmospheric
chemical mechanism for production of molecular constituents. J. Geophys. Res., 107(D17).
Guerova, G., Bey, I., Attié, J., Cui, J., and Sprenger, M. (2006). Impact of transatlatic transport
episodes on summertime ozone in Europe. Atmos. Chem. Phys., 6(2057-2072).
Hanna, S. R., Lu, Z., Frey, H. C., Wheeler, N., Vukovich, J., Arunachalam, S., Fernau, M., and
Hansen, D. A. (2001). Uncertainties in predicted ozone concentrations due to input uncer-
tainties for the UAM-V photochemical grid model applied to the July 1995 OTAG domain.
Atmos. Env., 35(5):891–903.
He, S. and al (2000). Application of ADIFOR for air pollution model sensitivity studies. Envi-
ronmental Modeling and Software, (15):475–489.
Heard, A., Pilling, M., and Tomlin, A. (1998). Mechanism reduction techniques applied to
tropospheric chemistry. Atmos. Env., 32(6):1059–1073.
BIBLIOGRAPHY 67
Hertel, O., Berkowicz, R., and Christensen, J. (1993). Test of two numerical schemes for use in
atmospheric transport-chemistry models. Atmos. Env., 27A(16):2591–2611.
Higdon, R. and Bennett, A. (1996). Stability analysis of operator splitting for large scale ocean
modeling. J. Comp. Phys., 123:311–329.
Hourdin, F. and Issartel, J.-P. (2000). Sub-surface nuclear tests monitoring through the CTBT
133 Xe network. Geophys. Res. Lett., 27:2245–2248.
Hu, Y., Odman, M., and Russel, A. (2006). Mass conservation in the CMAQ model. Atmos.
Env., 40:1199–1204.
Hudman, R. and al (2004). Ozone production in transpacific Asian pollution plumes and impli-
cations for ozone in air quality in California. J. Geophys. Res., 109(D23S10).
Isnard, O., Krysta, M., Bocquet, M., Dubiau, P., and Sportisse, B. (2005). Data assimilation
of radionuclides atmospheric dispersion at small scale: a tool to assess the consequences of
radiological emergencies. In Proceedings of the IAEA Conference. Rio Conference.
Jacob, D. J., Logan, J., and Murti, P. (1999). Effect of rising Asian emissions on surface ozone
in the United States. Geophys. Res. Lett., 26:2175–2178.
Jacobson, M. and Turco, R. (1995). Simulating condensational growth, evaporation and coagu-
lation of aerosols using a combined moving and stationary size grid. Aerosol Sci. and Technol.,
22:73:92.
Jacobson, M., Turco, R., Jensen, E., and Toon, O. (1994). Modeling coagulation among particles
of different composition and size. Atmos. Env., 28(7):1327–1338.
Jaubertie, A., Plion, P., Sportisse, B., Aumont, B., and Toupance, G. (1998). Reduction of
kinetic schemes for atmospheric chemistry. In Proceedings of the Conference APMS98. INRIA,
France.
Jay, L., Sandu, A., Potra, F., and Carmichael, G. (1997). Improved QSSA methods for atmo-
spheric chemistry integration. SIAM J. Sci. Comp.
Kalachev, L. and Field, R. (2001). Reduction of a model describing ozone oscillations in the tro-
posphere: example of an algorithmic approach to model reduction in atmospheric chemistry.
J. Atmos. Chem., 39:65–93.
Karamchandani, P., Santos, L., Sykes, I., Zhang, Y., Tonne, C., and Seigneur, C. (2000). Devel-
opment and application of a state-of-the-science reactive plume model. Environ. Sci. Tech.,
34:870–880.
Krishnamurti, T., Kishtawal, T., LaRow, T., Bachiocchi, D., Zhang, Z., Williford, C., Gadhil, S.,
and Surendran, S. (1999). Improved weather and seasonal climate forecasts from multi-model
superensemble. Science, 285:1548–1550.
Krysta, M., Bocquet, M., Isnard, O., Issartel, J.-P., and Sportisse, B. (2004). Data assimilation
of radionuclides at short and regional scales. In Farago, I., Georgiev, K., and Havasi, A.,
editors, Proceedings of the NATO Advanced Research Workshop on Advances in Air Pollution
Modeling for Environmental security, pages 275–285. Kluiwer.
Krysta, M., Bocquet, M., Sportisse, B., and Isnard, O. (2006). Data assimilation for short-range
dispersion of radionuclides: an application to wind tunnel data. Atmos. Env., 40(38):7267–
7279.
Lacour, S., Megueulle, C., Carlotti, P., Soulhac, L., and Sportisse, B. (2003). Unsteady effects on
pollutant dispersion around a tunnel portal. In Proceedings of Urban Air Quality Conference
UAQ4, Prague.
Lacour, S. and Sportisse, B. (2007). Estimation of indoor deposition velocity for ozone with a
simplified reactive box model. Atmos. Env. In revision.
Lagache, R., Declercq, C., Quélo, D., Sportisse, B., Palmier, P., Quetelard, B., and Haziak, F.
(2006). Evaluation du PDU de Lille-Métropole sur le trafic, les concentrations de polluants
atmosphériques et la mortalité. In Actes de la 15ième conférence sur les transports et la
pollution de l’air.
Lamb, R. (1973). Note on the application of K-Theory to diffusion problems involving nonlinear
chemical reactions. Atmos. Env., 7:257–263.
Larrouturou, B. and Sportisse, B. (1997). Some mathematical and numerical aspects of reduction
in chemical kinetics. In Bristeau, M. and al, editors, Computational Science for the 21st
century, conference in honor of Professor R.Glowinski. John Wiley and Sons.
Larson, V., Golaz, J.-C., Jiang, H., and Cotton, W. (2005). Supplying local microphyisics
parameterizations with information about subgrid variability: Latin Hypercube sampling. J.
Atmos. Sci., 62:4010–4026.
Le Dimet, F. and Sportisse, B., editors (2002). Data Assimilation for Geophysical flows, INRIA-
CEA-EDF School on Applied Non-linear Problems. INRIA.
Lee, S., Yoon, S., and Byun, D. (2004). The effect of mass inconsistency of the meteorological
field generated by a common meteorological model on air quality modeling. Atmos. Env.,
38:2917–2926.
Li, Q., Jacob, D., Bey, I., Palmer, P., Duncan, B., Field, B., Martin, R., Fiore, A., Yantosca, R.,
Simmonds, P., and Oltmans, S. (2002). Transatlantic transport of pollution and its effects on
surface ozone in Europe and North America. J. Geophys. Res., 107.
BIBLIOGRAPHY 69
Louis, J. (1979). A parametric model of vertical eddy fluxes in the atmosphere. Boundary-Layer
Meteor., 17:197:202.
Lovas, T., Mastorakos, E., and Goussis, D. (2006). Reduction of the RACM scheme using
Computational Singular Perturbation analysis. J. Geophys. Res., 111(D13302).
Lowe, R. and Tomlin, A. (2000). Low-dimensional manifolds and reduced chemical models for
tropospheric chemistry simulations. Atmos. Env., 34:2425–2436.
Mallet, V. (2007). Pitfalls arising in the numerical integration of the aerosol condensational
growth. Technical report, CWI.
Mallet, V., Korssakissok, I., Quélo, D., and Sportisse, B. (2007a). Polyphemus : un système
modulaire de modélisation pour la dispersion atmosphérique et l’évaluation des risques. Pol-
lution Atmosphérique.
Mallet, V., Pourchet, A., Quélo, D., and Sportisse, B. (2007b). Investigation of some nu-
merical issues in a Chemistry-Transport Model: gas-phase simulations. J. Geophys. Res.,
112(D15301). doi:10.1029/2006JD008373.
Mallet, V., Quélo, D., Sportisse, B., Ahmed de Biasi, M., Debry, É., Korsakissok, I., Wu, L.,
Roustan, Y., Sartelet, K., Tombette, M., and Foudhil, H. (2007c). Technical Note: The air
quality modeling system Polyphemus. Atmos. Chem. Phys. Discuss., 7(3):6,459–6,486.
Mallet, V. and Sportisse, B. (2005). A comprehensive study of ozone sensitivity with respect to
emissions over Europe with a chemistry-transport model. J. Geophys. Res., 110(D22).
Mallet, V. and Sportisse, B. (2006a). Air quality modeling: from deterministic to stochastic
modeling. Computers and Mathematics with Application. In press.
Mallet, V. and Sportisse, B. (2006c). Peut-on modéliser la qualité de l’air de manière déterministe
? In Actes 38ièmes Journées de Statistiques. Société Française de Statistique.
Marchuk, G. (1986). Mathematical models in environmental problems, volume 16. North Holland.
Mott, D., Oran, E., and Van Leer, B. (2000). A Quasi-Steady-State solver for the stiff Ordinary
Differential Equations of reaction kinetics. J. Comp. Phys., 164:407–428.
Naumann, K. (2003). COSIMA: a computer program simulating the dynamics of fractal aerosols.
J. Aerosol Sci., 34:1371–1397.
Neophytou, M., Goussis, D., Van Loon, M., and Mastorakos, E. (2004). Reduced chemical
mechanisms for atmospheric pollution using Computational Singular Perturbation analysis.
Atmos. Env., 38:3661–3673.
Nguyen, K. and Dabdub, D. (2003). Development and analysis of a non-splitting solution for
three-dimensional air quality models. Atmos. Env., 37:3741–3748.
Pagowski, M., Grell, G. A., McKeen, S. A., Dévényi, D., Wilczak, J. M., Bouchet, V., Gong, W.,
McHenry, J., Peckham, S., McQueen, J., Moffet, R., and Tang, Y. (2005). A simple method
to improve ensemble-based ozone forecasts. Geophys. Res. Lett., 32.
Pandis, S. and Seinfeld, J. (1989). Sensitivity analysis of a chemical mechanism for aqueous-
phase atmospheric chemistry. J. Geophys. Res., 94(D1):1105–1126.
Petersen, A., Spee, E., Van Dop, H., and Hundsdorfer, W. (1998). An evaluation and intercom-
parison of four new advection schemes for use in global chemistry models. J. Geophys. Res.,
103.
Pham, D. (1996). A singular evolutive interpolated Kalman filter for data assimilation in
oceanography. Technical Report RT 163, LMC/IMAG.
Pham, D., Verron, J., and Roubaud, M. C. (1998). A singular evolutive extended Kalman filter
for data assimilation in oceanography. J. Marine Systems, 16:323–340.
Pilinis, C. (1990). Derivation and numerical solution of the species mass distribution equations
for multicomponent particulate systems. Atmos. Env., 24A(7):1923–1928.
Pilinis, C., Capaldo, K., Nenes, A., and Pandis, S. (2000). MADM - a new Multi-component
Aerosol Dynamic Model. Aerosol Sci. and Technol., 32:482–502.
Preussner, P. and Brand, K. (1981). Application of a semi-implicit Euler method to mass action
kinetics. Chem. Eng. Sci., 36(10):1633–1641.
Pun, B., Griffin, R., Seigneur, C., and Seinfeld, J. (2002). Secondary Organic Aerosol 2. Ther-
modynamic model for gas/particle partitioning of molecular constituents. J. Geophys. Res.,
107(D17).
Putaud, J., Dingenen, R., Baltensperger, U., Brüggemann, E., Charron, A., Facchini, M., Dece-
sari, S., Fuzzi, S., Gehrig, R., Hansson, H., Harrison, R., Jones, A., Laj, P., Lorbeer, G.,
Maenhaut, W., Mihalopoulos, N., Müller, K., Palmgren, F., Querol, X., Rodriguez, S., Schnei-
der, J., Spindler, G., Brink, H., Tunved, P., Torseth, K., Weingartner, E., Wiedensohler, A.,
Wahlin, P., and Raes, F. (2003). A european aerosol phenomenology. Technical report, Joint
Research Centre, Institute for Environment and Sustainability.
Quélo, D., Krysta, M., Bocquet, M., Isnard, O., Minier, Y., and Sportisse, B. (2007). Validation
of the polyphemus system: the ETEX, Chernobyl and Algeciras cases. Atmos. Env., 41:5300–
5315.
BIBLIOGRAPHY 71
Quélo, D., Mallet, V., and Sportisse, B. (2005a). Inverse modeling of NO x emissions at re-
gional scale over Northern France. Preliminary investigation of the second-order sensitivity.
J. Geophys. Res., 110(D24310).
Quélo, D., Sportisse, B., Berroir, J., and Charpentier, I. (2002). Some remarks concerning inverse
modeling and data assimilation for slow-fast atmospheric chemical kinetics. In Sportisse, B.,
editor, APMS 2001, Geosciences, pages 499–513. Springer.
Quélo, D., Sportisse, B., and Isnard, O. (2005b). Data assimilation for short-range dispersion of
radionuclides: a case study for second-order sensitivity. J. Environ. Radioactivity, 84:393–408.
Quélo, D., Lagache, R., and Sportisse, B. (2004). Etude de l’impact qualité de l’air des scénarios
PDUs sur Lille à l’aide du modèle de Chimie-Transport Polair3D. Technical Report 2004-28,
CEREA. Rapport PREDIT.
Ramabhadran, T., Peterson, T., and Seinfeld, J. (1976). Dynamics of aerosol coagulation and
condensation. AIChe Journal, 22(5):840–851.
Rodriguez, M. and Dabdub, D. (2003). Monte Carlo uncertainty and sensitivity analysis of the
CACM chemical mechanism. J. Geophys. Res., 108(D15).
Rouchon, P. and Sportisse, B. (September 1997). Slow and fast kinetic scheme with slow diffu-
sion. In Proceedings of the Workshop Numerical Aspects of Reduction in Chemical Kinetics,
ENPC-CERMICS.
Roustan, Y. and Bocquet, M. (2006). Sensitivity analysis for mercury over Europe. J. Geophys.
Res., 111(D14304).
Roustan, Y., Bocquet, M., Musson-Genon, L., and Sportisse, B. (2005). Modeling of mercury at
the European scale with the chemistry-transport model Polair3D. In Proceedings GLOREAM
2004, pages 121–130. Also as Technical Report CEREA 2005-04.
Roustan, Y., Bocquet, M., Musson-Genon, L., and Sportisse, B. (2006). Modélisation du mer-
cure, du plomb et du cadmium à l’échelle du continent européen. Pollution Atmosphérique,
(191):317–327.
Roustan, Y., Lasry, F., and Sportisse, B. (2007). Modélisation de l’impact des émissions des cen-
trales EDF et SNET en France sur le transport atmosphérique des PM 10 et PM2.5 . Technical
Report 2007-9, CEREA.
Russell, A. and Dennis, R. (2000). NARSTO critical review of photochemical models and
modeling. Atmos. Env., 34:2,283–2,234.
Sandu, A. (2001a). Discretizing aerosol dynamics with B-splines. Technical report, Michigan
Technological University, Houghton, MI 49931.
Sandu, A. (2001b). Positive numerical integration methods for chemical kinetic systems. J.
Comp. Phys., 170:1–14.
Sandu, A. (2006). Targeted observations for atmospheric chemistry and transport models. In
International Conference for Computational Science (ICCS-2006).
72 BIBLIOGRAPHY
Sandu, A. and Borden, C. (2003). A framework for the numerical treatment of aerosol dynamics.
Appl. Numer. Math., 45:475:497.
Sandu, A., Daescu, D., and Carmichael, G. (2003). Direct and adjoint sensitivity analysis of
chemical kinetics with KPP: I Theory and software tools. Atmos. Env.
Sandu, A., Daescu, D., Carmichael, G., and Chai, T. (2005). Adjoint sensitivity analysis of
regional air quality models. J. Comp. Phys., 204:222–252.
Sandu, A., Potra, F., Carmichael, G., and Damian, V. (1996). Efficient implementation of fully
implicit methods for atmospheric chemical kinetics. J. Comp. Phys., 129:101–110.
Sandu, A., Verwer, J., Blom, J., Spee, E., and Carmichael, G. (1997a). Benchmarking stiff ODEs
solvers for atmospheric chemistry problems II: Rosenbrock solvers. Atmos. Env., 31:3459–3472.
Sandu, A., Verwer, J., Van Loon, M., Carmichael, G., Potra, F., Dabdub, D., and Seinfeld,
J. (1997b). Benchmarking stiff ODEs solvers for atmospheric chemistry problems I: implicit
versus explicit. Atmos. Env., 31:3151–3166.
Sartelet, K., Boutahar, J., Quélo, D., Coll, I., Plion, P., and Sportisse, B. (2002). Develop-
ment and validation of a 3d chemistry-transport model, Polair3D, by comparison with data
from ESQUIF campaign. In Proceedings of the 6th Gloream workshop: Global and regional
atmospheric modelling, page 140:146. Un. Aveiro, Portugal.
Sartelet, K., Debry, E., Fahey, K., Tombette, M., Roustan, Y., and Sportisse, B. (2007a).
Simulation of aerosols and gas phase species over Europe with the Polyphemus sys-
tem. Part I: model-to-data comparison for year 2001. Atmos. Env., 41:6116–6131.
doi:10.1016/j.atmosenv.2007.04.024.
Sartelet, K., Hayami, H., and Sportisse, B. (2007b). Dominant aerosol processes during high-
pollution episodes over Greater Tokyo. J. Geophys. Res. Accepted for publication.
Sartelet, K., Hayami, H., and Sportisse, B. (2007c). MICS-Asia Phase II: sensitivity to the
aerosol module. Atmos. Env. doi:10.1016/j.atmosenv.2007.03.005.
Sartelet, K. N., Hayami, H., Albriet, B., and Sportisse, B. (2005). Development and preliminary
validation of a Modal Aerosol Model for tropospheric chemistry: MAM. Aerosol Sci. and
Technol., 40(2):118–127.
Schmidt, H., Derognat, C., Vautard, R., and Beekmann, M. (2001). A comparison of simulated
and observed ozone mixing ratios for the summer of 1998 in Western Europe. Atmos. Env.,
35:6277:6297.
Segers, A. (2002). Data assimilation in atmospheric chemistry models using Kalman filtering.
PhD thesis, Delft University.
Seigneur, C., Stephanopoulos, G., and Carr Jr, R. (1982). Dynamic sensitivity analysis of
chemical reaction systems. Chem. Eng. Sci., 37(6):845–853.
Shorter, J. and Precila, C. (1999). An efficient chemical kinetics solver using high dimensional
model representation. Journal of Physical Chemistry (1999), 103:7192–7198.
BIBLIOGRAPHY 73
Sirkes, Z. and Tziperman, E. (1997). Finite difference of adjoint or adjoint of finite difference?
Mon. Wea. Rev., 125:3373–3378.
Sommeijer, B., Van der Houwen, P., and Verwer, J. (1981). On the treatment of time-dependent
boundary conditions in splitting methods for parabolic differential equations. Int.J.for Num.
Met. in Eng., 17:335–346.
Sportisse, B., editor (1998). Air Pollution Modeling and Simulation APMS98, INRIA/ENPC
APMS Conference. INRIA.
Sportisse, B. (2000). An analysis of operator splitting techniques in the stiff case. J. Comp.
Phys., 161:140–168.
Sportisse, B., editor (2001a). Air Pollution Modeling and Simulation APMS01. Springer Verlag.
Sportisse, B. (2001b). Box models versus Eulerian models in air pollution modeling. Atmos.
Env., 35:173–178.
Sportisse, B. (2003a). Ozone des villes, ozone des champs: qualité de l’air et pollution pho-
tochimique. PCM Le Pont.
Sportisse, B. (2007c). A review of current issues in air pollution modeling and simulation.
Computational Geosciences, 11(2):159–181.
Sportisse, B. (2007d). A review of parameterizations for modeling dry deposition and scavenging
of radionuclides. Atmos. Env., (41):2683–2698.
Sportisse, B., Bencteux, G., and Plion, P. (2000). Method of Lines versus Operator Splitting
methods for reaction-diffusion systems in air pollution modelling. Env. Mod. Soft., 15(6-
7):673–679.
Sportisse, B., Boutahar, J., Debry, E., Quélo, D., and Sartelet, K. (2002). Some tracks in air
pollution modeling. RACSAM Journal of the Spanish Science Academy, Real Academia de
Ciencias de Espana, 96(3).
74 BIBLIOGRAPHY
Sportisse, B., Debry, E., Fahey, K., Roustan, Y., Sartelet, K., and Tombette, M. (2006a). PAM
project (Multiphase Air Pollution): description of the aerosol models SIREAM and MAM.
Technical Report 2006-08, CEREA. Available at http://www.enpc.fr/cerea/polyphemus/.
Sportisse, B., Debry, E., Fahey, K., Roustan, Y., Sartelet, K., and Tombette, M.
(2006b). Rapport du projet Primequal-Predit PAM (Pollution Atmosphérique Multiphasique).
Modélisation. Technical Report 2006-1, CEREA. 171 pages.
Sportisse, B. and Djouad, R. (2000). Reduction of chemical kinetics in air pollution modelling.
J.Comp.Phys., 164:354–376.
Sportisse, B. and Djouad, R. (2002). Some aspects of multi-timescales issues for the numerical
modeling of atmospheric chemistry, volume 130, chapter 4, page 39:61. IMA Volume in
Mathematics and its applications. Springer New York.
Sportisse, B. and Djouad, R. (2003). Mathematical investigation of mass transfer for atmospheric
pollutants into a fixed droplet with aqueous chemistry. J. Geophys. Res., 108(D2):4073.
Sportisse, B. and Djouad, R. (2007). Use of Proper Orthogonal Decompositions for the reduction
of atmospheric chemistry. J. Geophys. Res., 112(D06303).
Sportisse, B., Jaubertie, A., and Plion, P. (1997). Reducing mechanisms in chemical kinetics for
the simulation of reactive transport; an application to air pollution modelling. In Proceedings
of the St Venant Symposium.
Sportisse, B. and Mallet, V. (2005). Calcul Scientifique pour l’Environnement. Cours ENSTA.
86 pages.
Sportisse, B. and Mallet, V. (2007). Data assimilation and inverse modelling of slow/fast atmo-
spheric chemistry. Journal of Computational Geosciences. Submitted.
Sportisse, B. and Quélo, D. (2002). Observational network design, adaptive observations and
targeting strategies in atmospheric data assimilation: sketch of a methodological review. In
INRIA-CEA-EDF Summer School on Data Assimilation for Geophysical flows.
Sportisse, B. and Quélo, D. (2003). Data assimilation and inverse modeling of atmospheric
chemistry. Proc. of Indian National Science Academy. Part A Physical Sciences, 69.
Sportisse, B., Quélo, D., and Mallet, V. (2007). Impact of mass consistency errors for atmo-
spheric dispersion. Atmos. Env., 41:6132–6142.
Sportisse, B., Sartelet, K., Debry, E., Fahey, K., Roustan, Y., and Tombette, M. (2006c). The
SIze REsolved Aerosol Model (SIREAM) and the Modal Aerosol Model (MAM). Technical
documentation. Technical Report 2006-8, CEREA. 115 pages.
Srivastava, R., Mc Rae, D., and Odman, M. (2001). Simulation of dispersion of a power plant
plume using an adaptive grid algorithm. Atmos. Env., 35:4801–4818.
BIBLIOGRAPHY 75
Stockwell, W., Kirchner, F., Kuhn, M., and Seefeld, S. (1997). A new mechanism for Regional
Atmospheric Chemistry Modeling. J. Geophys. Res., 95(D10):16343:16367.
Strang, G. (1968). On the construction and comparison of difference schemes. SIAM J. Numer.
Anal., 5:506–517.
Stratton, D., Gans, J., and Williams, E. (1994). Coagulation algorithms with size binning. J.
Comp. Phys., 112:364–369.
Sun, P. (1996). A pseudo Non-Time Splitting Method in Air Quality modeling. J. Comp. Phys.,
(127):152–157.
Taghavi, M., Musson-Genon, L., and Sportisse, B. (2004). Impact des rejets de la centrale
thermique de Martigues sur la qualité de l’air. Modélisation avec Polair3D. Technical Report
2004-27, CEREA. Rapport pour la Mission Thermique EDF.
Tatang, M., Pan, W., Prinn, R., and McRae, G. (1997). An efficient implementation for paramet-
ric uncertainty analysis of numerical geophysical models. J. Geophys. Res., 102:21,925–21,932.
Tikhonov, A., Vasileva, A., and Sveshnikov, A. (1985). Differential equations. Springer Verlag.
Tombette, M., Chazette, P., and Sportisse, B. (2007). Simulation of the Aerosol Optical Thick-
ness over Europe with a 3D size-resolved aerosol model: Comparisons with AERONET data.
Atmos. Chem. Phys. Submitted.
Tombette, M., Fahey, K., Sartelet, K., and Sportisse, B. (2005). Aerosol modelling at regional
scale: a sensitivity study with the Polyphemus platform. In Proceedings of GLOREAM 2005.
Also as Report CEREA 2005-50.
Tombette, M. and Sportisse, B. (2006). Rapport de contrat préliminaire: évaluation des impacts
du CPT Porcheville. Technical Report 2006-55, CEREA. Rapport pour la Mission Thermique
EDF.
Tomlin, A., Berzins, M., Ware, J., Smith, J., and Pilling, M. (1997). On the use of adaptative
gridding methods for modelling chemical transport from multi-scales sources. Atmos. Env.,
31(18):2945–2959.
Troen, I. B. and Mahrt, L. (1986). A simple model of the atmospheric boundary layer; sensitivity
to surface evaporation. Boundary-Layer Meteor., 37:129–148.
Tzivion, S., Reisin, T., and Levin, Z. (1999). A numerical solution of the kinetic collection
equation using high spectral grid resolution : A proposed reference. J. Comp. Phys., 148:527–
544.
Van Loon, M. (1996). Numerical methods in smog prediction. PhD thesis, Univ. Amsterdam.
Vasileva, A., Butuzov, V., and Kalachev, L. (1994). The boundary function method for singular
perturbation problems. SIAMs monographs.
76 BIBLIOGRAPHY
Verlaan, M. (1998). Efficient Kalman filtering algorithms for hydrodynamic models. PhD thesis,
TU Delft.
Verwer, J. (1994). Gauss-Seidel iteration for stiff ODEs from chemical kinetics. SIAM J. Sci.
Comp., 15:1243–1250.
Verwer, J., Blom, J., and Hundsdorfer, W. (1996a). An implicit-explicit approach for atmo-
spheric transport-chemistry problems. Appl. Numer. Math., 20:191–209.
Verwer, J., Blom, J., Van Loon, M., and Spee, E. (1996b). A comparison of stiff ODEs solvers
for atmospheric chemistry problems. Atmos. Env., 30:49–58.
Verwer, J., Hundsdorfer, W., and Blom, J. (1998). Numerical time integration for air pollution
models. In Proceedings of the Conference APMS’98. ENPC-INRIA.
Verwer, J. and Simpson, D. (1995). Explicit methods for stiff ODEs from atmospheric chemistry.
Appl. Num. Math., 18:413–430.
Verwer, J., Spee, E., Blom, J., and Hundsdorfer, W. (1999). A second order Rosenbrock method
applied to photochemical dispersion problem. SIAM J. Sci. Comp., 20(4):1456–1480.
Verwer, J. and Sportisse, B. (1998). A note on operator splitting analysis in the stiff linear case.
Technical Report MAS-R9830, CWI.
Verwer, J. and Van Loon, M. (1994). An evaluation of explicit pseudo-steady state approximation
schemes for stiff ODE systems from chemical kinetics. J. Comp. Phys., 113:347–352.
Vila-Guerau de Arellano, J., Dosio, A., Vinuesa, J., Holtslag, A., and Galmarini, S. (2004).
The dispersion of chemically reactive species in the atmospheric boundary layer. Me-
teor.Atmos.Phys., 87(23-38).
Vila-Guerau de Arellano, J., Duynkerke, P., Jonker, P., and Builtjes, P. (1993). An observational
study on the effects of time and space averaging in photochemical models. Atmos. Env.,
27A(3):353–362.
Vinuesa, J.-F. and Vila-Guerau de Arellano, J. (2003). Fluxes and covariances of reacting scalars
in the convective boundary layer. Tellus, B 55(4):935–949.
Vinuesa, J.-F. and Vila-Guerau de Arellano, J. (2005). Introducing effective reaction rates to
account for the inefficient mixing of the convective boundary layer. Atmos. Env., 39:445–461.
Wang, S., Levy II, H., Li, G., and Rabitz, H. (1999). Fully equivalent operational models for
atmospheric chemical kinetics within global chemistry-transport models. J. Geophys. Res.,
104(D23):30417–30426.
Wang, Z., Navon, I., Le Dimet, F., and Zou, X. (1992). The second order adjoint analysis:
theory and applications. Meteor.Atmos.Phys., (50):3:20.
Whitby, E. and McMurry, P. (1997). Modal aerosol dynamics modeling. Aerosol Sci. and
Technol., 27:673–688.
BIBLIOGRAPHY 77
Whitby, E., Stratmann, F., and Wilck, M. (2002). Merging and remapping modes in a modal
aerosol dynamics models: a ”dynamic mode manager”. J. Aerosol Sci., 33:623–645.
Yang, B. and Pope, S. (1998). An investigation of the accuracy of manifold methods and
splitting schemes in the computational implementation of combustion chemistry. Combustion
and Flame, 112:16–32.
Young, T. and Boris, J. (1977). A numerical technique for solving stiff ODE associated with the
chemical kinetics of reaction flow problems. J.Phys.Chem., 81:2424.
Zhang, Y., Bischof, C., Easter, R., and Wu, P. (1998a). Sensitivity analysis of a mixed-phase
chemical mechanism using automatic differentiation. J. Geophys. Res., 103(D15):18,953–
18,979.
Zhang, Y., Bischof, C., Easter, R., and Wu, P. (2005). Sensitivity analysis of photochemical
indicators for O3 chemistry using automatic differentiation. J. Atmos. Chem., 51:1–41.
Zhang, Y., Seigneur, C., Seinfeld, J., Jacobson, M., and Binkowski, F. (1999). Simulation of
Aerosol dynamic: A comparative Review of Algorithms used in Air Quality Models. Aerosol
Sci. and Technol., 31:487–514.
Zhang, Y., Seinfeld, J., Jacobson, M., Clegg, S., and Binkowski, F. (1998b). A comparative
review of inorganic aerosol thermodynamic equilibrium modules: similarities, differences and
their likely causes. Atmos. Env.