code/day89-tool-map/fused_axpy.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 89: one fused axpy-and-clamp, written twice.
//
// This file is the CUDA C++ half. triton_axpy.py is the Triton half, with
// the same decomposition (one program per tile of 1024, a mask on the
// tail) so the decision map's "what you still own" column is something
// you read off two sources rather than take on trust.
//
// The program times nothing, on purpose. Day 89 publishes no
// tool-versus-tool performance number: no comparison this project could
// verify exists after CUDA 13.3 (research/FACT-SHEET.md section 8). What
// it does produce is a correctness gate, a clamp-coverage gate, and one
// printed line about the hardware floor Triton documents.
//
// The CUDA half passed on a Tesla T4 with CUDA 12.6 on 2026-09-02. The first
// attempt exposed an FMA-boundary input; both transcripts are retained.
//
// Build and run (needs a GPU at compute capability 7.5 or above):
// nvcc -std=c++17 -O3 -arch=sm_75 -o fused_axpy fused_axpy.cu
// ./fused_axpy
//
// PTX only, which needs no GPU at all and is the artifact the lesson
// leans on:
// nvcc -std=c++17 -O3 -arch=compute_75 -ptx -o fused_axpy.ptx \
// fused_axpy.cu
#include <cfloat>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <vector>
#include <cuda_runtime.h>
// The one error macro, byte-identical to day 5's.
#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)
// 1,000,003 is prime, so no block size divides it and the guard runs on
// every launch. kTile is 1024 to match the Triton file's BLOCK_SIZE,
// which is the whole point of the pairing.
constexpr size_t kElems = 1000003;
constexpr int kTile = 1024;
constexpr float kAlpha = 2.5f;
// Depth 1: one multiply and one add per element, no reduction, so the
// tolerance does not scale with n the way day 66's K-scaled bound does.
// The host does the multiply and the add separately; the device may
// contract them into a single FMA (day 47), which changes the result by
// at most one rounding of the product. 4 * FLT_EPSILON leaves room for
// that plus the add's own rounding. The absolute floor exists because
// the clamp produces exact zeros, and a relative tolerance says nothing
// at zero.
constexpr float kRelTolerance = 4.0f * FLT_EPSILON;
constexpr float kAbsTolerance = 1e-6f;
static_assert(kTile % 32 == 0, "tile must be a whole number of warps");
// snippet: cuda-kernel
// You own the decomposition here: which thread reads which element, the
// guard that keeps the tail in bounds, and the launch shape that puts
// blockIdx where you want it. The arithmetic is two lines in the middle.
__global__ void fusedAxpyClamp(const float* __restrict__ x,
float* __restrict__ y, float a, size_t n) {
const size_t i = static_cast<size_t>(blockIdx.x) * blockDim.x + threadIdx.x;
if (i < n) {
const float v = a * x[i] + y[i];
y[i] = v > 0.0f ? v : 0.0f;
}
}
// end snippet: cuda-kernel
// Inputs that are not exactly representable in float, so the gate below
// is testing agreement between two roundings rather than comparing two
// copies of the same integer.
static float inputX(size_t i) {
return 1.0f / static_cast<float>(i % 13 + 1);
}
static float inputY(size_t i) {
// Avoid exact cancellation: at -0.5, x=1/5 rounds to zero with separate
// host operations but to 2^-27 with a device FMA, changing gate 2 alone.
return -0.51f + static_cast<float>(i % 17) / 17.0f;
}
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);
// The map's hardware column, asked of the card in front of you rather
// than read off a blog. Triton documents NVIDIA GPUs at compute
// capability 8.0 and above
// (https://github.com/triton-lang/triton/blob/main/README.md#compatibility
// , checked 2026-09-01); cuTile documents 8.x and above
// (https://docs.nvidia.com/cuda/cutile-python/quickstart.html ,
// checked 2026-08-29). This kernel needs neither.
const int cc = prop.major * 10 + prop.minor;
std::printf("tile-DSL floor is compute capability 8.0: this device %s\n",
(cc >= 80) ? "clears it" : "does not clear it");
const size_t bytes = kElems * sizeof(float);
std::vector<float> h_x(kElems);
std::vector<float> h_y(kElems);
std::vector<float> h_want(kElems);
size_t clampedInReference = 0;
for (size_t i = 0; i < kElems; ++i) {
h_x[i] = inputX(i);
h_y[i] = inputY(i);
const float v = kAlpha * h_x[i] + h_y[i];
h_want[i] = v > 0.0f ? v : 0.0f;
if (v <= 0.0f) {
clampedInReference += 1;
}
}
float* d_x = nullptr;
float* d_y = nullptr;
CUDA_CHECK(cudaMalloc(&d_x, bytes));
CUDA_CHECK(cudaMalloc(&d_y, bytes));
CUDA_CHECK(cudaMemcpy(d_x, h_x.data(), bytes, cudaMemcpyHostToDevice));
CUDA_CHECK(cudaMemcpy(d_y, h_y.data(), bytes, cudaMemcpyHostToDevice));
const size_t blocks = (kElems + kTile - 1) / kTile;
std::printf("n=%zu, tile=%d, grid=%zu, tail lanes=%zu\n", kElems, kTile,
blocks, blocks * kTile - kElems);
fusedAxpyClamp<<<static_cast<unsigned int>(blocks), kTile>>>(
d_x, d_y, kAlpha, kElems);
CUDA_CHECK(cudaGetLastError());
CUDA_CHECK(cudaDeviceSynchronize());
std::vector<float> h_got(kElems);
CUDA_CHECK(cudaMemcpy(h_got.data(), d_y, bytes, cudaMemcpyDeviceToHost));
// Free before gating, so every failure path frees too. Status
// fall-through, one return at the end, no early exit past this line.
CUDA_CHECK(cudaFree(d_x));
CUDA_CHECK(cudaFree(d_y));
size_t mismatches = 0;
size_t firstBad = kElems;
size_t clampedOnDevice = 0;
for (size_t i = 0; i < kElems; ++i) {
if (h_got[i] == 0.0f) {
clampedOnDevice += 1;
}
const float tol = kAbsTolerance + kRelTolerance * std::fabs(h_want[i]);
if (std::fabs(h_got[i] - h_want[i]) > tol) {
mismatches += 1;
if (firstBad == kElems) {
firstBad = i;
}
}
}
int failures = 0;
// Gate 1: the two implementations of the same arithmetic agree.
if (mismatches != 0) {
std::fprintf(stderr,
"gate 1 FAIL: %zu mismatches, first at %zu: gpu %.9g "
"vs host %.9g\n",
mismatches, firstBad, static_cast<double>(h_got[firstBad]),
static_cast<double>(h_want[firstBad]));
failures += 1;
}
// Gate 2: the clamp was exercised. A relu gate on inputs that are
// never negative passes without testing the branch, which is the
// tautological test day 66 rules out. Both counts must agree and
// both must be nonzero.
if (clampedInReference == 0 || clampedOnDevice != clampedInReference) {
std::fprintf(stderr,
"gate 2 FAIL: clamped on device %zu, in reference "
"%zu (both must be equal and nonzero)\n",
clampedOnDevice, clampedInReference);
failures += 1;
}
if (failures != 0) {
return EXIT_FAILURE;
}
std::printf("all gates pass: %zu elements, 0 mismatches, %zu clamped\n",
kElems, clampedOnDevice);
return EXIT_SUCCESS;
}