CUDA Сокращение на разделяемой памяти с несколькими массивами

Question

Я в настоящее время использую следующую функцию Reduction просуммировать все элементы в массиве с CUDA:

__global__ void reduceSum(int *input, int *input2, int *input3, int *outdata, int size){ 
    extern __shared__ int sdata[]; 

    unsigned int tID = threadIdx.x; 
    unsigned int i = tID + blockIdx.x * (blockDim.x * 2); 
    sdata[tID] = input[i] + input[i + blockDim.x]; 
    __syncthreads(); 

    for (unsigned int stride = blockDim.x/2; stride > 32; stride >>= 1) 
    { 
     if (tID < stride) 
     { 
      sdata[tID] += sdata[tID + stride]; 
     } 
     __syncthreads(); 
    } 

    if (tID < 32){ warpReduce(sdata, tID); } 

    if (tID == 0) 
    { 
     outdata[blockIdx.x] = sdata[0]; 
    } 
} 

Однако, как вы можете видеть из параметров функции я бы чтобы иметь возможность суммировать три отдельных массива внутри одной функции сокращения. Теперь, очевидно, простой способ сделать это - это запустить Kernel три раза и каждый раз передавать разные массивы, и это будет нормально работать. Я просто пишу это как тестовое ядро ​​только сейчас, но реальное ядро ​​закончит тем, что возьмет массив структур, и мне нужно будет выполнить добавление для всех значений X, Y и Z каждой структуры, поэтому Мне нужно суммировать их все в одном ядре.

Я initalised и выделяется память для всех трех массивов

int test[1000]; 
    std::fill_n(test, 1000, 1); 
    int *d_test; 

    int test2[1000]; 
    std::fill_n(test2, 1000, 2); 
    int *d_test2; 

    int test3[1000]; 
    std::fill_n(test3, 1000, 3); 
    int *d_test3; 

    cudaMalloc((void**)&d_test, 1000 * sizeof(int)); 
    cudaMalloc((void**)&d_test2, 1000 * sizeof(int)); 
    cudaMalloc((void**)&d_test3, 1000 * sizeof(int)); 

Я не уверен, что сетки и блоков измерения следует использовать для такого рода ядро, и я не совсем уверен, как изменить цикл сокращения в разместить данные, как я хочу, то есть выходного массива:

Block 1 Result|Block 2 Result|Block 3 Result|Block 4 Result|Block 5 Result|Block 6 Result| 

     Test Array 1 Sums    Test Array 2 Sums   Test Array 3 Sums   

Я надеюсь, что имеет смысл. Или есть лучший способ иметь только одну функцию сокращения, но иметь возможность вернуть суммирование Struct.X, Struct.Y или struct.Z?

Вот структура:

template <typename T> 
struct planet { 
    T x, y, z; 
    T vx, vy, vz; 
    T mass; 
}; 

Мне нужно сложить все VX и хранить его, все VY и хранить его и все VZ и хранить его.

Answer 1

Или есть лучший способ иметь только одну функцию сокращения, но иметь возможность возвращать суммирование Struct.X, Struct.Y или struct.Z?

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

Сокращения, подобные этому, не имеют большой арифметики/интенсивности работы, поэтому наш фокус на производительность будет в основном вращаться вокруг хранилища данных для эффективного доступа. При доступе к глобальной памяти графические процессоры обычно делают это в больших кусках - 32 байта или 128 байтовых блоков. Чтобы эффективно использовать подсистему памяти, мы хотим использовать все 32 или 128 запрошенных байтов по каждому запросу.

Но подразумеваемая модель хранения данных вашей структуры:

template <typename T> 
struct planet { 
    T x, y, z; 
    T vx, vy, vz; 
    T mass; 
}; 

в значительной степени исключает это. По этой проблеме вы заботитесь о vx, vy и vz.Эти 3 элементы должны быть смежными в пределах данной структуры (элемент), но в массиве из этих структур, они будут отделены друг от друга, необходимого хранения для других элементов структуры, по меньшей мере:

planet0:  T x 
       T y 
       T z    --------------- 
       T vx  <--   ^
       T vy  <--   | 
       T vz  <--  32-byte read 
       T mass     | 
planet1:  T x      | 
       T y      v 
       T z    --------------- 
       T vx  <-- 
       T vy  <-- 
       T vz  <-- 
       T mass 
planet2:  T x 
       T y 
       T z 
       T vx  <-- 
       T vy  <-- 
       T vz  <-- 
       T mass 

(ради Например, если предположить T является float)

Это указывает на то, ключевой недостаток массив структур (форматов AOS) хранения данных в GPU. Доступ к одному и тому же элементу из последовательных структур является неэффективным из-за гранулярности доступа (32 байта) графического процессора. Обычное предложение для работы в таких случаях заключается в преобразовании хранения AOS в SoA (структура массивов):

template <typename T> 
struct planets { 
    T x[N], y[N], z[N]; 
    T vx[N], vy[N], vz[N]; 
    T mass[N]; 
}; 

выше только один из возможных примеров, вероятно, не то, что вы на самом деле использовать, так как структура будет служить мало, потому что у нас будет только одна структура для планет N. Дело в том, что теперь, когда я получаю доступ к vx для последовательных планет, отдельные элементы vx смежны в памяти, поэтому 32-байтовое чтение дает мне 32 байта vx данных без каких-либо неиспользуемых или неиспользуемых элементов.

