In this tutorial, we explore TileGym GPU programming by building a practical Colab workflow that runs across different hardware conditions. We begin by probing the available CUDA environment, checking whether NVIDIA cuTile runs directly, and falling back to Triton when standard Colab GPUs lack the required cuTile stack. Through this setup, we learn the core tile-programming idea: instead of writing code for one thread at a time, we operate on entire data tiles, load them into the kernel, compute on them efficiently, and store the results back. We use this model to implement vector addition, fused GELU, row-wise softmax, tiled matrix multiplication, and flash attention, while comparing each result against PyTorch for correctness and benchmarking.

CUDA Environment Probe

import os, sys, math, time, textwrap
def rule(t=""):
   print("\n" + "=" * 78)
   if t: print(t)
   print("=" * 78)
rule("0. ENVIRONMENT PROBE")
try:
   import torch
except ImportError:
   print("Installing torch ...")
   os.system(f"{sys.executable} -m pip install -q torch")
   import torch
HAS_CUDA = torch.cuda.is_available()
DEV = "cuda" if HAS_CUDA else "cpu"
cc = (0, 0)
if HAS_CUDA:
   cc = torch.cuda.get_device_capability()
   print(f"GPU                 : {torch.cuda.get_device_name(0)}")
   print(f"Compute capability  : {cc[0]}.{cc[1]}")
   print(f"Torch CUDA runtime  : {torch.version.cuda}")
   print(f"Driver / torch      : {torch.__version__}")
else:
   print("No CUDA GPU found. In Colab: Runtime -> Change runtime type -> GPU (T4).")
   print("The tutorial will still run its correctness math on CPU where possible.")
CUTILE_HW_OK = HAS_CUDA and cc[0] >= 8
CUDA_MAJOR = int((torch.version.cuda or "0").split(".")[0]) if HAS_CUDA else 0
CUTILE_TOOLKIT_OK = CUDA_MAJOR >= 13
rule("1. ATTEMPTING REAL cuTile (NVIDIA CUDA Tile) BACKEND")
ct = None
CUTILE_READY = False
if CUTILE_HW_OK and CUTILE_TOOLKIT_OK:
   try:
       import cuda.tile as ct
       CUTILE_READY = True
       print("cuda.tile is already importable.")
   except Exception:
       print("Installing cuda-tile[tileiras] (this can take a while)...")
       os.system(f"{sys.executable} -m pip install -q 'cuda-tile[tileiras]' cupy-cuda13x")
       try:
           import cuda.tile as ct
           CUTILE_READY = True
       except Exception as e:
           print("cuTile import still failed:", repr(e))
else:
   reasons = []
   if not HAS_CUDA:            reasons.append("no CUDA GPU")
   if HAS_CUDA and cc[0] < 8:  reasons.append(f"compute capability {cc[0]}.{cc[1]} < 8.0 (Turing/T4 unsupported)")
   if not CUTILE_TOOLKIT_OK:   reasons.append(f"CUDA {torch.version.cuda} < 13.1 required by tileiras")
   print("Skipping real cuTile install because:", "; ".join(reasons) + ".")
   print("This is expected on a standard Colab T4 — we fall back to Triton below,")
   print("which teaches the exact same tile-based programming model.")
if CUTILE_READY:
   BACKEND = "cutile"
else:
   try:
       import triton, triton.language as tl
       BACKEND = "triton" if HAS_CUDA else "torch"
   except ImportError:
       if HAS_CUDA:
           print("Installing triton ...")
           os.system(f"{sys.executable} -m pip install -q triton")
           try:
               import triton, triton.language as tl
               BACKEND = "triton"
           except Exception:
               BACKEND = "torch"
       else:
           BACKEND = "torch"
