Verified Math: Exact Decimals, Certified Series, and Declared Laws

By Fabian Schindler (@Cypoe) — this whole lane is his work, ported from his own isar-proofs (Lean 4) and plex lineages.

Koru has numbers now — not machine floats with better marketing, but a surface where “number” means something a proof can hold onto. Four modules carry it: std/num (arbitrary-precision decimals, exact or externally tainted), std/quantity (SI dimensions and uncertainty with covariance), std/sym (canonical symbolic expressions), and the certificate lane on top of it — whose semantics are machine-checked in Lean 4 (Holonomic.lean, HolonomicClosure.lean, HolonomicCompose.lean, QuantityKernel.lean, and friends). The newest surface, sym.rel, is where relationships between abstract generators live.

std/num — exact means exact

A std/num value is an arbitrary-precision decimal — mantissa as a big integer plus a base-10 exponent — and it knows whether it is exact or external. A literal is exact; text ingested from the outside world is tainted until ingest discharges it, and arithmetic on tainted values refuses rather than laundering provenance:

Worth saying plainly what “decimal” means here, because the word is doing dishonest work elsewhere. .NET’s decimal — the type everyone is told to use for money — is a 128-bit floating decimal: 28–29 significant digits, then it rounds. It fixes the binary-representation lie (0.1 is representable), not the precision lie; “recommended for money” means the losses are small, not that they’re zero. A std/num addition, subtraction, multiplication, or power is exact with no ceiling — there is no bit-width where rounding silently begins. Where a decimal expansion cannot terminate — division, square root — the callsite has to name the precision: std/num:div(a, b, prec: 12). The rounding lives in your code where a reviewer can see it, not inside the type where you can’t. For money that difference is not philosophical; it is the difference between an audit trail and a shrug.

import std/io
import std/num

std/num:lit(123456789012345678901234567890)
| exact a |> std/num:lit(987654321098765432109876543210)
    | exact b |> std/num:add(a, b): r |> std/io:print.ln("{{ r:s }}") |> std/num:free(n: r)
    | bad x |> std/io:print.ln("bad")
| bad y |> std/io:print.ln("bad")

The operation surface is the one a calculator wants — add/sub/mul/div, pow, sqrt, exp, ln, sin/cos/atan, powf, mod, base conversion — with fixed decimal precision requested per call, never machine precision smuggled in. The multiplication underneath is Koru’s own limb-level Karatsuba (Karatsuba & Ofman, 1962): the race against bc and Python’s decimal (tests/benchmarks/016_num_precision_race, replayable) showed Zig’s stock schoolbook mul was the bottleneck, so the limb layer got owned — transcendentals dropped roughly an order of magnitude, and ReleaseFast beats bc on eight of ten workloads, including a Machin-formula π (1706, still the fast classroom recipe). The C bigint implementations still win the pure exact-op races; the honest remaining gap is there in the results table, not smoothed over.

std/quantity — dimensions are a group, not a string

std/quantity packs a magnitude, a canonical Dim, and optional uncertainty-with-covariance into one linear string. Dim is [7]i32 — the free abelian group on the seven SI base generators, the same object QuantityKernel.lean proves things about — so m/s · s → m is group law, not string surgery, and a unit mismatch on add refuses dim-mismatch while equiv still compares canonical text. Uncertainty propagates through the covariance-carrying lanes; lanes the sigma math can’t honestly cover refuse sym-unc rather than approximate silently.

The holonomic lane — certificates, not witnesses

std/sym gives canonical symbolic expressions; the certificate parts on top are where the semantics get heavy. A holonomic function is one annihilated by a linear ODE with polynomial coefficients — Zeilberger’s name (1990) for what Stanley called differentiably finite (1980); the friendly modern reference is Kauers & Paule’s The Concrete Tetrahedron (2011). The annihilating ODE is a certificate: finite data that proves a property of an infinite object.

• cert-seed(body) registers a declared law for an opaque generator — p0=[-1]|p1=[1] is f' − f = 0, exp’s law, given as truth.
• derive-ode(expr, max-order, max-deg) discovers a certificate by elimination over the truncated series — the “guessing” half of gfun’s guess-and-prove (Salvy & Zimmermann, 1994) — witness-grade, verified=N honestly.
• cert-sum/cert-prod compose certificates through companion-basis elimination over ℚ(x) — the standard linear-algebra proof of the classical closure theorems, so the output is theorem-grade; HolonomicClosure.lean checks the same statement. cert-compose refuses by name (not-generally-holonomic) because composition of D-finite functions isn’t, in general.
• cert-bridge projects an ODE certificate into its P-recursive coefficient recurrence — the D-finite ⇔ P-recursive equivalence (Stanley 1980; Kauers & Paule §7), implemented in sym.certcomp with its Lean-side proof pending in isar-proofs; det=1 flags exactly the hypothesis that makes the recurrence determine the sequence.
• cert-of/cert-series are the theorem upgrade: derive the ODE structurally, then let the recurrence generate the tail — every coefficient forced by the annihilating equation rather than merely witnessed.
• decide-integral/decide-ode/decide-field are the decision tower — polynomial integrands get constructed antiderivatives, holonomic integrands get shifted certificates (antiderivatives of D-finite functions stay D-finite — a closure corollary, same sources). decide-ode applies Kovacic’s pole-order necessary conditions (Kovacic, J. Symbolic Computation 1986) to classify Liouvillian solvability, so Airy’s equation comes back non-Liouvillian — correct, and Risch’s integration-in-finite-terms framing (1969) is why that answer is honest rather than a shrug. What the tower can’t reach says unresolved with its search bounds named.

