J.-F. Müller and T. Stavrakou IASB-BIRA Avenue Circulaire 3, 1180 Brussels jfm@aeronomie.be

1. Inverse modelling of emissions based on the adjoint model technique J.-F. Müller and T. Stavrakou IASB-BIRA Avenue Circulaire 3, 1180 Brussels jfm@aeronomie.be Seminar at Harvard University, June 2nd, 2006

2. Outline • Short introduction on carbon monoxide • Adjoint-based inverse modeling: methodology • The IMAGES model used in two inversion exercises constrained by: • A) 1997 CMDL data & GOME NO2 columns • B) the 2000-2001 MOPITT CO columns • Big-region vs. grid-based inversion approach • Related work at IASB-BIRA: satellite retrievals of tropospheric gases, chemistry of terpenes • Conclusions and perspectives

3. Carbon monoxide: sources and sinks (units: Tg C/year) deposition deposition 85 30 OH OH, hv OH CH2O CO2 CO CH4 1100 570 360 CO2 340 410 deposition 100 ??? OH,O3 80 NMVOC(non-methane volatile organic compounds) 250 SOA= Secondary Organic Aerosols 200 50 700 100

4. Inversion methodology and setup The a priori emission distributions for a given species can be expressed as : where j runs over the base functions. The posterior flux estimate is given by where f is a vector of dimensionless control parameters to be optimized, so that the posterior fluxes are close enough to the prior bottom-up fluxes and the resulting abundances exhibit minimal deviation from the observed concentrations. The solution of this problem corresponds to the minimum of the cost function.

5. Inversion methodology and setup Cost function: measures the bias between the model and the observations J(f)=½Σi (Hi(f)-yi)TE-1(Hi(f)-yi) + ½ (f-fB)TB-1(f-fB) Vector of the control parameters Matrix of errors on the observations Matrix of errors on the control parameters observations Model operator acting on the control parameters 1st guess values of the control parameters For what values of f is the cost function minimal?

6. Minimizing the cost Control variables f Forward CTM Integration from t0 to t Adjoint model Integration from t to t0 Transport Adjoint transport Observations Chemistry Checkpointing Adjoint chemistry Cost function J(f) Adjoint cost function Gradient of the cost function Calculation of new parameters f with a descent algorithm Minimum of J(f) ? no yes Optimized control parameters

7. Adjoint modelling: pros and cons Calculated derivatives are exact Non-linearities (chemical feedbacks) are taken into account The emissions of different compounds can be optimized simultaneously, their chemical interactions being taken into account Computational time not dependent on the number of control variables  grid-based inversions can be addressed High computational cost: calculation of derivatives requires 3 times more CPU time than a forward model run, and on the order of 20-50 iterations are needed to attain convergence (reduction of gradient by a factor >1000) The exact estimation of posterior error is not possible within this framework; instead, iterative approximations of the inverse Hessian can be used

8. The IMAGES model • 60 chemical compounds, 5°x5° resolution, 25 σ levels (Müller and Brasseur, 1995) • Use monthly averaged meteorological fields from ECMWF analyses, impact of wind variability represented as horizontal diffusion • Semi-lagrangian transport • Anthropogenic emissions : 1997 EDGAR v3 • Biomass burning emissions : GFED (Van der Werf et al., 2003) or the POET inventory (Olivier et al., 2003) • Biogenic emissions for isoprene and monoterpenes from Guenther et al., 1995, and for CO from Müller and Brasseur, 1995 • Two main modes: (A) with or (B) without diurnal cycle calculations • Mode B (Δt=1 day) uses info. on diurnal profiles of chemical species calculated in mode A (Δt=20 min) to correct the kinetic rate constants and photorates • Inverse modeling: only in mode B (emission updates not expected to affect the diurnal behavior of chemical compounds) • 16 months simulations, including spin-up of 4 months

9. A. Big-region inversion of the 1997 CO emissions The inversion is constrained by: • NOAA/CMDL CO mixing ratios • Ground-based FTIR CO vertical column abundances • GOME tropospheric NO2 columns • Simultaneous optimization of the • total annual CO & NOx emissions • over large regions (39 flux parameters) • chemical feedbacks via the adjoint • constant seasonality of the sources • B is assumed diagonal Müller and Stavrakou,ACP, 2005

11. Impact of emission changes on OH

12. Comparison to aircraft observations

13. Estimation of posterior errors • Direct calculation of the Hessian matrix using finite differences on the adjoint model • Use of the inverse BFGS formula and the output of the minimization algorithm at each iteration • Use of the DFP update formula

14. B. Big-region vs. Grid-based inversion for optimizing CO&VOC emissions • The inversion is constrained by the MOPITT daytime CO columns from May 2000 to April 2001 • The columns and averaging kernels are binned onto the IMAGES grid and monthly averaged  total : ~ 6000 observations • Error on the column is assumed 50% of the observed value « Big-region approach »: optimize the global CO fluxes over large regions as in case A (18 variables) « Grid-based » inversion: optimize the fluxes emitted from every model grid cell by month ( ~30000 param.) seasonality and geographical distribution varied source-specific correlations among prior errors on the flux parameters  B non-diagonal • In both cases, • distinguish between anthropogenic, biomass burning and biogenic emissions Stavrakou and Müller,JGR,in press

