diff --git a/doc/NamesNotInDictionary.txt b/doc/NamesNotInDictionary.txt index 57566b4d..99c2e48b 100644 --- a/doc/NamesNotInDictionary.txt +++ b/doc/NamesNotInDictionary.txt @@ -1,7 +1,7 @@ ####################### Date/time of when script was run: -2024-11-18 15:44:54.993446 +2024-12-16 11:19:41.130567 ####################### Non-dictionary standard names found in the following metadata files: @@ -108,12 +108,6 @@ atmospheric_physics/schemes/dry_adiabatic_adjust/dadadj.meta -------------------------- -atmospheric_physics/schemes/dry_adiabatic_adjust/dadadj_apply_qv_tendency.meta - - - tendency_of_water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water - --------------------------- - atmospheric_physics/schemes/zhang_mcfarlane/zm_conv_momtran.meta - atmosphere_detrainment_convective_mass_flux_for_deep_convection_for_convective_columns @@ -252,6 +246,12 @@ atmospheric_physics/schemes/utilities/geopotential_temp.meta -------------------------- +atmospheric_physics/schemes/utilities/physics_tendency_updaters.meta + + - ccpp_constituent_tendencies + +-------------------------- + atmospheric_physics/schemes/check_energy/check_energy_save_teout.meta - vertically_integrated_total_energy_at_end_of_physics_timestep @@ -342,17 +342,17 @@ atmospheric_physics/schemes/tropopause_find/tropopause_find.meta - tropopause_air_pressure - tropopause_air_pressure_from_chemical_method - tropopause_air_pressure_from_climatological_method - - tropopause_air_pressure_from_climatology_dataset - tropopause_air_pressure_from_cold_point_method - tropopause_air_pressure_from_hybrid_stobie_linoz_with_climatological_backup_method - tropopause_air_pressure_from_lapse_rate_method + - tropopause_air_pressure_from_tropopause_climatology_dataset - tropopause_air_temperature - tropopause_air_temperature_from_chemical_method - tropopause_air_temperature_from_climatological_method - tropopause_air_temperature_from_cold_point_method - tropopause_air_temperature_from_hybrid_stobie_linoz_with_climatological_backup_method - tropopause_air_temperature_from_lapse_rate_method - - tropopause_calendar_days_from_climatology + - tropopause_calendar_days_from_tropopause_climatology - tropopause_geopotential_height_wrt_surface - tropopause_geopotential_height_wrt_surface_from_chemical_method - tropopause_geopotential_height_wrt_surface_from_climatological_method @@ -374,10 +374,18 @@ atmospheric_physics/schemes/tropopause_find/tropopause_find.meta atmospheric_physics/schemes/musica/musica_ccpp.meta - blackbody_temperature_at_surface + - cloud_area_fraction - dynamic_constituents_for_musica_ccpp - - micm_solver_type - - number_of_grid_cells + - earth_sun_distance + - extraterrestrial_radiation_flux - photolysis_wavelength_grid_interfaces + - solar_zenith_angle - surface_albedo_due_to_UV_and_VIS_direct +-------------------------- + +atmospheric_physics/test/test_schemes/initialize_constituents.meta + + - dynamic_constituents_for_initialize_constituents + ####################### diff --git a/schemes/check_energy/check_energy_chng.F90 b/schemes/check_energy/check_energy_chng.F90 index 91af5151..7f7780f8 100644 --- a/schemes/check_energy/check_energy_chng.F90 +++ b/schemes/check_energy/check_energy_chng.F90 @@ -89,6 +89,9 @@ subroutine check_energy_chng_timestep_init( & character(len=512), intent(out) :: errmsg ! error message integer, intent(out) :: errflg ! error flag + errmsg = '' + errflg = 0 + !------------------------------------------------ ! Physics total energy. !------------------------------------------------ @@ -276,6 +279,9 @@ subroutine check_energy_chng_run( & integer :: i + errmsg = '' + errflg = 0 + !------------------------------------------------ ! Physics total energy. !------------------------------------------------ diff --git a/schemes/check_energy/check_energy_fix.F90 b/schemes/check_energy/check_energy_fix.F90 index 1bdc7ae8..5a3163c8 100644 --- a/schemes/check_energy/check_energy_fix.F90 +++ b/schemes/check_energy/check_energy_fix.F90 @@ -11,7 +11,7 @@ module check_energy_fix ! Add heating rate required for global mean total energy conservation !> \section arg_table_check_energy_fix_run Argument Table !! \htmlinclude arg_table_check_energy_fix_run.html - subroutine check_energy_fix_run(ncol, pver, pint, gravit, heat_glob, ptend_s, eshflx, scheme_name) + subroutine check_energy_fix_run(ncol, pver, pint, gravit, heat_glob, ptend_s, eshflx, scheme_name, errmsg, errflg) ! Input arguments integer, intent(in) :: ncol ! number of atmospheric columns integer, intent(in) :: pver ! number of vertical layers @@ -22,11 +22,17 @@ subroutine check_energy_fix_run(ncol, pver, pint, gravit, heat_glob, ptend_s, es real(kind_phys), intent(out) :: eshflx(:) ! effective sensible heat flux [W m-2] ! for check_energy_chng + ! Output arguments character(len=64), intent(out) :: scheme_name ! scheme name + character(len=512), intent(out) :: errmsg ! error message + integer, intent(out) :: errflg ! error flag ! Local variables integer :: i + errmsg = '' + errflg = 0 + ! Set scheme name for check_energy_chng scheme_name = "check_energy_fix" diff --git a/schemes/check_energy/check_energy_fix.meta b/schemes/check_energy/check_energy_fix.meta index b5f935fb..176a57d7 100644 --- a/schemes/check_energy/check_energy_fix.meta +++ b/schemes/check_energy/check_energy_fix.meta @@ -53,3 +53,15 @@ type = character | kind = len=64 dimensions = () intent = out +[ errmsg ] + standard_name = ccpp_error_message + units = none + type = character | kind = len=512 + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + units = 1 + type = integer + dimensions = () + intent = out diff --git a/schemes/check_energy/check_energy_gmean/check_energy_gmean.F90 b/schemes/check_energy/check_energy_gmean/check_energy_gmean.F90 index 4d50c428..7aa37f98 100644 --- a/schemes/check_energy/check_energy_gmean/check_energy_gmean.F90 +++ b/schemes/check_energy/check_energy_gmean/check_energy_gmean.F90 @@ -19,7 +19,8 @@ subroutine check_energy_gmean_run( & pint, & te_ini_dyn, teout, & tedif_glob, heat_glob, & - teinp_glob, teout_glob, psurf_glob, ptopb_glob) + teinp_glob, teout_glob, psurf_glob, ptopb_glob, & + errmsg, errflg) ! This scheme is non-portable due to dependency: Global mean module gmean from src/utils use gmean_mod, only: gmean @@ -43,10 +44,16 @@ subroutine check_energy_gmean_run( & real(kind_phys), intent(out) :: psurf_glob ! global mean surface pressure [Pa] real(kind_phys), intent(out) :: ptopb_glob ! global mean top boundary pressure [Pa] + character(len=512), intent(out) :: errmsg ! error message + integer, intent(out) :: errflg ! error flag + ! Local variables real(kind_phys) :: te(ncol, 4) ! total energy of input/output states (copy) real(kind_phys) :: te_glob(4) ! global means of total energy + errmsg = '' + errflg = 0 + ! Copy total energy out of input and output states. ! These four fields will have their global means calculated respectively te(:ncol, 1) = te_ini_dyn(:ncol) ! Input energy using dycore energy formula [J m-2] diff --git a/schemes/check_energy/check_energy_gmean/check_energy_gmean.meta b/schemes/check_energy/check_energy_gmean/check_energy_gmean.meta index 8ebf6f83..9aa7b57d 100644 --- a/schemes/check_energy/check_energy_gmean/check_energy_gmean.meta +++ b/schemes/check_energy/check_energy_gmean/check_energy_gmean.meta @@ -84,3 +84,15 @@ type = real | kind = kind_phys dimensions = () intent = out +[ errmsg ] + standard_name = ccpp_error_message + units = none + type = character | kind = len=512 + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + units = 1 + type = integer + dimensions = () + intent = out diff --git a/schemes/check_energy/check_energy_save_teout.F90 b/schemes/check_energy/check_energy_save_teout.F90 index 6ccdc15c..6f60eb3a 100644 --- a/schemes/check_energy/check_energy_save_teout.F90 +++ b/schemes/check_energy/check_energy_save_teout.F90 @@ -13,7 +13,7 @@ module check_energy_save_teout !> \section arg_table_check_energy_save_teout_run Argument Table !! \htmlinclude arg_table_check_energy_save_teout_run.html - subroutine check_energy_save_teout_run(ncol, te_cur_dyn, teout) + subroutine check_energy_save_teout_run(ncol, te_cur_dyn, teout, errmsg, errflg) ! Input arguments integer, intent(in) :: ncol ! number of atmospheric columns @@ -21,6 +21,11 @@ subroutine check_energy_save_teout_run(ncol, te_cur_dyn, teout) ! Output arguments real(kind_phys), intent(out) :: teout(:) ! total energy for global fixer in next timestep [J m-2] + character(len=512), intent(out) :: errmsg ! error message + integer, intent(out) :: errflg ! error flag + + errmsg = '' + errflg = 0 ! nb hplin: note that in physpkg.F90, the pbuf is updated to the previous dyn timestep ! through itim_old. Need to check if we need to replicate such pbuf functionality diff --git a/schemes/check_energy/check_energy_save_teout.meta b/schemes/check_energy/check_energy_save_teout.meta index 57a712a0..d587b6f7 100644 --- a/schemes/check_energy/check_energy_save_teout.meta +++ b/schemes/check_energy/check_energy_save_teout.meta @@ -23,3 +23,15 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = out +[ errmsg ] + standard_name = ccpp_error_message + units = none + type = character | kind = len=512 + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + units = 1 + type = integer + dimensions = () + intent = out diff --git a/schemes/check_energy/check_energy_zero_fluxes.F90 b/schemes/check_energy/check_energy_zero_fluxes.F90 index bda0d74c..52a444e2 100644 --- a/schemes/check_energy/check_energy_zero_fluxes.F90 +++ b/schemes/check_energy/check_energy_zero_fluxes.F90 @@ -12,7 +12,7 @@ module check_energy_zero_fluxes !> \section arg_table_check_energy_zero_fluxes_run Argument Table !! \htmlinclude arg_table_check_energy_zero_fluxes_run.html - subroutine check_energy_zero_fluxes_run(ncol, name, flx_vap, flx_cnd, flx_ice, flx_sen) + subroutine check_energy_zero_fluxes_run(ncol, name, flx_vap, flx_cnd, flx_ice, flx_sen, errmsg, errflg) ! Input arguments integer, intent(in) :: ncol ! number of atmospheric columns @@ -22,6 +22,11 @@ subroutine check_energy_zero_fluxes_run(ncol, name, flx_vap, flx_cnd, flx_ice, f real(kind_phys), intent(out) :: flx_cnd(:) ! boundary flux of liquid+ice (precip?) [m s-1] real(kind_phys), intent(out) :: flx_ice(:) ! boundary flux of ice (snow?) [m s-1] real(kind_phys), intent(out) :: flx_sen(:) ! boundary flux of sensible heat [W m-2] + character(len=512), intent(out) :: errmsg ! error message + integer, intent(out) :: errflg ! error flag + + errmsg = '' + errflg = 0 ! reset values to zero name = '' diff --git a/schemes/check_energy/check_energy_zero_fluxes.meta b/schemes/check_energy/check_energy_zero_fluxes.meta index 65782599..2c40e524 100644 --- a/schemes/check_energy/check_energy_zero_fluxes.meta +++ b/schemes/check_energy/check_energy_zero_fluxes.meta @@ -41,3 +41,15 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = out +[ errmsg ] + standard_name = ccpp_error_message + units = none + type = character | kind = len=512 + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + units = 1 + type = integer + dimensions = () + intent = out diff --git a/schemes/check_energy/dycore_energy_consistency_adjust.F90 b/schemes/check_energy/dycore_energy_consistency_adjust.F90 index 14fce382..4eaa6d97 100644 --- a/schemes/check_energy/dycore_energy_consistency_adjust.F90 +++ b/schemes/check_energy/dycore_energy_consistency_adjust.F90 @@ -16,7 +16,8 @@ subroutine dycore_energy_consistency_adjust_run( & do_consistency_adjust, & scaling_dycore, & tend_dTdt, & - tend_dTdt_local) + tend_dTdt_local, & + errmsg, errflg) ! Input arguments integer, intent(in) :: ncol ! number of atmospheric columns @@ -26,7 +27,12 @@ subroutine dycore_energy_consistency_adjust_run( & real(kind_phys), intent(in) :: tend_dTdt(:,:) ! model physics temperature tendency [K s-1] ! Output arguments - real(kind_phys), intent(out) :: tend_dTdt_local(:,:) ! (scheme) temperature tendency [K s-1] + real(kind_phys), intent(out) :: tend_dTdt_local(:,:) ! (scheme) temperature tendency [K s-1] + character(len=512), intent(out) :: errmsg ! error message + integer, intent(out) :: errflg ! error flag + + errmsg = '' + errflg = 0 if (do_consistency_adjust) then ! original formula for scaling of temperature: diff --git a/schemes/check_energy/dycore_energy_consistency_adjust.meta b/schemes/check_energy/dycore_energy_consistency_adjust.meta index e1afb294..05593fd3 100644 --- a/schemes/check_energy/dycore_energy_consistency_adjust.meta +++ b/schemes/check_energy/dycore_energy_consistency_adjust.meta @@ -41,3 +41,15 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = out +[ errmsg ] + standard_name = ccpp_error_message + units = none + type = character | kind = len=512 + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + units = 1 + type = integer + dimensions = () + intent = out diff --git a/schemes/dry_adiabatic_adjust/dadadj.meta b/schemes/dry_adiabatic_adjust/dadadj.meta index ee423db0..0818d3e4 100644 --- a/schemes/dry_adiabatic_adjust/dadadj.meta +++ b/schemes/dry_adiabatic_adjust/dadadj.meta @@ -119,6 +119,7 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = out + constituent = True [ dadpdf ] standard_name = binary_indicator_for_dry_adiabatic_adjusted_grid_cell units = fraction diff --git a/schemes/dry_adiabatic_adjust/dadadj_apply_qv_tendency.F90 b/schemes/dry_adiabatic_adjust/dadadj_apply_qv_tendency.F90 deleted file mode 100644 index 5987e313..00000000 --- a/schemes/dry_adiabatic_adjust/dadadj_apply_qv_tendency.F90 +++ /dev/null @@ -1,33 +0,0 @@ -module dadadj_apply_qv_tendency - - use ccpp_kinds, only: kind_phys - - implicit none - private - - public :: dadadj_apply_qv_tendency_run - -CONTAINS - - !> \section arg_table_dadadj_apply_qv_tendency_run Argument Table - !! \htmlinclude dadadj_apply_qv_tendency_run.html - subroutine dadadj_apply_qv_tendency_run(q_tend, state_q, dt, errmsg, errcode) - - ! update the constituent state. - ! Replace this with standard constitutent update function. - - ! Dummy arguments - real(kind_phys), intent(in) :: q_tend(:,:) ! water vapor tendency - real(kind_phys), intent(inout) :: state_q(:,:) ! water vapor - real(kind_phys), intent(in) :: dt ! physics time step - character(len=512), intent(out) :: errmsg - integer, intent(out) :: errcode - - errcode = 0 - errmsg = '' - - state_q = state_q + (q_tend * dt) - - end subroutine dadadj_apply_qv_tendency_run - -end module dadadj_apply_qv_tendency diff --git a/schemes/dry_adiabatic_adjust/dadadj_apply_qv_tendency.meta b/schemes/dry_adiabatic_adjust/dadadj_apply_qv_tendency.meta deleted file mode 100644 index 7b2b2f59..00000000 --- a/schemes/dry_adiabatic_adjust/dadadj_apply_qv_tendency.meta +++ /dev/null @@ -1,42 +0,0 @@ -[ccpp-table-properties] - name = dadadj_apply_qv_tendency - type = scheme -######################################################### -[ccpp-arg-table] - name = dadadj_apply_qv_tendency_run - type = scheme -[ q_tend ] - standard_name = tendency_of_water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water - units = kg kg-1 s-1 - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = in -[ state_q ] - standard_name = water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water - long_name = mass mixing ratio of water vapor / moist air and condensed water - advected = True - units = kg kg-1 - dimensions = (horizontal_loop_extent, vertical_layer_dimension) - type = real | kind = kind_phys - intent = inout -[ dt ] - standard_name = timestep_for_physics - long_name = time step - units = s - dimensions = () - type = real | kind = kind_phys - intent = in -[ errmsg ] - standard_name = ccpp_error_message - long_name = Error message for error handling in CCPP - units = none - type = character | kind = len=512 - dimensions = () - intent = out -[ errcode ] - standard_name = ccpp_error_code - long_name = Error flag for error handling in CCPP - units = 1 - type = integer - dimensions = () - intent = out diff --git a/schemes/held_suarez/held_suarez_1994.meta b/schemes/held_suarez/held_suarez_1994.meta index c43017c6..59b8c48b 100644 --- a/schemes/held_suarez/held_suarez_1994.meta +++ b/schemes/held_suarez/held_suarez_1994.meta @@ -50,7 +50,7 @@ intent = in [ clat ] standard_name = latitude - units = radians + units = rad type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = in diff --git a/schemes/musica/micm/musica_ccpp_micm.F90 b/schemes/musica/micm/musica_ccpp_micm.F90 index 7b9fbae9..c2e40d95 100644 --- a/schemes/musica/micm/musica_ccpp_micm.F90 +++ b/schemes/musica/micm/musica_ccpp_micm.F90 @@ -18,14 +18,15 @@ module musica_ccpp_micm contains !> Registers MICM constituent properties with the CCPP - subroutine micm_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode) + subroutine micm_register(solver_type, number_of_grid_cells, constituent_props, & + errmsg, errcode) use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t use musica_micm, only: Rosenbrock, RosenbrockStandardOrder use musica_util, only: error_t use iso_c_binding, only: c_int integer(c_int), intent(in) :: solver_type - integer(c_int), intent(in) :: num_grid_cells + integer(c_int), intent(in) :: number_of_grid_cells type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituent_props(:) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode @@ -37,7 +38,12 @@ subroutine micm_register(solver_type, num_grid_cells, constituent_props, errmsg, logical :: is_advected integer :: i, species_index - micm => micm_t(trim(filename_of_micm_configuration), solver_type, num_grid_cells, error) + if (associated( micm )) then + deallocate( micm ) + micm => null() + end if + 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) @@ -108,7 +114,6 @@ subroutine micm_run(time_step, temperature, pressure, dry_air_density, & type(string_t) :: solver_state type(solver_stats_t) :: solver_stats type(error_t) :: error - integer :: i_elem call micm%solve(real(time_step, kind=c_double), & c_loc(temperature), & diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90 index d3512f4a..4c79511f 100644 --- a/schemes/musica/musica_ccpp.F90 +++ b/schemes/musica/musica_ccpp.F90 @@ -2,7 +2,7 @@ module musica_ccpp use musica_ccpp_micm, only: micm_register, micm_init, micm_run, micm_final use musica_ccpp_namelist, only: filename_of_tuvx_micm_mapping_configuration - use musica_ccpp_tuvx, only: tuvx_init, tuvx_run, tuvx_final + use musica_ccpp_tuvx, only: tuvx_register, tuvx_init, tuvx_run, tuvx_final use musica_util, only: index_mappings_t implicit none @@ -14,37 +14,68 @@ module musica_ccpp !> \section arg_table_musica_ccpp_register Argument Table !! \htmlinclude musica_ccpp_register.html - subroutine musica_ccpp_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode) + subroutine musica_ccpp_register(constituent_props, errmsg, & + errcode) use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t + use musica_ccpp_namelist, only: micm_solver_type - integer, intent(in) :: solver_type - integer, intent(in) :: num_grid_cells type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituent_props(:) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode - call micm_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode) + type(ccpp_constituent_properties_t), allocatable :: constituent_props_subset(:) + integer :: number_of_grid_cells + + ! Temporary fix until the number of grid cells is only needed to create a MICM state + ! instead of when the solver is created. + ! The number of grid cells is not known at this point, so we set it to 1 and recreate + ! the solver when the number of grid cells is known at the init stage. + number_of_grid_cells = 1 + call micm_register(micm_solver_type, number_of_grid_cells, constituent_props_subset, & + errmsg, errcode) + if (errcode /= 0) return + constituent_props = constituent_props_subset + deallocate(constituent_props_subset) + + call tuvx_register(constituent_props_subset, errmsg, errcode) + if (errcode /= 0) return + constituent_props = [ constituent_props, constituent_props_subset ] 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, & - photolysis_wavelength_grid_interfaces, errmsg, errcode) + subroutine musica_ccpp_init(horizontal_dimension, vertical_layer_dimension, & + vertical_interface_dimension, & + photolysis_wavelength_grid_interfaces, & + constituent_props, errmsg, errcode) + use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t, ccpp_constituent_prop_ptr_t use ccpp_kinds, only : kind_phys use musica_ccpp_micm, only: micm + use musica_ccpp_namelist, only: micm_solver_type 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 - character(len=512), intent(out) :: errmsg - integer, intent(out) :: errcode - + integer, intent(in) :: horizontal_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 + type(ccpp_constituent_prop_ptr_t), intent(in) :: constituent_props(:) + character(len=512), intent(out) :: errmsg + integer, intent(out) :: errcode + + integer :: number_of_grid_cells + type(ccpp_constituent_properties_t), allocatable :: micm_species_props(:) + + ! Temporary fix until the number of grid cells is only needed to create a MICM state + ! instead of when the solver is created. + ! Re-create the MICM solver with the correct number of grid cells + number_of_grid_cells = horizontal_dimension * vertical_layer_dimension + call micm_register(micm_solver_type, number_of_grid_cells, micm_species_props, errmsg, errcode) call micm_init(errmsg, errcode) if (errcode /= 0) return call tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & photolysis_wavelength_grid_interfaces, & - micm%user_defined_reaction_rates, errmsg, errcode) + micm%user_defined_reaction_rates, & + constituent_props, errmsg, errcode) if (errcode /= 0) return end subroutine musica_ccpp_init @@ -55,11 +86,14 @@ end subroutine musica_ccpp_init !! 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, & - constituents, geopotential_height_wrt_surface_at_midpoint, & - geopotential_height_wrt_surface_at_interface, surface_temperature, & - surface_geopotential, surface_albedo, & - standard_gravitational_acceleration, errmsg, errcode) + 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_geopotential, & + surface_temperature, surface_albedo, & + photolysis_wavelength_grid_interfaces, extraterrestrial_flux, & + standard_gravitational_acceleration, cloud_area_fraction, & + air_pressure_thickness, solar_zenith_angle, & + earth_sun_distance, errmsg, errcode) use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t use ccpp_kinds, only: kind_phys use musica_ccpp_micm, only: number_of_rate_parameters @@ -74,10 +108,16 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co real(kind_phys), target, intent(inout) :: constituents(:,:,:) ! kg kg-1 real(kind_phys), intent(in) :: geopotential_height_wrt_surface_at_midpoint(:,:) ! m (column, layer) 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_temperature(:) ! K real(kind_phys), intent(in) :: surface_albedo ! unitless + real(kind_phys), intent(in) :: photolysis_wavelength_grid_interfaces(:) ! nm + real(kind_phys), intent(in) :: extraterrestrial_flux(:) ! photons cm-2 s-1 nm-1 real(kind_phys), intent(in) :: standard_gravitational_acceleration ! m s-2 + real(kind_phys), intent(in) :: cloud_area_fraction(:,:) ! unitless (column, level) + real(kind_phys), intent(in) :: air_pressure_thickness(:,:) ! Pa (column, level) + real(kind_phys), intent(in) :: solar_zenith_angle(:) ! radians (column) + real(kind_phys), intent(in) :: earth_sun_distance ! AU character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode @@ -89,13 +129,20 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co integer :: i_elem ! Calculate photolysis rate constants using TUV-x - call 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, & - rate_parameters, & + call tuvx_run(temperature, dry_air_density, & + constituents, & + geopotential_height_wrt_surface_at_midpoint, & + geopotential_height_wrt_surface_at_interface, & + surface_geopotential, surface_temperature, & + surface_albedo, & + photolysis_wavelength_grid_interfaces, & + extraterrestrial_flux, & + standard_gravitational_acceleration, & + cloud_area_fraction, & + air_pressure_thickness, & + solar_zenith_angle, & + earth_sun_distance, & + rate_parameters, & errmsg, errcode) ! Get the molar mass that is set in the call to instantiate() diff --git a/schemes/musica/musica_ccpp.meta b/schemes/musica/musica_ccpp.meta index 85e01f3e..801c3c79 100644 --- a/schemes/musica/musica_ccpp.meta +++ b/schemes/musica/musica_ccpp.meta @@ -6,18 +6,6 @@ [ccpp-arg-table] name = musica_ccpp_register type = scheme -[ solver_type ] - standard_name = micm_solver_type - units = none - type = integer - dimensions = () - intent = in -[ num_grid_cells ] - standard_name = number_of_grid_cells - units = count - type = integer - dimensions = () - intent = in [ constituent_props ] standard_name = dynamic_constituents_for_musica_ccpp units = none @@ -41,15 +29,21 @@ [ccpp-arg-table] name = musica_ccpp_init type = scheme +[ horizontal_dimension ] + standard_name = horizontal_dimension + units = count + type = integer + dimensions = () + intent = in [ vertical_layer_dimension ] standard_name = vertical_layer_dimension - units = none + units = count type = integer dimensions = () intent = in [ vertical_interface_dimension ] standard_name = vertical_interface_dimension - units = none + units = count type = integer dimensions = () intent = in @@ -59,6 +53,12 @@ type = real | kind = kind_phys dimensions = (photolysis_wavelength_grid_interface_dimension) intent = in +[ constituent_props ] + standard_name = ccpp_constituent_properties + units = None + type = ccpp_constituent_prop_ptr_t + dimensions = (number_of_ccpp_constituents) + intent = in [ errmsg ] standard_name = ccpp_error_message units = none @@ -123,30 +123,66 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent,vertical_interface_dimension) intent = in -[ surface_temperature ] - standard_name = blackbody_temperature_at_surface - type = real | kind = kind_phys - units = K - dimensions = (horizontal_loop_extent) - intent = in [ surface_geopotential ] standard_name = surface_geopotential type = real | kind = kind_phys units = m2 s-2 dimensions = (horizontal_loop_extent) intent = in +[ surface_temperature ] + standard_name = blackbody_temperature_at_surface + type = real | kind = kind_phys + units = K + 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 + units = fraction dimensions = () intent = in +[ photolysis_wavelength_grid_interfaces ] + standard_name = photolysis_wavelength_grid_interfaces + type = real | kind = kind_phys + units = m + dimensions = (photolysis_wavelength_grid_interface_dimension) + intent = in +[ extraterrestrial_flux ] + standard_name = extraterrestrial_radiation_flux + type = real | kind = kind_phys + units = photons cm-2 s-1 nm-1 + dimensions = (photolysis_wavelength_grid_section_dimension) + intent = in [ standard_gravitational_acceleration ] standard_name = standard_gravitational_acceleration units = m s-2 type = real | kind = kind_phys dimensions = () intent = in +[ cloud_area_fraction ] + standard_name = cloud_area_fraction + units = fraction + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent,vertical_layer_dimension) + intent = in +[ air_pressure_thickness ] + standard_name = air_pressure_thickness + units = Pa + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent,vertical_layer_dimension) + intent = in +[ solar_zenith_angle ] + standard_name = solar_zenith_angle + units = rad + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ earth_sun_distance ] + standard_name = earth_sun_distance + units = AU + type = real | kind = kind_phys + dimensions = () + intent = in [ errmsg ] standard_name = ccpp_error_message units = none diff --git a/schemes/musica/musica_ccpp_namelist.xml b/schemes/musica/musica_ccpp_namelist.xml index 5015f2d3..c36d1b81 100644 --- a/schemes/musica/musica_ccpp_namelist.xml +++ b/schemes/musica/musica_ccpp_namelist.xml @@ -76,7 +76,19 @@ This is the CCPP unit specification of the variable (e.g., m s-1). --> - + + integer + musica_ccpp + musica_ccpp + micm_solver_type + none + + The type of MICM solver to use. + + + 1 + + char*512 musica_ccpp @@ -115,5 +127,5 @@ UNSET_PATH - - + + \ No newline at end of file diff --git a/schemes/musica/musica_ccpp_util.F90 b/schemes/musica/musica_ccpp_util.F90 index 5a4eda9c..c87cce2f 100644 --- a/schemes/musica/musica_ccpp_util.F90 +++ b/schemes/musica/musica_ccpp_util.F90 @@ -2,11 +2,16 @@ ! SPDX-License-Identifier: Apache-2.0 module musica_ccpp_util + use ccpp_kinds, only: kind_phys + implicit none private public :: has_error_occurred + real(kind_phys), parameter, public :: PI = 3.14159265358979323846_kind_phys + real(kind_phys), parameter, public :: DEGREE_TO_RADIAN = PI / 180.0_kind_phys + contains !> @brief Evaluate a MUSICA error for failure and convert to CCPP error data diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 index 5a12a3ee..6e95bde0 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90 @@ -5,33 +5,134 @@ module musica_ccpp_tuvx use ccpp_kinds, only: kind_phys use musica_ccpp_namelist, only: filename_of_tuvx_configuration use musica_ccpp_util, only: has_error_occurred - use musica_tuvx, only: tuvx_t, grid_t, profile_t + use musica_tuvx, only: tuvx_t, grid_t, profile_t, radiator_t use musica_util, only: mappings_t, index_mappings_t implicit none private - public :: tuvx_init, tuvx_run, tuvx_final + public :: tuvx_register, tuvx_init, tuvx_run, tuvx_final + + real(kind_phys), parameter :: MAX_SOLAR_ZENITH_ANGLE = 110.0_kind_phys ! degrees + real(kind_phys), parameter :: MIN_SOLAR_ZENITH_ANGLE = 0.0_kind_phys ! degrees + + 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 = & + 'cloud_liquid_water_mixing_ratio_wrt_moist_air_and_condensed_water' + 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 - 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() +contains - type(index_mappings_t), pointer :: photolysis_rate_constants_mapping => null( ) - integer :: number_of_photolysis_rate_constants = 0 + !> Deallocates TUV-x resources + subroutine reset_tuvx_map_state( grids, profiles, radiators ) + use musica_tuvx, only: grid_map_t, profile_map_t, radiator_map_t -contains + type(grid_map_t), pointer :: grids + type(profile_map_t), pointer :: profiles + type(radiator_map_t), pointer :: radiators + + if (associated( grids )) deallocate( grids ) + if (associated( profiles )) deallocate( profiles ) + if (associated( radiators )) deallocate( radiators ) + + end subroutine reset_tuvx_map_state + + !> This is a helper subroutine created to deallocate objects associated with TUV-x + subroutine cleanup_tuvx_resources() + + 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 + + if (associated( surface_albedo_profile )) then + deallocate( surface_albedo_profile ) + surface_albedo_profile => null() + end if + + if (associated( extraterrestrial_flux_profile )) then + deallocate( extraterrestrial_flux_profile ) + extraterrestrial_flux_profile => null() + end if + + if (associated( cloud_optics )) then + deallocate( cloud_optics ) + cloud_optics => null() + end if + + if (associated( photolysis_rate_constants_mapping )) then + deallocate( photolysis_rate_constants_mapping ) + photolysis_rate_constants_mapping => null() + end if + + end subroutine cleanup_tuvx_resources + + !> Registers constituent properties with the CCPP needed by TUV-x + subroutine tuvx_register(constituent_props, errmsg, errcode) + use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t + use musica_util, only: error_t + + type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituent_props(:) + character(len=512), intent(out) :: errmsg + integer, intent(out) :: errcode + + allocate(constituent_props(1), stat=errcode) + if (errcode /= 0) then + errmsg = "[MUSICA Error] Failed to allocate memory for constituent properties." + return + end if + + ! Register cloud liquid water content needed for cloud optics calculations + call constituent_props(1)%instantiate( & + std_name = CLOUD_LIQUID_WATER_CONTENT_LABEL, & + long_name = CLOUD_LIQUID_WATER_CONTENT_LONG_NAME, & + units = CLOUD_LIQUID_WATER_CONTENT_UNITS, & + vertical_dim = "vertical_layer_dimension", & + default_value = 0.0_kind_phys, & + min_value = 0.0_kind_phys, & + molar_mass = CLOUD_LIQUID_WATER_CONTENT_MOLAR_MASS, & + advected = .true., & + errcode = errcode, & + errmsg = errmsg & + ) + if (errcode /= 0) return + + end subroutine tuvx_register !> Initializes TUV-x subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & wavelength_grid_interfaces, micm_rate_parameter_ordering, & - errmsg, errcode) + constituent_props, errmsg, errcode) + use ccpp_const_utils, only: ccpp_const_get_idx + use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t use musica_tuvx, only: grid_map_t, profile_map_t, radiator_map_t use musica_util, only: error_t, configuration_t use musica_ccpp_namelist, only: filename_of_tuvx_micm_mapping_configuration - use musica_ccpp_tuvx_util, only: tuvx_deallocate + use musica_ccpp_util, only: PI use musica_ccpp_tuvx_height_grid, & only: create_height_grid, height_grid_label, height_grid_unit use musica_ccpp_tuvx_wavelength_grid, & @@ -40,13 +141,19 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & 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) - real(kind_phys), intent(in) :: wavelength_grid_interfaces(:) ! m - type(mappings_t), intent(in) :: micm_rate_parameter_ordering ! index mappings for MICM rate parameters - character(len=512), intent(out) :: errmsg - integer, intent(out) :: errcode + use musica_ccpp_tuvx_extraterrestrial_flux, & + only: create_extraterrestrial_flux_profile, extraterrestrial_flux_label, & + extraterrestrial_flux_unit + use musica_ccpp_tuvx_cloud_optics, & + only: create_cloud_optics_radiator, cloud_optics_label + + integer, intent(in) :: vertical_layer_dimension ! (count) + integer, intent(in) :: vertical_interface_dimension ! (count) + real(kind_phys), intent(in) :: wavelength_grid_interfaces(:) ! m + type(mappings_t), intent(in) :: micm_rate_parameter_ordering ! index mappings for MICM rate parameters + type(ccpp_constituent_prop_ptr_t), intent(in) :: constituent_props(:) + character(len=512), intent(out) :: errmsg + integer, intent(out) :: errcode ! local variables type(grid_map_t), pointer :: grids @@ -56,177 +163,291 @@ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, & type(mappings_t), pointer :: photolysis_rate_constants_ordering type(error_t) :: error + ! Get needed indices in constituents array + call ccpp_const_get_idx(constituent_props, CLOUD_LIQUID_WATER_CONTENT_LABEL, & + index_cloud_liquid_water_content, errmsg, errcode) + if (errcode /= 0) return + if (index_cloud_liquid_water_content == DEFAULT_INDEX_NOT_FOUND) then + errmsg = "[MUSICA Error] Unable to find index for cloud liquid water content." + errcode = 1 + return + end if + grids => grid_map_t( error ) if (has_error_occurred( error, errmsg, errcode )) return height_grid => create_height_grid( vertical_layer_dimension, & vertical_interface_dimension, errmsg, errcode ) if (errcode /= 0) then - deallocate( grids ) + call reset_tuvx_map_state( grids, null(), null() ) return endif 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(), null() ) + call reset_tuvx_map_state( grids, null(), null() ) + call cleanup_tuvx_resources() 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(), null() ) + call reset_tuvx_map_state( grids, null(), null() ) + call cleanup_tuvx_resources() 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(), null() ) + call reset_tuvx_map_state( grids, null(), null() ) + call cleanup_tuvx_resources() 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(), null() ) + call reset_tuvx_map_state( grids, null(), null() ) + call cleanup_tuvx_resources() 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(), null() ) + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() 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, null() ) + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() 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() ) + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() 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 ) + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + end if + + extraterrestrial_flux_profile => create_extraterrestrial_flux_profile( & + wavelength_grid, wavelength_grid_interfaces, errmsg, errcode ) + if (errcode /= 0) then + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + endif + + call profiles%add( extraterrestrial_flux_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 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, surface_albedo_profile ) + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + end if + + cloud_optics => create_cloud_optics_radiator( height_grid, wavelength_grid, & + errmsg, errcode ) + if (errcode /= 0) then + call reset_tuvx_map_state( grids, profiles, radiators ) + call cleanup_tuvx_resources() + return + endif + + call radiators%add( cloud_optics, error ) + if (has_error_occurred( error, errmsg, errcode )) then + call reset_tuvx_map_state( grids, profiles, radiators ) + call cleanup_tuvx_resources() return end if tuvx => tuvx_t( trim(filename_of_tuvx_configuration), grids, profiles, & radiators, error ) if (has_error_occurred( error, errmsg, errcode )) then - call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, & - wavelength_grid, temperature_profile, surface_albedo_profile ) + call reset_tuvx_map_state( grids, profiles, radiators ) + call cleanup_tuvx_resources() return end if - call tuvx_deallocate( grids, profiles, radiators, null(), height_grid, & - wavelength_grid, temperature_profile, surface_albedo_profile ) + call reset_tuvx_map_state( grids, profiles, radiators ) + call cleanup_tuvx_resources() + ! Gets resources associated with TUV-x from 'tuvx' pointer grids => tuvx%get_grids( error ) - if (has_error_occurred( error, errmsg, errcode )) return + if (has_error_occurred( error, errmsg, errcode )) then + deallocate( tuvx ) + tuvx => null() + call cleanup_tuvx_resources() + return + end if 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(), null() ) + deallocate( tuvx ) + tuvx => null() + call reset_tuvx_map_state( grids, null(), null() ) + call cleanup_tuvx_resources() 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(), null() ) + deallocate( tuvx ) + tuvx => null() + call reset_tuvx_map_state( grids, null(), null() ) + call cleanup_tuvx_resources() 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(), null() ) + deallocate( tuvx ) + tuvx => null() + call reset_tuvx_map_state( grids, null(), null() ) + call cleanup_tuvx_resources() 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(), null() ) + deallocate( tuvx ) + tuvx => null() + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() return end if 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() ) + deallocate( tuvx ) + tuvx => null() + call reset_tuvx_map_state( grids, profiles, null() ) + call cleanup_tuvx_resources() + return + end if + + extraterrestrial_flux_profile => & + profiles%get( extraterrestrial_flux_label, extraterrestrial_flux_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 + + radiators => tuvx%get_radiators( 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 + + cloud_optics => radiators%get( cloud_optics_label, error ) + if (has_error_occurred( error, errmsg, errcode )) then + deallocate( tuvx ) + tuvx => null() + call reset_tuvx_map_state( grids, profiles, radiators ) + call cleanup_tuvx_resources() return end if - call tuvx_deallocate( grids, profiles, null(), null(), null(), null(), & - null(), null() ) + call reset_tuvx_map_state( grids, profiles, radiators ) + ! 'photolysis_rate_constants_ordering' is a local variable photolysis_rate_constants_ordering => & tuvx%get_photolysis_rate_constants_ordering( error ) - if (has_error_occurred( error, errmsg, errcode )) return + if (has_error_occurred( error, errmsg, errcode )) then + deallocate( tuvx ) + tuvx => null() + call cleanup_tuvx_resources() + return + end if number_of_photolysis_rate_constants = photolysis_rate_constants_ordering%size() call config%load_from_file( trim(filename_of_tuvx_micm_mapping_configuration), error ) if (has_error_occurred( error, errmsg, errcode )) then + deallocate( tuvx ) + tuvx => null() + call cleanup_tuvx_resources() deallocate( photolysis_rate_constants_ordering ) - photolysis_rate_constants_ordering => null() return end if photolysis_rate_constants_mapping => & index_mappings_t( config, photolysis_rate_constants_ordering, & micm_rate_parameter_ordering, error ) + if (has_error_occurred( error, errmsg, errcode )) then + deallocate( tuvx ) + tuvx => null() + call cleanup_tuvx_resources() + deallocate( photolysis_rate_constants_ordering ) + return + end if + deallocate( photolysis_rate_constants_ordering ) - photolysis_rate_constants_ordering => null() - if (has_error_occurred( error, errmsg, errcode )) return end subroutine tuvx_init !> Calculates photolysis rate constants for the current model conditions - 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, & - 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 + subroutine tuvx_run(temperature, dry_air_density, & + constituents, & + geopotential_height_wrt_surface_at_midpoint, & + geopotential_height_wrt_surface_at_interface, & + surface_geopotential, surface_temperature, & + surface_albedo, & + photolysis_wavelength_grid_interfaces, & + extraterrestrial_flux, & + standard_gravitational_acceleration, & + cloud_area_fraction, & + air_pressure_thickness, & + solar_zenith_angle, & + earth_sun_distance, & + 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, PI + 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 real(kind_phys), intent(in) :: temperature(:,:) ! K (column, layer) real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3 (column, layer) + real(kind_phys), intent(in) :: constituents(:,:,:) ! various (column, layer, constituent) real(kind_phys), intent(in) :: geopotential_height_wrt_surface_at_midpoint(:,:) ! m (column, layer) 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_temperature(:) ! K real(kind_phys), intent(in) :: surface_albedo ! unitless + real(kind_phys), intent(in) :: photolysis_wavelength_grid_interfaces(:) ! nm + real(kind_phys), intent(in) :: extraterrestrial_flux(:) ! photons cm-2 s-1 nm-1 real(kind_phys), intent(in) :: standard_gravitational_acceleration ! m s-2 + real(kind_phys), intent(in) :: cloud_area_fraction(:,:) ! unitless (column, layer) + real(kind_phys), intent(in) :: air_pressure_thickness(:,:) ! Pa (column, layer) + real(kind_phys), intent(in) :: solar_zenith_angle(:) ! radians + real(kind_phys), intent(in) :: earth_sun_distance ! m real(kind_phys), intent(inout) :: rate_parameters(:,:,:) ! various units (column, layer, reaction) character(len=512), intent(out) :: errmsg integer, intent(out) :: errcode @@ -234,12 +455,11 @@ 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) :: reciprocal_of_gravitational_acceleration ! s2 m-1 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) - real(kind_phys) :: solar_zenith_angle ! degrees - real(kind_phys) :: earth_sun_distance ! AU + real(kind_phys) :: reciprocal_of_gravitational_acceleration ! s2 m-1 + real(kind_phys) :: solar_zenith_angle_degrees type(error_t) :: error integer :: i_col, i_level @@ -249,33 +469,49 @@ subroutine tuvx_run(temperature, dry_air_density, & call set_surface_albedo_values( surface_albedo_profile, surface_albedo, errmsg, errcode ) if (errcode /= 0) return + call set_extraterrestrial_flux_values( extraterrestrial_flux_profile, & + photolysis_wavelength_grid_interfaces, & + extraterrestrial_flux, errmsg, errcode ) + 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, & + + ! check if solar zenith angle is within the range to calculate photolysis rate constants + solar_zenith_angle_degrees = solar_zenith_angle(i_col) * 180.0_kind_phys / PI + if (solar_zenith_angle_degrees > MAX_SOLAR_ZENITH_ANGLE .or. & + solar_zenith_angle_degrees < MIN_SOLAR_ZENITH_ANGLE) then + photolysis_rate_constants(:,:) = 0.0_kind_phys + else + 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 ) - if (errcode /= 0) return + if (errcode /= 0) return - call set_temperature_values( temperature_profile, temperature(i_col,:), & + call set_temperature_values( temperature_profile, temperature(i_col,:), & surface_temperature(i_col), errmsg, errcode ) - if (errcode /= 0) return - - ! temporary values until these are available from the host model - solar_zenith_angle = 0.0_kind_phys - earth_sun_distance = 1.0_kind_phys - - ! calculate photolysis rate constants and heating rates - call tuvx%run( solar_zenith_angle, earth_sun_distance, & - photolysis_rate_constants(:,:), heating_rates(:,:), & - error ) - if (has_error_occurred( error, errmsg, errcode )) return - - ! filter out negative photolysis rate constants - photolysis_rate_constants(:,:) = & - max( photolysis_rate_constants(:,:), 0.0_kind_phys ) + if (errcode /= 0) return + + call set_cloud_optics_values( cloud_optics, cloud_area_fraction(i_col,:), & + air_pressure_thickness(i_col,:), & + constituents(i_col,:,index_cloud_liquid_water_content), & + reciprocal_of_gravitational_acceleration, & + errmsg, errcode ) + if (errcode /= 0) return + + ! calculate photolysis rate constants and heating rates + call tuvx%run( solar_zenith_angle(i_col), earth_sun_distance, & + photolysis_rate_constants(:,:), heating_rates(:,:), & + error ) + if (has_error_occurred( error, errmsg, errcode )) return + + ! filter out negative photolysis rate constants + photolysis_rate_constants(:,:) = & + max( photolysis_rate_constants(:,:), 0.0_kind_phys ) + end if ! solar zenith angle check ! map photolysis rate constants to the host model's rate parameters and vertical grid do i_level = 1, size(rate_parameters, dim=2) @@ -295,36 +531,13 @@ subroutine tuvx_final(errmsg, errcode) errmsg = '' errcode = 0 - 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 - - if (associated( surface_albedo_profile )) then - deallocate( surface_albedo_profile ) - surface_albedo_profile => null() - end if + call cleanup_tuvx_resources() if (associated( tuvx )) then deallocate( tuvx ) tuvx => null() end if - if (associated( photolysis_rate_constants_mapping )) then - deallocate( photolysis_rate_constants_mapping ) - photolysis_rate_constants_mapping => null() - end if - end subroutine tuvx_final end module musica_ccpp_tuvx \ No newline at end of file diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_cloud_optics.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_cloud_optics.F90 new file mode 100644 index 00000000..8632176f --- /dev/null +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_cloud_optics.F90 @@ -0,0 +1,127 @@ +! Copyright (C) 2024 National Science Foundation-National Center for Atmospheric Research +! SPDX-License-Identifier: Apache-2.0 +module musica_ccpp_tuvx_cloud_optics + implicit none + + private + public :: create_cloud_optics_radiator, set_cloud_optics_values + + ! This module is used to set the optical properties of clouds in TUV-x. + ! Optical properties are defined as a function of wavelength and height, + ! and include the cloud optical depth, single scattering albedo, + ! and asymmetry parameter. + ! + ! See musica_ccpp_tuvx_height_grid for the definition of the height grid + ! and its mapping to the CAM-SIMA vertical grid. + + !> Label for cloud optical properties in TUV-x + character(len=*), parameter, public :: cloud_optics_label = "clouds" + !> Default value of number of vertical levels + integer, parameter :: DEFAULT_NUM_VERTICAL_LEVELS = 0 + !> Number of vertical levels + integer, protected :: num_vertical_levels = DEFAULT_NUM_VERTICAL_LEVELS + !> Default value of number of wavelength bins + integer, parameter :: DEFAULT_NUM_WAVELENGTH_BINS = 0 + !> Number of wavelength bins + integer, protected :: num_wavelength_bins = DEFAULT_NUM_WAVELENGTH_BINS + +contains + + !> Creates a TUV-x cloud optics radiator from the host-model wavelength grid + function create_cloud_optics_radiator( height_grid, wavelength_grid, & + errmsg, errcode ) result( radiator ) + use musica_ccpp_util, only: has_error_occurred + use musica_tuvx_grid, only: grid_t + use musica_tuvx_radiator, only: radiator_t + use musica_util, only: error_t + + type(grid_t), intent(inout) :: height_grid + type(grid_t), intent(inout) :: wavelength_grid + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + type(radiator_t), pointer :: radiator + + ! local variables + type(error_t) :: error + + num_vertical_levels = height_grid%number_of_sections( error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + num_wavelength_bins = wavelength_grid%number_of_sections( error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + radiator => radiator_t( cloud_optics_label, height_grid, wavelength_grid, & + error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + end function create_cloud_optics_radiator + + !> Sets TUV-x cloud optics values + subroutine set_cloud_optics_values( radiator, cloud_fraction, delta_pressure, & + cloud_liquid_water_content, & + reciprocal_of_gravitational_acceleration, & + errmsg, errcode ) + use ccpp_kinds, only: kind_phys + use musica_ccpp_util, only: has_error_occurred + use musica_tuvx_radiator, only: radiator_t + use musica_util, only: error_t + + type(radiator_t), intent(inout) :: radiator + real(kind_phys), intent(in) :: cloud_fraction(:) ! (unitless) + real(kind_phys), intent(in) :: delta_pressure(:) ! pressure delta about vertical level midpoints (Pa) + real(kind_phys), intent(in) :: cloud_liquid_water_content(:) ! (kg/kg) + real(kind_phys), intent(in) :: reciprocal_of_gravitational_acceleration ! (s^2/m) + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + + ! local variables + type(error_t) :: error + real(kind_phys) :: optical_depth(num_vertical_levels) ! working array for cloud optical depth + real(kind_phys) :: cloud_optical_depth(num_vertical_levels, num_wavelength_bins) + integer :: i_level, size_cloud_fraction + + size_cloud_fraction = size(cloud_fraction) + if ( size_cloud_fraction + 1 /= num_vertical_levels ) then + errmsg = "[MUSICA Error] Invalid size of cloud fraction for TUV-x." + errcode = 1 + return + end if + if ( size(delta_pressure) /= size_cloud_fraction ) then + errmsg = "[MUSICA Error] Invalid size of cloud pressure delta for TUV-x." + errcode = 1 + return + end if + if ( size(cloud_liquid_water_content) /= size_cloud_fraction ) then + errmsg = "[MUSICA Error] Invalid size of cloud liquid water content for TUV-x." + errcode = 1 + return + end if + + ! Estimate cloud optical depth (od) [unitless] from cloud fraction (cf) + ! [unitless] and liquid water content (lwc) [kg kg-1] by first calculating + ! the cloud liquid water path (lwp) [kg m-2]: + ! lwp = 1/g * lwc * dP / cf + ! where g is the gravitational acceleration [m s-2] and dP is the change in + ! pressure across the vertical level [Pa]. + ! The cloud optical depth is then estimated as: + ! od = lwp * 155 * cf^1.5 + ! A constant cloud optical depth is used for all wavelengths. + do i_level = 1, size_cloud_fraction + if ( cloud_fraction(i_level) > 0.0_kind_phys ) then + optical_depth(i_level) = ( reciprocal_of_gravitational_acceleration & + * cloud_liquid_water_content(i_level) * delta_pressure(i_level) & + / cloud_fraction(i_level) ) * 155.0_kind_phys * cloud_fraction(i_level)**1.5_kind_phys + else + optical_depth(i_level) = 0.0_kind_phys + end if + end do + do i_level = 1, size_cloud_fraction + cloud_optical_depth(i_level, :) = optical_depth(size_cloud_fraction-i_level+1) + end do + cloud_optical_depth(num_vertical_levels, :) = 0.0_kind_phys + call radiator%set_optical_depths( cloud_optical_depth, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + end subroutine set_cloud_optics_values + +end module musica_ccpp_tuvx_cloud_optics diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_extraterrestrial_flux.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_extraterrestrial_flux.F90 new file mode 100644 index 00000000..568954e9 --- /dev/null +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_extraterrestrial_flux.F90 @@ -0,0 +1,110 @@ +module musica_ccpp_tuvx_extraterrestrial_flux + use ccpp_kinds, only: kind_phys + + implicit none + + private + public :: create_extraterrestrial_flux_profile, set_extraterrestrial_flux_values + + !> Label for extraterrestrial_flux in TUV-x + character(len=*), parameter, public :: extraterrestrial_flux_label = "extraterrestrial flux" + !> Unit for extraterrestrial_flux in TUV-x + character(len=*), parameter, public :: extraterrestrial_flux_unit = "photon cm-2 s-1" + !> Wavelength grid interface values + real(kind_phys), protected, allocatable :: wavelength_grid_interfaces_(:) ! nm + !> Default value of number of wavelength grid bins + integer, parameter :: DEFAULT_NUM_WAVELENGTH_BINS = 0 + !> Number of wavelength grid bins + integer, protected :: num_wavelength_bins_ = DEFAULT_NUM_WAVELENGTH_BINS + +contains + + !> Creates a TUV-x extraterrestrial flux profile from the host-model wavelength grid + function create_extraterrestrial_flux_profile(wavelength_grid, & + wavelength_grid_interfaces, errmsg, errcode) result( profile ) + use musica_util, only: error_t + use musica_ccpp_util, only: has_error_occurred + use musica_ccpp_tuvx_wavelength_grid, only: m_to_nm + use musica_tuvx_grid, only: grid_t + use musica_tuvx_profile, only: profile_t + + type(grid_t), intent(inout) :: wavelength_grid + real(kind_phys), intent(in) :: wavelength_grid_interfaces(:) ! m + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + type(profile_t), pointer :: profile + + ! local variables + type(error_t) :: error + + profile => profile_t( extraterrestrial_flux_label, extraterrestrial_flux_unit, & + wavelength_grid, error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + num_wavelength_bins_ = wavelength_grid%number_of_sections( error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + + allocate(wavelength_grid_interfaces_( size( wavelength_grid_interfaces ) )) + + wavelength_grid_interfaces_(:) = wavelength_grid_interfaces(:) * m_to_nm + + end function create_extraterrestrial_flux_profile + + !> Sets TUV-x extraterrestrial flux midpoints + ! + ! Extraterrestrial flux is read from data files and interpolated to the + ! TUV-x wavelength grid. CAM extraterrestrial flux values are multiplied by the + ! width of the wavelength bins to get the TUV-x units of photon cm-2 s-1 + ! + ! TUV-x only uses mid-point values for extraterrestrial flux + subroutine set_extraterrestrial_flux_values(profile, photolysis_wavelength_grid_interfaces, & + extraterrestrial_flux, errmsg, errcode) + use musica_ccpp_util, only: has_error_occurred + use musica_tuvx_profile, only: profile_t + use musica_util, only: error_t + use ccpp_kinds, only: kind_phys + use ccpp_tuvx_utils, only: rebin + + type(profile_t), intent(inout) :: profile + real(kind_phys), intent(in) :: photolysis_wavelength_grid_interfaces(:) ! nm + real(kind_phys), intent(in) :: extraterrestrial_flux(:) ! photons cm-2 s-1 nm-1 + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errcode + + ! local variables + type(error_t) :: error + real(kind_phys) :: midpoints(num_wavelength_bins_) + + if (.not. allocated(wavelength_grid_interfaces_)) then + errmsg = "[MUSICA Error] Failed to allocate the TUV-x wavelength grid interface array" + errcode = 1 + return + end if + + if (num_wavelength_bins_ <= DEFAULT_NUM_WAVELENGTH_BINS) then + errmsg = "[MUSICA Error] Invalid size of TUV-x wavelength bins." + errcode = 1 + deallocate( wavelength_grid_interfaces_ ) + return + end if + + ! Regrid normalized flux to TUV-x wavelength grid + call rebin( size(photolysis_wavelength_grid_interfaces) - 1, num_wavelength_bins_, & + photolysis_wavelength_grid_interfaces, wavelength_grid_interfaces_, & + extraterrestrial_flux, midpoints ) + + ! Convert normalized flux to flux on TUV-x wavelength grid + midpoints = midpoints * ( wavelength_grid_interfaces_(2 : num_wavelength_bins_ + 1) & + - wavelength_grid_interfaces_(1 :num_wavelength_bins_) ) + + call profile%set_midpoint_values( midpoints, error) + if ( has_error_occurred( error, errmsg, errcode ) ) then + deallocate( wavelength_grid_interfaces_ ) + return + end if + + deallocate( wavelength_grid_interfaces_ ) + + end subroutine set_extraterrestrial_flux_values + +end module musica_ccpp_tuvx_extraterrestrial_flux \ No newline at end of file diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 index 751d18c6..d2b119b4 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_surface_albedo.F90 @@ -2,17 +2,16 @@ module musica_ccpp_tuvx_surface_albedo implicit none private - public :: create_surface_albedo_profile, set_surface_albedo_values, & - surface_albedo_label, surface_albedo_unit + public :: create_surface_albedo_profile, set_surface_albedo_values !> Label for surface albedo in TUV-x - character(len=*), parameter :: surface_albedo_label = "surface albedo" + character(len=*), parameter, public :: surface_albedo_label = "surface albedo" !> Unit for surface albedo in TUV-x - character(len=*), parameter :: surface_albedo_unit = "none" + character(len=*), parameter, public :: surface_albedo_unit = "none" !> 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 + integer, protected :: num_wavelength_bins_ = DEFAULT_NUM_WAVELENGTH_BINS contains @@ -32,13 +31,13 @@ function create_surface_albedo_profile( wavelength_grid, errmsg, errcode ) & ! local variables type(error_t) :: error - 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 + num_wavelength_bins_ = wavelength_grid%number_of_sections( error ) + if ( has_error_occurred( error, errmsg, errcode ) ) return + end function create_surface_albedo_profile !> Sets TUV-x surface albedo values @@ -58,10 +57,10 @@ subroutine set_surface_albedo_values( profile, host_surface_albedo, & ! local variables type(error_t) :: error - real(kind_phys) :: surface_albedo_interfaces(num_wavelength_bins + 1) + real(kind_phys) :: surface_albedo_interfaces(num_wavelength_bins_ + 1) - if (size(surface_albedo_interfaces) <= DEFAULT_NUM_WAVELENGTH_BINS + 1) then - errmsg = "[MUSICA Error] Invalid size of TUV-x wavelength interfaces." + if (num_wavelength_bins_ <= DEFAULT_NUM_WAVELENGTH_BINS) then + errmsg = "[MUSICA Error] Invalid size of TUV-x wavelength bins." errcode = 1 return end if diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_temperature.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_temperature.F90 index c8ac9cc7..7e9961b7 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_temperature.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_temperature.F90 @@ -2,13 +2,12 @@ module musica_ccpp_tuvx_temperature implicit none private - public :: create_temperature_profile, set_temperature_values, & - temperature_label, temperature_unit + public :: create_temperature_profile, set_temperature_values !> Label for temperature in TUV-x - character(len=*), parameter :: temperature_label = "temperature" + character(len=*), parameter, public :: temperature_label = "temperature" !> Unit for temperature in TUV-x - character(len=*), parameter :: temperature_unit = "K" + character(len=*), parameter, public :: temperature_unit = "K" contains diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 deleted file mode 100644 index 25ac1a28..00000000 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_util.F90 +++ /dev/null @@ -1,55 +0,0 @@ -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, surface_albedo_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 - type(profile_t), pointer :: surface_albedo_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 - - 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 96528d5d..20e34f5c 100644 --- a/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 +++ b/schemes/musica/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 @@ -1,4 +1,5 @@ module musica_ccpp_tuvx_wavelength_grid + use ccpp_kinds, only: kind_phys implicit none @@ -20,6 +21,8 @@ module musica_ccpp_tuvx_wavelength_grid character(len=*), parameter, public :: wavelength_grid_label = "wavelength" !> Unit for wavelength grid in TUV-x character(len=*), parameter, public :: wavelength_grid_unit = "nm" + !> Conversion factor from meters to nanometers (CAM-SIMA -> TUV-x) + real(kind_phys), parameter, public :: m_to_nm = 1.0e9 contains @@ -42,7 +45,7 @@ function create_wavelength_grid( wavelength_grid_interfaces, errmsg, errcode ) & reaL(kind_phys) :: midpoints( size( wavelength_grid_interfaces ) - 1 ) ! nm type(error_t) :: error - interfaces(:) = wavelength_grid_interfaces(:) * 1.0e9 + interfaces(:) = wavelength_grid_interfaces(:) * m_to_nm midpoints(:) = & 0.5 * ( interfaces( 1: size( interfaces ) - 1 ) & + interfaces( 2: size( interfaces ) ) ) diff --git a/schemes/tj2016/tj2016_sfc_pbl_hs.meta b/schemes/tj2016/tj2016_sfc_pbl_hs.meta index 9fdffe53..0cac5ac0 100644 --- a/schemes/tj2016/tj2016_sfc_pbl_hs.meta +++ b/schemes/tj2016/tj2016_sfc_pbl_hs.meta @@ -121,7 +121,7 @@ intent = in [ clat ] standard_name = latitude - units = radians + units = rad type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = in @@ -169,7 +169,7 @@ intent = in [ dudt ] standard_name = tendency_of_eastward_wind - units = m s-1 + units = m s-2 type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = out @@ -181,7 +181,7 @@ intent = in [ dvdt ] standard_name = tendency_of_northward_wind - units = m s-1 + units = m s-2 type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = out diff --git a/schemes/tropopause_find/tropopause_find.meta b/schemes/tropopause_find/tropopause_find.meta index bc218575..bc0978e1 100644 --- a/schemes/tropopause_find/tropopause_find.meta +++ b/schemes/tropopause_find/tropopause_find.meta @@ -65,7 +65,7 @@ intent = in [ lat ] standard_name = latitude - units = radians + units = rad type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = in @@ -112,16 +112,16 @@ dimensions = () intent = in [ tropp_p_loc ] - standard_name = tropopause_air_pressure_from_climatology_dataset + standard_name = tropopause_air_pressure_from_tropopause_climatology_dataset units = Pa type = real | kind = kind_phys - dimensions = (horizontal_dimension, number_of_months_in_year) + dimensions = (horizontal_dimension, number_of_time_slices_in_tropopause_climatology_dataset) intent = in [ tropp_days ] - standard_name = tropopause_calendar_days_from_climatology + standard_name = tropopause_calendar_days_from_tropopause_climatology units = 1 type = real | kind = kind_phys - dimensions = (number_of_months_in_year) + dimensions = (number_of_time_slices_in_tropopause_climatology_dataset) intent = in [ tropLev ] standard_name = tropopause_vertical_layer_index diff --git a/schemes/utilities/physics_tendency_updaters.F90 b/schemes/utilities/physics_tendency_updaters.F90 index 936e7bc8..e58ac366 100644 --- a/schemes/utilities/physics_tendency_updaters.F90 +++ b/schemes/utilities/physics_tendency_updaters.F90 @@ -9,6 +9,7 @@ module physics_tendency_updaters public :: apply_tendency_of_northward_wind_run public :: apply_heating_rate_run public :: apply_tendency_of_air_temperature_run + public :: apply_constituent_tendencies_run CONTAINS @@ -18,7 +19,7 @@ subroutine apply_tendency_of_eastward_wind_run(nz, dudt, u, dudt_total, dt, errcode, errmsg) ! Dummy arguments integer, intent(in) :: nz ! Num vertical layers - real(kind_phys), intent(in) :: dudt(:,:) ! tendency of eastward wind + real(kind_phys), intent(inout) :: dudt(:,:) ! tendency of eastward wind real(kind_phys), intent(inout) :: u(:,:) ! eastward wind real(kind_phys), intent(inout) :: dudt_total(:,:) ! total tendency of eastward wind real(kind_phys), intent(in) :: dt ! physics time step @@ -36,6 +37,8 @@ subroutine apply_tendency_of_eastward_wind_run(nz, dudt, u, dudt_total, dt, dudt_total(:, klev) = dudt_total(:, klev) + dudt(:, klev) end do + dudt = 0.0_kind_phys + end subroutine apply_tendency_of_eastward_wind_run !> \section arg_table_apply_tendency_of_northward_wind_run Argument Table @@ -44,7 +47,7 @@ subroutine apply_tendency_of_northward_wind_run(nz, dvdt, v, dvdt_total, dt, errcode, errmsg) ! Dummy arguments integer, intent(in) :: nz ! Num vertical layers - real(kind_phys), intent(in) :: dvdt(:,:) ! tendency of northward wind + real(kind_phys), intent(inout) :: dvdt(:,:) ! tendency of northward wind real(kind_phys), intent(inout) :: v(:,:) ! northward wind real(kind_phys), intent(inout) :: dvdt_total(:,:) ! total tendency of northward wind real(kind_phys), intent(in) :: dt ! physics time step @@ -62,6 +65,8 @@ subroutine apply_tendency_of_northward_wind_run(nz, dvdt, v, dvdt_total, dt, dvdt_total(:, klev) = dvdt_total(:, klev) + dvdt(:, klev) end do + dvdt = 0.0_kind_phys + end subroutine apply_tendency_of_northward_wind_run !> \section arg_table_apply_heating_rate_run Argument Table @@ -70,7 +75,7 @@ subroutine apply_heating_rate_run(nz, heating_rate, temp, dTdt_total, dt, cpair, errcode, errmsg) ! Dummy arguments integer, intent(in) :: nz ! Num vertical layers - real(kind_phys), intent(in) :: heating_rate(:,:) ! heating rate + real(kind_phys), intent(inout) :: heating_rate(:,:) ! heating rate real(kind_phys), intent(inout) :: temp(:,:) ! air temperature real(kind_phys), intent(inout) :: dTdt_total(:,:) ! total temperature tend. real(kind_phys), intent(in) :: dt ! physics time step @@ -89,6 +94,8 @@ subroutine apply_heating_rate_run(nz, heating_rate, temp, dTdt_total, dt, cpair, dTdt_total(:, klev) = dTdt_total(:, klev) + (heating_rate(:, klev) / cpair(:,klev)) end do + heating_rate = 0.0_kind_phys + end subroutine apply_heating_rate_run !> \section arg_table_apply_tendency_of_air_temperature_run Argument Table @@ -97,7 +104,7 @@ subroutine apply_tendency_of_air_temperature_run(nz, t_tend, temp, dtdT_total, dt, errcode, errmsg) ! Dummy arguments integer, intent(in) :: nz ! Num vertical layers - real(kind_phys), intent(in) :: t_tend(:,:) ! temperature tendency + real(kind_phys), intent(inout) :: t_tend(:,:) ! temperature tendency real(kind_phys), intent(inout) :: temp(:,:) ! air temperature real(kind_phys), intent(inout) :: dTdt_total(:,:) ! total temp. tendency real(kind_phys), intent(in) :: dt ! physics time step @@ -115,6 +122,33 @@ subroutine apply_tendency_of_air_temperature_run(nz, t_tend, temp, dtdT_total, dTdt_total(:, klev) = dTdt_total(:, klev) + t_tend(:, klev) end do + t_tend = 0.0_kind_phys + end subroutine apply_tendency_of_air_temperature_run + !> \section arg_table_apply_constituent_tendencies_run Argument Table + !!! \htmlinclude apply_constituent_tendencies_run.html + subroutine apply_constituent_tendencies_run(nz, const_tend, const, dt, errcode, errmsg) + ! Dummy arguments + integer, intent(in) :: nz ! Num vertical layers + real(kind_phys), intent(inout) :: const_tend(:,:,:) ! constituent tendency array + real(kind_phys), intent(inout) :: const(:,:,:) ! constituent state array + real(kind_phys), intent(in) :: dt ! physics time step + integer, intent(out) :: errcode + character(len=512), intent(out) :: errmsg + + ! Local variables + integer :: klev + + errcode = 0 + errmsg = '' + + do klev = 1, nz + const(:, klev, :) = const(:, klev, :) + (const_tend(:, klev, :) * dt) + end do + + const_tend = 0._kind_phys + + end subroutine apply_constituent_tendencies_run + end module physics_tendency_updaters diff --git a/schemes/utilities/physics_tendency_updaters.meta b/schemes/utilities/physics_tendency_updaters.meta index 8bbe77ef..c2272705 100644 --- a/schemes/utilities/physics_tendency_updaters.meta +++ b/schemes/utilities/physics_tendency_updaters.meta @@ -16,7 +16,7 @@ units = m s-2 type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = in + intent = inout [ u ] standard_name = eastward_wind units = m s-1 @@ -69,7 +69,7 @@ units = m s-2 type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = in + intent = inout [ v ] standard_name = northward_wind units = m s-1 @@ -122,7 +122,7 @@ units = J kg-1 s-1 type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = in + intent = inout [ temp ] standard_name = air_temperature units = K @@ -181,7 +181,7 @@ units = K s-1 type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = in + intent = inout [ temp ] standard_name = air_temperature units = K @@ -215,4 +215,53 @@ type = character | kind = len=512 dimensions = () intent = out +##################################################################### +[ccpp-table-properties] + name = apply_constituent_tendencies + type = scheme +[ccpp-arg-table] + name = apply_constituent_tendencies_run + type = scheme +[ nz ] + standard_name = vertical_layer_dimension + long_name = Number of vertical layers + units = count + type = integer + dimensions = () + intent = in +[ const_tend ] + standard_name = ccpp_constituent_tendencies + long_name = ccpp constituent tendencies + units = none + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension, number_of_ccpp_constituents) + intent = inout +[ const ] + standard_name = ccpp_constituents + long_name = ccpp constituents + units = none + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension, number_of_ccpp_constituents) + intent = inout +[ dt ] + standard_name = timestep_for_physics + long_name = time step + units = s + dimensions = () + type = real | kind = kind_phys + intent = in +[ errcode ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=512 + dimensions = () + intent = out ######################################################### diff --git a/schemes/utilities/state_converters.meta b/schemes/utilities/state_converters.meta index c0e57b15..170cc9e3 100644 --- a/schemes/utilities/state_converters.meta +++ b/schemes/utilities/state_converters.meta @@ -223,7 +223,7 @@ [ exner ] standard_name = dimensionless_exner_function type = real | kind = kind_phys - units = count + units = 1 dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = out [ errmsg ] diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 6a7d1e30..9ac2d8c0 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -10,13 +10,13 @@ set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH};${CMAKE_CURRENT_LIST_DIR}/cmake) set(CMAKE_USER_MAKE_RULES_OVERRIDE ${CMAKE_MODULE_PATH}/SetDefaults.cmake) set(CMAKE_RUNTIME_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}) -# -------------------------------------------------------------------------------- -# NOTE: If 'CCPP_ENABLE_MUSICA_TESTS' on, this is not a stand-alone cmake project anymore. -# MUSICA CCPP wrapper needs MUSICA library and ccpp-framework/src. -# To 'CCPP_ENABLE_MUSICA_TESTS', you either build a cmake project through -# 'docker/Dockerfile.musica' or follow the build instructions in the file. -# The following '$ENV' variables are set by the docker file. -# -------------------------------------------------------------------------------- +# --------------------------------------------------------------------------------------------- +# NOTE: If 'CCPP_ENABLE_MUSICA_TESTS' is enabled, this is no longer a stand-alone CMake project. +# The MUSICA CCPP wrapper requires both the MUSICA library and ccpp-framework/src. +# To enable 'CCPP_ENABLE_MUSICA_TESTS', you can either build a CMake project using +# 'docker/Dockerfile.musica' or follow the build instructions outlined in that file. +# The following '$ENV' variables are configured by the Dockerfile. +# --------------------------------------------------------------------------------------------- option(CCPP_ENABLE_MUSICA_TESTS "Build the MUSICA tests" OFF) option(CCPP_ENABLE_MEMCHECK "Enable memory checks in tests" OFF) @@ -24,12 +24,14 @@ set(INSTALL_GTEST OFF CACHE BOOL "" FORCE) set(BUILD_GMOCK OFF CACHE BOOL "" FORCE) if (CCPP_ENABLE_MUSICA_TESTS) - set(MUSICA_VERSION $ENV{MUSICA_VERSION}) set(MUSICA_SRC_PATH ${CMAKE_SOURCE_DIR}/../schemes/musica) + set(TO_BE_CCPPIZED_SRC_PATH ${CMAKE_SOURCE_DIR}/../to_be_ccppized) set(CCPP_SRC_PATH ${CMAKE_SOURCE_DIR}/$ENV{CCPP_SRC_PATH}) set(CCPP_TEST_SRC_PATH ${CMAKE_SOURCE_DIR}/include) - add_subdirectory(musica) - + include(TestUtils) + include(CTest) enable_testing() + + add_subdirectory(musica) endif() diff --git a/test/docker/Dockerfile.musica b/test/docker/Dockerfile.musica index 83dea485..f83ccdfb 100644 --- a/test/docker/Dockerfile.musica +++ b/test/docker/Dockerfile.musica @@ -6,10 +6,12 @@ FROM ubuntu:22.04 ARG MUSICA_GIT_TAG=326b5119768d5be9654baf96ae3bd6a1b757fdc8 +ARG CAM_SIMA_CHEMISTRY_DATA_TAG=abc7cacbec3d33d5c0ed5bb79a157e93b42c45c0 +ARG BUILD_TYPE=Debug RUN apt update \ && apt install -y sudo \ - && adduser test_user \ + && useradd -m test_user \ && echo "test_user ALL=(root) NOPASSWD: ALL" >> /etc/sudoers.d/test_user \ && chmod 0440 /etc/sudoers.d/test_user @@ -34,7 +36,6 @@ RUN sudo apt update \ libopenmpi-dev \ m4 \ make \ - nlohmann-json3-dev \ openmpi-bin \ python3 \ tree \ @@ -43,25 +44,22 @@ RUN sudo apt update \ zlib1g-dev \ && sudo apt clean -ENV PATH="${PATH}:/usr/lib/openmpi/bin" - ENV FC=mpif90 ENV FFLAGS="-I/usr/include/" # Install MUSICA (MUSICA-C) -RUN git clone https://github.com/NCAR/musica.git - -RUN cd musica \ - && git fetch \ +RUN git clone https://github.com/NCAR/musica.git \ + && cd musica \ && git checkout ${MUSICA_GIT_TAG} \ - && cmake \ - -S . \ - -B build \ - -D CMAKE_BUILD_TYPE=Release \ - -D MUSICA_ENABLE_TESTS=OFF \ - -D MUSICA_BUILD_FORTRAN_INTERFACE=OFF \ - -D MUSICA_ENABLE_MICM=ON \ - -D MUSICA_ENABLE_TUVX=ON \ + && cmake -S . -B build \ + -D CMAKE_BUILD_TYPE=${BUILD_TYPE} \ + -D MUSICA_ENABLE_TESTS=OFF \ + -D MUSICA_BUILD_FORTRAN_INTERFACE=OFF \ + -D MUSICA_ENABLE_MICM=ON \ + -D MUSICA_ENABLE_TUVX=ON \ + -D CMAKE_Fortran_COMPILER=mpif90 \ + -D CMAKE_C_COMPILER=mpicc \ + -D CMAKE_CXX_COMPILER=mpicxx \ && cd build \ && sudo make install @@ -71,10 +69,10 @@ RUN sudo chown -R test_user:test_user atmospheric_physics # Clone the MUSICA chemistry data set repository RUN git clone https://github.com/NCAR/cam-sima-chemistry-data.git \ && cd cam-sima-chemistry-data \ - && git checkout 144d982715f5bd6fa61e76f255a52db1843c55f2 \ + && git checkout ${CAM_SIMA_CHEMISTRY_DATA_TAG} \ && mv mechanisms /home/test_user/atmospheric_physics/schemes/musica/configurations -# Must make ccpp-framework available before building test +# Make CCPP-framework available before building test RUN cd atmospheric_physics/test \ && mkdir lib \ && cd lib \ @@ -82,10 +80,14 @@ RUN cd atmospheric_physics/test \ ENV CCPP_SRC_PATH="lib/ccpp-framework/src" ENV CCPP_FORTRAN_TOOLS_PATH="lib/ccpp-framework/scripts/fortran_tools" +# Make the CCPPStandardNames available +RUN cd atmospheric_physics/test/lib \ + && git clone --depth 1 https://github.com/ESCOMP/CCPPStandardNames.git +ENV CCPP_STD_NAMES_PATH="lib/CCPPStandardNames" + RUN cd atmospheric_physics/test \ - && cmake -S . \ - -B build \ - -D CMAKE_BUILD_TYPE=Debug \ + && cmake -S . -B build \ + -D CMAKE_BUILD_TYPE=${BUILD_TYPE} \ -D CCPP_ENABLE_MUSICA_TESTS=ON \ -D CCPP_ENABLE_MEMCHECK=ON \ && cmake --build ./build diff --git a/test/docker/Dockerfile.musica.no_install b/test/docker/Dockerfile.musica.no_install index e8f70fb0..5baec757 100644 --- a/test/docker/Dockerfile.musica.no_install +++ b/test/docker/Dockerfile.musica.no_install @@ -9,10 +9,12 @@ FROM ubuntu:22.04 ARG MUSICA_GIT_TAG=326b5119768d5be9654baf96ae3bd6a1b757fdc8 +ARG CAM_SIMA_CHEMISTRY_DATA_TAG=abc7cacbec3d33d5c0ed5bb79a157e93b42c45c0 +ARG BUILD_TYPE=Debug RUN apt update \ && apt install -y sudo \ - && adduser test_user \ + && useradd -m test_user \ && echo "test_user ALL=(root) NOPASSWD: ALL" >> /etc/sudoers.d/test_user \ && chmod 0440 /etc/sudoers.d/test_user @@ -37,7 +39,6 @@ RUN sudo apt update \ libopenmpi-dev \ m4 \ make \ - nlohmann-json3-dev \ openmpi-bin \ python3 \ tree \ @@ -46,8 +47,6 @@ RUN sudo apt update \ zlib1g-dev \ && sudo apt clean -ENV PATH="${PATH}:/usr/lib/openmpi/bin" - ENV FC=mpif90 ENV FFLAGS="-I/usr/include/" ENV MUSICA_GIT_TAG=${MUSICA_GIT_TAG} @@ -58,10 +57,10 @@ RUN sudo chown -R test_user:test_user atmospheric_physics # Clone the MUSICA chemistry data set repository RUN git clone https://github.com/NCAR/cam-sima-chemistry-data.git \ && cd cam-sima-chemistry-data \ - && git checkout 144d982715f5bd6fa61e76f255a52db1843c55f2 \ + && git checkout ${CAM_SIMA_CHEMISTRY_DATA_TAG} \ && mv mechanisms /home/test_user/atmospheric_physics/schemes/musica/configurations -# Must make ccpp-framework available before building test +# Make ccpp-framework available before building test RUN cd atmospheric_physics/test \ && mkdir lib \ && cd lib \ @@ -69,10 +68,14 @@ RUN cd atmospheric_physics/test \ ENV CCPP_SRC_PATH="lib/ccpp-framework/src" ENV CCPP_FORTRAN_TOOLS_PATH="lib/ccpp-framework/scripts/fortran_tools" +# Make the CCPPStandardNames available +RUN cd atmospheric_physics/test/lib \ + && git clone --depth 1 https://github.com/ESCOMP/CCPPStandardNames.git +ENV CCPP_STD_NAMES_PATH="lib/CCPPStandardNames" + RUN cd atmospheric_physics/test \ - && cmake -S . \ - -B build \ - -D CMAKE_BUILD_TYPE=Debug \ + && cmake -S . -B build \ + -D CMAKE_BUILD_TYPE=${BUILD_TYPE} \ -D CCPP_ENABLE_MUSICA_TESTS=ON \ -D CCPP_ENABLE_MEMCHECK=ON \ && cmake --build ./build diff --git a/test/musica/CMakeLists.txt b/test/musica/CMakeLists.txt index 74b29693..fbedfba3 100644 --- a/test/musica/CMakeLists.txt +++ b/test/musica/CMakeLists.txt @@ -1,5 +1,4 @@ include(FetchContent) -include(TestUtils) FetchContent_Declare(musica GIT_REPOSITORY https://github.com/NCAR/musica.git @@ -29,6 +28,8 @@ file(GLOB MUSICA_CCPP_SOURCES target_sources(test_musica_api PUBLIC ${MUSICA_CCPP_SOURCES} + ${TO_BE_CCPPIZED_SRC_PATH}/ccpp_tuvx_utils.F90 + ${TO_BE_CCPPIZED_SRC_PATH}/ccpp_const_utils.F90 ${CCPP_SRC_PATH}/ccpp_constituent_prop_mod.F90 ${CCPP_SRC_PATH}/ccpp_hash_table.F90 ${CCPP_SRC_PATH}/ccpp_hashable.F90 @@ -46,8 +47,6 @@ set_target_properties(test_musica_api LINKER_LANGUAGE Fortran ) -include(CTest) - add_test( NAME test_musica_api COMMAND $ @@ -65,7 +64,7 @@ add_custom_target( add_subdirectory(micm) add_subdirectory(tuvx) -# Test metdadata +# Test metadata against the source code find_package(Python REQUIRED) file(MAKE_DIRECTORY ${CMAKE_BINARY_DIR}/metadata_test) @@ -77,7 +76,15 @@ add_custom_target( ) add_test( - NAME test_metadata + NAME test_metadata_against_source_code COMMAND ${Python_EXECUTABLE} ${CMAKE_BINARY_DIR}/../$ENV{CCPP_FORTRAN_TOOLS_PATH}/offline_check_fortran_vs_metadata.py --directory ${CMAKE_BINARY_DIR}/metadata_test +) + +# Test musica scheme metadata against the CCPP standard names +add_test( + NAME test_musica_metadata_against_ccpp_standard_names + COMMAND ${Python_EXECUTABLE} ${CMAKE_BINARY_DIR}/../$ENV{CCPP_STD_NAMES_PATH}/tools/meta_stdname_check.py + --metafile-loc ${CMAKE_BINARY_DIR}/metadata_test/musica_ccpp.meta + --stdname-dict ${CMAKE_BINARY_DIR}/../$ENV{CCPP_STD_NAMES_PATH}/standard_names.xml ) \ No newline at end of file diff --git a/test/musica/musica_ccpp_namelist.F90 b/test/musica/musica_ccpp_namelist.F90 index d27cc754..0dcf3170 100644 --- a/test/musica/musica_ccpp_namelist.F90 +++ b/test/musica/musica_ccpp_namelist.F90 @@ -4,7 +4,8 @@ module musica_ccpp_namelist implicit none private - + + integer, public :: micm_solver_type = 1 character(len=250), public :: filename_of_micm_configuration = 'musica_configurations/chapman/micm/config.json' character(len=250), public :: filename_of_tuvx_configuration = 'musica_configurations/chapman/tuvx/config.json' character(len=250), public :: filename_of_tuvx_micm_mapping_configuration = 'musica_configurations/chapman/tuvx_micm_mapping.json' diff --git a/test/musica/test_musica_api.F90 b/test/musica/test_musica_api.F90 index 61cca2f7..6ed8b0e8 100644 --- a/test/musica/test_musica_api.F90 +++ b/test/musica/test_musica_api.F90 @@ -1,5 +1,6 @@ program run_test_musica_ccpp + use ccpp_kinds, only: kind_phys use musica_ccpp implicit none @@ -7,6 +8,8 @@ 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 + real(kind_phys), parameter :: DEGREE_TO_RADIAN = 3.14159265358979323846_kind_phys / 180.0_kind_phys + call test_chapman() call test_terminator() @@ -133,8 +136,6 @@ end subroutine get_wavelength_edges !> Tests the Chapman chemistry scheme subroutine test_chapman() - use musica_micm, only: Rosenbrock, RosenbrockStandardOrder - 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_micm, only: micm @@ -145,42 +146,50 @@ subroutine test_chapman() implicit none 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. - integer, parameter :: NUM_COLUMNS = 2 - integer, parameter :: NUM_LAYERS = 2 - integer, parameter :: NUM_WAVELENGTH_BINS = 102 - integer :: NUM_GRID_CELLS = NUM_COLUMNS * NUM_LAYERS - integer :: solver_type = Rosenbrock - integer :: errcode - character(len=512) :: errmsg - 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 - 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 - 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 - real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: initial_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, base_conc - character(len=512) :: species_name, units - character(len=:), allocatable :: micm_species_name - logical :: tmp_bool, is_advected - integer :: i, j, k - integer :: N2_index, O2_index, O_index, O1D_index, O3_index - real(kind_phys) :: total_O, total_O_init + integer, parameter :: NUM_COLUMNS = 2 + integer, parameter :: NUM_LAYERS = 2 + integer, parameter :: NUM_WAVELENGTH_BINS = 102 + integer :: NUM_GRID_CELLS = NUM_COLUMNS * NUM_LAYERS + integer :: errcode + character(len=512) :: errmsg + 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 + real(kind_phys), dimension(NUM_COLUMNS) :: surface_geopotential ! m2 s-2 + real(kind_phys), dimension(NUM_COLUMNS) :: surface_temperature ! K + real(kind_phys) :: surface_albedo ! unitless + integer, parameter :: num_photolysis_wavelength_grid_sections = 8 ! (count) + real(kind_phys), dimension(num_photolysis_wavelength_grid_sections+1) :: flux_data_photolysis_wavelength_interfaces ! nm + real(kind_phys), dimension(num_photolysis_wavelength_grid_sections) :: extraterrestrial_flux ! photons cm-2 s-1 nm-1 + 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 + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: dry_air_density ! kg m-3 + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: cloud_area_fraction ! unitless + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: air_pressure_thickness ! Pa + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS, & + NUM_SPECIES+NUM_TUVX_CONSTITUENTS) :: constituents ! kg kg-1 + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS, & + NUM_SPECIES+NUM_TUVX_CONSTITUENTS) :: initial_constituents ! kg kg-1 + real(kind_phys), dimension(NUM_COLUMNS) :: solar_zenith_angle ! radians + real(kind_phys) :: earth_sun_distance ! AU + 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, base_conc + character(len=512) :: species_name, units + character(len=:), allocatable :: micm_species_name + logical :: tmp_bool, is_advected + integer :: i, j, k + integer :: N2_index, O2_index, O_index, O1D_index, O3_index + real(kind_phys) :: total_O, total_O_init call get_wavelength_edges(photolysis_wavelength_grid_interfaces) - solver_type = Rosenbrock time_step = 60._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 /) @@ -196,18 +205,30 @@ subroutine test_chapman() pressure(:,2) = (/ 8000.04_kind_phys, 9000.04_kind_phys /) dry_air_density(:,1) = (/ 3.5_kind_phys, 4.5_kind_phys /) dry_air_density(:,2) = (/ 5.5_kind_phys, 6.5_kind_phys /) + flux_data_photolysis_wavelength_interfaces(:) = & + (/ 200.0_kind_phys, 210.0_kind_phys, 220.0_kind_phys, 230.0_kind_phys, & + 240.0_kind_phys, 250.0_kind_phys, 260.0_kind_phys, 270.0_kind_phys, 280.0_kind_phys /) + extraterrestrial_flux(:) = & + (/ 1.5e13_kind_phys, 1.5e13_kind_phys, 1.4e13_kind_phys, 1.4e13_kind_phys, & + 1.3e13_kind_phys, 1.2e13_kind_phys, 1.1e13_kind_phys, 1.0e13_kind_phys /) + cloud_area_fraction(:,1) = (/ 0.1_kind_phys, 0.2_kind_phys /) + cloud_area_fraction(:,2) = (/ 0.3_kind_phys, 0.4_kind_phys /) + air_pressure_thickness(:,1) = (/ 900.0_kind_phys, 905.0_kind_phys /) + air_pressure_thickness(:,2) = (/ 910.0_kind_phys, 915.0_kind_phys /) + solar_zenith_angle = (/ 0.0_kind_phys, 2.1_kind_phys /) + earth_sun_distance = 1.04_kind_phys 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' - call musica_ccpp_register(solver_type, NUM_GRID_CELLS, constituent_props, errmsg, errcode) + call musica_ccpp_register(constituent_props, errmsg, errcode) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 endif ASSERT(allocated(constituent_props)) - ASSERT(size(constituent_props) == NUM_SPECIES) + ASSERT(size(constituent_props) == NUM_SPECIES+NUM_TUVX_CONSTITUENTS) do i = 1, size(constituent_props) ASSERT(constituent_props(i)%is_instantiated(errcode, errmsg)) ASSERT(errcode == 0) @@ -221,7 +242,9 @@ subroutine test_chapman() (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) + (trim(species_name) == "N2" .and. molar_mass == 0.0280134_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) ASSERT(tmp_bool) call constituent_props(i)%units(units, errcode, errmsg) if (errcode /= 0) then @@ -241,8 +264,8 @@ subroutine test_chapman() call constituent_props_ptr(i)%set(const_prop, errcode, errmsg) end do - call musica_ccpp_init(NUM_LAYERS, NUM_LAYERS+1, photolysis_wavelength_grid_interfaces, & - errmsg, errcode) + call musica_ccpp_init(NUM_COLUMNS, NUM_LAYERS, NUM_LAYERS+1, photolysis_wavelength_grid_interfaces, & + constituent_props_ptr, errmsg, errcode) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 @@ -275,6 +298,12 @@ subroutine test_chapman() end do end do end do + ! set initial cloud liquid water mixing ratio to ~1e-3 kg kg-1 + do j = 1, NUM_COLUMNS + do k = 1, NUM_LAYERS + constituents(j,k,NUM_SPECIES+1) = 1.0e-3_kind_phys * (1.0 + 0.1 * (j-1) + 0.01 * (k-1)) + end do + end do initial_constituents(:,:,:) = constituents(:,:,:) write(*,*) "[MUSICA INFO] Initial Time Step" @@ -286,11 +315,14 @@ subroutine test_chapman() 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, surface_albedo, 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_geopotential, & + surface_temperature, surface_albedo, & + flux_data_photolysis_wavelength_interfaces, extraterrestrial_flux, & + standard_gravitational_acceleration, cloud_area_fraction, & + air_pressure_thickness, solar_zenith_angle, earth_sun_distance, errmsg, & + errcode ) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 @@ -322,9 +354,16 @@ subroutine test_chapman() constituents(i,j,O2_index) + constituents(i,j,O3_index) total_O_init = initial_constituents(i,j,O_index) + initial_constituents(i,j,O1D_index) + & initial_constituents(i,j,O2_index) + initial_constituents(i,j,O3_index) + ! cloud liquid water mixing ratio should be unchanged + ASSERT_NEAR(constituents(i,j,NUM_SPECIES+1), initial_constituents(i,j,NUM_SPECIES+1), 1.0e-13) ASSERT_NEAR(total_O, total_O_init, 1.0e-13) end do end do + do j = 1, NUM_LAYERS + ! O and O1D should be lower in the nighttime column + ASSERT(constituents(2,j,O_index) < constituents(1,j,O_index)) + ASSERT(constituents(2,j,O1D_index) < constituents(1,j,O1D_index)) + end do deallocate(constituent_props_ptr) @@ -332,8 +371,6 @@ end subroutine test_chapman !> Tests the simple Terminator chemistry scheme subroutine test_terminator() - use musica_micm, only: Rosenbrock, RosenbrockStandardOrder - 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_micm, only: micm @@ -344,42 +381,50 @@ subroutine test_terminator() implicit none integer, parameter :: NUM_SPECIES = 2 + 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. - integer, parameter :: NUM_COLUMNS = 2 - integer, parameter :: NUM_LAYERS = 2 - integer, parameter :: NUM_WAVELENGTH_BINS = 102 - integer :: NUM_GRID_CELLS = NUM_COLUMNS * NUM_LAYERS - integer :: solver_type = Rosenbrock - integer :: errcode - character(len=512) :: errmsg - 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 - 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 - 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 - real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: initial_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, base_conc - character(len=512) :: species_name, units - character(len=:), allocatable :: micm_species_name - logical :: tmp_bool, is_advected - integer :: i, j, k - integer :: Cl_index, Cl2_index - real(kind_phys) :: total_Cl, total_Cl_init + integer, parameter :: NUM_COLUMNS = 2 + integer, parameter :: NUM_LAYERS = 2 + integer, parameter :: NUM_WAVELENGTH_BINS = 102 + integer :: NUM_GRID_CELLS = NUM_COLUMNS * NUM_LAYERS + integer :: errcode + character(len=512) :: errmsg + 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 + real(kind_phys), dimension(NUM_COLUMNS) :: surface_geopotential ! m2 s-2 + real(kind_phys), dimension(NUM_COLUMNS) :: surface_temperature ! K + real(kind_phys) :: surface_albedo ! unitless + integer, parameter :: num_photolysis_wavelength_grid_sections = 8 ! (count) + real(kind_phys), dimension(num_photolysis_wavelength_grid_sections+1) :: flux_data_photolysis_wavelength_interfaces ! nm + real(kind_phys), dimension(num_photolysis_wavelength_grid_sections) :: extraterrestrial_flux ! photons cm-2 s-1 nm-1 + 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 + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: dry_air_density ! kg m-3 + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: cloud_area_fraction ! unitless + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: air_pressure_thickness ! Pa + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS, & + NUM_SPECIES+NUM_TUVX_CONSTITUENTS) :: constituents ! kg kg-1 + real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS, & + NUM_SPECIES+NUM_TUVX_CONSTITUENTS) :: initial_constituents ! kg kg-1 + real(kind_phys), dimension(NUM_COLUMNS) :: solar_zenith_angle ! radians + real(kind_phys) :: earth_sun_distance ! AU + 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, base_conc + character(len=512) :: species_name, units + character(len=:), allocatable :: micm_species_name + logical :: tmp_bool, is_advected + integer :: i, j, k + integer :: Cl_index, Cl2_index + real(kind_phys) :: total_Cl, total_Cl_init call get_wavelength_edges(photolysis_wavelength_grid_interfaces) - solver_type = Rosenbrock time_step = 60._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 /) @@ -395,18 +440,30 @@ subroutine test_terminator() pressure(:,2) = (/ 8000.04_kind_phys, 9000.04_kind_phys /) dry_air_density(:,1) = (/ 3.5_kind_phys, 4.5_kind_phys /) dry_air_density(:,2) = (/ 5.5_kind_phys, 6.5_kind_phys /) + flux_data_photolysis_wavelength_interfaces(:) = & + (/ 200.0_kind_phys, 210.0_kind_phys, 220.0_kind_phys, 230.0_kind_phys, & + 240.0_kind_phys, 250.0_kind_phys, 260.0_kind_phys, 270.0_kind_phys, 280.0_kind_phys /) + extraterrestrial_flux(:) = & + (/ 1.5e13_kind_phys, 1.5e13_kind_phys, 1.4e13_kind_phys, 1.4e13_kind_phys, & + 1.3e13_kind_phys, 1.2e13_kind_phys, 1.1e13_kind_phys, 1.0e13_kind_phys /) + cloud_area_fraction(:,1) = (/ 0.1_kind_phys, 0.2_kind_phys /) + cloud_area_fraction(:,2) = (/ 0.3_kind_phys, 0.4_kind_phys /) + air_pressure_thickness(:,1) = (/ 900.0_kind_phys, 905.0_kind_phys /) + air_pressure_thickness(:,2) = (/ 910.0_kind_phys, 915.0_kind_phys /) + solar_zenith_angle = (/ 0.0_kind_phys, 2.1_kind_phys /) + earth_sun_distance = 1.04_kind_phys filename_of_micm_configuration = 'musica_configurations/terminator/micm/config.json' filename_of_tuvx_configuration = 'musica_configurations/terminator/tuvx/config.json' filename_of_tuvx_micm_mapping_configuration = 'musica_configurations/terminator/tuvx_micm_mapping.json' - call musica_ccpp_register(solver_type, NUM_GRID_CELLS, constituent_props, errmsg, errcode) + call musica_ccpp_register(constituent_props, errmsg, errcode) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 endif ASSERT(allocated(constituent_props)) - ASSERT(size(constituent_props) == NUM_SPECIES) + ASSERT(size(constituent_props) == NUM_SPECIES+NUM_TUVX_CONSTITUENTS) do i = 1, size(constituent_props) ASSERT(constituent_props(i)%is_instantiated(errcode, errmsg)) ASSERT(errcode == 0) @@ -417,7 +474,9 @@ subroutine test_terminator() call constituent_props(i)%is_advected(is_advected, errcode, errmsg) 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) + (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) ASSERT(tmp_bool) call constituent_props(i)%units(units, errcode, errmsg) if (errcode /= 0) then @@ -437,8 +496,8 @@ subroutine test_terminator() call constituent_props_ptr(i)%set(const_prop, errcode, errmsg) end do - call musica_ccpp_init(NUM_LAYERS, NUM_LAYERS+1, photolysis_wavelength_grid_interfaces, & - errmsg, errcode) + call musica_ccpp_init(NUM_COLUMNS, NUM_LAYERS, NUM_LAYERS+1, photolysis_wavelength_grid_interfaces, & + constituent_props_ptr, errmsg, errcode) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 @@ -462,6 +521,12 @@ subroutine test_terminator() end do end do end do + ! set initial cloud liquid water mixing ratio to ~1e-3 kg kg-1 + do j = 1, NUM_COLUMNS + do k = 1, NUM_LAYERS + constituents(j,k,NUM_SPECIES+1) = 1.0e-3_kind_phys * (1.0 + 0.1 * (j-1) + 0.01 * (k-1)) + end do + end do initial_constituents(:,:,:) = constituents(:,:,:) write(*,*) "[MUSICA INFO] Initial Time Step" @@ -473,11 +538,14 @@ subroutine test_terminator() 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, surface_albedo, 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_geopotential, & + surface_temperature, surface_albedo, & + flux_data_photolysis_wavelength_interfaces, extraterrestrial_flux, & + standard_gravitational_acceleration, cloud_area_fraction, & + air_pressure_thickness, solar_zenith_angle, earth_sun_distance, errmsg, & + errcode ) if (errcode /= 0) then write(*,*) trim(errmsg) stop 3 @@ -504,8 +572,16 @@ subroutine test_terminator() total_Cl = constituents(i,j,Cl_index) + constituents(i,j,Cl2_index) total_Cl_init = initial_constituents(i,j,Cl_index) + initial_constituents(i,j,Cl2_index) ASSERT_NEAR(total_Cl, total_Cl_init, 1.0e-13) + ! cloud liquid water should be unchanged + ASSERT_NEAR(constituents(i,j,NUM_SPECIES+1), initial_constituents(i,j,NUM_SPECIES+1), 1.0e-13) end do end do + do j = 1, NUM_LAYERS + ! Cl should be lower in the nighttime column + ASSERT(constituents(2,j,Cl_index) < constituents(1,j,Cl_index)) + ! Cl2 should be higher in the nighttime column + ASSERT(constituents(2,j,Cl2_index) > constituents(1,j,Cl2_index)) + end do deallocate(constituent_props_ptr) diff --git a/test/musica/tuvx/CMakeLists.txt b/test/musica/tuvx/CMakeLists.txt index 2c8623fd..10024759 100644 --- a/test/musica/tuvx/CMakeLists.txt +++ b/test/musica/tuvx/CMakeLists.txt @@ -110,4 +110,65 @@ add_test( WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY} ) -add_memory_check_test(test_tuvx_surface_albedo $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) \ No newline at end of file +add_memory_check_test(test_tuvx_surface_albedo $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) + +# Extraterrestrial flux +add_executable(test_tuvx_extraterrestrial_flux test_tuvx_extraterrestrial_flux.F90) + +target_sources(test_tuvx_extraterrestrial_flux + PUBLIC + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_extraterrestrial_flux.F90 + ${MUSICA_SRC_PATH}/musica_ccpp_util.F90 + ${TO_BE_CCPPIZED_SRC_PATH}/ccpp_tuvx_utils.F90 + ${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90 +) + +target_link_libraries(test_tuvx_extraterrestrial_flux + PRIVATE + musica-fortran +) + +set_target_properties(test_tuvx_extraterrestrial_flux + PROPERTIES + LINKER_LANGUAGE Fortran +) + +add_test( + NAME test_tuvx_extraterrestrial_flux + COMMAND $ + WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY} +) + +add_memory_check_test(test_tuvx_extraterrestrial_flux $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) + +# Cloud optics +add_executable(test_tuvx_cloud_optics test_tuvx_cloud_optics.F90) + +target_sources(test_tuvx_cloud_optics + PUBLIC + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_height_grid.F90 + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_wavelength_grid.F90 + ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_cloud_optics.F90 + ${MUSICA_SRC_PATH}/musica_ccpp_util.F90 + ${TO_BE_CCPPIZED_SRC_PATH}/ccpp_tuvx_utils.F90 + ${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90 +) + +target_link_libraries(test_tuvx_cloud_optics + PRIVATE + musica::musica-fortran +) + +set_target_properties(test_tuvx_cloud_optics + PROPERTIES + LINKER_LANGUAGE Fortran +) + +add_test( + NAME test_tuvx_cloud_optics + COMMAND $ + WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY} +) + +add_memory_check_test(test_tuvx_cloud_optics $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}) diff --git a/test/musica/tuvx/configs/ts1_tsmlt.json b/test/musica/tuvx/configs/ts1_tsmlt.json deleted file mode 100644 index b2a2f986..00000000 --- a/test/musica/tuvx/configs/ts1_tsmlt.json +++ /dev/null @@ -1,2061 +0,0 @@ -{ - "__description": "TUV-x configuration for the MOZART-TS1 and MOZART-TSMLT chemical mechanisms", - "O2 absorption" : { - "cross section parameters file": "data/cross_sections/O2_parameters.txt" - }, - "grids": [ - { - "name": "time", - "type": "from config file", - "units": "hours", - "values": [ 12.0, 14.0 ] - } - ], - "profiles": [ - { - "name": "O3", - "type": "O3", - "units": "molecule cm-3", - "file path": "data/profiles/atmosphere/ussa.ozone" - }, - { - "name": "air", - "type": "air", - "units": "molecule cm-3", - "file path": "data/profiles/atmosphere/ussa.dens" - }, - { - "name": "O2", - "type": "O2", - "units": "molecule cm-3", - "file path": "data/profiles/atmosphere/ussa.dens" - }, - { - "name": "solar zenith angle", - "type": "solar zenith angle", - "units": "degrees", - "year" : 2002, - "month": 3, - "day": 21, - "longitude": 0.0, - "latitude": 0.0 - }, - { - "name": "Earth-Sun distance", - "type": "Earth-Sun distance", - "units": "AU", - "year" : 2002, - "month": 3, - "day": 21 - }, - { - "name": "extraterrestrial flux", - "enable diagnostics" : true, - "type": "extraterrestrial flux", - "units": "photon cm-2 s-1", - "file path": ["data/profiles/solar/susim_hi.flx", - "data/profiles/solar/atlas3_1994_317_a.dat", - "data/profiles/solar/sao2010.solref.converted", - "data/profiles/solar/neckel.flx"], - "interpolator": ["","","","fractional target"] - } - ], - "radiative transfer": { - "__output": true, - "solver" : { - "type" : "delta eddington" - }, - "cross sections": [ - { - "name": "air", - "type": "air" - }, - { - "name": "O3", - "netcdf files": [ - { "file path": "data/cross_sections/O3_1.nc" }, - { "file path": "data/cross_sections/O3_2.nc" }, - { "file path": "data/cross_sections/O3_3.nc" }, - { "file path": "data/cross_sections/O3_4.nc" } - ], - "type": "O3" - }, - { - "name": "O2", - "netcdf files": [ - { - "file path": "data/cross_sections/O2_1.nc", - "lower extrapolation": { "type": "boundary" } - } - ], - "type": "base" - } - ], - "radiators": [ - { - "enable diagnostics" : true, - "name": "air", - "type": "base", - "treat as air": true, - "cross section": "air", - "vertical profile": "air", - "vertical profile units": "molecule cm-3" - }, - { - "enable diagnostics" : true, - "name": "O2", - "type": "base", - "cross section": "O2", - "vertical profile": "O2", - "vertical profile units": "molecule cm-3" - }, - { - "enable diagnostics" : true, - "name": "O3", - "type": "base", - "cross section": "O3", - "vertical profile": "O3", - "vertical profile units": "molecule cm-3" - }, - { - "enable diagnostics" : true, - "name": "aerosols", - "type": "aerosol", - "optical depths": [2.40e-01, 1.06e-01, 4.56e-02, 1.91e-02, 1.01e-02, 7.63e-03, - 5.38e-03, 5.00e-03, 5.15e-03, 4.94e-03, 4.82e-03, 4.51e-03, - 4.74e-03, 4.37e-03, 4.28e-03, 4.03e-03, 3.83e-03, 3.78e-03, - 3.88e-03, 3.08e-03, 2.26e-03, 1.64e-03, 1.23e-03, 9.45e-04, - 7.49e-04, 6.30e-04, 5.50e-04, 4.21e-04, 3.22e-04, 2.48e-04, - 1.90e-04, 1.45e-04, 1.11e-04, 8.51e-05, 6.52e-05, 5.00e-05, - 3.83e-05, 2.93e-05, 2.25e-05, 1.72e-05, 1.32e-05, 1.01e-05, - 7.72e-06, 5.91e-06, 4.53e-06, 3.46e-06, 2.66e-06, 2.04e-06, - 1.56e-06, 1.19e-06, 9.14e-07], - "single scattering albedo": 0.99, - "asymmetry factor": 0.61, - "550 nm optical depth": 0.235 - } - ] - }, - "photolysis": { - "reactions": [ - { - "name": "jo2_a", - "__reaction": "O2 + hv -> O + O1D", - "cross section": { - "apply O2 bands": true, - "netcdf files": [ - { - "file path": "data/cross_sections/O2_1.nc", - "lower extrapolation": { "type": "boundary" }, - "interpolator": { "type": "fractional target" } - } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 0, - "override bands": [ - { - "band": "lyman-alpha", - "value": 0.53 - }, - { - "band": "schumann-runge continuum", - "value": 1.0 - } - ] - }, - "heating" : { - "energy term": 175.05 - } - }, - { - "name": "jo2_b", - "__reaction": "O2 + hv -> O + O", - "cross section": { - "apply O2 bands": true, - "netcdf files": [ - { - "file path": "data/cross_sections/O2_1.nc", - "lower extrapolation": { "type": "boundary" }, - "interpolator": { "type": "fractional target" } - } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0, - "override bands": [ - { - "band": "lyman-alpha", - "value": 0.47 - }, - { - "band": "schumann-runge continuum", - "value": 0.0 - } - ] - }, - "heating" : { - "energy term": 242.37 - } - }, - { - "name": "jo3_a", - "__reaction": "O3 + hv -> O2 + O(1D)", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/O3_1.nc" }, - { "file path": "data/cross_sections/O3_2.nc" }, - { "file path": "data/cross_sections/O3_3.nc" }, - { "file path": "data/cross_sections/O3_4.nc" } - ], - "type": "O3" - }, - "quantum yield": { - "type": "O3+hv->O2+O(1D)" - }, - "heating" : { - "energy term": 310.32 - } - }, - { - "name": "jo3_b", - "__reaction": "O3 + hv -> O2 + O(3P)", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/O3_1.nc" }, - { "file path": "data/cross_sections/O3_2.nc" }, - { "file path": "data/cross_sections/O3_3.nc" }, - { "file path": "data/cross_sections/O3_4.nc" } - ], - "type": "O3" - }, - "quantum yield": { - "type": "O3+hv->O2+O(3P)" - }, - "heating" : { - "energy term": 1179.87 - } - }, - { - "name": "jn2o", - "__reaction": "N2O + hv -> N2 + O(1D)", - "cross section": { - "type": "N2O+hv->N2+O(1D)" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jno2", - "__reaction": "NO2 + hv -> NO + O(3P)", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/NO2_1.nc" } - ], - "type": "NO2 tint" - }, - "quantum yield": { - "netcdf files": ["data/quantum_yields/NO2_1.nc"], - "type": "NO2 tint", - "lower extrapolation": { "type": "boundary" } - } - }, - { - "name": "jn2o5_a", - "__reaction": "N2O5 + hv -> NO2 + NO3", - "cross section": { - "type":"temperature based", - "netcdf file": "data/cross_sections/N2O5_JPL06.nc", - "parameterization": { - "type": "HARWOOD", - "aa": [ -18.27, -18.42, -18.59, -18.72, -18.84, - -18.90, -18.93, -18.87, -18.77, -18.71, - -18.31, -18.14, -18.01, -18.42, -18.59, - -18.13 ], - "bb": [ -91.0, -104.0, -112.0, -135.0, -170.0, - -226.0, -294.0, -388.0, -492.0, -583.0, - -770.0, -885.0, -992.0, -949.0, -966.0, - -1160.0 ], - "base temperature": 0.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "minimum wavelength": 260.0, - "maximum wavelength": 410.0, - "temperature ranges": [ - { - "maximum": 199.999999999999, - "fixed value": 200 - }, - { - "minimum": 200, - "maximum": 295 - }, - { - "minimum": 295.00000000001, - "fixed value": 295.0 - } - ] - }, - "parameterization wavelength grid": { - "name": "custom wavelengths", - "type": "from config file", - "units": "nm", - "values": [ - 255.0, 265.0, 275.0, 285.0, 295.0, 305.0, - 315.0, 325.0, 335.0, 345.0, 355.0, 365.0, - 375.0, 385.0, 395.0, 405.0, 415.0 - ] - } - }, - "quantum yield": { - "type": "Taylor series", - "constant value": 0.0, - "coefficients": [ -2.832441, 0.012809638 ], - "override bands": [ - { - "band": "range", - "minimum wavelength": 300.0, - "value": 1.0 - } - ] - } - }, - { - "name": "jn2o5_b", - "__reaction": "N2O5 + hv -> NO + O + NO3", - "cross section": { - "type":"temperature based", - "netcdf file": "data/cross_sections/N2O5_JPL06.nc", - "parameterization": { - "type": "HARWOOD", - "aa": [ -18.27, -18.42, -18.59, -18.72, -18.84, - -18.90, -18.93, -18.87, -18.77, -18.71, - -18.31, -18.14, -18.01, -18.42, -18.59, - -18.13 ], - "bb": [ -91.0, -104.0, -112.0, -135.0, -170.0, - -226.0, -294.0, -388.0, -492.0, -583.0, - -770.0, -885.0, -992.0, -949.0, -966.0, - -1160.0 ], - "base temperature": 0.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "minimum wavelength": 260.0, - "maximum wavelength": 410.0, - "temperature ranges": [ - { - "maximum": 199.999999999999, - "fixed value": 200 - }, - { - "minimum": 200, - "maximum": 295 - }, - { - "minimum": 295.00000000001, - "fixed value": 295.0 - } - ] - }, - "parameterization wavelength grid": { - "name": "custom wavelengths", - "type": "from config file", - "units": "nm", - "values": [ - 255.0, 265.0, 275.0, 285.0, 295.0, 305.0, - 315.0, 325.0, 335.0, 345.0, 355.0, 365.0, - 375.0, 385.0, 395.0, 405.0, 415.0 - ] - } - }, - "quantum yield": { - "type": "Taylor series", - "constant value": 0.0, - "coefficients": [ 3.832441, -0.012809638 ], - "override bands": [ - { - "band": "range", - "minimum wavelength": 300.0, - "value": 0.0 - } - ] - } - }, - { - "name": "jhno3", - "__reaction": "HNO3 + hv -> OH + NO2", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/HNO3_JPL06.nc" } - ], - "type": "HNO3+hv->OH+NO2" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jno3_a", - "__reaction": "NO3 + hv -> NO2 + O(3P)", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/NO3_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "netcdf files": [ - "data/quantum_yields/NO3-NO2+O(3P)_1.nc" - ], - "type": "tint", - "lower extrapolation": { - "type": "constant", - "value": 1.0 - } - } - }, - { - "name": "jno3_b", - "__reaction": "NO3 + hv -> NO + O2", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/NO3_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "netcdf files": [ - "data/quantum_yields/NO3-NO+O2_1.nc" - ], - "type": "tint" - } - }, - { - "name": "jch3ooh", - "__reaction": "CH3OOH + hv -> CH3O + OH", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CH3OOH_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jch2o_a", - "__reaction": "CH2O + hv -> H + HCO", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CH2O_1.nc" } - ], - "type": "CH2O" - }, - "quantum yield": { - "netcdf files": [ - "data/quantum_yields/CH2O_1.nc" - ], - "type": "base", - "lower extrapolation": { - "type": "boundary" - } - } - }, - { - "name": "jch2o_b", - "__reaction": "CH2O + hv -> H2 + CO", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CH2O_1.nc" } - ], - "type": "CH2O" - }, - "quantum yield": { - "netcdf files": [ - "data/quantum_yields/CH2O_1.nc" - ], - "type": "CH2O", - "lower extrapolation": { - "type": "boundary" - } - } - }, - { - "name": "jh2o2", - "__reaction": "H2O2 + hv -> OH + OH", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/H2O2_1.nc" } - ], - "type": "H2O2+hv->OH+OH" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jch3cho", - "__reaction": "CH3CHO + hv -> CH3 + HCO", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CH3CHO_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "netcdf files": [ - "data/quantum_yields/CH3CHO_1.nc" - ], - "type": "CH3CHO+hv->CH3+HCO" - } - }, - { - "name": "jpan", - "__reaction": "PAN + hv -> 0.6*CH3CO3 + 0.6*NO2 + 0.4*CH3O2 + 0.4*NO3 + 0.4*CO2", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/PAN_1.nc" } - ], - "type": "CH3ONO2+hv->CH3O+NO2" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jmvk", - "__reaction": "MVK + hv -> 0.7*C3H6 + 0.7*CO + 0.3*CH3O2 + 0.3*CH3CO3", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/MVK_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "MVK+hv->Products" - } - }, - { - "name": "jacet", - "__reaction": "CH3COCH3 + hv -> CH3CO + CH3", - "cross section": { - "type": "temperature based", - "parameterization": { - "type": "TAYLOR_SERIES", - "netcdf file": { - "file path": "data/cross_sections/ACETONE_JPL06.nc" - }, - "base temperature": 0.0, - "temperature ranges": [ - { - "maximum": 234.999999999999, - "fixed value": 235.0 - }, - { - "minimum": 235.0, - "maximum": 298.0 - }, - { - "minimum": 298.00000000001, - "fixed value": 298.0 - } - ] - } - }, - "quantum yield": { - "type": "CH3COCH3+hv->CH3CO+CH3", - "branch": "CO+CH3CO", - "low wavelength value": 1, - "minimum temperature": 218, - "maximum temperature": 295 - } - }, - { - "name": "jmgly", - "__reaction": "CH3COCHO + hv -> CH3CO3 + CO + HO2", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CH3COCHO_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "CH3COCHO+hv->CH3CO+HCO" - } - }, - { - "name": "jglyald", - "__reaction": "GLYALD + hv -> 2*HO2 + CO + CH2O", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/HOCH2CHO_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 0.5 - } - }, - { - "name": "jbrcl", - "__reaction": "BrCl + hv -> Br + Cl", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/BrCl_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jbro", - "__reaction": "BrO + hv -> Br + O", - "cross section": { - "netcdf files": [ - { - "file path": "data/cross_sections/BRO_JPL06.nc" - } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jbrono2_a", - "__reaction": "BrONO2 + hv -> Br + NO3", - "cross section": { - "type": "temperature based", - "parameterization": { - "type": "TAYLOR_SERIES", - "netcdf file": { - "file path": "data/cross_sections/BRONO2_JPL06.nc" - }, - "base temperature": 296.0, - "temperature ranges": [ - { - "maximum": 199.999999999999, - "fixed value": 200.0 - }, - { - "minimum": 200.0, - "maximum": 296.0 - }, - { - "minimum": 296.00000000001, - "fixed value": 296.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 0.85 - } - }, - { - "name": "jbrono2_b", - "__reaction": "BrONO2 + hv -> BrO + NO2", - "cross section": { - "type": "temperature based", - "parameterization": { - "type": "TAYLOR_SERIES", - "netcdf file": { - "file path": "data/cross_sections/BRONO2_JPL06.nc" - }, - "base temperature": 296.0, - "temperature ranges": [ - { - "maximum": 199.999999999999, - "fixed value": 200.0 - }, - { - "minimum": 200.0, - "maximum": 296.0 - }, - { - "minimum": 296.00000000001, - "fixed value": 296.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 0.15 - } - }, - { - "name": "jccl4", - "__reaction": "CCl4 + hv -> Products", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CCl4_1.nc" } - ], - "type": "CCl4+hv->Products" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcf2clbr", - "__reaction": "CF2BrCl + hv -> Products", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CF2BrCl_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcf3br", - "__reaction": "CF3Br + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/H1301_JPL06.nc", - "parameterization": { - "AA": [ 62.563, -2.0068, 1.6592e-2, -5.6465e-5, 6.7459e-8 ], - "BB": [ -9.1755e-1, 1.8575e-2, -1.3857e-4, 4.5066e-7, -5.3803e-10 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 178.0, - "maximum wavelength": 280.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210.0, - "maximum": 300.0 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcfcl3", - "__reaction": "CCl3F + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CFCL3_JPL06.nc", - "parameterization": { - "AA": [ -84.611, 7.9551e-1, -2.0550e-3, -4.4812e-6, 1.5838e-8 ], - "BB": [ -5.7912, 1.1689e-1, -8.8069e-4, 2.9335e-6, -3.6421e-9 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 174.1, - "maximum wavelength": 230.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210, - "maximum": 300 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcfc113", - "__reaction": "CFC-113 + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CFC113_JPL06.nc", - "parameterization": { - "AA": [ -1087.9, 20.004, -1.3920e-1, 4.2828e-4, -4.9384e-7 ], - "BB": [ 12.493, -2.3937e-1, 1.7142e-3, -5.4393e-6, 6.4548e-9 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 182.0, - "maximum wavelength": 230.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210, - "maximum": 300 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcfc114", - "__reaction": "CFC-114 + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CFC114_JPL10.nc", - "parameterization": { - "AA": [ -160.50, 2.4807, -1.5202e-2, 3.8412e-5, -3.4373e-8 ], - "BB": [ -1.5296, 3.5248e-2, -2.9951e-4, 1.1129e-6, -1.5259e-9 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 172.0, - "maximum wavelength": 220.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210, - "maximum": 300 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcfc115", - "__reaction": "CFC-115 + hv -> Products", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/CFC115_JPL10.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcf2cl2", - "__reaction": "CCl2F2 + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CF2CL2_JPL06.nc", - "parameterization": { - "AA": [ -43.8954569, -2.403597e-1, -4.2619e-4, 9.8743e-6, 0.0 ], - "BB": [ 4.8438e-3, 4.96145e-4, -5.6953e-6, 0.0, 0.0 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 200.0, - "maximum wavelength": 231.0, - "base temperature": 296.0, - "base wavelength": 200.0, - "logarithm": "natural", - "temperature ranges": [ - { - "maximum": 219.999999999999, - "fixed value": 220.0 - }, - { - "minimum": 220, - "maximum": 296 - }, - { - "minimum": 296.00000000001, - "fixed value": 296.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jch2br2", - "__reaction": "CH2BR2 + hv -> 2*BR", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CH2BR2_1.nc", - "parameterization": { - "AA": [ -70.211776, 1.940326e-1, 2.726152e-3, -1.695472e-5, 2.500066e-8 ], - "BB": [ 2.899280, -4.327724e-2, 2.391599e-4, -5.807506e-7, 5.244883e-10 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 210.0, - "maximum wavelength": 290.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210, - "maximum": 300 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jch3br", - "__reaction": "CH3Br + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CH3BR_JPL06.nc", - "parameterization": { - "AA": [ 46.520, -1.4580, 1.1469e-2, -3.7627e-5, 4.3264e-8 ], - "BB": [ 9.3408e-1, -1.6887e-2, 1.1487e-4, -3.4881e-7, 3.9945e-10 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 200.0, - "maximum wavelength": 280.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210, - "maximum": 300 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jch3ccl3", - "__reaction": "CH3CCl3+hv->Products", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CH3CCl3_1.nc" } - ], - "type": "tint" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jch3cl", - "__reaction": "CH3Cl + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CH3CL_JPL06.nc", - "parameterization": { - "AA": [ -299.80, 5.1047, -3.3630e-2, 9.5805e-5, -1.0135e-7 ], - "BB": [ -7.1727, 1.4837e-1, -1.1463e-3, 3.9188e-6, -4.9994e-9 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 174.1, - "maximum wavelength": 216.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210, - "maximum": 300 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jchbr3", - "__reaction": "CHBr3 + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/CHBR3_JPL10.nc", - "parameterization": { - "AA": [ -32.6067, 0.10308, 6.39e-5, -7.7392e-7, -2.2513e-9, 6.1376e-12 ], - "BB": [ 0.1582, -0.0014758, 3.8058e-6, 9.187e-10, -1.0772e-11, 0.0 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0, 5.0 ], - "minimum wavelength": 260.0, - "maximum wavelength": 362.0, - "base temperature": 296.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "invert temperature offset": true, - "temperature ranges": [ - { - "maximum": 259.999999999999, - "fixed value": 260.0 - }, - { - "minimum": 260.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcl2", - "__reaction": "Cl2 + hv -> Cl + Cl", - "cross section": { - "type": "Cl2+hv->Cl+Cl" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcl2o2", - "__reaction": "ClOOCl + hv -> Cl + ClOO", - "__comments": "TODO - this doesn't exactly match the products in TS1", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CL2O2_JPL10.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jclo", - "__reaction": "ClO + hv -> Cl + O", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CLO_JPL06.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jclono2_a", - "__reaction": "ClONO2 + hv -> Cl + NO3", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/ClONO2_1.nc" } - ], - "type": "ClONO2" - }, - "quantum yield": { - "type": "ClONO2+hv->Cl+NO3" - } - }, - { - "name": "jclono2_b", - "__reaction": "ClONO2 + hv -> ClO + NO2", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/ClONO2_1.nc" } - ], - "type": "ClONO2" - }, - "quantum yield": { - "type": "ClONO2+hv->ClO+NO2" - } - }, - { - "name": "jcof2", - "__reaction": "CF2O + hv -> Products", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CF2O_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jcofcl", - "__reaction": "CClFO + hv -> Products", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/CClFO_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jh2402", - "__reaction": "H2402 + hv -> 2*BR + 2*COF2", - "__comments": "TUV data set name CF2BrCF2Br", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/H2402_JPL06.nc", - "parameterization": { - "AA": [ 34.026, -1.152616, 8.959798e-3, -2.9089e-5, 3.307212e-8 ], - "BB": [ 4.010664e-1, -8.358968e-3, 6.415741e-5, -2.157554e-7, 2.691871e-10 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 190.0, - "maximum wavelength": 290.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210.0, - "maximum": 300.0 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhcfc141b", - "__reaction": "HCFC-141b + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/HCFC141b_JPL10.nc", - "parameterization": { - "AA": [ -682.913042, 12.122290, -8.187699e-2, 2.437244e-4, -2.719103e-7 ], - "BB": [ 4.074747, -8.053899e-2, 5.946552e-4, -1.945048e-6, 2.380143e-9 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 172.0, - "maximum wavelength": 240.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210.0, - "maximum": 300.0 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhcfc142b", - "__reaction": "HCFC-142b + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/HCFC142b_JPL10.nc", - "parameterization": { - "AA": [ -328.092008, 6.342799, -4.810362e-2, 1.611991e-4, -2.042613e-7 ], - "BB": [ 4.289533e-1, -9.042817e-3, 7.018009e-5, -2.389064e-7, 3.039799e-10 ], - "lp": [ 0.0, 1.0, 2.0, 3.0, 4.0 ], - "minimum wavelength": 172.0, - "maximum wavelength": 230.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210.0, - "maximum": 300.0 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhcfc22", - "__reaction": "HCFC-22 + hv -> Products", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/HCFC22_JPL06.nc", - "parameterization wavelength grid": { - "name": "custom wavelengths", - "type": "from config file", - "units": "nm", - "values": [ - 169.0, 171.0, 173.0, 175.0, 177.0, 179.0, 181.0, 183.0, 185.0, - 187.0, 189.0, 191.0, 193.0, 195.0, 197.0, 199.0, 201.0, 203.0, - 205.0, 207.0, 209.0, 211.0, 213.0, 215.0, 217.0, 219.0, 221.0 - ] - }, - "parameterization": { - "AA": [ -106.029, 1.5038, -8.2476e-3, 1.4206e-5 ], - "BB": [ -1.3399e-1, 2.7405e-3, -1.8028e-5, 3.8504e-8 ], - "lp": [ 0.0, 1.0, 2.0, 3.0 ], - "minimum wavelength": 174.0, - "maximum wavelength": 204.0, - "base temperature": 273.0, - "base wavelength": 0.0, - "logarithm": "base 10", - "temperature ranges": [ - { - "maximum": 209.999999999999, - "fixed value": 210.0 - }, - { - "minimum": 210.0, - "maximum": 300.0 - }, - { - "minimum": 300.00000000001, - "fixed value": 300.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhcl", - "__reaction": "HCl + hv -> H + Cl", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/HCl_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhobr", - "__reaction": "HOBr + hv -> OH + Br", - "cross section": { - "type": "HOBr+hv->OH+Br" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhocl", - "__reaction": "HOCl + hv -> HO + Cl", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/HOCl_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "joclo", - "__reaction": "OClO + hv -> Products", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/OClO_1.nc" }, - { "file path": "data/cross_sections/OClO_2.nc" }, - { "file path": "data/cross_sections/OClO_3.nc" } - ], - "type": "OClO+hv->Products" - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jho2no2_a", - "__reaction": "HNO4 + hv -> OH + NO3", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/HO2NO2_JPL06.nc", - "parameterization": { - "type": "BURKHOLDER", - "netcdf file": { - "file path": "data/cross_sections/HO2NO2_temp_JPL06.nc" - }, - "A": -988.0, - "B": 0.69, - "temperature ranges": [ - { - "maximum": 279.999999999999, - "fixed value": 280.0 - }, - { - "minimum": 280.0, - "maximum": 350.0 - }, - { - "minimum": 350.00000000001, - "fixed value": 350.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 0.30, - "override bands": [ - { - "band": "range", - "minimum wavelength": 200.0, - "value": 0.20 - } - ] - } - }, - { - "name": "jho2no2_b", - "__reaction": "HNO4 + hv -> HO2 + NO2", - "cross section": { - "type": "temperature based", - "netcdf file": "data/cross_sections/HO2NO2_JPL06.nc", - "parameterization": { - "type": "BURKHOLDER", - "netcdf file": { - "file path": "data/cross_sections/HO2NO2_temp_JPL06.nc" - }, - "A": -988.0, - "B": 0.69, - "temperature ranges": [ - { - "maximum": 279.999999999999, - "fixed value": 280.0 - }, - { - "minimum": 280.0, - "maximum": 350.0 - }, - { - "minimum": 350.00000000001, - "fixed value": 350.0 - } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 0.70, - "override bands": [ - { - "band": "range", - "minimum wavelength": 200.0, - "value": 0.80 - } - ] - } - }, - { - "name": "jmacr_a", - "__reaction": "CH2=C(CH3)CHO->1.34HO2+0.66MCO3+1.34CH2O+CH3CO3", - "__comments": "Methacrolein photolysis channel 1", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/Methacrolein_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 0.005 - } - }, - { - "name": "jmacr_b", - "__reaction": "CH2=C(CH3)CHO->0.66OH+1.34CO", - "__comments": "Methacrolein photolysis channel 2", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/Methacrolein_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 0.005 - } - }, - { - "name": "jhyac", - "__reaction": "CH2(OH)COCH3->CH3CO3+HO2+CH2O", - "__comments": "hydroxy acetone TODO: the products of this reaction differ from standalone TUV-x", - "cross section": { - "netcdf files": [ - { "file path": "data/cross_sections/Hydroxyacetone_1.nc" } - ], - "type": "base" - }, - "quantum yield": { - "type": "base", - "constant value": 0.65 - } - }, - { - "name": "jh2o_a", - "__reaction": "H2O + hv -> OH + H", - "cross section": { - "type": "base", - "merge data": true, - "netcdf files": [ - { - "file path": "data/cross_sections/H2O_1.nc", - "zero above": 183.0 - }, - { - "file path": "data/cross_sections/H2O_2.nc", - "zero below": 183.00000000001, - "zero above": 190.0 - }, - { - "file path": "data/cross_sections/H2O_3.nc", - "zero below": 190.00000000001 - } - ] - }, - "quantum yield" : { - "type": "base", - "netcdf files": [ "data/quantum_yields/H2O_H_OH.nc" ] - } - }, - { - "name": "jh2o_b", - "__reaction": "H2O + hv -> H2 + O1D", - "cross section": { - "type": "base", - "merge data": true, - "netcdf files": [ - { - "file path": "data/cross_sections/H2O_1.nc", - "zero above": 183.0 - }, - { - "file path": "data/cross_sections/H2O_2.nc", - "zero below": 183.00000000001, - "zero above": 190.0 - }, - { - "file path": "data/cross_sections/H2O_3.nc", - "zero below": 190.00000000001 - } - ] - }, - "quantum yield" : { - "type": "base", - "netcdf files": [ "data/quantum_yields/H2O_H2_O1D.nc" ] - } - }, - { - "name": "jh2o_c", - "__reaction": "H2O + hv -> 2*H + O", - "cross section": { - "type": "base", - "merge data": true, - "netcdf files": [ - { - "file path": "data/cross_sections/H2O_1.nc", - "zero above": 183.0 - }, - { - "file path": "data/cross_sections/H2O_2.nc", - "zero below": 183.00000000001, - "zero above": 190.0 - }, - { - "file path": "data/cross_sections/H2O_3.nc", - "zero below": 190.00000000001 - } - ] - }, - "quantum yield" : { - "type": "base", - "netcdf files": [ "data/quantum_yields/H2O_2H_O3P.nc" ] - } - }, - { - "name": "jch4_a", - "__reaction": "CH4 + hv -> H + CH3O2", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/CH4_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 0.45 - } - }, - { - "name": "jch4_b", - "__reaction": "CH4 + hv -> 1.44*H2 + 0.18*CH2O + 0.18*O + 0.33*OH + 0.33*H + 0.44*CO2 + 0.38*CO + 0.05*H2O", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/CH4_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 0.55 - } - }, - { - "name": "jco2", - "__reaction": "CO2 + hv -> CO + O", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/CO2_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhbr", - "__reaction": "HBR + hv -> BR + H", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/HBr_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jhf", - "__reaction": "HF + hv -> H + F", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/HF_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jsf6", - "__reaction": "SF6 + hv -> sink", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/SF6_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jh2so4", - "__reaction": "H2SO4 + hv -> SO3 + H2O", - "cross section": { - "type": "base", - "data": { - "default value": 0.0, - "point values": [ - { "wavelength": 121.65, "value": 6.3e-17 }, - { "wavelength": 525.0, "value": 1.43e-26 }, - { "wavelength": 625.0, "value": 1.8564e-25 }, - { "wavelength": 725.0, "value": 3.086999e-24 } - ] - } - }, - "quantum yield": { - "type": "H2SO4 Mills", - "netcdf files": [ - "data/quantum_yields/H2SO4_mills.nc" - ], - "parameterized wavelengths": [ - 525, - 625, - 725 - ], - "collision interval s": [ - 1.1e-9, - 8.9e-9, - 1.7e-7 - ], - "molecular diameter m": 4.18e-10, - "molecular weight kg mol-1": 98.078479e-3 - } - }, - { - "name": "jocs", - "__reaction": "OCS + hv -> S + CO", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/OCS_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jso", - "__reaction": "SO + hv -> S + O", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/SO_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jso2", - "__reaction": "SO2 + hv -> SO + O", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/SO2_Mills.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jso3", - "__reaction": "SO3 + hv -> SO2 + O", - "cross section": { - "type": "base", - "netcdf files": [ - { "file path": "data/cross_sections/SO3_1.nc" } - ] - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - }, - { - "name": "jno_i", - "__reaction": "NO + hv -> NOp + e", - "cross section": { - "type": "base", - "data": { - "default value": 0.0, - "point values": [ - { "wavelength": 121.65, "value": 2.0e-18 } - ] - } - }, - "quantum yield": { - "type": "base", - "constant value": 1.0 - } - } - ] - }, - "__CAM options": { - "aliasing": { - "default matching": "backup", - "pairs": [ - { - "to": "jalknit", - "__reaction": "ALKNIT + hv -> NO2 + 0.4*CH3CHO + 0.1*CH2O + 0.25*CH3COCH3 + HO2 + 0.8*MEK", - "from": "jch3ooh" - }, - { - "to": "jpooh", - "__reaction": "POOH (C3H6OHOOH) + hv -> CH3CHO + CH2O + HO2 + OH", - "from": "jch3ooh" - }, - { - "to": "jch3co3h", - "__reaction": "CH3COOOH + hv -> CH3O2 + OH + CO2", - "from": "jh2o2", - "scale by": 0.28 - }, - { - "to": "jmpan", - "__reaction": "MPAN + hv -> MCO3 + NO2", - "from": "jpan" - }, - { - "to": "jc2h5ooh", - "__reaction": "C2H5OOH + hv -> CH3CHO + HO2 + OH", - "from": "jch3ooh" - }, - { - "to": "jc3h7ooh", - "__reaction": "C3H7OOH + hv -> 0.82*CH3COCH3 + OH + HO2", - "from": "jch3ooh" - }, - { - "to": "jc6h5ooh", - "__reaction": "C6H5OOH + hv -> PHENO + OH", - "from": "jch3ooh" - }, - { - "to": "jeooh", - "__reaction": "EOOH + hv -> EO + OH", - "from": "jch3ooh" - }, - { - "to": "jrooh", - "__reaction": "ROOH + hv -> CH3CO3 + CH2O + OH", - "from": "jch3ooh" - }, - { - "to": "jxooh", - "__reaction": "XOOH + hv -> OH", - "from": "jch3ooh" - }, - { - "to": "jonitr", - "__reaction": "ONITR + hv -> NO2", - "from": "jch3cho" - }, - { - "to": "jisopooh", - "__reaction": "ISOPOOH + hv -> 0.402*MVK + 0.288*MACR + 0.69*CH2O + HO2", - "from": "jch3ooh" - }, - { - "to": "jmek", - "__reaction": "MEK + hv -> CH3CO3 + C2H5O2", - "from": "jacet" - }, - { - "to": "jalkooh", - "__reaction": "ALKOOH + hv -> .4*CH3CHO + .1*CH2O + .25*CH3COCH3 + .9*HO2 + .8*MEK + OH", - "from": "jch3ooh" - }, - { - "to": "jbenzooh", - "__reaction": "BENZOOH + hv -> OH + GLYOXAL + 0.5*BIGALD1 + HO2", - "from": "jch3ooh" - }, - { - "to": "jbepomuc", - "__reaction": "BEPOMUC + hv -> BIGALD1 + 1.5*HO2 + 1.5*CO", - "from": "jno2", - "scale by": 0.1 - }, - { - "to": "jbigald", - "__reaction": "BIGALD + hv -> 0.45*CO + 0.13*GLYOXAL + 0.56*HO2 + 0.13*CH3CO3 + 0.18*CH3COCHO", - "from": "jno2", - "scale by": 0.2 - }, - { - "to": "jbigald1", - "__reaction": "BIGALD1 + hv -> 0.6*MALO2 + HO2", - "from": "jno2", - "scale by": 0.14 - }, - { - "to": "jbigald2", - "__reaction": "BIGALD2 + hv -> 0.6*HO2 + 0.6*DICARBO2", - "from": "jno2", - "scale by": 0.2 - }, - { - "to": "jbigald3", - "__reaction": "BIGALD3 + hv -> 0.6*HO2 + 0.6*CO + 0.6*MDIALO2", - "from": "jno2", - "scale by": 0.2 - }, - { - "to": "jbigald4", - "__reaction": "BIGALD4 + hv -> HO2 + CO + CH3COCHO + CH3CO3", - "from": "jno2", - "scale by": 0.006 - }, - { - "to": "jbzooh", - "__reaction": "BZOOH + hv -> BZALD + OH + HO2", - "from": "jch3ooh" - }, - { - "to": "jmekooh", - "__reaction": "MEKOOH + hv -> OH + CH3CO3 + CH3CHO", - "from": "jch3ooh" - }, - { - "to": "jtolooh", - "__reaction": "TOLOOH + hv -> OH + .45*GLYOXAL + .45*CH3COCHO + .9*BIGALD", - "from": "jch3ooh" - }, - { - "to": "jterpooh", - "__reaction": "TERPOOH + hv -> OH + .1*CH3COCH3 + HO2 + MVK + MACR", - "from": "jch3ooh" - }, - { - "to": "jhonitr", - "__reaction": "HONITR + hv -> NO2 + 0.67*HO2 + 0.33*CH3CHO + 0.33*CH2O + 0.33*CO + 0.33*GLYALD + 0.33*CH3CO3 + 0.17*HYAC + 0.17*CH3COCH3", - "from": "jch2o_a" - }, - { - "to": "jhpald", - "__reaction": "HPALD + hv -> BIGALD3 + OH + HO2", - "from": "jno2", - "scale by": 0.006 - }, - { - "to": "jisopnooh", - "__reaction": "ISOPNOOH + hv -> NO2 + HO2 + ISOPOOH", - "from": "jch3ooh" - }, - { - "to": "jnc4cho", - "__reaction": "NC4CHO + hv -> BIGALD3 + NO2 + HO2", - "from": "jch2o_a" - }, - { - "to": "jnoa", - "__reaction": "NOA + hv -> NO2 + CH2O + CH3CO3", - "from": "jch2o_a" - }, - { - "to": "jnterpooh", - "__reaction": "NTERPOOH + hv -> TERPROD1 + NO2 + OH", - "from": "jch3ooh" - }, - { - "to": "jphenooh", - "__reaction": "PHENOOH + hv -> OH + HO2 + 0.7*GLYOXAL", - "from": "jch3ooh" - }, - { - "to": "jtepomuc", - "__reaction": "TEPOMUC + hv -> 0.5*CH3CO3 + HO2 + 1.5*CO", - "from": "jno2", - "scale by": 0.1 - }, - { - "to": "jterp2ooh", - "__reaction": "TERP2OOH + hv -> OH + 0.375*CH2O + 0.3*CH3COCH3 + 0.25*CO + CO2 + TERPROD2 + HO2 + 0.25*GLYALD", - "from": "jch3ooh" - }, - { - "to": "jterpnit", - "__reaction": "TERPNIT + hv -> TERPROD1 + NO2 + HO2", - "from": "jch3ooh" - }, - { - "to": "jterprd1", - "__reaction": "TERPROD1 + hv -> HO2 + CO + TERPROD2", - "from": "jch3cho" - }, - { - "to": "jterprd2", - "__reaction": "TERPROD2 + hv -> 0.15*RO2 + 0.68*CH2O + 0.8*CO2 + 0.5*CH3COCH3 + 0.65*CH3CO3 + 1.2*HO2 + 1.7*CO", - "from": "jch3cho" - }, - { - "to": "jxylenooh", - "__reaction": "XYLENOOH + hv -> OH + HO2 + 0.34*GLYOXAL + 0.54*CH3COCHO + 0.06*BIGALD1 + 0.2*BIGALD2 + 0.15*BIGALD3 + 0.21*BIGALD4", - "from": "jch3ooh" - }, - { - "to": "jxylolooh", - "__reaction": "XYLOLOOH + hv -> OH + 0.17*GLYOXAL + 0.51*CH3COCHO + HO2", - "from": "jch3ooh" - }, - { - "to": "jsoa1_a1", - "__reaction": "soa1_a1 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa1_a2", - "__reaction": "soa1_a2 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa2_a1", - "__reaction": "soa2_a1 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa2_a2", - "__reaction": "soa2_a2 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa3_a1", - "__reaction": "soa3_a1 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa3_a2", - "__reaction": "soa3_a2 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa4_a1", - "__reaction": "soa4_a1 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa4_a2", - "__reaction": "soa4_a2 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa5_a1", - "__reaction": "soa5_a1 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jsoa5_a2", - "__reaction": "soa5_a2 + hv -> Products", - "from": "jno2", - "scale by": 0.0004 - }, - { - "to": "jglyoxal", - "__reaction": "GLYOXAL + hv -> 2*CO + 2*HO2", - "from": "jmgly" - } - ] - } - } -} diff --git a/test/musica/tuvx/test_tuvx_cloud_optics.F90 b/test/musica/tuvx/test_tuvx_cloud_optics.F90 new file mode 100644 index 00000000..b01e75a8 --- /dev/null +++ b/test/musica/tuvx/test_tuvx_cloud_optics.F90 @@ -0,0 +1,115 @@ +! Copyright (C) 2024 National Science Foundation-National Center for Atmospheric Research +! SPDX-License-Identifier: Apache-2.0 +program test_tuvx_cloud_optics + + use musica_ccpp_tuvx_cloud_optics + +#define ASSERT(x) if (.not.(x)) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: x"; stop 1; endif +#define ASSERT_NEAR( a, b, abs_error ) if( (abs(a - b) >= abs_error) .and. (abs(a - b) /= 0.0) ) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: a, b"; stop 1; endif + + real, parameter :: ABS_ERROR = 1e-5 + + call test_create_cloud_optics_radiator() + +contains + + subroutine test_create_cloud_optics_radiator() + + use musica_util, only: error_t + use musica_ccpp_tuvx_height_grid, only: create_height_grid + use musica_ccpp_tuvx_wavelength_grid, only: create_wavelength_grid + use musica_tuvx_grid, only: grid_t + use musica_tuvx_radiator, only: radiator_t + use ccpp_kinds, only: kind_phys + + integer, parameter :: NUM_HOST_HEIGHT_MIDPOINTS = 2 + integer, parameter :: NUM_HOST_HEIGHT_INTERFACES = 3 + integer, parameter :: NUM_WAVELENGTH_MIDPOINTS = 3 + integer, parameter :: NUM_WAVELENGTH_INTERFACES = 4 + real(kind_phys) :: host_wavelength_interfaces(NUM_WAVELENGTH_INTERFACES) = [180.0e-9_kind_phys, 200.0e-9_kind_phys, 240.0e-9_kind_phys, 300.0e-9_kind_phys] + real(kind_phys) :: delta_pressure(NUM_HOST_HEIGHT_MIDPOINTS) = [100.0_kind_phys, 200.0_kind_phys] + real(kind_phys) :: cloud_fraction(NUM_HOST_HEIGHT_MIDPOINTS) = [0.1_kind_phys, 0.0_kind_phys] + real(kind_phys) :: liquid_water_content(NUM_HOST_HEIGHT_MIDPOINTS) = [0.0003_kind_phys, 0.0004_kind_phys] + real(kind_phys) :: reciprocal_of_gravitational_acceleration = 0.1_kind_phys + real(kind_phys) :: cloud_optical_depth(NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS) + real(kind_phys) :: expected_cloud_optical_depth(NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS) = & + reshape([ 0.0_kind_phys, 0.14704591_kind_phys, 0.0_kind_phys, 0.0_kind_phys, 0.14704591_kind_phys, 0.0_kind_phys, 0.0_kind_phys, 0.14704591_kind_phys, 0.0_kind_phys, 0.0_kind_phys, 0.14704591_kind_phys, 0.0_kind_phys ], & + [ NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS ]) + real(kind_phys) :: single_scattering_albedo(NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS) + real(kind_phys) :: asymmetry_parameter(NUM_HOST_HEIGHT_MIDPOINTS+1, NUM_WAVELENGTH_MIDPOINTS,1) + type(grid_t), pointer :: height_grid => null() + type(grid_t), pointer :: wavelength_grid => null() + type(radiator_t), pointer :: clouds => null() + type(error_t) :: error + character(len=512) :: errmsg + integer :: errcode + integer :: i + + height_grid => create_height_grid(NUM_HOST_HEIGHT_MIDPOINTS, NUM_HOST_HEIGHT_INTERFACES, & + errmsg, errcode) + ASSERT(errcode == 0) + ASSERT(associated(height_grid)) + + wavelength_grid => create_wavelength_grid(host_wavelength_interfaces, errmsg, errcode) + ASSERT(errcode == 0) + ASSERT(associated(wavelength_grid)) + + clouds => create_cloud_optics_radiator(height_grid, wavelength_grid, errmsg, errcode) + ASSERT(errcode == 0) + ASSERT(associated(clouds)) + + call set_cloud_optics_values(clouds, cloud_fraction, [ 0.0_kind_phys ], & + liquid_water_content, & + reciprocal_of_gravitational_acceleration, & + errmsg, errcode) + ASSERT(errcode == 1) + + call set_cloud_optics_values(clouds, cloud_fraction, delta_pressure, & + [ 1.0_kind_phys ], & + reciprocal_of_gravitational_acceleration, & + errmsg, errcode) + ASSERT(errcode == 1) + + call set_cloud_optics_values(clouds, [ 1.0_kind_phys ], delta_pressure, & + liquid_water_content, & + reciprocal_of_gravitational_acceleration, & + errmsg, errcode) + ASSERT(errcode == 1) + + call set_cloud_optics_values(clouds, cloud_fraction, delta_pressure, & + liquid_water_content, & + reciprocal_of_gravitational_acceleration, & + errmsg, errcode) + ASSERT(errcode == 0) + + call clouds%get_optical_depths(cloud_optical_depth, error) + ASSERT(error%is_success()) + do i = 1, size(cloud_optical_depth, dim=1) + do j = 1, size(cloud_optical_depth, dim=2) + ASSERT_NEAR(cloud_optical_depth(i,j), expected_cloud_optical_depth(i,j), ABS_ERROR) + end do + end do + + call clouds%get_single_scattering_albedos(single_scattering_albedo, error) + ASSERT(error%is_success()) + do i = 1, size(single_scattering_albedo, dim=1) + do j = 1, size(single_scattering_albedo, dim=2) + ASSERT_NEAR(single_scattering_albedo(i,j), 0.0_kind_phys, ABS_ERROR) + end do + end do + + call clouds%get_asymmetry_factors(asymmetry_parameter, error) + ASSERT(error%is_success()) + do i = 1, size(asymmetry_parameter, dim=1) + do j = 1, size(asymmetry_parameter, dim=2) + ASSERT_NEAR(asymmetry_parameter(i,j,1), 0.0_kind_phys, ABS_ERROR) + end do + end do + + deallocate( height_grid ) + deallocate( wavelength_grid ) + deallocate( clouds ) + + end subroutine test_create_cloud_optics_radiator + +end program test_tuvx_cloud_optics diff --git a/test/musica/tuvx/test_tuvx_extraterrestrial_flux.F90 b/test/musica/tuvx/test_tuvx_extraterrestrial_flux.F90 new file mode 100644 index 00000000..dcb94696 --- /dev/null +++ b/test/musica/tuvx/test_tuvx_extraterrestrial_flux.F90 @@ -0,0 +1,67 @@ +program test_tuvx_extraterrestrial_flux + + use musica_ccpp_tuvx_extraterrestrial_flux + + 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_extraterrestrial_flux() + +contains + + subroutine test_update_extraterrestrial_flux() + use musica_util, only: error_t + use musica_tuvx_grid, only: grid_t + use musica_tuvx_profile, only: profile_t + use musica_ccpp_tuvx_wavelength_grid, only: create_wavelength_grid + use ccpp_kinds, only: kind_phys + + integer, parameter :: NUM_WAVELENGTH_BINS = 4 + integer, parameter :: NUM_PHOTOLYSIS_WAVELENGTH_GRID_SECTIONS = 8 + real, parameter :: ABS_ERROR = 1e-6 + real, parameter :: MAGNITUDE_REDUCER = 1e-15 ! Its purpose is to adjust the magnitude of the values to reduce absolute errors + real(kind_phys) :: wavelength_grid_interfaces(NUM_WAVELENGTH_BINS + 1) = & + [200.0e-9_kind_phys, 220.0e-9_kind_phys, 240.0e-9_kind_phys, 260.0e-9_kind_phys, 280.0e-9_kind_phys] ! m + real(kind_phys) :: photolysis_wavelength_grid_interfaces(NUM_PHOTOLYSIS_WAVELENGTH_GRID_SECTIONS + 1) = & + [200.0_kind_phys, 210.0_kind_phys, 220.0_kind_phys, 230.0_kind_phys, 240.0_kind_phys, & + 250.0_kind_phys, 260.0_kind_phys, 270.0_kind_phys, 280.0_kind_phys] ! nm + real(kind_phys) :: extraterrestrial_flux(NUM_PHOTOLYSIS_WAVELENGTH_GRID_SECTIONS) = & + [1.5e13_kind_phys, 1.5e13_kind_phys, 1.4e13_kind_phys, 1.4e13_kind_phys, & + 1.3e13_kind_phys, 1.2e13_kind_phys, 1.1e13_kind_phys, 1.0e13_kind_phys] + real(kind_phys) :: expected_extraterrestrial_flux_midpoints(NUM_WAVELENGTH_BINS) = & + [3.0e14_kind_phys, 2.8e14_kind_phys, 2.5e14_kind_phys, 2.1e14_kind_phys] + real(kind_phys) :: extraterrestrial_flux_midpoints(NUM_WAVELENGTH_BINS) + type(grid_t), pointer :: wavelength_grid + type(profile_t), pointer :: profile + type(error_t) :: error + character(len=512) :: errmsg + integer :: errcode + integer :: i + + wavelength_grid => create_wavelength_grid (wavelength_grid_interfaces, errmsg, errcode ) + ASSERT(errcode == 0) + ASSERT(associated(wavelength_grid)) + + profile => create_extraterrestrial_flux_profile( wavelength_grid, wavelength_grid_interfaces, errmsg, errcode ) + ASSERT(errcode == 0) + ASSERT(associated(profile)) + + call set_extraterrestrial_flux_values( profile, photolysis_wavelength_grid_interfaces, & + extraterrestrial_flux, errmsg, errcode ) + ASSERT(errcode == 0) + + call profile%get_midpoint_values( extraterrestrial_flux_midpoints, error ) + ASSERT(error%is_success()) + + do i = 1, size(extraterrestrial_flux_midpoints) + ASSERT_NEAR(extraterrestrial_flux_midpoints(i) * MAGNITUDE_REDUCER, expected_extraterrestrial_flux_midpoints(i) * MAGNITUDE_REDUCER, ABS_ERROR) + end do + + deallocate( profile ) + deallocate( wavelength_grid ) + + end subroutine test_update_extraterrestrial_flux + +end program test_tuvx_extraterrestrial_flux \ No newline at end of file diff --git a/test/test_schemes/initialize_constituents.F90 b/test/test_schemes/initialize_constituents.F90 new file mode 100644 index 00000000..e9cc3fae --- /dev/null +++ b/test/test_schemes/initialize_constituents.F90 @@ -0,0 +1,155 @@ +module initialize_constituents + +implicit none +private + +public :: initialize_constituents_register +public :: initialize_constituents_run + +contains + +!> \section arg_table_initialize_constituents_register Argument Table +!! \htmlinclude initialize_constituents_register.html +subroutine initialize_constituents_register(constituents, errmsg, errcode) + use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t + use cam_initfiles, only: initial_file_get_id + use pio, only: pio_inquire, file_desc_t, pio_inq_varname + use ccpp_kinds, only: kind_phys + ! Dummy variables + type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituents(:) + character(len=512), intent(out) :: errmsg + integer, intent(out) :: errcode + ! Local variables + type(file_desc_t), pointer :: ncdata + integer :: num_variables + integer :: ierr + integer :: var_index + integer :: constituent_index + integer :: known_const_index + integer :: found_const_count + logical :: known_constituent + character(len=256) :: variable_name + character(len=512) :: alloc_err_msg + character(len=256), allocatable :: constituent_names(:) + character(len=65), parameter :: water_species_std_names(6) = & + (/'water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water ', & + 'cloud_liquid_water_mixing_ratio_wrt_moist_air_and_condensed_water', & + 'rain_mixing_ratio_wrt_moist_air_and_condensed_water ', & + 'cloud_ice_mixing_ratio_wrt_moist_air_and_condensed_water ', & + 'snow_mixing_ratio_wrt_moist_air_and_condensed_water ', & + 'graupel_water_mixing_ratio_wrt_moist_air_and_condensed_water '/) + + character(len=11), parameter :: const_file_names(6) = (/'cnst_Q ', & + 'cnst_CLDLIQ', & + 'cnst_RAINQM', & + 'cnst_CLDICE', & + 'cnst_SNOWQM', & + 'cnst_GRAUQM'/) + character(len=11), parameter :: water_species_number_concentrations(5) = & + (/'cnst_NUMLIQ', & + 'cnst_NUMRAI', & + 'cnst_NUMICE', & + 'cnst_NUMSNO', & + 'cnst_NUMGRA'/) + + errcode = 0 + errmsg = '' + constituent_index = 0 + found_const_count = 0 + + ncdata => initial_file_get_id() + ! See how many variables are present on the file + ierr = pio_inquire(ncdata, nVariables=num_variables) + allocate(constituent_names(num_variables), stat=ierr, errmsg=alloc_err_msg) + if (ierr /= 0) then + errcode = 1 + write(errmsg,*) 'Failed to allocate "constituent_names" in initialize_constituents_register: ', trim(alloc_err_msg) + return + end if + + ! Loop over all variables in the file and add each constituent to the + ! dynamic constituent array + do var_index = 1, num_variables + ierr = pio_inq_varname(ncdata, var_index, variable_name) + known_constituent = .false. + if (index(variable_name, 'cnst_') > 0) then + constituent_index = constituent_index + 1 + ! Replace with standard name if known, to avoid duplicates + if (found_const_count < size(water_species_std_names)) then + do known_const_index = 1, size(const_file_names) + if (trim(const_file_names(known_const_index)) == trim(variable_name)) then + constituent_names(constituent_index) = water_species_std_names(known_const_index) + found_const_count = found_const_count + 1 + known_constituent = .true. + exit + end if + end do + end if + if (.not. known_constituent) then + constituent_names(constituent_index) = variable_name + end if + end if + end do + + allocate(constituents(constituent_index), stat=ierr, errmsg=alloc_err_msg) + if (ierr /= 0) then + errcode = 1 + write(errmsg,*) 'Failed to allocate "constituents" in initialize_constituents_register: ', trim(alloc_err_msg) + return + end if + + do var_index = 1, size(constituents) + if (any(water_species_number_concentrations == trim(constituent_names(var_index)))) then + call constituents(var_index)%instantiate( & + std_name = constituent_names(var_index), & + long_name = constituent_names(var_index), & + units = 'kg-1', & + vertical_dim = 'vertical_layer_dimension', & + min_value = 0.0_kind_phys, & + advected = .true., & + water_species = .true., & + mixing_ratio_type = 'wet', & + errcode = errcode, & + errmsg = errmsg) + else if (any(water_species_std_names == trim(constituent_names(var_index)))) then + ! Do not set water_species = .true. for water species mixing ratios + ! Avoiding mismatch in properties vs. metadata-specified constituents + ! Water species properties are set in air_composition.F90 in CAM-SIMA + call constituents(var_index)%instantiate( & + std_name = constituent_names(var_index), & + long_name = constituent_names(var_index), & + units = 'kg kg-1', & + vertical_dim = 'vertical_layer_dimension', & + min_value = 0.0_kind_phys, & + advected = .true., & + mixing_ratio_type = 'wet', & + errcode = errcode, & + errmsg = errmsg) + else + call constituents(var_index)%instantiate( & + std_name = constituent_names(var_index), & + long_name = constituent_names(var_index), & + units = 'kg kg-1', & + vertical_dim = 'vertical_layer_dimension', & + min_value = 0.0_kind_phys, & + advected = .true., & + errcode = errcode, & + errmsg = errmsg) + end if + if (errcode /= 0) then + exit + end if + end do + + end subroutine initialize_constituents_register + +!> \section arg_table_initialize_constituents_run Argument Table +!! \htmlinclude initialize_constituents_run.html + subroutine initialize_constituents_run(errmsg, errcode) + character(len=512), intent(out) :: errmsg + integer, intent(out) :: errcode + + errcode = 0 + errmsg = '' + end subroutine initialize_constituents_run +end module initialize_constituents diff --git a/test/test_schemes/initialize_constituents.meta b/test/test_schemes/initialize_constituents.meta new file mode 100644 index 00000000..0618c16d --- /dev/null +++ b/test/test_schemes/initialize_constituents.meta @@ -0,0 +1,44 @@ +[ccpp-table-properties] + name = initialize_constituents + type = scheme +[ccpp-arg-table] + name = initialize_constituents_register + type = scheme +[ constituents ] + standard_name = dynamic_constituents_for_initialize_constituents + units = none + dimensions = (:) + allocatable = True + type = ccpp_constituent_properties_t + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=512 + dimensions = () + intent = out +[ errcode ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out +[ccpp-arg-table] + name = initialize_constituents_run + type = scheme +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=512 + dimensions = () + intent = out +[ errcode ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out diff --git a/test/test_suites/suite_dry_adiabatic_adjust.xml b/test/test_suites/suite_dry_adiabatic_adjust.xml index 5bcb11a2..12fe91ca 100644 --- a/test/test_suites/suite_dry_adiabatic_adjust.xml +++ b/test/test_suites/suite_dry_adiabatic_adjust.xml @@ -4,7 +4,7 @@ dadadj - dadadj_apply_qv_tendency + apply_constituent_tendencies apply_heating_rate qneg geopotential_temp diff --git a/test/test_suites/suite_tj2016_precip.xml b/test/test_suites/suite_tj2016_precip.xml index 7f6dd40b..381767bc 100644 --- a/test/test_suites/suite_tj2016_precip.xml +++ b/test/test_suites/suite_tj2016_precip.xml @@ -5,6 +5,8 @@ tj2016_precip apply_heating_rate qneg + sima_state_diagnostics + sima_tend_diagnostics diff --git a/test/test_suites/suite_tj2016_sfc_pbl_hs.xml b/test/test_suites/suite_tj2016_sfc_pbl_hs.xml index 2da13e64..267d8ed1 100644 --- a/test/test_suites/suite_tj2016_sfc_pbl_hs.xml +++ b/test/test_suites/suite_tj2016_sfc_pbl_hs.xml @@ -7,6 +7,8 @@ apply_tendency_of_eastward_wind apply_tendency_of_northward_wind qneg + sima_state_diagnostics + sima_tend_diagnostics diff --git a/test/valgrind.supp b/test/valgrind.supp index ed4ee497..ee1ba85e 100644 --- a/test/valgrind.supp +++ b/test/valgrind.supp @@ -35,4 +35,47 @@ fun:__tuvx_core_MOD_constructor fun:InternalCreateTuvx ... +} +{ + Suppress_MUSICA_TUV-x_CreateRadiator + Memcheck:Leak + match-leak-kinds: definite + fun:malloc + fun:__musica_string_MOD_string_assign_char + fun:__tuvx_radiator_from_host_MOD_constructor_char + fun:__tuvx_radiator_from_host_MOD_constructor_string + fun:InternalCreateRadiator + ... +} +{ + Suppress_MUSICA_TUV-x_AddRadiator + Memcheck:Leak + match-leak-kinds: definite + fun:malloc + fun:__tuvx_radiator_from_host_MOD___copy_tuvx_radiator_from_host_Radiator_from_host_t + fun:__tuvx_radiator_warehouse_MOD_add_radiator + fun:InternalAddRadiator + ... +} +{ + Suppress_MUSICA_TUV-x_GetRadiator + Memcheck:Leak + match-leak-kinds: definite + fun:malloc + fun:__tuvx_radiator_from_host_MOD___copy_tuvx_radiator_from_host_Radiator_from_host_t + fun:InternalGetRadiator + ... +} +{ + Suppress_MUSICA_TUV-x_CreateTuvx-RadiatorFromHost + Memcheck:Leak + match-leak-kinds: definite + fun:malloc + fun:__tuvx_radiator_from_host_MOD___copy_tuvx_radiator_from_host_Radiator_from_host_t + fun:__tuvx_radiator_warehouse_MOD_add_radiator + fun:__tuvx_radiator_warehouse_MOD_add_radiators + fun:__tuvx_radiative_transfer_MOD_constructor + fun:__tuvx_core_MOD_constructor + fun:InternalCreateTuvx + ... } \ No newline at end of file diff --git a/to_be_ccppized/ccpp_const_utils.F90 b/to_be_ccppized/ccpp_const_utils.F90 index 902d6059..e20aed31 100644 --- a/to_be_ccppized/ccpp_const_utils.F90 +++ b/to_be_ccppized/ccpp_const_utils.F90 @@ -13,13 +13,13 @@ subroutine ccpp_const_get_idx(constituent_props, name, cindex, errmsg, errflg) use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t ! Input arguments - type(ccpp_constituent_prop_ptr_t), pointer, intent(in) :: constituent_props(:) - character(len=*), intent(in) :: name ! constituent name + type(ccpp_constituent_prop_ptr_t), intent(in) :: constituent_props(:) + character(len=*), intent(in) :: name ! constituent name ! Output arguments - integer, intent(out) :: cindex ! global constituent index - character(len=512), intent(out) :: errmsg ! error message - integer, intent(out) :: errflg ! error flag + integer, intent(out) :: cindex ! global constituent index + character(len=512), intent(out) :: errmsg ! error message + integer, intent(out) :: errflg ! error flag ! Local variables integer :: t_cindex diff --git a/to_be_ccppized/ccpp_tuvx_utils.F90 b/to_be_ccppized/ccpp_tuvx_utils.F90 new file mode 100644 index 00000000..837801bc --- /dev/null +++ b/to_be_ccppized/ccpp_tuvx_utils.F90 @@ -0,0 +1,91 @@ +module ccpp_tuvx_utils + + implicit none + + private + public :: rebin, read_extraterrestrial_flux + + character(len=50), dimension(4), public :: filepath_of_extraterrestrial_flux + +contains + + !> Regrids normalized flux data to match a specified wavelength grid + !! 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 + + integer, intent(in) :: source_dimension + integer, intent(in) :: target_dimension + real(kind_phys), intent(in) :: source_coordinates(source_dimension+1) + real(kind_phys), intent(in) :: target_coordinates(target_dimension+1) + real(kind_phys), intent(in) :: source(source_dimension) + real(kind_phys), intent(out) :: target(target_dimension) + + ! local variables + integer :: i, si, si1, sil, siu + real(kind_phys) :: y, sl, su, tl, tu + + do i = 1, target_dimension + tl = target_coordinates(i) + if( tl < source_coordinates( source_dimension + 1) ) then + do sil = 1, source_dimension + 1 + if( tl <= source_coordinates( sil ) ) then + exit + end if + end do + tu = target_coordinates( i + 1 ) + do siu = 1, source_dimension + 1 + if( tu <= source_coordinates( siu ) ) then + exit + end if + end do + y = 0._kind_phys + sil = max( sil, 2 ) + siu = min( siu, source_dimension + 1 ) + do si = sil, siu + si1 = si - 1 + sl = max( tl, source_coordinates( si1 ) ) + su = min( tu, source_coordinates( si ) ) + y = y + ( su - sl ) * source( si1 ) + end do + target(i) = y / (target_coordinates( i + 1 ) - target_coordinates( i ) ) + else + target(i) = 0._kind_phys + end if + end do + + end subroutine rebin + + !> Reads a data file to retrieve the extraterrestrial radiation flux values. + ! This function is a temporary implementation and will be replaced in + ! future versions of the code. + subroutine read_extraterrestrial_flux(num_wavelength_grid_sections, & + wavelength_grid_interfaces, extraterrestrial_flux) + use ccpp_kinds, only: kind_phys + + integer, intent(out) :: num_wavelength_grid_sections ! (count) + real(kind_phys), allocatable, intent(out) :: wavelength_grid_interfaces(:) ! nm + real(kind_phys), allocatable, intent(out) :: extraterrestrial_flux(:) ! photons cm-2 s-1 nm-1 + + filepath_of_extraterrestrial_flux(1) = 'musica_configurations/chapman/tuvx/data/profiles/solar/susim_hi.flx' + filepath_of_extraterrestrial_flux(2) = 'musica_configurations/chapman/tuvx/data/profiles/solar/atlas3_1994_317_a.dat' + filepath_of_extraterrestrial_flux(3) = 'musica_configurations/chapman/tuvx/data/profiles/solar/sao2010.solref.converted' + filepath_of_extraterrestrial_flux(4) = 'musica_configurations/chapman/tuvx/data/profiles/solar/neckel.flx' + + num_wavelength_grid_sections = 8 + + allocate(wavelength_grid_interfaces(num_wavelength_grid_sections + 1)) + allocate(extraterrestrial_flux(num_wavelength_grid_sections)) + + wavelength_grid_interfaces(:) = & + [200.0_kind_phys, 210.0_kind_phys, 220.0_kind_phys, 230.0_kind_phys, & + 240.0_kind_phys, 250.0_kind_phys, 260.0_kind_phys, 270.0_kind_phys, 280.0_kind_phys] + + extraterrestrial_flux(:) = & + [1.5e13_kind_phys, 1.5e13_kind_phys, 1.4e13_kind_phys, 1.4e13_kind_phys, & + 1.3e13_kind_phys, 1.2e13_kind_phys, 1.1e13_kind_phys, 1.0e13_kind_phys] + + end subroutine read_extraterrestrial_flux + +end module ccpp_tuvx_utils \ No newline at end of file