Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

162 euler option #197

Closed
wants to merge 20 commits into from
Closed

162 euler option #197

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
20 commits
Select commit Hold shift + click to select a range
a0da3c9
Added MICM::CreateBackwardEuler and MICM::CreateBackwardEulerStandard…
dwfncar Jul 17, 2024
8cf0a04
Added BackwardEuler and BackwardEulerStandard solver types.
dwfncar Jul 17, 2024
8cf94f4
Added MICMSolver types BackwardEuler and BackwardEulerStandardOrder.
dwfncar Jul 17, 2024
f9239bb
Added MICM::CreateBackwardEuler.
dwfncar Jul 17, 2024
d99b8df
Added MICM::CreateBackwardEuler.
dwfncar Jul 17, 2024
d5b8206
Merge branch 'main' into 162-euler-option
dwfncar Aug 5, 2024
8be7180
Start on test for BackwardEulerStandardOrder.
dwfncar Aug 5, 2024
993d7af
Test set up for BackwardEulerStandardOrder.
dwfncar Aug 5, 2024
a5765f4
Start on CreateBackwardEulerStandardOrder.
dwfncar Aug 7, 2024
05ca29b
Merge branch 'main' into 162-euler-option
dwfncar Aug 10, 2024
2e96ee9
Merged main.
dwfncar Aug 10, 2024
f7de2ef
Fixed micm::CreateBackwardEuler.
dwfncar Aug 10, 2024
715ccc2
Added MICM::CreateBackwardEulerStandardOrder.
dwfncar Aug 10, 2024
0cd931d
Start on unit test SolveUsingVectorOrderedBackwardEuler.
dwfncar Aug 10, 2024
4ac136e
Updated micm.cpp for Backward Euler.
dwfncar Aug 10, 2024
6e89178
Added unit test check on number of user defined rates for Backward Eu…
dwfncar Aug 10, 2024
0a18361
Set custom rate parameters in Backward Euler C API unit tests.
dwfncar Aug 10, 2024
6ac37a4
Added Backward Euler options to MicmSolve.
dwfncar Aug 11, 2024
d4d4ed1
Added MicmSolve to Backward Euler unit tests.
dwfncar Aug 11, 2024
c82938b
Added MicmSolve to Backward Euler unit tests.
dwfncar Aug 11, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 29 additions & 5 deletions include/musica/micm.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -39,8 +39,10 @@ namespace musica
/// @brief Types of MICM solver
enum MICMSolver
{
Rosenbrock = 1, // Vector-ordered Rosenbrock solver
RosenbrockStandardOrder, // Standard-ordered Rosenbrock solver
Rosenbrock = 1, // Vector-ordered Rosenbrock solver
BackwardEuler, // Vector-ordered BackwardEuler solver
RosenbrockStandardOrder, // Standard-ordered Rosenbrock solver
BackwardEulerStandardOrder, // Standard-ordered BackwardEuler solver
};
Copy link
Collaborator

Choose a reason for hiding this comment

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

The change of the enum order will fail the fortran micm APIs. I guess the ordering based on the matrix type is fine as it is because we only have two types of solver, but if we think about adding more solvers, it'd probably be easier for maintainers to order them based on the solver type, then the enum value doesn't have to change every time we add a new solver. Also for the users, if a new release has to change the ordering, it could be confusing.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I agree it will better to order by solver type / matrix type ... will also rearrange source with this ordering.


struct SolverResultStats
Expand Down Expand Up @@ -140,11 +142,21 @@ namespace musica
/// @param error Error struct to indicate success or failure
void CreateRosenbrock(const std::string &config_path, Error *error);

/// @brief Create a BackwardEuler solver of vector-ordered matrix type by reading and parsing configuration file
Copy link
Collaborator

Choose a reason for hiding this comment

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

I have the same thoughts on the ordering of the functions related to the above comments.

/// @param config_path Path to configuration file or directory containing configuration file
/// @param error Error struct to indicate success or failure
void CreateBackwardEuler(const std::string &config_path, Error *error);

/// @brief Create a Rosenbrock solver of standard-ordered matrix type by reading and parsing configuration file
/// @param config_path Path to configuration file or directory containing configuration file
/// @param error Error struct to indicate success or failure
void CreateRosenbrockStandardOrder(const std::string &config_path, Error *error);

/// @brief Create a BackwardEuler solver of standard-ordered matrix type by reading and parsing configuration file
/// @param config_path Path to configuration file or directory containing configuration file
/// @param error Error struct to indicate success or failure
void CreateBackwardEulerStandardOrder(const std::string &config_path, Error *error);

/// @brief Solve the system
/// @param solver Pointer to solver
/// @param time_step Time [s] to advance the state by
Expand Down Expand Up @@ -210,24 +222,36 @@ namespace musica
public:
MICMSolver solver_type_;

