Specifying activity models#

Activity models are needed to accurately represent the non-ideal thermodynamic behavior of aqueous, liquid, gaseous, solid, and plasma phases. They are required when computing thermodynamic properties of phases and their species, such as molar enthalpies, molar heat capacities, molar volume, molar entropies, molar internal energies, species activities, species chemical potentials, etc.. These properties can only be accurate if proper activity models are assigned to the phases. For example, for an aqueous solution at extreme saline conditions, an activity model designed for low salinity conditions will not perform well.

In Reaktoro, an activity model assigned to a phase \(\pi\) is a function for which the following inputs are given:

  • temperature (in K), \(T\)

  • pressure (in Pa), \(P\)

  • mole fractions of the species in phase \(\pi\), \(x=(x_1,\ldots,x_\mathrm{N_{\pi}})\)

The evaluation of this activity model produces the following primary excess/residual thermodynamic properties for the corresponding phase and its species:

  • the excess/residual molar Gibbs energy of the phase (in units of J/mol), \(G_\mathrm{ex}\)

  • the excess/residual molar enthalpy of the phase (in units of J/mol), \(H_\mathrm{ex}\)

  • the excess/residual molar volume of the phase (in m3/mol), \(V_\mathrm{ex}\)

  • the temperature derivative of the excess/residual molar volume at constant pressure (in m3/(mol·K)), \((\partial V_\mathrm{ex}/\partial T)_P\)

  • the pressure derivative of the excess/residual molar volume at constant temperature (in m3/(mol·Pa)), \((\partial V_\mathrm{ex}/\partial P)_T\)

  • the excess/residual molar isobaric heat capacity of the phase (in units of J/(mol·K)), \(C_{P,\mathrm{ex}}\)

  • the activity coefficient \(\gamma_i\) of every species in the phase, for \(i=1,\ldots,\mathrm{N_{\pi}}\)

  • the activity \(a_i\) of every species in the phase, for \(i=1,\ldots,\mathrm{N_{\pi}}\)

Once we have computed these primary excess/residual thermodynamic properties, we can determine all others if needed (i.e., the secondary excess/residual thermodynamic properties). For example, the excess internal energy can be computed as \(U_\mathrm{ex} = H_\mathrm{ex} - P V_\mathrm{ex}\).

By default, all phases in Reaktoro are created with ideal activity models if not explicitly specified. The ideal models produce zero excess/residual properties for the phase and unit activity coefficients (thus, no correction for the non-ideal thermodynamic behavior of phases).

The code block that we present below:

from reaktoro import *

db = PhreeqcDatabase("phreeqc.dat")

solution = AqueousPhase(speciate("H O Na Cl C Ca Mg"))
gases = GaseousPhase("CO2(g) H2O(g)")
solidsolution = MineralPhase("Siderite Rhodochrosite")
ionexchange = IonExchangePhase("NaX CaX2")

system = ChemicalSystem(db, solution, gases, solidsolution, ionexchange)

is thus equivalent to this one:

from reaktoro import *

db = PhreeqcDatabase("phreeqc.dat")

solution = AqueousPhase(speciate("H O Na Cl C Ca Mg"))
solution.set(ActivityModelIdealAqueous())

gases = GaseousPhase("CO2(g) H2O(g)")
gases.set(ActivityModelIdealGas())

solidsolution = MineralPhase("Siderite Rhodochrosite")
solidsolution.set(ActivityModelIdealSolution(StateOfMatter.Solid))

ionexchange = IonExchangePhase("NaX CaX2")
ionexchange.set(ActivityModelIdealIonExchange())

system = ChemicalSystem(db, solution, gases, solidsolution, ionexchange)

Non-ideal thermodynamic behavior is often what most phases present in most conditions. Therefore, it is imperative to properly select an activity model for a phase suitable for the conditions of interest. Because there is no model that works best for all conditions, we leave the choice for the user.

We show in the following sections the available options for activity models for different phases and how to set and configure them.

Specifying activity models for aqueous phases#

The following are additional activity models available for aqueous phases in Reaktoro:

  • ActivityModelDavies

  • ActivityModelDebyeHuckel

  • ActivityModelDebyeHuckelKielland

  • ActivityModelDebyeHuckelLimitingLaw

  • ActivityModelDebyeHuckelParams

  • ActivityModelDebyeHuckelPHREEQC

  • ActivityModelDebyeHuckelWATEQ4F

  • ActivityModelHKF

  • ActivityModelPhreeqc

  • ActivityModelPitzer

For example, the code block below will construct a chemical system with an aqueous phase that uses PHREEQC’s activity model to compute the activities of the aqueous species:

solution = AqueousPhase(speciate("H O Na Cl"))
solution.set(ActivityModelPhreeqc(db))  # The PhreeqcDatabase object db contains the activity model parameters.

system = ChemicalSystem(db, solution)

Note

Ensuring consistency between PHREEQC’s thermodynamic databases and activity models is of paramount importance when modeling geochemical reactions. By opting for ActivityModelPhreeqc in Reaktoro when defining an aqueous phase whose species come from a PHREEQC database, the computation of the aqueous species activities mirrors that of PHREEQC, employing either the WATEQ-Debye-Hückel model [Truesdell and Jones, 1974] or the LLNL model [Wolery, 1992]. For a more detailed understanding of these models, including certain decisions like employing Davies for charged species in the absence of parameters, we refer to Parkhurst and Appelo [1999] and Parkhurst and Appelo [2013].

