Chapter 2#
An MICM Box Model Example in Fortran#
In this next MUSICA Fortran example, we will setup a MICM solver, starting with a set of MICM configuration files, and run the solver for a single integration time step.
The following three configuration files (config.json, species.json, and reactions.json)
should be saved in a subdirectory named configs/analytical relative directory
from which you plan to call the box model executable.
(You can find a copy of these files in the MUSICA repository at configs/analytical.)
The top-level MICM configuration config.json file
simply lists the configuration files to parse. In this case, these are the
chemical species configuration file species.json and
the reactions configuration file reactions.json.
The contents of the config.json file for this example are:
{ "camp-files": [ "species.json", "reactions.json" ] }
For this example, we will have a system of three chemical species
A, B, and C, defined in the JSON file species.json as follows:
{ "camp-data": [ { "name": "A", "type": "CHEM_SPEC" }, { "name": "B", "type": "CHEM_SPEC" }, { "name": "C", "type": "CHEM_SPEC" }, { "name": "D", "type": "CHEM_SPEC" }, { "name": "E", "type": "CHEM_SPEC" }, { "name": "F", "type": "CHEM_SPEC" } ] }
The reactions.json specifies a mechanism, or a set of reactions for the system.
Here, we will introduce two Arrhenius type reactions, the first
with B evolving to C, and specifying all five reaction parameters,
and the second reaction with A evolving to B and using only two reaction parameters.
We also include two reactions with rate constants that are provided by the host application at runtime. These types of reactions are useful when outside calculations are needed to determine the rate constants, such as in the case of photolysis reactions.
The reactions.json file for this example should look like this:
{ "camp-data": [ { "type": "MECHANISM", "name": "music box interactive configuration", "reactions": [ { "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 } } }, { "type": "USER_DEFINED", "MUSICA name": "reaction 1", "reactants": { "D": { "qty": 1 } }, "products": { "E": { "yield": 1 } } }, { "type": "USER_DEFINED", "MUSICA name": "reaction 2", "reactants": { "E": { "qty": 1 } }, "products": { "F": { "yield": 1 } } } ] } ] }
More information on MICM configurations and reactions can be found in the MICM documentation
To create a simple box model, save the following Fortran code to a file named micm_box_model.F90:
program test_micm_box_model use, intrinsic :: iso_c_binding use, intrinsic :: ieee_arithmetic use iso_fortran_env, only : real64 use musica_util, only: error_t, string_t, mapping_t use musica_micm, only: micm_t, solver_stats_t use musica_micm, only: Rosenbrock, RosenbrockStandardOrder implicit none call box_model_arrays() call box_model_c_ptrs() contains !> Runs a simple box model using the MICM solver and passing in fortran arrays subroutine box_model_arrays() character(len=256) :: config_path integer :: solver_type integer :: num_grid_cells real(real64), parameter :: GAS_CONSTANT = 8.31446261815324_real64 ! J mol-1 K-1 real(real64) :: time_step real(real64), target :: temperature(1) real(real64), target :: pressure(1) real(real64), target :: air_density(1) real(real64), target :: concentrations(1,6) real(real64), target :: user_defined_reaction_rates(1,2) type(string_t) :: solver_state type(solver_stats_t) :: solver_stats type(error_t) :: error type(micm_t), pointer :: micm integer :: i config_path = "configs/analytical" solver_type = RosenbrockStandardOrder num_grid_cells = 1 time_step = 200 temperature(1) = 273.0 pressure(1) = 1.0e5 air_density(:) = pressure(:) / (GAS_CONSTANT * temperature(:)) concentrations(1,:) = (/ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 /) user_defined_reaction_rates(1,:) = (/ 0.001, 0.002 /) write(*,*) "Creating MICM solver..." micm => micm_t(config_path, solver_type, num_grid_cells, error) do i = 1, micm%species_ordering%size() print *, "Species Name:", micm%species_ordering%name(i), & ", Index:", micm%species_ordering%index(i) end do write(*,*) "Solving starts..." call micm%solve(time_step, temperature, pressure, air_density, concentrations, & user_defined_reaction_rates, solver_state, solver_stats, error) write(*,*) "After solving, concentrations", concentrations deallocate( micm ) end subroutine box_model_arrays !> Runs a simple box model using the MICM solver and passing in C pointers subroutine box_model_c_ptrs() character(len=256) :: config_path integer :: solver_type integer :: num_grid_cells real(real64), parameter :: GAS_CONSTANT = 8.31446261815324_real64 ! J mol-1 K-1 real(real64) :: time_step real(real64), target :: temperature(1) real(real64), target :: pressure(1) real(real64), target :: air_density(1) real(real64), target :: concentrations(6) real(real64), target :: user_defined_reaction_rates(2) type(string_t) :: solver_state type(solver_stats_t) :: solver_stats type(error_t) :: error type(micm_t), pointer :: micm integer :: i config_path = "configs/analytical" solver_type = RosenbrockStandardOrder num_grid_cells = 1 time_step = 200 temperature(1) = 273.0 pressure(1) = 1.0e5 air_density(:) = pressure(:) / (GAS_CONSTANT * temperature(:)) concentrations = (/ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 /) user_defined_reaction_rates = (/ 0.001, 0.002 /) write(*,*) "Creating MICM solver..." micm => micm_t(config_path, solver_type, num_grid_cells, error) do i = 1, micm%species_ordering%size() print *, "Species Name:", micm%species_ordering%name(i), & ", Index:", micm%species_ordering%index(i) end do write(*,*) "Solving starts..." call micm%solve(time_step, c_loc(temperature), c_loc(pressure), & c_loc(air_density), c_loc(concentrations), & c_loc(user_defined_reaction_rates), & solver_state, solver_stats, error) write(*,*) "After solving, concentrations", concentrations deallocate( micm ) end subroutine box_model_c_ptrs end program test_micm_box_model
From the musica_util module we need the Fortran types
error_t, string_t, and mapping_t.
A pointer to a musica_micm::micm_t will serve as the interface to the MICM solver
(in the example the pointer name is micm).
Note that the config_path in the code sample has been set to configs/analytical,
so that subdirectory should be created relative to the main program and contain
the MICM JSON configuration files,
otherwise the config_path should be modified appropriately.
The initial species concentrations are initialized in the concentrations array,
which is an argument to the MICM solver.
Finally, a single time step solution is obtained through a call to micm%solve,
after which the updated concentrations may be displayed.
To build the example, follow the instructions in Chapter 1.
Assuming you name the executable micm_box_model, you can run the program as follows:
$ ./micm_box_model
Creating MICM solver...
Species Name:A, Index: 1
Species Name:B, Index: 2
Species Name:C, Index: 5
Species Name:D, Index: 3
Species Name:E, Index: 4
Species Name:F, Index: 6
Solving starts...
After solving, concentrations 0.38236272259073301 1.4676807523204496 0.67030703490468713 1.1155750798779909 1.1499565250888166 1.2141178852173222
$