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}}:
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{Species{Int64}, Int64, Species{Int64}, Int64, Int64}}:
H₂O@ + Ca²⁺ + 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@
H₂O@ + AlO₂⁻ + 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@
FeO₂⁻ + 4H⁺ + SO₄²⁻ = Fe(SO₄)⁺ + 2H₂O@
FeO₂⁻ + 4H⁺ + 2SO₄²⁻ = Fe(SO₄)₂⁻ + 2H₂O@
CO₃²⁻ + FeO₂⁻ + 4H⁺ + e⁻ = FeCO₃@ + 2H₂O@
CO₃²⁻ + Na⁺ = NaCO₃⁻
AlO₂⁻ + 4H⁺ + SO₄²⁻ = Al(SO₄)⁺ + 2H₂O@
CO₃²⁻ + FeO₂⁻ + 5H⁺ + e⁻ = FeHCO₃⁺ + 2H₂O@
Mg²⁺ + SO₄²⁻ = Mg(SO₄)@
Na⁺ + SO₄²⁻ = Na(SO₄)⁻
SO₄²⁻ + Sr²⁺ = Sr(SO₄)@
2Cl⁻ + FeO₂⁻ + 4H⁺ = FeCl₂⁺ + 2H₂O@
CO₃²⁻ + H⁺ + Sr²⁺ = SrHCO₃⁺
FeO₂⁻ + 5H⁺ + SO₄²⁻ = FeHSO₄²⁺ + 2H₂O@
CO₃²⁻ + H⁺ + Na⁺ = NaHCO₃@
FeO₂⁻ + 5H⁺ + SO₄²⁻ + e⁻ = FeHSO₄⁺ + 2H₂O@
AlO₂⁻ + 4H⁺ = Al³⁺ + 2H₂O@
4H₂O@ + Cl⁻ = ClO₄⁻ + 8H⁺ + 8e⁻
FeO₂⁻ + 4H⁺ + e⁻ = Fe²⁺ + 2H₂O@
CO₃²⁻ + H⁺ = HCO₃⁻
9H⁺ + SO₄²⁻ + 8e⁻ = HS⁻ + 4H₂O@
AlO₂⁻ + 3H⁺ = Al(OH)²⁺ + H₂O@
AlO₂⁻ + 2H⁺ = AlO⁺ + H₂O@
AlO₂⁻ + H⁺ = AlO₂H@
H₂O@ + Ca²⁺ = Ca(OH)⁺ + H⁺
FeO₂⁻ + 3H⁺ = Fe(OH)²⁺ + H₂O@
FeO₂⁻ + 2H⁺ = FeO⁺ + H₂O@
FeO₂⁻ + H⁺ = FeO₂H@
FeO₂⁻ + 3H⁺ + e⁻ = FeOH⁺ + H₂O@
3H⁺ + SO₄²⁻ + 2e⁻ = HSO₃⁻ + H₂O@
H⁺ + SO₄²⁻ = HSO₄⁻
H₂O@ + K⁺ = KOH@ + H⁺
H₂O@ + Mg²⁺ = Mg(OH)⁺ + H⁺
10H⁺ + NO₃⁻ + 8e⁻ = NH₄⁺ + 3H₂O@
H₂O@ + Na⁺ = NaOH@ + H⁺
H₂O@ = OH⁻ + H⁺
10H⁺ + 2SO₄²⁻ + 8e⁻ = S₂O₃²⁻ + 5H₂O@
2H⁺ + SO₄²⁻ + 2e⁻ = SO₃²⁻ + H₂O@
H₂O@ + Sr²⁺ = Sr(OH)⁺ + H⁺
FeO₂⁻ + 4H⁺ = Fe³⁺ + 2H₂O@
CO₃²⁻ + 10H⁺ + 8e⁻ = CH₄@ + 3H₂O@
Ca²⁺ + CO₃²⁻ + H⁺ = Ca(HCO₃)⁺
H₂O@ + Ca²⁺ + SiO₂@ = Ca(HSiO₃)⁺ + H⁺
10H⁺ + SO₄²⁻ + 8e⁻ = H₂S@ + 4H₂O@
CO₃²⁻ + Mg²⁺ = Mg(CO₃)@
CO₃²⁻ + H⁺ + Mg²⁺ = Mg(HCO₃)⁺
H₂O@ + Mg²⁺ + SiO₂@ = Mg(HSiO₃)⁺ + H⁺
CO₃²⁻ + 2H⁺ = CO₂@ + H₂O@
12H⁺ + 2NO₃⁻ + 10e⁻ = N₂@ + 6H₂O@
2H⁺ + 2e⁻ = H₂@
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₃)@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)