Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Implement diurnal ozone downscaling #2496

Open
wants to merge 45 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from 43 commits
Commits
Show all changes
45 commits
Select commit Hold shift + click to select a range
f00dcc1
add new diurnal ozone mod
Sep 1, 2022
fb2ba21
add do3 streams namelist info
Sep 2, 2022
37fed71
move to all lower case
Sep 2, 2022
898cead
Merge branch 'master' into diurnal_ozone
Sep 8, 2022
9079404
fix namelist variables
adrifoster Sep 8, 2022
743bfa0
update fieldlist names
Sep 8, 2022
40110ba
update bld
Sep 8, 2022
01c3efa
fix do3 streams
Sep 8, 2022
ea81976
create diurnalozonetype
Sep 13, 2022
c04e1fd
remove notes
Sep 13, 2022
1abfed3
add call to shr_ozone_readnl
Sep 13, 2022
eb989a4
move frequency flag to OzoneModType
Sep 13, 2022
53b6008
typo
Sep 13, 2022
46cd241
fix iulog error
Sep 13, 2022
57ab73a
fix private issue
Sep 13, 2022
56792c2
fix use statement
Sep 13, 2022
67d01db
add init call
Sep 13, 2022
1e1a540
remove ntimes
Sep 13, 2022
ce71f95
checking includes
Sep 13, 2022
ab88f6a
typo fix?
Sep 13, 2022
734ea15
update ozone type name
Sep 13, 2022
5617195
fix circular dependency
Sep 16, 2022
09d3c2f
adding introplation
Sep 27, 2022
f651f99
Merge branch 'master' into diurnal_ozone
Sep 27, 2022
f214f12
lots of updates, design doc
Sep 28, 2022
06a6271
add design doc
Sep 28, 2022
0c3ad39
previous updates
Jan 9, 2023
51f644f
merge in latest master
Apr 24, 2024
f85d84f
remove weird files?
Apr 24, 2024
673b90a
remove more weird files
Apr 24, 2024
f2d58ca
updating Interp function
Apr 24, 2024
0ecd18b
bug fix
Apr 24, 2024
f62bf19
various updates to actual read namelist args
adrifoster Apr 25, 2024
b017013
final updates
adrifoster Apr 26, 2024
526385b
remove ws diffs from buildnamelist
adrifoster Apr 26, 2024
541d806
remove ws diffs
adrifoster Apr 26, 2024
0586a7a
more ws
adrifoster Apr 26, 2024
dfe6489
even more ws
adrifoster Apr 26, 2024
314e66d
Merge branch 'master' into diurnal_ozone
adrifoster Apr 26, 2024
87c2183
Merge branch 'master' into diurnal_ozone
adrifoster May 6, 2024
f1778be
change logic for use_do3_streams
adrifoster May 6, 2024
8875fba
fix assert
adrifoster May 7, 2024
752870a
Merge branch 'master' into diurnal_ozone
adrifoster May 7, 2024
ebcdeda
update indentation
May 8, 2024
bfbe350
fixing indentation inconsistencies
May 17, 2024
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
33 changes: 32 additions & 1 deletion bld/CLMBuildNamelist.pm
Original file line number Diff line number Diff line change
Expand Up @@ -1581,6 +1581,7 @@ sub process_namelist_inline_logic {
setup_logic_luna($opts, $nl_flags, $definition, $defaults, $nl, $physv);
setup_logic_hillslope($opts, $nl_flags, $definition, $defaults, $nl);
setup_logic_o3_veg_stress_method($opts, $nl_flags, $definition, $defaults, $nl,$physv);
setup_logic_do3_streams($opts, $nl_flags, $definition, $defaults, $nl, $physv);
setup_logic_hydrstress($opts, $nl_flags, $definition, $defaults, $nl);
setup_logic_dynamic_roots($opts, $nl_flags, $definition, $defaults, $nl, $physv);
setup_logic_params_file($opts, $nl_flags, $definition, $defaults, $nl);
Expand Down Expand Up @@ -1718,6 +1719,11 @@ sub process_namelist_inline_logic {
##################################
setup_logic_cropcal_streams($opts, $nl_flags, $definition, $defaults, $nl);

##################################
# namelist group: dO3_streams #
##################################
setup_logic_do3_streams($opts, $nl_flags, $definition, $defaults, $nl);

##########################################
# namelist group: soil_moisture_streams #
##########################################
Expand Down Expand Up @@ -4223,6 +4229,31 @@ sub setup_logic_cropcal_streams {

#-------------------------------------------------------------------------------

sub setup_logic_do3_streams {
# input file for creating diurnal ozone from >daily data
my ($opts, $nl_flags, $definition, $defaults, $nl, $physv) = @_;

add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'use_do3_streams');
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • As you note in the design document, I'm not sure there's a use case where a user would want to change use_do3_streams. I'd suggest removing that and just basing the logic on atm_ozone_frequency. (You could potentially still avoid adding the stream definition stuff if ozone is off... not sure if it's worth doing that or not.)

if ( &value_is_true( $nl->get_value('use_do3_streams') ) ) {
if ($opts->{'driver'} ne "nuopc") {
$log->fatal_error("Cannot use do3_streams=.true. with MCT.");
}
Comment on lines +4238 to +4240
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • This check can probably be removed now that mct is no longer supported.

add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'stream_fldfilename_do3',
'hgrid'=>"360x720cru" );
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'stream_meshfile_do3',
'hgrid'=>"360x720cru" );
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'do3_mapalgo',
'hgrid'=>$nl_flags->{'res'} );
Comment on lines +4241 to +4246
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Is hgrid actually needed as an attribute to these add_default calls? It is confusing to me because the hgrid you're specifying doesn't match the grid of the files (based on the files given in namelist_defaults). And since there is only a single option in the namelist_defaults, I'm thinking that you probably don't need to specify hgrid (which I would think would be to choose between multiple options) - but I could be misunderstanding.

} else {
if ( defined($nl->get_value('stream_fldfilename_do3'))) {
$log->fatal_error("One of the do3 streams namelist items (stream_fldfilename_do3, " .
" is defined, but use_do3_streams option set to false");
}
}
}

