Homec4science
Diffusion SpecMiCP / ReactMiCP (reactmicp/unsaturated) Browse

specmicp/src/database/species/species.hpp
fbbaf4c5659ereactmicp/unsaturated

Last Change
Raw File

species.hpp
View Options

/*-------------------------------------------------------------------------------
Copyright (c) 2014,2015 F. 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:
1. Redistributions of source code must retain the above copyright notice, this
list of conditions and the following disclaimer.
2. 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.
3. Neither the name of the copyright holder 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 HOLDER 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_DATABASE_SPECIES_HPP
#define SPECMICP_DATABASE_SPECIES_HPP
#include "../../types.hpp"
#include "base_wrapper.hpp"
//! \file species.hpp
//!
//! \brief The base wrappers describing a species list
//!
//! The two classes of this file are the basis for any list of species
//! and are at the cores of the database and its algorithm.
//!
namespace specmicp {
namespace database {
//! \class SpeciesList
//! \brief Base class for a species list
//!
//! Manage the labels, id and the size of a species list.
//!
//! \ingroup database_species
class SpeciesList
{
public:
//! \brief An empty list
SpeciesList():
m_is_valid(false)
{}
//! \brief Initalise a list with 'size' species
SpeciesList(index_t size):
m_is_valid(false),
m_size(size),
m_labels(size),
m_hashes(size)
{}
virtual ~SpeciesList() = default;
//! \name Size
//! \brief Manage the size of the list
// ---------------
//! @{
//! \brief Return the number of species in the list
index_t size() const {return m_size;}
//! \brief Reset the number of species in the list
virtual void resize(index_t size) {
m_size = size;
m_labels.resize(size);
m_hashes.resize(size);
set_nonvalid();
}
//! @}
//! \name Getter
//! \brief Return values
// -------
//! @{
//! \brief Return the label of species
const std::string& get_label(index_t ind) const {return m_labels[ind];}
//! \brief Return the index of a species from the label
//! \return the index of the species or -1 if the species does not exist
index_t get_id(const std::string& label) const;
//! \brief to iterate over the elements
range_t range() const {
return boost::irange((index_t) 0, m_size);
}
//! \brief Return true if the list is valid
//!
//! The validity flag is automatically turned off when rsize, move_erase or move_erase_to is use.
//! The user must turn it on again by using the set_valid member
//!
//! \sa set_valid
//! \sa set_nonvalid
bool is_valid() const noexcept {
return m_is_valid;
}
//! \brief Return a hash of the list
size_t get_hash() const;
//! @}
//! \name Setter
//! \brief Set values
// ------
//! @{
//! \brief Set the label of species 'ind'
void set_label(index_t ind, const std::string& label);
//! \brief Set the label of species 'ind'
void set_label(index_t ind, std::string&& label);
//! \brief Set the list to be nonvalid
//! \sa is_valid
//! \sa set_valid
void set_nonvalid() noexcept {m_is_valid = false;}
//! \brief Set the list to be valid
//! \sa is_valid
//! \sa set_nonvalid
void set_valid() {m_is_valid = true;}
//! \brief Erase values at 'ind'
void reset(index_t ind) {
m_labels[ind] = "";
m_hashes[ind] = 0;
}
//! @}
//! \name Move species
//! \brief Move species in the list or to another list
// -------------
//! @{
//! \brief move species 'old_ind' to 'new_ind' and erase the previous species at 'new_ind'
virtual void move_erase(index_t old_ind, index_t new_ind)
{
set_nonvalid();
m_labels[new_ind] = m_labels[old_ind];
m_hashes[new_ind] = m_hashes[old_ind];
}
//! \brief move 'ind' to 'other' at 'other_ind'
virtual void move_erase_to(index_t ind, SpeciesList& other, index_t other_ind)
{
set_nonvalid();
other.set_nonvalid();
other.set_label(other_ind, std::move(get_label(ind)));
reset(ind);
}
//! @}
private:
using hash_result_type = std::hash<std::string>::result_type; //!< Result of the hash function
bool m_is_valid;
index_t m_size; //!< Size of the list
std::vector<std::string> m_labels; //!< Labels
std::vector<hash_result_type> m_hashes; //!< List of hashes
};
//! \class ReactiveSpeciesList
//! \brief A list of reactive species
//!
//! This list extends 'SpeciesList' and add the managements of
//! stoichiometric coefficients and equilibrium constant.
//! \ingroup database_species
class ReactiveSpeciesList: public SpeciesList
{
public:
//! \brief Initialize an empty list
ReactiveSpeciesList() {}
//! \brief Initialize a list of size 'size' with 'nb_components' components
ReactiveSpeciesList(index_t size, index_t nb_components):
SpeciesList(size),
m_nu_coeffs(size, nb_components),
m_log_k(size)
{}
virtual ~ReactiveSpeciesList() = default;
//! \name Size
//! \brief Manage the size of the list
// -----------
//! @{
//! \brief Resize
void resize(index_t size) override
{
m_nu_coeffs.resize(size);
m_log_k.resize(size);
SpeciesList::resize(size);
}
//! \brief Resize the list and the number of components
virtual void resize(index_t size, index_t nb_components)
{
m_nu_coeffs.resize(size, nb_components);
m_log_k.resize(size);
SpeciesList::resize(size);
}
//! \brief Return the number of components
index_t nb_components() {
return m_nu_coeffs.cols();
}
//! @}
//! \name Getter
//! \brief Return values
// ------
//! @{
//! \brief Return a const reference to the logk vector
const Vector& get_logk_vector() const {
return m_log_k.get_vector();
}
//! \brief Return a const reference to the stoichiometric matrix
const Matrix& get_nu_matrix() const {
return m_nu_coeffs.get_matrix();
}
//! \brief Return a const reference to a row of the stoichiometric coefficient matrix
Eigen::MatrixBase<Matrix>::ConstRowXpr get_nu_row(index_t k) const {
return m_nu_coeffs.get_row(k);
}
//! \brief Return the logK of species 'j'
const scalar_t& logk(index_t j) const {
return m_log_k(j);
}
//! \brief Return the stoichiometric coefficient of species j and component i
const scalar_t& nu_ji(index_t j, index_t i) const {
return m_nu_coeffs(j, i);
}
//! @}
//! \name Setter
//! \brief Set the values of the parameters
// ----------
//! @{
//! \brief set a logk
void set_logk(index_t j, scalar_t value) {
m_log_k.set_value(j, value);
}
//! \brief Set the logk vector
template <typename derived>
void set_logk_vector(const Eigen::MatrixBase<derived>& vector) {
specmicp_assert(vector.rows() == m_log_k.size);
m_log_k.set_vector(vector);
}
//! \brief Set the stoehiometric coefficient of species 'j' and component 'i' to 'value'
void set_nu_ji(index_t j, index_t i, scalar_t value) {
m_nu_coeffs.set_value(j, i, value);
}
//! \brief Set the stoichiometric coefficient
template <typename derived>
void set_nu_matrix(const Eigen::MatrixBase<derived>& matrix) {
specmicp_assert(matrix.rows() == size());
m_nu_coeffs.set_matrix(matrix);
}
//! \brief Switch the basis
//!
//! Reference : Bethke2008
template <typename derived1, typename derived2>
void switch_basis(
const Eigen::MatrixBase<derived1>& transform,
const Eigen::MatrixBase<derived2>& kswap)
{
m_nu_coeffs.multiply_by(transform);
m_log_k.add(-get_nu_matrix()*kswap);
}
//! \brief Reset a row of the stoichiometric coefficient matrix
void reset_nu_row(index_t k) {
m_nu_coeffs.reset_row(k);
}
//! @}
//! \name Adding species
//! \brief 'Add' one species to the other
//! Replace an aqueous species by its decomposition into components.
//! This is typically use during the canonicalization of the database
// --------------
//! @{
//! \brief Add the species 'species_to_add' to 'k'
void add_species_to(
index_t k,
index_t species_to_add,
scalar_t coeff
)
{
add_nu_species_to(k, species_to_add, coeff);
add_logk_species_to(k, species_to_add, coeff);
}
//! \brief Add the species 'species_to_add' from 'other' to 'k'
void add_alien_species_to(
index_t k,
const ReactiveSpeciesList& other,
index_t species_to_add,
scalar_t coeff
)
{
add_nu_alien_species_to(k, other, species_to_add, coeff);
add_logk_alien_species_to(k, other, species_to_add, coeff);
}
//! @}
//! \name Move operations
//! \brief Move species in the list or to other list
// -----
//! @{
//! \brief move species 'old_ind' to 'new_ind' and erase the previous species at 'new_ind'
void move_erase(index_t old_ind, index_t new_ind) override
{
m_nu_coeffs.move_erase(old_ind, new_ind);
m_log_k.move_erase(old_ind, new_ind);
SpeciesList::move_erase(old_ind, new_ind);
}
//! \brief move species 'old_ind' to 'new_ind' and erase the previous species at 'new_ind'
//!
//! This is the method that need to be used to remove components in the database.
virtual void move_erase(index_t old_ind, index_t new_ind, const std::vector<bool>& is_reactants_to_remove)
{
m_nu_coeffs.move_erase(old_ind, new_ind, is_reactants_to_remove);
m_log_k.move_erase(old_ind, new_ind);
SpeciesList::move_erase(old_ind, new_ind);
}
//! \brief Move species 'ind' to the list 'other' at 'other_ind'
void move_erase_to(index_t ind, ReactiveSpeciesList& other, index_t other_ind);
//! @}
private: // methods
//! \brief Add the logk of 'species_to_add' in 'other' at species 'k'
void add_logk_alien_species_to(
index_t k,
const ReactiveSpeciesList& other,
index_t species_to_add,
scalar_t coeff
)
{
auto tmp = m_log_k(k) + coeff*other.m_log_k(species_to_add);
m_log_k.set_value(k, tmp);
}
//! \brief Add the logk of species_to_add at species 'k'
void add_logk_species_to(
index_t k,
index_t species_to_add,
scalar_t coeff
)
{
auto tmp = m_log_k(k) + coeff*m_log_k(species_to_add);
m_log_k.set_value(k, tmp);
}
//! \brief Add the stoichiometric coefficient of species_to_add at species 'k'
void add_nu_alien_species_to(
index_t k,
const ReactiveSpeciesList& other,
index_t species_to_add,
scalar_t coeff
)
{
m_nu_coeffs.add_to_row(k, coeff*other.m_nu_coeffs.get_row(species_to_add));
}
//! \brief Add the stoichiometric coefficient of species_to_add at species 'k'
void add_nu_species_to(
index_t k,
index_t species_to_add,
scalar_t coeff
)
{
m_nu_coeffs.add_to_row(k, coeff*m_nu_coeffs.get_row(species_to_add));
}
private: // attributes
MatrixSpeciesWrapper m_nu_coeffs; //!< The stoichiometric coefficients
VectorSpeciesWrapper m_log_k; //!< The logk values
};
} // end namespace database
} // end namespace specmicp
namespace std {
template <>
struct hash<specmicp::database::SpeciesList> {
size_t operator ()(const specmicp::database::SpeciesList& alist) {return alist.get_hash();}
};
} //end namespace std
#endif // SPECMICP_DATABASE_SPECIES_HPP
c4science · Help