Rate laws, and every parameter in them
Before this page
Proving that an answer is the answer, for the meaning of partial equilibrium on which a kinetic run relies.
A Gibbs minimization answers what is stable. It says nothing about when, and a cement paste is a material whose whole engineering behavior lives in the when. This page is the catalog of what supplies that: every rate law the package ships, the physics each one encodes, the parameters each one carries, their units, and where the numbers come from.
It is written to be read beside the kinetics tutorial rather than instead of it: the tutorial says how to call these, this page says what they mean and what they are worth.
using ChemistryLab
using DynamicQuantities
using Printf
using Plots
default(framestyle = :box, grid = false)0. What a rate law is, here
Every rate in this package is a KineticFunc with one signature:
r = f(T, P, t, n::StateView, lna::StateView, n_initial::StateView) # mol/sSix arguments, and the choice of six is the point. A rate may depend on the temperature and pressure, on the elapsed time, on the current composition, on the activities of any species in the system, and on the initial composition — which is what lets a law be written in terms of a degree of reaction, α = 1 − n/n₀, without the caller having to track it.
Two consequences worth stating plainly:
A rate law may read any species, not only the ones in its own reaction. A dissolution rate catalyzed by
H+, an SCM whose reaction is gated by the portlandite left, a mechanism inhibited by the sulfate in solution — all of these are written by readingnorlnaat the species concerned. The reaction a rate is attached to fixes what the rate consumes; it does not restrict what the rate may depend on. One consequence has its own section in the tutorial, because it bites: a rate that reads a species the reaction itself consumes to exhaustion — see rate laws that depend on a consumed reactant.The evaluation path is differentiable. No
Float64casts, so aForwardDiff.Dualpropagates throughT, throught, through the composition and through the parameters. That is what makes a rate constant fittable (the calibration page) rather than only adjustable by hand.
There are two families below, and they answer different questions. The first is mechanistic: a rate proportional to how far the solution is from equilibrium with the mineral. The second is empirical: a degree of reaction as a function of time, fitted to calorimetry. Clinker and supplementary materials are described by the second, no mechanistic law of the first kind being established for a multiphase material in a paste.
1. Distance from equilibrium: the Palandri–Kharaka form
For a mineral dissolving into or precipitating from a solution, the rate is driven by the saturation ratio
computed by saturation_ratio from the stoichiometry, the log-activities and the standard Gibbs energies the database supplies. Ω < 1 is undersaturation and the mineral dissolves; Ω > 1 and it precipitates; Ω = 1 is equilibrium and the rate is exactly zero — which is the property that makes this form composable with an equilibrium solver rather than in competition with it.
The saturation ratio is the Gibbs energy of the reaction in another form. For the dissolution written with the mineral as reactant, at unit activity, the reaction quotient is the ion activity product, and Reactions and equilibrium constants gives
One mechanism (RateMechanism) contributes, per unit of reactive surface,
and a mineral's rate is the sum over its mechanisms, times its reactive area
the area being evaluated at every step from the current and initial amounts of the mineral by the surface model given to transition_state. The mechanisms are classically an acid one, a neutral one and a base one, which is why the products of activities are called catalysts (RateModelCatalyst): each contributes a_j^{n_j}, with n_j = 0.5 on H+ a typical acid mechanism. A species that slows a mechanism down is an inhibitor (RateModelInhibitor): its factor
| symbol | meaning | unit |
|---|---|---|
| reactive surface area of the mineral | m² | |
rate constant, usually arrhenius_rate_constant | mol m⁻² s⁻¹ | |
| catalyst activity and its exponent | — | |
| inhibitor activity, the inverse of the activity at which its factor is | —, —, — | |
| saturation exponent inside the bracket | — | |
| outer exponent | — | |
| saturation ratio | — |
The rate constant itself is the Arrhenius law referred to a reference temperature rather than to an absolute prefactor,
so that k₀ is the measured constant at T_ref and carries its own unit, instead of being an extrapolation to infinite temperature. The shipped parameter sets for this family come from Palandri and Kharaka (2004).
2. Clinker hydration: Parrott & Killoh
A clinker phase does not dissolve into a solution at a rate set by its own undersaturation — it is consumed behind a growing layer of hydrate, and the controlling step changes as the layer thickens. The 1984 description, as reported by Lothenbach et al. (2008) and used by Lavergne et al. (2018), writes three competing mechanisms in the degree of hydration and lets the slowest one limit:
Each has a physical reading.
the three multipliers being the corrections of section 4.
sets = ["C₃S" => PK84_PARAMS_C3S, "C₂S" => PK84_PARAMS_C2S,
"C₃A" => PK84_PARAMS_C3A, "C₄AF" => PK84_PARAMS_C4AF]
@printf("%-6s %8s %6s %10s %8s %6s %10s\n",
"phase", "k1 [1/d]", "n1", "k2 [1/d]", "k3 [1/d]", "n3", "Ea [kJ/mol]")
for (nm, p) in sets
@printf("%-6s %8.2f %6.2f %10.3f %8.2f %6.1f %10.0f\n",
nm, ustrip(us"1/d", p.k₁), p.n₁, ustrip(us"1/d", p.k₂),
ustrip(us"1/d", p.k₃), p.n₃, ustrip(us"J/mol", p.Ea) / 1000)
endphase k1 [1/d] n1 k2 [1/d] k3 [1/d] n3 Ea [kJ/mol]
C₃S 1.50 0.70 0.050 1.10 3.3 42
C₂S 0.50 1.00 0.006 0.20 5.0 21
C₃A 1.00 0.85 0.040 1.00 3.2 54
C₄AF 0.37 0.70 0.015 0.40 3.7 32Two features of that table are artifacts of the 1984 fit rather than chemistry, and both are useful as checks on an implementation. With
Two Parrott–Killoh variants ship, and only one is attributed
parrott_killoh is a different, smoothed variant that predates this one in the package (the two are set side by side in the tutorial) — min(max(r_{NG}, r_I), r_D) with a damped nucleation term — and its PK_PARAMS_* are not transferable to the canonical law. It is deprecated and the attribution to Parrott and Killoh (1984) was withdrawn rather than repaired: its nucleation term carries no Avrami logarithm, its shell coefficient sits in the diffusion expression, and no published set matches its parameters. The primary source, Parrott and Killoh (1984), is a conference proceedings without a DOI that could not be consulted here, so the attribution was withdrawn. Measured consequence: with PK_PARAMS_* all four phases land on the diffusion branch almost immediately, and a CEM I at w/c 0.40 reaches α ≈ 0.234 at seven days whatever the phase, the signature of a rate set by the diffusion branch alone. Use parrott_killoh_avrami.
The dormant period, which this law does not have
The law starts hydrating at once. Its Avrami branch is seeded at a small degree so that the integration can leave scripts/ionic_hydration.jl multiplies the rate of the two silicates by
A published law with a dormant period is that of Nguyen (2009) (Section 3.3.2): a constant rate
3. Supplementary materials: the Waller sigmoid
A slag or a fly ash has no phases, no formula and no single dissolution mechanism, so its reaction is described by the shape it is observed to have — a sigmoid in log time:
The second form is what the package evaluates, because writing the rate as a function of the current degree rather than of the clock is what lets the same temperature, fineness and humidity corrections multiply it as they multiply the clinker rate.
waller_sets = ["fly ash" => WALLER_PARAMS_FLY_ASH, "silica fume" => WALLER_PARAMS_SILICA_FUME]
α_waller(p, t_days) = 1 / (1 + (ustrip(us"d", p.τ) / t_days)^p.n)
@printf("%-12s %8s %6s %12s %10s %10s\n",
"material", "τ [d]", "n", "Ea [kJ/mol]", "α(28 d)", "α(90 d)")
for (nm, p) in waller_sets
@printf("%-12s %8.0f %6.2f %12.2f %10.3f %10.3f\n",
nm, ustrip(us"d", p.τ), p.n, ustrip(us"J/mol", p.Ea) / 1000,
α_waller(p, 28.0), α_waller(p, 90.0))
endmaterial τ [d] n Ea [kJ/mol] α(28 d) α(90 d)
fly ash 80 0.70 83.14 0.324 0.521
silica fume 80 0.70 83.14 0.324 0.521Pozzolanic and latent-hydraulic reactions are markedly more temperature sensitive than the clinker's: 83.14 kJ/mol against 21–54 kJ/mol above. That is why a blended cement gains so much from a warm cure and loses so much in a cold one, and it is a prediction of the parameter and not an extra rule.
These are not the same numbers as a direct measurement of reacted glass
The α(28 d) column above is ≈ 0.32 for fly ash. The RILEM TC 238-SCM round robin (Durdziński et al., 2017), which measured the reacted glass directly on its own materials rather than inferring it from heat, reports about 20 % for a siliceous fly ash at the same age, and 38–49 % for two slags. No slag set is shipped: the 100-day time this package once attributed to Waller (1999) is not in the thesis, which names a slag only among the additions an adiabatic test may contain. A slag's time is the caller's, from a source, merge(WALLER_PARAMS_FLY_ASH, (τ = τ_slag,)).
Neither is wrong. The parameters of Waller (1999) are a fit to particular materials, and "a slag" is not a substance — its reactivity depends on its glass content, its basicity and its fineness. The conclusion of the round robin (Durdziński et al., 2017) is the one to keep: the precision of any determination of an SCM's degree of reaction is "rather low, at best ± 4-5 %". So a reacted fraction is an input to be stated and swept, never a constant to be trusted to two digits — which is exactly how the blended binder pages treat it.
Silica fume carries the same τ and n as fly ash: its much higher reactivity is represented through the fineness, at an effective Blaine of 2000 m²/kg recommended by Lavergne et al. (2018). Its BET surface, about 20 000 m²/kg, is a different measurement of a different thing and must not be substituted — a factor of ten on the rate.
That is now enforced rather than advised. BETSurfaceArea and BlaineSurfaceArea are distinct types, and the ratio that produces a fineness factor is defined only between two of the same kind, so handing a BET area to blaine_factor raises instead of returning a plausible number. A bare quantity is still read as a Blaine fineness, which is the contract every existing call relies on; see Surface areas.
The laws of a blend that read the portlandite, and why they are not here
The Waller sigmoid ignores the clinker it reacts with. Two laws written for a blend do not: they make the rate of the addition a product of factors, one of which reads the portlandite. For a fly ash, Buffo-Lacarrière et al. (2007) multiply it by the volume fraction of portlandite in the paste (their Eq. 12); for a slag, Kolani et al. (2012) by that fraction plus the volume fraction of the lime still held in the unreacted slag (their Eq. 13). The other factors of the slag law are an activation proportional to the degree of hydration of the clinker, taken as the measure of the alkalis it has released (Eq. 10), an access of water through the hydrate layer,
with
Neither law is implemented here: three quantities a rate law needs are not among the published values. The constant
A filler effect on the clinker has no published formulation to implement either. Berodier and Scrivener (2014) relate the slope of the acceleration period of a cement paste to the distance between the particles, and find that a slag and a fly ash act as quartz does over the first day; the relation is shown as measured points, not as an expression a rate law can take.
4. The three multiplicative corrections
Each of the empirical laws above is multiplied by three dimensionless factors. They are separated because they answer separate questions, and because each can be switched off by leaving its keyword at nothing.
Temperature, the same Arrhenius factor as section 1:
Fineness (blaine_factor): the rate scales linearly with the specific surface, relative to the fineness the parameters were fitted at — 385 m²/kg for the clinker phases, 400 m²/kg for the Waller sets.
Internal humidity (humidity_factor): hydration stops when the pore water is no longer available, which is taken to happen below 80 % RH.
@printf("β_B at 462 m²/kg (clinker ref) : %.3f\n", blaine_factor(462u"m^2/kg"))
@printf("β_B at 2000 m²/kg (Waller ref) : %.3f\n",
blaine_factor(2000u"m^2/kg"; blaine_ref = 400u"m^2/kg"))
for h in (0.99, 0.90, 0.81, 0.801, 0.80, 0.75)
@printf("β_h(%.3f) = %.4f\n", h, humidity_factor(h))
endβ_B at 462 m²/kg (clinker ref) : 1.200
β_B at 2000 m²/kg (Waller ref) : 5.000
β_h(0.990) = 0.9140
β_h(0.900) = 0.3660
β_h(0.810) = 0.1114
β_h(0.801) = 0.0968
β_h(0.800) = 0.0000
β_h(0.750) = 0.0000That cut is a genuine discontinuity, and the page says so rather than smoothing it: the one-sided limit from above is 0.0953 and the value at 0.80 is exactly 0. It is mild in practice — the rate has already fallen by an order of magnitude from β_h(0.99) = 0.914 — but a solver stepping across it will feel it, which is why PoreHumidity exists to supply h from the current saturation rather than from a schedule.
When the area stops being constant
In every law shipped here
A grain that dissolves does not keep its area. ShrinkingCoreArea says so, and passing one where a fineness is expected makes the factor follow the amount left,
with
And that freedom is not new everywhere. Multiply the shell-formation branch by
so wherever that branch is the active one,
Which of these holds over a given dataset is a measurement and not an argument. identifiability is what makes it, and Where the numbers come from is where a parameter that the data turn out not to constrain gets labeled as such instead of quoted.
An evolving area needs its own calibration
The published constants were fitted with scripts/hydration_calibration.jl is the machinery — not to keep the old constants and add a new factor on top.
5. How far the reaction can go: α_max
Every law above is written in powers_alpha_max:
@printf("%6s %12s %12s\n", "w/c", "sealed", "under water")
for wc in (0.25, 0.30, 0.36, 0.40, 0.42, 0.50)
@printf("%6.2f %12.3f %12.3f\n", wc,
powers_alpha_max(wc), powers_alpha_max(wc; curing = :saturated))
end w/c sealed under water
0.25 0.595 0.694
0.30 0.714 0.833
0.36 0.857 1.000
0.40 0.952 1.000
0.42 1.000 1.000
0.50 1.000 1.000Read the last two rows: above 0.42 neither convention does anything, because water has stopped being what limits the reaction. Read the first: at w/c = 0.25 a sealed paste cannot pass 60 % hydration however long it is left, and a cured one cannot pass 69 %. The difference between the columns is the chemical shrinkage, refilled from the bath in one case and not in the other.
6. Which numbers are published and which are fitted
| quantity | where it comes from |
|---|---|
PK84_PARAMS_* | Parrott and Killoh (1984) as reported by Lothenbach et al. (2008), Lavergne et al. (2018) |
WALLER_PARAMS_* | Waller (1999), as used by Lavergne et al. (2018) |
| the dissolution rate constants and their exponents | Palandri and Kharaka (2004) |
blaine_ref 385 / 400 m²/kg | the finenesses those fits were made at |
| the 0.80 humidity cut and its exponent | (Lavergne et al., 2018; Eq. 10), who take it from earlier work on Portland cement |
powers_alpha_max 0.42 / 0.36 | Powers (1948) |
CALIBRATED_THETA | fitted here, on one record — the calibration page |
| an SCM's reacted fraction at a given age | an input, measured or assumed — never a constant of the code |
The last two rows are the ones to keep in view. Everything above them is somebody's published fit to somebody's materials; the calibration page shows what happens when those published parameters meet a calorimetry record they were not fitted to, and how much of the gap a five-parameter fit closes.
7. What a rate law may be
Section 0 says that a rate law may read anything. This section says what it must do with what it reads, so that a run does not leave the side of thermodynamics: where the rate vanishes, which sign it has, and how two laws closing a cycle must agree.
One reaction
For a single reaction the second law makes the entropy it produces non-negative,
Any law of the form
meets both requirements whatever transition_state are therefore parameters of a law, not consequences of thermodynamics.
What is not a parameter is where the rate vanishes. Written as a forward rate minus a backward one, a law that must vanish where the reaction quotient equals the equilibrium constant has the ratio of its two rate constants fixed by that constant (Lasaga, 1984; Eqs. (A6)–(A7), p. 4024). The constant is the one of the database the equilibrium solver uses, so that the kinetic reaction and the equilibrium partition agree on the state they lead to. A backward rate chosen on its own stops the reaction at some transition_state, first_order_rate and sorption_rate, compute
The laws that do not read the composition
The laws of sections 2 and 3, Parrott–Killoh and Waller, are functions of a degree of reaction and of time. They hold no
A law that runs one way
A phase that dissolves and never reprecipitates, or grows on its seeds and is never redissolved, has a rate that vanishes on one side of saturation while the affinity does not. That is allowed (the converse noted above), and dissolution_only and growth_only build it from any law:
The factor FixedSurfaceArea for instance; SHRINK_FLOOR. A restricted law keeps the sign of the affinity where it runs, so dissipation never lists it.
A phase that grows from its own seeds
A phase that precipitates kinetically is written as two species: a reservoir of seeds, a copy of the phase (with_symbol) that a nucleation rate fills, and the phase itself, growing on the surface those seeds offer. The number of seeds is the volume of the reservoir over the volume of one seed, SphericalSeeds), it grows on
the area of the seeds themselves before anything has grown on them. As needles of a fixed aspect ratio (NeedleSeeds), the length is SeededGrowthArea is that area, as the surface of a transition_state law, usually restricted to growth.
The reservoir fills at the frequency of nucleation
the classical form written through its value ClassicalNucleation), every derivative of which vanishes at saturation, or a frequency proportional to the excess over a critical ratio (ThresholdNucleation). nucleation_rate builds that law. What these tools do not do is share the substrate between generations of nuclei or let the seeds of one cohort stop nucleating others: those remain choices of a model, written in its own rate laws.
A layer that covers a surface
A thin layer that forms on a mineral while its phase is stable blocks part of the mineral's sites, and opens again, from the sites already free, when that phase is no longer stable (CoveringLayer). Counted in the sites it blocks, covered multiplies into the rate of the substrate, and it follows (covering_rate)
with
Several reactions
For several reactions the second law bounds the sum
Written as
whatever the composition. A rest state, every
nothing in total, so that unless every term vanishes at least one reaction runs against its affinity.
Checking a run
reaction_affinities gives speciated_states) and with the activities the rate laws read. dissipation adds the rate of each law evaluated there, the power test/kinetics/test_sorption_kinetics.jl).
See also
The water budget of a hydrating paste — where
comes from, and why it is not thermodynamicsCalibrating a hydration model — the fit, the holdout, and the sensitivity
A CEM I from its clinker phases — all of this running, over ninety days