diff --git a/src/specmicp/io/configuration.cpp b/src/specmicp/io/configuration.cpp index 28c2b16..cee6d51 100644 --- a/src/specmicp/io/configuration.cpp +++ b/src/specmicp/io/configuration.cpp @@ -1,617 +1,863 @@ /* ============================================================================= Copyright (c) 2014 - 2016 F. Georget 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. * ============================================================================= */ #include "configuration.hpp" -#include "../../utils/io/yaml.hpp" -#include "../adimensional/adimensional_system_solver_structs.hpp" -#include "../adimensional/config_default_options_solver.h" +#include "utils/io/yaml.hpp" +#include "utils/io/safe_config.hpp" -#include "../problem_solver/formulation.hpp" +#include "specmicp/adimensional/adimensional_system_solver_structs.hpp" +#include "specmicp/adimensional/config_default_options_solver.h" -#include "../../database.hpp" +#include "specmicp/problem_solver/formulation.hpp" -#include "../../physics/units.hpp" -#include "../../physics/io/configuration.hpp" -#include "../../physics/laws.hpp" +#include "database.hpp" + +#include "physics/units.hpp" +#include "physics/io/configuration.hpp" +#include "physics/laws.hpp" #define S_UNITS "units" #define S_SPECMICP "specmicp_options" #define S_SPECMICP_A_MAX_ITER "maximum_iteration" #define S_SPECMICP_A_RES_TOL "residual_tolerance" #define S_SPECMICP_A_STEP_TOL "step_tolerance" #define S_SPECMICP_A_MAX_STEP_LENGTH "maximum_step_length" #define S_SPECMICP_A_MAX_STEP_MAX_ITER "maximum_step_max_iteration" #define S_SPECMICP_A_SCALING "enable_scaling" #define S_SPECMICP_A_NONMONOTONE "enable_nonmonotone_linesearch" #define S_SPECMICP_A_DESCENT_DIRECTION "factor_descent_direction" #define S_SPECMICP_A_COND_CHECK "threshold_condition_check" #define S_SPECMICP_A_TRSHOLD_CYCLING_LSEARCH "threshold_cycling_linesearch" #define S_SPECMICP_A_NONIDEAL_TOL "non_ideality_tolerance" #define S_SPECMICP_A_CUTOFF_TOT_CONC "cutoff_total_concentration" #define S_FORMULATION "formulation" #define S_FORMULATION_A_SOLUTION "solution" #define S_FORMULATION_A_AQUEOUS "aqueous" #define S_FORMULATION_A_MINERALS "solid_phases" #define S_FORMULATION_A_AMOUNT "amount" #define S_FORMULATION_A_UNIT "unit" #define S_FORMULATION_A_LABEL "label" #define S_CONSTRAINTS "constraints" #define S_CONSTRAINTS_A_CHARGEKEEPER "charge_keeper" #define S_CONSTRAINTS_A_FIXEDACTIVITY "fixed_activity" #define S_CONSTRAINTS_A_FIXEDMOLALITY "fixed_molality" #define S_CONSTRAINTS_A_FIXEDFUGACITY "fixed_fugacity" #define S_CONSTRAINTS_A_SATURATED_SYSTEM "saturated_system" #define S_CONSTRAINTS_A_CONSERVATION_WATER "conservation_water" #define S_CONSTRAINTS_A_SURFACE_MODEL "surface_model" #define S_CONSTRAINTS_A_SURFACE_CONCENTRATION "surface_site_concentration" #define S_CONSTRAINTS_A_LABEL "label" #define S_CONSTRAINTS_A_AMOUNT "amount" #define S_CONSTRAINTS_A_AMOUNTLOG "amount_log" #define S_CONSTRAINTS_A_LABEL_COMPONENT "label_component" #define S_CONSTRAINTS_A_LABEL_GAS "label_gas" namespace specmicp { namespace io { // Additional declarations // ======================= void configure_specmicp_constraints_fixed_activity( AdimensionalSystemConstraints& constraints, const YAML::Node& conf_constraints_fixed_activity, const RawDatabasePtr& raw_db ); void configure_specmicp_constraints_fixed_fugacity( AdimensionalSystemConstraints& constraints, const YAML::Node& conf_constraints_fixed_fugacity, const RawDatabasePtr& raw_db ); void configure_specmicp_constraints_fixed_molality( AdimensionalSystemConstraints& constraints, const YAML::Node& conf_constraints_fixed_molality, const RawDatabasePtr& raw_db ); +void configure_specmicp_constraints_fixed_activity( + AdimensionalSystemConstraints& constraints, + const database::DataContainer * const raw_db, + YAMLConfigHandle&& conf_constraints_fixed_activity + ); + +void configure_specmicp_constraints_fixed_fugacity( + AdimensionalSystemConstraints& constraints, + const database::DataContainer * const raw_db, + YAMLConfigHandle&& conf_constraints_fixed_fugacity + ); + +void configure_specmicp_constraints_fixed_molality( + AdimensionalSystemConstraints& constraints, + const database::DataContainer * const raw_db, + YAMLConfigHandle&& conf_constraints_fixed_molality + ); + scalar_t molar_mass_aqueous_unknown_class(const std::string& label, const RawDatabasePtr& raw_db); scalar_t molar_mass_solid_unknown_class(const std::string& label, const RawDatabasePtr& raw_db ); scalar_t molar_volume_solid_unknown_class(const std::string& label, const RawDatabasePtr& raw_db ); // Implementation // ============== void configure_specmicp_options( AdimensionalSystemSolverOptions& options, const YAML::Node& configuration ) { check_mandatory_yaml_node(configuration, S_SPECMICP, "__main__"); const YAML::Node& conf = configuration[S_SPECMICP]; // fixme units //options.units_set = the_units_set; options.solver_options.fvectol = get_yaml_optional(conf, S_SPECMICP_A_RES_TOL, S_SPECMICP, SPECMICP_DEFAULT_RES_TOL); options.solver_options.steptol = get_yaml_optional(conf, S_SPECMICP_A_STEP_TOL, S_SPECMICP, SPECMICP_DEFAULT_STEP_TOL); options.solver_options.max_iter = get_yaml_optional(conf, S_SPECMICP_A_MAX_ITER, S_SPECMICP, SPECMICP_DEFAULT_MAX_ITER); options.solver_options.maxstep = get_yaml_optional(conf, S_SPECMICP_A_MAX_STEP_LENGTH, S_SPECMICP, SPECMICP_DEFAULT_MAX_STEP); options.solver_options.maxiter_maxstep = get_yaml_optional(conf, S_SPECMICP_A_MAX_STEP_MAX_ITER, S_SPECMICP, std::min(SPECMICP_DEFAULT_MAX_ITER_MAXSTEP, options.solver_options.max_iter) ); options.solver_options.use_scaling = get_yaml_optional(conf, S_SPECMICP_A_SCALING, S_SPECMICP, SPECMICP_DEFAULT_USE_SCALING); options.solver_options.non_monotone_linesearch = get_yaml_optional(conf, S_SPECMICP_A_NONMONOTONE, S_SPECMICP, SPECMICP_DEFAULT_USE_NONMONOTONE_LSEARCH); options.solver_options.factor_descent_condition = get_yaml_optional(conf, S_SPECMICP_A_DESCENT_DIRECTION, S_SPECMICP, SPECMICP_DEFAULT_FACTOR_DESC_COND ); options.solver_options.condition_limit = get_yaml_optional(conf, S_SPECMICP_A_COND_CHECK, S_SPECMICP, SPECMICP_DEFAULT_COND_THRESHOLD); options.solver_options.threshold_cycling_linesearch = get_yaml_optional(conf, S_SPECMICP_A_TRSHOLD_CYCLING_LSEARCH, S_SPECMICP, SPECMICP_DEFAULT_TRHSOLD_CYCLING_LSEARCH ); options.system_options.non_ideality_tolerance = get_yaml_optional(conf, S_SPECMICP_A_NONIDEAL_TOL, S_SPECMICP, SPECMICP_DEFAULT_NONIDEAL_TOL); options.system_options.cutoff_total_concentration = get_yaml_optional(conf, S_SPECMICP_A_CUTOFF_TOT_CONC, S_SPECMICP, SPECMICP_DEFAULT_CUTOFF_TOT_CONC); if (configuration[S_UNITS]) { options.units_set = configure_units(configuration[S_UNITS]); } } #define get_attribute(variable, attribute, type) \ if (configuration[attribute]) { \ options. variable = \ get_yaml_attribute( \ configuration, attribute, "AdimensionalSystemSolverOptions"); \ } //! \brief Configure specmicp solver options //! //! Conserve previous values if not given void SPECMICP_DLL_PUBLIC configure_specmicp_options_conservative( AdimensionalSystemSolverOptions& options, const YAML::Node& configuration ) { get_attribute(solver_options.fvectol, S_SPECMICP_A_RES_TOL, double) get_attribute(solver_options.steptol, S_SPECMICP_A_STEP_TOL, double) get_attribute(solver_options.max_iter, S_SPECMICP_A_MAX_ITER, int) get_attribute(solver_options.maxstep, S_SPECMICP_A_MAX_STEP_LENGTH, double) get_attribute(solver_options.maxiter_maxstep, S_SPECMICP_A_MAX_STEP_MAX_ITER, int) get_attribute(solver_options.use_scaling, S_SPECMICP_A_SCALING, bool) get_attribute(solver_options.non_monotone_linesearch, S_SPECMICP_A_NONMONOTONE, bool) get_attribute(solver_options.factor_descent_condition, S_SPECMICP_A_DESCENT_DIRECTION, double) get_attribute(solver_options.condition_limit, S_SPECMICP_A_COND_CHECK, double); get_attribute(solver_options.threshold_cycling_linesearch, S_SPECMICP_A_TRSHOLD_CYCLING_LSEARCH, double) get_attribute(system_options.non_ideality_tolerance, S_SPECMICP_A_NONIDEAL_TOL, double) get_attribute(system_options.cutoff_total_concentration, S_SPECMICP_A_CUTOFF_TOT_CONC, double) } #undef get_attribute void configure_specmicp_formulation( Formulation& formulation, const YAML::Node& conf_formulation, const RawDatabasePtr& raw_db, const units::UnitsSet& units_sets ) { configure_specmicp_formulation_solution(formulation, conf_formulation, units_sets); configure_specmicp_formulation_aqueous(formulation, conf_formulation, raw_db, units_sets); configure_specmicp_formulation_solid(formulation, conf_formulation, raw_db, units_sets); } void configure_specmicp_formulation_solution( Formulation& formulation, const YAML::Node& conf_formulation, const units::UnitsSet& _) { check_mandatory_yaml_node(conf_formulation, S_FORMULATION_A_SOLUTION, S_FORMULATION); const YAML::Node& sol_conf = conf_formulation[S_FORMULATION_A_SOLUTION]; check_mandatory_yaml_node(sol_conf, S_FORMULATION_A_AMOUNT, S_FORMULATION_A_SOLUTION); scalar_t factor=1.0; if (sol_conf[S_FORMULATION_A_UNIT]) { units::AmountUnit unit; units::parse_amount_unit(sol_conf[S_FORMULATION_A_UNIT].as(), unit); if (unit.type != units::AmountUnitType::Mass ) { throw std::invalid_argument("'" + sol_conf[S_FORMULATION_A_UNIT].as() + "'' is not a mass unit (in the solution configuration) !"); } factor = unit.factor_si; } formulation.set_mass_solution(sol_conf[S_FORMULATION_A_AMOUNT].as()*factor); } void configure_specmicp_formulation_aqueous( Formulation& formulation, const YAML::Node& conf_formulation, const RawDatabasePtr& raw_db, const units::UnitsSet& _) { check_mandatory_yaml_node(conf_formulation, S_FORMULATION_A_AQUEOUS, S_FORMULATION); const YAML::Node& aq_conf = conf_formulation[S_FORMULATION_A_AQUEOUS]; for (auto aqueous_conf: aq_conf) { check_mandatory_yaml_node(aqueous_conf, S_FORMULATION_A_LABEL, S_FORMULATION_A_AQUEOUS); check_mandatory_yaml_node(aqueous_conf, S_FORMULATION_A_AMOUNT, S_FORMULATION_A_AQUEOUS); check_mandatory_yaml_node(aqueous_conf, S_FORMULATION_A_UNIT, S_FORMULATION_A_AQUEOUS); std::string label = aqueous_conf[S_FORMULATION_A_LABEL].as(); scalar_t factor = 1.0; units::AmountUnit unit; units::parse_amount_unit(aqueous_conf[S_FORMULATION_A_UNIT].as(), unit); if (unit.type == units::AmountUnitType::Molality) { factor = unit.factor_si; } else if (unit.type == units::AmountUnitType::MoleConcentration) { factor = unit.factor_si/laws::density_water(units::celsius(25.0), units::LengthUnit::meter, units::MassUnit::kilogram); } else if (unit.type == units::AmountUnitType::MassConcentration) { const scalar_t molar_mass = molar_mass_aqueous_unknown_class(label, raw_db); factor = unit.factor_si / molar_mass / laws::density_water(units::celsius(25.0), units::LengthUnit::meter, units::MassUnit::kilogram); } else if (unit.type == units::AmountUnitType::Mass) { const scalar_t molar_mass = molar_mass_aqueous_unknown_class(label, raw_db); factor = unit.factor_si / molar_mass / formulation.mass_solution; } else if (unit.type == units::AmountUnitType::NumberOfMoles) { factor = unit.factor_si / formulation.mass_solution; } else { throw std::invalid_argument("'"+aqueous_conf[S_FORMULATION_A_UNIT].as() + "' is not a valid unit for an aqueous species"); } const scalar_t molality = factor*aqueous_conf[S_FORMULATION_A_AMOUNT].as(); formulation.add_aqueous_species(label, molality); } } void configure_specmicp_formulation_solid( Formulation& formulation, const YAML::Node& conf_formulation, const RawDatabasePtr& raw_db, const units::UnitsSet& unit_set) { check_mandatory_yaml_node(conf_formulation, S_FORMULATION_A_MINERALS, S_FORMULATION); const YAML::Node& minerals_conf = conf_formulation[S_FORMULATION_A_MINERALS]; for (auto solid_conf : minerals_conf) { check_mandatory_yaml_node(solid_conf, S_FORMULATION_A_LABEL, S_FORMULATION_A_AQUEOUS); check_mandatory_yaml_node(solid_conf, S_FORMULATION_A_AMOUNT, S_FORMULATION_A_AQUEOUS); check_mandatory_yaml_node(solid_conf, S_FORMULATION_A_UNIT, S_FORMULATION_A_AQUEOUS); const std::string label = solid_conf[S_FORMULATION_A_LABEL].as(); scalar_t factor = 1.0; units::AmountUnit unit; units::parse_amount_unit(solid_conf[S_FORMULATION_A_UNIT].as(), unit); // convert to mole concentration is SI if (unit.type == units::AmountUnitType::MoleConcentration or unit.type == units::AmountUnitType::NumberOfMoles) { factor = unit.factor_si; } else if (unit.type == units::AmountUnitType::MassConcentration or unit.type == units::AmountUnitType::Mass ) { const scalar_t molar_mass = molar_mass_solid_unknown_class(label, raw_db); factor = unit.factor_si / molar_mass; } else if (unit.type == units::AmountUnitType::Volume) { const scalar_t molar_volume = molar_volume_solid_unknown_class(label, raw_db); // may raise an error factor = unit.factor_si / molar_volume; } else { throw std::invalid_argument("'"+solid_conf[S_FORMULATION_A_UNIT].as() + "' is not a valid unit for a solid phase"); } // set to the correct units switch (unit_set.length) { case units::LengthUnit::decimeter: factor *= 1e-3; break; case units::LengthUnit::centimeter: factor *= 1e-6; break; case units::LengthUnit::meter: break; // just to remove warning } const scalar_t mole_conc = factor*solid_conf[S_FORMULATION_A_AMOUNT].as(); formulation.add_mineral(label, mole_conc); } } +// Constraints +// =========== + + +// new interface +// ------------- + +void configure_specmicp_constraints( + AdimensionalSystemConstraints& constraints, + const database::DataContainer * const raw_db, + YAMLConfigHandle&& conf_constraints + ) +{ + // equation for water + // ================== + // enable conservation of water + if (conf_constraints.has_attribute(S_CONSTRAINTS_A_CONSERVATION_WATER)) + { + + auto tmp_water = conf_constraints.get_attribute( + S_CONSTRAINTS_A_CONSERVATION_WATER); + if (tmp_water) + { + // only one constraints + if (conf_constraints.get_optional_attribute(S_CONSTRAINTS_A_SATURATED_SYSTEM, false)) + { + conf_constraints.report_error(YAMLConfigError::InvalidArgument, + "Attributes " S_CONSTRAINTS_A_CONSERVATION_WATER + "and " S_CONSTRAINTS_A_SATURATED_SYSTEM + " cannot be set to true at the same time"); + } + constraints.enable_conservation_water(); + } + else + constraints.disable_conservation_water(); + } + // or enable saturated system + if (conf_constraints.has_attribute(S_CONSTRAINTS_A_SATURATED_SYSTEM)) + { + if (conf_constraints.get_attribute(S_CONSTRAINTS_A_SATURATED_SYSTEM)) + { + constraints.set_saturated_system(); + } + } + // Aqueous components + // ================== + // charge keeper + // ------------- + if (conf_constraints.has_attribute(S_CONSTRAINTS_A_CHARGEKEEPER)) + { + const auto label = conf_constraints.get_attribute( + S_CONSTRAINTS_A_CHARGEKEEPER); + const index_t id = raw_db->get_id_component(label); + if (id == no_species) + { + conf_constraints.report_error( + YAMLConfigError::InvalidArgument, + "Species " + label + " does not exist in the database." + " It can't be the charge keeper." + ); + } + constraints.set_charge_keeper(id); + } + // Fixed activity + // --------------- + if (conf_constraints.has_section(S_CONSTRAINTS_A_FIXEDACTIVITY)) + { + configure_specmicp_constraints_fixed_activity( + constraints, + raw_db, + conf_constraints.get_section(S_CONSTRAINTS_A_FIXEDACTIVITY) + ); + } + // Fixed fugacity + // -------------- + if (conf_constraints.has_section(S_CONSTRAINTS_A_FIXEDFUGACITY)) + { + configure_specmicp_constraints_fixed_fugacity( + constraints, + raw_db, + conf_constraints.get_section(S_CONSTRAINTS_A_FIXEDFUGACITY) + ); + } + // Fixed activity + // --------------- + if (conf_constraints.has_section(S_CONSTRAINTS_A_FIXEDMOLALITY)) + { + configure_specmicp_constraints_fixed_molality( + constraints, + raw_db, + conf_constraints.get_section(S_CONSTRAINTS_A_FIXEDMOLALITY) + ); + } + // Surface model + // ============= + if (conf_constraints.has_attribute(S_CONSTRAINTS_A_SURFACE_MODEL)) + { + const bool is_enabled = conf_constraints.get_attribute( + S_CONSTRAINTS_A_SURFACE_MODEL + ); + if (is_enabled) + { + const auto value = conf_constraints.get_required_attribute( + S_CONSTRAINTS_A_SURFACE_CONCENTRATION); + constraints.enable_surface_model(value); + } + else + { + constraints.disable_surface_model(); + } + } +} + +// helpers functions + +index_t get_id_component( + const database::DataContainer* const raw_db, + YAMLConfigHandle& conf, + const std::string& attr_key=S_CONSTRAINTS_A_LABEL + ) +{ + const auto label = conf.get_required_attribute(attr_key); + const index_t id = raw_db->get_id_component(label); + if (id == no_species) + { + conf.report_error( + YAMLConfigError::InvalidArgument, + "'" + label + "' is not a valid component." + ); + } + return id; +} + +index_t get_id_gas( + const database::DataContainer* const raw_db, + YAMLConfigHandle& conf, + const std::string& attr_key=S_CONSTRAINTS_A_LABEL_GAS + ) +{ + const auto label = conf.get_required_attribute(attr_key); + const index_t id = raw_db->get_id_gas(label); + if (id == no_species) + { + conf.report_error( + YAMLConfigError::InvalidArgument, + "'" + label + "' is not a valid gas." + ); + } + return id; +} +scalar_t get_value(YAMLConfigHandle& conf) +{ + scalar_t value; + if (conf.has_attribute(S_CONSTRAINTS_A_AMOUNT)) { + value = std::log10(conf.get_attribute(S_CONSTRAINTS_A_AMOUNT)); + } + else if (conf.has_attribute(S_CONSTRAINTS_A_AMOUNTLOG)) { + value = conf.get_attribute(S_CONSTRAINTS_A_AMOUNTLOG); + } + else { + conf.report_error( + YAMLConfigError::InvalidArgument, + "Expected one of '" S_CONSTRAINTS_A_AMOUNT "' or '" + S_CONSTRAINTS_A_AMOUNTLOG "'." + ); + } + return value; +} + +void configure_specmicp_constraints_fixed_activity( + AdimensionalSystemConstraints& constraints, + const database::DataContainer * const raw_db, + YAMLConfigHandle&& conf_constraints_fixed_activity + ) +{ + for (uindex_t ind=0; ind(); const index_t id = raw_db->get_id_component(label); if (id == no_species) { throw std::invalid_argument("'"+label+"' is not a valid component (charge_keeper)."); } constraints.set_charge_keeper(id); } // Fixed activity // -------------- if (conf_constraints[S_CONSTRAINTS_A_FIXEDACTIVITY]) { configure_specmicp_constraints_fixed_activity(constraints, conf_constraints[S_CONSTRAINTS_A_FIXEDACTIVITY], raw_db); } // Fixed fugacity // -------------- if (conf_constraints[S_CONSTRAINTS_A_FIXEDFUGACITY]) { configure_specmicp_constraints_fixed_fugacity(constraints, conf_constraints[S_CONSTRAINTS_A_FIXEDFUGACITY], raw_db); } // Fixed molality // -------------- if (conf_constraints[S_CONSTRAINTS_A_FIXEDMOLALITY]) { configure_specmicp_constraints_fixed_molality(constraints, conf_constraints[S_CONSTRAINTS_A_FIXEDMOLALITY], raw_db); } // Water // ===== if (conf_constraints[S_CONSTRAINTS_A_CONSERVATION_WATER]) { if (conf_constraints[S_CONSTRAINTS_A_SATURATED_SYSTEM] and conf_constraints[S_CONSTRAINTS_A_SATURATED_SYSTEM].as()) { throw std::invalid_argument("Options '" S_CONSTRAINTS_A_CONSERVATION_WATER "' and '" S_CONSTRAINTS_A_SATURATED_SYSTEM "' cannot be used at the same time"); } const bool value = conf_constraints[S_CONSTRAINTS_A_CONSERVATION_WATER].as(); if (value) constraints.enable_conservation_water(); else constraints.disable_conservation_water(); } if (conf_constraints[S_CONSTRAINTS_A_SATURATED_SYSTEM]) { const bool value = conf_constraints[S_CONSTRAINTS_A_SATURATED_SYSTEM].as(); if (value) constraints.set_saturated_system(); } // Surface model // ============= if (conf_constraints[S_CONSTRAINTS_A_SURFACE_MODEL]) { const bool is_enabled = conf_constraints[S_CONSTRAINTS_A_SURFACE_MODEL].as(); if (is_enabled) { check_mandatory_yaml_node(conf_constraints, S_CONSTRAINTS_A_SURFACE_CONCENTRATION, S_CONSTRAINTS); const scalar_t value = conf_constraints[S_CONSTRAINTS_A_SURFACE_CONCENTRATION].as(); constraints.enable_surface_model(value); } else { constraints.disable_surface_model(); } } } + + void configure_specmicp_constraints_fixed_activity( AdimensionalSystemConstraints& constraints, const YAML::Node& conf_constraints_fixed_activity, const RawDatabasePtr& raw_db ) { for (auto fixact_conf: conf_constraints_fixed_activity) { check_mandatory_yaml_node(fixact_conf, S_CONSTRAINTS_A_LABEL, S_CONSTRAINTS_A_FIXEDACTIVITY); const std::string label = fixact_conf[S_CONSTRAINTS_A_LABEL].as(); const index_t id = raw_db->get_id_component(label); if (id == no_species) { throw std::invalid_argument("'"+label+"' is not a valid component (fixed_activity)."); } scalar_t value; if (fixact_conf[S_CONSTRAINTS_A_AMOUNT]) { value = std::log10(fixact_conf[S_CONSTRAINTS_A_AMOUNT].as()); } else if (fixact_conf[S_CONSTRAINTS_A_AMOUNTLOG]) { value = fixact_conf[S_CONSTRAINTS_A_AMOUNTLOG].as(); } else { throw std::invalid_argument("Either '" S_CONSTRAINTS_A_AMOUNT "' or '" S_CONSTRAINTS_A_AMOUNTLOG "' are required in section '" S_CONSTRAINTS_A_FIXEDACTIVITY "'."); } constraints.add_fixed_activity_component(id, value); } } + void configure_specmicp_constraints_fixed_fugacity( AdimensionalSystemConstraints& constraints, const YAML::Node& conf_constraints_fixed_fugacity, const RawDatabasePtr& raw_db ) { for(auto fixfug_conf: conf_constraints_fixed_fugacity) { check_mandatory_yaml_node(fixfug_conf, S_CONSTRAINTS_A_LABEL_COMPONENT, S_CONSTRAINTS_A_FIXEDFUGACITY); check_mandatory_yaml_node(fixfug_conf, S_CONSTRAINTS_A_LABEL_GAS, S_CONSTRAINTS_A_FIXEDFUGACITY); const std::string& label_c = fixfug_conf[S_CONSTRAINTS_A_LABEL_COMPONENT].as(); const index_t id_c = raw_db->get_id_component(label_c); if (id_c == no_species) { throw std::invalid_argument("'"+label_c+"' is not a valid component (fixed_fugacity)"); } const std::string& label_g = fixfug_conf[S_CONSTRAINTS_A_LABEL_GAS].as(); const index_t id_g = raw_db->get_id_gas(label_g); if (id_g == no_species) { throw std::invalid_argument("'"+label_g+"' is not a valid gas (fixed_fugacity)"); } scalar_t value; if (fixfug_conf[S_CONSTRAINTS_A_AMOUNT]) { value = std::log10(fixfug_conf[S_CONSTRAINTS_A_AMOUNT].as()); } else if (fixfug_conf[S_CONSTRAINTS_A_AMOUNTLOG]) { value = fixfug_conf[S_CONSTRAINTS_A_AMOUNTLOG].as(); } else { throw std::invalid_argument("Either '" S_CONSTRAINTS_A_AMOUNT "' or '" S_CONSTRAINTS_A_AMOUNTLOG "' are required in section '" S_CONSTRAINTS_A_FIXEDFUGACITY "'."); } constraints.add_fixed_fugacity_gas(id_g, id_c, value); } } + void configure_specmicp_constraints_fixed_molality( AdimensionalSystemConstraints& constraints, const YAML::Node& conf_constraints_fixed_molality, const RawDatabasePtr& raw_db ) { for (auto fixmol_conf: conf_constraints_fixed_molality) { check_mandatory_yaml_node(fixmol_conf, S_CONSTRAINTS_A_LABEL, S_CONSTRAINTS_A_FIXEDMOLALITY); const std::string& label = fixmol_conf[S_CONSTRAINTS_A_LABEL].as(); const index_t id = raw_db->get_id_component(label); if (id == no_species) { throw std::invalid_argument("'"+label+"' is not a valid component (fixed_molality)"); } scalar_t value; if (fixmol_conf[S_CONSTRAINTS_A_AMOUNT]) { value = std::log10(fixmol_conf[S_CONSTRAINTS_A_AMOUNT].as()); } else if (fixmol_conf[S_CONSTRAINTS_A_AMOUNTLOG]) { value = fixmol_conf[S_CONSTRAINTS_A_AMOUNTLOG].as(); } else { throw std::invalid_argument("Either '" S_CONSTRAINTS_A_AMOUNT "' or '" S_CONSTRAINTS_A_AMOUNTLOG "' are required in section '" S_CONSTRAINTS_A_FIXEDMOLALITY "'."); } constraints.add_fixed_molality_component(id, value); } } scalar_t molar_mass_aqueous_unknown_class(const std::string& label, const RawDatabasePtr& raw_db) { scalar_t molar_mass; if (raw_db->get_id_component(label) != no_species) { molar_mass = raw_db->molar_mass_basis(raw_db->get_id_component(label), units::MassUnit::kilogram); } else if (raw_db->get_id_aqueous(label) != no_species) { molar_mass = raw_db->molar_mass_aqueous(raw_db->get_id_aqueous(label), units::MassUnit::kilogram); } else if (raw_db->get_id_compound(label) != no_species) { molar_mass = raw_db->molar_mass_compound(raw_db->get_id_compound(label), units::MassUnit::kilogram); } else { throw std::invalid_argument("'"+label+"is not a valid species"); } return molar_mass; } scalar_t molar_mass_solid_unknown_class(const std::string& label, const RawDatabasePtr& raw_db ) { scalar_t molar_mass; if (raw_db->get_id_mineral(label) != no_species) { molar_mass = raw_db->molar_mass_mineral(raw_db->get_id_mineral(label), units::MassUnit::kilogram); } else if (raw_db->get_id_mineral_kinetic(label) != no_species) { molar_mass = raw_db->molar_mass_mineral_kinetic(raw_db->get_id_mineral_kinetic(label), units::MassUnit::kilogram); } else { throw std::invalid_argument("'"+label+"' is not a valid solid phase."); } return molar_mass; } scalar_t molar_volume_solid_unknown_class(const std::string& label, const RawDatabasePtr& raw_db ) { scalar_t molar_volume = -1; if (raw_db->get_id_mineral(label) != no_species) { molar_volume = raw_db->molar_volume_mineral(raw_db->get_id_mineral(label), units::LengthUnit::meter); } else if (raw_db->get_id_mineral_kinetic(label) != no_species) { molar_volume = raw_db->molar_volume_mineral_kinetic(raw_db->get_id_mineral_kinetic(label), units::LengthUnit::meter); } else { throw std::invalid_argument("'"+label+"' is not a valid solid phase."); } if (molar_volume < 0) { throw std::runtime_error("The molar volume for solid phase '"+label+"' is not defined"); } return molar_volume; } } //end namespace io } //end namespace specmicp diff --git a/src/specmicp/io/configuration.hpp b/src/specmicp/io/configuration.hpp index d08d0ad..ce23e0d 100644 --- a/src/specmicp/io/configuration.hpp +++ b/src/specmicp/io/configuration.hpp @@ -1,144 +1,179 @@ /* ============================================================================= Copyright (c) 2014 - 2016 F. Georget 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_SPECMICP_IO_HPP #define SPECMICP_SPECMICP_IO_HPP #include "types.hpp" #include //! \file specmicp/io/configuration.hpp //! \brief configuration for specmicp adim system namespace YAML { class Node; } //end namespace YAML namespace specmicp { struct AdimensionalSystemSolverOptions; struct Formulation; struct AdimensionalSystemConstraints; namespace database { struct DataContainer; } using RawDatabasePtr = std::shared_ptr; namespace units { struct UnitsSet; } //end namespace units namespace io { +class YAMLConfigHandle; + /*! \brief Configure specmicp solver options Erase values with default if not given Configuration must be of the form : \code specmicp_options: opts1: val1 opts2: val2 .... [units]: ... \endcode */ void SPECMICP_DLL_PUBLIC configure_specmicp_options( AdimensionalSystemSolverOptions& options, const YAML::Node& configuration - ); + ) SPECMICP_DEPRECATED("Use new safer config"); + + /*! \brief Configure specmicp solver options Conserve previous values if not given. Configuration must be of the form : \code opts1: val1 opts2: val2 .... \endcode \warning units must be set elsewhere ! */ void SPECMICP_DLL_PUBLIC configure_specmicp_options_conservative( AdimensionalSystemSolverOptions& options, const YAML::Node& configuration - ); + ) SPECMICP_DEPRECATED("Use new safer config"); + + +/*! \brief Configure specmicp solver options + +Conserve previous values if not given. + +Configuration must be of the form : +\code + opts1: val1 + opts2: val2 + .... +\endcode +\warning units must be set elsewhere ! +*/ +void SPECMICP_DLL_PUBLIC configure_specmicp_options_conservative( + AdimensionalSystemSolverOptions& options, + const YAML::Node& configuration + ) SPECMICP_DEPRECATED("Use new safer config"); void SPECMICP_DLL_PUBLIC configure_specmicp_formulation( Formulation& formulation, const YAML::Node& conf_formulation, const RawDatabasePtr& raw_db, const units::UnitsSet& units_sets ); void SPECMICP_DLL_PUBLIC configure_specmicp_formulation_solution( Formulation& formulation, const YAML::Node& conf_formulation, const units::UnitsSet& units_set ); void SPECMICP_DLL_PUBLIC configure_specmicp_formulation_aqueous( Formulation& formulation, const YAML::Node& conf_formulation, const RawDatabasePtr& raw_db, const units::UnitsSet& _ ); void SPECMICP_DLL_PUBLIC configure_specmicp_formulation_solid( Formulation& formulation, const YAML::Node& conf_formulation, const RawDatabasePtr& raw_db, const units::UnitsSet& _ ); +//! \brief Configure specmicp constraints +void SPECMICP_DLL_PUBLIC configure_specmicp_constraints( + AdimensionalSystemConstraints& constraints, + const database::DataContainer* const raw_db, + YAMLConfigHandle&& conf_constraints + ); + +//! \brief Configure specmicp constraints +//! +//! \deprecated void SPECMICP_DLL_PUBLIC configure_specmicp_constraints( AdimensionalSystemConstraints& constraints, const YAML::Node& conf_constraints, const RawDatabasePtr& raw_db - ); + ) SPECMICP_DEPRECATED("Use new safer config"); + + + } //end namespace io } //end namespace specmicp #endif // SPECMICP_SPECMICP_IO_HPP