Кэш и порядок по столбцам: как устроена скорость

Главный секрет быстрого численного кода — не процессор, а память: правильный порядок обхода массивов с учётом column-major хранения и кэша ускоряет программу в разы.

Column-major — порядок хранения многомерных массивов Fortran, при котором элементы одного столбца лежат в памяти подряд; локальность доступа — обращение к близко расположенным в памяти данным, что критично для эффективности кэша.

Почему память важнее арифметики

Современный процессор выполняет арифметику чудовищно быстро — миллиарды операций в секунду. Но он простаивает, ожидая данные из памяти: обращение к оперативной памяти занимает сотни тактов, тогда как умножение — единицы. Чтобы скрыть эту пропасть, между процессором и памятью стоит иерархия кэшей — маленькая, но очень быстрая память (L1, L2, L3). Когда процессор читает данные, он подтягивает в кэш не один элемент, а целую кэш-линию (обычно 64 байта — это 8 чисел real(8)). Если следующие нужные данные уже в этой линии, доступ почти мгновенный (попадание, cache hit); если нет — снова сотни тактов ожидания (промах, cache miss). Поэтому скорость численного кода определяется в первую очередь тем, насколько последовательно он обходит память. Алгоритм с той же арифметикой, но дружелюбным к кэшу обходом, бывает в 5-10 раз быстрее.

Column-major: фундаментальное свойство Fortran

Многомерный массив физически хранится в одномерной памяти, и язык должен выбрать порядок развёртки. Fortran использует column-major: быстрее всего меняется первый индекс. Для матрицы A(3,2) элементы лежат в памяти так: A(1,1), A(2,1), A(3,1), A(1,2), A(2,2), A(3,2) — сначала весь первый столбец, потом второй. Это противоположно языку C, где быстрее меняется последний индекс (row-major) и подряд лежат строки. Различие — не каприз, а историческое решение Fortran, и его нужно держать в голове постоянно, потому что оно диктует, как эффективно обходить массивы.

! Память для A(3,2) при column-major:
! адрес:   0      1      2      3      4      5
! элемент: A(1,1) A(2,1) A(3,1) A(1,2) A(2,2) A(3,2)
!          \___ столбец 1 ___/  \___ столбец 2 ___/

Правильный порядок вложенных циклов

Отсюда — важнейшее практическое правило. Когда вы обходите матрицу вложенными циклами, внутренний цикл должен идти по первому индексу (по строкам внутри столбца), чтобы доступ к памяти был последовательным. Сравним два варианта обнуления/обработки матрицы.

real :: a(n, n)
integer :: i, j

! ПЛОХО для Fortran: внутренний цикл по второму индексу j
! прыгает по памяти с шагом n — много кэш-промахов
do i = 1, n
  do j = 1, n
    a(i, j) = compute(i, j)
  end do
end do

! ХОРОШО для Fortran: внутренний цикл по первому индексу i
! идёт по памяти подряд — кэш-линии используются полностью
do j = 1, n
  do i = 1, n
    a(i, j) = compute(i, j)
  end do
end do

Разница огромна. В «плохом» варианте, перебирая j при фиксированном i, мы скачем по памяти с шагом n элементов — каждый доступ почти наверняка в новой кэш-линии, и для больших n кэш постоянно промахивается. В «хорошем» варианте, перебирая i при фиксированном j, мы движемся по памяти подряд: загрузив одну кэш-линию, используем все 8 её чисел, прежде чем грузить следующую. На матрице несколько тысяч на несколько тысяч «хороший» порядок даёт ускорение в разы — без единого изменения в арифметике. Это, пожалуй, самая важная оптимизация в численном Fortran, и ровно здесь чаще всего ошибаются программисты, пришедшие из C, где правило обратное.

Связь с операциями над массивами

