Cache-aware алгоритмы

Предпосылки: массив (непрерывный блок памяти, формула адреса); иерархия памяти и кеш (кеш-линия, локальность, три типа промахов, latency vs bandwidth).

Возьмём перемножение двух квадратных матриц A и B размера N×N с записью результата в C. Все элементы — double (8 байт). Наивная реализация — тройной цикл:

for (int i = 0; i < N; i++)
    for (int j = 0; j < N; j++)
        for (int k = 0; k < N; k++)
            C[i][j] += A[i][k] * B[k][j];

Сложность O(N³) — для 1024×1024 примерно миллиард умножений и сложений. У программиста, только что прошедшего курс алгоритмов, рассуждение простое: сложность определяет время, операций миллиард, процессор делает несколько миллиардов операций в секунду, программа отработает за доли секунды.

Запустим на двух размерах. 256×256 — порядок десятков миллисекунд, IPC (инструкций за такт) ~2.5. 1024×1024 — несколько секунд, IPC ~0.3–0.5. Объём работы вырос в 64 раза (куб отношения размеров), а время — почти в 500 раз. Процессор тот же, код тот же. Где теряется производительность?

Почему одинаковая сложность даёт разное время

Сложность считает арифметические операции. Но процессор работает не с абстрактными числами — с байтами, которые ему доставили из памяти. Сколько обращений попадёт в кеш, а сколько уйдёт в RAM, сложность не описывает, и именно это определяет время.

У наивного матричного умножения внутренний цикл по k читает строку A[i][*] — элементы подряд в памяти, кеш-линии работают. И столбец B[*][j] — между соседними элементами расстояние N×8 байт (шаг строки матрицы).

При N = 256 шаг — 2048 байт, элементы столбца разбросаны, но суммарно немного: строка A (2 КБ) плюс столбец B (256 линий × 64 байта = 16 КБ) — 18 КБ рабочего набора, с запасом помещается в L1d (32 КБ). Горячие данные остаются в L1, процессор кормится с задержкой 1 нс, наблюдаемое время совпадает с ожиданиями по сложности.

При N = 1024 шаг — 8192 байта, каждый элемент столбца сидит в своей кеш-линии. Столбец — 1024 линии × 64 байта = 64 КБ, L1d не удерживает даже один столбец внутреннего цикла. Хуже: за один проход внешнего цикла по j через кеш проходит вся матрица B — 8 МБ. На десктопном L3 в 25 МБ B занимает треть, и с учётом кода, других потоков и ядер может не дожить до следующего прохода внешнего цикла по i. Когда не доживает — это capacity miss, и каждый элемент столбца стоит уже 100 нс вместо 1 нс.

Один и тот же алгоритм одинаковой асимптотики, но один проход использует локальность, другой — нет. Сложность не различает эти случаи. Влиять на ситуацию можно с двух сторон: менять порядок обхода (какие ячейки трогаются в какой момент) и менять layout данных в памяти (как поля упакованы). Ниже — обе техники.

Tiling: умещаем рабочий набор в L1

Проблема наивного умножения — внутренний цикл за раз трогает столбец высотой N, а рядом с кешем нужно было трогать что-то, что в него помещается. Решение — блочное умножение (tiling, loop tiling). Матрица разбивается на квадратные блоки размера B×B, и перемножение идёт блоками:

for (ii = 0; ii < N; ii += B)
    for (jj = 0; jj < N; jj += B)
        for (kk = 0; kk < N; kk += B)
            // перемножаем блок A[ii..ii+B][kk..kk+B]
            // на блок B[kk..kk+B][jj..jj+B],
            // добавляем в блок C[ii..ii+B][jj..jj+B]
            for (i = ii; i < ii+B; i++)
                for (j = jj; j < jj+B; j++)
                    for (k = kk; k < kk+B; k++)
                        C[i][j] += A[i][k] * B[k][j];

Размер B подбирается так, чтобы три блока (из A, B и C) умещались в L1. При L1 = 32 КБ и double: B = 32 даёт блок 32×32×8 = 8 КБ, три блока — 24 КБ. Места хватает.

Эффект двойной. Первое: внутри блока столбец B — всего 32 элемента, 32 линии × 64 байта = 2 КБ. Помещается в L1 вместе с блоками A и C. Второе: блок B переиспользуется подряд всеми итерациями одного блока C и отпускается, как только этот блок C посчитан. Удерживать всю матрицу B в L3 между внешними итерациями больше не нужно. Предвыборка снова работает: внутри блока обход строк и столбцов идёт по коротким последовательным отрезкам.

