code/day84-cutlass/cute_layouts.cuThis is the source used by the lesson and its recorded evidence. Compile commands and expected output live in the directory README.
// SPDX-License-Identifier: MIT
//
// Day 84 part 1: CuTe's layout algebra, checked against the index arithmetic
// day 16's tiled matmul wrote out by hand.
//
// A cute::Layout is a (Shape, Stride) pair and nothing else. It is a function
// from a coordinate to an index, so composing two layouts composes two index
// maps, and dividing one layout by another splits it into "the tile" and
// "which tile". Every tile address day 16, day 43 and day 44 computed with
// multiplications and offsets is one of those functions, written by hand.
//
// Four gates, every one a real branch returning EXIT_FAILURE:
//
// 1. composition reproduces the worked example in NVIDIA's own CuTe layout
// algebra document, index for index. That document is the external
// answer key: the twelve values are theirs, not ours.
// 2. zipped_divide reproduces day 16's tile addressing at all 32,768
// elements of a 256 x 128 row-major matrix.
// 3. layout<0>(zipped_divide(A, B)) equals composition(A, B), the identity
// the same document states as always true.
// 4. a kernel gathers tiles through the same layout object on the device,
// and its output must equal the permutation the host built with day 16's
// arithmetic. Host and device evaluate the same Layout type.
//
// Gate 2 is the one worth reading. It is not a tautology: the left side is
// four lines of CuTe and the right side is the multiply-and-add a reader has
// been writing since day 16, and they were written from different places.
//
// Needs CUTLASS's headers, which are header only and build nothing:
// git clone --depth 1 --branch v4.7.1 https://github.com/NVIDIA/cutlass.git
//
// Build: nvcc -std=c++17 -O3 -arch=sm_75 --expt-relaxed-constexpr \
// -I cutlass/include -o cute_layouts cute_layouts.cu
// Run: ./cute_layouts
//
// Runs on any GPU this course supports; the layout algebra itself needs no
// GPU at all, and only gate 4 launches anything.
//
// VERIFIED 2026-09-02: built with CUDA 12.6 and all four gates passed on a
// Tesla T4. Transcript: evidence/run-2026-09-02.txt.
#include <cstdio>
#include <cstdlib>
#include <vector>
#include <cuda_runtime.h>
#include <cute/tensor.hpp>
using namespace cute;
// The one error macro. This file is standalone, the way a Compiler Explorer
// embed is, so it carries its own verbatim copy. `err_` carries a trailing
// underscore so it cannot collide with a variable at the call site, and the
// do/while makes the macro one statement so it survives a braceless `if`.
#define CUDA_CHECK(call) \
do { \
cudaError_t err_ = (call); \
if (err_ != cudaSuccess) { \
std::fprintf(stderr, "CUDA error %s:%d: %s: %s\n", __FILE__, \
__LINE__, #call, cudaGetErrorString(err_)); \
std::exit(EXIT_FAILURE); \
} \
} while (0)
// Day 16's matrix and day 16's tile. Static integers, so every shape and
// stride below is a compile-time value and the Layout types are empty: a
// kernel can take one by value and it costs no bytes of parameter space.
constexpr int kM = 256;
constexpr int kK = 128;
constexpr int kTileDim = 16;
constexpr int kElemsPerTile = kTileDim * kTileDim;
constexpr int kTilesM = kM / kTileDim;
constexpr int kTilesK = kK / kTileDim;
constexpr size_t kElems = static_cast<size_t>(kM) * kK;
constexpr int kThreadsPerBlock = 256;
static_assert(kM % kTileDim == 0 && kK % kTileDim == 0,
"the tiler must divide the matrix; a ragged edge is a "
"different lesson (CUTLASS calls it the epilogue path)");
static_assert(kThreadsPerBlock % 32 == 0,
"block size must be a whole number of warps");
// The twelve indices NVIDIA's CuTe layout algebra document prints for
// R = A o B with A = (6,2):(8,2) and B = (4,3):(3,1), copied from that
// document and used here as an answer key we did not write.
// Source: media/docs/cpp/cute/02_layout_algebra.md at CUTLASS v4.7.1,
// section "Computing Composition". The README carries the full URL.
constexpr int kDocComposition[12] = {0, 24, 2, 26, 8, 32,
10, 34, 16, 40, 18, 42};
// Gathers the matrix into tile-major order using the layout it is handed.
//
// One thread owns one element. Thread i writes out[i] and reads in[zd(i)],
// so a warp's 32 reads are 32 consecutive elements of one tile row pair,
// which is contiguous in `in` for the first 16 and jumps one matrix row for
// the second 16: the tile load day 16 wrote, expressed as one call.
//
// Launch assumption: gridDim.x * blockDim.x >= n. The layout is passed by
// value because its shape and stride are static, so the object is empty.
template <class ZippedLayout>
__global__ void gatherTiles(const float* __restrict__ in,
float* __restrict__ out, ZippedLayout zd,
size_t n) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
out[i] = in[zd(i)];
}
}
// Day 16's arithmetic, written the way day 16 wrote it, for one element of
// one tile. Nothing in this function knows that CuTe exists.
//
// snippet: hand-index
static int day16Index(int tileRow, int tileCol, int rowInTile, int colInTile) {
const int row = tileRow * kTileDim + rowInTile;
const int col = tileCol * kTileDim + colInTile;
return row * kK + col;
}
// end snippet
// Decomposes a tile-major linear position into the four coordinates
// day16Index takes. Mode 0 of a zipped divide is the tile and mode 1 is which
// tile, and inside each mode the leftmost coordinate moves fastest.
static void decomposeTileMajor(size_t pos, int* tileRow, int* tileCol,
int* rowInTile, int* colInTile) {
const int within = static_cast<int>(pos % kElemsPerTile);
const int which = static_cast<int>(pos / kElemsPerTile);
*rowInTile = within % kTileDim;
*colInTile = within / kTileDim;
*tileRow = which % kTilesM;
*tileCol = which / kTilesM;
}
int main() {
cudaDeviceProp prop;
CUDA_CHECK(cudaGetDeviceProperties(&prop, 0));
std::printf("GPU: %s (compute capability %d.%d), %d SMs\n", prop.name,
prop.major, prop.minor, prop.multiProcessorCount);
std::printf("CuTe from CUTLASS headers on the include path\n\n");
int status = EXIT_SUCCESS;
// Gate 1: composition against the document's own worked example.
//
// snippet: composition
auto docA = make_layout(make_shape(Int<6>{}, Int<2>{}),
make_stride(Int<8>{}, Int<2>{}));
auto docB = make_layout(make_shape(Int<4>{}, Int<3>{}),
make_stride(Int<3>{}, Int<1>{}));
auto docR = composition(docA, docB);
// end snippet
std::printf("Gate 1: composition, R = A o B\n");
std::printf(" A = ");
cute::print(docA);
std::printf("\n B = ");
cute::print(docB);
std::printf("\n R = ");
cute::print(docR);
std::printf("\n");
std::printf(" %4s %8s %8s %8s\n", "i", "B(i)", "A(B(i))", "R(i)");
for (int i = 0; i < 12; ++i) {
const int viaR = static_cast<int>(docR(i));
const int viaAB = static_cast<int>(docA(docB(i)));
std::printf(" %4d %8d %8d %8d\n", i, static_cast<int>(docB(i)), viaAB,
viaR);
if (viaR != viaAB || viaR != kDocComposition[i]) {
std::fprintf(stderr,
"composition disagrees at i = %d: R(i) = %d, "
"A(B(i)) = %d, the CuTe document says %d\n",
i, viaR, viaAB, kDocComposition[i]);
status = EXIT_FAILURE;
break;
}
}
if (status == EXIT_SUCCESS) {
std::printf(" 12 of 12 match the published table\n\n");
}
// Gates 2 and 3: the day 16 tiling, as algebra.
//
// mA is the row-major matrix: coordinate (row, col) maps to row * K + col,
// which is the stride pair (K, 1). The tiler is the tile's shape. The
// zipped divide splits mA into ((rowInTile, colInTile), (tileRow,
// tileCol)), so the first coordinate addresses inside a tile and the
// second picks the tile.
//
// snippet: tiling
auto mA = make_layout(make_shape(Int<kM>{}, Int<kK>{}),
make_stride(Int<kK>{}, _1{}));
auto tiler = make_shape(Int<kTileDim>{}, Int<kTileDim>{});
auto zd = zipped_divide(mA, tiler);
// end snippet
std::printf("Gate 2: zipped_divide against day 16's index arithmetic\n");
std::printf(" mA = ");
cute::print(mA);
std::printf("\n zd = ");
cute::print(zd);
std::printf("\n");
if (status == EXIT_SUCCESS) {
size_t checked = 0;
for (int tileCol = 0; tileCol < kTilesK && status == EXIT_SUCCESS;
++tileCol) {
for (int tileRow = 0; tileRow < kTilesM; ++tileRow) {
for (int colInTile = 0; colInTile < kTileDim; ++colInTile) {
for (int rowInTile = 0; rowInTile < kTileDim; ++rowInTile) {
const int fromCute = static_cast<int>(
zd(make_coord(rowInTile, colInTile),
make_coord(tileRow, tileCol)));
const int byHand =
day16Index(tileRow, tileCol, rowInTile, colInTile);
++checked;
if (fromCute != byHand) {
std::fprintf(
stderr,
"tile (%d,%d) element (%d,%d): CuTe says %d, "
"day 16's arithmetic says %d\n",
tileRow, tileCol, rowInTile, colInTile,
fromCute, byHand);
status = EXIT_FAILURE;
break;
}
}
if (status == EXIT_FAILURE) {
break;
}
}
if (status == EXIT_FAILURE) {
break;
}
}
}
if (status == EXIT_SUCCESS) {
std::printf(" %zu of %zu indices agree\n\n", checked, kElems);
}
}
// Gate 3. The document states this as an invariant, so a build where it
// does not hold is a build where the tile you compute on is not the tile
// the divide handed you.
if (status == EXIT_SUCCESS) {
auto tileOnly = layout<0>(zd);
auto composed = composition(mA, tiler);
std::printf(
"Gate 3: layout<0>(zipped_divide(mA, tiler)) "
"== composition(mA, tiler)\n");
std::printf(" layout<0>(zd) = ");
cute::print(tileOnly);
std::printf("\n mA o tiler = ");
cute::print(composed);
std::printf("\n");
for (int e = 0; e < kElemsPerTile; ++e) {
if (static_cast<int>(tileOnly(e)) !=
static_cast<int>(composed(e))) {
std::fprintf(stderr,
"the tile layout and the composition disagree "
"at element %d: %d against %d\n",
e, static_cast<int>(tileOnly(e)),
static_cast<int>(composed(e)));
status = EXIT_FAILURE;
break;
}
}
if (status == EXIT_SUCCESS) {
std::printf(" %d of %d agree\n\n", kElemsPerTile, kElemsPerTile);
}
}
// Gate 4: the same Layout object, evaluated on the device.
std::vector<float> h_in(kElems);
std::vector<float> h_want(kElems);
std::vector<float> h_got(kElems);
for (size_t e = 0; e < kElems; ++e) {
h_in[e] = static_cast<float>(e);
}
for (size_t pos = 0; pos < kElems; ++pos) {
int tileRow = 0;
int tileCol = 0;
int rowInTile = 0;
int colInTile = 0;
decomposeTileMajor(pos, &tileRow, &tileCol, &rowInTile, &colInTile);
h_want[pos] = h_in[static_cast<size_t>(
day16Index(tileRow, tileCol, rowInTile, colInTile))];
}
float* d_in = nullptr;
float* d_out = nullptr;
const size_t bytes = kElems * sizeof(float);
CUDA_CHECK(cudaMalloc(&d_in, bytes));
CUDA_CHECK(cudaMalloc(&d_out, bytes));
CUDA_CHECK(cudaMemcpy(d_in, h_in.data(), bytes, cudaMemcpyHostToDevice));
if (status == EXIT_SUCCESS) {
const int blocks = static_cast<int>((kElems + kThreadsPerBlock - 1) /
kThreadsPerBlock);
gatherTiles<<<blocks, kThreadsPerBlock>>>(d_in, d_out, zd, kElems);
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaDeviceSynchronize());
CUDA_CHECK(
cudaMemcpy(h_got.data(), d_out, bytes, cudaMemcpyDeviceToHost));
std::printf("Gate 4: the device gather through the same layout\n");
for (size_t pos = 0; pos < kElems; ++pos) {
if (h_got[pos] != h_want[pos]) {
std::fprintf(stderr,
"tile-major position %zu: device gathered %.1f, "
"the host permutation wants %.1f\n",
pos, static_cast<double>(h_got[pos]),
static_cast<double>(h_want[pos]));
status = EXIT_FAILURE;
break;
}
}
if (status == EXIT_SUCCESS) {
std::printf(" %zu of %zu elements match, exactly\n\n", kElems,
kElems);
}
}
CUDA_CHECK(cudaFree(d_in));
CUDA_CHECK(cudaFree(d_out));
if (status == EXIT_SUCCESS) {
std::printf("all four gates pass\n");
}
return status;
}