Getting started
This quickstart shows a few common, minimal examples to get you productive with ChemistryLab. It demonstrates loading species from a database, building reactions and solving a thermodynamic equilibrium problem.
Simplified example
Let us start with a minimal example in which we compute the thermodynamic properties of a reaction. As a first illustration, we consider the equilibrium of calcite in water. This equilibrium can be written as:
\[\text{CaCO}_3 \rightleftharpoons \text{Ca}^{2+} + {\text{CO}_3}^{2-}\]
It is possible to calculate the thermodynamic properties of the reaction, in particular the solubility constant of the reaction ($\ln K$) which is related to the Gibbs free energy of the reaction ($\Delta_r G^°$). This solubility constant is a function of temperature and the calculation is performed at a reference temperature of 298 K and at a pressure of 1 Atm using the following equation:
\[\Delta_r G^° = - RT\ln K\]
where $\Delta_r G^°$ is deduced from the Gibbs energies of formation ($\Delta_f {G_i}^°$) of the other chemical species involved in the reaction:
\[\Delta_r G^° = \sum_i \nu_i \Delta_f {G_i}^°\]
To do this, we load the species from one of the databases integrated into ChemistryLab, filter those relevant to the calcite–water system, then build the stoichiometric matrix and derive the reactions.
In this example, the database is cemdata. The .json file is included in ChemistryLab but is a copy of a file which can be found in ThermoHub.
build_species reads the database file and returns a Vector{Species} with compiled thermodynamic functions. speciation then filters this list to the species whose atomic composition is a subset of the seed atoms (here Ca, C, H and O from Cal, H2O@ and CO2):
using ChemistryLab
using DynamicQuantities # for unit management
all_species = build_species("../../data/cemdata18-thermofun.json")
species_calcite = speciation(all_species, split("Cal H2O@ CO2");
aggregate_state=[AS_AQUEOUS],
exclude_species=split("H2@ O2@ CH4@"))
dict_species_calcite = Dict(symbol(s) => s for s in species_calcite)Dict{String, Species{Int64}} with 12 entries:
"H+" => H+ {H+} [H+ ◆ H⁺]
"OH-" => OH- {OH-} [OH- ◆ OH⁻]
"Ca(HCO3)+" => Ca(HCO3)+ {CaHCO3+} [Ca(HCO3)+ ◆ Ca(HCO₃)⁺]
"Cal" => Cal {Calcite} [CaCO3 ◆ CaCO₃]
"Ca+2" => Ca+2 {Ca+2} [Ca+2 ◆ Ca²⁺]
"CO2@" => CO2@ {CO2 aq} [CO2@ ◆ CO₂@]
"HCO3-" => HCO3- {HCO3-} [HCO3- ◆ HCO₃⁻]
"CO2" => CO2 {CO2 g} [CO2 ◆ CO₂]
"CaOH+" => CaOH+ {CaOH+} [Ca(OH)+ ◆ Ca(OH)⁺]
"H2O@" => H2O@ {H2O l} [H2O@ ◆ H₂O@]
"CO3-2" => CO3-2 {CO3-2} [CO3-2 ◆ CO₃²⁻]
"Ca(CO3)@" => Ca(CO3)@ {CaCO3 aq} [CaCO3@ ◆ CaCO₃@]During species creation, ChemistryLab calculates the molar mass of the species. It also constructs thermodynamic functions (heat capacity, entropy, enthalpy, and Gibbs free energy of formation) as a function of temperature. The evolution of thermodynamic properties as a function of temperature, such as heat capacity, can thus be easily plotted.
dict_species_calcite["Cal"]Species{Int64}
name: Calcite
symbol: Cal
formula: CaCO3 ◆ CaCO₃
atoms: Ca => 1, C => 1, O => 3
charge: 0
aggregate_state: AS_CRYSTAL
class: SC_COMPONENT
properties: M = 0.10008599996541243 kg mol⁻¹
Tref = 298.15 K
Pref = 100000.0 m⁻¹ kg s⁻²
Cp⁰ = 104.5163192749 + 0.02192415855825T + -2.59408e6 / (T^2) [m² kg s⁻² K⁻¹ mol⁻¹] ◆ vars=(T) ◆ T=298.15 K
ΔₐH⁰ = -1.2482415842895252e6 + 104.5163192749T + 2.59408e6 / T + 0.010962079279125(T^2) [m² kg s⁻² mol⁻¹] ◆ vars=(T) ◆ T=298.15 K
S⁰ = -523.9438829693111 + 0.02192415855825T + 104.5163192749log(T) + 1.29704e6 / (T^2) [m² kg s⁻² K⁻¹ mol⁻¹] ◆ vars=(T) ◆ T=298.15 K
ΔₐG⁰ = -1.1423813547027335e6 + 628.460202244211T + 1.29704e6 / T - 0.010962079279125(T^2) - 104.5163192749T*log(T) [m² kg s⁻² mol⁻¹] ◆ vars=(T) ◆ T=298.15 K
V⁰ = 3.6933999061584004e-5 [m³ mol⁻¹]
Cp⁰_Tref = 81.87109375 m² kg s⁻² K⁻¹ mol⁻¹
ΔₐH⁰_Tref = -1.207405e6 m² kg s⁻² mol⁻¹
S⁰_Tref = 92.675598144531 m² kg s⁻² K⁻¹ mol⁻¹
ΔₐG⁰_Tref = -1.129176e6 m² kg s⁻² mol⁻¹
V⁰_Tref = 3.6933999061584004e-5 m³ mol⁻¹using Plots
p1 = plot(xlabel="Temperature [°C]", ylabel="Cp⁰ [J/mol/K]", title="Heat capacity of calcite \nas a function of temperature")
plot!(p1, θ -> dict_species_calcite["Cal"].Cp⁰(T = θ*ua"degC"), 0:0.1:100, label="Cp⁰")
Obtaining stoichiometric matrices requires the choice of a species-independent basis.
primaries = [dict_species_calcite[s] for s in split("H2O@ H+ CO3-2 Ca+2")]
SM = StoichMatrix(collect(values(dict_species_calcite)), primaries)
pprint(SM)┌───────┬────┬─────┬───────────┬─────┬──────┬──────┬───────┬─────┬───────┬──────┬───────┬──────────┐
│ │ H+ │ OH- │ Ca(HCO3)+ │ Cal │ Ca+2 │ CO2@ │ HCO3- │ CO2 │ CaOH+ │ H2O@ │ CO3-2 │ Ca(CO3)@ │
├───────┼────┼─────┼───────────┼─────┼──────┼──────┼───────┼─────┼───────┼──────┼───────┼──────────┤
│ H2O@ │ │ 1 │ │ │ │ -1 │ │ -1 │ 1 │ 1 │ │ │
│ H+ │ 1 │ -1 │ 1 │ │ │ 2 │ 1 │ 2 │ -1 │ │ │ │
│ CO3-2 │ │ │ 1 │ 1 │ │ 1 │ 1 │ 1 │ │ │ 1 │ 1 │
│ Ca+2 │ │ │ 1 │ 1 │ 1 │ │ │ │ 1 │ │ │ 1 │
└───────┴────┴─────┴───────────┴─────┴──────┴──────┴───────┴─────┴───────┴──────┴───────┴──────────┘These stoichiometric matrices thus allow us to write the chemical reactions at work.
list_reactions = reactions(SM)
dict_reactions_calcite = Dict(r.symbol => r for r in list_reactions)Dict{String, Reaction{Species{Int64}, Int64, Species{Int64}, Int64, Int64}} with 8 entries:
"Ca(CO3)@" => CO₃²⁻ + Ca²⁺ = CaCO₃@
"Cal" => CO₃²⁻ + Ca²⁺ = CaCO₃
"OH-" => H₂O@ = OH⁻ + H⁺
"CaOH+" => H₂O@ + Ca²⁺ = Ca(OH)⁺ + H⁺
"CO2@" => 2H⁺ + CO₃²⁻ = CO₂@ + H₂O@
"HCO3-" => H⁺ + CO₃²⁻ = HCO₃⁻
"CO2" => 2H⁺ + CO₃²⁻ = CO₂ + H₂O@
"Ca(HCO3)+" => H⁺ + CO₃²⁻ + Ca²⁺ = Ca(HCO₃)⁺Again, when constructing the reactions, the thermodynamic properties of the reactions as a function of temperature are deduced. It is thus possible to see, for example, the expression for the solubility product of calcite for the reaction under study and to plot its evolution.
dict_reactions_calcite["Cal"].logK⁰NumericFunc:
Unit: [m² kg s⁻² mol⁻¹]
Variables: T, P
References: T=298.15 K, P=100000.0 m⁻¹ kg s⁻²p2 = plot(xlabel="Temperature [°C]", ylabel="pKs", title="Solubility product (pKs) of calcite \nas a function of temperature")
plot!(p2, θ -> dict_reactions_calcite["Cal"].logK⁰(T = θ*ua"degC"), 0:0.1:100, label="pKs")
Equilibrium solving
The previous section computed thermodynamic properties of reactions analytically. ChemistryLab can go further and solve the full thermodynamic equilibrium — that is, find the species amounts that minimise the Gibbs free energy of the system given initial conditions.
Three objects are needed:
| Object | Role |
|---|---|
ChemicalSystem | Immutable description of the system: species list, primary species, stoichiometric matrices, and derived index maps. Built once and reused. |
ChemicalState | Mutable thermodynamic state: amounts n (mol), temperature T and pressure P. Modified in-place before and after solving. |
equilibrate | Convenience function that wraps a ChemicalSystem + ChemicalState into an optimisation problem and solves it. Returns a new equilibrated ChemicalState. |
using DynamicQuantities
# ChemicalSystem: declare the species and which four are the independent basis
primaries_eq = [dict_species_calcite[s] for s in split("H2O@ H+ CO3-2 Ca+2")]
cs = ChemicalSystem(collect(values(dict_species_calcite)), primaries_eq)
# ChemicalState: set the initial amounts
state = ChemicalState(cs)
set_quantity!(state, "Cal", 1e-3u"mol") # 1 mmol of calcite
set_quantity!(state, "H2O@", 1.0u"kg") # 1 kg of water
# Seed H⁺ and OH⁻ at pH 4 (trace amounts to help convergence)
V = volume(state)
set_quantity!(state, "H+", 1e-4u"mol/L" * V.liquid)
set_quantity!(state, "OH-", 1e-10u"mol/L" * V.liquid)
# Solve: find the Gibbs-energy minimum
state_eq = equilibrate(state)ChemicalState{Species{Int64}, AbstractReaction, DynamicQuantities.Quantity{Float64, DynamicQuantities.Dimensions{DynamicQuantities.FRInt32}}}
┌─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┐
│ T : 298.15 K │
│ P : 1.0 bar │
╞═════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╡
│ # liquid #│ n [mol]│ m [g]│ V [cm³]│ c [mol/L]│ │
├┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┤
│ tot. liquid│ 55.5096│ 1000.02│ 1002.96│ │ │
├┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┤
│ H2O@│ 55.5093│ 999.999│ 1002.96│ 55.3453│ │
│ Ca+2│ 0.000155243│ 0.00622184│ -0.00286249│ 0.000154785│ │
│ HCO3-│ 0.000133904│ 0.00817031│ 0.00324194│ 0.000133509│ │
│ OH-│ 3.39749e-5│ 0.000577811│ -0.000159948│ 3.38745e-5│ │
│ CO3-2│ 2.13374e-5│ 0.00128041│ -0.000129256│ 2.12743e-5│ │
│ Ca(CO3)@│ 5.55362e-6│ 0.00055584│ -8.69055e-5│ 5.53721e-6│ │
│ Ca(HCO3)+│ 2.69587e-7│ 2.72536e-5│ 3.59354e-6│ 2.6879e-7│ │
│ CO2@│ 9.30748e-8│ 4.09613e-6│ 3.05347e-6│ 9.27998e-8│ │
│ CaOH+│ 9.19469e-8│ 5.24879e-6│ 5.29841e-7│ 9.16752e-8│ │
│ H+│ 2.23643e-9│ 2.25432e-9│ 0.0│ 2.22982e-9│ │
╞═════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╡
│ # solid #│ n [mol]│ m [g]│ V [cm³]│ │ │
├┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┤
│ tot. solid│ 0.000838842│ 0.0839563│ 0.0309818│ │ │
├┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┤
│ Cal│ 0.000838842│ 0.0839563│ 0.0309818│ │ │
╞═════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╡
│ # gas #│ n [mol]│ m [g]│ V [cm³]│ │ p [bar]│
├┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┤
│ tot. gas│ 3.51884e-10│ 1.54861e-8│ 8.7231e-6│ │ │
├┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┤
│ CO2│ 3.51884e-10│ 1.54861e-8│ 8.7231e-6│ │ 1.0│
╞═════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╡
│ # TOTAL #│ n [mol]│ m [g]│ V [cm³]│ │ │
├┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┤
│ │ 55.5105│ 1000.1│ 1002.99│ │ │
╞═════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╡
│ pH : 9.53 │
│ pOH : 4.4701 │
│ porosity : 0.999969 │
│ saturation : 1.0 │
└─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘
println("pH = ", round(pH(state_eq), digits = 2))pH = 9.53Derived quantities such as pH, pOH, phase volumes and individual species amounts are all accessible on the returned ChemicalState. For a detailed description of the solver options, activity models, and temperature sweeps, see the equilibrium tutorial (Tutorial → Chemical Equilibrium).
Notes and next steps
- The
Formula,Species,ReactionandStoichMatrixAPIs are intentionally small and composable — explore thedocs/src/pages for detailed examples. - For equilibrium calculations, see
docs/src/man/equilibrium.mdand the worked examplesco2_carbonate_systemandcement_carbonation. - For cement-specific workflows, use
CemSpeciesand thedatabasesutilities to convert between oxide- and atom-based representations.
Now try the quickstart examples interactively in the REPL and then follow the next pages of the tutorial for deeper coverage.