MUSICA TUVX scheme: create aerosol radiator, set_aerosol_optics_values (#182)

dwfncar · davidfillmore · David Fillmore · web-flow · commit d4bd2025a221 · 2025-01-09T13:40:35.000-07:00
Originator(s): @dwfncar Summary (include the keyword ['closes', 'fixes', 'resolves'] and issue number): - Closes #99 Describe any changes made to the namelist: N/A List all files eliminated and why: N/A List all files added and what they do: ``` A schemes/musica/tuvx/musica_ccpp_tuvx_aerosol_optics.F90 A test/musica/tuvx/test_tuvx_aerosol_optics.F90 ``` List all existing files that have been modified, and describe the changes: ``` M schemes/musica/tuvx/musica_ccpp_tuvx.F90 M schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 M test/docker/Dockerfile.musica M test/docker/Dockerfile.musica.no_install M test/musica/tuvx/CMakeLists.txt ``` List any test failures: N/A Is this a science-changing update? New physics package, algorithm change, tuning changes, etc? No --------- Co-authored-by: davidfillmore <fillmore.david.winslow@gmail.com> Co-authored-by: David Fillmore <fillmore@Andromeda.local> Co-authored-by: Jiwon Gim <jiwongim@ucar.edu>
diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90
@@ -23,6 +23,7 @@ module musica_ccpp_tuvx
   type(profile_t),         pointer :: surface_albedo_profile => null()
   type(profile_t),         pointer :: extraterrestrial_flux_profile => null()
   type(radiator_t),        pointer :: cloud_optics => null()
+  type(radiator_t),        pointer :: aerosol_optics => null()
   type(index_mappings_t),  pointer :: photolysis_rate_constants_mapping => null( )
   integer,               parameter :: DEFAULT_NUM_PHOTOLYSIS_RATE_CONSTANTS = 0
   integer                          :: number_of_photolysis_rate_constants = DEFAULT_NUM_PHOTOLYSIS_RATE_CONSTANTS
@@ -84,6 +85,11 @@ subroutine cleanup_tuvx_resources()
       cloud_optics => null()
     end if
 
+    if (associated( aerosol_optics )) then
+      deallocate( aerosol_optics )
+      aerosol_optics => null()
+    end if
+
     if (associated( photolysis_rate_constants_mapping )) then
       deallocate( photolysis_rate_constants_mapping )
       photolysis_rate_constants_mapping => null()
@@ -146,6 +152,8 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, &
             extraterrestrial_flux_unit
     use musica_ccpp_tuvx_cloud_optics, &
       only: create_cloud_optics_radiator, cloud_optics_label
+    use musica_ccpp_tuvx_aerosol_optics, &
+      only: create_aerosol_optics_radiator, aerosol_optics_label
 
     integer,                           intent(in)  :: vertical_layer_dimension      ! (count)
     integer,                           intent(in)  :: vertical_interface_dimension  ! (count)
@@ -278,6 +286,21 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, &
       return
     end if
 
+    aerosol_optics => create_aerosol_optics_radiator( height_grid, wavelength_grid, &
+                                                      errmsg, errcode )
+    if (errcode /= 0) then
+      call reset_tuvx_map_state( grids, profiles, radiators )
+      call cleanup_tuvx_resources()
+      return
+    endif
+
+    call radiators%add( aerosol_optics, error )
+    if (has_error_occurred( error, errmsg, errcode )) then
+      call reset_tuvx_map_state( grids, profiles, radiators )
+      call cleanup_tuvx_resources()
+      return
+    end if
+
     tuvx => tuvx_t( trim(filename_of_tuvx_configuration), grids, profiles, &
                     radiators, error )
     if (has_error_occurred( error, errmsg, errcode )) then
@@ -372,6 +395,15 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, &
       return
     end if
 
+    aerosol_optics => radiators%get( aerosol_optics_label, error )
+    if (has_error_occurred( error, errmsg, errcode )) then
+      deallocate( tuvx )
+      tuvx => null()
+      call reset_tuvx_map_state( grids, profiles, radiators )
+      call cleanup_tuvx_resources()
+      return
+    end if
+
     call reset_tuvx_map_state( grids, profiles, radiators )
 
     ! 'photolysis_rate_constants_ordering' is a local variable
@@ -432,6 +464,7 @@ subroutine tuvx_run(temperature, dry_air_density,                  &
     use musica_ccpp_tuvx_surface_albedo,           only: set_surface_albedo_values
     use musica_ccpp_tuvx_extraterrestrial_flux,    only: set_extraterrestrial_flux_values
     use musica_ccpp_tuvx_cloud_optics,             only: set_cloud_optics_values
+    use musica_ccpp_tuvx_aerosol_optics,           only: set_aerosol_optics_values
 
     real(kind_phys),    intent(in)    :: temperature(:,:)                                  ! K (column, layer)
     real(kind_phys),    intent(in)    :: dry_air_density(:,:)                              ! kg m-3 (column, layer)
@@ -502,6 +535,9 @@ subroutine tuvx_run(temperature, dry_air_density,                  &
                                       errmsg, errcode )
         if (errcode /= 0) return
 
+        call set_aerosol_optics_values( aerosol_optics, errmsg, errcode )
+        if (errcode /= 0) return
+
         ! calculate photolysis rate constants and heating rates
         call tuvx%run( solar_zenith_angle(i_col), earth_sun_distance, &
                        photolysis_rate_constants(:,:), heating_rates(:,:), &
@@ -540,4 +576,4 @@ subroutine tuvx_final(errmsg, errcode)
 
   end subroutine tuvx_final
 
-end module musica_ccpp_tuvx
+end module musica_ccpp_tuvx
diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_aerosol_optics.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_aerosol_optics.F90
@@ -0,0 +1,103 @@
+! Copyright (C) 2024 National Science Foundation-National Center for Atmospheric Research
+! SPDX-License-Identifier: Apache-2.0
+module musica_ccpp_tuvx_aerosol_optics
+  implicit none
+
+  private
+  public :: create_aerosol_optics_radiator, set_aerosol_optics_values
+
+  !> Label for aerosol optical properties in TUV-x
+  character(len=*), parameter, public :: aerosol_optics_label = "aerosols"
+  !> Label
+  character(len=*), parameter, public :: \
+    aerosol_optical_depth_label = "optical depths"
+  character(len=*), parameter, public :: \
+    aerosol_single_scattering_albedo_label = "single scattering albedos"
+  character(len=*), parameter, public :: \
+    aerosol_asymmetry_factor_label = "asymmetry factor"
+  !> Unit
+  character(len=*), parameter, public :: aerosol_optical_depth_unit = "none"
+  character(len=*), parameter, public :: aerosol_single_scattering_albedo_unit = "none"
+  character(len=*), parameter, public :: aerosol_asymmetry_factor_unit = "none"
+  !> Default value of number of vertical levels
+  integer, parameter :: DEFAULT_NUM_VERTICAL_LEVELS = 0
+  !> Number of vertical levels
+  integer, protected :: num_vertical_levels = DEFAULT_NUM_VERTICAL_LEVELS
+  !> Default value of number of wavelength bins
+  integer, parameter :: DEFAULT_NUM_WAVELENGTH_BINS = 0
+  !> Number of wavelength bins
+  integer, protected :: num_wavelength_bins = DEFAULT_NUM_WAVELENGTH_BINS
+  !> Default value of number of streams
+  integer, parameter :: DEFAULT_NUM_STREAMS = 1
+  !> Number of streams
+  integer, protected :: num_streams = DEFAULT_NUM_STREAMS
+
+contains
+
+  !> Creates a TUV-x aerosol optics radiator
+  function create_aerosol_optics_radiator( height_grid, wavelength_grid, &
+      errmsg, errcode ) result( radiator )
+    use musica_ccpp_util,     only: has_error_occurred
+    use musica_tuvx_grid,     only: grid_t
+    use musica_tuvx_radiator, only: radiator_t
+    use musica_util,          only: error_t
+
+    type(grid_t),     intent(inout) :: height_grid
+    type(grid_t),     intent(inout) :: wavelength_grid
+    character(len=*), intent(out)   :: errmsg
+    integer,          intent(out)   :: errcode
+    type(radiator_t), pointer       :: radiator
+
+    ! local variables
+    type(error_t) :: error
+
+    num_vertical_levels = height_grid%number_of_sections( error )
+    if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+    num_wavelength_bins = wavelength_grid%number_of_sections( error )
+    if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+    radiator => radiator_t( aerosol_optics_label, height_grid, wavelength_grid, &
+                            error )
+    if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+  end function create_aerosol_optics_radiator
+
+  !> Sets TUV-x aerosol optics values
+  ! Temporarily setting optical properties to zero until aerosol optical
+  ! property calculations are ported to CAM-SIMA.
+  subroutine set_aerosol_optics_values( radiator, errmsg, errcode )
+    use ccpp_kinds,           only: kind_phys
+    use musica_ccpp_util,     only: has_error_occurred
+    use musica_tuvx_radiator, only: radiator_t
+    use musica_util,          only: error_t
+
+    type(radiator_t), intent(inout) :: radiator
+    character(len=*), intent(out)   :: errmsg
+    integer,          intent(out)   :: errcode
+
+    ! local variables
+    type(error_t)   :: error
+    real(kind_phys) :: \
+      aerosol_optical_depth(num_vertical_levels, num_wavelength_bins)
+    real(kind_phys) :: \
+      aerosol_single_scattering_albedo(num_vertical_levels, num_wavelength_bins)
+    real(kind_phys) :: \
+      aerosol_asymmetry_factor(num_vertical_levels, num_wavelength_bins, num_streams)
+
+    aerosol_optical_depth(:,:) = 0.0_kind_phys
+    aerosol_single_scattering_albedo(:,:) = 0.0_kind_phys
+    aerosol_asymmetry_factor(:,:,:) = 0.0_kind_phys
+
+    call radiator%set_optical_depths( aerosol_optical_depth, error )
+    if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+    call radiator%set_single_scattering_albedos( aerosol_single_scattering_albedo, error )
+    if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+    call radiator%set_asymmetry_factors( aerosol_asymmetry_factor, error )
+    if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+  end subroutine set_aerosol_optics_values
+
+end module musica_ccpp_tuvx_aerosol_optics
diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90
@@ -72,4 +72,4 @@ subroutine set_surface_albedo_values( profile, host_surface_albedo, &
 
   end subroutine set_surface_albedo_values
 
-end module musica_ccpp_tuvx_surface_albedo
+end module musica_ccpp_tuvx_surface_albedo
diff --git a/test/docker/Dockerfile.musica b/test/docker/Dockerfile.musica
@@ -6,7 +6,7 @@
 FROM ubuntu:22.04
 
 ARG MUSICA_GIT_TAG=326b5119768d5be9654baf96ae3bd6a1b757fdc8
-ARG CAM_SIMA_CHEMISTRY_DATA_TAG=abc7cacbec3d33d5c0ed5bb79a157e93b42c45c0
+ARG CAM_SIMA_CHEMISTRY_DATA_TAG=be87bc14822aa50b1afda0059ab6f5b5bd7397e6
 ARG BUILD_TYPE=Debug
 
 RUN apt update \
@@ -92,4 +92,4 @@ RUN cd atmospheric_physics/test \
              -D CCPP_ENABLE_MEMCHECK=ON \
     && cmake --build ./build
 
-WORKDIR /home/test_user/atmospheric_physics/test/build
+WORKDIR /home/test_user/atmospheric_physics/test/build
diff --git a/test/docker/Dockerfile.musica.no_install b/test/docker/Dockerfile.musica.no_install
@@ -9,7 +9,7 @@
 FROM ubuntu:22.04
 
 ARG MUSICA_GIT_TAG=326b5119768d5be9654baf96ae3bd6a1b757fdc8
-ARG CAM_SIMA_CHEMISTRY_DATA_TAG=abc7cacbec3d33d5c0ed5bb79a157e93b42c45c0
+ARG CAM_SIMA_CHEMISTRY_DATA_TAG=be87bc14822aa50b1afda0059ab6f5b5bd7397e6
 ARG BUILD_TYPE=Debug
 
 RUN apt update \
@@ -80,4 +80,4 @@ RUN cd atmospheric_physics/test \
              -D CCPP_ENABLE_MEMCHECK=ON \
     && cmake --build ./build
 
-WORKDIR /home/test_user/atmospheric_physics/test/build
+WORKDIR /home/test_user/atmospheric_physics/test/build
diff --git a/test/musica/tuvx/CMakeLists.txt b/test/musica/tuvx/CMakeLists.txt
@@ -172,3 +172,33 @@ add_test(
 )
 
 add_memory_check_test(test_tuvx_cloud_optics $<TARGET_FILE:test_tuvx_cloud_optics> "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY})
+
+# Aerosol optics
+add_executable(test_tuvx_aerosol_optics test_tuvx_aerosol_optics.F90)
+
+target_sources(test_tuvx_aerosol_optics
+  PUBLIC
+    ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_height_grid.F90
+    ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_wavelength_grid.F90
+    ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_aerosol_optics.F90
+    ${MUSICA_SRC_PATH}/musica_ccpp_util.F90
+    ${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90
+)
+
+target_link_libraries(test_tuvx_aerosol_optics
+  PRIVATE
+    musica::musica-fortran
+)
+
+set_target_properties(test_tuvx_aerosol_optics
+  PROPERTIES
+  LINKER_LANGUAGE Fortran
+)
+
+add_test(
+  NAME test_tuvx_aerosol_optics
+  COMMAND $<TARGET_FILE:test_tuvx_aerosol_optics>
+  WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}
+)
+
+add_memory_check_test(test_tuvx_aerosol_optics $<TARGET_FILE:test_tuvx_aerosol_optics> "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY})
diff --git a/test/musica/tuvx/test_tuvx_aerosol_optics.F90 b/test/musica/tuvx/test_tuvx_aerosol_optics.F90
@@ -0,0 +1,85 @@
+! Copyright (C) 2024 National Science Foundation-National Center for Atmospheric Research
+! SPDX-License-Identifier: Apache-2.0
+program test_tuvx_aerosol_optics
+
+  use musica_ccpp_tuvx_aerosol_optics
+
+#define ASSERT(x) if (.not.(x)) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: x"; stop 1; endif
+#define ASSERT_NEAR( a, b, abs_error ) if( (abs(a - b) >= abs_error) .and. (abs(a - b) /= 0.0) ) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: a, b"; stop 1; endif
+
+  call test_create_aerosol_optics_radiator()
+
+contains
+
+  subroutine test_create_aerosol_optics_radiator()
+
+    use musica_util,                      only: error_t
+    use musica_ccpp_tuvx_height_grid,     only: create_height_grid
+    use musica_ccpp_tuvx_wavelength_grid, only: create_wavelength_grid
+    use musica_tuvx_grid,                 only: grid_t
+    use musica_tuvx_radiator,             only: radiator_t
+    use ccpp_kinds,                       only: kind_phys
+
+    integer, parameter        :: NUM_HOST_HEIGHT_MIDPOINTS = 2
+    integer, parameter        :: NUM_HOST_HEIGHT_INTERFACES = 3
+    integer, parameter        :: NUM_WAVELENGTH_MIDPOINTS = 3
+    integer, parameter        :: NUM_WAVELENGTH_INTERFACES = 4
+    real(kind_phys)           :: host_wavelength_interfaces(NUM_WAVELENGTH_INTERFACES) = [180.0e-9_kind_phys, 200.0e-9_kind_phys, 240.0e-9_kind_phys, 300.0e-9_kind_phys]
+    real(kind_phys)           :: aerosol_optical_depth(NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS)
+    real(kind_phys)           :: single_scattering_albedo(NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS)
+    real(kind_phys)           :: asymmetry_parameter(NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS,1)
+    type(grid_t),     pointer :: height_grid => null()
+    type(grid_t),     pointer :: wavelength_grid => null()
+    type(radiator_t), pointer :: aerosols => null()
+    type(error_t)             :: error
+    character(len=512)        :: errmsg
+    integer                   :: errcode
+    integer                   :: i
+
+    height_grid => create_height_grid(NUM_HOST_HEIGHT_MIDPOINTS, NUM_HOST_HEIGHT_INTERFACES, &
+                                      errmsg, errcode)
+    ASSERT(errcode == 0)
+    ASSERT(associated(height_grid))
+
+    wavelength_grid => create_wavelength_grid(host_wavelength_interfaces, errmsg, errcode)
+    ASSERT(errcode == 0)
+    ASSERT(associated(wavelength_grid))
+
+    aerosols => create_aerosol_optics_radiator(height_grid, wavelength_grid, errmsg, errcode)
+    ASSERT(errcode == 0)
+    ASSERT(associated(aerosols))
+
+    call set_aerosol_optics_values(aerosols, errmsg, errcode)
+    ASSERT(errcode == 0)
+
+    call aerosols%get_optical_depths(aerosol_optical_depth, error)
+    ASSERT(error%is_success())
+    do i = 1, size(aerosol_optical_depth, dim=1)
+      do j = 1, size(aerosol_optical_depth, dim=2)
+        ASSERT_NEAR(aerosol_optical_depth(i,j), 0.0_kind_phys, ABS_ERROR)
+      end do
+    end do
+
+    call aerosols%get_single_scattering_albedos(single_scattering_albedo, error)
+    ASSERT(error%is_success())
+    do i = 1, size(single_scattering_albedo, dim=1)
+      do j = 1, size(single_scattering_albedo, dim=2)
+        ASSERT_NEAR(single_scattering_albedo(i,j), 0.0_kind_phys, ABS_ERROR)
+      end do
+    end do
+
+    call aerosols%get_asymmetry_factors(asymmetry_parameter, error)
+    ASSERT(error%is_success())
+    do i = 1, size(asymmetry_parameter, dim=1)
+      do j = 1, size(asymmetry_parameter, dim=2)
+        ASSERT_NEAR(asymmetry_parameter(i,j,1), 0.0_kind_phys, ABS_ERROR)
+      end do
+    end do
+
+    deallocate( height_grid )
+    deallocate( wavelength_grid )
+    deallocate( aerosols )
+
+  end subroutine test_create_aerosol_optics_radiator
+
+end program test_tuvx_aerosol_optics

Original file line number	Diff line number	Diff line change
`@@ -72,4 +72,4 @@ subroutine set_surface_albedo_values( profile, host_surface_albedo, &`
`72`	`72`
`73`	`73`	`end subroutine set_surface_albedo_values`
`74`	`74`
`75`		`-end module musica_ccpp_tuvx_surface_albedo`
	`75`	`+end module musica_ccpp_tuvx_surface_albedo`