Lab 1.3 — The Reduction Ladder: Naive to Warp Shuffles

← Course 4 syllabus · Module 1 · Prev: « Lab 1.2 · Next: Lab 1.4 »

Goal

The canonical GPU optimization exercise: a parallel sum reduction, climbed rung by rung — interleaved addressing → sequential addressing → first-add-during-load → a warp-shuffle finish with __shfl_down_sync — each rung measured, each speedup explained by a named mechanism (divergence, bank conflicts, idle threads, shared-memory traffic). Then its sibling primitive, the exclusive scan (prefix sum), built once with the standard two-phase work-efficient structure. Reductions and scans are the building blocks of GPU algorithms — compaction, sorting, culling — and this lab’s scan is reused directly by Lab 5.5’s splat depth-sort, so build it as a component, not a demo. The ladder ends where honest engineering does: compared against CUB/Thrust and cupy.sum, the libraries you should usually call instead — the point of climbing by hand is to read kernels like theirs, not to beat them. A numerical thread runs through the lab too: floating-point summation order changes the answer (Course 1 §3), and a parallel reduction is a different — usually better — summation order than the sequential loop.

Recommended reading

• Motta — the material on parallel reduction / cooperative algorithms and warp-level programming (title-level reference; confirm against the copy in hand).
• CUDA C++ Programming Guide — Warp Shuffle Functions (the _sync variants, mask semantics) and the warp-level-primitives discussion; also the section on synchronization within a block.
• CUDA C++ Best Practices Guide — the reduction-flavored parts of the execution-configuration and instruction-optimization sections; NVIDIA’s classic “Optimizing Parallel Reduction in CUDA” slide deck is the free ur-text of this exact ladder and worth an hour.
• Numba CUDA documentation — the warp intrinsics (cuda.shfl_down_sync) and the @cuda.reduce decorator (the convenience path you’ll compare against).
• CuPy documentation — cupy.sum, cupy.cumsum, and the ReductionKernel user-defined-kernel interface.
• Course 1 §3 — conditioning of summation; recall pairwise vs. recursive error growth before predicting the accuracy row.

Prerequisites

• Lab 1.1 (event timing, error-check discipline, effective-bandwidth arithmetic) and Lab 1.2 (shared memory, __syncthreads(), bank-conflict eyes) — this ladder uses all of it.
• The Lab 0.1 CMake conventions; repo building on the Linux desktop (RTX 4090).

Project & environment setup

• New C++ target reduce in cuda/ (rungs selectable by flag, so one binary sweeps the ladder), plus scan. Thrust ships with the CUDA toolkit and CUB ships inside Thrust — #include <cub/cub.cuh> needs no new dependency; note the include in notes.md when you add it.
• Keep --ptxas-options=-v on; per-rung register counts feed Lab 1.5.
• Python: same venv; a reduce_numba.py covering the shuffle-based kernel, @cuda.reduce, cupy.sum, and a cp.ReductionKernel variant.

Where results go:

Artifact	Path
Notes, per-rung table, accuracy experiment, reconciliation	labs/lab-1-3/notes.md
Nsight Compute per-rung reports	labs/lab-1-3/captures/
Ladder timing CSV (C++ rungs, CUB, Thrust, Python variants)	labs/lab-1-3/benchmarks/

Background

• Why reduction is the teaching kernel. A sum touches each element exactly once — \(4N\) bytes moved, \(N\) FLOPs — so the ideal reduction is purely memory-bound and its ceiling is the same bandwidth number SAXPY found. Every rung of the ladder removes an obstacle standing between the naive kernel and that ceiling; when a rung stops helping, you’ve arrived. That is the whole optimization-loop lesson in miniature.
• The rungs, and what each one fixes. Interleaved addressing: tree reduction in shared memory with a modulo-indexed stride — heavy warp divergence (every other thread idles within active warps) and bank conflicts. Sequential addressing: reverse the stride direction so active threads stay contiguous — divergence gone, conflicts gone, but half the threads do nothing after their first add. First add during load: each thread adds two (or more) global elements while loading — the block does real work with every thread and global traffic per block halves. Warp shuffle: the last 32-wide stage swaps __syncthreads()+shared memory for __shfl_down_sync, moving values register-to-register within the warp — this is where Course 3 Module 5’s horizontal NEON reductions reappear wearing a warp.
• Scan. The exclusive prefix sum is reduction’s non-obvious sibling: the standard work-efficient form is an up-sweep (a reduction tree that leaves partial sums in place) then a down-sweep that distributes prefixes back — \(O(N)\) work, two passes, plus a block-partials pass to stitch blocks together. The structure is specified here; every line is yours.
• Order changes the answer. Sequential summation has error growth like \(O(N\varepsilon)\) in the worst case; a tree reduction is effectively pairwise summation, with much slower error growth. So the GPU result at large \(N\) will generally be closer to the double reference than the naive CPU float loop is — a prediction worth savoring because it surprises people.

