Lab 1.1 — First Kernels: SAXPY in CUDA C++ and Python
← Course 4 syllabus · Module 1 · Prev: « Lab 0.4 · Next: Lab 1.2 »
Goal
First contact with the GPU as a throughput machine, on the Linux desktop (RTX 4090) where the CUDA tooling is deepest. The lab builds the kernel-launch mental model — the grid/block/thread hierarchy the API speaks, and the warp (32 threads in lockstep) that the hardware actually executes, the same lanes-first thinking as Course 3 Module 5’s NEON work with the lane count turned up and the scheduling inverted — then makes it concrete with the “hello world” of GPU computing: SAXPY, \(y_i \leftarrow \alpha x_i + y_i\), written and verified twice. Once in CUDA C++ under nvcc and CMake, once in Python with a Numba kernel and CuPy arrays — because prototype-in-Python, ship-in-C++ is the working rhythm of this whole module. Two habits are installed here and never relaxed: every CUDA call is error-checked, and nothing is timed without a synchronization story — the lab deliberately times a kernel wrong before timing it right.
Recommended reading
- Motta — the opening material on the CUDA programming model and writing/launching first kernels (title-level reference; confirm against the copy in hand).
- CUDA C++ Programming Guide (NVIDIA, free) — the Programming Model chapter: kernels, thread hierarchy, and the host/device split. This is the primary source the whole module leans on.
- CUDA C++ Best Practices Guide (NVIDIA, free) — the sections on timing with CUDA events and on error handling; short, and exactly this lab’s two habits.
- Numba CUDA documentation — Writing CUDA Kernels and Kernel Invocation: the decorator, the launch syntax, and what Numba compiles your Python into.
- CuPy documentation — CuPy Basics and the Elementwise kernels section; note where CuPy is “NumPy on the GPU” and where it hands you raw CUDA machinery.
Prerequisites
- Lab 0.1: the repo skeleton and its CMake conventions — and the repo cloned and building on the Linux box, since all of Module 1 runs there.
- The CUDA Toolkit and NVIDIA driver installed on the Linux box, along with Nsight Systems and Nsight Compute;
nvcc --versionandnvidia-smiboth answer. - Comfort with the Course 3 Module 2 measurement discipline — predicted bandwidth math before the benchmark, not after.
Project & environment setup
Module 1 lives in cuda/ of the labs repo’s course4/ workspace, C++ and Python side by side. The C++ side is a CUDA-enabled CMake subtree:
# cuda/CMakeLists.txt
enable_language(CUDA)
set(CMAKE_CUDA_ARCHITECTURES 89) # RTX 4090: Ada Lovelace generation
set(CMAKE_CUDA_STANDARD 17)
add_executable(saxpy saxpy.cu)The Python side is a uv-managed venv in the same directory:
cd cuda && uv init && uv add numpy numba matplotlib
# CuPy: install the wheel matching the installed CUDA Toolkit's major
# version — e.g. `uv add cupy-cuda12x`; check the CuPy installation docs.
uv run python -c "from numba import cuda; cuda.detect()"Where results go:
| Artifact | Path |
|---|---|
| Notes, device-properties table, predicted-vs-measured, wrong-vs-right timing postmortem | labs/lab-1-1/notes.md |
Nsight Systems timeline (.nsys-rep), sanitizer output |
labs/lab-1-1/captures/ |
| Timing CSVs / benchmark output from both languages | labs/lab-1-1/benchmarks/ |
Background
- The launch hierarchy vs. the execution reality. You launch a grid of blocks of threads; the hardware runs warps of 32 threads in lockstep on a streaming multiprocessor. Blocks are the unit of shared-memory cooperation (Lab 1.2); warps are the unit of execution and divergence. The GPU hides memory latency not with the caches and speculation of Course 3’s CPU, but by having enough warps in flight to switch to — the philosophical hinge of this whole module.
- SAXPY is memory-bound by construction. Per element: read \(x_i\) and \(y_i\), write \(y_i\) — 12 bytes for 2 FLOPs, an arithmetic intensity of \(2/12 = 1/6\) FLOP/byte. No kernel cleverness changes that ratio; SAXPY is a bandwidth meter wearing a math costume, which is why it’s the right first kernel.
- Error checking is not optional. Every runtime call returns a
cudaError_t; kernel launches return errors asynchronously, so a launch is only known-good aftercudaGetLastError()plus a synchronizing call. The course rule from here on: every call wrapped in a check macro, and every new kernel run once undercompute-sanitizer(the successor tocuda-memcheck) before it’s believed. - Launches are asynchronous — so naive timing lies. A host timer around an unsynchronized launch measures launch overhead (order microseconds), not the kernel. Correct GPU timing brackets work with CUDA events (
cudaEvent_tin C++,cp.cuda.Eventin Python) or synchronizes explicitly. This lab commits the crime on purpose, then fixes it. - The Python stack. Numba JIT-compiles a decorated Python function to PTX for the same GPU — same grid/block launch, same warps underneath. CuPy supplies the array layer (allocation, transfer, NumPy-mirroring ops). Together they are the prototyping bench; the C++ target is what ships.
Tasks
CUDA C++
- Interrogate the device. Write a small
device_querytarget that prints thecudaDevicePropfields that matter for this module — name, compute capability, warp size, SM count, global memory size, memory clock/bus (as reported), shared memory per block, registers per block, max threads per block. Record the table innotes.md; every predicted number in Module 1 traces back to this query, not to a spec sheet. - Bring up the CMake CUDA target. Get the
saxpyexecutable configuring and building with the fragment above, from the same presets as Lab 0.1. Confirm a deliberately wrongCMAKE_CUDA_ARCHITECTURESvalue fails in a recognizable way, then set it back to89. - Install the error-check discipline. Write the
CUDA_CHECK(...)macro (file/line, error string, abort) and the post-launch check idiom. Force one real error — e.g. an absurd block size — and record what the checked vs. unchecked failure mode looks like. - SAXPY with a grid-stride loop. Implement the kernel (your code, not the write-up’s) using the grid-stride loop idiom so one launch configuration handles any \(N\); host side allocates, copies, launches, copies back — every call checked. Run once under
compute-sanitizerand save the clean output. - Correctness vs. a CPU reference. Compare against a straightforward CPU SAXPY at \(N = 2^{20}\) and \(2^{24}\): max absolute error should be explainable by float rounding alone (Course 1 §3 — same-order operations, so exact agreement is in fact common; say why in your notes).
- Time it wrong, then right. First: host timer around the launch, no sync — record the absurd number. Then: CUDA events bracketing the kernel, warm-up launch excluded, median of ≥10 runs. Convert to effective bandwidth via \(12N/t\) and keep the number for Lab 1.5’s roofline.
CUDA Python
- Numba SAXPY. Same kernel in Numba (
@cuda.jit,cuda.grid(1)or an explicit grid-stride loop), launched over CuPy or Numba device arrays. First call includes JIT compilation — measure it, then exclude it. - CuPy, two ways. SAXPY as a plain CuPy expression (
a * x + y— note the temporaries it implies) and as a single fusedcp.ElementwiseKernel. This is the array-layer lesson: what the library does per expression, and what fusing buys. - Time it wrong, then right — again. Show that a wall-clock timer around a CuPy expression measures launch enqueue, not execution; redo with
cp.cuda.Event(or asynchronize()bracket). Tabulate C++ vs. Numba vs. fused-CuPy effective bandwidth at \(N = 2^{24}\).
Deliverable & expected results
saxpy(C++) and the Python scripts running clean on the Linux box: sanitizer-clean, error-checked, correctness-verified, with the device-properties table and both timing postmortems innotes.md.
| Quantity | Predicted | Measured |
|---|---|---|
| Bytes moved per SAXPY element | 12 B (two 4 B reads + one 4 B write) → \(12N\) total | … |
| What bounds the kernel | memory bandwidth, not FLOPs — \(I = 1/6\) FLOP/byte is far below any machine balance | … |
| Effective bandwidth at \(N = 2^{24}\) | a large fraction of the bandwidth implied by the device query — compute the ceiling from your device_query output |
… |
| “Wrong” timing (no sync) vs. event timing | wrong ≈ launch overhead, µs-scale, nearly independent of \(N\); right grows ∝ \(N\) | … |
| Numba / fused CuPy vs. C++ throughput | near parity once JIT-compiled — same GPU, same memory bus; unfused CuPy worse (temporaries) | … |
Profiling & performance
First Nsight Systems capture of the course: nsys profile -o labs/lab-1-1/captures/saxpy ./saxpy, opened in the GUI. Find the kernel on the timeline — then notice what dwarfs it: the H2D/D2H cudaMemcpy bars — genuine PCIe transfers into and out of the 4090’s dedicated VRAM, the classic discrete-GPU picture — plus context setup. Measure kernel time vs. total copy time off the timeline and write both down; this exact picture is why Lab 1.4 exists. Repeat once for the Numba script (nsys profile uv run python ...) and confirm the same kernel-shaped bar appears with a Numba-mangled name.
Analysis & reconciliation
In notes.md: derive the bandwidth ceiling from the device query, compare each of the three implementations against it, and explain every gap with a named mechanism — launch overhead amortization, JIT warm-up, CuPy temporaries, copy time wrongly included — not hand-waving. Reconcile the wrong-timing number against launch-overhead expectations (µs order). Close with the module’s opening thesis in your own words: what, exactly, did the GPU spend its time on, and how do you know?
Going further
- Sweep \(N\) from \(2^{10}\) to \(2^{26}\) and plot effective bandwidth vs. size — find where launch overhead stops mattering and where bandwidth saturates.
- Read the PTX Numba generated (
inspect_asm()) next tonvcc -ptxoutput for the C++ kernel — same machine underneath, different front doors. - Re-run the event-timed SAXPY while watching
nvidia-smi dmon— first sighting of the power axis that Lab 1.5 makes first-class.