MICM API

namespace micm

Typedefs

using DenseMatrixVector = VectorMatrix<double, MICM_DEFAULT_VECTOR_SIZE>

using SparseMatrixVector = SparseMatrix<double, SparseMatrixVectorOrdering<MICM_DEFAULT_VECTOR_SIZE>>

using DenseMatrixStandard = Matrix<double>

using SparseMatrixStandard = SparseMatrix<double, SparseMatrixStandardOrdering>

using VectorState = State<DenseMatrixVector, SparseMatrixVector>

using StandardState = State<DenseMatrixStandard, SparseMatrixStandard>

using RosenbrockVectorType = typename RosenbrockSolverParameters::template SolverType<ProcessSet<DenseMatrixVector, SparseMatrixVector>, LinearSolver<SparseMatrixVector, LuDecomposition>, ConstraintSet<DenseMatrixVector, SparseMatrixVector>>

using Rosenbrock = Solver<RosenbrockVectorType, State<DenseMatrixVector, SparseMatrixVector>>

using RosenbrockStandardType = typename RosenbrockSolverParameters::template SolverType<ProcessSet<DenseMatrixStandard, SparseMatrixStandard>, LinearSolver<SparseMatrixStandard, LuDecomposition>, ConstraintSet<DenseMatrixStandard, SparseMatrixStandard>>

using RosenbrockStandard = Solver<RosenbrockStandardType, State<DenseMatrixStandard, SparseMatrixStandard>>

using BackwardEulerVectorType = typename BackwardEulerSolverParameters::template SolverType<ProcessSet<DenseMatrixVector, SparseMatrixVector>, LinearSolver<SparseMatrixVector, LuDecomposition>, ConstraintSet<DenseMatrixVector, SparseMatrixVector>>

using BackwardEuler = Solver<BackwardEulerVectorType, State<DenseMatrixVector, SparseMatrixVector>>

using BackwardEulerStandardType = typename BackwardEulerSolverParameters::template SolverType<ProcessSet<DenseMatrixStandard, SparseMatrixStandard>, LinearSolver<SparseMatrixStandard, LuDecomposition>, ConstraintSet<DenseMatrixStandard, SparseMatrixStandard>>

using BackwardEulerStandard = Solver<BackwardEulerStandardType, State<DenseMatrixStandard, SparseMatrixStandard>>

using RosenbrockThreeStageBuilder = CpuSolverBuilder<RosenbrockSolverParameters, DenseMatrixVector, SparseMatrixVector>

using BackwardEulerBuilder = CpuSolverBuilder<BackwardEulerSolverParameters, DenseMatrixVector, SparseMatrixVector, LuDecompositionDoolittle>

using CudaLuDecomposition = CudaLuDecompositionMozartInPlace

Alias for the default CUDA LU decomposition algorithm.

template<class SolverParametersPolicy, std::size_t L = MICM_DEFAULT_VECTOR_SIZE>
using CudaSolverBuilderInPlace = SolverBuilder<SolverParametersPolicy, CudaDenseMatrix<double, L>, CudaSparseMatrix<double, SparseMatrixVectorOrdering<L>>, CudaProcessSet<CudaDenseMatrix<double, L>, CudaSparseMatrix<double, SparseMatrixVectorOrdering<L>>>, CudaLuDecompositionMozartInPlace, CudaLinearSolverInPlace<CudaSparseMatrix<double, SparseMatrixVectorOrdering<L>>, CudaLuDecompositionMozartInPlace>, CudaState<CudaDenseMatrix<double, L>, CudaSparseMatrix<double, SparseMatrixVectorOrdering<L>>, CudaLuDecompositionMozartInPlace>>

Builder of CUDA-based general solvers.

GPU solvers only work with vector-ordered matrices

Template Parameters:

• SolverParametersPolicy – Policy for the ODE solver

• L – Vector size

using CudaDenseMatrixVector = CudaDenseMatrix<double, MICM_DEFAULT_VECTOR_SIZE>

using CudaSparseMatrixVector = CudaSparseMatrix<double, SparseMatrixVectorOrdering<MICM_DEFAULT_VECTOR_SIZE>>

using GpuState = CudaState<CudaDenseMatrixVector, CudaSparseMatrixVector, CudaLuDecompositionMozartInPlace>

using CudaRosenbrockVectorType = typename CudaRosenbrockSolverParameters::template SolverType<CudaProcessSet<CudaDenseMatrixVector, CudaSparseMatrixVector>, CudaLinearSolverInPlace<CudaSparseMatrixVector>, ConstraintSet<CudaDenseMatrixVector, CudaSparseMatrixVector>>

using CudaRosenbrock = Solver<CudaRosenbrockVectorType, GpuState>

using GpuRosenbrockThreeStageBuilder = CudaSolverBuilderInPlace<CudaRosenbrockSolverParameters>

using RateConstantVariant = std::variant<ArrheniusRateConstantParameters, TroeRateConstantParameters, TernaryChemicalActivationRateConstantParameters, BranchedRateConstantParameters, TunnelingRateConstantParameters, TaylorSeriesRateConstantParameters, ReversibleRateConstantParameters, UserDefinedRateConstantParameters, SurfaceRateConstantParameters, LambdaRateConstantParameters>

Value-typed union of all supported rate constant parameter types. Stored by value in ChemicalReaction; consumed at store-build time by ReactionRateConstantStore::BuildFrom. Never sent to GPU.

using LuDecomposition = LuDecompositionDoolittle

Alias for the default LU decomposition algorithm.

using LuDecompositionInPlace = LuDecompositionMozartInPlace

Alias for the default in-place LU decomposition algorithm.

template<class SolverParametersPolicy, class DenseMatrixPolicy = Matrix<double>, class SparseMatrixPolicy = SparseMatrix<double, SparseMatrixStandardOrdering>, class LuDecompositionPolicy = LuDecomposition, class LMatrixPolicy = SparseMatrixPolicy, class UMatrixPolicy = SparseMatrixPolicy>
using CpuSolverBuilder = SolverBuilder<SolverParametersPolicy, DenseMatrixPolicy, SparseMatrixPolicy, ProcessSet<DenseMatrixPolicy, SparseMatrixPolicy>, LuDecompositionPolicy, LinearSolver<SparseMatrixPolicy, LuDecompositionPolicy, LMatrixPolicy, UMatrixPolicy>, State<DenseMatrixPolicy, SparseMatrixPolicy, LuDecompositionPolicy, LMatrixPolicy, UMatrixPolicy>>

Builder of CPU-based general solvers.

Template Parameters:

• SolverParametersPolicy – Parameters for the ODE solver

• DenseMatrixPolicy – Policy for dense matrices

• SparseMatrixPolicy – Policy for sparse matrices

• LuDecompositionPolicy – Policy for the LU decomposition

• LMatrixPolicy – Policy for the Lower matrix

• UMatrixPolicy – Policy for the Upper matrix

template<class SolverParametersPolicy, class DenseMatrix = Matrix<double>, class SparseMatrixPolicy = SparseMatrix<double, SparseMatrixStandardOrdering>, class LuDecompositionPolicy = LuDecompositionInPlace>
using CpuSolverBuilderInPlace = SolverBuilder<SolverParametersPolicy, DenseMatrix, SparseMatrixPolicy, ProcessSet<DenseMatrix, SparseMatrixPolicy>, LuDecompositionPolicy, LinearSolverInPlace<SparseMatrixPolicy, LuDecompositionPolicy>, State<DenseMatrix, SparseMatrixPolicy, LuDecompositionPolicy>>

Builder of CPU-based general solvers with in-place LU decomposition.

Template Parameters:

• SolverParametersPolicy – Parameters for the ODE solver

• DenseMatrixPolicy – Policy for dense matrices

• SparseMatrixPolicy – Policy for sparse matrices

• LuDecompositionPolicy – Policy for the LU decomposition

using StandardDenseMatrix = Matrix<double>

using StandardSparseMatrix = SparseMatrix<double, SparseMatrixStandardOrdering>

using SparseMatrixStandardOrdering = SparseMatrixStandardOrderingCompressedSparseRow

Alias for the default sparse matrix standard ordering.

template<std::size_t L = MICM_DEFAULT_VECTOR_SIZE>
using SparseMatrixVectorOrdering = SparseMatrixVectorOrderingCompressedSparseRow<L>

Alias for the default sparse matrix vector ordering.

using DefaultVectorSparseMatrix = SparseMatrix<double, SparseMatrixVectorOrdering<MICM_DEFAULT_VECTOR_SIZE>>

template<typename T>
using ViewCategory_t = typename ViewCategory<std::remove_cvref_t<T>>::type

Helper alias.

template<typename T>
using GroupingStrategy_t = typename GroupingStrategy<std::remove_cvref_t<T>>::type

Helper alias.

Enums

enum class SolverState

The final state the solver was in after the Solve function finishes.

Values:

enumerator NotYetCalled

This is the initial value at the start of the Solve function.

enumerator Running

This is only used for control flow in the Solve function.

enumerator Converged

A successful integration will have this value.

enumerator ConvergenceExceededMaxSteps

If the number of steps exceeds the maximum value on the solver parameter, this value will be returned.

enumerator StepSizeTooSmall

Very stiff systems will likely result in a step size that is not useable for the solver.

enumerator RepeatedlySingularMatrix

Matrices that are singular more than once will set this value. At present, this should never be returned.

enumerator NaNDetected

Mostly this value is returned by systems that tend toward chemical explosions.

enumerator InfDetected

Can happen when unititialized memory is used in the solver.

enumerator AcceptingUnconvergedIntegration

Used for backward euler. This allows us to “succeed” in the same way that cam-chem does.

enumerator ConstraintInitializationFailed

Newton iteration to initialize algebraic constraint variables failed to converge.

Functions

template<class FalloffParams> inline MICM_CONSTEXPR double FalloffKernel (const FalloffParams &p, double temperature, double air_density, double numerator_scale)

Shared falloff kernel for Troe and TernaryChemicalActivation. result = k0 * numerator_scale / (1 + ratio) * Fc^(N/(N + log10(ratio)^2)) Troe passes air_density as numerator_scale; Ternary passes 1.0.

inline MICM_CONSTEXPR double CalculateArrhenius (const ArrheniusRateConstantParameters &p, double temperature, double pressure)

Calculate Arrhenius rate constant. k = A * exp(C/T) * (T/D)^B * (1 + E*P)

inline MICM_CONSTEXPR double CalculateTroe (const TroeRateConstantParameters &p, double temperature, double air_density)

Calculate Troe rate constant.

inline MICM_CONSTEXPR double CalculateTernaryChemicalActivation (const TernaryChemicalActivationRateConstantParameters &p, double temperature, double air_density)

Calculate Ternary Chemical Activation rate constant.

inline MICM_CONSTEXPR double CalculateTunneling (const TunnelingRateConstantParameters &p, double temperature)

Calculate Tunneling rate constant. k = A * exp(-B/T + C/T^3)

inline MICM_CONSTEXPR double CalculateBranched (const BranchedRateConstantParameters &p, double temperature, double air_density)

Calculate Branched rate constant. Requires p.k0_ and p.z_ to be precomputed by ReactionRateConstantStore::BuildFrom.

inline MICM_CONSTEXPR double CalculateTaylorSeries (const TaylorSeriesRateConstantParameters &p, double temperature, double pressure)

Calculate Taylor Series rate constant. k = (sum_{j=0}^{n-1} c_j * T^j) * A * exp(C/T) * (T/D)^B * (1 + E*P)

inline MICM_CONSTEXPR double CalculateReversible (const ReversibleRateConstantParameters &p, double temperature)

Calculate Reversible rate constant. k = A * exp(C/T) * k_r.

inline MICM_CONSTEXPR double CalculateUserDefined (const UserDefinedRateConstantData &p, double custom_param_value)

Calculate user-defined rate constant. k = custom_param_value * scaling_factor.

inline MICM_CONSTEXPR double CalculateSurfaceOne (const SurfaceRateConstantData &p, double temperature, double radius, double num_conc)

Calculate one surface rate constant given pre-fetched aerosol parameters.

Parameters:

• radius – Aerosol effective radius [m]

