diff --git a/.github/workflows/test.yaml b/.github/workflows/test.yaml
index 5d7bf4aa..238e0c1d 100644
--- a/.github/workflows/test.yaml
+++ b/.github/workflows/test.yaml
@@ -12,4 +12,14 @@ jobs:
- name: build Docker image
run: docker build -t musica -f test/docker/Dockerfile.musica .
- name: run tests in container
- run: docker run --name test-container -t musica bash -c 'make test ARGS="--rerun-failed --output-on-failure -j8"'
\ No newline at end of file
+ run: docker run --name test-container -t musica bash -c 'make test ARGS="--rerun-failed --output-on-failure -j8"'
+ test_musica_api_no_install:
+ runs-on: ubuntu-latest
+ steps:
+ - uses: actions/checkout@v4
+ with:
+ submodules: recursive
+ - name: build Docker image
+ run: docker build -t musica-no-install -f test/docker/Dockerfile.musica.no_install .
+ - name: run tests in container
+ run: docker run --name test-container -t musica-no-install bash -c 'make test ARGS="--rerun-failed --output-on-failure -j8"'
\ No newline at end of file
diff --git a/schemes/musica/micm/musica_ccpp_micm.F90 b/schemes/musica/micm/musica_ccpp_micm.F90
index d11e86d2..2419d40e 100644
--- a/schemes/musica/micm/musica_ccpp_micm.F90
+++ b/schemes/musica/micm/musica_ccpp_micm.F90
@@ -1,5 +1,4 @@
module musica_ccpp_micm
- use iso_c_binding
! Note: "micm_t" is included in an external pre-built MICM library that the host
! model is responsible for linking to during compilation
@@ -17,50 +16,52 @@ module musica_ccpp_micm
contains
- !> Register MICM constituents with the CCPP
- subroutine micm_register(constituents, solver_type, num_grid_cells, errmsg, errcode)
+ !> Register MICM constituent properties with the CCPP
+ subroutine micm_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode)
use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t
- use musica_micm, only: Rosenbrock, RosenbrockStandardOrder
- use musica_util, only: error_t
+ use musica_micm, only: Rosenbrock, RosenbrockStandardOrder
+ use musica_util, only: error_t
+ use iso_c_binding, only: c_int
- type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituents(:)
integer(c_int), intent(in) :: solver_type
integer(c_int), 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
! local variables
- type(error_t) :: error
- real(kind=kind_phys) :: molar_mass
- logical :: is_advected
- integer :: i
-
- errcode = 0
- errmsg = ''
+ type(error_t) :: error
+ real(kind=kind_phys) :: molar_mass
+ character(len=:), allocatable :: species_name
+ logical :: is_advected
+ integer :: i, species_index
micm => micm_t(filename_of_micm_configuration, solver_type, num_grid_cells, error)
if (has_error_occurred(error, errmsg, errcode)) return
- allocate(constituents(size(micm%species_ordering)), stat=errcode)
+ allocate(constituent_props(micm%species_ordering%size()), stat=errcode)
if (errcode /= 0) then
- errmsg = "[MUSICA Error] Failed to allocate memory for constituents."
+ errmsg = "[MUSICA Error] Failed to allocate memory for constituent properties."
return
end if
- do i = 1, size(micm%species_ordering)
- associate( map => micm%species_ordering(i) )
- molar_mass = micm%get_species_property_double(map%name(), &
+ do i = 1, micm%species_ordering%size()
+ associate( map => micm%species_ordering )
+ species_name = map%name(i)
+ species_index = map%index(i)
+
+ molar_mass = micm%get_species_property_double(species_name, &
"molecular weight [kg mol-1]", &
error)
if (has_error_occurred(error, errmsg, errcode)) return
- is_advected = micm%get_species_property_bool(map%name(), &
- "__is advected", &
- error)
+ is_advected = micm%get_species_property_bool(species_name, &
+ "__is advected", &
+ error)
if (has_error_occurred(error, errmsg, errcode)) return
- call constituents(map%index())%instantiate( &
- std_name = map%name(), &
- long_name = map%name(), &
+ call constituent_props(species_index)%instantiate( &
+ std_name = species_name, &
+ long_name = species_name, &
units = 'kg kg-1', &
vertical_dim = 'vertical_layer_dimension', &
default_value = 0.0_kind_phys, &
@@ -78,27 +79,28 @@ end subroutine micm_register
!> Intitialize MICM
subroutine micm_init(errmsg, errcode)
character(len=512), intent(out) :: errmsg
- integer, intent(out) :: errcode
+ integer, intent(out) :: errcode
- errcode = 0
errmsg = ''
+ errcode = 0
end subroutine micm_init
!> Solve chemistry at the current time step
- subroutine micm_run(time_step, temperature, pressure, dry_air_density, constituents, &
- rate_params, errmsg, errcode)
- use musica_micm, only: solver_stats_t
- use musica_util, only: string_t, error_t
-
- real(kind_phys), intent(in) :: time_step ! s
- real(c_double), target, intent(in) :: temperature(:) ! K
- real(c_double), target, intent(in) :: pressure(:) ! Pa
- real(c_double), target, intent(in) :: dry_air_density(:) ! kg m-3
- real(c_double), target, intent(inout) :: constituents(:) ! mol m-3
- real(c_double), target, intent(inout) :: rate_params(:)
- character(len=512), intent(out) :: errmsg
- integer, intent(out) :: errcode
+ subroutine micm_run(time_step, temperature, pressure, dry_air_density, &
+ user_defined_rate_parameters, constituents, errmsg, errcode)
+ use musica_micm, only: solver_stats_t
+ use musica_util, only: string_t, error_t
+ use iso_c_binding, only: c_double
+
+ real(kind_phys), intent(in) :: time_step ! s
+ real(c_double), intent(in) :: temperature(:) ! K
+ real(c_double), intent(in) :: pressure(:) ! Pa
+ real(c_double), intent(in) :: dry_air_density(:) ! kg m-3
+ real(c_double), intent(in) :: user_defined_rate_parameters(:) ! various units
+ real(c_double), intent(inout) :: constituents(:) ! mol m-3
+ character(len=512), intent(out) :: errmsg
+ integer, intent(out) :: errcode
! local variables
type(string_t) :: solver_state
@@ -107,18 +109,16 @@ subroutine micm_run(time_step, temperature, pressure, dry_air_density, constitue
real(c_double) :: c_time_step
integer :: i_elem
- errcode = 0
- errmsg = ''
c_time_step = real(time_step, c_double)
- call micm%solve(c_time_step, &
- temperature, &
- pressure, &
- dry_air_density, &
- constituents, &
- rate_params, &
- solver_state, &
- solver_stats, &
+ call micm%solve(c_time_step, &
+ temperature, &
+ pressure, &
+ dry_air_density, &
+ constituents, &
+ user_defined_rate_parameters, &
+ solver_state, &
+ solver_stats, &
error)
if (has_error_occurred(error, errmsg, errcode)) return
@@ -127,10 +127,15 @@ end subroutine micm_run
!> Finalize MICM
subroutine micm_final(errmsg, errcode)
character(len=512), intent(out) :: errmsg
- integer, intent(out) :: errcode
+ integer, intent(out) :: errcode
- errcode = 0
errmsg = ''
+ errcode = 0
+
+ if (associated( micm )) then
+ deallocate( micm )
+ micm => null()
+ end if
end subroutine micm_final
diff --git a/schemes/musica/micm/micm_util.F90 b/schemes/musica/micm/musica_ccpp_micm_util.F90
similarity index 50%
rename from schemes/musica/micm/micm_util.F90
rename to schemes/musica/micm/musica_ccpp_micm_util.F90
index 3272ee91..47b92c81 100644
--- a/schemes/musica/micm/micm_util.F90
+++ b/schemes/musica/micm/musica_ccpp_micm_util.F90
@@ -1,4 +1,4 @@
-module micm_util
+module musica_ccpp_micm_util
implicit none
private
@@ -7,94 +7,70 @@ module micm_util
contains
- subroutine reshape_into_micm_arr(temperature, pressure, dry_air_density, constituents, &
- rate_params, m_temperature, m_pressure, m_dry_air_density, m_constituents, m_rate_params)
+ !> Reshape array (2D/3D -> 1D) and convert type (kind_phys -> c_double)
+ subroutine reshape_into_micm_arr(temperature, pressure, dry_air_density, constituents, &
+ micm_temperature, micm_pressure, micm_dry_air_density, &
+ micm_constituents)
use iso_c_binding, only: c_double
use ccpp_kinds, only: kind_phys
- real(kind_phys), target, intent(in) :: temperature(:,:) ! K
- real(kind_phys), target, intent(in) :: pressure(:,:) ! Pa
- real(kind_phys), target, intent(in) :: dry_air_density(:,:) ! kg m-3
- real(kind_phys), target, intent(in) :: constituents(:,:,:) ! kg kg-1
- real(kind_phys), target, intent(in) :: rate_params(:,:,:)
- real(c_double), target, intent(out) :: m_temperature(:) ! K
- real(c_double), target, intent(out) :: m_pressure(:) ! Pa
- real(c_double), target, intent(out) :: m_dry_air_density(:) ! kg m-3
- real(c_double), target, intent(out) :: m_constituents(:) ! kg kg-1
- real(c_double), target, intent(out) :: m_rate_params(:)
+ real(kind_phys), intent(in) :: temperature(:,:) ! K
+ real(kind_phys), intent(in) :: pressure(:,:) ! Pa
+ real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3
+ real(kind_phys), intent(in) :: constituents(:,:,:) ! kg kg-1
+ real(c_double), intent(out) :: micm_temperature(:) ! K
+ real(c_double), intent(out) :: micm_pressure(:) ! Pa
+ real(c_double), intent(out) :: micm_dry_air_density(:) ! kg m-3
+ real(c_double), intent(out) :: micm_constituents(:) ! kg kg-1
! local variables
- integer :: num_columns, num_layers
- integer :: num_constituents, num_rate_params
- integer :: i_column, i_layer, i_elem, i_constituents, i_rate_params
+ integer :: num_columns, num_layers, num_constituents
+ integer :: i_column, i_layer, i_elem, i_constituents
num_columns = size(constituents, dim=1)
num_layers = size(constituents, dim=2)
num_constituents = size(constituents, dim=3)
- num_rate_params = size(rate_params, dim=3)
- ! Reshape into 1-D arry in species-column first order
- ! refers to: state.variables_[i_cell][i_species] = concentrations[i_species_elem++]
+ ! Reshape into 1-D arry in species-column first order, referring to
+ ! state.variables_[i_cell][i_species] = concentrations[i_species_elem++]
i_elem = 1
i_constituents = 1
- i_rate_params = 1
do i_layer = 1, num_layers
do i_column = 1, num_columns
- m_temperature(i_elem) = real(temperature(i_column, i_layer), c_double)
- m_pressure(i_elem) = real(pressure(i_column, i_layer), c_double)
- m_dry_air_density(i_elem) = real(dry_air_density(i_column, i_layer), c_double)
- m_constituents(i_constituents : i_constituents + num_constituents - 1) &
+ micm_temperature(i_elem) = real(temperature(i_column, i_layer), c_double)
+ micm_pressure(i_elem) = real(pressure(i_column, i_layer), c_double)
+ micm_dry_air_density(i_elem) = real(dry_air_density(i_column, i_layer), c_double)
+ micm_constituents(i_constituents : i_constituents + num_constituents - 1) &
= real(constituents(i_column, i_layer, :), c_double)
- m_rate_params(i_rate_params : i_rate_params + num_rate_params - 1) &
- = real(rate_params(i_column, i_layer, :), c_double)
i_elem = i_elem + 1
i_constituents = i_constituents + num_constituents
- i_rate_params = i_rate_params + num_rate_params
end do
end do
end subroutine reshape_into_micm_arr
- subroutine reshape_into_ccpp_arr(temperature, pressure, dry_air_density, constituents, &
- rate_params, m_temperature, m_pressure, m_dry_air_density, m_constituents, m_rate_params)
+ !> Reshape array (1D -> 3D) and convert type (c_double -> kind_phys)
+ subroutine reshape_into_ccpp_arr(micm_constituents, constituents)
use iso_c_binding, only: c_double
use ccpp_kinds, only: kind_phys
- real(kind_phys), intent(out) :: temperature(:,:) ! K
- real(kind_phys), intent(out) :: pressure(:,:) ! Pa
- real(kind_phys), intent(out) :: dry_air_density(:,:) ! kg m-3
- real(kind_phys), intent(out) :: constituents(:,:,:) ! kg kg-1
- real(kind_phys), intent(out) :: rate_params(:,:,:)
- real(c_double), intent(in) :: m_temperature(:) ! K
- real(c_double), intent(in) :: m_pressure(:) ! Pa
- real(c_double), intent(in) :: m_dry_air_density(:) ! kg m-3
- real(c_double), intent(in) :: m_constituents(:) ! kg kg-1
- real(c_double), intent(in) :: m_rate_params(:)
+
+ real(c_double), intent(in) :: micm_constituents(:) ! kg kg-1
+ real(kind_phys), intent(inout) :: constituents(:,:,:) ! kg kg-1
! local variables
- integer :: num_columns, num_layers
- integer :: num_constituents, num_rate_params
- integer :: i_column, i_layer, i_elem, i_constituents, i_rate_params
+ integer :: num_columns, num_layers, num_constituents
+ integer :: i_column, i_layer, i_constituents
num_columns = size(constituents, dim=1)
num_layers = size(constituents, dim=2)
num_constituents = size(constituents, dim=3)
- num_rate_params = size(rate_params, dim=3)
- i_elem = 1
i_constituents = 1
- i_rate_params = 1
do i_layer = 1, num_layers
do i_column = 1, num_columns
- temperature(i_column, i_layer) = real(m_temperature(i_elem), kind_phys)
- pressure(i_column, i_layer) = real(m_pressure(i_elem), kind_phys)
- dry_air_density(i_column, i_layer) = real(m_dry_air_density(i_elem), kind_phys)
constituents(i_column, i_layer, :) &
- = real(m_constituents(i_constituents : i_constituents + num_constituents - 1), kind_phys)
- rate_params(i_column, i_layer, :) &
- = real(m_rate_params(i_rate_params : i_rate_params + num_rate_params - 1), kind_phys)
- i_elem = i_elem + 1
+ = real(micm_constituents(i_constituents : i_constituents + num_constituents - 1), kind_phys)
i_constituents = i_constituents + num_constituents
- i_rate_params = i_rate_params + num_rate_params
end do
end do
@@ -108,6 +84,7 @@ subroutine convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, const
real(kind_phys), intent(in) :: molar_mass_arr(:) ! kg mol-1
real(kind_phys), intent(inout) :: constituents(:,:,:) ! in: kg kg-1 | out: mol m-3
+ ! local variables
integer :: num_columns, num_layers, num_constituents
integer :: i_column, i_layer, i_elem
real(kind_phys) :: val
@@ -156,4 +133,4 @@ subroutine convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constit
end subroutine convert_to_mass_mixing_ratio
-end module micm_util
\ No newline at end of file
+end module musica_ccpp_micm_util
\ No newline at end of file
diff --git a/schemes/musica/musica_ccpp.F90 b/schemes/musica/musica_ccpp.F90
index bf4fae37..c486bee5 100644
--- a/schemes/musica/musica_ccpp.F90
+++ b/schemes/musica/musica_ccpp.F90
@@ -1,7 +1,7 @@
!> Top-level wrapper for MUSICA chemistry components
module musica_ccpp
- use musica_ccpp_micm, only : micm_register, micm_init, micm_run, micm_final
- use musica_ccpp_tuvx, only : tuvx_init, tuvx_run, tuvx_final
+ use musica_ccpp_micm, only: micm_register, micm_init, micm_run, micm_final
+ use musica_ccpp_tuvx, only: tuvx_init, tuvx_run, tuvx_final
implicit none
private
@@ -10,69 +10,89 @@ module musica_ccpp
contains
- subroutine musica_ccpp_register(constituents, solver_type, num_grid_cells, errmsg, errcode)
+ subroutine musica_ccpp_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode)
use ccpp_constituent_prop_mod, only: ccpp_constituent_properties_t
- type(ccpp_constituent_properties_t), allocatable, intent(out) :: constituents(:)
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(constituents, solver_type, num_grid_cells, errmsg, errcode)
+ call micm_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode)
end subroutine musica_ccpp_register
!> \section arg_table_musica_ccpp_init Argument Table
!! \htmlinclude musica_ccpp_init.html
- subroutine musica_ccpp_init(errmsg, errcode)
+ subroutine musica_ccpp_init(vertical_layer_dimension, vertical_interface_dimension, &
+ errmsg, errcode)
+ integer, intent(in) :: vertical_layer_dimension ! (count)
+ integer, intent(in) :: vertical_interface_dimension ! (count)
character(len=512), intent(out) :: errmsg
integer, intent(out) :: errcode
- call tuvx_init(errmsg, errcode)
+ call tuvx_init(vertical_layer_dimension, vertical_interface_dimension, &
+ errmsg, errcode)
+ if (errcode /= 0) return
call micm_init(errmsg, errcode)
+ if (errcode /= 0) return
end subroutine musica_ccpp_init
!> \section arg_table_musica_ccpp_run Argument Table
!! \htmlinclude musica_ccpp_run.html
subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, constituent_props, &
- constituents, rate_params, height, errmsg, errcode)
- use micm_util, only: reshape_into_micm_arr, reshape_into_ccpp_arr
- use micm_util, only: convert_to_mol_per_cubic_meter, convert_to_mass_mixing_ratio
+ constituents, geopotential_height_wrt_surface_at_midpoint, &
+ geopotential_height_wrt_surface_at_interface, surface_geopotential, &
+ standard_gravitational_acceleration, errmsg, errcode)
+ use musica_ccpp_micm_util, only: reshape_into_micm_arr, reshape_into_ccpp_arr
+ use musica_ccpp_micm_util, only: convert_to_mol_per_cubic_meter, convert_to_mass_mixing_ratio
use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t
use ccpp_kinds, only: kind_phys
use iso_c_binding, only: c_double
- real(kind_phys), intent(inout) :: time_step ! s
- real(kind_phys), target, intent(inout) :: temperature(:,:) ! K
- real(kind_phys), target, intent(inout) :: pressure(:,:) ! Pa
- real(kind_phys), target, intent(inout) :: dry_air_density(:,:) ! kg m-3
- type(ccpp_constituent_prop_ptr_t), intent(in) :: constituent_props(:)
- real(kind_phys), target, intent(inout) :: constituents(:,:,:) ! kg kg-1
- real(kind_phys), target, intent(inout) :: rate_params(:,:,:)
- real(kind_phys), target, intent(in) :: height(:,:) ! km
- character(len=512), intent(out) :: errmsg
- integer, intent(out) :: errcode
+
+ real(kind_phys), intent(in) :: time_step ! s
+ real(kind_phys), intent(in) :: temperature(:,:) ! K
+ real(kind_phys), intent(in) :: pressure(:,:) ! Pa
+ real(kind_phys), intent(in) :: dry_air_density(:,:) ! kg m-3
+ type(ccpp_constituent_prop_ptr_t), &
+ intent(in) :: constituent_props(:)
+ real(kind_phys), 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_geopotential(:) ! m2 s-2
+ real(kind_phys), intent(in) :: standard_gravitational_acceleration ! m s-2
+ character(len=512), intent(out) :: errmsg
+ integer, intent(out) :: errcode
! local variables
- real(c_double), target, dimension(size(temperature, dim=1) &
- * size(temperature, dim=2)) :: m_temperature
- real(c_double), target, dimension(size(pressure, dim=1) &
- * size(pressure, dim=2)) :: m_pressure
- real(c_double), target, dimension(size(dry_air_density, dim=1) &
- * size(dry_air_density, dim=2)) :: m_dry_air_density
- real(c_double), target, dimension(size(constituents, dim=1) &
- * size(constituents, dim=2) &
- * size(constituents, dim=3)) :: m_constituents ! mol m-3
- real(c_double), target, dimension(size(rate_params, dim=1) &
- * size(rate_params, dim=2) &
- * size(rate_params, dim=3)) :: m_rate_params
- real(kind_phys), target, dimension(size(constituents, dim=3)) :: molar_mass_arr ! kg mol-1
+ real(c_double), dimension(size(temperature, dim=1) &
+ * size(temperature, dim=2)) :: micm_temperature
+ real(c_double), dimension(size(pressure, dim=1) &
+ * size(pressure, dim=2)) :: micm_pressure
+ real(c_double), dimension(size(dry_air_density, dim=1) &
+ * size(dry_air_density, dim=2)) :: micm_dry_air_density
+ real(c_double), dimension(size(constituents, dim=1) &
+ * size(constituents, dim=2) &
+ * size(constituents, dim=3)) :: micm_constituents ! mol m-3
+ real(kind_phys), dimension(size(constituents, dim=3)) :: molar_mass_arr ! kg mol-1
+
+ ! temporarily dimensioned to Chapman mechanism until mapping between MICM and TUV-x is implemented
+ real(c_double), dimension(size(constituents, dim=1) &
+ * size(constituents, dim=2) &
+ * 3) :: photolysis_rate_constants ! s-1
integer :: i_elem
- call tuvx_run(height, temperature, dry_air_density, errmsg, errcode)
+ call tuvx_run(temperature, dry_air_density, &
+ geopotential_height_wrt_surface_at_midpoint, &
+ geopotential_height_wrt_surface_at_interface, &
+ surface_geopotential, &
+ standard_gravitational_acceleration, &
+ photolysis_rate_constants, &
+ errmsg, errcode)
- ! Get the molar_mass that is set in the call to instantiate()
+ ! Get the molar mass that is set in the call to instantiate()
do i_elem = 1, size(molar_mass_arr)
call constituent_props(i_elem)%molar_mass(molar_mass_arr(i_elem), errcode, errmsg)
if (errcode /= 0) then
@@ -82,7 +102,7 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co
end do
! TODO(jiwon) Check molar mass is non zero as it becomes a denominator for unit converison
- ! this code needs to go when ccpp framework does the check
+ ! this code will be deleted when the framework does the check
do i_elem = 1, size(molar_mass_arr)
if (molar_mass_arr(i_elem) <= 0) then
errcode = 1
@@ -94,16 +114,16 @@ subroutine musica_ccpp_run(time_step, temperature, pressure, dry_air_density, co
! Convert CAM-SIMA unit to MICM unit (kg kg-1 -> mol m-3)
call convert_to_mol_per_cubic_meter(dry_air_density, molar_mass_arr, constituents)
- ! Reshape array (3D -> 1D) and convert type (kind_phys -> c_double)
- call reshape_into_micm_arr(temperature, pressure, dry_air_density, constituents, rate_params, &
- m_temperature, m_pressure, m_dry_air_density, m_constituents, m_rate_params)
+ ! Reshape array (2D/3D -> 1D) and convert type (kind_phys -> c_double)
+ call reshape_into_micm_arr(temperature, pressure, dry_air_density, constituents, &
+ micm_temperature, micm_pressure, micm_dry_air_density, micm_constituents)
- call micm_run(time_step, m_temperature, m_pressure, m_dry_air_density, m_constituents, &
- m_rate_params, errmsg, errcode)
+ ! temporarily pass in unmapped photolysis rate constants until mapping between MICM and TUV-x is implemented
+ call micm_run(time_step, micm_temperature, micm_pressure, micm_dry_air_density, &
+ photolysis_rate_constants, micm_constituents, errmsg, errcode)
! Reshape array (1D -> 3D) and convert type (c_double -> kind_phys)
- call reshape_into_ccpp_arr(temperature, pressure, dry_air_density, constituents, rate_params, &
- m_temperature, m_pressure, m_dry_air_density, m_constituents, m_rate_params)
+ call reshape_into_ccpp_arr(micm_constituents, constituents)
! Convert MICM unit back to CAM-SIMA unit (mol m-3 -> kg kg-1)
call convert_to_mass_mixing_ratio(dry_air_density, molar_mass_arr, constituents)
@@ -113,8 +133,8 @@ end subroutine musica_ccpp_run
!> \section arg_table_musica_ccpp_final Argument Table
!! \htmlinclude musica_ccpp_final.html
subroutine musica_ccpp_final(errmsg, errcode)
- integer, intent(out) :: errcode
character(len=512), intent(out) :: errmsg
+ integer, intent(out) :: errcode
call tuvx_final(errmsg, errcode)
call micm_final(errmsg, errcode)
diff --git a/schemes/musica/musica_ccpp.meta b/schemes/musica/musica_ccpp.meta
index 30b7c7c0..f62966f9 100644
--- a/schemes/musica/musica_ccpp.meta
+++ b/schemes/musica/musica_ccpp.meta
@@ -1,12 +1,24 @@
[ccpp-table-properties]
name = musica_ccpp
type = scheme
- dependencies = micm/musica_ccpp_micm.F90,micm/musica_ccpp_tuvx.F90,musica_ccpp_util.F90
+ dependencies = micm/musica_ccpp_micm.F90,micm/musica_ccpp_micm_util.F90,tuvx/musica_ccpp_tuvx.F90,tuvx/musica_ccpp_tuvx_height_grid.F90,musica_ccpp_util.F90
dynamic_constituent_routine = musica_ccpp_register
[ccpp-arg-table]
name = musica_ccpp_init
type = scheme
+[ vertical_layer_dimension ]
+ standard_name = vertical_layer_dimension
+ units = none
+ type = integer
+ dimensions = ()
+ intent = in
+[ vertical_interface_dimension ]
+ standard_name = vertical_interface_dimension
+ units = none
+ type = integer
+ dimensions = ()
+ intent = in
[ errmsg ]
standard_name = ccpp_error_message
units = none
@@ -59,6 +71,30 @@
type = real | kind = kind_phys
dimensions = (horizontal_loop_extent,vertical_layer_dimension,number_of_ccpp_constituents)
intent = inout
+[ geopotential_height_wrt_surface_at_midpoint ]
+ standard_name = geopotential_height_wrt_surface
+ units = m
+ type = real | kind = kind_phys
+ dimensions = (horizontal_loop_extent,vertical_layer_dimension)
+ intent = in
+[ geopotential_height_wrt_surface_at_interface ]
+ standard_name = geopotential_height_wrt_surface_at_interface
+ units = m
+ type = real | kind = kind_phys
+ dimensions = (horizontal_loop_extent,vertical_interface_dimension)
+ intent = in
+[ surface_geopotential ]
+ standard_name = surface_geopotential
+ type = real | kind = kind_phys
+ units = m2 s-2
+ dimensions = (horizontal_loop_extent)
+ intent = in
+[ standard_gravitational_acceleration ]
+ standard_name = standard_gravitational_acceleration
+ units = m s-2
+ 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 9ab6b133..a8417876 100644
--- a/schemes/musica/musica_ccpp_namelist.xml
+++ b/schemes/musica/musica_ccpp_namelist.xml
@@ -97,7 +97,7 @@
filename_of_tuvx_configuration
none
- A configuration file for the TUVX photolysis rate calculator
+ A configuration file for the TUV-x photolysis rate calculator
UNSET_PATH
diff --git a/schemes/musica/musica_ccpp_util.F90 b/schemes/musica/musica_ccpp_util.F90
index 06c24964..5a4eda9c 100644
--- a/schemes/musica/musica_ccpp_util.F90
+++ b/schemes/musica/musica_ccpp_util.F90
@@ -15,10 +15,11 @@ module musica_ccpp_util
!> @param[out] error_code The CCPP error code.
!> @return True for an error, false for success.
logical function has_error_occurred(error, error_message, error_code)
- use musica_util, only : error_t
- type(error_t), intent(in) :: error
+ use musica_util, only: error_t
+
+ type(error_t), intent(in) :: error
character(len=512), intent(out) :: error_message
- integer, intent(out) :: error_code
+ integer, intent(out) :: error_code
character(len=30) :: error_code_str
@@ -33,6 +34,7 @@ logical function has_error_occurred(error, error_message, error_code)
error_message = '[MUSICA Error]: ' // error%category( ) // '[' // &
trim( adjustl( error_code_str ) ) // ']: ' // error%message( )
has_error_occurred = .true.
+
end function has_error_occurred
end module musica_ccpp_util
diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx.F90
index 2413a476..d39eb5f4 100644
--- a/schemes/musica/tuvx/musica_ccpp_tuvx.F90
+++ b/schemes/musica/tuvx/musica_ccpp_tuvx.F90
@@ -1,11 +1,10 @@
module musica_ccpp_tuvx
- use iso_c_binding
! Note: "tuvx_t" is included in an external pre-built tuvx library that the host
! model is responsible for linking to during compilation
- use musica_tuvx, only: tuvx_t
- use musica_ccpp_util, only: has_error_occurred
- use ccpp_kinds, only: kind_phys
+ use musica_tuvx, only: tuvx_t, grid_t
+ use musica_ccpp_util, only: has_error_occurred
+ use ccpp_kinds, only: kind_phys
use musica_ccpp_namelist, only: filename_of_tuvx_configuration
implicit none
@@ -13,15 +12,21 @@ module musica_ccpp_tuvx
public :: tuvx_init, tuvx_run, tuvx_final
- type(tuvx_t), pointer :: tuvx => null( )
+ type(tuvx_t), pointer :: tuvx => null( )
+ type(grid_t), pointer :: height_grid => null( )
contains
- !> Intitialize TUVX
- subroutine tuvx_init(errmsg, errcode)
+ !> Intitialize TUV-x
+ subroutine tuvx_init(vertical_layer_dimension, vertical_interface_dimension, &
+ errmsg, errcode)
use musica_tuvx, only: grid_map_t, profile_map_t, radiator_map_t
- use musica_util, only: error_t, mapping_t
+ use musica_util, only: error_t
+ use musica_ccpp_tuvx_height_grid, only: create_height_grid, &
+ height_grid_label, height_grid_unit
+ integer, intent(in) :: vertical_layer_dimension ! (count)
+ integer, intent(in) :: vertical_interface_dimension ! (count)
character(len=512), intent(out) :: errmsg
integer, intent(out) :: errcode
@@ -31,43 +36,120 @@ subroutine tuvx_init(errmsg, errcode)
type(radiator_map_t), pointer :: radiators
type(error_t) :: error
- errcode = 0
- errmsg = ''
-
grids => grid_map_t( error )
- if (has_error_occurred(error, errmsg, errcode)) return
+ 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 )
+ return
+ endif
+
+ call grids%add( height_grid, error )
+ if (has_error_occurred( error, errmsg, errcode )) then
+ deallocate( grids )
+ deallocate( height_grid )
+ height_grid => null()
+ return
+ end if
profiles => profile_map_t( error )
- if (has_error_occurred(error, errmsg, errcode)) return
-
- radiators =>radiator_map_t( error )
- if (has_error_occurred(error, errmsg, errcode)) return
-
- tuvx => tuvx_t( filename_of_tuvx_configuration, error )
- if (has_error_occurred(error, errmsg, errcode)) return
+ if (has_error_occurred( error, errmsg, errcode )) then
+ deallocate( grids )
+ deallocate( height_grid )
+ height_grid => null()
+ return
+ end if
+
+ radiators => radiator_map_t( error )
+ if (has_error_occurred( error, errmsg, errcode )) then
+ deallocate( grids )
+ deallocate( profiles )
+ deallocate( height_grid )
+ height_grid => null()
+ return
+ end if
+
+ tuvx => tuvx_t( filename_of_tuvx_configuration, grids, profiles, &
+ radiators, error )
+ if (has_error_occurred( error, errmsg, errcode )) then
+ deallocate( grids )
+ deallocate( profiles )
+ deallocate( radiators )
+ deallocate( height_grid )
+ height_grid => null()
+ return
+ end if
deallocate( grids )
deallocate( profiles )
deallocate( radiators )
- deallocate( tuvx )
+ deallocate( height_grid )
+ height_grid => null()
+
+ grids => tuvx%get_grids( error )
+ if (has_error_occurred( error, errmsg, errcode )) then
+ deallocate( tuvx )
+ tuvx => null()
+ return
+ end if
+
+ height_grid => grids%get( height_grid_label, height_grid_unit, error )
+ if (has_error_occurred( error, errmsg, errcode )) then
+ deallocate( tuvx )
+ tuvx => null()
+ deallocate( grids )
+ return
+ end if
+
+ deallocate( grids )
end subroutine tuvx_init
!> Calculates photolysis rate constants for the current model conditions
- subroutine tuvx_run( height, temperature, dry_air_density, errmsg, errcode )
- use musica_util, only: error_t
-
- real(kind_phys), intent(in) :: height(:,:) ! km (layer, column)
- real(kind_phys), intent(in) :: temperature(:,:) ! K (layer, column)
- real(kind_phys), intent(in) :: dry_air_density(:,:) ! molecule cm-3 (layer, column)
+ subroutine tuvx_run(temperature, dry_air_density, &
+ geopotential_height_wrt_surface_at_midpoint, &
+ geopotential_height_wrt_surface_at_interface, &
+ surface_geopotential, &
+ standard_gravitational_acceleration, &
+ photolysis_rate_constants, errmsg, errcode)
+ use musica_util, only: error_t
+ use musica_ccpp_tuvx_height_grid, only: set_height_grid_values, calculate_heights
+
+ 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) :: 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_geopotential(:) ! m2 s-2
+ real(kind_phys), intent(in) :: standard_gravitational_acceleration ! m s-2
+ ! temporarily set to Chapman mechanism and 1 dimension
+ ! until mapping between MICM and TUV-x is implemented
+ real(kind_phys), intent(out) :: photolysis_rate_constants(:) ! s-1 (column, reaction)
character(len=512), intent(out) :: errmsg
integer, intent(out) :: errcode
! local variables
- type(error_t) :: error
-
- errcode = 0
- errmsg = ''
+ 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
+ integer :: i_col
+
+ reciprocal_of_gravitational_acceleration = 1.0_kind_phys / standard_gravitational_acceleration
+
+ do i_col = 1, size(temperature, dim=1)
+ call calculate_heights( geopotential_height_wrt_surface_at_midpoint(i_col,:), &
+ geopotential_height_wrt_surface_at_interface(i_col,:), &
+ surface_geopotential(i_col), &
+ reciprocal_of_gravitational_acceleration, &
+ height_midpoints, height_interfaces )
+ call set_height_grid_values( height_grid, height_midpoints, height_interfaces, &
+ errmsg, errcode )
+ if (errcode /= 0) return
+ end do
+
+ ! stand-in until actual photolysis rate constants are calculated
+ photolysis_rate_constants(:) = 1.0e-6_kind_phys
end subroutine tuvx_run
@@ -76,8 +158,18 @@ subroutine tuvx_final(errmsg, errcode)
character(len=512), intent(out) :: errmsg
integer, intent(out) :: errcode
- errcode = 0
errmsg = ''
+ errcode = 0
+
+ if (associated( height_grid )) then
+ deallocate( height_grid )
+ height_grid => null()
+ end if
+
+ if (associated( tuvx )) then
+ deallocate( tuvx )
+ tuvx => null()
+ end if
end subroutine tuvx_final
diff --git a/schemes/musica/tuvx/musica_ccpp_tuvx_height_grid.F90 b/schemes/musica/tuvx/musica_ccpp_tuvx_height_grid.F90
new file mode 100644
index 00000000..cff0d2a2
--- /dev/null
+++ b/schemes/musica/tuvx/musica_ccpp_tuvx_height_grid.F90
@@ -0,0 +1,171 @@
+module musica_ccpp_tuvx_height_grid
+
+ implicit none
+
+ private
+ public :: create_height_grid, set_height_grid_values, calculate_heights
+
+ ! Conversions between the CAM-SIMA height grid and the TUV-x height grid
+ !
+ !-----------------------------------------------------------------------
+ ! Notes on the conversion between the host-model height grid and the TUV-x
+ !
+ ! TUV-x heights are "bottom-up" and require atmospheric constituent
+ ! concentrations at interfaces. Therefore, CAM-SIMA mid-points are used
+ ! as TUV-x grid interfaces, with an additional layer introduced between
+ ! the surface and the lowest CAM-SIMA mid-point, and a layer at the
+ ! top of the TUV-x grid to hold species densities above the top CAM-SIMA
+ ! mid-point.
+ !
+ ! Here,
+ ! - i_int is the index of an interface
+ ! - i_mid is the index of a mid-point
+ ! - pver is the CCPP vertical_layer_dimension
+ ! - pver+1 is the CCPP vertical_interface_dimension
+
+ ! ---- (interface) ===== (mid-point)
+ !
+ ! CAM TUV-x
+ ! ************************ (exo values) *****************************
+ ! ------(top)------ i_int = 1 -------(top)------ i_int = pver + 2
+ ! ================== i_mid = pver + 1
+ ! ================= i_mid = 1 ------------------ i_int = pver + 1
+ ! ----------------- i_int = 2 ================== i_mid = pver
+ ! ------------------ i_int = pver
+ ! ||
+ ! || ||
+ ! ||
+ ! ----------------- i_int = pver
+ ! ================= i_imd = pver ------------------ i_int = 2
+ ! ================== i_mid = 1
+ ! -----(ground)---- i_int = pver+1 -----(ground)----- i_int = 1
+ !
+ !-----------------------------------------------------------------------
+
+ !> Label for height grid in TUV-x
+ character(len=*), parameter, public :: height_grid_label = "height"
+ !> Unit for height grid in TUV-x
+ character(len=*), parameter, public :: height_grid_unit = "km"
+
+contains
+
+ !> Creates a TUV-x height grid
+ function create_height_grid(vertical_layer_dimension, vertical_interface_dimension, &
+ errmsg, errcode) result(height_grid)
+ use musica_ccpp_util, only: has_error_occurred
+ use musica_tuvx_grid, only: grid_t
+ use musica_util, only: error_t
+
+ integer, intent(in) :: vertical_layer_dimension ! (count)
+ integer, intent(in) :: vertical_interface_dimension ! (count)
+ character(len=*), intent(out) :: errmsg
+ integer, intent(out) :: errcode
+ type(grid_t), pointer :: height_grid
+
+ ! local variable
+ type(error_t) :: error
+
+ height_grid => null()
+ if ( vertical_layer_dimension < 1 ) then
+ errmsg = "[MUSICA Error] Invalid vertical_layer_dimension."
+ errcode = 1
+ return
+ end if
+ if ( vertical_interface_dimension - vertical_layer_dimension /= 1 ) then
+ errmsg = "[MUSICA Error] Invalid vertical_interface_dimension."
+ errcode = 1
+ return
+ end if
+ height_grid => grid_t( height_grid_label, height_grid_unit, &
+ vertical_interface_dimension, error )
+ if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+ end function create_height_grid
+
+ !> Sets TUV-x height grid values from the host-model height grid
+ subroutine set_height_grid_values(height_grid, host_midpoints, &
+ host_interfaces, errmsg, errcode)
+ use ccpp_kinds, only: kind_phys
+ use musica_ccpp_util, only: has_error_occurred
+ use musica_tuvx_grid, only: grid_t
+ use musica_util, only: error_t
+
+ type(grid_t), intent(inout) :: height_grid
+ real(kind_phys), intent(in) :: host_midpoints(:) ! km
+ real(kind_phys), intent(in) :: host_interfaces(:) ! km
+ character(len=*), intent(out) :: errmsg
+ integer, intent(out) :: errcode
+
+ ! local variables
+ type(error_t) :: error
+ real(kind_phys) :: midpoints(size(host_midpoints)+1)
+ real(kind_phys) :: interfaces(size(host_interfaces)+1)
+ integer :: n_host_midpoints, n_host_interfaces
+
+ if ( size(midpoints) /= height_grid%number_of_sections( error ) ) then
+ errmsg = "[MUSICA Error] Invalid size of TUV-x mid-point heights."
+ errcode = 1
+ return
+ end if
+
+ if ( size(interfaces) /= height_grid%number_of_sections( error ) + 1 ) then
+ errmsg = "[MUSICA Error] Invalid size of TUV-x interface heights."
+ errcode = 1
+ return
+ end if
+
+ n_host_midpoints = size(host_midpoints)
+ n_host_interfaces = size(host_interfaces)
+
+ interfaces(1) = host_interfaces(n_host_interfaces)
+ interfaces(2:n_host_interfaces) = host_midpoints(n_host_midpoints:1:-1)
+ interfaces(n_host_interfaces+1) = host_interfaces(1)
+
+ midpoints(1) = 0.5_kind_phys * &
+ ( host_midpoints(n_host_midpoints) + host_interfaces(n_host_interfaces) )
+ midpoints(2:n_host_midpoints) = host_interfaces(n_host_midpoints:2:-1)
+ midpoints(n_host_midpoints+1) = 0.5_kind_phys * &
+ ( interfaces(n_host_interfaces) + interfaces(n_host_interfaces+1) )
+
+ call height_grid%set_edges( interfaces, error )
+ if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+ call height_grid%set_midpoints( midpoints, error )
+ if ( has_error_occurred( error, errmsg, errcode ) ) return
+
+ end subroutine set_height_grid_values
+
+ !> Calculates the heights needed for the TUV-x grid based on available data
+ !!
+ !! Uses the reciprocal of gravitational acceleration, the surface geopotential,
+ !! and the geopotential height wrt the surface and midpoints/interfaces to calculate
+ !! the midpoint/interface height above sea level in km needed for the TUV-x grid.
+ !!
+ !! The equation used is taked from CAMChem
+ !! (see https://github.com/ESCOMP/CAM/blob/f0e489e9708ce7b91635f6d4997fbf1e390b0dbb/src/chemistry/mozart/mo_gas_phase_chemdr.F90#L514-L526)
+ subroutine calculate_heights(geopotential_height_wrt_surface_at_midpoint, &
+ geopotential_height_wrt_surface_at_interface, &
+ surface_geopotential, &
+ reciprocal_of_gravitational_acceleration, &
+ height_midpoints, height_interfaces)
+ use ccpp_kinds, only: kind_phys
+
+ real(kind_phys), intent(in) :: geopotential_height_wrt_surface_at_midpoint(:) ! m
+ real(kind_phys), intent(in) :: geopotential_height_wrt_surface_at_interface(:) ! m
+ real(kind_phys), intent(in) :: surface_geopotential ! m2 s-2
+ real(kind_phys), intent(in) :: reciprocal_of_gravitational_acceleration ! s2 m-1
+ real(kind_phys), intent(out) :: height_midpoints(:) ! km
+ real(kind_phys), intent(out) :: height_interfaces(:) ! km
+
+ ! local variable
+ real(kind_phys) :: surface_height ! m
+
+ surface_height = surface_geopotential * reciprocal_of_gravitational_acceleration
+ height_midpoints(:) = &
+ 0.001_kind_phys * ( geopotential_height_wrt_surface_at_midpoint(:) + surface_height )
+ height_interfaces(:) = &
+ 0.001_kind_phys * ( geopotential_height_wrt_surface_at_interface(:) + surface_height )
+
+ end subroutine calculate_heights
+
+end module musica_ccpp_tuvx_height_grid
\ No newline at end of file
diff --git a/test/cmake/TestUtils.cmake b/test/cmake/TestUtils.cmake
index 3f887db6..df31b9cd 100644
--- a/test/cmake/TestUtils.cmake
+++ b/test/cmake/TestUtils.cmake
@@ -3,6 +3,9 @@
if(CCPP_ENABLE_MEMCHECK)
find_program(MEMORYCHECK_COMMAND "valgrind")
+
+ # Set the Valgrind suppressions file for tests
+ set(MEMCHECK_SUPPRESS "--suppressions=${CMAKE_SOURCE_DIR}/valgrind.supp")
endif()
################################################################################
@@ -11,10 +14,10 @@ endif()
function(add_memory_check_test test_name test_binary test_args working_dir)
if(CCPP_ENABLE_MEMCHECK)
add_test(NAME memcheck_${test_name}
- COMMAND mpirun -v -np 1 ${MEMORYCHECK_COMMAND} --leak-check=full --error-exitcode=1 --trace-children=yes
+ COMMAND mpirun -v -np 1 ${MEMORYCHECK_COMMAND} --leak-check=full --error-exitcode=1 --trace-children=yes --gen-suppressions=all ${MEMCHECK_SUPPRESS}
${test_binary} ${test_args}
WORKING_DIRECTORY ${working_dir})
endif()
endfunction(add_memory_check_test)
-################################################################################
+################################################################################
\ No newline at end of file
diff --git a/test/docker/Dockerfile.musica b/test/docker/Dockerfile.musica
index 2ad24b71..65940f1a 100644
--- a/test/docker/Dockerfile.musica
+++ b/test/docker/Dockerfile.musica
@@ -1,6 +1,11 @@
+# This Dockerfile is designed for testing MUSICA CCPP functionality.
+# It includes:
+# - Unit tests for MUSICA utility functions
+# - Integration tests for MUSICA CCPP APIs
+
FROM ubuntu:22.04
-ARG MUSICA_GIT_TAG=aa0854ecee54bd7a5aeb7ea1ba0eebb2cf656146
+ARG MUSICA_GIT_TAG=dbbdb22f5f2807e27c2695db85291951ba178634
RUN apt update \
&& apt install -y sudo \
@@ -63,15 +68,18 @@ RUN cd musica \
COPY . atmospheric_physics
RUN sudo chown -R test_user:test_user atmospheric_physics
+# Must make ccpp-framework available before building test
RUN cd atmospheric_physics/test \
&& mkdir lib \
&& cd lib \
- && git clone -b add_const_interface --depth 1 https://github.com/peverwhee/ccpp-framework.git
+ && git clone -b develop --depth 1 https://github.com/NCAR/ccpp-framework.git
ENV CCPP_SRC_PATH="lib/ccpp-framework/src"
+ENV CCPP_FORTRAN_TOOLS_PATH="lib/ccpp-framework/scripts/fortran_tools"
RUN cd atmospheric_physics/test \
&& cmake -S . \
-B build \
+ -D CMAKE_BUILD_TYPE=Debug \
-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
new file mode 100644
index 00000000..b59489af
--- /dev/null
+++ b/test/docker/Dockerfile.musica.no_install
@@ -0,0 +1,80 @@
+# This Dockerfile is designed for testing MUSICA CCPP functionality.
+# It includes:
+# - Unit tests for MUSICA utility functions
+# - Integration tests for MUSICA CCPP APIs
+#
+# No MUSICA library installation is required, as the CMake FetchContent module
+# retrieves the MUSICA library automatically.
+
+FROM ubuntu:22.04
+
+ARG MUSICA_GIT_TAG=dbbdb22f5f2807e27c2695db85291951ba178634
+
+RUN apt update \
+ && apt install -y sudo \
+ && adduser test_user \
+ && echo "test_user ALL=(root) NOPASSWD: ALL" >> /etc/sudoers.d/test_user \
+ && chmod 0440 /etc/sudoers.d/test_user
+
+USER test_user
+WORKDIR /home/test_user
+
+RUN sudo apt update \
+ && sudo apt -y install \
+ cmake \
+ cmake-curses-gui \
+ curl \
+ g++ \
+ gcc \
+ gfortran \
+ git \
+ libblas-dev \
+ liblapack-dev \
+ lcov \
+ libcurl4-openssl-dev \
+ libhdf5-dev \
+ libnetcdff-dev \
+ libopenmpi-dev \
+ m4 \
+ make \
+ nlohmann-json3-dev \
+ openmpi-bin \
+ python3 \
+ tree \
+ valgrind \
+ vim \
+ 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}
+
+COPY . atmospheric_physics
+RUN sudo chown -R test_user:test_user atmospheric_physics
+
+# Must make ccpp-framework available before building test
+RUN cd atmospheric_physics/test \
+ && mkdir lib \
+ && cd lib \
+ && git clone -b develop --depth 1 https://github.com/NCAR/ccpp-framework.git
+ENV CCPP_SRC_PATH="lib/ccpp-framework/src"
+ENV CCPP_FORTRAN_TOOLS_PATH="lib/ccpp-framework/scripts/fortran_tools"
+
+RUN cd atmospheric_physics/test \
+ && cmake -S . \
+ -B build \
+ -D CMAKE_BUILD_TYPE=Debug \
+ -D CCPP_ENABLE_MUSICA_TESTS=ON \
+ -D CCPP_ENABLE_MEMCHECK=ON \
+ && cmake --build ./build
+
+RUN cd atmospheric_physics/test \
+ && mkdir third_party \
+ && cd third_party \
+ && git clone --depth 1 https://github.com/NCAR/tuv-x.git \
+ && cp -r tuv-x/data ../build/data
+
+WORKDIR /home/test_user/atmospheric_physics/test/build
\ No newline at end of file
diff --git a/test/musica/CMakeLists.txt b/test/musica/CMakeLists.txt
index 3fcdd843..487d181d 100644
--- a/test/musica/CMakeLists.txt
+++ b/test/musica/CMakeLists.txt
@@ -5,8 +5,9 @@ FetchContent_Declare(musica
GIT_REPOSITORY https://github.com/NCAR/musica.git
GIT_TAG $ENV{MUSICA_GIT_TAG} # Set by docker
)
+message(STATUS "Using MUSICA commit: $ENV{MUSICA_GIT_TAG}")
-set(MUSICA_BUILD_C_CXX_INTERFACE OFF)
+set(MUSICA_BUILD_C_CXX_INTERFACE ON)
set(MUSICA_BUILD_FORTRAN_INTERFACE ON)
set(MUSICA_ENABLE_TESTS OFF)
set(MUSICA_ENABLE_INSTALL OFF)
@@ -22,8 +23,9 @@ add_executable(test_musica_api test_musica_api.F90 musica_ccpp_namelist.F90)
target_sources(test_musica_api
PUBLIC
${MUSICA_SRC_PATH}/micm/musica_ccpp_micm.F90
- ${MUSICA_SRC_PATH}/micm/micm_util.F90
+ ${MUSICA_SRC_PATH}/micm/musica_ccpp_micm_util.F90
${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx.F90
+ ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_height_grid.F90
${MUSICA_SRC_PATH}/musica_ccpp.F90
${MUSICA_SRC_PATH}/musica_ccpp_util.F90
${CCPP_SRC_PATH}/ccpp_constituent_prop_mod.F90
@@ -63,4 +65,22 @@ add_custom_target(
${CMAKE_CURRENT_SOURCE_DIR}/tuvx/configs ${CMAKE_BINARY_DIR}/configs
)
-add_subdirectory(micm)
\ No newline at end of file
+add_subdirectory(micm)
+add_subdirectory(tuvx)
+
+# Test metdadata
+find_package(Python REQUIRED)
+
+file(MAKE_DIRECTORY ${CMAKE_BINARY_DIR}/metadata_test)
+add_custom_target(
+ copy_metadata_test_files ALL ${CMAKE_COMMAND} -E copy
+ ${CMAKE_CURRENT_SOURCE_DIR}/../../schemes/musica/musica_ccpp.meta
+ ${CMAKE_CURRENT_SOURCE_DIR}/../../schemes/musica/musica_ccpp.F90
+ ${CMAKE_BINARY_DIR}/metadata_test
+)
+
+add_test(
+ NAME test_metadata
+ COMMAND ${Python_EXECUTABLE} ${CMAKE_BINARY_DIR}/../$ENV{CCPP_FORTRAN_TOOLS_PATH}/offline_check_fortran_vs_metadata.py
+ --directory ${CMAKE_BINARY_DIR}/metadata_test
+)
\ No newline at end of file
diff --git a/test/musica/micm/CMakeLists.txt b/test/musica/micm/CMakeLists.txt
index d239da23..0c1c7853 100644
--- a/test/musica/micm/CMakeLists.txt
+++ b/test/musica/micm/CMakeLists.txt
@@ -1,9 +1,8 @@
-
add_executable(test_micm_util test_micm_util.F90)
target_sources(test_micm_util
PUBLIC
- ${MUSICA_SRC_PATH}/micm/micm_util.F90
+ ${MUSICA_SRC_PATH}/micm/musica_ccpp_micm_util.F90
${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90
)
diff --git a/test/musica/micm/test_micm_util.F90 b/test/musica/micm/test_micm_util.F90
index 76f2c139..5581d5d8 100644
--- a/test/musica/micm/test_micm_util.F90
+++ b/test/musica/micm/test_micm_util.F90
@@ -1,8 +1,6 @@
program test_micm_util
- use iso_c_binding
- use micm_util
- use ccpp_kinds, only: kind_phys
+ use musica_ccpp_micm_util
implicit none
@@ -12,33 +10,28 @@ program test_micm_util
call test_reshape()
call test_unit_conversion()
-
contains
subroutine test_reshape()
use iso_c_binding, only: c_double
use ccpp_kinds, only: kind_phys
- integer, parameter :: NUM_SPECIES = 4
- integer, parameter :: NUM_RATES = 3
- integer, parameter :: NUM_COLUMNS = 2
- integer, parameter :: NUM_LAYERS = 2
- integer, parameter :: NUM_GRID_CELLS = 4
- real(kind_phys), target :: temperature(NUM_COLUMNS,NUM_LAYERS)
- real(kind_phys), target :: pressure(NUM_COLUMNS,NUM_LAYERS)
- real(kind_phys), target :: dry_air_density(NUM_COLUMNS,NUM_LAYERS)
- real(kind_phys), target :: constituents(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES)
- real(kind_phys), target :: rate_params(NUM_COLUMNS,NUM_LAYERS,NUM_RATES)
- real(c_double), target :: m_temperature(NUM_GRID_CELLS)
- real(c_double), target :: m_pressure(NUM_GRID_CELLS)
- real(c_double), target :: m_dry_air_density(NUM_GRID_CELLS)
- real(c_double), target :: m_constituents(NUM_GRID_CELLS*NUM_SPECIES)
- real(c_double), target :: m_rate_params(NUM_GRID_CELLS*NUM_RATES)
+ integer, parameter :: NUM_SPECIES = 4
+ integer, parameter :: NUM_COLUMNS = 2
+ integer, parameter :: NUM_LAYERS = 2
+ integer, parameter :: NUM_GRID_CELLS = 4
+ real(kind_phys) :: temperature(NUM_COLUMNS,NUM_LAYERS)
+ real(kind_phys) :: pressure(NUM_COLUMNS,NUM_LAYERS)
+ real(kind_phys) :: dry_air_density(NUM_COLUMNS,NUM_LAYERS)
+ real(kind_phys) :: constituents(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES)
+ real(c_double) :: micm_temperature(NUM_GRID_CELLS)
+ real(c_double) :: micm_pressure(NUM_GRID_CELLS)
+ real(c_double) :: micm_dry_air_density(NUM_GRID_CELLS)
+ real(c_double) :: micm_constituents(NUM_GRID_CELLS*NUM_SPECIES)
! local variables
real(c_double), dimension(NUM_GRID_CELLS) :: arr_conditions
real(c_double), dimension(NUM_GRID_CELLS*NUM_SPECIES) :: arr_constituents
- real(c_double), dimension(NUM_GRID_CELLS*NUM_RATES) :: arr_rates
integer :: i_column, i_layer, i_elem, i_arr
real :: abs_error = 1e-7
@@ -52,44 +45,23 @@ subroutine test_reshape()
constituents(1,2,:) = (/ 0.41_kind_phys, 0.42_kind_phys, 0.43_kind_phys, 0.44_kind_phys /)
constituents(2,1,:) = (/ 0.21_kind_phys, 0.22_kind_phys, 0.23_kind_phys, 0.24_kind_phys /)
constituents(2,2,:) = (/ 0.31_kind_phys, 0.32_kind_phys, 0.33_kind_phys, 0.34_kind_phys /)
- rate_params(1,1,:) = (/ 100._kind_phys, 200._kind_phys, 300._kind_phys /)
- rate_params(1,2,:) = (/ 400._kind_phys, 500._kind_phys, 600._kind_phys /)
- rate_params(2,1,:) = (/ 700._kind_phys, 800._kind_phys, 900._kind_phys /)
- rate_params(2,2,:) = (/ 1000._kind_phys, 1100._kind_phys, 1200._kind_phys /)
-
arr_conditions = (/ 100.0, 200.0, 300.0, 400.0 /)
arr_constituents = (/ 0.1, 0.2, 0.3, 0.4, 0.21, 0.22, 0.23, 0.24, 0.41, 0.42, 0.43, 0.44, 0.31, 0.32, 0.33, 0.34 /)
- arr_rates = (/ 100., 200., 300., 700., 800., 900., 400., 500., 600., 1000., 1100., 1200. /)
call reshape_into_micm_arr(temperature, pressure, dry_air_density, constituents, &
- rate_params, m_temperature, m_pressure, m_dry_air_density, m_constituents, m_rate_params)
+ micm_temperature, micm_pressure, micm_dry_air_density, micm_constituents)
do i_elem = 1, NUM_GRID_CELLS
- ASSERT(m_temperature(i_elem) == arr_conditions(i_elem))
- ASSERT(m_pressure(i_elem) == arr_conditions(i_elem))
- ASSERT(m_dry_air_density(i_elem) == arr_conditions(i_elem))
+ ASSERT(micm_temperature(i_elem) == arr_conditions(i_elem))
+ ASSERT(micm_pressure(i_elem) == arr_conditions(i_elem))
+ ASSERT(micm_dry_air_density(i_elem) == arr_conditions(i_elem))
end do
- do i_elem = 1, size(m_constituents)
- ASSERT_NEAR(m_constituents(i_elem), arr_constituents(i_elem), abs_error)
+ do i_elem = 1, size(micm_constituents)
+ ASSERT_NEAR(micm_constituents(i_elem), arr_constituents(i_elem), abs_error)
end do
- do i_elem = 1, size(m_rate_params)
- ASSERT_NEAR(m_rate_params(i_elem), arr_rates(i_elem), abs_error)
- end do
-
- call reshape_into_ccpp_arr(temperature, pressure, dry_air_density, constituents, &
- rate_params, m_temperature, m_pressure, m_dry_air_density, m_constituents, m_rate_params)
-
- i_elem = 1
- do i_layer = 1, NUM_LAYERS
- do i_column = 1, NUM_COLUMNS
- ASSERT(temperature(i_column, i_layer) == arr_conditions(i_elem))
- ASSERT(pressure(i_column, i_layer) == arr_conditions(i_elem))
- ASSERT(dry_air_density(i_column, i_layer) == arr_conditions(i_elem))
- i_elem = i_elem + 1
- end do
- end do
+ call reshape_into_ccpp_arr(micm_constituents, constituents)
i_arr = 1
do i_layer = 1, NUM_LAYERS
@@ -101,31 +73,22 @@ subroutine test_reshape()
end do
end do
- i_arr = 1
- do i_layer = 1, NUM_LAYERS
- do i_column = 1, NUM_COLUMNS
- do i_elem = 1, NUM_RATES
- ASSERT_NEAR(rate_params(i_column, i_layer, i_elem), arr_rates(i_arr), abs_error)
- i_arr = i_arr + 1
- end do
- end do
- end do
-
end subroutine test_reshape
subroutine test_unit_conversion()
- use ccpp_kinds, only: kind_phys
+ use ccpp_kinds, only: kind_phys
+
+ integer, parameter :: NUM_COLUMNS = 2
+ integer, parameter :: NUM_LAYERS = 2
+ integer, parameter :: NUM_SPECIES = 4
+ real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: dry_air_density ! kg m-3
+ real(kind_phys), dimension(NUM_SPECIES) :: molar_mass_arr
+ real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: constituents
+ real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: ccpp_constituents ! kg kg-1
+ real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: micm_constituents ! mol m-3
+ integer :: i_column, i_layer, i_elem
+ real :: abs_error = 1e-3
- integer, parameter :: NUM_COLUMNS = 2
- integer, parameter :: NUM_LAYERS = 2
- integer, parameter :: NUM_SPECIES = 4
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS) :: dry_air_density ! kg m-3
- real(kind_phys), target, dimension(NUM_SPECIES) :: molar_mass_arr
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: constituents ! kg kg-1
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: ccpp_constituents ! mol m-3
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: micm_constituents ! mol m-3
- integer :: i_column, i_layer, i_elem
- real :: abs_error = 1e-3
dry_air_density(:,1) = (/ 3.5_kind_phys, 4.5_kind_phys /)
dry_air_density(:,2) = (/ 5.5_kind_phys, 6.5_kind_phys /)
molar_mass_arr(:) = (/ 200._kind_phys, 200._kind_phys, 200._kind_phys, 200._kind_phys /)
diff --git a/test/musica/musica_ccpp_namelist.F90 b/test/musica/musica_ccpp_namelist.F90
index a3e507c9..69c803b3 100644
--- a/test/musica/musica_ccpp_namelist.F90
+++ b/test/musica/musica_ccpp_namelist.F90
@@ -3,9 +3,9 @@ module musica_ccpp_namelist
implicit none
- public :: filename_of_micm_configuration, filename_of_tuvx_configuration
private
-
+ public :: filename_of_micm_configuration, filename_of_tuvx_configuration
+
character(len=*), parameter :: filename_of_micm_configuration = 'chapman'
character(len=*), parameter :: filename_of_tuvx_configuration = 'configs/ts1_tsmlt.json'
diff --git a/test/musica/test_musica_api.F90 b/test/musica/test_musica_api.F90
index e69f514f..9fdcfd1f 100644
--- a/test/musica/test_musica_api.F90
+++ b/test/musica/test_musica_api.F90
@@ -1,130 +1,135 @@
-subroutine test_musica_ccpp_api()
-
-#define ASSERT(x) if (.not.(x)) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: x"; stop 1; endif
+program run_test_musica_ccpp
- use iso_c_binding
use musica_ccpp
- 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
implicit none
- integer, parameter :: NUM_SPECIES = 4
- integer, parameter :: NUM_RATES = 3
- integer, parameter :: NUM_COLUMNS = 2
- integer, parameter :: NUM_LAYERS = 2
- integer :: solver_type
- integer :: errcode
- character(len=512) :: errmsg
- real(kind_phys) :: time_step ! s
- real(kind_phys), dimension(NUM_COLUMNS,NUM_LAYERS) :: height ! km
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS) :: temperature ! K
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS) :: pressure ! Pa
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS) :: dry_air_density ! kg m-3
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS,NUM_SPECIES) :: constituents ! kg kg-1
- real(kind_phys), target, dimension(NUM_COLUMNS,NUM_LAYERS,NUM_RATES) :: user_defined_reaction_rates
- type(ccpp_constituent_prop_ptr_t), allocatable :: constituent_props_ptr(:)
-
- ! local variables
- type(ccpp_constituent_properties_t), allocatable, target :: constituent_props(:)
- type(ccpp_constituent_properties_t), pointer :: const_prop
- real(kind_phys) :: molar_mass
- character(len=512) :: species_name, units
- logical :: tmp_bool, is_advected
- integer :: num_grid_cells
- integer :: i
-
- solver_type = Rosenbrock
- num_grid_cells = NUM_COLUMNS * NUM_LAYERS
- time_step = 60._kind_phys
- temperature(:,1) = (/ 100._kind_phys, 200._kind_phys /)
- temperature(:,2) = (/ 300._kind_phys, 400._kind_phys /)
- pressure(:,1) = (/ 6000.04_kind_phys, 7000.04_kind_phys /)
- 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 /)
- constituents(1,1,:) = (/ 0.1_kind_phys, 0.2_kind_phys, 0.3_kind_phys, 0.4_kind_phys /)
- constituents(1,2,:) = (/ 0.41_kind_phys, 0.42_kind_phys, 0.43_kind_phys, 0.44_kind_phys /)
- constituents(2,1,:) = (/ 0.21_kind_phys, 0.22_kind_phys, 0.23_kind_phys, 0.24_kind_phys /)
- constituents(2,2,:) = (/ 0.31_kind_phys, 0.32_kind_phys, 0.33_kind_phys, 0.34_kind_phys /)
- user_defined_reaction_rates(1,1,:) = (/2.7e-19_kind_phys, 1.13e-9_kind_phys, 5.8e-8_kind_phys/)
- user_defined_reaction_rates(1,2,:) = (/2.7e-19_kind_phys, 1.13e-9_kind_phys, 5.8e-8_kind_phys/)
- user_defined_reaction_rates(2,1,:) = (/2.7e-19_kind_phys, 1.13e-9_kind_phys, 5.8e-8_kind_phys/)
- user_defined_reaction_rates(2,2,:) = (/2.7e-19_kind_phys, 1.13e-9_kind_phys, 5.8e-8_kind_phys/)
-
- call musica_ccpp_register(constituent_props, solver_type, num_grid_cells, errmsg, errcode)
- ASSERT(allocated(constituent_props))
- ASSERT(size(constituent_props) == NUM_SPECIES)
- do i = 1, size(constituent_props)
- ASSERT(constituent_props(i)%is_instantiated(errcode, errmsg))
- ASSERT(errcode == 0)
- call constituent_props(i)%standard_name(species_name, errcode, errmsg)
- ASSERT(errcode == 0)
- call constituent_props(i)%molar_mass(molar_mass, errcode, errmsg)
- ASSERT(errcode == 0)
- call constituent_props(i)%is_advected(is_advected, errcode, errmsg)
- ASSERT(errcode == 0)
- tmp_booL = (trim(species_name) == "O2" .and. molar_mass == 0.032_kind_phys .and. .not. is_advected) .or. &
- (trim(species_name) == "O" .and. molar_mass == 0.016_kind_phys .and. .not. is_advected) .or. &
- (trim(species_name) == "O1D" .and. molar_mass == 0.016_kind_phys .and. .not. is_advected) .or. &
- (trim(species_name) == "O3" .and. molar_mass == 0.048_kind_phys .and. is_advected)
- ASSERT(tmp_bool)
- call constituent_props(i)%units(units, errcode, errmsg)
- ASSERT(errcode == 0)
- ASSERT(trim(units) == 'kg kg-1')
- end do
- if (errcode /= 0) then
- write(*,*) errcode, trim(errmsg)
- stop 3
- end if
-
- allocate(constituent_props_ptr(size(constituent_props)))
- do i = 1, size(constituent_props)
- const_prop => constituent_props(i)
- call constituent_props_ptr(i)%set(const_prop, errcode, errmsg)
- end do
-
- call musica_ccpp_init(errmsg, errcode)
- if (errcode /= 0) then
- write(*,*) trim(errmsg)
- stop 3
- endif
-
- write(*,*) "[MUSICA INFO] Initial Time Step"
- write(*,fmt="(1x,f10.2)") time_step
- write(*,*) "[MUSICA INFO] Initial Temperature"
- write(*,fmt="(4(1x,f10.4))") temperature
- write(*,*) "[MUSICA INFO] Initial Pressure"
- write(*,fmt="(4(1x,f10.4))") pressure
- write(*,*) "[MUSICA INFO] Initial Concentrations"
- write(*,fmt="(4(3x,e13.6))") constituents
- write(*,*) "[MUSICA INFO] Initial User-defined Reaction Rates"
- write(*,fmt="(3(3x,e13.6))") user_defined_reaction_rates
-
- call musica_ccpp_run(time_step, temperature, pressure, dry_air_density, constituent_props_ptr, &
- constituents, user_defined_reaction_rates, height, errmsg, errcode)
- if (errcode /= 0) then
- write(*,*) trim(errmsg)
- stop 3
- endif
-
- write(*,*) "[MUSICA INFO] Solved Concentrations"
- write(*,fmt="(4(3x,e13.6))") constituents
- write(*,*) "[MUSICA INFO] Solved User-defined Reaction Rates"
- write(*,fmt="(3(3x,e13.6))") user_defined_reaction_rates
-
- call musica_ccpp_final(errmsg, errcode)
-
- if (errcode /= 0) then
- write(*,*) trim(errmsg)
- stop 3
- endif
-
-end subroutine test_musica_ccpp_api
+#define ASSERT(x) if (.not.(x)) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: x"; stop 1; endif
-program run_test_musica_ccpp
-implicit none
call test_musica_ccpp_api()
-end program
\ No newline at end of file
+
+contains
+
+ subroutine test_musica_ccpp_api()
+ 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
+
+ implicit none
+
+ integer, parameter :: NUM_SPECIES = 4
+ integer, parameter :: NUM_COLUMNS = 2
+ integer, parameter :: NUM_LAYERS = 2
+ integer :: solver_type
+ integer :: errcode
+ character(len=512) :: errmsg
+ real(kind_phys) :: time_step ! s
+ 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) :: 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
+ type(ccpp_constituent_prop_ptr_t), allocatable :: constituent_props_ptr(:)
+
+ ! local variables
+ type(ccpp_constituent_properties_t), allocatable, target :: constituent_props(:)
+ type(ccpp_constituent_properties_t), pointer :: const_prop
+ real(kind_phys) :: molar_mass
+ character(len=512) :: species_name, units
+ logical :: tmp_bool, is_advected
+ integer :: num_grid_cells
+ integer :: i
+
+ solver_type = Rosenbrock
+ num_grid_cells = NUM_COLUMNS * NUM_LAYERS
+ time_step = 60._kind_phys
+ geopotential_height_wrt_surface_at_midpoint(1,:) = (/ 2000.0_kind_phys, 500.0_kind_phys /)
+ geopotential_height_wrt_surface_at_midpoint(2,:) = (/ 2000.0_kind_phys, -500.0_kind_phys /)
+ geopotential_height_wrt_surface_at_interface(1,:) = (/ 3000.0_kind_phys, 1000.0_kind_phys, 0.0_kind_phys /)
+ geopotential_height_wrt_surface_at_interface(2,:) = (/ 3000.0_kind_phys, 500.0_kind_phys, -1500.0_kind_phys /)
+ surface_geopotential = (/ 100.0_kind_phys, 200.0_kind_phys /)
+ standard_gravitational_acceleration = 10.0_kind_phys
+ temperature(:,1) = (/ 100._kind_phys, 200._kind_phys /)
+ temperature(:,2) = (/ 300._kind_phys, 400._kind_phys /)
+ pressure(:,1) = (/ 6000.04_kind_phys, 7000.04_kind_phys /)
+ 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 /)
+ constituents(1,1,:) = (/ 0.1_kind_phys, 0.2_kind_phys, 0.3_kind_phys, 0.4_kind_phys /)
+ constituents(1,2,:) = (/ 0.41_kind_phys, 0.42_kind_phys, 0.43_kind_phys, 0.44_kind_phys /)
+ constituents(2,1,:) = (/ 0.21_kind_phys, 0.22_kind_phys, 0.23_kind_phys, 0.24_kind_phys /)
+ constituents(2,2,:) = (/ 0.31_kind_phys, 0.32_kind_phys, 0.33_kind_phys, 0.34_kind_phys /)
+
+ call musica_ccpp_register(solver_type, num_grid_cells, constituent_props, errmsg, errcode)
+ ASSERT(allocated(constituent_props))
+ ASSERT(size(constituent_props) == NUM_SPECIES)
+ do i = 1, size(constituent_props)
+ ASSERT(constituent_props(i)%is_instantiated(errcode, errmsg))
+ ASSERT(errcode == 0)
+ call constituent_props(i)%standard_name(species_name, errcode, errmsg)
+ ASSERT(errcode == 0)
+ call constituent_props(i)%molar_mass(molar_mass, errcode, errmsg)
+ ASSERT(errcode == 0)
+ call constituent_props(i)%is_advected(is_advected, errcode, errmsg)
+ ASSERT(errcode == 0)
+ tmp_bool = (trim(species_name) == "O2" .and. molar_mass == 0.032_kind_phys .and. .not. is_advected) .or. &
+ (trim(species_name) == "O" .and. molar_mass == 0.016_kind_phys .and. .not. is_advected) .or. &
+ (trim(species_name) == "O1D" .and. molar_mass == 0.016_kind_phys .and. .not. is_advected) .or. &
+ (trim(species_name) == "O3" .and. molar_mass == 0.048_kind_phys .and. is_advected)
+ ASSERT(tmp_bool)
+ call constituent_props(i)%units(units, errcode, errmsg)
+ ASSERT(errcode == 0)
+ ASSERT(trim(units) == 'kg kg-1')
+ end do
+ if (errcode /= 0) then
+ write(*,*) errcode, trim(errmsg)
+ stop 3
+ end if
+
+ allocate(constituent_props_ptr(size(constituent_props)))
+ do i = 1, size(constituent_props)
+ const_prop => constituent_props(i)
+ call constituent_props_ptr(i)%set(const_prop, errcode, errmsg)
+ end do
+
+ call musica_ccpp_init(NUM_LAYERS, NUM_LAYERS+1, errmsg, errcode)
+ if (errcode /= 0) then
+ write(*,*) trim(errmsg)
+ stop 3
+ endif
+
+ write(*,*) "[MUSICA INFO] Initial Time Step"
+ write(*,fmt="(1x,f10.2)") time_step
+ write(*,*) "[MUSICA INFO] Initial Temperature"
+ write(*,fmt="(4(1x,f10.4))") temperature
+ write(*,*) "[MUSICA INFO] Initial Pressure"
+ write(*,fmt="(4(1x,f10.4))") pressure
+ 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_geopotential, &
+ standard_gravitational_acceleration, errmsg, errcode)
+ if (errcode /= 0) then
+ write(*,*) trim(errmsg)
+ stop 3
+ endif
+
+ write(*,*) "[MUSICA INFO] Solved Concentrations"
+ write(*,fmt="(4(3x,e13.6))") constituents
+
+ call musica_ccpp_final(errmsg, errcode)
+
+ if (errcode /= 0) then
+ write(*,*) trim(errmsg)
+ stop 3
+ endif
+
+ end subroutine test_musica_ccpp_api
+
+end program run_test_musica_ccpp
\ No newline at end of file
diff --git a/test/musica/tuvx/CMakeLists.txt b/test/musica/tuvx/CMakeLists.txt
new file mode 100644
index 00000000..7b7daacd
--- /dev/null
+++ b/test/musica/tuvx/CMakeLists.txt
@@ -0,0 +1,26 @@
+add_executable(test_tuvx_height_grid test_tuvx_height_grid.F90)
+
+target_sources(test_tuvx_height_grid
+ PUBLIC
+ ${MUSICA_SRC_PATH}/tuvx/musica_ccpp_tuvx_height_grid.F90
+ ${MUSICA_SRC_PATH}/musica_ccpp_util.F90
+ ${CCPP_TEST_SRC_PATH}/ccpp_kinds.F90
+)
+
+target_link_libraries(test_tuvx_height_grid
+ PRIVATE
+ musica::musica-fortran
+)
+
+set_target_properties(test_tuvx_height_grid
+ PROPERTIES
+ LINKER_LANGUAGE Fortran
+)
+
+add_test(
+ NAME test_tuvx_height_grid
+ COMMAND $
+ WORKING_DIRECTORY ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}
+)
+
+add_memory_check_test(test_tuvx_height_grid $ "" ${CMAKE_RUNTIME_OUTPUT_DIRECTORY})
diff --git a/test/musica/tuvx/configs/ts1_tsmlt.json b/test/musica/tuvx/configs/ts1_tsmlt.json
index 7ea1806e..817ae96b 100644
--- a/test/musica/tuvx/configs/ts1_tsmlt.json
+++ b/test/musica/tuvx/configs/ts1_tsmlt.json
@@ -4,14 +4,6 @@
"cross section parameters file": "data/cross_sections/O2_parameters.txt"
},
"grids": [
- {
- "name": "height",
- "type": "equal interval",
- "units": "km",
- "begins at" : 0.0,
- "ends at" : 120.0,
- "cell delta" : 1.0
- },
{
"name": "wavelength",
"type": "from csv file",
diff --git a/test/musica/tuvx/test_tuvx_height_grid.F90 b/test/musica/tuvx/test_tuvx_height_grid.F90
new file mode 100644
index 00000000..c0f5c6ed
--- /dev/null
+++ b/test/musica/tuvx/test_tuvx_height_grid.F90
@@ -0,0 +1,103 @@
+program test_tuvx_height_grid
+
+ use musica_ccpp_tuvx_height_grid
+
+ implicit none
+
+#define ASSERT(x) if (.not.(x)) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: x"; stop 1; endif
+#define ASSERT_NEAR( a, b, abs_error ) if( (abs(a - b) >= abs_error) .and. (abs(a - b) /= 0.0) ) then; write(*,*) "Assertion failed[", __FILE__, ":", __LINE__, "]: a, b"; stop 1; endif
+
+ call test_create_height_grid()
+ call test_calculate_height_grid_values()
+
+contains
+
+ subroutine test_create_height_grid()
+ use musica_util, only: error_t
+ use musica_tuvx_grid, only: grid_t
+ use ccpp_kinds, only: kind_phys
+
+ integer, parameter :: NUM_HOST_MIDPOINTS = 2
+ integer, parameter :: NUM_HOST_INTERFACES = 3
+ real(kind_phys), target :: host_midpoints(NUM_HOST_MIDPOINTS) = [200.8_kind_phys, 100.6_kind_phys]
+ real(kind_phys), target :: host_interfaces(NUM_HOST_INTERFACES) = [250.3_kind_phys, 150.2_kind_phys, 50.1_kind_phys]
+ type(grid_t), pointer :: height_grid
+ character(len=512) :: errmsg
+ integer :: errcode
+ real(kind_phys) :: abs_error = 1e-5
+ integer :: i
+
+ ! local variables
+ real(kind_phys), dimension(NUM_HOST_MIDPOINTS+1) :: midpoints
+ real(kind_phys), dimension(NUM_HOST_INTERFACES+1) :: interfaces
+ type(error_t) :: error
+
+ height_grid => create_height_grid(-1, 0, errmsg, errcode)
+ ASSERT(errcode == 1)
+ ASSERT(.not. associated(height_grid))
+
+ height_grid => create_height_grid(3, 3, errmsg, errcode)
+ ASSERT(errcode == 1)
+ ASSERT(.not. associated(height_grid))
+
+ height_grid => create_height_grid(NUM_HOST_MIDPOINTS, NUM_HOST_INTERFACES, &
+ errmsg, errcode)
+ ASSERT(errcode == 0)
+ ASSERT(associated(height_grid))
+
+ call set_height_grid_values(height_grid, host_midpoints, host_interfaces, &
+ errmsg, errcode)
+ ASSERT(errcode == 0)
+
+ ASSERT(height_grid%number_of_sections(error) == NUM_HOST_MIDPOINTS + 1)
+ ASSERT(error%is_success())
+
+ call height_grid%get_midpoints(midpoints, error)
+ ASSERT(error%is_success())
+
+ call height_grid%get_edges(interfaces, error)
+ ASSERT(error%is_success())
+
+ ASSERT_NEAR(midpoints(1), (100.6 + 50.1) * 0.5, abs_error)
+ ASSERT_NEAR(midpoints(2), 150.2, abs_error)
+ ASSERT_NEAR(midpoints(3), (250.3 + 200.8) * 0.5, abs_error)
+ ASSERT_NEAR(interfaces(1), 50.1, abs_error)
+ ASSERT_NEAR(interfaces(2), 100.6, abs_error)
+ ASSERT_NEAR(interfaces(3), 200.8, abs_error)
+ ASSERT_NEAR(interfaces(4), 250.3, abs_error)
+
+ deallocate( height_grid )
+
+ end subroutine test_create_height_grid
+
+ subroutine test_calculate_height_grid_values()
+
+ use ccpp_kinds, only: kind_phys
+
+ integer, parameter :: NUM_LAYERS = 2
+ real(kind_phys), dimension(NUM_LAYERS) :: geopotential_height_wrt_surface_at_midpoint ! m
+ real(kind_phys), dimension(NUM_LAYERS+1) :: geopotential_height_wrt_surface_at_interface ! m
+ real(kind_phys) :: surface_geopotential ! m2 s-2
+ real(kind_phys) :: reciprocal_of_gravitational_acceleration ! s2 m-1
+ real(kind_phys), dimension(NUM_LAYERS) :: height_midpoints ! km
+ real(kind_phys), dimension(NUM_LAYERS+1) :: height_interfaces ! km
+
+ geopotential_height_wrt_surface_at_midpoint(:) = (/ 2000.0_kind_phys, 500.0_kind_phys /)
+ geopotential_height_wrt_surface_at_interface(:) = (/ 3000.0_kind_phys, 1000.0_kind_phys, 0.0_kind_phys /)
+ surface_geopotential = 100.0_kind_phys
+ reciprocal_of_gravitational_acceleration = 10.0_kind_phys
+
+ call calculate_heights(geopotential_height_wrt_surface_at_midpoint, &
+ geopotential_height_wrt_surface_at_interface, &
+ surface_geopotential, reciprocal_of_gravitational_acceleration, &
+ height_midpoints, height_interfaces)
+
+ ASSERT_NEAR(height_midpoints(1), 3.0, 1e-5)
+ ASSERT_NEAR(height_midpoints(2), 1.5, 1e-5)
+ ASSERT_NEAR(height_interfaces(1), 4.0, 1e-5)
+ ASSERT_NEAR(height_interfaces(2), 2.0, 1e-5)
+ ASSERT_NEAR(height_interfaces(3), 1.0, 1e-5)
+
+ end subroutine test_calculate_height_grid_values
+
+end program test_tuvx_height_grid
diff --git a/test/valgrind.supp b/test/valgrind.supp
new file mode 100644
index 00000000..2910823e
--- /dev/null
+++ b/test/valgrind.supp
@@ -0,0 +1,27 @@
+##############################################################
+#
+# MUSICA TUV-x suppressions
+#
+# TODO(jiwon) We are experiencing memory leak issues in certain
+# functions of TUV-x. It appears that these leaks occur only
+# occasionally during initialization. We believe it’s acceptable
+# to add a Valgrind suppression for now, and we will investigate
+# further if it becomes a significant concern.
+#
+##############################################################
+{
+ Suppress_MUSICA_TUV-x_Leak
+ Memcheck:Leak
+ obj:/usr/libexec/valgrind/vgpreload_memcheck-amd64-linux.so
+ fun:__musica_config_MOD_get_string
+ fun:__tuvx_radiator_aerosol_MOD_constructor
+ fun:__tuvx_radiator_factory_MOD_radiator_builder
+ fun:__tuvx_radiator_warehouse_MOD_constructor
+ fun:__tuvx_radiative_transfer_MOD_constructor
+ fun:__tuvx_core_MOD_constructor
+ fun:InternalCreateTuvx
+ fun:musica::TUVX::Create
+ fun:CreateTuvx
+ fun:__musica_tuvx_MOD_constructor
+ fun:__musica_ccpp_tuvx_MOD_tuvx_init
+}
\ No newline at end of file