Saroja Polavarapu

1
Greenhouse gas simulation with GEM The story of
mass conservation

• Saroja Polavarapu
• Climate Research Division, Environment Canada

CRD Michael Neish, Douglas Chan, Shuzhan
Ren MRD Monique Tanguay, Claude Girard, Michel
Roch AQRD Jean de Grandpré, Sylvie Gravel
RPN seminar, 6 Feb. 2015
2
OUTLINE

• What are we trying to do and why?
• The mass conservation story of GEM for CO2
• Ensemble Kalman Filtering

3
1. The FLUX ESTIMATION PROBLEM
4
The Global Carbon Cycle
http//www.scidacreview.org/0703/html/biopilot.htm
l
1 Pg 1 Gt 1015 g
Net surface to atmosphere flux for biosphere or
ocean is a small difference between two very
large numbers
Earths crust 100,000

• The natural carbon cycle involves CO2 exchange
  between the terrestrial biosphere, oceans/lakes
  and the atmosphere.
• Fossil fuel combustion and anthropogenic land use
  are additional sources of CO2 to the atmosphere.

5
Net perturbations to global carbon budget
LeQuere et al. (2013, ESSDD)

• Based on 2002-2012
• 50 of emissions remain in atmosphere
• 25 is taken up by terrestrial biosphere
• 25 is taken up by oceans

6
Interannual variability
http//www.carboscope.eu/?qco2_budget
25 uptake by land
50 to atmosphere
25 uptake by ocean

• The uncertainty and interannual variability in
  the global CO2 uptake is mainly attributed to the
  terrestrial biosphere
• Thus, we must first learn more about biospheric
  sources/sinks

7
But what is the spatial distribution of the
fluxes and how is it changing?
Mean XCO2 Aug. 2009 GOSAT Greenhouse Gas
Observing Satellite
v2.0 averaged at 0.9x0.9
Figure courtesy of Ray Nassar, CCMR
Figure courtesy of Elton Chan, CCMR
With the increased coverage from new satellite
data, can we get flux estimates at higher spatial
resolution?
OCO-2 launch July 2014
http//oco.jpl.nasa.gov/
8
Atmospheric observations give feedback on model
forecasts

• If forecast does not match observation,
  difference could be due to errors in CO2 initial
  conditions, meteorological analyses, prescribed
  fluxes, model formulation, representativeness, or
  observation errors.

CO2 forecast
Fluxes
Forecast model
Model error
Meteorology analysis
CO2 analysis
9
The inverse problem for Carbon Flux estimation

• In flux inversions, if one solves for fluxes
  only, the transport model is needed to relate the
  flux to the observation model is a strong
  constraint
• Exact mass conservation in transport model over
  years of simulation is needed to attribute
  model-data mismatch to fluxes.
• Techniques used to solve inverse problem 4D-Var,
  EnsKF, Bayesian Inversion, Markov Chain Monte
  Carlo (MCMC)
• Perfect model assumption since forecast model is
  used as a strong constraint
• No allowance for imperfect meteorological
  analyses
• Extension for imperfect tracer initial conditions
  is not hard

flux
conc obs
Prior flux
Forecast model
Spatial interpolation
10
Conventional inverse problem setup
11
Inversions using surface network
Peylin et al. (2013)

• Inversion methods differ in
• Methodology
• Observations
• Sfc 100 flask continuous
• A priori fluxes
• Transport models
• Interannual variability is similar and due to
  land

1 5-6 2 7-8 3 9 4 10 11
12
The changing observing system
Mean XCO2 Aug. 2009 GOSAT Greenhouse Gas
Observing Satellite
v2.0 averaged at 0.9x0.9
GOSAT figure courtesy of Ray Nassar, EC

• 100 highly accurate surface stations with weekly
  or hourly data
• Regular aircraft obs over Pacific
• Satellites GOSAT (2009), OCO-2 (2014)

13
Environment Canada Carbon Assimilation System
(EC-CAS)

• An advanced, state-of-the-art assimilation system
  which will use ensemble forecasts to directly
  simulate various sources of model error.
  Comparable to systems in development for Japan,
  US, etc.
• Will be run routinely but behind real time since
  it takes time for flux to reach measurement
  locations
• Forward model global operational air quality
  model 0.9 x 0.9
• Statistical method Ensemble Kalman Smoother
  (ext. of oper.)