• num_conc – Particle number concentration [# m-3]

template<template<class> class MatrixPolicy>
std::vector<std::size_t> DiagonalMarkowitzReorder(const MatrixPolicy<int> &matrix)

Reorders a set of state variables using Diagonal Markowitz algorithm.

Parameters:

matrix – Original matrix non-zero elements

Returns:

Reordered mapping vector (reordered[i] = original[map[i]])

template<class MatrixPolicy>
inline std::vector<std::size_t> DiagonalMarkowitzReorder(const MatrixPolicy &matrix)

inline std::string SolverStateToString(const SolverState &state)

template<class SparseMatrixPolicy>
SparseMatrixPolicy BuildJacobian(const std::set<std::pair<std::size_t, std::size_t>> &nonzero_jacobian_elements, std::size_t number_of_grid_cells, std::size_t state_size, bool indexing_only)

template<class DenseMatrixPolicy, class ForcingFunc>
DenseMatrixPolicy FiniteDifferenceJacobian(ForcingFunc forcing_func, const DenseMatrixPolicy &base_variables, std::size_t num_species, double perturbation = 1.0e-8)

Compute a dense finite-difference Jacobian approximation using central differences.

The forcing callable should have the signature: void(const DenseMatrixPolicy& variables, DenseMatrixPolicy& forcing) where variables has shape [num_blocks x num_species] and forcing has the same shape. Callers should bind any additional arguments (rate constants, state parameters, etc.) into the callable via a lambda capture.

Returns a DenseMatrixPolicy of shape [num_blocks x (num_species * num_species)] where element [block][row * num_species + col] = df_row/dx_col.

When a perturbation would push a variable below zero, one-sided differences are used instead.

template<class DenseMatrixPolicy, class SparseMatrixPolicy>
JacobianComparisonResult CompareJacobianToFiniteDifference(const SparseMatrixPolicy &analytical_jacobian, const DenseMatrixPolicy &fd_jacobian, std::size_t num_species, double atol = 1.0e-7, double rtol = 1.0e-7)

Compare an analytical sparse Jacobian (which stores -df/dx per MICM convention) against a finite-difference dense Jacobian (which stores +df/dx).

Uses a combined tolerance: |analytical - fd| < atol + rtol * max(|analytical|, |fd|)

Only non-zero elements in the sparse Jacobian are compared.

template<class DenseMatrixPolicy, class SparseMatrixPolicy>
JacobianComparisonResult CheckJacobianSparsityCompleteness(const SparseMatrixPolicy &analytical_jacobian, const DenseMatrixPolicy &fd_jacobian, std::size_t num_species, double threshold = 1.0e-6)

Check that no significant Jacobian entry exists outside the declared sparsity pattern.

This catches missing NonZeroJacobianElements declarations. Any FD entry outside the sparsity pattern that exceeds the threshold indicates an undeclared dependency.

inline std::string GenerateRandomString()

template<class MatrixPolicy>
void CheckCopyToDevice(MatrixPolicy &matrix)

template<class MatrixPolicy>
void CheckCopyToHost(MatrixPolicy &matrix)

Variables

template<typename T>
std::size_t GROUP_VECTOR_SIZE_V = GroupVectorSize<T>::value

Helper variable template.

template<class RatesPolicy, class LinearSolverPolicy, class ConstraintSetPolicy>
class AbstractBackwardEuler

#include <micm/solver/backward_euler.hpp>

An implementation of the fully implicit backward euler method.

Public Types

using ParametersType = BackwardEulerSolverParameters

Solver parameters typename.

Public Functions

inline AbstractBackwardEuler(LinearSolverPolicy &&linear_solver, RatesPolicy &&rates, ConstraintSetPolicy &&constraints)

Default constructor.

Parameters:

• linear_solver – Linear solver

• rates – Rates calculator

• constraints – Algebraic constraints (not used by BackwardEuler, for API compatibility)

inline SolverResult Solve(double time_step, auto &state, const BackwardEulerSolverParameters &parameters) const

Advances the given step over the specified time step.

Parameters:

• time_step – Time [s] to advance the state by

• state – The state to advance

Returns:

result of the solver (success or failure, and statistics)

Public Static Functions

template<class DenseMatrixPolicy>
static inline bool IsConverged(const BackwardEulerSolverParameters &parameters, const DenseMatrixPolicy &residual, const DenseMatrixPolicy &Yn1, const std::vector<double> &absolute_tolerance, double relative_tolerance)

Determines whether the residual is small enough to stop the internal solver iteration.

Parameters:

• residual – The residual to check

• state – The current state being solved for

Returns:

true if the residual is small enough to stop the iteration

template<class RatesPolicy, class LinearSolverPolicy, class ConstraintSetPolicy, class Derived>
class AbstractRosenbrockSolver

#include <micm/solver/rosenbrock.hpp>

An implementation of the Rosenbrock ODE solver.

This implements the Curiously Recurring Template Pattern to allow the AlphaMinusJacobian and NormalizedError functions to be implemented in extending classes and called from the base class Solve() function. https://en.cppreference.com/w/cpp/language/crtp

Template Parameters:

• RatesPolicy – Calculator of forcing and Jacobian terms

• LinearSolverPolicy – Linear solver

• ConstraintSetPolicy – Constraint set for algebraic constraints

• Derived – Implementation of the Rosenbock solver

Public Types

using ParametersType = RosenbrockSolverParameters

Solver parameters typename.

Public Functions

inline AbstractRosenbrockSolver(LinearSolverPolicy &&linear_solver, RatesPolicy &&rates, ConstraintSetPolicy &&constraints)

Default constructor.

Parameters:

• linear_solver – Linear solver

• rates – Rates calculator

• constraints – Algebraic constraints Note: This constructor is not intended to be used directly. Instead, use the SolverBuilder to create a solver

inline SolverResult Solve(double time_step, auto &state, const RosenbrockSolverParameters &parameters) const noexcept

Advances the given step over the specified time step.

Parameters:

time_step – Time [s] to advance the state by

Returns:

A struct containing results and a status code

inline SolverState InitializeConstraints(auto &state, const RosenbrockSolverParameters &parameters, SolverStats &stats) const noexcept

Newton-iterate algebraic variables to satisfy G(y) = 0 before time integration.

Parameters:

• state – The solver state

• parameters – Solver parameters (provides max iterations and tolerance)

• stats – Solver stats to update with iteration counts

Returns:

SolverState::Converged on success, or an error state on failure

template<class SparseMatrixPolicy>
inline void AlphaMinusJacobian(auto &state, const double &alpha) const

compute [alpha * I - dforce_dy]

Parameters:

• jacobian – Jacobian matrix (dforce_dy)

• alpha –

inline void LinearFactor(const double alpha, SolverStats &stats, auto &state) const

Perform the LU decomposition of the matrix.

Parameters:

• alpha – The alpha value

• number_densities – The number densities

• stats – The solver stats

• state – The state

template<class DenseMatrixPolicy>
inline double NormalizedError(const DenseMatrixPolicy &y, const DenseMatrixPolicy &y_new, const DenseMatrixPolicy &errors, auto &state) const

Computes the scaled norm of the vector errors.

Parameters:

• y – the original vector

• y_new – the new vector

• errors – The computed errors

Returns:

struct ArrheniusRateConstantParameters

Public Members

double A_ = {1}

Pre-exponential factor [(mol m−3)^(−(𝑛−1)) s−1].

double B_ = {0}

Unitless exponential factor.

double C_ = {0}

Activation threshold, expected to be the negative activation energy divided by the boltzman constant [-E_a / k_b), K].

double D_ = {300}

A factor that determines temperature dependence [K].

double E_ = {0}

A factor that determines pressure dependence [Pa-1].

struct BackwardEulerSolverParameters

#include <micm/solver/backward_euler_solver_parameters.hpp>

Backward Euler solver parameters.

template<class DenseMatrixPolicy>
class BackwardEulerTemporaryVariables : public micm::TemporaryVariables

Public Functions

inline virtual std::unique_ptr<TemporaryVariables> Clone() const override

Clone this object, preserving the derived type.

struct BlockVariableTag

#include <micm/util/view_category.hpp>

Tag for block variables (vector-like data holders)

struct BranchedRateConstantParameters

Public Members

Branch branch_

reaction branch

double X_

pre-exponential factor

double Y_

exponential factor

double a0_

branching factor

int n_

number of heavy atoms in the RO2 reacting species (excluding the peroxy moiety)

double k0_ = {0.0}

Precomputed low-pressure rate factor: 2e-22 * N_A * 1e-6 * exp(n_) Set by ReactionRateConstantStore::BuildFrom; do not set manually.

double z_ = {0.0}

Precomputed branching ratio factor: A(293, [M]_ref) * (1 - a0_) / a0_ Set by ReactionRateConstantStore::BuildFrom; do not set manually.

class ChemicalReaction

#include <micm/process/chemical_reaction.hpp>

Represents a chemical reaction with reactants, products, rate constant and phase.

class ChemicalReactionBuilder

Public Functions

inline ChemicalReactionBuilder &SetReactants(const std::vector<Species> &reactants)

Sets the list of reactant species involved in the chemical reaction.

Parameters:

reactants – A list of Species objects representing the reactants

Returns:

Reference to the builder

inline ChemicalReactionBuilder &SetProducts(const std::vector<StoichSpecies> &products)

Sets the list of product species and their yields for the chemical reaction.

Parameters:

products – A list of StoichSpecies objects representing the products

Returns:

Reference to the builder

template<class T>
inline ChemicalReactionBuilder &SetRateConstant(T &&rate_constant)

Sets the rate constant from any supported parameter struct. Accepts any type that is a member of RateConstantVariant.

Parameters:

rate_constant – Parameter struct (e.g. ArrheniusRateConstantParameters)

Returns:

Reference to the builder

inline ChemicalReactionBuilder &SetPhase(const Phase &phase)

Sets the phase in which the chemical reaction occurs (e.g., gas, aqueous)

Parameters:

phase – Phase object representing the reaction phase

Returns:

Reference to the builder

inline Process Build()

Transfers ownership of all internally stored data into a ChemicalReaction, then wraps it into a Process using std::variant.

Throws:

MicmException – if the rate constant has not been set

Returns:

A Process containing the constructed ChemicalReaction

struct Conditions

#include <micm/system/conditions.hpp>

Environemental conditions.

class Constraint

#include <micm/constraint/constraint.hpp>

This class uses std::variant to hold different constraint types. Each constraint provides:

• A residual function G(y) that should equal zero when the constraint is satisfied

• Jacobian entries dG/dy for each species the constraint depends on

Public Functions

inline std::string GetName() const

Get the constraint name.

Returns:

Constraint name

inline std::vector<std::string> GetParameterNames() const

Get the custom parameter names.

Returns:

A set of parameter names

inline const std::string &AlgebraicSpecies() const

Returns the species whose state row should be replaced by this algebraic constraint.

Returns:

Algebraic species name

inline const std::vector<std::string> &SpeciesDependencies() const

Get species dependencies.

Returns:

Vector of species names this constraint depends on

inline std::size_t NumberOfDependencies() const

Get the number of species this constraint depends on.

Returns:

Number of dependent species

template<typename DenseMatrixPolicy>
inline auto ConstraintParameterFunction(const ConstraintInfo &info) const

Get a function object to update constraint-specific parameters Returns a function that updates parameters (e.g., temperature-dependent K_eq) based on current conditions Called during solver build to create the update function, which is then invoked by UpdateStateParameters before each solve to recompute parameters from the latest temperature and other conditions.

Parameters:

info – Constraint information including state parameter indices

Returns:

Function object that takes (conditions, state_param) and updates constraint parameters

template<typename DenseMatrixPolicy>
inline auto ResidualFunction(const ConstraintInfo &info, const auto &state_variable_indices, const auto &state_parameter_indices) const

Get a function object to compute the constraint residual This returns a reusable function that can be invoked multiple times.

Parameters:

• info – Constraint information including species indices and row index

• state_variable_indices – Map from species names to state variable indices

• state_parameter_indices – Map from parameter names to state parameter indices

Returns:

Function object that takes (state_variables, state_parameters, forcing) and computes the residual

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
inline auto JacobianFunction(const ConstraintInfo &info, const auto &state_variable_indices, const auto &state_parameter_indices, auto jacobian_flat_ids, SparseMatrixPolicy &jacobian) const

Get a function object to compute the constraint Jacobian This returns a reusable function that can be invoked multiple times.

Parameters:

• info – Constraint information including species indices and Jacobian flat IDs

• state_variable_indices – Map from species names to state variable indices

• state_parameter_indices – Map from parameter names to state parameter indices

• jacobian_flat_ids – Iterator to the jacobian flat IDs for this constraint

• jacobian – Sparse matrix to store Jacobian values

Returns:

Function object that takes (state_variables, state_parameters, jacobian) and computes partials

struct ConstraintInfo

#include <micm/constraint/constraint_info.hpp>

Information for each constraint (built during ConstraintSet construction)

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
class ConstraintSet

#include <micm/constraint/constraint_set.hpp>

Manages a collection of algebraic constraints for DAE solvers ConstraintSet handles the computation of constraint residuals (forcing terms) and Jacobian contributions for a set of constraints. It follows the same pattern as ProcessSet for integration with the Rosenbrock solver.

Public Functions

ConstraintSet() = default

Default constructor.

inline ConstraintSet(std::vector<Constraint> &&constraints, const std::unordered_map<std::string, std::size_t> &variable_map)

Construct a ConstraintSet from constraints and variable mapping Constraints replace selected species rows in the state/Jacobian (DAE formulation)

Parameters:

• constraints – Vector of constraints

• variable_map – Map from species names to state variable indices

ConstraintSet(ConstraintSet &&other) noexcept = default

Move constructor - default implementation.

ConstraintSet &operator=(ConstraintSet &&other) noexcept = default

Move assignment operator.

ConstraintSet(const ConstraintSet&) = default

Copy constructor.

ConstraintSet &operator=(const ConstraintSet&) = default

Copy assignment.

inline std::size_t Size() const

Get the number of constraints (built-in + external model)

inline const std::set<std::size_t> &AlgebraicVariableIds() const

Returns species ids whose rows are algebraic when constraints replace state rows.

Returns:

Set of variable ids for algebraic rows

inline void SetUniqueParameterNames()

Deduplicates parameter names across all constraints in the set Ensures all constraint parameters have globally unique names by appending numeric suffixes (_1, _2, etc.) to duplicates. This should be called immediately after construction so that parameter names are finalized before the solver builder creates the parameter map. This logic is not part of the constructor because it mutates the constraint parameters, which is considered beyond the scope of construction.

inline std::unordered_set<std::string> GetParameterNames() const

Returns all unique parameter names from all constraints in the set.

Returns:

Set of parameter names

inline void AddForcingTerms(const DenseMatrixPolicy &state_variables, const DenseMatrixPolicy &state_parameters, DenseMatrixPolicy &forcing) const

Add constraint residuals to forcing vector (constraint rows) For each constraint G_i, writes or adds G_i(x) to forcing[constraint_row].

Parameters:

• state_variables – Current species concentrations (grid cell, species)

• state_parameters – Current state parameters (grid cell, parameter) - e.g., temperature-dependent K_eq values

• forcing – Forcing terms (grid cell, state variable) - constraint rows will be modified

inline void SubtractJacobianTerms(const DenseMatrixPolicy &state_variables, const DenseMatrixPolicy &state_parameters, SparseMatrixPolicy &jacobian) const

Subtract constraint Jacobian terms from Jacobian matrix For each constraint G_i, subtracts dG_i/dx_j from jacobian[constraint_row, j] (Subtraction matches the convention used by ProcessSet)

Parameters:

• state_variables – Current species concentrations (grid cell, species)

• state_parameters – Current state parameters (grid cell, parameter) - e.g., temperature-dependent K_eq values

• jacobian – Sparse Jacobian matrix (grid cell, row, column)

inline void SetAlgebraicErrors(DenseMatrixPolicy &Yerror, const DenseMatrixPolicy &Y, const DenseMatrixPolicy &Ynew) const

Set algebraic variable error estimates using step changes For each algebraic variable a: Yerror[a] = Ynew[a] - Y[a].

Parameters:

• Yerror – Error vector — algebraic entries are overwritten with step changes

• Y – State at beginning of step

• Ynew – Proposed state at end of step (after constraint enforcement)

inline std::set<std::pair<std::size_t, std::size_t>> NonZeroJacobianElements() const

Returns positions of all non-zero Jacobian elements for constraint rows.

Returns:

Set of (row, column) index pairs

template<typename OrderingPolicy>
inline void SetJacobianFlatIds(const SparseMatrix<double, OrderingPolicy> &matrix)

Computes and stores flat indices for Jacobian elements.

Parameters:

matrix – The sparse Jacobian matrix

inline void SetConstraintFunctions(const auto &state_variable_indices, const auto &state_parameter_indices, SparseMatrixPolicy &jacobian)

Pre-compiles constraint residual and Jacobian functions for efficient evaluation Creates reusable function objects from each constraint’s ResidualFunction and JacobianFunction. Must be called after SetJacobianFlatIds and before solver execution.

Parameters:

• state_variable_indices – Map from species names to state variable indices

• jacobian – The sparse Jacobian matrix (used for function template instantiation)

inline auto GetUpdateStateParamFunctions()

Returns pre-compiled constraint parameter update functions These functions compute temperature-dependent parameters (e.g., K_eq) for each constraint Called by solver builder to retrieve functions for UpdateStateParameters pipeline.

Returns:

Vector of function objects that take (conditions, state_param) and update constraint parameters

inline void SetExternalConstraintModels(std::vector<ExternalModelConstraintSet<DenseMatrixPolicy, SparseMatrixPolicy>> &&models)

Set external model constraint wrappers.

Parameters:

models – Vector of type-erased external model constraint wrappers

inline void ResolveExternalConstraints(const std::unordered_map<std::string, std::size_t> &variable_map)

Resolve external model constraints at runtime.

Calls each external model’s algebraic_variable_names_func_() to determine which (if any) algebraic variables it contributes. Models returning empty sets are skipped. Populates external_constraint_count_ and adds to algebraic_variable_ids_.

Parameters:

variable_map – Map from species names to state variable indices

inline std::set<std::pair<std::size_t, std::size_t>> ExternalNonZeroJacobianElements(const std::unordered_map<std::string, std::size_t> &variable_map) const

Returns non-zero Jacobian elements contributed by external model constraints.

Parameters:

variable_map – Map from species names to state variable indices

Returns:

Set of (row, column) index pairs

inline std::vector<std::string> ExternalConstraintParameterNames() const

Returns all unique state parameter names from external constraint models.

