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

Add multi-grid cell solving to MUSICA API for MICM #192

Merged
merged 9 commits into from
Aug 13, 2024
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}}
Comment on lines +10 to +19
Copy link
Collaborator

Choose a reason for hiding this comment

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

My only concern with this is that USER_DEFINED isn't technically part of the camp configuration. However, I do think this is incredibly useful. We could remove photolysis, first order loss, emission, from the technical specification of the configuration because they are all just logical versions of user defined and it allows us to not have to do funny things with music box for integrated reaction rates.

I guess there's nothing actionable with that paragraph; I'm just putting thoughts out there

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Do you think we should add USER_DEFINED to the CAMP configuration? I agree, it might make sense to replace the photolysis, emissions and first-order loss with this one type.

},
{
"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"
boulderdaze marked this conversation as resolved.
Show resolved Hide resolved
},
{
"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
boulderdaze marked this conversation as resolved.
Show resolved Hide resolved
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
Loading