• Observations GAW global surface-based in-situ
  and remote sensing stations, satellite, aircraft,
  Total Carbon Column Observing Network (TCCON)
• Emissions biosphere (Canadian Terrestrial
  Ecosystem Model from CCCma), ocean, fossil fuel,
  biomass burning

13
14
EC-CAS Carbon Assimilation System
Flask, continuous, aircraft, satellite
Perturb initial conc., met fields, fluxes
Perturb obs
15
The future vision Comprehensive Carbon Data
Assimilation System
GEO Carbon Strategy Report (2010)

• Comprehensive carbon assimilation systems are
  being built by NASA, NOAA, agency-consortiums in
  Europe, Japan and EC.

16
EC-CAS team

• Saroja Polavarapu (lead) data assimilation
• Mike Neish system development
• Ray Nassar satellite observations, modeling
  anthropogenic and other emissions
• Douglas Chan carbon cycle science and modeling
• Bakr Badawy biospheric modeling
• University collaborators
• Prof. Dylan Jones (U Toronto), Feng Deng
• Prof. John Lin (U Waterloo), Myung Kim (U
  Waterloo)

17
2. MASS CONSERVATION WITH Gem
18
Operational model forecasts

• 10-day forecasts from the current operational
  model (experiment K4H1RA2B Strato2B final cycle)
  also show a steady drop in global mean surface
  pressure
• About 0.1 hPa is lost in 10 days

19
Why do we need mass conservation?

• When we assimilate greenhouse gas observations,
  mass conservation will not be possible. So no
  need for this if we want to estimate the CO2 or
  CH4 state only.
• However, we need mass conservation for
• Accurate forward simulations of CO2/CH4. With a
  good initial state and good source/sink inputs,
  can match observations.
• Estimation of sources/sinks of CO2/CH4 with
  inverse methods. Such methods will serve as a
  benchmark for the non-traditional EnKF state/flux
  estimation scheme.

20
Greenhouse gas simulation with GEM-MACH

• Need good simulation of GHG with exact mass
  conservation
• Priorities (1) CO2 global (2) CH4 global (3)
  CO2, CH4 regional
• CH4 chemistry (D.Chan, CCMR)
• stratosphere LINOZ (de Grandpré, McLinden,
  AQRD),
• troposphere
• OH climatology from CMAM (D. Plummer, CCCma)
• Start with global CO2 simulation
  

21
Model validation experiment setup

• How well can GEM-MACH simulate Carbon?
• Simulation for January 1 December 31, 2009
• Initial condition from CarbonTracker for Jan. 1,
  2009
• Meteorology surface fields (archived surface
  analyses), 3D winds, archived analyses from
  4D-Var
• Set up series of 24h forecasts, no CO2
  assimilation
• Emissions
• Every 3 hours (area type) though GEM-MACH set up
  for monthly fields with diurnal variation
• biosphere (CarbonTracker a posteriori)
• ocean (CarbonTracker a posteriori)
• Fossil Fuel (CarbonTracker but based on CDIAC)
• Biomass burning (GFED v3)
• Idea With CarbonTracker emissions and initial
  conditions, simulation should match CarbonTracker
  if transport is similar

22
Lack of global mass conservation

• Because of emissions, large gradients near the
  surface are created. The semi-Lagrangian
  advection scheme does not conserve mass.
• The poor vertical mixing of CO2 from the surface
  exacerbates the non-conservation issues.

10
20
50
Increasing KTmin reduces increase in global CO2
Petagrams C
Time
23
Hypothesis

• Hypothesis Adding fluxes at the surface creates
  large horizontal gradients. Without sufficient
  vertical mixing in the boundary layer, these
  unrealistically large gradients are smoothed by
  the semi-Lagrangian advection scheme
  (nonconservative flavour). With fast vertical
  mixing, horizontal gradients are reduced before
  the SL scheme can act.
• To prove this Check global mass before and after
  advection and before and after diffusion
• Mass change due to diffusion is machine precision
• Mass changes due to advection!

24
Further proof Turn off advection
control
4.5 Pg C 2.1 ppm too much
Annual growth is 7.7 Pg but should be 4
Petagrams C
No advection
Time
25
How did we accumulate so much CO2?
Mass change due to advection over 24 h

• The mass change due to advection over 24h is
  shown
• 7.7 Pg C per year is 0.0004 Pg C per time step
  for monotonic changes.
• For global CO2 of 818 Pg C only 0.00006 Pg C can
  be represented with 32-bits.
