merge main

NCAR · Dec 19, 2024 · 8c2c785 · 8c2c785
2 parents 583a307 + 1debffb
commit 8c2c785
Show file tree

Hide file tree

Showing 39 changed files with 1,690 additions and 598 deletions.
diff --git a/.github/workflows/pip.yml b/.github/workflows/pip.yml
@@ -13,7 +13,7 @@ jobs:
       fail-fast: false
       matrix:
         platform: [windows-latest, macos-13, ubuntu-latest]
-        python-version: ["3.7", "3.8", "3.9", "3.10", "3.11", "3.12"]
+        python-version: ["3.8", "3.9", "3.10", "3.11", "3.12"]
         # python versions: https://devguide.python.org/versions/
 
     runs-on: ${{ matrix.platform }}

diff --git a/CITATION.cff b/CITATION.cff
@@ -64,4 +64,4 @@ number: 10
 page: "E1743 - E1760"
 doi: "10.1175/BAMS-D-19-0331.1"
 url: "https://journals.ametsoc.org/view/journals/bams/101/10/bamsD190331.xml"
-version: 0.8.0
+version: 0.9.0
diff --git a/CMakeLists.txt b/CMakeLists.txt
@@ -1,7 +1,7 @@
 cmake_minimum_required(VERSION 3.21)
 
 # must be on the same line so that pyproject.toml can correctly identify the version
-project(musica-distribution VERSION 0.8.1)
+project(musica-distribution VERSION 0.9.0)
 
 set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH};${PROJECT_SOURCE_DIR}/cmake)
 set(CMAKE_USER_MAKE_RULES_OVERRIDE ${CMAKE_MODULE_PATH}/SetDefaults.cmake)
@@ -30,7 +30,7 @@ option(MUSICA_BUILD_DOCS "Build the documentation" OFF)
 option(MUSICA_ENABLE_MICM "Enable MICM" ON)
 option(MUSICA_ENABLE_TUVX "Enable TUV-x" ON)
 
-set(MUSICA_SET_MICM_VECTOR_MATRIX_SIZE "1" CACHE STRING "Set MICM vector-ordered matrix dimension")
+set(MUSICA_SET_MICM_VECTOR_MATRIX_SIZE "4" CACHE STRING "Set MICM vector-ordered matrix dimension")
 
 cmake_dependent_option(
   MUSICA_ENABLE_PYTHON_LIBRARY "Adds pybind11, a lightweight header-only library that exposes C++ types in Python and vice versa" OFF "MUSICA_BUILD_C_CXX_INTERFACE" OFF)

