Thermochemistry: the quantities and how they connect
This page is the shortest path from "what is a chemical potential" to the expressions this package actually evaluates, with the notation of the code rather than of a textbook. Everything here is either a definition or an identity; nothing is fitted.
1. One potential, and the two pieces it is made of
The Gibbs energy of a system holding nᵢ moles of each species is , and the chemical potential of a species is what one more mole of it costs:
Two pieces, and they come from two different places in the package:
, the standard potential, is thermodynamic data — it belongs to the species alone and comes from the database (§3);
, the activity, is where the mixture enters, and it comes from the activity model.
Divide by and the code appears verbatim. build_potentials returns the closure
Working in is not cosmetic: it is dimensionless, it is what the minimizer differentiates, and it makes and the standard potential directly comparable in size.
2. Activities, standard states, and the conventions in use
An activity is always a ratio to a standard state, and the choice of standard state is a convention that must be stated or the numbers are meaningless. The three in use here:
class
activity
standard state
convention
aqueous solute
mol/kg of solvent
molality
aqueous solute (ideal model)
mol/L
molarity
solvent (water)
pure water
mole fraction
pure solid, pure phase
⇒
the pure substance
—
gas in an ideal mixture
the pure gas at
mole fraction
solid-solution end-member
the pure end-member
mole fraction (§5)
concentration_scale is how a model declares which of the first two it uses, and the aqueous accessors read it rather than guessing:
julia
using ChemistryLabconcentration_scale(DiluteSolutionModel()), concentration_scale(HKFActivityModel())
(:molarity, :molality)
An activity coefficient is therefore not a property of a species but of a species in a solution, and in the standard state's own limit — infinite dilution for a solute, purity for a solid.
The database stores, per species, the coefficients of a heat-capacity polynomial, and the package integrates it. That is the whole of the temperature dependence, and it is exact rather than a series truncation: THERMO_MODELS (src/thermodynamics/thermo_models.jl) holds and the analytic integrals of the same coefficients for , and , so the three identities
hold by construction rather than numerically. Aqueous solutes use a second model, :solute_hkf88_reaktoro, the Helgeson-Kirkham-Flowers equation of state, whose electrostatic part is the Born term evaluated with the dielectric constant of water — the same water model that supplies and to the activity models.
Two consequences worth carrying:
a reaction's temperature dependence is not an input; it follows from the of its reactants and products, which is the Gibbs-Helmholtz relation (van 't Hoff, when is taken constant);
a species whose database entry has no model cannot be moved off at all, and the package says so rather than extrapolating.
4. Equilibrium is a constrained minimization, and its dual is the useful part
At fixed and , equilibrium is the composition that minimizes subject to the matter available. "The matter available" is a linear constraint, because atoms and charge are conserved whatever the reactions do:
where is the conservation matrix cs.SM.A — one row per component, one column per species — and the element (and charge) budget. This is the formulation of (Leal et al., 2017), and it is why no reaction list is needed: reactions are whatever moves inside the null space of .
The Lagrangian multipliers of the equality constraints are the interesting output. Writing for the multiplier of row , the first-order conditions are
The are component potentials (element potentials, when the components are elements): one number per component, from which the chemical potential of any species is a dot product. The two lines above are the complementarity conditions of the KKT system, and they are exactly what optimality_certificate checks — see Proving that an answer is the answer for the certificate itself and why convexity makes it sufficient rather than merely necessary.
5. Reaction quotients, , and the saturation index
Take a species and write its formation from the system's primary species, — the column of is that reaction. For any reaction,
so : the sign of says which way the reaction runs, and is equilibrium. For a solid dissolving into its ions, is the ion activity product IAP and is the solubility product, and the logarithm of their ratio is the saturation index.
Because the rows of are labeled by primary species, a component potential is that primary's , and the index becomes a difference of potentials with no equilibrium constant to look up:
Negative means undersaturated, zero means in equilibrium with the solution, positive means the phase should have precipitated. This is what saturation_indices returns.
Two consequences of the identity are worth carrying, because they are what make the index trustworthy rather than merely conventional:
no equilibrium constant is looked up. never appears in the computation; it is implied by the standard potentials, so an index cannot disagree with the the rest of the calculation used;
every phase present at an equilibrium must come out at exactly zero. That is not a property of the answer, it is the first-order condition of §4 restated, so it is a check on the solve: measured on a CEM I paste, the twelve present solids land within . If they do not, the state is not an equilibrium and no other index in the result means anything.
The identity itself — that the difference of component potentials above equals what saturation_indices computes — is asserted in the test suite (test/aqueous_properties.jl), where it is recomputed by hand from this formula and compared to the function.
Volumes are treated as ideal: the volume of a phase is the sum of its species' standard molar volumes,
with no excess volume of mixing. That is an assumption, and it is the one behind every porosity this package reports. volume returns the split by aggregate state, porosity the void fraction relative to a reference state, and chemical_shrinkage the volume the reaction itself consumes — hydrates occupying less than the water and clinker they were made from, which is the mechanism behind self-desiccation (Self-desiccation).
A species carrying no contributes zero and would corrupt a porosity in silence, which is why the volume machinery reports missing_molar_volumes(state) and why CapillaryWater refuses to be constructed when that list is not empty.
that the answer is a local minimum: the certificate proves global (see below)
that is the whole of the conservation
which phases appear: they are a result
Convexity is what makes the last cell true. is convex in for the built-in activity models, so a KKT point is the minimum and the certificate is sufficient rather than necessary — the one documented exception being CapillaryWater, where a composition-dependent shift of breaks the convexity argument and the certificate falls back to proving a KKT point.