Разница между наивной и блочной реализацией для 1024×1024 — 5–10 раз, IPC возвращается к 2–3. Количество операций то же — N³ умножений и сложений. Но доставленных из RAM байт — на порядки меньше.

Библиотеки линейной алгебры (BLAS, Basic Linear Algebra Subprograms — набор базовых операций над векторами и матрицами) именно этим и занимаются: подбирают размер блока под L1/L2 текущего процессора и генерируют код с оптимальным обходом. Матричное умножение из Intel MKL достигает 90–95% теоретического пика вычислительной производительности — и это в первую очередь результат правильной работы с кешем, а не арифметических оптимизаций.

Timsort: кеш-осознанная сортировка

Тот же принцип работает и вне линейной алгебры. Классическая merge sort на массивах больше L3 обходит два больших подмассива параллельно, создавая промахи на каждом шаге слияния: хвосты подмассивов лежат далеко друг от друга, и каждое сравнение приходит из RAM.

Timsort (используется в Python и Java как стандартная сортировка) снимает проблему переформатированием: сначала делит массив на блоки-runs размером 32–64 элемента и сортирует каждый методом вставок (insertion sort — эффективна на почти отсортированных данных), затем сливает готовые блоки. Выбор размера run — не случайность: 32 элемента × 8 байт = 256 байт, четыре кеш-линии, полностью умещаются в L1. Пока блок сортируется, ни одна линия не вытесняется. На слиянии работает другая часть: сливающиеся куски уже компактные, головы обоих — в кеше, а не в RAM.

Асимптотика Timsort — та же O(n log n), что у обычного merge sort. Но константа перед ней на современных процессорах — в 1.5–3 раза меньше именно за счёт того, что фаза сортировки блоков не трогает RAM.

AoS vs SoA: раскладка полей

Tiling меняет порядок обхода. Вторая рукоятка — как вообще разложены поля в памяти. На типичном цикле по массиву структур:

struct Point { float x, y, z; int flags; };  // 16 байт
struct Point points[1000000];                 // ~16 МБ
 
float sum = 0;
for (int i = 0; i < 1000000; i++)
    sum += points[i].x;                       // читаем только x

Структуры лежат подряд (Array of Structures, AoS). Кеш-линия 64 байта содержит 4 структуры. Обращение к points[0].x подтягивает в кеш заодно points[0].y, z, flags, а также points[1..3] полностью. Но цикл использует только x — 4 байта из 16. 75% загруженного из RAM не нужно.

Альтернатива — разложить поля отдельно (Structure of Arrays, SoA):

float x[1000000], y[1000000], z[1000000];
int flags[1000000];
 
float sum = 0;
for (int i = 0; i < 1000000; i++)
    sum += x[i];    // кеш-линия = 16 подряд идущих x

Теперь одна кеш-линия несёт 16 значений x вместо 4. Пропускная способность памяти тратится только на нужное поле. Разница на такой задаче — 2–4 раза, та же O(n) сложность.

Но SoA не бесплатна. Если цикл читает два или три поля подряд (sum += points[i].x + points[i].y + points[i].z), AoS уже выигрывает: все три поля приходят одной линией, а предвыборка прогревает следующую структуру целиком. SoA в этом случае требует трёх отдельных потоков обращений, каждый со своим давлением на bandwidth.

Правило выбора:

  • SoA — когда циклы обрабатывают одно-два поля у всех элементов (физические симуляции, базы данных в колоночном формате, векторизованные численные расчёты).
  • AoS — когда каждый элемент обычно используется целиком (графические объекты, игровые сущности, записи бизнес-логики).

Решение зависит не от «что красивее выглядит в коде», а от преобладающего паттерна доступа. Если в горячем цикле половина полей не трогается, AoS впустую прокачивает память; если трогаются все — SoA добавляет давления на bandwidth без выигрыша.

В базах данных это различие отражается в выборе между row-oriented хранением (строки целиком, удобно OLTP — пиши и читай запись полностью) и column-oriented (колонки отдельно, удобно OLAP — сканируй одну колонку из миллиардов строк). Тот же компромисс, что AoS vs SoA, только на уровне страниц на диске.

B-tree: cache-aware структура данных

