Stratospheric chemistry

From Geos-chem
Jump to navigation Jump to search

On this page we list information about the stratospheric chemistry mechanism in GEOS-Chem.

For the influence of stratospheric composition on tropospheric photolysis rates, see the Fast-J page.

Overview

Though GEOS-Chem is predominantly used for simulating tropospheric chemistry, it is necessary to deal with species that may have net stratospheric sources and transport to the troposphere (e.g., NOx), or net stratospheric sinks of gases transported from the troposphere (e.g., CO). An updated linearized stratospheric chemistry scheme has been prepared for implementation into GEOS-Chem v9-01-03 to complement the linearized ozone chemistry afforded by Linoz.

Authors and collaborators:

  • Lee Murray (Harvard)

New implementation of linearized stratospheric chemistry

An updated linearized stratospheric chemistry scheme has been implementated into GEOS-Chem v9-01-03o. It currently works with GCAP/GEOS4/GEOS5/MERRA/GEOS5.7 simulations at all resolutions, for full chemistry, SOA, dicarbonyl, isoprene, H2/HD, tagged Ox, and tagged CO.

One may turn on the linearized stratospheric chemistry scheme in the chemistry menu of the input.geos file

        %%% CHEMISTRY MENU %%%  : 
        Turn on Chemistry?      : T
        Use linear. strat. chem?: T
         => Use Linoz for O3?   : T
        Chemistry Timestep [min]: 60
        Read and save CSPEC_FULL: T
        USE solver coded by KPP : F

At the beginning of a new month, the model reads in the archived 3D monthly mean production rate, P (v/v/s), and loss frequency, k (s-1), for each species. The new concentration at the end of each chemistry time step is then determined by solving

        ! Simple analytic solution to dM/dt = P - kM over [0,t]
        if ( k .gt. 0d0 ) then
           STT(I,J,L,N) = M0 * exp(-k*t) + (P/k)*(1d0-exp(-k*t))
        else
           STT(I,J,L,N) = M0 + P*t
        endif

This is applied for each pertinent species to every grid box for which the tropospheric chemistry solver is not applied in a given chemistry time step.

Also, monthly mean stratospheric OH fields are read in and used to generated decay for certain species (currently, CH3Br, CHBr3, CH2Br2).

The bromine species are hard-wired to be ignored for now. Bry in the stratosphere is prescribed using monthly mean nighttime and daytime concentrations from the GEOS CCM as in its original (v9-01-03m) implementation. For more information, see the bromine chemistry mechanism wiki page.

Note: Linoz is still the recommended option for stratospheric ozone, and ozone prod/loss rates were only implemented here for completeness. If one turns off Linoz, then Synoz is applied in the code.

Preparation of monthly mean 3D prod/loss rates

The production rates (v/v/s), loss frequencies (s-1) and concentrations (v/v) have been archived from the GMI "Combo" model simulations, which one can think of as having the GEOS-Chem troposphere + the GMI full stratospheric chemistry model. For GEOS4/GEOS-Chem simulations, long-term monthly means were calculated from the "Aura4" GMI monthly concentrations, which used GEOS-4 assimilated met fields for 2004-2006. For GEOS5/ or MERRA/GEOS-Chem simulations, long-term monthly means were determined from the "MERRA" GMI rerun (containing the Synoz bug fix in the GMI model) using MERRA met fields for 2004-2010.

The GMI model outputs as diagnostic the monthly mean rate of each kinetic or photolysis reaction. That output was aggregated to get the total production rate (v/v/s). Photolysis frequencies are also available, but kinetic loss frequencies had to be estimated from the burden/loss rate, to determine the mean loss frequencies (s-1).

Species with linearized stratospheric chemistry

Numbers indicate number of GMI Combo model reactions aggregated, see [reaction details document] for more information