Returns:

Vector of parameter names (duplicates across models are detected and rejected)

inline auto GetExternalUpdateStateParamFunctions() const

Returns pre-compiled external constraint parameter update functions.

inline std::vector<std::string> ExternalInitializeConstraintParameterNames() const

Returns all unique state parameter names that need initialization from state variables.

Returns:

Vector of parameter names for state-diagnosed constraint parameters

inline auto GetExternalInitializeConstraintParamFunctions() const

Returns pre-compiled external constraint parameter initialization functions.

inline void InitializeConstraintParameters(const DenseMatrixPolicy &state_variables, DenseMatrixPolicy &state_parameters) const

Initializes constraint parameters from current state variables Called at the beginning of each Solve() to diagnose state-dependent constraint constants.

Parameters:

• state_variables – Current species concentrations

• state_parameters – State parameters to be updated with diagnosed values

inline void SetExternalModelConstraintFunctions(const std::unordered_map<std::string, std::size_t> &state_parameter_indices, const std::unordered_map<std::string, std::size_t> &state_variable_indices, const SparseMatrixPolicy &jacobian)

Pre-compiles external constraint residual, Jacobian, and parameter update functions Must be called after ResolveExternalConstraints and after Jacobian is built.

Parameters:

• state_parameter_indices – Map from parameter names to state parameter indices

• state_variable_indices – Map from species names to state variable indices

• jacobian – The sparse Jacobian matrix

template<class T, std::size_t L = MICM_DEFAULT_VECTOR_SIZE>
class CudaDenseMatrix : public micm::VectorMatrix<T, MICM_DEFAULT_VECTOR_SIZE>

Public Functions

inline void Axpy(const double alpha, const CudaDenseMatrix<T, L> &x)

For each element in the VectorMatrix x and y, perform y = alpha * x + y, where alpha is a scalar constant.

Parameters:

• alpha – The scaling scalar to apply to the VectorMatrix x

• x – The input VectorMatrix

Returns:

0 if successful, otherwise an error code

inline void Max(const T x)

For each element of the VectorMatrix, perform y = max(y, x), where x is a scalar constant.

Parameters:

x – The scalar constant to compare against

inline void Min(const T x)

For each element of the VectorMatrix, perform y = min(y, x), where x is a scalar constant.

Parameters:

x – The scalar constant to compare against

inline void Fill(T val)

Set every matrix element to a given value on the GPU.

Parameters:

val – Value to set each element to

template<class SparseMatrixPolicy, class LuDecompositionPolicy = CudaLuDecompositionMozartInPlace>
class CudaLinearSolverInPlace : public micm::LinearSolverInPlace<SparseMatrixPolicy, CudaLuDecompositionMozartInPlace>

Public Functions

inline CudaLinearSolverInPlace()

This is the default constructor, taking no arguments;.

inline CudaLinearSolverInPlace(const SparseMatrixPolicy &matrix, typename SparseMatrixPolicy::value_type initial_value)

This constructor takes two arguments: a sparse matrix and its values The base class here takes three arguments: the third argument is a lamda function that creates an instance of LuDecompositionPolicy; in this case, we will use the CudaLuDecompositionInPlace specified at line 13; See line 17 of “linear_solver_in_place.inl” for more details about how this lamda function works;

inline ~CudaLinearSolverInPlace()

This is the destructor that will free the device memory of the constant data from the class “CudaLinearSolverInPlace”

Public Members

LinearSolverInPlaceParam devstruct_

This is an instance of struct “LinearSolverInPlaceParam” that holds the constant data of “CudaLinearSolverInPlace” class on the device

class CudaLuDecompositionMozartInPlace : public micm::LuDecompositionMozartInPlace

#include <micm/cuda/solver/cuda_lu_decomposition_mozart_in_place.hpp>

This CudaLuDecompositionMozartInPlace class inherits everything from the base class “LuDecompositionMozartInPlace”.

Public Functions

inline CudaLuDecompositionMozartInPlace()

This is the default constructor, taking no arguments;.

template<class SparseMatrixPolicy>
inline CudaLuDecompositionMozartInPlace(const SparseMatrixPolicy &matrix)

This is the overloaded constructor that takes one argument called “matrix”; We need to specify the type (e.g., double, int, etc) and ordering (e.g., vector-stored, non-vector-stored, etc) of the “matrix”;

inline ~CudaLuDecompositionMozartInPlace()

This is destructor that will free the device memory of the constant data from the class “CudaLuDecompositionMozartInPlace”

template<class SparseMatrixPolicy>
void Decompose(SparseMatrixPolicy &ALU) const

This is the function to perform an LU decomposition on a given A matrix on the GPU.

Parameters:

ALU – Sparse matrix to decompose (will be overwritten with L and U matrices)

Public Members

LuDecomposeMozartInPlaceParam devstruct_

This is an instance of struct “LuDecomposeMozartInPlaceParam” that holds the constant data of “CudaLuDecompositionMozartInPlace” class on the device

Public Static Functions

template<class SparseMatrixPolicy>
static inline CudaLuDecompositionMozartInPlace Create(const SparseMatrixPolicy &matrix)

Create an LU decomposition algorithm for a given sparse matrix policy.

Parameters:

matrix – Sparse matrix

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
class CudaProcessSet : public micm::ProcessSet<DenseMatrixPolicy, SparseMatrixPolicy>

#include <micm/cuda/process/cuda_process_set.hpp>

A GPU-based implementation of ProcessSet.

Template Parameters:

• DenseMatrixPolicy – Policy for dense matrices (must satisfy CudaMatrix concept)

• SparseMatrixPolicy – Policy for sparse matrices (must satisfy CudaMatrix concept)

Public Functions

inline CudaProcessSet(const std::vector<Process> &processes, const std::unordered_map<std::string, std::size_t> &variable_map)

Create a process set calculator for a given set of processes.

Parameters:

• processes – Processes to create calculator for

• variable_map – A mapping of species names to concentration index

inline CudaProcessSet(const std::vector<Process> &processes, const std::unordered_map<std::string, std::size_t> &variable_map, const std::vector<ExternalModelProcessSet<DenseMatrixPolicy, SparseMatrixPolicy>> &external_models)

Create a process set calculator for a given set of processes with external models.

Parameters:

• processes – Processes to create calculator for

• variable_map – A mapping of species names to concentration index

• external_models – External models to include

inline void BuildCudaStore(const ReactionRateConstantStore &cpu_store)

Upload all analytic parameter arrays from cpu_store to device memory. Called once by Solver after ReactionRateConstantStore is built.

template<class StatePolicy>
inline void GpuCalculateRateConstants(const ReactionRateConstantStore &cpu_store, StatePolicy &state)

GPU-accelerated rate constant calculation.

1. Evaluate any lambda entries on CPU; upload rate_constants_ to device.

2. Upload conditions and custom_rate_parameters_ to device.

3. Evaluate parameterized multipliers on CPU; pack and upload to device.

4. Launch CalculateRateConstantsKernel to fill analytic slots and apply multipliers.

After this call, device rate_constants_ is fully populated for the current step.

inline void SetJacobianFlatIds(const SparseMatrixPolicy &matrix)

Set the indexes for the elements of Jacobian matrix before we could copy it to the device;.

this will override the “SetJacobianFlatIds” function from the “ProcessSet” class

Parameters:

matrix –

inline void SetAlgebraicVariableIds(const std::set<std::size_t> &variable_ids)

Marks species rows that should be treated as algebraic (constraints replace ODE rows). Updates algebraic variable IDs after ProcessSetParam construction. If algebraic variable IDs are not set post-construction, then this function may not be necessary.

Parameters:

variable_ids – Set of variable ids whose forcing/Jacobian rows should not receive kinetic contributions

Public Members

ProcessSetParam devstruct_

This is an instance of struct “ProcessSetParam” that holds the constant data of “ProcessSet” class on the device

CudaReactionRateStore cuda_rate_store_

GPU-resident analytic rate constant parameter store (built once per solver build)

class CudaReactionRateStore

#include <micm/cuda/process/cuda_reaction_rate_store.hpp>

GPU-resident mirror of ReactionRateConstantStore analytic data.

Constructed once per solver build; never modified during a run. The device conditions buffer grows on demand (amortised allocation).

Public Functions

inline void BuildFrom(const ReactionRateConstantStore &cpu_store)

Upload all analytic parameter arrays from cpu_store to device memory.

   Called once after the ReactionRateConstantStore is built in Solver's constructor.
   Any previous device allocations are freed before re-uploading.

inline const double *UploadMultiplierValues(const ReactionRateConstantStore &cpu_store, const std::vector<Conditions> &conditions, std::size_t L)

Evaluate parameterized multipliers on CPU, pack into interleaved layout, and upload. Layout: [group * n_mults * L + mult * L + lane].

Returns:

Device pointer to multiplier values, or nullptr if there are no multipliers.

inline const Conditions *UploadConditions(const std::vector<Conditions> &conditions)

Upload the current conditions array to device, growing the buffer if needed.

Returns:

Device pointer valid until the next call to UploadConditions.

template<class RatesPolicy, class LinearSolverPolicy, class ConstraintSetPolicy>
class CudaRosenbrockSolver : public micm::AbstractRosenbrockSolver<RatesPolicy, LinearSolverPolicy, ConstraintSetPolicy, CudaRosenbrockSolver<RatesPolicy, LinearSolverPolicy, ConstraintSetPolicy>>

Public Types

using ParametersType = CudaRosenbrockSolverParameters

Default constructor.

Solver parameters typename

Public Functions

inline CudaRosenbrockSolver()

Default constructor.

inline CudaRosenbrockSolver(LinearSolverPolicy &&linear_solver, RatesPolicy &&rates, ConstraintSetPolicy &&constraints)

Builds a CUDA Rosenbrock solver for the given system and solver parameters.

Parameters:

• linear_solver – Linear solver

• rates – Rates calculator

• constraints – Algebraic constraints

inline ~CudaRosenbrockSolver()

This is the destructor that will free the device memory of the constant data from the class “CudaRosenbrockSolver”

template<class SparseMatrixPolicy>
inline void AlphaMinusJacobian(auto &state, const double &alpha) const

Computes [alpha * I - jacobian] on the GPU.

Template Parameters:

SparseMatrixPolicy –

Parameters:

• jacobian – Jacobian matrix

• jacobian_diagonal_elements – Diagonal elements of the Jacobian matrix, not used

• alpha –

template<class DenseMatrixPolicy>
inline double NormalizedError(const DenseMatrixPolicy &y_old, const DenseMatrixPolicy &y_new, const DenseMatrixPolicy &errors, auto &state) const

Computes the scaled norm of the vector errors on the GPU; assume all the data are GPU resident already.

Template Parameters:

DenseMatrixPolicy –

Parameters:

• y_old – the original vector

• y_new – the new vector

• errors – The computed errors

Returns:

The scaled norm of the errors

struct CudaRosenbrockSolverParameters : public micm::RosenbrockSolverParameters

#include <micm/cuda/solver/cuda_solver_parameters.hpp>

Parameters for the CUDA Rosenbrock solver.

Public Functions

inline CudaRosenbrockSolverParameters(const RosenbrockSolverParameters &base)

Constructor from base class.

Parameters:

base –

template<class T, class OrderingPolicy>
class CudaSparseMatrix : public micm::SparseMatrix<T, OrderingPolicy>

Public Functions

inline void Fill(T val)

Set every matrix element to a given value on the GPU.

Parameters:

val – Value to set each element to

template<class DenseMatrixPolicy, class SparseMatrixPolicy, class LuDecompositionPolicy>
struct CudaState : public micm::State<DenseMatrixPolicy, SparseMatrixPolicy, LuDecompositionPolicy>

#include <micm/cuda/solver/cuda_state.hpp>

Construct a state variable for CUDA tests.

Public Functions

inline CudaState(const StateParameters &parameters, const std::size_t number_of_grid_cells)

Constructor which takes the state dimension information as input.

Parameters:

• parameters – State dimension information

• number_of_grid_cells – Number of grid cells

inline CudaState(CudaState &&other)

Move constructor.

inline CudaState &operator=(CudaState &&other)

Move assignment operator.

inline virtual void SetAbsoluteTolerances(const std::vector<double> &absoluteTolerance) override

Set the absolute tolerances per species.

Parameters:

absoluteTolerance – absolute tolerance

inline void SyncInputsToDevice()

Copy input variables to the device.

   Rate constants are NOT copied here; they are computed directly on the
   GPU by Solver::UpdateStateParameters (via CudaProcessSet::GpuCalculateRateConstants).

inline void SyncOutputsToHost()

Copy output variables to the host.

struct DenseMatrixColumnViewTag

#include <micm/util/view_category.hpp>

Tag for dense matrix column views (have ColumnIndex + GetMatrix)

class EquilibriumConstraint

#include <micm/constraint/types/equilibrium_constraint.hpp>

Constraint for chemical equilibrium with temperature-dependent K_eq using Van’t Hoff equation For a reversible reaction: aA + bB <-> cC + dD The equilibrium constraint is: G = K_eq(T) * [A]^a * [B]^b - [C]^c * [D]^d = 0 where K_eq(T) = K_HLC_ref * exp((delta_H / R) * (1/T - 1/T_ref))

Public Functions

EquilibriumConstraint() = default

Default constructor.

inline EquilibriumConstraint(const std::string &name, const Species &algebraic_species, std::vector<StoichSpecies> reactants, std::vector<StoichSpecies> products, VantHoffParam vant_hoff_param)

Construct an equilibrium constraint. Validates that equilibrium constraint > 0. Builds species_dependencies_ by concatenating reactants then products. Stores index mappings for efficient Jacobian computation. Stores a temperature-dependent equilibrium constant function.

Parameters:

• name – Constraint identifier

• algebraic_species – Species whose row is replaced by this algebraic constraint

• reactants – Vector of StoichSpecies (species, stoichiometry) for reactants

• products – Vector of StoichSpecies (species, stoichiometry) for products

• vant_hoff_param – Parameters for Van’t Hoff equation

inline const std::string &AlgebraicSpecies() const

Returns the species whose row should be replaced by this algebraic constraint.

Returns:

Species name of the explicitly set algebraic variable

template<typename DenseMatrixPolicy>
inline std::function<void(const std::vector<Conditions>&, DenseMatrixPolicy&)> ConstraintParameterFunction(const ConstraintInfo &info) const

Create function object to update temperature-dependent K_eq parameter Returns a function that computes K_eq(T) for each grid cell using Van’t Hoff equation Called during solver’s UpdateStateParameters phase before each solve.

Parameters:

info – Constraint information including state parameter indices

Returns:

Function object that takes (conditions, state_param) and writes K_eq(T) to state_param[K_eq_idx]

template<typename DenseMatrixPolicy>
inline std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, DenseMatrixPolicy&)> ResidualFunction(const ConstraintInfo &info, const auto &state_variable_indices, const auto &state_parameter_indices) const

Create function object to compute equilibrium constraint residual for all grid cells Computes G = K_eq(T) * prod([reactants]^stoich) - prod([products]^stoich) for the algebraic constraint Called during solver build (SetConstraintFunctions) to pre-compile residual computation.

Parameters:

• info – Constraint information including row index, species indices, and parameter indices

• state_variable_indices – Mapping of state variable names to column indices in state matrix

• state_parameter_indices – Mapping of state parameter names to column indices in state_param matrix

Returns:

