Stoichiometric Matrix

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

Stoichiometric matrix for a subset of species

The recommended workflow is:

Load all species from the database with build_species.
Filter to the relevant chemical space with speciation.
Build a ChemicalSystem, which automatically computes the stoichiometric matrices.

using ChemistryLab

# 1. Load all species from the database
all_species = build_species("../../../data/cemdata18-merged.json")

Identify the secondary species likely to appear during the reactions of interest (here: C₃S, Portlandite, Jennite and water):

# 2. Filter species compatible with the chosen seeds
species = speciation(all_species, split("C3S Portlandite Jennite H2O@");
              aggregate_state=[AS_AQUEOUS], exclude_species=split("H2@ O2@"))
dict_species = Dict(symbol(s) => s for s in species)

Deduce the primary species present in this subset:

# 3. Select primaries available in the subset
candidate_primaries = [dict_species[s] for s in CEMDATA_PRIMARIES if haskey(dict_species, s)]

Build the ChemicalSystem. Both stoichiometric matrices (CSM and SM) are computed once at construction time:

cs = ChemicalSystem(species, candidate_primaries)

14-element ChemicalSystem{Species, AbstractReaction, StoichMatrix{Real, Symbol, Vector{Symbol}, Matrix{Real}, Species}, StoichMatrix{Real, Species, Vector{Species}, Matrix{Real}, Species}, Nothing}:
 H2O@ {H2O  l} [H2O@ ◆ H₂O@]
 CaSiO3@ {CaSiO3  aq } [CaSiO3@ ◆ CaSiO₃@]
 SiO3-2 {SiO3-2  aq } [SiO3-2 ◆ SiO₃²⁻]
 Si4O10-4 {Si4O10-4  aq } [Si4O10-4 ◆ Si₄O₁₀⁴⁻]
 SiO2@ {SiO2  aq ( + 2 H2O = Si(OH)4  aq )} [SiO2@ ◆ SiO₂@]
 Ca+2 {Ca+2} [Ca+2 ◆ Ca²⁺]
 H+ {H+} [H+ ◆ H⁺]
 CaOH+ {CaOH+} [Ca(OH)+ ◆ Ca(OH)⁺]
 OH- {OH-} [OH- ◆ OH⁻]
 Ca(HSiO3)+ {CaHSiO3+  ( + H2O = CaSiO(OH)3+ )} [Ca(HSiO3)+ ◆ Ca(HSiO₃)⁺]
 HSiO3- {HSiO3-  ( + H2O = SiO(OH)3- )} [HSiO3- ◆ HSiO₃⁻]
 Jennite {Jennite   10/6 end-member of  id CSH SS (norm per 1 Si)  } [(SiO2)1(CaO)1.666667(H2O)2.1 ◆ (SiO₂)₁(CaO)₅//₃(H₂O)₂.₁]
 C3S {C3S} [(CaO)3SiO2 ◆ (CaO)₃SiO₂]
 Portlandite {Portlandite} [Ca(OH)2 ◆ Ca(OH)₂]

All independent reactions are then reconstructed from cs.SM:

lr = reactions(cs.SM)

10-element Vector{Reaction{Species, Real, Species, Real, Int64}}:
 H₂O@ + Ca²⁺ + SiO₂@ = CaSiO₃@ + 2H⁺
 H₂O@ + SiO₂@ = SiO₃²⁻ + 2H⁺
 2H₂O@ + 4SiO₂@ = Si₄O₁₀⁴⁻ + 4H⁺
 H₂O@ + Ca²⁺ = Ca(OH)⁺ + H⁺
 H₂O@ = OH⁻ + H⁺
 H₂O@ + Ca²⁺ + SiO₂@ = Ca(HSiO₃)⁺ + H⁺
 H₂O@ + SiO₂@ = HSiO₃⁻ + H⁺
 3.76667H₂O@ + 5//3Ca²⁺ + SiO₂@ = (SiO₂)₁(CaO)₅//₃(H₂O)₂.₁ + 10//3H⁺
 3H₂O@ + 3Ca²⁺ + SiO₂@ = (CaO)₃SiO₂ + 6H⁺
 2H₂O@ + Ca²⁺ = Ca(OH)₂ + 2H⁺

Stoichiometric matrix for all aqueous species in the database

To work with the full set of aqueous species, filter by aggregate state directly in speciation:

# Keep only aqueous species
aqueous_species = speciation(all_species, collect(keys(first(all_species).atoms));
                      aggregate_state=[AS_AQUEOUS])

A more direct alternative is to filter using all atom symbols present in the database:

all_atoms = union_atoms(all_species)
aqueous_species = speciation(all_species, all_atoms; aggregate_state=[AS_AQUEOUS])
dict_aqueous_species = Dict(symbol(s) => s for s in aqueous_species)
candidate_primaries = [dict_aqueous_species[s] for s in CEMDATA_PRIMARIES if haskey(dict_aqueous_species, s)]

cs_aq = ChemicalSystem(aqueous_species, candidate_primaries_aq)
lr_aq = reactions(cs_aq.SM)

67-element Vector{Reaction}:
 Si₄O₁₀⁴⁻ = 2SiO₃²⁻ + 2SiO₂@
 4HCN@ + 8AlSiO₅³⁻ + 15Si₄O₁₀⁴⁻ + 20Fe(SO₄)@ + 60e⁻ = 2H₂O@ + 68SiO₃²⁻ + 4SCN⁻  + 20FeCl₃@ + 8Al(SO₄)₂⁻
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 2Fe(SO₄)⁺ + 6e⁻ = 5SiO₃²⁻ + 2FeCl₃@ + Al(SO₄)₂⁻
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + Fe(SO₄)₂⁻ + 3e⁻ = 5SiO₃²⁻ + FeCl₃@ + Al(SO₄)₂⁻
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + Al(SO₄)₂⁻ = 5SiO₃²⁻ + 2Al(SO₄)⁺
 4SiO₃²⁻ + 4FeHCO₃⁺ = 2H₂O@ + Si₄O₁₀⁴⁻ + 4FeCO₃@
 4HCN@ + 8AlSiO₅³⁻ + 15Si₄O₁₀⁴⁻ + 20FeCO₃@ + 20Na(SO₄)⁻ + 60e⁻ = 2H₂O@ + 68SiO₃²⁻ + 4SCN⁻  + 20FeCl₃@ + 8Al(SO₄)₂⁻ + 20NaCO₃⁻
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 2Sr(SO₄)@ = 3SiO₃²⁻ + 2SrSiO₃@ + Al(SO₄)₂⁻
 FeCl₂⁺ + e⁻ = FeCl₃@
 29SiO₃²⁻ + 2SCN⁻  + AlSiO₅³⁻ + 10FeCl₃@ + 10SrHCO₃⁺ = 4H₂O@ + 10SrSiO₃@ + 2HCN@ + 5Si₄O₁₀⁴⁻ + Al(SO₄)₂⁻ + 10FeCO₃@ + 30e⁻
 2AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 4FeHSO₄²⁺ + 12e⁻ = 2H₂O@ + 6SiO₃²⁻ + 4FeCl₃@ + 2Al(SO₄)₂⁻
 4SiO₃²⁻ + 4NaHCO₃@ = 2H₂O@ + Si₄O₁₀⁴⁻ + 4NaCO₃⁻
 2HCN@ + 4AlSiO₅³⁻ + 5Si₄O₁₀⁴⁻ + 10FeHSO₄⁺ + 30e⁻ = 6H₂O@ + 24SiO₃²⁻ + 2SCN⁻  + 10FeCl₃@ + 4Al(SO₄)₂⁻
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ = 5SiO₃²⁻ + Al³⁺
 SiO₃²⁻ + Ca²⁺ = CaSiO₃@
 ∅ + e⁻ = Cl⁻
 4SiO₃²⁻ + 8HCN@ + 4Al(SO₄)₂⁻ + 5e⁻ = 4H₂O@ + 8SCN⁻  + 4AlSiO₅³⁻ + 5ClO₄⁻
 4HCN@ + 5Si₄O₁₀⁴⁻ + 2Al(SO₄)₂⁻ + 20Fe²⁺ + 60e⁻ = 2H₂O@ + 18SiO₃²⁻ + 4SCN⁻  + 2AlSiO₅³⁻ + 20FeCl₃@
 4SiO₃²⁻ + 4H⁺ = 2H₂O@ + Si₄O₁₀⁴⁻
 19SiO₃²⁻ + 2SCN⁻  + AlSiO₅³⁻ + 10FeCl₃@ + 10HCO₃⁻ = 4H₂O@ + 2HCN@ + 5Si₄O₁₀⁴⁻ + Al(SO₄)₂⁻ + 10FeCO₃@ + 30e⁻
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 2Mg(SO₄)@ = 5SiO₃²⁻ + Al(SO₄)₂⁻ + 2Mg²⁺
 4HCN@ + 5Si₄O₁₀⁴⁻ + 2Al(SO₄)₂⁻ + 20FeCO₃@ + 20Na⁺ + 60e⁻ = 2H₂O@ + 18SiO₃²⁻ + 4SCN⁻  + 2AlSiO₅³⁻ + 20FeCl₃@ + 20NaCO₃⁻
 49SiO₃²⁻ + 28HCN@ + 10FeCl₃@ + 9Al(SO₄)₂⁻ = 14H₂O@ + 18SCN⁻  + 9AlSiO₅³⁻ + 10Si₄O₁₀⁴⁻ + 10FeCO₃@ + 10NO₃⁻ + 30e⁻
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 2SO₄²⁻ = 5SiO₃²⁻ + Al(SO₄)₂⁻
 SrSiO₃@ = SiO₃²⁻ + Sr²⁺
 AlSiO₅³⁻ = SiO₃²⁻ + AlO₂⁻
 4HCN@ + 5Si₄O₁₀⁴⁻ + 2Al(SO₄)₂⁻ + 20FeCO₃@ + 60e⁻ = 2H₂O@ + 18SiO₃²⁻ + 4SCN⁻  + 2AlSiO₅³⁻ + 20FeCl₃@ + 20CO₃²⁻
 Si₄O₁₀⁴⁻ + FeO₂⁻ + 3e⁻ = 4SiO₃²⁻ + FeCl₃@
 32HCN@ + 5Si₄O₁₀⁴⁻ + 6Al(SO₄)₂⁻ + 20HS⁻ = 26H₂O@ + 14SiO₃²⁻ + 32SCN⁻  + 6AlSiO₅³⁻
 16SiO₃²⁻ + 4Al(OH)²⁺ = 2H₂O@ + 4AlSiO₅³⁻ + 3Si₄O₁₀⁴⁻
 2AlSiO₅³⁻ + Si₄O₁₀⁴⁻ = 6SiO₃²⁻ + 2AlO⁺
 8SiO₃²⁻ + 4AlO₂H@ = 2H₂O@ + 4AlSiO₅³⁻ + Si₄O₁₀⁴⁻
 Si₄O₁₀⁴⁻ + 4Ca(OH)⁺ = 2H₂O@ + 4CaSiO₃@
 Si₄O₁₀⁴⁻ + 4Fe(OH)²⁺ + 12e⁻ = 2H₂O@ + 4SiO₃²⁻ + 4FeCl₃@
 Si₄O₁₀⁴⁻ + 2FeO⁺ + 6e⁻ = 4SiO₃²⁻ + 2FeCl₃@
 3Si₄O₁₀⁴⁻ + 4FeO₂H@ + 12e⁻ = 2H₂O@ + 12SiO₃²⁻ + 4FeCl₃@
 2HCN@ + 5Si₄O₁₀⁴⁻ + Al(SO₄)₂⁻ + 10FeOH⁺ + 30e⁻ = 6H₂O@ + 19SiO₃²⁻ + 2SCN⁻  + AlSiO₅³⁻ + 10FeCl₃@
 8HCN@ + 6AlSiO₅³⁻ + 5Si₄O₁₀⁴⁻ + 20HSO₃⁻ = 14H₂O@ + 26SiO₃²⁻ + 8SCN⁻  + 6Al(SO₄)₂⁻
 2AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 4HSO₄⁻ = 2H₂O@ + 6SiO₃²⁻ + 2Al(SO₄)₂⁻
 Si₄O₁₀⁴⁻ + 4KOH@ = 2H₂O@ + 4SiO₃²⁻ + 4K⁺
 6SiO₃²⁻ + 2Al(SO₄)₂⁻ + 4Mg(OH)⁺ = 2H₂O@ + 2AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 4Mg(SO₄)@
 2SCN⁻  + AlSiO₅³⁻ + 5Si₄O₁₀⁴⁻ + 10FeCO₃@ + 10NH₄⁺ + 30e⁻ = 14H₂O@ + 21SiO₃²⁻ + 12HCN@ + 10FeCl₃@ + Al(SO₄)₂⁻
 2HCN@ + 5Si₄O₁₀⁴⁻ + Al(SO₄)₂⁻ + 10FeCO₃@ + 10NaOH@ + 30e⁻ = 6H₂O@ + 19SiO₃²⁻ + 2SCN⁻  + AlSiO₅³⁻ + 10FeCl₃@ + 10NaCO₃⁻
 Si₄O₁₀⁴⁻ + 4OH⁻ = 2H₂O@ + 4SiO₃²⁻
 16HCN@ + 2AlSiO₅³⁻ + 5Si₄O₁₀⁴⁻ + 10S₂O₃²⁻ = 8H₂O@ + 22SiO₃²⁻ + 16SCN⁻  + 2Al(SO₄)₂⁻
 4HCN@ + 3AlSiO₅³⁻ + 5Si₄O₁₀⁴⁻ + 10SO₃²⁻ = 2H₂O@ + 23SiO₃²⁻ + 4SCN⁻  + 3Al(SO₄)₂⁻
 Si₄O₁₀⁴⁻ + 4Sr(OH)⁺ = 2H₂O@ + 4SrSiO₃@
 Fe³⁺ + 3e⁻ = FeCl₃@
 74SiO₃²⁻ + 28HCN@ + 20FeCl₃@ + 14Al(SO₄)₂⁻ + 20CH₄@ = 54H₂O@ + 28SCN⁻  + 14AlSiO₅³⁻ + 15Si₄O₁₀⁴⁻ + 20FeCO₃@ + 60e⁻
 29SiO₃²⁻ + 2SCN⁻  + AlSiO₅³⁻ + 10FeCl₃@ + 10Ca(HCO₃)⁺ = 4H₂O@ + 10CaSiO₃@ + 2HCN@ + 5Si₄O₁₀⁴⁻ + Al(SO₄)₂⁻ + 10FeCO₃@ + 30e⁻
 4SiO₃²⁻ + 4Ca(HSiO₃)⁺ = 2H₂O@ + 4CaSiO₃@ + Si₄O₁₀⁴⁻
 3SiO₃²⁻ + 16HCN@ + 3Al(SO₄)₂⁻ + 10H₂S@ = 18H₂O@ + 16SCN⁻  + 3AlSiO₅³⁻
 4HCN@ + 8AlSiO₅³⁻ + 15Si₄O₁₀⁴⁻ + 20FeCO₃@ + 20Mg(SO₄)@ + 60e⁻ = 2H₂O@ + 68SiO₃²⁻ + 4SCN⁻  + 20FeCl₃@ + 8Al(SO₄)₂⁻ + 20Mg(CO₃)@
 22SiO₃²⁻ + SCN⁻  + 5FeCl₃@ + 2Al(SO₄)₂⁻ + 5Mg(HCO₃)⁺ = 2H₂O@ + HCN@ + 2AlSiO₅³⁻ + 5Si₄O₁₀⁴⁻ + 5FeCO₃@ + 5Mg(SO₄)@ + 15e⁻
 10SiO₃²⁻ + 2Al(SO₄)₂⁻ + 4Mg(HSiO₃)⁺ = 2H₂O@ + 2AlSiO₅³⁻ + 3Si₄O₁₀⁴⁻ + 4Mg(SO₄)@
 4HCN@ + 15Si₄O₁₀⁴⁻ + 2Al(SO₄)₂⁻ + 20FeCO₃@ + 60e⁻ = 2H₂O@ + 58SiO₃²⁻ + 4SCN⁻  + 2AlSiO₅³⁻ + 20FeCl₃@ + 20CO₂@
 68SiO₃²⁻ + 36HCN@ + 20FeCl₃@ + 8Al(SO₄)₂⁻ = 18H₂O@ + 16SCN⁻  + 8AlSiO₅³⁻ + 15Si₄O₁₀⁴⁻ + 20FeCO₃@ + 10N₂@ + 60e⁻
 SiO₃²⁻ + 2HCN@ + Al(SO₄)₂⁻ + 5H₂@ = 6H₂O@ + 2SCN⁻  + AlSiO₅³⁻
 2SiO₃²⁻ + 4HCN@ + 2Al(SO₄)₂⁻ = 2H₂O@ + 4SCN⁻  + 2AlSiO₅³⁻ + 5O₂@
 4SCN⁻  + 2AlSiO₅³⁻ + 15Si₄O₁₀⁴⁻ + 20FeCO₃@ + 20NH₃@ + 60e⁻ = 18H₂O@ + 62SiO₃²⁻ + 24HCN@ + 20FeCl₃@ + 2Al(SO₄)₂⁻
 20CaSiO₃@ + 4HCN@ + 5Si₄O₁₀⁴⁻ + 2Al(SO₄)₂⁻ + 20FeCO₃@ + 60e⁻ = 2H₂O@ + 38SiO₃²⁻ + 4SCN⁻  + 2AlSiO₅³⁻ + 20FeCl₃@ + 20CaCO₃@
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 2CaSO₄@ = 2CaSiO₃@ + 3SiO₃²⁻ + Al(SO₄)₂⁻
 4HSiO₃⁻ = 2H₂O@ + Si₄O₁₀⁴⁻
 FeCl²⁺ + 2e⁻ = FeCl₃@
 4HCN@ + 5Si₄O₁₀⁴⁻ + 2Al(SO₄)₂⁻ + 20FeCl⁺ + 40e⁻ = 2H₂O@ + 18SiO₃²⁻ + 4SCN⁻  + 2AlSiO₅³⁻ + 20FeCl₃@
 AlSiO₅³⁻ + Si₄O₁₀⁴⁻ + 2KSO₄⁻ = 5SiO₃²⁻ + Al(SO₄)₂⁻ + 2K⁺
 20SrSiO₃@ + 4HCN@ + 5Si₄O₁₀⁴⁻ + 2Al(SO₄)₂⁻ + 20FeCO₃@ + 60e⁻ = 2H₂O@ + 38SiO₃²⁻ + 4SCN⁻  + 2AlSiO₅³⁻ + 20FeCl₃@ + 20Sr(CO₃)@

Here only ionic species appear, given the AS_AQUEOUS filter.

The exercise can also be done on solid species (AS_CRYSTAL) or gases (AS_GAS). Primaries are still drawn from the aqueous subset in those cases:

# Solid species
solid_species = speciation(all_species, union_atoms(all_species); aggregate_state=[AS_CRYSTAL])
cs_solid = ChemicalSystem(solid_species, candidate_primaries)
lr_solid = reactions(cs_solid.SM)

# Gas species
gas_species = speciation(all_species, union_atoms(all_species); aggregate_state=[AS_GAS])
cs_gas = ChemicalSystem(gas_species, candidate_primaries)
lr_gas = reactions(cs_gas.SM)