Общение тредов: скалярное произведение и редукция

1. Задача: перемножили — а сложить кто будет?

Скалярное произведение двух векторов — сумма попарных произведений: a·b = a₀b₀ + a₁b₁ + … + aₙ₋₁bₙ₋₁. Первая половина задачи — чистый map из C04: тред i умножает a[i] на b[i], итерации независимы. А вот сумма — паттерн reduce: результат зависит от всех элементов сразу. Наивная попытка «пусть каждый прибавит к общей переменной» — это гонка данных: треды будут одновременно читать-и-перезаписывать одно место, теряя чужие слагаемые (в P02 соседней линии эта авария воспроизводится вживую на Python).

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

#define N 512

__global__ void dot(const float* a, const float* b, float* result)
{
    __shared__ float temp[N];              // общий стол блока (C03)
    int i = threadIdx.x;

    temp[i] = a[i] * b[i];                 // map-часть: каждый — своё произведение

    __syncthreads();                       // планёрка: ждём, пока стол заполнится

    if (i == 0) {                          // сложение — в одни руки
        float sum = 0;
        for (int k = 0; k < N; k++)
            sum += temp[k];
        *result = sum;
    }
}

Без __syncthreads() программа неверна: варпы блока идут не в ногу (C02), и тред 0 мог бы начать суммировать стол, на который варп 7 ещё ничего не положил. Синхронизация барьером даёт гарантию: ни один тред не пройдёт дальше этой строки, пока до неё не дошли все треды блока.

Решение работает — и разбазаривает параллелизм: пока тред 0 крутит цикл из 512 сложений, 511 тредов смотрят. Время работы по-прежнему O(n), как на CPU.

2. Дерево: сумма за логарифм

Сложение ассоциативно: слагаемые можно группировать как угодно. Сгруппируем так, чтобы работа делилась пополам. Шаг 1: каждый тред из первой половины прибавляет к своей ячейке ячейку из второй половины — 256 сложений одновременно. Шаг 2: то же с половиной оставшегося — 128 сложений. И так далее: 512 → 256 → 128 → … → 1 за log₂ 512 = 9 шагов вместо 512.

__global__ void dot_tree(const float* a, const float* b, float* result)
{
    __shared__ float temp[N];
    int i = threadIdx.x;

    temp[i] = a[i] * b[i];
    __syncthreads();

    // дерево: на каждом шаге работает первая половина оставшихся тредов
    for (int stride = blockDim.x / 2; stride > 0; stride /= 2) {
        if (i < stride)
            temp[i] += temp[i + stride];
        __syncthreads();               // барьер ПОСЛЕ каждого шага — вне if!
    }

    if (i == 0) *result = temp[0];     // ответ скопился в нулевой ячейке
}

Присмотритесь к позиции барьера: __syncthreads() стоит после if, на пути каждого треда. Так задумано, и это не мелочь:

Почему деление именно «первая половина + вторая половина», а не «соседние пары» (temp[2i] += temp[2i+1])? Оба варианта математически верны, но соседние пары на каждом шаге оставляют работающие треды через одного — варп из 32 тредов, где активен каждый второй, всё равно занимает такт целиком (SIMT, C02). Схема со stride держит работающие треды плотной пачкой в младших варпах: остальные варпы отключаются целиком и не тратят тактов. Плюс, как бонус, подряд идущие треды читают подряд идущие ячейки — привет банкам из C03: конфликтов нет.

3. Больше одного блока: atomicAdd

Пока вектор помещался в один блок (N ≤ 1024), хватало разделяемой памяти. Настоящие векторы длиннее, значит, блоков много, а __syncthreads между блоками не существует — блоки принципиально независимы (C02). Схема такая: каждый блок деревом сводит свой кусок к одной частичной сумме, а частичные суммы блоков собираются в глобальной памяти атомарным сложением:

__global__ void dot_big(const float* a, const float* b, float* result, int n)
{
    __shared__ float temp[256];
    int tid = threadIdx.x;
    int i = blockIdx.x * blockDim.x + threadIdx.x;    // глобальный индекс (C04)

    temp[tid] = (i < n) ? a[i] * b[i] : 0.0f;         // хвостовые треды кладут ноль
    __syncthreads();

    for (int stride = blockDim.x / 2; stride > 0; stride /= 2) {
        if (tid < stride) temp[tid] += temp[tid + stride];
        __syncthreads();
    }

    if (tid == 0)
        atomicAdd(result, temp[0]);    // частичную сумму блока — в общий результат
}