Function object that takes (state, state_param, forcing) and writes residual G to forcing[constraint_row]

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
inline std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, SparseMatrixPolicy&)> JacobianFunction(const ConstraintInfo &info, const auto &state_variable_indices, const auto &state_parameter_indices, auto jacobian_flat_ids, SparseMatrixPolicy &jacobian) const

Create function object to compute Jacobian partial derivatives dG/d[species] for all grid cells For reactant R with stoichiometry n: dG/d[R] = K_eq(T) * n * [R]^(n-1) * prod([other_reactants]^stoich) For product P with stoichiometry m: dG/d[P] = -m * [P]^(m-1) * prod([other_products]^stoich) Called during solver build (SetConstraintFunctions) to pre-compile Jacobian computation.

Parameters:

• info – Constraint information including row index, species indices, and parameter indices

• state_variable_indices – Mapping of state variable names to column indices in state matrix

• state_parameter_indices – Mapping of state parameter names to column indices in state_param matrix

• jacobian_flat_ids – Iterator to this constraint’s flat Jacobian indices in sparse matrix storage

• jacobian – Sparse matrix reference (used for type information)

Returns:

Function object that takes (state, state_param, jacobian_values) and writes partials to sparse Jacobian

Public Members

std::string name_

Name of the constraint, used when generating state parameter name.

Species algebraic_species_

Algebraic species.

std::vector<std::string> species_dependencies_

Names of species this constraint depends on.

std::vector<StoichSpecies> reactants_

Reactant species and their stoichiometric coefficients.

std::vector<StoichSpecies> products_

Product species and their stoichiometric coefficients.

std::vector<std::string> parameters_

For equilibrium constraints, this contains a single parameter K_eq.

std::function<double(const Conditions&)> equilibrium_constant_function_

Temperature-dependent Henry’s Law Constant.

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
struct ExternalModelConstraintSet

#include <micm/external_model.hpp>

Wrapper for external model constraint information.

This struct encapsulates an external model’s constraint definitions (algebraic variables, residual functions, Jacobian functions) and provides a type-erased interface that MICM can use to incorporate the model’s constraints into the DAE solver. Instances are constructed when AddExternalModel() is called on a solver builder for a model satisfying HasConstraints.

Template Parameters:

• DenseMatrixPolicy – Policy for dense matrices (state variables, forcing)

• SparseMatrixPolicy – Policy for sparse matrices (Jacobian)

Public Functions

template<typename ModelType, typename = std::enable_if_t<!std::is_same_v<std::decay_t<ModelType>, ExternalModelConstraintSet>>>
inline ExternalModelConstraintSet(ModelType &&model)

Constructs a type-erased wrapper from an external model instance.

Template Parameters:

ModelType – Type of the external model (must satisfy HasConstraints)

Parameters:

model – External model instance

Public Members

std::function<std::set<std::string>()> algebraic_variable_names_func_

Type-erased function returning names of algebraic variables (may be empty at runtime)

std::function<std::set<std::string>()> species_dependencies_func_

Type-erased function returning species dependencies for constraints.

std::function<std::set<std::pair<std::size_t, std::size_t>>(const std::unordered_map<std::string, std::size_t>&)> non_zero_jacobian_elements_func_

Type-erased function returning non-zero constraint Jacobian element positions.

std::function<std::set<std::string>()> state_parameter_names_func_

Type-erased function returning constraint state parameter names.

std::function<std::function<void(const std::vector<micm::Conditions>&, DenseMatrixPolicy&)>(const std::unordered_map<std::string, std::size_t> &state_parameter_indices)> update_state_parameters_function_

Type-erased function factory for constraint state parameter updates Returns a function that updates constraint parameters based on environmental conditions.

std::function<std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, DenseMatrixPolicy&)>(const std::unordered_map<std::string, std::size_t> &state_parameter_indices, const std::unordered_map<std::string, std::size_t> &state_variable_indices)> get_residual_function_

Type-erased function factory for constraint residual computation Returns a function that computes G(y) for the constraint rows.

std::function<std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, SparseMatrixPolicy&)>(const std::unordered_map<std::string, std::size_t> &state_parameter_indices, const std::unordered_map<std::string, std::size_t> &state_variable_indices, const SparseMatrixPolicy &jacobian)> get_jacobian_function_

Type-erased function factory for constraint Jacobian computation Returns a function that computes ∂G/∂y for the constraint rows.

std::function<std::set<std::string>()> initialize_constraint_parameter_names_func_

Type-erased function returning parameter names that need initialization from state Returns an empty set if the model does not need state-diagnosed parameters.

std::function<std::function<void(const DenseMatrixPolicy&, DenseMatrixPolicy&)>(const std::unordered_map<std::string, std::size_t> &state_parameter_indices, const std::unordered_map<std::string, std::size_t> &state_variable_indices)> get_initialize_constraint_parameters_function_

Type-erased function factory for constraint parameter initialization from state Returns a function that diagnoses constraint parameters from current state variables at the beginning of each Solve() call.

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
struct ExternalModelProcessSet

#include <micm/external_model.hpp>

Wrapper for external model process information.

This struct encapsulates an external model’s process definitions (forcing functions, Jacobian functions, and state parameter updates) and provides a type-erased interface that MICM can use to incorporate the model’s processes into the ODE solver. Instances are constructed when AddExternalModel() is called on a solver builder.

The wrapped functions are used during the solve to:

• Update state parameters based on environmental conditions (temperature, pressure, etc.)

• Calculate forcing terms (tendencies) for the model’s processes

• Calculate Jacobian contributions for implicit time integration

Note

Users typically do not construct this directly; instead, pass external model instances to solver builder’s AddExternalModel() method.

Template Parameters:

• DenseMatrixPolicy – Policy for dense matrices (state variables, parameters, forcing)

• SparseMatrixPolicy – Policy for sparse matrices (Jacobian)

Public Functions

template<typename ModelType, typename = std::enable_if_t<!std::is_same_v<std::decay_t<ModelType>, ExternalModelProcessSet>>>
inline ExternalModelProcessSet(ModelType &&model)

Constructs a type-erased wrapper from an external model instance.

This constructor captures the model instance and wraps its process definition methods in std::function objects. The wrapped functions can then be called by MICM’s solver without knowledge of the original model type, enabling seamless integration of external processes.

Template Parameters:

ModelType – Type of the external model (must implement the external model interface)

Parameters:

model – External model instance (will be moved or copied into shared ownership)

Public Members

std::function<std::set<std::pair<std::size_t, std::size_t>>(const std::unordered_map<std::string, std::size_t>&)> non_zero_jacobian_elements_func_

Type-erased function returning non-zero Jacobian element positions.

std::function<std::set<std::string>()> species_used_func_

Type-erased function returning the set of species used by the model’s processes.

std::function<std::function<void(const std::vector<micm::Conditions>&, DenseMatrixPolicy&)>(const std::unordered_map<std::string, std::size_t> &state_parameter_indices)> update_state_parameters_function_

Type-erased function factory for state parameter updates Returns a function that updates state parameters based on environmental conditions.

std::function<std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, DenseMatrixPolicy&)>(const std::unordered_map<std::string, std::size_t> &state_parameter_indices, const std::unordered_map<std::string, std::size_t> &state_variable_indices)> get_forcing_function_

Type-erased function factory for forcing term calculation Returns a function that computes forcing terms (tendencies) for the model’s processes.

std::function<std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, SparseMatrixPolicy&)>(const std::unordered_map<std::string, std::size_t> &state_parameter_indices, const std::unordered_map<std::string, std::size_t> &state_variable_indices, const SparseMatrixPolicy &jacobian)> get_jacobian_function_

Type-erased function factory for Jacobian calculation Returns a function that computes Jacobian contributions for the model’s processes.

struct ExternalModelSystem

#include <micm/external_model.hpp>

Wrapper for external model state information.

This struct encapsulates an external model’s state definition (variables and parameters) and provides a type-erased interface that MICM can use to query the model’s state requirements. Instances are constructed automatically when an external model is passed to SolverBuilder::AddExternalModel().

Note

Users typically do not construct this directly; instead, pass external model instances to SolverBuilder::AddExternalModel() and they will be wrapped automatically.

Public Functions

ExternalModelSystem() = delete

Default constructor is deleted (must construct from an external model instance)

template<typename ModelType, typename = std::enable_if_t<!std::is_same_v<std::decay_t<ModelType>, ExternalModelSystem>>>
inline ExternalModelSystem(ModelType &&model)

Constructs a type-erased wrapper from an external model instance.

This constructor captures the model instance and wraps its state definition methods in std::function objects that can be called without knowledge of the original model type.

Template Parameters:

ModelType – Type of the external model (must implement the external model interface)

Parameters:

model – External model instance (will be moved or copied into shared ownership)

Public Members

std::function<std::tuple<std::size_t, std::size_t>()> state_size_func_

Type-erased function returning the state size (number of variables, number of parameters)

std::function<std::set<std::string>()> variable_names_func_

Type-erased function returning the set of state variable names.

std::function<std::set<std::string>()> parameter_names_func_

Type-erased function returning the set of state parameter names.

template<typename T>
struct GroupingStrategy

Determines the grouping strategy of a matrix type (no default - must be specialized)

template<typename T>
struct GroupingStrategy<Matrix<T>>

#include <micm/util/matrix.hpp>

Matrix always uses simple grouping (L==1)

template<typename T, std::size_t L>
struct GroupingStrategy<VectorMatrix<T, L>>

#include <micm/util/vector_matrix.hpp>

VectorMatrix uses simple grouping when L==1, tiered grouping when L>1.

template<typename T>
struct GroupVectorSize : public std::integral_constant<std::size_t, 1>

#include <micm/util/sparse_matrix.hpp>

Type trait to extract GroupVectorSize (L) from matrix types at compile-time Default: L=1 for types without GroupVectorSize.

template<typename T>
struct GroupVectorSize<T> : public std::integral_constant<std::size_t, 1>, public std::integral_constant<std::size_t, T::GroupVectorSize()>

#include <micm/util/sparse_matrix.hpp>

Specialization for types with static GroupVectorSize method.

template<typename T, typename = void>
struct HasCategory : public std::false_type

#include <micm/util/view_category.hpp>

Helper to check if a type has a nested ‘category’ type.

template<typename T>
struct HasCategory<T, std::void_t<typename T::category>> : public std::false_type, public std::true_type

struct JacobianComparisonResult

#include <micm/util/jacobian_verification.hpp>

Result of comparing an analytical Jacobian against a finite-difference approximation.

struct LambdaRateConstantParameters

Public Members

std::string label_

Label for the reaction used to identify user-defined parameters.

std::function<double(const Conditions&)> lambda_function_

Lambda function for calculating the rate constant.

class LinearConstraint

#include <micm/constraint/types/linear_constraint.hpp>

Constraint for linear relationships: sum(coeff[i] * [species[i]]) = constant For example: A + B + C = 1.0 represents a conservation law The linear constraint is: G = c1*[A] + c2*[B] + c3*[C] - constant = 0.

Public Functions

LinearConstraint() = default

Default constructor.

inline LinearConstraint(const std::string &name, const Species &algebraic_species, const std::vector<StoichSpecies> &terms, double constant)

Construct a linear constraint Validates that terms are non-empty Builds species_dependencies_ from terms.

Parameters:

• name – Constraint identifier

• algebraic_species – Species whose row is replaced by this algebraic constraint

• terms – Vector of StoichSpecies (species, coefficient) in the linear sum

• constant – The value that sum(coeff[i] * [species[i]]) should equal

inline const std::string &AlgebraicSpecies() const

Returns the species whose row should be replaced by this algebraic constraint.

Returns:

Species name of the explicitly set algebraic variable

template<typename DenseMatrixPolicy>
inline std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, DenseMatrixPolicy&)> ResidualFunction(const ConstraintInfo &info, const auto &state_variable_indices, const auto &state_parameter_indices) const

Create a function to compute the linear constraint residual Returns a reusable function object that evaluates: G = sum(coeff[i] * [species[i]]) - constant.

Parameters:

• info – Constraint information including species indices and row index

• state_variable_indices – Mapping of state variable names to indices

• state_parameter_indices – Mapping of parameter names to indices (unused for LinearConstraint)

Returns:

Function object that takes (state_variables, state_parameters, forcing) and computes residual

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
inline std::function<void(const DenseMatrixPolicy&, const DenseMatrixPolicy&, SparseMatrixPolicy&)> JacobianFunction(const ConstraintInfo &info, const auto &state_variable_indices, const auto &state_parameter_indices, auto jacobian_flat_ids, SparseMatrixPolicy &jacobian) const

Create a function to compute Jacobian entries dG/d[species] For a linear constraint, the Jacobian is simply the coefficients: dG/d[species[i]] = coeff[i].

Parameters:

• info – Constraint information including species indices

• state_variable_indices – Mapping of state variable names to indices

• state_parameter_indices – Mapping of parameter names to indices (unused for LinearConstraint)

• jacobian_flat_ids – Iterator to this constraint’s flat Jacobian indices in shared storage

• jacobian – Sparse matrix to store Jacobian values

Returns:

Function object that takes (state_variables, state_parameters, jacobian) and computes partials

Public Members

std::string name_

Name of the constraint.

Species algebraic_species_

Algebraic species.

std::vector<std::string> species_dependencies_

Names of species this constraint depends on.

std::vector<StoichSpecies> terms_

Species and their coefficients in the linear sum.

double constant_

The constant value the linear sum should equal.

std::vector<std::string> parameters_

Parameter set (unused for this class, always empty).

template<class SparseMatrixPolicy, class LuDecompositionPolicy = LuDecomposition, class LMatrixPolicy = SparseMatrixPolicy, class UMatrixPolicy = SparseMatrixPolicy>
class LinearSolver

#include <micm/solver/linear_solver.hpp>

A general-use block-diagonal sparse-matrix linear solver.

The sparsity pattern of each block in the block diagonal matrix is the same.

Public Functions

inline LinearSolver()

default constructor

inline LinearSolver(const SparseMatrixPolicy &matrix, typename SparseMatrixPolicy::value_type initial_value)

Constructs a linear solver for the sparsity structure of the given matrix.

Parameters:

• matrix – Sparse matrix

• initial_value – Initial value for matrix elements

inline LinearSolver(const SparseMatrixPolicy &matrix, typename SparseMatrixPolicy::value_type initial_value, const std::function<LuDecompositionPolicy(const SparseMatrixPolicy&)> create_lu_decomp)

Constructs a linear solver for the sparsity structure of the given matrix.

Parameters:

• matrix – Sparse matrix

• initial_value – Initial value for matrix elements

• create_lu_decomp – Function to create an LU Decomposition object that adheres to LuDecompositionPolicy

inline void Factor(const SparseMatrixPolicy &matrix, LMatrixPolicy &lower_matrix, UMatrixPolicy &upper_matrix) const

Decompose the matrix into upper and lower triangular matrices.

Parameters:

matrix – Matrix to decompose into lower and upper triangular matrices

template<class MatrixPolicy>
inline void Solve(MatrixPolicy &x, const LMatrixPolicy &lower_matrix, const UMatrixPolicy &upper_matrix) const

Solve for x in Ax = b. x should be a copy of b and after Solve finishes x will contain the result.

template<class SparseMatrixPolicy, class LuDecompositionPolicy = LuDecompositionInPlace>
class LinearSolverInPlace

#include <micm/solver/linear_solver_in_place.hpp>

