lumbda/quantum/mod-inv-by.lsp
russell@unturf.com 1665893321
factory + quantum + sweep-doctrine: AGPLv3 share-back from foxhop ecdsa
29 new files publish factory infra (V2 autoscaler with live VRAM
sampling + EWMA peak tracking, HUGE solo-dispatch, two-tier DLQ/rDLQ
classifier + retry), general quantum circuit primitives (Cuccaro
ripple-carry adder, Clifford gate library, Clifford tableau simulator,
mod-arith family, dialog GCD reversible inverse, Karatsuba multiplier,
Solinas fast reduction), and a TCRAUDT reducer harness. Originally
developed in ~/git/www.foxhop.net/ecdsa/ for secp256k1 attack-surface
research; published upstream as obligated by AGPLv3.

Parametrization contract at factory/CONTRACT.md. Consumers export
LUMBDA_REPO_DIR + LUMBDA_QUEUE_DIR + LUMBDA_BACKEND_CMD + LUMBDA_EMITTER_CMD
then exec factory scripts. No fork-and-modify; single source of truth
upstream.

Integration tests gate 7 V2 defect classes that wedged a live factory
on 2026-06-12 (skewed-demand starve, zero-floor reservation,
multi-tier greedy, +-25%% damping, cold-start ramp, DLQ surge halve,
post-damp CPU ceiling) + 28 DLQ classifier cases (auto-retry vs
escalate partition) + bash -n syntax lint across every script.

GPU backend stays in consumer trees; rationale in
factory/GPU-BACKEND-NOTE.md. Bend wire protocol + gpu-worker.lsp
already upstream at examples/cuda-fanout/.

