CUDA — почему параллельное сокращение на основе деформации происходит медленнее?

У меня возникла идея о параллельном сокращении на основе варпа, поскольку все потоки варпа по определению синхронизированы.

Итак, идея заключалась в том, что входные данные можно уменьшить в 64 раза (каждый поток уменьшает два элемента) без необходимости синхронизации.

Как и в исходной реализации Марка Харриса, сокращение применяется на уровне блоков, а данные находятся в общей памяти. http://gpgpu.org/static/sc2007/SC07_CUDA_5_Optimization_Harris.pdf

Я создал ядро, чтобы протестировать его версию и мою версию, основанную на варпе.
Само ядро ​​полностью идентично хранит элементы BLOCK_SIZE в разделяемой памяти и выводит результат по своему уникальному блочному индексу в выходном массиве.

Сам алгоритм работает нормально. Протестировано с полным набором для проверки «подсчета».

Тело функций реализаций:

/**
 * Performs a parallel reduction with operator add 
 * on the given array and writes the result with the thread 0
 * to the given target value
 *
 * @param inValues T* Input float array, length must be a multiple of 2 and equal to blockDim.x
 * @param targetValue float 
 */
__device__ void reductionAddBlockThread_f(float* inValues,
    float &outTargetVar)
{
    // code of the below functions
}

<сильный>1. Реализация его версии:

if (blockDim.x >= 1024 && threadIdx.x < 512)
    inValues[threadIdx.x] += inValues[threadIdx.x + 512];
__syncthreads();
if (blockDim.x >= 512 && threadIdx.x < 256)
    inValues[threadIdx.x] += inValues[threadIdx.x + 256];
__syncthreads();
if (blockDim.x >= 256 && threadIdx.x < 128)
    inValues[threadIdx.x] += inValues[threadIdx.x + 128];
__syncthreads();
if (blockDim.x >= 128 && threadIdx.x < 64)
    inValues[threadIdx.x] += inValues[threadIdx.x + 64];
__syncthreads();

//unroll last warp no sync needed
if (threadIdx.x < 32)
{
    if (blockDim.x >= 64) inValues[threadIdx.x] += inValues[threadIdx.x + 32];
    if (blockDim.x >= 32) inValues[threadIdx.x] += inValues[threadIdx.x + 16];
    if (blockDim.x >= 16) inValues[threadIdx.x] += inValues[threadIdx.x + 8];
    if (blockDim.x >= 8) inValues[threadIdx.x] += inValues[threadIdx.x + 4];
    if (blockDim.x >= 4) inValues[threadIdx.x] += inValues[threadIdx.x + 2];
    if (blockDim.x >= 2) inValues[threadIdx.x] += inValues[threadIdx.x + 1];

    //set final value
    if (threadIdx.x == 0)
        outTargetVar = inValues[0];
}

Ресурсы:

Используются 4 синхреда
Используются 12 операторов if
11 операций чтения + добавления + записи
1 окончательная операция записи
5 использование регистров

Производительность:

среднее время выполнения пяти тестов: ~ 19,54 мс

<сильный>2. Подход на основе деформации: (То же тело функции, что и выше)

/*
 * Perform first warp based reduction by factor of 64
 *
 * 32 Threads per Warp -> LOG2(32) = 5
 *
 * 1024 Threads / 32 Threads per Warp = 32 warps
 * 2 elements compared per thread -> 32 * 2 = 64 elements per warp
 *
 * 1024 Threads/elements divided by 64 = 16
 * 
 * Only half the warps/threads are active
 */
if (threadIdx.x < blockDim.x >> 1)
{
    const unsigned int warpId = threadIdx.x >> 5;
    // alternative threadIdx.x & 31
    const unsigned int threadWarpId = threadIdx.x - (warpId << 5);
    const unsigned int threadWarpOffset = (warpId << 6) + threadWarpId;

    inValues[threadWarpOffset] += inValues[threadWarpOffset + 32];
    inValues[threadWarpOffset] += inValues[threadWarpOffset + 16];
    inValues[threadWarpOffset] += inValues[threadWarpOffset + 8];
    inValues[threadWarpOffset] += inValues[threadWarpOffset + 4];
    inValues[threadWarpOffset] += inValues[threadWarpOffset + 2];
    inValues[threadWarpOffset] += inValues[threadWarpOffset + 1];
}

// synchronize all warps - the local warp result is stored
// at the index of the warp equals the first thread of the warp
__syncthreads();

// use first warp to reduce the 16 warp results to the final one
if (threadIdx.x < 8)
{
    // get first element of a warp
    const unsigned int warpIdx = threadIdx.x << 6;

    if (blockDim.x >= 1024) inValues[warpIdx] += inValues[warpIdx + 512];
    if (blockDim.x >= 512) inValues[warpIdx] += inValues[warpIdx + 256];
    if (blockDim.x >= 256) inValues[warpIdx] += inValues[warpIdx + 128];
    if (blockDim.x >= 128) inValues[warpIdx] += inValues[warpIdx + 64];

    //set final value
    if (threadIdx.x == 0)
        outTargetVar = inValues[0];
}