Species Name Photolysis Prod Kinetic Prod Photolysis Loss Kinetic Loss
NOx 14 39 33
Ox 9 8 37
PAN 1 1 1
CO 16 20
ALK4 2
ISOP 3
HNO3 38 1 2
H2O2 4 1 2
ACET 1 7 1 1
MEK 1 18 1 2
ALD2 4 26 2 2
RCHO 9 33 1 2
MVK 1 7 3 3
MACR 1 7 2 4
PMN 1 3
PPN 1 1
R4N2 4 1 1
PRPE 1 1 3
C3H8 2
CH2O 10 70 2 5
C2H6 2
N2O5 1 2 4
HNO4 1 2 2
MP 1 1 2
GLYX 5 3 2
MGLY 4 23 2 2
GLYC 2 10 1 1
HAC 2 12 1 1
GPAN 1 1
ACTA 20 1
RIP 2 1 1
MAP 1 1 1
CH4 1 5

Additionally available species

The following species from the GMI model have also had their monthly mean production rates (v/v/s), loss frequencies (s-1), and concentrations (v/v) archived for GEOS-Chem. The linearized chemistry will automatically have prod/loss applied if a tracer sharing one of the following names is introduced into the model's STT array of transported tracers and DO_STRAT_CHEM is called.

H, H2, HCOOH, H2O, HO2, MO2, MOH, N, N2O, O, O1D, OH, Br, BrCl, BrO, BrONO2, HBr, HOBr, Cl, Cl2, ClO, Cl2O2, ClONO2, HCl, HOCl, OClO, CH3Br, CH3Cl, CH3CCl3, CCl4, CFCl3, CF2Cl2, CFC113, CFC114, CFC115, HCFC22, HCFC141b, HCFC142b, CF2Br2, CF2ClBr, CF3Br, H24O2, A3O2, ATO2, B3O2, EOH, ETO2, ETP, GCO3, GP, IALD, IAO2, IAP, INO2, INPN, ISN1, ISNP, KO2, MAN2, MAO3, MAOP, MCO3, MRO2, MRP, MVN2, PO2, PP, PRN1, PRPN, R4N1, R4O2, R4P, RA3P, RB3P, RCO3, RCOOH, RIO1, RIO2, ROH, RP, VRO2, VRP.

The bromine species are hard-wired to be ignored for now. Bry in the stratosphere is prescribed using monthly mean nighttime and daytime concentrations from the GEOS CCM, located in the bromine_201205/CCM_stratosphere_Bry data directory.

Also note, the input file may be accessed from anywhere to obtain monthly mean concentrations (e.g., 4D stratospheric OH) for use within the model.

Code structure

All stratospheric-relevant chemistry routines (including Linoz) have been consolidated into a single driving routine DO_STRAT_CHEM in GeosCore/strat_chem_mod.F90, called from DO_CHEMISTRY in GeosCore/chemistry_mod.F for full chemistry, tagged Ox and H2/HD simulations. The tagged CO simulation has been updated in GeosCore/tagged_co_mod.F.

  • GeosCore/strat_chem_mod.F: Source code file with stratospheric chemistry. The Linoz routines are called from here.
  • GEOS_NATIVE/strat_chem_201206/gmi.clim.*.nc: Species-specific vertical and horizontal resolution dependent files containing mean production rates and loss frequencies from the GMI combo model in the primary data directory

Note, the input file is in the netCDF file format, and therefore also requires GEOS-Chem be compiled with netCDF.

Original implementation

The original stratospheric chemistry scheme through GEOS-Chem v9-01-02 (v9-01-03n) consisted of:

  • A fixed stratosphere-to-tropophere flux of ozone from Synoz (an alternative and recommended linearized ozone chemistry Linoz option introduced in v8-02-04)

