Kernels

RecursiveAI.OneGpu.Kernels holds the compute kernels that machine-learning inference needs, written in HLSL and compiled by the shader build for SPIR-V (Vulkan), Metal Shading Language and WGSL (browser). This page describes the element-wise, reduction, softmax, normalisation and layout families; each is a public class created from a GpuDevice.

Using a family

using var elementwise = new ElementwiseKernels(device);
var encoder = device.BeginCommands();
elementwise.EncodeBinary(encoder, GpuBinaryOp.Add, x, 0, [rows, cols], bias, 0, [cols], y, 0, new GpuActivation(GpuUnaryOp.GeluErf));
device.Submit(encoder);

• A family is created once; pipelines are built the first time a variant is used and kept. Calling Encode… records one compute pass into the caller's encoder. With the same buffers, a call allocates nothing (bind groups are cached by buffer identity; offsets and sizes travel in push constants).
• Operands are a GpuBuffer<T>, an element offset and a shape. Tensors are dense and row-major. T is float or Half for arithmetic; int, uint, long where an operation only moves or converts bits. Comparison results and conditions are bool tensors (see below).
• Every size, shape, axis, permutation and offset is checked in safe code before anything is recorded. A bad argument throws GpuKernelShapeException, which names the operation and the argument or dimension (Reduce: axes: axis 3 is outside [-3, 3)).
• Families are not thread-safe. Disposing a family releases its pipelines once the GPU has finished with them.
• GpuKernelMode.Fast is the default. GpuKernelMode.Deterministic is accepted by every family; for the element-wise and layout families there is no reduction order, so both modes give the same bits.

Families

| GatherKernels | ONNX Gather (any axis, negative indices, int32 or int64 indices; the embedding lookup), GatherElements, GatherND (with batch dimensions), over 8-, 16-, 32- and 64-bit payloads |

| ScatterKernels | ScatterND and ScatterElements with reduction none, add, mul, min or max |

| LayoutKernels | Slice (including negative steps), Expand, Concat (up to four inputs per pass, more by chaining; inputs may live in the output's buffer when their ranges do not overlap the output's, read in place, or through a scratch copy in a browser); Expand may run in place (below) |

| Int8Gemm | MatMulInteger (uint8/int8 in any mix, scalar or per-row/column zero points, int32 result), the float-scaled form and the requantized form of QLinearMatMul |

| ConvIntegerKernels | ConvInteger (uint8/int8 input and weights in any mix, scalar x zero point, scalar or per-channel w zero point, groups, dilation, padding, rank 3 and 4, int32 result) |

| DynamicQuantizeKernels | DynamicQuantizeLinear (float32 to uint8 with scale and zero point) |

| DequantizeKernels | DequantizeLinear: uint8, int8 or int32 to float32 or Half, per tensor or per axis, optional zero point; QuantizeLinear: float32 or Half to uint8 or int8, per tensor or per axis |

| RecurrentKernels | ONNX RNN, GRU and LSTM over float32 or Half, forward, reverse and bidirectional (see recurrent.md) |

Adding an operation to the element-wise family is a case in Shaders/Functions.hlsl and an enum member: the operation is a push constant, not a pipeline.

Integer and bool tensors

