Skip to content

Commit

Permalink
Emission, loss rate constants (#166)
Browse files Browse the repository at this point in the history
* create general user-defined rate constant class

* fix missing include

* fix cuda test

* add photolysis config test

* add tests for emission and loss configs
  • Loading branch information
mattldawson committed Jul 31, 2023
1 parent 371d575 commit d20e289
Show file tree
Hide file tree
Showing 47 changed files with 813 additions and 235 deletions.
63 changes: 35 additions & 28 deletions include/micm/configure/solver_config.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@
#include <iostream>
#include <micm/process/arrhenius_rate_constant.hpp>
#include <micm/process/branched_rate_constant.hpp>
#include <micm/process/photolysis_rate_constant.hpp>
#include <micm/process/user_defined_rate_constant.hpp>
#include <micm/process/process.hpp>
#include <micm/process/ternary_chemical_activation_rate_constant.hpp>
#include <micm/process/troe_rate_constant.hpp>
Expand Down Expand Up @@ -90,7 +90,7 @@ namespace micm
std::vector<Species> species_arr_;

// Read from reaction configure
std::vector<PhotolysisRateConstant> photolysis_rate_arr_;
std::vector<UserDefinedRateConstant> user_defined_rate_arr_;
std::vector<ArrheniusRateConstant> arrhenius_rate_arr_;
std::vector<BranchedRateConstant> branched_rate_arr_;
std::vector<TroeRateConstant> troe_rate_arr_;
Expand Down Expand Up @@ -404,11 +404,11 @@ namespace micm
auto reactants = ParseReactants(object[REACTANTS]);
auto products = ParseProducts(object[PRODUCTS]);

std::string name = object[MUSICA_NAME].get<std::string>();
std::string name = "PHOTO." + object[MUSICA_NAME].get<std::string>();

photolysis_rate_arr_.push_back(PhotolysisRateConstant(name));
user_defined_rate_arr_.push_back(UserDefinedRateConstant(name));

std::unique_ptr<PhotolysisRateConstant> rate_ptr = std::make_unique<PhotolysisRateConstant>(name);
std::unique_ptr<UserDefinedRateConstant> rate_ptr = std::make_unique<UserDefinedRateConstant>(name);
processes_.push_back(Process(reactants, products, std::move(rate_ptr), gas_phase_));

return ConfigParseStatus::Success;
Expand Down Expand Up @@ -664,35 +664,55 @@ namespace micm

ConfigParseStatus ParseEmission(const json& object)
{
std::vector<std::string> required_keys = { "species" };
std::vector<std::string> optional_keys = { "MUSICA name" };
for (const auto& key : required_keys)
const std::string SPECIES = "species";
const std::string MUSICA_NAME = "MUSICA name";
for (const auto& key : { SPECIES, MUSICA_NAME })
{
if (!ValidateJsonWithKey(object, key))
return ConfigParseStatus::RequiredKeyNotFound;
}

std::string name = object["species"].get<std::string>();
std::string species = object["species"].get<std::string>();
json reactants_object{};
json products_object{};
products_object[species] = { { "YIELD", 1.0 } };
auto reactants = ParseReactants(reactants_object);
auto products = ParseProducts(products_object);

std::string name = "EMIS." + object[MUSICA_NAME].get<std::string>();

emission_arr_.push_back(Species(name));
user_defined_rate_arr_.push_back(UserDefinedRateConstant(name));

std::unique_ptr<UserDefinedRateConstant> rate_ptr = std::make_unique<UserDefinedRateConstant>(name);
processes_.push_back(Process(reactants, products, std::move(rate_ptr), gas_phase_));

return ConfigParseStatus::Success;
}

ConfigParseStatus ParseFirstOrderLoss(const json& object)
{
std::vector<std::string> required_keys = { "species" };
std::vector<std::string> optional_keys = { "MUSICA name" };
for (const auto& key : required_keys)
const std::string SPECIES = "species";
const std::string MUSICA_NAME = "MUSICA name";
for (const auto& key : { SPECIES, MUSICA_NAME })
{
if (!ValidateJsonWithKey(object, key))
return ConfigParseStatus::RequiredKeyNotFound;
}

std::string name = object["species"].get<std::string>();
std::string species = object["species"].get<std::string>();
json reactants_object{};
json products_object{};
reactants_object[species] = { { } };
auto reactants = ParseReactants(reactants_object);
auto products = ParseProducts(products_object);

first_order_loss_arr_.push_back(Species(name));
std::string name = "LOSS." + object[MUSICA_NAME].get<std::string>();

user_defined_rate_arr_.push_back(UserDefinedRateConstant(name));

std::unique_ptr<UserDefinedRateConstant> rate_ptr = std::make_unique<UserDefinedRateConstant>(name);
processes_.push_back(Process(reactants, products, std::move(rate_ptr), gas_phase_));

return ConfigParseStatus::Success;
}
};
Expand Down Expand Up @@ -729,19 +749,6 @@ namespace micm
std::move(System(std::move(this->gas_phase_), std::move(this->phases_))), std::move(this->processes_));
}