And performed a simple linearized chemistry scheme from archived monthly mean production rates, J-values, and OH fields from the 2D stratospheric chemistry model of Hans Schneider and Dylan Jones using GEOS-1 and/or GEOS-STRAT met fields (36 vertical levels) that have since been regridded to subsequent GEOS grid vertical resolutions.

  • Fixed monthly zonal mean production rates of NOy in pnoy_200106
  • Fixed monthly zonal mean production rates and loss frequencies of CO in pco_lco_200203
  • Loss of selected species by reaction with archived monthly zonal mean stratospheric OH fields in stratOH_200203
  • Loss of selected species by photolysis using archived monthly zonal mean J-values archived in stratjv_200203

These were variously performed by the subroutines

  • upbdflx_mod.f
    • DO_UPBDFLX - Called from main.f before transport, calls Synoz/Linoz ozone.
    • UPBDFLX_NOY - Called first from DO_UPBDFLX for NOy production, and then again after transport in main.f to repartition NOy into NOx and HNO3.
  • schem.f - Called just after SMVGEAR/KPP solves tropospheric chemistry in chemdr.f. Performs
    • Photolysis loss
    • CO prod/loss
    • Reaction with OH

Original Species

Species Photolysis Loss Loss by OH Special Treatment
NOx pnoy_200106 (P&L)
Ox Synoz (P) or Linoz (P&L)
PAN X
CO pco_lco_200203 (P&L)
ALK4 X
ISOP X
HNO3 pnoy_200106 (P&L)
H2O2 X X
ACET X X
MEK X
ALD2 X X
RCHO X X
MVK X X
MACR X X
PMN X
R4N2 X X
PRPE X
C3H8 X
CH2O X X
C2H6 X
N2O5 X
HNO4 X X
MP X X

Additional Documentation

Validation

Coming soon, with STE fluxes.

References

Murray, L. T., J. A. Logan, and D. J. Jacob, Interannual variability in tropical tropospheric ozone and OH: the role of lightning, J. Geophys. Res., 118, 1-13, doi:10.1002/jgrd.50857, 2013.

Murray, L. T., D. J. Jacob, J. A. Logan, R. C. Hudman, and W. J. Koshak, Optimized regional and interannual variability of lightning in a global chemical transport model constrained by LIS/OTD satellite data, J. Geophys. Res., 117, D20307, doi:10.1029/2012JD017934, 2012.

Previous issues that are now resolved

Reduce memory footprint of the stratospheric chemistry module

We have made a modification that should reduce the amount of memory that the stratospheric chemistry module requires.

Bob Yantosca wrote:

When I ran the 2° x 2.5° standard code (compiled with MET=geosfp GRID=2x25) in our regular 8-CPU computational queue, the run died with a segmentation fault error. This error is indicative of the code using more memory than is available. The seg fault occurred in the transport routine TPCORE_FVDAS, so for a while I thought the issue was there. But when I turned off the strat chem, the run proceeded to completion. So this meant that the problem was in the strat_chem_mod.F90 module.
Long story short: what was happening is that the code was allocating BIG arrays in strat_chem_mod.F90 at the start of the run, and then not enough memory was left by the time you got to transport routines. To fix this problem, I declared the following arrays as REAL*4 instead of REAL*8:
   !------------------------------------------------------------------------------
   ! Prior to 4/3/14:
   ! Change arrays from REAL*8 to REAL*4 to reduce memory footprint (bmy, 4/3/14)
   !  REAL*8,  ALLOCATABLE, TARGET :: PROD(:,:,:,:)   ! Production rate [v/v/s]
   !  REAL*8,  ALLOCATABLE, TARGET :: LOSS(:,:,:,:)   ! Loss frequency [s-1]
   !  REAL*8,  ALLOCATABLE, TARGET :: STRAT_OH(:,:,:) ! Monthly mean OH [v/v]
   !------------------------------------------------------------------------------
     REAL*4,  ALLOCATABLE, TARGET :: PROD(:,:,:,:)   ! Production rate [v/v/s]
     REAL*4,  ALLOCATABLE, TARGET :: LOSS(:,:,:,:)   ! Loss frequency [s-1]
     REAL*4,  ALLOCATABLE, TARGET :: STRAT_OH(:,:,:) ! Monthly mean OH [v/v]

     . . .

   !----------------------------------------------------------------------------
   ! Prior to 4/3/14:
   ! Change arrays from REAL*8 to REAL*4 to reduce memory footprint (bmy, 4/3/14)
   !  REAL*8, ALLOCATABLE  :: MInit(:,:,:,:)      ! Init. atm. state for STE period
   !  REAL*8, ALLOCATABLE  :: SChem_Tend(:,:,:,:) ! Stratospheric chemical tendency
   !----------------------------------------------------------------------------
     REAL*4, ALLOCATABLE  :: MInit(:,:,:,:)      ! Init. atm. state for STE period
     REAL*4, ALLOCATABLE  :: SChem_Tend(:,:,:,:) ! Stratospheric chemical tendency
                                                 !   (total P - L) [kg period-1]