rule(f"ACTIVE EXECUTION BACKEND:  {BACKEND.upper()}")
print({
   "cutile": "Running NVIDIA cuTile kernels on your Ampere+/CUDA13 GPU. Nice hardware!",
   "triton": "Running Triton tile kernels on your GPU (the standard Colab path).",
   "torch":  "No usable GPU kernel backend; showing reference math on CPU only.",
}[BACKEND])
print(textwrap.dedent("""
   ------------------------------------------------------------------
   SIMT (classic CUDA)          |   TILE model (cuTile / Triton)
   ------------------------------------------------------------------
   You write code for ONE       |   You write code for ONE BLOCK that
   thread. You compute a global |   owns a whole TILE (e.g. 1024 elems
   index, bounds-check it, and   |   or a 128x128 sub-matrix). You load
   touch a single element.      |   the tile, do math on the WHOLE tile,
                                |   store it. The compiler maps the tile
   C[i] = A[i] + B[i]           |   onto threads / tensor cores for you.
   ------------------------------------------------------------------
   cuTile primitives:  ct.bid(0), ct.load(...), ct.store(...), a @ b, ct.launch
   Triton primitives:  tl.program_id, tl.load, tl.store, tl.dot, grid[...]
   Same idea, two spellings. Below, every kernel is shown in BOTH.
   """))
CUTILE_SOURCE = {
"vector_add": '''
import cuda.tile as ct
@ct.kernel
def vector_add(a, b, c, tile_size: ct.Constant[int]):
   pid    = ct.bid(0)
   a_tile = ct.load(a, index=(pid,), shape=(tile_size,))
   b_tile = ct.load(b, index=(pid,), shape=(tile_size,))
   ct.store(c, index=(pid,), tile=a_tile + b_tile)
''',
"matmul": '''
import cuda.tile as ct
@ct.kernel
def matmul(A, B, C, K: ct.Constant[int],
          BM: ct.Constant[int], BN: ct.Constant[int], BK: ct.Constant[int]):
   m, n = ct.bid(0), ct.bid(1)
   acc  = ct.zeros((BM, BN), dtype=ct.float32)
   for k in range(ct.cdiv(K, BK)):
       a = ct.load(A, index=(m, k), shape=(BM, BK))
       b = ct.load(B, index=(k, n), shape=(BK, BN))
       acc = a @ b + acc
   ct.store(C, index=(m, n), tile=acc)
'''}

We begin by setting up the environment, importing the required libraries, and checking whether CUDA is available on the current runtime. We inspect the GPU capabilities, CUDA version, and PyTorch setup to determine whether the real cuTile backend is usable. We then select the active execution backend, explain the tile programming model, and store reference cuTile kernel source strings for comparison.

Defining Triton Kernels

