Atmospheric Chemistry
In this article it is described how to perform different kinds of atmospheric chemistry simulations. This includes the description of simulations with a simplified chemistry and MECCA-based (full) chemistry, their nameless settings, possible modules to make use of and information about initialization data. Further, there are given some examples of typical simulation you can do with ICON-ART including atmospheric chemistry.
Simplified Chemistry
When we talk about simplified calculated chemistry in ICON-ART, we mean that the concentration of the gases we want to simulate is calculated with a parametrization. Here production and depletion rates are calculated to solve the differential equation
numerically. Here, describes the number concentration of a certain tracer, describes the chemical production and is the belonging life time of tracer . For the namelist settings you are able to use for atmospheric chemistry, check out the ART-namelist parameters (see ART namelists). The procedure of creating an ICON-ART simulation in Atmospheric Chemistry always comes back to switching on a namelist parameter and providing the path of the respective XML-file. How to create these for several cases, please check the examples below in the Configurations part.
To learn more about technical details of simplified chemistry, see also Weimer et. al. (2017)
Note: When enabling simplified chemistry with the switch lart_chemtracer = .TRUE., you can improve your runtime but the simulated concentration values are less exact compared to MECCA-based chemistry.
MECCA-based Chemistry
General Information
The MECCA(=Module Efficiently Calculating the Chemistry of the Atmosphere) based chemistry describes a full gas phase chemistry that can be applied as an extension to the parametrized Simplified Chemistry (see above). MECCA based chemistry is generally more exact in the concentration values but the overall runtime is longer compared to purely simplified chemistry simulations. MECCA itself is originally a submodule of the CAABA box model where an air parcel is described as a box and outgoing from this model all exchange processes in- and outward of the box are calculated. As MECCA is part of this model, it contains a wide collection of the most important reactions, including Ozone-, Methane-, HOx-, NOx-, Carbonhydrogen-, Halogene- and Sulfur chemistry. MECCA is available in a supplement, available to download for free and containing all auxiliaries to perform MECCA-simulations.
Including MECCA-based Chemistry in a ICON-ART Simulation
(Note: It is recommended to perform all the following steps in the shell environment.)
The above mentioned collection of the gase phase chemistry reactions can be found in the supplement in the gas.eqn (path: caaba3.0/Mecca/gas.eqn).
Additionally it is also possible to edit existing reactions as well as creating new reactions with the help of "Replacement-files" (see an example in the Configurations part). Inside the gas.eqn every reaction is marked with a certain code. To select the specific reactions for the machanism labels can be set to your belonging reactions or, more easily, a new Gas-Equation-file gas_Mechanism1.eqn can be created, containing only the wanted reactions. (Note: Never edit the original gas.eqn! Better copy it in the first place and then rename and edit it, depending on the respective scientific goal.)
After that the following steps have to be fulfilled to create the code of your specific mechanism and to be able to execute an ICON-ART simulation with MECCA-based chemistry:
- set up a batch file: all previously set information about the mechanism can be selected and stated here (an example can be found below or also inside the supplement in
/caaba3.0/mecca/batch/example.bat). - execute
./meccainside the folder/caaba3.0/mecca. Here the previously created batch file has to be selected and the Fortran files with the mechanism are created. - since the created Fortran code is only located inside Mecca and not in ICON-ART so far, a transfer has to be carried out. A script that performs this transfer can be obtained via
git clone https://gitlab.dkrz.de/art/mecca preproc.git. - in a new directory
Mecca_preprochas been generated and the scriptcreate_icon_code4.shcan be found inside of it. By executing/.create_icon_code4.sh -hpaths to the Mecca- and ICON home directories can be provided as well as a name for the XML-file that is going to be linked in the unscript later. - the Mecca-XML-file is now generated and can be found in ICON in
/icon home>/runctrl examples/xml ctrl.
Now, in the respective runscript the namelist parameter lart_mecca has be set to .TRUE and for cart_mecca_xml the path to the Mecca file can be provided.
Important: As a final step, the ICON code has to be recompiled with the command ./config/dkrz/levante.intel --enable-art --enable-ecrad and after executed make -j 8.
Configurations
Soon some examples for Atmospheric Chemistry simulations will be shown here.
- work in progress -
Example: How to create a simulation with simplified chemistry?
In this first example it is shown how to perform a simulation with simplified chemistry in ICON-ART. Emission data will be applied on the simulation as well. The depicted case is about simulating the tropospheric hydroxyl radical (OH), one of the most important oxidants of the atmosphere. It's main source in the lower troposphere is the photolysis of ozone and its consequent reaction of an excited oxygen atom with the surrounding water vapor:
Additionally the excited Oxygen atom reacts further with Nitrogen and Oxygen:
The main sink of OH in the Troposphere is methane and carbon monooxide:
Now the OH concentrations are calculated with the respective kinetic and photolysis constants, based on chemical kinetic laws: Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle [\ce{OH}]=\frac{\mathrm{2[\ce{O(^1D)}]}k_{\ce{H2O}}[\ce{H2O}]}{k_{\ce{CH4}}[\ce{CH4}]+(k_{\ce{CO,1}}+k_{\ce{CO,2}})[\ce{CO}]}} with Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://wikimedia.org/api/rest_v1/":): {\displaystyle [\ce{O(^1D)}]=\frac{J_{\ce{O3}}[\ce{O3}]}{k_{\ce{O2}}[\ce{O2}]+k_{\ce{N2}}[\ce{N2}]+k_{\ce{H2O}}[\ce{H2O}]}} .
Additionally emission data of the main sinks of OH are implemented. Since the simulation is performed on a R2B04-grid the following emission data are the most suitable ones for the respective trace gases:
- : anthropogenic (EDGAR-432 monthly), biomass-burning (GFED3), biogenic (MEGAN-MACC)
- : anthropogenic (EDGAR-432 monthly)
For more information on recommended emission data see the abstract, dealing with Emission Data. Since emission data are relatively large, they can also just be left out.
Let's start with the runscript that has to be prepared. Here it is named ohsim_simple_icon.run but you can call it differently as well.
Inside of that, check that all your directories are correct, probably they have to be adjusted. Abbreviations used here are the following:
- CENTER: Your organization
- EXPNAME: name of your ICON-Simulation
- OUTDIR: Directory where the simulation output will be stored
- ARTFOLDER: Directory where the ICON-ART code is stored
- INDIR: Directory where the necessary Input data are stored
- EXP:
- lart: For ICON-ART Simulation that has to be switched to
Oheim_simple_icon.run.
#!/bin/bash
CENTER=IMK
workspace=/hkfs/work/workspace/scratch/hp8526-ws_icon_oh
basedir=${workspace}/icon-kit-testsuite
icon_data_poolFolder=/hkfs/work/workspace/scratch/fb4738-dwd_ozone/icon/INPUT/AMIP/amip_input
EXPNAME=ohsim_icon_simple_atom1
OUTDIR=${workspace}/output/${EXPNAME}
ICONFOLDER=/home/hk-project-iconart/hp8526/icon-kit
ARTFOLDER=${ICONFOLDER}/externals/art
INDIR=/hkfs/work/workspace/scratch/fb4738-dwd_ozone/icon/INPUT
EXP=ECHAM_AMIP_LIFETIME
lart=.True.
FILETYPE=4
COMPILER=intel
restart=.False.
read_restart_namelists=.False.
# Remove folder from OUTDIR for postprocessing output
OUTDIR_PREFIX=${workspace}
# Create output directory and go to this directory
if [ ! -d $OUTDIR ]; then
mkdir -p $OUTDIR
fi
cd $OUTDIR
#input for global domain
ln -sf ${INDIR}/../INPUT/GRID/icon_grid_0014_R02B05_G.nc iconR2B05_DOM01.nc
ln -sf ${INDIR}/../INPUT/EXTPAR/icon_extpar_0014_R02B05_G.nc extpar_iconR2B05_DOM01.nc
ln -sf /hkfs/work/workspace/scratch/fb4738-dwd_ozone/icon/output/0014_R02B05/uc1_ifs_t1279_grb2_remap_rev832_0014_R02B05_2018010100.nc ifs2icon_R2B05_DOM01.nc
ln -sf $ICONFOLDER/data/rrtmg_lw.nc rrtmg_lw.nc
ln -sf $ICONFOLDER/data/ECHAM6_CldOptProps.nc ECHAM6_CldOptProps.nc
ln -sf ${ARTFOLDER}/runctrl_examples/init_ctrl/mozart_coord.nc ${OUTDIR}/mozart_coord.nc
ln -sf ${ARTFOLDER}/runctrl_examples/init_ctrl/Linoz2004Br.dat ${OUTDIR}/Linoz2004Br.dat
ln -sf ${ARTFOLDER}/runctrl_examples/init_ctrl/Simnoy2002.dat ${OUTDIR}/Simnoy2002.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_scat-aer.dat ${OUTDIR}/FJX_scat-aer.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_j2j.dat ${OUTDIR}/FJX_j2j.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_scat-cld.dat ${OUTDIR}/FJX_scat-cld.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_scat-ssa.dat ${OUTDIR}/FJX_scat-ssa.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_scat-UMa.dat ${OUTDIR}/FJX_scat-UMa.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_spec_extended.dat ${OUTDIR}/FJX_spec_extended.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_spec.dat ${OUTDIR}/FJX_spec.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/atmos_std.dat ${OUTDIR}/atmos_std.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/atmos_h2och4.dat ${OUTDIR}/atmos_h2och4.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_j2j_extended.dat ${OUTDIR}/FJX_j2j_extended.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_spec_extended.dat ${OUTDIR}/FJX_spec_extended.dat
ln -sf ${ARTFOLDER}/runctrl_examples/photo_ctrl/FJX_spec_extended_lyman.dat ${OUTDIR}/FJX_spec_extended_lyman.dat
# this if condition is necessary because otherwise
# a new link in ${INDIR}/${EXP}/emiss_minimal is generated
# linking to itself
if [ ! -L ${OUTDIR}/emissions ]; then
ln -sd /home/hk-project-iconart/hp8526/emissions ${OUTDIR}/emissions
fi
# this if condition is necessary because otherwise
# a new link in ${INDIR}/${EXP}/emiss_minimal is generated
# linking to itself
if [ ! -L ${OUTDIR}/emissions ]; then
ln -sd /home/hk-project-iconart/hp8526/emissions ${OUTDIR}/emissions
fi
# the namelist filename
atmo_namelist=NAMELIST_${EXPNAME}
#
#-----------------------------------------------------------------------------
Additionally in the next lines of code you set the timing. In this simulation we only simulate a few days. Because OH is dependent from the solar radiation, the output interval is set to 10 hours to calculate OH to every time of the day.
! run_nml: general switches ----------
&# model timing
start_date=${start_date:="2016-07-29T00:00:00Z"}
end_date=${end_date:="2016-08-23T00:00:00Z"}
output_start=${start_date:="2016-07-29T00:00:00Z"}
output_end=${end_date:="2016-08-23T00:00:00Z"}
output_interval="PT10H"
modelTimeStep="PT6M"
leadtime="P10D"
checkpoint_interval="P30D"
Further, all the ICON-namelist parameters have to be set. For a regular ICON-ART-Simulation the following settings are recommended - if not stated
differently.
# model parameters
model_equations=3 # equation system
# 1=hydrost. atm. T
# 1=hydrost. atm. theta dp
# 3=non-hydrost. atm.,
# 0=shallow water model
# -1=hydrost. ocean
#-----------------------------------------------------------------------------
# the grid parameters
declare -a atmo_dyn_grids=("iconR2B04_DOM01.nc" "iconR2B05_DOM02.nc" "iconR2B06_DOM03.nc")
# "iconR2B08_DOM03.nc"
#atmo_dyn_grids="iconR2B07_DOM02.nc","iconR2B08_DOM03.nc"
#atmo_rad_grids="iconR2B06_DOM01.nc"
#-----------------------------------------------------------------------------
# create ICON master namelist
# ------------------------
# For a complete list see Namelist_overview and Namelist_overview.pdf
no_of_models=1
cat > icon_master.namelist << EOF
&master_nml
lRestart = .false.
/
&master_time_control_nml
experimentStartDate = "$start_date"
experimentStopDate = "$end_date"
forecastLeadTime = "$leadtime"
checkpointTimeIntval = "$checkpoint_interval"
/
&master_model_nml
model_type=1
model_name="ATMO"
model_namelist_filename="$atmo_namelist"
model_min_rank=1
model_max_rank=65536
model_inc_rank=1
/
EOF
#-----------------------------------------------------------------------------
#
#-----------------------------------------------------------------------------
#
# write ICON namelist parameters
# ------------------------
# For a complete list see Namelist_overview and Namelist_overview.pdf
#
# ------------------------
# reconstrcuct the grid parameters in namelist form
dynamics_grid_filename=""
for gridfile in ${atmo_dyn_grids}; do
dynamics_grid_filename="${dynamics_grid_filename} '${gridfile}',"
done
radiation_grid_filename=""
for gridfile in ${atmo_rad_grids}; do
radiation_grid_filename="${radiation_grid_filename} '${gridfile}',"
done
# ------------------------
cat > ${atmo_namelist} << EOF
¶llel_nml
nproma = 8 ! optimal setting 8 for CRAY; use 16 or 24 for IBM
p_test_run = .false.
l_test_openmp = .false.
l_log_checks = .false.
num_io_procs = 0 ! up to one PE per output stream is possible
itype_comm = 1
iorder_sendrecv = 3 ! best value for CRAY (slightly faster than option 1)
/
&grid_nml
dynamics_grid_filename = 'iconR2B05_DOM01.nc'
dynamics_parent_grid_id = 0
!radiation_grid_filename = ${radiation_grid_filename}
lredgrid_phys = .false.
lfeedback = .true.
ifeedback_type = 2
/
&initicon_nml
lconsistency_checks = .false.
init_mode = 2 ! operation mode 2: IFS
zpbl1 = 500.
zpbl2 = 1000.
! l_sst_in = .true.
/
&run_nml
num_lev = 90
lvert_nest = .true. ! use vertical nesting if a nest is active
! nsteps = ${nsteps} ! 50 ! 1200 ! 7200 !
! dtime = ${dtime} ! timestep in seconds
modelTimeStep = "${modelTimeStep}"
ldynamics = .TRUE. ! dynamics
ltransport = .true.
iforcing = 3 ! NWP forcing
ltestcase = .false. ! false: run with real data
msg_level = 7 ! print maximum wind speeds every 5 time steps
ltimer = .true. ! set .TRUE. for timer output
timers_level = 10 ! can be increased up to 10 for detailed timer output
output = "nml"
lart = ${lart}
/
&nwp_phy_nml
inwp_gscp = 1
inwp_convection = 1
inwp_radiation = 1
inwp_cldcover = 1
inwp_turb = 1
inwp_satad = 1
inwp_sso = 1
inwp_gwd = 1
inwp_surface = 1
icapdcycl = 3 ! apply CAPE modification to improve diurnalcycle over tropical land (optimizes NWP scores)
latm_above_top = .false., .true. ! the second entry refers to the nested domain (if present)
efdt_min_raylfric = 7200.
itype_z0 = 2
icpl_aero_conv = 1
icpl_aero_gscp = 1
! resolution-dependent settings - please choose the appropriate one
dt_rad = 2160.
dt_conv = 720.
dt_sso = 1440.
dt_gwd = 1440.
/
&nwp_tuning_nml
tune_zceff_min = 0.075 ! ** default value to be used for R3B7; use 0.05 for R2B6 in order to get similar temperature biases in upper troposphere **
itune_albedo = 1 ! somewhat reduced albedo (w.r.t. MODIS data) over Sahara in order to reduce cold bias
/
&turbdiff_nml
tkhmin = 0.75 ! new default since rev. 16527
tkmmin = 0.75 ! ""
pat_len = 100.
c_diff = 0.2
rat_sea = 8.5 ! ** new since r20191: 8.5 for R3B7, 8.0 for R2B6 in order to get similar temperature biases in the tropics **
ltkesso = .true.
frcsmot = 0.2 ! these 2 switches together apply vertical smoothing of the TKE source terms
imode_frcsmot = 2 ! in the tropics (only), which reduces the moist bias in the tropical lower troposphere
! use horizontal shear production terms with 1/SQRT(Ri) scaling to prevent unwanted side effects:
itype_sher = 3
ltkeshs = .true.
a_hshr = 2.0
/
&lnd_nml
ntiles = 3 !!! 1 for assimilation cycle and forecast
nlev_snow = 3 !!! 1 for assimilation cycle and forecast
lmulti_snow = .true. !!! .false. for assimilation cycle and forecast
itype_heatcond = 2
idiag_snowfrac = 2
lsnowtile = .false. !! later on .true. if GRIB encoding issues are solved
lseaice = .true.
llake = .false.
itype_lndtbl = 3 ! minimizes moist/cold bias in lower tropical troposphere
itype_root = 2
/
&radiation_nml
irad_o3 = 10
irad_aero = 6
albedo_type = 2 ! Modis albedo
vmr_co2 = 390.e-06 ! values representative for 2012
vmr_ch4 = 1800.e-09
vmr_n2o = 322.0e-09
vmr_o2 = 0.20946
vmr_cfc11 = 240.e-12
vmr_cfc12 = 532.e-12
/
&nonhydrostatic_nml
iadv_rhotheta = 2
ivctype = 2
itime_scheme = 4
exner_expol = 0.333
vwind_offctr = 0.2
damp_height = 50000.
rayleigh_coeff = 0.10
lhdiff_rcf = .true.
divdamp_order = 24 ! for data assimilation runs, '2' provides extra-strong filtering of gravity waves
divdamp_type = 32 !!! optional: 2 for assimilation cycle if very strong gravity-wave filtering is desired
divdamp_fac = 0.004
l_open_ubc = .false.
igradp_method = 3
l_zdiffu_t = .true.
thslp_zdiffu = 0.02
thhgtd_zdiffu = 125.
htop_moist_proc= 22500.
hbot_qvsubstep = 22500. ! use 19000. with R3B7
/
&sleve_nml
min_lay_thckn = 20.
max_lay_thckn = 400. ! maximum layer thickness below htop_thcknlimit
htop_thcknlimit = 14000. ! this implies that the upcoming COSMO-EU nest will have 60 levels
top_height = 75000.
stretch_fac = 0.9
decay_scale_1 = 4000.
decay_scale_2 = 2500.
decay_exp = 1.2
flat_height = 16000.
/
&dynamics_nml
iequations = 3
idiv_method = 1
divavg_cntrwgt = 0.50
lcoriolis = .TRUE.
/
&transport_nml
! qv, qc, qi, qr, qs
itype_vlimit = 1,1,1,1,1
ivadv_tracer = 3, 3, 3, 3, 3
itype_hlimit = 3, 4, 4, 4 , 4
ihadv_tracer = 52, 2,2,2,2
iadv_tke = 0
/
&diffusion_nml
hdiff_order = 5
itype_vn_diffu = 1
itype_t_diffu = 2
hdiff_efdt_ratio = 24.0
hdiff_smag_fac = 0.025
lhdiff_vn = .TRUE.
lhdiff_temp = .TRUE.
/
&interpol_nml
nudge_zone_width = 8
lsq_high_ord = 3
l_intp_c2l = .true.
l_mono_c2l = .true.
support_baryctr_intp = .false.
/
&extpar_nml
itopo = 1
n_iter_smooth_topo = 1
heightdiff_threshold = 3000.
/
&io_nml
itype_pres_msl = 4 ! IFS method with bug fix for self-consistency between SLP and geopotential
itype_rh = 1 ! RH w.r.t. water
/
&output_nml
filetype = ${FILETYPE} ! output format: 2=GRIB2, 4=NETCDFv2
dom = -1 ! write all domains
output_start = "${output_start}"
output_end = "${output_end}"
output_interval = "${output_interval}"
steps_per_file = 1
include_last = .TRUE.
output_filename = 'icon-art-${EXPNAME}-chem' ! file name base
ml_varlist = 'temp','pres','group:ART_CHEMISTRY','u','v','OH_Nconc'
output_grid = .TRUE.
remap = 1
reg_lon_def = -180.,1,180.
reg_lat_def = -90.,1,90.
/
Now, we're getting to the ICON-ART settings. To enable chemistry in an ICON-ART Simulation inn general, the switch lart_chem has to be set to .TRUE.. With lart_diag_out output of the diagnostic fields can be enabled. Due to setting lart_chem=.TRUE. either lart_chemtracer or lart_mecca have to be set to .TRUE.. Because wee want to perform a simulation with simplified chemistry, we have to switch on lart_chemtracer. If this namelist parameter is set to .TRUE., also cart_chemtracer_xml has to be fulfilled.
&art_nml
lart_chem = .TRUE.
lart_diag_out = .TRUE.
lart_aerosol = .FALSE.
lart_mecca = .FALSE.
lart_chemtracer = .TRUE.
cart_emiss_xml_file = '${ARTFOLDER}/runctrl_examples/emiss_ctrl/emissions_R2B05_0014_cs.xml'
cart_chemtracer_xml = '/home/hk-project-iconart/hp8526/icon-kit/externals/art/runctrl_examples/xml_ctrl/chemtracer_reimus.xml'
cart_input_folder = '${OUTDIR}'
cart_io_suffix = '0014'
iart_init_gas = 0
/
EOF
Example: How to create a simulation with MECCA-based chemistry?
-Work in progress-