Recent Developments in Computational Tools for Chemical Data Assimilation

1
Recent Developments in Computational Tools for
Chemical Data Assimilation

• A. Sandu Virginia Tech
• G.R. Carmichael University of Iowa
• J.H. Seinfeld Caltech
• D. Daescu Portland State University
• NASA, NOAA, EPA for Trace-P and ICARTT data
• NSF CAREER ACI-0093139
• NSF ITR APIM-0205198

2
Prediction is difficult

• There is no reason anyone would want a computer
  in their home.
• Ken Olson, President, Chairman, and Founder of
  Digital Equipment Corp. (DEC), 1977.

3
Information feedback loops between CTMs and
observations data assimilation and targeted meas.
Transport Meteorology
Optimal analysis state
Chemical kinetics
Observations
CTM
Data Assimilation
Aerosols
Targeted Observ.

• Improved
• forecasts
• science
• field experiment design
• models
• emission estimates

Emissions
4
In the 4D-Var approach D.A. is formulated as an
optimization problem (gradient based algorithms)
5
Adjoints of stiff chemical kinetics formulation,
challenges, and automatic implementation
Transport Meteorology
Optimal analysis state
Chemical kinetics
Observations
CTM
Data Assimilation
Aerosols
Targeted Observ.

• Improved
• forecasts
• science
• field experiment design
• models
• emission estimates

Emissions
6
KPP automatically generates simulation and
direct/adjoint sensitivity code for chemistry
Simulation code
Chemical mechanism
SUBROUTINE FunVar ( V, F, RCT, DV ) INCLUDE
'small.h' REAL8 V(NVAR), F(NFIX)
REAL8 RCT(NREACT), DV(NVAR)C A - rate for each
equation
REAL8 A(NREACT)C Computation of equation
rates
A(1) RCT(1)F(2) A(2) RCT(2)V(2)F(2)
A(3) RCT(3)V(3) A(4)
RCT(4)V(2)V(3) A(5) RCT(5)V(3)
A(6) RCT(6)V(1)F(1) A(7)
RCT(7)V(1)V(3) A(8) RCT(8)V(3)V(4)
A(9) RCT(9)V(2)V(5) A(10)
RCT(10)V(5)C Aggregate function
DV(1)
A(5)-A(6)-A(7) DV(2) 2A(1)-A(2)A(3)-A(4)
A(6)-A(9)A(10) DV(3)
A(2)-A(3)-A(4)-A(5)-A(7)-A(8) DV(4)
-A(8)A(9)A(10) DV(5) A(8)-A(9)-A(10)
END
INCLUDE atoms DEFVAR O O O1D O O3
O O O NO N O NO2 N O
O DEFFIX O2 O O M ignore
EQUATIONS Small Stratospheric O2 hv
2O 2.6E-10S O O2 O3
8.0E-17 O3 hv O O2
6.1E-04S O O3 2O2
1.5E-15 O3 hv O1D O2 1.0E-03S
O1D M O M 7.1E-11 O1D O3
2O2 1.2E-10 NO O3 NO2
O2 6.0E-15 NO2 O NO O2
1.0E-11 NO2 hv NO O 1.2E-02S
Damian et.al., 1996 Sandu et.al., 2002
7
Sparse Jacobians , Hessians, and sparse linear
algebra routines are automatically generated by
KPP
JACOBIAN ON OFF SPARSE JacVar(),
JacVarTR_SP_Vec() KppDecomp(), KppSolveTR()
HESSIAN ON OFF HessVar(),
HessVarTR_Vec()
SAPRC-99. NZ 848x2 (0.2)
SAPRC-99 79 spc./211 react. NZ839, NZLU920
8
Runge-Kutta methods and their adjoints are well
suited for inverse chemical kinetic problems
RK Method
Continuous Adjoint
Discrete Adjoint Hager, 2000
Consistency The discrete adjoint of RK method of
order p is an order p discretization of the
adjoint equation. (Proof using elementary
differentials of transfer functions). Note BDF
adjoints with variable step are not consistent
with continuous adjoint equation. Sandu, Daescu,
2003
9
Singular perturbation analysis is relevant for
studying the method behavior on stiff systems
Singular perturbation test problem
Distinguish between derivatives w.r.t.
stiff/non-stiff variables

• If RK with invertible coefficient matrix A
  and R(8) 0
• and the cost function depends only on
  the non-stiff variable y
• Then ?z 0 and ?y are solved with the same
  accuracy as the original
• method, within
  O(e).
• A similar conclusion holds for continuous
  RK adjoints.
• Sandu, Daescu, 2003

10
Rosenbrock methods and their adjoints are
efficiently implemented by KPP
Rosenbrock Method (Tfwd)
Continuous Adjoint (T1.2Tfwd)
Discrete Adjoint (T2.3Tfwd)
Sandu et.al., 2002
11
Second order adjoints provide Hessian-vector
products useful in optimization and analysis
RK TLM Methods (KPP)
First and Second Order RK Discrete
Adjoints (KPP)
Sandu et. al., 2005
12
Adjoints for Integral-PDE aerosol dynamic
equations formulation and challenges
Transport Meteorology
Optimal analysis state
Chemical kinetics
Observations
CTM
Data Assimilation
Aerosols
Targeted Observ.

• Improved
• forecasts
• science
• field experiment design
• models
• emission estimates

Emissions
13
Populations of aerosols (particles in the
atmosphere) are described by their mass density
Aerosol dynamic equation - IPDE
14
Adjoint aerosol dynamic models are needed to
solve inverse problems
Continuous adjoint equation
Observations of density in each bin allow the
recovery of initial distribution and of parameters
Sandu et. al., 2005 Henze et. al., 2004
15
Discrete adjoint models for numerical advection
formulation and challenges
Optimal analysis state
Transport Meteorology
Chemical kinetics
Observations
CTM
Data Assimilation
Aerosols
Targeted Observ.

