Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

196 backward euler fortran #210

Merged
merged 12 commits into from
Sep 12, 2024
8 changes: 5 additions & 3 deletions fortran/micm.F90
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -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
Expand Down
82 changes: 69 additions & 13 deletions fortran/test/fetch_content_integration/test_micm_api.F90
Original file line number Diff line number Diff line change
Expand Up @@ -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"
Expand All @@ -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

Expand Down Expand Up @@ -170,14 +172,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)
Expand All @@ -199,8 +202,6 @@ subroutine test_multiple_grid_cells(micm, NUM_GRID_CELLS)
real :: temp
type(ArrheniusReaction) :: r1, r2

time_step = 200

A_index = find_mapping_index( micm%species_ordering, "A", found )
ASSERT( found )
B_index = find_mapping_index( micm%species_ordering, "B", found )
Expand Down Expand Up @@ -271,12 +272,19 @@ 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)
print *, 'A', concentrations((i_cell-1)*NUM_SPECIES+A_index), A
print *, 'B', concentrations((i_cell-1)*NUM_SPECIES+B_index), B
print *, 'C', concentrations((i_cell-1)*NUM_SPECIES+C_index), C
print *, 'D', concentrations((i_cell-1)*NUM_SPECIES+D_index), D
print *, 'E', concentrations((i_cell-1)*NUM_SPECIES+E_index), E
print *, 'F', concentrations((i_cell-1)*NUM_SPECIES+F_index), F
print *
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

please remove diagnostic output

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
Expand All @@ -287,11 +295,15 @@ 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 )
print *, 'test_multiple_grid_cells vector Rosenbrock'
call test_multiple_grid_cells( micm, num_grid_cells, time_step, test_accuracy )

deallocate( micm )

Expand All @@ -303,14 +315,58 @@ 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 )
print *, 'test_multiple_grid_cells standard Rosenbrock'
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

please remove diagnostic output

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() )

print *, 'test_multiple_grid_cells vector Backward Euler'
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() )

print *, 'test_multiple_grid_cells standard Backward Euler'
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
Loading