Scaling Biomass Burning (BB) emission flux regionally

Hello,

I have been attempting to debug a problem for some time and now realise I need specialist expertise.

My suite u-ce487 is set up to perturb multiuple model parameter values simultaneously. Each model variant is allocated an index, so there are multiple atmos_main tasks. In this case, atmos_main001 and atmos_main002 have scale factors for BB emissions (0.5 and 2.0 respectively). But, output from these tasks is identical.

Note, that this style of perturbation to BB emission flux (and other parameters) worked without error on ARCHER, but has not yet worked on ARCHER2.

Debug commands are included in the ukca_add_emiss_mod.F90 routine
of my branch
branches/dev/leightonregayre/vn11.1_AEROCOM_MMPPEs@102825

A few pertinent debug WRITE statements can be found by searching for ‘ACURE_ERROR’.

The first of these writes out lat/lon indices and the value of an mask we use to apply scale factors regionally (acure_bb_emiss_region). In this simple example, the scale factors are the same for each region (0.5, or 2.0 for atmos_main001 and atmos_main002). I wrote these out to validate the scale factors are being applied to all lat/lon pairs, so kept all pe_output following advice in ticket #195 Retaining all pe output by avoiding housekeeping Scale factors are applied as intended, which confirmes our regional scaling ancilliary is read and used correctly as a mask.

The second debug statement writes out the values of a) the scale factor to be applied to BB emissions after they’re read in (aero_scaling_2d_carb_bb) and the value of the scale factor from our perturbation wrapper for one of the regions (acure_carb_bb_ems_rsh). These values are exactly as intended, 0.5 and 2.0 for all lat/lon pairs in atmos_main001 and atmos_main002 respectively.

The values of ‘aero_scaling_2d_carb_bb’ are then applied in the same routine as a scale factor to BC from BB. This worked as intended on ARCHER, but on ARCHER2 it is as if all scale factors are 1 (ignored - no change). There is a logical (l_acure_carb_bb_ems) to indicate if this emission source should be scaled, which logs show is true.

Note, I have used BB as an example here, but apply the same scale factors to fossil fuel (ff) and residential (res, or biofuel) emissions, with no effect on ARCHER2.

Summary: regional scale factors and regional ancillaries are combined to define a scale factor for BB emissions in each lat/lon correctly. The code for these scale factors is unchanged from ARCHER to ARCHER2, but on ARCHER2 there is no effect.

Any thought or ideas on how I can debug this problem further are welcome.

Thanks,

Leighton

Hi Leighton

We’re not ignoring this - waiting for inspiration more like. Is this still a problem?

Grenville

Hi Grenville,

Thanks for getting back to me. Yes, it’s still a problem. I’ve puzzled over this for some time before submitting this ticket. I added some further debug write statements to my branch, which may provide some inspiration.

Additional write statements are included in ukca_add_emiss.F90 for the maximum and minimum values of ‘em_field_2d’, within code segements for each tracer, in both the 2D and 3D code blocks (search for ‘DEBUG 2D’ or ‘DEBUG 3D’).

From the pe_output,
/work/n02/n02/lre/cylc-run/u-ce487/work/20170101T0000Z/atmos_main001/pe_output
it looks like emissions are added for SO2 (2D) and SO2_nat (3D), but none of the other tracers. For example, I can’t find any output with the text ‘BC_FF’, which should be written out whenever the tracer ‘BC_fossil’ is read in.

Does that help narrow down the possible causes of this problem?

Thanks,

Leighton

Hi Grenville,

Do you have any further thoughts on this problem? It’s preventing me from working on a problem with a once-remote deadline that is now not so remote.

Thanks,

Leighton

Hi Leighton

I’ve not had a chance to test this, but the following is very dodgy code. Testing against NINT(of a real number) is a dangerous thing to do – I suggest you check that the IF does what you think it should, and then stop testing the value of (acure_bb_emiss_region(i,j)) this way - check its value by bracketing it with .lt. and .gt. statements.

Grenville

 IF (NINT(acure_bb_emiss_region(i,j)) == 1) THEN
            aero_scaling_2d_carb_bb(i,j) = acure_carb_bb_ems_sam
          ELSE IF (NINT(acure_bb_emiss_region(i,j)) == 2) THEN
            aero_scaling_2d_carb_bb(i,j) = acure_carb_bb_ems_naf
          ELSE IF (NINT(acure_bb_emiss_region(i,j)) == 3) THEN
            aero_scaling_2d_carb_bb(i,j) = acure_carb_bb_ems_saf
          ELSE IF (NINT(acure_bb_emiss_region(i,j)) == 4) THEN
            aero_scaling_2d_carb_bb(i,j) = acure_carb_bb_ems_bnh
          ELSE IF (NINT(acure_bb_emiss_region(i,j)) == 5) THEN
            aero_scaling_2d_carb_bb(i,j) = acure_carb_bb_ems_rnh
          ELSE IF (NINT(acure_bb_emiss_region(i,j)) == 6) THEN
            aero_scaling_2d_carb_bb(i,j) = acure_carb_bb_ems_rsh
          ELSE
            aero_scaling_2d_carb_bb(i,j) = 1.0