• int32 is one 32-bit word per element. int64 is two words (low, then high) handled with 32-bit arithmetic, so a device needs no 64-bit integer support; every value from long.MinValue to long.MaxValue is exact, and so is division (a shift-and-subtract loop runs only when an operand does not fit 32 bits).
• bool is a GpuBuffer<byte>: one byte per element, 0 or 1, four to a 32-bit word. It is the same layout as a CPU bool[], so a model's bool tensor uploads as it is. Allocate GpuKernelBuffers.PaddedBoolLength(n) elements; a buffer that ends inside a word throws GpuKernelShapeException. A thread owns whole output words, so a kernel that writes a bool range at an odd offset keeps the neighbouring bytes of the words it shares, and no atomics are needed. Inputs treat any nonzero byte as true; outputs are 0 or 1.
• Arithmetic wraps (two's complement). Div truncates toward zero, is 0 for a zero divisor, and the minimum value divided by -1 is itself. Mod takes the divisor's sign, Fmod the dividend's; both give 0 for a divisor of 0 or -1. These are onemodel's CPU semantics, which the tests check bit for bit, as they do every operation of the family: integer, bool and cast results have no tolerance.
• Range takes its start and delta as host values (a graph's Range inputs are known when a model's shapes are planned) and its length from IntegerKernels.RangeLength, which follows the ONNX formula in exact integer arithmetic. The float range is start + i * delta as two separately rounded float32 operations.
• Cast to bool is "not equal to zero" (NaN counts, negative zero does not); from bool it gives 0 or 1.
• One variant serves Fast and Deterministic: integer arithmetic has no order to fix and no atomics are used.
• An int64 Where moves two words per element and so needs one more dimension than the other types: its shapes may have one fewer dimension after merging.

Binding caches

A family builds a bind group the first time it sees a set of buffers and reuses it on every later call with the same set, so a warm call allocates nothing. Every family constructor takes maxCachedBindGroups (default GpuKernelBindings.DefaultMaxCachedBindGroups, 1024). When a family has seen more sets than that it drops all of them and starts again, so the bound should exceed the number of distinct buffer sets one run records (a smaller bound is correct, only slower: a dropped group is destroyed after the command buffers that refer to it are submitted, in the browser too): a 72-weight transformer touches 72 sets in Gemm per pass, and with a bound of 64 every pass rebuilt them all (72, then 144 bind groups over two passes) where a bound of 72 or more builds 72 once. ReleaseBindings() on every family (and Gemm) drops the cached groups immediately; call it after disposing buffers that will not be used again, or after replacing a model's buffers, so stale entries do not keep them reachable. The pipelines stay, and the next call with a set of buffers builds its group again. BindGroupsCreated on every family (read-only, never reset) counts the groups built so far, so a consumer can record a run twice and check that the second run adds none; for RecurrentKernels it includes the input projection's Gemm.

In place

Slice and Concat refuse an output that is an input buffer. Expand accepts one buffer as both input and output (the way a memory plan places two tensors in one arena) when the element ranges are disjoint: the input's offset .. offset + count and the output's offset .. offset + total must not overlap, and for Half neither may share a 32-bit word with the other. Overlapping ranges throw GpuKernelShapeException (operation Expand, argument output) before anything is recorded. On Vulkan and Metal the same buffer is bound twice in the pass (read-only and read-write); the ranges being disjoint, no thread reads what another writes. A browser refuses one buffer that is both in one dispatch, so there the call expands into a scratch buffer the family keeps (grown on demand) and copies it to the output: two passes instead of one.

Half precision

Half tensors are stored as pairs in 32-bit words and converted in the shader with integer arithmetic, so they need no 16-bit storage feature and run on every backend, including a browser without shader-f16. Arithmetic is in single precision and the result is rounded to nearest even once; the rounding is written out in the shader because the drivers' own float-to-half conversion does not promise a rounding mode.

• A buffer the kernels touch must end on a word boundary: allocate GpuKernelBuffers.PaddedHalfLength(n) halves. Otherwise GpuKernelShapeException.
• An output may overlap an input only with the same offset and extent.
• Softmax and normalisation write whole rows as words, so half tensors need the softmax axis last, an even row length and even element offsets.
• Index arithmetic avoids integer division where it can (same-shape and bias-style broadcasts, a single run of reduced dimensions): on Apple GPUs a chain of divisions per element made binary add and reductions two to three times slower, so those cases are compiled as separate variants.
• Reduction outputs are written with two integer atomic AND/OR operations on disjoint bits, so neighbouring outputs of different workgroups never clash. There are no floating-point atomics anywhere.

Tolerances

The error measure is max|gpu − reference| / max|reference| over an output, with the reference computed on the CPU in double precision (inputs first rounded to half for half rows). The limits below are what the tests assert; the measured values are far inside them (see the parity table the tests print).

The tanh form of GELU computes its exponential with its own range reduction and polynomial (relative error below 1e-7) rather than the driver's exp, whose accuracy is implementation-defined. Erf is the Abramowitz–Stegun 7.1.26 approximation (absolute error below 1.5e-7); the Gelu erf form uses its complementary function so the negative tail keeps its relative accuracy. A float sum of n values has an error of about ε·√n·(Σ|x|/|Σx|) in the worst conditioning: a sum that cancels to almost nothing cannot be accurate in relative terms in any float implementation. The tests use mostly positive data for sums.

Float32 accumulation order

Float32 products and convolutions accumulate k in blocks with block sums added in increasing order (see gemm.md), so their error against a double-precision sum grows far more slowly than one chain over k would. Measured scaled error at K = 3072: 1.5e-7 to 3.5e-7 on every device (a single chain gave 1.8e-6 to 4.5e-6). The tests that check it are AccumulationAccuracyTests.

Fast and deterministic

The reductions, softmax and normalisation have both modes.

• Fast: each output (a reduction element, a softmax row, a norm row) is computed by one 256-thread workgroup. Every thread folds a strided slice in index order (reductions add short runs of 32 first, which keeps the rounding error near that of a pairwise sum), and the 256 partial results are combined with the device's subgroup arithmetic when GpuCapabilities.Compute reports basic, vote and arithmetic subgroup operations and no subgroup is smaller than 4 lanes; otherwise with the shared-memory tree. Variant on each family says which. Constructing a family with allowSubgroups: false forces the shared-memory path.
• Deterministic: the same slices, then a fixed halving tree in shared memory (thread i combines with i + s for s = 128, 64, … 1). No subgroup operation, no floating-point atomic, no dependence on occupancy, scheduling or subgroup size. Two runs on the same device and driver give identical bits, and because nothing depends on the subgroup size the result is also the same for the fast variant without subgroups.
• A reduction with fewer than 32 elements per output, or between 32 and 255 along a non-contiguous axis with thousands of outputs, and a softmax row shorter than 16, is done by one thread per output in index order; that path is deterministic by construction and used by both modes.
• A softmax row of 16 to 1024 elements gets a power-of-two group of threads (about four elements each), so one workgroup handles 256 / group rows. The statistics are merged by the fixed halving tree in both modes.
• A reduction with fewer than 128 outputs and at least 16384 reduced elements (a sum over a whole tensor) runs in two passes: the reduced elements are cut into at most 256 contiguous chunks of at least 4096, a workgroup reduces each chunk to a partial result in a scratch buffer of 128 x 256 floats (ReduceKernels.ScratchLength; the family keeps its own unless the caller passes one), and one workgroup per output combines the partials in index order. The cut is a function of the shape, so the deterministic mode stays bit-identical run to run.
• Otherwise a workgroup takes each output. The choice depends on the shape only, never on the device.
• A workgroup beyond the rounded-up grid (the grid is two-dimensional, so it can hold more workgroups than there are outputs) returns before its first barrier. The workgroup id is uniform, so WGSL's uniformity analysis accepts the barriers that follow. The multi-row softmax keeps a per-thread bounds flag instead, because its rows differ between the threads of one workgroup.

Gather, scatter, layout and int8

• 8-bit payloads. GatherKernels also copies byte (uint8) and sbyte (int8) elements, for int8 embedding tables. They are packed four to a 32-bit word like bool tensors (allocate GpuKernelBuffers.PaddedBoolLength), one thread owns an output word, and the result is bit-exact. Scatter, concat and the other layout kernels do not take them.

• DequantizeLinear (DequantizeKernels) computes (float)(x - zeroPoint) * scale as onemodel's CPU kernel does: the subtraction wraps in 32 bits, the conversion rounds to nearest even, one float32 multiply follows. float32 output is bit-exact with the CPU on every device tested (including int32 inputs near the extremes); Half output rounds that value to nearest even. Per-axis scales and zero points follow the ONNX axis attribute; the scale is float32; a zero point has the element type of the input. Blocked quantization and 16-bit inputs are not covered.

• QuantizeLinear (DequantizeKernels.EncodeQuantize) computes saturate(round_half_even(x / scale) + zeroPoint) as onemodel's CPU kernel does, bit for bit: the division is the correctly rounded float32 division (done in integer arithmetic, as DynamicQuantizeKernels does, because a GPU's own division is not promised to be correctly rounded), the rounding is to nearest even (ONNX operator definition), the zero point is added in float32, the result is clamped to the output type's range, and NaN becomes the lowest value. The output is uint8 (byte) or int8 (sbyte) and so is the zero point; Half input is widened exactly first. Both modes give the same bits. Tested against the CPU formula on ties, saturation, NaN, infinities, subnormals and five scale sets: 0 mismatching bytes on every device.