A general-use block-diagonal sparse-matrix linear solver.

The sparsity pattern of each block in the block diagonal matrix is the same. The L and U matrices are decomposed in-place over the original A matrix.

Public Functions

inline LinearSolverInPlace()

default constructor

inline LinearSolverInPlace(const SparseMatrixPolicy &matrix, typename SparseMatrixPolicy::value_type initial_value)

Constructs a linear solver for the sparsity structure of the given matrix.

Parameters:

• matrix – Sparse matrix

• initial_value – Initial value for matrix elements

inline LinearSolverInPlace(const SparseMatrixPolicy &matrix, typename SparseMatrixPolicy::value_type initial_value, const std::function<LuDecompositionPolicy(const SparseMatrixPolicy&)> create_lu_decomp)

Constructs a linear solver for the sparsity structure of the given matrix.

Parameters:

• matrix – Sparse matrix

• initial_value – Initial value for matrix elements

• create_lu_decomp – Function to create an LU Decomposition object that adheres to LuDecompositionPolicy

inline void Factor(SparseMatrixPolicy &matrix) const

Decompose the matrix into upper and lower triangular matrices (matrix will be overwritten)

Parameters:

matrix – Matrix to decompose in-place into lower and upper triangular matrices

template<class MatrixPolicy>
inline void Solve(MatrixPolicy &x, const SparseMatrixPolicy &lu_matrix) const

Solve for x in Ax = b. x should be a copy of b and after Solve finishes x will contain the result.

Parameters:

• x – The solution vector

• LU – The LU decomposition of the matrix as a square sparse matrix

class LuDecompositionDoolittle

#include <micm/solver/lu_decomposition_doolittle.hpp>

LU decomposer for SparseMatrix following the Doolittle algorithm.

The LU decomposition uses the Doolittle algorithm following the naming used here: https://www.geeksforgeeks.org/doolittle-algorithm-lu-decomposition/

The sudo-code for the corresponding dense matrix algorithm for matrix A and lower (upper) triangular matrix L(U) would be:

for i = 0…n-1 // Outer loop over rows (columns) for upper (lower) triangular matrix for k = i…n-1 // Middle loop over columns for upper triangular matrix sum = 0 for j = 0…i-1 // Inner loop over columns (rows) for lower (upper) triangular matrix sum += L[i][j] * U[j][k] U[i][k] = A[i][k] - sum L[i][i] = 1 // Lower triangular matrix is 1 along the diagonal for k = i+1…n-1 // Middle loop over rows for lower triangular matrix sum = 0 for j = 0…i-1 // Inner loop over columns (rows) for lower (upper) triangular matrix sum += L[k][j] * U[j][i]; L[k][i] = (A[k][i] - sum) / U[i][i]

For the sparse matrix algorithm, the indices of non-zero terms are stored in several arrays during construction. These arrays are iterated through during calls to Decompose to do the actual decomposition. Our LU Decomposition only assigns the values of the jacobian to the LU matrices when the jacobian is nonzero. However, the sparsity pattern of the jacobian doesn’t necessarily match that of the LU matrices. There can be more nonzero elements in the LU matrices than in the jacobian. When this happens, we still need to assign the value of the jacobian matrix to the LU matrix. This value is implicitly zero when the sparsity pattern differs. The Fill values here do this implicit assignment More detail in this issue: NCAR/micm#625

Public Functions

inline LuDecompositionDoolittle()

default constructor

template<class SparseMatrixPolicy, class LMatrixPolicy = SparseMatrixPolicy, class UMatrixPolicy = SparseMatrixPolicy>
inline LuDecompositionDoolittle(const SparseMatrixPolicy &matrix)

Construct an LU decomposition algorithm for a given sparse matrix.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy>
inline void Decompose(const SparseMatrixPolicy &A, auto &L, auto &U) const

Perform an LU decomposition on a given A matrix.

Parameters:

• A – Sparse matrix to decompose

• L – The lower triangular matrix created by decomposition

• U – The upper triangular matrix created by decomposition

Public Static Functions

template<class SparseMatrixPolicy, class LMatrixPolicy = SparseMatrixPolicy, class UMatrixPolicy = SparseMatrixPolicy>
static inline LuDecompositionDoolittle Create(const SparseMatrixPolicy &matrix)

Create an LU decomposition algorithm for a given sparse matrix policy.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy, class LMatrixPolicy = SparseMatrixPolicy, class UMatrixPolicy = SparseMatrixPolicy>
static inline std::pair<LMatrixPolicy, UMatrixPolicy> GetLUMatrices(const SparseMatrixPolicy &A, typename SparseMatrixPolicy::value_type initial_value, bool indexing_only = false)

Create sparse L and U matrices for a given A matrix.

Parameters:

A – Sparse matrix that will be decomposed

Returns:

L and U Sparse matrices

class LuDecompositionDoolittleInPlace

#include <micm/solver/lu_decomposition_doolittle_in_place.hpp>

LU decomposer for SparseMatrix following the Doolittle algorithm.

The LU decomposition uses the Doolittle algorithm following the naming used here: https://www.geeksforgeeks.org/doolittle-algorithm-lu-decomposition/

The sudo-code for the corresponding dense matrix algorithm for matrix A (in-line) would be:

for i = 0…n-1 // Outer loop over rows (columns) for upper (lower) triangular matrix for k = i…n-1 // Middle loop over columns for upper triangular matrix for j = 0…i-1 // Inner loop over columns (rows) for lower (upper) triangular matrix A[i][k] -= A[i][j] * A[j][k] for k = i+1…n-1 // Middle loop over rows for lower triangular matrix for j = 0…i-1 // Inner loop over columns (rows) for lower (upper) triangular matrix A[k][i] -= A[k][j] * A[j][i]; A[k][i] /= A[i][i]

For the sparse matrix algorithm, the indices of non-zero terms are stored in several arrays during construction. These arrays are iterated through during calls to Decompose to do the actual decomposition. Our LU Decomposition only assigns the values of the jacobian to the LU matrices when the jacobian is nonzero. However, the sparsity pattern of the jacobian doesn’t necessarily match that of the LU matrices. There can be more nonzero elements in the LU matrices than in the jacobian. It is expected that the elements of the L and U matrices that are zero in the A matrix will be set to zero before the combined matrix is passed to the decomposition function.

Public Functions

inline LuDecompositionDoolittleInPlace()

default constructor

template<class SparseMatrixPolicy>
inline LuDecompositionDoolittleInPlace(const SparseMatrixPolicy &matrix)

Construct an LU decomposition algorithm for a given sparse matrix.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy>
inline void Decompose(SparseMatrixPolicy &ALU) const

Perform an LU decomposition on a given A matrix.

Parameters:

• A – Sparse matrix to decompose

• L – The lower triangular matrix created by decomposition

• U – The upper triangular matrix created by decomposition

Public Static Functions

template<class SparseMatrixPolicy>
static inline LuDecompositionDoolittleInPlace Create(const SparseMatrixPolicy &matrix)

Create an LU decomposition algorithm for a given sparse matrix policy.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy>
static inline SparseMatrixPolicy GetLUMatrix(const SparseMatrixPolicy &A, typename SparseMatrixPolicy::value_type initial_value, bool indexing_only = false)

Create sparse L and U matrices for a given A matrix.

Parameters:

A – Sparse matrix that will be decomposed

Returns:

L and U Sparse matrices

class LuDecompositionMozart

#include <micm/solver/lu_decomposition_mozart.hpp>

LU decomposer for SparseMatrix following the algorithm from the MOZART model.

This LU decomposition uses the algorithm from the MOZART chemistry preprocessor at: ESCOMP/CHEM_PREPROCESSOR

The MOZART function overwrote the A matrix with the L and U matrices. The pseudo-code in C++ for the corresponding dense matrix algorithm for matrix A (inline change) would be:

for i = 0…n-1 // Outer loop over columns of the sparse matrix A for j = i+1…n-1 // Multiply column below diagonal A[j][i] = A[j][i] / A[i][i] for k = i+1…n-1 // Modify sub-matrix for j = i+1…n-1 A[j][k] = A[j][k] – A[j][i] * A[i][k]

The pseudo-code in C++ for the corresponding dense matrix algorithm for matrix A and separate lower (upper) triangular matrix L(U) would be:

for i = 0…n-1 // Initialize U and L matrices to the A values for j = 0…i // Initialize U matrix including diagonal U[j][i] = A[j][i] L[i][i] = 1 // Lower triangular matrix is 1 along the diagonal for j = i+1…n-1 // Initialize L matrix excluding diagonal L[j][i] = A[j][i] for i = 0…n-1 for j = i+1…n-1 // Multiply column below diagonal L[j][i] = L[j][i] / U[i][i] for k = i+1…n-1 // Modify sub-matrix for j = i+1…k U[j][k] = U[j][k] - L[j][i] * U[i][k] for j = k+1…n-1 L[j][k] = L[j][k] - L[j][i] * U[i][k]

For the sparse matrix algorithm, the indices of non-zero terms are stored in several arrays during construction. These arrays are iterated through during calls to Decompose to do the actual decomposition. Our LU Decomposition only assigns the values of the jacobian to the LU matrices when the jacobian is nonzero. However, the sparsity pattern of the jacobian doesn’t necessarily match that of the LU matrices. There can be more nonzero elements in the LU matrices than in the jacobian. When this happens, we still need to assign the value of the jacobian matrix to the LU matrix. This value is implicitly zero when the sparsity pattern differs. The Fill values here do this implicit assignment More detail in this issue: NCAR/micm#625

Public Functions

inline LuDecompositionMozart()

default constructor

template<class SparseMatrixPolicy>
LuDecompositionMozart(const SparseMatrixPolicy &matrix)

Construct an LU decomposition algorithm for a given sparse matrix.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy, class LMatrixPolicy, class UMatrixPolicy>
inline LuDecompositionMozart(const SparseMatrixPolicy &matrix)

Construct an LU decomposition algorithm for a given sparse matrix.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy>
inline void Decompose(const SparseMatrixPolicy &A, auto &L, auto &U) const

Perform an LU decomposition on a given A matrix.

Parameters:

• A – Sparse matrix to decompose

• L – The lower triangular matrix created by decomposition

• U – The upper triangular matrix created by decomposition

Public Static Functions

template<class SparseMatrixPolicy, class LMatrixPolicy, class UMatrixPolicy>
static inline LuDecompositionMozart Create(const SparseMatrixPolicy &matrix)

Create an LU decomposition algorithm for a given sparse matrix policy.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy, class LMatrixPolicy, class UMatrixPolicy>
static inline std::pair<LMatrixPolicy, UMatrixPolicy> GetLUMatrices(const SparseMatrixPolicy &A, typename SparseMatrixPolicy::value_type initial_value, bool indexing_only = false)

Create sparse L and U matrices for a given A matrix.

Parameters:

A – Sparse matrix that will be decomposed

Returns:

L and U Sparse matrices

class LuDecompositionMozartInPlace

#include <micm/solver/lu_decomposition_mozart_in_place.hpp>

LU decomposer for SparseMatrix following the algorithm from the MOZART model.

This LU decomposition uses the algorithm from the MOZART chemistry preprocessor at: ESCOMP/CHEM_PREPROCESSOR

The MOZART function overwrote the A matrix with the L and U matrices. The pseudo-code in C++ for the corresponding dense matrix algorithm for matrix A (inline change) would be:

for i = 0…n-1 // Outer loop over columns of the sparse matrix A for j = i+1…n-1 // Multiply column below diagonal A[j][i] = A[j][i] / A[i][i] for k = i+1…n-1 // Modify sub-matrix for j = i+1…n-1 A[j][k] = A[j][k] – A[j][i] * A[i][k]

For the sparse matrix algorithm, the indices of non-zero terms are stored in several arrays during construction. These arrays are iterated through during calls to Decompose to do the actual decomposition.

The GetLUMatrices function creates a new sparse matrix that includes the superset of the non-zero elements in the L and U matrices. It is expected that the elements of the L and U matrices that are zero in the A matrix will be set to zero before the combined matrix is passed to the decomposition function.

Subclassed by micm::CudaLuDecompositionMozartInPlace

Public Functions

inline LuDecompositionMozartInPlace()

default constructor

template<class SparseMatrixPolicy>
inline LuDecompositionMozartInPlace(const SparseMatrixPolicy &matrix)

Construct an LU decomposition algorithm for a given sparse matrix.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy>
inline void Decompose(SparseMatrixPolicy &ALU) const

Perform an LU decomposition on a given A matrix. All elements of L and U that are zero in A should be set to zero before calling this function.

Parameters:

ALU – Sparse matrix to decompose (will be overwritten with L and U matrices)

Public Static Functions

template<class SparseMatrixPolicy>
static inline LuDecompositionMozartInPlace Create(const SparseMatrixPolicy &matrix)

Create an LU decomposition algorithm for a given sparse matrix policy.

Parameters:

matrix – Sparse matrix

template<class SparseMatrixPolicy>
static inline SparseMatrixPolicy GetLUMatrix(const SparseMatrixPolicy &A, typename SparseMatrixPolicy::value_type initial_value, bool indexing_only = false)

Create a combined sparse L and U matrix for a given A matrix.

Parameters:

A – Sparse matrix that will be decomposed

Returns:

combined L and U Sparse matrices

template<class T = double>
class Matrix

#include <micm/util/matrix.hpp>

A 2D array class with contiguous memory.

Public Functions

inline std::size_t RowStride() const

Get the number of elements in the underlying vector between adjacent rows for the same column.

Returns:

The number of elements in the underlying vector between adjacent rows for the same column

inline std::size_t ColumnStride() const

Get the number of elements in the underlying vector between adjacent columns for the same row.

Returns:

The number of elements in the underlying vector between adjacent columns for the same row

inline void Fill(T val)

Set every matrix element to a given value.

Parameters:

val – Value to set each element to

inline void Axpy(const double &alpha, const Matrix &x)

For each element in the Matrix x and y, perform y = alpha * x + y, where alpha is a scalar constant.

Parameters:

• alpha – The scaling scalar to apply to the Matrix x

• x – The input Matrix

inline void Max(const T &x)

For each element of the matrix, perform y = max(y, x), where x is a scalar constant.

Parameters:

x – The scalar constant to compare against

inline void Min(const T &x)

For each element of the matrix, perform y = min(y, x), where x is a scalar constant.

Parameters:

x – The scalar constant to compare against

inline ConstColumnView GetConstColumnView(std::size_t column_index) const

Create a const column view for accessing a column.

Parameters:

column_index – The index of the column

Returns:

A ConstColumnView descriptor

inline ColumnView GetColumnView(std::size_t column_index)

Create a mutable column view for accessing a column.

Parameters:

column_index – The index of the column

Returns:

A ColumnView descriptor

inline RowVariable GetRowVariable()

Get a row variable with persistent storage for temporary values.

Returns:

A RowVariable with stack-allocated storage

inline RowVariable GetRowVariable() const

Get a row variable with persistent storage for temporary values (const version)

Returns:

A RowVariable with stack-allocated storage

template<typename Func, typename ...Args>
inline void ForEachRow(Func &&func, Args&&... args)

Apply a function to each row of the matrix.

Template Parameters:

• Func – The lambda/function type

• Args – The types of the column view arguments

Parameters:

• func – The function to apply to each row

• args – Column views or row variables

template<typename Func, typename ...Args>
inline void ForEachRow(Func &&func, Args&&... args) const

