code/day56-graphs/graphs.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 56: CUDA graphs, priced against the launches they replace.
//
// Day 48 wrote a four-stage elementwise chain twice and found that at 16 Mi
// elements the chain's cost is traffic, not launches. This program takes the
// same four kernels, captures the four-launch chain as one CUDA graph with
// stream capture, and measures exactly what the graph gives back: the cost
// of submitting work, per iteration, at five sizes.
//
// Part 1 what the capture recorded: node count, one instantiation.
// Part 2 N iterations of the chain, four launches each, against N
// launches of the captured graph, from 1 Ki to 16 Mi elements.
//
// A graph moves no bytes and deletes no round trips through global memory.
// Both paths below run the same four kernels in the same order over the
// same buffers, so the difference between the two columns is submission
// cost and nothing else. The prediction: that difference is microseconds
// per iteration at the small end and noise at the large end.
//
// Build: nvcc -std=c++17 -O3 -arch=sm_75 -lineinfo -o graphs graphs.cu
// Run: ./graphs
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <vector>
#include <cuda_runtime.h>
#include <nvtx3/nvToolsExt.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)
// The timed sizes and the correctness size are day 48's, so the two pages
// describe the same chain. 611 is 13 x 47: no power-of-two launch shape
// covers the correctness size exactly, so every bounds check fires on every
// launch. The largest timed size is a power of two so no guard fires inside
// a timed kernel.
constexpr size_t kElems = 16777216;
constexpr size_t kCheckElems = 1048576 + 611;
constexpr int kThreadsPerBlock = 256; // 8 warps
constexpr int kWarmupRuns = 3;
constexpr int kTimedRuns = 10;
// Each timed run is a batch of whole-chain iterations, so the event pair
// brackets something far larger than its own half-microsecond resolution
// and a per-iteration difference of a few microseconds is resolvable.
constexpr int kItersPerBatch = 20;
// Part 2. Every value must fit the buffers allocated for kElems.
constexpr size_t kSweep[] = {1024, 16384, 262144, 4194304, 16777216};
constexpr int kSweepCount =
static_cast<int>(sizeof(kSweep) / sizeof(kSweep[0]));
// The chain the graph replays has four kernels, so a capture that records
// any other number of nodes did not record the chain.
constexpr size_t kChainKernels = 4;
// From code/day48-fusion/fusion.cu (day 48): the chain's parameters,
// unchanged, so this page's chain is byte for byte the one day 48 priced.
constexpr float kGamma = 1.5f;
constexpr float kBeta = -0.25f;
constexpr float kLo = -6.0f;
constexpr float kHi = 3.0f;
constexpr float kGeluA = 0.7978845608f;
constexpr float kGeluB = 0.044715f;
constexpr double kRelTolerance = 1e-5;
constexpr double kAbsTolerance = 1e-6;
static_assert(kThreadsPerBlock % 32 == 0,
"the block is a whole number of warps");
static_assert(kThreadsPerBlock <= 1024,
"1024 threads is the per-block maximum on every compute "
"capability this course targets");
static_assert(kElems % 32 == 0,
"the timed size must be a whole number of blocks at every "
"block size this file invites, down to one warp, so that no "
"bounds check runs inside a timed kernel");
static_assert(kCheckElems % 2 == 1,
"the correctness size must be odd, so that no power-of-two "
"launch shape covers it exactly and every bounds check fires "
"on every launch");
static_assert(kSweep[kSweepCount - 1] <= kElems,
"the sweep reuses the buffers allocated for kElems, so its "
"largest size cannot exceed them");
static_assert(kHi > kLo, "the clamp range must be a range");
// From code/day48-fusion/fusion.cu (day 48): the four stages as ordinary
// device functions, unchanged.
__device__ __forceinline__ float scaleBias(float v) {
return v * kGamma + kBeta;
}
__device__ __forceinline__ float gelu(float v) {
const float inner = kGeluA * (v + kGeluB * v * v * v);
return 0.5f * v * (1.0f + tanhf(inner));
}
__device__ __forceinline__ float addResidual(float v, float r) {
return v + r;
}
__device__ __forceinline__ float clampRange(float v) {
return fminf(fmaxf(v, kLo), kHi);
}
// From code/day48-fusion/fusion.cu (day 48): the four staged kernels,
// unchanged. One thread owns one element in each. 32 consecutive lanes read
// and write 32 consecutive floats, so every access coalesces into four
// 32-byte sectors per warp. Launch assumption for all four:
// gridDim.x * blockDim.x >= n.
__global__ void applyScaleBias(const float* __restrict__ in,
float* __restrict__ out, size_t n) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
out[i] = scaleBias(in[i]);
}
}
__global__ void applyGelu(const float* __restrict__ in, float* __restrict__ out,
size_t n) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
out[i] = gelu(in[i]);
}
}
__global__ void applyAddResidual(const float* __restrict__ in,
const float* __restrict__ r,
float* __restrict__ out, size_t n) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
out[i] = addResidual(in[i], r[i]);
}
}
__global__ void applyClampRange(const float* __restrict__ in,
float* __restrict__ out, size_t n) {
const size_t i = blockIdx.x * static_cast<size_t>(blockDim.x) + threadIdx.x;
if (i < n) {
out[i] = clampRange(in[i]);
}
}
// From code/day48-fusion/fusion.cu (day 48): the chain in double on the
// host, unchanged. Written for obvious correctness, not speed. It never
// allocates; the caller owns every buffer.
static void chainCpu(const float* x, const float* r, double* out, size_t n) {
for (size_t i = 0; i < n; ++i) {
double v = static_cast<double>(x[i]) * kGamma + kBeta;
const double inner = kGeluA * (v + kGeluB * v * v * v);
v = 0.5 * v * (1.0 + std::tanh(inner));
v = v + static_cast<double>(r[i]);
if (v < kLo) {
v = kLo;
}
if (v > kHi) {
v = kHi;
}
out[i] = v;
}
}
// Returns the first index where got and want differ by more than the
// tolerance, or n if they agree everywhere.
static size_t firstMismatch(const float* got, const double* want, size_t n) {
for (size_t i = 0; i < n; ++i) {
const double scale = std::fabs(want[i]);
const double allowed = kAbsTolerance + kRelTolerance * scale;
if (std::fabs(static_cast<double>(got[i]) - want[i]) > allowed) {
return i;
}
}
return n;
}
static int blocksFor(size_t n) {
return static_cast<int>((n + kThreadsPerBlock - 1) / kThreadsPerBlock);
}
// The four launches from day 48, in order, on one stream. This one function
// is both paths: part 2 times it as it stands, and captureChain records it.
// The graph and the launch loop cannot drift, because they are the same
// code.
//
// The third chevron argument is the dynamic shared-memory size, zero here;
// it has to be spelled out to reach the fourth, the stream. Day 51
// introduced that slot.
static void launchChain(cudaStream_t stream, const float* d_x, const float* d_r,
float* d_t1, float* d_t2, float* d_t3, float* d_y,
size_t n) {
const int blocks = blocksFor(n);
applyScaleBias<<<blocks, kThreadsPerBlock, 0, stream>>>(d_x, d_t1, n);
applyGelu<<<blocks, kThreadsPerBlock, 0, stream>>>(d_t1, d_t2, n);
applyAddResidual<<<blocks, kThreadsPerBlock, 0, stream>>>(d_t2, d_r, d_t3,
n);
applyClampRange<<<blocks, kThreadsPerBlock, 0, stream>>>(d_t3, d_y, n);
}
// Records the chain instead of running it, then bakes the recording into an
// executable graph. While a stream is capturing, "all operations pushed
// into the stream will not be executed, but will instead be captured into a
// graph" (cudaStreamBeginCapture, CUDA 12.6 runtime API). The capture
// freezes the grid dimensions along with everything else, so a different n
// needs a different graph; editing one in place is day 57.
// snippet: capture
static cudaGraphExec_t captureChain(cudaStream_t stream, const float* d_x,
const float* d_r, float* d_t1, float* d_t2,
float* d_t3, float* d_y, size_t n,
size_t* nodeCount) {
cudaGraph_t graph = nullptr;
CUDA_CHECK(cudaStreamBeginCapture(stream, cudaStreamCaptureModeGlobal));
launchChain(stream, d_x, d_r, d_t1, d_t2, d_t3, d_y, n);
CUDA_CHECK(cudaStreamEndCapture(stream, &graph));
CUDA_CHECK(cudaGraphGetNodes(graph, nullptr, nodeCount));
cudaGraphExec_t exec = nullptr;
CUDA_CHECK(cudaGraphInstantiate(&exec, graph, 0));
// The executable graph is self-contained, so the recording can go.
CUDA_CHECK(cudaGraphDestroy(graph));
return exec;
}
// end snippet: capture
// The course's timeKernel, with one change and a comment on why: the events
// are recorded on the stream that owns the work. Everything this program
// launches sits on its own stream, and an event pair recorded on the
// default stream would bracket work that is not there.
template <typename LaunchFn>
static float timeStreamMs(cudaStream_t stream, LaunchFn launch) {
cudaEvent_t start, stop;
CUDA_CHECK(cudaEventCreate(&start));
CUDA_CHECK(cudaEventCreate(&stop));
// Warm up the exact thing being timed, 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, and the
// first launch of a graph pays its upload.
for (int i = 0; i < kWarmupRuns; ++i) {
launch();
}
CUDA_CHECK(cudaStreamSynchronize(stream));
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaEventRecord(start, stream));
for (int i = 0; i < kTimedRuns; ++i) {
launch();
}
CUDA_CHECK(cudaEventRecord(stop, stream));
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;
}
int main() {
const int device = 0;
CUDA_CHECK(cudaSetDevice(device));
cudaDeviceProp prop;
CUDA_CHECK(cudaGetDeviceProperties(&prop, device));
std::printf("GPU: %s (compute capability %d.%d, %d SMs)\n", prop.name,
prop.major, prop.minor, prop.multiProcessorCount);
std::printf("Chain: 4 elementwise kernels from day 48, %d threads/block\n",
kThreadsPerBlock);
std::printf(
"Each timed run is %d whole-chain iterations; mean of %d runs "
"after %d warm-ups\n\n",
kItersPerBatch, kTimedRuns, kWarmupRuns);
const size_t bytes = kElems * sizeof(float);
std::vector<float> h_x(kElems);
std::vector<float> h_r(kElems);
std::vector<float> h_launch(kCheckElems);
std::vector<float> h_graph(kCheckElems);
std::vector<double> h_want(kCheckElems);
// Day 48's input: x sweeps [-4, 4) so the clamp at kHi = 3 catches a
// real share of the elements rather than none.
for (size_t i = 0; i < kElems; ++i) {
h_x[i] = (static_cast<float>(i % 2048) / 1024.0f - 1.0f) * 4.0f;
h_r[i] = static_cast<float>(i % 97) / 97.0f - 0.5f;
}
chainCpu(h_x.data(), h_r.data(), h_want.data(), kCheckElems);
float* d_x = nullptr;
float* d_r = nullptr;
float* d_t1 = nullptr;
float* d_t2 = nullptr;
float* d_t3 = nullptr;
float* d_y = nullptr;
CUDA_CHECK(cudaMalloc(&d_x, bytes));
CUDA_CHECK(cudaMalloc(&d_r, bytes));
CUDA_CHECK(cudaMalloc(&d_t1, bytes));
CUDA_CHECK(cudaMalloc(&d_t2, bytes));
CUDA_CHECK(cudaMalloc(&d_t3, bytes));
CUDA_CHECK(cudaMalloc(&d_y, bytes));
CUDA_CHECK(cudaMemcpy(d_x, h_x.data(), bytes, cudaMemcpyHostToDevice));
CUDA_CHECK(cudaMemcpy(d_r, h_r.data(), bytes, cudaMemcpyHostToDevice));
cudaStream_t stream = nullptr;
CUDA_CHECK(cudaStreamCreate(&stream));
int status = EXIT_SUCCESS;
// Correctness first, at a size no launch shape covers exactly. The
// launched chain and the graph-replayed chain run the same kernels in
// the same order over the same data, so beyond matching the CPU
// reference they must match each other bit for bit.
nvtxRangePushA("check");
launchChain(stream, d_x, d_r, d_t1, d_t2, d_t3, d_y, kCheckElems);
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaStreamSynchronize(stream));
CUDA_CHECK(cudaMemcpy(h_launch.data(), d_y, kCheckElems * sizeof(float),
cudaMemcpyDeviceToHost));
size_t checkNodes = 0;
cudaGraphExec_t checkExec = captureChain(stream, d_x, d_r, d_t1, d_t2, d_t3,
d_y, kCheckElems, &checkNodes);
CUDA_CHECK(cudaGraphLaunch(checkExec, stream));
CUDA_CHECK(cudaStreamSynchronize(stream));
CUDA_CHECK(cudaMemcpy(h_graph.data(), d_y, kCheckElems * sizeof(float),
cudaMemcpyDeviceToHost));
CUDA_CHECK(cudaGraphExecDestroy(checkExec));
nvtxRangePop();
// Real branches returning EXIT_FAILURE, not assert(). CI builds
// Release, Release defines NDEBUG, and NDEBUG deletes assert().
if (checkNodes != kChainKernels) {
std::fprintf(stderr,
"capture recorded %zu nodes, expected %zu: the graph is "
"not the chain\n",
checkNodes, kChainKernels);
status = EXIT_FAILURE;
}
const size_t badLaunch =
firstMismatch(h_launch.data(), h_want.data(), kCheckElems);
if (badLaunch != kCheckElems) {
std::fprintf(stderr,
"launched chain wrong at %zu: got %.9g, want %.9g\n",
badLaunch, static_cast<double>(h_launch[badLaunch]),
h_want[badLaunch]);
status = EXIT_FAILURE;
}
const size_t badGraph =
firstMismatch(h_graph.data(), h_want.data(), kCheckElems);
if (badGraph != kCheckElems) {
std::fprintf(stderr, "graph chain wrong at %zu: got %.9g, want %.9g\n",
badGraph, static_cast<double>(h_graph[badGraph]),
h_want[badGraph]);
status = EXIT_FAILURE;
}
if (std::memcmp(h_launch.data(), h_graph.data(),
kCheckElems * sizeof(float)) != 0) {
std::fprintf(stderr,
"the graph and the launched chain disagree: the capture "
"did not replay the same four kernels\n");
status = EXIT_FAILURE;
}
float sweepLaunchMs[kSweepCount] = {0.0f};
float sweepGraphMs[kSweepCount] = {0.0f};
if (status == EXIT_SUCCESS) {
std::printf("Part 1: what the capture recorded\n");
std::printf(
" %zu kernels launched between BeginCapture and EndCapture,\n"
" %zu nodes in the graph, instantiated once per size below.\n\n",
kChainKernels, checkNodes);
// Part 2. One graph per size, because the capture froze the grid
// dimensions. Capture and instantiation sit outside the timed
// region on purpose: they are a one-off cost, and the timeline in
// profile/ is where to read what that one-off costs.
for (int s = 0; s < kSweepCount; ++s) {
const size_t n = kSweep[s];
nvtxRangePushA("individual");
sweepLaunchMs[s] = timeStreamMs(stream, [&] {
for (int it = 0; it < kItersPerBatch; ++it) {
launchChain(stream, d_x, d_r, d_t1, d_t2, d_t3, d_y, n);
}
});
nvtxRangePop();
nvtxRangePushA("capture");
size_t nodes = 0;
cudaGraphExec_t exec = captureChain(stream, d_x, d_r, d_t1, d_t2,
d_t3, d_y, n, &nodes);
nvtxRangePop();
if (nodes != kChainKernels) {
std::fprintf(stderr,
"capture at n = %zu recorded %zu nodes, "
"expected %zu\n",
n, nodes, kChainKernels);
status = EXIT_FAILURE;
}
nvtxRangePushA("graph");
sweepGraphMs[s] = timeStreamMs(stream, [&] {
for (int it = 0; it < kItersPerBatch; ++it) {
CUDA_CHECK(cudaGraphLaunch(exec, stream));
}
});
nvtxRangePop();
CUDA_CHECK(cudaGraphExecDestroy(exec));
}
}
// Arithmetic, not a claim about the card: a zero or negative elapsed
// time means an event pair never separated, and every figure derived
// from it would be nonsense.
if (status == EXIT_SUCCESS) {
for (int s = 0; s < kSweepCount; ++s) {
if (sweepLaunchMs[s] <= 0.0f || sweepGraphMs[s] <= 0.0f) {
std::fprintf(stderr, "timing broken at n = %zu\n", kSweep[s]);
status = EXIT_FAILURE;
}
}
}
if (status == EXIT_SUCCESS) {
std::printf("Part 2: %d iterations of the chain, submitted two ways\n",
kItersPerBatch);
std::printf("%11s %16s %16s %12s %8s\n", "n", "4 launches (us)",
"1 graph (us)", "saved (us)", "ratio");
std::printf("%11s %16s %16s %12s %8s\n", "----------",
"---------------", "---------------", "-----------",
"-------");
for (int s = 0; s < kSweepCount; ++s) {
const double launchUs =
static_cast<double>(sweepLaunchMs[s]) * 1000.0 / kItersPerBatch;
const double graphUs =
static_cast<double>(sweepGraphMs[s]) * 1000.0 / kItersPerBatch;
std::printf("%11zu %16.3f %16.3f %12.3f %8.2f\n", kSweep[s],
launchUs, graphUs, launchUs - graphUs,
launchUs / graphUs);
}
std::printf(
"\nEvery column is per whole-chain iteration: four kernels, the\n"
"same bytes, the same arithmetic on both sides. The saved column\n"
"is what submitting four launches as one graph gives back, and\n"
"it cannot exceed what submission was costing in the first "
"place.\n");
}
CUDA_CHECK(cudaStreamDestroy(stream));
CUDA_CHECK(cudaFree(d_x));
CUDA_CHECK(cudaFree(d_r));
CUDA_CHECK(cudaFree(d_t1));
CUDA_CHECK(cudaFree(d_t2));
CUDA_CHECK(cudaFree(d_t3));
CUDA_CHECK(cudaFree(d_y));
return status;
}