SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

arXiv:2610.04277v1 [cs.DC] 3 Oct 2026

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs Yehong Jiang∗

Sheng Chen

Fangwen Fu

yehong28@mit.edu Massachusetts Institute of Technology Cambridge, Massachusetts, USA

sheng.s.chen@intel.com Intel Corporation Shanghai, China

fangwen.fu@intel.com Intel Corporation Santa Clara, California, USA

Yen-Kuang Chen

Xinmin Tian

Stuart H. Sul

y.k.chen@intel.com Intel Corporation Santa Clara, California, USA

xinmin.tian@intel.com Intel Corporation Santa Clara, California, USA

ssul@stanford.edu Stanford University Stanford, California, USA

Simran Arora simarora@caltech.edu California Institute of Technology Pasadena, California, USA

Abstract

CCS Concepts: • Computing methodologies → Parallel programming languages; • Computer systems organization → Single instruction, multiple data.

New AI accelerators arrive before the kernels that make them fast, because peak kernel performance requires architecturespecific expertise in operand pipelines and data-movement techniques. Coding agents can now write, compile, and tune kernels on their own, so they could greatly accelerate kernel development and optimization. What agents produce depends on the interfaces they are given. These interfaces may expose a machine’s raw capabilities or encode the knowngood methods for exploiting its hardware features efficiently. We test the performance impact of these two interfaces on Intel GPUs through controlled experiments with three coding models on three kernels. With raw SYCL and execution feedback, Opus 4.8, the strongest of the three models tested, writes a single-GPU GEMM kernel reaching only 57.5% of Intel’s tuned oneDNN library. We present SyclKittens, a hardware-aware tile programming model for Intel GPUs that encodes these known-good methods, so its operations place operands in matrix-engine layouts, prefetch through the L1 cache, and communicate over the fabric that links GPU stacks into a node. With SyclKittens under the same agent, task, and feedback budget, the agent reaches 82.1%, showing that encoded methods turn hardware capabilities into performance. Workload-specific policies stay programmable, and engineers and agents jointly refine these schedules in SyclKittens to build a kernel suite that reaches ∼96% of oneDNN in geometric mean across GEMM shapes. The suite runs Llama-3.1-8B inference 1.59× faster than torch.compile on one GPU and up to 2.91× faster than a matched multi-GPU decode path on Intel’s oneCCL. SyclKittens is open source and available at https://github.com/intel/SyclKittens.

Keywords: SyclKittens, tile-based programming, GPU kernel optimization, performance portability, SYCL, coding agents, Intel Max GPU, LLM inference

1

Introduction

Supporting diverse AI hardware is limited by the difficulty of building high-performance kernel ecosystems. Tile programming models such as ThunderKittens and HipKittens package architecture-specific methods for NVIDIA and AMD GPUs, while ParallelKittens extends this approach for multiple GPUs [8, 30, 32]. Coding agents can improve programs through execution feedback, and kernel benchmarks show that they can write, compile, benchmark, and revise singleand multi-GPU kernels [5, 27, 37]. We ask whether a tile programming model can make Intel GPU methods reusable by programmers and coding agents while preserving workloadspecific control. In our study, three agentic-coding models produce correct Intel GPU kernels that are substantially slower than the tuned baselines. They use raw SYCL, a C++ programming model for heterogeneous systems [19]. The tasks are single-GPU matrix multiplication (GEMM), causal groupedquery attention (GQA) [1], and distributed GEMM. Each model receives Intel documentation and execution feedback and works autonomously under a fixed search budget. On single-GPU GEMM, Opus 4.8 produces the fastest raw-SYCL kernel among the three models, but it reaches only 57.5% of the throughput of oneDNN, Intel’s tuned deep-learning library [24]. Its kernel already uses the GPU’s matrix engines and fetches operands ahead of computation. However, its

∗ Work performed while the author was an intern at Intel Corporation.

1

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Raw SYCL SyclKittens + examples Policy partial-tile ownership, prefetch distance, shape-dependent schedules

Capability matrix instructions, 2-D block prefetch, peer memory

100 Engineers + agents, 96.2% 80

% of tuned oneDNN

Method matrix-engine tile layouts, cache-centric operand delivery, cooperative prefetch

SyclKittens

60 40 20

exposed by raw SYCL 0

encoded by SyclKittens

1 searched by engineers and agents (a) Capabilities, methods, and policies.

4

8 12 16 Evaluator feedback calls

20

(b) Single-GPU GEMM tuning progress.

Figure 1. Programming interfaces and agent tuning progress. (a) SyclKittens supplies methods; programmers and agents choose policies. (b) For single-GPU GEMM, SyclKittens helps Opus 4.8 reach Raw SYCL’s final performance with fewer feedback calls. Curves show median best-so-far throughput relative to oneDNN on Intel Max GPU; the dashed line is the engineer-agent kernel. Examples come from NVIDIA ThunderKittens and AMD HipKittens. thread groups do not coordinate prefetching and repeatedly request the same operand tiles. We introduce SyclKittens, a hardware-aware tile programming model for Intel GPUs. We target the Intel Data Center GPU Max 1550 (Intel Max GPU) [4, 10], deployed in Aurora [2]. The tile interface must match the packed register layouts of Xe Matrix Extensions (XMX), Intel’s matrix engines, and coordinate operand delivery through cache prefetch and 2-D block loads. Communication operations must account for the bandwidth difference between MDFI (Multi-Die Fabric Interconnect), the faster link between the two GPU stacks within a package, and XeLink, the peer load/store network between packages. Built with Intel hardware experts, SyclKittens extends the tile abstractions of ThunderKittens and dataflow algorithms from ThunderKittens and HipKittens with three classes of Intel-specific operations.

tuned GEMM, replacing heuristic prefetch coordination with disjoint ownership raises geometric-mean throughput from 91.0% to 96.2% of oneDNN (§3). 3. Topology-aware GPU communication. To avoid repeated XeLink transfers, SyclKittens gathers each shard once per remote package and replicates it over MDFI. This yields 1.42–2.15 times the end-toend distributed-GEMM throughput of the flat gather pipeline (§3). Raw SYCL exposes hardware capabilities, including matrix instructions, prefetch, and peer memory. SyclKittens encodes methods that coordinate these capabilities. Programmers and coding agents choose workload policies, including tile shape, prefetch ownership and distance, and work distribution (Figure 1a). An accelerator interface should encode the architecture’s methods as first-class operations while leaving workload policy programmable. We perform three evaluations of SyclKittens. Controlled coding-agent experiments test whether encoding Intelspecific methods in the programming model helps agents produce faster kernels. The kernel-suite evaluation tests whether these abstractions support high performance across a broad range of operations. The inference evaluation measures whether those kernels deliver end-to-end gains for Llama-3.1-8B on single and multiple GPUs.

1. Matrix operand layouts and packing. SyclKittens maps register tiles to XMX’s native execution width and fuses 2-D block loads with hardware packing for the column-layout right-hand matrix operand B, avoiding software packing shuffles. In the cumulative GEMM ablation, changing from 16 × 16 to 32 × 32 register tiles raises geometric-mean throughput from 37.4% to 68.2% of oneDNN across 27 shapes (§3). 2. Cooperative cache prefetch. Tile-level prefetch lets groups of threads prefetch disjoint regions of shared operands into the L1 cache while computing. In our

• For controlled coding-agent experiments, we vary programming context under a fixed feedback protocol. 2

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

3. We present a co-designed Intel GPU kernel suite and evaluate its performance in individual kernels and Llama-3.1-8B inference (§5.2, §6).

The first context, Raw SYCL (Raw), provides Intel documentation. The second context, SyclKittens, adds its headers and interface documentation. The third context, SyclKittens + examples, provides kernel source code for NVIDIA and AMD GPUs, or for different operations on Intel GPUs. For Opus, the best-performing context on each task is either SyclKittens or SyclKittens + examples. Single-GPU GEMM performance rises from 57.5% of oneDNN under Raw to 82.1% with SyclKittens. In cold iteration, where agents start from scratch, SyclKittens helps Opus reach Raw’s final performance with fewer tuning steps (Figure 1b). For GQA, Opus reaches 79.5% of PyTorch’s native scaled dot-product attention (SDPA) with SyclKittens + examples, up from 57.7% with Raw SYCL. All selected Opus GQA kernels use XMX across the three contexts. The examples help the agent find more effective execution schedules. • We build a high-performance kernel suite with SyclKittens through engineer-agent co-design. The suite reaches ∼96% of oneDNN in geometric mean across all 27 GEMM shapes. On causal multi-head attention (MHA) [35], SyclKittens delivers 1.21× the geometricmean throughput of fused SDPA with Intel Extension for PyTorch (IPEX). The suite also includes normalization, rotary position embedding (RoPE) [31], matrixvector multiplication (GEMV), multi-GPU collectives, and fused distributed kernels. • We evaluate Llama-3.1-8B [7] inference on one GPU and on a full Intel Max node with six two-stack packages (12 GPUs). (i) On one GPU, SyclKittens works alongside torch.compile, improving full-request prefill-plus-decode throughput by 1.59× in geometric mean over torch.compile alone. SyclKittens also runs configurations where torch.compile alone runs out of memory. (ii) At node scale, SyclKittens combines tensor-parallel projections [29] with context-parallel attention. Its decode implementation outperforms the tuned inference path built on Intel’s oneCCL collective-communication library, by up to 2.91× on four GPUs, and serves longer contexts than the baseline can fit in memory (§6).

2

Background and Motivation

2.1

Tile Programming and Coding Agents

Tile programming models expose fixed-shape register and shared-memory tiles, layouts that address larger tensors in memory, and bulk movement, matrix, reduction, and synchronization operations. They sit between tensor compilers, which choose most schedules, and low-level GPU programming, which exposes individual lanes and instructions. This middle layer separates tile-level algorithms from the architecture-specific layouts, movement, and scheduling that realize them. ThunderKittens [30] and HipKittens [8] provide hardware-aware interfaces on NVIDIA and AMD, respectively, and their shared primitives provide the conceptual starting point for SyclKittens. A prior IWOCL study ported ThunderKittens directly to SYCL for Intel kernels on a single GPU [18]. SyclKittens extends that work with the Intelspecific operations in §3, including partial-tile cooperative prefetch and topology-aware communication, and evaluates their use in full-model inference. Other systems divide control between the programmer and the compiler differently. CuTe/CUTLASS provides finegrained layout and thread-participation control through C++ templates [22, 23], while TileLang separates tile dataflow from schedule selection [36]. For Intel GPUs, programmers can use portable SYCL [19], explicit SIMD (ESIMD) vectors [9], compiler-selected Triton schedules through Intel’s XPU backend [15, 33], or template libraries with native hardware control, such as XeTLA [11] and SYCL-TLA [14]. Intel’s oneAPI toolchain also provides joint_matrix, a portable SYCL extension that drives matrix engines while hiding the operand layouts SyclKittens makes explicit. At node scale, ParallelKittens [32] extends ThunderKittens to multiple GPUs and composes communication with computation through NVIDIA interconnects. Intel oneCCL provides tuned collective operations [34]. SyclKittens exposes peer-memory operations and synchronization within kernels so that programmers can coordinate data movement with computation across Intel GPU stacks and packages. These operations let kernels overlap communication with attention in the manner of ring attention [20]. AI coding agents can generate and iteratively improve programs, but for GPU kernels, correctness and performance are separate objectives. KernelBench [27] evaluates fast, correct single-GPU kernel generation over 250 PyTorch workloads. ParallelKernelBench [5] extends the setting to multi-GPU PyTorch+NCCL baselines and fine-grained CUDA or DSL kernels. Both establish that execution and profiling feedback help. Our controlled study on Intel GPUs separates the effect of the tile interface from that of example kernels. Raw SYCL

This paper makes three contributions. 1. We build and open-source SyclKittens, an Intel GPU tile programming model that encodes architectural methods while preserving workload control, and quantify their performance contributions through GEMM ablations (§3). 2. We provide a controlled coding-agent study that isolates the effects of encoded methods, example kernels, and execution feedback on kernel performance and workload-policy search (§4). 3

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

and SyclKittens test the interface change; SyclKittens and SyclKittens + examples test the additional source examples under the matched protocol in §4. 2.2

while cross-plane traffic traverses an additional MDFI link. In distributed GEMM, a topology-aware gather that crosses XeLink once per remote package and replicates over MDFI delivers 1.42–2.15× the end-to-end throughput of the flat gather pipeline (§3.4). Appendix A reports measurement methods, full values, and a cross-vendor terminology map. The tile interface separates algorithmic operations from their architecture-specific lowering. Intel’s upcoming Crescent Island data center GPU illustrates which hardware properties a new lowering must consider. It is based on the Intel Xe3p architecture, which retains matrix engines and an L1 operand path [21]. Its scale-up network uses PCIe switches instead of XeLink but remains a load/store fabric. Porting SyclKittens to it would require validating the operand layouts and prefetch parameters and adapting the communication schedules to the new topology.

Intel GPU Model

We study the Intel Data Center GPU Max 1550 (Intel Max GPU) [2, 4, 10, 17]. Each package is one PCIe device built from two stacks, Intel’s term for its physical GPU tiles. We run the flat device hierarchy, so each stack is exposed as its own SYCL device with 56 Xe cores and 64 GB of HBM2e, and a node presents 12 stacks across 6 packages. In workload and performance descriptions, we call each exposed stack a GPU ; we reserve stack and package for the physical topology. Throughout this paper, tile refers only to a register or shared tile in the language, and not to a partition of the hardware. Figure 2 summarizes the stack and node hierarchy, and the three hardware properties discussed below shape the SyclKittens design in §3. Execution and registers. Intel software addresses register operands as linear vectors along the SIMD direction, one element per lane, whereas NVIDIA tensor-core interfaces use matrix fragments distributed across a warp. The shared tile abstractions therefore require different physical register layouts on the two architectures. A subgroup is the SIMD execution group within a workgroup. The Intel Max GPU supports 16- or 32-lane subgroups, while its Xe Matrix Extensions (XMX) matrix engines only consume 16-lane subgroups. When all 16 lanes move one four-byte register element, the subgroup access spans exactly one 64-byte cache line. These execution and transaction widths determine the register-tile mapping in §3.2.

2.3

What Raw-SYCL Agents Implement

We inspect the cold-iteration GEMM kernels in §4 to identify which matrix and communication methods raw-SYCL agents implement. In single-GPU bfloat16 (BF16) GEMM, Sonnet and Haiku produce correct kernels that stage operands in SLM without using XMX. Their best correct kernels reach only 5.7% and 1.6% of Intel’s tuned oneDNN baseline. Opus builds an XMX pipeline from SPIR-V operations, including a Dot Product Accumulate Systolic (DPAS) matrix instruction, 2-D block loads, a 𝐵 operand in the interleaved VNNI layout that DPAS requires, and prefetch. However, its thread groups do not coordinate prefetching and repeatedly request the same operand tiles. Its best correct kernel reaches 57.5% of that baseline. oneDNN generates Intel GPU assembly through nGEN [25] and selects GEMM strategies from a tuned kernel database [26]. The comparison is therefore against oneDNN’s assembly-generating implementation. In the distributed GEMM task, input matrix 𝐴 is rowsharded, 𝐵 is replicated, and every GPU must return the full 𝐶 = 𝐴𝐵. A kernel can gather 𝐴 before computing the full output, or compute its local rows of 𝐶 and then gather the output. All nine raw-SYCL cold-iteration trials produce correct kernels, but the models combine these choices differently. Haiku computes local output rows but uses vector fused multiply-add (FMA) operations. Sonnet uses XMX but gathers the input and computes the full output on every GPU. Opus uses XMX and overlaps copies across compute and copy queues, but keeps the same replicated-compute decomposition and ignores the MDFI/XeLink topology. Across these trials, no raw-SYCL candidate combines local-output computation, XMX, and topology-aware communication. The best reaches 42.3% of our frozen, co-measured, topologyaware SyclKittens hand reference.

Operand delivery. The Intel Max GPU supports direct matrix-operand delivery through the cache hierarchy. A hardware 2-D prefetch warms an L1 cache of up to 512 KB per Xe core, and a 2-D demand load then moves the cached operand into subgroup registers. Shared local memory (SLM) remains available. The SyclKittens matrix path studied here uses L1 prefetch and 2-D demand loads without intermediate SLM staging (Figure 2(a)). This differs from the asynchronous global-to-shared-memory staging used by ThunderKittens. Node fabric. The Multi-Die Fabric Interconnect (MDFI) links the two GPU stacks within a package. XeLink, Intel’s peer load/store network between packages, connects packages all-to-all through independent even- and odd-stack planes. Together, these fabrics allow any pair of stacks to communicate directly (Figure 2(b)). Multicast and reduction are composed in software over these peer paths. XeLink also provides per-link copy engines for explicit transfers. However, peer-to-peer communications over these fabrics differ in bandwidth and the number of transfers they can sustain concurrently. MDFI provides substantially higher bandwidth than XeLink. Same-plane traffic uses a direct XeLink path, 4

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

(b) Node topology

(a) Stack hierarchy HBM2e 64 GB

0 10

Shared L2 cache 204 MB 56 Xe-cores execute workgroups in parallel Xe

Xe

Xe

···

Xe

Xe

2

8

Xe

4 6

Each Xe-core: 512 KB L1 + SLM 8 vector engines + 8 XMX matrix engines

1 3

11 9

5 7

0

1

2

3

4

5

6

7

8

9

10 11

Figure 2. Intel Max GPU at stack and node scale. (a) HBM2e memory and shared L2 cache feed the Xe-cores within a stack. (b) Six two-stack packages form two independent all-to-all XeLink planes, one for even stacks and one for odd stacks. MDFI connects the two stacks within each package and is the only link between the planes. T. The layout parameter selects row- or column-major organization, while shape selects a static physical base tile. The GEMM ablation compares 16×16 and 32×32 base tiles; rectangular shapes are also supported. Sub-tile coordinates address rectangular regions within a tile and support cooperative prefetch in §3.3. Bulk rectangular movement and matrix or reduction operations are shared with ThunderKittens. Kernels choose the ordering of load, compute, and store stages over workgroup tiles. A Global Layout (GL) pairs a device-memory pointer with logical batch, depth, row, and column dimensions. Tile movement operations use its coordinates rather than computing byte offsets in each kernel. A Parallel Global Layout (PGL) [32] extends this coordinate space across devices. Additional API and lowering details appear in Appendix B.

These kernels show that access to matrix instructions and peer memory leaves operand coordination and work distribution to the kernel author. SyclKittens encodes operand layouts, cooperative cache prefetch, and topology-aware communication as reusable operations. Programmers and coding agents use those operations while choosing tile shapes, prefetch ownership and distance, and the distribution of computation. Section 3 defines these operations, and §4 evaluates their effect on agent-written kernels.

3

SyclKittens

SyclKittens is a header-only C++ tile programming model for Intel GPUs. Kernels declare typed register or sharedmemory tiles, address tensors through global layouts, and express dataflow through tile movement, matrix and reduction operations, and synchronization. The backend lowers these operations to SYCL or SPIR-V and selects their physical lane mappings. Programmers retain C++ control flow and direct access to SYCL or SPIR-V intrinsics. Distributed layouts, peer memory, and GPU-resident synchronization extend the same programming model across a node. The Intel-specific operations encode matrix operand layouts and packing, cooperative cache prefetch, and topologyaware communication. Tile shape, prefetch ownership and distance, and work distribution remain programmable. In the cumulative single-GPU GEMM ablation, changes to operand layouts, prefetch, and workload-specific schedules raise the geometric mean from 37.4% to 96.3% of tuned oneDNN. The distributed movement ablation reaches 82– 90% of communication-free GEMM throughput across twelve GPUs. §3.1 defines the interface. §3.2 and §3.3 explain operand layouts, packing, and cooperative prefetch, and §3.4 describes communication within a node. 3.1

3.2

Matrix Operand Layouts and Packing

SyclKittens maps logical row- and column-major tiles onto linear vectors along the SIMD direction, using the 16-lane subgroups that XMX consumes. This differs from matrix fragments distributed across a 32-lane warp in ThunderKittens’ NVIDIA mapping. Fixing the subgroup width at 16 avoids a 32-to-16 operand reorganization inherited from that mapping. For BF16 data, a 32-element tile row matches one 64-byte L1 transaction. A 16-element row consumes only half the cache line, leaving the other half for a later access. The 32×32 base tile therefore halves L1 transactions per 32 values and consumes each line immediately. In the cumulative GEMM ablation, this change raises the geometric mean from 37.4% to 68.2% of oneDNN (Table 1).

Direct B-operand load and packing. The SyclKittens API supports two paths from global memory into registers. One stages through a shared tile in SLM, and the other moves directly from a Global Layout into a register tile.

Tile Interface

The register-tile type rt<T,M,N,layout,shape> represents a logical 𝑀×𝑁 two-dimensional array of elements of type 5

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 1. Six-node, 27-shape cumulative BF16-input/BF16output GEMM ablation on the Intel Max GPU. Values are geometric means of each shape’s median SyclKittens/oneDNN ratio across nodes; each row adds one change to the row above.

ThunderKittens and HipKittens use shared-memory staging followed by ldmatrix or ds_read_b64_tr_b16, respectively. The Intel Max GPU has no equivalent SLM matrixoperand load, so loading a packed B operand from SLM requires multiple ordinary reads. The SyclKittens matrix path instead loads operands directly through L1 cache, avoiding this load sequence and the additional staging and synchronization. For a column-layout right-hand matrix operand B, load_transform lowers to LoadTransformINTEL, a hardware operation that fuses a global 2-D block load with VNNI packing for DPAS. An ordinary load followed by software shuffles is measurably slower; Appendix B reports the packing implementation and its effect. 3.3

Configuration Direct 16×16 register tiles, no prefetch + 32×32 register tiles + selective subgroup-masked L1 prefetch + explicit unroll-32 + snake traversal + partial-tile de-duplication + shape-specific prefetch distance

Cooperative Cache Prefetch

SyclKittens separates operand delivery into a prefetch phase and a demand phase. prefetch_load issues a 2-D prefetch hint that warms L1 cache without writing the tile’s registers. A later 2-D demand load transfers the cached operands into subgroup registers. This split allows kernels to prefetch future operands while computing with the current tiles, without allocating a shared-memory staging buffer. prefetch_load_coop takes sub-tile offsets that identify the operand region to prefetch. Kernels assign these offsets to computing subgroups and choose the prefetch distance. The operation provides partial-tile prefetch; the schedule determines ownership and lookahead. In the cumulative GEMM ablation (Table 1), selective subgroup-masked prefetch raises the geometric mean to 88.8% of oneDNN. Loop unrolling and snake traversal, which keeps neighboring workgroups within a compact L2 reuse region, raise it to 91.0%. The masked schedule still warms shared operand regions more than once. With partial-tile de-duplication, each subgroup prefetches a disjoint region of the shared operands using the sub-tile coordinates from §3.1. This ownership assignment warms each shared cache region once and raises the mean to 96.2%. A final per-shape selection of prefetch distance closes the sequence at 96.3% of oneDNN. Appendix D reports the protocol, errors, and per-shape results for this GEMM ablation. 3.4

% oneDNN vs. prior 37.4% 68.2% 88.8%

– 1.822× 1.301×

89.1% 91.0% 96.2% 96.3%

1.003× 1.022× 1.056× 1.001×

Table 2. Distributed GEMM movement ablation on 12 GPUs, reported as end-to-end throughput relative to the same GEMM core measured without communication in the same run. Movement method

61443 122883

Flat kernel or copy-engine gather Pipelined gather and GEMM Fused gather and GEMM Topology-aware gather and GEMM

38% 40% 47% 82%

63% 62% 62% 90%

The topology API groups adjacent stacks into GPU packages with topo::group_pvc and assigns transfer ownership with topo::plan_all_gather. A hierarchical gather sends each shard over XeLink once per remote package, then replicates it to the partner stack over MDFI. This reduces repeated XeLink crossings relative to a flat gather. PGL permits both kernel-driven peer access and explicit copy-engine transfers. Kernels choose how to divide computation and schedule these transfers. Table 2 holds the GEMM core fixed and varies input-gather schedules to isolate communication changes. Each GPU gathers the missing rows of 𝐴 and computes the full 𝐶. The agent study in §4 also explores computing local rows of 𝐶 and then gathering the output. The topology-aware pipeline achieves 2.15× and 1.42× the end-to-end distributed-GEMM throughput of the flat gather pipeline at the smaller and larger shapes, respectively. The same primitives support collectives and ring attention. Appendices H, I, and J report full sweeps.

Topology-Aware GPU Communication

make_pgl combines per-device GLs and peer pointers, and pgl.dev_ptrs[peer] exposes remote high-bandwidth memory (HBM) to ordinary kernel loads and stores. The GPUs share a SYCL context, allowing kernels to dereference peer pointers directly. GPUs synchronize through peer-visible counters without returning to the host. After writing peer memory, the writing GPU publishes completion with a system-scope release atomic (sycl::atomic_ref). The reading GPU waits with a system-scope acquire atomic. Observing the matching release makes the preceding writes visible to that GPU. signal, wait, and barrier_all implement this protocol.

4

Controlled Coding-Agent Evaluation

We ask how encoded methods affect agent-written kernel performance, when example kernels help, and how execution feedback improves kernels started from scratch or from a correct implementation. We evaluate Opus 4.8, Sonnet 5, and 6

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

evaluator calls. Warm iterative starts from an eligible singleshot kernel, retains that source as checkpoint 0, and receives the same feedback budget. Under feedback, the harness retains the fastest candidate. After development, we grade the selected candidate on five unseen seeds, with inputs and results withheld from the agent. Every task-specific correctness criterion must hold on every seed.

Haiku 4.5 on single-GPU GEMM, distributed GEMM, and causal grouped-query attention (GQA), varying programming context and feedback under a fixed call budget. The SyclKittens contexts supply its complete interface, including cooperative prefetch, but not the co-designed kernels or tuned schedules evaluated in §3. Tables 3 and 4 report controlled endpoints. Figure 3 separately compares selected controlled and co-designed kernels in absolute throughput. 4.1

