Option to specify upper boundary conditions from a file

[dan800]

It would be advantageous to be able to specify an upper boundary condition in CAM, CAM-CHEM and WACCM from a file. CAM-CHEM currently has a zero flux condition and so misses upper atmospheric sources of NOy important for ozone in the stratosphere. Similarly, CAM and CAM-CHEM in "low top" configurations likely miss sources of water vapor that come from methane oxidation in the uppermost stratosphere and mesosphere. In the "high top" CAM version (lid ~80km) specification of mesopause temperature and water vapor, CO2, etc. is needed rather than a hard coded value that does not change with climate change/emissions. In WACCM some values at the UBC are taken from a TIEGCM file or specified from empirical models (MSIS, NOEM) but could easily be set to use an WACCM-X UBC.

I am not sure what the best way to implement this is, but since WACCM has the ability to specify it's UBC from a TIEGCM file then perhaps that code can be generalized.

I see the main source of input UBCs would likely be CMIP6 timeseries files from simulations done by WACCM, mostly likely zonal-mean monthly data. Just as in vertically resolved emission datasets, I think that it should not be a requirement that these inputs be in the native grid, but would be interpolated onto the model grid (as long as they span the boundary).

[fvitt]

This is related to issue #216.
I will see if I can somewhat unify the specification of the upper boundary conditions for the low-top and high-top configurations. Explicit namelist specifications of UB seems like a good place to start.

[islasimpson]

This is great to hear. Do you have any estimate of when this will be ready? We'll be deciding on what to do at the upper boundary for an upcoming seasonal prediction project and some estimate of the time frame would help us to decide whether to proceed as we currently are or wait for this option. I don't mean this to cause any pressure to get it done more quickly - more just an inquiry to make an informed decision as to how to proceed. Thanks!

[fvitt]

@islasimpson
At this time I don't have a good estimate on how long it might take me to code up. I think I can start looking at the details tomorrow. I may have a better idea next week.

[fvitt]

@dan800 @drmikemills Do you have sample CMIP6 timeseries files from simulations done by WACCM that could be used as UBC inputs for development / testing?

[drmikemills]

All runs are archived here:

/glade/campaign/collections/cmip/CMIP6/timeseries-cmip6

Look for runs with compset names, such as BWHIST or FWHIST or BW1850.

Timeseries are under subdirectory atm/proc/tseries. Monthly output is under month_1.

[islasimpson]

Thanks for the update @fvitt

[dan800]

@drmikemills Is it possilbe for you to generate a single file from the individual timeseries that would have zonal means? The most pressing need is for something that spans the HIST period that Simone is running (1970 on?). We could also probably use an annual climatology e.g., for 2005-2015, which would be cycled over for near present day. The latter would be used for S2D simulations for present day.

