Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F87405359
utils.cpp
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Sat, Oct 12, 12:19
Size
4 KB
Mime Type
text/x-c
Expires
Mon, Oct 14, 12:19 (2 d)
Engine
blob
Format
Raw Data
Handle
21545257
Attached To
rSPECMICP SpecMiCP / ReactMiCP
utils.cpp
View Options
#include "utils.hpp"
EquilibriumState blank_composition(std::shared_ptr<specmicp::database::DataContainer> database)
{
std::shared_ptr<specmicp::ReactionPathModel> model = std::make_shared<specmicp::ReactionPathModel>();
double m_c3s = 1e-7;
double m_c2s = 1e-7;
double m_water = 1;
double c_hco3 = 1e-4;
model->amount_aqueous = {
{"H2O", specmicp::reaction_amount_t(m_water/database->molar_mass_basis_si(0)+c_hco3, 0)},
{"HCO3[-]", specmicp::reaction_amount_t(c_hco3, 0)},
{"HO[-]", specmicp::reaction_amount_t(-c_hco3, 0)},
};
model->amount_minerals = {
{"C3S", specmicp::reaction_amount_t(m_c3s, 0.0)},
{"C2S", specmicp::reaction_amount_t(m_c2s, 0.0)}
};
model->database_path = "data/cemdata_specmicp.js";
Eigen::VectorXd x;
specmicp::ReactionPathDriver driver(model, database, x);
//driver.get_options().solver_options.penalization_factor =1.0;
driver.get_options().ncp_function = specmicp::micpsolver::NCPfunction::penalizedFB;
//driver.get_options().ncp_function = specmicp::micpsolver::NCPfunction::min;
driver.get_options().solver_options.use_scaling = false;
driver.get_options().solver_options.max_factorization_step = 1;
driver.get_options().solver_options.factor_gradient_search_direction = 50;
driver.get_options().solver_options.maxstep = 50;
driver.get_options().solver_options.maxiter_maxstep = 50;
driver.get_options().solver_options.max_iter = 100;
driver.get_options().allow_restart = true;
specmicp::micpsolver::MiCPPerformance perf = driver.one_step(x);
if (perf.return_code != specmicp::micpsolver::MiCPSolverReturnCode::ResidualMinimized)
throw std::runtime_error("Error - system cannot be solved");
return driver.get_current_solution();
}
EquilibriumState sample_carbo_composition(std::shared_ptr<specmicp::database::DataContainer> database, double mult)
{
std::shared_ptr<specmicp::ReactionPathModel> model = std::make_shared<specmicp::ReactionPathModel>();
const double m_c3s = mult*0.7;
const double m_c2s = mult*0.3;
const double wc = 1.0;
const double m_water = wc*((3*56.08+60.08)*m_c3s + (2*56.08+60.08)*m_c2s)*1e-3;
const double delta_h2co3 = 0.005;
model->amount_aqueous = {
{"H2O", specmicp::reaction_amount_t(m_water/database->molar_mass_basis_si(0), delta_h2co3)},
{"HCO3[-]", specmicp::reaction_amount_t(delta_h2co3, 0.0)},
{"HO[-]", specmicp::reaction_amount_t(-delta_h2co3, 0.0)},
};
model->amount_minerals = {
{"C3S", specmicp::reaction_amount_t(m_c3s, 0.0)},
{"C2S", specmicp::reaction_amount_t(m_c2s, 0.0)}
};
model->database_path = "data/cemdata_specmicp.js";
Eigen::VectorXd x;
specmicp::ReactionPathDriver driver(model, database, x);
//driver.get_options().solver_options.penalization_factor =1.0;
driver.get_options().ncp_function = specmicp::micpsolver::NCPfunction::penalizedFB;
//driver.get_options().ncp_function = specmicp::micpsolver::NCPfunction::min;
driver.get_options().solver_options.use_scaling = false;
driver.get_options().solver_options.max_factorization_step = 1;
driver.get_options().solver_options.factor_gradient_search_direction = 50;
driver.get_options().solver_options.maxstep = 50;
driver.get_options().solver_options.maxiter_maxstep = 50;
driver.get_options().solver_options.max_iter = 100;
driver.get_options().allow_restart = true;
specmicp::micpsolver::MiCPPerformance perf = driver.one_step(x);
if (perf.return_code != specmicp::micpsolver::MiCPSolverReturnCode::ResidualMinimized)
throw std::runtime_error("Error - system cannot be solved");
return driver.get_current_solution();
}
//! \brief Return a simple fully processed database for testing purposes
//!
//! Correspond to the carbonation problem
std::shared_ptr<specmicp::database::DataContainer> get_test_carbo_database()
{
specmicp::database::Database database("data/cemdata_specmicp.js");
std::map<std::string, std::string> swapping ({
{"H[+]","HO[-]"},
{"Si(OH)4", "SiO(OH)3[-]"}
});
database.swap_components(swapping);
std::vector<std::string> to_keep = {"H2O", "HO[-]", "Ca[2+]", "SiO(OH)3[-]", "HCO3[-]"};
database.keep_only_components(to_keep);
return database.get_database();
}
Event Timeline
Log In to Comment