Performance reporting. Only correct kernels receive performance scores, and peak@3 is the fastest held-outcorrect result among the three starts. Single-GPU single-shot and cold-iterative throughput uses a fixed 197.5-TFLOP/s oneDNN baseline reached through the torch.mm dispatch on Intel XPU. Warm-iterative results use oneDNN co-measured with checkpoint 0 and the final kernel on the assigned GPU. Distributed throughput uses our frozen, co-measured SyclKittens hand reference, which gathers the input with a topology-aware schedule and computes the full output on every GPU. GQA uses native PyTorch XPU SDPA with enable_gqa=True and is_causal=True, co-measured on the same GPU and trial without importing IPEX. Figure 1(b) instead reports median best-so-far trajectories against comeasured oneDNN. Table 3 reports endpoints from correct kernels under the baseline rules above, so the values differ. Appendix C reports exact correctness gates, distributions, and protocol deviations.

Experimental Design

Automated harness and tasks. Our fully automated Claude Code harness evaluates Opus 4.8, Sonnet 5, and Haiku 4.5 on 40963 BF16 GEMM on one GPU, row-sharded 61443 all-gather GEMM on 12 GPUs, and causal GQA at 𝐵=1, 𝐻 =32, 𝐻𝐾𝑉 =8, 𝑁 =2048, and 𝐷=128. Programming contexts. Raw supplies Intel SYCL/SPIRV documentation. SyclKittens adds headers and API documentation. SyclKittens + examples additionally supplies task-specific kernel source code. For single-GPU GEMM, the examples are ThunderKittens and HipKittens GEMMs. For distributed GEMM, the context also adds a two-GPU SyclKittens all-reduce that shows Intel peer access but not an all-gather GEMM or 12-GPU schedule. For GQA, it adds two HipKittens GQA kernels. None provides an Intel implementation of the target kernel.

4.2

Hardware access and launch control. All contexts permit direct SYCL/SPIR-V intrinsics and forbid oneAPI joint_matrix. All GQA conditions use a byte-identical rawpointer harness and candidate function signature. Candidates control launch geometry and submissions in every context. Every GQA condition also forbids ESIMD, a separate Intel-specific SIMD programming model that neither raw SYCL nor SyclKittens builds on. We exclude Triton-XPU as a context because it changes both the code generator and programming interface, confounding the comparison. It appears instead as an autotuned performance baseline in §5.2.

Results

Encoded methods improve single-GPU GEMM. In cold iteration, all nine Opus single-GPU trials are correct, and peak@3 rises from 57.5% with Raw to 82.1% with SyclKittens and 82.6% with SyclKittens + examples. Sonnet also improves from 5.7% with Raw to 85.2% with SyclKittens. The tile programming model supplies XMX computation, operand layout, and data movement as operations that agents can compose while tuning schedules. Without feedback, SyclKittens raises Opus correctness from 1/3 to 3/3 and peak@3 from 6.6% to 63.4% on one GPU, while raising distributed GEMM correctness from 0/3 to 2/3. The median Opus SyclKittens trajectory in Figure 1(b) passes Raw’s final level on the second evaluator call. This measures search efficiency separately from the peak@3 endpoints.

Replicates and scope. Each model-context pair receives three independent starts in both single-shot and cold-iterative evaluation. Opus, Sonnet, and Haiku are one vendor’s large, mid-size, and small tiers, chosen to span capability. We measure how programming context affects generated kernels under a fixed harness and call budget.

Context effects vary across models and tasks. All 27 cold-iterative GQA trials are correct, but their peak@3 context ordering changes by model. SyclKittens + examples leads for Opus at 79.5%, whereas Raw leads for Sonnet at 18.6%, above SyclKittens + examples at 8.8% and SyclKittens at 4.6%. Haiku remains between 0.7% and 1.3% across all contexts. The best Sonnet result is below the weakest Opus context at 41.6%. Sonnet’s reversal from GEMM to GQA shows that a context’s benefit also depends on the task. Every selected Opus GQA kernel uses XMX through direct subgroup DPAS, so its gain from examples reflects schedule search rather

Feedback and selection. We use a custom Claude Code harness rather than the mini-SWE-agent implementation used by ParallelKernelBench [5, 37]. Both protocols use compilation, correctness, and performance feedback, while ours fixes condition-specific staged files and one evaluator. In single-shot runs, agents can inspect and edit staged files but cannot compile or execute candidates. Cold iterative starts from the context’s blank task stub and receives at most 20 7

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 3. Controlled BF16 GEMM on Intel Max GPUs, with 40963 on one GPU and 61443 on 12 GPUs. Single-shot and cold-iterative cells report held-out-correct trials / 3 and peak@3. Warm-iterative cells report held-out-correct finals / eligible sources and the fastest correct final. Warm rows here report Opus. The single-GPU denominator is oneDNN. The distributed denominator is our frozen, co-measured topology-aware input-gather SyclKittens hand reference. Programming contexts are defined in §4.1. Single-GPU paired details are in Appendix C.4, and the distributed all-model study is in Appendix C.5. A dash means no correct or applicable cell. Protocol

Single-GPU GEMM (% oneDNN)

Model

Distributed GEMM (% SK ref.)

Raw

SyclKittens

+ examples

Raw

SyclKittens

+ examples

Single-shot

Opus Sonnet Haiku

1/3; 6.6 0/3; – 2/3; 0.8

3/3; 63.4 3/3; 53.4 2/3; 47.5

3/3; 81.2 3/3; 31.0 0/3; –

0/3; – 3/3; 5.7 2/3; 0.5

2/3; 36.7 3/3; 5.4 1/3; 0.5

3/3; 59.4 3/3; 38.8 1/3; 0.5

Cold iterative

Opus Sonnet Haiku

3/3; 57.5 3/3; 5.7 3/3; 1.6

3/3; 82.1 3/3; 85.2 3/3; 0.6

3/3; 82.6 3/3; 84.8 3/3; 50.3

3/3; 42.3 3/3; 25.8 3/3; 17.2

3/3; 46.8 3/3; 50.9 3/3; 0.7

3/3; 44.0 3/3; 41.1 2/3; 1.5

Warm iterative

Opus

1/1; 14.2

3/3; 87.8

3/3; 87.7

–

2/2; 106.6

3/3; 61.5

Table 4. Controlled cold-iterative causal GQA on one Intel Max GPU, using the setup in §4.1. Each cell reports heldout-correct trials / 3 and peak@3, the fastest held-out-correct result among the three cold-iterative starts, as a percentage of each trial’s co-measured native PyTorch XPU SDPA baseline. All 27 trials are correct. Replicate-level and single-shot results appear in Appendix C.7. Model

Raw

Opus 4.8 3/3; 57.7 Sonnet 5 3/3; 18.6 Haiku 4.5 3/3; 0.7

output-row shards before gathering the output, avoiding redundant full-GEMM work on every GPU (Appendix C.5). Source inspection identifies implementation choices. As §2.3 shows, selected Raw kernels either bypass XMX or fail to coordinate the full operand pipeline. Every context exposes the same SPIR-V and XMX instructions, but the strongest no-feedback Opus kernel in SyclKittens + examples combines mma_AB, 2-D loads, VNNI packing, a 32×32 register base shape, and three-tile prefetch, while the correct singleshot Raw kernel only uses vector FMA. The SyclKittens interface supplies the matrix and movement operations, while the generated kernel chooses the tile shape and prefetch distance. Full source inspection appears in Appendix C.

SyclKittens + examples 3/3; 41.6 3/3; 4.6 3/3; 1.0

3/3; 79.5 3/3; 8.8 3/3; 1.3

5

than access to a different low-level mechanism. Full GQA results appear in Appendix C.7.

SyclKittens Inference Kernel Suite

Engineers and coding agents co-designed the SyclKittens inference kernel suite. Unlike the autonomous trials in §4, they jointly refine tile shapes, prefetch ownership and distance, and work distribution under correctness and performance checks set by engineers. This section reports the suite’s throughput against tuned vendor paths.

Example kernels help in some settings. Without feedback, examples raise Opus peak@3 from 63.4% to 81.2% on one GPU and from 36.7% to 59.4% on distributed GEMM, where correctness also rises from 2/3 to 3/3. In cold iteration, Opus reaches about 82% on single-GPU GEMM in both SyclKittens contexts, and Sonnet and Haiku vary by task and protocol. The examples’ measured benefit therefore depends on the model, task, and starting condition.

5.1

Experimental Setup

PyTorch 2.10 with IPEX 2.10.10 [12] dispatches operator baselines to oneDNN [24, 28], fused SDPA, or oneCCL. Our standalone decode GEMV baseline calls Intel’s optimized oneMKL math library directly through oneapi::mkl::blas::gemm [13]. Fusion and end-to-end results use torch.compile. For the Triton comparison, we use identical shapes, correctness gates, and measurement windows with Triton-XPU 3.6’s standard GEMM implementation, autotuned. Every sweep warms up before timing and co-measures each SyclKittens candidate with its baseline in the same

Agents consistently refine correct kernels. Warm iterative refinement retains the original correct kernel while feedback improves its schedule. All seven eligible Opus singleGPU GEMM sources and all six Opus GQA sources remain correct and improve, reaching 87.8% and 83.4% of their respective baselines. All thirteen eligible distributed sources across Opus, Sonnet, and Haiku in the two SyclKittens contexts also remain correct and improve. The strongest reaches 106.6% of our co-measured reference by computing local 8

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

125.7 co-design

123.2 remeasured

123.8 remeasured

Engineers + agents

49.8 cold

80

77.3 co-design

71.1 warm

43.1 cold

120

SyclKittens + examples

60.2 cold

160

168.3 warm

113.6 cold

TFLOP/s

200

163.1 warm

240

SyclKittens

187.6 co-design

Raw SYCL

280

40 0 Single-GPU GEMM (40963 )

Causal GQA (𝑁 =2048)

Distributed GEMM (61443 , 12 GPUs)

Figure 3. Selected implementations from different development settings. Bars compare selected correct Opus kernels with engineer-agent kernels in TFLOP/s. Annotations identify cold or warm iteration, remeasurement, and co-design. Raw single-GPU GEMM uses peak@3; the two SyclKittens contexts use remeasured kernels from warm iteration. GQA reports peak@3 after remeasuring all kernels with the same harness. Distributed GEMM uses a single-node Raw measurement and node-aggregated measurements for the other kernels. Programming contexts are defined in §4.1. These bars are separate from the controlled endpoint summaries in Tables 3 and 4. Decode-critical operators. The remaining SyclKittens kernels outperform their vendor baselines on the operations that dominate autoregressive decode. These operations either run briefly, so launch overhead matters, or stream weights from memory. Fused root mean square layer normalization (RMSNorm) [38], layer normalization (LayerNorm) [3], and RoPE each use one launch. Split-K GEMV matters most because the five 𝑀=1 Llama projections account for about 64% of per-token decode latency. SyclKittens maps each outputrow reduction across cooperating subgroups so the GPU stays occupied as tensor-parallel shards shrink. Figure 5(b) reports these gains before §6.1 measures them inside a full decode loop. Appendices F and G report the full grids.

thermal window. We report geometric means of per-shape speedups. SyclKittens GEMM requires relative Frobenius error below 0.02. Appendix D reports per-shape data and additional correctness gates.

5.2

Kernel Results

Dense GEMM across shapes and nodes. The GEMM of §3.3 reaches 96.3% of oneDNN in geometric mean across the 27 shapes in Figure 5(a). Across six independently measured nodes at 81923 , the per-node SyclKittens/oneDNN relative performance ranges from 95.5% to 96.3%. An interleaved six-node sweep finds that the final cooperative revision exceeds autotuned Triton-XPU 3.6 on all 27 shapes, with 1.087× geometric-mean throughput. Appendix D reports the direct comparison, exact tables, Stream-K results, and ablations.

6

End-to-End Inference

We measure how the kernel gains translate to Llama-3.1-8B inference [7] on one GPU and across a node. Composed with torch.compile, SyclKittens improves full-request throughput by 1.59× in geometric mean over torch.compile alone on the single-GPU configurations. Across twelve GPUs, SyclKittens also serves a 256K context. That memorycapacity result is separate from the matched-GPU-count speedups below.

Attention reuses the HipKittens algorithm. The attention kernel reuses the complete HipKittens dataflow algorithm through the tile interface. For each query tile, the GQA kernel retains HipKittens’ sequence of tiled 𝑄𝐾 ⊤ , runningmaximum and normalizer updates, output rescaling, and 𝑃𝑉 accumulation. The dataflow therefore needs no algorithmic change. SyclKittens re-implements its tile operations with Intel’s matrix, movement, and prefetch instructions. Figure 4 reports the resulting MHA and GQA throughput. The GQA baseline dispatches through its MHA path and repeats eight KV heads to match the 32 query heads, so that margin is against a fallback rather than a tuned vendor GQA kernel.

6.1

Framework Composition on One GPU

Operator ownership. We compare SyclKittens composed with torch.compile with torch.compile alone. torch.compile generates Triton for the pointwise, standalone normalization, and epilogue regions it fuses, but 9

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

(b) Causal GQA (𝐵=1, 𝐻 =32, 𝐻 KV =8)

(a) Causal MHA (𝐵=4, 𝐻 =𝐻 KV =16) PyTorch

80

SyclKittens

TFLOP/s

TFLOP/s

80 60 40 20

60 40 20

0

0 512

1024

2048

4096

6144

8192

512

Sequence Length 𝑁

1024

2048

4096

6144

8192

Sequence Length 𝑁

Figure 4. Forward causal attention on one Intel Max GPU. (a) MHA reaches 1.21× geometric mean over fused SDPA. (b) GQA reaches 4.0× in geometric mean over the baseline’s repeat-KV path, which dispatches GQA through MHA. Both panels use the same scale.

(b) Decode-critical operators

1.15

Speedup over baseline

SyclKittens / oneDNN

(a) BF16 GEMM

1.1 1.05 1 0.95 0.9 square

rectangular

0.85 1

5

9

13

17

21

25 27

Canonical shape index

HBM-scale

Few-token

GEMV

2

1

0 RMSNorm

LayerNorm

RoPE

GEMV

Figure 5. SyclKittens kernel suite. (a) The SyclKittens GEMM reaches 96.3% of the co-measured tuned oneDNN library path in geometric mean across 27 GEMMs, exceeds it on three shapes, and peaks at 1.17×. Shape indices 1–13 are square and 14–27 rectangular. (b) Fused normalization and RoPE are near the torch.compile HBM-scale path but lead it by 1.98–2.49× for few-token decode. GEMV has one bar because weight streaming is bandwidth-bound even in decode. SyclKittens’s wide-load split-K GEMV reaches 1.36× oneMKL in geometric mean and 1.54× on 𝑊𝑜 . graph-breaks at attention on XPU. Across this boundary, SyclKittens supplies a FlashAttention-style kernel [6], fused Q/K rotary, fused QKV projection, and residual-plusRMSNorm fusion. The kernels are exposed as PyTorch custom operations and installed by replacing the Llama decoder layer at runtime; HuggingFace source and model.generate() remain unchanged.

grow with context because the SyclKittens FlashAttentionstyle kernel avoids the baseline’s quadratic score materialization; decode gains come from eliminating launches around that attention boundary and from the projection path in §5.2. Neither system uses SYCL Graph, Intel’s equivalent of CUDA graph capture, in these measurements.

6.2

Prefill and decode. Figure 6 reports Llama-3.1-8B inference throughput over four batch sizes and six context lengths. Each point times one synchronized, deterministic model.generate() request with up to 128 new tokens. The reported full-request rate counts batch-wide prompt tokens plus the generated sequence length of one stream, as defined in Appendix K.4. Geometric means use the 21 configurations in which both systems fit. The composed path improves prefill by 1.22×, decode by 1.65×, and full-request prefill-plusdecode throughput by 1.59× (range 1.29–1.98×). Prefill gains

Node-Scale Composition

At node scale, decode performance depends on communication latency. In the evaluated batch-one decode, each tensorparallel all-reduce carries one 4096-element BF16 hidden vector (8 KiB), so launch, process, and synchronization latency matter more than peak bulk bandwidth. SyclKittens keeps all GPUs in one process and one shared SYCL context; a host thread launches each GPU’s stream, and inline peer movement with system-scope device barriers keeps Python and cross-process dispatch off the per-token path. 10

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

(a)

(b)

Prefill throughput (tok/s) N=256

N=512

6.1k

8.7k

9.5k

9.7k

9.3k

7.6k

6.7k

8.4k

8.9k

8.3k

6.8k

3.7k

8.3k

9.2k

9.5k

9.5k

9.5k

8.1k

8.1k

8.8k

8.6k

7.5k

6.1k

OOM

9.7k

10.0k

10.4k

10.0k

9.4k

8.2k

9.5k

9.8k

9.1k

7.8k

4.9k

OOM

10.3k

10.5k

10.5k

10.1k

9.5k

7.7k

10.1k

10.1k

9.3k

7.7k

5.0k

OOM

BS=1

BS=2

BS=4

BS=8

Decode throughput (tok/s per stream)

N=1024 N=2048 N=4096 N=8192

N=256

N=512

49

48

46

44

39

32

25

25

25

26

29

25

47

46

44

40

35

27

25

25

26

25

25

OOM

44

43

40

36

29

24

26

26

25

25

17

OOM

42

40

36

31

22

16

26

25

25

21

10

OOM

BS=1

BS=2

BS=4

BS=8

N=1024 N=2048 N=4096 N=8192

Figure 6. Llama-3.1-8B composition on one Intel Max GPU. (a) Prefill and (b) decode throughput across batch size and prompt length. Each cell reports SyclKittens composed with torch.compile (bold, top) and torch.compile alone (bottom); color encodes their ratio. Ratios use the 21 cells where both systems fit; in the three OOM cells the baseline runs out of memory while SyclKittens serves the request. (a) Six-package ring attention

60

400

TPOT (ms/token)

TFLOP/s per GPU

(b) Llama-3.1-8B decode

40

20

0 1.5K

6K

12K

24K

49K

oneCCL 𝑤𝑠=4 SyclKittens 𝑃 =4 SyclKittens best layout

300 200 100

98K 0

Total context length PyTorch

SyclKittens sequential

1K

4K

16K

64K

256K

Context length

SyclKittens overlap

Figure 7. Node-scale composition on Intel Max GPUs. (a) Non-causal BF16 ring attention uses six packages with one active stack each, comparing copy-engine overlap with sequential SyclKittens and blocking PyTorch schedules. (b) Multi-GPU decode time per output token; lower is better. The orange and dashed-blue curves compare matched four-GPU paths, where SyclKittens is 2.91× faster at 32K. The green curve uses four GPUs through 32K, six at 64K and 128K, and twelve at 256K. Its longer-context points demonstrate the best measured layouts, not matched-resource speedups. Inference-scale collectives. At 24 KiB, SyclKittens outperforms oneCCL on all four collectives by 1.81–3.31×. Reducescatter reaches 2.20× at 0.19 MiB and remains 1.10× faster at 96 MiB. Destination-split all-to-all reaches 3.49× at 1.5 MiB and remains 1.20× faster at 96 MiB. At that largest payload, oneCCL leads all-reduce and all-gather, where SyclKittens reaches 0.94× and 0.92× of oneCCL.

per step and merges the online softmax state. Copy-engine overlap leads the blocking PyTorch sendrecv ring at every context we measure, by 1.21–2.46× (Figure 7(a)). At 49K total tokens, the overlapped non-causal ring reaches 60.7 TFLOP/s per GPU, versus 42.0 for the blocking baseline. All-gather exceeds the sequential ring across the measured context sweep when the gathered KV cache fits in HBM; the ring provides the memory-scaling policy beyond that point.

Context-parallel attention. Long contexts also require splitting attention across GPUs. When the full KV cache fits in each GPU’s HBM, SyclKittens gathers it once and executes one flash pass. When it does not, a ring rotates one KV shard 11

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

6.3

Acknowledgments

Multi-GPU Decode

The Llama decoder composes tensor-parallel projections and inline all-reduces with context-parallel attention. Figure 7(b) shows that, at matched GPU count, the SyclKittens decode path outperforms the oneCCL-based path at every context, by up to 2.91× at 32K on four GPUs. This comparison includes both the source-partitioned collective and SyclKittens’s single-process execution model. The standalone results above measure complete collective calls without the surrounding decoder, while retaining each implementation’s process model. The best layout changes with context. Flat four-way tensor parallelism wins through 32K. At 64K and 128K, TP2×CP3 uses six GPUs to shard the growing KV scan, and at 256K, TP4×CP3 uses all twelve GPUs and serves at 73.3 ms/token. The last point is a capacity result, since the two-GPU tuned oneCCL path is out of memory there and the SyclKittens configuration uses more GPUs. Appendices I, J, and K provide per-shape tables, integration and correctness details, perlayer collective costs, the ring-attention step breakdown, timing-tail summaries, and layout sweeps.

7

Generative AI tools were used to assist with kernel development and optimization (documented in the appendices) and with editing the manuscript. The authors reviewed and verified all resulting content, including text, code, data, and figures, and take full responsibility for it.

References [1] Joshua Ainslie, James Lee-Thorp, Michiel de Jong, Yury Zemlyanskiy, Federico Lebrón, and Sumit Sanghai. 2023. GQA: Training Generalized Multi-Query Transformer Models from Multi-Head Checkpoints. In Proceedings of the 2023 Conference on Empirical Methods in Natural Language Processing. 4895–4901. doi:10.18653/v1/2023.emnlp-main. 298 [2] William E. Allcock, Benjamin S. Allen, James Anchell, Victor Anisimov, Thomas Applencourt, Abhishek Bagusetty, Ramesh Balakrishnan, Riccardo Balin, Solomon Bekele, Colleen Bertoni, Cyrus Blackworth, Renzo Bustamante, Kevin Canada, John Carrier, Christopher Channui, Lance C. Cheney, Taylor Childers, Paul Coffman, Susan Coghlan, Tanima Dey, Michael D’Mello, Ashok Emani, Murali Emani, Kyle G. Felker, Sam Foreman, Olivier Franza, Longfei Gao, Marta García, María Garzarán, Balazs Gerofi, Yasaman Ghadar, Subrata Goswami, Neha Gupta, Kevin Harms, Väinö Hatanpää, Brian Holland, Carissa Holohan, Brian Homerding, Khalid Hossain, Xue Hu, Louise Huot, Huda Ibeid, Joseph A. Insley, Sai Jayanthi, Hong Jiang, Wei Jiang, Xiao-Yong Jin, Jeongnim Kim, Christopher Knight, Panagiotis Kourdis, Kalyan Kumaran, JaeHyuk Kwack, Janghaeng Lee, Ti Leggett, Ben Lenard, Chris Lewis, Nevin Liber, Johann Lombardi, Raymond M. Loy, Ye Luo, Bethany Lusch, Nilakantan Mahadevan, Beth Markey, Victor A. Mateevitsi, Gordon McPheeters, Ryan Milner, Jerome Mitchell, Vitali A. Morozov, Servesh Muralidharan, Tom Musta, Mrigendra Nagar, Vikram Narayana, Marieme Ngom, Anthony-Trung Nguyen, Nathan Nichols, Aditya Nishtala, James C. Osborn, Michael E. Papka, Scott Parker, Saumil S. Patel, Julia Piotrowska, Adrian C. Pope, Sucheta Raghunanda, Esteban Rangel, Paul M. Rich, Katherine M. Riley, Silvio Rizzi, Kris Rowe, Varuni Sastry, Adam Scovel, Filippo Simini, Haritha Siddabathuni Som, Patrick Steinbrecher, Rick Stevens, Xinmin Tian, Peter Upton, Thomas Uram, Archit K. Vasan, Álvaro Vázquez-Mayagoitia, Kaushik Velusamy, Brice Videau, Venkatram Vishwanath, Brian Whitney, Timothy J. Williams, Michael Woodacre, Sam Zeltner, Chuanjun Zhang, Gengbin Zheng, and Huihuo Zheng. 2025. Aurora: Architecting Argonne’s First Exascale Supercomputer for Accelerated Scientific Discovery. arXiv:2509.08207 [cs.DC] https://arxiv.org/abs/2509.08207 [3] Jimmy Lei Ba, Jamie Ryan Kiros, and Geoffrey E. Hinton. 2016. Layer Normalization. arXiv preprint arXiv:1607.06450 (2016). https://arxiv. org/abs/1607.06450 [4] David Blythe. 2021. Intel’s Ponte Vecchio GPU Architecture. Hot Chips 33 Symposium on High-Performance Chips. https://hc33.hotchips. org/assets/program/conference/day2/hc2021_pvc_final.pdf [5] Willy Chan, Nathan Paek, Simon Guo, Simran Arora, and Daniel Fu. 2026. ParallelKernelBench: Can LLMs Write Fast Multi-GPU Kernels? https://github.com/togethercomputer/ParallelKernelBench Benchmark and evaluation harness for multi-GPU kernel generation. [6] Tri Dao, Daniel Y. Fu, Stefano Ermon, Atri Rudra, and Christopher Ré. 2022. FlashAttention: Fast and Memory-Efficient Exact Attention with IO-Awareness. In Advances in Neural Information Processing Systems, Vol. 35. 16344–16359. doi:10.52202/068431-1189 [7] Aaron Grattafiori et al. 2024. The Llama 3 Herd of Models. arXiv preprint arXiv:2407.21783 (2024). https://arxiv.org/abs/2407.21783 [8] William Hu, Drew Wadsworth, Sean Siddens, Stanley Winata, Daniel Y. Fu, Ryann Swann, Muhammad Osama, Christopher Ré, and Simran Arora. 2025. HipKittens: Fast and Furious AMD Kernels. arXiv preprint

Conclusion

SyclKittens is a tile programming model for Intel GPUs that provides matrix operand layouts, cooperative cache prefetch, and topology-aware GPU communication as reusable operations. With SyclKittens, the strongest coding agent in our study writes a faster single-GPU GEMM kernel than in raw SYCL under the same task and feedback budget. Through engineer-agent co-design, SyclKittens also supports a kernel suite that approaches oneDNN on GEMM and outperforms SDPA with IPEX on causal MHA. In Llama-3.1-8B inference, it outperforms torch.compile on one GPU and the oneCCLbased decode path on multiple GPUs. A study of human experts in raw SYCL and SyclKittens would test whether the same interface also reduces expert development effort. Once Crescent Island is available, we can test whether SyclKittens’s encoded compute methods transfer while its node schedules are retuned for a PCIebased fabric [21]. As coding agents take on kernel development for new accelerators, our results suggest that the programming interface should carry each architecture’s methods. Agents then start from known-good building blocks and spend their search on workload-specific policies. Contributions. Jiang led implementation, experiments, and writing. S. Chen contributed implementation and debugging. Fu, Y.-K. Chen, and Tian contributed Intel Xe-HPC architecture expertise and baselines. Sul contributed discussion, review, and writing. Arora provided direction, supervision, and review. 12

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