Cache-aware подход меняет не только алгоритмы, но и структуры данных. B-tree в базах данных — канонический пример. Каждый узел B-tree содержит десятки или сотни ключей в отсортированном массиве. Поиск внутри узла — бинарный поиск по этому массиву: данные лежат подряд, кеш-линии работают, предвыборка работает.

Если бы вместо B-tree использовалось обычное бинарное дерево поиска (BST), каждый узел содержал бы один ключ и два указателя, разбросанных по памяти — каждый шаг поиска был бы промахом кеша. Для дерева из 10 миллионов элементов:

Структура            Глубина   Доступ на шаг   Время полного поиска
-----------------------------------------------------------------
BST (1 ключ)         ~23       ~100 нс (RAM)   ~2.3 мкс
B-tree (256 ключей)  ~3        ~100 нс + binsearch   ~0.5 мкс

Разница — 4–5 раз, и почти вся она объясняется тем, что B-tree упаковывает работу в несколько кеш-линий на узел вместо одной линии на один ключ. То же самое происходит и на диске: узел B-tree обычно совпадает по размеру со страницей (4–16 КБ), и одно обращение к диску приносит все ключи узла сразу.

Размер узла — не произвольный параметр. Для дерева в оперативной памяти он часто подбирается под линии L2 или L3: типичный узел — 512 байт или 8 кеш-линий. На диске — под размер страницы файловой системы (4–16 КБ). В обоих случаях принцип один: сделать шаг обхода таким, чтобы он приносил как можно больше полезной работы с одной единицей доставки данных.

Когда cache-aware не нужен

Блочный алгоритм сложнее наивного (шесть циклов вместо трёх, параметры блоков зависят от железа, отладка труднее). SoA разносит связанные поля и усложняет работу с элементом как с целым. Применять эти приёмы имеет смысл не всегда.

Малые данные. Если рабочий набор и так помещается в L1 (десятки килобайт), кеш уже делает свою работу — tiling ничего не добавит. Сложность кода растёт без выигрыша.

Compute-bound нагрузки. Если программа тратит большую часть времени на арифметику, а не на ожидание памяти, локальность уже не узкое место. Криптография, хеширование, плотные численные расчёты с коротким состоянием — примеры, где cache-aware даёт 1–2% в лучшем случае.

Код исполняется редко. Оптимизация инициализации, однократного разбора конфига или редкого пути не окупается сложностью поддержки.

Нет стабильного паттерна обхода. Если программа обращается к данным по решению пользователя (поиск в БД по произвольному запросу), фиксированный layout под один паттерн не работает. Здесь нужны другие приёмы — индексы, хеш-таблицы, кеширование результатов.

Перед началом работы полезно показать, что программа действительно упирается в память. Промахи L3 и низкий IPC в профилировщике (perf stat, Intel VTune) — прямые индикаторы. Если эти метрики в норме, cache-aware переработка, скорее всего, бесполезна — узкое место в другом месте.

Граница: cache-oblivious

Cache-aware-алгоритмы требуют параметров: размер блока B в tiling, размер узла в B-tree. Для каждого процессора свои оптимальные значения, и неправильная настройка даёт часть выигрыша обратно.

Параллельная идея — cache-oblivious algorithms (кеш-безразличные). Они построены так, что работают почти оптимально на любом размере кеша, не зная его заранее: вместо одного уровня блокирования используется рекурсивное разбиение задачи, которое автоматически попадает в каждый уровень иерархии. Рекурсивное умножение матриц (Frigo, Leiserson, Prokop, 1999) — классический пример. Граница с тем, что разобрано выше: cache-aware — настраиваемый под конкретный процессор; cache-oblivious — автоматически адаптивный, но требует другого устройства алгоритма.

Для подбора параметров конкретного блока нужны размеры L1/L2/L3 и как адреса раскладываются по наборам линий — это в устройстве кеша. Для матричного умножения вопрос «32 vs 64 vs 128» решается по размеру L1d и ассоциативности конкретной микроархитектуры.

Sources

  • Matteo Frigo, Charles E. Leiserson, Harald Prokop, Sridhar Ramachandran, 1999, Cache-Oblivious Algorithms — FOCS 1999
  • Ulrich Drepper, 2007, What Every Programmer Should Know About Memory — Red Hat, Inc.
  • Goetz Graefe, 2011, Modern B-Tree Techniques — Foundations and Trends in Databases
  • Tim Peters, 2002, Timsort description — Python source (Objects/listsort.txt)