feat(fp16): accumulation-hardness PoUW for A100 (header-bound ZK consensus certificate + sm_80 miner pipeline)

Taskfile.yml

@@ -6,15 +6,25 @@ env:
 
 tasks:
   build:zk-cache:
-    desc: Generate verifier circuit caches (main + v1 + fp8, one-time)
+    desc: Generate verifier circuit caches (main + v1 + fp8 + fp16, one-time)
     run: once
     dir: zk-pow
+    env:
+      # The FP16 wrapper verifier cache is embedded via include_bytes! by the Go
+      # binding's `embedded_cache` feature, so fp16_cache.bin must exist for the
+      # binding to compile. FP16_CACHE_SAMPLE=1 builds the single smallest legal
+      # profile — a cheap bootstrap blob that makes the crate compile and exercises
+      # the load path. The full per-profile FP16 envelope is a heavy (multi-GB,
+      # minutes-per-profile) offline build required only before V5 is activated on a
+      # network; drop this env (or run build_cache with the explicit fp16 path) for it.
+      FP16_CACHE_SAMPLE: "1"
     status:
       - test -s src/v2/circuit/v2_cache.bin
       - test -s src/v1/v1_cache.bin
       - test -s src/api/fp8/fp8_cache.bin
+      - test -s src/api/fp16/fp16_cache.bin
     cmds:
-      - cargo run --release --no-default-features --bin build_cache src/api/fp8/fp8_cache.bin src/v2/circuit/v2_cache.bin src/v1/v1_cache.bin
+      - cargo run --release --no-default-features --bin build_cache src/api/fp8/fp8_cache.bin src/v2/circuit/v2_cache.bin src/v1/v1_cache.bin src/api/fp16/fp16_cache.bin
 
   build:zk-gobind:
     desc: Build the Rust ZK verification library for Go FFI
@@ -92,7 +102,7 @@ tasks:
     cmds:
       - rm -rf bin
       - rm -f xmss/libxmss.a
-      - rm -f zk-pow/src/api/fp8/fp8_cache.bin zk-pow/src/v2/circuit/v2_cache.bin zk-pow/src/v1/v1_cache.bin
+      - rm -f zk-pow/src/api/fp8/fp8_cache.bin zk-pow/src/v2/circuit/v2_cache.bin zk-pow/src/v1/v1_cache.bin zk-pow/src/api/fp16/fp16_cache.bin
       - rm -rf pearl-blake3/target py-pearl-mining/target zk-pow/target zk-pow/bindings/go/target
       - find . -name coverage.txt -type f -delete
 

docs/fp16_scheme/.gitignore

@@ -0,0 +1,12 @@
+*.aux
+*.log
+*.out
+*.toc
+*.fls
+*.fdb_latexmk
+*.synctex.gz
+
+# hardware capture build artifact
+validation/a100_hmma_capture
+validation/attack_cost_benchmark
+validation/__pycache__/

docs/fp16_scheme/README.md

@@ -0,0 +1,32 @@
+# Pearl FP16 scheme specification
+
+`fp16_scheme.tex` specifies an alternative Pearl proof-of-useful-work instantiation whose unit
+of work is **FP16** matrix multiplication on NVIDIA A100 (GA100, `sm_80`) tensor cores. It
+extends the FP8 certificate-v4 specification rather than replacing it: commitments, seed chain,
+state window, ticket, target, MoE extension, and the ZK split are shared, and FP8 remains a
+separate `Quant` value.
+
+The new idea: hardness comes from the **nonlinearity of the device's accumulation**, not from a
+coarse rounding grid. The A100 tensor core truncates each product onto a per-group alignment
+grid *before* summing (groups of 8, 24-bit window, round-toward-zero per group). That per-product
+truncation does not commute with the reduction, so it cannot be expressed as a matrix
+multiplication — the same obstruction that makes the integer transcript scheme hard, applied for
+free on every MAC. This lets recovered products carry FP16-level accuracy instead of
+FP8-residual accuracy.
+
+The experiments the spec's appendices summarize are reproduced by the committed harness in
+[`validation/`](validation/) (A100 HMMA accumulation-model capture + model-vs-silicon cross-check,
+RZ-vs-RNE resolution, policy `f_bp`/`rho` calibration, and the truncation-attack cost benchmark).
+See [`validation/README.md`](validation/README.md) for how to run it on `sm_80` silicon and the
+recorded results.
+
+## Building the PDF
+
+The checked-in `fp16_scheme.pdf` is the authoritative rendering. To rebuild:
+
+```sh
+# Any LaTeX engine works; the doc uses only amsmath/amssymb/booktabs/hyperref.
+tectonic fp16_scheme.tex        # self-contained, recommended
+# or
+latexmk -pdf fp16_scheme.tex
+```