Hi Grenville,

I appreciate the code is less than perfect. I will edit the subroutine logical tests as you suggest.

But, the problem is not caused by this poor coding. The values of the scale factors are assigned correctly. For example, I write out the value of ‘aero_scaling_2d_carb_bb’ and ‘acure_carb_bb_ems_rsh’ (search for “bb_ems_rsh” in the logs). Both have the same value 0.5, yet the emissions are effectively unscaled (value 1).

Any other ideas would be welcome.

Hi Leighton

Yes, my mistake - I should have actually tested it.

The only way text like BC_FF is not written out, is if

ELSE IF (l_acure_carb_ff_ems .AND.                        & ! LR apply 2D FF scaling to BC only
             (TRIM(emissions(l)%tracer_name) == 'BC_fossil')) THEN

is never satisfied. I added

      write(6,*) '^^^^^^^^^^^^^^^^^^^^^^^^^^^^^'
      write(6,*) 'l_acure_carb_bb_ems ', l_acure_carb_bb_ems
      write(6,*) 'TRIM(emissions(l)%tracer_name) ', l, TRIM(emissions(l)%tracer_name)
      write(6,*) '^^^^^^^^^^^^^^^^^^^^^^^^^^^^^'

at the beginning of the DO ilev = 1, model_levels loop.

the only values written out for TRIM(emissions(l)%tracer_name) are SO2_nat, mode_emiss, DMS, Monoterp, SO2_low, SO2_high. So, that explains why you don’t see BC_FF – I don’t know where emissions(l)%tracer_name is set. It seems unlikely that this is an ARCHER2 problem. Are you sure the UM configuration didn’t change between ARCHER and ARCHER2?

Grenville

Hi Grenville,

Thanks for working that out. So, all BC and OC sources are omitted.

The model version is identical. Version 11.1. I copied one of my ARCHER suites and adapted for ARCHER2. I also added some code changes for unrelated parameters, so branched from 11.1 then merged with an ealier branch of 11.1 as used on ARCHER. There are no important structural changes to the model version between platforms, that I’m aware of. Certainly, the carbonaceous aerosol should be read in from file.

I’ve needed to update the paths for my ‘ukca_em_files’. But, I’ve checked these on ARCHER2 and they seem sensible. For example, /work/y07/shared/umshared/CMIP6_ANCIL/data/ancils/n96e/timeslice_2014/AerosolChemistryEmissions/v3/BC_fossil_anthropogenic_2014_time_slice.nc exists and has a field labelled ‘BC_fossil’ which is expected in ukca_add_emiss.F90’.

Tracers are defined in ukca_setup_indices.F90, where my suite uses ukca_indices_sussbcocdu_5mode (as did my ARCHER suite; defined by i_mode_setup=2). This subroutine forces BC to be emitted into aitken insoluble mode only and the code has not changed.

I’m still using the i_ukca_chem=“Offline Oxidants(BE) (14)”, [with l_ukca_offline_be=TRUE], which has multiple BC emission sources defined in the list ‘em_chem_spec’, including BC_fossil. ‘em_chem_spec’ is defined within ‘ukca_setd1defs.F90’, which is used in ‘ukca_emiss_ctl_mod.F90’, which in turn calls ‘ukca_add_emiss.F90’. None of these have changed substantially (mostly debug write statements) between my ARCHER and ARCHER2 suites.

To add to my confusion, there are clearly BC emissions in my suite. The diagnostic ‘Primary Black Carbon to Aitken (Ins)’ is populated.

I guess I need to seek further help from someone with UKCA model domain expertise related to the handling of emissions.

Thanks for your help getting this far.

Leighton

Leighton

tracer_name is correctly set. The reason BC_biofuel etc are not processed in ukca_add_emiss_mod.F90 results from these lines of code:

 IF (k < 0) THEN
    IF (ANY(aero_ems_species == emissions(l)%tracer_name)) THEN
      CYCLE

this is true for BC_biofuel, BC_biomass, BC_fossil etc. CYCLE skips the emissions type. Maybe you could consult a UKCA expert to understand what is going on at the beginning of the DO l = 1, num_em_flds ! loop over emission fields loop in ukca_add_emiss_mod.F90.

Grenville

Hi Grenville,

I checked with CEMAC staff and the code section you pointed out in comment #8 is for ‘non-interactive’ or ‘offline’ emission inventories. My job.out log file shows SO2 files are treated in the same way, so I think the CYCLE code is a red herring.

No closer to understanding why the BC emissions are not read in by ukca_add_emiss.F90, yet the diagnostic for primary BC emiss is populated.

The BC ancillaries I’m using here are cyclical ‘timeslice’ .nc files and many come from the same CMIP6 shared area where SO2 emissions are sourced (these are read in fine), so shouldn’t be treated differently in my suite. I’ve checked the global attributes of these ancillaries and they’re identical between SO2 and BC sources in terms of ‘update_type’ and ‘update_freq_in_hours’.