• Payloads are 32-bit (float, int, uint for gather) or 64-bit (long) words, and Half (as 16-bit words). Results are bit-exact against the CPU for every payload, since these kernels only move values (scatter reductions on floats do one operation per update, in a fixed order).

• Invalid indices. A host-side Validate helper checks indices before encoding. On the device, an out-of-range index writes zero (gather) or is skipped (scatter), and, when the caller passes a GpuIndexStatus, the smallest offending position is recorded; ThrowIfFailedAsync raises GpuKernelIndexException naming it. Without a status object the zeros are silent.

• Scatter reductions are deterministic in both modes. Duplicate indices combine strictly in index order: every output element has one owner thread that scans the updates in order (a bucketed scheme, so large scatters still run in parallel). There are no floating-point atomics. indicesAreDistinct: true selects a plain parallel path when the caller knows no two updates meet. Ties in min/max keep the existing value, NaN propagates, integer add and mul wrap. Half scatter uses AND/OR on whole words, so it needs the same padding as other Half buffers. uint payloads are not supported for scatter.

• Int8 GEMM computes sum (a - zeroA) * (b - zeroB) in int32 (wrapping), so the result is exact and one variant serves both modes; it matches onemodel's CPU int8 operators bit for bit. The float form is float(sum + intBias) * scale followed by a float bias and the int32 form adds its bias to the sum; both are applied as each sum is stored, so a product is one dispatch. The float steps keep the two separate roundings of the unfused chain (the product of the multiply passes through an OR with a bit the host keeps zero, so no target can fuse it with the bias add); the tests compare against the chain on the CPU bit for bit. The requantized form is clamp(round_half_even(float(sum + bias) * multiplier) + zeroPointY); it writes bytes that neighbouring outputs share, so it keeps the sums in a scratch buffer for a second pass. That scratch is the instance's own and grows on demand (a larger product then builds a new bind group), unless the caller places it: Int8Scratch (buffer plus element offset) on EncodeRequantized, sized with Int8Gemm.ScratchElements, or ReserveScratch once at session creation. Int8Gemm.Variant reports the variant chosen for the device, and Int8GemmOptions.AllowPackedDot = false forces the portable one; both give the same bits. A product with more tiles along one axis than a dispatch grid holds (a very wide or tall one: more than 65535 tiles, that is 4.19 million columns at 64-wide tiles) numbers its tiles row-major over x and y; the tile a workgroup owns still depends only on its id, so the bits do not change.

  • Packed dot (packed-dot-tile64x64x32): chosen when the device reports an 8-bit dot product (GpuComputeCapabilities.Int8DotProduct: Vulkan shaderIntegerDotProduct, the browser's packed_4x8_integer_dot_product). Four byte pairs are multiplied and summed per instruction (dot4add_i8packed). That instruction takes signed bytes only, so unsigned operands are staged as u - 128 and the zero points move by the same amount; the zero points are then applied after the loop from byte sums collected while staging: sum (s - z')(t - w') = sum s*t - w' * sum s - z' * sum t + K * z' * w', exact modulo 2^32 like the sums themselves (Jacob et al., 2018, section 2.3 rearranges the same product). A byte sum is only collected when the other operand's shifted zero point can be non-zero. On a discrete GPU, products with fewer than 32 tiles of 64x64 run on 32x32 tiles so that the 128-row products of small models still fill the device (an integrated GPU took less time on the big tile everywhere).
  • Portable (portable-tile64x64x16): runs everywhere. Rows are staged a word at a time when the operand's offset, stride and leading dimension are multiples of four, and the products are formed in float, which is exact here (two values in -255..255 multiply to under 65025, so 256 of them sum to under 2^24); the running sums fold into int32 every 128 values of k. Metal has no packed dot product, so Apple GPUs run this one.

• DynamicQuantizeLinear runs in two passes (partial minimum and maximum per workgroup; then every workgroup combines the partials, the first stores scale and zero point, each quantizes its slice) and divides in integer arithmetic (a long division of the mantissas with round-to-nearest-even), because the targets' float division is not correctly rounded (Vulkan allows a few ulp, Metal compiles with fast math) and the scale and zero point must equal the CPU's bit for bit. An int8 linear layer is thereby three dispatches (two to quantize, one product) instead of six.

• Long submissions on Metal. The system aborts a command buffer that runs too long (kIOGPUCommandBufferCallbackErrorImpactingInteractivity). Encoders accept ReportWork(operations) and SplitSubmission(); Metal splits the command buffer at pass boundaries once a budget (GpuComputeCapabilities.SubmissionWorkBudget, 2^33 operations) is reached. Gemm cuts a product larger than the budget into row or batch pieces with the same accumulation order, so deterministic results are unchanged (GemmOptions.AllowSplitting). A command buffer that still fails surfaces as GpuCommandFailedException (Failure says watchdog, out of memory or invalid resource) from the next submit, wait or readback; never as a silent wrong result. The browser backend has no budget.

Special values

NaN and infinity follow IEEE 754 on Vulkan and Metal: Min, Max and the reductions Max/Min propagate NaN; ArgMax/ArgMin treat NaN as the extreme value (first NaN wins); softmax of a row with a NaN, a +infinity, or only −infinity gives NaN for the row, while −infinity beside finite values gives 0 (so additive masks work). The kernels test NaN with bit patterns, because Metal compiles with fast-math. In a browser, WGSL does not require NaN or infinity to behave (implementations may assume they never occur), so those cases are not tested on WebGPU and masks should use a large negative finite value there.

Limits

Element counts are below 2³¹ and every index is 32-bit. Broadcasting merges neighbouring dimensions that every operand walks contiguously, then needs at most 6; a reduction needs at most 4 runs of kept and 4 runs of reduced dimensions after merging; a transpose at most 6 unmergeable dimensions. Workgroups are 256 threads, which every Vulkan, Metal and WebGPU device offers; dispatches are two-dimensional so more than 65535 workgroups work everywhere.