Skip to content

Commit

Permalink
Add multi-grid cell solving to MUSICA API for MICM (#192)
Browse files Browse the repository at this point in the history
* add multi-cell solving option

* draft fortran API updates for multiple grid cells

* add multi-cell solving test for fortran API

* test both solvers with multiple grid cells in fortran tests

* update chapman test conditions and configuration

* add multicell solving to python API

* add multi-cell solver test for python API

* address msvc errors
  • Loading branch information
mattldawson authored Aug 13, 2024
1 parent 92b3365 commit a030011
Show file tree
Hide file tree
Showing 18 changed files with 1,020 additions and 278 deletions.
18 changes: 15 additions & 3 deletions configs/analytical/reactions.json
Original file line number Diff line number Diff line change
Expand Up @@ -7,14 +7,26 @@
"reactions":
[
{
"type": "ARRHENIUS", "A": 0.00012, "B": 7, "C" : 75, "D": 50, "E": 0.5,
"reactants": {"B": {"qty": 1}},
"products": {"C": {"yield": 1}}
"type": "USER_DEFINED",
"MUSICA name" : "reaction 2",
"reactants": {"E": {"qty": 1}},
"products": {"F": {"yield": 1}}
},
{
"type": "USER_DEFINED",
"MUSICA name" : "reaction 1",
"reactants": {"D": {"qty": 1}},
"products": {"E": {"yield": 1}}
},
{
"type": "ARRHENIUS", "A": 0.004, "C" : 50,
"reactants": {"A": {"qty": 1}},
"products": {"B": {"yield": 1}}
},
{
"type": "ARRHENIUS", "A": 0.012, "B": -2, "C" : 75, "D": 50, "E": 1.0e-6,
"reactants": {"B": {"qty": 1}},
"products": {"C": {"yield": 1}}
}
]
}
Expand Down
9 changes: 9 additions & 0 deletions configs/analytical/species.json
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,15 @@
},
{
"name": "C", "type": "CHEM_SPEC"
},
{
"name": "D", "type": "CHEM_SPEC"
},
{
"name": "E", "type": "CHEM_SPEC"
},
{
"name": "F", "type": "CHEM_SPEC"
}
]
}
6 changes: 3 additions & 3 deletions configs/chapman/reactions.json
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@
"yield": 2.0
}
},
"MUSICA name": "R1"
"MUSICA name": "jO2"
},
{
"type": "ARRHENIUS",
Expand All @@ -37,7 +37,7 @@
"O": {},
"O2": {}
},
"MUSICA name": "R3"
"MUSICA name": "jO3->O"
},
{
"type": "ARRHENIUS",
Expand All @@ -62,7 +62,7 @@
"O1D": {},
"O2": {}
},
"MUSICA name": "R5"
"MUSICA name": "jO3->O1D"
},
{
"type": "ARRHENIUS",
Expand Down
75 changes: 40 additions & 35 deletions fortran/micm.F90
Original file line number Diff line number Diff line change
Expand Up @@ -51,23 +51,21 @@ subroutine delete_micm_c(micm, error) bind(C, name="DeleteMicm")
type(error_t_c), intent(inout) :: error
end subroutine delete_micm_c

subroutine micm_solve_c(micm, time_step, temperature, pressure, air_density, num_concentrations, concentrations, &
num_user_defined_reaction_rates, user_defined_reaction_rates, solver_state, solver_stats, error) &
subroutine micm_solve_c(micm, time_step, temperature, pressure, air_density, concentrations, &
user_defined_reaction_rates, solver_state, solver_stats, error) &
bind(C, name="MicmSolve")
use musica_util, only: string_t_c, error_t_c
import c_ptr, c_double, c_int, solver_stats_t_c
type(c_ptr), value, intent(in) :: micm
real(kind=c_double), value, intent(in) :: time_step
real(kind=c_double), value, intent(in) :: temperature
real(kind=c_double), value, intent(in) :: pressure
real(kind=c_double), value, intent(in) :: air_density
integer(kind=c_int), value, intent(in) :: num_concentrations
real(kind=c_double), intent(inout) :: concentrations(num_concentrations)
integer(kind=c_int), value, intent(in) :: num_user_defined_reaction_rates
real(kind=c_double), intent(inout) :: user_defined_reaction_rates(num_user_defined_reaction_rates)
type(string_t_c), intent(out) :: solver_state
type(solver_stats_t_c), intent(out) :: solver_stats
type(error_t_c), intent(inout) :: error
type(c_ptr), value, intent(in) :: micm
real(kind=c_double), value, intent(in) :: time_step
type(c_ptr), value, intent(in) :: temperature
type(c_ptr), value, intent(in) :: pressure
type(c_ptr), value, intent(in) :: air_density
type(c_ptr), value, intent(in) :: concentrations
type(c_ptr), value, intent(in) :: user_defined_reaction_rates
type(string_t_c), intent(out) :: solver_state
type(solver_stats_t_c), intent(out) :: solver_stats
type(error_t_c), intent(inout) :: error
end subroutine micm_solve_c

function get_micm_version_c() bind(C, name="MicmVersion")
Expand Down Expand Up @@ -244,28 +242,35 @@ function constructor(config_path, solver_type, num_grid_cells, error) result( t

end function constructor

subroutine solve(this, time_step, temperature, pressure, air_density, num_concentrations, concentrations, &
num_user_defined_reaction_rates, user_defined_reaction_rates, solver_state, solver_stats, error)
subroutine solve(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 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), intent(in) :: temperature
real(c_double), intent(in) :: pressure
real(c_double), intent(in) :: air_density
integer(c_int), intent(in) :: num_concentrations
real(c_double), intent(inout) :: concentrations(*)
integer(c_int), intent(in) :: num_user_defined_reaction_rates
real(c_double), intent(inout) :: 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

call micm_solve_c(this%ptr, time_step, temperature, pressure, air_density, num_concentrations, concentrations, &
num_user_defined_reaction_rates, user_defined_reaction_rates, solver_state_c, solver_stats_c, error_c)
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

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)

solver_state = string_t(solver_state_c)
solver_stats = solver_stats_t(solver_stats_c)
Expand Down
Loading

0 comments on commit a030011

Please sign in to comment.