15. Spatial correlations for anthrop. emissions En= total emission of countryn, σ En= standard error din= fraction emitted by the countrynin the ith grid cell φi= total flux emitted by the cell i = fraction of the flux emitted by the cell i and country n σEn / En = 0.6, 0.35 for industrialized countries Anm= 1, whenn=m, 0.3 if n,mbelong to the same big region, 0 otherwise Cijnm = 0.7,0.85 whenn,mbelong to the industrialized countries, 1 when i=j

16. Correlation setup for pyrogenic and biogenic emissions Spatial correlations : • Based on the geographical distancedijbetween the grid cells i and j • Relative error on the flux : 0.7for pyrogenic / 0.6for biogenic • Decorrelation length : 2000 kmfor pyrogenic / 6000 kmfor biogenic • ein: fraction of the flux emitted by the cell i and ecosystemn (n=2for pyrogenic, 40for biogenic emissions) • Cnm : 1 or 0.5depending on whether the same or different ecosystems occupy the grid cells i and j Temporal correlations :linearly varying between0and 0.5for pyrogenic emissions, between0.7and0.9for biogenic emissions

17. MOPITT column Grid-based setup Big-region setup Optimization results • Both solutions succeed in reducing the model/MOPITT bias over most regions • Larger cost reduction in the grid-based case (4.6) as compared to the big-region setup (2.2)

18. Evolution of the cost and its gradient throughout the minimization The gradient is 10 times smaller than its initial value after 6 iterations

19. The gradient is 100 times smaller than its initial value after 24 iterations

20. The gradient is 1000 times smaller than its initial value after 42 iterations

21. Anthropogenic emissions by region

22. Remarkable convergence of optimizations using either GFED or POET prior emissions • Important changes in seasonality of biomass burning emissions • Increased S. African emissions in September, reduction in June when using GFED prior GFED prior POET big-region GFED grid-based GFED grid-based POET Vegetation fire emission updates Seasonal variation Big-region setup Grid-based setup

23. prior big-region grid-based Biogenic emission updates Seasonal variation grid-based inversion • Global enhancement of biogenic VOC emissions (~ +15%) • Higher NMVOCs oxidation source by 10%

24. prior big-region grid-based Comparison to independent data (CMDL, FTIR, aircraft campaigns) prior

25. Sensitivity inversions

26. Comparison of our results to past inverse modelling studies

27. East Asian anthropogenic emissions

28. Biogenic emissions error reduction After 6 iterations (grad./10) After 24 iterations (grad./100) After 42 iterations (grad./1000)

29. Anthropogenic emissions error reduction After 6 iterations (grad./10) After 24 iterations (grad./100) After 42 iterations (grad./1000)

30. Error reduction factors over large regions (estimated using the DFP-based update)

31. IMAGES model updates (in progress) • Sigma-pressure coordinate system, 40 levels • Use of ECMWF analyses for convective fluxes, PBL diffusion • clouds • washout/rainout • KPP as alternative chemical solver (not in adjoint model calculations - well for diurnal cycle calculations) • MEGAN model for BVOC emissions • Treatment of diurnal cycle • NMVOC chemical mechanisms • Optimize horizontal diffusion coefficients using adjoint technique and output using varying winds OR get rid of these coefficients and use varying winds done in progress future

32. Related work at IASB-BIRA : satellite retrievals (M. Van Roozendael et al.) • In collaboration with KNMI, determination of NO2 tropospheric columns from satellites (AMFs, stratosphere from KNMI model) • Retrieval of CH2O columns from GOME using IMAGES profiles

33. GOME-IMAGES CH2O : 1997-2001 (Courtesy of I. De Smet & M. Van Roozendael)

34. Related work at IASB-BIRA: chemistry of terpenes • State-of-the-art mechanism development for α-pinene, based on theoretical work of J. Peeters and co-workers (Uni. Leuven) • Mechanism validation by simulations of laboratory experiments using a box model • SOA parameterization based on original vapor pressure prediction method • Reduced mechanism (~30 compounds) (work in progress) • Future: ozonolysis of α-pinene and sesquiterpenes

35. Alpha-pinene + OH quasi-explicit mechanism : Peeters et al. (2001), Fantechi et al. (2002), Vereecken and Peeters (2004), Capouet et al. (2005) Very exotic chemistry (ring closure, isomeri-sations, peroxy radical decomposition, etc.) 800 species, 2400 reactions (ozonolysis included)

36. Model simulation of laboratory experiments Lamp spectra Capouet et al., 2005

37. J. Phys Chem. A (May 2005) Related work at IASB-BIRA : unexpected reaction sequences in the UT/LS Hermans et al., 2004; 2005 : +OH CH2O CO +hv +HO2 CH2OHO2 +NO HCOOH +HO2 Also for acetone and other carbonyls!

38. Conclusions and perspectives Feasibility of multi-compound and grid-based inversions Comparable results of big-region and grid-based approaches when averaged over large regions Importance of the error correlation setup for grid-based inversions -- further work needed to better quantify the correlations Posterior uncertainty analysis made possible by the DFP approximation, shows important error reductions for large-scale fluxes (e.g. Chinese anthropogenic emissions, African biomass burning), small error reductions for individual grid cells Synergetic use of different datasets is required to better quantify emissions, in particular the CO production from the NMVOCs CH2O from satellites promising in that perspective, but large differences between retrievals by different groups intercomparisons are mandatory Also, large differences between inversion studies based on same data but different models

39. THANK YOU!