Ресурсы:

1 используется синхропоток
7 операторов if
10 операций чтения и записи
1 окончательная операция записи
5 использование регистра

5 битовых сдвигов
1 добавить
1 подпрограмму

Производительность:

среднее время пяти тестовых прогонов: ~ 20,82 мс

Многократное тестирование обоих ядер на Geforce 8800 GT 512 МБ с 256 МБ значений с плавающей запятой. И работающее ядро ​​с 256 потоками на блок (100 % загрузка).

Версия на основе деформации примерно на 1,28 миллисекунды медленнее.

Если будущие карты позволят увеличить размер блока, подход, основанный на деформации, по-прежнему не будет нуждаться в дальнейшем операторе синхронизации, поскольку максимальное значение составляет 4096, которое уменьшается до 64, а при окончательной деформации уменьшается до 1

Почему не быстрее?, или где ошибка в идее, ядро?

От использования ресурсов подход warp должен быть впереди?

Редактирование 1: ядро ​​исправлено так, что активна только половина потоков, что не приводит к чтению за пределами границ, добавлены новые данные о производительности


person djmj    schedule 04.10.2012    source источник


Ответы (2)


Я думаю, причина того, что ваш код медленнее моего, заключается в том, что в моем коде вдвое меньше деформаций активны для каждого ADD на первой фазе. В вашем коде все деформации активны для всей первой фазы. Таким образом, в целом ваш код выполняет больше инструкций деформации. В CUDA важно учитывать общее количество выполненных «инструкций деформации», а не только количество инструкций, выполняемых одной деформацией.

Кроме того, нет смысла использовать только половину ваших варпов. Есть накладные расходы на запуск варпов только для того, чтобы они оценили две ветки и вышли.

Другая мысль заключается в том, что использование unsigned char и short может на самом деле стоить вам производительности. Я не уверен, но это точно не экономит ваши регистры, так как они не упакованы в отдельные 32-битные переменные.

Кроме того, в своем исходном коде я заменил blockDim.x параметром шаблона BLOCKDIM, что означает, что он использовал только 5 операторов if во время выполнения (если на втором этапе исключаются компилятором).

Кстати, более дешевый способ вычислить threadWarpId это

const int threadWarpId = threadIdx.x & 31;

Дополнительные идеи можно найти в этой статье.

РЕДАКТИРОВАТЬ: Вот альтернативное сокращение блоков на основе деформации.

template <typename T, int level>
__device__
void sumReduceWarp(volatile T *sdata, const unsigned int tid)
{
  T t = sdata[tid];
  if (level > 5) sdata[tid] = t = t + sdata[tid + 32];
  if (level > 4) sdata[tid] = t = t + sdata[tid + 16];
  if (level > 3) sdata[tid] = t = t + sdata[tid +  8];
  if (level > 2) sdata[tid] = t = t + sdata[tid +  4];
  if (level > 1) sdata[tid] = t = t + sdata[tid +  2];
  if (level > 0) sdata[tid] = t = t + sdata[tid +  1];
}

template <typename T>
__device__
void sumReduceBlock(T *output, volatile T *sdata)
{
  // sdata is a shared array of length 2 * blockDim.x

  const unsigned int warp = threadIdx.x >> 5;
  const unsigned int lane = threadIdx.x & 31;
  const unsigned int tid  = (warp << 6) + lane;

  sumReduceWarp<T, 5>(sdata, tid);
  __syncthreads();

  // lane 0 of each warp now contains the sum of two warp's values
  if (lane == 0) sdata[warp] = sdata[tid];

  __syncthreads();

  if (warp == 0) {
    sumReduceWarp<T, 4>(sdata, threadIdx.x);
    if (lane == 0) *output = sdata[0];
  }
}

Это должно быть немного быстрее, потому что он использует все варпы, запущенные на первом этапе, и не имеет ветвления на последнем этапе за счет дополнительной ветви, общей загрузки/хранения и __syncthreads() на новом промежуточном этапе. Я не тестировал этот код. Если вы запустите его, дайте мне знать, как он работает. Если вы используете шаблон для blockDim в своем исходном коде, он снова может быть быстрее, но я думаю, что этот код более лаконичен.

Обратите внимание, что временная переменная t используется, потому что Fermi и более поздние архитектуры используют чистую архитектуру загрузки/сохранения, поэтому += из общей памяти в общую память приводит к дополнительной загрузке (поскольку указатель sdata должен быть энергозависимым). Явная загрузка во временный файл позволяет избежать этого. На G80 это не повлияет на производительность.