arXiv:2511.08083 (2025). https://arxiv.org/abs/2511.08083 [9] Intel Corporation. 2023. Explicit SIMD Programming Extension for DPC++. https://github.com/intel/llvm/blob/sycl/sycl/doc/extensions/ supported/sycl_ext_intel_esimd/sycl_ext_intel_esimd.md [10] Intel Corporation. 2023. Intel Data Center GPU Max 1550 Specifications. https://www.intel.com/content/www/us/en/products/ sku/232873/intel- data- center- gpu- max- 1550/specifications.html Launched Q1 2023; codename Ponte Vecchio. [11] Intel Corporation. 2023. Intel Xe Templates for Linear Algebra (XeTLA) v0.3.7. https://github.com/intel/xetla/releases/tag/v0.3.7 [12] Intel Corporation. 2025. Intel Extension for PyTorch v2.10.10+xpu. https://github.com/intel/intel-extension-for-pytorch/releases/tag/ v2.10.10%2Bxpu aurora_frameworks 2025.3.1 (PyTorch 2.10). [13] Intel Corporation. 2025. Intel oneAPI Math Kernel Library – Data Parallel C++ Developer Reference: gemm. https://www.intel.com/ content/www/us/en/docs/onemkl/developer-reference-dpcpp/20250/gemm.html [14] Intel Corporation. 2025. SYCL* Templates for Linear Algebra (SYCL*TLA) v0.6. https://github.com/intel/sycl-tla/releases/tag/v0.6 [15] Intel Corporation. 2026. Intel XPU Backend for Triton v3.6.0. https: //github.com/intel/intel-xpu-backend-for-triton/releases/tag/v3.6.0 [16] Albert Q. Jiang, Alexandre Sablayrolles, Arthur Mensch, Chris Bamford, Devendra Singh Chaplot, Diego de las Casas, Florian Bressand, Gianna Lengyel, Guillaume Lample, Lucile Saulnier, Lélio Renard Lavaud, Marie-Anne Lachaux, Pierre Stock, Teven Le Scao, Thibaut Lavril, Thomas Wang, Timothée Lacroix, and William El Sayed. 2023. Mistral 7B. arXiv preprint arXiv:2310.06825 (2023). https: //arxiv.org/abs/2310.06825 [17] Hong Jiang. 2022. Intel’s Ponte Vecchio GPU: Architecture, System and Software. Hot Chips 34 Symposium on High-Performance Chips. https://www.hc34.hotchips.org/assets/program/conference/ day1 / GPU % 20HPC / Intel _ s % 20Ponte % 20Vecchio % 20GPU % 20 %20Architecture%20Systems%20and%20Software%20FINAL.pdf [18] Yehong Jiang, Sheng Chen, Fangwen Fu, Vincent Lu, Yen-Kuang Chen, Hong Wang, and Xinmin Tian. 2026. Porting ThunderKittens from CUDA to SYCL for Intel GPU: Process, Challenges, and Lessons Learned. In Proceedings of the International Workshop on OpenCL and SYCL. ACM, Heilbronn, Germany, 1–5. doi:10.1145/3811257.3811269 [19] Khronos Group. 2026. SYCL 2020 Specification. Technical Report. Khronos Group, Inc. https://registry.khronos.org/SYCL/specs/sycl2020/html/sycl-2020.html Revision 12, August 6, 2026. [20] Hao Liu, Matei Zaharia, and Pieter Abbeel. 2024. RingAttention with Blockwise Transformers for Near-Infinite Context. In The Twelfth International Conference on Learning Representations (ICLR). https: //openreview.net/forum?id=WsRHpHH4s0 [21] Sumit Mohan and Hong Jiang. 2026. Crescent Island: GPU Designed for Agentic AI Inference. Hot Chips 38 Symposium on High-Performance Chips. https://hc2026.hotchips.org/program/ [22] NVIDIA. 2024. CuTe: C++ CUDA Template Abstractions for Tensors and Layouts. https://github.com/NVIDIA/cutlass/blob/main/media/ docs/cpp/cute/00_quickstart.md [23] NVIDIA. 2025. CUTLASS: CUDA Templates for Linear Algebra Subroutines and Solvers. https://github.com/NVIDIA/cutlass [24] oneDNN Project. 2024. oneAPI Deep Neural Network Library (oneDNN). https://github.com/uxlfoundation/oneDNN Open-source project hosted by the UXL Foundation. [25] oneDNN Project. 2025. oneDNN v3.10 GPU JIT-compiled Kernels. https://github.com/uxlfoundation/oneDNN/blob/v3.10/src/gpu/ intel/jit/README.md [26] oneDNN Project. 2025. oneDNN v3.10 Intel GPU GEMM JIT Strategy Catalog. https://github.com/uxlfoundation/oneDNN/blob/v3.10/src/ gpu/intel/gemm/jit/selector/db/kernel.db [27] Anne Ouyang, Simon Guo, Simran Arora, Alex L. Zhang, William Hu, Christopher Ré, and Azalia Mirhoseini. 2025. KernelBench: Can LLMs

Write Efficient GPU Kernels?. In Proceedings of the 42nd International Conference on Machine Learning (Proceedings of Machine Learning Research, Vol. 267). PMLR, 47356–47415. https://proceedings.mlr. press/v267/ouyang25a.html [28] PyTorch Contributors. 2026. PyTorch v2.10.0 XPU ATen Matrix Multiplication Implementation. https://github.com/pytorch/pytorch/ blob/449b1768410104d3ed79d3bcfe4ba1d65c7f22c0/aten/src/ATen/ native/mkldnn/xpu/Blas.cpp#L182-L206 Release commit 449b176; xpu::mm_out calls onednn::matmul. [29] Mohammad Shoeybi, Mostofa Patwary, Raul Puri, Patrick LeGresley, Jared Casper, and Bryan Catanzaro. 2019. Megatron-LM: Training Multi-Billion Parameter Language Models Using Model Parallelism. arXiv preprint arXiv:1909.08053 (2019). https://arxiv.org/abs/1909. 08053 [30] Benjamin F. Spector, Simran Arora, Aaryan Singhal, Arjun Parthasarathy, Daniel Y. Fu, and Christopher Ré. 2025. ThunderKittens: Simple, Fast, and Adorable Kernels. In Proceedings of the 13th International Conference on Learning Representations. https://iclr.cc/ virtual/2025/poster/31243 [31] Jianlin Su, Murtadha Ahmed, Yu Lu, Shengfeng Pan, Wen Bo, and Yunfeng Liu. 2024. RoFormer: Enhanced Transformer with Rotary Position Embedding. Neurocomputing 568 (2024), 127063. doi:10.1016/ j.neucom.2023.127063 [32] Stuart H. Sul, Simran Arora, Benjamin F. Spector, and Christopher Ré. 2025. ParallelKittens: Systematic and Practical Simplification of Multi-GPU AI Kernels. arXiv preprint arXiv:2511.13940 (2025). https: //arxiv.org/abs/2511.13940 [33] Philippe Tillet, H. T. Kung, and David Cox. 2019. Triton: An Intermediate Language and Compiler for Tiled Neural Network Computations. In Proceedings of the 3rd ACM SIGPLAN International Workshop on Machine Learning and Programming Languages. ACM, 10–19. doi:10.1145/3315508.3329973 [34] UXL Foundation. 2025. Intel oneAPI Collective Communications Library (oneCCL) 2021.17. https://github.com/uxlfoundation/oneCCL/ releases/tag/2021.17 [35] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Lukasz Kaiser, and Illia Polosukhin. 2017. Attention Is All You Need. In Advances in Neural Information Processing Systems, Vol. 30. 5998–6008. https://proceedings.neurips.cc/paper/ 2017/hash/3f5ee243547dee91fbd053c1c4a845aa-Abstract.html [36] Lei Wang, Yu Cheng, Yining Shi, Zhengju Tang, Zhiwen Mo, Wenhao Xie, Lingxiao Ma, Yuqing Xia, Jilong Xue, Fan Yang, and Zhi Yang. 2025. TileLang: A Composable Tiled Programming Model for AI Systems. arXiv:2504.17577 https://arxiv.org/abs/2504.17577 [37] John Yang, Carlos E. Jimenez, Alexander Wettig, Kilian Lieret, Shunyu Yao, Karthik R. Narasimhan, and Ofir Press. 2024. SWE-agent: AgentComputer Interfaces Enable Automated Software Engineering. In Advances in Neural Information Processing Systems, Vol. 37. Neural Information Processing Systems Foundation, Inc., San Diego, California, USA, 50528–50652. doi:10.52202/079017-1601 [38] Biao Zhang and Rico Sennrich. 2019. Root Mean Square Layer Normalization. In Advances in Neural Information Processing Systems, Vol. 32. 12376–12385. https://proceedings.neurips.cc/paper/2019/ hash/1e8a19426224ca89e83cef47f1e7f53b-Abstract.html

Appendix overview These appendices separate outputs of the controlled-agent studies from the co-designed SyclKittens suite. Each experiment identifies its library baseline or author-built reference, along with its measurement and aggregation protocol. 13

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 5. Measured Intel Max GPU fabric behavior. XeLink rows report direct large-transfer peer-kernel bandwidth per link. MDFI reports delivered two-stack all-reduce bus bandwidth at a 96-MB payload. The last two rows report endto-end distributed-GEMM speedups with topology-aware XeLink-once plus MDFI replication over the flat gather pipeline on 12 stacks. Measurement XeLink peer write XeLink peer read MDFI two-stack all-reduce Topology-aware gather, 61443 Topology-aware gather, 122883

stores over both fabrics. For explicit copies, the Intel Max GPU also exposes per-link copy engines; simply offloading a naive all-pairs copy does not reproduce the bandwidth of a topology-aware schedule. A.2

The interfaces available on the same hardware differ in how much schedule structure they expose. SYCL joint_matrix is a portable matrix API. It leaves operand layout and lowering to the implementation, but does not expose the complete 2-D movement and prefetch schedule used here. ESIMD provides Intel-specific explicit-SIMD control and is the substrate for larger template libraries such as XeTLA. Direct SPIR-V builtins expose individual XMX and 2-D block operations but leave their legal and profitable composition to each kernel author. SyclKittens application kernels use a compact tile programming model whose backend lowers to the same Intel SPIR-V operations. The backend defines the SIMD16/XMX contract, physical register layout, VNNI operand transform, cachebased 2-D movement, and partial-tile prefetch operations. Kernels choose prefetch ownership and distance. The controlled study therefore evaluates a stable tile programming model with an explicitly architecture-specific lowering.

Result 18.8 GB/s/link 15.7 GB/s/link 189.0 GB/s 2.15× 1.42×

Shared measurement protocol. Unless stated otherwise, single-GPU experiments activate one stack per GPU package and leave its sibling idle. Every reported configuration passes the section-specific reference check before timing. Warmup, repetition, aggregation, and byte/FLOP conventions are stated with each experiment. Displayed measurements are rounded; ratios use the unrounded values.

A

Hardware and Programming-Interface Map

This section defines the terminology and programming interfaces used in the main paper’s Background and Motivation and SyclKittens sections. The cross-vendor table aligns functions; implementation costs remain platform-specific. A node contains six GPU packages, each with two stacks connected by MDFI, the load/store stack-to-stack interconnect. XeLink connects the GPU packages all-to-all through two independent stack planes, even and odd. A.1

Intel Programming-Interface Taxonomy

B

Tile Language Mechanisms

B.1

Realization Details

Prior work describes the initial SYCL port of ThunderKittens [18]. Our baseline for this port compiled and linked on Intel Max GPU but passed only 570 of 1032 unit tests (55%). Six root causes account for every failure. With these fixes, the library used by the measured kernels passed all 1032 unit tests. We describe these fixes, the added operations and performance primitives, the prefetch API behind §3.3, and the measured SLM limits that motivate the L1 operand path.

Measured Fabric Behavior

We characterize the two communication tiers on the same 12-stack Intel Max GPU node used in the main evaluation. Direct peer kernels issue contiguous remote loads or stores through pointers in one shared SYCL context. The collective measurements use five warmups followed by 15 timed iterations and report the median; the two-stack all-reduce bandwidth uses 2(𝑁 − 1)𝑆/(𝑁 𝑡), while distributed-GEMM speedups use raw elapsed-time ratios. These conventions are reported separately because they are not interchangeable estimates of a physical link rate. The topology-aware gather sends each source stripe over XeLink once to one landing stack per remote GPU package, then replicates it to the partner stack over MDFI. This halves XeLink egress relative to every stack fetching every remote stripe. Same-plane transfers can occupy both XeLink planes concurrently; cross-plane transfers additionally consume the six MDFI links, which lowers the number of independent flows the node can sustain. Kernels can issue peer loads and

B.1.1 Critical Correctness Fixes. We identified six distinct root causes accounting for all 462 test failures (Table 8). The largest groups are shared-to-register layout, incomplete vector operations, and row reduction. The failures include NVIDIA warp-32 and CUDA-shuffle assumptions that do not survive the move to Intel Max GPU’s SIMD16 execution model, XMX packing errors, and incomplete conversion or reduction paths. Each failure mode requires a separate fix. Fix 4 also resolves an interaction with the VNNI hardware transform load (§B.1.2). Incorrect hmma81616 operand marshalling made transform-loaded operands appear corrupt, prompting a temporary switch from LoadTransformINTEL to LoadINTEL. Once contiguous operand packing and Ksplit accumulator chaining corrected hmma81616, restoring LoadTransformINTEL retained 1032/1032 correctness and removed software VNNI repacking. 14

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 6. Cross-vendor terminology used in the paper. A dash denotes that the vendor does not expose a direct counterpart with the same semantics. Concept

NVIDIA

AMD

Intel

SIMT execution group Cooperative thread group Compute unit

warp thread block / CTA streaming multiprocessor (SM) shared memory Tensor Core / MMA ldmatrix global-to-shared cp.async NVLink / NVSwitch NVLS NCCL

wavefront workgroup compute unit (CU)

subgroup workgroup Xe core

LDS Matrix Core / MFMA ds_read_b64_tr_b16 global-to-LDS pipeline Infinity Fabric / xGMI – RCCL

shared local memory (SLM) XMX / DPAS LoadTransformINTEL 2-D global-to-L1 prefetch XeLink – oneCCL

Programmer-managed local memory Matrix engine Matrix operand reformat/load Asynchronous movement used by tile kernels Inter-package scale-up network Hardware multicast collective Collective library

NVIDIA ThunderKittens Tile interface

AMD HipKittens

Intel SyclKittens

Measured BF16 GEMM

rt / st / gl · mma, load/store, reductions · warp / group / gang

shared

Capability what the silicon exposes

Method how to use it well

Policy left programmable

Tensor Core, shared memory, NVLink

MFMA, shared memory

XMX (DPAS), 2-D block load, L1, MDFI and XeLink

raw SYCL agent 57.5% of oneDNN

ldmatrix reformat, shared-memory staging, warp specialization

ds_read_b64_tr, shared-memory staging, ping-pong waves

LoadTransformINTEL, L1 cache prefetch, cooperative partial tiles

SyclKittens agent 82.1%

partial-tile ownership · prefetch distance · shapedependent schedule · stack and package placement

engineers and agents ≈96% 27-shape geometric mean

Figure 8. Capabilities, methods, and policy across three vendors. The right column reports Opus 4.8’s controlled cold-iteration BF16 GEMM on Intel Max GPU (57.5% of oneDNN with Raw and 82.1% with SyclKittens) and the separate engineer-agent suite’s roughly 96% geometric mean across 27 GEMM shapes. The tile abstractions are shared with ThunderKittens and HipKittens; only the Intel column is contributed by this work. A new hardware target requires validating the lowering and adapting workload policies to its execution and memory system. HipKittens targets a single GPU, so its column excludes scale-up interconnect. B.1.2 Performance Primitives and New Operations. Three PVC-specific primitives affect the measured kernels. move_local<T> annotates SLM pointers with __attribute__((opencl_local)) (§B.1.4). For rt_16x32 float tiles, the 2D block-load path splits a BlockCount=2 request into two BlockCount=1 loads to remain within Intel Max GPU’s 64-byte hardware load width. The __spirv_Subgroup2DBlockLoadTransformINTEL path performs the GEMM B-operand VNNI repack during the load. Combined with the corrected hmma81616, restoring this transform raises single-GPU 81923 GEMM throughput from 106 to 147 TFLOP/s, a 39% gain over the regular-load path with software repacking.

The kernels also use three groups of tile-language operations. load_part and store_part with a position parameter (the split-D strategy that processes 𝐷=128 as two 𝐷=64 passes within register budget); make_causal with partial-tile support for causal masking at non-aligned tile boundaries; and mul::reduce and min::reduce completions, since upstream ThunderKittens only exposed sum and max reductions. B.1.3 Prefetch API. SyclKittens separates cache warming from the register load that consumes a tile. prefetch_load(tile, global, next); prefetch_load_coop(tile, global, next, {row, col}); load(tile, global, next); 15

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 7. Architectural facts that determine tile movement and scheduling. NVIDIA Ampere

Intel Max GPU

Execution width Matrix width L1 capacity Local memory

32 warp-level 28–192 KB/SM up to 164 KB/SM

Preferred async path L2 capacity Matrix operand representation

global → shared 40 MB/GPU warp-distributed fragments

16 or 32 fixed SIMD16 64–512 KB/Xe core up to 448 KB/Xe core global → L1 204 MB/stack SIMD-direction linear vectors

Cache-based staging. This interface is distinct from two shared memory mechanisms often also called prefetching. SLM staging copies global data into programmer-managed local memory before a synchronized load, while SLM fetchahead combines a load, local-memory store, and namedbarrier handoff. The measured Intel Max GPU path uses neither. Hardware 2-D prefetch places data in the configurable 512 KB L1, and each consumer later loads its own register tile. This cache-centric split preserves ordinary load semantics while allowing the cooperative ownership policy in §3.3 to remove redundant warming. B.1.4 Measured Limits of SLM Operand Staging. Two measured costs keep shared local memory (SLM) off the operand-staging path on Intel Max GPU, namely narrow matrix-operand loads and compiler overhead from generic address-space qualification.

Table 8. Six distinct tile-language correctness fixes that brought SyclKittens from 55% to 100% unit-test pass. Counts in the tests-fixed column sum to 462, the gap between the 570/1032 initial-port baseline and the final 1032/1032. Fix 6 carries no count of its own because it enables the wider tile shapes the other fixes are tested on. Fix

Root cause

1. NVIDIA warp-32 register map applied to shared→register SIMD16; only 4 of 8 required packed elements written per lane 2. row reduction Within-lane stub computed column sums; rewritten as XOR butterfly across the 16 SIMD lanes 3. vector Vector ops (align copy, naive↔align inner-dim conversion, shared↔register, reductions) touched only [pack_idx][0]; seven of eight inner elements never visited 4. XMX Per-element bit_cast scattered the hmma81616 A-operand of the 8×16×16 BF16 XMX variant; replaced by contiguous reinterpret_cast and K-split accumulator chaining 5. BF16→float Stray continue skipped the type conversion load after the 2D block load on rt_16x32 6. row-map lane data_iter*2 produced lane indices up to 30 overflow for rt_16x32; wrapped modularly (unlocks 16×32 tiles)

A- and B-operand costs. Intel Max GPU has no SLM matrix load (ldmatrix equivalent). Loading an XMX B-operand (16×16 bf16) from SLM requires ∼64 individual reads vs. a single LoadTransformINTEL from HBM. Multiple SLM-based attention-kernel attempts ran an order of magnitude below the HBM-native path (the measured SLM B-operand variant reached 2.3–3.0 TFLOPS against tens of TFLOPS for the HBM path) because that 64:1 read penalty compounds with the compiler overhead below.

Tests 192

76

120

Address-space qualification matters. ThunderKittens’ st<> shared-tile type stores SLM base pointers as uintptr_t, which Intel’s IGC compiler treats as generic address space. IGC emits approximately 1,400 runtime disambiguation IR operations, adding more cost than the SLM staging removes. Typed local_accessor<T> or multi_ptr<T,local_space> fixes this, and SyclKittens’s move_local<T> of §B.1.2 applies the equivalent annotation. Address-space annotations do not remove the 64:1 readcount difference above, which motivates the direct operand path used by the measured kernels.

56

18 n/a

prefetch_load emits Intel 2-D block-prefetch instructions for the base regions implied by the register-tile shape. It warms L1 but does not write the destination registers; the later load performs the XMX-compatible 2-D block load. The prefetch is therefore a performance hint, not an asynchronous copy whose completion is required for correctness. prefetch_load_coop exposes ownership within a workgroup. One call warms a single 16-column, 8-row, two-block region at the supplied raw row and column offset. The GEMM schedule assigns disjoint offsets to its computing subgroups, so every shared operand region is warmed once while all subgroups continue to issue matrix work. The destination argument supplies tile shape and coordinate types only; it is not modified.

C

Controlled Coding-Agent Studies

The controlled studies cover single-GPU BF16 GEMM, causal grouped-query attention (GQA) [1], and 12-GPU all-gather BF16 GEMM under single-shot and independent cold-iterative feedback protocols. All three tasks use Opus, Sonnet, and Haiku in the single-shot and cold-iterative arms. Warm-iterative follow-ups use Opus on single-GPU GEMM and GQA and all eligible model–context sources in the two SyclKittens contexts on distributed GEMM. They seed feedback with held-out-correct single-shot sources, separating blank-start search from warm refinement. The task, programming context, evaluator, and budget were fixed before measurement. 16

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

C.1

Table 9. Single-GPU coding-agent protocol.

Shared Controls

Across all three tasks, each model–context pair receives three independent starts for both single-shot and cold-iterative evaluation. The evaluated models are Opus 4.8, Sonnet 5, and Haiku 4.5. Each task uses three nested programming contexts. Raw provides Intel SYCL/SPIR-V documentation but no SyclKittens. The SyclKittens context adds its headers and API documentation without a target-kernel implementation. That documentation covers the whole interface, including prefetch_load_coop and its region-offset argument, so every mechanism the co-designed kernels use is available to the agent. The example context is SyclKittens + examples, whose task-specific contents are specified below. Every condition permits direct SYCL/SPIR-V intrinsics and forbids oneAPI joint_matrix, its matrix headers, and its namespace. The agent runtime is Claude Code 2.1.214, not mini-SWEagent. Our custom harness shares ParallelKernelBench’s edit–evaluate–refine structure [5, 37]. ParallelKernelBench uses mini-SWE-agent, a repository-rooted shell, and a budget of 20 model calls. Our harness fixes condition-specific staged files and programming context and limits feedback to 20 evaluator calls. After development, the selected candidate must pass five unseen correctness seeds; no held-out input or result is exposed to the agent. Single-shot agents can read, search, and edit staged task files but have no shell, compiler, evaluator, benchmark, network, or external repository access. Cold-iterative agents additionally receive at most 20 calls to one fixed compile– correctness–performance evaluator. Warm-iterative agents begin from an eligible held-out-correct single-shot source as checkpoint 0. Replicates start independently, no source transfers between trials, and the agent never sees the five heldout seeds. Under feedback, the harness retains the fastest development-correct candidate. C.2

Item

Setting

Models and runtime

claude-opus-4-8, claude-sonnet-5, and claude-haiku-4-5-20251001; Claude Code 2.1.214 20-minute timeout; $10 cost backstop; no evaluator calls 60-minute timeout; $25 cost backstop; at most 20 evaluator calls Eligible held-out-correct single-shot source retained as checkpoint 0; at most 20 evaluator calls Development: 42; held out: 101, 202, 303, 404, 505 Relative Frobenius error < 0.02 on all five seeds Development-correct candidate with the highest co-measured tuned oneDNN ratio

Single-shot budget Cold-iterative budget Warm-iterative initialization Seeds Held-out correctness Selection under feedback

C.3

Single-GPU Complete Controlled Results

Table 10 gives the ranges omitted from Table 3 in the main paper. A dash means that no trial passed all five held-out seeds. C.4

Single-GPU Warm-Iterative Refinement

The retained Raw Opus source rises from 6.47% to 14.24%. Raw has one eligible source, so its final median and maximum in Table 3 are both 14.2%. The matched follow-up begins from the six held-out-correct Opus sources in the two SyclKittens contexts. Together with Raw, all seven Opus sources remain correct and improve under paired grading. Warm-iterative feedback improves every source. The two SyclKittens contexts reach effectively tied frontiers (87.83% in SyclKittens and 87.71% in SyclKittens + examples), while examples raise the final median from 70.46% to 86.51%. The Raw source improves from 6.47% to 14.24%.

Single-GPU Task

Every trial implements 𝑀=𝑁 =𝐾=4096 BF16 GEMM with FP32 accumulation. SyclKittens + examples adds fixed ThunderKittens and HipKittens GEMM implementations. All conditions target the same GPU, compiler toolchain, output contract, and evaluator.

C.5

Distributed All-Gather GEMM Studies

The distributed task asks each agent to implement a complete 12-GPU, row-sharded BF16 all-gather GEMM at 𝑀=𝑁 =𝐾=6144. Each GPU initially owns 512 rows of dense random 𝐴; 𝐵 is replicated; and every GPU must produce the full replicated BF16 𝐶=𝐴𝐵 with FP32 accumulation. Candidates implement the complete communication-and-compute path. The timed interval includes both components, and throughput is reported per GPU as 2𝑀𝐾𝑁 /𝑡.

Protocol. Final correctness requires all five held-out seeds to pass. Performance is throughput divided by the fixed 197.5-TFLOP/s tuned oneDNN baseline, reached through the torch.mm dispatch on Intel XPU and measured on the same GPU. The main paper reports every model separately for single-shot and cold-iterative trials and summarizes Opus warm-iterative refinement; the full multi-model distributed warm-iterative campaign appears below. Both performance summaries use only correct trials, because the median and range describe reliability among held-out-correct trials, whereas peak@3 is the maximum held-out-correct result from the three independent trials.

Task-specific context. SyclKittens + examples adds a ThunderKittens H100 BF16 GEMM, a HipKittens AMD BF16 GEMM, and a two-GPU SyclKittens all-reduce. The GEMM examples target NVIDIA and AMD GPUs. The all-reduce 17

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 10. Complete controlled BF16 GEMM results. Each cell reports held-out-correct trials / 3; median (full range) / peak@3, in percent of the fixed tuned oneDNN baseline. A singleton is shown without a range, and a dash means no held-out-correct trial. Cold-iterative replicates are independent starts. Protocol

Model

Raw SYCL

SyclKittens

+ examples

Single-shot

Opus Sonnet Haiku

