Stoichiometric Matrix

Calculating stoichiometric matrices is a prerequisite for equilibrium calculations by minimizing Gibbs energy. The examples below show how they can be constructed.

Get Stoichiometric Matrix from a list of species

Let's imagine that we want to form the stoichiometric matrix of a list of solid and water species. For that, we need to read the database from which these species originate and retrieve the list of primary species from that database.

using ChemistryLab
using PrettyTables
df_elements, df_substances, df_reactions = read_thermofun("../../../data/cemdata18-merged.json")
df_primaries = extract_primary_species("../../../data/CEMDATA18-31-03-2022-phaseVol.dat")

It is then necessary to identify the list of secondary species likely to appear during the reactions.

given_species = filter(row -> row.symbol ∈ split("C3S Portlandite Jennite H2O@"), df_substances)
secondaries = filter(row->row.aggregate_state == "AS_AQUEOUS" 
                          && all(k->first(k) ∈ union_atoms(atoms.(given_species.species)), atoms(row.species))
                          && row.symbol ∉ split("H2@ O2@"),
                          df_substances)

We can then deduce the primary species concerned by the reaction.

all_species = unique(vcat(given_species, secondaries), :symbol)
species = [Species(f; symbol = phreeqc_to_unicode(n)) for (f, n) in zip(all_species.formula, all_species.symbol)]
candidate_primaries = [Species(f; symbol = phreeqc_to_unicode(n)) for (f, n) in zip(df_primaries.formula, df_primaries.symbol)]

And construct the stoichiometric matrix

SM = StoichMatrix(species, candidate_primaries)

Get Stoichiometric Matrix from a database file

using ChemistryLab
using PrettyTables
df_elements, df_substances, df_reactions = read_thermofun("../../../data/cemdata18-merged.json")
df_primaries = extract_primary_species("../../../data/CEMDATA18-31-03-2022-phaseVol.dat")

aqueous_species = filter(row->row.aggregate_state == "AS_AQUEOUS", df_substances)
species = [Species(f; symbol=phreeqc_to_unicode(n)) for (f,n) in zip(aqueous_species.formula, aqueous_species.symbol)]
candidate_primaries = [Species(f; symbol=phreeqc_to_unicode(n)) for (f,n) in zip(df_primaries.formula, df_primaries.symbol)]
SM = StoichMatrix(species, candidate_primaries) ;

All the independent reactions of the species contained in the database can thus be reconstructed. Here, only ionic species are listed given the choice to only read ionic species in the database ("AS_AQUEOUS").

67-element Vector{Reaction{Species{Int64}, Int64, Species{Int64}, Int64, Int64}}:
 Ca²⁺ + H₂O@ + SiO₂@ = CaSiO₃@ + 2H⁺
 H₂O@ + SiO₂@ = SiO₃²⁻ + 2H⁺
 H₂O@ + SiO₂@ + Sr²⁺ = SrSiO₃@ + 2H⁺
 CO₃²⁻ + 13H⁺ + NO₃⁻ + 10e⁻ = HCN@ + 6H₂O@
 CO₃²⁻ + 20H⁺ + NO₃⁻ + SO₄²⁻ + 15e⁻ = SCN⁻  + 10H₂O@
 AlO₂⁻ + H₂O@ + SiO₂@ = AlSiO₅³⁻ + 2H⁺
 2H₂O@ + 4SiO₂@ = Si₄O₁₀⁴⁻ + 4H⁺
 3Cl⁻ + FeO₂⁻ + 4H⁺ = FeCl₃@ + 2H₂O@
 AlO₂⁻ + 4H⁺ + 2SO₄²⁻ = Al(SO₄)₂⁻ + 2H₂O@
 FeO₂⁻ + 4H⁺ + SO₄²⁻ + e⁻ = Fe(SO₄)@ + 2H₂O@
 ⋮
 2H₂O@ = O₂@ + 4H⁺ + 4e⁻
 9H⁺ + NO₃⁻ + 8e⁻ = NH₃@ + 3H₂O@
 Ca²⁺ + CO₃²⁻ = CaCO₃@
 Ca²⁺ + SO₄²⁻ = CaSO₄@
 H₂O@ + SiO₂@ = HSiO₃⁻ + H⁺
 Cl⁻ + FeO₂⁻ + 4H⁺ = FeCl²⁺ + 2H₂O@
 Cl⁻ + FeO₂⁻ + 4H⁺ + e⁻ = FeCl⁺ + 2H₂O@
 K⁺ + SO₄²⁻ = KSO₄⁻
 CO₃²⁻ + Sr²⁺ = Sr(CO₃)@

The exercise can also be done on solid species. In this case, the data filter is carried out using the keyword "AS_CRYSTAL", in accordance with the terminology adopted in Thermofun.

140-element Vector{Reaction{Species, TR, Species, TP, Int64} where {TR<:Number, TP<:Number}}:
 9//2H₂O@ + 2AlO₂⁻ + Ca²⁺ + 2SiO₂@ = Ca(Al₂Si₂)O₈(H₂O)₉//₂
 6H₂O@ + 2AlO₂⁻ + Ca²⁺ + 4SiO₂@ = Ca(Al₂Si₄)O₁₂(H₂O)₆
 4H₂O@ + 3//2Mg²⁺ + 2SiO₂@ = Mg₃//₂Si₂O₁₁//₂(H₂O)₅//₂ + 3H⁺
 4H₂O@ + 3//2Mg²⁺ + SiO₂@ = Mg₃//₂SiO₇//₂(H₂O)₅//₂ + 3H⁺
 31//5H₂O@ + 2AlO₂⁻ + 2Na⁺ + 5//2SiO₂@ = Na₂(Al₂Si₅//₂)O₉(H₂O)₃₁//₅
 2H₂O@ + 2AlO₂⁻ + 2Na⁺ + 3SiO₂@ = Na₂(Al₂Si₃)O₁₀(H₂O)₂
 8H₂O@ + 2AlO₂⁻ + 2Na⁺ + 4SiO₂@ = Na₂(Al₂Si₄)O₁₂(H₂O)₈
 5H₂O@ + 14AlO₂⁻ + 12Ca²⁺ = (CaO)₁₂(Al₂O₃)₇ + 10H⁺
 4AlO₂⁻ + Ca²⁺ + 2H⁺ = CaO(Al₂O₃)₂ + H₂O@
 2AlO₂⁻ + Ca²⁺ = CaOAl₂O₃
 ⋮
 CO₃²⁻ + Mg²⁺ = MgCO₃
 SiO₂@ = SiO₂
 CO₃²⁻ + Sr²⁺ = SrCO₃
 SO₄²⁻ + Sr²⁺ = SrSO₄
 H₂O@ + SiO₂@ = (SiO₂H₂O)₁
 SiO₂@ = SiO₂
 H₂O@ + Ca²⁺ = CaO + 2H⁺
 H₂O@ + 2K⁺ = K₂O + 2H⁺
 H₂O@ + 2Na⁺ = Na₂O + 2H⁺

Or with gases ("AS_GAS")

6-element Vector{Reaction{Species{Int64}, Int64, Species{Int64}, Int64, Int64}}:
 2H₂O = O₂ + 4H⁺ + 4e⁻
 12H⁺ + 2NO₃⁻ + 10e⁻ = N₂ + 6H₂O
 CO₃²⁻ + 10H⁺ + 8e⁻ = CH₄ + 3H₂O
 CO₃²⁻ + 2H⁺ = CO₂ + H₂O
 2H⁺ + 2e⁻ = H₂
 10H⁺ + SO₄²⁻ + 8e⁻ = H₂S + 4H₂O