Fields needed (please everyone speak up if there's ones missing): T, H2O, NO, NO2, CFC11, CFC12, CO2, CH4.

@lkemmons does NO, NO2 cover enough of NOy at the CAM-CHEM low top UBC or are other species needed?

[lkemmons]

For NOy species: HNO3, NO3, N2O5 should also be included, maybe also HO2NO2, N?

[dan800]

I doubt there's much N in NOx or NOy at 80 or 35 km. Probably wouldn't hurt to add the others in CAM-CHEM runs. For CAM, we just need to make sure we have the water budget right and radiatively active gases. What about N2O?

[lkemmons]

N2O drops off rapidly in the stratosphere, as its source is the troposphere.

On Thu, Feb 17, 2022 at 6:57 PM Dan Marsh @.***> wrote:

I doubt there's much N in NOx or NOy at 80 or 35 km. Probably wouldn't hurt to add the others in CAM-CHEM runs. For CAM, we just need to make sure we have the water budget right and radiatively active gases. What about N2O?

— Reply to this email directly, view it on GitHub https://github.com/ESCOMP/CAM/issues/533#issuecomment-1043727840, or unsubscribe https://github.com/notifications/unsubscribe-auth/AH5BH7NA4QOV45GL45MWW3TU3WRRRANCNFSM5OSLCNIQ . Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

You are receiving this because you were mentioned.Message ID: @.***>

[drmikemills]

@dan800 It should be pretty easy for me to create such a zonal mean file with all of the constituents. Should it be for all model levels, or just one that is representative of an upper boundary? Let me know the complete list of constituents that should be included.

[dan800]

Yes, all model levels - then one file can cover L58 and L93. Full list should be a discussion point at today's CAM7 meeting.

[dan800]

@lkemmons re. N2O, agreed - I wondered if any got transported from the tropics and re-entered the high latitudes around the model top before being converted to NOy. Or if it's all long gone by then.

[drmikemills]

Here is a zonal mean UBC file from the ensemble average of the 3 WACCM AMIP runs:

/glade/p/acom/acom-climate/cesm2/postprocess/f.e21.FWHISTBgcCrop.f09_f09_mg17.CMIP6-AMIP-WACCM.003/atm/proc/tseries/month_1_zm/f.e21.FWHISTBgcCrop.f09_f09_mg17.CMIP6-AMIP-WACCM.ensAvg123.cam.h0zm.UBC.195001-201412.nc

[fvitt]

If I understand the vertical diffusion code correctly, the specified mixing ratios at the upper boundary are applied as a boundary conditions to the molecular diffusion solver. In the low top configuration molecular diffusion is not active and the upper boundary mixing ratios are not applied.

[lkemmons]

I am not the best person to decide which species should be included. Doug, Mike, Simone, Rolando, etc., should be included in this discussion.

[fvitt]

Currently in WACCM these species have specified upper boundary mixing ratios:

H2O (2.e-8 vmr)
CH4 (2.e-10 vmr)
F (1.e-15 mmr)
HF (1.e-15 mmr)
H (MSIS)
N (MSIS)
O (MSIS)
O2 (MSIS)
H2 (TGCM)
NO (SNOE)

[dan800]

@fvitt I'm surprised by this. I understand that species dependent molecular diffusion would be only turned on in a high top configuration, but eddy diffusion (Kzz) should operate throughout the vertical domain and I believe calls the same vertical diffusion solver. So I would think that setting the UBC there would be an option. Kzz comes out of all gravity wave modules so perhaps @JulioTBacmeister can comment.

[fvitt]

@dan800 @JulioTBacmeister et al: Take a look at diffusion_solver.F90. The ubc mmrs are used only inside a do_molec_diff conditional:

    866            if( do_molec_diff .and. diffuse(fieldlistm,'q',m)) then
    867 
    868               decomp = vd_lu_qdecomp( pcols , pver   , ncol              , cnst_fixed_ubc(m), cnst_mw(m), &
    869                    kvq   , kq_scal, mw_fac(:,:,m) ,dpidz_sq          , p_molec, &
    870                    interface_boundary, molec_boundary, &
    871                    tint  , ztodt  , nbot_molec       , &
    872                    lchnk , t                , m         , no_molec_decomp)
    873 
    874               ! This to calculate the upper boundary flux of H.    -Hanli Liu
    875               if ((cnst_fixed_ubflx(m))) then
    876 
    877                  ! ubc_flux is a flux of mass density through space, i.e.:
    878                  ! ubc_flux = rho_i * dz/dt = q_i * rho * dz/dt
    879                  ! For flux of mmr through pressure level, multiply by g:
    880                  ! q_i * rho * gravit * dz/dt = q_i * dp/dt
    881 
    882                  call decomp%left_div(q(:ncol,:,m), &
    883                       l_cond=BoundaryFlux( &
    884                       -gravit*ubc_flux(:ncol,m), ztodt, &
    885                       p%del(:,1)))
    886 
    887               else
    888                  call decomp%left_div(q(:ncol,:,m), &
    889                       l_cond=BoundaryData(ubc_mmr(:ncol,m)))
    890               end if
    891 
    892               call decomp%finalize()
    893 
    894            else
    895 
    896               if (.not. use_spcam) then
    897                  ! Currently, no ubc for constituents without molecular
    898                  ! diffusion (they cannot diffuse out the top of the model).
    899                  call no_molec_decomp%left_div(q(:ncol,:,m))
    900               end if
    901 
    902            end if

[dan800]

@fvitt @JulioTBacmeister I tried to read through this again - it's been a few years since I worked through the serpentine code. Frankly it needs a complete rewrite and a ton more commenting. Lots of comments that are in there indicate short term fixes were applied. Not sure what to do at this point. We could turn on molecular diffusion everywhere - I think all it's doing is calculating an additional diffusion coefficient to add to kzz (and Di would be negligible in comparison to kzz in the lower atmosphere). I think all these are collected in the tri-diagonal solver routine. To save computation the diffusion coefficient calc could be skipped and Dis set to zero if pressure at top greater than 0.01mb. That's if you can find where the Di coefficients are actually calculated.

[fvitt]

If I understand the code correctly, when molecular diffusion is active (as it is in waccm) all of the tracers have fixed mixing ratio at the upper boundary. The majority of them are set to zero at the upper boundary. The tracers with non-zero fixed mixing ratios at the upper boundary are listed above for waccm.

[JulioTBacmeister]

I'll look at diffusion_solver.F90 soon, but I believe tracer diffusion is also being done in CLUBB. There is also a vertical diffusion calculation done for EACH separate gravity wave source. I don't know where UBCs for these are coming from. CLUBB is probably not an issue since any diffusion estimated by CLUBB is likely to be zero at the model top even for low-top CAM. However, there will be GW breaking and diffusion near model top. It looks to me like every tracer passed in to the GW schemes is being diffused. Perhaps there is only a limited subset being passed in.

[dan800]

Maybe now is the time to bite the bullet and combine the various Kzzs and Dis and call the solver once, returning one tendency. At that time we can set the boundary condition, which should be a namelist setting for each transported chemical species, temperature and various states of water. It would also be good to separate out the horizontal diffusion operators (@PeterHjortLauritzen may have an opinion on how that can be done).

[dan800]

" It looks to me like every tracer passed in to the GW schemes is being diffused. Perhaps there is only a limited subset being passed in." I think that if it's a tracer it should be diffused - makes little sense to advect it and not diffuse it. We save a lot in chemistry by not transporting all the short lived species.

[islasimpson]

@fvitt - from what I heard at the CAM meeting today, it sounds like this is a bigger task than originally anticipated. From your perspective, if we want to get going on some simulations in the next couple of weeks, do you think we should go ahead without this fix, or do you think it's something that will become available in that time frame? If it's a big task that's going to , I'd be inclined to continue without this upper boundary fix as we need to run a control of order 100 years before this group starts their simulations in April. Let me know what you think. Thanks!

[fvitt]

@islasimpson This turned into a can of worms. I would like some guidance from @dan800 and others on how I should proceed. One option is to read in the new upper boundary conditions and turn on molecular diffusion for the lower top models.

[fvitt]

Would it be reasonable for configurations without a high-top upper boundary (without molecular diffusion) to explicitly overwrite the mixing ratios in the top layer with the upper boundary values, rather than specifying upper boundary conditions to the molecular diffusion solver? The fixed concentration lower boundary conditions are handled this way in the surface layer.

[fvitt]

Do the following namelist options for specifying upper boundary conditions sound reasonable?

ubc_specifier

List of fields that have fixed upper boundary conditions (temperature and mass mixing ratios)
and corresponding data sources. Each entry in this list is of the form:
  FLD->source
where source may be of the following:
  MSIS
  TGCM
  SNOE
  UBC_FILE
  value (mmr|vmr)

These fields can be set to constant values with units of "mmr" (mass mixing ratio)
or "vmr" (volume mixing ratio).  The upper boundary values of source UBC_FILE are prescribed
in the UBC_FILE_path file. Fields of source UBC_FILE may may include a CAM constituent name
and a file field name separated by ":".

Example:
  ubc_specifier = 'Q:H2O->UBC_FILE', 'CH4->2.D-10vmr', 'F->1.0D-15mmr', 'HF->1.0D-15mmr',
                  'T->MSIS', 'H->MSIS', 'N->MSIS','O->MSIS', 'O2->MSIS', 'H2->TGCM', 'NO->SNOE'

ubc_file_path

Full path of the file containing prescribed upper boundary conditions.

ubc_file_input_type

Type of time interpolation for data in ubc_file_path file.
Can be set to 'CYCLICAL', 'SERIAL', 'INTERP_MISSING_MONTHS', or 'FIXED'.

ubc_file_cycle_yr

The cycle year of the prescribed upper boundary data
if ubc_file_input_type is 'CYCLICAL'.
Format: YYYY

ubc_file_fixed_ymd

The date at which the upper boundary data is fixed
if ubc_file_input_type is 'FIXED'.
Format: YYYYMMDD

[dan800]

@fvitt

Thanks for thinking this through and making a plan. I like the suggestions so far. I think we can do away with the TGCM option. We have WACCM-X now which should be at least as good. What I don't understand is why this doesn't or couldn't fall under file based UBC specification. We can include H2 in the same UBC file that includes H2O (both using a WACCM-X climatology).

For WACCM-X would this section be blank - i.e., the default is zero flux condition for everything?

Dan

[fvitt]

@dan800 We could certainly remove the TGCM option or leave it in for backwards compatibility. Any and all these fields could be included in the new ubc data file.

For configurations that have zero flux upper boundary conditions for all constituents (including waccmx) would have a blank ubc_* namelist settings (or not set).

[dan800]

Any and all these fields could be included in the new ubc data file.

I think that is preferable to fixing the value or TGCM for most of what we do (assuming we can take it from a WACCM-X run). Even taking it from a FX2000 case is going to be more self consistent. If you want to drastically scale the UBC you still have a fixed single number option or you could apply a scale factor within the UBC file.

[fvitt]

@dan800 @islasimpson I have a draft implementation ready to test. The github repo for my branch is: https://github.com/fvitt/cam/tree/ubc_spcr

This reproduces bit-for-bit the current upper boundary functionality of WACCM.

I have done some very limited testing of specifying upper boundary concentrations using Mike's zonal mean UBC file in a low-top configuration. It seems to be working as expected.

Should I move ahead with pushing the changes onto the cam_development branch? Or do we need to do more testing of the new functionality first?

[dan800]

Hi Francis, That's great. I think we would want a few tests first. Have you tested in L32 version, where it would replace a zero flux condition. I think running a few years with the only change in CAM L32 is specified UBC only on H2O (q(1)?) and the standard zero flux config would be useful and not too expensive. @cecilehannay @JulioTBacmeister Thoughts on that?

• Dan

[islasimpson]

That's great. Given their timeframe, Scripps is going without this addition, but it might be good to test with L83 in any case. Thanks.

[fvitt]

I ran a 2-year test of F2000climo (32 levels) with H2O UBC from Mike's zonal mean file.

History files of the test run can be found under: /glade/scratch/fvitt/archive/f2000_ubc_f09_32lev_test01

This can be compared with a control run: /glade/scratch/fvitt/archive/f2000_ubc_f09_32lev_cntrl01

Here are some plots of H2O in the final month of the 2-year run:

Test run H2O:

Control run H2O:

UBC file H2O:

[dan800]

@fvitt Seems to me to be a huge improvement. Worth running a few more years to see if the high latitude water vapor fills in even more and allows looking at an annual mean.

On Mar 29, 2022, at 19:55, Francis Vitt @.***> wrote:

 I ran a 2-year test of F2000climo (32 levels) with H2O UBC from Mike's zonal mean file.

History files of the test run can be found under: /glade/scratch/fvitt/archive/f2000_ubc_f09_32lev_test01

This can be compared with a control run: /glade/scratch/fvitt/archive/f2000_ubc_f09_32lev_cntrl01

Here are some plots of H2O in the final month of the 2-year run:

Test run H2O:

Control run H2O:

UBC file H2O:

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were mentioned.

[fvitt]

Last month of 6 years run:

[dan800]

Even better! Seems that it might have improved the tape recorder as well - possible mixing in from mid-latitudes or a simple radiative impact on the cold point? @fvitt can you please point me to the h0s for test/control runs and also the UBC file used. @cecilehannay @JulioTBacmeister if we go forward with this change we probably need to do a set of CAM diagnositics. FWIW, I don't expect RESTOM to be the same given these differences, but better to have the right strat. GHG forcing when trying to balance the model.

[fvitt]

@dan800

History files for the control run can be found in: /glade/scratch/fvitt/archive/f2000_ubc_f09_32lev_cntrl01/atm/hist

For the test run: /glade/scratch/fvitt/archive/f2000_ubc_f09_32lev_test01/atm/hist

The UBC namelist settings:

ubc_file_path='/glade/p/cesmdata/cseg/inputdata/atm/cam/chem/ubc/f.e21.FWHISTBgcCrop.f09_f09_mg17.CMIP6-AMIP-WACCM.ensAvg123.cam.h0zm.UBC.195001-201412_c220322.nc'
ubc_file_input_type='CYCLICAL'
ubc_file_cycle_yr  = 2000
ubc_specifier = 'Q:H2O->ubc_file'

[dan800]

Annual mean H2O averaged over years 2-6 from zero flux UBC (left), WACCM UBC (middle) and WACCM-zero flux (right) zoom in of right figure: I checked RESTOM and got -2.0187 W/m2 for the zero flux condition. The WACCM UBC run further decreased RESTOM by 0.031 W/m2, which I think makes sense - strat water would decrease LW at TOM.

[dan800]

For reference, the WACCM H2O average 1995-2004 (the source of the UBC for the L32 WACCM UBC run).

[JulioTBacmeister]

I'd like some help interpreting the figures above. Should I regard the "WACCM H2O average 1995-2004" immediately above as the target?

Is the WACCM run with zero flux BCs close to this target? It seems to be, but the color scale difference makes it a little difficult to judge.

So then the two runs labeled "zero flux UBC" and "WACCM UBC" are with the L92 model, L58 model? Either way the take-home seems to be that specifying UBCs from WACCM is working well to bring the low/mid top model close to the full WACCM simulation. Is this the correct interpretation?

[dan800]

@JulioTBacmeister I hope this is easier to interpret. Presumably the WACCM H2O is pretty good. Control and Test are L32. Kind of surprised the L32 doesn't look too bad in the tropics at the lowest level - hmmmm. What am I missing? It could be the annual mean is washing out the phase in the tape recorder where we see real issues. For these plots I picked levels where the WACCM levels match the L32 levels (73 & 35 hPa) or was pretty close 7hPa.

[JulioTBacmeister]

So as you approach the tropopause the answer gets worse in L32, presumably due to a bad stratospheric circulation. I thought our concern would be how this looks in an L92 run.

[dan800]

@fvitt - is there some representation of the CH4 oxidation in the L32 model? Either a fixed CH4 oxidation H2O source or CH4 is a transported species that has a simple stratospheric loss rate?

@JulioTBacmeister my concern in that in the low top versions the strat water vapor was low due to a low entry value and missing CH4 oxidation source. I think at least half of the water near the stratopause is from CH4 (two H2Os per CH4 molecule). I would like to get the radiative forcing from stratospheric GHGs in the L32/58/92 models as good as we can get it.

I think the next step is to try this approach in the latest L58 tag - my guess is that it wont look as good as in the L32 due to the large dry bias and we may see a bigger RESTOM difference.

In the L92 I hope that if we can put in CH4 or a source of H2O from CH4 we can do a pretty good job and getting the UBC right is going to just make it a little better.

Either way, we want to minimize the stratospheric radiative forcing differences between the L92 and L58. My guess is the urge will be to run coupled with the L58 due to cost (no one seems itching to switch to the L92 configuration yet) and we want to minimize climate drift when we switch between the two. So getting H2O close is a priority. The UBC looks to be getting us part of the way there - certainly better than the current L32 config - and is low hanging fruit with almost zero cost.

[fvitt]

@dan800 The low top model includes this forcing of water vapor: 'H2O -> /glade/p/cesmdata/cseg/inputdata/atm/cam/chem/emis/elev/H2O_emission_CH4_oxidationx2_elev_3Dmonthly_L70_2000climo_c180511.nc'

[fvitt]

I would be happy to run tests for other configurations -- L58, L92, or whatever. Do these higher-top configurations require using someone's branch code, or is it as simple as changing the number of layers and using an appropriate IC file?

[dan800]

@fvitt Thanks Francis. I think Cecile has a particular candidate tag that is being evaluated, so best to run that. Both L58 and L92 would be good but start with L58. The L58 has the same model top as the L32, just higher resolution (AFAIK). Interesting that in spite of the CH4 oxidation H2O source, the stratosphere is still dry, indicating a large contribution of water from above ~3 hPa.

On Mar 31, 2022, at 19:56, Francis Vitt @.***> wrote:

 I would be happy to run tests for other configurations -- L58, L92, or whatever. Do these higher-top configurtions require using someone's branch code, or is it as simple as changing the number of layers and using an appropriate IC file?

— Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were mentioned.

[fvitt]

@dan800 This is the L58 with water vapor upper boundary test case:

/glade/work/fvitt/cesm/cases/f.e21.FWscHIST.ne30_L48_BL10_cam6_3_041_control.hf.001_ubc01

The source code root: /glade/work/fvitt/camdev/cam6_3_041_ubc/

UBC namelist settings:

ubc_specifier = 'Q:H2O->ubc_file'
ubc_file_path='/glade/p/cesmdata/cseg/inputdata/atm/cam/chem/ubc/f.e21.FWHISTBgcCrop.f09_f09_mg17.CMIP6-AMIP-WACCM.ensAvg123.cam.h0zm.UBC.195001-201412_c220322.nc'
ubc_file_input_type='SERIAL'

10 years run. History files are in /glade/scratch/fvitt/archive/f.e21.FWscHIST.ne30_L48_BL10_cam6_3_041_control.hf.001_ubc01/atm/hist

58L model H2O:

UBC file H2O:

[dan800]

@cecilehannay Can you please point me to the h0s for the baseline case to compare with this new L58 run? thanks, Dan