#-------------------------------------------------------------------------------

sub setup_logic_soilwater_movement {
my ($opts, $nl_flags, $definition, $defaults, $nl) = @_;

Expand Down Expand Up @@ -4703,7 +4734,7 @@ sub write_output_files {
soil_resis_inparm bgc_shared canopyfluxes_inparm aerosol
clmu_inparm clm_soilstate_inparm clm_nitrogen clm_snowhydrology_inparm hillslope_hydrology_inparm hillslope_properties_inparm
cnprecision_inparm clm_glacier_behavior crop_inparm irrigation_inparm
surfacealbedo_inparm water_tracers_inparm tillage_inparm);
surfacealbedo_inparm water_tracers_inparm tillage_inparm do3_streams);

#@groups = qw(clm_inparm clm_canopyhydrology_inparm clm_soilhydrology_inparm
# finidat_consistency_checks dynpft_consistency_checks);
Expand Down
6 changes: 6 additions & 0 deletions bld/namelist_files/namelist_defaults_ctsm.xml
Original file line number Diff line number Diff line change
Expand Up @@ -1732,6 +1732,12 @@ lnd/clm2/surfdata_esmf/NEON/surfdata_1x1_NEON_TOOL_hist_78pfts_CMIP6_simyr2000_c
<lai_mapalgo hgrid="1x1_asphaltjungleNJ" >nn</lai_mapalgo>
<lai_mapalgo hgrid="5x5_amazon" >nn</lai_mapalgo>

<!-- do3 streams namelist defaults -->
<use_do3_streams>.true.</use_do3_streams>
<stream_fldfilename_do3>/glade/work/afoster/ozone_update/ozone_damage_files/converted/diurnal_factor_O3surface_ssp370_2015-2024_c20220502.nc</stream_fldfilename_do3>
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • This file should be moved into inputdata.

<stream_meshfile_do3>share/meshes/fv0.9x1.25_141008_polemod_ESMFmesh.nc</stream_meshfile_do3>
<do3_mapalgo>bilinear</do3_mapalgo>

<!-- crop calendar streams namelist defaults -->
<stream_year_first_cropcal >1850</stream_year_first_cropcal>
<stream_year_last_cropcal >2100</stream_year_last_cropcal>
Expand Down
32 changes: 32 additions & 0 deletions bld/namelist_files/namelist_definition_ctsm.xml
Original file line number Diff line number Diff line change
Expand Up @@ -1768,6 +1768,38 @@ Mapping method from LAI input file to the model resolution
copy = copy using the same indices
</entry>

<!-- dO3_streams streams Namelist -->
<!-- ======================================================================================== -->

<!-- diurnal ozone anomaly -->
<entry id="use_do3_streams" type="logical" category="physics"
group="clm_inparm" valid_values="" value=".false.">
Toggle to turn on reading in of diurnal ozone anomaly file. Used to convert >daily ozone to sub-daily
</entry>

