Module df64 

Double-single (df64) emulated double-precision GEMM in raw WGSL.

WGSL has no f64 — only f32/f16/i32/u32. The native exact-double GPU path is therefore CUDA-only (crate::wgsl_forge::emit::cuda_c::GEMM_F64_SRC), which covers NVIDIA. This module fills the other half of the “best f64 path on every machine” story: on any wgpu-capable GPU (AMD, Intel, Apple, Adreno, Mali, llvmpipe, …) it emulates each f64 as a hi/lo pair of f32 — a double-single number carrying ~44–48 effective mantissa bits — and does the GEMM accumulation with error-free transforms (Dekker two_prod, Knuth two_sum, quick_two_sum). This is the df64/double-single technique.

Why this is a RAW WGSL string, not the portable IR

The df64 arithmetic (the two_sum/two_prod/df_add/df_mul helpers operating on vec2<f32> hi/lo pairs, with a Veltkamp-split Dekker two_prod) does not map onto the forge’s portable scalar-op IR. So — exactly like the cooperative-matrix kernel (crate::wgsl_forge::oracle::evaluate_coopmat_loadstore) — it is shipped as a hand-written WGSL source string, compiled via crate::wgsl_forge::execute::WgpuPipeline::compile and dispatched directly, rather than emitted from a crate::wgsl_forge::KernelSpec.

Buffer ABI (matches the CUDA-f64 / forge GEMM binding layout)

binding	name	usage	layout
0	a	storage, read	M*K df64 = 2*M*K f32, [hi,lo,hi,lo…]
1	b	storage, read	K*N df64 = 2*K*N f32, [hi,lo,…]
2	c	storage, read_write	M*N df64 = 2*M*N f32, [hi,lo,…]
3	dims	storage, read	[m, n, k] as u32

Row-major C[M×N] = A[M×K] · B[K×N]. One invocation computes one output element (@workgroup_size(64), gid.x over m*n). The inner kk accumulation order matches the f64 CPU reference crate::wgsl_forge::dispatch::gemm_cpu_f64, so the two agree to df64 precision (~1e-12 for K≈64 O(1) data) on adapters whose WGSL float arithmetic is not reassociated by the driver.

Correctness depends on per-op IEEE rounding (probed at runtime)

Every df64 algorithm (Veltkamp split, Dekker two_prod, Knuth two_sum) relies on each f32 +/-/* rounding exactly once, with no algebraic reassociation. Some GPU shader toolchains break that: on the naga→SPIR-V→NVIDIA-Vulkan path on this hardware, the compiler simplifies c - (c - a) to a (and an fma-based residual fma(x,y,-(x*y)) to 0), which collapses the residual (lo) terms — df64 then silently degrades to f32 (~2e-7 error) instead of ~double precision. Switching two_prod from fma to the Veltkamp split gave a byte-identical wrong result, confirming reassociation (not a missing fma) as the cause; WGSL has no portable pragma to disable it. The crate::wgsl_forge::dispatch::gemm_f64 chain therefore probes this kernel at runtime (df64_usable) and uses it only on adapters where it actually delivers f64 precision; elsewhere it falls to native CUDA-f64 or the exact CPU floor. The kernel below is correct on a faithful-IEEE adapter and is kept as the portable f64-GPU path for those.

Constants

GEMM_DF64_ENTRY

Entry-point name of GEMM_DF64_WGSL.

GEMM_DF64_WGSL

Raw WGSL source for the df64 (double-single) emulated-f64 GEMM. See the module docs for the binding ABI and precision contract. The error-free-transform helpers are transcribed verbatim — quick_two_sum/two_sum/two_prod and the df_add/df_mul pair-arithmetic are subtle and must not be “simplified”.