person harrism    schedule 04.10.2012
comment
Я должен извиниться за то, что сказал, что длина массива разделяемой памяти равна BLOCK_SIZE, но она должна быть либо равна размеру двойного блока (что не имеет смысла при инициализации из глобальной памяти), либо только половина потоков должна быть активной. И поскольку все потоки активны, чтение и запись выходят за рамки. Я исправлю это и проверю еще раз. Исправленное ядро ​​теперь точно правильное. Как я уже сказал, он всегда дает правильные математические результаты. Странно то, почему эти обращения к общей памяти из-за границы не привели к неправильным результатам. - person djmj; 05.10.2012
comment
О параметре шаблона BLOCKDIM: Вот почему я опубликовал свою реализацию вашего алгоритма. И мне все еще нужен более дорогой способ вычисления warpId, так как мне нужно, чтобы он вычислял threadWarpOffset. - person djmj; 05.10.2012
comment
Извините, я пропустил 6, а не 5 в третьей строке. Отредактировано. Я думаю, что мои рассуждения о том, почему ваш код работает медленнее, скорее всего, по крайней мере часть причины. - person harrism; 05.10.2012
comment
Из первого неправильного ядра, в котором участвовали все потоки, total warp instructions было определенно выше. Но какое это имеет значение, если активна только половина потоков, а остальные ждут следующего оператора синхронизации? Значит, вместо того, чтобы ждать, они могли бы еще и что-то сделать? Особенно, если потоки в варпе синхронизированы? - person djmj; 05.10.2012
comment
Бездействующие варпы не используют циклы обработки. Активные деформации делают. Ваш код в целом имеет более активные деформации для большего количества инструкций, поэтому общее количество циклов выше. Мой код выполняет те же вычисления за меньшее количество циклов. Я думаю, что деформационное сокращение выполнимо (сам сделал это), но я бы сделал это по-другому. Попробую раскопать код... - person harrism; 05.10.2012
comment
Хм, я думаю, что я понимаю это медленно. Поэтому я должен рассчитать реальное количество необходимых инструкций. О, чувак, это казалось таким простым. Черт, никогда не бывает реального способа сделать это с cuda. О unsigned short и unsigned char. Я заменил оба на int, и это привело к тому, что на один регистр меньше??? О, чувак... ты думаешь, что правильно программируешь все свои ядра и экономишь ресурсы, а потом. Упакован ли vec2 Short в 32-битный регистр? - person djmj; 05.10.2012
comment
Все еще хотелось бы увидеть ваше уменьшение варпа. Просто надеюсь, что размер блока 4096 в будущих выпусках cuda превзойдет вашу версию;) - person djmj; 05.10.2012
comment
Большие блоки не обязательно помогут. Вам нужен достаточный параллелизм, чтобы заполнить машину и скрыть задержку памяти, и не более того. Вот почему быстрый алгоритм выполняет итеративное суммирование в каждом потоке до тех пор, пока общее количество частичных сумм не станет меньше, чем блок потока может суммироваться с деревом. Вы не хотите, чтобы дерево было больше, чем нужно, потому что это неэффективно. - person harrism; 05.10.2012
comment
Добавлен пример уменьшения деформации. Это очень похоже на сканирование, используемое в CUDPP. - person harrism; 05.10.2012
comment
Спасибо, также была идея о if (lane == 0) sdata[warp] = sdata[tid]; для хранения результатов деформации в индексе деформации. (Актуальна ли моя пошаговая индексация в общей памяти? Нет ничего лучше объединенного доступа к общей памяти?) Но я подумал, что использование warpIdx = threadIdx.x << 6; спасет это, если условие и синхронизация. В вашем последнем сокращении участвуют все 32 потока, что не обязательно, поскольку максимально необходимо 8 потоков. Но это и удаленная вами проверка BLOCKDIM делают его более масштабируемым для будущих устройств. Но должна стоить производительность по сравнению с моей. Проверит вашу работоспособность. - person djmj; 05.10.2012
comment
У вас не может быть менее 32 активных потоков, иначе вы оставляете оборудование бездействующим, поэтому нет причин использовать только 8 потоков. Хорошая мысль об индексации шагов. В вашей версии могут быть 8-сторонние конфликты общего банка памяти, что довольно плохо. Никаких конфликтов с банками у меня нет. - person harrism; 06.10.2012
comment
Томас, не обманываете ли вы себя, позволяя себе предположить, что blockDim.x равно 1024? Вместо того, чтобы использовать его как переменную, я имею в виду? Во всяком случае, я отредактировал код, чтобы более четко отразить это предположение. - person einpoklum; 02.12.2013
comment
Кроме того, вы должны использовать эти операторы «если» во время операций +32, +16, +8, +4, +2, +1. В противном случае вы можете получить неправильный результат. Например. на шаге +16 он может иногда добавлять 32-й элемент к 16-му элементу до того, как 16-й элемент будет добавлен к 0-му. - person Olex; 04.03.2017

Вы также должны проверить примеры в SDK. Помню один очень красивый пример с реализациями нескольких способов редукций. По крайней мере, один из них также использует сокращение на основе деформации.

(Я не могу сейчас найти имя, потому что оно установлено только на моей другой машине)

person Michael    schedule 05.10.2012