1/3; 6.6 / 6.6 0/3; – 2/3; 0.7 (0.6–0.8) / 0.8

3/3; 49.1 (42.1–63.4) / 63.4 3/3; 31.9 (31.7–53.4) / 53.4 2/3; 24.0 (0.5–47.5) / 47.5

3/3; 49.4 (49.1–81.2) / 81.2 3/3; 27.5 (20.2–31.0) / 31.0 0/3; –

3/3; 49.5 (14.0–57.5) / 57.5 3/3; 5.4 (5.2–5.7) / 5.7 3/3; 1.1 (0.8–1.6) / 1.6

3/3; 81.0 (71.5–82.1) / 82.1 3/3; 84.7 (84.6–85.2) / 85.2 3/3; 0.4 (0.4–0.6) / 0.6

3/3; 82.5 (81.1–82.6) / 82.6 3/3; 71.0 (70.9–84.8) / 84.8 3/3; 3.3 (0.1–50.3) / 50.3

Opus Cold iterative Sonnet Haiku

Table 11. Matched single-GPU Opus warm-iterative refinement. Performance is relative to the co-measured oneDNN baseline; gains are percentage points. Condition

Replicate

Checkpoint 0, % vendor

Final, % vendor

Paired gain

Calls

Raw SYCL

r1

6.47

14.24

+7.77

20

SyclKittens

r1 r2 r3

49.80 66.84 44.67

87.83 70.46 56.64

+38.03 +3.62 +11.97

11 6 14

+ examples

r1 r2 r3

51.41 85.21 51.23

87.71 86.51 82.91

+36.30 +1.30 +31.68

10 5 7

demonstrates Intel peer access, but none of the three programs implements an all-gather GEMM, a 12-GPU schedule, or the topology-aware hand baseline, so the agent must compose the communication and compute path. All other treatment and tool controls follow the shared protocol above.

all-gather with PyTorch’s oneDNN GEMM. Timed end to end, it spans 128.9–144.3 TFLOPS, with a 134.9-TFLOPS median. That value is unrelated to the numerically identical SyclKittens GEMM-only value in Table 49. The matched factor sweep in Table 50 reproduces both SyclKittens operating points and attributes the difference primarily to the evaluator-versus-random-normal operand distribution, not to agent tuning. Table 12 reports only held-out-correct trials and uses the node-matched hand baseline as its denominator.

Correctness, timing, and anti-amortization. Correctness requires finite sampled output, relative L2 error below 0.03 on all five held-out seeds, and agreement on every GPU. Timing uses a 2-second power warmup, five explicit warmups, and the median of 15 complete operations. Before every warmup and scored invocation, the harness regenerates a distinct deterministic 𝐴 epoch outside the timed interval, then verifies the final epoch on all 12 outputs. A negative control that deliberately reuses gathered data fails this gate with post-timing relative L2 error 1.38, whereas valid raw-SYCL and SyclKittens candidates pass it; generated code therefore cannot remove communication by retaining gathered scratch across calls. Performance is normalized to the co-measured topologyaware hand baseline on the assigned node. This campaignlocal reference gathers 𝐴 and then computes the full 𝐶 on every GPU. It uses the evaluator’s deterministic input distribution and the frozen masked-prefetch build, rather than the random-normal cooperative-prefetch configuration in Table 49. In the single-shot campaign, the hand baseline spans 126.7–131.0 TFLOPS across six nodes, with a 128.9-TFLOPS median. The tuned vendor library path combines oneCCL

Model-stratified result. Single-shot effects depend on the model. Examples raise the Opus frontier from 36.7% to 59.4% and the Sonnet median from 0.3% to 34.2%, whereas Haiku remains below 1% in both SyclKittens contexts. Four of Raw’s five correct kernels also remain below 0.5%. No generated kernel exceeds the hand baseline or tuned oneCCL library path. Cold-iterative protocol. The feedback study keeps the task, programming contexts, evaluator, and three replicates per model and condition fixed. Each replicate begins independently from its treatment’s initial task template, with no candidate code transferred from the single-shot study or another replicate. Agents receive at most 20 calls to the fixed compilecorrectness-performance evaluator. The fastest developmentcorrect candidate is then graded on the five held-out seeds. Cold feedback changes Sonnet most, whose two SyclKittens contexts reach 50.9% and 41.1%, compared with its 25.8% 18

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 12. Complete controlled 12-GPU all-gather GEMM single-shot result. Performance is relative to the co-measured topology-aware hand baseline. peak@3 is the maximum held-out-correct result from the three trials; a dash means no trial passed all five held-out seeds. Model

Condition

Correct

Median (range), % hand

peak@3, % hand

Opus Opus Opus

Raw SYCL SyclKittens + examples

0/3 2/3 3/3

– 33.1 (29.4–36.7) 33.2 (24.3–59.4)

– 36.7 59.4

Sonnet Sonnet Sonnet

Raw SYCL SyclKittens + examples

3/3 3/3 3/3

0.3 (0.3–5.7) 0.3 (0.1–5.4) 34.2 (16.7–38.8)

5.7 5.4 38.8

Haiku Haiku Haiku

Raw SYCL SyclKittens + examples

2/3 1/3 1/3

0.3 (0.03–0.5) 0.5 0.5

0.5 0.5 0.5

Table 13. 12-GPU all-gather GEMM cold-iterative results. The median, range, and peak@3 include only held-out-correct trials and are percentages of the co-measured topology-aware hand baseline. Model

Condition

Correct

Median (range), % hand

peak@3, % hand

Opus

Raw SYCL SyclKittens + examples

3/3 3/3 3/3

42.0 (39.8–42.3) 41.7 (38.4–46.8) 42.0 (41.2–44.0)

42.3 46.8 44.0

Sonnet

Raw SYCL SyclKittens + examples

3/3 3/3 3/3

12.9 (12.5–25.8) 46.1 (32.1–50.9) 39.6 (36.3–41.1)

25.8 50.9 41.1

Haiku

Raw SYCL SyclKittens + examples

3/3 3/3 2/3

1.6 (1.1–17.2) 0.5 (0.1–0.7) 0.8 (0.1–1.5)

17.2 0.7 1.5

Raw frontier. Opus changes modestly from Raw, and Haiku remains low-throughput in both SyclKittens contexts.

Every eligible source improves. Median paired gains are 5.4 and 5.5 percentage points in SyclKittens and SyclKittens + examples, respectively, and one SyclKittens Opus source reaches 106.6% of hand (136.311 versus 127.91 TFLOP/s). Because eligibility produces different source populations, these results establish per-source refinement rather than a balanced comparison between contexts.

Distributed warm-iterative refinement. Warmiterative refinement starts from every eligible held-outcorrect single-shot source in the two SyclKittens contexts, namely six in SyclKittens and seven in SyclKittens + examples, spanning five Opus, six Sonnet, and two Haiku trials. All 13 sources remain held-out correct and improve. Table 14 reports each paired result relative to the node-matched hand baseline. The 106.6% source improves the distributed algorithm rather than the dense GEMM microkernel. It computes each GPU’s 𝐶 row shard from the local 𝐴 shard and replicated 𝐵, then performs a topology-aware all-gather of 𝐶. Because this task has 𝑁 =𝐾, gathering 𝐶 (𝑀×𝑁 ) moves the same BF16 tensor volume as gathering 𝐴 (𝑀×𝐾), but the transformation reduces aggregate GEMM work from 12 full products to one full product partitioned across 12 GPUs. It remains a valid result, because the fixed harness regenerates 𝐴 before every invocation, waits both queue sets, and verifies the full replicated 𝐶 on every GPU for all five held-out seeds.

Engineer-agent co-design. The controlled sources above remain frozen for attribution. A matched Box–Muller comparison found that the SyclKittens + examples kernel reached 131.616 TFLOP/s, exceeding the 128.545 of our initial refinement of a SyclKittens kernel. We restarted from the unchanged SyclKittens + examples kernel and retained it as the incumbent under a monotonic promotion rule. The selected refinement retains its DPAS core and cross-package route, then applies cyclic MDFI band order with a two/four split over the compute/copy queues. It reaches 125.724 TFLOP/s versus 123.839 for the incumbent in a new four-node randomnormal study, with a 1.007× node-first paired ratio that sits within run-to-run variation (Table 48). In this workflow, agents expand the algorithm family, engineers retain the strongest validated source as the incumbent, and joint tuning 19

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 14. 12-GPU all-gather GEMM warm-iterative refinement. Each eligible source remains correct. Gains are paired percentage-point changes from the fixed checkpoint 0 to the final selected checkpoint. Condition

Model

Replicate

Checkpoint 0, % hand

Final, % hand

Paired gain

Calls

SyclKittens

Opus Sonnet Haiku Sonnet Opus Sonnet

r1 r1 r1 r2 r3 r3

29.49 0.28 0.47 5.43 35.75 0.15

106.57 6.19 0.72 8.86 40.64 36.72

+77.08 +5.91 +0.24 +3.43 +4.89 +36.58

15 18 20 11 13 15

+ examples

Opus Sonnet Opus Sonnet Opus Sonnet Haiku

r1 r1 r2 r2 r3 r3 r3

24.13 35.18 59.96 39.76 33.42 16.60 0.49

38.51 66.57 61.51 41.17 38.90 38.85 0.70

+14.38 +31.39 +1.55 +1.41 +5.48 +22.24 +0.21

16 10 10 14 9 13 20

is accepted only under a common correctness and performance gate. C.6

two HipKittens causal-GQA kernels for AMD GPUs rather than an Intel implementation. Each of the 27 cold-iterative trials independently starts from a no-op condition stub and receives at most 20 evaluator calls. The 27 single-shot trials use the same frozen inputs without evaluator feedback. Correctness requires finite output, relative L2 error below 0.05, and elementwise error rate below 0.01 on all five held-out seeds. Performance is normalized by the native PyTorch XPU SDPA baseline co-measured in each trial, with enable_gqa=True and is_causal=True. The baseline process does not import IPEX. The observed cold-iterative baselines range from 73.0 to 85.3 TFLOP/s. The kernel-suite sweep in §E imports IPEX and measures its grouped-head fallback.

Single-GPU Mechanism Inspection

The single-shot Opus frontier in SyclKittens + examples is replicate r2. Its selected source declares a 32×32 register base shape, presents B in the VNNI-compatible column layout, sets the K-tile prefetch distance to three, uses 2-D loads, and reaches DPAS through mma_AB. It reaches 160.4 TFLOP/s, or 81.2% of the fixed oneDNN baseline. The cold-iterative frontiers in the remaining columns also use mma_AB, 2-D loads, and VNNI-compatible B, but select different schedules. Raw cold frontiers are 57.5% for Opus, 5.7% for Sonnet, and 1.6% for Haiku. Direct SPIR-V/XMX access is available in every condition. The generated sources differ in their tiling and operand-delivery implementations. C.7

Context effects differ by model. Examples lead for Opus at 79.5%, while Raw leads for Sonnet at 18.6%, followed by examples at 8.8% and SyclKittens at 4.6%. Haiku remains between 0.7% and 1.3% across all contexts. The best Sonnet cell falls below the weakest Opus cell at 41.6%. All 27 kernels are correct, so correctness alone would hide this model effect. A pooled context average would also hide the reversal between Opus and Sonnet and is therefore not reported. The single-shot arms contain 27 infrastructure-valid cells and 23 held-out-correct endpoints. The only incorrect Sonnet or Haiku endpoint is Haiku SyclKittens replicate 1, whose generated candidate fails the build gate. This is a correctness outcome rather than a harness failure.

Three-Model Grouped-Query Attention

We evaluate Opus, Sonnet, and Haiku on causal GQA at 𝐵=1, 𝐻 =32, 𝐻𝐾𝑉 =8, 𝑁 =2048, and 𝐷=128 under a corrected equal-ownership protocol. The Sonnet and Haiku extension changes only the model alias. Its five frozen task inputs are byte-identical to the Opus campaign, including the raw-pointer harness with SHA-256 bdd6bd8ff88d83c2 9b0f6a4655268407906a2fd63e9179cce173935db143b9c. Raw, SyclKittens, and SyclKittens + examples use that harness and the same candidate function signature. Every candidate owns its launch geometry, work-group and subgroup mapping, and submissions. ESIMD is forbidden in every condition because it is a separate Intel-specific SIMD programming model that neither raw SYCL nor SyclKittens builds on. Direct SPIR-V subgroup matrix and 2-D block-I/O intrinsics remain available. Raw provides Intel SYCL/SPIR-V documentation, SyclKittens adds its headers, and SyclKittens + examples additionally provides

Opus warm refinement. The endpoint-seeded warm campaign is separate from the balanced cold grid. It refines the six held-out-correct Opus single-shot sources and does not include Sonnet or Haiku. All six remain correct and improve. Paired gains range from 9.5 to 81.7 percentage points, and the strongest endpoint reaches 71.0 TFLOP/s or 83.4% of its co-measured baseline. Unequal eligibility makes these 20

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 15. Mechanisms in the single-shot Opus frontier in SyclKittens + examples and the fastest single-GPU cold kernel for each model. Percentages use the fixed tuned oneDNN baseline.

Technique

Single-shot Cold Opus + examples Opus + examples

Cold Sonnet SK

Performance DPAS through mma_AB 2-D loads and VNNI 𝐵 Register base shape K-tile prefetch distance K-loop unrolling

160.4 TF, 81.2% yes yes 32×32 3 none

168.3 TF, 85.2% 99.3 TF, 50.3% yes yes yes yes 32×32 16×16 4 none none 16×

163.1 TF, 82.6% yes yes 32×32 2 none

Cold Haiku + examples

Table 16. Complete equal-launch cold-iterative GQA results. Replicate columns report held-out candidate TFLOP/s and the parenthesized percentage of that trial’s co-measured native PyTorch XPU SDPA baseline, without IPEX imported. The final column reports median / peak@3 among correct trials. All 27 trials are correct. Model

Context

Correct

Opus 4.8

Raw SYCL SyclKittens + examples

3/3 3/3 3/3

42.30 (57.7) 13.81 (18.6) 26.38 (36.1) 24.84 (29.7) 34.41 (41.6) 15.89 (18.6) 59.95 (79.5) 38.04 (50.4) 32.77 (43.2)

r1

36.1 / 57.7 29.7 / 41.6 50.4 / 79.5

Sonnet 5

Raw SYCL SyclKittens + examples

3/3 3/3 3/3

14.47 (17.6) 15.36 (18.6) 14.29 (17.4) 3.48 (4.6) 1.38 (1.8) 3.33 (4.4) 6.82 (8.8) 2.96 (3.9) 3.41 (4.4)

17.6 / 18.6 4.4 / 4.6 4.4 / 8.8

Raw SYCL Haiku 4.5 SyclKittens + examples

3/3 3/3 3/3

0.56 (0.7) 0.29 (0.4) 1.00 (1.3)

r2

0.45 (0.6) 0.73 (1.0) 0.40 (0.5)

r3 Median / peak, % SDPA

0.36 (0.5) 0.21 (0.3) 0.24 (0.3)

0.6 / 0.7 0.4 / 1.0 0.5 / 1.3

Table 17. Complete equal-launch single-shot GQA results. Replicate columns report percentages of the co-measured native PyTorch XPU SDPA baseline, without IPEX imported. A dash denotes a trial without an admissible held-out-correct endpoint. The final column reports peak@3 among correct trials. Model

Context

Correct

r1

r2

r3

peak@3, % SDPA

Opus 4.8

Raw SYCL SyclKittens + examples

1/3 2/3 3/3

– – 1.3

– 1.7 8.3

1.4 11.3 12.7

1.4 11.3 12.7

Sonnet 5

Raw SYCL SyclKittens + examples

3/3 3/3 3/3

0.1 2.4 2.4

2.9 1.2 3.2

2.6 3.0 0.4

2.9 3.0 3.2

Haiku 4.5

Raw SYCL SyclKittens + examples

3/3 2/3 3/3

0.2 – 0.2

0.01 0.2 0.01

0.01 0.03 0.3

0.2 0.2 0.3

Table 18. Opus endpoint-seeded warm GQA results. Final values are percentages of the co-measured native PyTorch XPU SDPA baseline, without IPEX imported. Context

Eligible

Final values

Median / maximum

Raw SYCL SyclKittens + examples

1/1 2/2 3/3

50.1 83.4, 22.0 28.8, 48.2, 22.1

50.1 / 50.1 52.7 / 83.4 28.8 / 48.2

Opus mechanisms and variance. All nine selected Opus cold-iterative kernels reach XMX through direct SPIR-V subgroup DPAS, and eight also use direct 2-D block I/O. The strongest Raw kernel reaches 42.3 TFLOP/s with transformed 2-D loads and prefetch, while the strongest example-assisted kernel reaches 59.9 TFLOP/s with direct 2-D loads, transposed loads, and prefetch. The selected Opus kernels from SyclKittens and SyclKittens + examples use direct intrinsics rather than SyclKittens register-tile or MMA APIs. The treatment therefore changes the schedules Opus discovers even when the final implementations use the same low-level mechanism family.

paired refinements rather than a balanced context comparison. 21

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Key techniques.

The HipKittens examples expose tiled causal-attention algebra and pipelining, but the generated kernels are Intel implementations. The strongest Opus result translates that structure into PVC DPAS and cache-prefetch operations, while the other two example-assisted replicates do not reproduce the same schedule. Only the strongest of the three example-assisted Opus starts exceeds Raw’s cold-iteration peak. C.8

• Register-resident DPAS tiling. A 32×64 fp32 accumulator per subgroup (falling back to 16×16 register tiles for sub-32 output tiles), bf16 operand tiles, and up to 32 subgroups per work-group keep the systolic array fed. Profiling shows the kernel is DPAS-bound (XVE about 94.5% array-active, L3 bandwidth about 1.3% of peak), so large, well-populated work-groups matter more than additional cache optimization. • VNNI 2D block loads. Hardware block-load instructions with the VNNI operand transform deliver bf16 operands to the DPAS unit in their native packed layout, avoiding a local-memory staging round trip. • Masked prefetch baseline. Before de-duplication, full-tile A prefetches are issued when sg_id % 5 != 0, and full-tile B prefetches are issued when sg_id & 8 != 0. This selective schedule reaches 88.8% of oneDNN. Loop unrolling and snake traversal then raise the geometric mean to 91.0%. • Snake-curve group swizzle. Work-group tiles are traversed in an L2-friendly super-stack order, with the swizzle width fixed per shape (4 at 20483 , 16 at 30723 /81923 , 4 at 122883 /163843 , 8 for other square shapes, 8 for all rectangular shapes) from a full snakewidth sweep. A guard falls back to linear indexing for non-aligned grids, and grid padding with an earlyreturn bounds check (returning when a work-group’s row or column tile index exceeds the grid) prevents tile skips and output corruption on non-divisible snake grids. • Collaborative de-duplicated prefetch. The partialtile policy assigns every subgroup a disjoint slice of the shared 256×32 A and 32×256 B workgroup panels. Each panel byte is warmed once per K-step, rather than once per consuming subgroup; the prologue uses the same assignment. • Shape-adaptive lookahead. De-duplication reduces explicit prefetch traffic by about 3.6× and enables deeper lookahead, namely distance 2 for all squares and mostly distance 3 for rectangles, with measured distance-1 and distance-2 exceptions. The retained masked baseline instead uses distance 1 for 𝐾 > 4096 and distance 2 otherwise. • K-loop unroll (32-way). Exposing 32 K-iterations to the compiler overlaps loads with MMA, measured +7% at 81923 and +16% at 102403 over a 16-way unroll. • bf16-output store. The fp32 accumulator is downconverted in-register by a 2-wide round-to-nearesteven convertor, then written with a 2-byte block store, which halves the store traffic of the earlier fp32 path and removes the fp32 global round trip. This matches the oneDNN output dtype; the co-measured bf16/fp32 time ratio is 0.95 on a healthy node.

Data Quality and Protocol Deviations

The deviations below do not alter the planned balance of the single-shot or cold-iterative comparisons. Warm-iterative analyses remain paired because eligibility determines their source populations. The accompanying artifact retains prompts, programming contexts, generated sources, agent interactions, development evaluations, selected candidates, and held-out logs for independent inspection. C.9

AI Use in SyclKittens development

The autonomous controlled study is distinct from the engineer-agent suite. Intel hardware experts helped design the programming model. Engineers and coding agents collaborated on implementation, call-site integration, variant generation, and tuning. Engineers owned the experiments and correctness and performance acceptance. Generated code entered the suite only after compilation, operation-specific tests, interleaved measurement, and human review, as described by the main-paper AI disclosure.

D

BF16 GEMM

Across the 27-shape evaluation suite, the SyclKittens bf16 GEMM with collaborative de-duplicated prefetch achieves 0.963× oneDNN in geometric-mean throughput. Measurements use bf16 output, shape-adaptive launch policies, within-trial co-measurement, and six-node aggregation. D.1

Configuration

GEMM FLOPs are 2𝑀𝐾𝑁 for both implementations, so their TFLOPS figures are directly comparable. D.2

Kernel Implementation and Key Techniques

The GEMM kernel is a single compile-time-templated translation unit built into a callable per-shape module; each shape in the 27-shape evaluation suite is compiled separately. Data flow (per output tile). Each subgroup owns a BM × BN (32×64) fp32 accumulator register tile c_reg and streams the 𝐾 dimension in BK=32 chunks. For each K-step it loads a BM × BK bf16 𝐴 tile (row-major) and a BK × BN bf16 𝐵 tile (column-major) via 2D block loads, with the VNNI transform applied to B. It issues the DPAS contraction into c_reg and prefetches future 𝐴/𝐵 panels at the configured lookahead distance. After the K loop, c_reg is scaled and stored. 22

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 19. Protocol deviations and their effect on the reported studies. Deviation

Resolution

Effect

Earlier definitions of the two SyclKittens contexts Off-workspace edit in a single-GPU single-shot trial Missing trusted grading record in a distributed single-shot trial Two distributed pilot protocols

Excluded before analysis; the reported dataset uses the final context definitions. Discarded before evaluation; repeated from the same initial task and programming context. Repeated from the unchanged initial task before analysis.

All reported comparisons use the same final definitions. No execution feedback entered evaluation; the replacement completed the planned trial. Only the repeat contributes to the reported result.

Datasets that omitted treatment-specific guides or reused one input epoch were excluded.

These datasets do not enter the results; the reported protocol includes both guides and fresh inputs.

Table 20. GEMM configuration. Parameter

const int coop_a_row_offset = coop_h_idx * pref_a_row_offset; const int coop_b_row_offset = coop_b_row_idx * pref_b_row_offset; const int coop_b_col_offset = coop_b_col_idx * 32;

Value

Problem (square)

𝑀=𝐾=𝑁 , with 𝑁 ∈ {1024, 2048, 3072, 4096, 5120, 6144, 7168, 8192, 9216, 10240, 12288, 14336, 16384} Problem (rect- 14 shapes drawn from Llama-3/BERT-style angular) inference workloads (QKV, output projection, FFN up/down, decode with 𝑀 ∈ {256, 512}) Dtype bf16 𝐴/𝐵/𝐶, fp32 DPAS accumulation Hardware Intel Max GPU, one active stack per GPU package, sibling stack idle Nodes Six nodes; results are aggregated by crossnode median Baselines oneDNN, reached through the torch.mm dispatch used on this platform Tiling BM=32, BN=64, BK=32 with a 32×32 register tile and 32-way K-loop unrolling; per-shape snake width, shape-adaptive root-sync launch, and collaborative de-duplicated prefetch (distance 2 for squares; shape-selected distances 1–3 for rectangles)

int prefetch_a = GEMM_PF_A; int prefetch_b = GEMM_PF_B; // De-duplicated prologue: use the same disjoint slices as steady state. #pragma unroll for (int k = 0; k < prefetch_a; k++) prefetch_load_coop(a_reg, g_a, {0, 0, wg_m, k}, {coop_a_row_offset, 0}); #pragma unroll for (int k = 0; k < prefetch_b; k++) prefetch_load_coop(b_reg, g_b, {0, 0, k, wg_n}, {coop_b_row_offset, coop_b_col_offset}); #pragma unroll GEMM_UNROLL for (int k = 0; k < GK / BK; k++, prefetch_a++, prefetch_b++) { load(a_reg, g_a, {0, 0, wg_m, k}); load(b_reg, g_b, {0, 0, k, wg_n}); if (prefetch_a < GK / BK) prefetch_load_coop(a_reg, g_a, {0, 0, wg_m, prefetch_a}, {coop_a_row_offset, 0}); if (prefetch_b < GK / BK) prefetch_load_coop(b_reg, g_b, {0, 0, prefetch_b, wg_n}, {coop_b_row_offset, coop_b_col_offset}); mma_AB(c_reg, a_reg, b_reg, c_reg); }

D.2.1 Collaborative de-duplicated inner loop. The partial-tile policy keeps the same demand loads and DPAS contraction but replaces full-tile warming with partial-tile ownership. The 32 subgroups form an 8×4 grid. For A, the four subgroup columns partition each 32×32 tile into row slices; for B, the subgroup-row indices partition its K extent and select one 32-column half. These offsets cover the distinct workgroup panel exactly once. Unlike the earlier full-tile cooperative-prefetch experiment, this path has no split barrier. // Collaborative partial-tile ownership arithmetic. const int coop_h_idx = sg_id % (256 / BN); const int coop_v_idx = sg_id / (256 / BN); const int pref_a_row_offset = BM / (256 / BN); const int pref_b_sg_row_size = 256 / (BM * 2); const int pref_b_row_offset = BK / (256 / BN); const int coop_b_row_idx = coop_v_idx % pref_b_sg_row_size; const int coop_b_col_idx = coop_v_idx / pref_b_sg_row_size;

The loop needs no producer-consumer split, SLM staging, multi-stage buffer rotation, or barrier. ThunderKittens instead splits its pipeline across producer and consumer warpgroups and waits on a semaphore at every K step, while the HipKittens loop issues two workgroup barriers per matrix instruction. 23

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 21. Matched cumulative prefetch tail. The 81923 row reports cross-node medians; the 27-shape row reports geometric means.

Table 22. Shape-specific-PF square (SyclKittens/oneDNN geometric mean of 0.970×). Size

SyclKittens/oneDNN

Suite

Snake

De-duplicated

Shape-specific PF

10243

81923 All 27

94.3% 91.0%

95.8% 96.2%

95.9% 96.3%

20483 30723 40963 51203 61443 71683 81923 92163 102403 122883 143363 163843

1.165× 0.881× 0.999× 0.962× 0.940× 0.947× 0.944× 0.959× 0.989× 0.948× 0.96× 0.972× 0.964×