Хорошая новость: операции над массивами целиком (a = b + c) и встроенные функции (sum, matmul) уже знают про column-major и обходят память оптимально. Поэтому «массивный» стиль Fortran не только короче, но и автоматически кэш-дружелюбен — ещё один довод писать a = b + c вместо ручных циклов. Когда же цикл необходим (сложная логика внутри), помните о порядке индексов. Срез столбца a(:, j) — это непрерывный кусок памяти, и операции над ним эффективны; срез строки a(i, :) — разреженный (с шагом n), и работа с ним медленнее. Структурируйте алгоритмы так, чтобы «горячие» данные были столбцами, а не строками.

Локальность во времени и пространстве

Кэш вознаграждает два вида локальности. Пространственная локальность — обращение к соседним адресам (ровно её обеспечивает правильный порядок циклов). Временная локальность — повторное обращение к одним и тем же данным за короткое время (пока они ещё в кэше). Блочные алгоритмы линейной алгебры (из урока про BLAS) эксплуатируют именно временную локальность: разбивают задачу на блоки, помещающиеся в кэш, и многократно используют каждый блок, пока он горяч. Понимание обоих видов локальности позволяет проектировать быстрые алгоритмы: минимизировать «прыжки» по памяти и максимизировать переиспользование загруженных данных.

ПриёмКакую локальность улучшаетЭффект
Внутренний цикл по 1-му индексупространственнуюпоследовательный доступ, меньше промахов
Операции над срезами столбцов a(:,j)пространственнуюнепрерывная память
Блочные алгоритмывременнуюпереиспользование данных в кэше
Слияние циклов (fusion)временнуюодин проход вместо нескольких

Как работает под капотом

Когда процессор обращается к адресу, контроллер кэша проверяет, есть ли содержащая его кэш-линия в кэше. При промахе линия (64 байта) целиком загружается из памяти, вытесняя какую-то старую. Аппаратный предвыборщик (prefetcher) распознаёт последовательный доступ и начинает подгружать следующие линии заранее — но только если шаг доступа регулярный и небольшой. Последовательный обход (правильный порядок циклов) активирует предвыборщик, и данные приходят до того, как понадобятся; хаотичный обход с большим шагом обманывает предвыборщик, и процессор простаивает. Кроме того, при последовательном доступе каждая загруженная линия используется на все 8 чисел, а при шаге n — лишь на одно из восьми, то есть 7/8 пропускной способности памяти тратится впустую. Вот арифметика, стоящая за «ускорением в разы»: дело не в числе операций, а в том, сколько раз процессор ждёт память и насколько полно использует каждую загруженную линию.

Иерархия памяти в цифрах

Чтобы интуиция о кэше стала конкретной, полезно увидеть реальные масштабы задержек разных уровней памяти. Они различаются на порядки, и именно эта разница объясняет, почему расположение данных важнее числа арифметических операций.

УровеньТипичная задержкаОбъём
Регистры процессора< 1 тактдесятки значений
Кэш L1~4 такта~32-64 КБ
Кэш L2~12 тактов~256 КБ-1 МБ
Кэш L3~40 тактов~8-32 МБ
Оперативная память~200-300 тактовгигабайты

Разница между попаданием в L1 (4 такта) и обращением к оперативной памяти (300 тактов) — почти стократная. За время одного промаха в память процессор мог бы выполнить сотни арифметических операций. Поэтому вычислительное ядро, которое постоянно промахивается мимо кэша, простаивает большую часть времени, ожидая данные, — независимо от того, насколько мало в нём операций. И наоборот: код, удерживающий горячие данные в кэше, загружает процессор работой, а не ожиданием. Это переворачивает наивную модель «скорость = число операций»: на современном железе скорость памяти-ориентированного кода = число обращений к медленной памяти, а арифметика почти бесплатна на их фоне. Отсюда и важность column-major порядка обхода: он минимизирует промахи, превращая случайные дальние обращения в последовательные, которые кэш и предвыборщик обслуживают эффективно.

Блочные алгоритмы: проектирование под кэш