С таким преобразованием проблема сокращения становится относительно простой снова, с точки зрения организации кода. Вы можете использовать по существу то же самое, что и ваш код сокращения одного массива, который называется 3 раза подряд, либо с прямым расширением кода ядра, чтобы по существу обрабатывать все 3 массива независимо. Ядро «3-в-1» может выглядеть следующим образом:

template <typename T> 
__global__ void reduceSum(T *input_vx, T *input_vy, T *input_vz, T *outdata_vx, T *outdata_vy, T *outdata_vz, int size){ 
    extern __shared__ T sdata[]; 

    const int VX = 0; 
    const int VY = blockDim.x; 
    const int VZ = 2*blockDim.x; 

    unsigned int tID = threadIdx.x; 
    unsigned int i = tID + blockIdx.x * (blockDim.x * 2); 
    sdata[tID+VX] = input_vx[i] + input_vx[i + blockDim.x]; 
    sdata[tID+VY] = input_vy[i] + input_vy[i + blockDim.x]; 
    sdata[tID+VZ] = input_vz[i] + input_vz[i + blockDim.x]; 
    __syncthreads(); 

    for (unsigned int stride = blockDim.x/2; stride > 32; stride >>= 1) 
    { 
     if (tID < stride) 
     { 
      sdata[tID+VX] += sdata[tID+VX + stride]; 
      sdata[tID+VY] += sdata[tID+VY + stride]; 
      sdata[tID+VZ] += sdata[tID+VZ + stride]; 
     } 
     __syncthreads(); 
    } 

    if (tID < 32){ warpReduce(sdata+VX, tID); } 
    if (tID < 32){ warpReduce(sdata+VY, tID); } 
    if (tID < 32){ warpReduce(sdata+VZ, tID); } 

    if (tID == 0) 
    { 
     outdata_vx[blockIdx.x] = sdata[VX]; 
     outdata_vy[blockIdx.x] = sdata[VY]; 
     outdata_vz[blockIdx.x] = sdata[VZ]; 
    } 
} 

(закодировано в браузере - не тестировались - просто продолжение того, что вы показали, как «эталонное ядро»)

выше AoS -> Преобразование данных SoA, скорее всего, принесет производительность и в другом месте вашего кода. Поскольку предлагаемое ядро ​​будет обрабатывать 3 массивов сразу, размеры сетки и блока должны быть точно так же как то, что вы используете бы для справки ядра в случае одного массива. Для хранения общего хранилища необходимо увеличить (три раза) на каждый блок.

Answer 2

Robert Crovella дал отличный ответ, в котором подчеркивается важность трансформации AoS -> SoA, которая часто повышает производительность на графическом процессоре, я хотел бы предложить среднюю площадку, которая может быть более удобной. Язык CUDA предоставляет несколько типов векторов только для цели, которую вы описываете (см. this section of the CUDA programming guide).

Например, CUDA определяет int3, тип данных, который хранит 3 целых числа.

struct int3 
{ 
    int x; int y; int z; 
}; 

Подобные типы существуют для поплавков, гольцов, двойников и т.д. Что приятно об этих типах данных является то, что они могут быть загружены с одной командой, которая может дать вам небольшой прирост производительности. См. this NVIDIA blog post для обсуждения этого вопроса. Это также более «естественный» тип данных для этого случая, и это может облегчить работу с другими частями вашего кода. Можно определить, например:

struct planets { 
    float3 position[N]; 
    float3 velocity[N]; 
    int mass[N]; 
}; 

Ядро сокращение, которое использует этот тип данных может выглядеть примерно так (адаптированный Роберта).

__inline__ __device__ void SumInt3(int3 const & input1, int3 const & input2, int3 & result) 
{ 
    result.x = input1.x + input2.x; 
    result.y = input1.y + input2.y; 
    result.z = input1.z + input2.z; 
} 

__inline__ __device__ void WarpReduceInt3(int3 const & input, int3 & output, unsigned int const tID) 
{ 
    output.x = WarpReduce(input.x, tID); 
    output.y = WarpReduce(input.y, tID); 
    output.z = WarpReduce(input.z, tID);  
} 

__global__ void reduceSum(int3 * inputData, int3 * output, int size){ 
    extern __shared__ int3 sdata[]; 

    int3 temp; 

    unsigned int tID = threadIdx.x; 
    unsigned int i = tID + blockIdx.x * (blockDim.x * 2); 

    // Load and sum two integer triplets, store the answer in temp. 
    SumInt3(input[i], input[i + blockDim.x], temp); 

    // Write the temporary answer to shared memory. 
    sData[tID] = temp; 

    __syncthreads(); 

    for (unsigned int stride = blockDim.x/2; stride > 32; stride >>= 1) 
    { 
     if (tID < stride) 
     { 
      SumInt3(sdata[tID], sdata[tID + stride], temp); 
      sData[tID] = temp; 
     } 
     __syncthreads(); 
    } 

    // Sum the intermediate results accross a warp. 
    // No need to write the answer to shared memory, 
    // as only the contribution from tID == 0 will matter. 
    if (tID < 32) 
    { 
     WarpReduceInt3(sdata[tID], tID, temp); 
    } 

    if (tID == 0) 
    { 
     output[blockIdx.x] = temp; 
    } 
}