• The error we are looking for is only 6.7 times
  machine epsilon

Mass change (Pg C)
0.028 Pg C
Hour
26
Mass conservation Semi-Lagrangian advection
Houweling et al. (2010, ACP)

• ECMWF IFS (red) shows spurious increase of 1 ppm
  (Recall GEM spurious increase of 2.1 ppm in one
  year)
• South Pole (or Darwin) shows background CO2
  values best and better illustrates annual trend

Simulated and observed XCO2
1 ppm
IFS LMDZ TM3 TM5 obs
Pg C
27
EC-CAS version GEM v4.6.0-rc8

• Factors found to reduce spurious mass gain
• Reducing time step
• Adding horizontal diffusion to tracers consistent
  with meteorological fields
• Including convective transport of tracers
  (Zhang-McFarlane scheme)
• Adding tracer mass conservation scheme with
  global mass fixer (Bermejo-Conde)

Without tracer mass conservation scheme
With Bermejo-CondeILMC scheme
Petagrams C
Time
28
Global mass of water vapour in GEM analyses
during 2009
The total mass of the atmosphere varies mainly
due to water vapour loading.

• Trenberth and Smith (2005, J.Clim.)
• The mass of dry air is constant
• Water vapour cycle amplitude is 0.36 hPa or
  0.00037 rel. to dry air
• Here amplitude is 1900/5E60.00038

Petagrams
Time
29
Next puzzle

• Why dont we get exact mass conservation even
  after we (Monique Tanguay) implemented a tracer
  mass conservation scheme?
• Tracer variable in model
• Presently a pseudo moist mixing ratio (mass
  CO2/mass moist air). (Mixing ratio is NOT
  adjusted whenever water vapour changes, e.g.
  after physics step.)
• Tracer mixing ratios defined w.r.t. dry air is
  another way
• Observations are of mixing ratio w.r.t. dry air
  (mass CO2/mass dry air)
• Lets check global mass of CO2, air and water
  vapour at various points in the model time step.

30
Terminology used here
Mass change due to advection actually includes
any changes anywhere in the dynamics or physics
steps.
31
GEM does not conserve global moist air mass
without use of switch PSADJ
Global moist air mass

• About 2500 Pg of moist air is lost in 10 days.
• This is a relative loss of 0.0005

PSADJon
PSADJoff
Petagrams
Day in January 2009
32
Impact of switch PSADJ
Mass change due to dynamicsphysics steps

• Turning PSADJ on removed the H2O signal from
  moist air mass

2 0 -2 -4 -6 -8 -10
10 5 0 -5 -10 -15 -20
PSADJon PSADJoff
PSADJon PSADJoff
PSADJ seems to remove water! Air mass include H2O
signal
Petagrams
Ideally red curve should look like this
Moist air
Water vapour
1 3 5 7 9
1 3 5 7 9
January 2009
January 2009
33
A missing source of mass for Ps

• Claude Girard (RPN-A) determined that the changes
  in global water vapour due to physics impact
  GEMs thermodynamic equation but the flux of
  water across the Earths surface was not properly
  accounted for.
• We need to use the surface pressure to compute
  the total air mass and we assume that the water
  vapour mass is accounted for when computing
  tracer mass. If the surface pressure does not
  reflect the current mass of water vapour, global
  tracer mass calculations will not be accurate.
• Claude devised a means of adding this source of
  mass to surface pressure at the end of the
  physics time step (see his Note from August 19,
  2014). http//iweb.cmc.ec.gc.ca/armasmp/docs/mass
  -cons/Total_mass_variation_iin_GEM_girard.pdf
• Monique Tanguay (RPN-A) implemented this in GEM
  v4.7.0 and v4.6.0-rc8

34
Impact of new mass source of Ps
Mass change due to dynamicsphysics steps

• With adw_source_psTrue, moist air mass change
  from one time step to the next resembles water
  vapour mass change, as hoped

Adw_source_psTrue
Adw_source_psFalse
Moist air
Water vapour
Hours from May 1, 2009 0UTC
Hours from May 1, 2009 0UTC
35
Impact of ps_source on air mass
Air mass evolution over 10 days

• PSADJ acts on moist air
• We actually need PSADJ to act on dry air

PS_source on PS_source off
PS_source on PS_source off
Petagrams
Moist air mass
Dry air mass
Hours from Jan. 1, 2009 00Z
Hours from Jan. 1, 2009 00Z
36
New ps_source term will impact meteorological
forecasts
Expt by Monique Tanguay
2 day forecast valid 3 Jan 2009 Difference in
surface pressure due to new Ps mass source