Each prefetch_load_coop call lowers its assigned offset to one Intel 2-D block-prefetch instruction with a 16×8 element block and count two. Across the workgroup, the calls warm 32,768 bytes per K-step, equal to the distinct A+B footprint, versus approximately 118,784 bytes for the masked full-tile path. The same assignment is used in the prologue; de-duplicated priming adds a 0.8% geometric-mean gain on small and medium squares. D.3

Table 23. Shape-specific-PF rectangular (SyclKittens/oneDNN geometric mean of 0.953×).

Methodology

Measurements use the shared single-GPU isolation protocol. oneDNN and SyclKittens are loaded into one process and timed in an interleaved round robin. Each trial uses 10 warmup and 30 measured iterations and rotates the implementation that starts first. Interleaving reduces clock drift in the SyclKittens/oneDNN ratio. Within a node, each value is the per-shape median of per-trial ratios. The reported ratios are medians of the six node medians. Before a shape enters the results tables, its output is checked against torch.mm by relative Frobenius norm. All measured configurations pass, with maximum error 1.663 × 10−3 . D.4

De-Duplicated-Prefetch Performance

Matched cumulative tail. The final three cumulative stages share shapes, output type, launch policy, correctness gate, and six-node co-measurement. They isolate snake traversal, partial-tile ownership de-duplication, and shapespecific prefetch distance. Across all shapes, de-duplication improves the geometric mean by 1.056×. Shape-specific prefetch distances produce the final 96.3% configuration. Every value is a six-node cross-node median from comeasured trials using the loop in §D.2.1. Ratios use measured oneDNN throughput as the denominator. SyclKittens exceeds oneDNN on three of the 27 shapes. D.5

Shape (𝑀×𝐾×𝑁 )

Workload

2048×4096×12288 4096×4096×12288 8192×4096×12288 2048×4096×4096 8192×4096×4096 2048×4096×11008 4096×4096×11008 2048×11008×4096 4096×11008×4096 2048×1024×4096 2048×4096×1024 256×4096×4096 512×4096×4096 512×4096×11008

QKV 𝑀=2K QKV 𝑀=4K QKV 𝑀=8K Out 𝑀=2K Out 𝑀=8K FFN up 𝑀=2K FFN up 𝑀=4K FFN down 𝑀=2K FFN down 𝑀=4K BERT FFN up BERT FFN down Inference 𝑀=256 Inference 𝑀=512 Inference FFN

D.6

GEMM

GEMM

SyclKittens/oneDNN 0.977× 0.960× 0.953× 0.924× 0.926× 0.925× 0.935× 0.884× 0.941× 0.909× 1.042× 0.962× 1.039× 0.976×

Triton-XPU Comparison

One sweep co-measures oneDNN, the final cooperative SyclKittens revision, and Intel Triton-XPU 3.6 over the same 27 BF16 input/output shapes. For each shape on each node, we independently autotune Triton and record its best configuration for comparison. Each trial rotates their execution order, with 10 warmup iterations and 30 timed iterations per implementation; host timing includes an explicit XPU completion barrier. Each implementation is checked against the same fp32 torch.mm reference. Ratios are reduced first by the median across 10 trials within each node and then by the median across six nodes. The final cooperative SyclKittens reaches 1.087× geometricmean and 1.030× median throughput over Triton. The crossnode median exceeds Triton on all 27 shapes; the smallest margin is 1.0125×. Table 24 shows that the advantage holds for both square and rectangular workloads.

Correctness

All 27 shapes pass the relative-Frobenius-norm gate ∥𝐶 − ref ∥ 𝐹 /∥ref ∥ 𝐹 < 0.02. Every shape reports a uniform relative error of 1.66 × 10−3 , the expected bf16 rounding floor for these accumulation depths, with no per-shape outlier. The bf16-output store retains the fp32 output path’s accuracy. 24

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 24. Direct final-cooperative SyclKittens/Triton-XPU comparison. Geometric means and medians are over crossnode-median per-shape throughput ratios. Shape set

Shapes

Geomean

Median

SyclKittens > Triton

27 13 14

1.087× 1.095× 1.079×

1.030× 1.031× 1.019×

27/27 13/13 14/14

All Square Rectangular

Table 26. Attention-forward configuration (MHA and GQA). Parameter

Value

Problem (MHA)

𝐵 = 4, 𝐻 = 𝐻𝑘𝑣 = 16, 𝐷 = 128 (multi-head, no grouping) 𝐵 = 1, 𝐻 = 32, 𝐻𝑘𝑣 = 8, 𝐷 = 128 (group size 4; full-sequence queries) 𝑁 ∈ {512, 1024, 2048, 4096, 6144, 8192}

Problem (GQA) Sequence lengths Dtype Hardware Peak reference

Table 25. Split-K delivered performance, 𝑀=𝑁 =512, K sweep; throughput columns are TFLOPS.

Baselines 𝐾

oneDNN

Default

Auto (𝑆=8)

Auto/oneDNN

Auto/default

2048 4096 8192 16384 32768

46.8 65.0 100.5 130.2 169.9

20.8 22.4 23.5 24.0 23.8

50.0 78.4 110.5 143.3 163.7

1.07× 1.21× 1.10× 1.10× 0.96×

2.4× 3.5× 4.7× 6.0× 6.9×

Kernel (MHA) Kernel (GQA)

GPU. Across non-causal and causal modes, its geometricmean speedups over IPEX SDPA are 1.09–1.21× for MHA and 1.70–4.0× for GQA.

Across the six nodes, the 27-shape SyclKittens/Triton geometric mean ranges from 1.084× to 1.095×, so the advantage does not depend on which node measures it. The maximum SyclKittens and Triton relative Frobenius errors against the reference are both 1.663 × 10−3 , so both implementations clear the same accuracy bar. D.7

E.1

Configuration

MHA is characterized at 𝐵 = 4 (throughput regime); the GQA performance tables in Sec. E.4 measure batch-one fullsequence queries at 𝐵 = 1. The GQA runtime selection policy is a 2-D function of (𝑁 , 𝐵) (Sec. E.5), retuned across 𝐵 ∈ {1, 2, 4, 8} for full-sequence attention. GQA is compared against SDPA with IPEX imported and enable_gqa=True. This suite measurement uses its grouped-head fallback, which repeats KV heads, rather than the separately measured pre-expanded-KV fast path. The controlled agent study uses native PyTorch XPU SDPA without importing IPEX. There is no Triton grouped-attention baseline.

Split-K / Stream-K Design

The SyclKittens atomic-combine split-K path targets occupancy-starved decode, QKV-projection, and MoEexpert shapes with small 𝑀, 𝑁 , large 𝐾, and an L2-resident output. Eight compile-time K planes reuse the base mainloop and atomically combine fp32 accumulators before one bf16 output round. D.7.1 Delivered performance (𝑀=𝑁 =512, K sweep). The oneDNN baseline is torch.mm (bf16 in, bf16 out); “auto” denotes the auto-selected split (𝑆=8). All results pass relativeFrobenius < 2 × 10−6 . The default kernel remains at 21–24 TFLOPS because 128 workgroups under-fill the machine while the serial mainloop grows with 𝐾. Split-K is 2.4–6.9× faster than that default and exceeds oneDNN by up to 1.21× for 𝐾 ≤ 16384. At 𝐾 = 32768, oneDNN also uses K-parallelism and reaches 169.9 TFLOPS versus 163.7 TFLOPS for split-K.

E.2

Kernel Implementation and Key Techniques

The MHA implementation is a single-pass, register-resident FlashAttention-style kernel [6] built on ThunderKittens tile primitives and specialized for the Intel Max GPU (Xe) DPAS matrix engine. It is templated on head dimension, workgroup size, query-tile length, and causality, so one compiled module holds every tuned variant and runtime dispatch selects the measured work-group size and query-tile length (co-loading separately compiled variants is not possible, since their generated kernel symbols collide).

D.7.2 Correctness. The atomic-add epilogue is bit-exact against the plain-store epilogue in all 10 unit-test cases (relative-Frobenius error zero).

E

bf16 inputs/outputs, fp32 accumulation Intel Max GPU, one GPU 157.5 TFLOPS bf16 (oneDNN 512×4096×11008); MFU quoted against this IPEX SDPA (Flash backend), Triton Flash (fp16 operand) Fence-free adaptive MHA dispatch Fence-free adaptive dispatch with a 2-D (𝑁 , 𝐵) GQA branch

Data flow (per query tile). Each subgroup owns a seq_q × 𝐷 query tile and streams the KV sequence in kv_stride = 32 chunks. For every KV tile, the kernel computes 𝑆 = 𝑄𝐾 ⊤ with mma<N,T>, accumulating inner_loop = 𝐷/32 = 4 DPAS contractions into an 𝑆 [seq_q × 32] register tile. It then applies the optional causal mask, updates the online softmax’s running maximum and normalizer, rescales

Attention Forward (MHA and GQA)

The fence-free adaptive forward-attention kernel reaches 71.7 TFLOPS for MHA and 75.7 TFLOPS for GQA on one 25

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 27. MHA forward, non-causal (𝐵 = 4, 𝐻 = 𝐻𝑘𝑣 = 16, 𝐷 = 128). geomean 1.09× over SDPA; peak 71.7 TFLOPS (45% MFU) at 𝑁 = 6144. All points average two runs except 𝑁 = 6144, which uses one valid run.

the output accumulator, and accumulates 𝑂 += 𝑃𝑉 with mma<N,N>. After the KV loop, 𝑂 is normalized by the running sum and stored. Key techniques. • Online softmax and scaling. The 𝑁 × 𝑁 score matrix is never materialized. The running max (max_vec), normalizer (norm_vec), and output accumulator remain in registers, with the accumulator rescaled by √ 𝑒 (𝑚old −𝑚new ) at each KV step. Folding 1/ 𝐷 with log2 𝑒 = 1.44269504 maps 𝑝 = exp2(𝑆 · scale − 𝑚) directly to the hardware exp2 unit. • Register reuse and pipelining. The seq_q×𝐷 query tile is loaded once before the KV loop, avoiding a reload for every KV tile. Software prefetch issues the next K/V loads while DPAS consumes the current tile. Hardware dependency tokens order load-to-DPAS and DPAS-to-consumer hazards without software scoreboard drains; a 4-seed × 20-repeat ablation measured up to a ∼15% cost from those drains at large-𝑁 causal, with no correctness benefit. • DPAS accumulator and operand layouts. For 𝐷 > 64, split-D accumulation carries the output in two 𝐷/2-wide accumulators with independent V loads and DPAS, preserving the 32 × 32 DPAS shape and register budget. The block loads map row-major A, transposed B for mma<N,T>, and VNNI-packed column-major B for mma<N,N>; a mismatched role corrupts 𝑆 or 𝑂. The transposed-B load lets the matrix engine consume 𝐾 ⊤ directly without explicit lane-shuffle repacking. • Query tiling and causal scheduling. Setting seq_q = 8 halves rows per subgroup and doubles concurrent query tiles; it is selected at small 𝑁 , with the crossover to seq_q = 16 determined by the MHA or GQA dispatch table below. The causal loop omits fully masked upper-triangle KV tiles and masks only the diagonal tile across kv_stride/seq_q sub-blocks, namely two for seq_q = 16 and four for seq_q = 8. Reverse linear query ordering dispatches the heaviest causal blocks first to balance workgroups. • Head mapping and grid residency. Query head ℎ maps to KV head ⌊ℎ/(𝐻 /𝐻𝑘𝑣 )⌋; MHA is the groupsize-one specialization, with separate tuned dispatch parameters (Sec. E.5). Avoiding forced root synchronization lets the full grid run concurrently instead of limiting residency to approximately one workgroup per compute core, improving large-𝑁 throughput by up to ∼16×. • Split-K forward. When KV depth per GPU leaves execution resources idle, workers process separate KV partitions and write partial 𝑂/max/norm states. A small log-sum-exp kernel merges them, improving performance by 0.6–2.9% at 𝑁 ≥ 4096.

𝑁

SyclKittens (T)

IPEX (T)

SyclKittens/IPEX

Triton (T)

SyclKittens/Triton

MFU

512 1024 2048 4096 6144 8192

49.3 59.2 67.2 69.0 71.7 70.0

48.8 57.1 59.6 61.3 62.6 62.7

1.01× 1.04× 1.13× 1.12× 1.14× 1.11×

45.7 64.3 68.5 64.8 69.1 64.6

1.08× 0.92× 0.98× 1.06× 1.04× 1.08×

31% 38% 43% 44% 45% 44%

Table 28. MHA forward, causal (𝐵 = 4, 𝐻 = 𝐻𝑘𝑣 = 16, 𝐷 = 128). Geomean 1.21× over SDPA; peak 68.7 TFLOPS (44% MFU) at 𝑁 = 8192; 1.42× at 𝑁 = 512. All points average two runs except 𝑁 = 6144, which uses one valid run. 𝑁

SyclKittens (T)

IPEX (T)

SyclKittens/IPEX

Triton (T)

SyclKittens/Triton

MFU

512 1024 2048 4096 6144 8192

46.9 56.5 61.3 66.4 68.5 68.7

32.9 48.8 52.5 56.6 58.0 58.4

1.42× 1.16× 1.17× 1.17× 1.18× 1.18×

25.1 42.5 52.4 55.9 58.3 60.7

1.87× 1.33× 1.17× 1.19× 1.18× 1.13×

30% 36% 39% 42% 43% 44%

E.3

Methodology

Each (𝑁 , config) runs in a separate subprocess pinned to one GPU, with at most one heavy task per GPU. Baseline, SyclKittens, and Triton are interleaved within each trial; GQA compares only the baseline and SyclKittens. Each of 10 independent trials uses 10 warmup and 30 measured iterations. TFLOPS is the trial mean, with standard deviation shown where applicable; GQA values are additionally averaged over three independent runs. Non-causal forward attention counts 4 𝐵𝐻 𝑁 2 𝐷 FLOPs. The causal count uses the same coefficient over lowertriangular work, so each path is directly comparable with its own baseline. Every configuration is checked against a PyTorch reference before inclusion (Sec. E.6). E.4

Delivered Performance

Numbers use the validated adaptive runtime policy rather than separately compiled per-shape kernels. E.4.1 MHA (adaptive kernel, 𝐵 = 4, 𝐻 = 16, 𝐷 = 128). The model output contract and SDPA baseline require bf16 output, so the reported measurements include the fp32-tobf16 epilogue conversion. At 50% occupancy, this conversion spills the seq_q = 16 epilogue; omitting it and returning fp32 is 12% faster at small 𝑁 and 4% at large 𝑁 . Store-payload reshaping and a compiler register-pressure cap both increased runtime. At 𝑁 = 512, the selected causal variant reaches 1.42× IPEX and 1.87× Triton. Its seq_q = 8 tile doubles concurrent query tiles, skips fully masked upper-triangle KV blocks, and applies masking only to the diagonal block, reducing executed 𝑄𝐾 ⊤ , softmax, and 𝑃𝑉 work. 26

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 29. GQA forward, non-causal (𝐵 = 1, 𝐻 = 32, 𝐻𝑘𝑣 = 8, 𝐷 = 128). geomean 1.70× over SDPA; peak 75.7 TFLOPS (48% MFU) at 𝑁 = 6144. 𝑁

SyclKittens (T)

IPEX (T)

SyclKittens/IPEX

MFU

512 1024 2048 4096 6144 8192

46.3 ± 0.2 52.5 ± 0.0 67.8 ± 1.0 72.2 ± 0.6 75.7 ± 0.3 72.3 ± 0.1

27.2 45.9 42.7 41.8 38.4 31.1

1.70× 1.14× 1.59× 1.73× 1.97× 2.33×

29% 33% 43% 46% 48% 46%

Table 32. GQA 2-D (𝑁 , 𝐵) autotune selections, causal. Cells give wg/seq_q (median latency, µs).

SyclKittens (T)

IPEX (T)

SyclKittens/IPEX

MFU

512 1024 2048 4096 6144 8192

37.1 ± 0.2 66.0 ± 0.1 69.8 ± 1.0 69.8 ± 0.8 71.6 ± 0.3 72.1 ± 0.3

9.5 17.6 19.1 19.2 17.4 13.7

3.91× 3.74× 3.65× 3.64× 4.11× 5.27×

24% 42% 44% 44% 45% 46%

1

2

4

8

256 512 1024 2048

256/8 (48) 64/8 (72) 128/8 (145) 128/16 (473)

128/8 (53) 128/8 (96) 64/16 (266) 256/16 (1045)

256/8 (71) 64/16 (159) 128/16 (569) 256/16 (2100)

128/8 (111) 128/16 (332) 256/16 (1168) 256/16 (4249)

Table 33. GQA 2-D (𝑁 , 𝐵) autotune selections, non-causal. Cells give work-items/query rows (median latency, µs).

Table 30. GQA forward, causal (𝐵 = 1, 𝐻 = 32, 𝐻𝑘𝑣 = 8, 𝐷 = 128). Geomean 4.0× over SDPA; peak 72.1 TFLOPS (46% MFU) at 𝑁 = 8192. 𝑁

𝑁 \𝐵

𝑁 \𝐵

1

2

4

8

256 512 1024 2048

256/8 (53) 128/8 (107) 64/8 (296) 256/16 (1088)

128/8 (73) 128/8 (162) 256/16 (578) 256/16 (2007)

256/8 (101) 256/16 (307) 256/16 (1017) 256/16 (3988)

128/16 (171) 256/16 (536) 256/16 (2088) 256/16 (7828)

E.4.2 GQA (adaptive kernel, 𝐵 = 1, 𝐻 = 32, 𝐻𝑘𝑣 = 8, 𝐷 = 128). The tables below characterize full-sequence GQA at 𝐵 = 1, compared against SDPA with grouped-query attention enabled; numbers are the mean of three independent singleGPU runs. The runtime selection policy also covers 𝐵 ∈ {2, 4, 8} (Sec. E.5), where the selected variant differs from 𝐵 = 1. The causal separation follows from the baseline’s loss of parallelism at the single-batch grouped shape beyond 𝑁 = 4096. The adaptive GQA kernel retains parallelism, and its 𝑁 ≤ 1024 causal configurations select seq_q = 8 (Sec. E.5).

GQA 2-D (𝑁 , 𝐵) policy (compile-time else branch, group size 4). Unlike MHA, the selected GQA configuration depends on both sequence length 𝑁 and batch 𝐵. At low grid occupancy, seq_q = 8 doubles the number of concurrent query tiles. As 𝐵 fills the machine, the selected seq_q = 16 crossover moves to smaller 𝑁 , and larger workgroups are selected more often. The batch-aware sweep spans 𝐵 ∈ {1, 2, 4, 8}, 𝑁 ∈ {256, 512, 1024, 2048}, workgroup size ∈ {64, 128, 256}, and seq_q ∈ {8, 16} at decoder-model dimensions. Each cell has 24 repetitions and must pass a ≥ 200 TFLOPS calibration gate with coefficient of variation ≤ 7%. The tables report the selected workgroup size and seq_q with median latency in parentheses. Indexing the policy only by the 𝐵 = 1 measurements selects the wrong configuration after the crossover at 𝐵 ≥ 2. For causal 𝑁 = 1024, applying the 𝐵 = 1 selection (128 workitems and eight query rows) at 𝐵 = 4 takes 681 µs, 20% longer than 128 work-items and 16 query rows at 569 µs; at 𝐵 = 8, the difference is 38%. The 2-D dispatch function encodes these selections using pre-instantiated variants. Workgroup size depends on both batch and sequence length, even for a fixed query-tile size. The 𝐵 = 1 selections and the 𝑁 > 2048 long-context selection (256 work-items and 16 query rows) are preserved.

E.5

E.6

Table 31. Selected MHA configurations (compile-time MHA branch, group size 1). Regime

𝑁 ≤ 512

512 < 𝑁 ≤ 2048

𝑁 ≥ 4096

Non-causal Causal

seq_q = 8, wg=256 seq_q = 8, wg=128

seq_q = 16, wg=256 seq_q = 16, wg=128

seq_q = 16, wg=256 seq_q = 16, wg=256

Autotune Grid and Selection Policy

The grid sweep covers workgroup size ∈ {64, 128, 256} × seq_q ∈ {8, 16}, swept over all 𝑁 for both masks. The dispatch function encodes this table in a compile-time MHA branch, independent of the GQA policy. Handling the narrow seq_q = 8 register tile in the causal masking routine makes this configuration available and yields the 1.42× result at 𝑁 = 512.

Correctness

MHA (adaptive kernel, 𝐵 = 4, 𝐻 = 16, 𝐷 = 128, vs. a PyTorch SDPA reference). The table reports maximum absolute error with bf16 inputs and fp32 accumulation. All configurations pass. The causal seq_q = 8 tile passes at ≤ 0.011. GQA (adaptive kernel, 𝐵 = 1, 𝐻 = 32, 𝐻𝑘𝑣 = 8, vs. a KV-replicated SDPA reference). 27

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Causal GQA, 𝐵 = 1, 𝐻 /𝐻𝑘𝑣 = 32/8

Throughput (TFLOPS)

Causal MHA, 𝐵 = 4, 𝐻 = 16 80

80

60

60

40

40

20

20

SyclKittens IPEX SDPA

0

SyclKittens IPEX SDPA

0 512

1K

2K 4K 6K Sequence length 𝑁

512

8K

1K

2K 4K 6K Sequence length 𝑁

8K

Figure 9. Causal BF16 attention scaling on one Intel Max GPU with 𝐷 = 128. MHA uses 𝐵 = 4 and 𝐻 = 16; GQA uses 𝐵 = 1, 𝐻 = 32, and 𝐻𝑘𝑣 = 8. The baseline is IPEX SDPA, with its grouped-head fallback for GQA. SyclKittens sustains 66–72 TFLOPS after the shortest GQA point, while the single-batch grouped baseline peaks near 19 TFLOPS and declines at long sequence lengths. Table 34. MHA forward correctness, maximum absolute error vs. a PyTorch reference. Mask

𝑁 = 512

Causal Non-causal

0.0080

1024

2048

4096

8192

0.0078 0.0100 0.0109 ≤ 0.0019 across all 𝑁

0.0103

on one GPU. Plain norms read and write 2𝐵𝑁 𝐷 elements; fused residual-add, norm, and residual writeback read and write 4𝐵𝑁 𝐷 elements. The measured kernels assign eight cooperating subgroups to each row and use 128-bit load_wide/store_wide transactions. Reported useful traffic is 2 𝐵𝑁 𝐷 · 2 bytes for plain norm and 4 𝐵𝑁 𝐷 ·2 bytes for fused norm at 2 B/elem. This common accounting convention makes the GB/s columns comparable; it excludes extra intermediate traffic in implementations that launch multiple kernels.

Table 35. GQA forward correctness, maximum absolute error vs. a KV-replicated PyTorch reference. Mask Causal Non-causal

𝑁 = 512

1024

2048

4096

6144

8192

F.2

0.0134 0.0019

0.0078 0.0010

0.0079 0.0008

0.0091 0.0008

0.0098 0.0005

0.0097 0.0005

Each norm is a single translation unit templated on row width, subgroups per row, and whether to add and write back the residual. It is built into a per-width pybind11 module with separate plain and fused entry points. The kernel is written entirely in SyclKittens, with no hand-rolled intrinsics, on top of the load_wide/store_wide 128-bit vector-memory ops. Data flow (per row). The grid is (1, 𝐵, 𝑁 ) with WPR · 16 threads per block. One work-group owns each row, and its WPR subgroups cooperate. Subgroup 𝑤 owns the contiguous slice 𝑥 [row, 𝑤 ·𝐷𝐶 : (𝑤+1)·𝐷𝐶] in an rv_fl<DC, naive> register vector. This uses 𝐷𝐶/16 fp32 values per lane, 16 at 𝐷 = 2048 and 32 at 𝐷 = 4096, both within the register file. Each subgroup loads its chunk with load_wide. In the fused path it adds the residual chunk and writes the residual sum with store_wide. The subgroup then computes its partial moment or moments and publishes them to SLM. After one group_barrier, every subgroup reads all partials, forms the row statistic, applies the affine transform, and writes its output chunk with store_wide. Key techniques.

All GQA configurations pass (tolerance 0.10); the causal seq_q = 8 path (𝑁 ≤ 1024) passes at ≤ 0.0134.

F

Memory-Bound Norms, RoPE, and Decode GEMV

Root mean square layer normalization (RMSNorm) [38], layer normalization (LayerNorm) [3], rotary position embedding (RoPE) [31], and single-token GEMV all depend on wide coalesced memory access, with subgroup cooperation added where a row reduction is required. Norms and RoPE are compared with eager PyTorch and torch.compile (TorchInductor to Triton on Intel Max GPU); decode GEMV is compared with the tuned vendor library path (oneMKL). F.1

Configuration

RMSNorm uses 𝐷 = 2048 and LayerNorm 𝐷 = 4096, with bf16 inputs/outputs, fp32 accumulation, and 𝜖 = 10−5 . We test (𝐵, 𝑁 ) ∈ {(1, 1), (1, 64), (1, 512), (4, 1024), (16, 2048), (64, 2048)}

Kernel Implementation and Key Techniques

• Cooperative multi-subgroup row decomposition. Splitting the row across WPR = 8 subgroups reduces 28

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

per-lane accumulators from 𝐷/16 (single-subgroup layout with 128 at 𝐷=2048, 256 at 𝐷=4096, which spills) to 𝐷/(16 WPR) (16 / 32, no spill). This decomposition increases measured LayerNorm bandwidth from ∼570 to ∼970–1100 GB/s. • 128-bit coalesced load_wide/store_wide. Each lane issues a full uint4 (16 B) per transaction; 16 lanes touch 256 contiguous bytes in one cache-linealigned burst (interleaved-128 mapping), the widest transaction the Intel Max GPU memory pipe accepts. The old scalar rv_fl fallback issued 32-bit strided loads at ≤half the effective bandwidth. Requires 𝐷𝐶 a multiple of 128 and a 16 B-aligned base (satisfied for 𝐷 ∈ {1024, 2048, 4096}). • SLM moment reduction. RMSNorm reduces one moÍ ment ( 𝑥 2 ) through a WPR-float SLM scratch; LayerÍ Í Norm reduces two ( 𝑥, 𝑥 2 ) through 2 · WPR floats, Í Í then forms 𝜇 = 𝐷1 𝑥, 𝜎 2 = 𝐷1 𝑥 2 − 𝜇 2 . A single group_barrier separates publish from read-back. • Single-launch fused residual path. The kernel folds the 𝑥 +residual add into the load and emits the residual sum, so the pre-norm transformer’s reused hidden state ℎ = 𝑥 + sub(𝑥) is materialized without a second kernel or a re-read of ℎ. • Register-fused affine + rsqrt. inv_std = rsqrt(var+ eps) maps to the hardware reciprocal-sqrt unit; the normalize, scale (𝛾), and shift (𝛽, LayerNorm) are fused in registers before the store, with 𝛾/𝛽 themselves streamed via load_wide. F.3

