Reaktoro 
A unified framework for modeling chemically reactive systems
PhreeqcUtils.hpp
1 // Reaktoro is a unified framework for modeling chemically reactive systems.
2 //
3 // Copyright (C) 2014-2015 Allan Leal
4 //
5 // This program is free software: you can redistribute it and/or modify
6 // it under the terms of the GNU General Public License as published by
7 // the Free Software Foundation, either version 3 of the License, or
8 // (at your option) any later version.
9 //
10 // This program is distributed in the hope that it will be useful,
11 // but WITHOUT ANY WARRANTY; without even the implied warranty of
12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 // GNU General Public License for more details.
14 //
15 // You should have received a copy of the GNU General Public License
16 // along with this program. If not, see <http://www.gnu.org/licenses/>.
17 
18 #pragma once
19 
20 // C++ includes
21 #include <map>
22 #include <string>
23 #include <vector>
24 
25 // Reaktoro includes
26 #include <Reaktoro/Math/Matrix.hpp>
27 #include <Reaktoro/Common/ScalarTypes.hpp>
28 
29 // Phreeqc includes
30 #include <Reaktoro/Interfaces/PhreeqcLegacy.hpp>
31 
32 namespace Reaktoro {
33 
34 // Forward declarations
35 class ReactionEquation;
36 
37 namespace PhreeqcUtils {
38 
42 auto load(PHREEQC& phreeqc, std::string database) -> void;
43 
48 auto execute(PHREEQC& phreeqc, std::string input, std::string output) -> void;
49 
54 auto findElement(const PHREEQC& phreeqc, std::string name) -> PhreeqcElement*;
55 
60 auto findSpecies(const PHREEQC& phreeqc, std::string name) -> PhreeqcSpecies*;
61 
66 auto findPhase(const PHREEQC& phreeqc, std::string name) -> PhreeqcPhase*;
67 
70 auto elements(const PhreeqcSpecies* species) -> std::map<PhreeqcElement*, double>;
71 
74 auto elements(const PhreeqcPhase* phase) -> std::map<PhreeqcElement*, double>;
75 
80 auto stoichiometry(std::string element, const PhreeqcSpecies* species) -> double;
81 
86 auto stoichiometry(std::string element, const PhreeqcPhase* phase) -> double;
87 
89 auto name(const PhreeqcElement* element) -> std::string;
90 
92 auto name(const PhreeqcSpecies* species) -> std::string;
93 
95 auto name(const PhreeqcPhase* phase) -> std::string;
96 
102 auto reactionEquation(const PhreeqcSpecies* species) -> ReactionEquation;
103 
108 auto reactionEquation(const PhreeqcPhase* phase) -> ReactionEquation;
109 
111 auto isAqueousSpecies(const PhreeqcSpecies* species) -> bool;
112 
114 auto isExchangeSpecies(const PhreeqcSpecies* species) -> bool;
115 
117 auto isGaseousSpecies(const PhreeqcPhase* phase) -> bool;
118 
120 auto isMineralSpecies(const PhreeqcPhase* phase) -> bool;
121 
125 auto index(std::string name, const std::vector<PhreeqcSpecies*>& species) -> unsigned;
126 
130 auto index(std::string name, const std::vector<PhreeqcPhase*>& phases) -> unsigned;
131 
134 auto activeAqueousSpecies(const PHREEQC& phreeqc) -> std::vector<PhreeqcSpecies*>;
135 
137 auto activeExchangeSpecies(const PHREEQC& phreeqc) -> std::vector<PhreeqcSpecies*>;
138 
142 auto activeProductSpecies(const PHREEQC& phreeqc) -> std::vector<PhreeqcSpecies*>;
143 
146 auto activeGaseousSpecies(const PHREEQC& phreeqc) -> std::vector<PhreeqcPhase*>;
147 
150 auto activePhasesInEquilibriumPhases(const PHREEQC& phreeqc) -> std::vector<PhreeqcPhase*>;
151 
154 auto activePhasesInSaturationList(const PHREEQC& phreeqc) -> std::vector<PhreeqcPhase*>;
155 
159 auto speciesAmounts(const PHREEQC& phreeqc, const std::vector<PhreeqcSpecies*>& species) -> Vector;
160 
164 auto speciesAmounts(const PHREEQC& phreeqc, const std::vector<PhreeqcPhase*>& phases) -> Vector;
165 
170 auto lnEquilibriumConstant(const PhreeqcSpecies* species, double T, double P) -> ThermoScalar;
171 
176 auto lnEquilibriumConstant(const PhreeqcPhase* phase, double T, double P) -> ThermoScalar;
177 
178 } // namespace PhreeqcUtils
179 } // namespace Reaktoro
auto index(const T &value, const std::vector< T > &values) -> Index
Find the index of a value in a container of values.
Definition: SetUtils.hxx:21
Eigen::VectorXd Vector
Define an alias to the vector type of the Eigen library.
Definition: Matrix.hpp:384
ThermoScalarBase< double > ThermoScalar
A type that defines a scalar thermo property.
Definition: ScalarTypes.hpp:40
auto database(std::string name) -> std::string
Return the contents of a built-in database as a string.
Definition: DatabaseUtils.cpp:65
The namespace containing all components of the Reaktoro library.
Definition: ChemicalScalar.hpp:24
auto elements() -> std::vector< std::string >
Return a vector of all known 116 chemical elements.
Definition: ElementUtils.cpp:155