code/day90-library-or-kernel/library_or_kernel.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 90, the module 9 checkpoint: library or kernel?
//
// The audit written as a program. It launches nothing and needs no GPU.
// Every input is a row this course already measured and committed, and
// every output is arithmetic over those rows plus a four-gate decision.
//
// Inputs, by transcript:
// code/day20-convolution/evidence/run-2026-08-30.txt capstone 1
// code/day39-cccl/evidence/run-2026-08-30.txt the calibration
// code/day40-pagerank/evidence/run-2026-08-30.txt capstone 2
// code/day44-matmul-2/evidence/run-2026-09-01.txt capstone 4's ancestor
// code/day48-fusion/evidence/run-2026-09-01.txt the toll model
// code/day60-capstone-3/evidence/run-2026-09-01.txt capstone 3
//
// Each constant below names the row it came from. Nothing is estimated:
// where a transcript prints a rate, this file recomputes it from the
// bytes and the milliseconds and refuses to continue if the two
// disagree, so a mistyped constant cannot reach the audit table.
//
// VERIFIED: Tesla T4 and CUDA_VISIBLE_DEVICES="" paths, CUDA 12.6,
// driver 580.173.02, 2026-09-02.
//
// Build: nvcc -std=c++17 -O3 -arch=sm_75 -o library_or_kernel \
// library_or_kernel.cu
// Run: ./library_or_kernel
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cuda_runtime.h>
// The one error macro, copied verbatim. This file is standalone.
#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)
// A measured row, copied out of one transcript. `recordedGbs` is what
// that transcript printed in its own GB/s column. The program recomputes
// it from `bytes` and `ms`; the two disagreeing means this file was
// typed wrong, not that the card changed.
struct Row {
const char* label;
size_t bytes;
double ms;
double recordedGbs;
};
// Bytes are counted from each program's own constants, the way its
// transcript counts them. Day 20's image buffers are 4096 x 4096 floats,
// 67,108,864 bytes each; day 48's are 16,777,216 floats; day 39's reduce
// row counts input bytes only, which is the one figure all three of its
// implementations share; day 40's rows are one iteration times fifty.
constexpr Row kRows[] = {
// day 20, capstone 1: "timed on 4096 x 4096, mean of 10 runs"
{"d20 copy (the baseline)", 134217728, 0.566, 237.0},
{"d20 blur, tiled, pitch 36", 134217728, 3.033, 44.3},
{"d20 sobel magnitude", 134217728, 2.221, 60.4},
{"d20 threshold", 83886080, 0.344, 244.0},
{"d20 sobel + threshold, two kernels", 218103808, 2.565, 85.0},
{"d20 sobel + threshold, fused", 83886080, 1.323, 63.4},
// day 39, the calibration: parts 1 and 2
{"d39 reduce, hand, day 24 version 4", 67106420, 0.349, 192.2},
{"d39 reduce, cub::DeviceReduce::Sum", 67106420, 0.254, 264.4},
{"d39 scan, hand, two-level Hillis-Steele", 134212840, 1.331, 100.8},
{"d39 scan, cub::DeviceScan::ExclusiveSum", 134212840, 0.571, 235.1},
// day 40, capstone 2: the per-stage table and the 50-iteration table
{"d40 copy, 128 MiB (DRAM)", 134217728, 0.5483, 244.8},
{"d40 50 iterations, our kernels", 171873500, 1.7963, 95.7},
{"d40 50 iterations, Thrust + CUB", 308694300, 9.6324, 32.0},
// day 48, the toll model: part 2 at 16,777,216 elements
{"d48 copyFloor (the floor)", 134217728, 0.548, 244.8},
{"d48 staged chain, timed as one", 603979776, 2.419, 249.7},
{"d48 chainFused", 201326592, 0.785, 256.5},
};
constexpr int kRowCount = sizeof(kRows) / sizeof(kRows[0]);
// Row indices used below by name, so a reordering of the table cannot
// silently repoint a fusion case at the wrong pair.
enum RowId {
kD20Copy,
kD20Blur,
kD20Sobel,
kD20Threshold,
kD20SobelThresholdTwo,
kD20SobelThresholdFused,
kD39ReduceHand,
kD39ReduceCub,
kD39ScanHand,
kD39ScanCub,
kD40Copy,
kD40Ours,
kD40Library,
kD48Copy,
kD48Staged,
kD48Fused,
};
static_assert(kRowCount == kD48Fused + 1,
"the row table and the RowId enum have drifted apart");
// The transcripts print milliseconds to four significant figures, so a
// recomputed rate can differ from the printed one by about 0.15 percent
// from rounding alone: half a unit in the last place of day 40's 0.0359
// per-iteration figure is 0.14 percent of it. Half a percent is that,
// doubled, and covers nothing else.
constexpr double kRateTolerance = 0.005;
// A stage is "already at bandwidth" when its own achieved rate sits
// within a tenth of the copy row measured in the same run. Ten percent
// rather than five because a stage that re-reads a buffer small enough
// to sit in L2 can print above the copy row: day 48's fused chain reads
// a residual array of 64 MiB against a 4 MiB L2 and still lands 4.8
// percent over the floor.
constexpr double kAtBandwidth = 0.10;
// Gate 4's threshold: how much faster your kernel has to be before
// owning it is worth the maintenance. A kernel you keep is a kernel you
// re-tune per architecture, and day 44 spent six steps to buy 6x. This
// is a judgement, not a measurement, which is why it is a named constant
// and why the audit prints how close each verdict sits to it.
constexpr double kWorthOwning = 1.2;
static double gigabytesPerSecond(size_t bytes, double ms) {
return static_cast<double>(bytes) / (ms * 1.0e6);
}
static double relativeGap(double a, double b) {
return std::fabs(a - b) / b;
}
// A pair of chains that compute the same thing, one fused and one split
// across library-shaped calls, plus the copy row from the same run.
struct FusionCase {
const char* label;
int fused;
int staged;
int copy;
bool expectBytesExplainIt;
};
constexpr FusionCase kCases[] = {
{"d48 four elementwise kernels vs one", kD48Fused, kD48Staged, kD48Copy,
true},
{"d20 sobel + threshold", kD20SobelThresholdFused, kD20SobelThresholdTwo,
kD20Copy, false},
{"d40 one PageRank iteration", kD40Ours, kD40Library, kD40Copy, false},
};
constexpr int kCaseCount = sizeof(kCases) / sizeof(kCases[0]);
// snippet: toll
// The toll a library call pays when it cannot fuse with its neighbours:
// the extra bytes the split pushes through global memory, divided by the
// rate this card copies at. It is a prediction and not a measurement,
// and printing it beside the observed difference is the whole point,
// because where it is wrong it is wrong in one direction and for one
// reason.
static double predictedTollMs(size_t extraBytes, double copyGbs) {
return static_cast<double>(extraBytes) / (copyGbs * 1.0e6);
}
// end snippet
enum Semantics { kExact, kWrapper, kNoCall };
enum Layout { kAsIs, kNeedsTransform };
enum Neighbours { kStandalone, kFusedNow };
enum Verdict { kOffPath, kKeepKernel, kCallLibrary };
struct Candidate {
const char* where;
const char* op;
const char* call;
bool onCriticalPath;
Semantics semantics;
Layout layout;
Neighbours neighbours;
bool gapMeasured;
// library milliseconds over hand milliseconds for the same work.
// Above 1 the hand kernel won. Meaningless when gapMeasured is false.
double handWinRatio;
};
// snippet: gates
// The four gates, in the order that ends the audit soonest.
//
// Gate 0 is not one of the four. It is day 50's first checklist
// question, and an op that owns none of the run time is not worth
// replacing whatever the other answers are.
//
// Gate 2, layout, decides nothing here on purpose. A layout mismatch is
// a price, not a veto: cuBLAS wanting column-major costs a swap of the
// operands, not a rewrite. It becomes a veto only when the conversion
// costs more than the op, which is a COO matrix rebuilt into CSR on
// every call rather than once.
static Verdict verdict(const Candidate& c) {
if (!c.onCriticalPath) {
return kOffPath; // gate 0, day 50 step 1
}
if (c.semantics == kNoCall) {
return kKeepKernel; // gate 1: nothing computes what you compute
}
if (!c.gapMeasured) {
return kCallLibrary; // gates 2 and 3 only price it
}
if (c.handWinRatio >= kWorthOwning) {
return kKeepKernel; // gate 4: the win pays for the maintenance
}
return kCallLibrary;
}
// end snippet
constexpr Candidate kCandidates[] = {
{"day 39", "sum 16,776,605 floats", "cub::DeviceReduce::Sum", true, kExact,
kAsIs, kStandalone, true, 0.254 / 0.349},
{"day 39", "exclusive scan of 16,776,605 uints",
"cub::DeviceScan::ExclusiveSum", true, kExact, kAsIs, kStandalone, true,
0.571 / 1.331},
{"capstone 1", "5x5 Gaussian blur", "nppiFilterGaussBorder", true, kExact,
kAsIs, kStandalone, false, 0.0},
{"capstone 1", "sobel magnitude then threshold, fused",
"NPP filter plus NPP threshold, two calls", true, kNoCall, kAsIs,
kFusedNow, false, 0.0},
{"capstone 2", "SpMV inside the iteration",
"thrust::gather plus cub::DeviceSegmentedReduce", true, kWrapper,
kNeedsTransform, kFusedNow, true, 0.1524 / 0.0273},
{"capstone 2", "the same SpMV", "cusparseSpMV on the same CSR", true,
kExact, kAsIs, kStandalone, false, 0.0},
{"capstone 3", "the per-frame kernel chain", "any of them", false, kExact,
kAsIs, kStandalone, false, 0.0},
{"day 44", "step 6 FP32 matmul at N = 2048",
"cublasGemmEx, CUDA_R_32F, CUBLAS_COMPUTE_32F", true, kExact,
kNeedsTransform, kStandalone, true, 2.858 / 4.219},
{"capstone 4", "the WMMA GEMM",
"cublasGemmEx, CUDA_R_16F in, CUDA_R_32F out", true, kExact,
kNeedsTransform, kStandalone, false, 0.0},
};
constexpr int kCandidateCount = sizeof(kCandidates) / sizeof(kCandidates[0]);
static const char* semanticsText(Semantics s) {
if (s == kExact) {
return "exact";
}
if (s == kWrapper) {
return "composed";
}
return "no call";
}
static const char* layoutText(Layout l) {
return (l == kAsIs) ? "as is" : "transform";
}
static const char* neighboursText(Neighbours n) {
return (n == kStandalone) ? "standalone" : "fused now";
}
static const char* verdictText(Verdict v) {
if (v == kOffPath) {
return "off the path";
}
return (v == kKeepKernel) ? "keep the kernel" : "call the library";
}
// Prints the card this happens to run on, when there is one. The audit
// does not use it: every number below came from the Tesla T4 named in
// the transcripts, and a different card here would not change a single
// row. CUDA_CHECK does not wrap the count call, because zero devices is
// the expected case on the machine this day was written for and exiting
// would be wrong.
static void printHost() {
int devices = 0;
const cudaError_t status = cudaGetDeviceCount(&devices);
if (status != cudaSuccess || devices == 0) {
std::printf("host: no CUDA device visible (%s)\n",
status == cudaSuccess ? "the driver reports zero devices"
: cudaGetErrorString(status));
std::printf(
"that changes nothing here: every row below was\n"
"measured on the Tesla T4 named in the transcripts\n\n");
return;
}
cudaDeviceProp prop{};
CUDA_CHECK(cudaGetDeviceProperties(&prop, 0));
std::printf("host: %s (compute capability %d.%d), unused by the audit\n\n",
prop.name, prop.major, prop.minor);
}
int main() {
int failures = 0;
printHost();
std::printf("Part 1: every row recomputed from its own bytes and ms\n");
std::printf("%-42s %12s %8s %8s\n", "row", "bytes", "GB/s", "printed");
for (int i = 0; i < kRowCount; ++i) {
const double got = gigabytesPerSecond(kRows[i].bytes, kRows[i].ms);
std::printf("%-42s %12zu %8.1f %8.1f\n", kRows[i].label, kRows[i].bytes,
got, kRows[i].recordedGbs);
if (relativeGap(got, kRows[i].recordedGbs) > kRateTolerance) {
std::fprintf(stderr,
"FAIL: %s recomputes to %.2f GB/s, transcript "
"prints %.1f\n",
kRows[i].label, got, kRows[i].recordedGbs);
++failures;
}
}
std::printf(
"\nPart 2: what an unfused split costs, predicted and "
"observed\n");
std::printf("%-36s %12s %9s %9s %7s %6s\n", "case", "extra bytes",
"predicted", "observed", "obs/pre", "at bw");
for (int c = 0; c < kCaseCount; ++c) {
const Row& fused = kRows[kCases[c].fused];
const Row& staged = kRows[kCases[c].staged];
const Row& copy = kRows[kCases[c].copy];
if (staged.bytes <= fused.bytes || staged.ms <= fused.ms) {
std::fprintf(stderr, "FAIL: %s is not a split of the fused row\n",
kCases[c].label);
++failures;
continue;
}
const size_t extra = staged.bytes - fused.bytes;
const double predicted = predictedTollMs(extra, copy.recordedGbs);
const double observed = staged.ms - fused.ms;
const double ratio = observed / predicted;
// Do the two chains already run at the card's copy rate? That is
// a separate question from the one above, answered from each
// stage's own achieved bandwidth, and the whole model rests on
// the two answers agreeing.
const bool atBandwidth =
relativeGap(gigabytesPerSecond(fused.bytes, fused.ms),
copy.recordedGbs) <= kAtBandwidth &&
relativeGap(gigabytesPerSecond(staged.bytes, staged.ms),
copy.recordedGbs) <= kAtBandwidth;
std::printf("%-36s %12zu %9.4f %9.4f %7.2f %6s\n", kCases[c].label,
extra, predicted, observed, ratio,
atBandwidth ? "yes" : "no");
if (atBandwidth != kCases[c].expectBytesExplainIt) {
std::fprintf(stderr,
"FAIL: %s, stages at copy bandwidth is %s, the "
"table below assumes %s\n",
kCases[c].label, atBandwidth ? "yes" : "no",
kCases[c].expectBytesExplainIt ? "yes" : "no");
++failures;
}
// The claim, stated so a run can kill it: bytes predict the toll
// when both chains are already at bandwidth, and under-predict it
// otherwise. Never the other way round, because a split cannot
// move fewer bytes than the byte count says.
const bool agrees = relativeGap(observed, predicted) <= 0.02;
const bool underPredicts = ratio > 1.5;
if (kCases[c].expectBytesExplainIt ? !agrees : !underPredicts) {
std::fprintf(stderr,
"FAIL: %s, observed/predicted is %.2f, which is "
"not what a chain %s bandwidth should give\n",
kCases[c].label, ratio,
kCases[c].expectBytesExplainIt ? "at" : "below");
++failures;
}
}
std::printf("\nPart 3: the audit, gate by gate\n");
std::printf(
"gate 4 unmeasured plus \"call the library\" means nothing\n"
"blocks the call, not that the call is faster\n");
std::printf("%-11s %-38s %-9s %-9s %-10s %8s %s\n", "where", "op",
"gate 1", "gate 2", "gate 3", "gate 4", "verdict");
int library = 0;
int kernel = 0;
double closestToThreshold = 1.0e9;
for (int i = 0; i < kCandidateCount; ++i) {
const Candidate& c = kCandidates[i];
const Verdict v = verdict(c);
char gap[16];
if (c.gapMeasured) {
std::snprintf(gap, sizeof(gap), "%.2fx", c.handWinRatio);
const double margin = relativeGap(c.handWinRatio, kWorthOwning);
if (margin < closestToThreshold) {
closestToThreshold = margin;
}
} else {
std::snprintf(gap, sizeof(gap), "%s", "unmeasured");
}
// A candidate that fails gate 0 never reaches the four, so its
// row prints dashes rather than answers nobody worked out.
const bool asked = (v != kOffPath);
std::printf("%-11s %-38s %-9s %-9s %-10s %8s %s\n", c.where, c.op,
asked ? semanticsText(c.semantics) : "-",
asked ? layoutText(c.layout) : "-",
asked ? neighboursText(c.neighbours) : "-",
asked ? gap : "-", verdictText(v));
std::printf("%-11s -> %s\n", "", c.call);
library += (v == kCallLibrary) ? 1 : 0;
kernel += (v == kKeepKernel) ? 1 : 0;
}
std::printf(
"\n%d call the library, %d keep the kernel, %d off the "
"path\n",
library, kernel, kCandidateCount - library - kernel);
std::printf(
"gate 4's threshold is %.2fx; the nearest measured gap "
"sits %.0f percent away\n",
kWorthOwning, closestToThreshold * 100.0);
// A procedure that returns one answer is not a procedure. This is
// the cheapest check that the gates do any work at all.
if (library == 0 || kernel == 0) {
std::fprintf(stderr,
"FAIL: the audit reached one verdict for "
"every candidate\n");
++failures;
}
// And the audit's answers must not hang on a judgement call. If a
// measured gap sits near the threshold, the verdict is an opinion
// wearing a number, and the page has to say which one.
if (closestToThreshold < 0.20) {
std::fprintf(stderr,
"FAIL: a measured gap sits within 20 percent of "
"gate 4's threshold\n");
++failures;
}
if (failures != 0) {
std::fprintf(stderr, "%d check(s) failed\n", failures);
return EXIT_FAILURE;
}
std::printf("\nevery check passed\n");
return EXIT_SUCCESS;
}