the HBM roofline. Warm steady-state latency is recorded separately for decode. Traffic is 2 𝐵𝑁 𝐷 · 2 B for plain norm and 4 𝐵𝑁 𝐷 · 2 B for fused norm; GB/s is traffic divided by median time. Before timing, each configuration checks the plain norm, fused norm, and residual store against a PyTorch reference with relative error below 2e−2. F.4

Delivered Performance

The norm tables below report per-configuration cold-HBM medians. SyclKittens/eager is the speedup over eager PyTorch, and SyclKittens/compile is the speedup over torch.compile. F.4.1 RMSNorm (𝐷 = 2048). Plain, geomean 1.12× over torch.compile, 15.6× over eager. Fused (residual + norm, 4×), geomean 1.04× over torch.compile. F.4.2 LayerNorm (𝐷 = 4096). Plain, geomean 1.42× over torch.compile, 2.78× over eager. Fused (residual + norm, 4×), geomean 1.12× over torch.compile. For the HBM-bound shapes (𝐵 ≥ 16), eager PyTorch sustains ∼45–320 GB/s because its decomposition traverses the data in several kernels. Triton fusion raises torch.compile to 517–877 GB/s on LayerNorm, while SyclKittens sustains 623–997 GB/s (1.08–1.63×). At the two large RMSNorm cells, both paths approach the ∼1.3 TB/s per-GPU roofline and differ by 0.99–1.05×; at 𝐵 = 4, lower launch overhead gives SyclKittens a 1.13–1.33× advantage.

Methodology

Measurements follow the shared isolation protocol on one selected GPU. Eager F.rms_norm/F.layer_norm decomposes RMSNorm into multiple kernels; torch.compile lowers the same operations through Inductor to Triton on XPU. For the fused comparison, SyclKittens writes both the normalized output and residual sum (4× traffic). The expression torch.compile(F.rms_norm(x+res)) discards the sum (3×) and is approximately 30% faster at identical bandwidth. All three measured implementations therefore return (norm, x+res) and use the same useful-byte denominator. An unfused add-then-norm that retains the residual moves 5× traffic; the fused both-output kernel moves 4×, a 20% traffic reduction. The three implementations are interleaved round-robin within each of 10 trials, each with 10 warmup and 30 measured iterations; tables report the per-configuration median. A 2 s GEMM power-condition burn-in precedes measurement, and the Inductor graph-cache limit is 64 to prevent eager fallback. All cold-HBM tables use one cache protocol, in which a ∼256 MB buffer is zeroed before every timed call to evict L2, following Triton’s convention, so all reads miss to HBM. Without eviction, small-𝐵 working sets below the ∼204 MB per-GPU L2 report cache bandwidth above

F.5

Decode/Dispatch-Bound Latency

For few-token decode shapes, each norm touches at most 512 rows and fixed launch overhead dominates. SyclKittens uses one Level-Zero submission; torch.compile incurs a Dynamo guard and Triton launch, while eager uses multiple kernels. The preceding matmul or attention leaves the decode norm input L2-resident, so the table reports warm steadystate latency for that cache state. For fused residual-add and norm, SyclKittens holds a warm median ∼24 µs dispatch floor across the sweep while torch.compile sits at ∼55–61 µs. SyclKittens’s near-flat latency (24.2–24.6 µs for RMSNorm and 23.6–23.7 µs for LayerNorm) reflects one fused submission. The torch.compile floor is ∼2.3–2.6× higher, while eager latency grows with 𝑁 as its decomposed operations re-traverse each row. F.6

Correctness

Every configuration is validated before timing, namely the plain norm and fused norm against a PyTorch fp32-reference (relative error), and the residual sum 𝑥 + res byte-exactly. Representative maxima appear below. 29

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 36. RMSNorm, plain, 𝐷 = 2048. 𝐵

𝑁

SyclKittens (µs)

SyclKittens GB/s

eager (µs)

SyclKittens/eager

compile (µs)

compile GB/s / SyclKittens/compile

4 16 64

1024 2048 2048

61.8 378 1090

543 711 989

735 5810 22400

11.9× 15.4× 20.7×

82.6 395 1090

406 / 1.33× 680 / 1.05× 981 / 1.01×

Table 37. RMSNorm, fused, 𝐷 = 2048. 𝐵

𝑁

SyclKittens (µs)

SyclKittens GB/s

eager (µs)

SyclKittens/eager

compile (µs)

compile GB/s / SyclKittens/compile

4 16 64

1024 2048 2048

112 655 2190

600 820 980

815 6700 25900

7.28× 10.2× 11.8×

126 669 2170

534 / 1.13× 803 / 1.02× 992 / 0.99×

Table 38. LayerNorm, plain, 𝐷 = 4096. 𝐵

𝑁

SyclKittens (µs)

SyclKittens GB/s

eager (µs)

SyclKittens/eager

compile (µs)

compile GB/s / SyclKittens/compile

4 16 64

1024 2048 2048

108 641 2160

623 838 993

255 1790 7000

2.37× 2.80× 3.24×

130 921 3530

517 / 1.21× 583 / 1.44× 608 / 1.63×

Table 39. LayerNorm, fused, 𝐷 = 4096. 𝐵

𝑁

SyclKittens (µs)

SyclKittens GB/s

eager (µs)

SyclKittens/eager

compile (µs)

compile GB/s / SyclKittens/compile

4 16 64

1024 2048 2048

210 1160 4310

639 922 997

424 3570 13700

2.02× 3.07× 3.18×

228 1320 4900

589 / 1.08× 812 / 1.14× 877 / 1.14×

All configurations pass; the residual store is exact (both SyclKittens and Inductor write the same 𝑥 + res). The larger RMSNorm rel-err (2.8e-3 versus 1.4e-5 for LayerNorm) reflects bf16 rounding of the un-centered second moment. F.7

The kernel instead re-tiles the raw 128-bit packing from the earlier hand-tuned norm kernel for rotary pairs. • One work-item owns one 8-element chunk-pair. Work-item (row, 𝑐) takes the 𝑐-th 8-wide chunk of the first half, 𝑥 1 [𝑐·8 : 𝑐·8+8], and its partner chunk in the second half, 𝑥 2 [𝑐·8 : 𝑐·8+8], 64 bf16 elements away. Since 𝐷 = 128, there are 8 chunk-pairs per head; consecutive work-items read consecutive 16byte words, so each subgroup’s transaction is a fully coalesced 256-byte burst. • Four 16-byte loads, two 16-byte stores, all coalesced. Each item issues one uint4 load apiece for its 𝑥 1 , 𝑥 2 , cos, and sin chunks. Each load is 16 B, or eight bf16 elements. The item computes the eight rotary pairs in fp32 registers and stores the two output chunks as two uint4 writes. Every 𝑥/out transaction is the widest the Intel Max GPU pipe accepts. • Local pair computation. RoPE requires no SLM, barrier, or shuffle because each chunk-pair is local to one work-item. This mapping exposes independent work at high occupancy. It requires 𝐷 to be a multiple of 16 and the base to be 16 B aligned, as it is for contiguous [𝐵, 𝑁 , 𝐻, 𝐷] bf16 with 𝐷 = 128.

Rotary Position Embedding (RoPE) with 128-bit Chunk-Pair Mapping

RoPE is the transformer block’s other memory-bound elementwise op. The 128-bit coalesced kernel reaches 1020 GB/s. Relative to torch.compile, its geometric-mean speedups are 1.09× for cold bandwidth-bound shapes and 1.98× for warm decode shapes. F.7.1 Configuration. We evaluate halved rotary with 𝐷 = 128, 𝐻 = 32, bf16 inputs/outputs, and fp32 computation on the same six (𝐵, 𝑁 ) cells as the norms. Baselines are eager PyTorch and torch.compile; useful traffic is 2𝐵𝑁 𝐻 𝐷 · 2 bytes for one input read and output write. F.7.2 Wide Chunk-Pair Mapping. RoPE couples element 𝑑 with 𝑑 + 64 across the two halves of a head. A contiguous wide load groups eight adjacent bf16 elements, so it never places a rotary pair on the same lane. The SyclKittens load_wide also requires a 128-element multiple and is designed for reductions, so it cannot be used unchanged. 30

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 40. Decode-regime latency, fused norm, warm. op

𝐵

𝑁

SyclKittens (µs)

compile (µs)

eager (µs)

SyclKittens/compile / SyclKittens/eager

RMSNorm RMSNorm RMSNorm LayerNorm LayerNorm LayerNorm

1 1 1 1 1 1

1 64 512 1 64 512

24.3 24.6 24.2 23.6 23.7 23.6

58.0 61.3 61.3 54.8 60.9 61.3

87.6 105 169 46.7 46.5 66.2

2.30× / 3.60× 2.49× / 4.25× 2.53× / 7.00× 2.31× / 1.98× 2.57× / 1.96× 2.60× / 2.80×

Table 41. Correctness, RMSNorm/LayerNorm. Check SyclKittens plain rel-err torch.compile plain rel-err residual-store rel-err

RMSNorm (𝐷=2048)

LayerNorm (𝐷=4096)

Tolerance

2.8e-3 2.8e-3 0.0

1.4e-5 1.2e-5 0.0

2e-2 2e-2 2e-2

F.7.3 Delivered Performance. RoPE follows the shared interleaved-trial and power-conditioning methodology. The bandwidth table uses the cold-HBM protocol; the decode table uses warm, L2-resident inputs. SyclKittens/eager is the speedup over eager, and SyclKittens/compile is the speedup over torch.compile. Bandwidth (HBM-bound, L2 flushed before every call, cold), geomean 1.09× over torch.compile. At 𝐵=64, SyclKittens measures 1020 GB/s and torch.compile measures 882 GB/s; the large-𝐵 norm kernels measure 980–997 GB/s. Decode/dispatch-bound (warm, L2-resident input), latency (µs), geomean 1.98× over torch.compile. As with the norms, SyclKittens holds a near-flat ∼27–31 µs single-launch floor across the few-token sweep while torch.compile sits ∼1.9–2.1× higher. The difference is one Level-Zero submission versus a Dynamo guard followed by a Triton launch.

The DPAS systolic array requires 𝑀 ≥ 8, so 𝑀 = 1 uses at most one eighth of that matrix shape. For this GEMV comparison, we call oneMKL oneapi::mkl::blas::gemm directly. Eager F.linear and torch.compile instead reach the corresponding oneDNN GEMM path through PyTorch; on this platform the two entry points share the underlying GEMM implementation. F.8.1 Configuration. We evaluate the five Llama-3.18B [7] 𝑀 = 1 projections at TP= 1, namely𝑊𝑞𝑘𝑣 (6144×4096), 𝑊𝑜 (4096 × 4096), 𝑊𝑔𝑎𝑡𝑒 and 𝑊𝑢𝑝 (14336 × 4096), and 𝑊𝑑𝑜𝑤𝑛 (4096 × 14336). Inputs and outputs are bf16 with fp32 accumulation. The tuned vendor library path is oneMKL GEMM plus its fp32-to-bf16 cast; SyclKittens uses eight split-K subgroups per output row, 256-element software-prefetched chunks, and a direct bf16 store. F.8.2 Split-K Streaming Reduction. A GEMV row is a dot product, structurally the RMSNorm per-row reduction with the weight row in place of 𝑥. The kernel reuses the same SyclKittens wide-load reduction, restructured for the 𝑀 = 1 occupancy problem.

F.7.4 Correctness. Every config is validated before timing against an fp32 reference (relative-error tolerance 2e-2). SyclKittens computes the halved rotary in fp32 and rounds to bf16 identically to the reference, making the result bit-exact. The measured relative error is 0.0 at all six shapes, matching torch.compile. F.8

• Split-K occupancy. One workgroup owns each output row and WPR subgroups divide its 𝐾 reduction. The norms use the same decomposition to fit registers; here it also increases parallel work. Total launched work is 𝑁 ·WPR·16 threads, so it still scales with 𝑁 , but each row supplies multiple cooperating subgroups. This reduces sensitivity to the output-row count relative to oneMKL, as the tensor-parallel sweep below shows. • 128-bit coalesced load_wide over the 𝐾-slice. Each subgroup streams its contiguous slice 𝑊 [𝑛, 𝑤 ·𝐷𝐶 : (𝑤+1)·𝐷𝐶] (𝐷𝐶 = 𝐾/WPR) in CHUNK-wide pieces

Single-Token Decode GEMV (𝑀 = 1)

The decode-time weight projections (𝑊𝑞𝑘𝑣 , 𝑊𝑜 , 𝑊𝑔𝑎𝑡𝑒 , 𝑊𝑢𝑝 , 𝑊𝑑𝑜𝑤𝑛 ) are matrix-vector products out[𝑁 ] = 𝑊 [𝑁 , 𝐾]@𝑥 [𝐾]. They process one token, so 𝑀 = 1. Their arithmetic intensity is approximately one FLOP per byte of streamed weights. The kernel streams all of 𝑊 (𝑁 ·𝐾 ·2 B) once and does ∼2 FLOP/element. These five projections account for approximately 64% of per-token decode latency, and their throughput is determined by bf16 weight-streaming bandwidth. 31

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 42. RoPE bandwidth, cold. 𝐵

𝑁

SyclKittens (µs)

SyclKittens GB/s

eager (µs)

SyclKittens/eager

compile (µs)

compile GB/s / SyclKittens/compile

4 16 64

1024 2048 2048

115 646 2110

589 832 1020

370 3960 15800

3.24× 6.14× 7.46×

124 678 2440

544 / 1.08× 792 / 1.05× 882 / 1.15×

Table 43. RoPE decode-regime latency, warm. 𝐵

𝑁

SyclKittens (µs)

compile (µs)

eager (µs)

SyclKittens/compile / SyclKittens/eager

1 1 1

1 64 512

26.9 27.3 30.7

51.5 58.5 57.9

90.6 90.4 99.6

1.92× / 3.37× 2.15× / 3.32× 1.89× / 3.25×

oneMKL (42.2 ms, 0.86×), establishing that the wide-load and split-K structure produces the measured improvement. The kernel is latency-stable, because cold p99 is within ≤5% of p50 on every shape, so the p99 per-token GEMM is 28.4 ms (SyclKittens) vs 38.1 ms (oneMKL), the 1.33× holds at the tail. The L2-resident regime (per-op timed) gives 14.7 ms to 10.7 ms per token, 1.37×. Cold and L2-resident regimes. A single cold GEMV reads 33–117 MB from an idle memory pipeline and sustains 490–531 GB/s; the 0.5–1 GB norm streams approach the 990 GB/s long-stream asymptote. With weights resident in the 204 MB L2 and no flush, SyclKittens reaches 1060–1380 GB/s, 82–106% of the 1.3 TB/s HBM reference, and 10.7 ms/token. Values above 100% reflect L2 bandwidth. The approximately 16-GB model weights do not fit in L2, while the 160 GEMVs per token execute consecutively, placing decode between these two measurements. Across both regimes, SyclKittens exceeds oneMKL, with ratios of 1.24–1.53× cold and 1.02–1.90× resident, with an approximately 1.36× geometric mean in each. Split-K occupancy and removal of the cast and scratch path apply in both cache states.

via the same load_wide (uint4/lane, 256 B/subgrouptransaction) that took the norms to line rate. A narrow 16-bit, one-workgroup-per-row implementation reaches approximately 35% of that rate. • Software-prefetched chunks, fp32 accumulate. Chunk 𝑐+1’s load_wide is issued before the fp32 FMA of chunk 𝑐 consumes chunk 𝑐, following the GEMM pipeline in Appendix D. 𝑥 is reused across all rows and stays L2-resident. • Direct bf16 output. A single-round subgroup sum plus a tiny WPR-float SLM combine (one group_barrier) yields out[𝑛], written straight to bf16. The oneMKL path uses an fp32 scratch buffer and a separate cast kernel, adding one launch and 𝑁 ·(4+2) B of traffic. WPR = 8, CHUNK = 256 requires 𝐾 mod 2048 = 0 (holds for all five unsharded shapes); TP-sharded variants must also satisfy this constraint, and other 𝐾 fall back to oneMKL. The kernel is written on the SyclKittens load_wide/rv_fl/sum operations, no hand-rolled intrinsics. A hardware L1-prefetch variant of load_wide also exists in SyclKittens but is not used here. With only 2–7 chunks per slice, the register double-buffer already covers the latency; measurements showed no benefit and occasional regression from the extra prefetch.

F.8.4 TP Sharding and Output-Dimension Occupancy. Inference runs tensor-parallel, so the column-parallel projections (𝑊𝑞𝑘𝑣 , 𝑊𝑔𝑎𝑡𝑒 , 𝑊𝑢𝑝 ) shard the output dimension 𝑁 to 𝑁 /𝑛 at TP= 𝑛. We emulate TP= 𝑛 by running the per-shard (𝑁 /𝑛, 𝐾) GEMV on a single GPU with the same algorithm. The SyclKittens/oneMKL cold-HBM ratio rises monotonically with 𝑛. Values are per-shape cold p50, reproducible across two runs, with fixed WPR = 8. oneMKL’s 𝑀=1 GEMM parallelism is its output-tile count ⌈𝑁 /𝐵 𝑁 ⌉, so as 𝑁 shrinks it loses tiles and its bandwidth falls from 360 to 135 GB/s for 𝑊𝑞𝑘𝑣 as 𝑁 goes from 6144 to 768. SyclKittens’s split-K makes launched work 𝑁 · WPR (still ≥ 6144 subgroups at 𝑁 = 768), so it declines only gently (503 to 323). The ratio widens because oneMKL loses outputtile parallelism faster; both paths slow at smaller 𝑁 . Ratios separate this occupancy effect from the latency reduction caused by smaller shards. Row-parallel 𝑊𝑜 /𝑊𝑑𝑜𝑤𝑛 instead shard 𝐾 with 𝑁 fixed at 4096. Their parallelism and ratio remain unchanged, isolating output-dimension occupancy as the cause of the widening column-parallel ratio. The ratio reaches 2.4× at TP= 8. A fixed WPR = 8 captures this effect; a per-shape WPR(𝑁 ) sweep found no further gain, since even the smallest shard

F.8.3 Delivered Performance. Each shape reports the median (p50) over 300 individually event-timed iterations; the harness also records p95/p99. The table uses the shared cold-HBM cache protocol, so 𝑊 is read from HBM on every call. SyclKittens/oneMKL is the speedup over the oneMKL decode path; %roof uses the measured ∼1.3 TB/s per-GPU streaming-HBM reference. The specification is 1.6 TB/s, and the large-operand kernels above approach a ∼990 GB/s asymptote. Cold HBM, geomean 1.36× over oneMKL. Summed per-token GEMM latency (5 GEMVs × 32 layers, TP= 1) falls from oneMKL 36.3 ms to SyclKittens 27.2 ms (1.33×). A naive one-work-group-per-row bf16 stream (narrow 16-bit loads plus a 256-wide reduce) is slower than 32

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 44. Decode GEMV, cold HBM. shape

𝑁 ×𝐾

oneMKL (µs)

oneMKL GB/s

SyclKittens (µs)

SyclKittens GB/s

%roof

SyclKittens/oneMKL

𝑊𝑞𝑘𝑣 𝑊𝑜 𝑊𝑔𝑎𝑡𝑒 𝑊𝑢𝑝 𝑊𝑑𝑜𝑤𝑛

6144×4096 4096×4096 14336×4096 14336×4096 4096×14336

140 106 273 274 342

360 317 431 429 344

99.5 69.0 221 222 240

506 486 531 530 490

39% 37% 41% 41% 38%

1.40× 1.54× 1.23× 1.24× 1.42×

Table 45. TP-sharding SyclKittens/oneMKL cold-HBM ratio (cells are oneMKL / SyclKittens GB/s, then the ratio). proj 𝑊𝑞𝑘𝑣 (𝑁 =6144/𝑛) 𝑊𝑔𝑎𝑡𝑒 (𝑁 =14336/𝑛)

𝐾 4096 4096

TP= 1

TP= 2

TP= 4

TP= 8

359 / 503, 1.40× 429 / 530, 1.24×

279 / 496, 1.78× 365 / 513, 1.40×

183 / 410, 2.24× 273 / 480, 1.76×

135 / 323, 2.39× 186 / 411, 2.21×

already saturates at WPR = 8 and larger WPR only adds combine overhead. F.8.5 Correctness. Each shape is checked against an fp32 reference out_ref = (𝑊 .float()@𝑥 .float()) before timing. SyclKittens accumulates in fp32 and applies one bf16 round on store, so its max-abs √ error is identical to the oneMKL path (both round an 𝑂 ( 𝐾)-magnitude sum once). The error is in {0.50, 0.50, 0.80, 0.50, 0.99} for the five shapes, reflecting bf16 output quantization. All shapes match oneMKL to the last bf16 ULP band.

G

Value

Problem

Mistral-7B / Llama-3.1-8B Q/K geometry: 𝐻 = 32, 𝐻𝑘𝑣 = 8, 𝐷 = 128; shared cos/sin (𝑁 , 𝐷/2) decode/dispatch-bound (𝐵, 𝑁 ) ∈ { (1, 1), (1, 64), (1, 512) }; bandwidth/HBM-bound { (4, 1024), (16, 2048), (64, 2048) } bf16 in/out, fp32-internal rotary math; cos/sin bf16 (𝑁 , 𝐷/2) Halved rotary: 𝑜𝑑 = 𝑥 1,𝑑 𝑐𝑑 − 𝑥 2,𝑑 𝑠𝑑 , 𝑜𝑑+64 = 𝑥 2,𝑑 𝑐𝑑 + 𝑥 1,𝑑 𝑠𝑑 Intel Max GPU, one GPU, flat device hierarchy eager PyTorch (two halved rotaries on Q and K); torch.compile (Inductor to Triton, XPU) of the same op-sequence One fused Q/K launch with 𝐷 = 128 and four workers

Sweep

Dtype Convention Hardware Baselines

Op-Fusion Microbenchmarks Launch

Fusing Q/K rotary combines two launches and allows cached cos/sin data to be reused across both tensors. We compare this path with eager PyTorch and torch.compile (TorchInductor to Triton on Intel Max GPU) using the cold/warm methodology from Appendix F; Section K.4 measures the complete composed path. G.1

Parameter

G.1.1

Configuration.

G.1.2 Kernel implementation and key techniques. The kernel reuses the 128-bit coalesced chunk-pair rotary worker from Appendix F and fuses Q and K by laying their work items end to end in one nd_range, so the linear-id space is the concatenation of Q’s 𝐵 𝑁 𝐻 · CHUNKS items and K’s 𝐵 𝑁 𝐻𝑘𝑣 · CHUNKS items. A single tid < q_threads test routes each work-item to its tensor with the correct head count; both branches call the same per-row worker. One work-item owns one 8-element chunk-pair of one headrow; with 𝐷 = 128, (𝐷/2)/8 = 8 chunk-pairs per row. RoPE couples element 𝑑 with 𝑑+64, so each item loads its firsthalf chunk x1 and its second-half partner x2 (64 bf16 apart) plus the shared cos and sin chunks. That is four coalesced 16-byte loads, the rotary math in fp32 registers, then two 16-byte stores.

Fused Q/K Rotary

Every decoder layer, before attention, applies RoPE to both the projected query and key. For Mistral-7B [16] and Llama3.1-8B [7] GQA, 𝐻 = 32 query heads, 𝐻𝑘𝑣 = 8 key heads, and 𝐷 = 128. The baseline launches Q and K rotary separately over the shared (𝑁 , 𝐷/2) tables. The fused kernel assigns either a Q or K head-row to each work-item and uses identical halved-rotary math. Traffic ledger. Both implementations read Q, write 𝑄 ′ , read K, and write 𝐾 ′ exactly once. Launch overhead controls few-token decode; fusion reduces two launches to one. The shared (𝑁 , 𝐷/2) cos/sin tables are loaded by both Q and K work-items and can be reused through the cache within the same launch. Useful traffic excludes these shared table reads and counts the Q/K tensors. It is therefore 2 𝐵 𝑁 (𝐻 +𝐻𝑘𝑣 )𝐷 ·2 bytes (read and write Q and K); at 𝐵=64, 𝑁 =2048 that is 2.68 GB moved in a single launch at 973 GB/s.

// ONE nd_range whose linear-id space concatenates the Q and K work items: // q_threads = B*N*H*CHUNKS come first, then k_threads = B*N*HKV*CHUNKS. // CHUNKS = (D/2)/8 chunk-pairs per row (8 at D=128). 33

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

// grid dispatch: one branch selects the sub-grid, both call one worker const long tid = item.get_global_linear_id(); if (tid < q_threads) // Q sub-grid: heads = H rotary_row_chunk<D>(q, q_out, cos, sin, /*heads=*/ H, N, tid); else // K sub-grid: heads = HKV rotary_row_chunk<D>(k, k_out, cos, sin, /*heads=*/ HKV, N, tid - q_threads);

and 30 measured iterations. Each config is measured twice, once warm (operands stay L2-resident, the decode metric) and once cold (a scratch buffer is zeroed to flush L2 before every timed call, following Triton’s cold-cache convention). The cold panel reports HBM bandwidth for the bandwidth-bound shapes; the warm panel reports few-token decode latency. Before timing, rotated Q and K from both SyclKittens and torch.compile are checked against an fp32 reference with relative-error tolerance 2e-2.

// the shared per-row worker (identical math to the unfused rotary) const long row = local_tid / CHUNKS; const int chunk = local_tid % CHUNKS; const int coff = chunk * 8; // bf16 offset within a half const long row_off = row * D; // (.,N,heads,D) row base const int n = (row / heads) % N; // seq index for cos /sin const long cs = n * (D / 2) + coff; // shared cos/sin base

G.1.4 Delivered performance. Per-config medians come from 10 co-measured trials of eager PyTorch, torch.compile, and SyclKittens. Bandwidth-bound shapes are reported as cold, L2-flushed HBM bandwidth; the few-token decode shapes as warm dispatch latency. SyclKittens/eager is the speedup over eager, and SyclKittens/compile is the speedup over torch.compile. Table 46. Fused Q/K rotary, bandwidth-bound shapes, cold HBM bandwidth (L2 flushed). Geomean 1.41× over torch.compile, 9.03× over eager; SyclKittens peaks at 973 GB/s.