if BACKEND == "triton":
   @triton.jit
   def _vadd_kernel(a_ptr, b_ptr, c_ptr, n, BLOCK: tl.constexpr):
       pid  = tl.program_id(0)
       offs = pid * BLOCK + tl.arange(0, BLOCK)
       mask = offs < n
       a = tl.load(a_ptr + offs, mask=mask)
       b = tl.load(b_ptr + offs, mask=mask)
       tl.store(c_ptr + offs, a + b, mask=mask)
   @triton.jit
   def _fused_gelu_kernel(x_ptr, w_ptr, b_ptr, o_ptr, n, BLOCK: tl.constexpr):
       pid  = tl.program_id(0)
       offs = pid * BLOCK + tl.arange(0, BLOCK)
       mask = offs < n
       x = tl.load(x_ptr + offs, mask=mask)
       w = tl.load(w_ptr + offs, mask=mask)
       b = tl.load(b_ptr + offs, mask=mask)
       h = x * w + b
       c = 0.7978845608028654
       z = c * (h + 0.044715 * h * h * h)
       e = tl.exp(-2.0 * z)
       tanh = (1.0 - e) / (1.0 + e)
       g = 0.5 * h * (1.0 + tanh)
       tl.store(o_ptr + offs, g, mask=mask)
   @triton.jit
   def _softmax_kernel(x_ptr, o_ptr, stride, n_cols, BLOCK: tl.constexpr):
       row  = tl.program_id(0)
       cols = tl.arange(0, BLOCK)
       mask = cols < n_cols
       ptr  = x_ptr + row * stride + cols
       x    = tl.load(ptr, mask=mask, other=-float("inf"))
       x    = x - tl.max(x, axis=0)
       num  = tl.exp(x)
       den  = tl.sum(num, axis=0)
       tl.store(o_ptr + row * stride + cols, num / den, mask=mask)
   @triton.jit
   def _matmul_kernel(A, B, C, M, N, K,
                      sam, sak, sbk, sbn, scm, scn,
                      BM: tl.constexpr, BN: tl.constexpr, BK: tl.constexpr):
       pid_m = tl.program_id(0)
       pid_n = tl.program_id(1)
       offs_m = pid_m * BM + tl.arange(0, BM)
       offs_n = pid_n * BN + tl.arange(0, BN)
       offs_k = tl.arange(0, BK)
       a_ptr = A + offs_m[:, None] * sam + offs_k[None, :] * sak
       b_ptr = B + offs_k[:, None] * sbk + offs_n[None, :] * sbn
       acc = tl.zeros((BM, BN), dtype=tl.float32)
       for k in range(0, K, BK):
           a = tl.load(a_ptr, mask=offs_k[None, :] < K - k, other=0.0)
           b = tl.load(b_ptr, mask=offs_k[:, None] < K - k, other=0.0)
           acc += tl.dot(a, b)
           a_ptr += BK * sak
           b_ptr += BK * sbk
       c_ptr = C + offs_m[:, None] * scm + offs_n[None, :] * scn
       cmask = (offs_m[:, None] < M) & (offs_n[None, :] < N)
       tl.store(c_ptr, acc.to(C.dtype.element_ty), mask=cmask)
   @triton.jit
   def _flash_kernel(Q, K, V, O, sqz, skz, svz, soz,
                     L, D, scale,
                     BL: tl.constexpr, BD: tl.constexpr):
       pid_l = tl.program_id(0)
       z     = tl.program_id(1)
       offs_l = pid_l * BL + tl.arange(0, BL)
       offs_d = tl.arange(0, BD)
       q_ptr = Q + z * sqz + offs_l[:, None] * D + offs_d[None, :]
       q = tl.load(q_ptr, mask=offs_l[:, None] < L, other=0.0)
       m_i = tl.full((BL,), -float("inf"), dtype=tl.float32)
       l_i = tl.zeros((BL,), dtype=tl.float32)
       acc = tl.zeros((BL, BD), dtype=tl.float32)
       for start in range(0, L, BL):
           offs_k = start + tl.arange(0, BL)
           k_ptr = K + z * skz + offs_k[:, None] * D + offs_d[None, :]
           v_ptr = V + z * svz + offs_k[:, None] * D + offs_d[None, :]
           k = tl.load(k_ptr, mask=offs_k[:, None] < L, other=0.0)
           v = tl.load(v_ptr, mask=offs_k[:, None] < L, other=0.0)
           s = tl.dot(q, tl.trans(k)) * scale
           s = tl.where(offs_k[None, :] < L, s, -float("inf"))
           m_ij = tl.maximum(m_i, tl.max(s, axis=1))
           p    = tl.exp(s - m_ij[:, None])
           alpha = tl.exp(m_i - m_ij)
           l_i  = l_i * alpha + tl.sum(p, axis=1)
           acc  = acc * alpha[:, None] + tl.dot(p.to(v.dtype), v)
           m_i  = m_ij
       acc = acc / l_i[:, None]
       o_ptr = O + z * soz + offs_l[:, None] * D + offs_d[None, :]
       tl.store(o_ptr, acc.to(O.dtype.element_ty), mask=offs_l[:, None] < L)
   def run_vadd(a, b):
       c = torch.empty_like(a); n = a.numel()
       grid = (triton.cdiv(n, 1024),)
       _vadd_kernel[grid](a, b, c, n, BLOCK=1024)
       return c
   def run_fused_gelu(x, w, b):
       o = torch.empty_like(x); n = x.numel()
       grid = (triton.cdiv(n, 1024),)
       _fused_gelu_kernel[grid](x, w, b, o, n, BLOCK=1024)
       return o
   def run_softmax(x):
       m, ncols = x.shape
       o = torch.empty_like(x)
       BLOCK = triton.next_power_of_2(ncols)
       _softmax_kernel[(m,)](x, o, x.stride(0), ncols, BLOCK=BLOCK)
       return o
   def run_matmul(a, b):
       M, K = a.shape; K2, N = b.shape
       c = torch.empty((M, N), device=a.device, dtype=a.dtype)
       BM = BN = 64; BK = 32
       grid = (triton.cdiv(M, BM), triton.cdiv(N, BN))
       _matmul_kernel[grid](a, b, c, M, N, K,
                            a.stride(0), a.stride(1), b.stride(0), b.stride(1),
                            c.stride(0), c.stride(1), BM=BM, BN=BN, BK=BK)
       return c
   def run_flash(q, k, v):
       Z, L, D = q.shape
       o = torch.empty_like(q)
       scale = 1.0 / math.sqrt(D)
       BL = 64
       grid = (triton.cdiv(L, BL), Z)
       _flash_kernel[grid](q, k, v, o,
                           q.stride(0), k.stride(0), v.stride(0), o.stride(0),
                           L, D, scale, BL=BL, BD=D)
       return o
