Homec4science
Diffusion SpecMiCP / ReactMiCP (reactivetransport) Browse

specmicp/src/specmicp/reaction_path.hpp
b9a808db21fbreactivetransport

Last Change
Raw File

reaction_path.hpp
View Options

/*-------------------------------------------------------
- Module : specmicp
- File : reaction_path.hpp
- Author : Fabien Georget
Copyright (c) 2014, Fabien Georget <fabieng@princeton.edu>, Princeton University
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright
notice, this list of conditions and the following disclaimer in the
documentation and/or other materials provided with the distribution.
* Neither the name of the Princeton University nor the
names of its contributors may be used to endorse or promote products
derived from this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
---------------------------------------------------------*/
#ifndef SPECMICP_SPECMICP_REACTIONPATH_HPP
#define SPECMICP_SPECMICP_REACTIONPATH_HPP
//! \file reaction_path.hpp reaction path driver
#include <memory>
#include <map>
#include "database.hpp"
#include "micpsolver/micpsolver_structs.hpp"
#include "reduced_system_solver.hpp"
#include "equilibrium_data.hpp"
namespace specmicp {
//! \brief typedef to a reaction path amount
//!
//! The first element is the inital amount,
//! The second element is the increase to add at every step
using reaction_amount_t = std::pair<double, double>;
//! \brief A reaction path model
//!
//! Contains every data we need to run a reaction path model
struct ReactionPathModel
{
//! The amount of components
std::map<std::string, reaction_amount_t> amount_aqueous;
//! The amount of minerals
std::map<std::string, reaction_amount_t> amount_minerals;
//! List of minerals to take into account
std::vector<std::string> minerals_to_keep;
//! the number of step
int nb_step;
//! the path to the database - optional
std::string database_path;
};
//! \brief Driver to solve a reaction path model
class ReactionPathDriver
{
public:
//! \brief Initialise from a model
//!
//! \param model struct containing the simulation data
//! \param options the options of the equilibrium solver and the MiCP solver
//!
//! Build a database from the attribute database_path of the model
ReactionPathDriver(std::shared_ptr<ReactionPathModel> model,
ReducedSystemSolverOptions options=ReducedSystemSolverOptions()):
m_model(model),
m_solver_options(options),
m_must_init(false)
{
read_database();
}
//! \brief Initialise the driver from a model and a raw database
//!
//! \param model struct containing the simulation data
//! \param data shared_ptr to the raw database
//! \param options the options of the equilibrium solver and the MiCP solver
//!
ReactionPathDriver(std::shared_ptr<ReactionPathModel> model,
RawDatabasePtr data,
ReducedSystemSolverOptions options=ReducedSystemSolverOptions()):
m_model(model),
m_data(data),
m_solver_options(options),
m_must_init(false)
{}
//! \brief Initialise the driver from a model and a raw database
//!
//! \param model struct containing the simulation data
//! \param data shared_ptr to the raw database
//! \param x an Eigen vector, will be initialized by the solver
//! \param options the options of the equilibrium solver and the MiCP solver
//!
ReactionPathDriver(std::shared_ptr<ReactionPathModel> model,
RawDatabasePtr data,
Eigen::VectorXd& x,
ReducedSystemSolverOptions options=ReducedSystemSolverOptions()):
m_model(model),
m_data(data),
m_solver_options(options),
m_must_init(false)
{
init_solution(x);
}
//! Initialize the database
void read_database();
//! Dissolve everything into components
void dissolve_to_components();
//! \brief Perform one step
//!
//! \param[in,out] x the solution vector
//! \param init if true, the solution vector will be initialized
micpsolver::MiCPPerformance one_step(Vector& x, bool init=false);
//! \brief run all steps
void run_all();
//! \brief Return the current solution as an EquilibriumState object
EquilibriumState get_current_solution() const {return m_current_solution;}
//! \brief initialize the solution vector
//!
//! \param x an Eigen vector
//!
//! Resize the solution vector, prepare the total concentration vector and
//! set the flag to initialize the solution vector with a reasonable guess
//! on the next call of the solver
void init_solution(Vector& x);
//! Return a read-only reference to the options
const ReducedSystemSolverOptions& get_options() const {return m_solver_options;}
//! \brief Read-Write access to the options
//!
//! Note the options also contain the options of the MiCP solver
ReducedSystemSolverOptions& get_options() {return m_solver_options;}
ReducedSystemSolverOptions* get_ptr_options() {return &m_solver_options;}
private:
std::shared_ptr<ReactionPathModel> m_model; // the simulation to solve
RawDatabasePtr m_data; // raw container to the database
Eigen::VectorXd m_tot_conc; // The initial total concentration
Eigen::VectorXd m_tot_conc_increase; // the increase to add at each step
ReducedSystemSolverOptions m_solver_options; // the options
EquilibriumState m_current_solution; // the current solution for
// - temporary storage
// - check
// - return EquilibriumState
bool m_must_init; // true if the solution must be initialized at the next run of the solver
};
} // end namespace specmicp
#endif // SPECMICP_SPECMICP_REACTIONPATH_HPP
c4science · Help