and then I had to make the correesponding changes in module routines GET_RATES:
   !------------------------------------------------------------------------------
   ! Prior to 4/3/14:
   ! We no longer need this (bmy, 4/3/14)
   !    REAL*8             :: ARRAY2( IIPAR, JJPAR, LLPAR )  ! Actual vertical res
   !------------------------------------------------------------------------------
       . . .
 
       ! Intialize arrays
   !------------------------------------------------------------------------------
   ! Prior to 4/3/14:
   ! Initialize w/ 0e0 for REAL*4 (bmy, 4/3/14)
   !    LOSS = 0d0
   !    PROD = 0d0
   !------------------------------------------------------------------------------
       LOSS = 0e0
       PROD = 0e0

       . . .

   !-----------------------------------------------------------------------------
   ! Prior to 4/3/14:
   ! PROD is now REAL*4 (bmy, 4/13/14)
   !       ! Cast from REAL*4 to REAL*8 and resize to 1:LLPAR
   !       call transfer_3D( array, array2 )
   !
   !       PROD(:,:,:,N) = ARRAY2
   !-----------------------------------------------------------------------------
          ! Resize the data to LLPAR levels (if necessary)
          ptr_3D => PROD(:,:,:,N)
          CALL TRANSFER_3D( array, ptr_3D )
          NULLIFY( ptr_3D )

       . . .

   !-----------------------------------------------------------------------------
   ! Prior to 4/3/14:
   ! LOSS is now REAL*4 (bmy, 4/13/14)
   !       ! Cast from REAL*4 to REAL*8 and resize to 1:LLPAR
   !       call transfer_3D( array, array2 )
   !
   !       LOSS(:,:,:,N) = ARRAY2
   !-----------------------------------------------------------------------------
          ! Resize the data to LLPAR levels (if necessary)
          ptr_3D => LOSS(:,:,:,N)
          CALL TRANSFER_3D( array, ptr_3D )
          NULLIFY( ptr_3D )
with the corresponding modifications also in routine GET_RATES_INTERP.
Furthermore, I had to modify routine GeosUtil/transfer_mod.F in the following way:
(1) Rename routine TRANSFER_3D to TRANSFER_3D_R8
(2) Copy TRANSFER_3D_R8 to a new routine TRANSFER_3D_R4, anc declare the OUT argument as REAL*4 instead of REAL*8:
   !
   ! !OUTPUT PARAMETERS:
   !
         REAL*4,  INTENT(OUT) :: OUT(IIPAR,JJPAR,LLPAR)   ! Output data
(3) Add the following code at the top of the GeosUtil/transfer_mod.F:
     INTERFACE TRANSFER_3D
        MODULE PROCEDURE TRANSFER_3D_R4
        MODULE PROCEDURE TRANSFER_3D_R8
     END INTERFACE

       . . .

     PRIVATE :: TRANSFER_3D_R4
     PRIVATE :: TRANSFER_3D_R8
The modifications in GeosUtil/transfer_mod.F will allow the 72-level REAL*4 data to be resized and stored in the REAL*4 arrays in GeosCore/strat_chem_mod.F90.
These modifications allow the 2° x 2.5° simulations to proceed in our regular 8-CPU queue. But even with this modification, you may still hit