<entry id="stream_fldfilename_do3" type="char*256(30)" category="datasets"
input_pathname="abs" group="do3_streams" valid_values="" >
Filename of input stream data for diurnal ozone anomaly
</entry>

<entry id="stream_meshfile_do3" type="char*256" category="datasets"
input_pathname="abs" group="do3_streams" valid_values="" >
Filename of mesh file for diurnal ozone anomaly file
</entry>

<entry id="do3_mapalgo" type="char*256" category="datasets"
group="do3_streams" valid_values="bilinear,nn,nnoni,nnonj,spval,copy" >
Mapping method from diurnal ozone input file to the model resolution
bilinear = bilinear interpolation
nn = nearest neighbor
nnoni = nearest neighbor on the "i" (longitude) axis
nnonj = nearest neighbor on the "j" (latitude) axis
spval = set to special value
copy = copy using the same indices
</entry>

<!-- ======================================================================================== -->

<!-- ======================================================================================== -->
<!-- cropcal_streams streams Namelist -->
<!-- ======================================================================================== -->
Expand Down
116 changes: 116 additions & 0 deletions doc/design/design_doc_diurnal_ozone.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,116 @@
# Software Design Documentation
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for writing this design doc and adding it here! This is especially helpful when coming back to this a couple of years later like I'm doing now!!!


<span style="font-size:larger;"><b>Project Name</b></span>
<!-- Note you can add this to an issue or pull request in github. That is the normal usage of this template. -->

**Date:**: 09/27/22

**Written By**: Adrianna Foster (@adrifoster)

## Introduction

---------------------------------------

As laid out in CTSM Issue [#270](https://github.com/ESCOMP/CTSM/issues/270) we want to downscale input ozone partial pressure (mol/mol) temporally if we receive a multi-day average input ozone (usually in DATM mode).
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • I would change the parenthetical to say "(in DATM mode or when CAM is running without full chemistry)"


We have the infrastructure laid out for CAM and DATM to inform CTSM which type of input we are receiving (i.e., 'multiday_average' or 'subdaily'). We should only apply this downscaling when receiving multiday average ozone.

We have a gridded file provided by Louisa Emmons for a diurnal ozone variation factor. The data is additionally dimensioned by 'secs' (seconds of day), and depending on the model time of day, can be used as a multiplicative factor for converting multi-day average ozone to sub-daily ozone partial pressure.

## Solutions

---------------------------------------

Read in the diurnal anomaly file as a streams file, being careful to not assume too much so that other files may be provided in the future. Create a diurnal ozone type that has an `Interp` method to downscale and interpolate multi-day average ozone data to sub-daily based on a factor attribute. Inside the existing `CalcOzoneUptake` routine, if we are using multi-day average ozone, interpolate/downscale the `forc_o3` array using the `Interp` subroutine.

**Some alternate ideas**:

1. have the existing `ozone_type` type have this diurnal ozone anomaly type as an attribute that gets initiated only if we are reading in multi-day average ozone.
2. have the diurnal ozone type be fed into relevant `ozone_type` methods as arguments. We would only want to do this if other modules are going to use the diurnal ozone type. This has a possibility, but no clear plans are on the horizon for it.

Consider going with the first solution for now.

## Design Considerations

---------------------------------------

### Assumptions and Dependencies

Right now we are assuming that the units of the input diurnal file are going to be in seconds, and that the file covers the whole day.

There is no need to assume that the dimension will be 1-24 (i.e. on the hour) or that the dimensions be equal intervals.


## Design and Architecture

---------------------------------------

### System diagram or flowchart

#### Diurnal Factor Streams File

We need to read in the diurnal ozone factor file provided by Louisa Emmons as a streams file. Some modifications are required on the file so that it can be accurately read and parsed by ESMF and to facilitate our chosen implementation strategy:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Here and below: I would changed "parsed by ESMF" to "parsed by CDEPS"


1. Add a new 'time' variable (even though the data is not dimensioned by time, except for seconds of day) so that ESMF can parse it. The date is arbitrarily set to 2000-01-01, and the dimensions are set to 'UNLIMITED'.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is really helpful that you have laid out these changes and their reasons - thank you!


2. Shift the 'secs' array up 1800 seconds so that it provides the midpoint of the averaged time interval, rather than the beginning. This conversion is done to facilitate interpolation between the time-of-day dimensions.

These file modifications can be seen in the Jupyter notebook /glade/u/home/afoster/Diurnal_ozone.ipynb.

The file is read in using a new `src/share_esmf/diurnalOzoneStreamMod` module, with associated updates to `bld/CLMBuildNamelist.pm`, `bld/namelist_files/namelist_defaults_ctsm.xml` and `bld/namelist_files/namelist_definition_csvm.xml`
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • typo: csvm.xml -> ctsm.xml


Right now the new variable `use_do3_streams` is set in the `user_nl_clm` file. And the other namelist variables (`stream_fldfilename_do3`, `stream_meshfile_do3`, and `do3_mapalgo`) are inside the `lnd_in` file under a `do3_streams` namelist.

**I'm not sure setting up the `use_do3_streams` namelist variable is the best way to go about this, since we essentially want this feature on whenever the ozone frequency is 'multiday_average'**
Comment on lines +62 to +64
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • As I noted in a comment above, I agree with you here: we should probably remove use_do3_streams and just use the multiday_average value.


Right now we are hard-coding the `stream_var_name`, and `stream_lev_dimname`. I think this is okay. Otherwise, we would want them to be namelist variables that get read in.

Because the diurnal anomaly file is not dimensioned by time other than seconds of day, we just need to read it in and advance one time. This is all currently done in the `diurnalOzoneStreamMod`'s `read_O3_stream` subroutine. We set the date to an arbitrary date (the same as on the file, I'm not sure if this is necessary?).