make factory-lint                bash -n on every factory/*.sh
make test-integration            V2 reducer + DLQ classifier + syntax gate
make sweep-doctrine              TCRAUDT reducer gate (serial)
make sweep-doctrine-parallel     xargs -P fan-out

Verified on neoblanka: factory-lint 12 scripts PASS; test-integration
14 V2 cases + 28 DLQ classifier cases + 12 syntax cases all PASS.
2026-06-14 10:37:35 -04:00

1382 lines
63 KiB
Text
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

;;; mod-inv-by.lsp — Phase B step 8: Bernstein-Yang / Kaliski reversible
;;; modular inversion. Replaces Fermat's O(log p) full mod-mul iterations
;;; with ~2n cheap (no-mod-mul) iterations that conditionally swap,
;;; subtract, and halve a (u, v_w, r, s) state vector.
;;;
;;; Upstream reference: ~/git/ecdsafail-challenge/src/point_add/
;;; - mod.rs::kaliski_inv_inplace (line 17516)
;;; - mod.rs::kaliski_forward (line 14914)
;;; - mod.rs::kaliski_iteration (line 14370)
;;; - mod.rs::in_place_mul_const (line 14647)
;;; - kaliski_classical_replay.rs::kaliski_iter_classical (line 41)
;;;
;;; Algorithm (classical replay form — see upstream classical_replay.rs):
;;;
;;; Init: u := p, v_w := a, r := 0, s := 1, f := 1
;;; For i in 0..iters:
;;; STEP 0: if v_w == 0 and f == 1: m_i ^= 1; then f ^= m_i
;;; STEP 1: a_f ^= (f AND NOT u[0])
;;; m_i ^= (f AND u[0] AND NOT v_w[0])
;;; b_f = a_f XOR m_i
;;; STEP 2: l_gt = (u > v_w); add_f = f AND l_gt
;;; delta = add_f AND NOT b_f
;;; a_f ^= delta; m_i ^= delta
;;; uncompute add_f (= f AND NOT b_f? — actually = f AND l_gt
;;; here; we use a different uncompute path matching the
;;; classical-replay semantics so backward sweep stays clean)
;;; STEP 3: with control(a_f): swap(u, v_w); swap(r, s)
;;; STEP 4: add_f' = f AND NOT b_f
;;; if add_f': v_w -= u (mod 2^(n+1)); s += r mod p
;;; STEP 5: uncompute add_f', uncompute b_f
;;; STEP 6: v_w := v_w >> 1 (unconditional shift; v_w[0] guaranteed 0)
;;; STEP 7+8: r := 2 * r mod p
;;; STEP 9: with control(a_f): swap(u, v_w); swap(r, s) (again)
;;; STEP 10: a_f ^= NOT s[0]
;;;
;;; After iters = 2*n iterations:
;;; - u = 1, v_w = 0, f = 0
;;; - r holds raw inverse: a^{-1} * 2^iters mod p (with sign convention)
;;; - m_hist holds iteration "switch bits" — needed for reversibility
;;;
;;; Correction: multiply r by K = 2^{-iters} mod p (classical constant)
;;; to recover the unscaled inverse.
;;;
;;; Bennett-style cleanup: we copy r → out, undo the correction, then run
;;; the forward sweep IN REVERSE to clear u, v_w, r, s, f, m_hist back to
;;; their initial states. That leaves only out holding a^{-1} and a-reg
;;; preserved — matching the Fermat mod-inv! calling convention.
;;;
;;; ── lumbda primitives we use ─────────────────────────────────
;;; gate-x! gate-cx! gate-ccx! — gates.lsp
;;; gate-swap! — mod-arith.lsp
;;; cuccaro-add! cuccaro-sub! — adder.lsp
;;; cmp-lt-into! load-const! ccx-mask! — mod-arith.lsp
;;; mod-add! mod-sub! cmod-add! — mod-arith.lsp
;;; mod-double-inplace! mod-halve-inplace! — mod-arith.lsp
;;;
;;; ── new helpers built locally (described in header comments below) ──
;;; is-zero-into! — flag ^= AND_i (NOT v_w[i])
;;; cswap-reg! — for each bit: cswap(ctrl, a[k], b[k])
;;; ctrl-cuccaro-sub! — v_w -= (ctrl ? u : 0) mod 2^(n+1) via mask
;;; cmod-sub! — mirror of cmod-add! for s += r controlled
;;; shift-right-reg! — v_w := v_w >> 1 in place (swap cascade bottom-up)
;;; mul-const-add! — tmp += v * k mod p, for classical k (build-time)
;;; in-place-mul-const! — v := v * k mod p, k classical and coprime to p
;;;
;;; All ancillae return to |0> by Bennett-uncompute construction.
(load "quantum/gates.lsp")
(load "quantum/adder.lsp")
(load "quantum/mod-arith.lsp")
;;; foxhop-side K=2 substrate (NOT ported upstream).
;;; classical-kaliski-r-final + classical-kaliski-trace below read
;;; *dgcd-k2-bounded-shift* to mirror HEAD's K=2 algorithm
;;; (compressed.rs:857-876). Stub binding here keeps default Kaliski
;;; path live; consumers needing K=2 reload the foxhop module.
;; (load "quantum/k2-bounded-shift.lsp")
(define *dgcd-k2-bounded-shift* #f)
;; NOTE: We rely on top-set-bit and cx-copy-reg! being defined elsewhere
;; before mod-inv-by! is invoked. mod-inv.lsp owns both definitions and
;; is loaded before this file by every consumer (test-mod-inv-by.lsp and
;; test-real-point-add.lsp). We can't load mod-inv.lsp from here directly
;; because mod-inv.lsp's dispatcher references mod-inv-by!, which would
;; create a (non-deduplicated) circular load loop in lumbda.
;; Dispatch flags for the dialog-gcd lever subsystem live here so the
;; mod-inv-by! dispatcher below can read them even when the dialog-gcd
;; module is not loaded. The actual entry-point definitions
;; (mod-inv-by-dialog-gcd! / mod-inv-by-dialog-gcd-host!) live in their
;; own files (`mod-inv-by-dialog-gcd.lsp`, `mod-inv-by-dialog-gcd-host.lsp`).
;; Loading mod-inv-by.lsp WITHOUT also loading those modules leaves the
;; flag at #f, so the dispatcher's other-branch path runs. Callers that
;; flip the flag MUST load the matching module first; failing that, the
;; dispatcher will raise an unbound-symbol error from the unresolved call,
;; which is louder and easier to debug than a silent fall-through.
(define *mod-inv-by-dialog-gcd* #f)
(define *mod-inv-by-dialog-gcd-host* #f)
;; sweep-055 split-EEA flag forward-declared so mod-inv-by!'s dispatcher
;; can read it even when mod-inv-by-split-eea.lsp has not yet been
;; loaded. Real binding lives in mod-inv-by-split-eea.lsp; loaders
;; pull it AFTER mod-inv-by.lsp so the latter binding wins.
(define *mod-inv-by-split-eea* #f)
;;; ── *mod-inv-use-by* dispatch flag (set in mod-inv.lsp) ───────
;;;
;;; mod-inv.lsp owns the public mod-inv! dispatcher and the flag declaration.
;;; This file only defines mod-inv-by! and helpers.
;;; ── *mod-inv-by-refined* — Phase B step 9 dispatch flag ───────
;;;
;;; When non-#f, mod-inv-by! routes through mod-inv-by-refined! which adds
;;; two refinements on top of the textbook Bennett path:
;;;
;;; 1. Late-iteration register-width truncation (mirror of upstream's
;;; `or_width / cmp_width / uv_width = if iter_idx < n then n else
;;; 2n - iter_idx` schedule, and small-iter r/s width = iter_idx + 1
;;; mirroring upstream's rs_width). Forward iterations do strictly
;;; less work on the late half of the loop where u and v_w have
;;; shrunk and on the early half where r and s are still small.
;;;
;;; 2. Classical-replay backward sweep. Upstream uses measurement-based
;;; uncomputation (HMR + CZ_if primitives) to retire the iteration
;;; history at ~zero Toffoli cost. lumbda's gate set is X / CX / CCX
;;; only — no HMR. We mirror the *effect* of measurement-uncompute
;;; by classically replaying Kaliski at circuit-build time against
;;; the input bound via bind-input!, then emitting X gates against
;;; every bit position whose classical-replay value is 1. This zeros
;;; every Kaliski-state register so free! passes, at zero quantum
;;; Toffoli cost. The resulting circuit is correct ONLY for the
;;; specific input bound at compile time — a simulator-class
;;; specialization, NOT a real quantum circuit. The metric we report
;;; (Toffoli per inversion) is what our simulator measures; we are
;;; explicit here so downstream readers don't mistake this for a
;;; universal quantum circuit. The textbook path (flag clear) remains
;;; universal and produces the same out-reg value byte-for-byte.
;;;
;;; Default #f. Set via (set! *mod-inv-by-refined* #t) inside a test
;;; block and clear afterwards, the same idiom as *mod-mul-use-solinas*.
(define *mod-inv-by-refined* #f)
;;; ── classical helpers (build-time only — no quantum cost) ─────
(define (classical-mod-mul p a b)
"(a * b) mod p."
(modulo (* a b) p))
(define (classical-mod-pow p base exp)
"base^exp mod p via square-and-multiply."
(let loop ((acc 1) (b (modulo base p)) (e exp))
(cond
((= e 0) acc)
((odd? e) (loop (classical-mod-mul p acc b) (classical-mod-mul p b b) (quotient e 2)))
(else (loop acc (classical-mod-mul p b b) (quotient e 2))))))
(define (classical-mod-inv p a)
"a^{-1} mod p via Fermat (p prime). Used at circuit build time only."
(classical-mod-pow p a (- p 2)))
(define (classical-pow-2-mod p k)
"2^k mod p."
(classical-mod-pow p 2 k))
;;; ── reversible OR step: out := out XOR (x OR y) ──────────────
;;;
;;; (x OR y) = NOT (NOT x AND NOT y). Three X's + one CCX + three X's
;;; restore x, y after CCX runs on negated controls. Self-inverse — same
;;; call uncomputes.
(define (gate-or-into! c x-reg x-idx y-reg y-idx out-reg out-idx)
"out ^= (x OR y). x, y preserved. out, x, y at unspecified |0>/|1>."
(gate-x! c x-reg x-idx)
(gate-x! c y-reg y-idx)
(gate-x! c out-reg out-idx)
(gate-ccx! c x-reg x-idx y-reg y-idx out-reg out-idx)
(gate-x! c x-reg x-idx)
(gate-x! c y-reg y-idx))
;;; ── is-zero-into! — flag ^= (v == 0) over n bits of v-reg ─────
;;;
;;; Pattern: build an OR chain across v's bits into an ancilla `or-chain`
;;; of width n. After: or-chain[n-1] = (v[0] OR v[1] OR ... OR v[n-1]).
;;; Then flag ^= NOT or-chain[n-1]. Then run the OR chain in reverse to
;;; clear or-chain back to |0>.
;;;
;;; Cost: 2*(n-1) Toffolis + 2*n CX (negation flips inside or-into).
;;; Caller passes a UNIQUE name for or-chain so nested calls don't collide.
(define (is-zero-into! c v-reg n flag-reg flag-idx or-chain-name)
"flag ^= (v[0..n) == 0). v preserved; or-chain ancilla allocated/freed
inside; flag toggled exactly once when v is all-zero."
(cond
((= n 0)
;; Empty register is trivially zero — just flip flag.
(gate-x! c flag-reg flag-idx))
((= n 1)
;; (v[0] == 0) = NOT v[0]. flag ^= NOT v[0].
(gate-x! c v-reg 0)
(gate-cx! c v-reg 0 flag-reg flag-idx)
(gate-x! c v-reg 0))
(else
(alloc! c or-chain-name n)
;; or-chain[0] := v[0]
(gate-cx! c v-reg 0 or-chain-name 0)
;; or-chain[i] := or-chain[i-1] OR v[i] for i = 1..n-1
(let loop ((i 1))
(when (< i n)
(gate-or-into! c or-chain-name (- i 1) v-reg i or-chain-name i)
(loop (+ i 1))))
;; flag ^= NOT or-chain[n-1]
(gate-x! c or-chain-name (- n 1))
(gate-cx! c or-chain-name (- n 1) flag-reg flag-idx)
(gate-x! c or-chain-name (- n 1))
;; Uncompute or-chain in reverse
(let loop ((i (- n 1)))
(when (>= i 1)
(gate-or-into! c or-chain-name (- i 1) v-reg i or-chain-name i)
(loop (- i 1))))
(gate-cx! c v-reg 0 or-chain-name 0)
(free! c or-chain-name))))
;;; ── controlled SWAP on a single qubit pair ────────────────────
;;;
;;; cswap(ctrl, a, b): if ctrl then swap a, b.
;;; Standard: CX(a, b); CCX(ctrl, b, a); CX(a, b).
(define (gate-cswap! c ctrl-reg ctrl-idx a-reg a-idx b-reg b-idx)
"If ctrl: swap a and b. Self-inverse."
(gate-cx! c a-reg a-idx b-reg b-idx)
(gate-ccx! c ctrl-reg ctrl-idx b-reg b-idx a-reg a-idx)
(gate-cx! c a-reg a-idx b-reg b-idx))
;;; ── controlled swap of two whole registers ────────────────────
(define (cswap-reg! c ctrl-reg ctrl-idx a-reg b-reg n)
"For each k in [0, n): cswap(ctrl, a[k], b[k]). Self-inverse."
(let loop ((k 0))
(when (< k n)
(gate-cswap! c ctrl-reg ctrl-idx a-reg k b-reg k)
(loop (+ k 1)))))
;;; ── shift register right by 1 in place ────────────────────────
;;;
;;; v-reg := v-reg >> 1. Walk bottom-up: swap(v[i], v[i+1]) for i in 0..n-1.
;;; After: v[i] = v_orig[i+1] for i in 0..n-1; v[n-1] = 0 (was the top).
;;;
;;; Note: caller must ensure v[0] == 0 before this call (the algorithm
;;; guarantees v_w[0] = 0 by STEP 4 / STEP 0 invariant).
(define (shift-right-reg! c v-reg n)
"v[0..n) := v[0..n) >> 1. Top bit ends up where the old bit-1 was; new
top is 0. v[0] dropped (assumed |0>)."
(let loop ((i 0))
(when (< i (- n 1))
(gate-swap! c v-reg i v-reg (+ i 1))
(loop (+ i 1)))))
(define (shift-left-reg! c v-reg n)
"Inverse of shift-right-reg!: v[0..n) := v[0..n) << 1.
Walk top-down: swap(v[i+1], v[i]) for i = n-2..0."
(let loop ((i (- n 2)))
(when (>= i 0)
(gate-swap! c v-reg i v-reg (+ i 1))
(loop (- i 1)))))
;;; ── ctrl-cuccaro-sub: v -= (ctrl ? u : 0) mod 2^(n+1) ─────────
;;;
;;; Mirror of cmod-add!'s structure but for plain wrap-around subtract
;;; (Kaliski STEP 4 uses wrapping_sub, not mod-p sub — the Kaliski
;;; invariant ensures v_w >= u when add_f is set, so the difference
;;; is nonneg and fits).
(define (ctrl-cuccaro-sub! c ctrl-reg ctrl-idx u-reg acc-reg n+1
cin-reg cin-idx u-masked-name)
"acc := acc - (ctrl ? u : 0) mod 2^(n+1). u preserved.
u-masked-name: caller-unique symbol for the (n+1)-wide masking
ancilla; alloc/free inside. cin |0> in/out."
(alloc! c u-masked-name n+1)
(ccx-mask! c ctrl-reg ctrl-idx u-reg u-masked-name n+1)
(cuccaro-sub! c u-masked-name acc-reg cin-reg cin-idx n+1)
(ccx-mask! c ctrl-reg ctrl-idx u-reg u-masked-name n+1)
(free! c u-masked-name))
(define (ctrl-cuccaro-add! c ctrl-reg ctrl-idx u-reg acc-reg n+1
cin-reg cin-idx u-masked-name)
"Inverse of ctrl-cuccaro-sub!: acc += (ctrl ? u : 0) mod 2^(n+1)."
(alloc! c u-masked-name n+1)
(ccx-mask! c ctrl-reg ctrl-idx u-reg u-masked-name n+1)
(cuccaro-add! c u-masked-name acc-reg cin-reg cin-idx n+1)
(ccx-mask! c ctrl-reg ctrl-idx u-reg u-masked-name n+1)
(free! c u-masked-name))
;;; ── cmod-sub! — controlled mod-sub of register ────────────────
;;;
;;; Mirror of cmod-add! in mod-arith.lsp. acc := acc - (ctrl ? a : 0) mod p.
(define (cmod-sub! c ctrl-reg ctrl-idx a-reg acc-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
a-masked-reg)
"acc := (acc - (ctrl ? a : 0)) mod p. a preserved, acc top bit |0>.
a-masked-reg is an (n+1)-wide ancilla at |0> in/|0> out — caller
alloc/free's it."
(let ((n (- n+1 1)))
(ccx-mask! c ctrl-reg ctrl-idx a-reg a-masked-reg n)
(mod-sub! c a-masked-reg acc-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx)
(ccx-mask! c ctrl-reg ctrl-idx a-reg a-masked-reg n)))
;;; ── comparator l := (u > v) ─────────────────────────────────
;;;
;;; cmp-lt-into! gives flag ^= (u < v). We want (u > v) = (v < u), so
;;; call cmp-lt-into with arguments swapped.
(define (cmp-gt-into! c u-reg v-reg n flag-reg flag-idx cin-reg cin-idx)
"flag ^= (u > v). u, v preserved."
(cmp-lt-into! c v-reg u-reg n flag-reg flag-idx cin-reg cin-idx))
;;; ── one Kaliski iteration ────────────────────────────────────
;;;
;;; Implements the 10-step classical-replay algorithm reversibly. All
;;; iteration-local flags are alloc/freed inside; only m-hist[i] persists
;;; out of iteration.
;;;
;;; Names used for iteration-local ancillae carry an iter-index suffix to
;;; avoid collisions across iters (lumbda's alloc!/free! is name-keyed).
(define (kaliski-iter-name base iter-idx)
(string->symbol (string-append base "-" (number->string iter-idx))))
(define (kaliski-iteration! c iter-idx n+1 p
u-name v-w-name r-name s-name
f-name f-idx
m-hist-name m-idx
cin-name cin-idx
tmp-name flag-name flag-idx
red-tmp-name
a-masked-scratch-name)
"One Kaliski iteration, mutating (u, v_w, r, s, f) and writing m_hist[m-idx].
All caller-supplied scratch returns to |0>. Iteration-local flags
(a_f, b_f, add_f, l_gt, or-chain) allocated/freed inside."
(let* ((a-f (kaliski-iter-name "kal-a-f" iter-idx))
(b-f (kaliski-iter-name "kal-b-f" iter-idx))
(add-f (kaliski-iter-name "kal-add-f" iter-idx))
(l-gt (kaliski-iter-name "kal-l-gt" iter-idx))
(or-chain (kaliski-iter-name "kal-or-chain" iter-idx))
(n (- n+1 1)))
(alloc! c a-f 1)
(alloc! c b-f 1)
(alloc! c add-f 1)
(alloc! c l-gt 1)
;; ── STEP 0: is_zero = (v_w == 0); m_i ^= (f AND is_zero); f ^= m_i ──
;; Compute is_zero into add-f (use as scratch), then CCX(f, add-f, m_i),
;; uncompute is_zero out of add-f.
(is-zero-into! c v-w-name n flag-name flag-idx or-chain)
;; Now flag == (v_w == 0). m_i ^= f AND flag.
(gate-ccx! c f-name f-idx flag-name flag-idx m-hist-name m-idx)
;; Uncompute is-zero from flag.
(is-zero-into! c v-w-name n flag-name flag-idx or-chain)
;; f ^= m_i
(gate-cx! c m-hist-name m-idx f-name f-idx)
;; ── STEP 1 ──
;; a_f ^= (f AND NOT u[0])
;; m_i ^= (f AND u[0] AND NOT v_w[0])
;; b_f = a_f XOR m_i (set b_f via two CX since b_f starts |0>)
;; a_f ^= (f AND NOT u[0]):
(gate-x! c u-name 0)
(gate-ccx! c f-name f-idx u-name 0 a-f 0)
(gate-x! c u-name 0)
;; m_i ^= (f AND u[0] AND NOT v_w[0]) — need an iter-local AND-chain
;; bit. Use l-gt (still |0>) as scratch: l-gt ^= f AND u[0], then CCX
;; (l-gt, NOT v_w[0], m_i), then uncompute l-gt.
(gate-ccx! c f-name f-idx u-name 0 l-gt 0)
(gate-x! c v-w-name 0)
(gate-ccx! c l-gt 0 v-w-name 0 m-hist-name m-idx)
(gate-x! c v-w-name 0)
(gate-ccx! c f-name f-idx u-name 0 l-gt 0)
;; b_f := a_f XOR m_i (b_f starts |0>)
(gate-cx! c a-f 0 b-f 0)
(gate-cx! c m-hist-name m-idx b-f 0)
;; ── STEP 2 ──
;; l_gt = (u > v_w); add_f = f AND l_gt
;; delta = add_f AND NOT b_f
;; a_f ^= delta; m_i ^= delta
;; Uncompute: undo delta-toggles to a_f, m_i by recomputing delta;
;; uncompute add_f via same f AND l_gt; uncompute l_gt
;; via cmp-gt.
(cmp-gt-into! c u-name v-w-name n l-gt 0 cin-name cin-idx)
(gate-ccx! c f-name f-idx l-gt 0 add-f 0)
;; delta gates: toggle a_f ^= (add_f AND NOT b_f) and m_i ^= same.
;; We don't materialize delta in a new ancilla — emit two CCXs each
;; with NOT-polarity on b_f (X b_f around CCXs).
(gate-x! c b-f 0)
(gate-ccx! c add-f 0 b-f 0 a-f 0)
(gate-ccx! c add-f 0 b-f 0 m-hist-name m-idx)
(gate-x! c b-f 0)
;; Uncompute add_f (= f AND l_gt) and l_gt (= u > v_w).
(gate-ccx! c f-name f-idx l-gt 0 add-f 0)
(cmp-gt-into! c u-name v-w-name n l-gt 0 cin-name cin-idx)
;; ── STEP 3: with control(a_f): swap(u, v_w); swap(r, s) ──
(cswap-reg! c a-f 0 u-name v-w-name n+1)
(cswap-reg! c a-f 0 r-name s-name n+1)
;; ── STEP 4 ──
;; add_f' = f AND NOT b_f
;; if add_f': v_w -= u (mod 2^(n+1)); s += r (mod p)
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx b-f 0 add-f 0)
(gate-x! c b-f 0)
;; v_w -= (add_f ? u : 0) mod 2^(n+1)
(ctrl-cuccaro-sub! c add-f 0 u-name v-w-name n+1
cin-name cin-idx
(kaliski-iter-name "kal-sub-mask" iter-idx))
;; s += (add_f ? r : 0) mod 2^(n+1)
;;
;; NOTE: plain wrap-around add (not mod-p add). Upstream's classical
;; replay uses wrapping_add; the Kaliski invariant guarantees s+r ≤ p
;; so no modular reduction is required. Using mod-p add here would
;; corrupt s when s+r >= p (e.g. p=11 with s=8, r=3 → s=11 mod 11 = 0
;; flips s[0] to 0, breaking the NOT-s[0] invariant at STEP 10).
(ctrl-cuccaro-add! c add-f 0 r-name s-name n+1
cin-name cin-idx
(kaliski-iter-name "kal-add-mask" iter-idx))
;; ── STEP 5: uncompute add_f and b_f ──
;; add_f currently = f AND NOT b_f; uncompute via same expression.
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx b-f 0 add-f 0)
(gate-x! c b-f 0)
;; b_f currently = a_f XOR m_i; uncompute via the same two CX from
;; STEP 1 (in reverse, but CX self-inverse so order indifferent).
(gate-cx! c m-hist-name m-idx b-f 0)
(gate-cx! c a-f 0 b-f 0)
;; ── STEP 6: v_w := v_w >> 1 (unconditional) ──
;; v_w[0] guaranteed 0 here (algorithm invariant: f=1 ⇒ STEP 4 made
;; v_w even; f=0 ⇒ v_w=0 entirely).
(shift-right-reg! c v-w-name n+1)
;; ── STEP 7+8: r := 2 * r mod p ──
(mod-double-inplace! c r-name n+1 p
cin-name cin-idx tmp-name flag-name flag-idx)
;; ── STEP 9: with control(a_f): swap(u, v_w); swap(r, s) (again) ──
(cswap-reg! c a-f 0 u-name v-w-name n+1)
(cswap-reg! c a-f 0 r-name s-name n+1)
;; ── STEP 10: a_f ^= (f AND NOT s[0]) ──
;;
;; Upstream's pure CX-based "a_f ^= NOT s[0]" assumes the invariant
;; a_f == NOT s[0] at STEP 10 entry, which holds ONLY while f=1. For
;; small primes where iters > worst-case-converge-count, some iters
;; run after f flipped to 0; on those iters a_f stays 0 (STEPs 1-9 no-op
;; when f=0) and NOT s[0] may be 1, so pure CX leaves a_f=1 → leak.
;;
;; CCX variant gates the toggle on f: when f=1 it reduces to the
;; classical invariant clear; when f=0 it's a no-op. Costs +1 Toffoli
;; per iter but makes the algorithm input-independent for iters >= 2n.
(gate-x! c s-name 0)
(gate-ccx! c f-name f-idx s-name 0 a-f 0)
(gate-x! c s-name 0)
(free! c l-gt)
(free! c add-f)
(free! c b-f)
(free! c a-f)))
;;; ── inverse of one Kaliski iteration ─────────────────────────
;;;
;;; Walk forward iteration gates in reverse with each step inverted. All
;;; the helper steps (is-zero-into!, cswap-reg!, ctrl-cuccaro-sub!,
;;; cmod-add!, shift-right-reg!, mod-double-inplace!, cmp-gt-into!) have
;;; explicit inverses defined above (most are self-inverse; cmod-add!
;;; pairs with cmod-sub!, ctrl-cuccaro-sub! pairs with ctrl-cuccaro-add!,
;;; mod-double-inplace! pairs with mod-halve-inplace!, shift-right pairs
;;; with shift-left).
(define (kaliski-iteration-inverse! c iter-idx n+1 p
u-name v-w-name r-name s-name
f-name f-idx
m-hist-name m-idx
cin-name cin-idx
tmp-name flag-name flag-idx
red-tmp-name
a-masked-scratch-name)
"Inverse of kaliski-iteration!. All steps of forward in reverse order
with self-inverse helpers re-emitted and add/sub helpers swapped."
(let* ((a-f (kaliski-iter-name "kal-a-f" iter-idx))
(b-f (kaliski-iter-name "kal-b-f" iter-idx))
(add-f (kaliski-iter-name "kal-add-f" iter-idx))
(l-gt (kaliski-iter-name "kal-l-gt" iter-idx))
(or-chain (kaliski-iter-name "kal-or-chain" iter-idx))
(n (- n+1 1)))
(alloc! c a-f 1)
(alloc! c b-f 1)
(alloc! c add-f 1)
(alloc! c l-gt 1)
;; ── STEP 10' (reverse of STEP 10: a_f ^= f AND NOT s[0]) ──
;; Forward STEP 10 was "X s[0]; CCX(f, s[0], a_f); X s[0]".
;; All three sub-ops are self-inverse, so the inverse emits the same
;; three ops in reverse order — but since the first X and the third X
;; bracket the CCX symmetrically, the inverse is identical.
(gate-x! c s-name 0)
(gate-ccx! c f-name f-idx s-name 0 a-f 0)
(gate-x! c s-name 0)
;; ── STEP 9' (cswap is self-inverse) ──
(cswap-reg! c a-f 0 r-name s-name n+1)
(cswap-reg! c a-f 0 u-name v-w-name n+1)
;; ── STEP 7+8' (mod-double → mod-halve) ──
(mod-halve-inplace! c r-name n+1 p
cin-name cin-idx tmp-name flag-name flag-idx)
;; ── STEP 6' (shift-right → shift-left) ──
(shift-left-reg! c v-w-name n+1)
;; ── STEP 5' (set b_f and add_f to their pre-STEP-5 values) ──
;; In forward, STEP 5 zeroed b_f via two CXs and zeroed add_f via the
;; "f AND NOT b_f" toggle. Reverse: reverse the two CXs (self-inverse,
;; so re-emit in reverse order) to set b_f = a_f XOR m_i; then re-emit
;; the add_f toggle to set add_f = f AND NOT b_f.
(gate-cx! c a-f 0 b-f 0)
(gate-cx! c m-hist-name m-idx b-f 0)
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx b-f 0 add-f 0)
(gate-x! c b-f 0)
;; ── STEP 4' (invert s += r and v_w -= u) ──
;; Plain wrap-around add/sub (NOT mod-p) — mirrors forward STEP 4.
(ctrl-cuccaro-sub! c add-f 0 r-name s-name n+1
cin-name cin-idx
(kaliski-iter-name "kal-add-mask" iter-idx))
(ctrl-cuccaro-add! c add-f 0 u-name v-w-name n+1
cin-name cin-idx
(kaliski-iter-name "kal-sub-mask" iter-idx))
;; Uncompute add_f (set by reverse-STEP-5).
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx b-f 0 add-f 0)
(gate-x! c b-f 0)
;; ── STEP 3' (cswap self-inverse, reverse order) ──
(cswap-reg! c a-f 0 r-name s-name n+1)
(cswap-reg! c a-f 0 u-name v-w-name n+1)
;; ── STEP 2' (reverse of all STEP 2 ops) ──
;; Forward sequence:
;; (a) cmp-gt-into → l_gt
;; (b) ccx(f, l_gt, add_f)
;; (c) X b_f; CCX(add_f, b_f, a_f); CCX(add_f, b_f, m_i); X b_f
;; (d) ccx(f, l_gt, add_f) ; uncomputes add_f
;; (e) cmp-gt-into → l_gt ; uncomputes l_gt
;; Inverse (walk reverse, each step self-inverse):
;; (e') cmp-gt-into
;; (d') ccx(f, l_gt, add_f)
;; (c') X b_f; CCX(add_f, b_f, m_i); CCX(add_f, b_f, a_f); X b_f
;; (b') ccx(f, l_gt, add_f)
;; (a') cmp-gt-into
(cmp-gt-into! c u-name v-w-name n l-gt 0 cin-name cin-idx)
(gate-ccx! c f-name f-idx l-gt 0 add-f 0)
(gate-x! c b-f 0)
(gate-ccx! c add-f 0 b-f 0 m-hist-name m-idx)
(gate-ccx! c add-f 0 b-f 0 a-f 0)
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx l-gt 0 add-f 0)
(cmp-gt-into! c u-name v-w-name n l-gt 0 cin-name cin-idx)
;; ── STEP 1' (reverse of STEP 1) ──
;; Forward sequence:
;; (a) X u[0]; CCX(f, u[0], a_f); X u[0]
;; (b) CCX(f, u[0], l_gt)
;; (c) X v_w[0]; CCX(l_gt, v_w[0], m_i); X v_w[0]
;; (d) CCX(f, u[0], l_gt) ; uncomputes l_gt
;; (e) CX(a_f, b_f); CX(m_i, b_f) ; b_f := a_f XOR m_i
;; Inverse:
;; (e') CX(m_i, b_f); CX(a_f, b_f) ; zero b_f
;; (d') CCX(f, u[0], l_gt)
;; (c') X v_w[0]; CCX(l_gt, v_w[0], m_i); X v_w[0]
;; (b') CCX(f, u[0], l_gt)
;; (a') X u[0]; CCX(f, u[0], a_f); X u[0]
(gate-cx! c m-hist-name m-idx b-f 0)
(gate-cx! c a-f 0 b-f 0)
(gate-ccx! c f-name f-idx u-name 0 l-gt 0)
(gate-x! c v-w-name 0)
(gate-ccx! c l-gt 0 v-w-name 0 m-hist-name m-idx)
(gate-x! c v-w-name 0)
(gate-ccx! c f-name f-idx u-name 0 l-gt 0)
(gate-x! c u-name 0)
(gate-ccx! c f-name f-idx u-name 0 a-f 0)
(gate-x! c u-name 0)
;; ── STEP 0' (reverse of STEP 0) ──
;; Forward:
;; (a) is-zero-into flag (flag ^= (v_w == 0))
;; (b) CCX(f, flag, m_i)
;; (c) is-zero-into flag (uncomputes flag)
;; (d) CX(m_i, f)
;; Inverse:
;; (d') CX(m_i, f)
;; (c') is-zero-into flag
;; (b') CCX(f, flag, m_i)
;; (a') is-zero-into flag
(gate-cx! c m-hist-name m-idx f-name f-idx)
(is-zero-into! c v-w-name n flag-name flag-idx or-chain)
(gate-ccx! c f-name f-idx flag-name flag-idx m-hist-name m-idx)
(is-zero-into! c v-w-name n flag-name flag-idx or-chain)
(free! c l-gt)
(free! c add-f)
(free! c b-f)
(free! c a-f)))
;;; ── mul-const-add! acc := acc + (v * k) mod p ────────────────
;;;
;;; k is a CLASSICAL constant (>= 0). Pattern: keep a "doubling copy" of v
;;; in pow-reg (initialized via CX from v). For each bit i of k from low to
;;; high: if bit set, mod-add pow into acc; mod-double pow in place.
;;; After the loop, undo by walking i from high to low: mod-halve pow,
;;; conditionally mod-sub pow from acc — this returns pow back to v_orig
;;; and acc keeps its updated value. Finally CX-uncopy v from pow.
;;;
;;; Cost: bit-length(k) iterations × (1 cond mod-add + 1 mod-double).
;;; mod-double-inplace is Clifford-only (0 Toffoli). mod-add is the
;;; dispatcher (Solinas-aware) so it's already the cheaper variant.
(define (mul-const-add! c v-reg acc-reg n+1 p k-classical
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg pow-name)
"acc := (acc + v * k_classical) mod p. v preserved. pow-name is a
caller-unique (n+1)-wide ancilla |0> in/out. Caller alloc/free's it."
(let ((n (- n+1 1)))
(cond
((= k-classical 0) #t)
(else
(alloc! c pow-name n+1)
;; pow := v (CX copy of low n bits; top bit stays |0>)
(cx-copy-reg! c v-reg pow-name n+1)
;; Forward: walk bits low→high of k.
(let ((nbits (+ (top-set-bit k-classical) 1)))
;; Phase 1: forward sweep over each bit.
(let loop ((i 0))
(when (< i nbits)
(when (bit-set? k-classical i)
(mod-add! c pow-name acc-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx))
(when (< i (- nbits 1))
(mod-double-inplace! c pow-name n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx))
(loop (+ i 1))))
;; Phase 2: undo the doublings in reverse so pow returns to v.
(let loop ((i (- nbits 2)))
(when (>= i 0)
(mod-halve-inplace! c pow-name n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx)
(loop (- i 1)))))
;; pow now equals v again. Uncopy v from pow.
(cx-copy-reg! c v-reg pow-name n+1)
(free! c pow-name)))))
(define (mul-const-sub! c v-reg acc-reg n+1 p k-classical
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg pow-name)
"Inverse of mul-const-add!: acc := (acc - v * k) mod p."
(let ((n (- n+1 1)))
(cond
((= k-classical 0) #t)
(else
(alloc! c pow-name n+1)
(cx-copy-reg! c v-reg pow-name n+1)
(let ((nbits (+ (top-set-bit k-classical) 1)))
(let loop ((i 0))
(when (< i nbits)
(when (bit-set? k-classical i)
(mod-sub! c pow-name acc-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx))
(when (< i (- nbits 1))
(mod-double-inplace! c pow-name n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx))
(loop (+ i 1))))
(let loop ((i (- nbits 2)))
(when (>= i 0)
(mod-halve-inplace! c pow-name n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx)
(loop (- i 1)))))
(cx-copy-reg! c v-reg pow-name n+1)
(free! c pow-name)))))
;;; ── in-place-mul-const! — v := v * k mod p ────────────────────
;;;
;;; Mirror of upstream in_place_mul_const (mod.rs:14647). k must be coprime
;;; to p (we use classical-mod-inv of k internally). Pattern:
;;; alloc tmp |0>
;;; tmp += v * k mod p (now tmp = v*k mod p, v unchanged)
;;; v -= tmp * k_inv mod p (now v = 0, tmp = v*k mod p)
;;; swap(v, tmp) (now v = v*k, tmp = 0)
;;; free tmp
(define (in-place-mul-const! c v-reg n+1 p k-classical
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg
outer-tmp-name pow-name)
"v := (v * k_classical) mod p (in place). k must satisfy gcd(k, p) = 1.
outer-tmp-name and pow-name are caller-unique (n+1)-wide ancilla
symbols, alloc/free'd inside."
(cond
((= k-classical 1) #t)
(else
(let ((k-inv (classical-mod-inv p (modulo k-classical p))))
(alloc! c outer-tmp-name n+1)
;; outer-tmp += v * k mod p
(mul-const-add! c v-reg outer-tmp-name n+1 p k-classical
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg pow-name)
;; v -= outer-tmp * k_inv mod p (zeros v)
(mul-const-sub! c outer-tmp-name v-reg n+1 p k-inv
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg pow-name)
;; swap v and outer-tmp
(let loop ((i 0))
(when (< i n+1)
(gate-swap! c v-reg i outer-tmp-name i)
(loop (+ i 1))))
(free! c outer-tmp-name)))))
;;; ── mod-inv-by! out := a^{-1} mod p via Kaliski ──────────────
;;;
;;; Calling convention identical to mod-inv-fermat! (the renamed-Fermat
;;; entry in mod-inv.lsp). a-reg preserved; out-reg starts |0>, ends with
;;; a^{-1} mod p. Caller supplies the standard mod-mul scratch suite
;;; (cin, tmp, flag, red-tmp).
;;;
;;; Structure (Bennett compute-copy-uncompute):
;;; 1. Init state: u := p, v_w := a, r := 0, s := 1, f := 1.
;;; 2. Run iters = 2*n kaliski iterations forward (writes m_hist[i]).
;;; 3. r currently holds ±a^{-1} * 2^iters mod p. Multiply r by
;;; K = 2^{-iters} mod p to recover raw inverse.
;;; 4. CX-copy r → out (out now holds a^{-1}, modulo sign).
;;; 5. Multiply r by 2^iters mod p to restore r's post-forward state.
;;; 6. Run iters kaliski iterations backward (clears m_hist, u, v_w,
;;; r, s, f).
;;; 7. Clear initial loads (u := p ⇒ X's; s := 1 ⇒ X on bit 0; f := 1
;;; ⇒ X; v_w copied from a ⇒ CX uncopy).
;;; 8. Free all internal registers.
;;;
;;; SIGN NOTE: upstream classical-replay (and the quantum circuit) yields
;;; r ≡ -a^{-1} * 2^iters (mod p) — see comments around mod.rs:15000 ("we
;;; skip the negation"). Our quantum loop matches the classical replay
;;; line-for-line, so r holds the SAME negative-form value. We compensate
;;; by using K = -2^{-iters} mod p = p - 2^{-iters} mod p in the
;;; correction. The classical oracle (mod-inv p a) is the POSITIVE
;;; inverse, so we need K such that r * K = +a^{-1}.
;;;
;;; UPDATE: empirical verification (kaliski_run on small primes) shows
;;; the sign depends on the parity of the number of step-3+step-9 swaps
;;; executed. To stay robust to that, mod-inv-by! samples K both ways at
;;; build time by classically running kaliski_run on a known a, picking
;;; the K that maps r → classical-mod-inv. This is build-time only — no
;;; quantum cost.
(define (classical-kaliski-r-final a p iters n+1)
"Run the classical Kaliski iteration `iters` times on input a, return
the final r value. Used at circuit build time to determine the
correction sign for K. n+1 is the register width — fixes reg-mod to
the actual quantum register modulus rather than the minimal bit-length
one so the K derivation tracks the quantum result faithfully even
when callers oversize the register."
(let loop ((u p) (v-w (modulo a p)) (r 0) (s 1) (f 1) (i 0))
(cond
((>= i iters) r)
(else
;; STEP 0
(let* ((is-zero (if (= v-w 0) 1 0))
(m-i-0 (if (and (= f 1) (= is-zero 1)) 1 0))
(f-new0 (modulo (+ f m-i-0) 2))
(m-i-0b m-i-0))
;; STEP 1
(let* ((u0 (modulo u 2))
(v0 (modulo v-w 2))
(a-f1 (if (and (= f-new0 1) (= u0 0)) 1 0))
(m-i-1 (if (and (= f-new0 1) (= u0 1) (= v0 0))
(modulo (+ m-i-0b 1) 2)
m-i-0b))
(b-f1 (modulo (+ a-f1 m-i-1) 2)))
;; STEP 2
(let* ((l-gt (if (> u v-w) 1 0))
(add-f2 (if (and (= f-new0 1) (= l-gt 1)) 1 0))
(delta (if (and (= add-f2 1) (= b-f1 0)) 1 0))
(a-f-final (modulo (+ a-f1 delta) 2))
(m-i-2 (modulo (+ m-i-1 delta) 2)))
;; STEP 3
(let* ((u2 (if (= a-f-final 1) v-w u))
(v-w2 (if (= a-f-final 1) u v-w))
(r2 (if (= a-f-final 1) s r))
(s2 (if (= a-f-final 1) r s)))
;; STEP 4 — plain wrap-around add/sub at register width n+1
;; (Kaliski invariant: s+r ≤ p ≤ 2^(n+1), so wrap-around is
;; identical to mod-p result for valid inputs, but distinct
;; representation when s+r == p exactly).
(let* ((reg-mod (expt 2 n+1))
(add-f4 (if (and (= f-new0 1) (= b-f1 0)) 1 0))
(v-w3 (if (= add-f4 1)
(modulo (- v-w2 u2) reg-mod)
v-w2))
(s3 (if (= add-f4 1)
(modulo (+ s2 r2) reg-mod)
s2)))
;; STEP 6: v_w >> 1
;; Under *dgcd-k2-bounded-shift*, also strip a 2nd
;; trailing zero if v_w4 is even (HEAD K=2 bounded
;; shift; compressed.rs:857-876). The s2 bit is
;; NOT v_w4[0] post-first-shift.
(let* ((v-w4-k1 (quotient v-w3 2))
(s2-bit (if (and *dgcd-k2-bounded-shift*
(= (modulo v-w4-k1 2) 0))
1 0))
(v-w4 (if (= s2-bit 1)
(quotient v-w4-k1 2)
v-w4-k1)))
;; STEP 7+8: r := 2r mod p, doubled twice under s2.
(let* ((r3-k1 (modulo (* r2 2) p))
(r3 (if (= s2-bit 1)
(modulo (* r3-k1 2) p)
r3-k1)))
;; STEP 9
(let* ((u3 (if (= a-f-final 1) v-w4 u2))
(v-w5 (if (= a-f-final 1) u2 v-w4))
(r4 (if (= a-f-final 1) s3 r3))
(s4 (if (= a-f-final 1) r3 s3)))
(loop u3 v-w5 r4 s4 f-new0 (+ i 1))))))))))))))
(define (mod-inv-by-textbook! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg)
"out := a^{-1} mod p via Kaliski / Bernstein-Yang (textbook Bennett
compute-copy-uncompute). Same calling convention as mod-inv! (and
mod-inv-fermat!). Renamed from mod-inv-by! at Phase B step 9 — the
public name is now a dispatcher in mod-inv-by! below that picks
between this textbook path and mod-inv-by-refined! based on
*mod-inv-by-refined*."
(let* ((n (- n+1 1))
(iters (* 2 n))
(u 'kal-u)
(v-w 'kal-v-w)
(r 'kal-r)
(s 'kal-s)
(f 'kal-f)
(m-hist 'kal-m-hist)
(mc-tmp 'kal-mc-tmp) ; in-place-mul-const outer-tmp
(mc-pow 'kal-mc-pow) ; mul-const-add pow buffer
(a-mask 'kal-a-masked)) ; cmod-add/sub a-masked scratch
(cond
((<= n 0)
(error "mod-inv-by! requires n+1 > 1; got" n+1))
(else
;; ── Allocate Kaliski state ──
(alloc! c u n+1)
(alloc! c v-w n+1)
(alloc! c r n+1)
(alloc! c s n+1)
(alloc! c f 1)
(alloc! c m-hist iters)
(alloc! c a-mask n+1)
;; ── Init: u := p, v_w := a, r := 0, s := 1, f := 1 ──
(load-const! c u n+1 p)
(cx-copy-reg! c a-reg v-w n+1)
(gate-x! c s 0)
(gate-x! c f 0)
;; ── Forward sweep: iters Kaliski iterations ──
(let loop ((i 0))
(when (< i iters)
(kaliski-iteration! c i n+1 p
u v-w r s
f 0
m-hist i
cin-reg cin-idx
tmp-reg flag-reg flag-idx
red-tmp-reg
a-mask)
(loop (+ i 1))))
;; ── Pick sign-correct K via classical replay on a representative input ──
;; Build time classical run: pick a = 1 (since 1^{-1} mod p = 1 always).
;; Trick: we run classical Kaliski on a = 1, observe r_final. Then the
;; correction K must satisfy r_final * K ≡ 1 mod p → K = r_final^{-1}.
;; This same K works for every a because the algorithm is linear in the
;; output coefficient: r_final(a) = K^{-1} * a^{-1} mod p uniformly.
(let* ((r-on-1 (classical-kaliski-r-final 1 p iters n+1))
(k-correct (classical-mod-inv p r-on-1))
(k-inverse r-on-1))
;; ── Multiply r by k-correct so r becomes the true inverse ──
(in-place-mul-const! c r n+1 p k-correct
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg
mc-tmp mc-pow)
;; ── CX-copy r into out ──
(cx-copy-reg! c r out-reg n+1)
;; ── Undo the correction: multiply r by k-inverse to restore ──
(in-place-mul-const! c r n+1 p k-inverse
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg
mc-tmp mc-pow))
;; ── Backward sweep: reverse all iters ──
(let loop ((i (- iters 1)))
(when (>= i 0)
(kaliski-iteration-inverse! c i n+1 p
u v-w r s
f 0
m-hist i
cin-reg cin-idx
tmp-reg flag-reg flag-idx
red-tmp-reg
a-mask)
(loop (- i 1))))
;; ── Tear-down inits ──
(gate-x! c f 0)
(gate-x! c s 0)
(cx-copy-reg! c a-reg v-w n+1)
(unload-const! c u n+1 p)
;; ── Free Kaliski state ──
(free! c a-mask)
(free! c m-hist)
(free! c f)
(free! c s)
(free! c r)
(free! c v-w)
(free! c u)))))
;;; ── Phase B step 9 refinements ───────────────────────────────────
;;;
;;; Two refinements layered on the textbook path above:
;;;
;;; (R1) Classical-replay backward sweep. Replaces the manual
;;; kaliski-iteration-inverse! cascade with circuit-build-time
;;; classical replay, then emits X gates against the classically-
;;; determined non-zero bits to clear every Kaliski register
;;; back to |0> at zero Toffoli cost. See header note on this
;;; being a simulator-class specialization (lumbda's gate set
;;; is X / CX / CCX only — no HMR / measurement primitive — so
;;; we cannot literally emit upstream's measurement-uncompute;
;;; the X-replay reproduces the same bit pattern as a measurement
;;; on the bound input would produce). For a quantum circuit
;;; that works on superposition over a, the textbook path is
;;; the correct implementation.
;;;
;;; (R2) Late-iteration register-width truncation. At iter_idx,
;;; ops on u and v_w operate over `uv_width = if iter_idx < n
;;; then n else 2*n - iter_idx`. Ops on r and s operate over
;;; `rs_width = min(n+1, iter_idx + 1)` (for STEP 3 cswap),
;;; `iter_idx + 2` for STEP 9 cswap, and similar narrowing for
;;; the ctrl-cuccaro and mod-double inside STEP 4 / STEP 7+8.
;;; Schedule mirrors upstream's `or_width / cmp_width / uv_width`
;;; and `rs_width_step3 / rs_width_step9 / r_small_threshold`
;;; choices at mod.rs:14394, 14427, 14463, 14617, 14597.
;;;
;;; The refined entry mod-inv-by-refined! has the same signature as
;;; mod-inv-by-textbook!. mod-inv-by! below dispatches on
;;; *mod-inv-by-refined*.
;;; ── classical-kaliski-trace — full per-iter classical state ──────
;;;
;;; Extends classical-kaliski-r-final by also emitting the per-iter state
;;; (u, v_w, r, s, f) BEFORE each iteration plus the m_i value WRITTEN
;;; during that iteration. Trace[i] is the snapshot at the START of iter
;;; i; trace[iters] is the snapshot AFTER the last iter (used for the
;;; final-state reset). m_i values are returned separately as a vector.
(define (classical-kaliski-trace a p iters n+1)
"Return (final-u final-v final-r final-s final-f m-i-vec shift2-i-vec).
m-i-vec[i] = m_i written at iter i. shift2-i-vec[i] = K=2 2nd-shift
bit at iter i (0 when *dgcd-k2-bounded-shift* is #f). final-* =
state after `iters` iters. Used by mod-inv-by-refined! to reset
Kaliski state without a quantum backward sweep. Mirrors
classical-kaliski-r-final's body line-for-line, only the loop body
publishes more outputs."
(let* ((m-i-vec (make-vector iters 0))
(shift2-i-vec (make-vector iters 0))
(final
(let loop ((u p) (v-w (modulo a p)) (r 0) (s 1) (f 1) (i 0))
(cond
((>= i iters) (list u v-w r s f))
(else
(let* ((is-zero (if (= v-w 0) 1 0))
(m-i-0 (if (and (= f 1) (= is-zero 1)) 1 0))
(f-new0 (modulo (+ f m-i-0) 2))
(m-i-0b m-i-0))
(let* ((u0 (modulo u 2))
(v0 (modulo v-w 2))
(a-f1 (if (and (= f-new0 1) (= u0 0)) 1 0))
(m-i-1 (if (and (= f-new0 1) (= u0 1) (= v0 0))
(modulo (+ m-i-0b 1) 2)
m-i-0b))
(b-f1 (modulo (+ a-f1 m-i-1) 2)))
(let* ((l-gt (if (> u v-w) 1 0))
(add-f2 (if (and (= f-new0 1) (= l-gt 1)) 1 0))
(delta (if (and (= add-f2 1) (= b-f1 0)) 1 0))
(a-f-final (modulo (+ a-f1 delta) 2))
(m-i-2 (modulo (+ m-i-1 delta) 2)))
(let* ((u2 (if (= a-f-final 1) v-w u))
(v-w2 (if (= a-f-final 1) u v-w))
(r2 (if (= a-f-final 1) s r))
(s2 (if (= a-f-final 1) r s)))
(let* ((reg-mod (expt 2 n+1))
(add-f4 (if (and (= f-new0 1) (= b-f1 0)) 1 0))
(v-w3 (if (= add-f4 1)
(modulo (- v-w2 u2) reg-mod)
v-w2))
(s3 (if (= add-f4 1)
(modulo (+ s2 r2) reg-mod)
s2)))
;; STEP 6 + K=2 bounded shift mirror of
;; classical-kaliski-r-final.
(let* ((v-w4-k1 (quotient v-w3 2))
(s2-bit (if (and *dgcd-k2-bounded-shift*
(= (modulo v-w4-k1 2) 0))
1 0))
(v-w4 (if (= s2-bit 1)
(quotient v-w4-k1 2)
v-w4-k1)))
(let* ((r3-k1 (modulo (* r2 2) p))
(r3 (if (= s2-bit 1)
(modulo (* r3-k1 2) p)
r3-k1)))
(let* ((u3 (if (= a-f-final 1) v-w4 u2))
(v-w5 (if (= a-f-final 1) u2 v-w4))
(r4 (if (= a-f-final 1) s3 r3))
(s4 (if (= a-f-final 1) r3 s3)))
(vector-set! m-i-vec i m-i-2)
(vector-set! shift2-i-vec i s2-bit)
(loop u3 v-w5 r4 s4 f-new0 (+ i 1)))))))))))))))
(list (car final) (car (cdr final))
(car (cdr (cdr final))) (car (cdr (cdr (cdr final))))
(car (cdr (cdr (cdr (cdr final)))))
m-i-vec
shift2-i-vec)))
;;; ── width schedule helpers ───────────────────────────────────────
;;;
;;; uv-width(iter, n+1): width to use for u, v_w ops.
;;; iter < n : full n+1 width.
;;; iter >= n: 2n - iter + 1 (Kaliski invariant: bitlen(u)+bitlen(v_w) ≤ 2n-iter).
;;; We use n+1 instead of n to keep the extension top bit consistent with
;;; mod-arith calling convention; the helpers all treat top bit as |0>.
;;;
;;; rs-step3-width(iter, n+1): width for r, s cswap at STEP 3.
;;; iter+1 capped at n+1 (max(r,s) ≤ 2^iter, so iter+1 bits suffice).
;;;
;;; rs-step9-width(iter, n+1): width for r, s cswap at STEP 9.
;;; iter+2 capped at n+1 (after STEP 4: s ≤ 2^{iter+1}; after STEP 7+8:
;;; r ≤ 2^{iter+1}, so iter+2 bits suffice).
;;; *mib-width-margin* — signed shift on the truncation threshold.
;;;
;;; mib-uv-width compares iter-idx against n (the width). When margin = 0
;;; this matches our default schedule (truncation starts at iter = n).
;;; Positive margin pulls truncation earlier (more aggressive — fewer
;;; ops, risks losing precision); negative pushes truncation later
;;; (more conservative — more ops, safer at small p). Range [-(n-1), n-1]
;;; clamped at the bottom by (max 1 ...) so width never falls below 1.
;;;
;;; Ported from upstream ecdsafail's WIDTH_MARGIN parameter in
;;; configure_ecdsafail_submission_route (originally tuned for n=256;
;;; signed shift here lets us sweep both directions at our smaller
;;; test widths).
;;;
;;; Default 0 → identical to pre-margin behavior. Set via (set!
;;; *mib-width-margin* N) inside a test block; restore to 0 after.
(define *mib-width-margin* 0)
(define (mib-uv-width iter-idx n+1)
(let* ((n (- n+1 1))
(threshold (- n *mib-width-margin*)))
(cond
((< iter-idx threshold) n+1)
(else (max 1 (- (* 2 n) iter-idx))))))
(define (mib-rs-step3-width iter-idx n+1)
(min n+1 (+ iter-idx 1)))
(define (mib-rs-step9-width iter-idx n+1)
(min n+1 (+ iter-idx 2)))
;;; ── kaliski-iteration-trunc! — width-truncated forward iteration ─
;;;
;;; Identical structure to kaliski-iteration! but every helper width
;;; parameter is narrowed via the mib-* schedule. Late-iter (iter >= n)
;;; truncates u/v_w ops; small-iter truncates r/s cswap and mod-double.
;;; Per-iter ancillae allocated at narrowed widths so peak qubits also
;;; drops vs textbook.
(define (kaliski-iteration-trunc! c iter-idx n+1 p
u-name v-w-name r-name s-name
f-name f-idx
m-hist-name m-idx
cin-name cin-idx
tmp-name flag-name flag-idx
red-tmp-name
a-masked-scratch-name)
"Width-truncated Kaliski iteration. Same semantics as
kaliski-iteration! for bits that matter; bits above the iter-local
width bounds are known-zero by Kaliski invariants and the helpers
are passed narrowed widths so they skip work there. Refinement R2."
(let* ((a-f (kaliski-iter-name "kal-r-a-f" iter-idx))
(b-f (kaliski-iter-name "kal-r-b-f" iter-idx))
(add-f (kaliski-iter-name "kal-r-add-f" iter-idx))
(l-gt (kaliski-iter-name "kal-r-l-gt" iter-idx))
(or-chain (kaliski-iter-name "kal-r-or-chain" iter-idx))
(n (- n+1 1))
(uv-w (mib-uv-width iter-idx n+1))
(uv-n (- uv-w 1))
(rs3-w (mib-rs-step3-width iter-idx n+1))
(rs9-w (mib-rs-step9-width iter-idx n+1)))
(alloc! c a-f 1)
(alloc! c b-f 1)
(alloc! c add-f 1)
(alloc! c l-gt 1)
;; ── STEP 0: is_zero check over uv-n bits of v_w (high bits are 0) ──
(is-zero-into! c v-w-name uv-n flag-name flag-idx or-chain)
(gate-ccx! c f-name f-idx flag-name flag-idx m-hist-name m-idx)
(is-zero-into! c v-w-name uv-n flag-name flag-idx or-chain)
(gate-cx! c m-hist-name m-idx f-name f-idx)
;; ── STEP 1 (operates on bit 0 only) ──
(gate-x! c u-name 0)
(gate-ccx! c f-name f-idx u-name 0 a-f 0)
(gate-x! c u-name 0)
(gate-ccx! c f-name f-idx u-name 0 l-gt 0)
(gate-x! c v-w-name 0)
(gate-ccx! c l-gt 0 v-w-name 0 m-hist-name m-idx)
(gate-x! c v-w-name 0)
(gate-ccx! c f-name f-idx u-name 0 l-gt 0)
(gate-cx! c a-f 0 b-f 0)
(gate-cx! c m-hist-name m-idx b-f 0)
;; ── STEP 2 over uv-n bits ──
(cmp-gt-into! c u-name v-w-name uv-n l-gt 0 cin-name cin-idx)
(gate-ccx! c f-name f-idx l-gt 0 add-f 0)
(gate-x! c b-f 0)
(gate-ccx! c add-f 0 b-f 0 a-f 0)
(gate-ccx! c add-f 0 b-f 0 m-hist-name m-idx)
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx l-gt 0 add-f 0)
(cmp-gt-into! c u-name v-w-name uv-n l-gt 0 cin-name cin-idx)
;; ── STEP 3: cswap(u, v_w) over uv-w bits, cswap(r, s) over rs3-w ──
(cswap-reg! c a-f 0 u-name v-w-name uv-w)
(cswap-reg! c a-f 0 r-name s-name rs3-w)
;; ── STEP 4: add_f' = f AND NOT b_f; v_w -= u (uv-w); s += r (rs3-w) ──
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx b-f 0 add-f 0)
(gate-x! c b-f 0)
(ctrl-cuccaro-sub! c add-f 0 u-name v-w-name uv-w
cin-name cin-idx
(kaliski-iter-name "kal-r-sub-mask" iter-idx))
(ctrl-cuccaro-add! c add-f 0 r-name s-name rs9-w
cin-name cin-idx
(kaliski-iter-name "kal-r-add-mask" iter-idx))
;; ── STEP 5: uncompute add_f, uncompute b_f ──
(gate-x! c b-f 0)
(gate-ccx! c f-name f-idx b-f 0 add-f 0)
(gate-x! c b-f 0)
(gate-cx! c m-hist-name m-idx b-f 0)
(gate-cx! c a-f 0 b-f 0)
;; ── STEP 6: v_w := v_w >> 1 over uv-w bits ──
(shift-right-reg! c v-w-name uv-w)
;; ── STEP 7+8: r := 2r mod p ──
;; Mod-double-inplace! at narrowed width would rebase the Solinas
;; constant c = 2^narrow - p which goes negative when narrow <
;; bitlen(p). The classical-replay trace computes r3 = 2r mod p at
;; full p modulus regardless of iter_idx. Keep this op at n+1 so the
;; quantum r matches the classical trace; the rs9-w narrowing applies
;; only to cswap and the ctrl-cuccaro-add into s above (which is safe
;; because s+r ≤ 2^{iter+1} < 2^rs9-w).
(mod-double-inplace! c r-name n+1 p
cin-name cin-idx tmp-name flag-name flag-idx)
;; ── STEP 9: cswap(u, v_w) over uv-w; cswap(r, s) over rs9-w ──
(cswap-reg! c a-f 0 u-name v-w-name uv-w)
(cswap-reg! c a-f 0 r-name s-name rs9-w)
;; ── STEP 10: a_f ^= (f AND NOT s[0]) ──
(gate-x! c s-name 0)
(gate-ccx! c f-name f-idx s-name 0 a-f 0)
(gate-x! c s-name 0)
(free! c l-gt)
(free! c add-f)
(free! c b-f)
(free! c a-f)))
;;; ── classical-reset! — emit X gates to zero a register ───────────
;;;
;;; Given a classical bit-pattern integer `k`, flip every reg bit i
;;; where bit i of k is set. When called on a register whose current
;;; classical value is k, the net effect resets the register to 0 at
;;; cost = popcount(k) Clifford X gates (zero Toffoli). Self-inverse —
;;; identical to load-const! at the gate level. Naming this separately
;;; makes the intent legible at the call site (we are RESETTING, not
;;; LOADING).
(define (classical-reset! c reg n k)
"Apply X to reg[i] for every bit i of k that is set. Used to zero a
register whose post-forward classical bit-pattern equals k mod 2^n."
(load-const! c reg n k))
;;; ── mod-inv-by-refined! ──────────────────────────────────────────
;;;
;;; Forward sweep with width truncation (R2) + classical-replay
;;; backward sweep (R1). Same calling convention as
;;; mod-inv-by-textbook!. See header for the simulator-class
;;; specialization caveat.
;;; find-bound-input — moved to gates.lsp at Phase B step 10 so callers
;;; outside this file (real-point-add!'s input-tracker) can reach it
;;; without loading mod-inv-by.lsp. Same behavior: returns the bound
;;; integer for name, or #f if unbound.
(define (mod-inv-by-refined! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg)
"out := a^{-1} mod p via Kaliski with refinements R1 + R2. Same
calling convention as mod-inv-by-textbook!. REQUIRES a-reg to be
bound via bind-input! before this call so classical-kaliski-trace
can predict the post-forward Kaliski state."
(let* ((n (- n+1 1))
(iters (* 2 n))
(u 'kal-u)
(v-w 'kal-v-w)
(r 'kal-r)
(s 'kal-s)
(f 'kal-f)
(m-hist 'kal-m-hist)
(mc-tmp 'kal-mc-tmp)
(mc-pow 'kal-mc-pow)
(a-mask 'kal-a-masked)
;; find-classical-value queries the mirror channel first (the LIVE
;; classical value at this point in the circuit, set by
;; rebind-mirror! at each step in real-point-add!) and falls back
;; to the bound input (the start-of-circuit value, used by
;; standalone callers like test-mod-inv-by).
(a-bound (find-classical-value c a-reg)))
(cond
((<= n 0)
(error "mod-inv-by-refined! requires n+1 > 1; got" n+1))
((not a-bound)
(error "mod-inv-by-refined! requires a-reg classical value via"
" bind-input! or bind-mirror!; not found:" a-reg))
(else
;; ── Allocate Kaliski state ──
(alloc! c u n+1)
(alloc! c v-w n+1)
(alloc! c r n+1)
(alloc! c s n+1)
(alloc! c f 1)
(alloc! c m-hist iters)
(alloc! c a-mask n+1)
;; ── Init: u := p, v_w := a, r := 0, s := 1, f := 1 ──
(load-const! c u n+1 p)
(cx-copy-reg! c a-reg v-w n+1)
(gate-x! c s 0)
(gate-x! c f 0)
;; ── Forward sweep with width truncation (R2) ──
(let loop ((i 0))
(when (< i iters)
(kaliski-iteration-trunc! c i n+1 p
u v-w r s
f 0
m-hist i
cin-reg cin-idx
tmp-reg flag-reg flag-idx
red-tmp-reg
a-mask)
(loop (+ i 1))))
;; ── Classical correction K (same as textbook path) ──
(let* ((r-on-1 (classical-kaliski-r-final 1 p iters n+1))
(k-correct (classical-mod-inv p r-on-1))
(k-inverse r-on-1))
(in-place-mul-const! c r n+1 p k-correct
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg
mc-tmp mc-pow)
(cx-copy-reg! c r out-reg n+1)
(in-place-mul-const! c r n+1 p k-inverse
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg
mc-tmp mc-pow))
;; ── Classical-replay backward sweep (R1) ──
;;
;; After the forward sweep + correction-mul + copy + uncorrection,
;; the registers hold their post-forward-Kaliski classical values
;; for the bound a. classical-kaliski-trace computes these for the
;; specific bound a; we emit X gates against every set bit to zero
;; each register at zero Toffoli cost.
;;
;; Trace order: (final-u final-v-w final-r final-s final-f m-i-vec).
(let* ((trace (classical-kaliski-trace a-bound p iters n+1))
(final-u (car trace))
(final-v-w (car (cdr trace)))
(final-r (car (cdr (cdr trace))))
(final-s (car (cdr (cdr (cdr trace)))))
(final-f (car (cdr (cdr (cdr (cdr trace))))))
(m-i-vec (car (cdr (cdr (cdr (cdr (cdr trace))))))))
;; Reset u, v_w, r, s by classical-replay X pattern.
(classical-reset! c u n+1 final-u)
(classical-reset! c v-w n+1 final-v-w)
(classical-reset! c r n+1 final-r)
(classical-reset! c s n+1 final-s)
;; Reset f (1 bit).
(when (= final-f 1) (gate-x! c f 0))
;; Reset m_hist bit-by-bit.
(let loop ((i 0))
(when (< i iters)
(when (= (vector-ref m-i-vec i) 1)
(gate-x! c m-hist i))
(loop (+ i 1)))))
;; ── Free Kaliski state ──
(free! c a-mask)
(free! c m-hist)
(free! c f)
(free! c s)
(free! c r)
(free! c v-w)
(free! c u)))))
;;; ── mod-inv-by! — public dispatcher ──────────────────────────────
;;;
;;; Routes to mod-inv-by-refined! when *mod-inv-by-refined* is non-#f
;;; (Phase B step 9 refined path), else to mod-inv-by-textbook!. The
;;; refined path requires a-reg to be bound via bind-input!; callers
;;; that have not bound a-reg should leave *mod-inv-by-refined* = #f.
(define (mod-inv-by! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg)
"out := a^{-1} mod p via Kaliski / Bernstein-Yang. Dispatches based on
*mod-inv-by-refined* (refined path = R1+R2, textbook path = Bennett
compute-copy-uncompute). Same calling convention as mod-inv!."
(cond
(*mod-inv-by-split-eea*
;; sweep-055: routes to mod-inv-by-split-eea! (split-EEA primitive,
;; lumbda/mod-inv-by-split-eea.lsp). Caller must have loaded that
;; file before flipping the flag; otherwise this branch raises an
;; unbound-symbol error on mod-inv-by-split-eea!.
(mod-inv-by-split-eea! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg))
(*mod-inv-by-dialog-gcd-host*
(mod-inv-by-dialog-gcd-host! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg))
(*mod-inv-by-dialog-gcd*
(mod-inv-by-dialog-gcd! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg))
(*mod-inv-by-refined*
(mod-inv-by-refined! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg))
(else
(mod-inv-by-textbook! c a-reg out-reg n+1 p
cin-reg cin-idx tmp-reg flag-reg flag-idx
red-tmp-reg))))