const wide8_t ra = xw[(row_off + coff) / 8]; const wide8_t rb = xw[(row_off + D / 2 + coff) / 8]; const wide8_t rc = cw[cs / 8], rs = sw[cs / 8]; float a[8], b[8], c[8], s[8]; unpack8(ra, a); unpack8(rb, b); unpack8(rc, c); unpack8(rs, s);

𝐵

𝑁 SyclKittens µs SyclKittens GB/s eager µs SyclKittens/eager compile µs compile GB/s SyclKittens/compile

4 1024 16 2048 64 2048

228 878 2760

369 765 973

2810 6800 21500

12.2× 7.74× 7.77×

385 1110 3650

219 605 735

1.69× 1.26× 1.32×

Table 47. Fused Q/K rotary, decode/dispatch-bound shapes, warm latency. Geomean 2.06× over torch.compile, 17.2× over eager.

float o1[8], o2[8]; #pragma unroll for (int i = 0; i < 8; ++i) { o1[i] = a[i] * c[i] - b[i] * s[i]; // rotate first half o2[i] = b[i] * c[i] + a[i] * s[i]; // rotate second half } ow[(row_off + coff) / 8] = pack8(o1); // two coalesced stores ow[(row_off + D / 2 + coff) / 8] = pack8(o2);

𝐵

𝑁 SyclKittens µs compile µs eager µs SyclKittens/compile SyclKittens/eager

1 1 1 64 1 512

149 138 148

299 283 293

2620 2370 2390

2.13× 2.05× 1.98×

18.3× 17.3× 16.2×

Relative to torch.compile, fused Q/K rotary is 1.26– 1.69× faster for the cold bandwidth-bound shapes and 1.98– 2.13× faster for the warm decode shapes, where its latency remains 138–149 µs. The gain comes from one joint Q/K launch and reuse of the shared cosine/sine data; torch.compile emits separate guarded rotary launches, while eager decomposes each rotary into several elementwise kernels.

The concatenated id space reduces two Level-Zero submissions to one and permits cached reuse of the shared cos/sin data across Q and K rows. Each work-item uses 128-bit coalesced transactions, the widest accepted by the Intel Max GPU memory path, and requires no SLM, barrier, or shuffle because one work-item owns each rotary pair. Both paths call the same fp32 rotary_row_chunk worker.

G.1.5 Correctness. Against the fp32 reference rounded to bf16, fused SyclKittens is bit-exact (relative error 0.0) at all six benchmark shapes, while torch.compile matches within tolerance. Its maximum absolute difference is 0.0 from both the two-launch unfused rotary path and transformers.apply_rotary_pos_emb, with standalone coverage over 𝐵 ∈ {1, 2, 4, 8} × 𝑁 ∈ {256, 512, 1024, 2048}.

G.1.3 Methodology. One heavy task per GPU package, with one stack active, so there is no intra-package HBM contention. The baselines are eager PyTorch and torch.compile(backend="inductor") applied to the same function. Interleaved round-robin timing within each trial keeps implementations under the same clock conditions, controlling drift in the SyclKittens/eager and SyclKittens/compile ratios. Each of 10 trials uses 10 warmup

H

Distributed GEMM

The operation starts with row-sharded 𝐴 and replicated 𝐵 and must return the full replicated 𝐶 = 𝐴𝐵 on every GPU. 34

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Within the original hand-designed input-gather family, the topology-aware path combines an all-gather of 𝐴 with the SyclKittens collaborative de-duplicated-prefetch GEMM core. Each source stripe crosses XeLink once per remote GPU package, is replicated to the partner stack over MDFI, and is then consumed by one local GEMM. On 12 GPUs, this path reaches 84% and 93% of its co-measured GEMM-only baseline at 61443 and 122883 . The controlled agent study also explores computing local 𝐶 shards first and gathering the output. H.1

compute kernel and differ in how and when they place 𝐴 in local HBM.

Input-gather data flow (per GPU).. (1) each GPU owns one 𝑀/12×𝐾 row-shard of 𝐴 in local HBM; (2) the missing 11 shards are gathered from peers into a contiguous 𝑀×𝐾 landing buffer; (3) one register-resident GEMM computes the full 𝐶 = 𝐴full 𝐵 over all 𝑀 rows. The methods trade off how step (2) overlaps step (3).

Configuration Key techniques.

A node contains six GPU packages, each with two stacks connected by MDFI, the load/store stack-to-stack interconnect. XeLink connects the GPU packages all-to-all through two independent stack planes, even and odd. Parameter

Value

Problem

Row-sharded GEMM, replicated output: 𝐴local [𝑑] : bf16[𝑀/12, 𝐾], replicated 𝐵 [𝑑] : bf16[𝐾, 𝑁 ]; every GPU returns the full 𝐶 [𝑑] = 𝐴𝐵 : bf16[𝑀, 𝑁 ] 𝑀=𝑁 =𝐾 ∈ {6144, 12288} (square); shard rows per GPU rpd = 𝑀/12 ∈ {512, 1024} bf16 inputs/outputs, fp32 accumulation (inregister RNE down-convert on store) Six GPU packages × two stacks; intra-package MDFI and all-to-all inter-package XeLink over the even and odd stack planes Co-located SyclKittens GEMM-only reference (A already local, no communication, all 12 GPUs concurrent): 139 T @61443 , 135 T @122883 . Isolated GEMM on one GPU is approximately 158–160 T; the 12-GPU-concurrent baseline is lower because each GPU package’s two stacks share its power/clock budget. bf16 register-resident 32×64×32 DPAS tiling with snake swizzle, reused verbatim from the single-GPU kernel Hand input-gather family; agent local𝐶/output-gather family; engineer-agent refinement of the SyclKittens + examples kernel

Shapes Dtype Topology

Peak reference

Hand GEMM core Method

• Topology-aware gather. Cross XeLink once per stripe and replicate within each GPU package. The 12 stacks are organized as six GPU packages with two stacks each; each pair shares a fast intra-package MDFI link outside the XeLink fabric. Flat all-gather has every stack pull all 11 remote shards, roughly 10 fabric crossings per stack. The topology-aware path sends each source stripe across XeLink exactly once to one landing stack per remote GPU package (landing = 2*q + (p&1)), then the landing stack replicates the shard to its partner (partner(t) = tˆ1) over MDFI. This halves XeLink egress on the communicationbound workload, giving 1.42× (122883 ) to 2.15× (61443 ) end-to-end distributed-GEMM speedup over the flat gather pipeline in the main paper’s fixed-core movement ablation. • Within input-gather schedules, materialize-thenGEMM. The topology-aware path gathers all shards to local HBM and then runs one GEMM. Interleaving the gather into the GEMM mainloop or reading remote 𝐴 per DPAS inside the kernel is slower in the measured implementations. Remote operand reads stall the DPAS pipeline, while materialization lets the compute kernel use its local operand-delivery path. • Gather topology. DMA memcpy gather and EU vec128 copy gather land within 1% of each other at both shapes; both are fabric-bound. The topology-aware path changes the traffic topology. • Large shapes hide communication. At 122883 the GEMM is big enough that even flat gather reaches 63% of the baseline and the topology-aware path reaches 90%; at 61443 the communication-to-compute ratio is higher, so flat gather collapses to 38% while the topology-aware path reaches 82%. • Compute-core reuse. The topology-aware path uses the bf16 GEMM kernel (snake-swizzled 32×64×32 DPAS, per-shape snake width) verbatim; no distributed-specific compute changes, so the single-GPU tuning carries over.

TFLOP/s uses the semantic operation count 2𝑀𝑁 𝐾/𝑡 per GPU because every GPU returns the full replicated 𝐶. Inputgather methods physically recompute those FLOPs on every GPU; the local-𝐶 family computes one row band per GPU and gathers the completed bands. H.2

Kernel Implementation and Key Techniques

For the original hand-designed family, the GEMM kernel is DPAS-bound and reaches approximately 90% of its local roofline. Its distributed overhead therefore comes from moving sharded 𝐴 across the fabric. These methods use the same 35

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 48. Engineer-agent refinement of the SyclKittens + examples kernel at 61443 . Values are medians across four nodes of each node’s seven-repetition median; the ratio is paired within repetition.

H.3

Source

TFLOP/s

vs. incumbent

Agent kernel (incumbent) Refined kernel

123.839 125.724

1.000× 1.007×

confidence interval on that ratio spans 0.93× to 1.04×, so the refined kernel is best read as matching the incumbent rather than beating it. We report it as the certified engineer-agent source because it wins the promotion protocol. The study suggests enumerating collective placement over inputs, partial results, and outputs before tuning schedules, carrying topology routing across the chosen placement, and then retuning for the resulting local shard geometry. Our refinement protocol retains the strongest validated agent source as a selectable incumbent throughout tuning. An engineered variant is promoted only after it wins the same input, correctness, timing, and reduction protocol.

Methodology

All 12 GPUs run concurrently in a single process (12 SYCL queues, one per GPU), so the measurement includes intranode fabric and shared power contention. Each method is timed as the median of 30 iterations after 10 warmup iterations and a 2 s turbo burn-in, identical to the single-GPU GEMM reference protocol. Every method is co-measured in the same run against the GEMM-only baseline, so all ratios are within-run. Inputs use Box-Muller 𝑁 (0, 1) random initialisation (matching torch.randn) for the performance numbers; a separate constant-fill run gates correctness analytically (see the correctness discussion below). Per-GPU throughput is 2𝑀𝑁 𝐾/𝑡.

H.5

This path uses the same collaborative de-duplicated prefetch core as Appendix D. This cooperative-prefetch measurement is separate from the earlier fixed-core movement ablation’s 82% and 90% endpoints. Results are three-repetition medians on one quiet node; the distributed measurement includes the complete topology-aware all-gather and is normalized to the co-measured GEMM-only baseline. Table 49. Production distributed-GEMM run with Box– Muller 𝑁 (0, 1) inputs and the cooperative-prefetch build. “GEMM-only” is the local SyclKittens compute reference; “distributed” includes topology-aware all-gather of 𝐴 followed by the full SyclKittens GEMM on every GPU.

H.3.1 Controlled coding-agent task. Appendix C.5 reports the controlled coding-agent evaluation of the 61443 row-sharded operation, including its task definition, protocols, and correctness-gated results. H.4

Topology-Aware All-Gather with Collaborative Prefetch

Engineer-Agent Co-Design

Our first engineering pass refined a kernel from the SyclKittens context. A subsequent common-harness comparison with Box–Muller 𝑁 (0, 1) inputs measured the best SyclKittens and SyclKittens + examples kernels alongside that initial refinement. The SyclKittens + examples kernel reached 131.616 TFLOP/s, compared with 123.207 for the SyclKittens kernel and 128.545 for the initial refinement. We restarted from the unchanged SyclKittens + examples kernel and retained it as the incumbent. The selected refinement preserves its DPAS compute core and staggered even/odd-plane cross-package gather. It changes only the MDFI replication phase, in which package bands are visited cyclically, with the first two copies assigned to the compute queue and the remaining four to the copy queue. Candidate promotion requires both a paired ratio above one and an absolute median-of-node-medians no lower than the incumbent. Each repetition uses a 2-second burn-in, 10 warmups, and the median of 30 timed operations, with source order rotated across pairs. Both sources pass 4,096 sampled references on all 12 outputs before and after timing in all 56 executions; maximum post-timing relative L2 error is 8.67×10−7 . The margin is small next to run-to-run variation. Across the 28 paired repetitions the median ratio is 1.007× and the 95%

Shape

SyclKittens GEMM-only (T/GPU)

SyclKittens distributed (T/GPU)

Dist./GEMM

Snake width

61443 122883

134.9 134.7

113.1 124.9

84% 93%

2 4

The larger shape amortizes communication and reaches 93% of the local compute ceiling. The smaller shape is communication-bound but still retains 84%. Both shapes pass the analytic constant-fill correctness gate. The 113.1TFLOP/s result is not the campaign-local hand reference in Appendix C.5. That campaign uses its evaluator’s deterministic input distribution and a separately frozen build, while this table uses random-normal inputs and the cooperativeprefetch build. The coincident 134.9 values in the two appendices are also different quantities. Here it is SyclKittens GEMM-only throughput with no communication; there it is the median end-to-end oneCCL all-gather plus oneDNN GEMM throughput. Matched root-cause decomposition. We rebuilt both configurations from the same frozen source and compiler, then co-measured them on four healthy nodes using three repetitions per condition, a 2-second burn-in, 10 warmups, and 30 timed operations. Table 50 reports the cross-node median of each node’s repetition median at 61443 . 36

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 51. Context-parallel attention configuration.

Table 50. Controlled decomposition of the two distributedGEMM operating points. “Eval” draws from 32 values in [−1, 0.9375] at 1/16 spacing; “randn” is Box–Muller 𝑁 (0, 1). Build / snake width

Inputs

GEMM-only (T/GPU)

Distributed (T/GPU)

Cooperative / 2 Cooperative / 2 Masked / 8

randn Eval Eval

137.5 156.2 159.8

113.5 125.3 128.0

Problem

𝐻 = 𝐻 kv = 16, 𝐷 = 128, non-causal full attention over the ring; 𝐵 ∈ {1, 4} 𝑁 per ∈ {256, 1024, 2048, 4096, 8192, 16384} for the 𝐵 = 1 sweep (six GPU packages, one active stack each) 𝑁 total = 6 𝑁 per ∈ {1536, . . . , 98304} bf16 Q/K/V, fp32 output/LSE and accumulation Six GPU packages, one active stack each, inter-package XeLink ring 157.5 TFLOPS bf16 per active GPU (single-GPU attention roofline) O(𝑁 per )-memory ring attention and O(𝑁 total )-memory all-gather attention; both use the fence-free 256-work-item flash core with 16-query blocks PyTorch ring: sendrecv/aten, sendrecv/triton, async/aten; PyTorch all-gather: allgather/aten (tuned oneCCL library path), allgather/triton

Total context Dtype Topology Peak reference SyclKittens kernels

Baselines

Correctness

I.1

A separate constant-fill run (𝐴local [𝑑] = (𝑑+1)·0.01, 𝐵 = 0.01) gives an analytic reference. The cooperative-prefetch GEMMonly and topology-aware distributed paths pass the max-diff < 0.5 gate; the in-register bf16 RNE store is the only rounding source.

I

Value

Per-GPU shard

The paired campaign-to-appendix ratio is 1.130× across the four nodes. Holding the cooperative build and snake width fixed, changing only randn to Eval inputs contributes 1.109×; holding Eval inputs fixed, changing to the campaign’s masked-prefetch, snakewidth-8 configuration contributes another 1.022×. The cooperative build itself is 1.033× faster than the masked build at matched snake width 2 and randn inputs. The difference between the two operating points is primarily due to operand-dependent power behavior on Intel Max GPU, with a smaller contribution from the schedule difference. The controlled appendix-like result, 113.5 TFLOP/s, also reproduces the original 113.1 TFLOP/s measurement.

H.6

Parameter

Configuration

Per-GPU TFLOPS uses 4 𝐵𝐻 𝑁 per 𝑁 total 𝐷 because each GPU’s query shard attends the full context; this is directly comparable across methods and baselines. I.2

Kernel Implementation and Key Techniques

Ring attention shards context across six GPU packages, with one active GPU stack each. At every step, each GPU computes partial attention for its local queries against the resident KV shard, rotates that shard one XeLink hop, and updates the online-softmax state. The ring implementation pairs independently tuned collective and compute implementations. Both methods use ThunderKittens tile primitives.

Context-Parallel Attention (Ring and All-Gather)

Context-parallel attention combines the fence-free core from Appendix E with ring and all-gather KV distribution. Their memory footprints differ, so both paths are measured under the same performance and correctness gates on three clean nodes.

Data flow (per ring step, per active GPU).. Each step has three stages. First, the fence-free 256-work-item flash core with 16-query blocks from Appendix E computes partial 𝑂/LSE for the local queries and resident KV shard. Second, the kernel merges the max, normalizer, and output accumulator into the running online-softmax state. Third, a doublebuffered bulk memcpy moves the next KV shard from the ring neighbor through the shared SYCL context into the alternate local-HBM buffer. The in-order queue serializes compute and transfer, so the next step reads local HBM. The overlapped follow-up below uses copy engines to overlap rotation with compute.

• Ring attention. Each active GPU holds one KV shard and rotates shards around the XeLink ring, so perGPU memory is O(𝑁 per ). It supports contexts beyond one GPU’s HBM. • All-gather attention. Each active GPU gathers the full 𝑁 total KV once (fan-out push all-gather) and runs a single flash pass, so per-GPU memory is O(𝑁 total ). When the gathered context fits in HBM, it pays the inter-package cost once instead of at each ring step. Measurements use a clock-gated, three-node, ten-trial harness. The PyTorch paths pair ring (send/receive) and all-gather communication with both the autotuned 2-Dblock-pointer Triton flash kernel used in Appendix E and an IPEX/aten fused flash kernel. Every functional variant passes the same correctness gate.

Key techniques. • Bulk-memcpy KV rotation. KV shards are copied in bulk to local HBM and attended locally. A measured zero-copy variant instead reads K/V over XeLink 37

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

for each DPAS operation. Because these reads do not overlap asynchronously with the DPAS engine, that variant stalls compute and is 5–50× slower. • Double-buffered rotation. Two KV buffers (current/next) keep the active and incoming shards separate. They also support the overlapped follow-up, which hides rotation latency behind compute at large 𝑁 per . • Fence-free 256-work-item compute core. The same software-fence-free single-pass flash kernel used for single-GPU attention (Q preloaded once, IGC SWSB hazard ordering, split-D output, LPT launch order); about 15% faster at large 𝑁 than the earlier fence-heavy 128-work-item core it replaces. • In-order queue serialization (no explicit fences between phases). The per-GPU in-order SYCL queue orders compute then rotate implicitly; no host-side event fence is inserted between the flash kernel and the rotation memcpy.

by the cross-node median and min–max spread over three nodes. The PyTorch variants run on the same shards in the same harness. Ring baselines are sendrecv/aten, sendrecv/triton, and async/aten, the last using overlapped batched non-blocking send/receive. All-gather baselines are allgather/aten, the tuned vendor library path combining oneCCL all_gather with IPEX fused flash, and allgather/triton, which combines the oneCCL gather with autotuned Triton flash. These two are the algorithm-matched baselines for SyclKittens all-gather attention. Every functional method is checked against a full-context attention reference before inclusion. The SyclKittens allgather path is also verified in-binary against an onlinesoftmax reference (§I.5).

All-gather attention. Instead of rotating one shard per step, each active GPU gathers the entire 𝑁 total KV up front and runs one flash pass over the full context. Two techniques define this path.

Tables report cross-node median TFLOPS per active GPU over three clean nodes, with the node-level min–max spread in brackets; both methods use the same harness.

I.4

I.4.1 Ring attention. SyclKittens/PyTorch ratios compare SyclKittens ring attention with the blocking sendrecv ring baselines. The single-batch (𝐵 = 1) sweep covers 1.5K–98K total context (256–16384 tokens per GPU). Against sendrecv/aten, the SyclKittens ring margin narrows from 1.85× at 𝑁 per = 256 to 1.05× at 16384; against sendrecv/triton, it remains 1.76–2.89×. Table 52 reports the sequential schedule. Overlapping the shard rotation with compute on the copy engines raises the same six configurations to 1.7, 15.2, 31.7, 53.6, 60.7, and 69.0 TFLOP/s per GPU, which is 1.21–2.46× the aten ring and 1.45× at 49K total tokens (Figure 7(a)). The step breakdown shows how rotation and onlinesoftmax merge diminish relative to compute as the shard grows. From 𝑁 per = 256 to 8192, the flash-core share rises from 20% to 65%, while reduction falls from 41% to 6% and rotation from 40% to 29%. Ring throughput therefore approaches the single-GPU attention efficiency in Appendix E as the shard grows.

• Fan-out push all-gather. Each active GPU reads its own KV shard once and writes it into the full [𝐵, 𝑁 total, 𝐻 kv, 𝐷] buffer of all five peers at its own slot in one parallel loop over cached peer pointers. The self-shard uses a local memcpy. This is the fan-out push pattern from the general-purpose all-gather in Appendix J. It sustains approximately 400 GB/s aggregate across the node, compared with the approximately 16 GB/s non-adjacent link that limits strict ring rotation and approximately 108 GB/s for the earlier fan-in pull gather. • Single flash pass without online merge. The gathered KV is complete, so the validated 256-work-item core attends all 𝑁 total keys in one pass. The kernel already reads the KV buffer depth for the number of GPUs holding KV shards and uses 𝑁 per for the query grid, so 𝑁𝑞 ≠ 𝑁𝑘𝑣 needs only a thin dispatch wrapper, no softmax-merge reduction. Ring instead performs six rotations and online merges.

I.4.2 All-gather attention. SyclKittens all-gather attention gathers the full context once and runs a single flash pass. Ratios are SyclKittens over each PyTorch all-gather baseline (same algorithm), plus SyclKittens all-gather over SyclKittens ring on the identical shard. Single-batch (𝐵 = 1) query shard swept 256–16384 per GPU (gathered context 𝑁 total = 6 × 𝑁 per ). Relative to the sequential SyclKittens ring, the allgather advantage narrows from 4.26× at 𝑁 per = 256 to 1.17× at 16384, while remaining 1.65–2.79× faster than allgather/triton. The tuned oneCCL-plus-IPEX path overtakes it only at the two largest shards, where the

The trade-off is memory, because all-gather is O(𝑁 total ) per GPU and runs out of memory once the gathered context exceeds HBM, whereas ring stays O(𝑁 per ). Section I.4 quantifies the crossover. I.3

Delivered Performance

Methodology

The ring uses the shared isolation protocol across six GPU packages, avoiding intra-package HBM contention. Each configuration runs 10 independent trials of 15 iterations. PyTorch ring trials reduce timing to the slowest rank with a max all-reduce. Results use the per-node median followed 38

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 52. Sequential ring attention, 𝐵 = 1 context sweep (six GPU packages; one active stack each). The aten and Triton columns use blocking send/receive. 𝑁 per

𝑁 total SyclKittens ring (T) aten (T) SyclKittens/aten Triton (T) SyclKittens/Triton

256 1024 2048 4096 8192 16384

1536 6144 12288 24576 49152 98304

1.6 [2–2] 11.1 [11–11] 20.6 [21–21] 34.4 [34–34] 47.7 [47–48] 59.6 [59–60]

1.85× 1.71× 1.59× 1.40× 1.13× 1.05×

0.9 6.5 12.9 24.6 42.0 56.9

0.6 5.2 10.3 17.6 26.8 33.9

2.89× 2.13× 2.01× 1.95× 1.78× 1.76×

Table 53. Ring attention step breakdown. 𝑁 per

compute

reduce

P2P rotation

256 4096 8192

20% 48% 65%

41% 10% 6%

40% 43% 29%

Throughput per active GPU (TFLOPS)

Table 54. All-gather attention, 𝐵 = 1 context sweep (six GPU packages; one active stack each). Both baselines use oneCCL all-gather. 𝑁 per

𝑁 total SyclKittens (T) aten (T) SyclKittens/aten Triton (T) SyclKittens/Triton SyclKittens/ring

256 1024 2048 4096 8192 16384

1536 6144 12288 24576 49152 98304

7.0 [7–7] 32.4 [32–32] 41.5 [42–42] 61.3 [61–61] 61.2 [61–63] 69.6 [69–70]

3.1 20.6 34.9 50.2 66.3 75.3

2.24× 1.57× 1.19× 1.22× 0.92× 0.93×

Separated ring

2.5 14.4 21.8 29.2 37.1 39.7

2.79× 2.24× 1.91× 2.10× 1.65× 1.75×

4.26× 2.93× 2.02× 1.78× 1.28× 1.17×

All-gather attention

60

40

20

0 256

1K 2K 4K Context tokens per active GPU (𝑁 per )

8K

16K

Figure 10. Non-causal context-parallel attention on six Intel Max GPU packages, one active stack each, with 𝐵 = 1, 𝐻 = 𝐻𝑘𝑣 = 16, and 𝐷 = 128. All-gather exceeds the sequential SyclKittens ring across this sweep; ring attention retains 𝑂 (𝑁 per ) memory for longer contexts. Both curves are three-node medians. I.4.3 Memory-capacity trade-off. HBM capacity determines which path applies.

SyclKittens/aten ratios are 0.92× and 0.93×. At 𝑁 per = 8192, the fan-out gather takes 4.7 ms (9%) and flash takes 46.3 ms (91%), so the flash pass controls throughput at this shape. 39

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

Table 55. Ring attention baseline correctness, maximum relative error. 𝑁 per = 256

1024

2048

4096

8192

16384

0.0191 0.0204

0.0126 0.0126

0.0101 0.0102

0.0076 0.0076

0.0063 0.0063

0.0041 0.0041

sendrecv/aten sendrecv/triton

dispatches adaptively between a latency-oriented path and a topology-aware bandwidth path. J.2

Every SyclKittens collective builds on direct peer read/write, so GPU 𝑑’s kernel addresses GPU 𝑘’s HBM through the shared SYCL context, with a single host barrier per phase and no staging through host memory. Three measured techniques shape the per-collective kernels.

• Gathered context fits in HBM. All-gather pays the inter-package cost once and is 1.17–4.26× faster than the sequential ring across the measured 𝐵 = 1 sweep. • Context exceeds one GPU’s HBM. Ring retains an O(𝑁 per ) footprint and increases supported context with the active GPU count.

• Remote-write fan-out. On the measured peer route, a remote write sustains ∼18.8 GB/s versus ∼15.7 GB/s for a remote read. All-gather and large-payload allreduce use fan-out push. Each GPU reads its local slice once and writes it to all 𝑁 −1 remote destinations. Local HBM reads cost 𝑒𝑝𝑑 rather than (𝑁 −1) 𝑒𝑝𝑑, leaving bandwidth for incoming remote writes and avoiding the HBM collapse of large-payload pull kernels. • Hierarchical reduction. For all-reduce and reducescatter the hierarchical path first sums each GPU package’s two GPU stacks over MDFI, then reduces the partials from the six GPU packages over XeLink, halving the XeLink message count (the copy now carries already-summed data). • Destination-split fused EU copy. For all-to-all, one kernel per source GPU pushes all 𝑁 chunks with byvalue dispatch tables, while each workgroup owns one destination. This removes the per-message launch floor and exposes independent remote endpoints concurrently. The topology-aware all-gather uses the corresponding fused latency path.

