From a316b0d6c23f0fa42456bdac0609e0b8cde85680 Mon Sep 17 00:00:00 2001 From: David Fillmore <1524012+dwfncar@users.noreply.github.com> Date: Thu, 12 Sep 2024 08:07:24 -0600 Subject: [PATCH] 196 backward euler fortran (#210) * Added Backward Euler definitions to musics/mich.hpp. * Added Backward Euler definitions micm.cpp. * Added unit test BackwardEulerStandardOrder. * Added unit test BackwardEulerVectorOrder. * Save. * Added Backward Euler fortran API. * Added accuracy parameter for ASSERT_NEAR. * Reduced time_step for Backward Euler. * Added time_step and test_accuracy arguments to test_mulitple_grid_cells. * Removed some prints from test_micm_api.F90. --------- Co-authored-by: David Fillmore --- fortran/micm.F90 | 8 ++- .../test_micm_api.F90 | 71 +++++++++++++++---- 2 files changed, 63 insertions(+), 16 deletions(-) diff --git a/fortran/micm.F90 b/fortran/micm.F90 index 1a818b8e..5f54d396 100644 --- 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 + public :: Rosenbrock, RosenbrockStandardOrder, BackwardEuler, BackwardEulerStandardOrder private !> Wrapper for c solver stats @@ -29,8 +29,10 @@ 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 :: Rosenbrock = 1 - enumerator :: RosenbrockStandardOrder = 2 + enumerator :: Rosenbrock = 1 + enumerator :: RosenbrockStandardOrder = 2 + enumerator :: BackwardEuler = 3 + enumerator :: BackwardEulerStandardOrder = 4 end enum interface diff --git a/fortran/test/fetch_content_integration/test_micm_api.F90 b/fortran/test/fetch_content_integration/test_micm_api.F90 index f027c407..442a6639 100644 --- a/fortran/test/fetch_content_integration/test_micm_api.F90 +++ b/fortran/test/fetch_content_integration/test_micm_api.F90 @@ -6,7 +6,7 @@ program test_micm_api use, intrinsic :: iso_c_binding use, intrinsic :: ieee_arithmetic use musica_micm, only: micm_t, solver_stats_t, get_micm_version - use musica_micm, only: Rosenbrock, RosenbrockStandardOrder + use musica_micm, only: Rosenbrock, RosenbrockStandardOrder, BackwardEuler, BackwardEulerStandardOrder use musica_util, only: assert, error_t, mapping_t, string_t, find_mapping_index #include "micm/util/error.hpp" @@ -29,6 +29,8 @@ program test_micm_api call test_api() call test_multiple_grid_cell_vector_Rosenbrock() call test_multiple_grid_cell_standard_Rosenbrock() + call test_multiple_grid_cell_vector_BackwardEuler() + call test_multiple_grid_cell_standard_BackwardEuler() contains @@ -167,14 +169,15 @@ subroutine test_api() end subroutine test_api - subroutine test_multiple_grid_cells(micm, NUM_GRID_CELLS) + subroutine test_multiple_grid_cells(micm, NUM_GRID_CELLS, time_step, test_accuracy) type(micm_t), pointer, intent(inout) :: micm integer, intent(in) :: NUM_GRID_CELLS + real(c_double), intent(in) :: time_step + real, intent(in) :: test_accuracy integer, parameter :: NUM_SPECIES = 6 integer, parameter :: NUM_USER_DEFINED_REACTION_RATES = 2 - real(c_double) :: time_step real(c_double), target :: temperature(NUM_GRID_CELLS) real(c_double), target :: pressure(NUM_GRID_CELLS) real(c_double), target :: air_density(NUM_GRID_CELLS) @@ -195,8 +198,6 @@ subroutine test_multiple_grid_cells(micm, NUM_GRID_CELLS) real :: temp type(ArrheniusReaction) :: r1, r2 - time_step = 200 - A_index = micm%species_ordering%index( "A", error ) ASSERT( error%is_success() ) B_index = micm%species_ordering%index( "B", error ) @@ -267,12 +268,12 @@ subroutine test_multiple_grid_cells(micm, NUM_GRID_CELLS) D = initial_D * exp( -k1 * time_step ) E = initial_D * (k1 / (k2 - k1)) * (exp(-k1 * time_step) - exp(-k2 * time_step)) F = initial_F + initial_D * (1.0 + (k1 * exp(-k2 * time_step) - k2 * exp(-k1 * time_step)) / (k2 - k1)) - ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+A_index), A, 5.0e-3) - ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+B_index), B, 5.0e-3) - ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+C_index), C, 5.0e-3) - ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+D_index), D, 5.0e-3) - ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+E_index), E, 5.0e-3) - ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+F_index), F, 5.0e-3) + ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+A_index), A, test_accuracy) + ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+B_index), B, test_accuracy) + ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+C_index), C, test_accuracy) + ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+D_index), D, test_accuracy) + ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+E_index), E, test_accuracy) + ASSERT_NEAR(concentrations((i_cell-1)*NUM_SPECIES+F_index), F, test_accuracy) end do end subroutine test_multiple_grid_cells @@ -283,11 +284,14 @@ subroutine test_multiple_grid_cell_vector_Rosenbrock() integer :: num_grid_cells type(error_t) :: error + real(c_double), parameter :: time_step = 200 + real, parameter :: test_accuracy = 5.0e-3 + num_grid_cells = 3 micm => micm_t( "configs/analytical", Rosenbrock, num_grid_cells, error ) ASSERT( error%is_success() ) - call test_multiple_grid_cells( micm, num_grid_cells ) + call test_multiple_grid_cells( micm, num_grid_cells, time_step, test_accuracy ) deallocate( micm ) @@ -299,14 +303,55 @@ subroutine test_multiple_grid_cell_standard_Rosenbrock() integer :: num_grid_cells type(error_t) :: error + real(c_double), parameter :: time_step = 200 + real, parameter :: test_accuracy = 5.0e-3 + num_grid_cells = 3 micm => micm_t( "configs/analytical", RosenbrockStandardOrder, num_grid_cells, error ) ASSERT( error%is_success() ) - call test_multiple_grid_cells( micm, num_grid_cells ) + call test_multiple_grid_cells( micm, num_grid_cells, time_step, test_accuracy ) deallocate( micm ) end subroutine test_multiple_grid_cell_standard_Rosenbrock + subroutine test_multiple_grid_cell_vector_BackwardEuler() + + type(micm_t), pointer :: micm + integer :: num_grid_cells + type(error_t) :: error + + real(c_double), parameter :: time_step = 10 + real, parameter :: test_accuracy = 0.1 + + num_grid_cells = 3 + micm => micm_t( "configs/analytical", BackwardEuler, num_grid_cells, error ) + ASSERT( error%is_success() ) + + call test_multiple_grid_cells( micm, num_grid_cells, time_step, test_accuracy ) + + deallocate( micm ) + + end subroutine test_multiple_grid_cell_vector_BackwardEuler + + subroutine test_multiple_grid_cell_standard_BackwardEuler() + + type(micm_t), pointer :: micm + integer :: num_grid_cells + type(error_t) :: error + + real(c_double), parameter :: time_step = 10 + real, parameter :: test_accuracy = 0.1 + + num_grid_cells = 3 + micm => micm_t( "configs/analytical", BackwardEulerStandardOrder, num_grid_cells, error ) + ASSERT( error%is_success() ) + + call test_multiple_grid_cells( micm, num_grid_cells, time_step, test_accuracy ) + + deallocate( micm ) + + end subroutine test_multiple_grid_cell_standard_BackwardEuler + end program