Apply a function to each row of the matrix (const version)

Template Parameters:

• Func – The lambda/function type

• Args – The types of the column view arguments

Parameters:

• func – The function to apply to each row

• args – Column views or row variables

Public Static Functions

template<typename Func, typename ...Args>
static inline auto Function(Func &&func, Args&... args)

Create a function that can be applied to matrices and vectors.

Creates a reusable callable that validates matrix dimensions and applies a user function row-by-row. For standard Matrix (L=1), each row is processed individually.

Note

Validation occurs in two phases:

1. At function creation: Validates row counts match across all matrices and vector sizes

2. At invocation: Re-validates dimensions in case matrices/vectors were resized

Note

Column view creation happens inside user lambda and is validated at invocation time, not at function creation time. Ensure all column indices are within matrix bounds to avoid runtime errors.

Template Parameters:

• Func – The lambda/function type

• Args – The matrix and vector types

Parameters:

• func – The function to wrap - receives GroupView objects for matrices and vectors

• args – The matrices and vectors to validate and capture dimensions from

Throws:

std::system_error – if column counts don’t match at creation, vectors have wrong sizes at creation, or if at invocation time: matrices/vectors have mismatched row counts, column counts don’t match creation, or column indices are out of bounds

Returns:

A callable that validates dimensions and applies the function

class ColumnView

#include <micm/util/matrix.hpp>

A lightweight descriptor for a mutable column in a matrix.

class ConstColumnView

#include <micm/util/matrix.hpp>

A lightweight descriptor for a const column in a matrix.

class ConstGroupView

#include <micm/util/matrix.hpp>

ConstGroupView provides a const view of a single row (group of size 1) for iteration.

class GroupView

#include <micm/util/matrix.hpp>

GroupView provides a view of a single row (group of size 1) for iteration.

class RowVariable

#include <micm/util/matrix.hpp>

A row-local temporary variable with its own storage.

struct MicmException : public std::runtime_error

class Phase

#include <micm/system/phase.hpp>

Represents a chemical phase (e.g., gaseous, aqueous) Each phase defines a set of species that participate in chemical reactions within that phase.

Public Functions

Phase() = default

Defaulted constructors and assignment operators.

inline Phase(const std::string &name, const std::vector<PhaseSpecies> &phase_species)

Create a phase with a name and a set of species.

inline std::size_t StateSize() const

Returns the number of non-parameterized species.

inline std::vector<std::string> UniqueNames() const

Returns a set of unique names for each non-parameterized species.

inline std::vector<std::string> SpeciesNames() const

Returns a set of unique names for each non-parameterized species (excludes phase name prefix)

Public Members

std::vector<PhaseSpecies> phase_species_

The list of phase-specific species.

class PhaseSpecies

#include <micm/system/phase.hpp>

Represents a chemical species within a specific phase, storing the species information and its optional diffusion coefficient.

class Process

template<typename DenseMatrixPolicy, typename SparseMatrixPolicy>
class ProcessSet

#include <micm/process/process_set.hpp>

Solver function calculators for a collection of processes.

Template Parameters:

• DenseMatrixPolicy – Policy for dense matrices

• SparseMatrixPolicy – Policy for sparse matrices

Subclassed by micm::CudaProcessSet< DenseMatrixPolicy, SparseMatrixPolicy >

Public Functions

ProcessSet() = default

Default constructor.

inline ProcessSet(const std::vector<Process> &processes, const std::unordered_map<std::string, std::size_t> &variable_map)

Constructs a ProcessSet by mapping species in each process to their corresponding indices Initializes internal data structures related to a set of processes, mapping them to variable indices using a provided variable_map. Also prepares the data needed for computing Jacobian contributions.

Parameters:

• processes – A list of processes, each with reactants and products

• variable_map – A map from species names to their corresponding index in the solver’s state

Throws:

std::system_error – If a reactant or product name in a process is not found in variable_map

inline ProcessSet(const std::vector<Process> &processes, const std::unordered_map<std::string, std::size_t> &variable_map, const std::vector<ExternalModelProcessSet<DenseMatrixPolicy, SparseMatrixPolicy>> &external_process_sets)

Constructs a ProcessSet as above, but also includes contributions from external models.

Parameters:

• processes – A list of processes, each with reactants and products

• variable_map – A map from species names to their corresponding index in the solver’s state

• external_process_sets – A list of external process sets that provide additional processes and Jacobian contributions

Throws:

std::system_error – If a reactant or product name in a process is not found in variable_map

inline std::set<std::pair<std::size_t, std::size_t>> NonZeroJacobianElements() const

Returns the positions of all non-zero Jacobian elements.

Returns:

A set of (row, column) index pairs, each representing a non-zero entry

template<typename OrderingPolicy>
inline void SetJacobianFlatIds(const SparseMatrix<double, OrderingPolicy> &matrix)

Computes and stores flat (1D) indices for non-zero Jacobian elements Stores combination of process ids and reactant ids to support column-wise Jacobian updates.

Parameters:

matrix – The sparse Jacobian matrix used to compute flat indices.

inline void SetAlgebraicVariableIds(const std::set<std::size_t> &variable_ids)

Marks species rows that should be treated as algebraic (constraints replace ODE rows)

Parameters:

variable_ids – Set of variable ids whose forcing/Jacobian rows should not receive kinetic contributions

inline void SetExternalModelFunctions(const std::unordered_map<std::string, std::size_t> &state_parameter_indices, const std::unordered_map<std::string, std::size_t> &state_variable_indices, const SparseMatrixPolicy &jacobian)

Sets external model functions for forcing terms and Jacobian contributions.

Parameters:

• state_parameter_indices – Map of state parameter names to their indices

• state_variable_indices – Map of state variable names to their indices

• jacobian – The sparse Jacobian matrix used by the solver

inline void AddForcingTerms(const auto &state, const DenseMatrixPolicy &state_variables, DenseMatrixPolicy &forcing) const

Adds forcing terms for the set of processes for the current conditions.

Parameters:

• state – Current state containing rate constants and other relevant data

• state_variables – Current state variable values (grid cell, state variable)

• forcing – Forcing terms for each state variable (grid cell, state variable)

inline void SubtractJacobianTerms(const auto &state, const DenseMatrixPolicy &state_variables, SparseMatrixPolicy &jacobian) const

Subtracts Jacobian terms for the set of processes for the current conditions.

Parameters:

• state – Current state containing rate constants and other relevant data

• state_variables – Current state variable values (grid cell, state variable)

• jacobian – Jacobian matrix for the system (grid cell, dependent variable, independent variable)

inline std::set<std::string> SpeciesUsed(const std::vector<Process> &processes) const

Extracts all species involved in the given processes.

Parameters:

processes – A list of Process objects, each with reactants and products

Returns:

A set of species names

struct ReactionRateConstantStore

#include <micm/process/reaction_rate_store.hpp>

Structure-of-arrays store for all reaction rate constant parameters.

Processes must be sorted by RateConstantTypeOrder before BuildFrom is called.

Public Static Functions

static inline ReactionRateConstantStore BuildFrom(std::vector<Process> &processes)

Build a ReactionRateConstantStore from a sorted process list.

Parameters:

processes – Non-const ref so LambdaRateConstantParameters pointers remain mutable at runtime.

template<class StatePolicy>
static inline void EvaluateCpuRateConstants(const ReactionRateConstantStore &store, StatePolicy &state)

Evaluate all lambda rate constants into state.rate_constants_. Must be called each step before CpuCalculateRateConstants.

template<class StatePolicy>
static inline void CpuCalculateRateConstants(const ReactionRateConstantStore &store, StatePolicy &state)

Calculate all analytic rate constants into state.rate_constants_. Lambda entries are untouched; parameterized multipliers applied last. Uses DenseMatrixPolicy::Function so a single implementation handles both Matrix (scalar) and VectorMatrix (interleaved) layouts. Each reaction type is computed across all cells at once via ForEachRow, which is more SIMD-friendly than the previous per-cell loop.

struct LambdaEntry

Public Members

LambdaRateConstantParameters *source_

Non-owning; valid for the lifetime of the owning Solver.

std::size_t rc_index_

Column index in state.rate_constants_[cell].

struct ParameterizedMultiplier

#include <micm/process/reaction_rate_store.hpp>

One entry per reaction with at least one parameterized reactant.

struct ReversibleRateConstantParameters

Public Members

double A_ = {1}

Pre-exponential factor [(mol m−3)^(−(𝑛−1)) s−1].

double C_ = {0}

Activation threshold [K].

double k_r_ = {0}

Reverse rate constant [s−1], indicating how fast the species leaves the condensed phase to return to the gas phase.

template<class RatesPolicy, class LinearSolverPolicy, class ConstraintSetPolicy>
class RosenbrockSolver : public micm::AbstractRosenbrockSolver<RatesPolicy, LinearSolverPolicy, ConstraintSetPolicy, RosenbrockSolver<RatesPolicy, LinearSolverPolicy, ConstraintSetPolicy>>

Public Functions

inline RosenbrockSolver(LinearSolverPolicy &&linear_solver, RatesPolicy &&rates, ConstraintSetPolicy &&constraints)

Default constructor.

Note: This constructor is not intended to be used directly. Instead, use the SolverBuilder to create a solver

Parameters:

• linear_solver – Linear solver

• rates – Rates calculator

• constraints – Algebraic constraints

struct RosenbrockSolverParameters

#include <micm/solver/rosenbrock_solver_parameters.hpp>

Rosenbrock solver parameters.

Subclassed by micm::CudaRosenbrockSolverParameters

Public Static Functions

static inline RosenbrockSolverParameters TwoStageRosenbrockParameters()

an L-stable method, 2 stages, order 2

Returns:

static inline RosenbrockSolverParameters ThreeStageRosenbrockParameters()

an L-stable method, 3 stages, order 3, 2 function evaluations

Parameters:

reorder_state –

Returns:

static inline RosenbrockSolverParameters FourStageRosenbrockParameters()

L-stable rosenbrock method of order 4, with 4 stages.

Returns:

static inline RosenbrockSolverParameters FourStageDifferentialAlgebraicRosenbrockParameters()

A stiffly-stable method, 4 stages, order 3.

Returns:

static inline RosenbrockSolverParameters SixStageDifferentialAlgebraicRosenbrockParameters()

stiffly-stable rosenbrock method of order 4, with 6 stages

Returns:

template<class DenseMatrixPolicy>
class RosenbrockTemporaryVariables : public micm::TemporaryVariables

Public Functions

inline virtual std::unique_ptr<TemporaryVariables> Clone() const override

Clone this object, preserving the derived type.

struct SimpleGroupingTag

#include <micm/util/view_category.hpp>

Simple grouping: L==1, group index directly maps to element Used by: Matrix (always), VectorMatrix (when L==1), Standard ordering sparse (always), Vector ordering sparse (when L==1)

template<class SolverPolicy, class StatePolicy>
class Solver

Public Functions

inline std::size_t MaximumNumberOfGridCells() const

Returns the maximum number of grid cells per state.

This is the maximum number of grid cells that can fit within one group for vectorized solvers. For non-vectorized solvers, there is no limit other than the maximum size of a std::size_t.

Returns:

Number of grid cells

inline void UpdateStateParameters(StatePolicy &state)

Update state parameters based on current conditions (temperature, pressure, etc.) Invokes all registered parameter update functions for external models and constraints to recompute temperature-dependent values (e.g., aerosol rate constants, equilibrium constants) Should be called before solving if conditions have changed since the last solve.

Parameters:

state – State object containing conditions and custom_rate_parameters to be updated

template<class SolverParametersPolicy, class DenseMatrixPolicy, class SparseMatrixPolicy, class RatesPolicy, class LuDecompositionPolicy, class LinearSolverPolicy, class StatePolicy>
class SolverBuilder

#include <micm/solver/solver_builder.hpp>

Builder of general solvers.

Template Parameters:

• SolverParametersPolicy – Policy for the ODE solver

• DenseMatrixPolicy – Policy for dense matrices

• SparseMatrixPolicy – Policy for sparse matrices

• RatesPolicy – Calculator of forcing and Jacobian terms

• LinearSolverPolicy – Policy for the linear solver

Public Functions

inline SolverBuilder &SetSystem(const System &system)

Set the chemical system.

Parameters:

system – The chemical system

Returns:

Updated SolverBuilder

inline SolverBuilder &SetReactions(const std::vector<Process> &reactions)

Set the reactions.

Parameters:

reactions – The reactions

Returns:

Updated SolverBuilder

inline SolverBuilder &SetConstraints(std::vector<Constraint> &&constraints)

Set algebraic constraints for DAE solving.

Parameters:

constraints – Vector of constraints

Returns:

Updated SolverBuilder

inline SolverBuilder &SetIgnoreUnusedSpecies(bool ignore_unused_species)

Set whether to ignore unused species.

Parameters:

ignore_unused_species – True if unused species should be ignored

Returns:

Updated SolverBuilder

inline SolverBuilder &SetReorderState(bool reorder_state)

Set whether to reorder the state to optimize the LU decomposition.

Parameters:

reorder_state – True if the state should be reordered

Returns:

Updated SolverBuilder

template<class ExternalModel>
inline SolverBuilder &AddExternalModel(ExternalModel model)

Add an external model (state variables, processes, and/or constraints)

If the model satisfies HasState, its state variables and parameters are registered with the solver. The model must satisfy at least one of HasProcesses (process wrappers are created) or HasConstraints (constraint wrappers are created).

Parameters:

model – The external model (taken by value; caller decides whether to copy or move)

Returns:

Updated SolverBuilder

inline auto Build()

Creates an instance of Solver with a properly configured ODE solver.

Returns:

An instance of Solver

struct SolverResult

Public Members

SolverState state_ = SolverState::NotYetCalled

The final state the solver was in.

SolverStats stats_ = {}

A collection of runtime state for this call of the solver.

struct SolverStats

Public Members

uint64_t function_calls_ = {}

The number of forcing function calls.

uint64_t jacobian_updates_ = {}

The number of jacobian function calls.

uint64_t number_of_steps_ = {}

The total number of internal time steps taken.

uint64_t accepted_ = {}

The number of accepted integrations.

uint64_t rejected_ = {}

The number of rejected integrations.

uint64_t decompositions_ = {}

The number of LU decompositions.

uint64_t solves_ = {}

The number of linear solves.

uint64_t constraint_init_iterations_ = {}

The number of constraint initialization iterations performed.

double final_time_ = {}

The final time the solver iterated to.

template<class T = double, class OrderingPolicy>
class SparseMatrix : public OrderingPolicy

#include <micm/util/sparse_matrix.hpp>

A sparse block-diagonal 2D matrix class with contiguous memory.

Each block sub-matrix is square and has the same structure of non-zero elements

The template parameters are the type of the matrix elements and a class that defines the sizing and ordering of the data elements

Public Types

using BlockVariable = typename OrderingPolicy::template BlockVariable<T>

Alias for the ordering policy’s BlockVariable type.

using ConstGroupView = typename OrderingPolicy::template ConstGroupView<SparseMatrix>

Alias for the ordering policy’s ConstGroupView type.

using GroupView = typename OrderingPolicy::template GroupView<SparseMatrix>

Alias for the ordering policy’s GroupView type.

Public Functions

inline SparseMatrix(const SparseMatrixBuilder<T, OrderingPolicy> &builder, bool indexing_only = false)