The all-gather buffer holds six KV shards per GPU in this configuration. Ring stores only a constant number of shards, including its rotation buffers, so its KV memory per GPU scales with local rather than total context length. I.5

Correctness

Every functional method is validated against an aten fusedflash full-context reference (relative error, bf16 in, fp32 accumulate). Representative maxima across the three nodes appear below. All functional configurations pass (tolerance 0.05); the error falls with 𝑁 per as the softmax normalizer averages more terms. SyclKittens ring attention is validated by its own in-binary reference check on each shard, and SyclKittens all-gather attention is checked in-binary against an onlinesoftmax full-context reference (relative error 0.0012 at 𝑁 per = 256, well under tolerance) before any timing is recorded.

J

The reduction and gather kernels dispatch adaptively between latency- and bandwidth-oriented paths. All-to-all uses the destination-split path at every size because it is neutral at the launch floor and wins as endpoint traffic grows.

Multi-GPU Collectives

The SyclKittens collectives reduce per-call launch and handshake overhead. In a six-node matrix, all four win at the smallest standalone payload, 24 KiB, on every node. The evaluated batch-one decode uses two 8-KiB all-reduces per layer; this smaller payload is exercised in the end-to-end measurements in §K.5, not in the standalone sweep. Reduce-scatter also wins 0.19 MiB on every node, and destination-split allto-all wins every measured payload on every node. Most other medium cells are node-sensitive. At 96 MiB, all-to-all is 1.20× faster on the median and wins all six nodes; oneCCL retains the median lead for all-reduce and all-gather, while reduce-scatter favors SyclKittens on four of six nodes. These paths also support distributed GEMM (Appendix H). J.1

Kernel Implementation and Key Techniques

J.2.1 Multi-GPU API. The collective kernels use the PGL peer-pointer table, Intel Max GPU topology planner, and system-scope atomic synchronization introduced in the tilelanguage appendix. These replace ThunderKittens’ CUDA NVLS multicast lowering while preserving explicit control over slice ownership, transfer direction, and schedule. We validated the API conversion of all-gather, reducescatter, all-reduce, and all-to-all against immutable preconversion binaries with complementary raw–API–API–raw and API–raw–raw–API orderings. Equal-weighting the two per-node bracket estimates gives an API/raw geomean of 1.013× over all 28 operation–payload cells; the worst cell median is 0.974×. Every cell passes the predefined 0.95× guard and the aggregate passes the 0.99× guard, establishing no material regression under the declared protocol. A separate clean current-source recollection rebuilt all four production kernels and completed the full SyclKittens– oneCCL–SyclKittens matrix. All 336 SyclKittens rows passed

Configuration

We compare bf16 all-reduce, all-gather, reduce-scatter, and all-to-all on all 12 GPUs of one node containing six GPU packages. SyclKittens uses one process with shared-context peer pointers; the tuned vendor library path is oneCCL/xCCL with one process per GPU. Each SyclKittens collective 40

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

full-buffer validation with zero mismatches, confirming the production build and oneCCL integration after conversion.

Payload

Latency (µs)

Bandwidth (GB/s) Speedup

SyclKittens oneCCL SyclKittens oneCCL

J.3

Methodology

Each same-node bracket runs SyclKittens and oneCCL on the same contiguous 12 GPUs and traverses the same physical MDFI and XeLink links under the same topology. The formal matrix repeats that bracket on six nodes. Each pass reports the median of 15 timed iterations after five warmups. For oneCCL, all_reduce(MAX) selects the slowest rank as the collective completion time. All-reduce bandwidth uses 2(𝑁 −1)/𝑁 · 𝑆/𝑡. All-gather, reduce-scatter, and all-to-all use the p2p-bytes convention 𝑁 (𝑁 −1) · 𝑒𝑝𝑑 · 2/𝑡. The SyclKittens binary and its oneCCL counterpart print the same formula for aligned size sweeps, making their GB/s values directly comparable. Every node runs SyclKittens–oneCCL–SyclKittens; the SyclKittens value is the mean of its two bracketing medians, and speedup is paired on that node before taking the six-node median. Every SyclKittens point is full-buffer verified against a reference sum, gather, or scatter; maximum absolute error remains within bf16 rounding. J.4

Tables report six-node marginal median 𝜇s (12 GPUs), the matched-convention GB/s for both, and the median paired same-node speedup (oneCCL divided by SyclKittens; greater than 1 means SyclKittens is faster). Absolute performance tiers vary by node, so the paired speedup need not equal the ratio of the two marginal median times. All-Reduce (busbw 2(𝑁 −1)/𝑁 ). Latency (µs)

92.8 108.0 144.5 227.2 259.2 368.5 987.8

161.5 160.5 180.0 245.0 285.0 378.0 903.0

Payload

Latency (µs)

2.9 20.0 119.7 304.5 533.9 751.2 1 121.0

1.7 13.5 96.1 282.5 485.7 732.3 1 226.2

1.81× 1.51× 1.26× 1.10× 1.12× 1.03× 0.92×

Bandwidth (GB/s) Speedup

SyclKittens oneCCL SyclKittens oneCCL

Delivered Performance

Payload

24 KiB 0.19 MiB 1.5 MiB 6 MiB 12 MiB 24 MiB 96 MiB

24 KiB 0.19 MiB 1.5 MiB 6 MiB 12 MiB 24 MiB 96 MiB

80.8 80.8 110.5 223.2 249.5 315.0 713.5

238.5 238.5 265.5 265.5 283.5 357.5 781.0

Payload

Latency (µs)

3.3 26.8 156.6 310.0 554.8 878.8 1 551.9

1.1 9.1 65.2 260.7 488.2 774.3 1 417.8

2.19× 2.20× 1.97× 1.22× 1.21× 1.19× 1.10×

Bandwidth (GB/s) Speedup

SyclKittens oneCCL SyclKittens oneCCL 24 KiB 0.19 MiB 1.5 MiB 6 MiB 12 MiB 24 MiB 96 MiB

91.8 91.0 90.2 136.5 207.5 342.0 1 188.5

293.0 302.5 302.0 355.0 464.0 612.5 1 421.5

2.9 23.8 191.7 507.0 667.0 809.4 931.7

0.9 7.1 57.3 194.9 298.3 452.0 779.0

3.31× 3.46× 3.49× 2.62× 2.24× 1.79× 1.20×

Bandwidth (GB/s) Speedup

SyclKittens oneCCL SyclKittens oneCCL 24 KiB 0.19 MiB 1.5 MiB 6 MiB 12 MiB 24 MiB 96 MiB

81.2 131.8 208.2 264.5 352.0 516.8 1 511.8

173.5 172.5 274.5 319.0 417.5 596.0 1 460.5

0.6 2.7 13.8 43.6 65.5 89.3 122.1

0.3 2.1 10.5 36.2 55.3 77.4 126.4

2.14× 1.16× 1.06× 1.04× 1.04× 1.12× 0.94×

per-collective crossover from 24 KiB through 96 MiB; the smaller 8-KiB decode payload is evaluated in the complete decode path. The decode orchestrator also replaces 64–96 per-token host waits with a system-scope atomic peer-counter barrier on the GPUs. This removes host serialization and improves multi-GPU decode TPOT by 1.1–1.54×, with the largest change on the 2D TP×CP layout. The synchronization change does not alter arithmetic and passes the same correctness checks. It is used for the decode results in §K.5.

All-Gather (p2p bytes). Reduce-Scatter (p2p bytes). All-to-All (p2p bytes). The crossover follows the dominant cost. At 24 KiB, all four SyclKittens collectives win all six nodes; reduce-scatter extends that robust win to 0.19 MiB and all-to-all extends it through 96 MiB. The medium AG/RS/AR ratios favor SyclKittens in the six-node median but cross parity between nodes. At 96 MiB, oneCCL wins AG, and AR in the median, RS favors SyclKittens on four of six nodes, and A2A favors SyclKittens on all six nodes. Thus the standalone sweep establishes the

J.5

Correctness

Every SyclKittens collective is full-buffer verified on the participating GPUs against a reference on each run and each size, with all-reduce and reduce-scatter checked against an fp32 elementwise sum, all-gather and all-to-all against the exact gathered/permuted layout. Maximum absolute error is at most 0.01, within bf16 rounding. The oneCCL baseline uses the torch.distributed collective. 41

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

All-reduce

All-gather

Reduce-scatter

All-to-all

Speedup over oneCCL

4

3

2

1

iB

24 K

0.19

MiB

1.5

MiB

iB

6M

iB iB 12 M 24 M

iB 96 M

Payload per GPU Figure 11. BF16 collective speedup over oneCCL on 12 Intel Max GPUs (higher is better). SyclKittens uses one process with shared-context peer pointers; oneCCL uses one process per GPU. Curves are six-node medians of paired same-node brackets. SyclKittens wins all six nodes for every collective at 24 KiB and for all-to-all through 96 MiB; medium AG/RS/AR cells vary by node, and oneCCL leads AG and AR at 96 MiB.

K

End-to-End Llama-3.1-8B

K.2

On one GPU, composing SyclKittens operations with TorchInductor-generated regions improves Llama-3.1-8B full-request throughput by 1.59× in geometric mean over torch.compile alone on the 21 configurations that both paths run. The composed path also serves the complete 𝑁 = 8192 grid, including cells where the baselines run out of memory. The final subsection evaluates multi-GPU decode separately. K.1

Eager PyTorch executes the operator sequence directly. The torch.compile baseline fuses pointwise and reduction regions but graph-breaks at SDPA on this backend; all primary ratios compare the composed path with Inductor alone under identical compilation settings. K.3

Table 56 lists the implementations used for the measured end-to-end path. Each selection is the fastest correct implementation in its operation-specific study, and every run records the implementation that executes.

K.4

End-to-End Composition with TorchInductor

The integration boundary is the SDPA graph break. SyclKittens replaces regions that span this boundary, while TorchInductor continues to generate the remaining pointwise, reduction, and epilogue regions. Stand-alone norm and SiLU remain with Inductor because it matches or exceeds the corresponding stand-alone SyclKittens kernels.

Table 56. Evaluated kernel set for the end-to-end Llama-3.18B path. Implementation

Operations across the SDPA Graph Break

SyclKittens supplies four operations across the attention graph break, namely GQA flash attention, joint Q/K rotary, a concatenated QKV projection, and residual-add plus RMSNorm. Inductor retains RMSNorm, SwiGLU, and the remaining pointwise and epilogue regions.

Evaluated Kernel Set

Block op

Baselines

Notes

Full-request metric. Each repeat synchronizes the GPU, times one complete deterministic model.generate() call with max_new_tokens=128, and synchronizes again after return. For batch size 𝐵, prompt length 𝑁 , mean actual generated sequence length 𝐺¯ ≤ 128, and mean synchronized elapsed time 𝑡¯ in seconds over the timed repeats, the har¯ 𝑡¯. This convention combines batchness records (𝐵𝑁 + 𝐺)/ wide prompt tokens with the mean generated length of one stream. We retain it for the reported full-request rates and

Residual+RMSNorm Fused SyclKittens residual 1.04× Inductor and RMSNorm GQA (𝑔=4) Fence-free SyclKittens GQA shape-tuned MHA Fence-free SyclKittens MHA non-grouped fallback Fused SyclKittens Q/K rotary fusion Q/K rotary microbenchmark QKV projection concatenated 𝑊𝑞𝑘𝑣 one prefill GEMM

42

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

Table 58. Per-shape prefill throughput (tok/s).

ratios. Total batch-token throughput for equally long out¯ 𝑡¯. When compared runs puts would instead be 𝐵(𝑁 + 𝐺)/ generate the same lengths, the reported ratio equals their inverse elapsed-time ratio. The timed interval includes prefill and every generated decode step; compilation occurs during untimed warmup.

BS

1 256 1 512 1 1024 1 2048 1 4096 1 8192 2 256 2 512 2 1024 2 2048 2 4096 2 8192 4 256 4 512 4 1024 4 2048 4 4096 4 8192 8 256 8 512 8 1024 8 2048 8 4096 8 8192

K.4.1 Configuration. Both composed and Inductor-alone runs use identical compilation settings, eight untimed warmups, and five synchronized repeats. The only difference is the four SyclKittens operations above; Inductor owns RMSNorm, SwiGLU, and all other pointwise and epilogue code generation in both runs. K.4.2 Delivered performance. We evaluate bf16 Llama3.1-8B on one GPU over 24 shapes with 𝐵 ∈ {1, 2, 4, 8} and 𝑁 ∈ {256, 512, 1024, 2048, 4096, 8192}. Composed, Inductoralone, and eager runs are node-matched; the primary ratio is composed over Inductor alone. oom marks a baseline capacity failure. Table 57. Geometric-mean summary for the 24-shape singleGPU grid. Ratios include only cells where both compared paths run. comp/ind uses 21 cells; Inductor runs out of memory at 𝑁 = 8192 with 𝐵 = 2, 4, 8. The composed path runs all 24 cells, and eager also runs the 𝐵 = 2 cell. Regime Prefill Decode End-to-end

vs. inductor-alone

vs. eager

Range (comp/ind)

1.22× 1.65× 1.59×

1.40× 1.58× 1.55×

0.91×–2.03× 1.30×–2.20× 1.29×–1.98×

N eager inductor composed comp/ind comp/eager 5340 6890 7470 7160 5990 3510 6890 7540 7720 7060 5860 3700 7890 8220 7980 7090 5390 oom 8550 8640 8090 7090 5520 oom

6720 8450 8920 8340 6790 3730 8140 8780 8580 7510 6130 oom 9520 9800 9060 7810 4930 oom 10100 10100 9260 7740 5040 oom

6140 8700 9490 9750 9300 7590 8350 9210 9550 9550 9460 8110 9660 10000 10400 10000 9430 8250 10300 10500 10500 10100 9500 7740

0.91× 1.03× 1.06× 1.17× 1.37× 2.03× 1.03× 1.05× 1.11× 1.27× 1.54× – 1.02× 1.02× 1.14× 1.28× 1.91× – 1.02× 1.05× 1.14× 1.30× 1.89× –

1.15× 1.26× 1.27× 1.36× 1.55× 2.16× 1.21× 1.22× 1.24× 1.35× 1.61× 2.19× 1.22× 1.22× 1.30× 1.41× 1.75× – 1.21× 1.22× 1.30× 1.42× 1.72× –

Table 59. Decode throughput per stream (tok/s). BS

N eager inductor composed comp/ind comp/eager

1 256 1 512 1 1024 1 2048 1 4096 1 8192 2 256 2 512 2 1024 2 2048 2 4096 2 8192 4 256 4 512 4 1024 4 2048 4 4096 4 8192 8 256 8 512 8 1024 8 2048 8 4096 8 8192

For prefill at prompt length 256, the composed/Inductor ratio is 0.91×–1.03×. The projection GEMMs are too small to amortize the fused QKV, rotary, and residual launches, with 0.91× at the 1 × 256 cell. From prompt length 512 onward, the ratio increases with context, reaching 1.30× at 8 × 2048, 1.91× at 4 × 4096, and 2.03× at 1 × 8192 as the flash-attention path replaces the baseline attention path across the SDPA graph break. Decode improves by 1.65× geometrically and reaches 1.96× at 1 × 256. Decode is dispatch-bound; the SDPA graph break prevents Inductor from fusing the per-step attention prologue and epilogue, whereas the composed path uses fewer launches. The harness retains only last-token logits, preventing a full-vocabulary all-position logits tensor from determining memory capacity. At 𝑁 = 8192, torch.compile runs out of memory at 𝐵𝑆 ≥ 2 and eager at 𝐵𝑆 ≥ 4. The composed path serves BS 1–8 at the prefill throughputs in Table 58 because its flash attention has the smaller memory footprint. Across the 𝑁 ≤ 2048 subset, the composed path retains a 1.56× end-to-end geometric mean. The gain comes from combining Inductor’s pointwise/norm fusion with the crossop and attention fusions that span its SDPA graph break; the

29 30 30 30 28 22.7 29 29 29 28 23 15 30 30 28 23 15 oom 30 28 23 15 9 oom

25 25 25 26 29 24.9 25 25 26 25 25 oom 26 26 25 25 17 oom 26 25 25 21 10 oom

49 48 46 44 39 32.5 47 46 44 40 35 27 44 43 40 36 29 24 42 40 36 31 22 16

1.96× 1.90× 1.82× 1.71× 1.37× 1.30× 1.87× 1.81× 1.72× 1.58× 1.40× – 1.74× 1.67× 1.59× 1.44× 1.73× – 1.64× 1.58× 1.43× 1.48× 2.20× –

1.67× 1.62× 1.54× 1.49× 1.38× 1.43× 1.62× 1.56× 1.50× 1.43× 1.53× 1.82× 1.49× 1.43× 1.43× 1.57× 1.87× – 1.43× 1.43× 1.57× 2.00× 2.32× –

flash-attention memory footprint also preserves coverage at long context. 43

Yehong Jiang, Sheng Chen, Fangwen Fu, Yen-Kuang Chen, Xinmin Tian, Stuart H. Sul, and Simran Arora

¯ 𝑡¯, using Table 60. Per-shape full-request rate, (𝐵𝑁 + 𝐺)/ batch-wide prompt tokens and generated length per stream. All columns are node-matched. BS

N eager inductor composed comp/ind comp/eager

1 256 1 512 1 1024 1 2048 1 4096 1 8192 2 256 2 512 2 1024 2 2048 2 4096 2 8192 4 256 4 512 4 1024 4 2048 4 4096 4 8192 8 256 8 512 8 1024 8 2048 8 4096 8 8192

K.5

for the five projections, inline peer all-reduce, and a GPUside atomic rendezvous. The entire tile-by-layer loop is one C++ orchestrator call, so Python dispatch is removed from the inner loop. The C++ orchestrator still submits GPU work, and its cost is included in the measured decode time.

253 472 874 1520 2280 2230 465 863 1530 2400 2930 2480 880 1580 2500 3230 3200 oom 1570 2540 3410 3730 3510 oom

218 409 768 1390 2400 2360 411 766 1400 2290 3140 oom 775 1430 2400 3540 3130 oom 1440 2450 3810 4480 3370 oom

413 749 1320 2210 3290 4070 738 1310 2200 3360 4610 5100 1280 2180 3460 4830 5760 6190 2170 3450 4930 6140 6670 6270

1.89× 1.83× 1.72× 1.59× 1.37× 1.72× 1.80× 1.70× 1.58× 1.47× 1.47× – 1.65× 1.53× 1.44× 1.36× 1.84× – 1.51× 1.41× 1.29× 1.37× 1.98× –

K.5.4 TPOT versus context. Table 61 reports steadystate time-per-output-token (TPOT), the geometric mean over 100 decode steps after prefilling to length 𝑁 ; prefill is excluded. The p95 and p99 remain within 2 ms of the geometric mean throughout, and the medians coincide. Across matched two- and four-GPU configurations, SyclKittens is 1.52–2.91× faster than the oneCCL-based decode path over the contexts where both run.

1.64× 1.59× 1.51× 1.46× 1.44× 1.83× 1.59× 1.51× 1.44× 1.40× 1.57× 2.06× 1.45× 1.39× 1.39× 1.49× 1.80× – 1.38× 1.36× 1.45× 1.65× 1.90× –

Table 61. Multi-GPU decode TPOT (ms/token, geometric mean over 100 steps; lower is better). oom marks a baseline capacity failure. All SyclKittens configurations pass the first-new-token logit-vector check against the single-GPU reference with cosine similarity at least 0.99997. ctx SyclKittens 𝑃 = 1 SyclKittens 𝑃 = 2 SyclKittens 𝑃 = 4 xCCL 𝑤𝑠 = 2 xCCL 𝑤𝑠 = 4 1K 2K 4K 8K 16K 32K 64K 128K 256K

Distributed Multi-GPU Decode (Single-Process TP+CP vs. oneCCL)

28.4 30.0 30.6 31.3 43.4 54.1 117.1 209.7 –

22.8 24.6 25.7 26.1 27.4 36.2 73.3 128.0 370.7

19.2 20.3 21.7 22.4 23.1 25.3 54.5 95.6 271.9

42.5 43.0 43.3 49.2 66.6 104.6 173.9 310.8 oom

40.7 40.5 40.7 41.8 49.3 73.7 128.1 222.3 413.0

At 32K context SyclKittens on 2 GPUs is 2.89× faster than oneCCL on 2 (36.2 vs. 104.6 ms) and on 4 GPUs 2.91× (25.3 vs. 73.7 ms); equivalently SyclKittens sustains 39.5 tok/s on 4 GPUs at 32K versus oneCCL’s 13.6. At 128K the speedups are 2.43× (𝑃 = 2) and 2.33× (𝑃 = 4); at 256K, 𝑃 = 4 is 1.52× over xCCL 𝑤𝑠 = 4 and 𝑃 = 2 runs at all where xCCL 𝑤𝑠 = 2 runs out of memory.

We measure Llama-3.1-8B decode on one node with six GPU packages, each with two stacks connected by MDFI. XeLink connects the GPU packages all-to-all through independent even and odd stack planes. A flat device hierarchy exposes all 12 stacks. We compare two implementations of tensor parallelism [29]. SyclKittens uses a single-process C++ orchestrator. The tuned Intel oneCCL library path uses torch.distributed with the xccl backend and one process per stack.

K.5.5 Per-layer collective cost. Each layer performs two all-reduces as inline single-process peer pulls over the shared context. A separate category-level profile assigns all-reduces about 1.7 ms/token at 𝑃 = 2 and 3.1 ms/token at 𝑃 = 4. In that instrumented profile, they account for 3.4–3.7% at short context and less than 1.5% after the KV scan dominates. The profiler synchronizes around categories, so these percentages describe its instrumented step breakdown, not the final overlapped TPOT table. At fixed GPU count, this collective cost stays nearly constant as context grows, so its fraction of TPOT falls. The oneCCL path invokes a cross-process collective at each of the two reduction sites per layer; its two-GPU TPOT rises from 43 ms at 1K to 105 ms at 32K, while the SyclKittens path rises from 23 to 36 ms.

K.5.1 Configuration and compared implementations. The model has hidden size 4096, 32 query and 8 KV heads, head dimension 128, intermediate size 14336, and 32 layers. The node exposes six GPU packages and 12 stacks. Both systems use the same TP weight and KV sharding; SyclKittens uses inline peer collectives and can additionally shard KV by context. K.5.2 Parallelization strategy. Decode combines a context-invariant weight stream with a KV scan that grows linearly with context. Tensor parallelism shards the weight stream, whereas context parallelism shards the growing KV term. The 2-D configurations below combine these schemes.

K.5.6 GPU count and capacity. Flat TP=4 has the lowest TPOT across the measured TP-only range. The 2-D measurements below evaluate longer contexts. The 128-KB/token

K.5.3 What SyclKittens runs per token. Each layer uses split-K flash decode for the KV scan, 128-bit wide-load GEMV 44

SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs

The lowest-TPOT mapping is flat 𝑃 = 4 through 32K, TP = 2 × CP = 3 on six GPUs at 64K–128K, and TP = 4 × CP = 3 on 12 GPUs at 256K. The six-GPU path measures 29.4, 31.0, and 45.4 ms at 32K, 64K, and 128K. At 64K and 128K it is respectively 1.76× and 2.11× faster than flat 𝑃 = 4. The 12-GPU path remains at 59.0–60.0 ms through 128K, then reaches 73.3 ms at 256K, compared with 93.9 ms for six GPUs and 271.9 ms for flat 𝑃 = 4. Both 2-D layouts use CP = 3; the additional TP sharding affects the largest KV scan. With these additional GPUs, the measured ratio over the oneCCL-based decode path is 4.90× at 128K (6-GPU 45.4 vs. 𝑤𝑠 = 4’s 222.3) and 5.63× at 256K (12-GPU 73.3 vs. 413.0, where 𝑤𝑠 = 2 runs out of memory). These are scaling comparisons against four GPUs, not matched-resource speedups. Both selected layouts use the same single-process collective and GPU-side barrier.

KV cache exceeds one-GPU capacity at approximately 384K context, after which sharding is required for capacity. K.5.7 2D scaling to the full node (TP × CP). Plain TP is capped at 𝑃 = 8 and loses efficiency before that limit. Context parallelism (CP) instead shards the KV cache by sequence position, computes a local flash partial on each GPU, and merges the softmax statistics across GPUs without a headcount limit. A standalone CP-attention prototype accelerates the 128K KV scan by 6.6× across 4 GPUs and 8.5× across 8, with maximum absolute error 10−4 . The full 2-D orchestrator shards weights by TP and replicates them across CP, so CP leaves per-GPU GEMM work unchanged but replicates that work across CP groups. It shards the KV cache by sequence, merges the CP softmax partials within each TP group, and performs the two TP all-reduces per layer. One process drives all 12 GPUs.

K.5.8 Correctness methodology. The correctness gate uses first-new-token logit-vector cosine similarity against the single-GPU (𝑃 = 1) run. Greedy token matching varies across runs because fused SDPA reduction order changes near-tied argmax logits. In the broader validation sweep, which also includes eight GPUs, 𝑃 ∈ {2, 4, 8} configurations reproduce 𝑃 = 1 to cosine 0.99998 at dense context and 0.99973 at 128K–256K. The maximum absolute logit difference is 0.125, corresponding to bf16 rounding.

Table 62. Measured 2D end-to-end decode TPOT (ms/token, geomean over 20 steps). Full end-to-end TP × CP decode, all layers, one Intel Max GPU node. Every configuration passes the first-new-token logit-vector check against the single-GPU reference with cosine similarity at least 0.9999. Compared against the best flat-TP point and oneCCL. ctx 6-GPU TP2 × CP3 12-GPU TP4 × CP3 best flat TP (𝑃 = 4) xCCL 𝑤𝑠 = 4 32K 64K 128K 256K

29.4 31.0 45.4 93.9

59.8 59.0 60.0 73.3

25.3 54.5 95.6 271.9

73.7 128.1 222.3 413.0

45

Related documents

• SyclKittens: A Tile Programming Model for Programmers and Coding Agents on Intel GPUs#1199551
• Advanced Computer Architecture#1190332
• CodeSpec: Dual Executable Specifications for Agentic Long-Horizon Feature Development#411107
• Exo-GPU: Safe, Imperative, User-schedulable Programming for Tensor Cores#1110750
• Credo: Reusable Declarative Primitives for Agentic Workflows#781081