/// @brief Get a collection of 'PhotolysisRateConstant'
/// @return a collection of 'PhotolysisRateConstant'
std::vector<PhotolysisRateConstant>& GetPhotolysisRateConstants()
{
if (last_parse_status_ != ConfigParseStatus::Success)
{
std::string msg = "Parsing configuration files failed. The parsing failed with error: " +
configParseStatusToString(last_parse_status_);
throw std::runtime_error(msg);
}

return this->photolysis_rate_arr_;
}
};

} // namespace micm
13 changes: 11 additions & 2 deletions include/micm/process/rate_constant.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -18,10 +18,19 @@ namespace micm
public:
/// @brief Virtual destructor
virtual ~RateConstant(){};

/// @brief Deep copy
virtual std::unique_ptr<RateConstant> clone() const = 0;
/// @brief Returns the number of doubles needed to hold user-defined rate constant parameters
/// @return Number of user-defined rate constant parameters

/// @brief Returns a set of labels for user-defined rate constant parameters
/// @return Vector of custom parameter labels
virtual std::vector<std::string> CustomParameters() const
{
return std::vector<std::string>{};
}

/// @brief Returns the number of custom parameters
/// @return Number of custom parameters
virtual std::size_t SizeCustomParameters() const
{
return 0;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -11,28 +11,27 @@ namespace micm
class System;

/// @brief A photolysis rate constant
class PhotolysisRateConstant : public RateConstant
class UserDefinedRateConstant : public RateConstant
{
public:
std::string name_;

public:
/// @brief Default constructor.
PhotolysisRateConstant();
UserDefinedRateConstant();

/// @brief
/// @param name A name for this reaction
PhotolysisRateConstant(const std::string& name);
UserDefinedRateConstant(const std::string& name);

/// @brief Deep copy
std::unique_ptr<RateConstant> clone() const override;

/// @brief Returns the number of parameters (1) that can be set at runtime
/// for photolysis reactions
///
/// The single editable parameter is the unscaled rate constant for
/// the photolysis reaction
/// @return Number of custom rate constant parameters
/// @brief Returns a label for the user-defined rate constant parameter
/// @return Rate constant label
std::vector<std::string> CustomParameters() const override;

/// @brief Returns the number of custom parameters
/// @return Number of custom parameters
std::size_t SizeCustomParameters() const override;

/// @brief Calculate the rate constant
Expand All @@ -43,29 +42,34 @@ namespace micm
const override;
};

inline PhotolysisRateConstant::PhotolysisRateConstant()
inline UserDefinedRateConstant::UserDefinedRateConstant()
: name_()
{
}

inline PhotolysisRateConstant::PhotolysisRateConstant(const std::string& name)
inline UserDefinedRateConstant::UserDefinedRateConstant(const std::string& name)
: name_(name)
{
}

inline std::unique_ptr<RateConstant> PhotolysisRateConstant::clone() const
inline std::unique_ptr<RateConstant> UserDefinedRateConstant::clone() const
{
return std::unique_ptr<RateConstant>{ new PhotolysisRateConstant{ *this } };
return std::unique_ptr<RateConstant>{ new UserDefinedRateConstant{ *this } };
}

inline double PhotolysisRateConstant::calculate(
inline double UserDefinedRateConstant::calculate(
const Conditions& conditions,
const std::vector<double>::const_iterator& custom_parameters) const
{
return (double)*custom_parameters;
}

inline std::size_t PhotolysisRateConstant::SizeCustomParameters() const
inline std::vector<std::string> UserDefinedRateConstant::CustomParameters() const
{
return std::vector<std::string>{ name_ };
}

inline std::size_t UserDefinedRateConstant::SizeCustomParameters() const
{
return 1;
}
Expand Down
9 changes: 4 additions & 5 deletions include/micm/solver/rosenbrock.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -695,15 +695,14 @@ namespace micm
template<template<class> class MatrixPolicy, template<class> class SparseMatrixPolicy>
inline State<MatrixPolicy> RosenbrockSolver<MatrixPolicy, SparseMatrixPolicy>::GetState() const
{
std::size_t n_params = 0;
std::vector<std::string> param_labels{};
for (const auto& process : processes_)
{
if (process.rate_constant_)
n_params += process.rate_constant_->SizeCustomParameters();
}
for (auto& label : process.rate_constant_->CustomParameters())
param_labels.push_back(label);
return State<MatrixPolicy>{ micm::StateParameters{ .state_variable_names_ = system_.UniqueNames(state_reordering_),
.custom_rate_parameter_labels_ = param_labels,
.number_of_grid_cells_ = parameters_.number_of_grid_cells_,
.number_of_custom_parameters_ = n_params,
.number_of_rate_constants_ = processes_.size() } };
}

Expand Down
Loading

0 comments on commit d20e289

Please sign in to comment.