• Differences are largest in the tropics and in
  synoptic scale systems over the ocean
• Max differences are about 3 hPa!
• There may be an impact on meteo forecasts

37
Forecast impact of ps_source NH summer/winter
(neg/neutral)
Expts by Michel Roch
GYY15 96h
SH
NH
Tropic
Control Ps_source
Summer
Winter
38
Another new flag PSADJ-dry

• For both experiments PS_source is on
• Good conservation in first 24 hours

PSADJ on moist air PSADJ on dry air
Spurious increase in dry air mass
10 day change is 100/5.1E62E-5
Petagrams
10 day change is 10/5.1E62E-6
Better but not perfect conservation of dry air
Dry air mass
Hours from Jan. 1, 2009 00Z
39
Impact of ps_dry Neutral
Expts by Michel Roch
GYY15 96h
SH
NH
Tropic
Summer
Winter
40
Need a dry mixing ratio for tracers
Change in CO2 mass

• No tracer mass conservation. Keep PSADJ-dryon,
  PS_sourceon.
• Define CO2 mixing ratio as ????CO2/??dry-air
  instead of ????CO2/??air
• Dry mixing ratio removes water signal, and CO2 is
  more constant but still need to add tracer mass
  conservation scheme for this variable.
• This was done by Monique Tanguay

CO2 mixing ratio w.r.t. moist air CO2 mixing
ratio w.r.t. dry air

• Spin-up problem
• CO2 change is 0.0001 Pg C
• Diffusion error is lt0.0000005 PgC

Petagrams
Hours from Jan. 1, 2009 00Z
41
CO2 global mass with a dry tracer mixing ratio
Dry tracer mixing ratio

• With a tracer mixing ratio w.r.t. dry air, mass
  conservation is good until October then drifts
• Accounting for local changes in air mass (surface
  pressure) when analyses are inserted every 24h
  yields exact conservation! But CO2 fields are
  terribleso the story is not over...

Petagrams C
Local offset scheme Dry tracer mixing ratio
normalized Expected mass
Time
42
Good agreement with surface obs
Dry tracer mixing ratio, global offset scheme
Without assimilation, CT fluxes
Obs GEM
Alert
Sable Island
Toronto
43
CO2 transport with GEM compared to chemistry
transport models
Column mean CO2 for 2009
The model runs use the same initial state and
fluxes
GEM
GEOS-CHEM (US academia)
CarbonTracker (NOAA)
44
Relevance of greenhouse gas modeling work to GEM

• Feedback to GEM-MACH
• https//wiki.cmc.ec.gc.ca/wiki/EC-CAS_Technology_T
  ransfer
• Vertical diffusion equation
• Emissions coding error
• Assessing vertical diffusion equation in dry
  mixing ratio
• Feedback to GEM
• Helping to test tracer conservation schemes
• Illustrated a missing source of mass for surface
  pressure
• Involved in redesign of global mean surface
  pressure fix
• Helping to convert tracer equation to dry mixing
  ratio

45
3. Ensemble kalman filter
46
Ensemble Kalman Filter first look

• No tracer assimilation, only passive advection
• Testing with 64 ensemble members, 0.9 grid
  spacing
• Start on 28 Dec 2008. Run for 4 weeks to 23 Jan
  2009
• All members have same initial CO2 and same
  fluxes. Spread is due to spread in winds only.
• Winds differ among ensemble members due to
  differences in model parameters (convection
  scheme, parameters involved in PBL model,
  diffusion of potential temperature, etc. ),
  observation error perturbations
• How does uncertainty in winds affect CO2 spread?

47
Evolution of ensemble spread
Animation of column mean CO2
Dec. 28, 2008 to Jan. 23, 2009
Ensemble mean
Ensemble spread
48
Other coupled meteorology/tracer forecast systems

• ECMWF
• Real time operational 5-day CO2 forecasts since
  2013. No assimilation of CO2 obs. Updated
  initial conditions from flux inversions every
  Jan. 1. Plans Near-real time assimilation of
  surface obs of CO2 with coupled
  meteorological/tracer assimilation
• NASA/Goddard GEOS5
• Coupled CO and CO2 assimilation to meteorological
  assimilation. Weakly couple ocean and land data
  assimilation systems to atmospheric assimilation
  system.