Let’s create a chemical state for this system (using ChemicalState class), and set it so that we have a 1 molal NaCl saline solution. We then compute its chemical properties (using class ChemicalProps) and print the activity coefficients, \(\gamma_i\), of the species.

state = ChemicalState(system)
state.temperature(25, "celsius")
state.pressure(1, "bar")
state.set("H2O", 1.0, "kg")
state.set("Na+", 1.0, "mol")
state.set("Cl-", 1.0, "mol")

print(f"{'Species':<20}{'Activity Coefficient'}")
props = ChemicalProps(state)
for species in system.species():
    print(f"{species.name():<20}{props.speciesActivityCoefficient(species.name())}")
Species             Activity Coefficient
H+                  0.7431715908205375
H2O                 1.000806912
Cl-                 0.6086519103778879
H2                  1.2589254117941673
Na+                 0.7312452684857986
OH-                 0.5790923293396838
NaOH                1.2589254117941673
O2                  1.2589254117941673

Note above that the activity coefficients for neutral species other than H2O are identical. This is because most aqueous activity models do not properly compute activity coefficients for these species (or their default variants do not contain model parameters for them).

One way to overcome or improve this is to combine an aqueous activity model with Setschenow model, so that the following equation:

\[\log_{10}\gamma_i = b_i I\]

is used instead to compute the activity coefficient of a neutral species \(i\), where \(I\) is the ionic strength of the solution and \(b_i\) is a Setschenow constant.

We recreate the chemical system in the next code block in which a chained activity model is constructed so that the above Setschenow equation is applied for the activity coefficients of species O2, H2, and NaOH:

solution = AqueousPhase(speciate("H O Na Cl"))
solution.set(chain(
    ActivityModelHKF(),
    ActivityModelSetschenow("O2", 0.123),
    ActivityModelSetschenow("H2", 0.234),
    ActivityModelSetschenow("NaOH", 0.345),
))

system = ChemicalSystem(db, solution)

We can then repeat the process of computing the chemical properties of the new system and printing the activity coefficients of the species:

state = ChemicalState(system)
state.temperature(25, "celsius")
state.pressure(1, "bar")
state.set("H2O", 1.0, "kg")
state.set("Na+", 1.0, "mol")
state.set("Cl-", 1.0, "mol")

print(f"{'Species':<20}{'Activity Coefficient'}")
props = ChemicalProps(state)
for species in system.species():
    print(f"{species.name():<20}{props.speciesActivityCoefficient(species.name())}")
Species             Activity Coefficient
H+                  0.6144274878149923
H2O                 1.0014153316335355
Cl-                 0.6651569578616804
H2                  1.7139573075084253
Na+                 0.6519897227119619
OH-                 0.7188624482783554
NaOH                2.2130947096056377
O2                  1.3273944577297396

Note that the activity coefficients of O2, H2, and NaOH are now different from each other. No changes were applied to the activity coefficients of H2O and charged species.

CO2 is a common gas dissolved in aqueous solutions. Most of the activity models above also do not properly calculate the activity of this aqueous species (with exception of ActivityModelPitzer). Thus, when using such models, it is recommended to create a chained activity model in which one of the following models below are combined with one of the above:

  • ActivityModelDrummond

  • ActivityModelDuanSun

  • ActivityModelRumpf

Here is an example:

solution = AqueousPhase(speciate("H O Na Cl C Ca Mg"))
solution.set(chain(
    ActivityModelHKF(),
    ActivityModelDrummond("CO2")
))

system = ChemicalSystem(db, solution)

With the change above, the Drummond model [Drummond, 1981] will be used for the computation of the activity of CO2 aqueous species (which can be named differently across the databases supported by Reaktoro, such as CO2 in PHREEQC databases, CO2(aq) in SUPCRT/SUPCRTB databases, or CO2@ in ThermoFun databases).

Tip

When computing gas solubilities in aqueous saline solutions, it is important that the activity of the dissolved gas species (e.g., CO2(aq)) is adequately computed.

Specifying activity models for gaseous and liquid phases#

Reaktoro implements a general form of the cubic equation of state (EOS) as presented in Smith et al. [2005] from which all classic cubic models are derived. You can assign one of the following models to a GaseousPhase or LiquidPhase object:

  • ActivityModelVanDerWaals (van der Waals [1873])

  • ActivityModelRedlichKwong (Redlich and Kwong [1949])

  • ActivityModelSoaveRedlichKwong (Soave [1972])

  • ActivityModelPengRobinson76 (Peng and Robinson [1976])

  • ActivityModelPengRobinson78 (Robinson and Peng [1978])

  • ActivityModelPengRobinson (identical to ActivityModelPengRobinson78)

