diff --git a/.gitignore b/.gitignore index 83c5d6be..8639122f 100644 --- a/.gitignore +++ b/.gitignore @@ -4,4 +4,5 @@ **/CMakeCache.txt **/CMakeFiles/ .vscode -xcode/ \ No newline at end of file +xcode/ +*.mod \ No newline at end of file diff --git a/schemes/musica/micm/musica_ccpp_micm.F90 b/schemes/musica/micm/musica_ccpp_micm.F90 index d2c6044d..5ba9eaf8 100644 --- a/schemes/musica/micm/musica_ccpp_micm.F90 +++ b/schemes/musica/micm/musica_ccpp_micm.F90 @@ -24,10 +24,14 @@ subroutine micm_register(solver_type, number_of_grid_cells, constituent_props, & use musica_micm, only: Rosenbrock, RosenbrockStandardOrder use musica_util, only: error_t use iso_c_binding, only: c_int + use musica_ccpp_gas_species, only: species_t, species_group_t, & + species_profiled_t, species_profiled_groupt_t,INDEX_UNINITIALIZED integer(c_int), intent(in) :: solver_type integer(c_int), intent(in) :: number_of_grid_cells type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituent_props(:) + type(species_groupt_t), allocatable, intent(out) :: species_group(:) + type(species_profiled_group_t), allocatable, intent(out) :: species_profiled_group(:) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode @@ -35,20 +39,36 @@ subroutine micm_register(solver_type, number_of_grid_cells, constituent_props, & type(error_t) :: error real(kind=kind_phys) :: molar_mass character(len=:), allocatable :: species_name - logical :: is_advected - integer :: i, species_index + logical :: is_advected, is_profiled + integer :: i, i_species_profiled + integer :: species_index + integer :: number_of_species + real(kind=kind_phys) :: species_unit, species_scale_height + + type(species_profiled_group_t), allocatable, :: temp_species_profiled_group(:) micm => micm_t(trim(filename_of_micm_configuration), solver_type, & number_of_grid_cells, error) if (has_error_occurred(error, errmsg, errcode)) return - allocate(constituent_props(micm%species_ordering%size()), stat=errcode) + number_of_species = micm%species_ordering%size() + + allocate(constituent_props(number_of_species), stat=errcode) if (errcode /= 0) then errmsg = "[MUSICA Error] Failed to allocate memory for constituent properties." return end if - do i = 1, micm%species_ordering%size() + allocate(temp_species_profiled_group(number_of_species), stat=errcode) + if (errcode /= 0) then + errmsg = "[MUSICA Error] Failed to allocate memory for musica species group." + return + end if + + i_species_profiled = 0 + + ! Sets constituents, creates species_t and species_profiled_t + do i = 1, number_of_species associate( map => micm%species_ordering ) species_name = map%name(i) species_index = map%index(i) @@ -74,8 +94,40 @@ subroutine micm_register(solver_type, number_of_grid_cells, constituent_props, & errcode = errcode, & errmsg = errmsg) if (errcode /= 0) return + + ! set musica species + musica_species_group(i) = musica_species_t(species_name, molar_mass, & + species_index, INDEX_UNINITIALIZED) + + is_profiled = micm%get_species_property_bool(species_name, & + "__is profiled", error) + if (has_error_occurred(error, errmsg, errcode)) return + + if (is_profiled) then + i_species_profiled = i_species_profiled + 1 + + species_unit = micm%get_species_property_double(species_name, & + "__unit", error) + if (has_error_occurred(error, errmsg, errcode)) return + + species_scale_height = micm%get_species_property_double(species_name, & + "__scale height [km]", error) + if (has_error_occurred(error, errmsg, errcode)) return + + musica_species_profiled_group(i_species_profiled) = & + musica_species_profiled_t( species_name, species_unit, molar_mass, & + species_scale_height, species_index, INDEX_UNINITIALIZED) + end if + end associate ! map end do + + if (i_species_profiled /= 0) then + allocate( species_profiled_group(i_species_profiled) ) + species_profiled_group(:) = temp_species_profiled_group(1:i_species_profiled) + end if + deallocate(temp_species_profiled_group) + number_of_rate_parameters = micm%user_defined_reaction_rates%size() end subroutine micm_register diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90 index d5526f52..55ec4806 100644 --- a/schemes/musica/musica_ccpp.F90 +++ b/schemes/musica/musica_ccpp.F90 @@ -10,24 +10,28 @@ module musica_ccpp public :: musica_ccpp_register, musica_ccpp_init, musica_ccpp_run, musica_ccpp_final + type(musica_sequence_t) :: musica_sequence contains !> \section arg_table_musica_ccpp_register Argument Table !! \htmlinclude musica_ccpp_register.html subroutine musica_ccpp_register(micm_solver_type, number_of_grid_cells, & - constituent_props, errmsg, errcode) + constituent_props, species_group, species_profiled_group, errmsg, errcode) use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t - + use musica_ccpp_gas_species, only: species_groupt_t, species_profiled_group_t + integer, intent(in) :: micm_solver_type integer, intent(in) :: number_of_grid_cells type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituent_props(:) + type(species_groupt_t), allocatable, intent(out) :: species_group(:) + type(species_profiled_group_t), allocatable, intent(out) :: species_profiled_group(:) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode type(ccpp_constituent_properties_t), allocatable :: constituent_props_subset(:) call micm_register(micm_solver_type, number_of_grid_cells, constituent_props_subset, & - errmsg, errcode) + species_group, species_profiled_group, errmsg, errcode) if (errcode /= 0) return constituent_props = constituent_props_subset deallocate(constituent_props_subset) @@ -44,9 +48,11 @@ subroutine musica_ccpp_init(vertical_layer_dimension, vertical_interface_dimensi photolysis_wavelength_grid_interfaces, & constituent_props, errmsg, errcode) use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t - use ccpp_kinds, only : kind_phys - use musica_ccpp_micm, only: micm - use musica_ccpp_util, only: has_error_occurred + use ccpp_kinds, only: kind_phys + use ccpp_const_utils, only: ccpp_const_get_idx + use musica_ccpp_micm, only: micm + use musica_ccpp_util, only: has_error_occurred + integer, intent(in) :: vertical_layer_dimension ! (count) integer, intent(in) :: vertical_interface_dimension ! (count) real(kind_phys), intent(in) :: photolysis_wavelength_grid_interfaces(:) ! m @@ -54,6 +60,31 @@ subroutine musica_ccpp_init(vertical_layer_dimension, vertical_interface_dimensi character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode + !!!!!!!!! JIWON + ! Search for the index based on the species name + do i_species = 1, size(species) + name = species%name + call ccpp_const_get_idx(constituent_props, name, index_constituents, errmsg, errcode) + if (errcode /= 0) return + species%index_constituents = index_constituents + array(i_species) = index_constituents + end do + + ! Get the molar mass that is set in the call to instantiate() + do i_species = 1, size(species) + call constituent_props(array(i_species))%molar_mass(molar_mass_arr(i_elem), errcode, errmsg) + if (errcode /= 0) then + errmsg = "[MUSICA Error] Unable to get molar mass." + return + end if + end do + !!!!!!!!! + !!!!!!!!! + + ! set molar mass, set species from tuvx + call set_musica_species() + + call micm_init(errmsg, errcode) if (errcode /= 0) return call tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & @@ -111,6 +142,39 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co number_of_rate_parameters) :: rate_parameters ! various units integer :: i_elem + real(kind_phys), dimension(size(constituents, dim=1), & + size(constituents, dim=2), & + number_of_musica_speices) :: constituents_musica + + ! ! TODO(jiwon) this might not be correct because it doesn't know the index + ! ! Get the molar mass that is set in the call to instantiate() + ! do i_elem = 1, size(molar_mass_arr) + ! call constituent_props(i_elem)%molar_mass(molar_mass_arr(i_elem), errcode, errmsg) + ! if (errcode /= 0) then + ! errmsg = "[MUSICA Error] Unable to get molar mass." + ! return + ! end if + ! end do + + ! ! TODO(jiwon) Check molar mass is non zero as it becomes a denominator for unit converison + ! ! this code will be deleted when the framework does the check + ! do i_elem = 1, size(molar_mass_arr) + ! if (molar_mass_arr(i_elem) <= 0) then + ! errcode = 1 + ! errmsg = "[MUSICA Error] Molar mass must be greater than zero." + ! return + ! end if + ! end do + + ! ! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3) + ! call convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents) + + ! find musica constituents array + call get_musica_constituents(constituents, constituents_musica, sequence_musica_species) + + ! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3) + call convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents) + ! Calculate photolysis rate constants using TUV-x call tuvx_run(temperature, dry_air_density, & geopotential_height_wrt_surface_at_midpoint, & @@ -121,39 +185,20 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co photolysis_wavelength_grid_interfaces, & extraterrestrial_flux, & standard_gravitational_acceleration, & - cloud_area_fraction, constituents, & + cloud_area_fraction, constituents_musica, & ! need to reset the arrangement air_pressure_thickness, rate_parameters, & errmsg, errcode) - ! Get the molar mass that is set in the call to instantiate() - do i_elem = 1, size(molar_mass_arr) - call constituent_props(i_elem)%molar_mass(molar_mass_arr(i_elem), errcode, errmsg) - if (errcode /= 0) then - errmsg = "[MUSICA Error] Unable to get molar mass." - return - end if - end do - - ! TODO(jiwon) Check molar mass is non zero as it becomes a denominator for unit converison - ! this code will be deleted when the framework does the check - do i_elem = 1, size(molar_mass_arr) - if (molar_mass_arr(i_elem) <= 0) then - errcode = 1 - errmsg = "[MUSICA Error] Molar mass must be greater than zero." - return - end if - end do - - ! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3) - call convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents) - ! Solve chemistry at the current time step call micm_run(time_step, temperature, pressure, dry_air_density, rate_parameters, & - constituents, errmsg, errcode) + constituents_musica, errmsg, errcode) ! Convert MICM unit back to CAM-SIMA unit (mol m-3 -> kg kg-1) - call convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constituents) + ! return as CAM-SIMA constituents array + call convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constituents_musica) + call update_musica_constituents(constituents_musica, constituents, sequence_musica_species) + end subroutine musica_ccpp_run !> \section arg_table_musica_ccpp_final Argument Table diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index 700b894e..7045cf64 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -7,30 +7,35 @@ module musica_ccpp_tuvx use musica_ccpp_util, only: has_error_occurred use musica_tuvx, only: tuvx_t, grid_t, profile_t, radiator_t use musica_util, only: mappings_t, index_mappings_t + use musica_ccpp_tuvx_gas_species_profile, only: gas_species_t, profile_group_t implicit none private public :: tuvx_register, tuvx_init, tuvx_run, tuvx_final - type(tuvx_t), pointer :: tuvx => null() - 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() - type(profile_t), pointer :: extraterrestrial_flux_profile => null() - type(radiator_t), pointer :: cloud_optics => null() - type(index_mappings_t), pointer :: photolysis_rate_constants_mapping => null( ) - integer, parameter :: DEFAULT_NUM_PHOTOLYSIS_RATE_CONSTANTS = 0 - integer :: number_of_photolysis_rate_constants = DEFAULT_NUM_PHOTOLYSIS_RATE_CONSTANTS - integer, parameter :: DEFAULT_INDEX_NOT_FOUND = -1 - character(len=*), parameter :: CLOUD_LIQUID_WATER_CONTENT_LABEL = & + type(tuvx_t), pointer :: tuvx => null() + 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() + type(profile_t), pointer :: extraterrestrial_flux_profile => null() + type(radiator_t), pointer :: cloud_optics => null() + type(index_mappings_t), pointer :: photolysis_rate_constants_mapping => null( ) + + type(gas_species_t), allocatable :: gas_species_group(:) + type(profile_group_t), allocatable :: profile_gas_species_group(:) + + integer, parameter :: DEFAULT_NUM_PHOTOLYSIS_RATE_CONSTANTS = 0 + integer :: number_of_photolysis_rate_constants = DEFAULT_NUM_PHOTOLYSIS_RATE_CONSTANTS + integer, parameter :: DEFAULT_INDEX_NOT_FOUND = -1 + character(len=*), parameter :: CLOUD_LIQUID_WATER_CONTENT_LABEL = & 'cloud_liquid_water_mixing_ratio_wrt_moist_air_and_condensed_water' - character(len=*), parameter :: CLOUD_LIQUID_WATER_CONTENT_LONG_NAME = & + character(len=*), parameter :: CLOUD_LIQUID_WATER_CONTENT_LONG_NAME = & 'Cloud water mass mixing ratio with respect to moist air plus all airborne condensates' - character(len=*), parameter :: CLOUD_LIQUID_WATER_CONTENT_UNITS = 'kg kg-1' - real(kind_phys), parameter :: CLOUD_LIQUID_WATER_CONTENT_MOLAR_MASS = 0.018_kind_phys ! kg mol-1 - integer :: index_cloud_liquid_water_content = DEFAULT_INDEX_NOT_FOUND + character(len=*), parameter :: CLOUD_LIQUID_WATER_CONTENT_UNITS = 'kg kg-1' + real(kind_phys), parameter :: CLOUD_LIQUID_WATER_CONTENT_MOLAR_MASS = 0.018_kind_phys ! kg mol-1 + integer :: index_cloud_liquid_water_content = DEFAULT_INDEX_NOT_FOUND contains @@ -48,8 +53,13 @@ subroutine reset_tuvx_map_state( grids, profiles, radiators ) end subroutine reset_tuvx_map_state - !> This is a helper subroutine created to deallocate objects associated with TUV-x + !> Deallocates objects associated with TUV-x subroutine cleanup_tuvx_resources() + use musica_ccpp_tuvx_gas_species_profile, & + only: deallocate_gas_species_profile_group + + ! local variable + integer :: i_species if (associated( height_grid )) then deallocate( height_grid ) @@ -86,6 +96,8 @@ subroutine cleanup_tuvx_resources() photolysis_rate_constants_mapping => null() end if + call deallocate_gas_species_profile_group( profile_gas_species_group ) + end subroutine cleanup_tuvx_resources !> Registers constituent properties with the CCPP needed by TUV-x @@ -142,6 +154,8 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & extraterrestrial_flux_unit use musica_ccpp_tuvx_cloud_optics, & only: create_cloud_optics_radiator, cloud_optics_label + use musica_ccpp_tuvx_gas_species_profile, & + only: configure_gas_species, create_gas_species_profile_group integer, intent(in) :: vertical_layer_dimension ! (count) integer, intent(in) :: vertical_interface_dimension ! (count) @@ -152,12 +166,14 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & integer, intent(out) :: errcode ! local variables - type(grid_map_t), pointer :: grids - type(profile_map_t), pointer :: profiles - type(radiator_map_t), pointer :: radiators - type(configuration_t) :: config - type(mappings_t), pointer :: photolysis_rate_constants_ordering - type(error_t) :: error + type(error_t) :: error + type(grid_map_t), pointer :: grids + type(profile_map_t), pointer :: profiles + type(radiator_map_t), pointer :: radiators + type(configuration_t) :: config + type(mappings_t), pointer :: photolysis_rate_constants_ordering + integer :: number_of_gas_species + integer :: i_species ! Get needed indices in constituents array call ccpp_const_get_idx(constituent_props, CLOUD_LIQUID_WATER_CONTENT_LABEL, & @@ -252,6 +268,34 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & return end if + call configure_gas_species( constituent_props, gas_species_group, & + errmsg, errcode ) + if (errcode /= 0) then + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + end if + + number_of_gas_species = size(gas_species_group) + allocate( profile_gas_species_group( number_of_gas_species ) ) + + call create_gas_species_profile_group( height_grid, & + gas_species_group, profile_gas_species_group, errmsg, errcode) + if (errcode /= 0) then + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + endif + + do i_species = 1, number_of_gas_species + call profiles%add( profile_gas_species_group(i_species)%profile, error ) + if (has_error_occurred( error, errmsg, errcode )) then + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + end if + end do + radiators => radiator_map_t( error ) if (has_error_occurred( error, errmsg, errcode )) then call reset_tuvx_map_state( grids, profiles, null() ) @@ -350,6 +394,22 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & return end if + if (.not. allocated( profile_gas_species_group)) then + allocate( profile_gas_species_group( number_of_gas_species ) ) + end if + + do i_species = 1, number_of_gas_species + profile_gas_species_group(i_species)%profile => & + profiles%get( gas_species_group(i_species)%label, gas_species_group(i_species)%unit, error ) + if (has_error_occurred( error, errmsg, errcode )) then + deallocate( tuvx ) + tuvx => null() + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + end if + end do + radiators => tuvx%get_radiators( error ) if (has_error_occurred( error, errmsg, errcode )) then deallocate( tuvx ) @@ -418,13 +478,14 @@ subroutine tuvx_run(temperature, dry_air_density, & cloud_area_fraction, constituents, & air_pressure_thickness, rate_parameters, & 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_ccpp_util, only: has_error_occurred - use musica_ccpp_tuvx_surface_albedo, only: set_surface_albedo_values - use musica_ccpp_tuvx_extraterrestrial_flux, only: set_extraterrestrial_flux_values - use musica_ccpp_tuvx_cloud_optics, only: set_cloud_optics_values + use musica_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_util, only: has_error_occurred + use musica_ccpp_tuvx_surface_albedo, only: set_surface_albedo_values + use musica_ccpp_tuvx_extraterrestrial_flux, only: set_extraterrestrial_flux_values + use musica_ccpp_tuvx_cloud_optics, only: set_cloud_optics_values + use musica_ccpp_tuvx_gas_species_profile, only: set_gas_species_values, MOLAR_MASS_DRY_AIR real(kind_phys), intent(in) :: temperature(:,:) ! K (column, layer) real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 (column, layer) @@ -447,6 +508,8 @@ subroutine tuvx_run(temperature, dry_air_density, & ! local variables real(kind_phys), dimension(size(geopotential_height_wrt_surface_at_midpoint, dim = 2)) :: height_midpoints real(kind_phys), dimension(size(geopotential_height_wrt_surface_at_interface, dim = 2)) :: height_interfaces + real(kind_phys), dimension(size(height_interfaces)) :: height_deltas ! km + real(kind_phys), dimension(size(constituents, dim = 2)) :: gas_species_constituents real(kind_phys), dimension(size(rate_parameters, dim=2)+2, & number_of_photolysis_rate_constants) :: photolysis_rate_constants, & ! s-1 heating_rates ! K s-1 (TODO: check units) @@ -454,7 +517,7 @@ subroutine tuvx_run(temperature, dry_air_density, & real(kind_phys) :: solar_zenith_angle ! degrees real(kind_phys) :: earth_sun_distance ! AU type(error_t) :: error - integer :: i_col, i_level + integer :: i_col, i_level, i_gas_species reciprocal_of_gravitational_acceleration = 1.0_kind_phys / standard_gravitational_acceleration @@ -469,13 +532,14 @@ subroutine tuvx_run(temperature, dry_air_density, & if (errcode /= 0) return do i_col = 1, size(temperature, dim=1) + call calculate_heights( geopotential_height_wrt_surface_at_midpoint(i_col,:), & geopotential_height_wrt_surface_at_interface(i_col,:), & surface_geopotential(i_col), & reciprocal_of_gravitational_acceleration, & height_midpoints, height_interfaces ) call set_height_grid_values( height_grid, height_midpoints, height_interfaces, & - errmsg, errcode ) + height_deltas, errmsg, errcode ) if (errcode /= 0) return call set_temperature_values( temperature_profile, temperature(i_col,:), & @@ -489,6 +553,19 @@ subroutine tuvx_run(temperature, dry_air_density, & errmsg, errcode ) if (errcode /= 0) return + do i_gas_species = 1, size(gas_species_group) + if (gas_species_group(i_gas_species)%label == "air") then + gas_species_constituents = dry_air_density(i_col,:) * MOLAR_MASS_DRY_AIR + else + gas_species_constituents = constituents(i_col,:,gas_species_group(i_gas_species)%index_constituent_props) + end if + call set_gas_species_values(profile_gas_species_group(i_gas_species)%profile, & + gas_species_group(i_gas_species), & + gas_species_constituents, & ! mol m-3 + height_deltas, errmsg, errcode) + if (errcode /= 0) return + end do + ! temporary values until these are available from the host model solar_zenith_angle = 0.0_kind_phys earth_sun_distance = 1.0_kind_phys @@ -515,12 +592,15 @@ end subroutine tuvx_run !> Finalizes TUV-x subroutine tuvx_final(errmsg, errcode) + use musica_ccpp_tuvx_gas_species_profile, only: deallocate_gas_species + character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode errmsg = '' errcode = 0 + call deallocate_gas_species( gas_species_group ) call cleanup_tuvx_resources() if (associated( tuvx )) then diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_gas_species_profile.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_gas_species_profile.F90 new file mode 100644 index 00000000..8a13d230 --- /dev/null +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_gas_species_profile.F90 @@ -0,0 +1,195 @@ +module musica_ccpp_tuvx_gas_species_profile + use ccpp_kinds, only: kind_phys + use musica_tuvx_profile, only: profile_t + use musica_ccpp_gas_species, only: gas_species_t + implicit none + + private + public :: configure_gas_species, create_gas_species_profile_group, & + set_gas_species_values, deallocate_gas_species, deallocate_gas_species_profile_group + + !> Conversion factor from km to cm + real(kind_phys), parameter, public :: km_to_cm = 1.0e5 + !> Conversion factor from m3 to cm3 + real(kind_phys), parameter, public :: m_3_to_cm_3 = 1.0e6 + !> Molar mass value of dry air is obtained from 'CAM-SIMA/src/utils/std_atm_profile.F90' + ! TODO(jiwon) + real(kind_phys), parameter, public :: MOLAR_MASS_DRY_AIR = 0.0289644_kind_phys ! kg mol-1 + + !> Defines a type to contaion a pointer to 'profile_t' + type, public :: profile_group_t + type(profile_t), pointer :: profile + end type profile_group_t + +contains + + !> Configures gas species. + ! Note: The user can customize this function to configure gas species in the system + ! This function allocates memory for the gas species that is needed in tuvx_init and tuvx_run. + ! The user is responsible for deallocating the memory. + subroutine configure_gas_species(constituent_props, gas_species_group, errmsg, errcode) + use ccpp_const_utils, only: ccpp_const_get_idx + use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t + + type(ccpp_constituent_prop_ptr_t), intent(in) :: constituent_props(:) + type(gas_species_t), allocatable, intent(out) :: gas_species_group(:) + character(len=512), intent(out) :: errmsg + integer, intent(out) :: errcode + + ! local variable + integer :: num_gas_species = 3 + real(kind_phys), parameter :: SCALE_HEIGHT_DRY_AIR = 8.01_kind_phys ! km + real(kind_phys), parameter :: SCALE_HEIGHT_O2 = 7.0_kind_phys ! km + real(kind_phys), parameter :: SCALE_HEIGHT_O3 = 7.0_kind_phys ! km + real(kind_phys) :: molar_mass_O2 ! kg mol-1 + real(kind_phys) :: molar_mass_O3 ! kg mol-1 + integer, parameter :: INDEX_NOT_KNOWN = -9999 + integer :: index_air, index_O2, index_O3 + + allocate(gas_species_group(num_gas_species)) ! need to deallocate + + ! TODO(jiwon) - I commented out this block of code that searches for the molar mass + ! of air and instead hard-coded it using the value imported from CAM-SIMA. + ! Should we register air during the registration phase, similar to cloud liquid water content? + ! Or should we include it in the configuration file, or is it fine as it is? + + ! Air + ! Find the index of the species and create a 'gas_species_t' for it + ! call ccpp_const_get_idx(constituent_props, "air", index_air, errmsg, errcode) + ! write(*,*) "air index: ", index_air + ! if (errcode /= 0) return + + ! call constituent_props(index_air)%molar_mass(molar_mass_air, errcode, errmsg) + ! if (errcode /= 0) return + + is_advected = micm%get_species_property_bool(species_name, & + "__is advected", & + error) + + + index_air = INDEX_NOT_KNOWN + gas_species_group(1) = gas_species_t("air", "molecule cm-3", & + MOLAR_MASS_DRY_AIR, SCALE_HEIGHT_DRY_AIR, index_air) + ! O2 + call ccpp_const_get_idx(constituent_props, "O2", index_O2, errmsg, errcode) + if (errcode /= 0) return + + call constituent_props(index_O2)%molar_mass(molar_mass_O2, errcode, errmsg) + if (errcode /= 0) return + + gas_species_group(2) = gas_species_t("O2", "molecule cm-3", & + molar_mass_O2, SCALE_HEIGHT_O2, index_O2) + + ! O3 + call ccpp_const_get_idx(constituent_props, "O3", index_O3, errmsg, errcode) + if (errcode /= 0) return + + call constituent_props(index_O3)%molar_mass(molar_mass_O3, errcode, errmsg) + if (errcode /= 0) return + + gas_species_group(3) = gas_species_t("O3", "molecule cm-3", & + molar_mass_O3, SCALE_HEIGHT_O3, index_O3) + + end subroutine configure_gas_species + + !> Creates an array of profile group that contains a 'profile_t' for each gas species + subroutine create_gas_species_profile_group(height_grid, gas_species_group, & + profile_group, errmsg, errcode) + 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 + + type(grid_t), intent(in) :: height_grid + type(gas_species_t), intent(in) :: gas_species_group(:) + type(profile_group_t), intent(inout) :: profile_group(:) + character(len=512), intent(out) :: errmsg + integer, intent(out) :: errcode + + ! local variables + type(error_t) :: error + integer :: i_species + + do i_species = 1, size(gas_species_group) + profile_group(i_species)%profile => profile_t( gas_species_group(i_species)%label, & + gas_species_group(i_species)%unit, height_grid, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + end do + + end subroutine create_gas_species_profile_group + + !> Sets interface, density, and above-column density values for gas species + subroutine set_gas_species_values(profile, gas_species, constituent, & + height_deltas, errmsg, errcode) + use musica_ccpp_util, only: has_error_occurred + use musica_tuvx_profile, only: profile_t + use musica_util, only: error_t + + type(profile_t), intent(inout) :: profile + type(gas_species_t), intent(in) :: gas_species + real(kind_phys), intent(in) :: constituent(:) ! mol m-3 + real(kind_phys), intent(in) :: height_deltas(:) ! km, change in height in each vertical layer + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + + ! local variables + type(error_t) :: error + integer :: num_constituents + real(kind_phys) :: constituent_mol_per_cm_3(size(constituent)) ! mol cm-3 + real(kind_phys) :: interfaces(size(constituent) + 2) + real(kind_phys) :: densities(size(constituent) + 1) + + num_constituents = size(constituent) + constituent_mol_per_cm_3(:) = constituent(:) / m_3_to_cm_3 + interfaces(1) = constituent_mol_per_cm_3(num_constituents) + interfaces(2:num_constituents+1) = constituent_mol_per_cm_3(num_constituents:1:-1) + interfaces(num_constituents+2) = constituent_mol_per_cm_3(1) + + if (gas_species%label == "O3") then + densities(:) = height_deltas(:) * km_to_cm & + * ( interfaces(1:num_constituents+1) & + + interfaces(2:num_constituents+2) ) * 0.5_kind_phys + else + densities(:) = height_deltas(:) * km_to_cm & + * sqrt(interfaces(1:num_constituents+1)) & + * sqrt(interfaces(2:num_constituents+2)) + end if + call profile%set_edge_values( interfaces, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + call profile%set_layer_densities( densities, error) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + call profile%calculate_exo_layer_density( gas_species%scale_height, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + end subroutine set_gas_species_values + + !> Deallocates gas species object + subroutine deallocate_gas_species(gas_species_group) + type(gas_species_t), allocatable, intent(inout) :: gas_species_group(:) + + if (allocated( gas_species_group )) deallocate( gas_species_group ) + + end subroutine + + !> Deallocates profiles of each gas species + subroutine deallocate_gas_species_profile_group(profile_group) + type(profile_group_t), allocatable, intent(inout) :: profile_group(:) + + ! local variables + integer :: i_species + + if (allocated( profile_group )) then + do i_species = 1, size(profile_group) + if (associated( profile_group(i_species)%profile) ) then + deallocate( profile_group(i_species)%profile) + end if + end do + + deallocate( profile_group ) + end if + + end subroutine deallocate_gas_species_profile_group + +end module musica_ccpp_tuvx_gas_species_profile \ No newline at end of file diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_height_grid.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_height_grid.F90 index cff0d2a2..ec341196 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_height_grid.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_height_grid.F90 @@ -84,7 +84,7 @@ end function create_height_grid !> Sets TUV-x height grid values from the host-model height grid subroutine set_height_grid_values(height_grid, host_midpoints, & - host_interfaces, errmsg, errcode) + host_interfaces, height_deltas, errmsg, errcode) use ccpp_kinds, only: kind_phys use musica_ccpp_util, only: has_error_occurred use musica_tuvx_grid, only: grid_t @@ -93,6 +93,7 @@ subroutine set_height_grid_values(height_grid, host_midpoints, & type(grid_t), intent(inout) :: height_grid real(kind_phys), intent(in) :: host_midpoints(:) ! km real(kind_phys), intent(in) :: host_interfaces(:) ! km + real(kind_phys), intent(inout) :: height_deltas(:) ! km character(len=*), intent(out) :: errmsg integer, intent(out) :: errcode @@ -102,6 +103,9 @@ subroutine set_height_grid_values(height_grid, host_midpoints, & real(kind_phys) :: interfaces(size(host_interfaces)+1) integer :: n_host_midpoints, n_host_interfaces + n_host_midpoints = size(host_midpoints) + n_host_interfaces = size(host_interfaces) + if ( size(midpoints) /= height_grid%number_of_sections( error ) ) then errmsg = "[MUSICA Error] Invalid size of TUV-x mid-point heights." errcode = 1 @@ -114,8 +118,11 @@ subroutine set_height_grid_values(height_grid, host_midpoints, & return end if - n_host_midpoints = size(host_midpoints) - n_host_interfaces = size(host_interfaces) + if ( size(height_deltas) /= n_host_interfaces ) then + errmsg = "[MUSICA Error] Invalid size of TUV-x height deltas." + errcode = 1 + return + end if interfaces(1) = host_interfaces(n_host_interfaces) interfaces(2:n_host_interfaces) = host_midpoints(n_host_midpoints:1:-1) @@ -127,6 +134,8 @@ subroutine set_height_grid_values(height_grid, host_midpoints, & midpoints(n_host_midpoints+1) = 0.5_kind_phys * & ( interfaces(n_host_interfaces) + interfaces(n_host_interfaces+1) ) + height_deltas(:) = interfaces(2:size(interfaces)) - interfaces(1:size(interfaces)-1) + call height_grid%set_edges( interfaces, error ) if ( has_error_occurred( error, errmsg, errcode ) ) return diff --git a/test/docker/Dockerfile.musica b/test/docker/Dockerfile.musica index aa00a0da..0eed2488 100644 --- a/test/docker/Dockerfile.musica +++ b/test/docker/Dockerfile.musica @@ -6,7 +6,7 @@ FROM ubuntu:22.04 ARG MUSICA_GIT_TAG=326b5119768d5be9654baf96ae3bd6a1b757fdc8 -ARG CAM_SIMA_CHEMISTRY_DATA_TAG=abc7cacbec3d33d5c0ed5bb79a157e93b42c45c0 +ARG CAM_SIMA_CHEMISTRY_DATA_TAG=9e84370d917fa4357cb3e7d9a0f3b57869fb342c ARG BUILD_TYPE=Debug RUN apt update \ diff --git a/test/docker/Dockerfile.musica.no_install b/test/docker/Dockerfile.musica.no_install index 5baec757..be0b63e0 100644 --- a/test/docker/Dockerfile.musica.no_install +++ b/test/docker/Dockerfile.musica.no_install @@ -9,7 +9,7 @@ FROM ubuntu:22.04 ARG MUSICA_GIT_TAG=326b5119768d5be9654baf96ae3bd6a1b757fdc8 -ARG CAM_SIMA_CHEMISTRY_DATA_TAG=abc7cacbec3d33d5c0ed5bb79a157e93b42c45c0 +ARG CAM_SIMA_CHEMISTRY_DATA_TAG=9e84370d917fa4357cb3e7d9a0f3b57869fb342c ARG BUILD_TYPE=Debug RUN apt update \ diff --git a/test/musica/test_musica_api.F90 b/test/musica/test_musica_api.F90 index eb185fbc..2976c42a 100644 --- a/test/musica/test_musica_api.F90 +++ b/test/musica/test_musica_api.F90 @@ -7,8 +7,13 @@ program run_test_musica_ccpp #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 + write(*,*) "[MUSICA Test] Running the Chapman test" call test_chapman() + write(*,*) "[MUSICA Test] Ends the Chapman test" + + write(*,*) "[MUSICA Test] Running the Terminator test" call test_terminator() + write(*,*) "[MUSICA Test] Ends the Terminator test" contains @@ -144,8 +149,8 @@ subroutine test_chapman() implicit none - integer, parameter :: NUM_SPECIES = 5 - integer, parameter :: NUM_TUVX_CONSTITUENTS = 1 + integer, parameter :: NUM_SPECIES = 5 + integer, parameter :: NUM_TUVX_CONSTITUENTS = 1 ! This test requires that the number of grid cells = 4, which is the default ! vector dimension for MICM. This restriction will be removed once ! https://github.com/NCAR/musica/issues/217 is finished. @@ -235,11 +240,11 @@ subroutine test_chapman() ASSERT(errcode == 0) call constituent_props(i)%is_advected(is_advected, errcode, errmsg) ASSERT(errcode == 0) - tmp_bool = (trim(species_name) == "O2" .and. molar_mass == 0.0319988_kind_phys .and. is_advected) .or. & + tmp_bool = (trim(species_name) == "O2" .and. molar_mass == 0.0319988_kind_phys .and. .not. is_advected) .or. & (trim(species_name) == "O" .and. molar_mass == 0.0159994_kind_phys .and. .not. is_advected) .or. & (trim(species_name) == "O1D" .and. molar_mass == 0.0159994_kind_phys .and. .not. is_advected) .or. & (trim(species_name) == "O3" .and. molar_mass == 0.0479982_kind_phys .and. is_advected) .or. & - (trim(species_name) == "N2" .and. molar_mass == 0.0280134_kind_phys .and. is_advected) .or. & + (trim(species_name) == "N2" .and. molar_mass == 0.0280134_kind_phys .and. .not. is_advected) .or. & (trim(species_name) == "cloud_liquid_water_mixing_ratio_wrt_moist_air_and_condensed_water" .and. & molar_mass == 0.018_kind_phys .and. is_advected) ASSERT(tmp_bool) @@ -372,8 +377,8 @@ subroutine test_terminator() implicit none - integer, parameter :: NUM_SPECIES = 2 - integer, parameter :: NUM_TUVX_CONSTITUENTS = 1 + integer, parameter :: NUM_SPECIES = 4 + integer, parameter :: NUM_TUVX_CONSTITUENTS = 1 ! This test requires that the number of grid cells = 4, which is the default ! vector dimension for MICM. This restriction will be removed once ! https://github.com/NCAR/musica/issues/217 is finished. @@ -412,7 +417,7 @@ subroutine test_terminator() character(len=:), allocatable :: micm_species_name logical :: tmp_bool, is_advected integer :: i, j, k - integer :: Cl_index, Cl2_index + integer :: Cl_index, Cl2_index, O2_index, O3_index real(kind_phys) :: total_Cl, total_Cl_init call get_wavelength_edges(photolysis_wavelength_grid_interfaces) @@ -465,8 +470,10 @@ subroutine test_terminator() ASSERT(errcode == 0) tmp_bool = (trim(species_name) == "Cl" .and. molar_mass == 0.035453_kind_phys .and. is_advected) .or. & (trim(species_name) == "Cl2" .and. molar_mass == 0.070906_kind_phys .and. is_advected) .or. & - (trim(species_name) == "cloud_liquid_water_mixing_ratio_wrt_moist_air_and_condensed_water" & - .and. molar_mass == 0.018_kind_phys .and. is_advected) + (trim(species_name) == "O2" .and. molar_mass == 0.0319988_kind_phys .and. .not. is_advected) .or. & + (trim(species_name) == "O3" .and. molar_mass == 0.0479982_kind_phys .and. .not. is_advected) .or. & + (trim(species_name) == "cloud_liquid_water_mixing_ratio_wrt_moist_air_and_condensed_water" & + .and. molar_mass == 0.018_kind_phys .and. is_advected) ASSERT(tmp_bool) call constituent_props(i)%units(units, errcode, errmsg) if (errcode /= 0) then @@ -501,6 +508,12 @@ subroutine test_terminator() else if (micm_species_name == "Cl2") then Cl2_index = i base_conc = 1.0e-6_kind_phys + else if (micm_species_name == "O2") then + O2_index = i + base_conc = 0.21_kind_phys + else if (micm_species_name == "O3") then + O3_index = i + base_conc = 1.0e-4_kind_phys else write(*,*) "Unknown species: ", micm_species_name stop 3 diff --git a/test/musica/tuvx/CMakeLists.txt b/test/musica/tuvx/CMakeLists.txt index 21301b0b..14f26b46 100644 --- a/test/musica/tuvx/CMakeLists.txt +++ b/test/musica/tuvx/CMakeLists.txt @@ -171,3 +171,36 @@ add_test( ) add_memory_check_test(test_tuvx_cloud_optics $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) + +# Gas species profile +add_executable(test_tuvx_gas_species_profile test_tuvx_gas_species_profile.F90) + +target_sources(test_tuvx_gas_species_profile + PUBLIC + ${MUSICA_CCPP_SOURCES} + ${TO_BE_CCPPIZED_SRC_PATH}/ccpp_const_utils.F90 + ${TO_BE_CCPPIZED_SRC_PATH}/ccpp_tuvx_utils.F90 + ${CCPP_SRC_PATH}/ccpp_constituent_prop_mod.F90 + ${CCPP_SRC_PATH}/ccpp_hash_table.F90 + ${CCPP_SRC_PATH}/ccpp_hashable.F90 + ${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90 + ${CMAKE_CURRENT_SOURCE_DIR}/../musica_ccpp_namelist.F90 +) + +target_link_libraries(test_tuvx_gas_species_profile + PRIVATE + musica::musica-fortran +) + +set_target_properties(test_tuvx_gas_species_profile + PROPERTIES + LINKER_LANGUAGE Fortran +) + +add_test( + NAME test_tuvx_gas_species_profile + COMMAND $ + WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY} +) + +add_memory_check_test(test_tuvx_gas_species_profile $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) \ No newline at end of file diff --git a/test/musica/tuvx/test_tuvx_gas_species_profile.F90 b/test/musica/tuvx/test_tuvx_gas_species_profile.F90 new file mode 100644 index 00000000..14ce341e --- /dev/null +++ b/test/musica/tuvx/test_tuvx_gas_species_profile.F90 @@ -0,0 +1,125 @@ +program test_tuvx_gas_species_profile + + use musica_ccpp_tuvx_gas_species_profile + + 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_gas_species_profile() + +contains + + subroutine test_create_gas_species_profile() + use musica_micm, only: Rosenbrock, RosenbrockStandardOrder + use musica_tuvx_grid, only: grid_t + use ccpp_kinds, only: kind_phys + use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t + use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t + use musica_ccpp, only: musica_ccpp_register + use musica_ccpp_tuvx_height_grid, only: create_height_grid + use musica_ccpp_namelist, only: filename_of_micm_configuration, & + filename_of_tuvx_configuration, & + filename_of_tuvx_micm_mapping_configuration + + integer, parameter :: NUM_COLUMNS = 2 + integer, parameter :: NUM_LAYERS = 2 + integer :: NUM_GRID_CELLS = NUM_COLUMNS * NUM_LAYERS + integer :: solver_type = Rosenbrock + ! Molar mass for O2, O3 is read from 'cam-sima-chemistry-data/mechanisms/chapman/micm/species.json' + real(kind_phys), parameter :: MOLAR_MASS_O2 = 0.0319988_kind_phys + real(kind_phys), parameter :: MOLAR_MASS_O3 = 0.0479982_kind_phys + integer, parameter :: NUMBER_GAS_SPECIES = 3 + integer, parameter :: INDEX_NOT_KNOWN = -9999 + 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 + type(gas_species_t), allocatable :: gas_species_group(:) + type(profile_group_t), allocatable :: profile_gas_species_group(:) + type(grid_t), pointer :: height_grid => null() + real(kind_phys) :: dry_air_density(NUM_COLUMNS,NUM_LAYERS) ! kg m-3 + real(kind_phys) :: gas_species_constituents(NUM_COLUMNS,NUM_LAYERS) + real(kind_phys) :: height_deltas(NUM_LAYERS) ! km + integer :: errcode + character(len=512) :: errmsg + character(len=50) :: label + character(len=50) :: unit + real(kind_phys) :: molar_mass ! kg mol-1 + real(kind_phys) :: scale_height ! km + integer :: index + logical :: tmp_bool + integer :: i, i_gas_species, i_col + + filename_of_micm_configuration = 'musica_configurations/chapman/micm/config.json' + filename_of_tuvx_configuration = 'musica_configurations/chapman/tuvx/config.json' + filename_of_tuvx_micm_mapping_configuration = 'musica_configurations/chapman/tuvx_micm_mapping.json' + + dry_air_density(:,1) = (/ 3.5_kind_phys, 4.5_kind_phys /) + dry_air_density(:,2) = (/ 5.5_kind_phys, 6.5_kind_phys /) + gas_species_constituents(:,1) = (/ 0.35_kind_phys, 0.45_kind_phys /) + gas_species_constituents(:,2) = (/ 0.55_kind_phys, 0.65_kind_phys /) + height_deltas(:) = (/ 0.5_kind_phys, 1.5_kind_phys /) + + call musica_ccpp_register(solver_type, NUM_GRID_CELLS, constituent_props, errmsg, errcode) + if (errcode /= 0) then + write(*,*) trim(errmsg) + stop 3 + endif + + allocate(constituent_props_ptr(size(constituent_props))) + do i = 1, size(constituent_props) + const_prop => constituent_props(i) + call constituent_props_ptr(i)%set( const_prop, errcode, errmsg ) + end do + + height_grid => create_height_grid( NUM_LAYERS, NUM_LAYERS + 1 , errmsg, errcode ) + ASSERT(errcode == 0) + ASSERT(associated(height_grid)) + + call configure_gas_species( constituent_props_ptr, gas_species_group, errmsg, errcode ) + ASSERT(errcode == 0) + ASSERT(size(gas_species_group) == NUMBER_GAS_SPECIES) + + allocate( profile_gas_species_group( size(gas_species_group) ) ) + + call create_gas_species_profile_group( height_grid, gas_species_group, & + profile_gas_species_group, errmsg, errcode ) + ASSERT(errcode == 0) + ASSERT(size(profile_gas_species_group) == NUMBER_GAS_SPECIES) + + do i_gas_species = 1, size(gas_species_group) + label = gas_species_group(i_gas_species)%label + unit = gas_species_group(i_gas_species)%unit + molar_mass = gas_species_group(i_gas_species)%molar_mass + scale_height = gas_species_group(i_gas_species)%scale_height + index = gas_species_group(i_gas_species)%index_constituent_props + + tmp_bool = (trim(label) == "air" .and. trim(unit) == "molecule cm-3" .and. molar_mass == MOLAR_MASS_DRY_AIR .and. & + scale_height == 8.01_kind_phys .and. index == INDEX_NOT_KNOWN) .or. & + (trim(label) == "O2" .and. trim(unit) == "molecule cm-3" .and. molar_mass == MOLAR_MASS_O2 .and. & + scale_height == 7.0_kind_phys .and. index /= INDEX_NOT_KNOWN .and. index /= 0) .or. & + (trim(label) == "O3" .and. trim(unit) == "molecule cm-3" .and. molar_mass == MOLAR_MASS_O3 .and. & + scale_height == 7.0_kind_phys .and. index /= INDEX_NOT_KNOWN .and. index /= 0) + ASSERT(tmp_bool) + + ASSERT(associated(profile_gas_species_group(i_gas_species)%profile)) + end do + + do i_col = 1, NUM_COLUMNS + do i_gas_species = 1, size(gas_species_group) + call set_gas_species_values(profile_gas_species_group(i_gas_species)%profile, & + gas_species_group(i_gas_species), & + gas_species_constituents(i_col, :), & ! mol m-3 + height_deltas, errmsg, errcode) + if (errcode /= 0) return + end do + end do + + call deallocate_gas_species(gas_species_group) + call deallocate_gas_species_profile_group(profile_gas_species_group) + + end subroutine test_create_gas_species_profile + +end program test_tuvx_gas_species_profile \ 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 6da26960..c18ce338 100644 --- a/test/musica/tuvx/test_tuvx_height_grid.F90 +++ b/test/musica/tuvx/test_tuvx_height_grid.F90 @@ -25,14 +25,19 @@ subroutine test_create_height_grid() real(kind_phys) :: host_interfaces(NUM_HOST_INTERFACES) = [250.3_kind_phys, 150.2_kind_phys, 50.1_kind_phys] 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) :: expected_height_deltas(NUM_HOST_INTERFACES) real(kind_phys) :: midpoints(NUM_HOST_MIDPOINTS+1) real(kind_phys) :: interfaces(NUM_HOST_INTERFACES+1) + real(kind_phys) :: height_deltas(NUM_HOST_INTERFACES) type(grid_t), pointer :: height_grid => null() type(error_t) :: error character(len=512) :: errmsg integer :: errcode integer :: i + expected_height_deltas(:) = expected_interfaces(2:size(expected_interfaces)) & + - expected_interfaces(1:size(expected_interfaces)-1) + height_grid => create_height_grid(-1, 0, errmsg, errcode) ASSERT(errcode == 1) ASSERT(.not. associated(height_grid)) @@ -47,12 +52,16 @@ subroutine test_create_height_grid() ASSERT(associated(height_grid)) call set_height_grid_values(height_grid, host_midpoints, host_interfaces, & - errmsg, errcode) + height_deltas, errmsg, errcode) ASSERT(errcode == 0) ASSERT(height_grid%number_of_sections(error) == NUM_HOST_MIDPOINTS + 1) ASSERT(error%is_success()) + do i = 1, size(height_deltas) + ASSERT_NEAR(height_deltas(i), expected_height_deltas(i), ABS_ERROR) + end do + call height_grid%get_midpoints(midpoints, error) ASSERT(error%is_success()) do i = 1, size(midpoints) diff --git a/to_be_ccppized/ccpp_tuvx_utils.F90 b/to_be_ccppized/ccpp_tuvx_utils.F90 index 837801bc..e7329d56 100644 --- a/to_be_ccppized/ccpp_tuvx_utils.F90 +++ b/to_be_ccppized/ccpp_tuvx_utils.F90 @@ -10,7 +10,7 @@ module ccpp_tuvx_utils contains !> Regrids normalized flux data to match a specified wavelength grid - !! This function is copied from CAM/src/chemistry/utils/mo_util.F90 + ! This function is copied from CAM/src/chemistry/utils/mo_util.F90 subroutine rebin( source_dimension, target_dimension, source_coordinates, & target_coordinates, source, target ) use ccpp_kinds, only: kind_phys