docs/fp16_scheme/validation/README.md

@@ -0,0 +1,118 @@
+# FP16 / A100 scheme — reproducible hardware validation
+
+This directory is the **reproducible hardware oracle and calibration harness** for
+the FP16 accumulation-hardness scheme (`docs/fp16_scheme/fp16_scheme.tex`). It
+exists to close three gaps a review of the scheme flagged:
+
+1. the committed reference corpus
+   (`zk-pow/src/api/fp16/testdata/a100_dot_vectors.txt`) topped out at **k=256**
+   and carried no explicit subnormal / cancellation / overflow edges;
+2. the accumulation model was validated only **circularly** in-tree
+   (`accumulate.rs`, the Python reference, and the vector file are three
+   expressions of the *same* model) — nothing checked it against the **device**;
+3. the **RZ-vs-RNE** rounding question (does the A100 round toward zero, as the
+   model assumes, or to nearest-even?) had no in-tree silicon resolution, and the
+   Appendix "Empirical validation" numbers (f_bp / rho distributions, shortcut
+   rates) were prose with **no runnable harness**.
+
+Everything here runs on a real **sm_80** GPU (A100 / GA100 / CMP 170HX) with a
+standard CUDA toolkit (tested: nvcc 12.8, CMP 170HX 64 GB, driver 610.43.02). It
+deliberately links only the CUDA runtime — **no torch, no pearl-gemm** — so it
+builds and runs on an sm_80 box that cannot install the py3.12 miner stack.
+
+## Contents
+
+| file | what it does |
+|---|---|
+| `a100_hmma_capture.cu` | Standalone capture tool. Issues the real `mma.sync.m16n8k16.f32.f16` (asm + fragment layout copied verbatim from the validated production kernel `fp16_gemm/_kernel_sm80.cu`) and dumps device FP32 result bits for each input dot product. |
+| `generate_and_capture.py` | Generates an edge-case corpus (k up to 1024; subnormal operands / accumulators / outputs; cancellation; grid-boundary; max-magnitude), captures it on silicon, and **cross-checks the software model against the device** bit-for-bit. Writes `a100_dot_vectors_k1024.txt`. |
+| `rounding_mode_probe.py` | Builds four model variants (per-product ∈ {RZ,RNE} × per-group ∈ {RZ,RNE}), finds vectors where they disagree, captures them on silicon, and reports which rounding the device matches. **Resolves RZ-vs-RNE empirically.** |
+| `calibration.py` | Reproduces the Appendix policy-metric study: per-strategy `f_bp`, `rho`, and shortcut-predictor reproduction rates, checked against the gate (`f_bp>=0.30`, `rho>=1.2`) and the paper's claimed regime. CPU only. `--rust-noise` drives the **real consensus noise pipeline** (next row) instead of a Python reimplementation. |
+| `attack_cost_benchmark.cu` | Times the honest tensor-core GEMM against the per-product truncation-correction kernel (which must run on CUDA cores) on silicon, measuring the attacker's cost multiplier — the paper's "6-13x / memory-bound far more" claim. |
+| `../../../zk-pow/examples/fp16_noise_tool.rs` | Rust example that emits seed-exact noised operands through the real `fp16_noised_operands` pipeline (root-derived seeds -> BLAKE3 noise lines -> `noisy_quantize`); `calibration.py --rust-noise` shells out to it. Build: `cd zk-pow && cargo build --release --example fp16_noise_tool`. |
+
+## How to run
+
+```sh
+cd docs/fp16_scheme/validation
+nvcc -arch=sm_80 -O2 -o a100_hmma_capture a100_hmma_capture.cu   # build once
+python3 generate_and_capture.py        # edge corpus + model-vs-silicon cross-check
+python3 rounding_mode_probe.py          # RZ-vs-RNE resolution on silicon
+python3 calibration.py                  # policy-metric calibration (Python noise; no GPU)
+
+# seed-exact consensus noise + the attack-cost benchmark:
+(cd ../../../zk-pow && cargo build --release --example fp16_noise_tool)
+python3 calibration.py --rust-noise     # same calibration, real consensus noise bytes
+nvcc -arch=sm_80 -O3 -o attack_cost_benchmark attack_cost_benchmark.cu
+./attack_cost_benchmark 512 512 1024 50 # honest vs truncation-correction cost
+```
+
+The regenerated `a100_dot_vectors_k1024.txt` is copied to
+`zk-pow/src/api/fp16/testdata/` and consumed by the Rust cross-check test
+`api::fp16::accumulate::tests::matches_hardware_capture_k1024`, so CI enforces
+"model == captured silicon bits" (including a k=1024 and a subnormal-output
+assertion). The capture binary itself is a build artifact (git-ignored).
+
+## Recorded results (on CMP 170HX, sm_80, nvcc 12.8)
+
+**Model vs silicon — bit-exact.** The capture tool reproduces all **238**
+pre-existing committed vectors bit-for-bit, and the model reproduces all **272**
+newly captured edge vectors (k up to 1024; 16 of them subnormal-FP32 outputs;
+64 subnormal-FP32 carry-ins). Total: **910** independent silicon dot products
+agree with the model, **0 mismatches** (238 committed + 272 edge + 400
+rounding-discriminating).
+
+**RZ-vs-RNE — RZ confirmed at both stages.** On **400** vectors constructed so the
+four rounding variants disagree, the device matched the committed `(RZ, RZ)`
+model on **400/400**; the alternatives matched only where they happen to coincide
+with RZ and were excluded elsewhere (`RZ/RNE` 113/400, `RNE/RZ` 109/400,
+`RNE/RNE` 60/400). **Conclusion: the A100/GA100 `HMMA.16816.F32` rounds toward
+zero** both when aligning each product/accumulator onto the `2^(eta-24)` grid and
+when encoding the per-group sum to FP32 — exactly what `accumulate.rs` models.
+
+**Overflow is unreachable from a single FP16 dot product.** Max-magnitude
+(|operand| up to 65504) k=1024 tiles with ±max-finite FP32 carry-in all stayed
+finite on silicon; no generated vector overflowed. The accumulation guard against
+non-finite results is therefore defensive, consistent with the feasibility note
+that FP16 products align at `eta >= 99` and the carry model forbids a cell-start
+carry-in.
+
+**Subnormal-FP32 output is reachable only via a surviving subnormal carry-in.**
+With FP16 operands the smallest nonzero product is `2^-48`, so any nonzero product
+dominates a subnormal (`< 2^-126`) accumulator and pushes the output normal; a
+nonzero subnormal output arises only when a subnormal FP32 carry-in survives
+products that are zero or exactly cancel. The corpus includes 16 such cases
+(captured and model-matched), which exercise the device's subnormal encode — the
+`matmul_a100_stark` MA13 branch the ZK audit noted was otherwise silicon-untested.
+This confirms MA13 is a sound guard on the reachable sub-case.
+
+**Policy calibration (`calibration.py`, 16×16 tile, k=1024, rank-32 δ=1/2 noise).**
+Across 8 adversarial/realistic strategies, with the **shape-faithful Python**
+noise: `f_bp ∈ [0.453, 0.512]`, `rho ∈ [1.67, 1.94]`, shortcut reproduction
+≤ 1.56%. With the **seed-exact Rust consensus** noise (`--rust-noise`, via
+`fp16_noise_tool`): `f_bp ∈ [0.448, 0.518]`, `rho ∈ [1.674, 1.939]`, shortcut
+≤ 2.34%. Both match the paper's claimed regime (`f_bp ∈ [0.36,0.51]`,
+`rho ≥ 1.38`, cheapest shortcut ≤ ~1.6%), all well above the `0.30 / 1.2` gate —
+and the two noise sources agree closely, so the Python proxy was faithful.
+
+**Truncation-correction attack cost (`attack_cost_benchmark.cu`, on CMP 170HX).**
+Granting the attacker the exact sums and the breakpoint oracle for free, the
+unavoidable per-product correction (decode + multiply + per-group truncation on
+CUDA cores) measured **37.6×** the honest bit-exact tensor-core GEMM at
+`512×512×1024` and **43.0×** at `1024³`. The honest baseline is the miner's own
+kernel style (one warp per 16×8 subtile, no shared-memory pipelining — exactly
+`fp16_gemm/_kernel_sm80.cu`), so the ratio is representative, not pessimistic. It
+sits well above the whitepaper's "6–13× compute-bound-optimistic" floor and below
+its "220–470× memory-bound" ceiling, consistent with the hardness claim that the
+corrections cannot be batched onto tensor cores.
+
+## Scope / remaining
+
+- The accumulation model and both roundings are silicon-validated; the
+  policy-metric distributions are reproducible with **both** a Python proxy and
+  the **real consensus noise pipeline** (`--rust-noise`); and the attack-cost
+  multiplier is measured on silicon. The three gaps the review flagged are closed.
+- The attack benchmark measures the *lower-bound* elementwise correction cost; a
+  production attacker's kernel could be tuned, but so could the honest one (which
+  benefits from tensor cores — the asymmetry the scheme relies on). The strategy
+  set in `calibration.py` is extensible for broader adversarial sweeps.