Constructs a SparseMatrix from a given builder and optional indexing mode. Initializes the SparseMatrix using the provided SparseMatrixBuilder, which defines the matrix structure, block size, and non-zero elements. Optionally, the constructor can be used in “indexing only” mode, where the data storage is not allocated.

Template Parameters:

• T – The type of the matrix elements.

• OrderingPolicy – The policy class that defines the ordering and storage of elements.

Parameters:

• builder – The builder object containing matrix configuration and initial values.

• indexing_only – If true, only indexing structures are initialized and data storage is omitted.

inline void Fill(T val)

Set every matrix element to a given value.

Parameters:

val – Value to set each element to

inline void PrintNonZeroElements(std::ostream &os) const

Print the sparse matrix with row index, column index, and non-zero value; useful to test other linear algebra libraries.

Parameters:

os – Output stream to print to, defaults to std::cout

inline ConstBlockView GetConstBlockView(std::size_t vector_index) const

Create a const block view for accessing the nth non-zero element.

Parameters:

vector_index – The data array index from VectorIndex(0, row, col) for the element

Returns:

A ConstBlockView descriptor

inline ConstBlockView GetConstBlockView(std::size_t row, std::size_t col) const

Create a const block view for accessing a block element.

Parameters:

• row – The row index of the block element

• col – The column index of the block element

Returns:

A ConstBlockView descriptor

inline BlockView GetBlockView(std::size_t vector_index)

Create a mutable block view for accessing the nth non-zero element.

Parameters:

vector_index – The data array index from VectorIndex(0, row, col) for the element

Returns:

A BlockView descriptor

inline BlockView GetBlockView(std::size_t row, std::size_t col)

Create a mutable block view for accessing a block element.

Parameters:

• row – The row index of the block element

• col – The column index of the block element

Returns:

A BlockView descriptor

inline BlockVariable GetBlockVariable()

Get a block variable with persistent storage for temporary values.

Returns:

A BlockVariable with stack-allocated storage

template<typename Func, typename ...Args>
inline void ForEachBlock(Func &&func, Args&&... args)

Apply a function to each block of the matrix.

Template Parameters:

• Func – The lambda/function type

• Args – The types of the block view arguments

Parameters:

• func – The function to apply to each block

• args – Block views or block variables

Public Static Functions

template<typename Func, typename ...Args>
static inline auto Function(Func &&func, Args&... args)

Create a function that can be applied to sparse matrices and vectors.

Creates a reusable callable that validates dimensions and applies a user function across block groups. The function iterates over groups of L blocks at a time, where L is determined by the OrderingPolicy::GroupVectorSize(). Supports mixing sparse matrices, dense matrices, and vector-like types.

Note

Validation occurs in two phases:

1. At function creation: Validates matrix dimensions, vector sizes, and ordering compatibility

2. At invocation: Re-validates dimensions in case matrices/vectors were resized

Note

Column/Block view creation happens inside user lambda and is validated at invocation time, not at function creation time. Ensure all view indices are within matrix bounds to avoid runtime errors.

Template Parameters:

• Func – The lambda/function type

• Args – The matrix and vector types (can mix SparseMatrix, VectorMatrix, Matrix, and vectors)

Parameters:

• func – The function to wrap - receives GroupView objects for matrices and forwarded vectors

• args – The matrices and vectors to validate and capture dimensions from

Throws:

std::system_error – if matrices have incompatible orderings (different L values), mismatched block counts, or vectors have wrong sizes

Returns:

A callable that validates dimensions and applies the function

class BlockView

#include <micm/util/sparse_matrix.hpp>

A lightweight descriptor for a mutable block element in a sparse matrix.

class ConstBlockView

#include <micm/util/sparse_matrix.hpp>

A lightweight descriptor for a const block element in a sparse matrix.

struct SparseMatrixBlockViewTag

#include <micm/util/view_category.hpp>

Tag for sparse matrix block views (have RowIndex + ColumnIndex)

template<class T, class OrderingPolicy = SparseMatrixStandardOrdering>
class SparseMatrixBuilder

class SparseMatrixStandardOrderingCompressedSparseColumn

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_column.hpp>

Defines the ordering of SparseMatrix object data in Compressed Sparse Column format.

Data is stored with blocks in the block diagonal matrix as the highest level structure, then by column, then by non-zero rows in each column.

Public Functions

inline std::size_t GroupSize(std::size_t number_of_non_zero_elements) const

Returns the size of each group of blocks in the compressed data vector.

Returns:

Size of each group of blocks

inline std::size_t NumberOfGroups(std::size_t number_of_blocks) const

Returns the total number of groups of blocks in the compressed data.

Parameters:

number_of_blocks – Total number of block sub-matrices in the overall matrix

Returns:

Number of groups of blocks (equal to number_of_blocks for standard ordering)

inline bool IsZero(std::size_t row, std::size_t column) const

Returns whether a particular element is always zero.

Parameters:

• row – Row index

• column – Column index

Returns:

true if the element is always zero, false otherwise

Public Static Functions

static inline std::size_t GroupVectorSize()

Returns the number of blocks included in each group of blocks.

Returns:

Number of blocks in each group (1 for standard ordering)

template<typename T>
class BlockVariable

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_column.hpp>

A block-local temporary variable with its own storage For standard ordering: single value.

template<typename SparseMatrixType>
class ConstGroupView

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_column.hpp>

ConstGroupView provides a const view of a single group of blocks for iteration For standard ordering: L=1, so each group contains 1 block.

template<typename SparseMatrixType>
class GroupView

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_column.hpp>

GroupView provides a view of a single group of blocks for iteration For standard ordering: L=1, so each group contains 1 block.

class SparseMatrixStandardOrderingCompressedSparseRow

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_row.hpp>

Defines the ordering of SparseMatrix object data in Compressed Sparse Row format.

Data is stored with blocks in the block diagonal matrix as the highest level structure, then by row, then by non-zero columns in each row.

Subclassed by micm::SparseMatrix< double, SparseMatrixStandardOrdering >

Public Functions

inline std::size_t GroupSize(std::size_t number_of_non_zero_elements) const

Returns the size of each group of blocks in the compressed data vector.

Returns:

Size of each group of blocks

inline std::size_t NumberOfGroups(std::size_t number_of_blocks) const

Returns the total number of groups of blocks in the compressed data.

Parameters:

number_of_blocks – Total number of block sub-matrices in the overall matrix

Returns:

Number of groups of blocks (equal to number_of_blocks for standard ordering)

inline bool IsZero(std::size_t row, std::size_t column) const

Returns whether a particular element is always zero.

Parameters:

• row – Row index

• column – Column index

Returns:

true if the element is always zero, false otherwise

Public Static Functions

static inline std::size_t GroupVectorSize()

Returns the number of blocks included in each group of blocks.

Returns:

Number of blocks in each group (1 for standard ordering)

template<typename T>
class BlockVariable

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_row.hpp>

A block-local temporary variable with its own storage For standard ordering: single value.

template<typename SparseMatrixType>
class ConstGroupView

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_row.hpp>

ConstGroupView provides a const view of a single group of blocks for iteration For standard ordering: L=1, so each group contains 1 block.

template<typename SparseMatrixType>
class GroupView

#include <micm/util/sparse_matrix_standard_ordering_compressed_sparse_row.hpp>

GroupView provides a view of a single group of blocks for iteration For standard ordering: L=1, so each group contains 1 block.

template<std::size_t L = MICM_DEFAULT_VECTOR_SIZE>
class SparseMatrixVectorOrderingCompressedSparseColumn

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_column.hpp>

Defines the ordering of SparseMatrix object data in Compressed Sparse Column format into blocks of rows to encourage vectorization.

Data is stored with sets of blocks in the block diagonal matrix as the highest level structure, then by column, then by non-zero rows in each column, then by individual blocks in the set of blocks.

The template argument is the number of blocks per set of blocks and should be approximately the size of the vector register.

Public Functions

inline std::size_t GroupSize() const

Returns the size of each group of blocks in the compressed data vector.

Parameters:

number_of_non_zero_elements – Number of non-zero elements in the matrix

Returns:

Size of each group of blocks

inline std::size_t NumberOfGroups(std::size_t number_of_blocks) const

Returns the total number of groups of blocks in the compressed data vector, including any partial groups.

Parameters:

number_of_blocks – Total number of block sub-matrices in the overall matrix

Returns:

Number of groups of blocks

inline bool IsZero(const std::size_t row, const std::size_t column) const

Returns whether a given row and column index is a zero element.

Parameters:

• row – Index of the row

• column – Index of the column

Returns:

True if the element is zero, false otherwise

Public Static Functions

static inline std::size_t GroupVectorSize()

Returns the number of blocks included in each group of blocks.

Returns:

Number of blocks per group

template<typename T>
class BlockVariable

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_column.hpp>

A block-local temporary variable with its own storage For vector ordering: array of L values when L>1, single value when L=1.

template<typename SparseMatrixType>
class ConstGroupView

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_column.hpp>

ConstGroupView provides a const view of a single group of blocks for iteration For vector ordering: each group contains L blocks (except possibly the last group)

template<typename SparseMatrixType>
class GroupView

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_column.hpp>

GroupView provides a view of a single group of blocks for iteration For vector ordering: each group contains L blocks (except possibly the last group)

template<std::size_t L = MICM_DEFAULT_VECTOR_SIZE>
class SparseMatrixVectorOrderingCompressedSparseRow

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_row.hpp>

Defines the ordering of SparseMatrix object data in Compressed Sparse Row format into blocks of rows to encourage vectorization.

Data is stored with sets of blocks in the block diagonal matrix as the highest level structure, then by row, then by non-zero columns in each row, then by individual blocks in the set of blocks.

The template argument is the number of blocks per set of blocks and should be approximately the size of the vector register.

Public Functions

inline std::size_t GroupSize() const

Returns the size of each group of blocks in the compressed data vector.

Parameters:

number_of_non_zero_elements – Number of non-zero elements in the matrix

Returns:

Size of each group of blocks

inline std::size_t NumberOfGroups(std::size_t number_of_blocks) const

Returns the total number of groups of blocks in the compressed data vector, including any partial groups.

Parameters:

number_of_blocks – Total number of block sub-matrices in the overall matrix

Returns:

Number of groups of blocks

inline bool IsZero(std::size_t row, std::size_t column) const

Returns whether a particular element is always zero.

Parameters:

• row – Row index

• column – Column index

Returns:

true if the element is always zero, false otherwise

Public Static Functions

static inline std::size_t GroupVectorSize()

Returns the number of blocks included in each group of blocks.

Returns:

Number of blocks in each group

template<typename T>
class BlockVariable

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_row.hpp>

A block-local temporary variable with its own storage For vector ordering: array of L values when L>1, single value when L=1.

template<typename SparseMatrixType>
class ConstGroupView

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_row.hpp>

ConstGroupView provides a const view of a single group of blocks for iteration For vector ordering: each group contains L blocks (except possibly the last group)

template<typename SparseMatrixType>
class GroupView

#include <micm/util/sparse_matrix_vector_ordering_compressed_sparse_row.hpp>

GroupView provides a view of a single group of blocks for iteration For vector ordering: each group contains L blocks (except possibly the last group)

class Species

#include <micm/system/species.hpp>

A representation of a chemcial species.

Public Functions

Species() = default

Default constructor.

inline Species &operator=(const Species &other)

Copy assignment.

Parameters:

other – species to copy

inline Species(const Species &other)

Copy Constructor.

Parameters:

other –

inline Species(const std::string &name)

Construct a species by name only.

Parameters:

name – The name of the species

inline Species(const std::string &name, const std::map<std::string, double> &properties)

Construct a species by name and properties.

Parameters:

• name – The name of the species

• properties – The properties of the species

inline bool IsParameterized() const

Returns whether a species is parameterized.

inline void SetThirdBody()

Set a Species instance parameterized on air density.

template<class T>
inline T GetProperty(const std::string &key) const

Return the value of a species property.

template<class T>
inline void SetProperty(const std::string &key, T value)

Set the value of a species property.

Public Members

std::string name_

The name of this species.

std::map<std::string, std::string> properties_string_

A list of properties of this species.

std::function<double(const Conditions)> parameterize_ = {nullptr}

A function that if provided will be used to parameterize the concentration of this species during solving. Species with this function defined will be excluded from the solver state.

template<class DenseMatrixPolicy = StandardDenseMatrix, class SparseMatrixPolicy = StandardSparseMatrix, class LuDecompositionPolicy = LuDecomposition, class LMatrixPolicy = SparseMatrixPolicy, class UMatrixPolicy = SparseMatrixPolicy>
struct State

Subclassed by micm::CudaState< CudaDenseMatrixVector, CudaSparseMatrixVector, CudaLuDecompositionMozartInPlace >

Public Types

using DenseMatrixPolicyType = DenseMatrixPolicy

Type of the DenseMatrixPolicy.

Public Functions

inline State()

Default constructor Only defined to be used to create default values in types, but a default constructed state is not useable.

inline State(const StateParameters &parameters, const std::size_t number_of_grid_cells)

Constructor with parameters.

Parameters:

parameters – State dimension information

inline State(const State &other)

Copy constructor.

Parameters:

other – The state object to be copied

inline State &operator=(const State &other)

Assignment operator.

Parameters:

other – The state object to be assigned

Returns:

Reference to the assigned state object

inline State(State &&other) noexcept

Move constructor.

Parameters:

other – The state object to be moved

inline State &operator=(State &&other) noexcept

Move assignment operator.

Parameters:

other – The state object to be moved

Returns:

Reference to the moved state object

inline std::size_t NumberOfGridCells() const

Get the number of grid cells.

Returns:

The number of grid cells

inline VariableProxy operator[](std::size_t index)

Square-bracket access operator for state variable index.

Parameters:

index – The index of the variable to access

Returns:

Reference to the variable matrix column corresponding to the given index

inline ConstVariableProxy operator[](std::size_t index) const

Square-bracket access operator for state variable index (const version)

Parameters:

index – The index of the variable to access

Returns:

Const reference to the variable matrix column corresponding to the given index

inline VariableProxy operator[](const std::string &name)

Square-bracket access operator for unique variable name.

Parameters:

name – The unique name of the variable to access

Returns:

VariableProxy proxy object providing access to the values of the named variable (e.g., its concentration) across grid cells. This is a proxy, not a direct reference to an internal matrix column; see VariableProxy documentation for details on single- vs multi-cell access patterns.

inline ConstVariableProxy operator[](const std::string &name) const

Square-bracket access operator for unique variable name (const version)

Parameters:

name – The unique name of the variable to access

Returns:

ConstVariableProxy proxy object providing read-only access to the values of the named variable across grid cells. This is a proxy, not a direct reference to an internal matrix column; see ConstVariableProxy documentation for details on single- vs multi-cell access patterns.

inline VariableProxy operator[](const Species &species)

Square-bracket access operator for species object.

Parameters:

species – The species object corresponding to the variable to access

Returns:

VariableProxy proxy object providing access to the values of the given species across grid cells. This is a proxy, not a direct reference to an internal matrix column; see VariableProxy documentation for details on single- vs multi-cell access patterns.

inline ConstVariableProxy operator[](const Species &species) const

Square-bracket access operator for species object (const version)

Parameters:

species – The species object corresponding to the variable to access

Returns:

ConstVariableProxy proxy object providing read-only access to the values of the given species across grid cells. This is a proxy, not a direct reference to an internal matrix column; see ConstVariableProxy documentation for details on single- vs multi-cell access patterns.