std/sym:lit(inv((1e0+(-1e0*x))))
| exact a |> std/sym:decide-integral(a, max-order: 3, max-deg: 3)
    | elementary m |> std/io:print.ln("{{ m:s }}") |> std/io:free(content: m)
    | holonomic-special m |> std/io:print.ln("hs:{{ m:s }}") |> std/io:free(content: m)
    | unresolved m |> std/io:print.ln("un:{{ m:s }}") |> std/io:free(content: m)
    | not-holonomic m |> std/io:print.ln("nh:{{ m:s }}") |> std/io:free(content: m)
    | non-rational |> std/io:print.ln("nr")
    | non-analytic |> std/io:print.ln("na")
    | bad |> std/io:print.ln("bad")
| bad |> std/io:print.ln("bad-lit")

sym.rel — relationships declared, not discovered

The newest band generalizes “ODE law on a generator” to declared relationships in general. A rel:v1 registry names generators and attaches derivative laws in a second phase — rel-seed declares the name, rel-law attaches D(g) as an expression over x and other generators, so cyclic pairs resolve:

import std/io
import std/sym

// condensed from 615_049: cyclic laws, then the conservation check.
std/sym:rel-new(): r0
|> std/sym:rel-seed(name: sin, reg: r0)
    | exact r1 |> std/sym:rel-seed(name: cos, reg: r1)
        | exact r2 |> std/sym:rel-law(name: sin, reg: r2) { g(cos) }
            | exact r3 |> std/sym:rel-law(name: cos, reg: r3) { -1e0*g(sin) }
                | exact r4 |> std/sym:rel-derive(reg: r4) { g(sin)*g(sin)+g(cos)*g(cos) }
                    | exact d |> std/io:print.ln("{{ d:s }}") |> std/io:free(content: d)
                        |> std/io:free(content: r4)
                    | unknown-generator |> std/io:print.ln("unknown-generator") |> std/io:free(content: r4)
                    | lawless |> std/io:print.ln("lawless") |> std/io:free(content: r4)
                    | bad |> std/io:print.ln("bad") |> std/io:free(content: r4)
                | unknown-generator |> std/io:print.ln("unknown-generator")
                | bad |> std/io:print.ln("bad")
            | unknown-generator |> std/io:print.ln("unknown-generator")
            | bad |> std/io:print.ln("bad")
        | dup |> std/io:print.ln("dup")
        | full |> std/io:print.ln("full")
        | bad |> std/io:print.ln("bad")
    | dup |> std/io:print.ln("dup")
    | full |> std/io:print.ln("full")
    | bad |> std/io:print.ln("bad")

rel-derive prints 0e0 — the conservation of sin² + cos² proven mechanically, since every rewrite step is a derivative rule or a declared-law substitution. rel-norm canonicalizes for caller-side equality checks. Generator references form their own expression universe — g(name) parses only under a registry, so <sexpr!> values can’t leak in and rel text can’t leak out into sym:add; undeclared names refuse unknown-generator, deriving a lawless generator refuses lawless.

What refuses, on purpose

The surface’s shape is defined by its edges as much as its operations. cert-compose refuses because composition isn’t generally holonomic. decide-ode says non-Liouvillian rather than pretending. cert-series on x says under-determined — its recurrence genuinely leaves a seed free. Tainted operands, undeclared axes, lawless generators, dim mismatches, malformed external certificates — each is a named refusal, not a silent path. Discovery is witness-grade; declaration and derivation are theorem-grade; and the line between them is drawn where the Lean proofs draw it.

Open rungs, honestly: g(name, arg) application (composition and addition theorems like exp(x+y)=exp(x)exp(y)) needs substitution machinery — band 2. Nonlinear jet relations are honest-refusal territory. And the exact-op multiplication gap versus the C bigint implementations narrows, not closes, with Karatsuba.

Where the math comes from

• Karatsuba & Ofman, “Multiplication of Multidigit Numbers on Automata” (1962) — the limb multiplication.
• Machin (1706) — π = 16·atan(1/5) − 4·atan(1/239), the benchmark’s π workload.
• Stanley, “Differentiably Finite Power Series” (Europ. J. Combinatorics, 1980; also Enumerative Combinatorics vol. 2, §6.4) — the D-finite theory, the closure theorems, the P-recursive equivalence.
• Zeilberger, “A Holonomic Systems Approach to Special Functions Identities” (J. Comput. Appl. Math., 1990) — the holonomic framing the certificate lane uses.
• Kauers & Paule, The Concrete Tetrahedron (2011) — the working reference: closure proofs, series↔recurrence, guessing.
• Salvy & Zimmermann, “GFUN: a Maple package for the manipulation of generating and holonomic functions” (ACM TOMS, 1994) — the implementation lineage derive-ode’s discovery half descends from.
• Kovacic, “An Algorithm for Solving Second Order Linear Homogeneous Differential Equations” (J. Symbolic Computation, 1986) — the pole-structure conditions behind decide-ode.
• Risch, “The Problem of Integration in Finite Terms” (1969) — the frame decide-integral answers honestly inside its (much smaller) reach.

The formal semantics the certificates are checked against live in Lean 4, file by file in isar-proofs: HolonomicClosure.lean (sum/product), HolonomicCertificate.integralShift (integration), HolonomicCompose.lean (the composition obstruction), QuantityKernel.lean (the Dim group). The P-recursive projection’s proof file is the one gap — pending, honestly.