;;; 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))))