I updated my code as you suggested in comment 5 above. As expected that did not solve the problem with missing ancillaries, just improved my coding :wink:

I also have write statements in ‘emiss_io_mod.F90’ that show the BC ancillaries are being opened and the correct emission types (e.g. ‘BC_fossil’) are being identified. There are no instances of ‘EM_GET’ in my log files, which I think would indicate an error reading the ancillary information.

I’m at a loss as to how to proceed with debugging this error. Are you aware of anyone else using these ‘timeslice’ ancillaries successfully on ARCHER2?

Thanks,

Leighton

IF (k < 0) THEN
    IF (ANY(aero_ems_species == emissions(l)%tracer_name)) THEN
     write(6,*) '################'
     write(6,*)  'l = ', l, k
     write(6,*) 'tracer_name ', emissions(l)%tracer_name
     write(6,*) '################'
      CYCLE

gives

################
l = 7, -1
tracer_name BC_biofuel
################

################
l = 10, -1
tracer_name BC_fossil
################

################
l = 13, -1
tracer_name BC_biomass
################
################
l = 16, -1
tracer_name OC_biomass
################

################
l = 19, -1
tracer_name OC_biofuel
################

################
l = 22, -1
tracer_name OC_fossil
################

these are skipped:

Thanks Grenville,

Perhaps there has been some misunderstanding locally about this code.

Could you tell me what your write statement outputs for SO2 emissions please?

Leighton

Leighton

(I can’t claim any domain experience - I’m just adding write statements to follow the code)

Only tracers for which k < 0 are skipped. For SO2, DMS,

tracer_name SO2_nat, k =  3

tracer_name DMS, k =  2

tracer_name SO2_low, k =  3

tracer_name SO2_high, k =  3

presumably

    DO m = 1, jpctr          ! loop over tracers
      IF (mapped_tracer == advt (m)) THEN
        k = m
        EXIT
      END IF

doesn’t find BC_biomass in advt. Perhaps you or CEMAC know how advt is created?

Grenville

Hi Grenville,

I apprecite you can’t have domain experience in all aspects of the model. This really does need input from our CEMAC colleagues who are aware. Your suggestions here have helped progress our understanding, so really valuable.

Thanks,

Leighton

Hi Grenville,

Would it be possible to have Luke Abraham comment on asad funtionality please?

When we set up out ensemble, we were under the assumption that ukca_add_emiss_mod.F90 would carry BC as a tracer. But, in hindsight, it seems this subroutine only carries CHEMICAL tracers.

My log files show:

ASAD IS TREATING ADVECTED TRACERS IN THE ORDER:
1 H2O2
2 DMS
3 SO2
4 H2SO4
5 DMSO
6 Monoterp
7 Sec_Org

so BC and OC are not being advected as tracers. So, we may have attempted to perturb regional BC emissions in the wrong subroutine. The confusion remains in that BC emission types are recognised as ‘non-interactive emiss’ types in the ukca_add_emiss_mod.F90 subroutine. We tested regional perturbations to SO2 emissions on ARCHER, but of course those are advected in this routine, so this may not be an ARCHER/ARCHER2 problem as initially assumed.

The use of ‘tracers’ rather than ‘CHEMICAL tracers’ in ukca_add_emiss_mod.F90 is confusing, as aerosol are also reffered to as tracers, for example, ukca_setup_indices.F90

Thanks,

Leighton

Hi Leighton,

Sorry for my delay in replying - I was completely snowed-under last week.

ASAD is purely the chemical solver, and actually has nothing to do with emissions or the emissions routines. The tracers it considers are only those mentioned in the chemical mechanism as defined in ukca_chem_master.F90 or one of the other ukca_chem scheme definitions. You’re running offline oxidants so the scheme is rather basic.

The emissions routines consider both chemical and aerosol emissions and are not as well integrated as they should be. The logic is considered in ukca_emiss_mode_mod.F90, and ukca_emiss_ctl_mod.F90 is where ukca_add_emiss is called. In ukca_emiss_mod.F90 there is the routine ukca_emiss_update_mode which “Convert(s) total aerosol emission mass into mass and number for each mode.”

You might want to follow the logic from the call to ukca_emiss_ctl down from ukca_main to see how the tracers are treated. I’m sorry this is rather opaque - it’s a complicated set of routines as it deals with NetCDF formats, the associated metadata, and the fact that historically aerosol and chemical emissions were done very differently and this basically maintains that distinction and those differences. Essentially the aerosol tracers are considered there, but there are separate routines for these.

Best wishes,
Luke

Hi Luke,

Thanks for the detailed answer. You’ve confirmed our misunderstanding of how ukca_add_emiss_mod.F90 works. Hopefully, whenever these routines are updated, the naming and/or header comments can be adapted to make the purpose of this code a little clearer. I’ll follow the logic from ukca_emiss_ctl as you suggest.

Thanks again Grenville for reflecting on this tricky and deceptive problem.

All the best,

Leighton

This topic was automatically closed 2 days after the last reply. New replies are no longer allowed.