code/day17-registers/registers.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 17: registers, local memory and spills.
//
// The interesting output of this program is not its stdout. It is what ptxas
// prints while compiling it, which is why the build line carries -Xptxas -v.
// That report needs no GPU: it comes out of the compiler.
//
// Four kernels, one filter, four places the same 64 values end up living.
//
// decayNarrow 4 taps. A handful of live values.
// decayWide 64 taps. The whole window live at once, in registers.
// decayWideIndexed 64 taps, second pass indexed with a runtime value. The
// compiler cannot name the element at compile time, so the
// window moves to local memory, which is off-chip.
// decayWideBounded decayWide plus __launch_bounds__, which caps the
// register budget and buys occupancy with spills.
//
// decayWideIndexed is handed shift = 0, so it must produce the same numbers
// as decayWide, element for element. Same arithmetic, same order, different
// storage. The program checks that rather than promising it.
//
// The run-time table is the other half of the same story: cudaFuncGetAttributes
// reports registers and local memory per kernel, and it must agree with what
// ptxas printed at compile time.
//
// Build: nvcc -std=c++17 -O3 -arch=sm_75 -Xptxas -v -o registers registers.cu
// Run: ./registers
//
// Verified 2026-08-30 on a Tesla T4 (compute capability 7.5), driver
// 595.84, CUDA 12.6 (V12.6.85). Transcript: evidence/run-2026-08-30.txt
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <vector>
#include <cuda_runtime.h>
// 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)
// kElems is not a multiple of the block size, so the bounds check runs on
// every launch instead of never. The input buffer carries kWideTaps - 1 extra
// elements so the widest window reads in[i + 63] without wrapping, which keeps
// a 64-bit modulo out of the inner loop. Day 11 paid for one of those and it
// moved the numbers it was there to measure.
constexpr size_t kElems = 1024ull * 1024ull + 611ull;
constexpr int kNarrowTaps = 4;
constexpr int kWideTaps = 64;
constexpr size_t kPaddedElems = kElems + kWideTaps - 1;
constexpr int kThreadsPerBlock = 256; // 8 warps
constexpr int kMinBlocksPerSm = 4;
constexpr int kWarmupRuns = 3;
constexpr int kTimedRuns = 10;
constexpr float kDecay = 0.5f;
constexpr float kRelTolerance = 1e-5f;
// Everything a compile-time constant already settles is a static_assert, not
// a run-time branch, because the compiler can answer it and the GPU should
// not have to.
static_assert(kThreadsPerBlock % 32 == 0,
"block size must be a whole number of warps");
static_assert(kElems % kThreadsPerBlock != 0,
"pick an element count the block size does not divide, so the "
"bounds check runs on every launch");
static_assert(kNarrowTaps < kWideTaps, "the ladder has to go up");
static_assert((kWideTaps & (kWideTaps - 1)) == 0,
"decayWideIndexed masks its index, so the tap count must be a "
"power of two");
static_assert(kThreadsPerBlock * kMinBlocksPerSm <= 1024,
"sm_75 holds at most 1024 resident threads per SM (compute "
"capabilities appendix, Table 30), so a larger launch bound is "
"one the hardware cannot honour");
// out[i] is a two-pass decay filter over in[i .. i + kNarrowTaps - 1]. The
// forward pass builds a running value and keeps every step; the backward pass
// folds those steps up again in the opposite order. Reading the window
// backwards is what makes every tap stay live between the two loops, which is
// what puts them in registers rather than letting the compiler stream them.
//
// One thread owns one output element. On every tap a warp's 32 lanes read 32
// consecutive floats, so each tap is four 32-byte sectors: the coalesced
// pattern from day 11.
//
// Launch assumption: gridDim.x * blockDim.x >= n, and `in` holds at least
// n + kWideTaps - 1 elements so the widest window never runs off the end.
__global__ void decayNarrow(const float* in, float* out, size_t n) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
float window[kNarrowTaps];
float carry = 0.0f;
#pragma unroll
for (int k = 0; k < kNarrowTaps; ++k) {
carry = carry * kDecay + in[i + static_cast<size_t>(k)];
window[k] = carry;
}
float acc = 0.0f;
#pragma unroll
for (int k = kNarrowTaps - 1; k >= 0; --k) {
acc = acc * kDecay + window[k];
}
out[i] = acc;
}
}
// The same filter over sixteen times the window. Every index into `window` is
// a literal after the unroll, so the array never needs an address and ptxas
// can hold it in registers. This is the kernel whose register count is worth
// comparing against the register file, because it is the one that asks for a
// lot and gets it.
//
// Memory and launch assumption: identical to decayNarrow.
// snippet: wide-kernel
__global__ void decayWide(const float* in, float* out, size_t n) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
float window[kWideTaps];
float carry = 0.0f;
#pragma unroll
for (int k = 0; k < kWideTaps; ++k) {
carry = carry * kDecay + in[i + static_cast<size_t>(k)];
window[k] = carry;
}
float acc = 0.0f;
#pragma unroll
for (int k = kWideTaps - 1; k >= 0; --k) {
acc = acc * kDecay + window[k];
}
out[i] = acc;
}
}
// end snippet: wide-kernel
// One expression of difference from decayWide, in the backward pass: the
// index is computed from `shift`, a kernel argument. The compiler cannot know
// that the caller passes zero, so it cannot name the element at compile time,
// so `window` needs a real address. An array with a real address lives in local
// memory, and local memory is in device DRAM.
//
// Memory: the window traffic is now per thread and off-chip. The reads of
// `in` are unchanged and still coalesced.
//
// Launch assumption: identical to decayNarrow. The caller must pass shift = 0
// for the result to match decayWide, and main() does.
__global__ void decayWideIndexed(const float* in, float* out, size_t n,
int shift) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
float window[kWideTaps];
float carry = 0.0f;
#pragma unroll
for (int k = 0; k < kWideTaps; ++k) {
carry = carry * kDecay + in[i + static_cast<size_t>(k)];
window[k] = carry;
}
// snippet: indexed-loop
float acc = 0.0f;
#pragma unroll
for (int k = kWideTaps - 1; k >= 0; --k) {
acc = acc * kDecay + window[(k + shift) & (kWideTaps - 1)];
}
// end snippet: indexed-loop
out[i] = acc;
}
}
// decayWide's body with a promise attached. __launch_bounds__(T, B) tells
// ptxas that no launch will use more than T threads per block and that at
// least B blocks should fit on an SM, which is a register budget: the
// register file divided by T * B. If the kernel wanted more than that, ptxas
// has to fit it anyway, and the values that no longer fit go to local memory
// as spill stores and come back as spill loads.
//
// Memory and launch assumption: identical to decayWide, plus the bound. A
// launch with more than kThreadsPerBlock threads per block now fails with
// "too many resources requested for launch" instead of running slowly.
// snippet: launch-bounds
__global__ __launch_bounds__(kThreadsPerBlock, kMinBlocksPerSm) void
decayWideBounded(const float* in, float* out, size_t n) {
// end snippet: launch-bounds
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
float window[kWideTaps];
float carry = 0.0f;
#pragma unroll
for (int k = 0; k < kWideTaps; ++k) {
carry = carry * kDecay + in[i + static_cast<size_t>(k)];
window[k] = carry;
}
float acc = 0.0f;
#pragma unroll
for (int k = kWideTaps - 1; k >= 0; --k) {
acc = acc * kDecay + window[k];
}
out[i] = acc;
}
}
// CPU reference for all four kernels; the tap count is the only thing that
// differs between them. Written for obvious correctness, not speed: plain
// loops, no OpenMP, no intrinsics. It accumulates in double where the kernels
// accumulate in float, because a reference exists to be right rather than to
// match bit for bit, and kRelTolerance covers the difference. It never
// allocates; the caller owns every buffer.
static void decayFilterCpu(const float* in, float* out, size_t n, int taps) {
double window[kWideTaps];
for (size_t i = 0; i < n; ++i) {
double carry = 0.0;
for (int k = 0; k < taps; ++k) {
carry = carry * kDecay + static_cast<double>(in[i + k]);
window[k] = carry;
}
double acc = 0.0;
for (int k = taps - 1; k >= 0; --k) {
acc = acc * kDecay + window[k];
}
out[i] = static_cast<float>(acc);
}
}
// Returns the first index where got and want differ by more than the relative
// tolerance, or n if they agree everywhere. Returning the index rather than a
// bool is the point: "wrong at 512" names the block, "wrong" does not.
static size_t firstMismatch(const float* got, const float* want, size_t n,
float relTolerance) {
for (size_t i = 0; i < n; ++i) {
const float scale = (want[i] == 0.0f) ? 1.0f : std::fabs(want[i]);
if (std::fabs(got[i] - want[i]) > relTolerance * scale) {
return i;
}
}
return n;
}
// Times a launch with CUDA events and returns the mean milliseconds per run.
// Copy it verbatim; it is the same helper day 9 handed you. The warm-up lives
// inside so it cannot go missing from the fourth kernel someone adds later.
template <typename LaunchFn>
static float timeKernel(LaunchFn launch) {
cudaEvent_t start, stop;
CUDA_CHECK(cudaEventCreate(&start));
CUDA_CHECK(cudaEventCreate(&stop));
// Warm up this kernel, not just the first kernel in the program. Lazy
// module loading has been the default since CUDA 12.2 on Linux, so the
// first launch of each kernel pays its own load.
for (int i = 0; i < kWarmupRuns; ++i) {
launch();
}
CUDA_CHECK(cudaDeviceSynchronize());
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaEventRecord(start));
for (int i = 0; i < kTimedRuns; ++i) {
launch();
}
CUDA_CHECK(cudaEventRecord(stop));
CUDA_CHECK(cudaEventSynchronize(stop));
CUDA_CHECK(cudaGetLastError());
float ms = 0.0f;
CUDA_CHECK(cudaEventElapsedTime(&ms, start, stop));
CUDA_CHECK(cudaEventDestroy(start));
CUDA_CHECK(cudaEventDestroy(stop));
return ms / kTimedRuns;
}
// One row of the table. Everything except the name and the time comes from
// the driver's own view of the compiled kernel, so the page cannot quote a
// register count the binary does not have.
static void printRow(const char* name, const cudaFuncAttributes& attr,
int blocksPerSm, int maxThreadsPerSm, float ms) {
const double occupancy = 100.0 * blocksPerSm * kThreadsPerBlock /
static_cast<double>(maxThreadsPerSm);
std::printf("%-17s %5d %9d %8d %7d %8.1f %9.3f\n", name, attr.numRegs,
static_cast<int>(attr.localSizeBytes), attr.maxThreadsPerBlock,
blocksPerSm, occupancy, ms);
}
int main() {
const int device = 0;
CUDA_CHECK(cudaSetDevice(device));
cudaDeviceProp prop;
CUDA_CHECK(cudaGetDeviceProperties(&prop, device));
// The ceiling every register count on this page is measured against, read
// from the device rather than copied from a table. A card that reports no
// resident-thread limit would make the division meaningless, so this is a
// branch and not an assumption.
if (prop.maxThreadsPerMultiProcessor <= 0) {
std::fprintf(stderr, "device reports %d resident threads per SM\n",
prop.maxThreadsPerMultiProcessor);
return EXIT_FAILURE;
}
const int regsForFullOccupancy =
prop.regsPerMultiprocessor / prop.maxThreadsPerMultiProcessor;
std::printf("GPU: %s (compute capability %d.%d)\n", prop.name, prop.major,
prop.minor);
std::printf("32-bit registers per SM: %d\n", prop.regsPerMultiprocessor);
std::printf("resident threads per SM: %d\n",
prop.maxThreadsPerMultiProcessor);
std::printf(
"so %d registers per thread is the most a kernel can use and still "
"fill this card\n\n",
regsForFullOccupancy);
std::printf("n = %zu elements, %d threads per block, %d and %d taps\n\n",
kElems, kThreadsPerBlock, kNarrowTaps, kWideTaps);
const size_t inBytes = kPaddedElems * sizeof(float);
const size_t outBytes = kElems * sizeof(float);
const int blocks =
static_cast<int>((kElems + kThreadsPerBlock - 1) / kThreadsPerBlock);
// Every host pointer carries h_ and every device pointer d_, because the
// prefix is the only thing standing between you and passing one where the
// other belongs.
std::vector<float> h_in(kPaddedElems);
std::vector<float> h_narrow(kElems);
std::vector<float> h_wide(kElems);
std::vector<float> h_indexed(kElems);
std::vector<float> h_bounded(kElems);
std::vector<float> h_want(kElems);
for (size_t i = 0; i < kPaddedElems; ++i) {
h_in[i] = static_cast<float>(i % 17) * 0.25f;
}
float* d_in = nullptr;
float* d_out = nullptr;
CUDA_CHECK(cudaMalloc(&d_in, inBytes));
CUDA_CHECK(cudaMalloc(&d_out, outBytes));
CUDA_CHECK(cudaMemcpy(d_in, h_in.data(), inBytes, cudaMemcpyHostToDevice));
// One correctness launch per kernel, each followed by both checks in this
// order: cudaGetLastError() reports a launch the driver refused, and
// cudaDeviceSynchronize() reports what the kernel did.
decayNarrow<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems);
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaDeviceSynchronize());
CUDA_CHECK(
cudaMemcpy(h_narrow.data(), d_out, outBytes, cudaMemcpyDeviceToHost));
decayWide<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems);
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaDeviceSynchronize());
CUDA_CHECK(
cudaMemcpy(h_wide.data(), d_out, outBytes, cudaMemcpyDeviceToHost));
decayWideIndexed<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems, 0);
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaDeviceSynchronize());
CUDA_CHECK(
cudaMemcpy(h_indexed.data(), d_out, outBytes, cudaMemcpyDeviceToHost));
decayWideBounded<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems);
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaDeviceSynchronize());
CUDA_CHECK(
cudaMemcpy(h_bounded.data(), d_out, outBytes, cudaMemcpyDeviceToHost));
// What the driver thinks each compiled kernel needs. These four numbers
// are the ptxas report read from the other end, and they have to agree
// with it.
cudaFuncAttributes attrNarrow;
cudaFuncAttributes attrWide;
cudaFuncAttributes attrIndexed;
cudaFuncAttributes attrBounded;
CUDA_CHECK(cudaFuncGetAttributes(&attrNarrow, decayNarrow));
CUDA_CHECK(cudaFuncGetAttributes(&attrWide, decayWide));
CUDA_CHECK(cudaFuncGetAttributes(&attrIndexed, decayWideIndexed));
CUDA_CHECK(cudaFuncGetAttributes(&attrBounded, decayWideBounded));
int blocksNarrow = 0;
int blocksWide = 0;
int blocksIndexed = 0;
int blocksBounded = 0;
CUDA_CHECK(cudaOccupancyMaxActiveBlocksPerMultiprocessor(
&blocksNarrow, decayNarrow, kThreadsPerBlock, 0));
CUDA_CHECK(cudaOccupancyMaxActiveBlocksPerMultiprocessor(
&blocksWide, decayWide, kThreadsPerBlock, 0));
CUDA_CHECK(cudaOccupancyMaxActiveBlocksPerMultiprocessor(
&blocksIndexed, decayWideIndexed, kThreadsPerBlock, 0));
CUDA_CHECK(cudaOccupancyMaxActiveBlocksPerMultiprocessor(
&blocksBounded, decayWideBounded, kThreadsPerBlock, 0));
const float msNarrow = timeKernel([&] {
decayNarrow<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems);
});
const float msWide = timeKernel(
[&] { decayWide<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems); });
const float msIndexed = timeKernel([&] {
decayWideIndexed<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems, 0);
});
const float msBounded = timeKernel([&] {
decayWideBounded<<<blocks, kThreadsPerBlock>>>(d_in, d_out, kElems);
});
// Freed before anything is reported, so the failure paths below free too.
// Every cudaMalloc has a matching cudaFree before every return, including
// the return you take when the answer is wrong.
CUDA_CHECK(cudaFree(d_in));
CUDA_CHECK(cudaFree(d_out));
decayFilterCpu(h_in.data(), h_want.data(), kElems, kNarrowTaps);
size_t bad =
firstMismatch(h_narrow.data(), h_want.data(), kElems, kRelTolerance);
if (bad != kElems) {
std::fprintf(stderr, "decayNarrow wrong at %zu: got %.9g, want %.9g\n",
bad, h_narrow[bad], h_want[bad]);
return EXIT_FAILURE;
}
decayFilterCpu(h_in.data(), h_want.data(), kElems, kWideTaps);
const char* wideNames[3] = {"decayWide", "decayWideIndexed",
"decayWideBounded"};
const float* wideResults[3] = {h_wide.data(), h_indexed.data(),
h_bounded.data()};
for (int v = 0; v < 3; ++v) {
bad =
firstMismatch(wideResults[v], h_want.data(), kElems, kRelTolerance);
if (bad != kElems) {
std::fprintf(stderr, "%s wrong at %zu: got %.9g, want %.9g\n",
wideNames[v], bad, wideResults[v][bad], h_want[bad]);
return EXIT_FAILURE;
}
}
// The claim the whole page rests on: moving the window to local memory
// changed where the values live and nothing else. A real branch rather
// than an assert, because CI builds Release, Release defines NDEBUG, and
// NDEBUG deletes assert() from exactly the build that matters.
bad = firstMismatch(h_indexed.data(), h_wide.data(), kElems, kRelTolerance);
if (bad != kElems) {
std::fprintf(stderr,
"decayWideIndexed and decayWide disagree at %zu (%.9g "
"against %.9g); they are meant to differ only in where "
"the window is stored\n",
bad, h_indexed[bad], h_wide[bad]);
return EXIT_FAILURE;
}
// Reported, not gated. The two kernels issue the same floating-point
// operations in the same order, so every element should also be equal bit
// for bit; a handful of last-bit differences would mean ptxas contracted a
// multiply and an add in one of them and not the other, which is worth
// knowing and is not a failure.
size_t exactDiffs = 0;
for (size_t i = 0; i < kElems; ++i) {
if (h_indexed[i] != h_wide[i]) {
++exactDiffs;
}
}
std::printf("%-17s %5s %9s %8s %7s %8s %9s\n", "kernel", "regs", "local B",
"maxTPB", "blk/SM", "occ %", "ms");
std::printf("%-17s %5s %9s %8s %7s %8s %9s\n", "-----------------", "-----",
"---------", "--------", "-------", "--------", "---------");
printRow("decayNarrow", attrNarrow, blocksNarrow,
prop.maxThreadsPerMultiProcessor, msNarrow);
printRow("decayWide", attrWide, blocksWide,
prop.maxThreadsPerMultiProcessor, msWide);
printRow("decayWideIndexed", attrIndexed, blocksIndexed,
prop.maxThreadsPerMultiProcessor, msIndexed);
printRow("decayWideBounded", attrBounded, blocksBounded,
prop.maxThreadsPerMultiProcessor, msBounded);
std::printf(
"\nall four kernels match the CPU reference; decayWideIndexed and "
"decayWide\ndiffer at %zu of %zu elements, so the storage change did "
"not change the answer\n",
exactDiffs, kElems);
return EXIT_SUCCESS;
}