The following variants of the Peng-Robinson EOS are also supported:

  • ActivityModelPengRobinsonPhreeqc, which combines ActivityModelPengRobinson with binary interaction parameters used in PHREEQC [Parkhurst and Appelo, 2013]

  • ActivityModelPengRobinsonSoreideWhitson, which combines ActivityModelPengRobinson with binary interaction parameters from Søreide and Whitson [1992]

  • ActivityModelPengRobinsonPhreeqcOriginal which is a ported version of the implementation of the Peng-Robinson EOS in PHREEQC [Parkhurst and Appelo, 2013]

Reaktoro also implements some equations of state designed for specific gaseous phases:

  • ActivityModelSpycherPruessEnnis for H2O-CO2 gas mixtures, Spycher et al. [2003]

  • ActivityModelSpycherReed for H2O-CO2-CH4 gas mixtures, Spycher and Reed [1988]

The example below demonstrates the use of the Peng-Robinson EOS as the activity model for a gaseous phase:

gases = GaseousPhase("CO2(g) CH4(g) H2O(g) O2(g) H2(g)")
gases.set(ActivityModelPengRobinson())

system = ChemicalSystem(db, gases)

Let’s create a chemical state for this system and compute its thermodynamic and chemical properties, from which we will obtain the fugacity coefficients of the gases:

state = ChemicalState(system)
state.temperature(100.0, "celsius")
state.pressure(1.0, "MPa")
state.set("CO2(g)", 0.80, "mol")
state.set("CH4(g)", 0.10, "mol")
state.set("H2O(g)", 0.05, "mol")
state.set("O2(g)",  0.03, "mol")
state.set("H2(g)",  0.02, "mol")

props = ChemicalProps(state)

print(f"{'Gas':<10}{'Fugacity Coefficient'}")
for i in range(system.species().size()):
    print(f"{system.species(i).name():<10}{props.speciesActivityCoefficient(i)}")
Gas       Fugacity Coefficient
CO2(g)    0.9738114880582376
CH4(g)    0.9927282935115731
H2O(g)    0.9307234180901945
O2(g)     1.0048772470887164
H2(g)     1.0215793244862943

Note

Note above that the method ChemicalProps.speciesActivityCoefficient is used for retrieving the fugacity coefficients of the gases. The activity coefficient of a gas is identical to its fugacity coefficient, and the fugacity of a gas has identical value to its activity when reference pressure is 1 bar, which is the case adopted in Reaktoro.

We print below the full computed properties of the system (in this case containing only a gaseous phase):