atomicAdd(адрес, x) — «прочитать, прибавить, записать» как одна неделимая операция: два блока, добравшиеся до result одновременно, выстроятся в очередь, и ни одно слагаемое не потеряется. Это ровно то, чего не хватало наивному «пусть каждый прибавит к общей переменной» — но заметьте цену: атомарные операции сериализуются. Здесь их всего по одной на блок (десятки штук) — дёшево. Делать atomicAdd каждым тредом — значит выстроить в очередь миллионы операций и убить параллелизм; дерево в разделяемой памяти существует именно затем, чтобы атомарной осталась одна капля на блок. Перед запуском ядра не забудьте обнулить result (cudaMemset(result, 0, sizeof(float))).

Замечание для дотошных: сложение float неассоциативно на уровне округлений, поэтому редукция и последовательная сумма могут расходиться в последних битах — и порядок блоков от запуска к запуску меняет эти биты. Для проверок сравнивайте с допуском (fabs(x − y) < 1e-5 · масштаб), а не оператором ==.

4. Тренажёр: дерево на живых числах

Контрольные вопросы

• Почему «пусть каждый тред прибавит своё произведение к общей переменной» — не решение?
  Это гонка данных: операция «прочитать, прибавить, записать» не атомарна, и два треда, прочитавшие одно значение одновременно, перезапишут результаты друг друга — слагаемые потеряются. Итог неверен и меняется от запуска к запуску.
• Что гарантирует __syncthreads() и какова область его действия?
  Барьер: ни один тред блока не продолжит выполнение, пока к этой строке не пришли все треды блока. Действует только внутри блока — синхронизации между блоками не существует, блоки независимы.
• Почему __syncthreads() нельзя ставить внутри ветки if, куда заходит не весь блок?
  Барьер ждёт все треды блока. Если половина тредов в ветку не зашла, она до барьера не дойдёт — оставшиеся будут ждать вечно: блок зависает либо ведёт себя неопределённо. Барьеры ставят на общем пути всех тредов.
• За сколько шагов дерево сводит блок из 1024 тредов и почему stride-схема лучше «соседних пар»?
  За log₂ 1024 = 10 шагов. Stride-схема держит активные треды плотной пачкой в младших варпах — остальные варпы отключаются целиком и не занимают такты; «соседние пары» оставляют активным каждый второй тред, и все варпы продолжают работать вполсилы. Вдобавок stride-доступ идёт по подряд идущим ячейкам — без конфликтов банков.
• Зачем нужен atomicAdd в межблочной версии и почему им не заменяют всё дерево?
  Между блоками нет барьера, и частичные суммы блоков собирают атомарным сложением в глобальной памяти — очередь гарантирует, что слагаемые не потеряются. Но атомарные операции сериализуются: одна на блок — дёшево, по одной на тред — миллионная очередь вместо параллелизма. Дерево существует, чтобы атомарной осталась капля на блок.
• Почему результат редукции float может отличаться от последовательной суммы и как правильно сравнивать?
  Сложение float неассоциативно из-за округлений, а редукция меняет порядок сложений (и порядок блоков ещё и разный от запуска к запуску). Расхождение в последних битах — норма. Сравнивать нужно с допуском: |x − y| меньше малого эпсилона, умноженного на масштаб значений, а не оператором ==.

Источники

1. Optimizing Parallel Reduction in CUDA / M. Harris // NVIDIA Developer Technology : [презентация]. — URL: https://developer.download.nvidia.com/assets/cuda/files/reduction.pdf (дата обращения: 09.07.2026).
2. CUDA C++ Programming Guide. Synchronization Functions; Atomic Functions // NVIDIA Docs : [сайт]. — URL: https://docs.nvidia.com/cuda/cuda-c-programming-guide/ (дата обращения: 09.07.2026).
3. Using Shared Memory in CUDA C/C++ // NVIDIA Developer Blog : [сайт]. — URL: https://developer.nvidia.com/blog/using-shared-memory-cuda-cc/ (дата обращения: 09.07.2026).