Skip to main content

cortiq_engine/
linear_core.rs

1//! Linear-attention cores, selected by `arch.linear_core.kind`
2//! (descriptor-driven operators — Patent 15 claim 8).
3//!
4//! Two tracks (owner decision 2026-07-04):
5//!
6//! * `gated_delta_net` — the faithful vendor operator (Qwen3.5 /
7//!   Qwen3-Next). Default for models that ship GDN weights: conversion
8//!   carries the tensors 1:1 and needs no training. Port of the
9//!   validated `gated_delta_net` (vmfcore/rust/src/forward.rs) against
10//!   the numpy/torch oracle (vmfcore/gdn_layer.py).
11//!
12//! * `vmf_phase` — the canonical core: token carries a phase θ; kernel
13//!   φ(θ) = [cos θ; sin θ] gives a linear factorization; the condensate
14//!   is a recurrent state S[head][p2, dv] with decay exp(−exp(A_log)).
15//!   Noise-robust and simpler than vendor recurrences. Exotic operators
16//!   are folded onto it at CONVERT time (`--linear-core vmf_phase`) and
17//!   quality is restored by the offline heal — the research track and
18//!   the production mechanism for Patent-15 skills (mask→heal→compress).
19//!
20//! Both cores implement the same contract: `*_forward` (one position,
21//! advances the state) and `*_pair` (fused two positions; lane 1
22//! commits, lane 2 is tentative in `scratch` for speculative verify).
23//! State lives in the layer's `linear_state: Vec<f32>` and is resized
24//! lazily by the core itself.
25
26use crate::pool::Pool;
27use crate::qtensor::QTensor;
28
29/// Weights of one vmf_phase layer (`model.layers.{i}.vmf_attn.*`).
30pub struct VmfPhaseWeights {
31    /// [nh·nphase, hidden] — query phase projection
32    pub thq: QTensor,
33    /// [nh·nphase, hidden] — key phase projection
34    pub thk: QTensor,
35    /// [nh·dv, hidden]
36    pub v_proj: QTensor,
37    /// [hidden, nh·dv]
38    pub out_proj: QTensor,
39    /// Per-component decay exp(−exp(A_log)), len nh·2·nphase (precomputed).
40    pub decay: Vec<f64>,
41    /// Selective-write input gate κ (hybrid_k core, stage 71): weight
42    /// [nh, hidden] + bias [nh]; κ_h = σ(W_k·x + b)_h multiplies the
43    /// state WRITE (S = decay·S + κ·φk⊗v). None = classic phase core,
44    /// bit-identical to the pre-κ kernel. Measured at mechanism level:
45    /// knee ×2–6 earlier, restores correlated-noise robustness, LM
46    /// crossover vs softmax at SEQ 512 (experiments/lc_final_merged.json).
47    pub k_gate: Option<(QTensor, Vec<f32>)>,
48}
49
50#[derive(Clone, Copy)]
51pub struct VmfPhaseCfg {
52    pub num_heads: usize,
53    pub nphase: usize,
54    pub value_head_dim: usize,
55    pub hidden_size: usize,
56    /// θ-mass (η′ correction): a restoring potential pulling the phase
57    /// toward 0 — θ_eff = θ/(1+mass) — which WIDENS the phase kernel.
58    /// Measured (experiments/vmf_native_core*.py) to restore noise
59    /// robustness when the phase projection is FIXED (exactly CMF's
60    /// fold-before-heal regime: thq/thk are init, not trained) — recall
61    /// 3%→91% at moderate noise; redundant once the projection is
62    /// healed. 0.0 = massless Goldstone (bit-identical to prior kernel).
63    /// Set via CMF_PHASE_MASS. Validated at mechanism level, not yet LM.
64    pub phase_mass: f32,
65}
66
67impl VmfPhaseCfg {
68    pub fn state_len(&self) -> usize {
69        self.num_heads * 2 * self.nphase * self.value_head_dim
70    }
71}
72
73/// One recurrent step for one head-set given projected phases/values.
74/// `state` is S[nh][p2, dv] stored f32 (per-element math in f64 — the
75/// storage halves, each step's arithmetic keeps the old precision).
76fn phase_step(
77    thq: &[f32],
78    thk: &[f32],
79    v: &[f32],
80    decay: &[f64],
81    kap: Option<&[f32]>,
82    cfg: &VmfPhaseCfg,
83    state: &mut [f32],
84    out: &mut [f32],
85) {
86    let (nh, nph, dv) = (cfg.num_heads, cfg.nphase, cfg.value_head_dim);
87    // θ-mass (η′): θ_eff = θ/(1+mass). mass=0 → factor 1 → no-op.
88    let mscale = 1.0f64 / (1.0 + cfg.phase_mass as f64);
89    let p2 = 2 * nph;
90    for h in 0..nh {
91        let s = &mut state[h * p2 * dv..(h + 1) * p2 * dv];
92        let thk_h = &thk[h * nph..(h + 1) * nph];
93        let thq_h = &thq[h * nph..(h + 1) * nph];
94        let vt = &v[h * dv..(h + 1) * dv];
95        let ot = &mut out[h * dv..(h + 1) * dv];
96        let dec = &decay[h * p2..(h + 1) * p2];
97        // Selective write (hybrid_k): κ scales what enters the condensate.
98        let kh = kap.map_or(1.0f64, |k| k[h] as f64);
99        for f in 0..p2 {
100            // φ(θ) = [cos·nph, sin·nph], θ scaled by the mass factor.
101            let (fk, fq) = if f < nph {
102                (
103                    (thk_h[f] as f64 * mscale).cos(),
104                    (thq_h[f] as f64 * mscale).cos(),
105                )
106            } else {
107                (
108                    (thk_h[f - nph] as f64 * mscale).sin(),
109                    (thq_h[f - nph] as f64 * mscale).sin(),
110                )
111            };
112            let fkw = fk * kh;
113            let row = &mut s[f * dv..(f + 1) * dv];
114            let dcf = dec[f];
115            for d in 0..dv {
116                // S = decay·S + κ·φk⊗v (f64 math, f32 cell)
117                let cell = dcf * row[d] as f64 + fkw * vt[d] as f64;
118                row[d] = cell as f32;
119                ot[d] += (fq * cell) as f32; // o = Σ φq·S
120            }
121        }
122    }
123}
124
125/// κ_h = σ(W_k·x + b)_h — the per-head write gate (None when the layer
126/// has no k_gate tensors: classic phase core).
127fn kappa_of(x: &[f32], w: &VmfPhaseWeights, nh: usize, pool: Option<&Pool>) -> Option<Vec<f32>> {
128    let (kw, kb) = w.k_gate.as_ref()?;
129    let mut k = vec![0.0f32; nh];
130    kw.matvec(x, &mut k, pool);
131    for (v, b) in k.iter_mut().zip(kb) {
132        *v = 1.0 / (1.0 + (-(*v + b)).exp());
133    }
134    Some(k)
135}
136
137/// Forward one position through a vmf_phase layer, advancing `state`.
138pub fn vmf_phase_forward(
139    x: &[f32],
140    w: &VmfPhaseWeights,
141    cfg: &VmfPhaseCfg,
142    state: &mut Vec<f32>,
143    pool: Option<&Pool>,
144) -> Vec<f32> {
145    if state.len() != cfg.state_len() {
146        *state = vec![0f32; cfg.state_len()];
147    }
148    let (nh, nph, dv) = (cfg.num_heads, cfg.nphase, cfg.value_head_dim);
149
150    let mut thq = vec![0.0f32; nh * nph];
151    w.thq.matvec(x, &mut thq, pool);
152    let mut thk = vec![0.0f32; nh * nph];
153    w.thk.matvec(x, &mut thk, pool);
154    let mut v = vec![0.0f32; nh * dv];
155    w.v_proj.matvec(x, &mut v, pool);
156
157    let kap = kappa_of(x, w, nh, pool);
158    let mut o = vec![0.0f32; nh * dv];
159    phase_step(&thq, &thk, &v, &w.decay, kap.as_deref(), cfg, state, &mut o);
160
161    let mut out = vec![0.0f32; cfg.hidden_size];
162    w.out_proj.matvec(&o, &mut out, pool);
163    out
164}
165
166/// Fused two-position forward (speculative verify). Lane 1 commits into
167/// `state` (its token is always committed); lane 2's tentative state
168/// goes into `scratch` — the caller swaps it in on draft acceptance and
169/// simply drops it on rejection.
170#[allow(clippy::too_many_arguments)]
171pub fn vmf_phase_pair(
172    x1: &[f32],
173    x2: &[f32],
174    w: &VmfPhaseWeights,
175    cfg: &VmfPhaseCfg,
176    state: &mut Vec<f32>,
177    scratch: &mut Vec<f32>,
178    pool: Option<&Pool>,
179) -> (Vec<f32>, Vec<f32>) {
180    if state.len() != cfg.state_len() {
181        *state = vec![0f32; cfg.state_len()];
182    }
183    let (nh, nph, dv) = (cfg.num_heads, cfg.nphase, cfg.value_head_dim);
184
185    let mut thq1 = vec![0.0f32; nh * nph];
186    let mut thq2 = vec![0.0f32; nh * nph];
187    w.thq.matvec2(x1, x2, &mut thq1, &mut thq2, pool);
188    let mut thk1 = vec![0.0f32; nh * nph];
189    let mut thk2 = vec![0.0f32; nh * nph];
190    w.thk.matvec2(x1, x2, &mut thk1, &mut thk2, pool);
191    let mut v1 = vec![0.0f32; nh * dv];
192    let mut v2 = vec![0.0f32; nh * dv];
193    w.v_proj.matvec2(x1, x2, &mut v1, &mut v2, pool);
194
195    // Lane 1 commits into the real state.
196    let kap1 = kappa_of(x1, w, nh, pool);
197    let mut o1 = vec![0.0f32; nh * dv];
198    phase_step(
199        &thq1,
200        &thk1,
201        &v1,
202        &w.decay,
203        kap1.as_deref(),
204        cfg,
205        state,
206        &mut o1,
207    );
208
209    // Lane 2 runs on a copy — tentative until the draft is verified.
210    let kap2 = kappa_of(x2, w, nh, pool);
211    scratch.clear();
212    scratch.extend_from_slice(state);
213    let mut o2 = vec![0.0f32; nh * dv];
214    phase_step(
215        &thq2,
216        &thk2,
217        &v2,
218        &w.decay,
219        kap2.as_deref(),
220        cfg,
221        scratch,
222        &mut o2,
223    );
224
225    let mut out1 = vec![0.0f32; cfg.hidden_size];
226    let mut out2 = vec![0.0f32; cfg.hidden_size];
227    w.out_proj.matvec2(&o1, &o2, &mut out1, &mut out2, pool);
228    (out1, out2)
229}
230
231// ───────────────────────── GatedDeltaNet (faithful vendor operator) ─────────────────────────
232
233/// Weights of one GatedDeltaNet layer (`model.layers.{i}.linear_attn.*`,
234/// names 1:1 with the source model — no fold, no training).
235pub struct GdnWeights {
236    /// [2·nk·dk + nv·dv, hidden] — fused q/k/v projection
237    pub in_proj_qkv: QTensor,
238    /// [nv·dv, hidden] — output-gate projection z
239    pub in_proj_z: QTensor,
240    /// [nv, hidden] — decay modulation a
241    pub in_proj_a: QTensor,
242    /// [nv, hidden] — write-strength b (β = σ(b))
243    pub in_proj_b: QTensor,
244    /// [c_dim · kk] — depthwise causal conv taps, flattened [c][tap]
245    pub conv1d: Vec<f32>,
246    /// [nv]
247    pub a_log: Vec<f32>,
248    /// [nv]
249    pub dt_bias: Vec<f32>,
250    /// [dv] — gated RMSNorm weight (plain x̂·w, validated by the oracle)
251    pub norm: Vec<f32>,
252    /// [hidden, nv·dv]
253    pub out_proj: QTensor,
254}
255
256#[derive(Clone, Copy)]
257pub struct GdnCfg {
258    pub num_v_heads: usize,
259    pub num_k_heads: usize,
260    pub key_head_dim: usize,
261    pub value_head_dim: usize,
262    pub conv_kernel: usize,
263    pub hidden_size: usize,
264    pub rms_eps: f64,
265}
266
267impl GdnCfg {
268    pub fn conv_dim(&self) -> usize {
269        2 * self.num_k_heads * self.key_head_dim + self.num_v_heads * self.value_head_dim
270    }
271
272    /// Packed state: [conv ring (kk−1)·c_dim | S nv·dk·dv], one Vec<f64>
273    /// so the speculative scratch-swap moves ring and condensate together.
274    pub fn state_len(&self) -> usize {
275        (self.conv_kernel - 1) * self.conv_dim()
276            + self.num_v_heads * self.key_head_dim * self.value_head_dim
277    }
278}
279
280fn softplus(x: f64) -> f64 {
281    if x > 20.0 { x } else { x.exp().ln_1p() }
282}
283
284fn sigmoid(x: f64) -> f64 {
285    1.0 / (1.0 + (-x).exp())
286}
287
288fn silu(x: f64) -> f64 {
289    x / (1.0 + (-x).exp())
290}
291
292/// `*mut f32` that may cross worker threads; safety comes from the
293/// disjoint (head, element) ranges each worker writes.
294#[derive(Clone, Copy)]
295struct SendMutF32(*mut f32);
296unsafe impl Send for SendMutF32 {}
297unsafe impl Sync for SendMutF32 {}
298
299/// One recurrent step given the raw (pre-conv) projections of this
300/// position. Advances the packed state (conv ring + S) and writes the
301/// gated per-head output into `of` [nv·dv].
302///
303/// The condensate math runs in f32 (the vendor operator's own dtype —
304/// `mamba_ssm_dtype: float32` in the source configs; the old f64 was
305/// over-precision at 4× the traffic and no SIMD). The two S passes are
306/// element-wise in `dj` with no cross-lane reduction, so LLVM
307/// auto-vectorizes them (fmla on NEON, FMA on AVX2). Heads are
308/// independent given the conv output and run across the pool — on a
309/// Qwen3.5-27B this loop is 48 heads × 128×128 × 48 layers per token,
310/// the single biggest serial block in the hybrid's decode.
311#[allow(clippy::too_many_arguments)]
312fn gdn_step(
313    qkv: &[f32],
314    z: &[f32],
315    a: &[f32],
316    b: &[f32],
317    w: &GdnWeights,
318    cfg: &GdnCfg,
319    state: &mut [f32],
320    of: &mut [f32],
321    pool: Option<&Pool>,
322) {
323    let (nv, nk, dk, dv, kk) = (
324        cfg.num_v_heads,
325        cfg.num_k_heads,
326        cfg.key_head_dim,
327        cfg.value_head_dim,
328        cfg.conv_kernel,
329    );
330    let c_dim = cfg.conv_dim();
331    let (kd, rep) = (nk * dk, nv / nk);
332    let (ring, s_all) = state.split_at_mut((kk - 1) * c_dim);
333
334    // Depthwise causal conv over [ring…, current] + SiLU. Taps are
335    // ordered oldest→newest; tap kk−1 multiplies the current position.
336    // (Tiny: c_dim × kk — f64 accumulation kept.)
337    let mut cq = vec![0f32; c_dim];
338    for c in 0..c_dim {
339        let taps = &w.conv1d[c * kk..(c + 1) * kk];
340        let mut acc = qkv[c] as f64 * taps[kk - 1] as f64;
341        for j in 0..kk - 1 {
342            acc += ring[j * c_dim + c] as f64 * taps[j] as f64;
343        }
344        cq[c] = silu(acc) as f32;
345    }
346    // Ring shift: drop the oldest position, append the raw current one.
347    if kk > 1 {
348        ring.copy_within(c_dim.., 0);
349        let tail = (kk - 2) * c_dim;
350        ring[tail..tail + c_dim].copy_from_slice(&qkv[..c_dim]);
351    }
352
353    let cq = &cq;
354    let s_ptr = SendMutF32(s_all.as_mut_ptr());
355    let of_ptr = SendMutF32(of.as_mut_ptr());
356    let head_range = |h0: usize, h1: usize| {
357        // Rebind the Sync wrappers whole — edition-2021 disjoint capture
358        // would otherwise grab the raw `.0` fields and lose Send/Sync.
359        let (s_ptr, of_ptr) = (s_ptr, of_ptr);
360        // Per-worker scratch, recycled across calls (thread-local freelists).
361        let mut kv = crate::attention::take_buf(dv);
362        let mut delta = crate::attention::take_buf(dv);
363        let mut o = crate::attention::take_buf(dv);
364        let mut kf = crate::attention::take_buf(dk);
365        let mut qf = crate::attention::take_buf(dk);
366        for h in h0..h1 {
367            let ko = h / rep; // source q/k head (GQA)
368            let (qs, ks) = (ko * dk, kd + ko * dk);
369            // l2-normalize q and k; q additionally scaled by 1/√dk.
370            let (mut nq, mut nkn) = (0f64, 0f64);
371            for d in 0..dk {
372                nq += (cq[qs + d] as f64) * (cq[qs + d] as f64);
373                nkn += (cq[ks + d] as f64) * (cq[ks + d] as f64);
374            }
375            let invq = (1.0 / ((nq + 1e-6).sqrt() * (dk as f64).sqrt())) as f32;
376            let invk = (1.0 / (nkn + 1e-6).sqrt()) as f32;
377            for d in 0..dk {
378                qf[d] = cq[qs + d] * invq;
379                kf[d] = cq[ks + d] * invk;
380            }
381
382            let g = (-(w.a_log[h] as f64).exp() * softplus(a[h] as f64 + w.dt_bias[h] as f64)).exp()
383                as f32;
384            let beta = sigmoid(b[h] as f64) as f32;
385
386            // SAFETY: disjoint per-head S and output slices per worker.
387            let s = unsafe { std::slice::from_raw_parts_mut(s_ptr.0.add(h * dk * dv), dk * dv) };
388            let oh = unsafe { std::slice::from_raw_parts_mut(of_ptr.0.add(h * dv), dv) };
389            let vt = &cq[2 * kd + h * dv..2 * kd + (h + 1) * dv];
390
391            // S ← g·S;  kv = kᵀS;  S += k ⊗ β(v − kv);  o = qᵀS —
392            // algebraically regrouped so S is READ twice and WRITTEN
393            // once: kv over S_old (then ×g), one fused update+query pass.
394            kv[..dv].fill(0.0);
395            for di in 0..dk {
396                let kfd = kf[di];
397                let row = &s[di * dv..(di + 1) * dv];
398                for dj in 0..dv {
399                    kv[dj] += row[dj] * kfd; // elementwise in dj → SIMD
400                }
401            }
402            for dj in 0..dv {
403                delta[dj] = (vt[dj] - g * kv[dj]) * beta;
404            }
405            o[..dv].fill(0.0);
406            for di in 0..dk {
407                let kfd = kf[di];
408                let qfd = qf[di];
409                let row = &mut s[di * dv..(di + 1) * dv];
410                for dj in 0..dv {
411                    let cell = g * row[dj] + kfd * delta[dj];
412                    row[dj] = cell;
413                    o[dj] += qfd * cell; // elementwise in dj → SIMD
414                }
415            }
416            // Gated RMSNorm per head: x̂·w·silu(z) (oracle-validated form).
417            let ss: f64 = o[..dv].iter().map(|&v| (v as f64) * (v as f64)).sum();
418            let inv = 1.0 / (ss / dv as f64 + cfg.rms_eps).sqrt();
419            for dj in 0..dv {
420                oh[dj] =
421                    ((o[dj] as f64 * inv) * w.norm[dj] as f64 * silu(z[h * dv + dj] as f64)) as f32;
422            }
423        }
424        crate::attention::recycle_buf(&mut kv);
425        crate::attention::recycle_buf(&mut delta);
426        crate::attention::recycle_buf(&mut o);
427        crate::attention::recycle_buf(&mut kf);
428        crate::attention::recycle_buf(&mut qf);
429    };
430    match pool {
431        Some(pool) if nv >= 4 => pool.run(&|widx, n| {
432            let chunk = nv.div_ceil(n);
433            let h0 = (widx * chunk).min(nv);
434            let h1 = (h0 + chunk).min(nv);
435            if h0 < h1 {
436                head_range(h0, h1);
437            }
438        }),
439        _ => head_range(0, nv),
440    }
441}
442
443/// Forward one position through a GatedDeltaNet layer, advancing `state`.
444pub fn gdn_forward(
445    x: &[f32],
446    w: &GdnWeights,
447    cfg: &GdnCfg,
448    state: &mut Vec<f32>,
449    pool: Option<&Pool>,
450) -> Vec<f32> {
451    if state.len() != cfg.state_len() {
452        *state = vec![0f32; cfg.state_len()];
453    }
454    let (c_dim, vd) = (cfg.conv_dim(), cfg.num_v_heads * cfg.value_head_dim);
455
456    let mut qkv = vec![0.0f32; c_dim];
457    let mut z = vec![0.0f32; vd];
458    let mut a = vec![0.0f32; cfg.num_v_heads];
459    let mut b = vec![0.0f32; cfg.num_v_heads];
460    // D5: two heavy projections (the GDN mixer is ~half a hybrid layer's
461    // bytes) — one GPU submission; a/b are tiny and stay on CPU. The
462    // Batch probe arbitrates GPU vs the fused-CPU dispatch per machine.
463    let cpu_projs = |qkv: &mut Vec<f32>, z: &mut Vec<f32>, a: &mut Vec<f32>, b: &mut Vec<f32>| {
464        QTensor::matvec_many(
465            [&w.in_proj_qkv, &w.in_proj_z, &w.in_proj_a, &w.in_proj_b],
466            x,
467            [
468                qkv.as_mut_slice(),
469                z.as_mut_slice(),
470                a.as_mut_slice(),
471                b.as_mut_slice(),
472            ],
473            pool,
474        );
475    };
476    let mut done = false;
477    if crate::gpu::enabled_here() && gdn_projs_eligible(w) {
478        match crate::gpu::probe_arm(crate::gpu::OpClass::Batch) {
479            crate::gpu::ProbeArm::Gpu => {
480                let t0 = std::time::Instant::now();
481                if gdn_projs_gpu(w, x, &mut qkv, &mut z) {
482                    crate::gpu::probe_record(crate::gpu::OpClass::Batch, true, t0.elapsed());
483                    w.in_proj_a.matvec(x, &mut a, pool);
484                    w.in_proj_b.matvec(x, &mut b, pool);
485                    done = true;
486                }
487            }
488            crate::gpu::ProbeArm::CpuTimed => {
489                let t0 = std::time::Instant::now();
490                crate::gpu::cpu_scope(|| cpu_projs(&mut qkv, &mut z, &mut a, &mut b));
491                crate::gpu::probe_record(crate::gpu::OpClass::Batch, false, t0.elapsed());
492                done = true;
493            }
494            crate::gpu::ProbeArm::Cpu => {
495                crate::gpu::cpu_scope(|| cpu_projs(&mut qkv, &mut z, &mut a, &mut b));
496                done = true;
497            }
498        }
499    }
500    if !done {
501        cpu_projs(&mut qkv, &mut z, &mut a, &mut b);
502    }
503
504    let mut of = vec![0.0f32; vd];
505    gdn_step(&qkv, &z, &a, &b, w, cfg, state, &mut of, pool);
506
507    let mut out = vec![0.0f32; cfg.hidden_size];
508    w.out_proj.matvec(&of, &mut out, pool);
509    out
510}
511
512/// Batched GDN forward (prefill-GEMM): the qkv/z/a/b and out_proj
513/// projections are matmat over the batch (a weight row once per chunk),
514/// the gdn_step recurrence runs sequentially over positions (state is the
515/// same as the sequential path; the math is elementwise identical).
516pub fn gdn_forward_batch(
517    xs: &[f32],
518    b: usize,
519    w: &GdnWeights,
520    cfg: &GdnCfg,
521    state: &mut Vec<f32>,
522    pool: Option<&Pool>,
523) -> Vec<f32> {
524    if state.len() != cfg.state_len() {
525        *state = vec![0f32; cfg.state_len()];
526    }
527    let (c_dim, vd) = (cfg.conv_dim(), cfg.num_v_heads * cfg.value_head_dim);
528    let nv = cfg.num_v_heads;
529
530    let mut qkv = vec![0.0f32; b * c_dim];
531    w.in_proj_qkv.matmat(xs, b, &mut qkv, pool);
532    let mut z = vec![0.0f32; b * vd];
533    w.in_proj_z.matmat(xs, b, &mut z, pool);
534    let mut a = vec![0.0f32; b * nv];
535    w.in_proj_a.matmat(xs, b, &mut a, pool);
536    let mut bb = vec![0.0f32; b * nv];
537    w.in_proj_b.matmat(xs, b, &mut bb, pool);
538
539    let mut of = vec![0.0f32; b * vd];
540    for bi in 0..b {
541        gdn_step(
542            &qkv[bi * c_dim..(bi + 1) * c_dim],
543            &z[bi * vd..(bi + 1) * vd],
544            &a[bi * nv..(bi + 1) * nv],
545            &bb[bi * nv..(bi + 1) * nv],
546            w,
547            cfg,
548            state,
549            &mut of[bi * vd..(bi + 1) * vd],
550            pool,
551        );
552    }
553    let mut out = vec![0.0f32; b * cfg.hidden_size];
554    w.out_proj.matmat(&of, b, &mut out, pool);
555    out
556}
557
558/// GDN qkv+z GPU eligibility: q1 mixers offload by default (the CPU q1
559/// kernel is compute-bound); q8 stays opt-in via CMF_GPU_GDN=1 (measured
560/// neutral). The probe in `gdn_forward` still arbitrates either way.
561fn gdn_projs_eligible(w: &GdnWeights) -> bool {
562    w.in_proj_qkv.is_q1()
563        || std::env::var("CMF_GPU_GDN")
564            .map(|v| v == "1")
565            .unwrap_or(false)
566}
567
568/// GDN qkv+z on GPU in a single submission (independent matvecs of one input).
569fn gdn_projs_gpu(w: &GdnWeights, x: &[f32], qkv: &mut [f32], z: &mut [f32]) -> bool {
570    use crate::gpu::matvec_batch;
571    use crate::qtensor::QTensor;
572    if !crate::gpu::enabled_here() {
573        return false;
574    }
575    fn part<'a>(
576        t: &'a QTensor,
577        x: &[f32],
578    ) -> Option<(
579        std::sync::Arc<cortiq_core::CmfModel>,
580        crate::gpu::BatchJob<'a>,
581    )> {
582        use crate::gpu::BatchJob;
583        use crate::qtensor::prescale;
584        use cortiq_core::TensorDtype;
585        match t {
586            QTensor::Mapped {
587                model,
588                idx,
589                dtype: dt @ (TensorDtype::Q8Row | TensorDtype::Q8_2f),
590                rows,
591                cols,
592                row_scale,
593                col_field,
594                ..
595            } => Some((
596                model.clone(),
597                BatchJob {
598                    idx: *idx,
599                    rows: *rows,
600                    cols: *cols,
601                    row_scale,
602                    xs: prescale(x, col_field, *dt).into_owned(),
603                    q1: false,
604                },
605            )),
606            QTensor::Mapped {
607                model,
608                idx,
609                dtype: TensorDtype::Q1,
610                rows,
611                cols,
612                ..
613            } => Some((
614                model.clone(),
615                BatchJob {
616                    idx: *idx,
617                    rows: *rows,
618                    cols: *cols,
619                    row_scale: &[],
620                    xs: x.to_vec(),
621                    q1: true,
622                },
623            )),
624            _ => None,
625        }
626    }
627    let Some((model, jq)) = part(&w.in_proj_qkv, x) else {
628        return false;
629    };
630    let Some((_, jz)) = part(&w.in_proj_z, x) else {
631        return false;
632    };
633    matvec_batch(&model, &[jq, jz], &mut [qkv, z])
634}
635
636/// Fused two-position forward (speculative verify): lane 1 commits into
637/// `state`, lane 2 is tentative in `scratch` (ring + S move together).
638#[allow(clippy::too_many_arguments)]
639pub fn gdn_pair(
640    x1: &[f32],
641    x2: &[f32],
642    w: &GdnWeights,
643    cfg: &GdnCfg,
644    state: &mut Vec<f32>,
645    scratch: &mut Vec<f32>,
646    pool: Option<&Pool>,
647) -> (Vec<f32>, Vec<f32>) {
648    if state.len() != cfg.state_len() {
649        *state = vec![0f32; cfg.state_len()];
650    }
651    let (c_dim, vd, nv) = (
652        cfg.conv_dim(),
653        cfg.num_v_heads * cfg.value_head_dim,
654        cfg.num_v_heads,
655    );
656
657    let mut qkv1 = vec![0.0f32; c_dim];
658    let mut qkv2 = vec![0.0f32; c_dim];
659    w.in_proj_qkv.matvec2(x1, x2, &mut qkv1, &mut qkv2, pool);
660    let mut z1 = vec![0.0f32; vd];
661    let mut z2 = vec![0.0f32; vd];
662    w.in_proj_z.matvec2(x1, x2, &mut z1, &mut z2, pool);
663    let mut a1 = vec![0.0f32; nv];
664    let mut a2 = vec![0.0f32; nv];
665    w.in_proj_a.matvec2(x1, x2, &mut a1, &mut a2, pool);
666    let mut b1 = vec![0.0f32; nv];
667    let mut b2 = vec![0.0f32; nv];
668    w.in_proj_b.matvec2(x1, x2, &mut b1, &mut b2, pool);
669
670    let mut of1 = vec![0.0f32; vd];
671    gdn_step(&qkv1, &z1, &a1, &b1, w, cfg, state, &mut of1, pool);
672
673    scratch.clear();
674    scratch.extend_from_slice(state);
675    let mut of2 = vec![0.0f32; vd];
676    gdn_step(&qkv2, &z2, &a2, &b2, w, cfg, scratch, &mut of2, pool);
677
678    let mut out1 = vec![0.0f32; cfg.hidden_size];
679    let mut out2 = vec![0.0f32; cfg.hidden_size];
680    w.out_proj.matvec2(&of1, &of2, &mut out1, &mut out2, pool);
681    (out1, out2)
682}
683
684// ───────────────────────── ShortConv (LFM2 gated short convolution) ─────────────────────────
685
686/// Weights of one LFM2 short-convolution mixer
687/// (`model.layers.{i}.short_conv.*`, renamed from the vendor `conv.*` at
688/// convert time). No recurrent condensate — the only state is the causal
689/// conv ring (the last `kernel−1` gated inputs per channel).
690pub struct ShortConvWeights {
691    /// [3·hidden, hidden] — fused (B, C, x) projection.
692    pub in_proj: QTensor,
693    /// [hidden · kernel] depthwise conv taps, flattened `[channel][tap]`
694    /// (the source `[hidden, 1, kernel]` with the singleton group axis
695    /// dropped). Tap `kernel−1` multiplies the current position.
696    pub conv: Vec<f32>,
697    /// [hidden, hidden] — output projection.
698    pub out_proj: QTensor,
699}
700
701#[derive(Clone, Copy)]
702pub struct ShortConvCfg {
703    pub hidden_size: usize,
704    /// Conv kernel width `L` (`conv_L_cache`; LFM2 uses 3).
705    pub kernel: usize,
706}
707
708impl ShortConvCfg {
709    /// Conv ring: the last `kernel−1` gated inputs per channel.
710    pub fn state_len(&self) -> usize {
711        (self.kernel - 1) * self.hidden_size
712    }
713}
714
715/// One position through the gated conv, given the fused projection
716/// `bcx = in_proj·x` [3·hidden] = [B | C | x]. Advances the conv ring and
717/// writes the gated conv output `y = C ⊙ conv(B ⊙ x)` [hidden] into `y`.
718///
719/// The conv is PyTorch's causal depthwise `Conv1d(padding=kernel−1)`
720/// truncated to the current length: for tap `k`, weight `w[c][k]` pairs
721/// with the input `kernel−1−k` steps in the past, so `w[c][kernel−1]` is
722/// the current position. The ring holds `in[t−1] … in[t−(kernel−1)]` at
723/// slots `0 … kernel−2`.
724fn short_conv_step(
725    bcx: &[f32],
726    conv: &[f32],
727    cfg: &ShortConvCfg,
728    ring_state: &mut [f32],
729    y: &mut [f32],
730) {
731    let (h, k) = (cfg.hidden_size, cfg.kernel);
732    let ring = k - 1;
733    let (bg, cg, xg) = (&bcx[0..h], &bcx[h..2 * h], &bcx[2 * h..3 * h]);
734    for c in 0..h {
735        let bx = bg[c] * xg[c];
736        let wc = &conv[c * k..(c + 1) * k];
737        // Current tap, then the past taps read from the channel's ring.
738        let mut acc = wc[k - 1] * bx;
739        let rc = &mut ring_state[c * ring..c * ring + ring];
740        for s in 0..ring {
741            acc += wc[k - 2 - s] * rc[s];
742        }
743        y[c] = cg[c] * acc;
744        // Shift newest-in-front: slot 0 becomes the just-seen input.
745        for s in (1..ring).rev() {
746            rc[s] = rc[s - 1];
747        }
748        if ring > 0 {
749            rc[0] = bx;
750        }
751    }
752}
753
754/// Forward one position through a short-conv layer, advancing `state`.
755pub fn short_conv_forward(
756    x: &[f32],
757    w: &ShortConvWeights,
758    cfg: &ShortConvCfg,
759    state: &mut Vec<f32>,
760    pool: Option<&Pool>,
761) -> Vec<f32> {
762    if state.len() != cfg.state_len() {
763        *state = vec![0f32; cfg.state_len()];
764    }
765    let h = cfg.hidden_size;
766    let mut bcx = vec![0.0f32; 3 * h];
767    w.in_proj.matvec(x, &mut bcx, pool);
768    let mut y = vec![0.0f32; h];
769    short_conv_step(&bcx, &w.conv, cfg, state, &mut y);
770    let mut out = vec![0.0f32; h];
771    w.out_proj.matvec(&y, &mut out, pool);
772    out
773}
774
775/// Batched short-conv forward (prefill-GEMM): in_proj/out_proj are matmat
776/// over the chunk (a weight row streamed once), the conv walks the
777/// positions in order — the chunk is contiguous, so the ring state is
778/// exactly the sequential path's and the math is elementwise identical.
779pub fn short_conv_forward_batch(
780    xs: &[f32],
781    b: usize,
782    w: &ShortConvWeights,
783    cfg: &ShortConvCfg,
784    state: &mut Vec<f32>,
785    pool: Option<&Pool>,
786) -> Vec<f32> {
787    if state.len() != cfg.state_len() {
788        *state = vec![0f32; cfg.state_len()];
789    }
790    let h = cfg.hidden_size;
791    let mut bcx = vec![0.0f32; b * 3 * h];
792    w.in_proj.matmat(xs, b, &mut bcx, pool);
793    let mut y = vec![0.0f32; b * h];
794    for bi in 0..b {
795        short_conv_step(
796            &bcx[bi * 3 * h..(bi + 1) * 3 * h],
797            &w.conv,
798            cfg,
799            state,
800            &mut y[bi * h..(bi + 1) * h],
801        );
802    }
803    let mut out = vec![0.0f32; b * h];
804    w.out_proj.matmat(&y, b, &mut out, pool);
805    out
806}
807
808/// Fused two-position forward (speculative verify). Lane 1 commits into
809/// `state`; lane 2's tentative ring goes into `scratch` — swapped in on
810/// draft acceptance, dropped on rejection. LFM2 ships no MTP head, so this
811/// is exercised only by the pair-fusion micro-benchmark; kept correct.
812#[allow(clippy::too_many_arguments)]
813pub fn short_conv_pair(
814    x1: &[f32],
815    x2: &[f32],
816    w: &ShortConvWeights,
817    cfg: &ShortConvCfg,
818    state: &mut Vec<f32>,
819    scratch: &mut Vec<f32>,
820    pool: Option<&Pool>,
821) -> (Vec<f32>, Vec<f32>) {
822    if state.len() != cfg.state_len() {
823        *state = vec![0f32; cfg.state_len()];
824    }
825    let h = cfg.hidden_size;
826    let mut bcx1 = vec![0.0f32; 3 * h];
827    let mut bcx2 = vec![0.0f32; 3 * h];
828    w.in_proj.matvec2(x1, x2, &mut bcx1, &mut bcx2, pool);
829
830    let mut y1 = vec![0.0f32; h];
831    short_conv_step(&bcx1, &w.conv, cfg, state, &mut y1);
832    scratch.clear();
833    scratch.extend_from_slice(state);
834    let mut y2 = vec![0.0f32; h];
835    short_conv_step(&bcx2, &w.conv, cfg, scratch, &mut y2);
836
837    let mut out1 = vec![0.0f32; h];
838    let mut out2 = vec![0.0f32; h];
839    w.out_proj.matvec2(&y1, &y2, &mut out1, &mut out2, pool);
840    (out1, out2)
841}
842
843#[cfg(test)]
844mod tests {
845    use super::*;
846
847    fn tiny() -> (VmfPhaseWeights, VmfPhaseCfg) {
848        let cfg = VmfPhaseCfg {
849            num_heads: 2,
850            nphase: 3,
851            value_head_dim: 4,
852            hidden_size: 8,
853            phase_mass: 0.0,
854        };
855        let synth = |rows: usize, cols: usize, salt: usize| {
856            QTensor::from_f32(
857                (0..rows * cols)
858                    .map(|i| (((i * 13 + salt * 7) % 97) as f32 / 97.0 - 0.5) * 0.4)
859                    .collect(),
860                rows,
861                cols,
862            )
863        };
864        let w = VmfPhaseWeights {
865            thq: synth(cfg.num_heads * cfg.nphase, cfg.hidden_size, 1),
866            thk: synth(cfg.num_heads * cfg.nphase, cfg.hidden_size, 2),
867            v_proj: synth(cfg.num_heads * cfg.value_head_dim, cfg.hidden_size, 3),
868            out_proj: synth(cfg.hidden_size, cfg.num_heads * cfg.value_head_dim, 4),
869            decay: (0..cfg.num_heads * 2 * cfg.nphase)
870                .map(|i| 0.9 + 0.005 * (i % 10) as f64)
871                .collect(),
872            k_gate: None,
873        };
874        (w, cfg)
875    }
876
877    #[test]
878    fn state_persists_and_changes_output() {
879        let (w, cfg) = tiny();
880        let x: Vec<f32> = (0..8).map(|i| (i as f32 * 0.3).sin()).collect();
881        let mut state = Vec::new();
882        let o1 = vmf_phase_forward(&x, &w, &cfg, &mut state, None);
883        let o2 = vmf_phase_forward(&x, &w, &cfg, &mut state, None);
884        // Same input, evolved condensate → different output.
885        assert!(o1.iter().zip(&o2).any(|(a, b)| (a - b).abs() > 1e-6));
886        assert_eq!(state.len(), cfg.state_len());
887    }
888
889    /// θ-mass (η′): mass=0 is bit-identical to the massless kernel; mass>0
890    /// changes the output (phase narrowed → kernel widened). Guards the
891    /// no-op default and that the knob is actually wired.
892    #[test]
893    fn phase_mass_zero_is_noop_and_positive_shifts() {
894        let (w, cfg0) = tiny();
895        let mut cfg_m = cfg0.clone();
896        cfg_m.phase_mass = 1.0;
897        let x: Vec<f32> = (0..8).map(|i| (i as f32 * 0.4).sin()).collect();
898
899        let mut s0 = Vec::new();
900        let base = vmf_phase_forward(&x, &w, &cfg0, &mut s0, None);
901        // Re-run with mass=0 → must be bit-identical.
902        let mut s0b = Vec::new();
903        let base2 = vmf_phase_forward(&x, &w, &cfg0, &mut s0b, None);
904        assert_eq!(base, base2, "mass=0 must be deterministic/no-op");
905        // mass=1 → output differs (θ halved before cos/sin).
906        let mut sm = Vec::new();
907        let massed = vmf_phase_forward(&x, &w, &cfg_m, &mut sm, None);
908        assert!(
909            base.iter().zip(&massed).any(|(a, b)| (a - b).abs() > 1e-5),
910            "mass>0 must change the output"
911        );
912        assert!(massed.iter().all(|v| v.is_finite()));
913    }
914
915    /// κ write gate (hybrid_k): saturated-open gate (bias ≫ 0 → κ→1)
916    /// matches the gateless kernel within fp tolerance; a closed gate
917    /// (bias ≪ 0 → κ→0) writes nothing — the state stays zero and the
918    /// output collapses to the empty-condensate readout.
919    #[test]
920    fn kappa_gate_open_matches_none_and_closed_writes_nothing() {
921        let (mut w, cfg) = tiny();
922        let x: Vec<f32> = (0..8).map(|i| (i as f32 * 0.3).sin()).collect();
923
924        let mut s_none = Vec::new();
925        let base1 = vmf_phase_forward(&x, &w, &cfg, &mut s_none, None);
926        let base2 = vmf_phase_forward(&x, &w, &cfg, &mut s_none, None);
927
928        // Open gate: W=0, bias=+20 → κ = σ(20) ≈ 1 − 2e−9.
929        w.k_gate = Some((
930            QTensor::from_f32(
931                vec![0.0; cfg.num_heads * cfg.hidden_size],
932                cfg.num_heads,
933                cfg.hidden_size,
934            ),
935            vec![20.0; cfg.num_heads],
936        ));
937        let mut s_open = Vec::new();
938        let o1 = vmf_phase_forward(&x, &w, &cfg, &mut s_open, None);
939        let o2 = vmf_phase_forward(&x, &w, &cfg, &mut s_open, None);
940        for (a, b) in base1.iter().zip(&o1).chain(base2.iter().zip(&o2)) {
941            assert!(
942                (a - b).abs() < 1e-5,
943                "open κ must match gateless: {a} vs {b}"
944            );
945        }
946
947        // Closed gate: bias=−20 → κ ≈ 0 → nothing is written.
948        w.k_gate = Some((
949            QTensor::from_f32(
950                vec![0.0; cfg.num_heads * cfg.hidden_size],
951                cfg.num_heads,
952                cfg.hidden_size,
953            ),
954            vec![-20.0; cfg.num_heads],
955        ));
956        let mut s_closed = Vec::new();
957        let oc = vmf_phase_forward(&x, &w, &cfg, &mut s_closed, None);
958        assert!(
959            s_closed.iter().all(|&v| v.abs() < 1e-7),
960            "closed κ: state must stay empty"
961        );
962        assert!(
963            oc.iter().all(|&v| v.abs() < 1e-6),
964            "closed κ: empty-condensate readout"
965        );
966    }
967
968    #[test]
969    fn pair_matches_two_singles_bitexact() {
970        let (w, cfg) = tiny();
971        let x1: Vec<f32> = (0..8).map(|i| (i as f32 * 0.2).cos()).collect();
972        let x2: Vec<f32> = (0..8).map(|i| (i as f32 * 0.5).sin()).collect();
973
974        // Reference: two sequential singles.
975        let mut s_ref = Vec::new();
976        let r1 = vmf_phase_forward(&x1, &w, &cfg, &mut s_ref, None);
977        let r2 = vmf_phase_forward(&x2, &w, &cfg, &mut s_ref, None);
978
979        // Pair: lane1 commits, lane2 tentative in scratch.
980        let mut s = Vec::new();
981        let mut scratch = Vec::new();
982        let (p1, p2) = vmf_phase_pair(&x1, &x2, &w, &cfg, &mut s, &mut scratch, None);
983        assert_eq!(r1, p1, "lane 1 must be bit-identical");
984        assert_eq!(r2, p2, "lane 2 must be bit-identical");
985        // Accepting the draft = swapping scratch in → equals s_ref.
986        std::mem::swap(&mut s, &mut scratch);
987        assert_eq!(s, s_ref, "accepted state must equal sequential state");
988    }
989
990    #[test]
991    fn rejected_draft_leaves_state_at_lane1() {
992        let (w, cfg) = tiny();
993        let x1: Vec<f32> = (0..8).map(|i| (i as f32 * 0.7).sin()).collect();
994        let x2 = vec![0.5f32; 8];
995
996        let mut s_ref = Vec::new();
997        let _ = vmf_phase_forward(&x1, &w, &cfg, &mut s_ref, None);
998
999        let mut s = Vec::new();
1000        let mut scratch = Vec::new();
1001        let _ = vmf_phase_pair(&x1, &x2, &w, &cfg, &mut s, &mut scratch, None);
1002        // Reject: state must be exactly the post-lane1 state.
1003        assert_eq!(s, s_ref);
1004    }
1005
1006    // ───────────── GatedDeltaNet ─────────────
1007
1008    fn tiny_gdn() -> (GdnWeights, GdnCfg) {
1009        let cfg = GdnCfg {
1010            num_v_heads: 4,
1011            num_k_heads: 2,
1012            key_head_dim: 3,
1013            value_head_dim: 5,
1014            conv_kernel: 4,
1015            hidden_size: 8,
1016            rms_eps: 1e-6,
1017        };
1018        let c_dim = cfg.conv_dim();
1019        let vd = cfg.num_v_heads * cfg.value_head_dim;
1020        let synth = |rows: usize, cols: usize, salt: usize| {
1021            QTensor::from_f32(
1022                (0..rows * cols)
1023                    .map(|i| (((i * 13 + salt * 7) % 97) as f32 / 97.0 - 0.5) * 0.4)
1024                    .collect(),
1025                rows,
1026                cols,
1027            )
1028        };
1029        let vecf = |n: usize, salt: usize| -> Vec<f32> {
1030            (0..n)
1031                .map(|i| (((i * 11 + salt * 5) % 89) as f32 / 89.0 - 0.5) * 0.6)
1032                .collect()
1033        };
1034        let w = GdnWeights {
1035            in_proj_qkv: synth(c_dim, cfg.hidden_size, 1),
1036            in_proj_z: synth(vd, cfg.hidden_size, 2),
1037            in_proj_a: synth(cfg.num_v_heads, cfg.hidden_size, 3),
1038            in_proj_b: synth(cfg.num_v_heads, cfg.hidden_size, 4),
1039            conv1d: vecf(c_dim * cfg.conv_kernel, 5),
1040            a_log: (0..cfg.num_v_heads).map(|i| 0.2 + 0.3 * i as f32).collect(),
1041            dt_bias: vecf(cfg.num_v_heads, 6),
1042            norm: vec![1.0; cfg.value_head_dim],
1043            out_proj: synth(cfg.hidden_size, vd, 7),
1044        };
1045        (w, cfg)
1046    }
1047
1048    #[test]
1049    fn gdn_state_persists_and_changes_output() {
1050        let (w, cfg) = tiny_gdn();
1051        let x: Vec<f32> = (0..8).map(|i| (i as f32 * 0.3).sin()).collect();
1052        let mut state = Vec::new();
1053        let o1 = gdn_forward(&x, &w, &cfg, &mut state, None);
1054        let o2 = gdn_forward(&x, &w, &cfg, &mut state, None);
1055        assert!(o1.iter().zip(&o2).any(|(a, b)| (a - b).abs() > 1e-6));
1056        assert_eq!(state.len(), cfg.state_len());
1057    }
1058
1059    #[test]
1060    fn gdn_pair_matches_two_singles_bitexact() {
1061        let (w, cfg) = tiny_gdn();
1062        let x1: Vec<f32> = (0..8).map(|i| (i as f32 * 0.2).cos()).collect();
1063        let x2: Vec<f32> = (0..8).map(|i| (i as f32 * 0.5).sin()).collect();
1064
1065        let mut s_ref = Vec::new();
1066        let r1 = gdn_forward(&x1, &w, &cfg, &mut s_ref, None);
1067        let r2 = gdn_forward(&x2, &w, &cfg, &mut s_ref, None);
1068
1069        let mut s = Vec::new();
1070        let mut scratch = Vec::new();
1071        let (p1, p2) = gdn_pair(&x1, &x2, &w, &cfg, &mut s, &mut scratch, None);
1072        assert_eq!(r1, p1, "lane 1 must be bit-identical");
1073        assert_eq!(r2, p2, "lane 2 must be bit-identical");
1074        std::mem::swap(&mut s, &mut scratch);
1075        assert_eq!(s, s_ref, "accepted state must equal sequential state");
1076    }
1077
1078    #[test]
1079    fn gdn_rejected_draft_leaves_state_at_lane1() {
1080        let (w, cfg) = tiny_gdn();
1081        let x1: Vec<f32> = (0..8).map(|i| (i as f32 * 0.7).sin()).collect();
1082        let x2 = vec![0.5f32; 8];
1083
1084        let mut s_ref = Vec::new();
1085        let _ = gdn_forward(&x1, &w, &cfg, &mut s_ref, None);
1086
1087        let mut s = Vec::new();
1088        let mut scratch = Vec::new();
1089        let _ = gdn_pair(&x1, &x2, &w, &cfg, &mut s, &mut scratch, None);
1090        assert_eq!(s, s_ref);
1091    }
1092
1093    /// The conv ring must give the same result as an explicit causal
1094    /// conv over the whole sequence (oracle semantics: zero left-pad,
1095    /// tap kk−1 on the current position).
1096    #[test]
1097    fn gdn_conv_ring_matches_explicit_causal_conv() {
1098        let (w, cfg) = tiny_gdn();
1099        let seq: Vec<Vec<f32>> = (0..6)
1100            .map(|t| (0..8).map(|i| ((t * 8 + i) as f32 * 0.17).sin()).collect())
1101            .collect();
1102
1103        // Reference: recompute position t from scratch each time with a
1104        // fresh state built by replaying the prefix.
1105        let mut s_inc = Vec::new();
1106        for (t, x) in seq.iter().enumerate() {
1107            let inc = gdn_forward(x, &w, &cfg, &mut s_inc, None);
1108            let mut s_replay = Vec::new();
1109            let mut replay = Vec::new();
1110            for xr in &seq[..=t] {
1111                replay = gdn_forward(xr, &w, &cfg, &mut s_replay, None);
1112            }
1113            assert_eq!(inc, replay, "position {t}: ring must equal replay");
1114        }
1115    }
1116
1117    fn tiny_short_conv() -> (ShortConvWeights, ShortConvCfg) {
1118        let cfg = ShortConvCfg {
1119            hidden_size: 8,
1120            kernel: 3,
1121        };
1122        let synth = |rows: usize, cols: usize, salt: usize| {
1123            QTensor::from_f32(
1124                (0..rows * cols)
1125                    .map(|i| (((i * 11 + salt * 5) % 89) as f32 / 89.0 - 0.5) * 0.5)
1126                    .collect(),
1127                rows,
1128                cols,
1129            )
1130        };
1131        let w = ShortConvWeights {
1132            in_proj: synth(3 * cfg.hidden_size, cfg.hidden_size, 1),
1133            conv: (0..cfg.hidden_size * cfg.kernel)
1134                .map(|i| ((i * 7 % 13) as f32 / 13.0 - 0.5) * 0.8)
1135                .collect(),
1136            out_proj: synth(cfg.hidden_size, cfg.hidden_size, 2),
1137        };
1138        (w, cfg)
1139    }
1140
1141    /// The incremental conv ring must equal a from-scratch causal replay
1142    /// of the prefix at every position — the decode/prefill contract.
1143    #[test]
1144    fn short_conv_ring_matches_explicit_causal_conv() {
1145        let (w, cfg) = tiny_short_conv();
1146        let seq: Vec<Vec<f32>> = (0..6)
1147            .map(|t| (0..8).map(|i| ((t * 8 + i) as f32 * 0.19).cos()).collect())
1148            .collect();
1149        let mut s_inc = Vec::new();
1150        for (t, x) in seq.iter().enumerate() {
1151            let inc = short_conv_forward(x, &w, &cfg, &mut s_inc, None);
1152            let mut s_replay = Vec::new();
1153            let mut replay = Vec::new();
1154            for xr in &seq[..=t] {
1155                replay = short_conv_forward(xr, &w, &cfg, &mut s_replay, None);
1156            }
1157            assert_eq!(inc, replay, "position {t}: ring must equal replay");
1158            assert_eq!(s_inc.len(), cfg.state_len());
1159        }
1160    }
1161
1162    /// The batched prefill path (matmat + sequential conv over the chunk)
1163    /// must reproduce the position-by-position decode path exactly.
1164    #[test]
1165    fn short_conv_batch_matches_sequential() {
1166        let (w, cfg) = tiny_short_conv();
1167        let b = 5;
1168        let xs: Vec<f32> = (0..b * cfg.hidden_size)
1169            .map(|i| (i as f32 * 0.13).sin() * 0.6)
1170            .collect();
1171
1172        let mut s_seq = Vec::new();
1173        let mut seq_out = vec![0.0f32; b * cfg.hidden_size];
1174        for bi in 0..b {
1175            let o = short_conv_forward(
1176                &xs[bi * cfg.hidden_size..(bi + 1) * cfg.hidden_size],
1177                &w,
1178                &cfg,
1179                &mut s_seq,
1180                None,
1181            );
1182            seq_out[bi * cfg.hidden_size..(bi + 1) * cfg.hidden_size].copy_from_slice(&o);
1183        }
1184
1185        let mut s_batch = Vec::new();
1186        let batch_out = short_conv_forward_batch(&xs, b, &w, &cfg, &mut s_batch, None);
1187        assert_eq!(
1188            seq_out, batch_out,
1189            "batch conv must match sequential decode"
1190        );
1191        assert_eq!(s_seq, s_batch, "ring state must match after the chunk");
1192    }
1193}