From 22e8c459538d542eda07dd13c2f4637bbad29c90 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Thu, 26 Sep 2024 11:49:00 -0600 Subject: [PATCH 01/13] add surface albedo --- .../tuvx/musica_ccpp_tuvx_surface_albedo.F90 | 69 +++++++++++++++++++ 1 file changed, 69 insertions(+) create mode 100644 musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 diff --git a/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 new file mode 100644 index 00000000..33c0f2dc --- /dev/null +++ b/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 @@ -0,0 +1,69 @@ +module musica_ccpp_tuvx_surface_albedo + implicit none + + private + public :: create_surface_albedo_profile, set_surface_albedo_values + + !> Label for surface albedo in TUV-x + character(len=*), parameter :: surface_albedo_label = "surface_albedo" + !> Units for surface albedo in TUV-x + character(len=*), parameter :: surface_albedo_units = "" + +contains + + !> Creates a TUV-x surface albedo profile from the host-model wavelength grid + function create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) & + result( profile ) + + use musica_ccpp_util, only: has_error_occurred + use musica_tuvx_grid, only: grid_t + use musica_tuvx_profile, only: profile_t + use musica_util, only: error_t + + ! Arguments + type(grid_t), intent(in) :: wavelength_grid + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + + ! Return value + type(profile_t), pointer :: profile + + ! Local variables + type(error_t) :: error + + profile => profile_t( surface_albedo_label, surface_albedo_units, & + wavelength_grid, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + end function create_surface_albedo_profile + + !> Sets TUV-x temperature edges from the host-model temperature midpoints + subroutine set_surface_albedo_values( profile, host_surface_albedo, & + errmsg, errcode ) + + use ccpp_kinds, only: kind_phys + use musica_ccpp_util, only: has_error_occurred + use musica_tuvx_profile, only: profile_t + use musica_util, only: error_t + + ! Arguments + type(profile_t), intent(inout) :: profile + real(kind_phys), intent(in) :: host_surface_albedo ! unitless + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + + ! Local variables + type(error_t) :: error + + ! TODO(jiwon) - is this correct? + ! real(r8) :: albedos(this%n_wavelength_bins_ + 1) + real(kind_phys) :: surface_albedo(size(host_surface_albedo)) + + surface_albedo = host_surface_albedo + + call profile%set_edge_values( surface_albedo, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + end subroutine set_surface_albedo + +end module musica_ccpp_tuvx_surface_albedo \ No newline at end of file From d82773eb912876dc970824306f2a4dadc010eda1 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Fri, 18 Oct 2024 17:50:01 -0600 Subject: [PATCH 02/13] add updating surface albedo functions --- .../tuvx/musica_ccpp_tuvx_surface_albedo.F90 | 44 ++++++++------- test/musica/CMakeLists.txt | 13 ++--- test/musica/tuvx/CMakeLists.txt | 31 ++++++++++- test/musica/tuvx/test_tuvx_surface_albedo.F90 | 53 +++++++++++++++++++ 4 files changed, 116 insertions(+), 25 deletions(-) create mode 100644 test/musica/tuvx/test_tuvx_surface_albedo.F90 diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 index 33c0f2dc..b6045162 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 @@ -6,8 +6,12 @@ module musica_ccpp_tuvx_surface_albedo !> Label for surface albedo in TUV-x character(len=*), parameter :: surface_albedo_label = "surface_albedo" - !> Units for surface albedo in TUV-x - character(len=*), parameter :: surface_albedo_units = "" + !> Unit for surface albedo in TUV-x + character(len=*), parameter :: surface_albedo_unit = "" ! unitless + !> 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 contains @@ -20,24 +24,26 @@ function create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) & use musica_tuvx_profile, only: profile_t use musica_util, only: error_t - ! Arguments - type(grid_t), intent(in) :: wavelength_grid - character(len=*), intent(out) :: errmsg - integer, intent(out) :: errcode + type(grid_t), intent(inout) :: wavelength_grid + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode - ! Return value + ! return value type(profile_t), pointer :: profile - ! Local variables + ! local variables type(error_t) :: error - profile => profile_t( surface_albedo_label, surface_albedo_units, & + num_wavelength_bins = wavelength_grid%number_of_sections( error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + profile => profile_t( surface_albedo_label, surface_albedo_unit, & wavelength_grid, error ) if ( has_error_occurred( error, errmsg, errcode ) ) return end function create_surface_albedo_profile - !> Sets TUV-x temperature edges from the host-model temperature midpoints + !> Sets TUV-x surface albedo value subroutine set_surface_albedo_values( profile, host_surface_albedo, & errmsg, errcode ) @@ -46,24 +52,26 @@ subroutine set_surface_albedo_values( profile, host_surface_albedo, & use musica_tuvx_profile, only: profile_t use musica_util, only: error_t - ! Arguments type(profile_t), intent(inout) :: profile real(kind_phys), intent(in) :: host_surface_albedo ! unitless character(len=*), intent(out) :: errmsg integer, intent(out) :: errcode - ! Local variables + ! local variables type(error_t) :: error + real(kind_phys) :: surface_albedos(num_wavelength_bins + 1) - ! TODO(jiwon) - is this correct? - ! real(r8) :: albedos(this%n_wavelength_bins_ + 1) - real(kind_phys) :: surface_albedo(size(host_surface_albedo)) + if (size(surface_albedos) <= DEFAULT_NUM_WAVELENGTH_BINS + 1) then + errmsg = "[MUSICA Error] Invalid dimension for number of wavelngth bins" + errcode = 1 + return + end if - surface_albedo = host_surface_albedo + surface_albedos(:) = host_surface_albedo - call profile%set_edge_values( surface_albedo, error ) + call profile%set_edge_values( surface_albedos, error ) if ( has_error_occurred( error, errmsg, errcode ) ) return - end subroutine set_surface_albedo + end subroutine set_surface_albedo_values end module musica_ccpp_tuvx_surface_albedo \ No newline at end of file diff --git a/test/musica/CMakeLists.txt b/test/musica/CMakeLists.txt index 487d181d..0e815ba3 100644 --- a/test/musica/CMakeLists.txt +++ b/test/musica/CMakeLists.txt @@ -20,14 +20,15 @@ FetchContent_MakeAvailable(musica) add_executable(test_musica_api test_musica_api.F90 musica_ccpp_namelist.F90) +file(GLOB MUSICA_CCPP_SOURCES + ${MUSICA_SRC_PATH}/*.F90 + ${MUSICA_SRC_PATH}/micm/*.F90 + ${MUSICA_SRC_PATH}/tuvx/*.F90 +) + target_sources(test_musica_api PUBLIC - ${MUSICA_SRC_PATH}/micm/musica_ccpp_micm.F90 - ${MUSICA_SRC_PATH}/micm/musica_ccpp_micm_util.F90 - ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx.F90 - ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_height_grid.F90 - ${MUSICA_SRC_PATH}/musica_ccpp.F90 - ${MUSICA_SRC_PATH}/musica_ccpp_util.F90 + ${MUSICA_CCPP_SOURCES} ${CCPP_SRC_PATH}/ccpp_constituent_prop_mod.F90 ${CCPP_SRC_PATH}/ccpp_hash_table.F90 ${CCPP_SRC_PATH}/ccpp_hashable.F90 diff --git a/test/musica/tuvx/CMakeLists.txt b/test/musica/tuvx/CMakeLists.txt index b863be59..8eb6678b 100644 --- a/test/musica/tuvx/CMakeLists.txt +++ b/test/musica/tuvx/CMakeLists.txt @@ -53,4 +53,33 @@ add_test( WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY} ) -add_memory_check_test(test_tuvx_temperature $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) \ No newline at end of file +add_memory_check_test(test_tuvx_temperature $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) + +# Surface albedo +add_executable(test_tuvx_surface_albedo test_tuvx_surface_albedo.F90) + +target_sources(test_tuvx_surface_albedo + PUBLIC + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_surface_albedo.F90 + ${MUSICA_SRC_PATH}/musica_ccpp_util.F90 + ${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90 +) + +target_link_libraries(test_tuvx_surface_albedo + PRIVATE + musica::musica-fortran +) + +set_target_properties(test_tuvx_surface_albedo + PROPERTIES + LINKER_LANGUAGE Fortran +) + +add_test( + NAME test_tuvx_surface_albedo + COMMAND $ + WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY} +) + +add_memory_check_test(test_tuvx_surface_albedo $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) \ No newline at end of file diff --git a/test/musica/tuvx/test_tuvx_surface_albedo.F90 b/test/musica/tuvx/test_tuvx_surface_albedo.F90 new file mode 100644 index 00000000..cf1f6ac1 --- /dev/null +++ b/test/musica/tuvx/test_tuvx_surface_albedo.F90 @@ -0,0 +1,53 @@ +program test_tuvx_surface_albedo + + use musica_ccpp_tuvx_surface_albedo + + implicit none + +#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_update_surface_albedo() + +contains + + subroutine test_update_surface_albedo() + use musica_ccpp_tuvx_wavelength_grid, only: create_wavelength_grid + use musica_util, only: error_t + use musica_tuvx_grid, only: grid_t + use musica_tuvx_profile, only: profile_t + use ccpp_kinds, only: kind_phys + + integer, parameter :: NUM_WAVELENGTH_BIN = 4 + type(grid_t), pointer :: wavelength_grid + type(profile_t), pointer :: profile + real(kind_phys), target :: host_surface_albedo = 500.5_kind_phys + real(kind_phys) :: surface_albedos(NUM_WAVELENGTH_BIN + 1) + character(len=512) :: errmsg + integer :: errcode + type(error_t) :: error + real(kind_phys) :: abs_error = 1e-4 + integer :: i + + wavelength_grid => create_wavelength_grid(NUM_WAVELENGTH_BIN, errmsg, errcode) + profile => create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) + ASSERT(errcode == 0) + ASSERT(associated(profile)) + + call set_surface_albedo_values( profile, host_surface_albedo, errmsg, errcode ) + ASSERT(errcode == 0) + + call profile%get_edge_values( surface_albedos, error) + ASSERT(error%is_success()) + ASSERT_NEAR(surface_albedos(1), 500.5, abs_error) + ASSERT_NEAR(surface_albedos(2), 500.5, abs_error) + ASSERT_NEAR(surface_albedos(3), 500.5, abs_error) + ASSERT_NEAR(surface_albedos(4), 500.5, abs_error) + ASSERT_NEAR(surface_albedos(5), 500.5, abs_error) + + deallocate( profile ) + deallocate( wavelength_grid ) + + end subroutine test_update_surface_albedo + +end program test_tuvx_surface_albedo \ No newline at end of file From 9413340ed692d6df62de51043d89fdd3c680d565 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Fri, 18 Oct 2024 17:51:09 -0600 Subject: [PATCH 03/13] add wavelength grid --- .../tuvx/musica_ccpp_tuvx_wavelength_grid.F90 | 37 +++++++++++++++++++ 1 file changed, 37 insertions(+) create mode 100644 schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 new file mode 100644 index 00000000..28f0e37f --- /dev/null +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 @@ -0,0 +1,37 @@ +module musica_ccpp_tuvx_wavelength_grid + + implicit none + + private + public :: create_wavelength_grid + + !> Label for height grid in TUV-x + character(len=*), parameter, public :: wavelength_grid_label = "wavelength" + !> Units for height grid in TUV-x + character(len=*), parameter, public :: wavelength_grid_unit = "nm" + +contains + + !> Creates a TUV-x wavelength grid + function create_wavelength_grid( num_wavelength_bin, errmsg, errcode ) & + result( wavelength_grid ) + + use musica_ccpp_util, only: has_error_occurred + use musica_tuvx_grid, only: grid_t + use musica_util, only: error_t + + integer, intent(in) :: num_wavelength_bin + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + type(grid_t), pointer :: wavelength_grid + + ! local variable + type(error_t) :: error + + wavelength_grid => grid_t( wavelength_grid_label, wavelength_grid_unit, & + num_wavelength_bin, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + end function create_wavelength_grid + +end module musica_ccpp_tuvx_wavelength_grid \ No newline at end of file From e9ece543e02c49c3130d8c8fd7b5d1f1f5c5bbc2 Mon Sep 17 00:00:00 2001 From: Matt Dawson Date: Wed, 23 Oct 2024 17:15:17 -0700 Subject: [PATCH 04/13] set wavelength grid from host (#142) --- schemes/musica/musica_ccpp.F90 | 7 +- schemes/musica/musica_ccpp.meta | 6 + schemes/musica/tuvx/musica_ccpp_tuvx.F90 | 54 +++++++- .../tuvx/musica_ccpp_tuvx_wavelength_grid.F90 | 31 ++++- test/musica/test_musica_api.F90 | 115 +++++++++++++++++- test/musica/tuvx/CMakeLists.txt | 28 +++++ test/musica/tuvx/configs/ts1_tsmlt.json | 6 - test/musica/tuvx/test_tuvx_surface_albedo.F90 | 4 +- .../musica/tuvx/test_tuvx_wavelength_grid.F90 | 53 ++++++++ 9 files changed, 287 insertions(+), 17 deletions(-) create mode 100644 test/musica/tuvx/test_tuvx_wavelength_grid.F90 diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90 index c486bee5..f0e95873 100644 --- a/schemes/musica/musica_ccpp.F90 +++ b/schemes/musica/musica_ccpp.F90 @@ -26,14 +26,17 @@ end subroutine musica_ccpp_register !> \section arg_table_musica_ccpp_init Argument Table !! \htmlinclude musica_ccpp_init.html subroutine musica_ccpp_init(vertical_layer_dimension, vertical_interface_dimension, & - errmsg, errcode) + photolysis_wavelength_grid_interfaces, errmsg, errcode) + use ccpp_kinds, only : kind_phys + integer, intent(in) :: vertical_layer_dimension ! (count) integer, intent(in) :: vertical_interface_dimension ! (count) + real(kind_phys), intent(in) :: photolysis_wavelength_grid_interfaces(:) ! m character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode call tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & - errmsg, errcode) + photolysis_wavelength_grid_interfaces, errmsg, errcode) if (errcode /= 0) return call micm_init(errmsg, errcode) if (errcode /= 0) return diff --git a/schemes/musica/musica_ccpp.meta b/schemes/musica/musica_ccpp.meta index f62966f9..97a07721 100644 --- a/schemes/musica/musica_ccpp.meta +++ b/schemes/musica/musica_ccpp.meta @@ -19,6 +19,12 @@ type = integer dimensions = () intent = in +[ photolysis_wavelength_grid_interfaces ] + standard_name = photolysis_wavelength_grid_interfaces + units = m + type = real | kind = kind_phys + dimensions = (photolysis_wavelength_grid_interface_dimension) + intent = in [ errmsg ] standard_name = ccpp_error_message units = none diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index d39eb5f4..c6445927 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -14,19 +14,24 @@ module musica_ccpp_tuvx type(tuvx_t), pointer :: tuvx => null( ) type(grid_t), pointer :: height_grid => null( ) + type(grid_t), pointer :: wavelength_grid => null( ) contains !> Intitialize TUV-x subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & - errmsg, errcode) + wavelength_grid_interfaces, errmsg, errcode) use musica_tuvx, only: grid_map_t, profile_map_t, radiator_map_t use musica_util, only: error_t use musica_ccpp_tuvx_height_grid, only: create_height_grid, & height_grid_label, height_grid_unit + use musica_ccpp_tuvx_wavelength_grid, only: create_wavelength_grid, & + wavelength_grid_label, & + wavelength_grid_unit - integer, intent(in) :: vertical_layer_dimension ! (count) - integer, intent(in) :: vertical_interface_dimension ! (count) + integer, intent(in) :: vertical_layer_dimension ! (count) + integer, intent(in) :: vertical_interface_dimension ! (count) + real(kind_phys), intent(in) :: wavelength_grid_interfaces(:) ! m character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode @@ -54,11 +59,32 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & return end if + wavelength_grid => create_wavelength_grid( wavelength_grid_interfaces, & + errmsg, errcode ) + if (errcode /= 0) then + deallocate( grids ) + deallocate( height_grid ) + height_grid => null() + return + endif + + call grids%add( wavelength_grid, error ) + if (has_error_occurred( error, errmsg, errcode )) then + deallocate( grids ) + deallocate( height_grid ) + height_grid => null() + deallocate( wavelength_grid ) + wavelength_grid => null() + return + end if + profiles => profile_map_t( error ) if (has_error_occurred( error, errmsg, errcode )) then deallocate( grids ) deallocate( height_grid ) height_grid => null() + deallocate( wavelength_grid ) + wavelength_grid => null() return end if @@ -68,6 +94,8 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & deallocate( profiles ) deallocate( height_grid ) height_grid => null() + deallocate( wavelength_grid ) + wavelength_grid => null() return end if @@ -79,6 +107,8 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & deallocate( radiators ) deallocate( height_grid ) height_grid => null() + deallocate( wavelength_grid ) + wavelength_grid => null() return end if @@ -87,6 +117,8 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & deallocate( radiators ) deallocate( height_grid ) height_grid => null() + deallocate( wavelength_grid ) + wavelength_grid => null() grids => tuvx%get_grids( error ) if (has_error_occurred( error, errmsg, errcode )) then @@ -103,6 +135,17 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & return end if + wavelength_grid => grids%get( wavelength_grid_label, wavelength_grid_unit, & + error ) + if (has_error_occurred( error, errmsg, errcode )) then + deallocate( tuvx ) + tuvx => null() + deallocate( grids ) + deallocate( height_grid ) + height_grid => null() + return + end if + deallocate( grids ) end subroutine tuvx_init @@ -166,6 +209,11 @@ subroutine tuvx_final(errmsg, errcode) height_grid => null() end if + if (associated( wavelength_grid )) then + deallocate( wavelength_grid ) + wavelength_grid => null() + end if + if (associated( tuvx )) then deallocate( tuvx ) tuvx => null() diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 index 28f0e37f..7fd7c1ad 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 @@ -5,6 +5,17 @@ module musica_ccpp_tuvx_wavelength_grid private public :: create_wavelength_grid + ! TUV-x Wavelegnth grid notes + ! + !----------------------------------------------------------------------- + ! The wavelength grid used with TUV-x is based on the grid used in the + ! CAM-Chem photolysis rate constant lookup tables. Slight modifications + ! were made to the grid in the Shumann-Runge and Lyman-alpha regions to + ! work with the expectations of the TUV-x code. + ! + ! The wavelength grid is defined by the host model. Any wavelength- + ! resolved quantities passed to TUV-x must be on this grid. + !> Label for height grid in TUV-x character(len=*), parameter, public :: wavelength_grid_label = "wavelength" !> Units for height grid in TUV-x @@ -13,23 +24,35 @@ module musica_ccpp_tuvx_wavelength_grid contains !> Creates a TUV-x wavelength grid - function create_wavelength_grid( num_wavelength_bin, errmsg, errcode ) & + function create_wavelength_grid( wavelength_grid_interfaces, errmsg, errcode ) & result( wavelength_grid ) + use ccpp_kinds, only: kind_phys use musica_ccpp_util, only: has_error_occurred + use musica_config, only: config_t use musica_tuvx_grid, only: grid_t use musica_util, only: error_t - integer, intent(in) :: num_wavelength_bin + real(kind_phys), intent(in) :: wavelength_grid_interfaces(:) ! m character(len=*), intent(out) :: errmsg integer, intent(out) :: errcode type(grid_t), pointer :: wavelength_grid - ! local variable + ! local variables + real(kind_phys) :: interfaces( size( wavelength_grid_interfaces ) ) ! [nm] + reaL(kind_phys) :: midpoints( size( wavelength_grid_interfaces ) - 1 ) ! [nm] type(error_t) :: error + interfaces(:) = wavelength_grid_interfaces(:) * 1.0e9 + midpoints(:) = & + 0.5 * ( interfaces( 1: size( interfaces ) - 1 ) & + + interfaces( 2: size( interfaces ) ) ) wavelength_grid => grid_t( wavelength_grid_label, wavelength_grid_unit, & - num_wavelength_bin, error ) + size( midpoints ), error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + call wavelength_grid%set_edges( interfaces, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + call wavelength_grid%set_midpoints( midpoints, error ) if ( has_error_occurred( error, errmsg, errcode ) ) return end function create_wavelength_grid diff --git a/test/musica/test_musica_api.F90 b/test/musica/test_musica_api.F90 index 9fdcfd1f..3400de5b 100644 --- a/test/musica/test_musica_api.F90 +++ b/test/musica/test_musica_api.F90 @@ -21,10 +21,12 @@ subroutine test_musica_ccpp_api() integer, parameter :: NUM_SPECIES = 4 integer, parameter :: NUM_COLUMNS = 2 integer, parameter :: NUM_LAYERS = 2 + integer, parameter :: NUM_WAVELENGTH_BINS = 102 integer :: solver_type integer :: errcode character(len=512) :: errmsg real(kind_phys) :: time_step ! s + real(kind_phys), dimension(NUM_WAVELENGTH_BINS+1) :: photolysis_wavelength_grid_interfaces ! m real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: geopotential_height_wrt_surface_at_midpoint ! m real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS+1) :: geopotential_height_wrt_surface_at_interface ! m real(kind_phys), dimension(NUM_COLUMNS) :: surface_geopotential ! m2 s-2 @@ -47,6 +49,116 @@ subroutine test_musica_ccpp_api() solver_type = Rosenbrock num_grid_cells = NUM_COLUMNS * NUM_LAYERS time_step = 60._kind_phys + ! These are the values that will be used in CAM-SIMA and correspond to the wavelength + ! bins used in the CAM-Chem photolysis rate constant lookup table. + ! + ! We're using the actual values here because several of the TS1/TSMLT photolysis + ! rate constant configurations are sensitive to the wavelength grid. + photolysis_wavelength_grid_interfaces = (/ & + 120.0e-9_kind_phys, & + 121.4e-9_kind_phys, & + 121.9e-9_kind_phys, & + 123.5e-9_kind_phys, & + 124.3e-9_kind_phys, & + 125.5e-9_kind_phys, & + 126.3e-9_kind_phys, & + 127.1e-9_kind_phys, & + 130.1e-9_kind_phys, & + 131.1e-9_kind_phys, & + 135.0e-9_kind_phys, & + 140.0e-9_kind_phys, & + 145.0e-9_kind_phys, & + 150.0e-9_kind_phys, & + 155.0e-9_kind_phys, & + 160.0e-9_kind_phys, & + 165.0e-9_kind_phys, & + 168.0e-9_kind_phys, & + 171.0e-9_kind_phys, & + 173.0e-9_kind_phys, & + 174.4e-9_kind_phys, & + 175.4e-9_kind_phys, & + 177.0e-9_kind_phys, & + 178.6e-9_kind_phys, & + 180.2e-9_kind_phys, & + 181.8e-9_kind_phys, & + 183.5e-9_kind_phys, & + 185.2e-9_kind_phys, & + 186.9e-9_kind_phys, & + 188.7e-9_kind_phys, & + 190.5e-9_kind_phys, & + 192.3e-9_kind_phys, & + 194.2e-9_kind_phys, & + 196.1e-9_kind_phys, & + 198.0e-9_kind_phys, & + 200.0e-9_kind_phys, & + 202.0e-9_kind_phys, & + 204.1e-9_kind_phys, & + 206.2e-9_kind_phys, & + 208.0e-9_kind_phys, & + 211.0e-9_kind_phys, & + 214.0e-9_kind_phys, & + 217.0e-9_kind_phys, & + 220.0e-9_kind_phys, & + 223.0e-9_kind_phys, & + 226.0e-9_kind_phys, & + 229.0e-9_kind_phys, & + 232.0e-9_kind_phys, & + 235.0e-9_kind_phys, & + 238.0e-9_kind_phys, & + 241.0e-9_kind_phys, & + 244.0e-9_kind_phys, & + 247.0e-9_kind_phys, & + 250.0e-9_kind_phys, & + 253.0e-9_kind_phys, & + 256.0e-9_kind_phys, & + 259.0e-9_kind_phys, & + 263.0e-9_kind_phys, & + 267.0e-9_kind_phys, & + 271.0e-9_kind_phys, & + 275.0e-9_kind_phys, & + 279.0e-9_kind_phys, & + 283.0e-9_kind_phys, & + 287.0e-9_kind_phys, & + 291.0e-9_kind_phys, & + 295.0e-9_kind_phys, & + 298.5e-9_kind_phys, & + 302.5e-9_kind_phys, & + 305.5e-9_kind_phys, & + 308.5e-9_kind_phys, & + 311.5e-9_kind_phys, & + 314.5e-9_kind_phys, & + 317.5e-9_kind_phys, & + 322.5e-9_kind_phys, & + 327.5e-9_kind_phys, & + 332.5e-9_kind_phys, & + 337.5e-9_kind_phys, & + 342.5e-9_kind_phys, & + 347.5e-9_kind_phys, & + 350.0e-9_kind_phys, & + 355.0e-9_kind_phys, & + 360.0e-9_kind_phys, & + 365.0e-9_kind_phys, & + 370.0e-9_kind_phys, & + 375.0e-9_kind_phys, & + 380.0e-9_kind_phys, & + 385.0e-9_kind_phys, & + 390.0e-9_kind_phys, & + 395.0e-9_kind_phys, & + 400.0e-9_kind_phys, & + 405.0e-9_kind_phys, & + 410.0e-9_kind_phys, & + 415.0e-9_kind_phys, & + 420.0e-9_kind_phys, & + 430.0e-9_kind_phys, & + 440.0e-9_kind_phys, & + 450.0e-9_kind_phys, & + 500.0e-9_kind_phys, & + 550.0e-9_kind_phys, & + 600.0e-9_kind_phys, & + 650.0e-9_kind_phys, & + 700.0e-9_kind_phys, & + 750.0e-9_kind_phys & + /) geopotential_height_wrt_surface_at_midpoint(1,:) = (/ 2000.0_kind_phys, 500.0_kind_phys /) geopotential_height_wrt_surface_at_midpoint(2,:) = (/ 2000.0_kind_phys, -500.0_kind_phys /) geopotential_height_wrt_surface_at_interface(1,:) = (/ 3000.0_kind_phys, 1000.0_kind_phys, 0.0_kind_phys /) @@ -96,7 +208,8 @@ subroutine test_musica_ccpp_api() call constituent_props_ptr(i)%set(const_prop, errcode, errmsg) end do - call musica_ccpp_init(NUM_LAYERS, NUM_LAYERS+1, errmsg, errcode) + call musica_ccpp_init(NUM_LAYERS, NUM_LAYERS+1, photolysis_wavelength_grid_interfaces, & + errmsg, errcode) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 diff --git a/test/musica/tuvx/CMakeLists.txt b/test/musica/tuvx/CMakeLists.txt index 8eb6678b..2c8623fd 100644 --- a/test/musica/tuvx/CMakeLists.txt +++ b/test/musica/tuvx/CMakeLists.txt @@ -26,6 +26,34 @@ add_test( add_memory_check_test(test_tuvx_height_grid $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) +# Wavelength grid +add_executable(test_tuvx_wavelength_grid test_tuvx_wavelength_grid.F90) + +target_sources(test_tuvx_wavelength_grid + PUBLIC + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 + ${MUSICA_SRC_PATH}/musica_ccpp_util.F90 + ${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90 +) + +target_link_libraries(test_tuvx_wavelength_grid + PRIVATE + musica::musica-fortran +) + +set_target_properties(test_tuvx_wavelength_grid + PROPERTIES + LINKER_LANGUAGE Fortran +) + +add_test( + NAME test_tuvx_wavelength_grid + COMMAND $ + WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY} +) + +add_memory_check_test(test_tuvx_wavelength_grid $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) + # Temperature add_executable(test_tuvx_temperature test_tuvx_temperature.F90) diff --git a/test/musica/tuvx/configs/ts1_tsmlt.json b/test/musica/tuvx/configs/ts1_tsmlt.json index 817ae96b..82e8d985 100644 --- a/test/musica/tuvx/configs/ts1_tsmlt.json +++ b/test/musica/tuvx/configs/ts1_tsmlt.json @@ -4,12 +4,6 @@ "cross section parameters file": "data/cross_sections/O2_parameters.txt" }, "grids": [ - { - "name": "wavelength", - "type": "from csv file", - "units": "nm", - "file path": "data/grids/wavelength/cam.csv" - }, { "name": "time", "type": "from config file", diff --git a/test/musica/tuvx/test_tuvx_surface_albedo.F90 b/test/musica/tuvx/test_tuvx_surface_albedo.F90 index cf1f6ac1..d3ae0997 100644 --- a/test/musica/tuvx/test_tuvx_surface_albedo.F90 +++ b/test/musica/tuvx/test_tuvx_surface_albedo.F90 @@ -22,6 +22,8 @@ subroutine test_update_surface_albedo() type(grid_t), pointer :: wavelength_grid type(profile_t), pointer :: profile real(kind_phys), target :: host_surface_albedo = 500.5_kind_phys + real(kind_phys) :: wavelength_grid_interfaces(NUM_WAVELENGTH_BIN + 1) = & + [200.0e-9_kind_phys, 210.0e-9_kind_phys, 240.0e-9_kind_phys, 300.0e-9_kind_phys, 400.0e-9_kind_phys] real(kind_phys) :: surface_albedos(NUM_WAVELENGTH_BIN + 1) character(len=512) :: errmsg integer :: errcode @@ -29,7 +31,7 @@ subroutine test_update_surface_albedo() real(kind_phys) :: abs_error = 1e-4 integer :: i - wavelength_grid => create_wavelength_grid(NUM_WAVELENGTH_BIN, errmsg, errcode) + wavelength_grid => create_wavelength_grid(wavelength_grid_interfaces, errmsg, errcode) profile => create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) ASSERT(errcode == 0) ASSERT(associated(profile)) diff --git a/test/musica/tuvx/test_tuvx_wavelength_grid.F90 b/test/musica/tuvx/test_tuvx_wavelength_grid.F90 new file mode 100644 index 00000000..7e974de8 --- /dev/null +++ b/test/musica/tuvx/test_tuvx_wavelength_grid.F90 @@ -0,0 +1,53 @@ +program test_tuvx_wavelength_grid + + use musica_ccpp_tuvx_wavelength_grid + + implicit none + +#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_wavelength_grid() + +contains + + subroutine test_create_wavelength_grid() + use musica_util, only: error_t + use musica_tuvx_grid, only: grid_t + use ccpp_kinds, only: kind_phys + + integer, parameter :: NUM_WAVELENGTH_GRID_MIDPOINTS = 2 + integer, parameter :: NUM_WAVELENGTH_GRID_INTERFACES = 3 + real(kind_phys), target :: host_interfaces(NUM_WAVELENGTH_GRID_INTERFACES) = [180.0e-9_kind_phys, 200.0e-9_kind_phys, 240.0e-9_kind_phys] + real(kind_phys), target :: expected_interfaces(NUM_WAVELENGTH_GRID_INTERFACES) = [180.0_kind_phys, 200.0_kind_phys, 240.0_kind_phys] + real(kind_phys), target :: expected_midpoints(NUM_WAVELENGTH_GRID_MIDPOINTS) = [190.0_kind_phys, 220.0_kind_phys] + real(kind_phys) :: interfaces(NUM_WAVELENGTH_GRID_INTERFACES) + real(kind_phys) :: midpoints(NUM_WAVELENGTH_GRID_MIDPOINTS) + type(grid_t), pointer :: wavelength_grid => null() + character(len=512) :: errmsg + integer :: errcode + real(kind_phys) :: abs_error = 1e-5 + integer :: i + type(error_t) :: error + + wavelength_grid => create_wavelength_grid(host_interfaces, errmsg, errcode) + ASSERT(errcode == 0) + ASSERT(associated(wavelength_grid)) + + call wavelength_grid%get_edges(interfaces, error) + ASSERT(error%is_success()) + do i = 1, NUM_WAVELENGTH_GRID_INTERFACES + ASSERT_NEAR(interfaces(i), expected_interfaces(i), abs_error) + end do + + call wavelength_grid%get_midpoints(midpoints, error) + ASSERT(error%is_success()) + do i = 1, NUM_WAVELENGTH_GRID_MIDPOINTS + ASSERT_NEAR(midpoints(i), expected_midpoints(i), abs_error) + end do + + deallocate(wavelength_grid) + + end subroutine test_create_wavelength_grid + +end program test_tuvx_wavelength_grid \ No newline at end of file From 51c6f2fd1eedad337da0fe342e31d50d70cc8f43 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Thu, 24 Oct 2024 11:14:39 -0600 Subject: [PATCH 05/13] fix a syntax bug --- schemes/musica/tuvx/musica_ccpp_tuvx.F90 | 1 + 1 file changed, 1 insertion(+) diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index 31a1ea74..46737a5c 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -276,6 +276,7 @@ subroutine tuvx_final(errmsg, errcode) if (associated( wavelength_grid )) then deallocate( wavelength_grid ) wavelength_grid => null() + end if if (associated( temperature_profile )) then deallocate( temperature_profile ) From ceadfb4c401ea5a7c259bb26637a35413a092b53 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Thu, 24 Oct 2024 13:56:08 -0600 Subject: [PATCH 06/13] add helper function for deallocating objects associated with TUV-x --- schemes/musica/tuvx/musica_ccpp_tuvx.F90 | 117 +++++------------- schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 | 49 ++++++++ .../tuvx/musica_ccpp_tuvx_wavelength_grid.F90 | 7 +- test/musica/tuvx/test_tuvx_surface_albedo.F90 | 3 + 4 files changed, 85 insertions(+), 91 deletions(-) create mode 100644 schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index 46737a5c..ecf8d70e 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -24,13 +24,13 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & wavelength_grid_interfaces, errmsg, errcode) use musica_tuvx, only: grid_map_t, profile_map_t, radiator_map_t use musica_util, only: error_t - use musica_ccpp_tuvx_height_grid, only: create_height_grid, & - height_grid_label, height_grid_unit - use musica_ccpp_tuvx_wavelength_grid, only: create_wavelength_grid, & - wavelength_grid_label, & - wavelength_grid_unit - use musica_ccpp_tuvx_temperature, only: create_temperature_profile, & - temperature_label, temperature_unit + use musica_ccpp_tuvx_util, only: tuvx_deallocate + use musica_ccpp_tuvx_height_grid, & + only: create_height_grid, height_grid_label, height_grid_unit + use musica_ccpp_tuvx_wavelength_grid, & + only: create_wavelength_grid, wavelength_grid_label, wavelength_grid_unit + use musica_ccpp_tuvx_temperature, & + only: create_temperature_profile, temperature_label, temperature_unit integer, intent(in) :: vertical_layer_dimension ! (count) integer, intent(in) :: vertical_interface_dimension ! (count) @@ -56,102 +56,62 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & call grids%add( height_grid, error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() + call tuvx_deallocate( grids, null(), null(), null(), height_grid, null(), null() ) return end if wavelength_grid => create_wavelength_grid( wavelength_grid_interfaces, & errmsg, errcode ) if (errcode /= 0) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() + call tuvx_deallocate( grids, null(), null(), null(), height_grid, null(), null() ) return endif call grids%add( wavelength_grid, error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() + call tuvx_deallocate( grids, null(), null(), null(), height_grid, wavelength_grid, & + null() ) return end if profiles => profile_map_t( error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() + call tuvx_deallocate( grids, null(), null(), null(), height_grid, wavelength_grid, & + null() ) return end if temperature_profile => create_temperature_profile( height_grid, errmsg, errcode ) if (errcode /= 0) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() - deallocate( profiles ) + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, wavelength_grid, & + null() ) return endif call profiles%add( temperature_profile, error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() - deallocate( profiles ) - deallocate( temperature_profile ) - temperature_profile => null() + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, wavelength_grid, & + temperature_profile ) return end if radiators => radiator_map_t( error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() - deallocate( profiles ) - deallocate( temperature_profile ) - temperature_profile => null() + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, wavelength_grid, & + temperature_profile ) return end if tuvx => tuvx_t( filename_of_tuvx_configuration, grids, profiles, & radiators, error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() - deallocate( profiles ) - deallocate( temperature_profile ) - temperature_profile => null() - deallocate( radiators ) + call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, wavelength_grid, & + temperature_profile ) return end if - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() - deallocate( profiles ) - deallocate( temperature_profile ) - temperature_profile => null() - deallocate( radiators ) + call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, wavelength_grid, & + temperature_profile ) grids => tuvx%get_grids( error ) if (has_error_occurred( error, errmsg, errcode )) then @@ -162,49 +122,32 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & height_grid => grids%get( height_grid_label, height_grid_unit, error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( tuvx ) - tuvx => null() - deallocate( grids ) + call tuvx_deallocate( grids, null(), null(), tuvx, null(), null(), null() ) return end if wavelength_grid => grids%get( wavelength_grid_label, wavelength_grid_unit, & error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( tuvx ) - tuvx => null() - deallocate( grids ) - deallocate( height_grid ) - height_grid => null() + call tuvx_deallocate( grids, null(), null(), tuvx, height_grid, null(), null() ) return end if - deallocate( grids ) - profiles => tuvx%get_profiles( error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( tuvx ) - tuvx => null() - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() + call tuvx_deallocate( grids, null(), null(), tuvx, height_grid, wavelength_grid, & + null() ) return end if temperature_profile => profiles%get( temperature_label, temperature_unit, error ) if (has_error_occurred( error, errmsg, errcode )) then - deallocate( tuvx ) - tuvx => null() - deallocate( height_grid ) - height_grid => null() - deallocate( wavelength_grid ) - wavelength_grid => null() - deallocate( profiles ) + call tuvx_deallocate( grids, profiles, null(), tuvx, height_grid, wavelength_grid, & + null() ) return end if - deallocate( profiles ) + call tuvx_deallocate( grids, profiles, null(), null(), null(), null(), null() ) end subroutine tuvx_init diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 new file mode 100644 index 00000000..aa524001 --- /dev/null +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 @@ -0,0 +1,49 @@ +module musica_ccpp_tuvx_util + implicit none + + private + public :: tuvx_deallocate + +contains + + !> This is a helper subroutine created to deallocate objects associated with TUV-x + subroutine tuvx_deallocate(grids, profiles, radiators, tuvx, height_grid, & + wavelength_grid, temperature_profile) + use musica_tuvx, only: tuvx_t, grid_map_t, profile_map_t, radiator_map_t, & + grid_t, profile_t + + type(grid_map_t), pointer :: grids + type(profile_map_t), pointer :: profiles + type(radiator_map_t), pointer :: radiators + type(tuvx_t), pointer :: tuvx + type(grid_t), pointer :: height_grid + type(grid_t), pointer :: wavelength_grid + type(profile_t), pointer :: temperature_profile + + if (associated( grids )) deallocate( grids ) + if (associated( profiles )) deallocate( profiles ) + if (associated( radiators )) deallocate( radiators ) + + if (associated( tuvx )) then + deallocate( tuvx ) + tuvx => null() + end if + + if (associated( height_grid )) then + deallocate( height_grid ) + height_grid => null() + end if + + if (associated( wavelength_grid )) then + deallocate( wavelength_grid ) + wavelength_grid => null() + end if + + if (associated( temperature_profile )) then + deallocate( temperature_profile ) + temperature_profile => null() + end if + + end subroutine tuvx_deallocate + +end module musica_ccpp_tuvx_util \ No newline at end of file diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 index 7fd7c1ad..3c6db02c 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 @@ -16,9 +16,9 @@ module musica_ccpp_tuvx_wavelength_grid ! The wavelength grid is defined by the host model. Any wavelength- ! resolved quantities passed to TUV-x must be on this grid. - !> Label for height grid in TUV-x + !> Label for wavelength grid in TUV-x character(len=*), parameter, public :: wavelength_grid_label = "wavelength" - !> Units for height grid in TUV-x + !> Unit for wavelength grid in TUV-x character(len=*), parameter, public :: wavelength_grid_unit = "nm" contains @@ -29,7 +29,6 @@ function create_wavelength_grid( wavelength_grid_interfaces, errmsg, errcode ) & use ccpp_kinds, only: kind_phys use musica_ccpp_util, only: has_error_occurred - use musica_config, only: config_t use musica_tuvx_grid, only: grid_t use musica_util, only: error_t @@ -41,7 +40,7 @@ function create_wavelength_grid( wavelength_grid_interfaces, errmsg, errcode ) & ! local variables real(kind_phys) :: interfaces( size( wavelength_grid_interfaces ) ) ! [nm] reaL(kind_phys) :: midpoints( size( wavelength_grid_interfaces ) - 1 ) ! [nm] - type(error_t) :: error + type(error_t) :: error interfaces(:) = wavelength_grid_interfaces(:) * 1.0e9 midpoints(:) = & diff --git a/test/musica/tuvx/test_tuvx_surface_albedo.F90 b/test/musica/tuvx/test_tuvx_surface_albedo.F90 index d3ae0997..24712fbd 100644 --- a/test/musica/tuvx/test_tuvx_surface_albedo.F90 +++ b/test/musica/tuvx/test_tuvx_surface_albedo.F90 @@ -32,6 +32,9 @@ subroutine test_update_surface_albedo() integer :: i wavelength_grid => create_wavelength_grid(wavelength_grid_interfaces, errmsg, errcode) + ASSERT(errcode == 0) + ASSERT(associated(wavelength_grid)) + profile => create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) ASSERT(errcode == 0) ASSERT(associated(profile)) From a674929129d0d2ea0409ad7a20daf30febb62323 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Thu, 24 Oct 2024 18:54:06 -0600 Subject: [PATCH 07/13] add update function for surface albedo prior to calculating rate constants --- schemes/musica/musica_ccpp.F90 | 8 +- schemes/musica/musica_ccpp.meta | 6 ++ schemes/musica/tuvx/musica_ccpp_tuvx.F90 | 93 +++++++++++++------ .../tuvx/musica_ccpp_tuvx_surface_albedo.F90 | 13 +-- schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 | 8 +- .../tuvx/musica_ccpp_tuvx_wavelength_grid.F90 | 4 +- test/musica/test_musica_api.F90 | 11 ++- 7 files changed, 101 insertions(+), 42 deletions(-) diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90 index 4be069b1..11e589dc 100644 --- a/schemes/musica/musica_ccpp.F90 +++ b/schemes/musica/musica_ccpp.F90 @@ -48,7 +48,8 @@ end subroutine musica_ccpp_init subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, constituent_props, & constituents, geopotential_height_wrt_surface_at_midpoint, & geopotential_height_wrt_surface_at_interface, surface_temperature, & - surface_geopotential, standard_gravitational_acceleration, errmsg, errcode) + surface_geopotential, surface_albedo, & + standard_gravitational_acceleration, errmsg, errcode) use musica_ccpp_micm_util, only: reshape_into_micm_arr, reshape_into_ccpp_arr use musica_ccpp_micm_util, only: convert_to_mol_per_cubic_meter, convert_to_mass_mixing_ratio use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t @@ -66,6 +67,7 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co real(kind_phys), intent(in) :: geopotential_height_wrt_surface_at_interface(:,:) ! m (column, interface) real(kind_phys), intent(in) :: surface_temperature(:) ! K real(kind_phys), intent(in) :: surface_geopotential(:) ! m2 s-2 + real(kind_phys), intent(in) :: surface_albedo ! unitless real(kind_phys), intent(in) :: standard_gravitational_acceleration ! m s-2 character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode @@ -91,8 +93,8 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co call tuvx_run(temperature, dry_air_density, & geopotential_height_wrt_surface_at_midpoint, & geopotential_height_wrt_surface_at_interface, & - surface_temperature, & - surface_geopotential, & + surface_temperature, surface_geopotential, & + surface_albedo, & standard_gravitational_acceleration, & photolysis_rate_constants, & errmsg, errcode) diff --git a/schemes/musica/musica_ccpp.meta b/schemes/musica/musica_ccpp.meta index c719720a..30d19bad 100644 --- a/schemes/musica/musica_ccpp.meta +++ b/schemes/musica/musica_ccpp.meta @@ -101,6 +101,12 @@ units = m2 s-2 dimensions = (horizontal_loop_extent) intent = in +[ surface_albedo ] + standard_name = surface_albedo_due_to_UV_and_VIS_direct + type = real | kind = kind_phys + units = None + dimensions = () + intent = in [ standard_gravitational_acceleration ] standard_name = standard_gravitational_acceleration units = m s-2 diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index cbf0fc87..d72fe4ce 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -16,6 +16,7 @@ module musica_ccpp_tuvx type(grid_t), pointer :: height_grid => null() type(grid_t), pointer :: wavelength_grid => null() type(profile_t), pointer :: temperature_profile => null() + type(profile_t), pointer :: surface_albedo_profile => null() contains @@ -31,6 +32,8 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & only: create_wavelength_grid, wavelength_grid_label, wavelength_grid_unit use musica_ccpp_tuvx_temperature, & only: create_temperature_profile, temperature_label, temperature_unit + use musica_ccpp_tuvx_surface_albedo, & + only: create_surface_albedo_profile, surface_albedo_label, surface_albedo_unit integer, intent(in) :: vertical_layer_dimension ! (count) integer, intent(in) :: vertical_interface_dimension ! (count) @@ -56,62 +59,79 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & call grids%add( height_grid, error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, null(), null(), null(), height_grid, null(), null() ) + call tuvx_deallocate( grids, null(), null(), null(), height_grid, null(), & + null(), null() ) return end if wavelength_grid => create_wavelength_grid( wavelength_grid_interfaces, & errmsg, errcode ) if (errcode /= 0) then - call tuvx_deallocate( grids, null(), null(), null(), height_grid, null(), null() ) + call tuvx_deallocate( grids, null(), null(), null(), height_grid, null(), & + null(), null() ) return endif call grids%add( wavelength_grid, error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, null(), null(), null(), height_grid, wavelength_grid, & - null() ) + call tuvx_deallocate( grids, null(), null(), null(), height_grid, & + wavelength_grid, null(), null() ) return end if profiles => profile_map_t( error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, null(), null(), null(), height_grid, wavelength_grid, & - null() ) + call tuvx_deallocate( grids, null(), null(), null(), height_grid, & + wavelength_grid, null(), null() ) return end if temperature_profile => create_temperature_profile( height_grid, errmsg, errcode ) if (errcode /= 0) then - call tuvx_deallocate( grids, profiles, null(), null(), height_grid, wavelength_grid, & - null() ) + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, & + wavelength_grid, null(), null() ) return endif call profiles%add( temperature_profile, error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, profiles, null(), null(), height_grid, wavelength_grid, & - temperature_profile ) + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, & + wavelength_grid, temperature_profile, null() ) + return + end if + + surface_albedo_profile => create_surface_albedo_profile( wavelength_grid, & + errmsg, errcode ) + if (errcode /= 0) then + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, & + wavelength_grid, temperature_profile, null() ) + return + endif + + call profiles%add( surface_albedo_profile, error ) + if (has_error_occurred( error, errmsg, errcode )) then + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, & + wavelength_grid, temperature_profile, surface_albedo_profile ) return end if radiators => radiator_map_t( error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, profiles, null(), null(), height_grid, wavelength_grid, & - temperature_profile ) + call tuvx_deallocate( grids, profiles, null(), null(), height_grid, & + wavelength_grid, temperature_profile, surface_albedo_profile ) return end if tuvx => tuvx_t( filename_of_tuvx_configuration, grids, profiles, & radiators, error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, wavelength_grid, & - temperature_profile ) + call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, & + wavelength_grid, temperature_profile, surface_albedo_profile ) return end if - call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, wavelength_grid, & - temperature_profile ) + call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, & + wavelength_grid, temperature_profile, surface_albedo_profile ) grids => tuvx%get_grids( error ) if (has_error_occurred( error, errmsg, errcode )) then @@ -122,32 +142,42 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & height_grid => grids%get( height_grid_label, height_grid_unit, error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, null(), null(), tuvx, null(), null(), null() ) + call tuvx_deallocate( grids, null(), null(), tuvx, null(), null(), & + null(), null() ) return end if wavelength_grid => grids%get( wavelength_grid_label, wavelength_grid_unit, & error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, null(), null(), tuvx, height_grid, null(), null() ) + call tuvx_deallocate( grids, null(), null(), tuvx, height_grid, null(), & + null(), null() ) return end if profiles => tuvx%get_profiles( error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, null(), null(), tuvx, height_grid, wavelength_grid, & - null() ) + call tuvx_deallocate( grids, null(), null(), tuvx, height_grid, & + wavelength_grid, null(), null() ) return end if temperature_profile => profiles%get( temperature_label, temperature_unit, error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, profiles, null(), tuvx, height_grid, wavelength_grid, & - null() ) + call tuvx_deallocate( grids, profiles, null(), tuvx, height_grid, & + wavelength_grid, null(), null() ) return end if - call tuvx_deallocate( grids, profiles, null(), null(), null(), null(), null() ) + surface_albedo_profile => profiles%get( surface_albedo_label, surface_albedo_unit, error ) + if (has_error_occurred( error, errmsg, errcode )) then + call tuvx_deallocate( grids, profiles, null(), tuvx, height_grid, & + wavelength_grid, temperature_profile, null() ) + return + end if + + call tuvx_deallocate( grids, profiles, null(), null(), null(), null(), & + null(), null() ) end subroutine tuvx_init @@ -160,11 +190,13 @@ subroutine tuvx_run(temperature, dry_air_density, & geopotential_height_wrt_surface_at_midpoint, & geopotential_height_wrt_surface_at_interface, & surface_temperature, surface_geopotential, & + surface_albedo, & standard_gravitational_acceleration, & photolysis_rate_constants, errmsg, errcode) - use musica_util, only: error_t - use musica_ccpp_tuvx_height_grid, only: set_height_grid_values, calculate_heights - use musica_ccpp_tuvx_temperature, only: set_temperature_values + use musica_util, only: error_t + use musica_ccpp_tuvx_height_grid, only: set_height_grid_values, calculate_heights + use musica_ccpp_tuvx_temperature, only: set_temperature_values + use musica_ccpp_tuvx_surface_albedo, only: set_surface_albedo_values real(kind_phys), intent(in) :: temperature(:,:) ! K (column, layer) real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 (column, layer) @@ -172,6 +204,7 @@ subroutine tuvx_run(temperature, dry_air_density, & real(kind_phys), intent(in) :: geopotential_height_wrt_surface_at_interface(:,:) ! m (column, interface) real(kind_phys), intent(in) :: surface_temperature(:) ! K real(kind_phys), intent(in) :: surface_geopotential(:) ! m2 s-2 + real(kind_phys), intent(in) :: surface_albedo ! unitless real(kind_phys), intent(in) :: standard_gravitational_acceleration ! m s-2 ! temporarily set to Chapman mechanism and 1 dimension ! until mapping between MICM and TUV-x is implemented @@ -202,6 +235,9 @@ subroutine tuvx_run(temperature, dry_air_density, & if (errcode /= 0) return end do + call set_surface_albedo_values( surface_albedo_profile, surface_albedo, errmsg, errcode ) + if (errcode /= 0) return + ! stand-in until actual photolysis rate constants are calculated photolysis_rate_constants(:) = 1.0e-6_kind_phys @@ -230,6 +266,11 @@ subroutine tuvx_final(errmsg, errcode) temperature_profile => null() end if + if (associated( surface_albedo_profile )) then + deallocate( surface_albedo_profile ) + surface_albedo_profile => null() + end if + if (associated( tuvx )) then deallocate( tuvx ) tuvx => null() diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 index b6045162..0dc913f3 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 @@ -2,7 +2,8 @@ module musica_ccpp_tuvx_surface_albedo implicit none private - public :: create_surface_albedo_profile, set_surface_albedo_values + public :: create_surface_albedo_profile, set_surface_albedo_values, & + surface_albedo_label, surface_albedo_unit !> Label for surface albedo in TUV-x character(len=*), parameter :: surface_albedo_label = "surface_albedo" @@ -18,7 +19,6 @@ module musica_ccpp_tuvx_surface_albedo !> Creates a TUV-x surface albedo profile from the host-model wavelength grid function create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) & result( profile ) - use musica_ccpp_util, only: has_error_occurred use musica_tuvx_grid, only: grid_t use musica_tuvx_profile, only: profile_t @@ -43,10 +43,11 @@ function create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) & end function create_surface_albedo_profile - !> Sets TUV-x surface albedo value + !> Sets TUV-x surface albedo values + !! + !! CAM uses a single value for surface albedo at all wavelengths subroutine set_surface_albedo_values( profile, host_surface_albedo, & - errmsg, errcode ) - + errmsg, errcode ) use ccpp_kinds, only: kind_phys use musica_ccpp_util, only: has_error_occurred use musica_tuvx_profile, only: profile_t @@ -62,7 +63,7 @@ subroutine set_surface_albedo_values( profile, host_surface_albedo, & real(kind_phys) :: surface_albedos(num_wavelength_bins + 1) if (size(surface_albedos) <= DEFAULT_NUM_WAVELENGTH_BINS + 1) then - errmsg = "[MUSICA Error] Invalid dimension for number of wavelngth bins" + errmsg = "[MUSICA Error] Invalid size of TUV-x interface wavelengths. " errcode = 1 return end if diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 index aa524001..25ac1a28 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 @@ -8,7 +8,7 @@ module musica_ccpp_tuvx_util !> This is a helper subroutine created to deallocate objects associated with TUV-x subroutine tuvx_deallocate(grids, profiles, radiators, tuvx, height_grid, & - wavelength_grid, temperature_profile) + wavelength_grid, temperature_profile, surface_albedo_profile) use musica_tuvx, only: tuvx_t, grid_map_t, profile_map_t, radiator_map_t, & grid_t, profile_t @@ -19,6 +19,7 @@ subroutine tuvx_deallocate(grids, profiles, radiators, tuvx, height_grid, & type(grid_t), pointer :: height_grid type(grid_t), pointer :: wavelength_grid type(profile_t), pointer :: temperature_profile + type(profile_t), pointer :: surface_albedo_profile if (associated( grids )) deallocate( grids ) if (associated( profiles )) deallocate( profiles ) @@ -44,6 +45,11 @@ subroutine tuvx_deallocate(grids, profiles, radiators, tuvx, height_grid, & temperature_profile => null() end if + if (associated( surface_albedo_profile )) then + deallocate( surface_albedo_profile ) + surface_albedo_profile => null() + end if + end subroutine tuvx_deallocate end module musica_ccpp_tuvx_util \ No newline at end of file diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 index 3c6db02c..96528d5d 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 @@ -38,8 +38,8 @@ function create_wavelength_grid( wavelength_grid_interfaces, errmsg, errcode ) & type(grid_t), pointer :: wavelength_grid ! local variables - real(kind_phys) :: interfaces( size( wavelength_grid_interfaces ) ) ! [nm] - reaL(kind_phys) :: midpoints( size( wavelength_grid_interfaces ) - 1 ) ! [nm] + real(kind_phys) :: interfaces( size( wavelength_grid_interfaces ) ) ! nm + reaL(kind_phys) :: midpoints( size( wavelength_grid_interfaces ) - 1 ) ! nm type(error_t) :: error interfaces(:) = wavelength_grid_interfaces(:) * 1.0e9 diff --git a/test/musica/test_musica_api.F90 b/test/musica/test_musica_api.F90 index 30c0ee74..2aa6417c 100644 --- a/test/musica/test_musica_api.F90 +++ b/test/musica/test_musica_api.F90 @@ -31,6 +31,7 @@ subroutine test_musica_ccpp_api() real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS+1) :: geopotential_height_wrt_surface_at_interface ! m real(kind_phys), dimension(NUM_COLUMNS) :: surface_temperature ! K real(kind_phys), dimension(NUM_COLUMNS) :: surface_geopotential ! m2 s-2 + real(kind_phys) :: surface_albedo ! unitless real(kind_phys) :: standard_gravitational_acceleration ! s2 m-1 real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: temperature ! K real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: pressure ! Pa @@ -166,6 +167,7 @@ subroutine test_musica_ccpp_api() geopotential_height_wrt_surface_at_interface(2,:) = (/ 3000.0_kind_phys, 500.0_kind_phys, -1500.0_kind_phys /) surface_temperature = (/ 300.0_kind_phys, 300.0_kind_phys /) surface_geopotential = (/ 100.0_kind_phys, 200.0_kind_phys /) + surface_albedo = 20.0_kind_phys standard_gravitational_acceleration = 10.0_kind_phys temperature(:,1) = (/ 100._kind_phys, 200._kind_phys /) temperature(:,2) = (/ 300._kind_phys, 400._kind_phys /) @@ -226,10 +228,11 @@ subroutine test_musica_ccpp_api() write(*,*) "[MUSICA INFO] Initial Concentrations" write(*,fmt="(4(3x,e13.6))") constituents - call musica_ccpp_run(time_step, temperature, pressure, dry_air_density, constituent_props_ptr, & - constituents, geopotential_height_wrt_surface_at_midpoint, & - geopotential_height_wrt_surface_at_interface, surface_temperature, & - surface_geopotential, standard_gravitational_acceleration, errmsg, errcode) + call musica_ccpp_run(time_step, temperature, pressure, dry_air_density, constituent_props_ptr, & + constituents, geopotential_height_wrt_surface_at_midpoint, & + geopotential_height_wrt_surface_at_interface, surface_temperature, & + surface_geopotential, surface_albedo, standard_gravitational_acceleration, & + errmsg, errcode) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 From 8f00fc5bb98a76f90f9c700e4d7853e78f535313 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Fri, 25 Oct 2024 12:31:06 -0600 Subject: [PATCH 08/13] refactor the test code for minor improvements --- .../tuvx/musica_ccpp_tuvx_surface_albedo.F90 | 2 +- test/musica/tuvx/test_tuvx_height_grid.F90 | 18 ++++++------ test/musica/tuvx/test_tuvx_surface_albedo.F90 | 24 ++++++++-------- test/musica/tuvx/test_tuvx_temperature.F90 | 28 +++++++++---------- 4 files changed, 36 insertions(+), 36 deletions(-) diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 index 0dc913f3..2d7dac71 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 @@ -63,7 +63,7 @@ subroutine set_surface_albedo_values( profile, host_surface_albedo, & real(kind_phys) :: surface_albedos(num_wavelength_bins + 1) if (size(surface_albedos) <= DEFAULT_NUM_WAVELENGTH_BINS + 1) then - errmsg = "[MUSICA Error] Invalid size of TUV-x interface wavelengths. " + errmsg = "[MUSICA Error] Invalid size of TUV-x wavelength interfaces." errcode = 1 return end if diff --git a/test/musica/tuvx/test_tuvx_height_grid.F90 b/test/musica/tuvx/test_tuvx_height_grid.F90 index c0f5c6ed..012de868 100644 --- a/test/musica/tuvx/test_tuvx_height_grid.F90 +++ b/test/musica/tuvx/test_tuvx_height_grid.F90 @@ -17,15 +17,15 @@ subroutine test_create_height_grid() use musica_tuvx_grid, only: grid_t use ccpp_kinds, only: kind_phys - integer, parameter :: NUM_HOST_MIDPOINTS = 2 - integer, parameter :: NUM_HOST_INTERFACES = 3 - real(kind_phys), target :: host_midpoints(NUM_HOST_MIDPOINTS) = [200.8_kind_phys, 100.6_kind_phys] - real(kind_phys), target :: host_interfaces(NUM_HOST_INTERFACES) = [250.3_kind_phys, 150.2_kind_phys, 50.1_kind_phys] - type(grid_t), pointer :: height_grid - character(len=512) :: errmsg - integer :: errcode - real(kind_phys) :: abs_error = 1e-5 - integer :: i + integer, parameter :: NUM_HOST_MIDPOINTS = 2 + integer, parameter :: NUM_HOST_INTERFACES = 3 + real(kind_phys) :: host_midpoints(NUM_HOST_MIDPOINTS) = [200.8_kind_phys, 100.6_kind_phys] + real(kind_phys) :: host_interfaces(NUM_HOST_INTERFACES) = [250.3_kind_phys, 150.2_kind_phys, 50.1_kind_phys] + type(grid_t), pointer :: height_grid + character(len=512) :: errmsg + integer :: errcode + real(kind_phys) :: abs_error = 1e-5 + integer :: i ! local variables real(kind_phys), dimension(NUM_HOST_MIDPOINTS+1) :: midpoints diff --git a/test/musica/tuvx/test_tuvx_surface_albedo.F90 b/test/musica/tuvx/test_tuvx_surface_albedo.F90 index 24712fbd..2f4b1385 100644 --- a/test/musica/tuvx/test_tuvx_surface_albedo.F90 +++ b/test/musica/tuvx/test_tuvx_surface_albedo.F90 @@ -18,13 +18,14 @@ subroutine test_update_surface_albedo() use musica_tuvx_profile, only: profile_t use ccpp_kinds, only: kind_phys - integer, parameter :: NUM_WAVELENGTH_BIN = 4 + integer, parameter :: NUM_WAVELENGTH_BINS = 4 type(grid_t), pointer :: wavelength_grid type(profile_t), pointer :: profile - real(kind_phys), target :: host_surface_albedo = 500.5_kind_phys - real(kind_phys) :: wavelength_grid_interfaces(NUM_WAVELENGTH_BIN + 1) = & + real(kind_phys) :: wavelength_grid_interfaces(NUM_WAVELENGTH_BINS + 1) = & [200.0e-9_kind_phys, 210.0e-9_kind_phys, 240.0e-9_kind_phys, 300.0e-9_kind_phys, 400.0e-9_kind_phys] - real(kind_phys) :: surface_albedos(NUM_WAVELENGTH_BIN + 1) + real(kind_phys) :: host_surface_albedo = 500.5_kind_phys + real(kind_phys) :: surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) + real(kind_phys) :: expected_surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) = host_surface_albedo character(len=512) :: errmsg integer :: errcode type(error_t) :: error @@ -42,17 +43,16 @@ subroutine test_update_surface_albedo() call set_surface_albedo_values( profile, host_surface_albedo, errmsg, errcode ) ASSERT(errcode == 0) - call profile%get_edge_values( surface_albedos, error) + call profile%get_edge_values( surface_albedo_interfaces, error) ASSERT(error%is_success()) - ASSERT_NEAR(surface_albedos(1), 500.5, abs_error) - ASSERT_NEAR(surface_albedos(2), 500.5, abs_error) - ASSERT_NEAR(surface_albedos(3), 500.5, abs_error) - ASSERT_NEAR(surface_albedos(4), 500.5, abs_error) - ASSERT_NEAR(surface_albedos(5), 500.5, abs_error) + ASSERT_NEAR(surface_albedo_interfaces(1), expected_surface_albedo_interfaces(1), abs_error) + ASSERT_NEAR(surface_albedo_interfaces(2), expected_surface_albedo_interfaces(2), abs_error) + ASSERT_NEAR(surface_albedo_interfaces(3), expected_surface_albedo_interfaces(3), abs_error) + ASSERT_NEAR(surface_albedo_interfaces(4), expected_surface_albedo_interfaces(4), abs_error) + ASSERT_NEAR(surface_albedo_interfaces(5), expected_surface_albedo_interfaces(5), abs_error) deallocate( profile ) deallocate( wavelength_grid ) - end subroutine test_update_surface_albedo - + end subroutine test_invalid_size_wavelength_interfaces end program test_tuvx_surface_albedo \ No newline at end of file diff --git a/test/musica/tuvx/test_tuvx_temperature.F90 b/test/musica/tuvx/test_tuvx_temperature.F90 index 2f78eb56..38a22744 100644 --- a/test/musica/tuvx/test_tuvx_temperature.F90 +++ b/test/musica/tuvx/test_tuvx_temperature.F90 @@ -18,20 +18,20 @@ subroutine test_update_temperature() use musica_tuvx_profile, only: profile_t use ccpp_kinds, only: kind_phys - integer, parameter :: NUM_HOST_MIDPOINTS = 5 - integer, parameter :: NUM_HOST_INTERFACES = 6 - real(kind_phys) :: host_midpoint_temperature(NUM_HOST_MIDPOINTS) & - = [800.8_kind_phys, 700.7_kind_phys, 600.6_kind_phys, 500.5_kind_phys, 400.4_kind_phys] - real(kind_phys) :: host_surface_temperature = 300.3_kind_phys - real(kind_phys) :: interface_temperatures(NUM_HOST_MIDPOINTS+2) - real(kind_phys) :: expected_interface_temperatures(NUM_HOST_MIDPOINTS+2) & - = [300.3_kind_phys, 400.4_kind_phys, 500.5_kind_phys, 600.6_kind_phys, 700.7_kind_phys, 800.8_kind_phys, 800.8_kind_phys] - type(grid_t), pointer :: height_grid - type(profile_t), pointer :: profile - character(len=512) :: errmsg - integer :: errcode - type(error_t) :: error - real(kind_phys) :: abs_error = 1e-4 + integer, parameter :: NUM_HOST_MIDPOINTS = 5 + integer, parameter :: NUM_HOST_INTERFACES = 6 + real(kind_phys) :: host_midpoint_temperature(NUM_HOST_MIDPOINTS) = & + [800.8_kind_phys, 700.7_kind_phys, 600.6_kind_phys, 500.5_kind_phys, 400.4_kind_phys] + real(kind_phys) :: host_surface_temperature = 300.3_kind_phys + real(kind_phys) :: interface_temperatures(NUM_HOST_MIDPOINTS+2) + real(kind_phys) :: expected_interface_temperatures(NUM_HOST_MIDPOINTS+2) = & + [300.3_kind_phys, 400.4_kind_phys, 500.5_kind_phys, 600.6_kind_phys, 700.7_kind_phys, 800.8_kind_phys, 800.8_kind_phys] + type(grid_t), pointer :: height_grid + type(profile_t), pointer :: profile + character(len=512) :: errmsg + integer :: errcode + type(error_t) :: error + real(kind_phys) :: abs_error = 1e-4 height_grid => create_height_grid(NUM_HOST_MIDPOINTS, NUM_HOST_INTERFACES, errmsg, errcode) ASSERT(errcode == 0) From 851e9a8cd70fb8b156c786be2c781c838de99cad Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Fri, 25 Oct 2024 14:47:43 -0600 Subject: [PATCH 09/13] refactor the test code for minor improvements --- schemes/musica/musica_ccpp.F90 | 4 +- test/musica/micm/test_micm_util.F90 | 50 +++++++------- test/musica/test_musica_api.F90 | 27 ++++---- test/musica/tuvx/test_tuvx_height_grid.F90 | 66 +++++++++---------- test/musica/tuvx/test_tuvx_surface_albedo.F90 | 21 +++--- test/musica/tuvx/test_tuvx_temperature.F90 | 20 +++--- .../musica/tuvx/test_tuvx_wavelength_grid.F90 | 30 ++++----- 7 files changed, 102 insertions(+), 116 deletions(-) diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90 index fe0f226f..b7deb9a8 100644 --- a/schemes/musica/musica_ccpp.F90 +++ b/schemes/musica/musica_ccpp.F90 @@ -47,8 +47,8 @@ end subroutine musica_ccpp_init !! \htmlinclude musica_ccpp_run.html !! !! The standard name for the variable 'surface_tempearture' is - !! `blackbody_temperature_at_surface` because this is what we have as - !! the standard name for `cam_in%ts`, whcih represents the same quantity. + !! 'blackbody_temperature_at_surface' because this is what we have as + !! the standard name for 'cam_in%ts', whcih represents the same quantity. subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, constituent_props, & constituents, geopotential_height_wrt_surface_at_midpoint, & geopotential_height_wrt_surface_at_interface, surface_temperature, & diff --git a/test/musica/micm/test_micm_util.F90 b/test/musica/micm/test_micm_util.F90 index 5581d5d8..798f1408 100644 --- a/test/musica/micm/test_micm_util.F90 +++ b/test/musica/micm/test_micm_util.F90 @@ -19,7 +19,11 @@ subroutine test_reshape() integer, parameter :: NUM_SPECIES = 4 integer, parameter :: NUM_COLUMNS = 2 integer, parameter :: NUM_LAYERS = 2 - integer, parameter :: NUM_GRID_CELLS = 4 + integer, parameter :: NUM_GRID_CELLS = NUM_COLUMNS * NUM_LAYERS + real, parameter :: ABS_ERROR = 1e-7 + real(c_double) :: expected_conditions(NUM_GRID_CELLS) = [100.0, 200.0, 300.0, 400.0] + real(c_double) :: expected_constituents(NUM_GRID_CELLS*NUM_SPECIES) = & + [0.1, 0.2, 0.3, 0.4, 0.21, 0.22, 0.23, 0.24, 0.41, 0.42, 0.43, 0.44, 0.31, 0.32, 0.33, 0.34] real(kind_phys) :: temperature(NUM_COLUMNS,NUM_LAYERS) real(kind_phys) :: pressure(NUM_COLUMNS,NUM_LAYERS) real(kind_phys) :: dry_air_density(NUM_COLUMNS,NUM_LAYERS) @@ -28,12 +32,7 @@ subroutine test_reshape() real(c_double) :: micm_pressure(NUM_GRID_CELLS) real(c_double) :: micm_dry_air_density(NUM_GRID_CELLS) real(c_double) :: micm_constituents(NUM_GRID_CELLS*NUM_SPECIES) - - ! local variables - real(c_double), dimension(NUM_GRID_CELLS) :: arr_conditions - real(c_double), dimension(NUM_GRID_CELLS*NUM_SPECIES) :: arr_constituents - integer :: i_column, i_layer, i_elem, i_arr - real :: abs_error = 1e-7 + integer :: i_column, i_layer, i_elem, i_arr temperature(:,1) = (/ 100._kind_phys, 200._kind_phys /) temperature(:,2) = (/ 300._kind_phys, 400._kind_phys /) @@ -45,20 +44,18 @@ subroutine test_reshape() constituents(1,2,:) = (/ 0.41_kind_phys, 0.42_kind_phys, 0.43_kind_phys, 0.44_kind_phys /) constituents(2,1,:) = (/ 0.21_kind_phys, 0.22_kind_phys, 0.23_kind_phys, 0.24_kind_phys /) constituents(2,2,:) = (/ 0.31_kind_phys, 0.32_kind_phys, 0.33_kind_phys, 0.34_kind_phys /) - arr_conditions = (/ 100.0, 200.0, 300.0, 400.0 /) - arr_constituents = (/ 0.1, 0.2, 0.3, 0.4, 0.21, 0.22, 0.23, 0.24, 0.41, 0.42, 0.43, 0.44, 0.31, 0.32, 0.33, 0.34 /) call reshape_into_micm_arr(temperature, pressure, dry_air_density, constituents, & micm_temperature, micm_pressure, micm_dry_air_density, micm_constituents) do i_elem = 1, NUM_GRID_CELLS - ASSERT(micm_temperature(i_elem) == arr_conditions(i_elem)) - ASSERT(micm_pressure(i_elem) == arr_conditions(i_elem)) - ASSERT(micm_dry_air_density(i_elem) == arr_conditions(i_elem)) + ASSERT(micm_temperature(i_elem) == expected_conditions(i_elem)) + ASSERT(micm_pressure(i_elem) == expected_conditions(i_elem)) + ASSERT(micm_dry_air_density(i_elem) == expected_conditions(i_elem)) end do do i_elem = 1, size(micm_constituents) - ASSERT_NEAR(micm_constituents(i_elem), arr_constituents(i_elem), abs_error) + ASSERT_NEAR(micm_constituents(i_elem), expected_constituents(i_elem), ABS_ERROR) end do call reshape_into_ccpp_arr(micm_constituents, constituents) @@ -67,7 +64,7 @@ subroutine test_reshape() do i_layer = 1, NUM_LAYERS do i_column = 1, NUM_COLUMNS do i_elem = 1, NUM_SPECIES - ASSERT_NEAR(constituents(i_column, i_layer, i_elem), arr_constituents(i_arr), abs_error) + ASSERT_NEAR(constituents(i_column, i_layer, i_elem), expected_constituents(i_arr), ABS_ERROR) i_arr = i_arr + 1 end do end do @@ -78,16 +75,16 @@ end subroutine test_reshape subroutine test_unit_conversion() use ccpp_kinds, only: kind_phys - integer, parameter :: NUM_COLUMNS = 2 - integer, parameter :: NUM_LAYERS = 2 - integer, parameter :: NUM_SPECIES = 4 - real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: dry_air_density ! kg m-3 - real(kind_phys), dimension(NUM_SPECIES) :: molar_mass_arr - real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: constituents - real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: ccpp_constituents ! kg kg-1 - real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: micm_constituents ! mol m-3 - integer :: i_column, i_layer, i_elem - real :: abs_error = 1e-3 + integer, parameter :: NUM_COLUMNS = 2 + integer, parameter :: NUM_LAYERS = 2 + integer, parameter :: NUM_SPECIES = NUM_COLUMNS * NUM_LAYERS + real, parameter :: ABS_ERROR = 1e-3 + real(kind_phys) :: dry_air_density(NUM_COLUMNS,NUM_LAYERS) + real(kind_phys) :: molar_mass_arr(NUM_SPECIES) + real(kind_phys) :: constituents(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) + real(kind_phys) :: ccpp_constituents(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) ! kg kg-1 + real(kind_phys) :: micm_constituents(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) ! mol m-3 + integer :: i_column, i_layer, i_elem dry_air_density(:,1) = (/ 3.5_kind_phys, 4.5_kind_phys /) dry_air_density(:,2) = (/ 5.5_kind_phys, 6.5_kind_phys /) @@ -112,7 +109,7 @@ subroutine test_unit_conversion() do i_column = 1, NUM_COLUMNS do i_layer = 1, NUM_LAYERS do i_elem = 1, NUM_SPECIES - ASSERT_NEAR(constituents(i_column, i_layer, i_elem), micm_constituents(i_column, i_layer, i_elem), abs_error) + ASSERT_NEAR(constituents(i_column, i_layer, i_elem), micm_constituents(i_column, i_layer, i_elem), ABS_ERROR) end do end do end do @@ -122,10 +119,11 @@ subroutine test_unit_conversion() do i_column = 1, NUM_COLUMNS do i_layer = 1, NUM_LAYERS do i_elem = 1, NUM_SPECIES - ASSERT_NEAR(constituents(i_column, i_layer, i_elem), ccpp_constituents(i_column, i_layer, i_elem), abs_error) + ASSERT_NEAR(constituents(i_column, i_layer, i_elem), ccpp_constituents(i_column, i_layer, i_elem), ABS_ERROR) end do end do end do + end subroutine test_unit_conversion end program test_micm_util \ No newline at end of file diff --git a/test/musica/test_musica_api.F90 b/test/musica/test_musica_api.F90 index 2aa6417c..7386de86 100644 --- a/test/musica/test_musica_api.F90 +++ b/test/musica/test_musica_api.F90 @@ -16,16 +16,15 @@ subroutine test_musica_ccpp_api() use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t - implicit none - integer, parameter :: NUM_SPECIES = 4 integer, parameter :: NUM_COLUMNS = 2 integer, parameter :: NUM_LAYERS = 2 integer, parameter :: NUM_WAVELENGTH_BINS = 102 - integer :: solver_type + integer :: num_grid_cells = NUM_COLUMNS * NUM_LAYERS + integer :: solver_type = Rosenbrock integer :: errcode character(len=512) :: errmsg - real(kind_phys) :: time_step ! s + real(kind_phys) :: time_step = 60._kind_phys ! s real(kind_phys), dimension(NUM_WAVELENGTH_BINS+1) :: photolysis_wavelength_grid_interfaces ! m real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: geopotential_height_wrt_surface_at_midpoint ! m real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS+1) :: geopotential_height_wrt_surface_at_interface ! m @@ -38,19 +37,13 @@ subroutine test_musica_ccpp_api() real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: dry_air_density ! kg m-3 real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: constituents ! kg kg-1 type(ccpp_constituent_prop_ptr_t), allocatable :: constituent_props_ptr(:) + type(ccpp_constituent_properties_t), allocatable, target :: constituent_props(:) + type(ccpp_constituent_properties_t), pointer :: const_prop + real(kind_phys) :: molar_mass + character(len=512) :: species_name, units + logical :: tmp_bool, is_advected + integer :: i - ! local variables - type(ccpp_constituent_properties_t), allocatable, target :: constituent_props(:) - type(ccpp_constituent_properties_t), pointer :: const_prop - real(kind_phys) :: molar_mass - character(len=512) :: species_name, units - logical :: tmp_bool, is_advected - integer :: num_grid_cells - integer :: i - - solver_type = Rosenbrock - num_grid_cells = NUM_COLUMNS * NUM_LAYERS - time_step = 60._kind_phys ! These are the values that will be used in CAM-SIMA and correspond to the wavelength ! bins used in the CAM-Chem photolysis rate constant lookup table. ! @@ -248,6 +241,8 @@ subroutine test_musica_ccpp_api() stop 3 endif + deallocate(constituent_props_ptr) + end subroutine test_musica_ccpp_api end program run_test_musica_ccpp \ No newline at end of file diff --git a/test/musica/tuvx/test_tuvx_height_grid.F90 b/test/musica/tuvx/test_tuvx_height_grid.F90 index 012de868..6da26960 100644 --- a/test/musica/tuvx/test_tuvx_height_grid.F90 +++ b/test/musica/tuvx/test_tuvx_height_grid.F90 @@ -7,6 +7,8 @@ program test_tuvx_height_grid #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 + real, parameter :: ABS_ERROR = 1e-5 + call test_create_height_grid() call test_calculate_height_grid_values() @@ -21,17 +23,16 @@ subroutine test_create_height_grid() integer, parameter :: NUM_HOST_INTERFACES = 3 real(kind_phys) :: host_midpoints(NUM_HOST_MIDPOINTS) = [200.8_kind_phys, 100.6_kind_phys] real(kind_phys) :: host_interfaces(NUM_HOST_INTERFACES) = [250.3_kind_phys, 150.2_kind_phys, 50.1_kind_phys] - type(grid_t), pointer :: height_grid + real(kind_phys) :: expected_midpoints(NUM_HOST_MIDPOINTS+1) = [(100.6 + 50.1) * 0.5, 150.2, (250.3 + 200.8) * 0.5] + real(kind_phys) :: expected_interfaces(NUM_HOST_INTERFACES+1) = [50.1_kind_phys, 100.6_kind_phys, 200.8_kind_phys, 250.3_kind_phys] + real(kind_phys) :: midpoints(NUM_HOST_MIDPOINTS+1) + real(kind_phys) :: interfaces(NUM_HOST_INTERFACES+1) + type(grid_t), pointer :: height_grid => null() + type(error_t) :: error character(len=512) :: errmsg integer :: errcode - real(kind_phys) :: abs_error = 1e-5 integer :: i - ! local variables - real(kind_phys), dimension(NUM_HOST_MIDPOINTS+1) :: midpoints - real(kind_phys), dimension(NUM_HOST_INTERFACES+1) :: interfaces - type(error_t) :: error - height_grid => create_height_grid(-1, 0, errmsg, errcode) ASSERT(errcode == 1) ASSERT(.not. associated(height_grid)) @@ -54,50 +55,47 @@ subroutine test_create_height_grid() call height_grid%get_midpoints(midpoints, error) ASSERT(error%is_success()) + do i = 1, size(midpoints) + ASSERT_NEAR(midpoints(i), expected_midpoints(i), ABS_ERROR) + end do call height_grid%get_edges(interfaces, error) ASSERT(error%is_success()) - - ASSERT_NEAR(midpoints(1), (100.6 + 50.1) * 0.5, abs_error) - ASSERT_NEAR(midpoints(2), 150.2, abs_error) - ASSERT_NEAR(midpoints(3), (250.3 + 200.8) * 0.5, abs_error) - ASSERT_NEAR(interfaces(1), 50.1, abs_error) - ASSERT_NEAR(interfaces(2), 100.6, abs_error) - ASSERT_NEAR(interfaces(3), 200.8, abs_error) - ASSERT_NEAR(interfaces(4), 250.3, abs_error) + do i = 1, size(interfaces) + ASSERT_NEAR(interfaces(i), expected_interfaces(i), ABS_ERROR) + end do deallocate( height_grid ) end subroutine test_create_height_grid subroutine test_calculate_height_grid_values() - use ccpp_kinds, only: kind_phys - integer, parameter :: NUM_LAYERS = 2 - real(kind_phys), dimension(NUM_LAYERS) :: geopotential_height_wrt_surface_at_midpoint ! m - real(kind_phys), dimension(NUM_LAYERS+1) :: geopotential_height_wrt_surface_at_interface ! m - real(kind_phys) :: surface_geopotential ! m2 s-2 - real(kind_phys) :: reciprocal_of_gravitational_acceleration ! s2 m-1 - real(kind_phys), dimension(NUM_LAYERS) :: height_midpoints ! km - real(kind_phys), dimension(NUM_LAYERS+1) :: height_interfaces ! km - - geopotential_height_wrt_surface_at_midpoint(:) = (/ 2000.0_kind_phys, 500.0_kind_phys /) - geopotential_height_wrt_surface_at_interface(:) = (/ 3000.0_kind_phys, 1000.0_kind_phys, 0.0_kind_phys /) - surface_geopotential = 100.0_kind_phys - reciprocal_of_gravitational_acceleration = 10.0_kind_phys + integer, parameter :: NUM_LAYERS = 2 + real(kind_phys) :: geopotential_height_wrt_surface_at_midpoint(NUM_LAYERS) = [2000.0_kind_phys, 500.0_kind_phys] + real(kind_phys) :: geopotential_height_wrt_surface_at_interface(NUM_LAYERS+1) = [3000.0_kind_phys, 1000.0_kind_phys, 0.0_kind_phys] + real(kind_phys) :: surface_geopotential = 100.0_kind_phys + real(kind_phys) :: reciprocal_of_gravitational_acceleration = 10.0_kind_phys + real(kind_phys) :: expected_height_midpoints(NUM_LAYERS) = [3.0_kind_phys, 1.5_kind_phys] + real(kind_phys) :: expected_height_interfaces(NUM_LAYERS+1) = [4.0_kind_phys, 2.0_kind_phys, 1.0_kind_phys] + real(kind_phys) :: height_midpoints(NUM_LAYERS) + real(kind_phys) :: height_interfaces(NUM_LAYERS+1) + integer :: i call calculate_heights(geopotential_height_wrt_surface_at_midpoint, & geopotential_height_wrt_surface_at_interface, & surface_geopotential, reciprocal_of_gravitational_acceleration, & height_midpoints, height_interfaces) - ASSERT_NEAR(height_midpoints(1), 3.0, 1e-5) - ASSERT_NEAR(height_midpoints(2), 1.5, 1e-5) - ASSERT_NEAR(height_interfaces(1), 4.0, 1e-5) - ASSERT_NEAR(height_interfaces(2), 2.0, 1e-5) - ASSERT_NEAR(height_interfaces(3), 1.0, 1e-5) + do i = 1, size(height_midpoints) + ASSERT_NEAR(height_midpoints(i), expected_height_midpoints(i), ABS_ERROR) + end do + + do i = 1, size(height_interfaces) + ASSERT_NEAR(height_interfaces(i), expected_height_interfaces(i), ABS_ERROR) + end do end subroutine test_calculate_height_grid_values -end program test_tuvx_height_grid +end program test_tuvx_height_grid \ No newline at end of file diff --git a/test/musica/tuvx/test_tuvx_surface_albedo.F90 b/test/musica/tuvx/test_tuvx_surface_albedo.F90 index 2f4b1385..0b4ce5dd 100644 --- a/test/musica/tuvx/test_tuvx_surface_albedo.F90 +++ b/test/musica/tuvx/test_tuvx_surface_albedo.F90 @@ -19,17 +19,17 @@ subroutine test_update_surface_albedo() use ccpp_kinds, only: kind_phys integer, parameter :: NUM_WAVELENGTH_BINS = 4 - type(grid_t), pointer :: wavelength_grid - type(profile_t), pointer :: profile + real, parameter :: ABS_ERROR = 1e-4 real(kind_phys) :: wavelength_grid_interfaces(NUM_WAVELENGTH_BINS + 1) = & [200.0e-9_kind_phys, 210.0e-9_kind_phys, 240.0e-9_kind_phys, 300.0e-9_kind_phys, 400.0e-9_kind_phys] real(kind_phys) :: host_surface_albedo = 500.5_kind_phys + real(kind_phys) :: expected_surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) = 500.5_kind_phys real(kind_phys) :: surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) - real(kind_phys) :: expected_surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) = host_surface_albedo + type(grid_t), pointer :: wavelength_grid + type(profile_t), pointer :: profile + type(error_t) :: error character(len=512) :: errmsg integer :: errcode - type(error_t) :: error - real(kind_phys) :: abs_error = 1e-4 integer :: i wavelength_grid => create_wavelength_grid(wavelength_grid_interfaces, errmsg, errcode) @@ -45,14 +45,13 @@ subroutine test_update_surface_albedo() call profile%get_edge_values( surface_albedo_interfaces, error) ASSERT(error%is_success()) - ASSERT_NEAR(surface_albedo_interfaces(1), expected_surface_albedo_interfaces(1), abs_error) - ASSERT_NEAR(surface_albedo_interfaces(2), expected_surface_albedo_interfaces(2), abs_error) - ASSERT_NEAR(surface_albedo_interfaces(3), expected_surface_albedo_interfaces(3), abs_error) - ASSERT_NEAR(surface_albedo_interfaces(4), expected_surface_albedo_interfaces(4), abs_error) - ASSERT_NEAR(surface_albedo_interfaces(5), expected_surface_albedo_interfaces(5), abs_error) + do i = 1, size(surface_albedo_interfaces) + ASSERT_NEAR(surface_albedo_interfaces(i), expected_surface_albedo_interfaces(i), ABS_ERROR) + end do deallocate( profile ) deallocate( wavelength_grid ) - end subroutine test_invalid_size_wavelength_interfaces + end subroutine test_update_surface_albedo + end program test_tuvx_surface_albedo \ No newline at end of file diff --git a/test/musica/tuvx/test_tuvx_temperature.F90 b/test/musica/tuvx/test_tuvx_temperature.F90 index 38a22744..cecde8b4 100644 --- a/test/musica/tuvx/test_tuvx_temperature.F90 +++ b/test/musica/tuvx/test_tuvx_temperature.F90 @@ -20,18 +20,19 @@ subroutine test_update_temperature() integer, parameter :: NUM_HOST_MIDPOINTS = 5 integer, parameter :: NUM_HOST_INTERFACES = 6 + real, parameter :: ABS_ERROR = 1e-4 real(kind_phys) :: host_midpoint_temperature(NUM_HOST_MIDPOINTS) = & [800.8_kind_phys, 700.7_kind_phys, 600.6_kind_phys, 500.5_kind_phys, 400.4_kind_phys] real(kind_phys) :: host_surface_temperature = 300.3_kind_phys - real(kind_phys) :: interface_temperatures(NUM_HOST_MIDPOINTS+2) - real(kind_phys) :: expected_interface_temperatures(NUM_HOST_MIDPOINTS+2) = & + real(kind_phys) :: expected_temperature_interfaces(NUM_HOST_MIDPOINTS+2) = & [300.3_kind_phys, 400.4_kind_phys, 500.5_kind_phys, 600.6_kind_phys, 700.7_kind_phys, 800.8_kind_phys, 800.8_kind_phys] + real(kind_phys) :: temperature_interfaces(NUM_HOST_MIDPOINTS+2) type(grid_t), pointer :: height_grid type(profile_t), pointer :: profile character(len=512) :: errmsg integer :: errcode type(error_t) :: error - real(kind_phys) :: abs_error = 1e-4 + integer :: i height_grid => create_height_grid(NUM_HOST_MIDPOINTS, NUM_HOST_INTERFACES, errmsg, errcode) ASSERT(errcode == 0) @@ -45,16 +46,11 @@ subroutine test_update_temperature() host_surface_temperature, errmsg, errcode ) ASSERT(errcode == 0) - call profile%get_edge_values( interface_temperatures, error) + call profile%get_edge_values( temperature_interfaces, error) ASSERT(error%is_success()) - - ASSERT_NEAR(interface_temperatures(1), expected_interface_temperatures(1), abs_error) - ASSERT_NEAR(interface_temperatures(2), expected_interface_temperatures(2), abs_error) - ASSERT_NEAR(interface_temperatures(3), expected_interface_temperatures(3), abs_error) - ASSERT_NEAR(interface_temperatures(4), expected_interface_temperatures(4), abs_error) - ASSERT_NEAR(interface_temperatures(5), expected_interface_temperatures(5), abs_error) - ASSERT_NEAR(interface_temperatures(6), expected_interface_temperatures(6), abs_error) - ASSERT_NEAR(interface_temperatures(7), expected_interface_temperatures(7), abs_error) + do i = 1, size(temperature_interfaces) + ASSERT_NEAR(temperature_interfaces(i), expected_temperature_interfaces(i), ABS_ERROR) + end do deallocate( profile ) deallocate( height_grid ) diff --git a/test/musica/tuvx/test_tuvx_wavelength_grid.F90 b/test/musica/tuvx/test_tuvx_wavelength_grid.F90 index 7e974de8..e90d49c8 100644 --- a/test/musica/tuvx/test_tuvx_wavelength_grid.F90 +++ b/test/musica/tuvx/test_tuvx_wavelength_grid.F90 @@ -16,19 +16,19 @@ subroutine test_create_wavelength_grid() use musica_tuvx_grid, only: grid_t use ccpp_kinds, only: kind_phys - integer, parameter :: NUM_WAVELENGTH_GRID_MIDPOINTS = 2 - integer, parameter :: NUM_WAVELENGTH_GRID_INTERFACES = 3 - real(kind_phys), target :: host_interfaces(NUM_WAVELENGTH_GRID_INTERFACES) = [180.0e-9_kind_phys, 200.0e-9_kind_phys, 240.0e-9_kind_phys] - real(kind_phys), target :: expected_interfaces(NUM_WAVELENGTH_GRID_INTERFACES) = [180.0_kind_phys, 200.0_kind_phys, 240.0_kind_phys] - real(kind_phys), target :: expected_midpoints(NUM_WAVELENGTH_GRID_MIDPOINTS) = [190.0_kind_phys, 220.0_kind_phys] - real(kind_phys) :: interfaces(NUM_WAVELENGTH_GRID_INTERFACES) - real(kind_phys) :: midpoints(NUM_WAVELENGTH_GRID_MIDPOINTS) - type(grid_t), pointer :: wavelength_grid => null() - character(len=512) :: errmsg - integer :: errcode - real(kind_phys) :: abs_error = 1e-5 - integer :: i - type(error_t) :: error + integer, parameter :: NUM_WAVELENGTH_GRID_MIDPOINTS = 2 + integer, parameter :: NUM_WAVELENGTH_GRID_INTERFACES = 3 + real, parameter :: ABS_ERROR = 1e-5 + real(kind_phys) :: host_interfaces(NUM_WAVELENGTH_GRID_INTERFACES) = [180.0e-9_kind_phys, 200.0e-9_kind_phys, 240.0e-9_kind_phys] + real(kind_phys) :: expected_interfaces(NUM_WAVELENGTH_GRID_INTERFACES) = [180.0_kind_phys, 200.0_kind_phys, 240.0_kind_phys] + real(kind_phys) :: expected_midpoints(NUM_WAVELENGTH_GRID_MIDPOINTS) = [190.0_kind_phys, 220.0_kind_phys] + real(kind_phys) :: interfaces(NUM_WAVELENGTH_GRID_INTERFACES) + real(kind_phys) :: midpoints(NUM_WAVELENGTH_GRID_MIDPOINTS) + type(grid_t), pointer :: wavelength_grid => null() + character(len=512) :: errmsg + integer :: errcode + type(error_t) :: error + integer :: i wavelength_grid => create_wavelength_grid(host_interfaces, errmsg, errcode) ASSERT(errcode == 0) @@ -37,13 +37,13 @@ subroutine test_create_wavelength_grid() call wavelength_grid%get_edges(interfaces, error) ASSERT(error%is_success()) do i = 1, NUM_WAVELENGTH_GRID_INTERFACES - ASSERT_NEAR(interfaces(i), expected_interfaces(i), abs_error) + ASSERT_NEAR(interfaces(i), expected_interfaces(i), ABS_ERROR) end do call wavelength_grid%get_midpoints(midpoints, error) ASSERT(error%is_success()) do i = 1, NUM_WAVELENGTH_GRID_MIDPOINTS - ASSERT_NEAR(midpoints(i), expected_midpoints(i), abs_error) + ASSERT_NEAR(midpoints(i), expected_midpoints(i), ABS_ERROR) end do deallocate(wavelength_grid) From 0b01c77f8c0939d05bb7a3123f1b363414cac4c7 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Fri, 25 Oct 2024 15:07:00 -0600 Subject: [PATCH 10/13] small adjustments to the musica code --- schemes/musica/musica_ccpp.F90 | 4 ++-- .../musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 | 12 +++++------- 2 files changed, 7 insertions(+), 9 deletions(-) diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90 index b7deb9a8..6890b1a3 100644 --- a/schemes/musica/musica_ccpp.F90 +++ b/schemes/musica/musica_ccpp.F90 @@ -29,8 +29,8 @@ subroutine musica_ccpp_init(vertical_layer_dimension, vertical_interface_dimensi photolysis_wavelength_grid_interfaces, errmsg, errcode) use ccpp_kinds, only : kind_phys - integer, intent(in) :: vertical_layer_dimension ! (count) - integer, intent(in) :: vertical_interface_dimension ! (count) + integer, intent(in) :: vertical_layer_dimension ! (count) + integer, intent(in) :: vertical_interface_dimension ! (count) real(kind_phys), intent(in) :: photolysis_wavelength_grid_interfaces(:) ! m character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 index 2d7dac71..3761a36f 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 @@ -27,9 +27,7 @@ function create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) & type(grid_t), intent(inout) :: wavelength_grid character(len=*), intent(out) :: errmsg integer, intent(out) :: errcode - - ! return value - type(profile_t), pointer :: profile + type(profile_t), pointer :: profile ! local variables type(error_t) :: error @@ -60,17 +58,17 @@ subroutine set_surface_albedo_values( profile, host_surface_albedo, & ! local variables type(error_t) :: error - real(kind_phys) :: surface_albedos(num_wavelength_bins + 1) + real(kind_phys) :: surface_albedo_interfaces(num_wavelength_bins + 1) - if (size(surface_albedos) <= DEFAULT_NUM_WAVELENGTH_BINS + 1) then + if (size(surface_albedo_interfaces) <= DEFAULT_NUM_WAVELENGTH_BINS + 1) then errmsg = "[MUSICA Error] Invalid size of TUV-x wavelength interfaces." errcode = 1 return end if - surface_albedos(:) = host_surface_albedo + surface_albedo_interfaces(:) = host_surface_albedo - call profile%set_edge_values( surface_albedos, error ) + call profile%set_edge_values( surface_albedo_interfaces, error ) if ( has_error_occurred( error, errmsg, errcode ) ) return end subroutine set_surface_albedo_values From 272be0f154459e405a6c08fb58011b3adcb3632d Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Fri, 25 Oct 2024 15:35:33 -0600 Subject: [PATCH 11/13] remove surface albedo entry from json --- .../musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 | 4 ++-- test/musica/test_musica_api.F90 | 2 +- test/musica/tuvx/configs/ts1_tsmlt.json | 10 ---------- 3 files changed, 3 insertions(+), 13 deletions(-) diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 index 3761a36f..751d18c6 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 @@ -6,9 +6,9 @@ module musica_ccpp_tuvx_surface_albedo surface_albedo_label, surface_albedo_unit !> Label for surface albedo in TUV-x - character(len=*), parameter :: surface_albedo_label = "surface_albedo" + character(len=*), parameter :: surface_albedo_label = "surface albedo" !> Unit for surface albedo in TUV-x - character(len=*), parameter :: surface_albedo_unit = "" ! unitless + character(len=*), parameter :: surface_albedo_unit = "none" !> Default value of number of wavelength bins integer, parameter :: DEFAULT_NUM_WAVELENGTH_BINS = 0 !> Number of wavelength bins diff --git a/test/musica/test_musica_api.F90 b/test/musica/test_musica_api.F90 index 7386de86..6b4d4154 100644 --- a/test/musica/test_musica_api.F90 +++ b/test/musica/test_musica_api.F90 @@ -160,7 +160,7 @@ subroutine test_musica_ccpp_api() geopotential_height_wrt_surface_at_interface(2,:) = (/ 3000.0_kind_phys, 500.0_kind_phys, -1500.0_kind_phys /) surface_temperature = (/ 300.0_kind_phys, 300.0_kind_phys /) surface_geopotential = (/ 100.0_kind_phys, 200.0_kind_phys /) - surface_albedo = 20.0_kind_phys + surface_albedo = 0.10_kind_phys standard_gravitational_acceleration = 10.0_kind_phys temperature(:,1) = (/ 100._kind_phys, 200._kind_phys /) temperature(:,2) = (/ 300._kind_phys, 400._kind_phys /) diff --git a/test/musica/tuvx/configs/ts1_tsmlt.json b/test/musica/tuvx/configs/ts1_tsmlt.json index c85d085d..b2a2f986 100644 --- a/test/musica/tuvx/configs/ts1_tsmlt.json +++ b/test/musica/tuvx/configs/ts1_tsmlt.json @@ -48,16 +48,6 @@ "month": 3, "day": 21 }, - { - "name": "surface albedo", - "type": "from config file", - "units": "none", - "uniform value": 0.10, - "grid": { - "name": "wavelength", - "units": "nm" - } - }, { "name": "extraterrestrial flux", "enable diagnostics" : true, From 066c8472840ad06e7868022b726d62a5fd2a68a2 Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Mon, 28 Oct 2024 11:54:21 -0600 Subject: [PATCH 12/13] address code review --- schemes/musica/micm/musica_ccpp_micm.F90 | 8 ++++---- schemes/musica/musica_ccpp.F90 | 2 +- schemes/musica/tuvx/musica_ccpp_tuvx.F90 | 7 ++++--- test/musica/tuvx/test_tuvx_surface_albedo.F90 | 6 +++--- 4 files changed, 12 insertions(+), 11 deletions(-) diff --git a/schemes/musica/micm/musica_ccpp_micm.F90 b/schemes/musica/micm/musica_ccpp_micm.F90 index 2419d40e..1ef50d5a 100644 --- a/schemes/musica/micm/musica_ccpp_micm.F90 +++ b/schemes/musica/micm/musica_ccpp_micm.F90 @@ -16,7 +16,7 @@ module musica_ccpp_micm contains - !> Register MICM constituent properties with the CCPP + !> Registers MICM constituent properties with the CCPP subroutine micm_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode) use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t use musica_micm, only: Rosenbrock, RosenbrockStandardOrder @@ -76,7 +76,7 @@ subroutine micm_register(solver_type, num_grid_cells, constituent_props, errmsg, end subroutine micm_register - !> Intitialize MICM + !> Initializes MICM subroutine micm_init(errmsg, errcode) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode @@ -86,7 +86,7 @@ subroutine micm_init(errmsg, errcode) end subroutine micm_init - !> Solve chemistry at the current time step + !> Solves chemistry at the current time step subroutine micm_run(time_step, temperature, pressure, dry_air_density, & user_defined_rate_parameters, constituents, errmsg, errcode) use musica_micm, only: solver_stats_t @@ -124,7 +124,7 @@ subroutine micm_run(time_step, temperature, pressure, dry_air_density, & end subroutine micm_run - !> Finalize MICM + !> Finalizes MICM subroutine micm_final(errmsg, errcode) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90 index 6890b1a3..09364397 100644 --- a/schemes/musica/musica_ccpp.F90 +++ b/schemes/musica/musica_ccpp.F90 @@ -46,7 +46,7 @@ end subroutine musica_ccpp_init !> \section arg_table_musica_ccpp_run Argument Table !! \htmlinclude musica_ccpp_run.html !! - !! The standard name for the variable 'surface_tempearture' is + !! The standard name for the variable 'surface_temperature' is !! 'blackbody_temperature_at_surface' because this is what we have as !! the standard name for 'cam_in%ts', whcih represents the same quantity. subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, constituent_props, & diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index 0d845310..86385a78 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -1,6 +1,6 @@ module musica_ccpp_tuvx - ! Note: "tuvx_t" is included in an external pre-built tuvx library that the host + ! Note: "tuvx_t" is included in an external pre-built TUV-x library that the host ! model is responsible for linking to during compilation use musica_tuvx, only: tuvx_t, grid_t, profile_t use musica_ccpp_util, only: has_error_occurred @@ -20,7 +20,7 @@ module musica_ccpp_tuvx contains - !> Intitialize TUV-x + !> Initializes TUV-x subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & wavelength_grid_interfaces, errmsg, errcode) use musica_tuvx, only: grid_map_t, profile_map_t, radiator_map_t @@ -231,6 +231,7 @@ subroutine tuvx_run(temperature, dry_air_density, & if (errcode /= 0) return end do + ! surface albedo related to short direct radiation call set_surface_albedo_values( surface_albedo_profile, surface_albedo, errmsg, errcode ) if (errcode /= 0) return @@ -239,7 +240,7 @@ subroutine tuvx_run(temperature, dry_air_density, & end subroutine tuvx_run - !> Finalize tuvx + !> Finalizes TUV-x subroutine tuvx_final(errmsg, errcode) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode diff --git a/test/musica/tuvx/test_tuvx_surface_albedo.F90 b/test/musica/tuvx/test_tuvx_surface_albedo.F90 index 0b4ce5dd..c47ebc63 100644 --- a/test/musica/tuvx/test_tuvx_surface_albedo.F90 +++ b/test/musica/tuvx/test_tuvx_surface_albedo.F90 @@ -19,11 +19,11 @@ subroutine test_update_surface_albedo() use ccpp_kinds, only: kind_phys integer, parameter :: NUM_WAVELENGTH_BINS = 4 - real, parameter :: ABS_ERROR = 1e-4 + real, parameter :: ABS_ERROR = 1e-6 real(kind_phys) :: wavelength_grid_interfaces(NUM_WAVELENGTH_BINS + 1) = & [200.0e-9_kind_phys, 210.0e-9_kind_phys, 240.0e-9_kind_phys, 300.0e-9_kind_phys, 400.0e-9_kind_phys] - real(kind_phys) :: host_surface_albedo = 500.5_kind_phys - real(kind_phys) :: expected_surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) = 500.5_kind_phys + real(kind_phys) :: host_surface_albedo = 0.3_kind_phys + real(kind_phys) :: expected_surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) = 0.3_kind_phys real(kind_phys) :: surface_albedo_interfaces(NUM_WAVELENGTH_BINS + 1) type(grid_t), pointer :: wavelength_grid type(profile_t), pointer :: profile From f4b5c293c006f9b12474953352f2b9bcc1a93d8b Mon Sep 17 00:00:00 2001 From: Jiwon Gim Date: Mon, 28 Oct 2024 12:09:10 -0600 Subject: [PATCH 13/13] update the comment for surface albedo --- schemes/musica/tuvx/musica_ccpp_tuvx.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index 86385a78..bcc7be09 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -231,7 +231,7 @@ subroutine tuvx_run(temperature, dry_air_density, & if (errcode /= 0) return end do - ! surface albedo related to short direct radiation + ! surface albedo with respect to direct UV/visible radiation call set_surface_albedo_values( surface_albedo_profile, surface_albedo, errmsg, errcode ) if (errcode /= 0) return