Skip to content

Proving that an answer is the answer

An interior-point method minimizes by walking the interior of the feasible set, and on a cement equilibrium it stops on MaxIters — at any tolerance. Whether the point it returns is the minimum is then an open question. This page is why the question has an answer at all, what a certificate checks, and how a solver can aim at the optimality conditions directly instead of at the objective.

The formulation being certified is the one set out in Thermochemistry §4: minimize subject to   and  .

Why the question has an answer

Write the Gibbs energy in RT units as G(n) = Σᵢ nᵢ μᵢ(n). Its ideal part

has Hessian  , and for any

by Cauchy–Schwarz applied to  . A pure phase has unit activity, so it contributes a term linear in its amount. Hence G is convex, the feasible set    is a polyhedron, and — the constraints being affine, so that the linearity constraint qualification holds everywhere — the KKT conditions are necessary and sufficient.

Two consequences follow. The minimizer is unique, so a solver returning different answers from different starting points is not finding local minima but stopping short of stationarity. And optimality can be checked: a composition satisfying the KKT conditions is proved globally optimal.

The certificate

optimality_certificate checks the three conditions, on any composition and whatever produced it. Writing    for the element potentials:

conditionon which speciesmeaning
  interior (n > floor)stationarity
 conservation of matter
 a pure phase at its boundthat phase undersaturated
   a mixing phase held entirely absentthat solution cannot form

The last row is Michelsen's tangent-plane measure, and it is a separate test because a mixing phase needs one: its members are never exactly zero while it exists, so they are neither interior nor at a bound, and a solid solution left out of the assemblage used to pass the certificate unexamined. The trial composition is refined against the phase's own activity model, so the test is not the ideal approximation.

Two subtleties decide whether the check is meaningful.

A species at its bound obeys the inequality, not the equality. Imposing the equality on an amount held at 1e-16 whose mass-action value is e⁻³⁰⁰ misstates its log-activity by 263 RT units, and the check then reports a residual of 74 for a composition solved to 5e-12.

A species carrying a vanished component is absent by the constraint, not by thermodynamics, and its saturation index is meaningless — the element potential of a component nobody supplies is determined by nothing. The test for that is not bₖ ≈ 0 but bₖ ≈ 0 with the non-zero entries of row k sharing a sign: only then does   with   force each term to vanish. The H⁺ row carries +1 for H⁺ and −1 for OH⁻, so its zero total is the ordinary state of pure water; treating it as degenerate kills the entire acid–base system and returns pH 7.000 with the calcite undissolved.

The certifying solver

DualEquilibriumSolver solves the KKT system directly, in element potentials. From    an aqueous species obeys the mass-action law   , and a pure phase is present exactly when  , absent when undersaturated — the classical phase-stability criterion.

Two levels. The inner one inverts the solutes' mass-action laws at fixed potentials and fixed solvent amount; the outer is a Newton on    unknowns — the solvent, the m element potentials, and the amounts of the active phases. Parameterizing the solutes by makes their positivity automatic, which is what removes the fraction-to-boundary limit that caps the interior-point step at every iteration.

The solvent is deliberately not inverted through its own mass-action law: its activity is a mole fraction, so   always, and an arbitrary y can demand more, for which no finite composition exists. It belongs to the outer system, where the balance determines it.

What it buys, measured

On calcite in pure water the certified pH is 9.90 against an interior-point 6.96 — not an imprecision but a wrong answer, and one nothing in that solver's output reveals. On the Reaktoro reference (calcite, CO₂ and water) both routes now agree with Reaktoro on every species: above 10⁻⁵ mol to 10⁻³ relative, the trace ions to 5 %, the worst being CaOH⁺ at ×1.032 on 1.6 nmol. That reference used to carry a @test_broken for CaOH⁺ at ×2.47; what closed it was the convergence test moving to the true KKT error at μ = 0 (OptimaSolver 0.4.1), and test/equilibrium_reference.jl is now 26 plain assertions.

speciated_states certifies every instant it replays and names any it cannot. On a full ordinary Portland cement over 28 days, all forty replayed instants are certified, with element balances between 1e-11 and 1e-13 mol.

See also