Skip to content

Latest commit

 

History

History
1818 lines (1489 loc) · 104 KB

codebase_notes.md

File metadata and controls

1818 lines (1489 loc) · 104 KB

Codebase Notes {#codebase_notes}

These maintained notes preserve implementation rationale and historical investigations referenced by source. Experimental paths and research assets mentioned in older notes are retained on development; their mention here does not make them available in stable main. The branch policy and current stage catalog define the supported build scope.

Longform engineering notes that back up decisions in the source: measurement tables, before/after numbers, bug postmortems, and "why not the obvious thing" arguments.

Why this file exists

A tuning constant or an unusual kernel launch shape is only defensible with the evidence behind it, and that evidence is bulky — measurement tables across eight fields, a bug's full symptom trail, the three approaches that were tried first. Keeping all of it inline pushes the actual code off the screen: several stages had 40-110 line comment blocks where the contract worth reading at the call site was four lines of it.

So the split is by audience and lifetime:

Stays in the source Lives here
The invariant or contract a caller must not violate The measurements that justified the constant
The rule ("target pardeg 131072, floor at 768") The table showing what each value scored
A one-line warning ("do not change without re-measuring") Why it is not safe to extrapolate, with the curve
Why this line is not the obvious thing, in a sentence The full postmortem of the bug that made it non-obvious

A note is worth writing when someone would otherwise redo the work — repeat a benchmark sweep, re-derive a race, or "simplify" a shape back into a bug.

Conventions

Every note has a stable ID (CN-<AREA>-<n>) used as its anchor. The source carries a short comment ending in a pointer:

// Grid is (kGatherBlocksPerSegment, n_segs), not one block per segment.
// Rationale and measurements: docs/codebase_notes.md CN-CONCAT-1

Rules that keep this from rotting:

  • Never move a contract here. If violating it breaks the build or corrupts data, it belongs at the call site.
  • IDs are permanent. Superseded notes are marked, not deleted or renumbered, so an old pointer still resolves.
  • Date and pin measurements to hardware, input, and configuration. An unattributed number cannot be reproduced or refuted.
  • Deleting the code deletes the note. A note about a kernel that no longer exists is a trap.

CN-BUILD-1 — HIP toolchain and source-language boundaries

Source: CMakeLists.txt, include/backend/, modules/coders/ans/CMakeLists.txt

CMake enables the HIP language for GPU sources while keeping ordinary .cpp sources in C++. The runtime target is hip::host: hip::device carries an interface -x hip flag that would also reach plain C++ sources, where GCC rejects it as an unknown language. CMake supplies the HIP compilation flags to sources explicitly tagged LANGUAGE HIP. Existing .cu filenames are retained and retagged rather than renamed; this path was exercised with ROCm 6.4 on MI100.

hipCUB, rocPRIM and rocThrust are header-only, but their imported targets can transitively link hip::device. The build therefore takes their include directories without linking those targets and marks the headers as system dependencies. The ROCm root must be discoverable through CMAKE_PREFIX_PATH. Explicit GPU_TARGETS avoids reliance on GPU autodetection when configuring on a build node.

The backend define selects the shared API/type/algorithm facades. HIP's sync warp intrinsics require HIP_ENABLE_WARP_SYNC_BUILTINS; teardown paths also intentionally ignore some runtime status results, whose HIP declarations are [[nodiscard]].

The vendored dietgpu ANS implementation uses inline NVPTX and a fixed 32-lane mapping. Its sources are excluded from HIP builds. ANSStage::isSupportedOnBackend() and the compile-time guard in Pipeline::addStage() reject that unavailable stage without instantiating references to its excluded constructor. Earlier CMake comments describing this guard as unfinished were stale.

CN-BUILD-2 — NVRTC headers and driver-link propagation

Source: modules/fused/common/CMakeLists.txt, modules/fused/common/nvrtc_jit.cpp, scripts/gen_nvrtc_headers.py

The build embeds the shipped device harnesses and their transitive project headers for in-memory NVRTC compilation. This keeps installed specialization independent of its source checkout and pins the device source to the library build. The generator runs every build to pick up edits anywhere in that include closure, but leaves its output timestamp unchanged when the contents match. See the CN-NVRTC-FUSE note for the implementation history and runtime compiler constraints.

CUDA::cuda_driver is a public dependency of fzgmod_modules. Its shared library references driver functions such as cuModuleLoadData and cuLaunchKernel, and consumers must resolve those references even in GPU-less build containers using the driver stub. Keeping this dependency private previously caused executable link failures. NVRTC itself stays private because its library resolves at build time.

CN-BUILD-3 — stage fingerprints and build provenance

Source: CMakeLists.txt, scripts/gen_stage_fingerprints.py

Stage fingerprints hash each registered stage's sources and their transitive project includes. The registry's source directories define the stage population. Consumers can use these fingerprints to distinguish changed stage implementations between builds rather than relying only on pipeline configuration hashes.

Generation runs on every build, because configure-time generation would miss a kernel edit followed by an ordinary rebuild. The generator rewrites the header only when the fingerprints change, avoiding unnecessary recompilation. The CLI's version and Git revision identify the build separately; shared-library SOVERSION follows the major release for ABI compatibility, while VERSION records the full release. The private monolithic SZp reference implementation remains linked to retain its factory for decoding existing FZM archives.

CN-ARCHIVE-1 — logical archive lengths versus padded execution

Source: src/pipeline/compressor_exec.cpp, src/pipeline/compressor_io.cpp

2026-10-12. A partial PFPL field exposed a mismatch between in-memory and file decode: execution padded the source and serialization recorded that physical length, while only the live pipeline retained the caller's original length. The fix writes logical lengths into existing core/per-source header fields. Stage configurations and compressed streams still drive physical inverse sizing; the file reader trims the copied result rather than reducing decode storage. Header-primed decode restores logical lengths for every source and clears stale trimming state when explicitly prepared for a different input.

No layout/version change is required. Legacy headers retain their stored padded lengths; a generic reader cannot reliably recover information those writers omitted. Direct caller-owned decode still requires physical capacity and refuses smaller buffers before launching. This does not solve runtime-size propagation through interior graph edges, which is a separate execution-planning project.

Regression coverage includes Off/Auto partial chunks, aligned reuse after a partial input, old padded metadata, header priming, a second primary source, and graph replay below captured capacity. The original malformed fixture padding assumption was corrected by reading the physical count from serialized stage metadata.

CN-JIT-1 — supported ownership of loaded specialization modules

Source: modules/fused/common/nvrtc_jit.cpp, tests/pipeline/test_jit_cache.cpp

2026-10-12. Loaded CUDA modules/functions are context-bound, but the original cache key held only architecture and generated source. The bounded fix retains the first device's primary context for process lifetime and validates ownership before every lookup, including cache hits. Architecture is read from that driver's current context. Custom contexts and switching devices/contexts are rejected.

Modules are deliberately retained until process exit instead of unloading during static destruction after the driver may have shut down. Context reset/destruction is unsupported: raw-handle equality cannot reliably identify a recreated context. Supporting arbitrary context lifetimes would require a separate lifecycle design.

Concurrent misses can compile duplicate modules; the cache retains one winner and unloads the losers in the same supported context. This establishes correctness, not a guarantee of one compilation under simultaneous first use. Regressions cover cold concurrent lookups, independent specialized pipelines, repeated round-trips, initial/custom-context rejection and reuse after restoring the owning context.

CN-OWN-1 — steady-state owned-buffer reclamation test

Source: tests/pipeline/test_device_buffer.cpp

2026-10-12, H100. Two merged-suite runs failed the existing global free-memory assertion while focused runs passed. An observed 160 MiB shift was not reproduced in subsequent probes, so its cause is unproven. The decoder builds its inverse DAG lazily and retains pool allocations; the test now warms and synchronizes that state before measuring repeated owned output allocation/destruction.

Root Release probes reported 84,431,536,128 free bytes before decode, after a synchronized first decode, and after 128 cycles. A temporary negative control retaining 128 one-MiB outputs reduced free bytes by exactly 128 MiB and failed the unchanged eight-MiB tolerance. The probe and retained-output code were removed. The warmup improves baseline hygiene; it does not establish the cause of the earlier fluctuation or prove that global free-memory accounting is noise-free.

CN-CONCAT-1 — concat gather kernel grid shape

Source: src/pipeline/concat_kernel.cu, src/pipeline/concat_kernel.h Measured: 2026-08, H100 (sm_90), CESM-2D CLDHGH, gpu_zstd_lossless.toml

launch_gather_kernel originally launched <<<n_segs, 256>>> — one block per output segment, regardless of segment size. That collapses whenever segments are unevenly sized, which is precisely the split-mode GPU-Zstd case the kernel exists to serve: the literals port is ~99% of the archive, so one 256-thread block was copying ~23 MB.

before after
gather_kernel 2.27 ms 0.020 ms (112x)
share of pipeline GPU time 60.7% negligible
effective bandwidth ~20 GB/s (on a 3 TB/s part) —
CLDHGH split-mode compress, host wall 3.60 ms 1.33 ms
end-to-end 7.2 GB/s 19.5 GB/s

Output was bit-identical across the change.

Why it hid so well: concatOutputs() runs after the DagEventTimer bracket, so this never appeared in dag_elapsed_ms or in the per-stage table. The only visible symptom was a 2.9x host/device ratio on gpu_zstd_lossless against 1.02x for the single-stream control — and nothing was looking at that.

Why the block count is a compile-time constant: the launch configuration must not depend on input data, or a captured CUDA Graph stops being valid when replayed on a differently-sized input. Blocks with nothing to do exit after one descriptor read, which costs nothing next to the copy they enable.


CN-HF-1 — Huffman coarse-encode partition length (sublen)

Source: modules/coders/huffman/phf/hf_bk.cc, capi_phf_coarse_tune_sublen() Measured: 2026-08, H100 (sm_90), LorenzoQuant -> Huffman, bklen 1024, book_source = Adaptive

sublen (elements per coarse-encode partition, with pardeg = ceil(inlen/sublen)) drives cost in three directions at once:

  • encode phase2/phase4 and the decode kernel launch pardeg-sized grids, so too large a sublen starves the GPU of parallelism;
  • GPU_coarse_encode_phase3_sync D2Hs two pardeg-sized arrays, runs a serial host prefix-sum plus two accumulates over pardeg, and H2Ds the result behind two stream syncs — all O(pardeg), so too small a sublen makes that host barrier dominate;
  • two pardeg-sized arrays ship inside the encoded stream, so too small a sublen also costs compression ratio.

It was previously the constant 768 (original cuSZ v1.x tuning), ignoring inlen despite taking it as a parameter — leaving throughput on the table in both directions at once.

Compress GB/s / decompress GB/s / ratio:

input n sublen=768 (old) current rule change
NYX 134,217,728 216.9 / 139.1 / 29.03 246.3 / 151.2 / 29.67 +14% +9% +2.2%
CESM-ATM 168,480,000 143.7 / 131.9 / 7.97 202.8 / 137.9 / 8.02 +41% +5% +0.6%
HACC 280,953,867 80.8 / 80.8 / 3.62 155.4 / 130.1 / 3.64 +92% +61% +0.6%
HURR 25,000,000 169.8 / 94.8 / 16.13 unchanged (floor) —
CESM 6,480,000 76.3 / 30.4 / 6.19 unchanged (floor) —

The rule targets pardeg ~= 131072 and floors at the historical 768, which makes every change a strict Pareto improvement: compress, decompress and ratio all improve or stay equal on every field measured. The floor is deliberate — below ~64M elements the trade stops being free, because a smaller sublen buys throughput by spending ratio.

That smaller-sublen regime is real and sometimes worth taking, but it is a trade, so it is opt-in via FZ_HF_SUBLEN. At sublen=256 against the 768 default: CESM 26 MB +27% compress / +161% decompress for -3.9% ratio; EXAALT 11.5 MB +30% / +168% for -1.1% ratio.

Do not extrapolate these constants without re-measuring. The curve is not monotonic: parallelism starvation past ~4096 is a steep cliff (CESM-ATM decompress falls 138 -> 92 GB/s going 1024 -> 4096), which is why kMaxSublen sits at the top of the measured range rather than being open-ended.

sublen is recorded in phf_header, so streams stay self-describing: archives written before this change still decode with their own geometry. Not a format change.


CN-GPULZ-1 — match level and split mode

Source: modules/coders/gpulz/gpulz_stage.h Measured: H100, 24.7 MB of Lorenzo-quantized CESM CLDHGH residuals, chunk_size = 2048, word_size = 4

Match level (setMatchLevel) is encode-side only — the stream format is identical either way, so a stream produced at one level decodes the same at the other and the level is not serialized.

level throughput ratio
0 — exact longest match over the 32-element near window only 170 GB/s 4.36x
1 — additionally consults a hashed two-word-key table (offsets to 255) 126 GB/s 5.13x

Split mode (setSplitMode) emits literals / lengths / offsets / meta as four ports instead of one interleaved stream. This is the Zstandard split, for the same reason: the parts have very different symbol distributions, and interleaving them raises the entropy a downstream coder sees. Measured across six SDRB fields, coding the four ports separately beats the single-stream form by 23-43% compression ratio.

The literals port keeps the data's natural word alphabet, so it should feed a symbol-width-matched coder (HuffmanStage<uint16_t> for uint16 quant codes) rather than a byte coder — that alphabet effect is the larger half of the gain.

Every port must be entropy coded and all four re-merged. Unlike the single-stream form, which codes the whole payload by construction, a split leaks any byte left out. Both the raw-fallback chunks and the per-chunk size table are folded into the ports above for exactly that reason: leaving them out cost 28% and 67% respectively in testing.


CN-HF-3 — Adaptive codebook staleness and the refit triggers

Source: modules/coders/huffman/huffman_stage.h, setRefitThreshold() / setRefitInterval() Measured: LorenzoQuant -> Huffman, eb 1e-4 NOA, radius 512, via profiling/huffman_book_drift.cpp

A codebook pinned by HuffmanBookSource::Adaptive goes stale as the symbol distribution moves. Drift in ratio against PerBlock:

corpus mean worst
23 CESM CLOUD levels 4.7% 9.6%
5 different CESM-ATM fields 8.3% 27.6%
20 HURR slabs 9.2% 22.6%

Real but bounded — and well under a fixed model book's 12.6–53.2%.

The sharp edge is a degenerate first block. A constant (zero-range) fitting block poisons the whole run: the same CESM field measures 42.0% mean when the fit lands on a constant level versus 4.7% when it does not, and constant slabs are common (CESM CLOUD levels 0-2, HURR CLOUDf48 slabs 18-19). That is the case that made Fixed look competitive (19.5% vs 42.0%) — a generic shape degrades gracefully where an overfitted one collapses.

Bit-rate trigger (setRefitThreshold, default 1.2) is free: encode() already reports total_nbit, so bits-per-symbol needs no histogram. With the degeneracy guard, stepping 26 CESM CLOUD levels where the first is constant goes from 42.0% mean / 52.5% worst to 0.86% / 3.06%.

Periodic refit (setRefitInterval, default 0 = off) exists because the rate trigger only fires when the rate gets worse. Across 5 different CESM-ATM variables the default is unchanged at 8.3%, because LHFLX's absolute bit rate fell under a mismatched book even while losing 27.6% ratio — a rate trigger structurally cannot see that. --refit-interval 2 brings it to 6.9%.

Deciding that case properly needs the fresh book's rate, i.e. the histogram this mode exists to avoid. When every call carries a genuinely different variable, PerBlock remains the right choice.


CN-HF-4 — why an unencodable book falls back instead of throwing

Source: modules/coders/huffman/huffman_stage.cu Measured: HACC, the two fields that hard-failed

A histogram spanning a wide enough dynamic range drives the rarest symbol past the 27-bit code field, which the builder clamps silently and then emits a stream nothing can decode.

Falling back to an Adaptive book is not a relaxation of the error bound: the bound belongs to the quantizer, and flooring frequencies changes how symbols are spelled, not what they mean. On the two fields that hard-failed:

field eb CR PSNR max abs err
HACC vy 1e-2 15.27x 44.78 dB exactly the bound
HACC xx 1e-3 13.79x 64.77 dB exactly the bound

Throwing instead cost the cusz preset every wide-dynamic-range field in the corpus, which is a worse answer than a book the format can actually hold. buildAdaptiveBook halves the floor shift until the book fits, and shift 0 is uniform, so the fallback always terminates.


CN-HUFFMAN-1 — DeviceResident throughput, latency, and fit-required boundary

Source: commit a2e888d, modules/coders/huffman/huffman_stage.cu, profiling/huffman_throughput.cpp Measured: 2026-08-22, NVIDIA H100 80 GB HBM3 (sm_90), driver 595.84, CUDA 12.9.86/runtime 12.9, GCC 13.3.0, CMake 3.28.3, Release, BUILD_PROFILING=ON, BUILD_TESTING=ON (historical development configuration; research profiling targets are available on development, not stable main).

Method

fzgmod-profile-huffman-throughput invokes HuffmanStage<uint16_t> directly. Input generation, allocations, H2D upload, persistent phf::Buf construction, and correctness copies are outside the timed regions. Allocations issued by execute() itself remain timed. Each steady row uses 10 warmups and 31 samples; the table reports medians.

Terminology: every reported timing is a warmed-up repeated measurement; none of the benchmark tables compares the process's first CUDA invocation against later invocations. “Adaptive cold fit” is the profiler's historical label for a fit-required call: it creates a fresh, preallocated stage for every warmup and sample, so codebook construction is included but allocation/input setup is not. “Adaptive warm reuse” performs one untimed fit, asserts that the book pinned, disables automatic refitting, and measures reuse of that same book. PerBlock is also fully warmed up, but by definition it histograms and rebuilds the codebook on every call; there is no reusable-book hot state for PerBlock. Hereafter “cold” means codebook state, not GPU/process warmup state.

Two independent clocks are reported:

  • host-observed: execute, stream completion, and postStreamSync, including DeviceResident's required PHF header/book-status D2H readback;
  • GPU event: from the event immediately before execute through the last work it enqueues. The deferred metadata readback is intentionally outside this event interval and remains visible in host-observed latency.

Throughput uses uncompressed input bytes and decimal GB/s. Inputs were 4 KiB, 64 KiB, 1 MiB, and 64 MiB, with an extra 2/4/8/16/32 MiB crossover sweep. The distributions were uniform random over 1024 symbols, skewed, eight-symbol low-entropy, and predictive integer codes synthesized from the in-repository data/CLDHGH.f32 field. No external data was downloaded.

Every one of the 128 mode/book/size/distribution configurations decoded exactly. For every deterministic format case, HostCoordinated and DeviceResident archives were byte-identical. All combinations are supported. Only terminal forward DeviceResident+Fixed reports CUDA Graph compatible; PerBlock, both Adaptive states, every HostCoordinated case, and all inverse cases correctly report graph incompatible rather than failing the benchmark.

Host-observed results: representative scientific codes

Compression latency and logical throughput follow. Ratio belongs to the archive and is identical between execution modes.

input book state HostCoordinated ms / GB/s DeviceResident ms / GB/s H/D latency speedup ratio
4 KiB PerBlock 0.194 / 0.02 0.842 / 0.00 0.23x 1.04x
4 KiB Adaptive cold fit 0.288 / 0.01 3.111 / 0.00 0.09x 1.01x
4 KiB Adaptive warm reuse 0.145 / 0.03 0.168 / 0.02 0.86x 1.01x
4 KiB Fixed 0.146 / 0.03 0.170 / 0.02 0.86x 0.95x
1 MiB PerBlock 0.289 / 3.63 1.903 / 0.55 0.15x 2.63x
1 MiB Adaptive cold fit 0.438 / 2.40 3.768 / 0.28 0.12x 2.63x
1 MiB Adaptive warm reuse 0.230 / 4.56 0.252 / 4.16 0.91x 2.63x
1 MiB Fixed 0.247 / 4.25 0.260 / 4.04 0.95x 2.18x
64 MiB PerBlock 0.835 / 80.40 3.019 / 22.23 0.28x 2.68x
64 MiB Adaptive cold fit 0.981 / 68.38 4.220 / 15.90 0.23x 2.68x
64 MiB Adaptive warm reuse 0.699 / 96.07 0.500 / 134.17 1.40x 2.68x
64 MiB Fixed 0.735 / 91.32 0.531 / 126.32 1.38x 2.19x

The large-input warm result is distribution-independent in direction:

64 MiB distribution Adaptive warm, Host / Device GB/s speedup Fixed, Host / Device GB/s speedup
uniform 83.45 / 113.02 1.35x 79.00 / 106.95 1.35x
skewed 103.66 / 146.44 1.41x 93.36 / 130.94 1.40x
low entropy 101.99 / 143.66 1.41x 91.02 / 125.60 1.38x
scientific codes 96.07 / 134.17 1.40x 91.32 / 126.32 1.38x

Decompression does not benefit consistently. At 64 MiB, DeviceResident/host latency speedup ranges from 0.87x to 1.02x across distributions and archives. For the scientific archive it is 50.9 vs 57.9 GB/s (Adaptive, 0.88x) and 63.0 vs 69.5 GB/s (Fixed, 0.91x). Book “source” affects decode only through the bytes in the self-describing archive; no forward codebook is reused by inverse.

Crossover and fit/reuse conclusions

  • Warm Adaptive and Fixed compression: HostCoordinated wins through 2 MiB; DeviceResident wins at the next tested point, 4 MiB, for all four distributions. The 4 MiB speedup is 1.17–1.22x for Adaptive warm and 1.21–1.26x for Fixed. The measured crossover is therefore bracketed at 2–4 MiB on this H100.
  • PerBlock and Adaptive fit-required calls: no crossover through 64 MiB. At 64 MiB, DeviceResident is 0.21–0.88x for PerBlock and 0.23–0.27x for Adaptive cold.
  • Decompression: no robust crossover through 64 MiB. A few uniform or low-entropy points reach 1.00–1.02x, but DeviceResident is usually 5–13% slower.
  • Adaptive amortization: at 64 MiB the DeviceResident cold-fit penalty over HostCoordinated is 2.48–3.24 ms, while each warm reuse saves 0.189–0.210 ms. It therefore needs roughly 13–17 warm calls after each fit (14–18 total calls per book) to amortize the device cold fit. A refit interval shorter than that loses overall on this machine even though every warm row is faster.

Bottleneck attribution

Nsight Systems on the 1 MiB scientific DeviceResident PerBlock case attributes 82.4% of CUDA kernel time to huffmanBuildCanonicalBookDeviceKernel<uint16_t>: 1.584 ms per call, with encode_phase2_deflate next at 0.155 ms. This matches the source launch <<<1,1>>> and the nearly payload-independent 3–4 ms Adaptive cold latency. Adaptive flooring activates all 1024 leaves, which is why its cold builder is slower than a sparse PerBlock tree.

The transient device-book scratch is 40,964 bytes at bklen=1024 and was allocated and freed each build. In the trace, the steady cudaMallocFromPoolAsync API median was 0.025 ms; CUB scan kernels were about 0.003 ms combined. These are measurable but far below the 1.584 ms serial tree kernel, so caching scratch is not the first optimization. Allocation and pool initialization may still matter to true single-shot latency and should remain in a separate first-use benchmark.

For warm DeviceResident compression, host latency exceeds the CUDA-event interval by about 0.05–0.06 ms on average because exact output size and status must be read back. On the 4 KiB scientific case the GPU interval improves from 0.141 ms to 0.105 ms, but host-observed latency regresses from 0.145 ms to 0.168 ms: metadata completion plus launch overhead erases the device-side assembly gain. At 64 MiB the same fixed cost is only about 0.02 ms and encoding kernels dominate, producing the 1.38–1.41x warm win.

Optimization 1: shared-memory builder workspace

The first implementation change keeps the canonical builder's 40,964-byte workspace in dynamic shared memory for the normal bklen <= 1024 case. This removes the per-fit pool allocate/free and, more importantly, moves the serial heap's dependent loads and stores out of global memory. Books above 1024 symbols retain the old pool-backed path; DeviceResidentLargeBookGlobalScratchFallback covers a 2048-symbol byte-identical round trip.

The direct-stage benchmark below repeats the original 10-warmup/31-median method. Values are host-observed compression latency; decompression is unchanged because it does not invoke the builder.

distribution input baseline ms shared-workspace ms latency reduction
uniform 4 KiB 2.967 2.561 13.7%
uniform 1 MiB 3.758 3.257 13.3%
uniform 64 MiB 4.097 3.652 10.9%
skewed 4 KiB 0.687 0.587 14.6%
skewed 1 MiB 3.780 3.274 13.4%
skewed 64 MiB 4.036 3.586 11.1%
low entropy 4 KiB 0.696 0.449 35.5%
low entropy 1 MiB 0.591 0.504 14.7%
low entropy 64 MiB 0.855 0.828 3.2%
scientific codes 4 KiB 0.841 0.722 14.1%
scientific codes 1 MiB 1.902 1.638 13.9%
scientific codes 64 MiB 3.017 2.706 10.3%

The end-to-end cuSZ-style LorenzoQuant<float,uint16> -> Huffman<uint16> PerBlock rerun used 31 samples on the same three local fields as the baseline. Compression ratio, reconstruction error, and Huffman payload bytes remain unchanged. The table compares DeviceResident host-observed compression latency; the HostCoordinated column is included to show the remaining crossover gap.

field bytes HostCoordinated baseline ms DeviceResident before ms DeviceResident after ms D/R reduction remaining H/D speedup
CESM CLDHGH 25,920,000 0.500 2.105 1.885 10.5% 0.27x
Hurricane TC 100,000,000 0.765 2.421 2.184 9.8% 0.35x
NYX temperature 536,870,912 2.545 3.516 3.333 5.2% 0.76x

There is still no PerBlock crossover through 537 MB, although the largest field improves from 1.38x slower than HostCoordinated to 1.31x slower.

This is a useful first reduction, but it does not change the diagnosis. In a post-change Nsight Systems trace of the 1 MiB uniform case, the builder still accounts for 90.3% of kernel time and takes 2.958 ms. The remaining cost is the serial compatibility heap itself, not allocation. A first warp-cooperative heap prototype was rejected because it failed the host/device byte-equivalence test on tied frequencies; deterministic compatibility is not traded for benchmark speed.

Optimization 2: cap DeviceResident decode CTAs

The device-header decode path previously launched 8 * numSMs CTAs for every archive because pardeg was read only inside the kernel. Each CTA copies the reverse book to shared memory before discovering whether it owns a partition, so small and medium fields launched hundreds of needless book copies. The launch is now capped by ceil(estimated_pardeg / 256) while the embedded header remains authoritative. The kernel already uses a grid-stride partition loop, so archives with older or environment-selected partition geometry remain correct even when their embedded pardeg exceeds the estimate. Test DeviceResidentDecodeUsesEmbeddedPartitionGeometry encodes with sublen=64, decodes with the default estimator and its much smaller grid, and round-trips.

cuSZ field DeviceResident decompress before ms after ms reduction HostCoordinated ms after D/H speedup
CESM CLDHGH 0.858 0.865 -0.8% 0.781 0.90x
Hurricane TC 1.117 0.954 14.6% 0.886 0.93x
NYX temperature 4.149 3.605 13.1% 3.606 1.00x

The decoder reaches parity on the 537 MB NYX field and substantially closes the 100 MB gap, but CESM remains launch/serial-decode limited. Direct-stage 64 MiB GPU-event time improves by 4.9–5.7% for uniform/scientific archives and is nearly unchanged for sparse books; the end-to-end gain is larger where the old overlaunch was most exposed.

Decision

The allocation/global-memory part of the builder is now addressed. The next step must replace the serial compatibility heap with a deterministic parallel builder, with byte-equivalence tests that deliberately contain frequency ties. Phase2/4 tuning remains lower priority. Decode CTA overlaunch is addressed; further inverse work would need to parallelize the per-partition serial bit walker rather than adjust launch geometry again.

Handoff state for the next optimization

The working tree is intentionally uncommitted. The current validated implementation contains two optimizations: shared-memory builder workspace for bklen <= 1024 with a global-scratch fallback above it, and a capped DeviceResident decode grid. The dedicated profiler, tests, changelog entry, and this note are part of the same working-tree change. Start by reading AGENTS.md, this entire CN-HUFFMAN-1 entry, and the current diff; do not reset or overwrite these changes.

The next target is the 1-thread compatibility heap inside huffmanBuildCanonicalBookDeviceKernel. A warp-cooperative heap prototype was attempted and removed after tied frequencies produced a different/invalid tree. Do not resurrect it without first proving the exact qinsert/qremove tie order. The existing host/device byte-identity contract is intentional: equal-frequency cases must not merely round-trip, they must produce the same PHF payload.

Acceptance gates for the next builder optimization:

  1. preserve PerBlock, Adaptive-fit, Fixed-fit, overlong fallback, and symbol-range behavior for uint8_t, uint16_t, and uint32_t;
  2. retain byte-for-byte HostCoordinated/DeviceResident output, including deliberate tied-frequency distributions;
  3. retain the bklen > 1024 fallback or replace it with a tested general path;
  4. run test_huffman, test_graph_capture, and Compute Sanitizer;
  5. rerun the direct 10-warmup/31-median profiler before the three-field cuSZ test;
  6. record failed prototypes and before/after evidence here, not in source comments;
  7. do not commit unless explicitly requested.

Reproduction

Run this historical research harness on development. Stable main preserves production timing APIs and offers ENABLE_NVTX for annotations, but does not ship BUILD_PROFILING programs.

CUDACXX=/software/u24/nvhpc/25.7/Linux_x86_64/25.7/cuda/12.9/bin/nvcc \
  cmake --preset cuda-h200 -DBUILD_TESTING=ON -DBUILD_PROFILING=ON \
  -DCUDAToolkit_ROOT=/software/u24/nvhpc/25.7/Linux_x86_64/25.7/cuda/12.9
cmake --build --preset cuda-h200 -j8 \
  --target test_huffman fzgmod-profile-huffman-throughput
./build/cuda-h200/tests/stages/test_huffman --gtest_color=no
./build/cuda-h200/bin/profiling/fzgmod-profile-huffman-throughput \
  --sizes 4K,64K,1M,64M --warmups 10 --repetitions 31 \
  --distribution all --csv /tmp/huffman-throughput-full.csv
./build/cuda-h200/bin/profiling/fzgmod-profile-huffman-throughput \
  --sizes 2M,4M,8M,16M,32M --warmups 10 --repetitions 31 \
  --distribution all --book adaptive_warm --operation compress \
  --csv /tmp/huffman-crossover-warm.csv
./build/cuda-h200/bin/profiling/fzgmod-profile-huffman-throughput \
  --sizes 2M,4M,8M,16M,32M --warmups 10 --repetitions 31 \
  --distribution all --book fixed --operation compress \
  --csv /tmp/huffman-crossover-fixed.csv
/software/u24/nvhpc/25.7/Linux_x86_64/25.7/compilers/bin/nsys profile \
  --trace=cuda,nvtx --stats=true --force-overwrite=true \
  -o /tmp/huffman-device-perblock-1m \
  ./build/cuda-h200/bin/profiling/fzgmod-profile-huffman-throughput \
  --sizes 1M --warmups 0 --repetitions 3 --distribution scientific \
  --mode device --book per_block --operation compress

CN-LRZ-1 — Lorenzo block-mode inverse: CTA width vs reset period

Source: modules/predictors/lorenzo/lorenzo_stage.cu, lorenzo_segmented_scan_kernel

The original block-mode inverse launched blockDim == block_size, tying CTA width to the reset period: a 1024-element segment meant 1024-thread blocks running Hillis-Steele with ~2·log2(1024) barriers. ncu flagged it as barrier-bound, and decompression fell off a cliff exactly there — geometric mean 102 GB/s at block_size = 512 down to 69 at 1024.

The replacement gives each thread Seq consecutive elements: serial scan in registers, then a warp shuffle scan, then one pass over the warp totals — 2 barriers per scan pass regardless of segment length.

field before after
CESM Z3 decompress @ block_size 1024 86 GB/s 120 GB/s
CESM U decompress @ block_size 1024 68 GB/s 107 GB/s

Compression ratio and PSNR bit-identical. The point is not just the speedup: 1024 now matches or beats 256/512 instead of trailing them, which unlocks the highest-ratio configuration. Reset periods that are not a multiple of 32 fall back to the barrier-based scan, which handles any width.


CN-TLRZ-1 — tiled Lorenzo inverse: why per-row, not phased per-tile

Source: modules/predictors/tiled_lorenzo/tiled_lorenzo_stage.cu, tiled_lorenzo_scan_kernel_rows

The earlier phased one-block-per-tile inverse was barrier-bound: ncu showed sm__throughput 51% at 91% resident warps, because it left most of its 64 threads idle across two __syncthreads — even the busiest phase used only ty*tz of tile_elems threads.

The per-row kernel removes barriers and idle lanes entirely: every thread is self-contained and owns one x-row. Exploiting the separable structure, a row thread re-derives its own seed by walking the tile's tiny x=0 "spine" — the z-column prefix then the y-column prefix — and then runs its own tx-length x-chain.

The spine walks are a handful of adds over L1-resident bytes (re-read across the tile's rows, but tile_elems is tiny); the dominant traffic is the coalesced read of each row's contiguous tx deltas plus the inherent tile→natural scatter store. The correctness identity is stated at the kernel itself, since it is what makes the rewrite legitimate rather than merely faster.


CN-QUANT-1 — rounding must happen in the input's precision

Source: modules/quantizers/quantizer/quantizer.cu, quantRound() Measured: 2026-08-07, S3D N2; see the FZGM paper notes

The ABS/NOA kernels rounded with an unconditional __float2int_rn, capping quantization at float32 precision no matter what TInput was. That silently broke every f64 field whose values sit far from zero relative to their own range.

S3D N2 spans 1.1e-5 about 0.7369. The float32 spacing there is 5.96e-08 — 54x the abs_eb that a 1e-4 NOA bound asks for (1.103e-09). The float32 rounding error alone blew the bound, so round-trips came back at 48-53 dB reporting status: ok.

int remains the bin type: codes are at most 32-bit, and widening q would change the wrap semantics that linear mode documents.


CN-QUANT-2 — why outlier-capacity overflow throws instead of logging

Source: modules/quantizers/quantizer/quantizer.cu Measured: 2026-08-07, S3D N2; see the FZGM paper notes

Overflow means excess outliers were dropped, and a dropped outlier reconstructs to whatever the code path leaves behind — for a field whose whole range sits far outside the radius, effectively zero. The caller gets a stream that decodes cleanly into wrong values.

It used to be a DEBUG log line, and the corruption was silent on real data: S3D N2 (range 1.1e-5 on an offset of 0.74) drives every element outside a 32768 radius, overflows a 10% capacity, and round-trips at -96 dB PSNR with status: ok.

Policy matches LorenzoQuantStage: outlier_capacity == 0 is an explicit opt-in to the lossy trade-off and stays a quiet drop; any other capacity that overflows is a failure to honour the error bound and must say so.


CN-TIMER-1 — why DAG timing uses CUDA events, not a host clock

Source: src/pipeline/dag_event_timer.h

DagEventTimer brackets dag->execute() (or a cudaGraphLaunch) with a start/stop CUDA event pair, so elapsedMs() is the device wall time of the pipeline: from when the stream reaches the start marker until every kernel — including those joined back from internal streams at the end of execute() — has finished.

This replaced a host steady_clock measurement around execute(). That was measuring the wrong thing: execute() returns to the host after merely enqueuing kernels (it ends in an async cudaStreamWaitEvent barrier), so a host timer there captured launch latency, not GPU compute.

The event bracket deliberately excludes host setup and PCIe transfers issued outside it. Note the corollary that bit CN-CONCAT-1: work performed outside the bracket — such as concatOutputs() — is invisible to dag_elapsed_ms, so a regression there shows up only in the host/device wall-time ratio.

CN-AB-TR — transpose-based encode: the rate-dependent crossover and threshold 4

Source: modules/coders/adaptive_bitpack/adaptive_bitpack_kernels.cu (encode_pack_kernel_warp_tr, encode_pack_outlier_kernel_warp_tr, encodeTransposeThreshold)

The warp-cooperative pack kernel is ALU-bound, not memory-bound (ncu on CLDHGH, H100: pack at ~72% ALU pipe, ~20% DRAM, store efficiency already ~1.1 sectors/request — the old scalar path's 88%-excessive-sectors store problem is gone in the warp kernels). The ALU cost is the O(rate) per-plane __ballot_sync

  • bit-extract loop. encode_pack_*_warp_tr replaces that loop with a single 32×32 warpBitTranspose32 per 32-element half — the exact inverse of the transpose decode (decode_unpack_*_warp_tr), and since the transpose is an involution, feeding lane l its magnitude yields plane word l at lane l. The archive format is unchanged; the sign region still uses one ballot per half.

The win is entirely rate-dependent, because the gather is O(rate) and the transpose is O(1). Per-kernel gpu__time_duration.sum (ncu, CLDHGH 3600×1800, cuszp2/outlier, isolated pack kernel):

eb gather transpose (thr 4) speedup
1e-3 51.8µs 52.3µs ~1.00×
1e-4 64.3µs 57.4µs 1.12×
1e-5 78.3µs 58.5µs 1.34×

The gather climbs with tightening bounds (higher residual magnitudes → higher rate); the transpose stays flat at ~58µs. End-to-end compress DAG throughput (non-graph, matched clock): 1e-3 unchanged, 1e-4 +2.7%, 1e-5 +7.3% — diluted because the pack kernel is ~¼ of the compress DAG (Lorenzo + quant + rate + scan

  • pack + concat). Decompress is untouched.

Threshold 4, not the decode side's 6. Blocks below the threshold fall back to the gather (cheaper at trivial rate). The encode crossover measured lower than decode's: at eb=1e-4 threshold 4 (57.4µs) beat threshold 6 (60.6µs) and beat transpose-all (58.4µs) — a rate-4 block's four ballots already cost more than the fixed 5-shuffle transpose, but forcing the transpose onto rate-1..3 blocks loses. FZ_ENCODE_TR overrides it (=0 disables entirely for A/B).

Store coalescing. The transpose has each of lanes 0..r-1 write one 4-byte plane word at base + word_bytes*(1+lane) (consecutive words → coalesced). A naïve byte-wise store desequences that across lanes (byte 0 of each lane is strided by 4) and raised store sectors ~50%. storePlaneWord does a single aligned uint32 store when the block's payload base is 4-byte aligned — a warp-uniform property, since every plane slot is base + a multiple of 4 — and falls back to the byte-wise store32le otherwise (outlier ob_bytes of 1–3 can misalign the base). The alignment branch is warp-uniform, so it never diverges.

CN-AB-FUSE — why single-pass fused encode (decoupled look-back) was rejected

Prototyped and reverted (2026-08). The 3-kernel outlier encode path (encode_rate → CUB ExclusiveSum → encode_pack) reads the quantization codes twice and recomputes absU in both kernels. The obvious fix is to fuse all three into one kernel that reads the codes once (keeping them in registers between the rate reduce and the pack) and computes the per-block byte offsets inline with a Merrill/Garland decoupled look-back at CTA granularity.

It was implemented, validated bit-exact against the 3-kernel path (identical compressed size and PSNR, 27 stage + 6 cuszp tests, compute-sanitizer memcheck + synccheck clean at 25k CTAs), and it was 2.4x slower — a clean loss. Two real concurrency bugs surfaced en route and are worth recording: (1) the tile id must be blockIdx.x, not a dynamic atomic ticket — decode recomputes offsets by scanning costs in natural block order, so a ticket-scrambled payload layout is undecodable; (2) aggregate and inclusive prefixes must share a single 64-bit {ready,value} descriptor read/written atomically — separate flag/value words tear at scale and over-count the offset (→ OOB write).

Why it loses (ncu, CLDHGH, H100): the fused kernel runs at 97% occupancy but 17% SM throughput, with 82% of stall cycles at the CTA __syncthreads barrier. The pack cannot start until it knows its byte offset, so all 256 threads block on the look-back while only thread 0 does a serial, global-memory-latency-bound scan of predecessor descriptors. Embedding the prefix scan into the heavy pack kernel entangles the scan's serial dependency with the expensive work. The slowdown is a constant ~2.4x across 32/64/128 MB (not superlinear, so not O(n²) degeneration — just a fixed structural cost). CUB's standalone scan avoids all of this: it is a lightweight kernel where the prefix wavefront propagates over trivial per-tile work in ~5 µs, then the pack runs separately at full occupancy. The ~16% saving from eliminating the redundant read is dwarfed by the barrier/look-back overhead.

Conclusion: the 3-kernel design is the right architecture for this pipeline. A warp-cooperative look-back (32 lanes reducing predecessor descriptors in parallel) would cut the barrier stall but its best case is ~parity for a lot of added lock-free-concurrency risk — not pursued. Related: the L2-blocked tiling probe (profiling/l2tile_profile.cpp) is the other rejected traffic/locality lever; this pipeline is compute/latency-bound (see CN-AB-TR), not traffic-bound.

CN-FUSE-PROOF — full block-local fusion reaches native-class throughput (go decision)

Go/no-go experiment (profiling/fuse_proof.cu), 2026-08. Counterpart to the CN-AB-FUSE negative result. That one fused only rate+scan+pack and lost 2.4x to an in-kernel look-back barrier. This one fuses the entire block-local chain — Quantizer(linear) + Lorenzo(32) + AdaptiveBitpack(32,outlier) — into per-warp work: one warp owns one 32-element block and computes quant + Lorenzo delta + adaptive-bitpack in registers, so the int32 codes are never written to or read from DRAM.

Measured (H100, compress, best-of-10), fused vs the staged FZGM pipeline:

field staged fused speedup
CLDHGH 256 MB 114 GB/s 358 GB/s 3.14x
HACC/xx 512 MB 113 GB/s 280 GB/s 2.49x
NYX 512 MB 114 GB/s 225 GB/s 1.97x

That is native-cuSZp2-class throughput (~200-400 GB/s on H100). Correctness: byte-identical archive to the staged path on all three fields (exact compressed-size match to the byte over 134M elements — any bug in quant/Lorenzo/rate/pack would shift a block cost and diverge the size), and round-trips within the error bound where the linear quantizer itself is in range (CLDHGH; NYX/HACC exceed int32 code range at eb=1e-3 in both paths, a pre-existing linear-quantizer property, not a fusion bug).

Why it worked where CN-AB-FUSE failed. Two design choices: (1) fuse the whole chain so the front-end quant+Lorenzo compute is hidden — ncu shows the fused rate kernel at 74% SM / 74% ALU (healthy), 360 us, vs the staged rate-only kernel's 385 us: quant+Lorenzo folded in for ~free while deleting their two separate kernels (~482 us) and the code round-trips. (2) Keep the lightweight CUB offset scan and recompute quant+Lorenzo in the pack kernel rather than materialise deltas or run an in-kernel look-back — this sidesteps the barrier that sank CN-AB-FUSE, at the cost of reading the input floats twice (still a big net win). A true single-kernel variant with an amortised look-back is the remaining headroom toward the top of the native range.

Speedup scales inversely with rate: low-rate/high-CR data (CLDHGH) benefits most (cheap pack, so front-end elimination dominates); high-rate data (NYX) least (the pack's inherent bit-packing ALU dominates). Floor observed ~2x.

Decision: GO on the auto-fusion framework. The runtime-modularity tax is recoverable and large. This hand-fused kernel is the reference implementation and the first "block-local predict+quant+fixedcoder" driver skeleton for the planner; the staged path is its byte-exact validation oracle. Aligns with the FZ group's domain-specific-compiler-for-lossy-compression direction (WSE case study) — same codegen-of-fused-kernels idea, retargeted to CUDA. Next: the FusionSpec stage contract + a DAG fusion planner (see CN-FUSE-PLAN once landed).

CN-FUSE-PLAN — the fusion planner: what it identifies and why

Source: include/stage/fusion.h, src/pipeline/fusion_planner.cpp

Step 1 of the auto-fusion framework (motivation: CN-FUSE-PROOF). Pure analysis — it changes no execution; it finds the stage chains a fused kernel could collapse, so a later pass can substitute a fused (registered or NVRTC-generated) implementation and leave everything else staged.

Stage::getFusionSpec() classifies a stage's data-access pattern — Map (element-wise), BlockLocal (resettable fixed-size neighbourhood), Cooperative (a per-block fixed-length coder), or Unfusable (default; opaque or genuinely global). Access pattern, not computation, is what decides fusability, so the enum is deliberately tiny. Config-dependent fusability is reported honestly: QuantizerStage is Map only in linear mode (the outlier/in-place paths scatter to side buffers), LorenzoStage is BlockLocal only in 1-D block-reset mode (the N-D stencil needs a different driver), AdaptiveBitpackStage is Cooperative only at the warp block sizes 32/64. Inverse stages are always Unfusable.

planFusionGroups() returns maximal groups under four rules, each load-bearing: (1) strictly linear — a group edge needs the producer to feed exactly one stage and the consumer to be fed by exactly one; fan-in/out inside a group would break register-resident composition. (2) one block size — block-local and cooperative members must agree (Map imposes none); mixed periods can't share a tile. (3) coder terminates — nothing fuses past a Cooperative variable- length coder, so it's always the tail. (4) size ≥ 2 — a lone stage is not a fusion. The cuszp2 pipeline resolves to exactly one group {Quantizer, Lorenzo, AdaptiveBitpack}, block 32, coder-terminated — the chain the CN-FUSE-PROOF kernel hand-fuses. Tests: tests/pipeline/test_fusion_planner.cpp.

Not yet built (next steps): a fused-implementation registry keyed by the stage-chain fingerprint (scripts/gen_stage_fingerprints.py gives the key), DAG substitution of a FusedStage for a matched group, and NVRTC codegen that assembles a driver skeleton from the members' device functors for chains with no hand-written implementation.

CN-FUSE-EXEC — wiring fused execution into the DAG (Step 2)

Source: include/pipeline/fusion_registry.h, src/pipeline/fusion_registry.cpp, modules/fused/fused_block/, the fusion hook in CompressionDAG::execute(), Pipeline::planAndInstallFusion(), FusionPolicy.

Turns the planner (CN-FUSE-PLAN) into measured speedup. Measured end-to-end through the real Pipeline on CLDHGH (H100): staged 107→234 GB/s at eb=1e-3 (2.18x), 105→220 at 1e-4 (2.10x), identical PSNR — native-class.

Integration choice: swap execution, not DAG structure. buildHeader() and buildInverseDAG() read the DAG nodes, so keeping all three stage nodes intact (only replacing forward execution) leaves the archive byte-identical and decompress completely unchanged. The registry (findFusedImpl) matches a group's exact shape; CompressionDAG::execute() runs the matched FusedImpl at the group's head node — writing the tail node's output buffer — and records every member's completion event so downstream waits are satisfied, skipping the members' individual execute()s. FusionPolicy::Auto (default Off, FZ_FUSION=off|auto override) installs groups at finalize; unmatched groups stay staged, so Auto is always safe.

Two non-obvious correctness requirements, both from decompress reusing the forward stage objects (setInverse(true); the in-memory decompress path does not call deserializeHeader, so it relies on forward-execute side effects):

  1. Prime forward-computed state the inverse needs. The quantizer computes computed_abs_eb_ during forward execute(); fusion skips that, so the reused inverse quant reconstructed with the default bound (1e-4) — a silent 10x too small reconstruction with a byte-identical archive. The fused runner must call QuantizerStage::primeAbsEbForFusion() (ABS only). The tail coder likewise needs num_elements_ via setFusedResult(). General rule: a fused runner must establish whatever forward state the group's reused stages read on the inverse path.
  2. Disable buffer coloring. Coloring aliases buffers by staged liveness, where a group's input dies after its first stage. A fused kernel keeps that input live across the whole group while writing the group's output, so an aliased input/output region corrupts. planAndInstallFusion() disables coloring (and CUDA graph mode — the fused runner synchronises to read the archive length) when any group installs.

Validated bit-exact (byte-identical archive to staged, round-trip within bound, rs == rf) and compute-sanitizer clean: tests/pipeline/test_fusion_planner.cpp. Next: NVRTC codegen keyed by the stage-chain fingerprint (see CN-FUSE-DRIVER for the parametric-driver step that generalised the runner to cuSZp3).

CN-FUSE-DRIVER — parametric block-local fused driver (cuSZp2 + cuSZp3) (Step 3)

Source: modules/fused/fused_block/ (the templated kernels + two predictor policies), matchesCuszp3/runCuszp3 in src/pipeline/fusion_registry.cpp, TiledLorenzoStage::getFusionSpec().

Generalises the single hand-written cuSZp2 runner (CN-FUSE-EXEC) into one driver for the whole predict+quant+fixed-rate-outlier-coder family. Key realisation: the fused rate/pack kernels ARE the AdaptiveBitpack outlier warp kernels (encode_{rate,pack}_outlier_kernel_warp, CN-AB-TR area), with the materialised int-codes array in[start+idx] replaced by an inline pred.delta(lane,b,m) call. So the two kernels are templated on <int ElemsPerLane, class Pred> (block_size = 32*EPL, EPL≤2 keeps the per-lane register arrays from spilling) and a predictor policy struct passed by value:

  • Lorenzo1DPredictor (EPL=1, cuSZp2): 1-D Lorenzo via intra-warp __shfl_up.
  • TiledLorenzo2DPredictor (EPL=2, cuSZp3): 2-D separable tiled Lorenzo (tz==1), each element re-quantises its own left/up predecessor from the float field — a pure map, so the neighbour code equals what the staged quantizer produced. Because the pack path is byte-for-byte the shipped AB outlier warp encode, and EPL=1 reduces to the old cuSZp2 kernel exactly, both archives stay byte-identical (the cuSZp2 end-to-end test still passes unchanged).

cuSZp3 element bookkeeping (the non-obvious part): the AB stage's element count is the padded tile-major count num_tiles*tile_elems (TiledLorenzo emits zero-padded edge tiles), so every block is full (count==64) and padding elements are ab_active with delta 0 — the driver passes n_ab (not the natural field size) and setFusedResult(n_ab, …). The predictor, by contrast, keys off the field bounds (gx<dx && gy<dy) and returns 0 for padding, matching tiled_lorenzo_delta_kernel. Matcher requires 2-D (tz==1) tile_elems==64 and AB block 64 + outlier; other shapes return 0 and fall back to staged.

Measured end-to-end on CLDHGH 3600x1800 (H100, eb=1e-3), compress steady-state: staged 178.9 GB/s dag / 151.5 host → fused 222.7 dag / 196.0 host (~1.25x dag, 1.29x host), ratio 8.31x and PSNR 63.7997 identical. Smaller multiple than cuSZp2's 2.1x because the cuSZp3 staged baseline (efficient per-row TiledLorenzo) is already faster; the win is still removing two int32-width DRAM round-trips. Validated byte-identical + round-trip on non-tile-aligned dims (300x180, exercises edge-tile padding) and compute-sanitizer clean: tests/pipeline/test_fusion_planner.cpp (Cuszp3*). fzgpu is a different family (BitplaneRZE coder, radius/zigzag quant with its own outliers) — deferred; the natural next step is NVRTC codegen keyed by the stage-chain fingerprint.

CN-BSHUF-SMEM — shared-memory staging fixes the bitshuffle uncoalesced-store bug

Source: modules/shufflers/bitshuffle/bitshuffle_stage.cu — bitshuffle{Encode,Decode}Kernel32Smem, butterfly32, FZ_BITSHUF_SMEM.

The 32-bit butterfly bitshuffle (bitshuffleEncodeKernel32) writes each bit-plane word to global at stride npp = N_chunk/32 (512 B for the default 16 KB chunk): out[i/32 + sublane*npp]. The 32 lanes of a warp share i/32 and differ only in sublane, so a warp's stores hit 32 addresses 512 B apart — fully uncoalesced.

Measured (PFPL on CLDHGH 3600x1800, eb=1e-4 NOA, H100, ncu):

scattered smem-staged
store sectors/request 32 (max) 4 (coalesced)
store byte efficiency 12.5% —
Compute (SM) 10% —
DRAM throughput 9% 31%
L2 pipeline 79% (bound) —
kernel time (steady) 122 µs 33 µs (3.7x)

The kernel was neither compute-bound (10% SM) nor DRAM-bound (9%) — it was L2-transaction bound: the scattered stores generated 8x the necessary sectors. Because it is one CTA per 16 KB chunk, the fix stages the permuted chunk in shared memory and flushes it coalesced: coalesced global read → register butterfly → conflict-free smem scatter → __syncthreads() → coalesced global flush. The smem plane layout is padded to npp + 1 words per plane; a warp writes s[sublane*pstride + col] with col warp-uniform, so stride pstride = npp+1 (coprime to 32, since npp is a multiple of 32) puts all 32 lanes in distinct banks. Decode is the mirror (coalesced load into padded smem → butterfly → coalesced store); s[sublane*pstride + i/32] is exactly the value the scattered decode gathered, so the butterfly (self-inverse) is unchanged and the output byte-identical.

Whole PFPL DAG 0.302 → 0.214 ms (1.41x) from this one kernel; PSNR unchanged (84.7564 dB), 16 bitshuffle stage tests + 35 stage tests pass, compute-sanitizer clean. Default on; FZ_BITSHUF_SMEM=0 forces the scattered kernel for A/B. Dynamic smem is opt-in raised past 48 KB via cudaFuncAttributeMaxDynamicSharedMemorySize, with a fallback to the scattered kernel if a chunk's 32*(npp+1)*4 bytes exceed the device ceiling. Only the 4-byte (primary) path is staged; 1/2/8-byte are unchanged.

This is the "optimize the bottleneck before fusing" result from the PFPL roofline study: bitshuffle was 38% of PFPL runtime and looked compute-bound from the stage timing, but was actually a coalescing bug. It also validates the chunk-cooperative (CTA-per-chunk, smem-resident) kernel shape a flexible fusion framework would target.

CN-RZE-HELPERS — folding away the RZE strip/add-offset helper kernels

Source: modules/coders/rze/rze_stage.cu — RzeStripFlagOp, rzePackKernel.

The RZE forward path launched six device operations serially on one stream: rzeEncodeKernel → sizes memcpy → rzeStripFlagKernel (& 0x7FFFFFFF) → CUB ExclusiveSum → rzeAddOffsetKernel (+ header_size) → rzePackKernel. Per-kernel duration (PFPL on CLDHGH, 1583 chunks, H100, ncu):

kernel µs note
rzeEncodeKernel 42.5 compute-bound (72% SM) — the RZE compaction, hard to shrink
rzeStripFlagKernel 3.1 7 blocks — trivial, pure overhead
DeviceScanInit + Scan 6.6 CUB offset scan (1583 elems)
rzeAddOffsetKernel 2.9 7 blocks — trivial, pure overhead
rzePackKernel 12.7 copies the 4.5 MB packed payload — real memory work
Σ GPU 67.8 vs a measured 92 µs stage

The strip and add-offset kernels are folded away: the flag strip becomes a thrust::transform_iterator feeding ExclusiveSum (so the scan reads flag-stripped sizes with no separate kernel and no d_clean_dev_ array), and the header offset is added inside rzePackKernel (dst = header_off + offset[cid], sz = sizes[cid] & 0x7FFFFFFF). The tail-size readback reads d_sizes_dev_ and masks/adds header_size on the host. Six device ops → four; one device scratch array removed; byte-identical archive (22 RZE + 35 stage + 16 pipeline tests pass, compute-sanitizer clean, PFPL PSNR unchanged 84.7564 dB).

Honest payoff: small. RZE stage 92 → 87 µs, PFPL DAG 0.214 → 0.212 ms (~1%). The naive "stage 92 − Σ GPU 68 = 24 µs of launch latency" estimate was optimistic: kernels queue asynchronously on the stream, so most of that gap is event-timing granularity and encode-tail overlap, not exposed launch stalls — removing two tiny kernels only recovers their ~6 µs of GPU time. The change stands as a simplification (fewer kernels/arrays, cleaner offset handling) with a marginal perf bonus, not a bottleneck fix. The RZE stage remains dominated by the compute-bound encode (42 µs) and the memory-bound pack (13 µs); shrinking it further means a faster compaction algorithm or fusing encode+pack, not trimming helpers. Contrast CN-BSHUF-SMEM, where the bottleneck really was addressable.

CN-PFPL-FUSE — chunk-cooperative fused PFPL proof (the ~2.2× ceiling holds)

Source: profiling/pfpl_fuse_proof.cu (fzgmod-profile-pfplfuse). Go/no-go for the second fusion strategy — chunk-cooperative (one CTA per chunk, intermediates in shared memory), the LC-style counterpart to the warp-register driver (CN-FUSE-*).

PFPL = Quantizer(NOA,zigzag) → Difference(1-D chunk, negabinary) → Bitshuffle(ew4) → RZE(w1). The three post-quant stages all chunk at 16384 B = 4096 int32 = one CTA, so one kernel does quant → diff+negabinary → bitshuffle → RZE-encode per chunk with every intermediate in smem; the cross-chunk scan + pack tail is kept (RZE offsets depend on all chunks). Device functions are reused verbatim (Zigzag::encode, Negabinary::encode, butterfly32 copied, and lc_detail::d_RZE), and the quantizer's ebx2_r is matched via a forced value_base, so the fused archive is byte-identical to a manually-run staged 4-stage forward.

Measured (CLDHGH, n truncated to 1582 full chunks, eb=1e-4 NOA, H100):

staged (4-stage) fused
archive / ratio 4.458 MB / 5.54x 4.458 MB / 5.54x (byte-identical)
compress wall 0.208 ms 0.103 ms (2.01x)
compress DAG 0.172 ms — (fused path 0.103 wall)
fused encode kernel — 70.5 µs, 47% SM, 11% DRAM
achieved occupancy — 92.6% warps active
smem / block — 36 KB (= RZE-encode alone)

The occupancy worry was unfounded. The fused kernel reuses RZE's own two 16 KB buffers (codes→bitshuffled in one, negabinary→RZE-s_out in the other) plus its 4 KB temp, so its footprint equals rzeEncodeKernel alone — 4 blocks/SM, 92.6% achieved occupancy, no collapse. The roofline prediction (~2.2×, from CN reassessment after the bitshuffle fix) holds up at ~2× measured; the fused encode is now compute-bound (47% SM — the RZE compaction + quant, the irreducible floor), DRAM down to 11% (the three intermediate round-trips are gone). compute-sanitizer clean.

Why this matters for the framework: it is the first working chunk-cooperative fusion, and it confirms the general recipe — a fusable chunk-aligned chain predict/transform* → variable-length coder collapses into [CTA-per-chunk fused kernel] + [scan + pack tail], reusing the coder's smem budget for free. Prototype limits (not blockers): outlier-free quant fast path (verified outlier_count==0 for this data/config via the byte compare), whole-chunk only (n truncated to a multiple of 4096; the partial-tail chunk needs the staged path's guard), and the fused-path timing still pays a per-call CUB-temp cudaMalloc a persistent buffer would remove.

CN-CHUNK-FUSE — composable device-op harness for chunk-cooperative fusion

Source: modules/fused/chunk_fusion/chunk_fusion.cuh (harness + ops), chunk_fusion.cu (launcher), proof profiling/pfpl_fuse_proof.cu.

Generalizes the hand-coded PFPL fused kernel (CN-PFPL-FUSE) into a composable mechanism — the first step toward automatic, pipeline-agnostic fusion (goal: a future NVRTC path that generates the glue for any compatible chain).

The taxonomy it encodes. Fusion has two orthogonal axes:

  • Strategy (execution template): warp-register (cuSZp — warp=block, registers, shuffles, no barriers; fused_block/) vs chunk-cooperative (LC/PFPL — CTA=chunk, shared memory, __syncthreads; this file). Determines the harness.
  • Role (access pattern = FusionAccess): Map (quant), stencil/block-local (diff, lorenzo), fixed cooperative (bitshuffle), variable-length coder / sink (RZE, RRE, AdaptiveBitpack), unfusable. Determines how a stage's device-op slots in.

The mechanism. Each stage contributes a small __device__ OP — QuantInplaceZigzag (Map: global floats → codes), DiffNegabinary / Bitshuffle32 (transforms: smem→smem, ping-ponged), RZECoder / RRECoder (the swappable sink, uniform LC signature d_XXX(csize, in[CS], out[CS], temp[CS])). The harness chunk_fused_kernel <QuantOp, Coder, Transforms...> is stage-agnostic glue: load chunk → Chain<Transforms...> ping-pongs the ops through two 16 KB smem buffers with the right syncs → coder sink → emit per-chunk bytes+size. The cross-chunk scan+pack tail (RZEStage-format archive) lives in the launcher. Swapping the coder or transform set re-composes with no new glue — that is the whole point, and what an NVRTC path would later emit as source (the ops stay hand-written headers; only the harness is generated).

Generalization proven (CLDHGH, 1582 chunks, eb=1e-4 NOA, H100): the same harness, given RZECoder vs RRECoder, fuses PFPL-with-RZE and PFPL-with-RRE, each byte-identical to its staged pipeline (RZE 4.458 MB @ 243 GB/s; RRE 4.695 MB @ 239 GB/s — different ratio, correctly a different coder), compute-sanitizer clean. Requires inplace-outlier quant (single "codes" port) so the chain is strictly linear; the op encodes out-of-radius values as raw float bits inline. Partial-tail chunk is handled (Bitshuffle32 copies the sub-chunk through, matching the staged bitshuffle's tail memcpy). Not yet wired to the planner/registry — that (chunk-cooperative group detection + mapping a matched chain to the harness instantiation) is the next step, then NVRTC codegen.

CN-CHUNK-WIRE — wiring chunk-cooperative fusion into the planner/registry

Source: getFusionSpec() on Quantizer(inplace)/Difference/Bitshuffle/RZE/RRE, matchesPfpl/runPfpl in src/pipeline/fusion_registry.cpp, QuantizerStage::primeComputedAbsEb, RZEStage::setFusedResult, and RREStage::setFusedResult.

Turns the chunk-fusion harness (CN-CHUNK-FUSE) into automatic fusion via FusionPolicy::Auto — PFPL now fuses inside a real Pipeline, byte-identical and round-trip-correct, for both RZE and RRE (tests/pipeline/test_fusion_planner.cpp Pfpl*, 10/10 fusion tests, sanitizer clean, 16/16 pipeline + 35/35 stages).

Reuses the existing four-layer machinery (FusionSpec → planFusionGroups → registry → execution substitution; see CN-FUSE-EXEC). New pieces specific to the chunk-cooperative strategy:

  • FusionSpec roles, block_size in BYTES. Quantizer(inplace+zigzag, ABS/NOA) is Map; Difference(negabinary chunk) and Bitshuffle(ew4) are BlockLocal{16384}; RZE/RRE are Cooperative{16384} (the sink — bitshuffle must NOT report Cooperative or the planner terminates before the coder). The planner unifies block_size at 16384 (bytes), a different unit from cuSZp's element counts (32/64) — fine because the registry dynamic_cast disambiguates the strategy.
  • NOA priming. The fused runner bypasses Quantizer::execute(), so it calls primeComputedAbsEb() (runs the value-range scan, sets computed_abs_eb_/ computed_value_base_) — the fused kernel's ebx2_r and the reused inverse/header both need it. (ABS-only primeAbsEbForFusion was the cuSZp version.)
  • Coder axis in the runner. runPfpl detects RZE vs RRE by dynamic_cast and passes the matching ChunkCoderKind — the same registry entry fuses either.

The decompress-sizing gotcha (cost an hour; general to any fused variable-length coder). RZEStage::estimateOutputSizes() in the inverse falls back to the compressed input size when cached_orig_bytes_ == 0, and the fused runner had skipped the forward execute() that normally sets it — so the inverse RZE output buffer was sized to the compressed length (~393 KB) and the decode wrote the full 4 MB one byte past, a heap overflow only compute-sanitizer caught (the sticky CUDA error surfaced two kernels later as a bogus "DifferenceStage launch: invalid argument"). Fix: setFusedResult(archive_bytes, orig_bytes) also sets cached_orig_bytes_. General rule (extends CN-FUSE-EXEC's "prime forward state the inverse needs"): a fused coder runner must set every field the inverse stage reads to size its output, not just the forward output size. Next: NVRTC codegen to emit the harness for arbitrary compatible chains (done — CN-NVRTC-FUSE).

CN-NVRTC-FUSE — runtime codegen of the chunk-cooperative harness (NVRTC)

Source: modules/fused/chunk_fusion/nvrtc_chunk_fusion.{h,cpp}, chunk_geometry.h, the chunk_fused_body split in chunk_fusion.cuh, the FZ_FUSION_NVRTC toggle + shared packChunks tail in chunk_fusion.cu.

Replaces the compile-time template instantiation of the fused kernel (CN-CHUNK-FUSE picks chunk_fused_kernel<Quant, Coder, Transforms...> in C++) with a runtime one: a ChunkFusionSpec (the stage-chain fingerprint — quant op, transform ops, coder op, all by device-op type name) is turned into CUDA source that wraps chunk_fused_body<...>, compiled with NVRTC, JIT-loaded via the CUDA driver API, cached by (source, arch), and launched with cuLaunchKernel. Swapping the coder or adding/removing a transform is now a data change to the spec — no new template instantiation, no registry glue. The per-stage __device__ ops stay hand-written in chunk_fusion.cuh; only the composing glue is generated. Selected by env FZ_FUSION_NVRTC=1 inside launchFusedChunkPfpl() (off by default; the template path is the fallback). Both paths share packChunks (the cross-chunk scan/pack tail), so only the encode-kernel launch differs.

Validated: byte-identical archives to the template path for RZE and RRE (test_fusion_planner Pfpl* under FZ_FUSION_NVRTC=1, 11/11), compute-sanitizer clean on the JIT'd kernel, throughput parity (242 vs 248 GB/s RZE, 231 vs 240 GB/s RRE — within run-to-run noise), 52/52 ctest unregressed on the default path. The codegen contract (spec → template-arg list) is locked host-only by NvrtcCodegenComposesSpecOps.

What it took to make the op stack NVRTC-compilable (the real work):

  • backend/api.h no-ops under __CUDACC_RTC__ — NVRTC compiles device code only and provides the intrinsics itself; pulling in the host <cuda_runtime.h> would fail. atomics.h/warp.h then reduce to plain atomicAdd_block/ __shfl_*_sync/__ballot_sync, all NVRTC-native.
  • Split chunk_geometry.h out of chunk_fusion.h. The device path must not pull backend/types.h (it names cudaStream_t/cudaMemPool_t, host-only). The geometry constants (CHUNK_BYTES/NELEM/TPB…) now live in a dependency-free header shared by the device harness, the host launcher, and the NVRTC codegen.
  • Minimal in-memory stdlib stubs. NVRTC ships no <cstdint>/<type_traits>/ <cstddef> on an include path; the op headers use only a tiny surface (is_integral/is_signed/is_unsigned/make_unsigned + the int aliases + size_t), supplied as three nvrtcCreateProgram headers. assert is an NVRTC builtin, so <cassert> is an empty shim.
  • __builtin_memcpy → __float_as_uint in the quant op's outlier bit-cast — the former is not an NVRTC builtin; the latter is native to both nvcc and NVRTC and is bit-identical.
  • __builtin_clz{,ll} → guarded FZ_CLZ32/64 in d_PRencode (RARE/RAZE): under __CUDACC_RTC__ use the __clz/__clzll device intrinsics NVRTC provides, else keep __builtin_* (which nvcc accepts in host+device passes). Bit-identical for the non-zero inputs used, so a fused RARE/RAZE chain matches the staged kernel. Same class of fix — surfaced when the RARE/RAZE coders became fusable (Phase D).

Params ABI (Phase B, generalized). The generated kernel takes one packed const uint8_t* params blob — each op declares a POD Params (stateless ops use EmptyParams → 0 bytes), chunk_fused_body hands each op its slice at a compile-time offset (OpParamBytes/SumParamBytes), and the host packs the blob. So the codegen no longer bakes in a per-pipeline param shape; any op set with any params composes. The one canonical chunk_fused_body still backs both the template and NVRTC paths (no drift). Interim: launchFusedChunkPfpl hardcodes PFPL's quant Params when packing; the generic runner (Phase C) assembles the blob from each stage's getFusedOp().params.

Historical prototype constraints. The original implementation resolved op headers through source-tree -I paths. As of 2026-10-12, the build embeds the warp, TI and chunk harnesses plus their quoted project-header closure into the library; nvrtcCreateProgram receives them in memory without source include paths. scripts/gen_nvrtc_headers.py reads the original headers on every build and rewrites its generated header only when content changes, so there is no second maintained copy of the device code. The installed-package regression removes its isolated source/build trees and relocates the installation before compiling kernels.

At the prototype stage, the registry matchesPfpl still hardcoded the 4-stage shape — the codegen dissolves the kernel per-shape glue, not yet the matcher (Phase C: route by declared FusionStrategy, assemble the chain from getFusedOp()).

CN-QUANT-3 — why linear-mode bin overflow is refused, not wrapped

2026-08-23. Linear (cuSZp-style, no-outlier) quantization stores round(value / (2*abs_eb)) directly as a two's-complement TCode, with no radius clamp and no outlier path. Tight NOA bounds on offset-valued fields need bins far wider than int32: S3D N2 at rel_range=1e-6 needs roughly 3.3e10. The old kernel rounded straight to int and the result wrapped silently, so the archive decoded to values far outside the bound with no error.

The forward kernel now rounds into a long long intermediate (quantRoundWide), rejects non-finite or out-of-guard scaled values, and checks the result against the signed TCode range. Any overflow sets a device flag (d_linear_overflow_scratch_, one persistent uint32_t) and writes code 0; postStreamSync() reads the flag and throws with remediation guidance (wider codes, a value-centering transform, an outlier-capable quantizer, or a looser bound).

Why throw instead of logging or clamping: a wrapped or clamped bin cannot be reconstructed within the bound, so any archive produced after an overflow is silently wrong. Refusing the run is the only behaviour that keeps the error-bound contract. The per-element check costs one compare on an already-computed value; the host read-back is one 4-byte copy at the existing stream sync.

Regression test: QuantizerLinear.RefusesCodeOverflow (f64 input, offset 0.7369, range 1.1e-5, NOA 1e-6, uint32_t codes; must throw).

CN-FSZTI-1 — thread-independent AdaptiveLorenzo forward kernels (negative result)

2026-09-09, H100. Two attempts to make FusedQuantAdaptiveLorenzoStage's forward kernel thread-independent (TI), mirroring native FSZ's execution shape. Both were byte-identical to the shipped CTA-cooperative kernel and both were 1.5-1.7x slower, so neither was merged. The branches (ti_predictor_wp1 at f2e9919, ti_predictor_arrays at c1a8b26) were deleted on 2026-10-09; this entry is the record.

Motivation. Native FSZ's compress kernel uses one-warp CTAs (FSZ_TBLOCK_SZ=32) in which each thread serially owns 4 whole 256-element tiles in local registers, with no shuffles and no shared memory beyond one per-CTA offset counter. FZGM's kernel is CTA-cooperative: 256 threads per tile, warp shuffles, __syncthreads(), 9 shared arrays, and a single-threaded final 4-variant cost selection. ncu showed ours compute-bound (77% SM, 96% occupancy), so native wins by doing less total work. The remaining FZGM-vs-native FSZ gap was ~3.7x after M1 (FusedQuantAdaptiveLorenzo) shipped.

Common scaffolding. A new kernel sat alongside the CTA kernel and wrote the same dense per-tile buffers (modes_dense, means_dense, flags, residuals), so the downstream CUB scan, compaction and AdaptiveBitpack needed no change and the archive format was untouched. Dispatch was gated by FZ_AL_TI=1 (default off) with FZ_AL_TI_TPT for tiles per thread, parsed once in a magic-static env config. A non-serialized Config::TIDispatch {Auto, ForceCTA, ForceTI} and ti_tiles_per_thread let one gtest binary compare both kernels without racing the env cache. Tests (MatchesCTAByteIdentical, MatchesCTAAcrossTilesPerThread, TIForwardRoundTrip) swept blocks_per_tile in {1,2,4,8} x order-2 on/off x centering on/off x TPT {2,4,8}. blocks_per_tile = 16 is excluded because the CTA kernel's __launch_bounds__(256,8) already forbids 512-thread tiles. The kernels were derived from this file's own CTA kernel; no FSZ kernel body was copied (only the thread mapping idea).

WP1: scalar-state, three-pass TI kernel (fused_quant_adaptive_lorenzo_forward_kernel_ti<T, TILES_PER_THREAD>). Per-block CoderStats (all/rest/first for LZ1 and LZ2) were kept as scalars, costed with blockCost() and reset at each 32-element block, so register use did not scale with blocks_per_tile. The cost was three passes per tile: (1) a full-tile scan for per-block stats and the mean sum, (2) a 32-element reread of block 0 to cost the centered variants once the mean was known, (3) a full-tile pass to emit the chosen residuals.

NYX/temperature, ncu --set full CTA TI TPT=2 TI TPT=4 TI TPT=8
Duration 1.84 ms 3.67 ms 3.45 ms 3.43 ms
Compute (SM) 77.40% 9.86% 10.26% 10.26%
Memory throughput 64.45% 63.57% 66.07% 64.85%
DRAM throughput 17.23% 17.69% 18.03% 17.59%
Registers/thread 32 48 52 53
Achieved occupancy 95.99% 52.45% 46.30% 24.36%

End to end (fzgmod-cli -b --report-json, device_ms.min, 20 reps x 3 invocations, fsz_fused_quant.toml): NYX/temperature CTA 2.906-2.908 ms vs TI 4.585-4.604 ms (1.58x slower); HACC/vx 6.545-6.550 ms vs 9.925-9.994 ms (1.52x slower). Compressed sizes were identical (NYX 24,421,280 B, CR 21.98).

WP1b: register-array TI kernel (fused_quant_adaptive_lorenzo_forward_kernel_ti_arrays<T, TPT, BPT>, one-warp CTAs, BPT = blocks_per_tile in {1,2,4,8} as a compile-time argument). Both predictors' CoderStats were held in CoderStats[BPT] arrays (block loop fully unrolled so they scalarize; tile loop #pragma unroll 1 so the arrays are reused, not multiplied by TPT). WP1's block-0 reread was removed algebraically, because centering changes only residual 0 (LZ1) or residuals 0 and 1 (LZ2):

  • quantized value q0, mean mu, c0 = q0 - mu, cm0 = |c0|;
  • centered LZ1 block 0 = {rest1 | cm0, rest1, cm0};
  • keep d1[1] and the OR of LZ2 magnitudes over elements 2..31 (tail2); with cm1 = |d1[1] - c0|, centered LZ2 block 0 = {tail2 | cm1 | cm0, tail2 | cm1, cm0};
  • centered variant cost = uncentered total - block-0 cost + centered block-0 cost + sizeof(T) for the mean; the minimum of the four variants picks the mode.

Prediction and mode analysis then read each input once; a second traversal emits the selected residual stream (the modular stage must materialize it for the separate AdaptiveBitpack launch).

NYX, TPT=8, ncu WP1 scalar TI WP1b array TI CTA
Kernel duration 3.43 ms 3.24 ms 1.84 ms
Compute (SM) 10.26% 8.66% 77.40%
Memory throughput 64.85% 68.50% 64.45%
DRAM throughput 17.59% 17.85% 17.23%
Registers/thread 53 76 32
Local load/store n/a 0 / 0 n/a
Achieved occupancy 24.36% 23.83% 95.99%

End to end: NYX 2.904-2.910 ms (CTA) vs 4.420-4.440 ms (1.53x slower); HACC 6.546-6.549 ms vs 10.705 ms (1.64x slower). NYX TPT sweep: 4.670 ms (2), 4.975-5.020 ms (4), 4.420-4.440 ms (8). cuobjdump confirmed 76 registers and zero local bytes, which rules out spilling. Two variants were also tried and reverted, both byte-identical: float4 loads (NYX 4.548 ms, 78 registers) and a one-raw-input-pass hybrid that stored LZ1 residuals during analysis and transformed them in place for LZ2 (NYX 4.778 ms).

Root cause. Both TI kernels are latency-bound on uncoalesced per-thread global access: each thread walks its own contiguous tile, so a warp touches 32 different tiles per instruction. Loads issued 8,388,608 requests for 268,435,456 sectors, i.e. 32 sectors per request with 4 useful bytes per 32-byte sector. Native FSZ has the same uncoalesced pattern (19 sectors/request measured) but tolerates it, because its single-pass, 67-register kernel does enough arithmetic per load to hide the latency and fuses prediction, selection and final coding in one kernel. FZGM's modular boundary forces a dense residual stream between AdaptiveLorenzo and AdaptiveBitpack, which adds traffic without adding compute. Removing WP1's block-0 reread gained only ~5.5% of kernel time.

What a future attempt would need. More TPT or register tuning will not help. Either (a) a coalescing-aware thread-to-data mapping that is no longer "thread owns a contiguous tile", or (b) predictor+coder fusion that consumes thread-local predictor state directly instead of materializing residuals. Either must keep byte identity with the CTA kernel and pass the same ncu checks (sectors/request, local memory, SM%). The same lesson appears in AdaptiveBitpack's scalar-coder experiment (FZ_AB_FORCE_SCALAR, removed): matching native's thread independence without its coalescing and predictor-coder locality loses.

CN-WARPPROBE-1 — tiled plain-chain adaptive-probe threshold

2026-09-09, H100, large-data corpus at rel_range 1e-3. The TI-vs-warp-cooperative probe threshold adaptive_thresh = 16.0 was a 2-point fit on 1-D cuSZp2 HACC data (xx r14: TI wins 393 vs 236 GB/s; vx r17.6: warp-cooperative wins). The probe always estimates rate with a flattened serial Lorenzo1D kernel (ti_rate_probe_kernel), so for tiled 2-D/3-D plain (non-outlier) chains it sent CESMATM-3D/T and SCALE-LETKF/T to TI at about half the single-pass throughput (25-26% of native cuSZp3 instead of 48-50%).

adaptive_thresh_tiled = 1.4 was calibrated from FZ_DEBUG_PROBE=1 avg_r on the 10 tiled corpus fields. The data are not separable by one global cut: CESMATM-3D/CLOUD wants single-pass at avg_r = 1.749, below EXAFEL/data, which wants TI at avg_r = 3.641. 1.4 routes 9/10 fields correctly, costs EXAFEL/data ~5% (331.9 -> 315.6 GB/s), and fixes CESMATM-3D/{T,U,CLOUD} and SCALE-LETKF/{T,U,QV} (+15-96%). Corpus-wide (24 cells, cuSZp2 + cuSZp3 plain/outlier vs native): compress throughput range 25.0-84.2% -> 45.7-80.5% of native, geomean 0.558 -> 0.601; output byte-identical (dispatch only). Resolving the EXAFEL/CLOUD ambiguity needs a dimension-aware or tile-shape-aware probe, not a constant. Full analysis: paper_organizer projects/FZGM/investigations/fusion/fused_execution_paths_map.md.

CN-WARPTILE-1 — compile-time tile shape in the tiled warp predictors

2026-09-08, H100, cuSZp3 2-D/3-D warp-cooperative path. Native-vs-FZGM instruction profiling showed the fused path executing 10.8x more instructions per element than native (14.53 vs 1.35), ALU-bound (65.9% ALU pipe) where native is memory-bound (68% DRAM). lx = local % tx and ly = local / tx are per-element and were runtime divisions against a uint32_t field. Making the tile shape a template argument (both shipped shapes are powers of two: 8x8 2-D, 4x4x4 3-D) turns them into mask/shift: 1.95B -> 1.73B executed instructions (-11%), NYX/temperature 177.7 -> 191.6 GB/s (+7.8%), byte-identical. Corpus-wide (95 cells vs native, geomean): cuszp3_plain compress unchanged at 0.45x, decompress 0.28x -> 0.39x; cuszp3_outlier 0.64x/0.46x -> 0.66x/0.56x compress/decompress. A per-tile tile_ctx() hoist of the tile-index divisions was measured separately and had no effect (<0.1% instructions; the compiler already CSEs it), so it was not kept. These numbers predate f64 specialization (2f5fe00); the port onto the Real-templated predictors (2026-10-09) was re-verified for correctness, not re-timed.

CN-F64SPEC-1 — f64 and predictor-free fixed-mode diagnostics

2026-10-05, H100. These are local performance diagnostics, not publication rows. ABS bounds were precomputed from the full-field range at relative bound 1e-3. Each CLI invocation ran 20 in-process iterations; the diagnostic reducer discarded the first three, then compared device-time medians. The discarded inverse iterations include cold NVRTC compilation. This differs from Benchkit's five retained repetitions plus one warmup with independent phase executions.

Field / pipeline Auto compression speedup over Staged Auto decompression speedup over Staged
brown-p2-plain 1.382x 1.736x
brown-p2-outlier 1.737x 2.098x
miranda-p2-plain 1.369x 1.674x
miranda-p2-outlier 1.727x 2.044x
miranda-p3-plain 0.800x 1.025x
miranda-p3-outlier 0.884x 1.147x
miranda-p3-fixed 0.921x 1.504x

The new fixed 1-D route has the following warmed diagnostic results:

Field / pipeline Compression Decompression Selected forward schedule
brown-fixed-f64 0.606x 1.384x single_pass
exaal-fixed-f32 0.511x 1.096x two_pass

Both fixed 1-D compression cases regress. The registry admits legal implementations and does not choose by measured profitability. This extension therefore closes an implementation gap without establishing a compression performance benefit; include these regressions if reporting its coverage.

The cuSZp2 results show the intended recovery of staged boundary costs. Tiled cuSZp3 compression is slower on this MIRANDA field; removing boundaries does not guarantee a win. Extra double quantization work, the separate input validation scan/synchronization, and the selected warp schedule are plausible costs, but this diagnostic does not establish their individual contributions. Float64 thread-independent dispatch remains unavailable. No native-relative claim follows from these Auto/Staged ratios. Some cuSZp2 inverse sample CVs are 6–8%; run the publication timing gate before including any measurements.

Raw samples and the reduction script are in build/f64-auto/performance-check/. Fixed 1-D follow-up artifacts are in build/f64-auto/fixed-validation/ and fixed-performance-check/.

CN-F64SPEC-2 — targeted publication rerun scope

The historical Auto source has 330 successful, bound-satisfying, timing-reliable f64 cuSZp rows: MIRANDA 7, S3D 11, NWCHEM 1 and BROWN 3 fields, times three bounds and five variants. Add 36 f32 fixed 1-D cells for six HACC and six EXAALT fields. Each cell reports both compression and decompression.

Existing native/FZGM complete-case compression coverage is 278 coordinates: p2 plain/outlier 66 each, p3 fixed 33, plain 62, outlier 51. Native S3D p3 fixed severe bound failures exclude all 33 coordinates; native plain/outlier validity removes another 4/15. Extending FZGM does not repair native failures. Keep coverage and valid populations distinct, and do not cherry-pick faster new Auto rows. Twelve of the 278 previous fallbacks were fixed 1-D structural mismatches, rather than dtype mismatches.

Historical precomputed-range H100 cell cadence is 6.6–7.2 seconds, including process, artifact, and quality work. Estimate 36–40 minutes raw for 330 Auto cells, budget 45–75 minutes practical; the extra 36 f32 cells add about 5–10 minutes. Matched Off+Auto gives 732 cells including f32 fixed: budget roughly 1.5–2.5 hours, subject to current per-process JIT and timing retries. This is an estimate, not a reserved runtime or completed campaign.

Recommend an additive campaign pinned to a clean new source revision, starting with a 12–24-cell calibration overlap. Preserve dataset hashes, effective bounds, coder geometry, and original manifests. Reuse native only where protocol/provenance and quality/timing gates permit; refresh Off for the small affected matrix to establish current-source comparisons. Apply Benchkit verify and all validity/timing/complete-case gates before regenerating paper tables. Do not combine old NOA timings with new precomputed-range ABS timings without explicit protocol reconciliation. Historical staged/native populations remain useful but are not automatically interchangeable with this new build.

Correctness gate: 60/60 CTest executables passed. Fixed-mode memcheck and racecheck passed under ordinary and forced single-pass (FZ_SP_BPW=16) dispatch; the f64 suite also passed memcheck after the shared validation change. Full BROWN f64 and EXAALT f32 checks have identical codec payloads, normalized codec metadata, and reconstructions between Staged and Auto. Full-file hashes differ in capacity metadata, as with the original f64 checks. EXAALT has the same 12 ppm maximum-bound overshoot in both arms; this satisfies Benchkit's existing 0.1% relative acceptance tolerance and is recorded explicitly in the artifact. It is not a new error introduced by specialization.

Publication protocol source: compression_benchmarking/configs/experiments/ range_precomputed_spec.yaml; FZGM sessions are rq2-range-precomputed-merged-h100-{off,auto}-20260930b, as pinned by specialization_native_performance_full_corpus.json. Native cuSZp2/3 evidence comes from Sep. 9 with dimension-matched fixed-mode native additions on Sep. 28; those bounds were resolved outside timing and reconciled in the merged session.

CN-INIT-1 — Defined bytes in materialized outlier capacities

During forward execution, LorenzoQuant reports the maximum capacities of its outlier-error and outlier-index ports. The actual count becomes available at postStreamSync(), after downstream stages have run. A Merge stage can therefore copy the full capacities into its payload while only the outlier prefixes have been written. Compute Sanitizer initcheck exposed this in MergeStage.PipelineReversal during the October 2026 release preparation.

The producer now clears both outlier capacities on its execution stream before launching its forward kernels. The kernels overwrite the valid prefixes; their unused tails remain zero. Empty capacities and absent pointers skip the clears. This preserves stream ordering, graph capture, segment geometry, and serialized counts. The previously undefined tail bytes may change a merged archive's bytes or downstream compressed size; decoder values still depend only on the valid prefixes recorded by the outlier count.

The additional writes per execution total max_outliers * (sizeof(TInput) + sizeof(uint32_t)) bytes. This is a bounded initialization fix, not a measured performance improvement. Propagating actual interior stream lengths before downstream execution would avoid copying and coding these tails, but requires a separate runtime and graph-capture contract.

The staged Quantizer forward path has the same delayed outlier count and now applies the same initialization to its value and index ports. Inplace and linear modes have no separate outlier arrays and skip these writes. Fused split-outlier execution already reports the populated side-output sizes before its downstream consumers run; this change targets staged execution only. Quantizer sanitizer validation is tracked separately from the observed LorenzoQuant failure.

CN-HFINIT-1 — Huffman assembly copies populated prefixes

Compute Sanitizer initcheck exposed a device-to-device read of uninitialized bitstream capacity in CLI.BenchmarkGraphModeFallsBackOnIncompatibleStage during the October 2026 release preparation. The compatibility Huffman encoder concatenates only its populated codeword cells, but phf::Buf::memcpy_merge() copied bitstream_max_len cells into its archive scratch buffer. That upper bound is an allocation capacity, not the encoded length. In the observed call, the copy spanned 15,552 bytes and the first undefined source bytes were at offset 512. Memcheck passed because the source remained inside its allocation.

Assembly now copies header.total_ncell bitstream cells and header.pardeg entries from each partition metadata array. Counts are checked against their allocated capacities before copying. Header and reverse-codebook copying retain the existing offsets and lengths. This fixes a specific over-copy rather than zeroing the entire scratch allocation or changing runtime stream propagation. Archive round-trips, identity, and sanitizer results must be verified separately; no throughput improvement is claimed.

CN-SZXINIT-1 — Initialize packed bytes before modifying bits

SZx writes residual codes one bit at a time using byte-wide read-modify-write operations. A fresh payload allocation has no defined initial byte contents. Compute Sanitizer initcheck exposed these reads in SZxStage.SmoothRoundTrip_F32 and other SZx tests during the October 2026 release preparation. Reconstructing valid values can still succeed because every valid code bit is eventually overwritten; unused high bits in the final packed byte remained undefined.

The bit writer now assigns zero when entering a new byte, before setting or clearing its first code bit. The bit position is continuous across elements of a block, and each block owns a separate byte-aligned payload slice. This defines both the bytes read by the bit writer and the unused padding bits without clearing an allocation-sized upper bound or adding another kernel launch. Archive metadata, code widths, offsets, and decoder behavior are unchanged by this source change. Existing archive bytes can differ in previously undefined padding bits; byte identity and performance are not claimed. Sanitizer results are tracked separately from this diagnosis.

CN-ANSINIT-1 — Defined ANS coalesced padding

October 2026 release initcheck identified three undefined-byte paths in the vendored ANS coalescer: reserved header fields (and upper option bits preserved by setProbBits()), the extra metadata slot needed for an odd block count, and the last 16-byte payload vector of a block whose encoded word count is not a multiple of eight. The encoder writes only the valid 16-bit words in that final vector. Copying the whole vector read unwritten scratch bytes and propagated them into archive copies and checksums.

The coalescer now value-initializes its header, zeroes unused metadata slots, and assembles a partial final vector from guarded valid-word reads plus zero padding. Full vectors retain their vectorized copy. Header size, block offsets, compressed sizes, and valid encoded words are unchanged; previously undefined reserved and padding bytes become zero. The decoder consumes the existing valid word counts, so these changes do not alter reconstruction or require a new archive format.

The fix writes only the padding already present in the coalesced representation. Clearing all uncoalesced scratch instead would add writes proportional to the full encoder scratch capacity. This bounded fix avoids that extra bandwidth; its runtime impact has not been benchmarked. The ANS round-trip, partial-block, reuse, and file tests exercise the relevant layouts under Compute Sanitizer.

CN-GIRACE-1 — Mode-2 profiling orders reset before sample loading

Compute Sanitizer racecheck reported shared-memory write/write warnings in GInterp mode-2 profiling during the October 2026 release preparation: the reset helper wrote zero to sample centers while the following loader wrote actual sample values. The affected cases were GI36_AutoTuneMode2_3D, GI37_AutoTuneMode2_2D, and GI38_AutoTuneMode2_FileRoundTrip. All 62 functional GInterp cases passed in that run; racecheck reported three warning groups and zero error groups. This does not establish that the phase overlap is harmless.

An unconditional block barrier now separates shared-memory reset from sample loading. All threads enter this sequence, and the reset helper has no early return; threads that have no reset items still reach the barrier. Existing geometry, neighbor writers, and the loader-end barrier are preserved. This change targets the mode-2 probe rather than the interpolation kernels or other tuning modes. The profiling launcher uses one CUDA block.

The change must be checked with racecheck and synccheck along with the existing mode-2 correctness and file-round-trip tests. Reconstruction failures or a performance improvement were not demonstrated by the diagnostic run.

CN-FILEIO-1 — File decode result has one owner

buildInverseDAG() marks each inverse sink result persistent, so MINIMAL reset() releases transient intermediates but keeps the result until DAG destruction. The static Pipeline::decompressFromFile() helper previously freed that result directly while its DAG still recorded ownership. The DAG destructor then retried the release. MemoryPool::free() rejected the untracked pointer and emitted a warning instead of issuing a second CUDA free, explaining why memcheck could pass while the CLI still warned.

The helper now leaves the result under DAG ownership. Its separately allocated caller output is guarded by OwnedDeviceBuffer until all copying and performance bookkeeping succeeds, then ownership is published. The compressed input blob is external to the inverse DAG; the helper frees it and clears its pointer so later exception cleanup does not retry that release.

The CLI/file regression failed before the change with the exact warning on all three decodes (CLI with performance output, then two static calls without it). It checks reconstruction and warnings through function return, including DAG destruction. This catches the ownership error that device sanitizers alone miss.

CN-LIBPRESSIO-1 — libpressio adapter configuration and archive ownership {#cn-libpressio-1}

Source: the separate libpressio repository's src/plugins/compressors/fzgpumodules/ and tools/swig/libpressio.py. These notes describe the integration using FZGPUModules v2.0.1 and the FZGPUModulesChain builder, rather than the archived v2.0.0 helper interface. They preserve the adapter implementation notes outside the libpressio source tree.

Builder and constructor phases

Encoder configuration starts with empty stage and connection lists. Encoding requires an explicit stage list or a TOML config file; creating a bare plugin and decoding an archive do not require either. The default absolute and relative error bounds remain 1e-3.

add() stores stage tokens, ordered connections, and per-stage parameters. configure(**options) stores pipeline-wide options with the fzgpumodules: prefix and returns the builder. It performs no native calls or validation. Repeated keys overwrite earlier builder values. Subsequent changes to a builder affect compressors constructed afterward; they do not update existing objects.

early_config() returns stages and connections, all stored top-level options, and all stored stage parameters in one dictionary. It does not construct a compressor, finalize a pipeline, or allocate GPU buffers. The dictionary and stage/connection lists are new objects; arbitrary option values are not deep-copied.

The generic libpressio Python constructor first applies early_config with inferred option types, queries the resulting compressor's available options, and uses that template to convert and apply compressor_config. Selecting stages first exposes their option names and types to the second pass. The FZGM builder includes its stage parameters in the first pass too: its return value is not limited to topology. compressor(eb=..., metric=...) supplies this dictionary automatically, puts pressio:abs and pressio:metric in the second pass, and returns a PressioCompressor. The native encoder is constructed on first encode.

Execution and ownership

Use one GPU with the standard CUDA runtime primary context. Do not switch the current device, use custom CUDA contexts, or reset the device while FZGM objects are alive. pressio:thread_safe=serialized requires external serialization of calls across all FZGM instances, including clones; the plugin adds no lock.

A clone copies options and stage parameters, then builds its own encoder and CUDA stream on first execution. It does not share a live pipeline, buffers, or execution metrics. Stream creation failures are reported through libpressio's error API.

NumPy input is uploaded for compression, and the resulting archive is returned to host memory for ordinary Python save/load. GPU input keeps the archive on the GPU; a GPU output request also keeps it there. Decoding returns host output by default or GPU output when requested. Libpressio's domain manager owns these transfers; the adapter copies core pool-owned buffers into independently owned output buffers.

fzgpumodules:execution_time_us measures wall time around the core compression call, excluding pipeline construction, input upload, archive assembly, and host output transfer. fzgpumodules:peak_memory reports the core pipeline pool's peak usage, excluding external input/archive allocations and host memory. Host rusage metrics measure a different quantity.

Configuration updates

FZGM option updates are transactional: an error restores all top-level values, per-stage parameters (including inactive stage settings), and the pipeline's previous rebuild state. Successful changes invalidate the cached encoder; unchanged options retain it. Type errors and invalid top-level numeric values are reported rather than silently ignored. Dimensions require one to three integer values; only the first may be zero to select automatic sizing.

Stage-specific ranges, stage tokens, connection topology, and kernel constraints are still validated lazily when the pipeline is constructed. A successful set_options() does not guarantee that a new pipeline can execute.

Self-contained output

Each encoded buffer is an FZLP v1 envelope containing the original libpressio array dtype/dimensions, an FZM metadata header, and the compressed payload. A fresh compressor can decode it without a prior encode, chain options, TOML file, or caller-provided output buffer. The installed SDK must support the stages used by the archive. Decoder specialization policy and memory strategy are local execution choices; stage topology and parameters come from the archive.

from pathlib import Path

Path("array.fzlp").write_bytes(compressed.tobytes())
# In another Python process with the same compatible installation:
decoder = lp.PressioCompressor.from_config({"compressor_id": "fzgpumodules"})
restored = decoder.decode(Path("array.fzlp").read_bytes())

A supplied output array must match the archived dtype and shape. NumPy inputs produce a host-resident encoded array that can be saved directly; GPU inputs keep encoded output on the GPU. Decoded host output is returned by default; supplying a GPU output preserves GPU residency. Previously encoded buffers remain independently decodable after later encodes or option changes. Metadata bytes count toward libpressio's reported compressed size and compression ratio. Reconstruction currently creates a fresh pipeline per decode; a decoder cache is a separate performance improvement. Decoder reconstruction currently uses the core default pool multiplier and stream setup, without CUDA Graph capture. Execution time, peak memory, rebuild, outlier-count, and stage-output metrics describe the last encode; specialization metrics describe the active pipeline.

This changes the plugin's output format from raw buffers to versioned envelopes. Legacy raw buffers are rejected because they lack reconstruction metadata. FZLP is a libpressio array envelope, not byte-identical to the core .fzm file format. Neither a Python helper nor temporary files are involved in decoding.

FZLP v1 layout

All integers are little endian. The 40-byte prefix contains the eight bytes FZLP0001, a uint32 version (1), a uint32 dtype code, uint32 rank (1–32), uint32 reserved flags (0), uint64 FZM metadata length, and uint64 payload length. It is followed by rank uint64 dimensions in libpressio dimension order, zero padding to a 16-byte boundary, the FZM metadata bytes, and the payload. Dtype codes 0–9 respectively denote float32, float64, int8, int16, int32, int64, uint8, uint16, uint32, uint64. The dimensions must describe exactly the uncompressed byte count in the FZM header. The FZM header checksum validates its metadata; FZLP v1 does not add a payload checksum. Array dimensions and dtype are validated against byte counts.

Adapter layout

  • fzgpumodules.cc: libpressio options, pipeline lifecycle, execution, and metrics.
  • configuration.h: option values, validation, and stage parameter ownership. Updates are parsed on a candidate and committed together after validation.
  • pipeline_builder.h: manual/TOML graph construction, named-port wiring, output references, and graph capture. Returns a complete encoder before replacing live state.
  • stage_registry.h: SDK-gated handler catalog, token lookup, and parameter snapshots.
  • stages/stage_kind.h: common handler contract and value-copy support.
  • stages/stage_*.h: stage-specific options, documentation, and construction.
  • archive.h: the self-contained array envelope format.

To add a stage family, implement a handler in stages/, include it and register it in stage_registry.h, and add FZGM integration tests. Use the existing SDK capability guards for optional families. Registration order retains the previous first-match behavior; individual handlers own their token aliases. The main plugin requires no stage-specific edits.