print(props)
+----------------------------------------+--------------+-----------+
| Property                               |        Value |      Unit |
+----------------------------------------+--------------+-----------+
| Temperature                            |     373.1500 |         K |
| Pressure                               |      10.0000 |       bar |
| Volume                                 |   3.0251e-03 |        m3 |
| Gibbs Energy                           | -117619.4634 |         J |
| Enthalpy                               |   16254.8858 |         J |
| Entropy                                |     358.7682 |       J/K |
| Internal Energy                        |   13229.8226 |         J |
| Helmholtz Energy                       | -120644.5266 |         J |
| Charge                                 |   0.0000e+00 |       mol |
| Element Amount:                        |              |           |
| :: H                                   |   5.4000e-01 |       mol |
| :: C                                   |   9.0000e-01 |       mol |
| :: O                                   |   1.7100e+00 |       mol |
| Species Amount:                        |              |           |
| :: CO2(g)                              |   8.0000e-01 |       mol |
| :: CH4(g)                              |   1.0000e-01 |       mol |
| :: H2O(g)                              |   5.0000e-02 |       mol |
| :: O2(g)                               |   3.0000e-02 |       mol |
| :: H2(g)                               |   2.0000e-02 |       mol |
| Mole Fraction:                         |              |           |
| :: CO2(g)                              |   8.0000e-01 |   mol/mol |
| :: CH4(g)                              |   1.0000e-01 |   mol/mol |
| :: H2O(g)                              |   5.0000e-02 |   mol/mol |
| :: O2(g)                               |   3.0000e-02 |   mol/mol |
| :: H2(g)                               |   2.0000e-02 |   mol/mol |
| Activity Coefficient:                  |              |           |
| :: CO2(g)                              |       0.9738 |         - |
| :: CH4(g)                              |       0.9927 |         - |
| :: H2O(g)                              |       0.9307 |         - |
| :: O2(g)                               |       1.0049 |         - |
| :: H2(g)                               |       1.0216 |         - |
| Activity:                              |              |           |
| :: CO2(g)                              |   7.7905e+00 |         - |
| :: CH4(g)                              |   9.9273e-01 |         - |
| :: H2O(g)                              |   4.6536e-01 |         - |
| :: O2(g)                               |   3.0146e-01 |         - |
| :: H2(g)                               |   2.0432e-01 |         - |
| lg(Activity):                          |              |           |
| :: CO2(g)                              |       0.8916 |         - |
| :: CH4(g)                              |      -0.0032 |         - |
| :: H2O(g)                              |      -0.3322 |         - |
| :: O2(g)                               |      -0.5208 |         - |
| :: H2(g)                               |      -0.6897 |         - |
| ln(Activity):                          |              |           |
| :: CO2(g)                              |       2.0529 |         - |
| :: CH4(g)                              |      -0.0073 |         - |
| :: H2O(g)                              |      -0.7649 |         - |
| :: O2(g)                               |      -1.1991 |         - |
| :: H2(g)                               |      -1.5881 |         - |
| Chemical Potential:                    |              |           |
| :: CO2(g)                              |  -1.2620e+05 |     J/mol |
| :: CH4(g)                              |  -2.5059e+05 |     J/mol |
| :: H2O(g)                              |  -2.3544e+03 |     J/mol |
| :: O2(g)                               |   4.4707e+05 |     J/mol |
| :: H2(g)                               |  -2.7242e+03 |     J/mol |
| Standard Volume:                       |              |           |
| :: CO2(g)                              |   0.0000e+00 |    m3/mol |
| :: CH4(g)                              |   0.0000e+00 |    m3/mol |
| :: H2O(g)                              |   0.0000e+00 |    m3/mol |
| :: O2(g)                               |   0.0000e+00 |    m3/mol |
| :: H2(g)                               |   0.0000e+00 |    m3/mol |
| Standard Gibbs Energy (formation):     |              |           |
| :: CO2(g)                              | -132568.1837 |     J/mol |
| :: CH4(g)                              | -250562.6135 |     J/mol |
| :: H2O(g)                              |      18.8176 |     J/mol |
| :: O2(g)                               |  450794.5594 |     J/mol |
| :: H2(g)                               |    2202.9019 |     J/mol |
| Standard Enthalpy (formation):         |              |           |
| :: CO2(g)                              |   29576.4051 |     J/mol |
| :: CH4(g)                              | -258559.1955 |     J/mol |
| :: H2O(g)                              |   41077.5707 |     J/mol |
| :: O2(g)                               |  563423.3786 |     J/mol |
| :: H2(g)                               |  -11992.3063 |     J/mol |
| Standard Entropy (formation):          |              |           |
| :: CO2(g)                              |     434.5292 | J/(mol*K) |
| :: CH4(g)                              |     -21.4299 | J/(mol*K) |
| :: H2O(g)                              |     110.0328 | J/(mol*K) |
| :: O2(g)                               |     301.8326 | J/(mol*K) |
| :: H2(g)                               |     -38.0416 | J/(mol*K) |
| Standard Internal Energy (formation):  |              |           |
| :: CO2(g)                              |   29576.4051 |     J/mol |
| :: CH4(g)                              | -258559.1955 |     J/mol |
| :: H2O(g)                              |   41077.5707 |     J/mol |
| :: O2(g)                               |  563423.3786 |     J/mol |
| :: H2(g)                               |  -11992.3063 |     J/mol |
| Standard Helmholtz Energy (formation): |              |           |
| :: CO2(g)                              | -132568.1837 |     J/mol |
| :: CH4(g)                              | -250562.6135 |     J/mol |
| :: H2O(g)                              |      18.8176 |     J/mol |
| :: O2(g)                               |  450794.5594 |     J/mol |
| :: H2(g)                               |    2202.9019 |     J/mol |
| Standard Heat Capacity (constant P):   |              |           |
| :: CO2(g)                              |     443.1980 | J/(mol*K) |
| :: CH4(g)                              |    -169.5362 | J/(mol*K) |
| :: H2O(g)                              |     -38.5340 | J/(mol*K) |
| :: O2(g)                               |    -167.9177 | J/(mol*K) |
| :: H2(g)                               |    -113.5628 | J/(mol*K) |
| Standard Heat Capacity (constant V):   |              |           |
| :: CO2(g)                              |     443.1980 | J/(mol*K) |
| :: CH4(g)                              |    -169.5362 | J/(mol*K) |
| :: H2O(g)                              |     -38.5340 | J/(mol*K) |
| :: O2(g)                               |    -167.9177 | J/(mol*K) |
| :: H2(g)                               |    -113.5628 | J/(mol*K) |
+----------------------------------------+--------------+-----------+

Specifying activity models for pure mineral and condensed phases#

Pure minerals and condensed phases (substances in solid or liquid states) constitute phases with a single species. For these pure phases, there is no need to specify an activity model other than the ideal. For solid solution phases, check the next section.

Specifying activity models for solid solution phases#

Reaktoro supports currently only the following activity models for solid solutions:

  • ActivityModelRedlichKister

The Redlich-Kister model assumes a solution with two solid species, and their activity coefficients, \(\gamma_1\) and \(\gamma_2\), are calculated as follows:

