draft burkholder parameterization and tests

src/cross_sections/temperature_based.F90

@@ -20,6 +20,7 @@ module tuvx_cross_section_temperature_based
 
   integer, parameter :: PARAM_BASE          = 1
   integer, parameter :: PARAM_TAYLOR_SERIES = 2
+  integer, parameter :: PARAM_BURKHOLDER    = 3
 
   !> Calculator for temperature-based cross sections
   type, extends(cross_section_t) :: cross_section_temperature_based_t
@@ -63,6 +64,8 @@ function constructor( config, grid_warehouse, profile_warehouse )           &
     use tuvx_grid_warehouse,           only : grid_warehouse_t
     use tuvx_netcdf,                   only : netcdf_t
     use tuvx_profile_warehouse,        only : profile_warehouse_t
+    use tuvx_temperature_parameterization_burkholder,                         &
+        only : temperature_parameterization_burkholder_t
     use tuvx_temperature_parameterization_taylor_series,                      &
         only : temperature_parameterization_taylor_series_t
 
@@ -135,6 +138,9 @@ function constructor( config, grid_warehouse, profile_warehouse )           &
         if( param_type == "TAYLOR_SERIES" ) then
           allocate( this%parameterization_, source =                          &
             temperature_parameterization_taylor_series_t( param_config ) )
