Chapter 12 — Control, estimation, and the engineering of precision

Chapters 6 and 7 proved properties of digital circuits — counters, pipelines, protocol FSMs — where the specification is exact and the question is "does the RTL meet it". Control hardware is different: the specification lives in ℝ ("this loop is stable", "this filter attenuates disturbances"), the implementation lives in BitVec, and the gap between them is quantization — a gap that is not a rounding footnote but the central engineering decision. Pick too few fractional bits and a provably-stable design oscillates on the bench. Pick too many and you burn LUTs, closing timing gets harder, and nothing improves.

The chapter has one argument, and every section is a step in it:

Write the system down once, in ℝ. Prove it stable there. Then carry
that ℝ text down to hardware mechanically, and account for exactly what
the trip costs.

Three layers, and the interesting work is at the joints between them:

1. ℝ is the reference. V(f(x)) ≤ ρ·V(x), kernel-checked, no sorry (§12.1, §12.3). 2. The fixed-point equation is *derived*, not re-typed. retype transports the ℝ definition to Q15.16 fixed point (the format is introduced from scratch in §12.1.2), and that transported equation is checked to be the one the RTL datapath actually computes (§12.2). This is the joint people skip, and it is where a proved model quietly stops describing the shipped circuit. 3. The circuit's own stability. Quantization is not a rounding footnote — flooring is a nonlinearity that can sustain a limit cycle — so the ℝ theorem is re-earned as an ISS ultimate bound, and the residual gap is stated plainly (§12.4, §12.10).

And one payoff worth flagging up front: keeping the ℝ equation is not only an obligation to discharge. It is the oracle that lets you search for the implementation — "are these two circuits the same design or two different designs that happen to agree on my tests?" is decidable over ℝ and merely plausible over BitVec (§12.2.5).

The vehicles are real synthesizable designs from IP/Control/: a PID loop, an LQR regulator, IIR filters at five precisions, and two state estimators (Kalman and H∞), one of which computes its own gains on-chip with a multi-cycle divider. Along the way the chapter answers a practical question with measured data — which verification tool covers which claim (SMT, Monte-Carlo falsification, Lyapunov induction), and where each one stops (§12.5).

Everything below that is stated as a number was measured or proven in this repository; the chapter cites the file each time. Where a link in the chain is checked-on-cases rather than proved, it says so.

12.1 The worked example: one PID loop, from equations to theorem

Everything in this chapter is demonstrated on one concrete system, small enough to hold in your head and real enough to synthesize. This section walks it end to end: the equations, the RTL that implements them, and the theorem that certifies them — each line of one mapped to the others.

12.1.1 The equations

A first-order plant (think: the rate response of one drone axis — command a torque, the rate follows with a lag), sampled at dt = 1/16 s:

x[n+1] = 0.9·x[n] + 0.1·u[n]

and the textbook discrete PID regulating it to a setpoint r:

e[n]   = r − x[n]                          error
I[n+1] = I[n] + Ki·e[n]                    integrator      Ki = 0.25
u[n]   = Kp·e[n] + I[n+1] + Kd·(e[n]−p[n]) control         Kp = 2, Kd = 0.125
p[n+1] = e[n]                              previous error

Two state registers in the controller (I, p), one in the plant (x).

12.1.2 The RTL, line for line

Fixed point, and the Q notation, in one paragraph. Hardware has no real numbers, and floating point costs far more than a control loop can usually justify. So a fractional value is stored as a plain integer with an implied scale factor: the 32-bit integer n represents the value n / 2¹⁶. Nothing in the hardware marks where the point is — the convention lives in your head and in the code that reads the register.

Q15.16 names that convention: 16 fractional bits, 15 integer bits, plus one sign bit — 32 in total. So the representable range is about ±32768 (2¹⁵) and the resolution is 2⁻¹⁶ ≈ 0.0000153. That resolution — the value of one integer step — is the LSB (least significant bit), the unit every error bound in this chapter is quoted in. Two consequences follow immediately, and both drive the rest of the chapter:

Adding two Q15.16 values is just integer addition — the scales

already match.

Multiplying them is not. (a/2¹⁶)·(b/2¹⁶) = a·b/2³², so the

integer product carries twice the intended scale and must be shifted back down by 16 to return to Q15.16. That shift discards the low 16 bits, and discarding them is where quantization error enters (§12.4).

So mulQSig below is a 32×32→64 multiply followed by an arithmetic shift right by 16 — widen to avoid overflowing the product, then rescale.

IP/Control/PID.lean implements exactly the four lines above in Q15.16:

def pid (iLim uLim : BitVec 32)
    (r y kp ki kd : Signal dom (BitVec 32)) : Signal dom (BitVec 32) :=
  circuit do
    let integReg ← Signal.reg (0#32)                    -- I
    let ePrevReg ← Signal.reg (0#32)                    -- p

    let e := r - y                                      -- e[n]   = r − x[n]
    let integNext := clampSymC iLim                     -- I[n+1] = I[n] + Ki·e[n]
      (integ + mulQSig ki e)                            --          (clamped: anti-windup)
    let d := e - ePrev                                  -- e[n] − p[n]
    let u := clampSymC uLim                             -- u[n]   = Kp·e + I⁺ + Kd·d
      (mulQSig kp e + integNext + mulQSig kd d)         --          (clamped: actuator limit)

    integReg <~ integNext                               -- register writes
    ePrevReg <~ e                                       -- p[n+1] = e[n]
    return u

The correspondence is one-to-one — that is the point of writing hardware in the same language as the specification. Note the two clampSymCs: they are not in the textbook equations. They are the implementation's two safety nonlinearities, and they carry the first of two very different claims:

1. Bounded, unconditionally. |I| ≤ iLim and |u| ≤ uLim for ANY gains, ANY input, ANY noise — because the clamp is inside the update, the claim is a one-line case split on two comparators, no control theory involved. Tests/IP/Control/PIDTest.lean slams the loop with a huge constant error and watches the integrator saturate at exactly ±16.0:

   ✓ integrator stays within ±16.0 under a huge sustained error
   ✓ integrator actually saturates (the clamp is exercised)

This is why the datapath can never overflow — and it is also why claim 1 is NOT stability: a badly tuned loop happily bangs between the rails forever while satisfying every bound.

2. Converges — the actual control claim. That needs a theorem about the dynamics, which is §12.1.3.

12.1.3 The theorem

Close the loop symbolically (set r = 0; substitute u into the plant) and the three states evolve linearly:

x[n+1] = 0.6625·x + 0.1·I − 0.0125·p      (0.6625 = 0.9 − 0.1·(Kp+Ki+Kd))
I[n+1] = −0.25·x  + I
p[n+1] = −x

Eigenvalues: 0.8717, 0.8085, −0.0177 — all inside the unit circle, so the loop is stable. But "I computed eigenvalues" is a numerical remark, not a proof. The machine-checked version (proofs/SparkleProofs/Control/PIDDesign.lean) exhibits a quadratic certificate instead:

P = ⎡ 8.0999  −5.4050  −0.0850 ⎤        V(s) = sᵀPs
    ⎢−5.4050   9.7880   0.0574 ⎥
    ⎣−0.0850   0.0574   1.0013 ⎦

and proves, over ℝ with Mathlib, zero sorry:

theorem pid_lyapunov_decrease (x I p : ℝ) :
    V (nextX x I p) (nextI x I p) (nextP x I p) ≤ (39/40) * V x I p := by
  unfold V nextX nextI nextP p11 p12 p13 p22 p23 p33 Kp Ki Kd pa pb
  nlinarith [sq_nonneg (x + (1694/10000)*I + (27/10000)*p),
             sq_nonneg (I − (25/10000)*p), sq_nonneg p, …]

Every sample, the energy V shrinks by at least the factor 39/40 — for every state, not the trajectories you happened to simulate. Together with P ≻ 0 (Sylvester, also proven) and the sandwich 0.999·‖s‖² ≤ V(s) ≤ 15·‖s‖², the corollary pid_geometric_decay gives geometric convergence of the state itself, by a five-line induction.

Where did P and those strange nlinarith hints come from? Offline — and this recipe is used for every quadratic-form proof in the repo, so learn it once:

1. iterate the discrete Lyapunov equation P ← AᵀPA + I numerically to a fixed point; round to 4 decimals; 2. sweep for the true worst ratio V(As)/V(s) (here 0.9306) and pick a round certified rate above it (39/40 = 0.975 — slack is robustness); 3. compute ρP − AᵀPA over exact rationals and take its exact LDLᵀ; if all pivots are positive (here 0.798 / 0.732 / 0.975) the matrix is PSD and the LDLᵀ rows are literally a sum-of-squares witness; 4. hand those rows to nlinarith as sq_nonneg hints. The proof lands on the first try, because you are not asking the tactic to find the certificate — only to check it.

Guess-and-verify is the honest shape of every Lyapunov argument; the kernel checking step 4 is what turns the numerical remark into a theorem.

12.1.4 What just happened, and what did not

The theorem is about the ℝ model. The RTL runs Q15.16. That gap is the subject of the rest of the chapter, and it has two distinct halves that are easy to conflate:

Is the RTL computing the same equation? A question about

transport — answered by deriving the fixed-point equation from the ℝ one instead of hand-writing it, then checking the datapath against the derivation (§12.2).

Does the equation still converge once every product is floored? A

question about dynamics — answered by ISS: the quantization error of each multiply is one-sided and < 1 LSB (mulQ_error), and a Lyapunov contraction absorbs bounded per-step disturbances into a computable ultimate bound (§12.4). Plus the circuit held to the pure model cycle-by-cycle on three backends (§12.5d).

Answering only the second — the usual practice — leaves you with a proof about an equation you have not established the hardware implements.

What is deliberately not claimed: that the clamps never engage (they exist precisely for the transients where they do), and that the certificate covers the saturated regions — extending V piecewise across the clamp boundaries is the natural next theorem, and is open.

12.2 The spine: one equation, three layers

Everything from here to §12.10 is one argument, and it is worth stating before the details bury it. There is one system. It is written down once, in ℝ, and that ℝ text is the reference — the thing the design means. Everything else in the chapter is an attempt to carry that meaning down to something you can put on an FPGA, and to be honest about what survives the trip.

The system is the one you already have in your head: the PID closed loop of §12.1.1. Written out over ℝ, its three states evolve by the three functions §12.1.3 proved a certificate for —

nextX x I p  =  0.6625·x + 0.1·I − 0.0125·p      -- plant
nextI x I p  =  −0.25·x + I                      -- integrator
nextP x I p  =  −x                               -- previous error

— and those three names are what the rest of this section transports.

Three layers, three different kinds of claim:

Layer	Artifact	What is claimed	How
1. Design	`nextX, nextI, nextP : ℝ → ℝ`	`V(f(s)) ≤ ρ·V(s)`, ρ = 39/40	Lean proof, kernel-checked
2. Implementation	`nextXQ, nextIQ, nextPQ : FixQ → FixQ`	the Q15.16 equations are the ℝ equations, retyped	mechanical transport, then equality vs the RTL
3. Circuit	`IP/Control/PID.lean` (Verilog)	ℝ stability + bounded error ⟹ ultimate bound	ISS, plus co-sim

Layer 1 is §12.1.3 — already done, same three functions. Layer 3 is §12.4. This section is layer 2, which is the one that is easy to skip and is exactly where designs go wrong: the ℝ model is proved, the RTL is simulated, and nobody ever checks that the RTL implements that model rather than a slightly different one someone typed in by hand.

12.2.1 Transporting the equation, not re-typing it

The usual practice is to write the fixed-point version by hand. That creates a second source of truth, and the two drift: a gain rounded differently, a - that became a saturating - on one side only. No amount of simulating the RTL finds this, because the RTL is being compared against itself.

Instead, generate it. retype replaces a type throughout a definition and substitutes the corresponding operations, so the fixed-point equation is derived from the ℝ equation rather than written beside it:

structure FixQ where            -- Q15.16: the stored integer is value·2¹⁶
  n : Int

instance : Add FixQ := ⟨fun a b => ⟨a.n + b.n⟩⟩
instance : Mul FixQ := ⟨fun a b => ⟨(a.n * b.n) / 65536⟩⟩   -- floors; see §12.2.2

declare_retype RealToFixQ : Real => FixQ

-- the gains and plant of §12.1.1, unchanged
noncomputable def Kp : ℝ := 2 ; noncomputable def Ki : ℝ := 1 / 4
noncomputable def Kd : ℝ := 1 / 8
noncomputable def pa : ℝ := 9 / 10 ; noncomputable def pb : ℝ := 1 / 10

noncomputable def nextX (x I p : ℝ) : ℝ :=
  (pa - pb * (Kp + Ki + Kd)) * x + pb * I - pb * Kd * p

retype_def KpQ := Kp using Real => FixQ    -- …and Ki, Kd, pa, pb
attribute [retype RealToFixQ] Kp Ki Kd pa pb   -- so nextX's references follow
retype_def nextXQ := nextX using Real => FixQ

nextXQ is now a Q15.16 function nobody wrote. Evaluating it:

#eval KpQ                     -- { n := 131072 }  = 2        ✓
#eval KiQ                     -- { n := 16384 }   = 0.25     ✓
#eval paQ                     -- { n := 58982 }   ≈ 0.89999  ✓
#eval nextXQ ⟨65536⟩ ⟨0⟩ ⟨0⟩   -- { n := 43419 }   = 0.66250  ✓

The last line is the check worth pausing on. Set x = 1, I = 0, p = 0: the ℝ equation collapses to its leading coefficient, which §12.1.3 displays as 0.6625 — and 43419/65536 = 0.66250. The transported equation reproduces the number the ℝ proof is about, and no human transcribed a coefficient.

The attribute line is the one non-obvious step. Without it, nextX's references to Kp/pa stay at type ℝ and the elaborator reports pa has type ℝ but is expected to have type FixQ; the attribute is what makes a definition's dependencies follow the transport too.

12.2.2 Does the RTL implement this?

Now the question that layer 2 exists to answer. Sparkle's datapath uses mulQ on BitVec 32 (IP/Control/FixedPoint.lean):

def mulQ (a b : BitVec 32) : BitVec 32 :=
  BitVec.extractLsb' 16 32 ((a.signExtend 64) * (b.signExtend 64))

and the retyped model uses (a.n * b.n) / 65536 on Int. These are the same function, and the reason is not obvious: extractLsb' 16 on a sign-extended product is an arithmetic right shift, which rounds toward −∞, and Lean's Int./ also floors. Had the hardware used a truncating shifter, or had FixQ used Int.tdiv, the two would agree on positives and disagree on every negative product.

That is a claim, so it is checked rather than asserted:

def refMul (a b : Int) : Int := (a * b) / 65536
def chk (a b : Int) : Bool :=
  (mulQ (BitVec.ofInt 32 a) (BitVec.ofInt 32 b)).toInt == refMul a b

#eval cases.filter (fun (a,b) => !chk a b)     -- []

over sign-crossing and boundary cases ((-1, 65535), (3, -65537), (-58982, -65536), …) — all agree.

Both halves are live code, not listings: proofs/SparkleProofs/Retype/FixedPointTransport.lean does the transport and pins every constant with #guard (131072, 16384, 8192, 58982, and the step result 43419 — each checkable by hand against §12.1.1), and Tests/IP/Control/PrecisionSweepTest.lean's Transport agreement suite runs the comparison on the Sparkle side under lake test:

Transport agreement (ℝ →retype→ Q15.16 vs Sparkle mulQ):
  ✓ every case agrees
  ✓ negative product with remainder floors down (not toward zero)
  ✓ positive product with remainder floors to zero
  ✓ pa = 0.9 transports to 58982
  ✓ one observer step matches the transported equation

So the datapath's multiply is the retyped multiply, and the chain

ℝ equation → (retype) → Q15.16 equation → (this check) → RTL datapath

closes. §12.5d then holds the whole circuit to the pure model cycle-by-cycle on three backends, which extends the agreement from one operation to one design.

12.2.3 What layer 2 does and does not buy

It buys the elimination of a specific, common, silent failure: the implemented equation not being the designed equation. After this section, "the RTL computes something slightly different from the model" is not on the table.

It also does not mean the transported equation equals the ℝ answer. Worth seeing concretely: nextXQ ⟨65536⟩ ⟨0⟩ ⟨0⟩ gives 43419, while quantizing the already-simplified constant 0.6625 in one go gives 43417. Neither is wrong. The transport evaluates the expression (pa − pb·(Kp+Ki+Kd))·x in Q15.16, flooring at each product; the other rounds the ℝ answer once. Two LSB apart — and that gap is precisely the per-step quantization error that §12.4 has to bound rather than wish away. (Both numbers are pinned in Tests/IP/Control/PrecisionSweepTest.lean, so neither can drift silently.)

It does not make the fixed-point loop stable. nextXQ is a different dynamical system from nextX — it floors every product — and a floor is not a small perturbation of the identity, it is a nonlinearity that can sustain a limit cycle (§12.4 shows one doing exactly that, and §12.4's "more bits ≠ better" warning shows a finer grid behaving worse). Layer 2 says the two systems are related by transport; layer 3 has to say what that relation preserves.

This layer's last sampled link is now proved (Sparkle/Verification/FixedPointProps.lean):

theorem mulQ_toInt (a b : BitVec 32)
    (hlo : -(2^31 : Int) ≤ (a.toInt * b.toInt) / 2^16)
    (hhi : (a.toInt * b.toInt) / 2^16 < 2^31) :
    (mulQ a b).toInt = (a.toInt * b.toInt) / 2^16

So §12.2.2's fixture sweep is no longer load-bearing: the datapath multiply is the transported multiply, for all inputs whose result is representable.

Two things about it are worth keeping.

The hypothesis is not decoration. Without it the statement is false — the extract wraps. At the 8/4 analogue (small enough to enumerate) a = 17, b = 121 gives 17·121/16 = 128, one past the signed 8-bit maximum, and the extract returns −128. All three facts — holds under the guard, fails without it, and that exact counterexample — are pinned in Tests/IP/Control/PrecisionSweepTest.lean, so the side condition cannot be dropped later as pedantry.

`bv_decide` cannot do it. The obvious move fails outright —

error: None of the hypotheses are in the supported BitVec fragment

— because toInt and Int./ leave the bitvector fragment it decides; no width helps. The proof instead goes: the 64-bit product of two sign-extended 32-bit values is exact (|a·b| ≤ 2^62 kills the bmod), then extractLsb' 16 32 is / 2^16 on toInt given the range hypothesis (Int.bmod unfolded and split, since omega cannot see through it), then compose. Mathlib-free, so it is checked by the same lake build as the RTL.

The other gap this section used to carry is now closed. The transport lives in proofs/SparkleProofs/Retype/, alongside the theorems, and it imports nextX/nextI/nextP from PIDDesign.lean — the very terms pid_lyapunov_decrease is proved about. Earlier it had to restate them (retype pinned a newer Lean than Sparkle, and a Lake graph has one toolchain), so "the transported equation is the certified equation" was underwritten by #guards comparing two copies. It now holds because there is only one copy.

12.2.4 Rewriting the equation: ℝ-equal is not circuit-equal

A practical question that arrives the moment you care about area: can you rewrite the control law to use fewer multiplies?

u = Kp·e + Ki·s + Kd·(e − p)        -- as written
u = (Kp+Kd)·e + Ki·s − Kd·p         -- Kd folded into the e coefficient

Over ℝ these are the same function, and ring proves it in one line (AlgebraicRewrite.uA_eq_uB). It is tempting to stop there.

Don't. In Q15.16 they are different circuits, because each mulQ floors and the rewrite moved the floors around. Measured over 200 000 random states:


states where the two disagree	87 751 / 200 000 (44 %)
worst disagreement	1 lsb

Disagreement is the common case, not a corner case. So the honest form of the claim is not "the rewrite is safe" but

exact over ℝ, and at most N lsb apart in Q15.16, where N comes from
the floor count

and both halves need proving. uAq_uBq_gap does the second: each shape is within 3 lsb of the shared ℝ value, so they are within 6 lsb of each other. Note the gap between proved 6 and measured 1 — floor-counting is sound but loose, because the floors are strongly correlated. A tight bound needs a finer argument.

This matters for the certificate, and it is the same trap as §12.2.3's gains one level up: prove stability for uA, synthesize uB, and the Lyapunov argument covers a system the hardware does not implement — unless the difference is carried as an extra disturbance into the ISS bound.

Using this to search for an RTL shape. If you want a tool to pick the implementation — try rewrites, keep the cheapest that still meets the error budget — the useful question is not "how big is the error" but "which rewrites cost nothing", because those need no bookkeeping at all:

rewrite	fixed point	cost to a search
reassociate / commute `+`	exact	free
`Kp·e + Kd·e → (Kp+Kd)·e`	exact iff one gain is integral	free under a check
`K·(a+b) → K·a + K·b`	1 lsb	budgeted

The middle row is the trap. Folding two gains saves a multiplier and looks like plain constant folding, so it is tempting to treat as always-safe. It is not — and the criterion is not the one you would guess:

2.0  + 0.125  → exact      (2.0 is a whole number)
2.0  + 0.618  → exact      (the OTHER gain can be anything)
0.5  + 0.25   → 25 % differ  ← both powers of two, still wrong
0.618 + 1.26  → differs

Being dyadic is not what matters; having no fractional part is, because then a·e / 2¹⁶ is an integer and the floor discards nothing. AlgebraicRewrite.fold_exact proves it under scale ∣ a, and PrecisionSweepTest pins all four rows above — including the counter-intuitive 0.5 + 0.25.

So a search can apply the free rules without touching the certificate, and only pay the ISS-disturbance price for the rewrites that genuinely cross a floor.

12.2.5 Why keep the ℝ equation at all — searching for the RTL

Layers 1–3 make the ℝ text sound like an obligation: something you must write down and then discharge. This section is the argument that it pays for itself, and it is the strongest practical reason to keep an ℝ model around after the RTL exists.

Ask the design question directly: which of the many ways to write this control law should the hardware use? The candidates all compute the same thing on paper and differ by a multiplier here, a shift there. Choosing between them needs two judgements:

1. Are these actually the same function? 2. What does each one cost — in area, and in error?

Without an ℝ model you can only answer the second. You have Q15.16 expressions; you can simulate them and compare outputs; and when two of them differ by 1 lsb you cannot tell whether you are looking at

the same design, rewritten, with quantization accounting for the gap, or
a different design that happens to be close on the states you sampled.

Those need opposite responses — accept the first, reject the second — and the fixed-point text alone does not distinguish them. The ℝ equation does, because equality over ℝ is decidable by ring: no sampling, no threshold, no judgement call. It is the oracle that turns "these agree on my tests" into "these are the same function, and here is the price of the rewrite".

That is worth stating plainly, because the ℝ model is usually sold as documentation or as a proof artifact. It is also a search oracle, and that role survives long after the stability proof is finished.

The search itself. IP/Control/ShapeSearch.lean implements this: enumerate the rewritings, keep the ℝ-equivalent ones, price each by multiplier count and measured lsb gap. On 0.618·x + 1.26·x it reports

  muls=1 gap=1 lsb  (1.878 * x)                  ← fold: one multiplier less
  muls=2 gap=0 lsb  (0.618*x + 1.26*x)           ← original: exact

which is the shape of the answer you want: the cheaper option is offered with its price attached, not silently applied. On the dyadic-gain PID it finds only muls=0 gap=0 shapes — every gain is already a shift, so reassociation is all that is on the table and all of it is free.

Two limits, both deliberate. The gap is measured over a sample, so it is a lower bound on the worst case — use the search to find candidates, then discharge the chosen one against the proved bound (uAq_uBq_gap). And admission is by ℝ-equivalence on samples, not by ring; a candidate that survives is a proposal, and ring is what makes it a theorem. This is the same division of labour as the Float falsifier in SparkleProofs/Retype/: search cheaply, prove the survivor.

The loop, end to end. Putting the pieces together, choosing an implementation looks like this — and every step is a check the repo already has:

   ℝ equation  ──────────────────────────┐
        │                                │  ring        (is it the same function?)
        │ enumerate rewrites             │
        ▼                                │
   candidates ──── ℝ-equivalent? ────────┘
        │              (ShapeSearch)
        ▼
   price each:  multipliers          (cost model)
                lsb gap             (fold_exact / floor counting)
        │
        ▼
   still within budget?  Vbound_mono + Precision.Vbound
        │
        ▼
   synthesize the winner, discharge its gap against uAq_uBq_gap

The budget step is the one that would otherwise be a re-proof. Because Vbound_mono says the ultimate bound is monotone in the per-step error, a candidate that adds extra lsb needs only a numeric comparison — no fresh Lyapunov argument per shape. That is what makes searching over dozens of shapes practical rather than a proof-engineering project.

And a caution the search makes concrete. Run it on the dyadic-gain PID and the ranking is not the one intuition offers. Folding Kp + Kd looks like the obvious saving — one multiply instead of two — but 2.0, 0.25 and 0.125 are all powers of two and lower to shifts, so the original costs zero multipliers and the "optimized" 2.125·e costs one. The rewrite that looks cheaper on the page is more expensive in silicon. A cost model that knows about shifts catches this; reading the expression does not.

Divisions carry one extra obligation. (a/b)·x = a·(x/b) needs b ≠ 0 (div_reassoc states it), so the rewrite is unlicensed where the guard fails. And in fixed point a division is divQref, which truncates — its error interval is two-sided (−1,1) lsb, against mulQ's one-sided (−1,0]. Moving a division across a rewrite changes the shape of the error, not just its size.

12.3 Lyapunov stability in general

What §12.1.3 did for one system is the general recipe. A Lyapunov function for s⁺ = f(s) is any V ≥ 0 with

V(f(s)) ≤ ρ·V(s),      ρ < 1     (for all s)

— existence of such a V is geometric stability, and for linear systems a quadratic V(s) = sᵀPs always works when the system is stable (solve the discrete Lyapunov equation, as above). Two more instances live in this repo, proved with the same LDLᵀ recipe:

the LQR double integrator (LQRDesign.lean) — the chapter's vehicle

for the quantization story, ρ = 39/40, true ratio 0.97179;

both estimator error dynamics (EstimatorDesign.lean, §12.6) — where

the same recipe also proves the H∞ dissipation inequality, a 4-variable quadratic form.

Why insist on the contraction form V∘f ≤ ρV rather than mere decrease V∘f < V? Because only the former survives disturbances: ρ < 1 leaves room to absorb a bounded perturbation into a geometric series (§12.4's ISS argument), whereas strict-decrease-by-an-unquantified-amount absorbs nothing. Certificates should always be stated with slack.

12.4 Quantization as a bounded disturbance

The synthesized datapath does not compute over ℝ. It computes Q15.16 (§12.1.2), and the shift that rescales every product throws the low 16 bits away. Layer 2 established that the RTL computes the right equation; this section is the other half — what flooring every product does to the ℝ theorem. The bridge rests on one small fact and one standard theorem.

The small fact. In Lean 4, Int./ floors — and BitVec.sshiftRight also floors. They are the same function. So the RTL's >>> 16 and the spec's / 2^16 agree exactly, and one Q15.16 multiply errs by less than one LSB, always downward — error confined to (−1, 0] LSB, independent of sign (proofs/SparkleProofs/Control/Transport.lean, mulQ_error). Contrast the divider below, which truncates toward zero and pays a two-sided (−1, 1) LSB interval for it.

The step that used to be folklore. "A per-step error bounded by ε" is the hypothesis the whole ISS argument runs on, and until recently nothing connected it to the integers the circuit iterates — ultimate_bound was stated for a QuantTraj, a structure whose fields assume x⁺ = nextX x + d with |d| ≤ ε. It read "IF the hardware trajectory has this shape, THEN it is ultimately bounded", and the antecedent was asserted rather than proved.

StepError.lean discharges it. Counting floors per component:

channel	update	floors	proved bound
`x1`	`x1 + mulQ dtQ x2`	1	1 lsb (`stepX1_err`)
`x2`	`x2 + mulQ dtQ (-(mulQ k1Q x1 + mulQ k2Q x2))`	3	3 lsb = ε (`stepX2_err`)

mkTraj then packages an integer state sequence into a QuantTraj — the disturbances are not invented, d n is defined as the difference between what the circuit computed and what ℝ would have, so the step equations hold by ring and the only content is the bound above. The consequence is stated with no QuantTraj in sight:

theorem intTraj_ultimate_bound (f1 f2 : Nat → ℤ) (h1 h2 : …) (n : Nat) :
    V (toR (f1 n)) (toR (f2 n))
      ≤ σ ^ n * V (toR (f1 0)) (toR (f2 0)) + Vbound

And the gains, too. stepX2_err is stated against the ℝ update using the quantized gains toR k1Q, toR k2Q — not the exact k1 = 0.6180, k2 = 1.2600, which round 0.248 and 0.360 LSB away. That difference cannot be absorbed into ε: a gain error multiplies the state, so its contribution grows with |x|.

The fix is not a bigger ε — it is to certify the system the hardware implements. QuantizedGains.lean does that: same V, same P, same ρ = 39/40, only the gains changed. It works because the certificate was built with slack — the true contraction ratio moves from 0.97178926 to 0.97178930, a shift of 3.7·10⁻⁸, while ρ sits at 0.975. This is the concrete payoff of §12.1.3's step 2 ("pick a round rate above the measured ratio"): an implementation detail perturbs the plant and the proof still holds. A ρ pinned at 0.9717893 would have needed a new P.

Composing the two gives the statement with nothing left conditional:

theorem circuit_ultimate_bound (f1 f2 : Nat → ℤ)
    (hs1 : ∀ n, f1 (n+1) = stepX1 (f1 n) (f2 n))
    (hs2 : ∀ n, f2 (n+1) = stepX2 40501 82575 (f1 n) (f2 n))
    (n : Nat) :
    V (toR (f1 n)) (toR (f2 n))
      ≤ σ ^ n * V (toR (f1 0)) (toR (f2 0)) + Vbound

f1/f2 are the integers in the registers, 40501/82575 are the gains in the ROM, stepX1/stepX2 are the floors the datapath performs. (dt = 1/16 is dyadic, so toR dtQ = dt exactly — proved, not assumed.)

The standard theorem (ISS). A per-step error bounded by ε is a bounded disturbance, and a Lyapunov contraction survives bounded disturbances in degraded form:

V(x⁺) ≤ σ·V(x) + c·ε²,     σ = (1+ρ)/2 < 1

Iterating telescopes into an ultimate bound: the state enters and stays in a ball whose radius scales with ε (Transport.ultimate_bound — an unbounded-horizon induction, checked by the kernel). Quantization does not threaten stability; it buys a steady-state error floor. How big a floor is exactly the precision question:

Vbound(f) = 583200 / 4^f        (f = fractional bits)

proofs/SparkleProofs/Control/Precision.lean proves this closed form and its consequences: each extra fractional bit cuts the floor by 4×, and for a representative budget V ≤ 0.01,

f = 12 fails (Vbound ≈ 0.035),
f = 13 passes (Vbound ≈ 0.0087), and monotonically every

f ≥ 13 passes (min_fracBits_for_budget).

"How many fractional bits do I need?" is a theorem now, not a guess.

A measured warning: more bits ≠ better behaviour

Tests/IP/Control/PrecisionSweepTest.lean runs one impulse through the same marginal IIR resonator (ℝ poles at radius 0.999) at five formats. Residual ringing after 200 samples, in 1e-3 units:

f = 4  → 62        f = 8  → 0        f = 16 → 52        f = 24 → 52

Non-monotone. Two separate mechanisms, and the bound above covers only one of them:

Datapath quantization (covered): adds the ε-disturbance, floor

shrinks 4× per bit.

Coefficient quantization (not covered — it perturbs ρ, not ε):

rounding the coefficients moved the poles inward at every format (radius 0.968 at f=4, 0.998 at f=8, 0.99899 at f=16 — all stable), and at coarse f the deadband — products flooring to zero once the state is below one LSB — kills the ringing outright, while fine f faithfully reproduces the marginal design and keeps ringing for hundreds of cycles.

Fine precision gives you a more faithful copy of whatever you designed — including its flaws. Precision and design margin are separate budgets; the worksheet in §12.9 keeps them separate.

Also measured in the same sweep: Q7.8 (16-bit) and Q23.8 (32-bit) produce bit-identical output — they share f = 8. Width buys range (later saturation), never accuracy. If a review comment says "widen the datapath for accuracy", this test is the counterexample.

12.5 Three verifiers, three coverage zones — measured

You now have three ways to check a fixed-point control claim. They are not interchangeable; each has a hard edge, and we measured where.

(a) SMT / `bv_decide` (SAT over bit-vectors, in Std). A decision procedure: it proves or hands you a concrete counterexample assignment. Measured on this machine, on a two-state biquad-shaped recurrence (4 multiplies per step, state boxed, asking for an output bound), timeout 420 s:

width	steps unrolled	result
8, 16, 32	1	proved, ≈0.5 s (even with a 32×32 multiply)
32	8	proved, 124 s
32	16	proved, 8 s (SAT solvers are non-monotone)
16	20	proved, 187 s
16	24	timeout
32	24	timeout

Two lessons. First, the folklore "SMT dies at 32 bits" is wrong for single-step goals — a one-step invariant with a full-width multiply is instant at any width, so inductive step obligations are ideal SMT targets. Second, the wall is unrolling depth, around N≈20–24 here regardless of width. Bounded-model-checking a trajectory further than that is hopeless, and an unbounded claim ("never oscillates, ever") is structurally out of reach — BMC can only ever say "no violation in the first N steps". This is precisely the gap DSVerifier-style tools live in.

(b) Monte-Carlo falsification (proofs/SparkleProofs/Retype/). The ℝ design is transported to executable Float by retype (retype_def VF := V using Real => Float) and hammered with 10⁵ random states before anyone writes a proof. Measured output, from the build log:

contraction: counterexample = none, worst ratio = 0.971789  (ρ = 0.975)
ISS:         counterexample = none, worst overshoot = 0
negative control (ρ = 0.97 < true 0.97179): 4832 violations found

The negative control is the important line: with the rate set just below the true worst case the search does find thousands of violations, so "found nothing" is evidence, not vacuity. Cost: milliseconds per candidate certificate. What it can never give: a guarantee. Its role is to kill wrong constants cheaply — the ρ and the Young-split δ in LQRDesign.lean each went through one wrong guess that the Float harness would have (and later did) flag instantly, versus a slow failed nlinarith round-trip.

(c) Lyapunov induction in the kernel (proofs/). Unbounded horizon, all states, quantization included — the only tool of the three that can state "for every n" — at the price of needing a certificate to verify and nlinarith labour to check it. The division of labour that emerges, and that this repo actually uses:

Float search  →  find/refute the certificate candidate   (milliseconds)
bv_decide     →  discharge one-step BitVec side goals     (seconds)
Lyapunov      →  the unbounded claim itself               (kernel-checked)

(d) …and none of them see the backend. A war story from building this chapter, kept because it is the sharpest possible illustration of coverage zones. The width-generic divider below passed every Lean-side check — pure-model-vs-reference on 40 cases, all proofs, #synthesizeVerilog accepted it, iverilog parsed the output — and the emitted RTL still computed 1.0/3.0 = saturate. Cause: the Verilog backend's width inference hit the symbolic width w + f + 1, could not reduce Nat.succ 48 to a literal, and silently defaulted to 8 bits, so a 49-bit divisor register was fed through an 8-bit wire and latched zero. Nothing on the Lean side can catch that class of bug, because the bug is in the translation itself. What caught it: simulating the emitted Verilog against the pure model — which then agreed bit-for-bit once the backend was fixed (five full time-varying-Kalman samples, FSM + shared divider + covariance recursion: 0, 0, 1581, 8167, 22409 from both). A second lesson from the same week, one layer up: the interpreted Signal.val co-sim of a multi-register circuit do FSM hangs (issue #95 — the Circuit monad composes per-write state-update closures, so evaluation cost grows ~k^k in the register count). An attempted fix made simulation linear, but it changed expression sharing enough to alter the generated RTL of nested engines behind multi-output records — caught only by re-running the emitted RTL against a known-good Keccak digest, and therefore not merged. The estimators here are co-simulated through the CSim JIT (lake exe control-jit-test) and iverilog instead.

The moral for the toolbox: theorem, SAT and search all verify the model; only executing the artifact — the emitted RTL, the compiled simulator — verifies the compiler that produced it. Keep one artifact-level cross-sim per backend in the loop no matter how much you have proven, and treat a green #synthesizeVerilog as "well-formed", never as "correct".

12.6 Estimators: Kalman and H∞ are the same circuit

Real loops close on estimated state — gyro rates are noisy, and position comes through a filter. IP/Control/Observer.lean implements the standard predictor-form observer for one drone axis (angle + rate, y = angle + noise):

x̂₁⁺ = x̂₁ + dt·x̂₂ + k₁·(y − x̂₁)
x̂₂⁺ = x̂₂ + dt·u  + k₂·(y − x̂₁)

Here is the fact that surprises people: the steady-state Kalman filter and the H∞ filter are this same RTL with different values in the two gain constants:

Kalman:  K = [0.4636, 1.3960]     (Riccati, q = 1/32, r = 0.01)
H∞:      K = [0.4974, 1.5472]     (H∞ Riccati at γ = 1.964; γ_min ≈ 1.309)

Every meaningful difference lives offline — in how the gains were designed and in what can be proven about them (proofs/SparkleProofs/Control/EstimatorDesign.lean):

Both gains: the error dynamics contract a quadratic V at

ρ = 0.98 (kf_contraction, hinf_contraction). Convergence under zero disturbance — shared, and not the interesting part.

Only the H∞ gain: the dissipation inequality

  V(e⁺) − V(e) ≤ γ²·(w²/q + v²/r) − ‖e‖²        γ = 2, for ALL w, v

which telescopes (hinf_energy_bound) into: over any horizon, the estimation-error energy is at most γ² times the weighted disturbance energy, plus the initial storage. No Gaussian assumption, no smallness assumption — w may be a worst-case gust, a sensor glitch, an adversary. Kalman's optimality claim is an average under the noise model you assumed; this is a bound under no model at all.

The on-chip Riccati (`tvKalman`)

The full time-varying Kalman filter propagates its covariance in hardware and divides each sample to get the gain — the operation Signal has no operator for. This is what the width-generic restoring divider IP/Control/DividerQ.lean exists for (the RV32 integer divider can't do fractional Q15.16 division; its proof is pinned to 32 bits, so the generic core is a new module in the same loop shape). Gain step, 50 cycles per division at Q15.16:

s  = p₁₁ + r          -- r > 0 and p₁₁ clamped ≥ 0  ⇒  never divides by zero
k₁ = (p₁₁ + dt·p₁₂)/s
k₂ = p₁₂/s

A 5-phase FSM (IDLE → DIV1 → START2 → DIV2 → UPDATE) shares one divider engine between the two divisions. Measured cross-validation (Tests/IP/Control/ObserverTest.lean): iterated from P = 0, the on-chip fixed-point Riccati converges to gains within 49 and 12 LSB of the offline design constants — the hardware recursion and the offline script confirming each other to ~0.1%.

12.7 Use case A: the fast drone

A racing-drone rate loop closes at 1–8 kHz on a controller running at 27 MHz. Sample-period budget at 4 kHz: 6750 cycles. What the numbers in this chapter say about that design:

Latency: the fixed-gain observer and the PID/LQR update are

single-cycle datapaths — measurement to actuator in one clock, 0.015% of the budget. The tvKalman FSM needs ~115 cycles (two 50-cycle divisions): 1.7% of the budget. Both fit trivially — run the fixed-gain observer in the fast loop and, if you want adaptive gains, let tvKalman update them at a slower rate.

Precision: the inner loop's signals are small (rate errors,

±10 rad/s) and the budget is tight. Vbound(f) = 583200/4^f says f = 13 meets a 0.01 budget; Q15.16 gives 30× margin, and §12.4's sweep says the 16-bit Q7.8 datapath — half the multiplier area — fails the same budget by three orders of magnitude. The correct cheap choice is a 16-bit container only if you can spend 13+ bits on fraction, i.e. Q2.13 for a ±4 signal range: the theorem, not the bit-width folklore, decides.

Verification split: the per-step overflow-freedom obligations

are 1-step bv_decide goals (instant at 32 bits); the "never diverges under quantization" claim is the transported Lyapunov bound; and every gain retune gets the millisecond Float-falsifier pass before anyone re-proves anything.

12.8 Use case B: the disturbance-heavy environment

Now the other regime: an inspection drone next to a building in gusty wind, a vehicle with unmodelled vibration — the disturbance is large, structured, and not the white Gaussian noise the Kalman design assumed. Measured head-to-head (deterministic seeds, ObserverTest.lean), error energy over the run:

profile	Kalman	H∞	winner
A: LCG noise (the KF's design case)	1 081 071	1 087 559	KF by 0.6 %
B: square-wave gust, half-period 10 (worst of a swept range)	1 282 428	1 175 037	H∞ by 9 %

Read this honestly, because the honest reading is the lesson:

On its home turf the Kalman filter wins, by less than a percent.

H∞'s insurance premium is small.

Under the adversarial gust H∞ wins by single-digit percent. If 9%

were the whole story, you might not bother.

The story is the certificate: hinf_energy_bound holds for

every disturbance sequence — including the ones the bench never produced. The Kalman gain has no such theorem at this γ. When the disturbance model is the thing you don't trust, the filter with the worst-case guarantee is the engineering choice, and the guarantee — not the 9% — is the product.

This is also where the three-verifier split from §12.5 pays off: the dissipation inequality is a 4-variable quadratic form — far beyond BMC's unrolling wall as a trajectory claim, but as a one-shot algebraic inequality it is exactly what an SOS certificate plus nlinarith handles, and the exact rational LDLᵀ that seeds those hints was found in milliseconds offline.

12.9 The precision-selection worksheet

The procedure this chapter justifies, step by step, for a new control datapath:

1. Design in ℝ (any tool). Extract gains and a Lyapunov / dissipation certificate candidate. 2. Falsify first (SparkleProofs.Retype pattern): retype the model to Float, sweep 10⁵ states against the candidate, include a negative control. Fix constants until the search goes quiet. 3. Prove the ℝ certificate (proofs/ pattern): nlinarith with LDLᵀ-derived square hints; keep ρ deliberately loose. 4. Pick f from the budget: Vbound(f) = c/4^f, take the minimal f under budget (here: theorem min_fracBits_for_budget), then one bit of margin. Pick w from range (saturation headroom), never from accuracy. 5. Sweep the formats in simulation (PrecisionSweepTest pattern): confirm the residuals, and check the coefficient-quantization effect separately — it moves poles, and the bound from step 4 does not see it. 6. Discharge step obligations with `bv_decide` where they are 1-step; never ask it for an unbounded claim. 7. Synthesize every format you compared — a comparison in which only one column is real hardware is a spreadsheet, not a sweep.

Every step exists as working code in this repository; the chapter's role was only to put them in order.

import IP.Control.IIRBiquadGen

open Sparkle.IP.Control.IIRBiquadGen
open Sparkle.IP.Control.FixedPointGen

-- The §12.4 sweep, live: the same marginal resonator at three formats.
-- Residual ringing amplitude (×1e-3) after 200 quiet samples:
def tail (w f : Nat) : Nat :=
  let impulse := (q w f 1 1) :: List.replicate 300 (BitVec.zero w)
  let ys := run w f (quantize w f marginalCoeffs) (limOf w f)
    ⟨BitVec.zero w, BitVec.zero w⟩ impulse
  (((ys.drop 200).map (fun y => (y.toInt * 1000 / (2 ^ f : Int)).natAbs))).foldl
    Nat.max 0

#eval tail 16 8    -- 0   : coarse f — the deadband killed the ringing
#eval tail 32 16   -- 52  : fine f — faithfully still ringing
#eval tail 32 8    -- 0   : same f as Q7.8, twice the width — identical

12.10 Where this stops, honestly

The Signal-level equality (circuit do = the pure step

functions, every cycle) is checked by cycle-accurate co-simulation, not yet by the loop_iterate proof — that bridge exists in the repo (Sparkle/Verification/Divider/ proved it for the RV32 divider) and is the natural next proof.

All certificates are verified, not synthesized: the Riccati and

LDLᵀ computations run in offline scripts, Lean checks the inequalities. A DARE solver inside Lean would close that gap.

~~The H∞ fixed-point transport is not written~~ — closed for the

contraction half. HinfTransport.lean certifies the quantized gains (32598/65536, 101397/65536; true ratio moves 0.97224763 → 0.97223470 against ρ = 0.98) and carries them through ISS to hinf_ultimate_bound. Two constants had to be recomputed rather than copied: the Young split is 1/100, not the LQR's 1/80 — (81/80)·0.98 = 0.99225 > 0.99 does not close — and the disturbance constant is 56 rather than 10, because hiP is scaled differently (p11+2|p12|+p22 = 55.71). The dissipation inequality (hinf_energy_bound) is still ℝ-only; transporting an energy bound is a different argument from transporting a contraction, and is open.

~~The ISS argument assumes the hardware trajectory is an ε-perturbed

ℝ trajectory~~ — closed. StepError.lean proves the per-step bound (1 lsb on x1, 3 lsb on x2) and intTraj_ultimate_bound restates the envelope directly about integer state sequences, with no QuantTraj hypothesis.

~~Gain quantization is not covered by ε~~ — closed.

QuantizedGains.lean re-derives the certificate for the gains the circuit holds (40501/65536, 82575/65536). Same P, same ρ = 39/40: the true contraction ratio moves only 3.7·10⁻⁸, well inside the slack. StepError.circuit_ultimate_bound composes it with the per-step bound, so the theorem's gains are the ROM's gains.

~~The mulQ agreement is checked on fixtures, not proved~~ —

closed. Sparkle/Verification/FixedPointProps.lean proves mulQ_toInt: the datapath multiply equals (a·b)/2¹⁶ for all inputs whose result is representable. Mathlib-free, so the RTL and its correctness proof are checked by one lake build.

~~retype pins a newer toolchain than Sparkle~~ — closed. Every

package is on Lean v4.32.1, the retypelab/ sidecar is folded into proofs/SparkleProofs/Retype/, and both the Float falsifier and the Q15.16 transport now import their ℝ models from LQRDesign.lean / PIDDesign.lean instead of restating them. So §12.2's claim — the transported equation is the equation the certificate is about — holds by construction rather than by a drift check.