else:
   def run_vadd(a, b):            return a + b
   def run_fused_gelu(x, w, b):   return torch.nn.functional.gelu(x * w + b, approximate="tanh")
   def run_softmax(x):            return torch.softmax(x, dim=-1)
   def run_matmul(a, b):          return a @ b
   def run_flash(q, k, v):        return torch.nn.functional.scaled_dot_product_attention(q, k, v)

We define the Triton implementations for vector addition, fused GELU, row softmax, tiled matrix multiplication, and flash attention. We express each operation using tile-level loads, computations, reductions, dot products, and stores, so that the GPU can handle blocks of data efficiently. We also provide pure PyTorch fallback functions so the tutorial still runs when Triton or a supported GPU backend is unavailable.

Vector Add and GELU

def bench(fn, *a, iters=50, warmup=10):
   if HAS_CUDA:
       for _ in range(warmup): fn(*a)
       torch.cuda.synchronize()
       t0 = time.perf_counter()
       for _ in range(iters): fn(*a)
       torch.cuda.synchronize()
       return (time.perf_counter() - t0) / iters * 1e3
   else:
       t0 = time.perf_counter()
       for _ in range(iters): fn(*a)
       return (time.perf_counter() - t0) / iters * 1e3
def check(name, got, ref, atol=1e-2, rtol=1e-2):
   got = got.float().cpu(); ref = ref.float().cpu()
   ok = torch.allclose(got, ref, atol=atol, rtol=rtol)
   md = (got - ref).abs().max().item()
   print(f"  correctness [{name:12s}] : {'PASS' if ok else 'FAIL'}  (max abs diff {md:.2e})")
   return ok
dtype = torch.float16 if HAS_CUDA else torch.float32
rule("KERNEL 1 — VECTOR ADD (load tile -> add -> store tile)")
print("cuTile version of this kernel:\n" + CUTILE_SOURCE["vector_add"])
n = 1 << 20
a = torch.randn(n, device=DEV, dtype=torch.float32)
b = torch.randn(n, device=DEV, dtype=torch.float32)
check("vector_add", run_vadd(a, b), a + b)
if BACKEND != "torch":
   print(f"  {BACKEND} time: {bench(run_vadd, a, b):.4f} ms   torch: {bench(lambda x,y:x+y, a, b):.4f} ms")
rule("KERNEL 2 — FUSED GELU(x*w + b)  (three ops fused into one memory pass)")
x = torch.randn(n, device=DEV, dtype=torch.float32)
w = torch.randn(n, device=DEV, dtype=torch.float32)
bb = torch.randn(n, device=DEV, dtype=torch.float32)
ref = torch.nn.functional.gelu(x * w + bb, approximate="tanh")
check("fused_gelu", run_fused_gelu(x, w, bb), ref)
if BACKEND != "torch":
   torch_fn = lambda x,w,b: torch.nn.functional.gelu(x*w+b, approximate="tanh")
   print(f"  {BACKEND} (1 pass): {bench(run_fused_gelu, x, w, bb):.4f} ms   "
         f"torch (3 passes): {bench(torch_fn, x, w, bb):.4f} ms")