/// @brief Vector-ordered Rosenbrock solver type
/// @brief Vector-ordered Rosenbrock and BackwardEuler solver types
using DenseMatrixVector = micm::VectorMatrix<double, MICM_VECTOR_MATRIX_SIZE>;
using SparseMatrixVector = micm::SparseMatrix<double, micm::SparseMatrixVectorOrdering<MICM_VECTOR_MATRIX_SIZE>>;

using RosenbrockVectorType = typename micm::RosenbrockSolverParameters::
template SolverType<micm::ProcessSet, micm::LinearSolver<SparseMatrixVector, micm::LuDecomposition>>;
using Rosenbrock = micm::Solver<RosenbrockVectorType, micm::State<DenseMatrixVector, SparseMatrixVector>>;
using VectorState = micm::State<DenseMatrixVector, SparseMatrixVector>;
std::unique_ptr<Rosenbrock> rosenbrock_;

/// @brief Standard-ordered Rosenbrock solver type
using BackwardEulerVectorType = typename micm::BackwardEulerSolverParameters::
template SolverType<micm::ProcessSet, micm::LinearSolver<SparseMatrixVector, micm::LuDecomposition>>;
using BackwardEuler = micm::Solver<BackwardEulerVectorType, micm::State<DenseMatrixVector, SparseMatrixVector>>;
std::unique_ptr<BackwardEuler> backward_euler_;

/// @brief Standard-ordered Rosenbrock and BackwardEuler solver types
using DenseMatrixStandard = micm::Matrix<double>;
using SparseMatrixStandard = micm::SparseMatrix<double, micm::SparseMatrixStandardOrdering>;

using RosenbrockStandardType = typename micm::RosenbrockSolverParameters::
template SolverType<micm::ProcessSet, micm::LinearSolver<SparseMatrixStandard, micm::LuDecomposition>>;
using RosenbrockStandard = micm::Solver<RosenbrockStandardType, micm::State<DenseMatrixStandard, SparseMatrixStandard>>;
using StandardState = micm::State<DenseMatrixStandard, SparseMatrixStandard>;
std::unique_ptr<RosenbrockStandard> rosenbrock_standard_;

using BackwardEulerStandardType = typename micm::BackwardEulerSolverParameters::
template SolverType<micm::ProcessSet, micm::LinearSolver<SparseMatrixStandard, micm::LuDecomposition>>;
using BackwardEulerStandard = micm::Solver<BackwardEulerStandardType, micm::State<DenseMatrixStandard, SparseMatrixStandard>>;
std::unique_ptr<BackwardEulerStandard> backward_euler_standard_;

private:
int num_grid_cells_;
std::unique_ptr<micm::SolverParameters> solver_parameters_;
Expand Down Expand Up @@ -294,4 +318,4 @@ namespace musica
*error = ToError(MUSICA_ERROR_CATEGORY, MUSICA_ERROR_CODE_SPECIES_NOT_FOUND, msg.c_str());
return T();
}
} // namespace musica
} // namespace musica
120 changes: 120 additions & 0 deletions src/micm/micm.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@
#include <musica/util.hpp>

#include <micm/solver/rosenbrock_solver_parameters.hpp>
#include <micm/solver/backward_euler_solver_parameters.hpp>
#include <micm/solver/solver_builder.hpp>
#include <micm/system/species.hpp>
#include <micm/version.hpp>
Expand All @@ -32,11 +33,21 @@ namespace musica
micm->SetSolverType(MICMSolver::Rosenbrock);
micm->CreateRosenbrock(std::string(config_path), error);
}
else if (solver_type == MICMSolver::BackwardEuler)
{
micm->SetSolverType(MICMSolver::BackwardEuler);
micm->CreateBackwardEuler(std::string(config_path), error);
}
else if (solver_type == MICMSolver::RosenbrockStandardOrder)
{
micm->SetSolverType(MICMSolver::RosenbrockStandardOrder);
micm->CreateRosenbrockStandardOrder(std::string(config_path), error);
}
else if (solver_type == MICMSolver::BackwardEulerStandardOrder)
{
micm->SetSolverType(MICMSolver::BackwardEulerStandardOrder);
micm->CreateBackwardEulerStandardOrder(std::string(config_path), error);
}
else
{
std::string msg = "Solver type '" + std::to_string(solver_type) + "' not found";
Expand Down Expand Up @@ -104,6 +115,22 @@ namespace musica
solver_stats,
error);
}
else if (micm->solver_type_ == MICMSolver::BackwardEuler)
{
micm->Solve(
micm->backward_euler_,
time_step,
temperature,
pressure,
air_density,
num_concentrations,
concentrations,
num_custom_rate_parameters,
custom_rate_parameters,
solver_state,
solver_stats,
error);
}
else if (micm->solver_type_ == MICMSolver::RosenbrockStandardOrder)
{
micm->Solve(
Expand All @@ -120,6 +147,22 @@ namespace musica
solver_stats,
error);
}
else if (micm->solver_type_ == MICMSolver::BackwardEulerStandardOrder)
{
micm->Solve(
micm->backward_euler_standard_,
time_step,
temperature,
pressure,
air_density,
num_concentrations,
concentrations,
num_custom_rate_parameters,
custom_rate_parameters,
solver_state,
solver_stats,
error);
}
};