inline void SetConcentrations(const std::unordered_map<std::string, std::vector<double>> &species_to_concentration)

Set species’ concentrations.

Parameters:

species_to_concentration –

inline void SetConcentration(const Species &species, double concentration)

Set a single species concentration.

Deprecated:

This method is deprecated in favor of using the operator[] with species or name to set concentrations, e.g., state[species] = concentration or state[“species_name”] = concentration

Parameters:

• species – the species to set the concentration for

• concentration – concentration [mol m-3]

inline void SetConcentration(const Species &species, const std::vector<double> &concentration)

Set concentrations for a single species across multiple grid cells.

Deprecated:

This method is deprecated in favor of using the operator[] with species or name to set concentrations, e.g., state[species] = concentrations or state[“species_name”] = concentrations

Parameters:

• species – the species to set the concentrations for

• concentration – vector of concentrations [mol m-3], one per grid cell

inline void SetConcentration(const std::string &element, double concentration)

Set the concentration for a named element (species or other variable)

Deprecated:

This method is deprecated in favor of using the operator[] with species or name to set concentrations, e.g., state[species] = concentration or state[“species_name”] = concentration

Parameters:

• species – the name of the element (can be a non-species variable, e.g., number_concentration)

• concentration – concentration value [mol m-3]

inline void SetConcentration(const std::string &element, const std::vector<double> &concentration)

Set concentrations for a named element (species or other variable) across multiple grid cells.

Deprecated:

This method is deprecated in favor of using the operator[] with species or name to set concentrations, e.g., state[species] = concentrations or state[“species_name”] = concentrations

Parameters:

• species – the name of the element (can be a non-species variable, e.g., number_concentration)

• concentration – vector of concentrations [mol m-3], one per grid cell

inline void UnsafelySetCustomRateParameters(const std::vector<std::vector<double>> &parameters)

Set custom parameters assuming the values are properly ordered.

Parameters:

parameters – map of custom rate parameters

inline void SetCustomRateParameters(const std::unordered_map<std::string, std::vector<double>> &parameters)

Set custom parameters for rate constant calculations by label.

Parameters:

parameters – map of custom rate parameters

inline void SetCustomRateParameter(const std::string &label, double value)

Set a single custom rate constant parameter.

Parameters:

• label – parameter label

• value – new parameter value

inline void SetRelativeTolerance(double relativeTolerance)

Set the relative tolerances.

Parameters:

relativeTolerance – relative tolerance

inline virtual void SetAbsoluteTolerances(const std::vector<double> &absoluteTolerance)

Set the absolute tolerances per species.

Parameters:

absoluteTolerance – absolute tolerance

inline void PrintHeader()

Print a header of species to display concentrations with respect to time.

inline void PrintState(double time)

Print state (concentrations) at the given time.

Parameters:

time – solving time

Public Members

std::size_t number_of_grid_cells_ = {1}

The number of grid cells stored in the state.

DenseMatrixPolicy variables_

The concentration of chemicals, varies through time.

DenseMatrixPolicy custom_rate_parameters_

Rate parameters particular to user-defined rate constants, may vary in time.

DenseMatrixPolicy rate_constants_

The reaction rates, may vary in time.

std::vector<Conditions> conditions_

Atmospheric conditions, varies in time.

std::vector<double> upper_left_identity_diagonal_

The block matrix with an upper left identity, zeros elsewhere.

SparseMatrixPolicy jacobian_

The jacobian structure, varies for each solve.

std::unordered_map<std::string, std::size_t> variable_map_

Immutable data required for the state.

class ConstVariableProxy

class VariableProxy

struct StateParameters

#include <micm/solver/state.hpp>

Invariants that can be used to construct a state.

struct StoichSpecies

#include <micm/system/stoich_species.hpp>

Represents a species in a chemical reaction, defined by its stoichiometric coefficient.

struct SurfaceRateConstantData

#include <micm/process/rate_constant/surface_rate_constant.hpp>

GPU-safe calculation data for a surface reaction. Populated by ReactionRateConstantStore::BuildFrom from SurfaceRateConstantParameters; do not construct directly.

Public Members

double diffusion_coefficient_

Gas-phase diffusion coefficient for the reacting species [m2 s-1].

double mean_free_speed_factor_

Precomputed factor for mean free speed: 8 * R / (pi * Mw) [K-1 m2 s-2].

double reaction_probability_

Reaction probability (0-1) [unitless].

std::size_t custom_param_base_index_

Index into custom_rate_parameters_[cell] for aerosol effective radius [m]; particle number concentration [# m-3] is at custom_param_base_index_ + 1.

struct SurfaceRateConstantParameters

Public Members

std::string label_

Label for the reaction used to identify user-defined parameters.

PhaseSpecies phase_species_

Gas-phase species reacting on surface.

double reaction_probability_ = {1.0}

Reaction probability (0-1) [unitless].

class System

#include <micm/system/system.hpp>

Defines the gas-phase species available in the chemical system.

Public Functions

inline size_t StateSize() const

Returns the number of gas-phase state variables.

inline std::vector<std::string> UniqueNames() const

Returns the unique gas-phase species names.

Returns:

vector of unique state variable names

inline std::vector<std::string> UniqueNames(const std::function<std::string(const std::vector<std::string> &variables, const std::size_t i)> f) const

Returns the unique gas-phase species names, optionally reordered.

Parameters:

f – Function used to apply a specific order to unique names

Returns:

vector of unique state variable names

Public Members

Phase gas_phase_

The gas phase, defining a set of species present in the system.

struct TaylorSeriesRateConstantParameters

Public Members

double A_ = {1}

Pre-exponential factor [(mol m−3)^(−(𝑛−1)) s−1].

double B_ = {0}

Unitless exponential factor.

double C_ = {0}

Activation threshold, expected to be the negative activation energy divided by the boltzman constant [-E_a / k_b), K].

double D_ = {300}

A factor that determines temperature dependence [K].

double E_ = {0}

A factor that determines pressure dependence [Pa-1].

double coefficients_[MAX_COEFFICIENTS] = {1.0}

Taylor coefficients for the series expansion. Only the first n_coefficients_ entries are used.

std::size_t n_coefficients_ = {1}

Number of active Taylor coefficients [1, MAX_COEFFICIENTS].

Public Static Attributes

static std::size_t MAX_COEFFICIENTS = 16

Maximum number of Taylor series coefficients supported.

class TemporaryVariables

#include <micm/solver/temporary_variables.hpp>

This is the base class for temporary variables; currently it is empty and will be expanded by a specific solver later.

Subclassed by micm::BackwardEulerTemporaryVariables< DenseMatrixPolicy >, micm::RosenbrockTemporaryVariables< DenseMatrixPolicy >

Public Functions

virtual std::unique_ptr<TemporaryVariables> Clone() const = 0

Clone this object, preserving the derived type.

struct TernaryChemicalActivationRateConstantParameters

Public Members

double k0_A_ = 1.0

low-pressure pre-exponential factor

double k0_B_ = 0.0

low-pressure temperature-scaling parameter

double k0_C_ = 0.0

low-pressure exponential factor

double kinf_A_ = 1.0

high-pressure pre-exponential factor

double kinf_B_ = 0.0

high-pressure temperature-scaling parameter

double kinf_C_ = 0.0

high-pressure exponential factor

double Fc_ = 0.6

TernaryChemicalActivation F_c parameter.

double N_ = 1.0

TernaryChemicalActivation N parameter.

struct TieredGroupingTag

#include <micm/util/view_category.hpp>

Tiered grouping: L>1, groups contain L elements, need block_in_group offset Used by: VectorMatrix (when L>1), Vector ordering sparse (when L>1)

struct TroeRateConstantParameters

Public Members

double k0_A_ = 1.0

low-pressure pre-exponential factor

double k0_B_ = 0.0

low-pressure temperature-scaling parameter

double k0_C_ = 0.0

low-pressure exponential factor

double kinf_A_ = 1.0

high-pressure pre-exponential factor

double kinf_B_ = 0.0

high-pressure temperature-scaling parameter

double kinf_C_ = 0.0

high-pressure exponential factor

double Fc_ = 0.6

Troe F_c parameter.

double N_ = 1.0

Troe N parameter.

struct TunnelingRateConstantParameters

Public Members

double A_ = 1.0

Pre-exponential factor [(mol m−3)^(−(𝑛−1)) s−1].

double B_ = 0.0

Linear temperature-dependent parameter [K].

double C_ = 0.0

Cubed temperature-dependent parameter [K^3].

struct UserDefinedRateConstantData

#include <micm/process/rate_constant/user_defined_rate_constant.hpp>

GPU-safe calculation data for a user-defined rate constant. Populated by ReactionRateConstantStore::BuildFrom from UserDefinedRateConstantParameters; do not construct directly.

Public Members

double scaling_factor_ = {1.0}

Scaling factor applied to the user-provided rate constant value.

std::size_t custom_param_index_ = {0}

Index into custom_rate_parameters_[cell] holding the rate constant value.

struct UserDefinedRateConstantParameters

Public Members

std::string label_

Label for the reaction used to identify user-defined parameters.

double scaling_factor_ = {1.0}

Scaling factor to apply to user-provided rate constants.

struct VantHoffParam

#include <micm/constraint/types/equilibrium_constraint.hpp>

Define parameters for Van’t Hoff equation.

template<class T, std::size_t L = MICM_DEFAULT_VECTOR_SIZE>
class VectorMatrix

#include <micm/util/vector_matrix.hpp>

A 2D array class with contiguous memory structured to encourage vectorization.

The memory layout groups rows into groups whose size can be set such that for a single column, the group of rows can fit in the vector register.

The template arguments are the type of the matrix elements and the size of the number of rows per group.

Subclassed by micm::CudaDenseMatrix< double, MICM_DEFAULT_VECTOR_SIZE >

Public Functions

inline std::size_t RowStride() const

Get the number of elements in the underlying vector between adjacent rows for the same column.

Returns:

The number of elements in the underlying vector between adjacent rows for the same column

inline std::size_t ColumnStride() const

Get the number of elements in the underlying vector between adjacent columns for the same row.

Returns:

The number of elements in the underlying vector between adjacent columns for the same row

inline void Fill(T val)

Set every matrix element to a given value.

Parameters:

val – Value to set each element to

inline void Axpy(const double &alpha, const VectorMatrix &x)

For each element in the VectorMatrix x and y, perform y = alpha * x + y, where alpha is a scalar constant.

Parameters:

• alpha – The scaling scalar to apply to the VectorMatrix x

• x – The input VectorMatrix

inline void Max(const T &x)

For each element of the VectorMatrix, perform y = max(y, x), where x is a scalar constant.

Parameters:

x – The scalar constant to compare against

inline void Min(const T &x)

For each element of the VectorMatrix, perform y = min(y, x), where x is a scalar constant.

Parameters:

x – The scalar constant to compare against

inline ConstColumnView GetConstColumnView(std::size_t column_index) const

Create a const column view for accessing a column.

Parameters:

column_index – The index of the column

Returns:

A ConstColumnView descriptor

inline ColumnView GetColumnView(std::size_t column_index)

Create a mutable column view for accessing a column.

Parameters:

column_index – The index of the column

Returns:

A ColumnView descriptor

inline RowVariable GetRowVariable()

Get a row variable with persistent storage for temporary values.

Returns:

A RowVariable with stack-allocated storage

inline RowVariable GetRowVariable() const

Get a row variable with persistent storage for temporary values (const version)

Returns:

A RowVariable with stack-allocated storage

template<typename Func, typename ...Args>
inline void ForEachRow(Func &&func, Args&&... args)

Apply a function to each row of the matrix (processes L rows at a time)

Template Parameters:

• Func – The lambda/function type

• Args – The types of the column view arguments

Parameters:

• func – The function to apply to each row

• args – Column views or row variables

template<typename Func, typename ...Args>
inline void ForEachRow(Func &&func, Args&&... args) const

Apply a function to each row of the matrix (const version)

Template Parameters:

• Func – The lambda/function type

• Args – The types of the column view arguments

Parameters:

• func – The function to apply to each row

• args – Column views or row variables

Public Static Functions

template<typename Func, typename ...Args>
static inline auto Function(Func &&func, Args&... args)

Create a function that can be applied to vector matrices and vectors.

Creates a reusable callable that validates matrix dimensions and applies a user function across row groups. The function iterates over groups of L rows at a time for vectorization, where L is the compile-time template parameter.

Note

Validation occurs in two phases:

1. At function creation: Validates row counts match across all matrices and vector sizes

2. At invocation: Re-validates dimensions in case matrices/vectors were resized

Note

Column view creation happens inside user lambda and is validated at invocation time, not at function creation time. Ensure all column indices are within matrix bounds to avoid runtime errors.

Template Parameters:

• Func – The lambda/function type

• Args – The matrix and vector types

Parameters:

• func – The function to wrap - receives GroupView objects for matrices and vectors

• args – The matrices and vectors to validate and capture dimensions from

Throws:

std::system_error – if column counts don’t match at creation, or if at invocation time: matrices/vectors have mismatched row counts, column counts don’t match creation, or dimensions mismatch

Returns:

A callable that validates dimensions and applies the function

class ColumnView

#include <micm/util/vector_matrix.hpp>

A lightweight descriptor for a mutable column in a matrix.

class ConstColumnView

#include <micm/util/vector_matrix.hpp>

A lightweight descriptor for a const column in a matrix.

class ConstGroupView

#include <micm/util/vector_matrix.hpp>

ConstGroupView provides a const view of a single group of L rows for iteration.

Public Functions

inline ConstGroupView(const VectorMatrix &matrix, std::size_t group)

Constructor that calculates num_rows_in_group from matrix dimensions.

inline ConstGroupView(const VectorMatrix &matrix, std::size_t group, std::size_t num_rows_in_group)

Constructor with explicit num_rows_in_group.

class GroupView

#include <micm/util/vector_matrix.hpp>

GroupView provides a view of a single group of L rows for iteration.

Public Functions

inline GroupView(VectorMatrix &matrix, std::size_t group)

Constructor that calculates num_rows_in_group from matrix dimensions.

inline GroupView(VectorMatrix &matrix, std::size_t group, std::size_t num_rows_in_group)

Constructor with explicit num_rows_in_group.

class RowVariable

#include <micm/util/vector_matrix.hpp>

A row-local temporary variable with its own storage.

template<typename T, typename = void>
struct ViewCategory

#include <micm/util/view_category.hpp>

Determines the category of a view type (checks for nested ‘category’ type first) Primary template: fallback for types without a nested ‘category’ type (e.g., std::vector) This enables SFINAE to work correctly in concepts like VectorLike.

template<typename T>
struct ViewCategory<T, std::enable_if_t<HasCategory<std::remove_cvref_t<T>>::value>>

#include <micm/util/view_category.hpp>

If type has a nested ‘category’ type, use it.

namespace constants

Variables

static double BOLTZMANN_CONSTANT = 1.380649e-23

static double AVOGADRO_CONSTANT = 6.02214076e23

static double GAS_CONSTANT = BOLTZMANN_CONSTANT * AVOGADRO_CONSTANT

namespace property_keys

Variables

static const char *DIFFUSION_COEFFICIENT = "diffusion coefficient [m2 s-1]"

static const char *MOLECULAR_WEIGHT = "molecular weight [kg mol-1]"