Homec4science
Diffusion SpecMiCP / ReactMiCP (yaml/safe) Browse

specmicp/src/database/selector.cpp
392b50f45f16yaml/safe

Last Change
Raw File

selector.cpp
View Options

/*-------------------------------------------------------
- Module : database
- File : selector.cpp
- Author : Fabien Georget
Copyright (c) 2014, Fabien Georget, Princeton University
---------------------------------------------------------*/
#include "selector.hpp"
namespace specmicp {
namespace database {
void DatabaseSelector::remove_component(const std::vector<int>& id_component_to_remove)
{
analyse_component(id_component_to_remove);
select_aqueous(id_component_to_remove);
select_minerals(id_component_to_remove);
select_minerals_kinetic(id_component_to_remove);
if (data->nb_gas != 0) select_gas(id_component_to_remove);
select_components();
}
void DatabaseSelector::analyse_component(const std::vector<int> &id_component_to_remove)
{
for (auto it=id_component_to_remove.begin(); it!=id_component_to_remove.end(); ++it)
{
m_is_component_to_remove[*it] = true;
}
m_nb_component_to_keep = data->nb_component - id_component_to_remove.size();
int new_i = 0;
for (int i=0; i<data->nb_component; ++i)
{
if (is_component_to_remove(i)) continue;
data->param_aq.row(new_i) = data->param_aq.row(i);
data->molar_mass_basis(new_i) = data->molar_mass_basis(i);
++new_i;
}
assert(new_i == nb_component_to_keep());
data->molar_mass_basis.conservativeResize(nb_component_to_keep());
}
void DatabaseSelector::select_aqueous(const std::vector<int>& id_component_to_remove)
{
// first we select which aqueous species we should remove
std::vector<bool> aqueous_to_remove(data->nb_aqueous, false);
for (int j=0; j<data->nb_aqueous; ++j)
{
for (auto it=id_component_to_remove.begin(); it!=id_component_to_remove.end(); ++it)
{
if (data->nu_aqueous(j, *it) != 0)
{
aqueous_to_remove[j] = true;
break;
}
}
}
// then we remove them, in two steps
const int nb_aq_to_keep = std::count(aqueous_to_remove.begin(), aqueous_to_remove.end(), false);
int new_j = 0;
// 1) first we copy data in the beginning of the arrays
for (int j=0; j<data->nb_aqueous; ++j)
{
if (aqueous_to_remove[j]) continue;
data->logk_aqueous(new_j) = data->logk_aqueous(j);
data->param_aq.row(nb_component_to_keep()+new_j) = data->param_aq.row(data->nb_component+j);
data->labels_aqueous[new_j] = data->labels_aqueous[j];
int new_i = 0;
for (int i=0; i<data->nb_component; ++i)
{
if (is_component_to_remove(i)) continue;
data->nu_aqueous(new_j, new_i) = data->nu_aqueous(j, i);
++new_i;
}
assert(new_i == nb_component_to_keep());
++new_j;
}
assert(new_j == nb_aq_to_keep);
// 2) then we resize the arrays
data->nu_aqueous.conservativeResize(nb_aq_to_keep, nb_component_to_keep());
data->logk_aqueous.conservativeResize(nb_aq_to_keep);
data->labels_aqueous.resize(nb_aq_to_keep);
data->param_aq.conservativeResize(nb_component_to_keep()+nb_aq_to_keep, Eigen::NoChange);
data->nb_aqueous = nb_aq_to_keep;
// By doing in this order, we avoid bulk copy of arrays
}
void DatabaseSelector::select_minerals(const std::vector<int>& id_component_to_remove)
{
// first we select which minerals we should remove
std::vector<bool> minerals_to_remove(data->nb_mineral, false);
for (int m=0; m<data->nb_mineral; ++m)
{
for (auto it=id_component_to_remove.begin(); it!=id_component_to_remove.end(); ++it)
{
if (data->nu_mineral(m, *it) != 0)
{
minerals_to_remove[m] = true;
break;
}
}
}
// then we remove them, in two steps
const int nb_min_to_keep = std::count(minerals_to_remove.begin(), minerals_to_remove.end(), false);
int new_m = 0;
// 1) first we copy data in the beginning of the arrays
for (int m=0; m<data->nb_mineral; ++m)
{
if (minerals_to_remove[m]) continue;
data->logk_mineral(new_m) = data->logk_mineral(m);
data->labels_minerals[new_m] = data->labels_minerals[m];
data->_molar_volume_mineral(new_m) = data->_molar_volume_mineral(m);
int new_i = 0;
for (int i=0; i<data->nb_component; ++i)
{
if (is_component_to_remove(i)) continue;
data->nu_mineral(new_m, new_i) = data->nu_mineral(m, i);
++new_i;
}
assert(new_i == nb_component_to_keep());
++new_m;
}
assert(new_m == nb_min_to_keep);
// 2) then we resize the arrays
data->nu_mineral.conservativeResize(nb_min_to_keep, nb_component_to_keep());
data->logk_mineral.conservativeResize(nb_min_to_keep);
data->labels_minerals.resize(nb_min_to_keep);
data->_molar_volume_mineral.conservativeResize(nb_min_to_keep);
data->nb_mineral = nb_min_to_keep;
// By doing in this order, we avoid bulk copy of arrays
}
void DatabaseSelector::select_minerals_kinetic(const std::vector<int>& id_component_to_remove)
{
// first we select which minerals we should remove
std::vector<bool> minerals_to_remove(data->nb_mineral_kinetic, false);
for (int m=0; m<data->nb_mineral_kinetic; ++m)
{
for (auto it=id_component_to_remove.begin(); it!=id_component_to_remove.end(); ++it)
{
if (data->nu_mineral_kinetic(m, *it) != 0)
{
minerals_to_remove[m] = true;
break;
}
}
}
// then we remove them, in two steps
const int nb_min_to_keep = std::count(minerals_to_remove.begin(), minerals_to_remove.end(), false);
int new_m = 0;
// 1) first we copy data in the beginning of the arrays
for (int m=0; m<data->nb_mineral_kinetic; ++m)
{
if (minerals_to_remove[m]) continue;
data->logk_mineral_kinetic(new_m) = data->logk_mineral_kinetic(m);
data->labels_minerals_kinetic[new_m] = data->labels_minerals_kinetic[m];
data->_molar_volume_mineral_kinetic(new_m) = data->_molar_volume_mineral_kinetic(m);
int new_i = 0;
for (int i=0; i<data->nb_component; ++i)
{
if (is_component_to_remove(i)) continue;
data->nu_mineral_kinetic(new_m, new_i) = data->nu_mineral_kinetic(m, i);
++new_i;
}
assert(new_i == nb_component_to_keep());
++new_m;
}
assert(new_m == nb_min_to_keep);
// 2) then we resize the arrays
data->nu_mineral_kinetic.conservativeResize(nb_min_to_keep, nb_component_to_keep());
data->logk_mineral_kinetic.conservativeResize(nb_min_to_keep);
data->labels_minerals_kinetic.resize(nb_min_to_keep);
data->_molar_volume_mineral_kinetic.conservativeResize(nb_min_to_keep);
data->nb_mineral_kinetic = nb_min_to_keep;
// By doing in this order, we avoid bulk copy of arrays
}
void DatabaseSelector::select_gas(const std::vector<int>& id_component_to_remove)
{
// first we select which aqueous species we should remove
std::vector<bool> gas_to_remove(data->nb_gas, false);
for (int j=0; j<data->nb_gas; ++j)
{
for (auto it=id_component_to_remove.begin(); it!=id_component_to_remove.end(); ++it)
{
if (data->nu_gas(j, *it) != 0)
{
gas_to_remove[j] = true;
break;
}
}
}
// then we remove them, in two steps
const int nb_gas_to_keep = std::count(gas_to_remove.begin(), gas_to_remove.end(), false);
int new_j = 0;
// 1) first we copy data in the beginning of the arrays
for (int j=0; j<data->nb_gas; ++j)
{
if (gas_to_remove[j]) continue;
data->logk_gas(new_j) = data->logk_gas(j);
data->labels_gas[new_j] = data->labels_gas[j];
int new_i = 0;
for (int i=0; i<data->nb_component; ++i)
{
if (is_component_to_remove(i)) continue;
data->nu_gas(new_j, new_i) = data->nu_gas(j, i);
++new_i;
}
assert(new_i == nb_component_to_keep());
++new_j;
}
assert(new_j == nb_gas_to_keep);
// 2) then we resize the arrays
data->nu_gas.conservativeResize(nb_gas_to_keep, nb_component_to_keep());
data->logk_gas.conservativeResize(nb_gas_to_keep);
data->labels_gas.resize(nb_gas_to_keep);
data->nb_gas = nb_gas_to_keep;
// By doing in this order, we avoid bulk copy of arrays
}
void DatabaseSelector::select_components()
{
int new_i = 0;
for (int i=0; i<data->nb_component; ++i)
{
if (is_component_to_remove(i)) continue;
data->labels_basis[new_i] = data->labels_basis[i];
++new_i;
}
assert(new_i == nb_component_to_keep());
data->labels_basis.resize(nb_component_to_keep());
data->nb_component = nb_component_to_keep();
}
void DatabaseSelector::keep_only_component(const std::vector<int>& id_to_keep)
{
// First build the list of components to remove
std::vector<int> id_to_remove;
id_to_remove.reserve(data->nb_component - id_to_keep.size());
for (int id=0; id<data->nb_component; ++id)
{
auto search = std::find(id_to_keep.begin(), id_to_keep.end(), id);
if (search == id_to_keep.end())
id_to_remove.push_back(id);
}
// Then remove them
remove_component(id_to_remove);
}
} // end namespace database
} // end namespace specmicp
c4science · Help