String MicmVersion()
Expand All @@ -137,10 +180,18 @@ namespace musica
{
map = micm->GetSpeciesOrdering(micm->rosenbrock_, error);
}
else if (micm->solver_type_ == MICMSolver::BackwardEuler)
{
map = micm->GetSpeciesOrdering(micm->backward_euler_, error);
}
else if (micm->solver_type_ == MICMSolver::RosenbrockStandardOrder)
{
map = micm->GetSpeciesOrdering(micm->rosenbrock_standard_, error);
}
else if (micm->solver_type_ == MICMSolver::BackwardEulerStandardOrder)
{
map = micm->GetSpeciesOrdering(micm->backward_euler_standard_, error);
}
if (!IsSuccess(*error))
{
return nullptr;
Expand Down Expand Up @@ -171,10 +222,18 @@ namespace musica
{
map = micm->GetUserDefinedReactionRatesOrdering(micm->rosenbrock_, error);
}
if (micm->solver_type_ == MICMSolver::BackwardEuler)
{
map = micm->GetUserDefinedReactionRatesOrdering(micm->backward_euler_, error);
}
else if (micm->solver_type_ == MICMSolver::RosenbrockStandardOrder)
{
map = micm->GetUserDefinedReactionRatesOrdering(micm->rosenbrock_standard_, error);
}
else if (micm->solver_type_ == MICMSolver::BackwardEulerStandardOrder)
{
map = micm->GetUserDefinedReactionRatesOrdering(micm->backward_euler_standard_, error);
}
if (!IsSuccess(*error))
{
return nullptr;
Expand Down Expand Up @@ -266,6 +325,40 @@ namespace musica
}
}

void MICM::CreateBackwardEuler(const std::string &config_path, Error *error)
{
try
{
micm::SolverConfig<> solver_config;
solver_config.ReadAndParse(std::filesystem::path(config_path));
solver_parameters_ = std::make_unique<micm::SolverParameters>(solver_config.GetSolverParams());

backward_euler_ = std::make_unique<BackwardEuler>(
micm::SolverBuilder<
micm::BackwardEulerSolverParameters,
micm::VectorMatrix<double, MICM_VECTOR_MATRIX_SIZE>,
micm::SparseMatrix<double, micm::SparseMatrixVectorOrdering<MICM_VECTOR_MATRIX_SIZE>>,
micm::ProcessSet,
micm::LinearSolver<
micm::SparseMatrix<double, micm::SparseMatrixVectorOrdering<MICM_VECTOR_MATRIX_SIZE>>,
micm::LuDecomposition>,
VectorState>(micm::BackwardEulerSolverParameters())
.SetSystem(solver_parameters_->system_)
.SetReactions(solver_parameters_->processes_)
.SetNumberOfGridCells(num_grid_cells_)
.SetIgnoreUnusedSpecies(true)
.Build());

DeleteError(error);
*error = NoError();
}
catch (const std::system_error &e)
{
DeleteError(error);
*error = ToError(e);
}
}

void MICM::CreateRosenbrockStandardOrder(const std::string &config_path, Error *error)
{
try
Expand Down Expand Up @@ -293,6 +386,33 @@ namespace musica
}
}

void MICM::CreateBackwardEulerStandardOrder(const std::string &config_path, Error *error)
{
try
{
micm::SolverConfig<> solver_config;
solver_config.ReadAndParse(std::filesystem::path(config_path));
solver_parameters_ = std::make_unique<micm::SolverParameters>(solver_config.GetSolverParams());

backward_euler_standard_ =
std::make_unique<BackwardEulerStandard>(micm::CpuSolverBuilder<micm::BackwardEulerSolverParameters>(
micm::BackwardEulerSolverParameters())
.SetSystem(solver_parameters_->system_)
.SetReactions(solver_parameters_->processes_)
.SetNumberOfGridCells(num_grid_cells_)
.SetIgnoreUnusedSpecies(true)
.Build());

DeleteError(error);
*error = NoError();
}
catch (const std::system_error &e)
{
DeleteError(error);
*error = ToError(e);
}
}

void MICM::Solve(
auto &solver,
double time_step,
Expand Down
Loading
Loading