Высшая форма оптимизации под память — блочные (tiled) алгоритмы, которые сознательно перестраивают вычисление так, чтобы максимально переиспользовать данные, пока они в кэше. Каноничный пример — умножение матриц. Наивный тройной цикл для больших матриц катастрофически промахивается: пройдя одну строку результата, он успевает вытеснить из кэша данные, которые понадобятся для следующей. Блочный вариант разбивает матрицы на небольшие подматрицы-блоки, помещающиеся в кэш, и перемножает их целиком, прежде чем перейти к следующим. Каждый загруженный блок используется на полную, пока горяч.

! Идея блочного умножения C = A*B (схематично):
integer, parameter :: B = 64        ! размер блока под кэш
do jj = 1, n, B
  do kk = 1, n, B
    do ii = 1, n, B
      ! перемножаем блоки [ii:ii+B] x [kk:kk+B] x [jj:jj+B]
      do j = jj, min(jj+B-1, n)
        do k = kk, min(kk+B-1, n)
          do i = ii, min(ii+B-1, n)
            c(i,j) = c(i,j) + a(i,k) * b(k,j)
          end do
        end do
      end do
    end do
  end do
end do

Тот же объём арифметики, но радикально иной паттерн доступа: маленькие блоки переиспользуются многократно из кэша, а не подгружаются из памяти заново. На больших матрицах блочный вариант обгоняет наивный в разы. Именно так устроены оптимизированные dgemm в BLAS — с многоуровневой блочностью под L1, L2, L3 и регистры. Блочность эксплуатирует временную локальность (повторное использование), тогда как правильный порядок циклов — пространственную (последовательность). Освоив оба принципа, вы понимаете, почему профессиональные численные библиотеки так быстры, и можете применять те же идеи в собственных алгоритмах — например, при обходе сеток в методах конечных разностей блочность по подобластям резко снижает промахи кэша. Главный сдвиг мышления, который даёт этот урок: проектируя численный алгоритм, думайте не только о математике и числе операций, но и о движении данных — откуда и куда едут байты, и как сделать так, чтобы процессор поменьше ждал память. На современном железе именно это, а не количество умножений, чаще всего определяет, будет ваш код быстрым или медленным.

Частые ошибки

  • Перенести C-привычку (row-major) в Fortran. В C внутренний цикл идёт по последнему индексу; в Fortran — по первому. Обратный порядок убивает производительность.
  • Активно работать со срезами строк a(i,:). Они разреженны в памяти (шаг n). Предпочитайте срезы столбцов a(:,j) там, где это горячий путь.
  • Игнорировать порядок индексов «потому что компилятор оптимизирует». Компилятор иногда переставляет циклы, но полагаться на это нельзя; пишите правильный порядок сами.
  • Многократные проходы по большому массиву вместо одного. Несколько отдельных циклов по одному массиву тратят пропускную способность; слияние в один проход улучшает временную локальность.
  • Забывать про column-major при передаче данных в/из C. Матрица из C придёт «транспонированной» с точки зрения Fortran.

Итоги

  • Скорость численного кода определяется памятью и кэшем, а не только арифметикой.
  • Fortran хранит массивы column-major: быстрее всего меняется первый индекс, столбцы лежат подряд.
  • Внутренний цикл должен идти по первому индексу — это даёт последовательный доступ и активирует предвыборку.
  • Операции над массивами и встроенные функции уже кэш-оптимальны — ещё один повод писать массивный стиль.
  • Срезы столбцов a(:,j) непрерывны и быстры; срезы строк a(i,:) разреженны и медленны.
Проверьте себя
1. В Fortran при обходе матрицы a(n,n) вложенными циклами какой индекс должен быть во внутреннем цикле для максимальной скорости?
AВторой индекс (столбцы) — как в C
BПервый индекс (строки внутри столбца) — из-за column-major хранения
CЛюбой, компилятор всё равно оптимизирует
DОба одновременно
2. Почему срез столбца a(:,j) эффективнее среза строки a(i,:) в Fortran?
AСтолбцов всегда меньше, чем строк
BСрез столбца — непрерывный кусок памяти, а срез строки разреженный (с шагом n)
CСрез строки вообще не поддерживается
DРазницы в производительности нет