Tasks

CUDA C++

1. CPU references. A float sequential sum, a double sequential sum (the accuracy referee), and a Kahan-compensated float sum, over the same \(N = 2^{24}\) uniform-random array.
2. Climb the ladder. Implement the four rungs above (shared-memory tree first, shuffle finish last), each verified against the double reference within an argued tolerance, each event-timed at \(N = 2^{24}\), each converted to effective bandwidth via \(4N/t\). One table row per rung.
3. Exclusive scan. Up-sweep/down-sweep within a block, then the block-partials stitch (scan of block sums + uniform add). Verify against np.cumsum-style CPU reference at sizes that cross the multi-block boundary. This is the component Lab 5.5 imports — leave it with a callable interface, not a main().
4. The library bar. cub::DeviceReduce::Sum and thrust::reduce on the same data, same timing rules; add cub::DeviceScan::ExclusiveSum for the scan. Record where your best rung lands relative to CUB.
5. Accuracy experiment. For \(N = 2^{24}\): absolute error vs. the double reference for — sequential float, Kahan, your ladder’s result, CUB’s result. Explain the ordering before you run it.

CUDA Python

6. Shuffle reduction in Numba. Port the top rung using cuda.shfl_down_sync — confirming Numba exposes the same warp machinery — verify and time it.
7. The convenience tier. @cuda.reduce, cupy.sum, and a custom cp.ReductionKernel; time all three. The lesson to extract: what the one-liner costs (or doesn’t) relative to your hand-built rung, and when that trade is obviously right.

Deliverable & expected results

• The rung table (time, effective bandwidth, % of ladder-best) for four rungs + CUB + Thrust + the three Python variants; the accuracy table; the scan verified and packaged for reuse.

Quantity	Predicted	Measured
Ladder ceiling	approaches memory-bandwidth bound — each element read once, \(4N\) bytes, so best-rung effective bandwidth should approach the Lab 1.1 SAXPY figure	…
Biggest single rung-to-rung jump	fixing interleaved addressing (divergence + conflicts) — qualitative, argue it first	…
CUB / Thrust vs. best hand rung	within a few percent of the best hand version — same bandwidth wall, professionally tuned approach to it	…
Error, sequential float vs. tree reduction vs. Kahan	tree ≪ sequential (pairwise-order growth), Kahan best or tied — direction, not magnitudes	…
Numba shuffle rung vs. C++ shuffle rung	near parity — same intrinsics, same bus	…

Profiling & performance

Nsight Compute across the ladder: one report per rung (ncu --set full on rungs 1 and 4 at minimum, --section Occupancy --section MemoryWorkloadAnalysis for the middle rungs). Watch three quantities travel as you climb: warp-execution efficiency / divergence metrics (rung 1’s disease), shared-memory traffic (falls off a cliff at the shuffle rung), and achieved memory throughput (should rise toward the bus limit). Save all reports; Lab 1.5 re-reads them through the occupancy and roofline lenses.

Analysis & reconciliation

Per rung, one paragraph: the mechanism the rung was supposed to fix, the counter in the Nsight report that proves it fixed it, and the fraction of the remaining gap-to-bandwidth it closed. Reconcile the accuracy table against the summation-order argument — if Kahan and the tree land in an unexpected order, work out why with a small hand example. Finish with the judgment call, stated plainly: given CUB’s margin over your best rung, when is hand-writing a reduction ever the right engineering decision, and what did climbing the ladder buy you anyway?

Going further

• Add cooperative_groups’ tiled partitions as a cleaner spelling of the warp stage, and a single-pass grid-synchronized reduction where supported — compare against the two-kernel stitch.
• Extend the scan to segmented scan (per-key prefix sums) — the exact shape GPU sorting and Lab 5.5’s splat binning want.
• Reduce in half precision with a float accumulator and revisit the accuracy table — a preview of the mixed-precision trades Module 6 meets in rendering.