+        else if( param_type == "BURKHOLDER" ) then
+          allocate( this%parameterization_, source =                          &
+            temperature_parameterization_burkholder_t( param_config ) )
         else
           call die_msg( 370773773, "Invalid temperature-based "//             &
                         "parameterization type: '"//param_type//"'" )
@@ -274,6 +280,8 @@ subroutine mpi_pack( this, buffer, position, comm )
 
     use musica_assert,                 only : assert, die
     use musica_mpi,                    only : musica_mpi_pack
+    use tuvx_temperature_parameterization_burkholder,                         &
+        only : temperature_parameterization_burkholder_t
     use tuvx_temperature_parameterization_taylor_series,                      &
         only : temperature_parameterization_taylor_series_t
 
@@ -299,6 +307,8 @@ subroutine mpi_pack( this, buffer, position, comm )
           param_type = PARAM_BASE
         type is( temperature_parameterization_taylor_series_t )
           param_type = PARAM_TAYLOR_SERIES
+        type is( temperature_parameterization_burkholder_t )
+          param_type = PARAM_BURKHOLDER
         class default
           call die( 424852458 )
       end select
@@ -317,6 +327,8 @@ subroutine mpi_unpack( this, buffer, position, comm )
 
     use musica_assert,                 only : assert, die
     use musica_mpi,                    only : musica_mpi_unpack
+    use tuvx_temperature_parameterization_burkholder,                         &
+        only : temperature_parameterization_burkholder_t
     use tuvx_temperature_parameterization_taylor_series,                      &
         only : temperature_parameterization_taylor_series_t
 
@@ -342,6 +354,9 @@ subroutine mpi_unpack( this, buffer, position, comm )
         case( PARAM_TAYLOR_SERIES )
           allocate( temperature_parameterization_taylor_series_t ::          &
                     this%parameterization_ )
+        case( PARAM_BURKHOLDER )
+          allocate( temperature_parameterization_burkholder_t ::             &
+                    this%parameterization_ )
         case default
           call die( 324803089 )
       end select

src/cross_sections/util/CMakeLists.txt

@@ -4,6 +4,7 @@
 target_sources(tuvx_object
   PRIVATE
     temperature_parameterization.F90
+    temperature_parameterization_burkholder.F90
     temperature_parameterization_taylor_series.F90
     temperature_range.F90  
 )

src/cross_sections/util/temperature_parameterization_burkholder.F90

@@ -0,0 +1,244 @@
+! Copyright (C) 2024 National Center for Atmospheric Research
+! SPDX-License-Identifier: Apache-2.0
+
+module tuvx_temperature_parameterization_burkholder
+! Calculates cross-section elements using a temperature-based
+! parameterization from Burkholder et al. Phys. Chem. Chem. Phys. 4, 1432-1437 (2002).
+
+  ! Including musica_config at the module level to avoid an ICE
+  ! with Intel 2022.1 compiler
+  use musica_config,                   only : config_t
+  use musica_constants,                only : dk => musica_dk
+  use tuvx_temperature_parameterization,                                      &
+      only : temperature_parameterization_t
+  use tuvx_temperature_range,          only : temperature_range_t
+
+  implicit none
+
+  private
+  public :: temperature_parameterization_burkholder_t
+
+  !> Parameters for calculating cross section values based on
+  !! temperature using the algoritm in Burkholder et al.
+  !! Phys. Chem. Chem. Phys. 4, 1432-1437 (2002).
+  !!
+  !! Cross section elements are calculated as:
+  !!
+  !! \f[
+  !! Q(T) = 1 + e^{\frac{A}{B*T}}
+  !! \sigma(T,\lambda) = \frac{aa(\lambda)}{Q(T)} + bb(\lambda)*\[1-\frac{1}{Q(T)}\]
+  !! \f]
+  !!
+  !! where A, B, aa, and bb are constants, T is temperature [K] and \f$\lambda\f$ is
+  !! wavelength [nm].
+  type, extends(temperature_parameterization_t) :: temperature_parameterization_burkholder_t
+    real(kind=dk) :: A_
+    real(kind=dk) :: B_
+  contains
+    !> Calculate the cross section value for a specific temperature and wavelength
+    procedure :: calculate
+    !> Returns the number of bytes required to pack the parameterization
+    !! onto a character buffer
+    procedure :: pack_size => pack_size
+    !> Packs the parameterization onto a character buffer
+    procedure :: mpi_pack => mpi_pack
+    !> Unpacks the parameterization from a character buffer
+    procedure :: mpi_unpack => mpi_unpack
+  end type temperature_parameterization_burkholder_t
+
+  !> Constructor for temperature_parameterization_burkholder_t
+  interface temperature_parameterization_burkholder_t
+    module procedure :: constructor
+  end interface temperature_parameterization_burkholder_t
+
+contains
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+  !> Constructs a Burkholder (2002) temperature-based parameterization
+  function constructor( config ) result ( this )
+
+    use musica_assert,                 only : assert_msg
+    use musica_config,                 only : config_t
+    use musica_iterator,               only : iterator_t
+    use musica_string,                 only : string_t
+    use tuvx_grid,                     only : grid_t
+    use tuvx_netcdf,                   only : netcdf_t
+
+    type(temperature_parameterization_burkholder_t) :: this
+    type(config_t), intent(inout)                      :: config
+
+    character(len=*), parameter :: my_name =                                  &
+        "Burkholder (2002) temperature parameterization constructor"
+    type(string_t) :: required_keys(3), optional_keys(2), file_path
+    type(config_t) :: temp_ranges, temp_range, netcdf_file
+    class(iterator_t), pointer :: iter
+    type(netcdf_t) :: netcdf
+    integer :: i_range, i_param, n_param
+    logical :: found
+
+    required_keys(1) = "netcdf file"
+    required_keys(2) = "A"
+    required_keys(3) = "B"
+    optional_keys(1) = "type"
+    optional_keys(2) = "temperature ranges"
+    call assert_msg( 235183546,                                               &
+                     config%validate( required_keys, optional_keys ),         &
+                     "Bad configuration for Burkholder (2002) temperature "// &
+                     "parameterization." )
+
+    ! Load NetCDF file
+    call config%get( "netcdf file", netcdf_file, my_name )
+    call netcdf_file%get( "file path", file_path, my_name )
+    call netcdf%read_netcdf_file( file_path = file_path%to_char( ),           &
+                                  variable_name = "temperature_" )
+    n_param = size( netcdf%parameters, dim = 2 )
+    call assert_msg( 164185428, n_param >= 2, "Burkholder (2002) "//          &
+                     "parameterization must have at two sets of "//           &
+                     "coefficients" )
+
+    ! Load parameters
+    call config%get( "A", this%A_, my_name )
+    call config%get( "B", this%B_, my_name )
+    this%wavelengths_ = netcdf%wavelength(:)
+    this%AA_ = netcdf%parameters(:,1)
+    this%BB_ = netcdf%parameters(:,2)
+    call config%get( "temperature ranges", temp_ranges, my_name,              &
+                     found = found )
+    if( .not. found ) then
+      allocate( this%ranges_( 1 ) )
+      return
+    end if
+    allocate( this%ranges_( temp_ranges%number_of_children( ) ) )
+    iter => temp_ranges%get_iterator( )
+    i_range = 0
+    do while( iter%next( ) )
+      i_range = i_range + 1
+      call temp_ranges%get( iter, temp_range, my_name )
+      this%ranges_( i_range ) = temperature_range_t( temp_range )
+    end do
+    deallocate( iter )
+
+    ! initialize unused data members
+    allocate( this%lp_(0) )
+
+  end function constructor
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+  subroutine calculate( this, temperature, wavelengths, cross_section )
+
+    use tuvx_profile,                  only : profile_t
+
+    class(temperature_parameterization_burkholder_t), intent(in) :: this
+    real(kind=dk), intent(in)    :: temperature
+    real(kind=dk), intent(in)    :: wavelengths(:)
+    real(kind=dk), intent(inout) :: cross_section(:)
+
+    real(kind=dk) :: temp, Q
+    integer :: i_range, w_min, w_max
+
+    w_min = this%min_wavelength_index_
+    w_max = this%max_wavelength_index_
+    do i_range = 1, size( this%ranges_ )
+    associate( temp_range => this%ranges_( i_range ) )
+      if( temperature < temp_range%min_temperature_ .or.                      &
+          temperature > temp_range%max_temperature_ ) cycle
+      if( temp_range%is_fixed_ ) then
+        temp = temp_range%fixed_temperature_ - this%base_temperature_
+      else
+        temp = temperature - this%base_temperature_
+      end if
+      Q = 1.0 + exp( this%A_ / ( this%B_ * temp ) )
+      cross_section( w_min:w_max ) = ( this%AA_(:) / Q +                      &
+                                       this%BB_(:) * ( 1.0 - 1.0 / Q )        &
+                                     ) * 1.0e-20
+    end associate
+    end do
+
+  end subroutine calculate
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+  !> Returns the size of a character buffer required to pack the
+  !! parameterization
+  integer function pack_size( this, comm )
+
+    use musica_mpi,                    only : musica_mpi_pack_size
+
+    !> Parameterization to be packed
+    class(temperature_parameterization_burkholder_t), intent(in) :: this
+    !> MPI communicator
+    integer,                                             intent(in) :: comm
+
+#ifdef MUSICA_USE_MPI
+    pack_size = this%temperature_parameterization_t%pack_size( comm ) +       &
+                musica_mpi_pack_size( this%A_, comm ) +                       &
+                musica_mpi_pack_size( this%B_, comm )
+#else
+    pack_size = 0
+#endif
+
+  end function pack_size
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+  !> Packs the parameterization onto a character buffer
+  subroutine mpi_pack( this, buffer, position, comm )
+
+    use musica_assert,                 only : assert
+    use musica_mpi,                    only : musica_mpi_pack
+
+    !> Parameterization to be packed
+    class(temperature_parameterization_burkholder_t), intent(in) :: this
+    !> Memory buffer
+    character, intent(inout) :: buffer(:)
+    !> Current buffer position
+    integer,   intent(inout) :: position
+    !> MPI communicator
+    integer,   intent(in)    :: comm
+
+#ifdef MUSICA_USE_MPI
+    integer :: prev_pos
+
+    prev_pos = position
+    call this%temperature_parameterization_t%mpi_pack( buffer, position, comm )
+    call musica_mpi_pack( buffer, position, this%A_, comm )
+    call musica_mpi_pack( buffer, position, this%B_, comm )
+    call assert( 190816083, position - prev_pos <= this%pack_size( comm ) )
+#endif
+
+  end subroutine mpi_pack
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+  !> Unpacks a parameterization from a character buffer
+  subroutine mpi_unpack( this, buffer, position, comm )
+
+    use musica_assert,                 only : assert
+    use musica_mpi,                    only : musica_mpi_unpack
+
+    !> The parameterization to be unpacked
+    class(temperature_parameterization_burkholder_t), intent(out) :: this
+    !> Memory buffer
+    character, intent(inout) :: buffer(:)
+    !> Current buffer position
+    integer,   intent(inout) :: position
+    !> MPI communicator
+    integer,   intent(in)    :: comm
+
+#ifdef MUSICA_USE_MPI
+    integer :: prev_pos
+
+    prev_pos = position
+    call this%temperature_parameterization_t%mpi_unpack( buffer, position, comm )
+    call musica_mpi_unpack( buffer, position, this%A_, comm )
+    call musica_mpi_unpack( buffer, position, this%B_, comm )
+    call assert( 634825156, position - prev_pos <= this%pack_size( comm ) )
+#endif
+
+  end subroutine mpi_unpack
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+end module tuvx_temperature_parameterization_burkholder

src/cross_sections/util/temperature_parameterization_taylor_series.F90

@@ -263,6 +263,7 @@ subroutine mpi_unpack( this, buffer, position, comm )
     call this%temperature_parameterization_t%mpi_unpack( buffer, position, comm )
     call musica_mpi_unpack( buffer, position, this%sigma_,       comm )
     call musica_mpi_unpack( buffer, position, this%A_,           comm )
+    call assert( 966515884, position - prev_pos <= this%pack_size( comm ) )
 #endif
 
   end subroutine mpi_unpack

test/data/cross_sections/util/burkholder.config.json

@@ -0,0 +1,22 @@
+{
+  "type": "BURKHOLDER",
+  "netcdf file": {
+    "file path": "test/data/cross_sections/util/burkholder.nc"
+  },
+  "A": 12.5,
+  "B": 202.3,
+  "temperature ranges": [
+    {
+      "maximum": 209.999999999999,
+      "fixed value": 210.0
+    },
+    {
+      "minimum": 210.0,
+      "maximum": 300.0
+    },
+    {
+      "minimum": 300.00000000001,
+      "fixed value": 300.0
+    }
+  ]  
+}

test/unit/cross_section/util/CMakeLists.txt

@@ -6,6 +6,8 @@ include(test_util)
 ################################################################################
 # Cross section utility tests
 
+create_standard_test( NAME temperature_parameterization_burkholder
+                      SOURCES temperature_parameterization_burkholder.F90 )
 create_standard_test( NAME temperature_parameterization_taylor_series
                       SOURCES temperature_parameterization_taylor_series.F90 )