Also -- the more diagnostics you turn on, the closer you get pushed to the memory limit. So even w/ this fix, you may still encounter a seg fault if you are trying to save out too many quantities to the bpch or timeseries files.

--Bob Y. 15:24, 4 April 2014 (EDT)

Eliminate floating-point invalid in INIT_STRAT_CHEM

This fix was included GEOS-Chem v9-02o (Approved 03 Sep 2013).

In routine INIT_STRAT_CHEM (in module GeosCore/strat_chem_mod.F90), we need to initialize the REAL*8 variable TPAUSEL_CNT as follows:

    Old code:
    TpauseL_Cnt              = 0.

    New code:
    TpauseL_Cnt              = 0.d0

The old code causes a floating point invalid situation when compiling with FPE=yes.

--Bob Y. 12:11, 13 September 2013 (EDT)

Minor bug in parallel DO loop for Bry concentrations

A minor parallelization bug existed in the benchmark v9-01-03o. We have now fixed this in benchmark v9-01-03p.

The BEFORE variable was declared as a 3-D allocatable array at the top of module strat_chem_mod.F90:

 REAL*8, ALLOCATABLE :: Before(:,:,:)      ! Init atm. state each chem. dt

and was allocated at the start of the GEOS-Chem simulation by routine INIT_STRAT_CHEM. BEFORE was also used in this parallel loop in routine DO_STRAT_CHEM:

       !===============================
       ! Prescribe Br_y concentrations
       !===============================

       !$OMP PARALLEL DO &
       !$OMP DEFAULT( SHARED ) &
       !$OMP PRIVATE( NN, I, J, L, BryDay, BryNight, IJWINDOW )
       DO NN=1,6

          IF ( GC_Bry_TrID(NN) > 0 ) THEN

             ! Make note of inital state for determining tendency later
             BEFORE = STT(:,:,:,GC_Bry_TrID(NN))

             ... etc ...

As you can see, the BEFORE array has (I,J,L) scope, but is located within an OpenMP parallel DO loop with (I,J,L,N) scope. Therefore, BEFORE should have been declared in the !$OMP PRIVATE listing, but was not.

To correct this situation we changed BEFORE from an allocatable array to a local array within routine DO_STRAT_CHEM:

   ! Arrays
   REAL*8                    :: STT0(IIPAR,JJPAR,LLPAR,N_TRACERS)
   REAL*8                    :: BEFORE(IIPAR,JJPAR,LLPAR)

and also added it to the !OMP PRIVATE declaration in this parallel DO loop:

       !===============================
       ! Prescribe Br_y concentrations
       !===============================

       !$OMP PARALLEL DO &
       !$OMP DEFAULT( SHARED ) &
       !$OMP PRIVATE( NN, BEFORE, I, J, L, BryDay, BryNight, IJWINDOW )
       DO NN=1,6

          IF ( GC_Bry_TrID(NN) > 0 ) THEN

             ! Make note of inital state for determining tendency later
             ! NOTE: BEFORE has to be made PRIVATE to the DO loop since
             ! it only has IJL scope, but the loop is over IJLN!
             ! (bmy, 8/7/12)
             BEFORE = STT(:,:,:,GC_Bry_TrID(NN))

             ... etc ...

This ensures that the strat-trop exchange diagnostic will be computed properly.

--Bob Y. 17:29, 8 August 2012 (EDT)

Bug fix for Intel Fortran Compiler 12

Prasad Kasibhatla wrote:

I am experiencing the IFORT 12 compile failure of GEOS-Chem v9-01-03 with nested grid options turned on in define.h. The compile seems to be failing in the Stratospheric chemistry module GeosCore/strat_chem_mod.F90.

Please see this wiki post on our Intel Fortran Compiler wiki page for a full description of the issue.

--Bob Y. 15:27, 20 December 2012 (EST)