docs/fp16_scheme/validation/a100_hmma_capture.cu

@@ -0,0 +1,160 @@
+// Hardware capture of the A100/GA100 (sm_80) HMMA.16816.F32 accumulation, for
+// the FP16 proof-of-useful-work scheme (docs/fp16_scheme). This is the
+// *hardware* oracle: it issues the real tensor-core `mma.sync` on silicon and
+// dumps the FP32 result bits, so the software model in
+// `zk-pow/src/api/fp16/accumulate.rs` can be checked against the device rather
+// than against another copy of itself.
+//
+// Deliberately standalone: it links only the CUDA runtime (no torch, no
+// pearl-gemm), so it builds and runs on an sm_80 box that cannot install the
+// py3.12 miner stack. The mma asm and the distributed fragment layout are copied
+// verbatim from the validated production kernel
+// `miner/pearl-gemm/src/pearl_gemm/fp16_gemm/_kernel_sm80.cu`, so a capture here
+// exercises exactly the datapath the miner uses.
+//
+// Build:  nvcc -arch=sm_80 -O2 -o a100_hmma_capture a100_hmma_capture.cu
+//
+// Input (stdin), one dot product per line:
+//     k  a_bits[0..k)  b_bits[0..k)  c_bits
+// where a_bits/b_bits are FP16 bit patterns (u16) and c_bits is the FP32
+// carry-in bit pattern (u32). `k` must be a positive multiple of 16.
+//
+// Output (stdout), one line per input: the FP32 result bit pattern (u32),
+// captured from the device D[0][0] of a 16x8 tile whose A-row 0 is `a`, B-row 0
+// is `b`, and C[0][0] is the carry-in (all other tile entries zero).
+//
+// Exit non-zero on any CUDA error or malformed input.
+
+#include <cuda_fp16.h>
+#include <cstdint>
+#include <cstdio>
+#include <cstdlib>
+#include <vector>
+#include <string>
+
+#define CUDA_CHECK(expr)                                                        \
+    do {                                                                        \
+        cudaError_t _e = (expr);                                                \
+        if (_e != cudaSuccess) {                                                \
+            fprintf(stderr, "CUDA error %s at %s:%d\n",                         \
+                    cudaGetErrorString(_e), __FILE__, __LINE__);                \
+            std::exit(2);                                                       \
+        }                                                                       \
+    } while (0)
+
+__device__ __forceinline__ unsigned pack2(const __half* p, int i0, int i1) {
+    __half2 h = __halves2half2(p[i0], p[i1]);
+    return *reinterpret_cast<unsigned*>(&h);
+}
+
+// One warp computes one 16x8 output subtile D = A(16xK) . B(8xK)^T with FP32
+// accumulation chained in ascending k order (groups of 16 per mma.sync, no
+// split-k, no atomics), carrying the accumulator forward across the whole k
+// axis -- the pinned reduction order that defines the device result. Identical
+// to fp16_gemm_a100_kernel but single-tile (M=16, N=8).
+__global__ void capture_kernel(const __half* __restrict__ A,
+                               const __half* __restrict__ B,
+                               const float* __restrict__ C,
+                               float* __restrict__ D, int K) {
+    int lane = threadIdx.x & 31;
+    int gid = lane >> 2;  // 0..7
+    int t4 = lane & 3;    // 0..3
+
+    float c0 = 0.f, c1 = 0.f, c2 = 0.f, c3 = 0.f;
+    if (C != nullptr) {
+        c0 = C[(gid) * 8 + t4 * 2];
+        c1 = C[(gid) * 8 + t4 * 2 + 1];
+        c2 = C[(gid + 8) * 8 + t4 * 2];
+        c3 = C[(gid + 8) * 8 + t4 * 2 + 1];
+    }
+    for (int k0 = 0; k0 < K; k0 += 16) {
+        unsigned a0 = pack2(A, (gid) * K + k0 + t4 * 2, (gid) * K + k0 + t4 * 2 + 1);
+        unsigned a1 = pack2(A, (gid + 8) * K + k0 + t4 * 2, (gid + 8) * K + k0 + t4 * 2 + 1);
+        unsigned a2 = pack2(A, (gid) * K + k0 + t4 * 2 + 8, (gid) * K + k0 + t4 * 2 + 9);
+        unsigned a3 = pack2(A, (gid + 8) * K + k0 + t4 * 2 + 8, (gid + 8) * K + k0 + t4 * 2 + 9);
+        unsigned b0 = pack2(B, (gid) * K + k0 + t4 * 2, (gid) * K + k0 + t4 * 2 + 1);
+        unsigned b1 = pack2(B, (gid) * K + k0 + t4 * 2 + 8, (gid) * K + k0 + t4 * 2 + 9);
+        asm volatile(
+            "mma.sync.aligned.m16n8k16.row.col.f32.f16.f16.f32 "
+            "{%0,%1,%2,%3}, {%4,%5,%6,%7}, {%8,%9}, {%0,%1,%2,%3};\n"
+            : "+f"(c0), "+f"(c1), "+f"(c2), "+f"(c3)
+            : "r"(a0), "r"(a1), "r"(a2), "r"(a3), "r"(b0), "r"(b1));
+    }
+    D[(gid) * 8 + t4 * 2] = c0;
+    D[(gid) * 8 + t4 * 2 + 1] = c1;
+    D[(gid + 8) * 8 + t4 * 2] = c2;
+    D[(gid + 8) * 8 + t4 * 2 + 1] = c3;
+}
+
+int main() {
+    // Persistent device buffers sized to the largest tile we expect.
+    int cap_k = 0;
+    __half *dA = nullptr, *dB = nullptr;
+    float *dC = nullptr, *dD = nullptr;
+    std::vector<__half> hA, hB;
+    float hC[128], hD[128];
+
+    // Each record is whitespace-separated: k, then k a-bits, k b-bits, c-bits.
+    long k;
+    while (scanf("%ld", &k) == 1) {
+        if (k <= 0) {
+            fprintf(stderr, "bad k=%ld (must be positive)\n", k);
+            return 1;
+        }
+        // The mma tile runs in units of 16 along the contraction axis. Any k is
+        // zero-padded up to the next multiple of 16; within each group of 8 the
+        // padded entries are zero (skipped by the model), and any whole trailing
+        // all-zero group is a no-op, so the padded capture equals the model's
+        // result for the true k -- which also validates those properties on
+        // silicon (the committed corpus includes k not a multiple of 16).
+        long kpad = (k + 15) / 16 * 16;
+        std::vector<uint16_t> a(k), b(k);
+        uint32_t cbits;
+        for (long i = 0; i < k; i++) {
+            unsigned v;
+            if (scanf("%u", &v) != 1) { fprintf(stderr, "truncated a\n"); return 1; }
+            a[i] = (uint16_t)v;
+        }
+        for (long i = 0; i < k; i++) {
+            unsigned v;
+            if (scanf("%u", &v) != 1) { fprintf(stderr, "truncated b\n"); return 1; }
+            b[i] = (uint16_t)v;
+        }
+        if (scanf("%u", &cbits) != 1) { fprintf(stderr, "truncated c\n"); return 1; }
+
+        if (kpad > cap_k) {
+            if (dA) { cudaFree(dA); cudaFree(dB); }
+            CUDA_CHECK(cudaMalloc(&dA, sizeof(__half) * 16 * kpad));
+            CUDA_CHECK(cudaMalloc(&dB, sizeof(__half) * 8 * kpad));
+            if (!dC) {
+                CUDA_CHECK(cudaMalloc(&dC, sizeof(float) * 128));
+                CUDA_CHECK(cudaMalloc(&dD, sizeof(float) * 128));
+            }
+            hA.resize(16 * kpad);
+            hB.resize(8 * kpad);
+            cap_k = kpad;
+        }
+        // Zero the tile, place the vectors in row 0 of A and B, carry-in at [0][0].
+        for (auto& h : hA) h = __ushort_as_half((unsigned short)0);
+        for (auto& h : hB) h = __ushort_as_half((unsigned short)0);
+        for (long i = 0; i < k; i++) hA[i] = __ushort_as_half((unsigned short)a[i]);
+        for (long i = 0; i < k; i++) hB[i] = __ushort_as_half((unsigned short)b[i]);
+        for (int i = 0; i < 128; i++) hC[i] = 0.f;
+        float cf;
+        memcpy(&cf, &cbits, 4);
+        hC[0] = cf;
+
+        CUDA_CHECK(cudaMemcpy(dA, hA.data(), sizeof(__half) * 16 * kpad, cudaMemcpyHostToDevice));
+        CUDA_CHECK(cudaMemcpy(dB, hB.data(), sizeof(__half) * 8 * kpad, cudaMemcpyHostToDevice));
+        CUDA_CHECK(cudaMemcpy(dC, hC, sizeof(float) * 128, cudaMemcpyHostToDevice));
+        capture_kernel<<<1, 32>>>(dA, dB, dC, dD, (int)kpad);
+        CUDA_CHECK(cudaGetLastError());
+        CUDA_CHECK(cudaDeviceSynchronize());
+        CUDA_CHECK(cudaMemcpy(hD, dD, sizeof(float) * 128, cudaMemcpyDeviceToHost));
+
+        uint32_t dbits;
+        memcpy(&dbits, &hD[0], 4);
+        printf("%u\n", dbits);
+    }
+    return 0;
+}