We also need to grab the `secs` array to do interpolation/downscaling with.

Both are read from the file and then an instance of `diurnal_ozone_anom_type` is initialized and relevant arrays are set to the ozone factor and seconds data.

### Diurnal Ozone Anomaly Type

The module `DiurnalOzoneType` sets up a `diurnal_ozone_anom_type` which has only a few attributes and methods:

**Attributes**:

1. `ntimes` - size of time/seconds-of-day dimension (private)
2. `o3_anomaly_grc` - o3 anomaly data, 2d, gridcells x ntimes
3. `time_arr` - time/seconds of day array, 1d, ntimes

**Methods**:

1. `Init` - Initializes the anomaly data structures, calls `InitAllocate`
2. `InitAllocate` - allocates arrays and sets them to nan
3. `Interp` - Interpolates/downscales an input multi-day average ozone (`forc_o3`) using the o3 anomaly data and outputs a downscaled `forc_o3_down` array. See below for `Interp` algorithm/pseudo code
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • It looks like you never added the algorithm / pseudo-code. That's fine, but I'd just remove this reference to a non-existent section.


### Implementation within OzoneMod

We add an instance of the `diurnal_ozone_anom_type` as an attribute of the `ozone_base_type`.

Within the `Init` method of the `ozone_base_type`, we read the `drv_flds_in` file to get the `atom_ozone_frequency_val`. If this value is `atm_ozone_frequency_multiday_average` then we initialize and read in the o3 anomaly data by calling the `read_O3_stream` described above. *Do we need an else decision here?*
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Typo: atom -> atm

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good question about the else. I have added a comment on this below, in the relevant part of the code.


We also add an integer `atm_ozone_freq` as a new attribute, and set it in `Init`.

Within the `CalcOzoneUptake` method, we check for the `atm_ozone_freq` flag and if it is `atm_ozone_frequency_multiday_average` we call the `Interp` method.

### Algorithm and Pseudo code for Interp



## Rollout Plan

---------------------------------------
Define the roll-out phases and tests you plan to do


## Review Sign-off

---------------------------------------

* Reviewer(s):

* Sign-off Completed on YYYY-MM-DD*
156 changes: 156 additions & 0 deletions src/biogeophys/DiurnalOzoneType.F90
Original file line number Diff line number Diff line change
@@ -0,0 +1,156 @@
module DiurnalOzoneType

!-----------------------------------------------------------------------
! !DESCRIPTION:
! Sets up type for converting multi-day ozone input to sub-daily values using an input
! ozone anomaly file
!
! !USES:
#include "shr_assert.h"
use shr_kind_mod , only : r8 => shr_kind_r8
use decompMod , only : bounds_type
use clm_varcon , only : spval
use clm_varctl , only : iulog
use abortutils , only : endrun

implicit none
save
private

! !PUBLIC TYPES:
type, public :: diurnal_ozone_anom_type
private
! Private data members
integer :: ntimes ! size of time dimension
real(r8), public, allocatable :: o3_anomaly_grc(:,:) ! o3 anomaly data [grc, ntimes]
real(r8), public, allocatable :: time_arr(:) ! time dimension (units = seconds of day)

contains
! Public routines
procedure, public :: Init
procedure, private :: InitAllocate
procedure, public :: Interp

end type diurnal_ozone_anom_type

