Reaktoro 
A unified framework for modeling chemically reactive systems
Phreeqc.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 // Reaktoro includes
21 #include <Reaktoro/Interfaces/Interface.hpp>
22 
23 // Forward declarations
24 class PHREEQC;
25 
26 namespace Reaktoro {
27 
28 // Forward declarations
29 class ReactionEquation;
30 
31 class Phreeqc : public Interface
32 {
33 public:
35  Phreeqc();
36 
39  Phreeqc(std::string database);
40 
42  virtual ~Phreeqc();
43 
45  virtual auto temperature() const -> double;
46 
48  virtual auto pressure() const -> double;
49 
51  virtual auto speciesAmounts() const -> Vector;
52 
54  virtual auto numElements() const -> unsigned;
55 
57  virtual auto numSpecies() const -> unsigned;
58 
60  virtual auto numPhases() const -> unsigned;
61 
63  virtual auto numSpeciesInPhase(Index iphase) const -> unsigned;
64 
66  virtual auto elementName(Index ielement) const -> std::string;
67 
69  virtual auto elementMolarMass(Index ielement) const -> double;
70 
72  virtual auto elementStoichiometry(Index ispecies, Index ielement) const -> double;
73 
75  virtual auto speciesName(Index ispecies) const -> std::string;
76 
78  virtual auto phaseName(Index iphase) const -> std::string;
79 
84  virtual auto properties(Index iphase, double T, double P) -> PhaseThermoModelResult;
85 
91  virtual auto properties(Index iphase, double T, double P, const Vector& nphase) -> PhaseChemicalModelResult;
92 
94  virtual auto clone() const -> std::shared_ptr<Interface>;
95 
102  auto set(double T, double P) -> void;
103 
112  auto set(double T, double P, const Vector& n) -> void;
113 
120  auto load(std::string database) -> void;
121 
129  auto execute(std::string input, std::string output) -> void;
130 
133  auto execute(std::string input) -> void;
134 
139  auto reset() -> void;
140 
142  auto reactions() const -> std::vector<ReactionEquation>;
143 
145  auto stoichiometricMatrix() const -> Matrix;
146 
148  auto standardMolarGibbsEnergies() const -> Vector;
149 
151  auto standardMolarEnthalpies() const -> Vector;
152 
154  auto standardMolarVolumes() const -> Vector;
155 
158 
161 
163  auto lnActivityCoefficients() const -> Vector;
164 
166  auto lnActivityConstants() const -> Vector;
167 
169  auto lnActivities() const -> Vector;
170 
172  auto lnEquilibriumConstants() const -> Vector;
173 
175  auto phaseMolarVolumes() const -> Vector;
176 
178  auto phreeqc() -> PHREEQC&;
179 
181  auto phreeqc() const -> const PHREEQC&;
182 
183 private:
184  struct Impl;
185 
186  std::shared_ptr<Impl> pimpl;
187 };
188 
189 } // namespace Reaktoro
Eigen::MatrixXd Matrix
Define an alias to the matrix type of the Eigen library.
Definition: Matrix.hpp:387
auto phreeqc() -> PHREEQC &
Return a reference to the low-level Phreeqc instance.
Definition: Phreeqc.cpp:1532
auto load(std::string database) -> void
Load a PHREEQC database.
auto execute(std::string input, std::string output) -> void
Execute a PHREEQC input script either provided as a file or input string.
virtual ~Phreeqc()
Destroy this Phreeqc instance.
Definition: Phreeqc.cpp:1380
virtual auto pressure() const -> double
Return the pressure (in units of Pa)
Definition: Phreeqc.cpp:1391
virtual auto properties(Index iphase, double T, double P) -> PhaseThermoModelResult
Return the thermodynamic properties of a phase.
virtual auto elementStoichiometry(Index ispecies, Index ielement) const -> double
Return the stoichiometry of an element in a species.
Definition: Phreeqc.cpp:1439
virtual auto elementName(Index ielement) const -> std::string
Return the name of an element.
Definition: Phreeqc.cpp:1427
Definition: Phreeqc.hpp:31
virtual auto numPhases() const -> unsigned
Return the number of phases.
Definition: Phreeqc.cpp:1415
auto standardMolarVolumes() const -> Vector
Return the standard molar volumes of the species (in units of m3/mol)
Definition: Phreeqc.cpp:1496
virtual auto clone() const -> std::shared_ptr< Interface >
Return a clone of this Phreeqc instance.
Definition: Phreeqc.cpp:1469
virtual auto numSpecies() const -> unsigned
Return the number of species.
Definition: Phreeqc.cpp:1409
The result of the chemical model function that calculates the chemical properties of a phase...
Definition: PhaseChemicalModel.hpp:30
auto standardMolarHeatCapacitiesConstV() const -> Vector
Return the standard molar isochoric heat capacities of the species (in units of J/(mol*K)) ...
Definition: Phreeqc.cpp:1508
virtual auto elementMolarMass(Index ielement) const -> double
Return the molar mass of an element (in units of kg/mol)
Definition: Phreeqc.cpp:1433
auto standardMolarEnthalpies() const -> Vector
Return the standard molar enthalpies of the species (in units of J/mol)
Definition: Phreeqc.cpp:1490
virtual auto temperature() const -> double
Return the temperature (in units of K)
Definition: Phreeqc.cpp:1385
The result of the thermodynamic model function that calculates the standard thermodynamic properties ...
Definition: PhaseThermoModel.hpp:29
auto phaseMolarVolumes() const -> Vector
Return the molar volumes of the phases.
Definition: Phreeqc.cpp:1526
auto lnActivities() const -> Vector
Return the ln activities of the species.
Definition: Phreeqc.cpp:1520
virtual auto speciesAmounts() const -> Vector
Return the amounts of the species (in units of mol)
Definition: Phreeqc.cpp:1397
Phreeqc()
Construct a default Phreeqc instance.
Definition: Phreeqc.cpp:1370
auto stoichiometricMatrix() const -> Matrix
Return the stoichiometric matrix of the system of reactions.
Eigen::VectorXd Vector
Define an alias to the vector type of the Eigen library.
Definition: Matrix.hpp:384
A class used to interface other codes with Reaktoro.
Definition: Interface.hpp:36
auto reset() -> void
Reset this Phreeqc instance to a clean state.
virtual auto speciesName(Index ispecies) const -> std::string
Return the name of a species.
Definition: Phreeqc.cpp:1445
std::size_t Index
Define a type that represents an index.
Definition: Index.hpp:26
virtual auto phaseName(Index iphase) const -> std::string
Return the name of a phase.
Definition: Phreeqc.cpp:1451
virtual auto numSpeciesInPhase(Index iphase) const -> unsigned
Return the number of species in a phase.
Definition: Phreeqc.cpp:1421
auto standardMolarHeatCapacitiesConstP() const -> Vector
Return the standard molar isobaric heat capacities of the species (in units of J/(mol*K)) ...
Definition: Phreeqc.cpp:1502
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
virtual auto numElements() const -> unsigned
Return the number of elements.
Definition: Phreeqc.cpp:1403
auto standardMolarGibbsEnergies() const -> Vector
Return the standard molar Gibbs free energies of the species (in units of J/mol)
Definition: Phreeqc.cpp:1484
auto lnActivityCoefficients() const -> Vector
Return the ln activity coefficients of the species.
Definition: Phreeqc.cpp:1514
auto lnEquilibriumConstants() const -> Vector
Return the ln equilibrium constants of the reactions.
auto lnActivityConstants() const -> Vector
Return the ln activity contants of the species.
auto reactions() const -> std::vector< ReactionEquation >
Return the system of reactions.