• Provide boundary conditions for regional
  modelling and flux inversions.
• Improve modelling of radiative transfer,
  evapotranspiration
• Feedback on modeling of boundary layer,
  convection, advection
• Provide a prioris for satellite retrievals of CO2
  and CH4

49
Future work

• EnKF development (Polavarapu/Neish)
• compare to obs without CO2 assimilation
• Extend EnKF for tracer assimilation
• Global methane simulations (D.Chan)
• Tropospheric chemistry uses CMAM OH climatology
  (D.Plummer)
• Stratospheric chemistry from Jean deGrandpré
  (LINOZ)
• OCO-2 OSSE work (Ray Nassar)
• Regional greenhouse gas simulations to support
  inverted Lagrangian trajectory work
• supports measurement network interpretation work
  of Elton Chan and Douglas Chan
• Coupling with CTEM (CCCma ecosystem model) (Bakr
  Badawy)
• Evaluate CTEM with GEM meteorology

50

• EXTRA SLIDES

51
Interannual variability 1870-2013
LeQueré et al. (ESSD, 2014)
Atmospheric accumulation has strong variability
due to land uptake. This is due to climate
variability.
Ocean uptake is not as variable.
52
Wet and dry air components of GEM analyses during
2009
This peak is due to water vapour

• Removing the water vapour from the air mass
  reduces the peak variation from 3000 to 1500
• However, long time scale variations still exist

53
Kalman Filter for coupled system

• System equations (truth)
• Equations for all 3 uncertainties also evolved,
  but not shown here

Weather forecast model
Meteorology CO2 tracers Sources (Fluxes)
Transport model
Flux evolution model
54
Kalman Filter for coupled system

• System equations (truth)
• Model errors separated into meteorology,
  transport, flux
• Model errors can be explicitly modelled and
  accounted for in actual forecast step (from
  analysis)

Weather forecast model error
Meteorology CO2 tracers Sources (Fluxes)
Transport model error
Flux evolution model error
55
Kalman Filter for coupled system

• Analysis step
• Meteorology and tracer analysis steps decoupled,
  initially
• CO2 obs used to update CO2 and flux estimates
• Analysis step includes observation and
  representativeness errors
• Equations for uncertainties not shown

Meteorology step already done operationally CO2
tracers Sources (Fluxes)
GHG obs
Weight matrix based on error models
56
EC-CAS mass conservation progress in the past
year
Polavarapu, Neish, deGrandpre, Gravel
May 2013 July 2013 Nov 2013 Jan 2014 July 2014

• Latest run has
• Surface fields from archives with Liebman
  filtering
• Physics recycling
• Half time step 900s
• Hor diffusion of tracers
• Bermejo-CondeILMC

GEM 4.5.0-a10 (7.2)
GEM 4.5.0-b7 (5.5)
GEM 4.5.0 (4.8)
(3.7)
GEM 4.5.0mass cons
GEM 4.6.0-rc8 (3.3)
Expected mass (3.3)
By the end of one year, the mass gain is exactly
right!
57
Change in global mass due to advection
24h run of GEM-MACH

• Even with the flux boundary condition the global
  mass changes during 24 h
• The mass change due to advection is shown here

ExperimentfluxBCtry2
58
Change in global mass due to diffusion
24h run of GEM-MACH

• The mass change due to the diffusion step is
  isolated here
• The maximum error of 0.0008 Pg C is 13 times
  machine epsilon
• The typical error of 0.0002 Pg C is 3.3 times
  machine epsilon
• This plot reflects only the precision of the
  diagnostic calculation. No conclusions about
  model behaviour can be drawn.

Jump of machine epsilon
ExperimentfluxBCtry2
59
Impact of greenhouse gas modeling work on GEM-MACH

• Using a model in a different way often leads to
  insights/feedback on the model
• GEM-MACH was designed for air quality forecasts
  regional scale up to 2 days, or global scale up
  to 10 days.
• But we used it for greenhouse gas simulations.
  Unlike the air quality problem, (1) we have no
  reactive chemistry to hide behind, and (2) we are
  looking at long time scales.
• We saw problems because of this different focus
  and we were able to provide feedback to GEM-MACH
• https//wiki.cmc.ec.gc.ca/wiki/EC-CAS_Technology_T
  ransfer
• Vertical diffusion equation
• Emissions coding error
• Assessing vertical diffusion equation in dry
  mixing ratio (in progress)