character(len=*), parameter, private :: sourcefile = &
__FILE__

contains

! ========================================================================
! Infrastructure routines (initialization, etc.)
! ========================================================================

!-----------------------------------------------------------------------
subroutine Init(this, bounds, n)
!
! DESCRIPTION:
! Initialize ozone anomaly data structures
!
!
! ARGUMENTS:
class(diurnal_ozone_anom_type), intent(inout) :: this
type(bounds_type), intent(in) :: bounds
integer, intent(in) :: n
!-----------------------------------------------------------------------

! set ntimes from input
this%ntimes = n

! allocate arrays
call this%InitAllocate(bounds)

end subroutine Init

!-----------------------------------------------------------------------
subroutine InitAllocate(this, bounds)
!
! !DESCRIPTION:
! Allocate module variables and data structures
!
! !USES:
use shr_infnan_mod, only: nan => shr_infnan_nan, assignment(=)
!
! !ARGUMENTS:
class(diurnal_ozone_anom_type), intent(inout) :: this
type(bounds_type), intent(in) :: bounds
!
! LOCAL VARIABLES:
integer :: begg, endg
!---------------------------------------------------------------------

begg = bounds%begg; endg = bounds%endg

allocate(this%o3_anomaly_grc(begg:endg,1:this%ntimes)); this%o3_anomaly_grc(:,:) = nan
allocate(this%time_arr(1:this%ntimes)); this%time_arr(:) = nan

end subroutine InitAllocate

!-----------------------------------------------------------------------

subroutine Interp(this, bounds, forc_o3, forc_o3_down)
!
! !DESCRIPTION:
! Downscale/interpolate multi-day ozone data to subdaily
!
! !USES:
use clm_time_manager , only : get_curr_date
use clm_varcon , only : secspday
!
! !ARGUMENTS:
class(diurnal_ozone_anom_type), intent(in) :: this
type(bounds_type), intent(in) :: bounds ! bounds type
real(r8), intent(in) :: forc_o3(:) ! ozone partial pressure (mol/mol)
real(r8), intent(out) :: forc_o3_down(:) ! ozone partial pressure, downscaled (mol/mol)

!
! LOCAL VARIABLES:
integer :: t ! looping index
integer :: yr ! year
integer :: mon ! month
integer :: day ! day of month
integer :: tod ! time of day (seconds past 0Z)
integer :: begg, endg ! bounds
integer :: t_prev ! previous time index
adrifoster marked this conversation as resolved.
Show resolved Hide resolved
real(r8) :: tdiff_end
real(r8) :: tdiff_start
real(r8) :: tdiff
!-----------------------------------------------------------------------

begg = bounds%begg; endg = bounds%endg

! Get current date/time - we really only need seconds
call get_curr_date(yr, mon, day, tod)

! find the time interval we are in
do t = 1, this%ntimes
if (real(tod) <= this%time_arr(t)) then
exit
end if
end do
Comment on lines +127 to +131
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Is it possible to have real(tod) > this%time_arr(this%ntimes)? My sense is that this case isn't currently handled, but I'm not sure if it's possible for it to arise based on the times on the file.


! interpolate, checking for edge cases
if (t == 1) then
! wrap around back
t_prev = this%ntimes
tdiff_end = secspday - this%time_arr(t_prev) + real(tod)
tdiff = this%time_arr(t) + secspday - this%time_arr(t_prev)
else
t_prev = t - 1
tdiff_end = real(tod) - this%time_arr(t_prev)
tdiff = this%time_arr(t) - this%time_arr(t_prev)
end if

tdiff_start = this%time_arr(t) - real(tod)

! interpolate
forc_o3_down(begg:endg) = forc_o3(begg:endg)* &
((this%o3_anomaly_grc(begg:endg, t_prev)*tdiff_start + &
this%o3_anomaly_grc(begg:endg, t_prev)*tdiff_end)/tdiff)
Comment on lines +148 to +150
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  • Is there a typo here? You use t_prev as the index for both; should one of them be t? I haven't thought carefully about which one should be t other than thinking: in the edge case where tod == this%time_arr(t), then I think we want to use the anomaly from time t exactly, right? I guess one other thing I'm thinking about here as I start to think about this more is that I wonder if the names of the tdiff_start and tdiff_end variables should be swapped: I haven't thought about this carefully, but it seems like tdiff_start gives the difference from the end of the time interval, which seems a bit counter-intuitive... but maybe I'm seeing this wrong.


end subroutine Interp

!-----------------------------------------------------------------------

end module DiurnalOzoneType
Loading
Loading