• Improved
• forecasts
• science
• field experiment design
• models
• emission estimates

Emissions
16
Inconsistency, instability of discrete adjoint
advection schemes when the forward scheme pattern
changes
upwind direction change
discretization change
Forward flux-limiter Status change
Von Neumann stability regions for different
upwindings
CFLlt0.87
unstable
CFLlt0.81
17
Loss-of-information hinders the optimization
process

Instability and the inconsistency do not seem to
affect the optimization for small Courant numbers

• Loss of information hinders the optimization
  process
• solution collapses into a sink
• exits the domain through an outflow boundary
• propagates too slow to reach one of the
  observation sites

Sandu et. al., 2005
18
The 4D-Var tools have been implemented in the
parallel adjoint STEM and are being applied to
real data
Transport Meteorology
Optimal analysis state
Chemical kinetics
Observations
CTM
4D-Var Data Assimilation
Aerosols
Targeted Observ.

• Improved
• forecasts
• science
• field experiment design
• models
• emission estimates

Emissions
19
High performance computing parallel adjoint
STEM-3
Distributed, 2-level checkpointing scheme

• Chemistry KPP
• Forward sparse Rosenbrock, RK, LMM
• DDM sensitivity (Rosenbrock, LMM)
• Discrete adjoints Rosenbrock
• Continuous adjoints Rosenbrock, RK, LMM
• Aerosols 0-D, not yet 3-D
• Transport
• Forward upwind FV, FD, FE
• Adjoints for linear upwind FD
• Parallelization with PAQMSG

Sandu et.al., 2003, 2004 Carmichael et. al.,
2003, 2004
20
Assimilation of P3-B observations taken during
the Trace-P campaign on March 07, 2001
P3-B observations O3 (8) NO, NO2 (20)
HNO3, PAN, RNO3
(100) Control Init. conc. 50 species
(sensitivity) Optimization problem with 5
million variables Parallel adjoint STEM, 2-level
checkpointing
NO2
O3
Chai et.al., 2004
21
Assimilation of AIRNOW O3 surface observations
for July 20, 2004 adjusts predictions
considerably
Observations circles, color coded by O3 mixing
ratio
Surface O3 (forecast)
Surface O3 (analysis)
22
Ensemble-based chemical data assimilation can
complement variational techniques
Transport Meteorology
Optimal analysis state
Chemical kinetics
Observations
CTM
Ensemble Data Assimilation
Aerosols
Targeted Observ.

• Improved
• forecasts
• science
• field experiment design
• models
• emission estimates

Emissions
23
AR model of background errors accounts for
flow-dependent correlations and is inexpensive

• Background error representation considerably
  impacts assimilation results
• Typically estimated empirically from multiple
  model runs (NMC)
• Correct mathematical models of background
  errors are of great interest

McMillian Reservoir, DC

• AR model of background errors
• N ?t lifetime of the species
• B
• flow dependent
• inexpensive to compute
• has full rank

24
Singular vectors characterize the directions of
maximum error growth (under different measures)

• Measures ECI,0 ? TESV EI, CPa-1 ?
  HSV
• Understand areas of maximum sensitivity to
  uncertainty
• Construct the subspace of initial perturbations
  for EnKF

Cause stiff transient Solution TLM ADJ chem.
projections

• Liao and Sandu, 2005

25
TESVs are shaped by both meteorology and
chemistry, as seen for different sections
NO2
O3
s2 Trace-P, March 1-2, 2001
1
2
26
Ensemble-based chemical data assimilation
considerably improves the solution
NO2 original error
O3 original error
O3 error w/ analysis
NO2 error w/ analysis
24h, March 1-2, 2001 50 members Perturbed I.C.,
B.C., and emissions AR TESV O3 NO2 obs (red)

27
Results for species not directly observed also
show marked improvements after assimilation
PAN original error
HCHO original error
CO original error
PAN error w/ analysis
HCHO error w/ analysis
CO error w/ analysis
28
The model can be used to place the observations
in the locations of maximum informational benefit
Transport Meteorology
Optimal analysis state
Chemical kinetics
Observations
CTM
Data Assimilation
Targeted Observations
Aerosols

• Improved
• forecasts
• science
• field experiment design
• models
• emission estimates

Emissions
29
Targeted observations are deployed to areas where
they provide maximum of information
Verification Korea, ground O3 0 GMT, Mar/4/2001
(criterion based on SVs)
O3
HCHO
NO2

• Liao and Sandu, 2005

30
Dynamic integration of chemical data and
atmospheric models is an important, growing field

• Assimilation of chemical observations into CTMs
• improves reanalysis of fields
• model forecast skills
• provides top-down estimate of emission
  inventories
• Current state of the art
• the tools needed for 4d-Var chemical data
  assimilation are in place
• adjoints for stiff systems, aerosols, transport
  singular vectors
• parallelization and multi-level checkpointing
  schemes
• models of background errors
• their strengths demonstrated using real (field
  campaign) data
• ambitious science projects are ongoing
• Emerging needs and research directions
• hybrid methods (combining ensemble and
  variational approaches)
• second order adjoints and optimization
• reduced order models (POD, singular vectors)
• algorithms for targeted observations

31
Quote of the day

• Persons pretending to forecast the future shall
  be considered disorderly under section 901(3) of
  the criminal code and liable to a fine of 250
  and/or 6 months in prison.
• Section 889, New York State Code of Criminal
  Procedure (after M.D. Webster)