15. Cache-Oblivious Algorithms

What this lecture covers

The previous lecture introduced cache-oblivious matrix multiply. This lecture extends the idea to two more problems: simulating the heat equation with a stencil computation, and sorting. For heat diffusion, it derives a naive looping implementation, measures its poor cache behavior, and replaces it with a trapezoidal divide-and-conquer decomposition that is cache-oblivious in any number of dimensions. For sorting, it walks from ordinary two-way merge sort through cache-aware multiway merge sort (tuned with a "voodoo parameter") to funnel sort, a cache-oblivious algorithm with the same optimal cache complexity.

Along the way, live demos compare looping and trapezoidal code in serial and parallel settings, exposing a surprising result about hardware prefetching and memory bandwidth, and the lecture closes with a practical checklist for diagnosing why a parallel program isn't speeding up as expected.

Key ideas

• Stencil computation: an update rule, like the finite-difference approximation of the heat equation, that recomputes each grid point from a small fixed neighborhood of points at the previous time step.
• Trapezoidal decomposition: cache-oblivious stencil code divides the space-time domain into trapezoids using space cuts (slope -1/+1, when a trapezoid is too wide) and time cuts (horizontal, when too tall), recursing until a trapezoid's base fits in cache.
• D-dimensional cache-oblivious bound: the trapezoidal algorithm achieves Theta(NT / (M^(1/D) * B)) cache misses versus Theta(NT/B) for naive looping, an improvement by a factor of M^(1/D).
• Parallel space cut: a "V-shaped" cut splits a wide trapezoid into two independent trapezoids that run in parallel plus a dependent one that runs after, giving the stencil recursion some parallelism without breaking correctness.
• Hardware prefetching versus cache-obliviousness: in the serial, single-core demo the cache-oblivious code has far fewer misses but only a modest speedup over looping code, because the looping code's regular access pattern lets hardware prefetching hide much of its memory latency; in parallel, prefetching competes with other cores for memory bandwidth, so the cache-oblivious code's advantage grows.
• Speedup bottleneck checklist: check for insufficient parallelism, scheduler overhead (diagnosable via work/span and tools like Cilkscale), memory bandwidth contention (diagnosable by running P copies of the serial code in parallel), and lock or true/false-sharing contention (hardest to detect).
• Multiway merge sort: merging R sorted sub-arrays at once with a tournament tree, tuned so R = Theta(M/B), saves a factor of Theta(B log M) in cache misses over binary merge sort, at the cost of being cache-aware (not oblivious).
• Funnel sort: a cache-oblivious sorting algorithm (recursively sorting N^(1/3) groups of N^(2/3) elements, merged with a K-funnel) that matches multiway merge sort's optimal cache complexity without tuning any parameter to the machine.

Walkthrough

The heat equation and stencil computation (1:01)

The lecture derives the finite-difference approximation of the 1D heat equation, showing that each new value depends on three neighboring values at the previous time step (a three-point stencil), and demonstrates that a straightforward row-by-row loop needs only two rows of storage but performs poorly in cache once a row exceeds the cache size, since each new row re-evicts the previous one.

The trapezoidal cache-oblivious algorithm (22:31)

The lecture defines trapezoidal regions in space-time whose points can be computed independently of anything outside the trapezoid, then gives the recursive rule: space-cut a too-wide trapezoid along a slope of -1 through its center (left half first, then right, since the right depends on the left); time-cut a too-tall trapezoid horizontally (bottom half first, then top). The base case is height one. A recursion-tree analysis shows this yields Theta(NT/(MB)) cache misses in 1D and generalizes to Theta(NT/(M^(1/D)*B)) in D dimensions, versus Theta(NT/B) for the naive loop.

Live demos and the prefetching surprise (34:39)

A simulated cache trace shows the trapezoidal code finishing with far fewer misses than the looping code. A real-time 2D heat-diffusion demo, however, shows only a modest speedup (about 1450 to 1830 iterations per second), which the lecture attributes to hardware prefetching helping the looping code's regular access pattern in the single-core, memory-bandwidth-rich serial case.

Parallelizing the trapezoidal algorithm (46:58)

Reusing the earlier theorem bounding parallel cache misses by serial misses plus O(steals * M/B), the lecture introduces the parallel space cut (a "V" cut producing two independent trapezoids plus a dependent one) so the stencil computation can run in parallel. On a larger offline benchmark (3000x3000 grid, four cores), the parallel looping code gets only about a 2x speedup while the parallel trapezoidal code gets nearly linear (about 3.96x), because parallel execution creates memory-bandwidth contention that hurts prefetching-reliant looping code far more than the cache-efficient trapezoidal code.

Diagnosing speedup bottlenecks (55:06)

The lecture lists four causes of poor parallel speedup: insufficient parallelism, scheduling overhead, lack of memory bandwidth, and contention (locking or false/true sharing), and gives a practical way to test for bandwidth limits: run several identical copies of the serial program in parallel and see if they slow each other down.

Cache-efficient sorting (59:10)

Starting from ordinary two-way merging (Theta(N/B) cache misses) and merge sort (Theta((N/B) log(N/M)) misses), the lecture introduces multiway merging with a tournament tree, showing that setting the fan-in R = Theta(M/B) gives Theta((N log N)/(B log M)) cache misses, a large improvement over binary merge sort for large inputs, but this "voodoo parameter" R makes it cache-aware. Funnel sort, built from recursively structured K-funnels, matches this optimal bound while remaining cache-oblivious.

Before you watch

• Review the previous lecture's ideal-cache model, tall-cache assumption, and cache-oblivious matrix multiplication.
• Recall the recursion-tree method for solving work and cache-miss recurrences, and the master theorem.
• Be comfortable with the parallel cache-miss theorem relating Q_P to Q_1 and the number of successful steals.

Check your understanding

1. Why must the left trapezoid in a space cut be computed before the right one?
2. How does the cache-miss bound for the D-dimensional trapezoidal algorithm change as D increases, and why?
3. Why did the serial demo show only a modest speedup for the cache-oblivious code despite far fewer cache misses, while the parallel demo showed a much larger gap?
4. What four categories of causes should you check, in order, when a parallel program isn't speeding up as expected?
5. Why does setting the merge fan-in R = Theta(M/B) in multiway merge sort make it cache-aware rather than cache-oblivious, and how does funnel sort avoid this?