Homec4science
Diffusion SpecMiCP / ReactMiCP (database) Browse

specmicp/src/specmicp/equilibrium_data.hpp
e4061f7461eddatabase

Last Change
Raw File

equilibrium_data.hpp
View Options

/*-------------------------------------------------------
- Module : specmicp
- File : equilibrium_data.hpp
- Author : Fabien Georget
Copyright (c) 2014, Fabien Georget, Princeton University
---------------------------------------------------------*/
#ifndef SPECMICP_SPECMICP_EQUILIBRIUMDATA_HPP
#define SPECMICP_SPECMICP_EQUILIBRIUMDATA_HPP
//! \file equilibrium_data.hpp store results from equilibrium computation
#include "database/database.hpp"
#define POW10(m) std::pow(10.0, m)
namespace specmicp {
class EquilibriumState
{
EIGEN_MAKE_ALIGNED_OPERATOR_NEW
public:
EquilibriumState():
m_has_been_properly_initialised(false)
{}
EquilibriumState(const Eigen::VectorXd& main_vars,
const Eigen::VectorXd& secondary,
const Eigen::VectorXd& loggamma,
double ionic_strength,
std::shared_ptr<database::DataContainer> data):
m_has_been_properly_initialised(true),
m_main_variables(main_vars),
m_concentration_secondary(secondary),
m_ionic_strength(ionic_strength),
m_loggamma(loggamma),
m_data(data)
{}
//! Return true if it contains actual data
bool is_valid() {return m_has_been_properly_initialised;}
// Main variables
// ==============
//! Return the mass of free water in the system
double mass_water() const {return m_main_variables(0);}
//! Return the molality of component i
double molality_component(int i) const {return POW10(m_main_variables(i));}
//! Return the number of moles of the mineral at equilibrium 'm'
double moles_mineral(int m) const {return m_main_variables(m_data->nb_component+m);}
//! Return the mass of a mineral (kg)
double mass_mineral(int m) const;
// molality
// ========
// return the pH
double pH() const;
//! Return the molality of secondary species j
double molality_secondary(int j) const {return m_concentration_secondary(j);}
//! Return a const reference to the vector of molality of the secondary array
const Eigen::VectorXd& molality_secondary() const {return m_concentration_secondary;}
//! Return total concentration of aqueous species for component i
double total_aqueous_concentration_component(int i) const;
//! Compute total aqeuous concentration for each component
void total_aqueous_concentrations(Eigen::VectorXd& totconc) const;
// Activity
// ========
//! Return a const reference to the activity coefficient vector
const Eigen::VectorXd& activity_coefficient() const {return m_loggamma;}
//! Return activity coefficient of component i
double activity_coefficient_component(int i) const {
assert(i > 0);
return POW10(m_loggamma(i));
}
//! Return activity coefficient of secondary species j
double activity_coefficient_secondary(int j) const {return POW10(m_loggamma(j+m_data->nb_component));}
//! Return activity of component i (i > 0)
double activity_component(int i) const {
assert(i > 0);
return POW10(m_loggamma(i) + m_main_variables(i));
}
//! Return activity of secondary species j
double activity_secondary(int j) const {
return POW10(m_loggamma(j+m_data->nb_component))*m_concentration_secondary(j);
}
//! return log(IAP) for mineral m
double logIAP(int m) const;
//! return log(IAP) for mineral "label mineral"
double logIAP(const std::string& label_mineral) const;
//! return saturation index for mineral m
double saturation_index(int m) const {return logIAP(m) - m_data->logk_mineral(m);}
//! return saturation index for mineral "label mineral"
double saturation_index(const std::string& label_mineral) const;
//! return log(IAP) for mineral m
double logIAP_kinetic(int m) const;
//! return log(IAP) for mineral "label_mineral"
double logIAP_kinetic(const std::string& label_mineral) const;
//! return saturation index for mineral m
double saturation_index_kinetic(int m) const {return logIAP(m) - m_data->logk_mineral_kinetic(m);}
//! return saturation index for mineral m
double saturation_index_kinetic(const std::string& label_mineral) const;
private:
bool m_has_been_properly_initialised;
Eigen::VectorXd m_main_variables;
Eigen::VectorXd m_concentration_secondary;
double m_ionic_strength;
Eigen::VectorXd m_loggamma;
std::shared_ptr<database::DataContainer> m_data;
};
} // end namespace specmicp
#endif // SPECMICP_SPECMICP_EQUILIBRIUMDATA_HPP
c4science · Help