diff --git a/cmake/dependencies.cmake b/cmake/dependencies.cmake
@@ -62,7 +62,7 @@ endif()
 if (MUSICA_ENABLE_MICM AND MUSICA_BUILD_C_CXX_INTERFACE)
 
   set_git_default(MICM_GIT_REPOSITORY https://github.com/NCAR/micm.git)
-  set_git_default(MICM_GIT_TAG b3c462a)
+  set_git_default(MICM_GIT_TAG v.3.7.0)
 
   FetchContent_Declare(micm
       GIT_REPOSITORY ${MICM_GIT_REPOSITORY}
@@ -86,14 +86,17 @@ if (MUSICA_ENABLE_TUVX AND MUSICA_BUILD_C_CXX_INTERFACE)
   set(TUVX_INSTALL_INCLUDE_DIR ${MUSICA_INSTALL_INCLUDE_DIR} CACHE STRING "" FORCE)
 
   set_git_default(TUVX_GIT_REPOSITORY https://github.com/NCAR/tuv-x.git)
-  set_git_default(TUVX_GIT_TAG v0.10.0)
+  set_git_default(TUVX_GIT_TAG v0.10.1)
 
   FetchContent_Declare(tuvx
     GIT_REPOSITORY ${TUVX_GIT_REPOSITORY}
     GIT_TAG ${TUVX_GIT_TAG}
     GIT_PROGRESS NOT ${FETCHCONTENT_QUIET}
   )
 
+  set(TUVX_ENABLE_TESTS OFF)
+  set(TUVX_ENABLE_REGRESSION_TESTS OFF)
+
   FetchContent_MakeAvailable(tuvx)
 endif()
 

diff --git a/docker/Dockerfile.fortran-intel b/docker/Dockerfile.fortran-intel
@@ -1,5 +1,5 @@
 # versions and sizes from here: https://hub.docker.com/r/intel/oneapi-hpckit/tags
-FROM intel/oneapi-hpckit:latest
+FROM intel/oneapi-hpckit:2024.0.1-devel-ubuntu22.04
 
 # Based off of this: https://dgpu-docs.intel.com/driver/installation.html#repository-public-key-used-for-package-and-repository-signing
 # however those docs (at the time of this writing are incorrect) and this is the correct url

diff --git a/fortran/micm.F90 b/fortran/micm.F90
@@ -10,7 +10,7 @@ module musica_micm
   implicit none
 
   public :: micm_t, solver_stats_t, get_micm_version
-  public :: Rosenbrock, RosenbrockStandardOrder, BackwardEuler, BackwardEulerStandardOrder
+  public :: UndefinedSolver, Rosenbrock, RosenbrockStandardOrder, BackwardEuler, BackwardEulerStandardOrder
   private
 
   !> Wrapper for c solver stats
@@ -28,6 +28,7 @@ module musica_micm
   ! We could use Fortran 2023 enum type feature if Fortran 2023 is supported
   ! https://fortran-lang.discourse.group/t/enumerator-type-in-bind-c-derived-type-best-practice/5947/2
   enum, bind(c)
+    enumerator :: UndefinedSolver            = 0
     enumerator :: Rosenbrock                 = 1
     enumerator :: RosenbrockStandardOrder    = 2
     enumerator :: BackwardEuler              = 3
@@ -135,10 +136,14 @@ end function get_user_defined_reaction_rates_ordering_c
   type :: micm_t
     type(mappings_t), pointer :: species_ordering => null()
     type(mappings_t), pointer :: user_defined_reaction_rates => null()
-    type(c_ptr), private         :: ptr = c_null_ptr
+    type(c_ptr), private      :: ptr = c_null_ptr
+    integer,     private      :: number_of_grid_cells = 0
+    integer,     private      :: solver_type = UndefinedSolver
   contains
     ! Solve the chemical system
-    procedure :: solve
+    procedure, private :: solve_arrays
+    procedure, private :: solve_c_ptrs
+    generic :: solve => solve_arrays, solve_c_ptrs
     ! Get species properties
     procedure :: get_species_property_string
     procedure :: get_species_property_double
@@ -191,9 +196,11 @@ function constructor(config_path, solver_type, num_grid_cells, error)  result( t
     use musica_util, only: error_t_c, error_t, copy_mappings
     type(micm_t), pointer         :: this
     character(len=*), intent(in)  :: config_path
-    integer(c_int), intent(in)    :: solver_type
-    integer(c_int), intent(in)    :: num_grid_cells
+    integer, intent(in)           :: solver_type
+    integer, intent(in)           :: num_grid_cells
     type(error_t), intent(inout)  :: error
+
+    ! local variables
     character(len=1, kind=c_char) :: c_config_path(len_trim(config_path)+1)
     integer                       :: n, i
     type(error_t_c)               :: error_c
@@ -206,7 +213,10 @@ function constructor(config_path, solver_type, num_grid_cells, error)  result( t
     end do
     c_config_path(n+1) = c_null_char
 
-    this%ptr = create_micm_c(c_config_path, solver_type, num_grid_cells, error_c)
+    this%number_of_grid_cells = num_grid_cells
+    this%solver_type = solver_type
+    this%ptr = create_micm_c( c_config_path, int(solver_type, kind=c_int), &
+                              int(num_grid_cells, kind=c_int), error_c )
     error = error_t(error_c)
     if (.not. error%is_success()) then
         deallocate(this)
@@ -233,41 +243,121 @@ function constructor(config_path, solver_type, num_grid_cells, error)  result( t
 
   end function constructor
 
-  subroutine solve(this, time_step, temperature, pressure, air_density, concentrations, &
-              user_defined_reaction_rates, solver_state, solver_stats, error)
+  !> Solves the chemical system
+  !!
+  !! This function accepts fortran arrays and checks their sizes
+  !! against the number of grid cells and the species/rate parameter ordering.
+  subroutine solve_arrays(this, time_step, temperature, pressure, air_density, &
+      concentrations, user_defined_reaction_rates, solver_state, solver_stats, error)
     use iso_c_binding, only: c_loc
+    use iso_fortran_env, only: real64
+    use musica_util, only: string_t, string_t_c, error_t_c, error_t
+    class(micm_t),                    intent(in)    :: this
+    real(real64),                     intent(in)    :: time_step
+    real(real64), target, contiguous, intent(in)    :: temperature(:)
+    real(real64), target, contiguous, intent(in)    :: pressure(:)
+    real(real64), target, contiguous, intent(in)    :: air_density(:)
+    real(real64), target, contiguous, intent(inout) :: concentrations(:,:)
+    real(real64), target, contiguous, intent(in)    :: user_defined_reaction_rates(:,:)
+    type(string_t),                   intent(out)   :: solver_state
+    type(solver_stats_t),             intent(out)   :: solver_stats
+    type(error_t),                    intent(out)   :: error
+
+    type(string_t_c)       :: solver_state_c
+    type(solver_stats_t_c) :: solver_stats_c
+    type(error_t_c)        :: error_c
+
+    if (size(temperature) .ne. this%number_of_grid_cells) then
+        error = error_t(1, "MICM_SOLVE", "Temperature array size does not match number of grid cells")
+        return
+    end if
+    if (size(pressure) .ne. this%number_of_grid_cells) then
+        error = error_t(1, "MICM_SOLVE", "Pressure array size does not match number of grid cells")
+        return
+    end if
+    if (size(air_density) .ne. this%number_of_grid_cells) then
+        error = error_t(1, "MICM_SOLVE", "Air density array size does not match number of grid cells")
+        return
+    end if
+    if (this%solver_type .eq. Rosenbrock .or. this%solver_type .eq. BackwardEuler) then
+        if (size(concentrations, 1) .ne. this%number_of_grid_cells) then
+            error = error_t(1, "MICM_SOLVE", "Concentrations array dimension 1 does not match number of grid cells")
+            return
+        end if
+        if (size(concentrations, 2) .ne. this%species_ordering%size()) then
+            error = error_t(1, "MICM_SOLVE", "Concentrations array dimension 2 does not match species ordering")
+            return
+        end if
+        if (size(user_defined_reaction_rates, 1) .ne. this%number_of_grid_cells) then
+            error = error_t(1, "MICM_SOLVE", "User defined reaction rates array dimension 1 does not match number of grid cells")
+            return
+        end if
+        if (size(user_defined_reaction_rates, 2) .ne. this%user_defined_reaction_rates%size()) then
+            error = error_t(1, "MICM_SOLVE", "User defined reaction rates array dimension 2 does not match user defined reaction rates ordering")
+            return
+        end if
+    else
+        if (size(concentrations, 1) .ne. this%species_ordering%size()) then
+            error = error_t(1, "MICM_SOLVE", "Concentrations array dimension 1 does not match species ordering")
+            return
+        end if
+        if (size(concentrations, 2) .ne. this%number_of_grid_cells) then
+            error = error_t(1, "MICM_SOLVE", "Concentrations array dimension 2 does not match number of grid cells")
+            return
+        end if
+        if (size(user_defined_reaction_rates, 1) .ne. this%user_defined_reaction_rates%size()) then
+            error = error_t(1, "MICM_SOLVE", "User defined reaction rates array dimension 1 does not match user defined reaction rates ordering")
+            return
+        end if
+        if (size(user_defined_reaction_rates, 2) .ne. this%number_of_grid_cells) then
+            error = error_t(1, "MICM_SOLVE", "User defined reaction rates array dimension 2 does not match number of grid cells")
+            return
+        end if
+    end if
+
+    call micm_solve_c(this%ptr, real(time_step, kind=c_double), c_loc(temperature), &
+                      c_loc(pressure), c_loc(air_density), c_loc(concentrations), &
+                      c_loc(user_defined_reaction_rates), &
+                      solver_state_c, solver_stats_c, error_c)
+
+    solver_state = string_t(solver_state_c)
+    solver_stats = solver_stats_t(solver_stats_c)
+    error = error_t(error_c)
+
+  end subroutine solve_arrays
+
+  !> Solves the chemical system
+  !!
+  !! This function accepts c pointers and does not check their sizes.
+  !! The user is responsible for ensuring the sizes are correct.
+  subroutine solve_c_ptrs(this, time_step, temperature, pressure, air_density, &
+      concentrations, user_defined_reaction_rates, solver_state, solver_stats, error)
+    use iso_fortran_env, only: real64
     use musica_util, only: string_t, string_t_c, error_t_c, error_t
-    class(micm_t)                         :: this
-    real(c_double),         intent(in)    :: time_step
-    real(c_double), target, intent(in)    :: temperature(:)
-    real(c_double), target, intent(in)    :: pressure(:)
-    real(c_double), target, intent(in)    :: air_density(:)
-    real(c_double), target, intent(inout) :: concentrations(:)
-    real(c_double), target, intent(in)    :: user_defined_reaction_rates(:)
-    type(string_t),         intent(out)   :: solver_state
-    type(solver_stats_t),   intent(out)   :: solver_stats
-    type(error_t),          intent(out)   :: error
+    class(micm_t),        intent(in)    :: this
+    real(real64),         intent(in)    :: time_step
+    type(c_ptr),          intent(in)    :: temperature
+    type(c_ptr),          intent(in)    :: pressure
+    type(c_ptr),          intent(in)    :: air_density
+    type(c_ptr),          intent(in)    :: concentrations
+    type(c_ptr),          intent(in)    :: user_defined_reaction_rates
+    type(string_t),       intent(out)   :: solver_state
+    type(solver_stats_t), intent(out)   :: solver_stats
+    type(error_t),        intent(out)   :: error
 
     type(string_t_c)       :: solver_state_c
     type(solver_stats_t_c) :: solver_stats_c
     type(error_t_c)        :: error_c
-    type(c_ptr)            :: temperature_c, pressure_c, air_density_c, concentrations_c, &
-                              user_defined_reaction_rates_c
-
-    temperature_c    = c_loc(temperature)
-    pressure_c       = c_loc(pressure)
-    air_density_c    = c_loc(air_density)
-    concentrations_c = c_loc(concentrations)
-    user_defined_reaction_rates_c = c_loc(user_defined_reaction_rates)
-    call micm_solve_c(this%ptr, time_step, temperature_c, pressure_c, air_density_c, &
-                      concentrations_c, user_defined_reaction_rates_c, solver_state_c, &
-                      solver_stats_c, error_c)
+
+    call micm_solve_c(this%ptr, real(time_step, kind=c_double), temperature, pressure, &
+                      air_density, concentrations, user_defined_reaction_rates, &
+                      solver_state_c, solver_stats_c, error_c)
 
     solver_state = string_t(solver_state_c)
     solver_stats = solver_stats_t(solver_stats_c)
     error = error_t(error_c)
 
-  end subroutine solve
+  end subroutine solve_c_ptrs
 
   !> Constructor for solver_stats_t object that takes ownership of solver_stats_t_c
   function solver_stats_t_constructor( c_solver_stats ) result( new_solver_stats )
@@ -359,8 +449,9 @@ end function solver_stats_t_solves
 
   !> Get the final time the solver iterated to
   function solver_stats_t_final_time( this ) result( final_time )
+    use iso_fortran_env, only: real64
     class(solver_stats_t), intent(in) :: this
-    real                              :: final_time
+    real(real64)                      :: final_time
 
     final_time = this%final_time_
 
@@ -382,15 +473,16 @@ function get_species_property_string(this, species_name, property_name, error) r
   end function get_species_property_string
 
   function get_species_property_double(this, species_name, property_name, error) result(value)
+    use iso_fortran_env, only: real64
     use musica_util, only: error_t_c, error_t, to_c_string
     class(micm_t)                :: this
     character(len=*), intent(in) :: species_name, property_name
     type(error_t), intent(inout) :: error
-    real(c_double)               :: value
+    real(real64)                 :: value
 
     type(error_t_c)              :: error_c
-    value = get_species_property_double_c(this%ptr, &
-              to_c_string(species_name), to_c_string(property_name), error_c)
+    value = real( get_species_property_double_c( this%ptr, to_c_string(species_name), &
+                  to_c_string(property_name), error_c ), kind=real64 )
     error = error_t(error_c)
   end function get_species_property_double
 
@@ -399,11 +491,11 @@ function get_species_property_int(this, species_name, property_name, error) resu
     class(micm_t)                :: this
     character(len=*), intent(in) :: species_name, property_name
     type(error_t), intent(inout) :: error
-    integer(c_int)               :: value
+    integer                      :: value
 
     type(error_t_c)              :: error_c
-    value = get_species_property_int_c(this%ptr, &
-              to_c_string(species_name), to_c_string(property_name), error_c)
+    value = int( get_species_property_int_c(this%ptr, &
+                 to_c_string(species_name), to_c_string(property_name), error_c) )
     error = error_t(error_c)
   end function get_species_property_int
 

diff --git a/fortran/test/fetch_content_integration/CMakeLists.txt b/fortran/test/fetch_content_integration/CMakeLists.txt
@@ -48,6 +48,8 @@ if (MUSICA_ENABLE_MICM)
       LINKER_LANGUAGE Fortran
   )
 
+  target_compile_definitions(test_micm_fortran_api PUBLIC MICM_VECTOR_MATRIX_SIZE=${MUSICA_SET_MICM_VECTOR_MATRIX_SIZE})
+
   add_test(
     NAME test_micm_fortran_api
     COMMAND $<TARGET_FILE:test_micm_fortran_api>