\[\ln\gamma_{1}=x_{2}^{2}[a_{0}+a_{1}(3x_{1}-x_{2})+a_{2}(x_{1}-x_{2})(5x_{1}-x_{2}),\]
\[\ln\gamma_{2}=x_{1}^{2}[a_{0}-a_{1}(3x_{2}-x_{1})+a_{2}(x_{2}-x_{1})(5x_{2}-x_{1}),\]

where \(a_0\), \(a_1\), \(a_2\) are model parameters; \(x_i\) is the mole fraction of the species \(i\) in the solution.

We show below how to create a chemical system with a single solid solution phase formed with minerals K-feldspar (KAlSi3O8) and albite (NaAlSi3O8) with assiged Redlich-Kister activity model:

a0, a1, a2 = 1.0, 2.0, 3.0  # the Redlich-Kister parameters (demonstration values!)

solidsolution = SolidPhase("K-feldspar Albite")
solidsolution.set(ActivityModelRedlichKister(a0, a1, a2))

system = ChemicalSystem(db, solidsolution)

The following code block creates a chemical state for this system with a single solid solution phase and prints its chemical properties.

state = ChemicalState(system)
state.set("K-feldspar", 0.5, "mol")
state.set("Albite", 0.5, "mol")

props = ChemicalProps(state)
print(props)
+----------------------------------------+------------+-----------+
| Property                               |      Value |      Unit |
+----------------------------------------+------------+-----------+
| Temperature                            |   298.1500 |         K |
| Pressure                               |     1.0000 |       bar |
| Volume                                 | 1.0473e-04 |        m3 |
| Gibbs Energy                           | 19899.3040 |         J |
| Enthalpy                               | 58945.9651 |         J |
| Entropy                                |   130.9631 |       J/K |
| Internal Energy                        | 58935.4921 |         J |
| Helmholtz Energy                       | 19888.8310 |         J |
| Charge                                 | 0.0000e+00 |       mol |
| Element Amount:                        |            |           |
| :: O                                   | 8.0000e+00 |       mol |
| :: Na                                  | 5.0000e-01 |       mol |
| :: Al                                  | 1.0000e+00 |       mol |
| :: Si                                  | 3.0000e+00 |       mol |
| :: K                                   | 5.0000e-01 |       mol |
| Species Amount:                        |            |           |
| :: K-feldspar :: KAlSi3O8              | 5.0000e-01 |       mol |
| :: Albite :: NaAlSi3O8                 | 5.0000e-01 |       mol |
| Mole Fraction:                         |            |           |
| :: K-feldspar :: KAlSi3O8              | 5.0000e-01 |   mol/mol |
| :: Albite :: NaAlSi3O8                 | 5.0000e-01 |   mol/mol |
| Activity Coefficient:                  |            |           |
| :: K-feldspar :: KAlSi3O8              |     2.1170 |         - |
| :: Albite :: NaAlSi3O8                 |     0.7788 |         - |
| Activity:                              |            |           |
| :: K-feldspar :: KAlSi3O8              | 1.0585e+00 |         - |
| :: Albite :: NaAlSi3O8                 | 3.8940e-01 |         - |
| lg(Activity):                          |            |           |
| :: K-feldspar :: KAlSi3O8              |     0.0247 |         - |
| :: Albite :: NaAlSi3O8                 |    -0.4096 |         - |
| ln(Activity):                          |            |           |
| :: K-feldspar :: KAlSi3O8              |     0.0569 |         - |
| :: Albite :: NaAlSi3O8                 |    -0.9431 |         - |
| Chemical Potential:                    |            |           |
| :: K-feldspar :: KAlSi3O8              | 1.2083e+04 |     J/mol |
| :: Albite :: NaAlSi3O8                 | 2.4279e+04 |     J/mol |
| Standard Volume:                       |            |           |
| :: K-feldspar :: KAlSi3O8              | 1.0815e-04 |    m3/mol |
| :: Albite :: NaAlSi3O8                 | 1.0131e-04 |    m3/mol |
| Standard Gibbs Energy (formation):     |            |           |
| :: K-feldspar :: KAlSi3O8              | 11941.9215 |     J/mol |
| :: Albite :: NaAlSi3O8                 | 26617.2081 |     J/mol |
| Standard Enthalpy (formation):         |            |           |
| :: K-feldspar :: KAlSi3O8              | 48025.2203 |     J/mol |
| :: Albite :: NaAlSi3O8                 | 68627.2314 |     J/mol |
| Standard Entropy (formation):          |            |           |
| :: K-feldspar :: KAlSi3O8              |   121.0240 | J/(mol*K) |
| :: Albite :: NaAlSi3O8                 |   140.9023 | J/(mol*K) |
| Standard Internal Energy (formation):  |            |           |
| :: K-feldspar :: KAlSi3O8              | 48014.4053 |     J/mol |
| :: Albite :: NaAlSi3O8                 | 68617.1004 |     J/mol |
| Standard Helmholtz Energy (formation): |            |           |
| :: K-feldspar :: KAlSi3O8              | 11931.1065 |     J/mol |
| :: Albite :: NaAlSi3O8                 | 26607.0771 |     J/mol |
| Standard Heat Capacity (constant P):   |            |           |
| :: K-feldspar :: KAlSi3O8              |  -123.5945 | J/(mol*K) |
| :: Albite :: NaAlSi3O8                 |  -123.5945 | J/(mol*K) |
| Standard Heat Capacity (constant V):   |            |           |
| :: K-feldspar :: KAlSi3O8              |  -123.5945 | J/(mol*K) |
| :: Albite :: NaAlSi3O8                 |  -123.5945 | J/(mol*K) |
+----------------------------------------+------------+-----------+

Specifying activity models for ion exchange phases#

By default, ion exchange phases are created in Reaktoro with an ideal activity model assigned to it. This ideal model computes the activity of the \(i\)-th exchange species (e.g., NaX, KX, CaX2) according to:

\[a_i = \dfrac{x_i z_{\mathrm{e},i}}{\sum_k x_k z_{\mathrm{e},k} } \]

where \(z_{\mathrm{e},i}\) is the exchanger’s equivalent (e.g., 1 for NaX, 2 for CaX2) and \(x_i\) is the mole fraction of the exchange species.

For non-ideal activity models for ion exchange phases, the following is currently available in Reaktoro:

  • ActivityModelIonExchangeGainesThomas

The Gaines-Thomas model is implemented according to what is used PHREEQC [Parkhurst and Appelo, 2013]. We show below a demonstration in which a chemical system with aqueous and ion exchange phases are created followed by the computation of the system’s chemical properties at a given chemical state:

db = PhreeqcDatabase("phreeqc.dat")

solution = AqueousPhase("H2O Na+ Cl- H+ OH- K+ Ca+2 Mg+2")
solution.set(ActivityModelPhreeqc(db))  # Ensure consistency with PHREEQC by using the same activity model PHREEQC would use with the phreeqc.dat database!

exchange = IonExchangePhase("NaX KX CaX2 MgX2")
exchange.set(ActivityModelIonExchangeGainesThomas())

system = ChemicalSystem(db, solution, exchange)

state = ChemicalState(system)
state.set("H2O" , 1.00, "kg")
state.set("Na+" , 1.00, "mmol")
state.set("K+"  , 1.00, "mmol")
state.set("Mg+2", 1.00, "mmol")
state.set("Ca+2", 1.00, "mmol")
state.set("NaX" , 0.06, "umol")
state.set("KX" ,  0.02, "umol")
state.set("CaX2" ,0.01, "umol")
state.set("MgX2" ,0.01, "umol")

props = ChemicalProps(state)
print(props)
+----------------------------------------+-------------+-----------+
| Property                               |       Value |      Unit |
+----------------------------------------+-------------+-----------+
| Temperature                            |    298.1500 |         K |
| Pressure                               |      1.0000 |       bar |
| Volume                                 |  1.0029e-03 |        m3 |
| Gibbs Energy                           |     -0.0002 |         J |
| Enthalpy                               |      0.0001 |         J |
| Entropy                                |      0.0000 |       J/K |
| Internal Energy                        |   -100.2933 |         J |
| Helmholtz Energy                       |   -100.2935 |         J |
| Charge                                 |  6.0000e-03 |       mol |
| Element Amount:                        |             |           |
| :: X                                   |  1.2000e-07 |       mol |
| :: H                                   |  1.1101e+02 |       mol |
| :: O                                   |  5.5506e+01 |       mol |
| :: Na                                  |  1.0001e-03 |       mol |
| :: Mg                                  |  1.0000e-03 |       mol |
| :: Cl                                  |  1.0000e-16 |       mol |
| :: K                                   |  1.0000e-03 |       mol |
| :: Ca                                  |  1.0000e-03 |       mol |
| Species Amount:                        |             |           |
| :: H2O                                 |  5.5506e+01 |       mol |
| :: Na+                                 |  1.0000e-03 |       mol |
| :: Cl-                                 |  1.0000e-16 |       mol |
| :: H+                                  |  1.0000e-16 |       mol |
| :: OH-                                 |  1.0000e-16 |       mol |
| :: K+                                  |  1.0000e-03 |       mol |
| :: Ca+2                                |  1.0000e-03 |       mol |
| :: Mg+2                                |  1.0000e-03 |       mol |
| :: NaX                                 |  6.0000e-08 |       mol |
| :: KX                                  |  2.0000e-08 |       mol |
| :: CaX2                                |  1.0000e-08 |       mol |
| :: MgX2                                |  1.0000e-08 |       mol |
| Mole Fraction:                         |             |           |
| :: H2O                                 |  9.9993e-01 |   mol/mol |
| :: Na+                                 |  1.8015e-05 |   mol/mol |
| :: Cl-                                 |  1.8015e-18 |   mol/mol |
| :: H+                                  |  1.8015e-18 |   mol/mol |
| :: OH-                                 |  1.8015e-18 |   mol/mol |
| :: K+                                  |  1.8015e-05 |   mol/mol |
| :: Ca+2                                |  1.8015e-05 |   mol/mol |
| :: Mg+2                                |  1.8015e-05 |   mol/mol |
| :: NaX                                 |  6.0000e-01 |   mol/mol |
| :: KX                                  |  2.0000e-01 |   mol/mol |
| :: CaX2                                |  1.0000e-01 |   mol/mol |
| :: MgX2                                |  1.0000e-01 |   mol/mol |
| Activity Coefficient:                  |             |           |
| :: H2O                                 |      1.0000 |         - |
| :: Na+                                 |      0.9278 |         - |
| :: Cl-                                 |      0.9265 |         - |
| :: H+                                  |      0.9336 |         - |
| :: OH-                                 |      0.9261 |         - |
| :: K+                                  |      0.9262 |         - |
| :: Ca+2                                |      0.7440 |         - |
| :: Mg+2                                |      0.7466 |         - |
| :: NaX                                 |      0.9276 |         - |
| :: KX                                  |      0.9261 |         - |
| :: CaX2                                |      0.7434 |         - |
| :: MgX2                                |      0.7460 |         - |
| Activity:                              |             |           |
| :: H2O                                 |  9.9993e-01 |         - |
| :: Na+                                 |  9.2783e-04 |         - |
| :: Cl-                                 |  9.2646e-17 |         - |
| :: H+                                  |  9.3362e-17 |         - |
| :: OH-                                 |  9.2608e-17 |         - |
| :: K+                                  |  9.2624e-04 |         - |
| :: Ca+2                                |  7.4401e-04 |         - |
| :: Mg+2                                |  7.4659e-04 |         - |
| :: NaX                                 |  4.6382e-01 |         - |
| :: KX                                  |  1.5434e-01 |         - |
| :: CaX2                                |  1.2391e-01 |         - |
| :: MgX2                                |  1.2434e-01 |         - |
| lg(Activity):                          |             |           |
| :: H2O                                 |     -0.0000 |         - |
| :: Na+                                 |     -3.0325 |         - |
| :: Cl-                                 |    -16.0332 |         - |
| :: H+                                  |    -16.0298 |         - |
| :: OH-                                 |    -16.0334 |         - |
| :: K+                                  |     -3.0333 |         - |
| :: Ca+2                                |     -3.1284 |         - |
| :: Mg+2                                |     -3.1269 |         - |
| :: NaX                                 |     -0.3336 |         - |
| :: KX                                  |     -0.8115 |         - |
| :: CaX2                                |     -0.9069 |         - |
| :: MgX2                                |     -0.9054 |         - |
| ln(Activity):                          |             |           |
| :: H2O                                 |     -0.0001 |         - |
| :: Na+                                 |     -6.9827 |         - |
| :: Cl-                                 |    -36.9177 |         - |
| :: H+                                  |    -36.9100 |         - |
| :: OH-                                 |    -36.9182 |         - |
| :: K+                                  |     -6.9844 |         - |
| :: Ca+2                                |     -7.2035 |         - |
| :: Mg+2                                |     -7.2000 |         - |
| :: NaX                                 |     -0.7683 |         - |
| :: KX                                  |     -1.8686 |         - |
| :: CaX2                                |     -2.0882 |         - |
| :: MgX2                                |     -2.0847 |         - |
| Chemical Potential:                    |             |           |
| :: H2O                                 | -1.6857e-01 |     J/mol |
| :: Na+                                 | -1.7310e+04 |     J/mol |
| :: Cl-                                 | -9.1518e+04 |     J/mol |
| :: H+                                  | -9.1498e+04 |     J/mol |
| :: OH-                                 | -1.1636e+04 |     J/mol |
| :: K+                                  | -1.7314e+04 |     J/mol |
| :: Ca+2                                | -1.7857e+04 |     J/mol |
| :: Mg+2                                | -1.7848e+04 |     J/mol |
| :: NaX                                 | -1.9045e+03 |     J/mol |
| :: KX                                  | -8.6277e+03 |     J/mol |
| :: CaX2                                | -9.7430e+03 |     J/mol |
| :: MgX2                                | -8.5928e+03 |     J/mol |
| Standard Volume:                       |             |           |
| :: H2O                                 |  1.8069e-05 |    m3/mol |
| :: Na+                                 | -1.5225e-06 |    m3/mol |
| :: Cl-                                 |  1.8045e-05 |    m3/mol |
| :: H+                                  |  0.0000e+00 |    m3/mol |
| :: OH-                                 | -4.1384e-06 |    m3/mol |
| :: K+                                  |  8.9794e-06 |    m3/mol |
| :: Ca+2                                | -1.8255e-05 |    m3/mol |
| :: Mg+2                                | -2.1937e-05 |    m3/mol |
| :: NaX                                 |  0.0000e+00 |    m3/mol |
| :: KX                                  |  0.0000e+00 |    m3/mol |
| :: CaX2                                |  0.0000e+00 |    m3/mol |
| :: MgX2                                |  0.0000e+00 |    m3/mol |
| Standard Gibbs Energy (formation):     |             |           |
| :: H2O                                 |      0.0000 |     J/mol |
| :: Na+                                 |      0.0000 |     J/mol |
| :: Cl-                                 |      0.0000 |     J/mol |
| :: H+                                  |      0.0000 |     J/mol |
| :: OH-                                 |  79882.1992 |     J/mol |
| :: K+                                  |      0.0000 |     J/mol |
| :: Ca+2                                |      0.0000 |     J/mol |
| :: Mg+2                                |      0.0000 |     J/mol |
| :: NaX                                 |     -0.0020 |     J/mol |
| :: KX                                  |  -3995.5945 |     J/mol |
| :: CaX2                                |  -4566.4315 |     J/mol |
| :: MgX2                                |  -3424.8346 |     J/mol |
| Standard Enthalpy (formation):         |             |           |
| :: H2O                                 |      0.0000 |     J/mol |
| :: Na+                                 |      0.0000 |     J/mol |
| :: Cl-                                 |      0.0000 |     J/mol |
| :: H+                                  |      0.0000 |     J/mol |
| :: OH-                                 |  56358.9533 |     J/mol |
| :: K+                                  |      0.0000 |     J/mol |
| :: Ca+2                                |      0.0000 |     J/mol |
| :: Mg+2                                |      0.0000 |     J/mol |
| :: NaX                                 |     -0.0020 |     J/mol |
| :: KX                                  |  -4299.9881 |     J/mol |
| :: CaX2                                |   7199.9758 |     J/mol |
| :: MgX2                                |   7399.9709 |     J/mol |
| Standard Entropy (formation):          |             |           |
| :: H2O                                 |      0.0000 | J/(mol*K) |
| :: Na+                                 |      0.0000 | J/(mol*K) |
| :: Cl-                                 |      0.0000 | J/(mol*K) |
| :: H+                                  |      0.0000 | J/(mol*K) |
| :: OH-                                 |    -78.8974 | J/(mol*K) |
| :: K+                                  |      0.0000 | J/(mol*K) |
| :: Ca+2                                |      0.0000 | J/(mol*K) |
| :: Mg+2                                |      0.0000 | J/(mol*K) |
| :: NaX                                 |      0.0000 | J/(mol*K) |
| :: KX                                  |     -1.0209 | J/(mol*K) |
| :: CaX2                                |     39.4647 | J/(mol*K) |
| :: MgX2                                |     36.3066 | J/(mol*K) |
| Standard Internal Energy (formation):  |             |           |
| :: H2O                                 |     -1.8069 |     J/mol |
| :: Na+                                 |      0.1522 |     J/mol |
| :: Cl-                                 |     -1.8045 |     J/mol |
| :: H+                                  |      0.0000 |     J/mol |
| :: OH-                                 |  56359.3671 |     J/mol |
| :: K+                                  |     -0.8979 |     J/mol |
| :: Ca+2                                |      1.8255 |     J/mol |
| :: Mg+2                                |      2.1937 |     J/mol |
| :: NaX                                 |     -0.0020 |     J/mol |
| :: KX                                  |  -4299.9881 |     J/mol |
| :: CaX2                                |   7199.9758 |     J/mol |
| :: MgX2                                |   7399.9709 |     J/mol |
| Standard Helmholtz Energy (formation): |             |           |
| :: H2O                                 |     -1.8069 |     J/mol |
| :: Na+                                 |      0.1522 |     J/mol |
| :: Cl-                                 |     -1.8045 |     J/mol |
| :: H+                                  |      0.0000 |     J/mol |
| :: OH-                                 |  79882.6131 |     J/mol |
| :: K+                                  |     -0.8979 |     J/mol |
| :: Ca+2                                |      1.8255 |     J/mol |
| :: Mg+2                                |      2.1937 |     J/mol |
| :: NaX                                 |     -0.0020 |     J/mol |
| :: KX                                  |  -3995.5945 |     J/mol |
| :: CaX2                                |  -4566.4315 |     J/mol |
| :: MgX2                                |  -3424.8346 |     J/mol |
| Standard Heat Capacity (constant P):   |             |           |
| :: H2O                                 |      0.0000 | J/(mol*K) |
| :: Na+                                 |      0.0000 | J/(mol*K) |
| :: Cl-                                 |      0.0000 | J/(mol*K) |
| :: H+                                  |      0.0000 | J/(mol*K) |
| :: OH-                                 |   -189.6441 | J/(mol*K) |
| :: K+                                  |      0.0000 | J/(mol*K) |
| :: Ca+2                                |      0.0000 | J/(mol*K) |
| :: Mg+2                                |      0.0000 | J/(mol*K) |
| :: NaX                                 |      0.0000 | J/(mol*K) |
| :: KX                                  |      0.0000 | J/(mol*K) |
| :: CaX2                                |      0.0000 | J/(mol*K) |
| :: MgX2                                |      0.0000 | J/(mol*K) |
| Standard Heat Capacity (constant V):   |             |           |
| :: H2O                                 |      0.0000 | J/(mol*K) |
| :: Na+                                 |      0.0000 | J/(mol*K) |
| :: Cl-                                 |      0.0000 | J/(mol*K) |
| :: H+                                  |      0.0000 | J/(mol*K) |
| :: OH-                                 |   -189.6441 | J/(mol*K) |
| :: K+                                  |      0.0000 | J/(mol*K) |
| :: Ca+2                                |      0.0000 | J/(mol*K) |
| :: Mg+2                                |      0.0000 | J/(mol*K) |
| :: NaX                                 |      0.0000 | J/(mol*K) |
| :: KX                                  |      0.0000 | J/(mol*K) |
| :: CaX2                                |      0.0000 | J/(mol*K) |
| :: MgX2                                |      0.0000 | J/(mol*K) |
+----------------------------------------+-------------+-----------+

Recap#

This guide demonstrated how you can create chemical systems with activity models of your choice assigned to the phases. Reaktoro is constantly under development and new activity models will be available over time. If you would like to help with the development of a specific activity model, please get in touch.

Attention

Ensure you have correctly selected and attached an activity model to a phase when comparing Reaktoro’s computations with other codes. When benchmarking, it is important to choose similar models. Ideally, you want to select the same one and use identical parameters for the model.