We build benchmarking and correctness-checking utilities that compare each custom kernel against a PyTorch reference implementation. We then run the vector-addition kernel and verify that the tile-based output matches the standard PyTorch addition. We also test the fused GELU kernel, demonstrating how multiplication, bias addition, and GELU activation are combined into a single efficient pass.

Softmax and Tiled Matmul

rule("KERNEL 3 — ROW SOFTMAX (tile reductions: max then sum, numerically stable)")
rows, cols = 4096, 1024
x = torch.randn(rows, cols, device=DEV, dtype=torch.float32)
check("softmax", run_softmax(x), torch.softmax(x, dim=-1))
if BACKEND != "torch":
   print(f"  {BACKEND} time: {bench(run_softmax, x):.4f} ms   "
         f"torch: {bench(lambda z: torch.softmax(z,-1), x):.4f} ms")
rule("KERNEL 4 — TILED MATMUL (K-loop accumulation -> tensor cores)")
print("cuTile version of this kernel:\n" + CUTILE_SOURCE["matmul"])
M = N = K = 1024
a = torch.randn(M, K, device=DEV, dtype=dtype)
b = torch.randn(K, N, device=DEV, dtype=dtype)
check("matmul", run_matmul(a, b), a @ b, atol=1e-1, rtol=1e-1)
if BACKEND != "torch":
   t = bench(run_matmul, a, b)
   flops = 2 * M * N * K
   print(f"  {BACKEND}: {t:.4f} ms  ({flops/ (t*1e-3) / 1e12:.2f} TFLOP/s)   "
         f"torch: {bench(lambda x,y:x@y, a, b):.4f} ms")

We run the row-wise softmax kernel and compare it against PyTorch’s softmax to verify numerical correctness. We then perform tiled matrix multiplication, multiplying matrix blocks and accumulating along the K dimension. We benchmark these kernels against PyTorch to observe how tile-based execution performs on the active backend.

Flash Attention Kernel

rule("KERNEL 5 — FLASH ATTENTION (online softmax; the advanced capstone)")
Z, L, D = 8, 512, 64
q = torch.randn(Z, L, D, device=DEV, dtype=dtype)
k = torch.randn(Z, L, D, device=DEV, dtype=dtype)
v = torch.randn(Z, L, D, device=DEV, dtype=dtype)
ref = torch.nn.functional.scaled_dot_product_attention(q, k, v)
check("flash_attn", run_flash(q, k, v), ref, atol=2e-2, rtol=2e-2)
if BACKEND != "torch":
   sdpa = lambda q,k,v: torch.nn.functional.scaled_dot_product_attention(q,k,v)
   print(f"  {BACKEND}: {bench(run_flash, q, k, v):.4f} ms   "
         f"torch sdpa: {bench(sdpa, q, k, v):.4f} ms")
rule("DONE")
print(f"""
Summary
-------
Backend that ran : {BACKEND}
What you learned : the tile programming model (whole-tile load/compute/store),
                  fusion, tile reductions, K-loop matmul on tensor cores, and
                  an online-softmax flash-attention kernel.
To run the REAL cuTile kernels shown above you need CUDA 13.1+ and an
Ampere/Ada/Blackwell GPU. On such a machine:
   pip install 'cuda-tile[tileiras]' cupy-cuda13x
   pip install tilegym[tileiras]
Then the strings in CUTILE_SOURCE run as-is via ct.launch(...).
""")

We finish with the flash attention kernel, which applies online softmax to compute attention without materializing the full attention matrix. We compare its output to PyTorch’s scaled dot-product attention and benchmark runtime performance when a GPU backend is available. We close the tutorial by summarizing the backend we used and the main tile programming concepts we learned.

Conclusion

In conclusion, we understood how tile-based kernels map high-level mathematical operations onto efficient GPU execution patterns. We saw how fusion reduces memory traffic, how tile reductions stabilize and make softmax efficient, how tiled matrix multiplication accumulates over K-blocks, and how flash attention uses online softmax to avoid materializing the full attention matrix. We also gained a path for experimentation: we ran Triton kernels on common Colab GPUs today while still seeing how the same concepts translate to real cuTile kernels on newer CUDA 13.1+ Ampere, Ada, or Blackwell systems.