Сокращение CUDA для нахождения максимума массива

#c #arrays #cuda #parallel-processing

#c #массивы #cuda #параллельная обработка

Вопрос:

Я изучаю курс Udacity по параллельному программированию (домашнее задание 3) и не могу понять, почему я не могу получить максимум в массиве, используя параллельное сокращение (форумы Udacity еще не предоставили решение). Я почти уверен, что я правильно настроил массивы и что алгоритм верен. Я подозреваю, что у меня проблема с управлением памятью (доступ за пределы, неправильные размеры массива, копирование туда и обратно). Пожалуйста, помогите! Я запускаю это в среде Udacity, а не локально. Ниже приведен код, который я использую в настоящее время. По какой-то причине, когда я меняю fmaxf ‘s на fminf ‘s, он находит минимум.

 #include "reference_calc.cpp"
#include "utils.h"
#include "math.h"
#include <stdio.h>
#include <cmath>

__global__ void reduce_max_kernel(float *d_out, const float *d_logLum, int size) {

    // Reduce log Lum with Max Operator
    int myId = threadIdx.x   blockDim.x * blockIdx.x;
    int tid  = threadIdx.x;

    extern __shared__ float temp[];

    if (myId < size) {
        temp[tid] = d_logLum[myId];
    }
    else {
        temp[tid] = d_logLum[tid];
    }

    for (unsigned int s = blockDim.x / 2; s > 0; s >>= 1) {
        if (tid < s) {
            if (myId < size) {
                temp[tid] = fmaxf(d_logLum[myId   s], d_logLum[myId]);
            } else {
                temp[tid] = d_logLum[tid];
            }
        }
        __syncthreads(); 
    }

    if (tid == 0) {
        d_out[blockIdx.x] = temp[0];
    }
}

__global__ void reduce_max_kernel2(float *d_out, float *d_in) {

    // Reduce log Lum with Max Operator
    int myId = threadIdx.x   blockDim.x * blockIdx.x;
    int tid  = threadIdx.x;

    for (unsigned int s = blockDim.x >> 1; s > 0; s >>= 1) {
        if (tid < s) {
            d_in[myId] = fmaxf(d_in[myId   s], d_in[myId]);
        }
        __syncthreads();   
    }

    if (tid == 0) {
        d_out[0] = d_in[0];
    }

}


void your_histogram_and_prefixsum(const float* const d_logLuminance,
                                  unsigned int* const d_cdf,
                                  float amp;min_logLum,
                                  float amp;max_logLum,
                                  const size_t numRows,
                                  const size_t numCols,
                                  const size_t numBins)
{
  //TODO
  /*Here are the steps you need to implement
    1) find the minimum and maximum value in the input logLuminance channel
       store in min_logLum and max_logLum
    2) subtract them to find the range
    3) generate a histogram of all the values in the logLuminance channel using
       the formula: bin = (lum[i] - lumMin) / lumRange * numBins
    4) Perform an exclusive scan (prefix sum) on the histogram to get
       the cumulative distribution of luminance values (this should go in the
       incoming d_cdf pointer which already has been allocated for you)       */
    //int size = 1 << 18;
    int points = numRows * numCols;
    int logPoints = ceil(log(points)/log(2));
    int sizePow = logPoints;
    int size = pow(2, sizePow);
    int numThreads = 1024;
    int numBlocks = size / numThreads;

    float *d_out;
    float *d_max_out;

    checkCudaErrors(cudaMalloc((void **) amp;d_out, numBlocks * sizeof(float)));
    checkCudaErrors(cudaMalloc((void **) amp;d_max_out, sizeof(float)));

    cudaDeviceSynchronize();
    reduce_max_kernel<<<numBlocks, numThreads, sizeof(float)*numThreads>>>(d_out, d_logLuminance, points);

    cudaDeviceSynchronize();
    reduce_max_kernel2<<<1, numBlocks>>>(d_max_out, d_out);

    float h_out_max;
    checkCudaErrors(cudaMemcpy(amp;h_out_max, d_max_out, sizeof(float), cudaMemcpyDeviceToHost));

    printf("%fn", h_out_max);

    checkCudaErrors(cudaFree(d_max_out));
    checkCudaErrors(cudaFree(d_out));

}
  

Комментарии:

1. Если вам нужна помощь в отладке вашего кода, нам нужно увидеть кратчайший полный пример, который кто-то другой мог скомпилировать и запустить.

2. В вашем reduce_max_kernel случае обычно рекомендуется использовать a __syncthreads() после заполнения или изменения разделяемой памяти. У вас есть это в вашем цикле reduce for, но, похоже, у вас его нет после начальной загрузки разделяемой памяти (до цикла for). Некоторые потоки могут забежать вперед и войти в цикл for до того, как другие получат возможность выполнить начальную загрузку разделяемой памяти. Кроме того, цель использования общей памяти — фактически использовать ее в вашем цикле сокращения. Прямо сейчас ваш цикл for загружает значения из глобальной, а не из общей памяти. Это нарушено.

Ответ №1:

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

Ниже, пожалуйста, найдите полный рабочий пример, построенный вокруг вашей попытки. Я оставил неправильные инструкции вашего подхода.

 #include <thrustdevice_vector.h>

#define BLOCKSIZE 256

/********************/
/* CUDA ERROR CHECK */
/********************/
#define gpuErrchk(ans) { gpuAssert((ans), __FILE__, __LINE__); }
inline void gpuAssert(cudaError_t code, char *file, int line, bool abort=true)
{
    if (code != cudaSuccess) 
    {
        fprintf(stderr,"GPUassert: %s %s %dn", cudaGetErrorString(code), file, line);
        if (abort) { getchar(); exit(code); }
    }
}

/*******************************************************/
/* CALCULATING THE NEXT POWER OF 2 OF A CERTAIN NUMBER */
/*******************************************************/
unsigned int nextPow2(unsigned int x)
{
    --x;
    x |= x >> 1;
    x |= x >> 2;
    x |= x >> 4;
    x |= x >> 8;
    x |= x >> 16;
    return   x;
}

__global__ void reduce_max_kernel(float *d_out, const float *d_logLum, int size) {

    int tid         = threadIdx.x;                              // Local thread index
    int myId        = blockIdx.x * blockDim.x   threadIdx.x;    // Global thread index

    extern __shared__ float temp[];

    // --- Loading data to shared memory. All the threads contribute to loading the data to shared memory.
    temp[tid] = (myId < size) ? d_logLum[myId] : -FLT_MAX;

    // --- Your solution
    // if (myId < size) { temp[tid] = d_logLum[myId]; } else { temp[tid] = d_logLum[tid]; }

    // --- Before going further, we have to make sure that all the shared memory loads have been completed
    __syncthreads();

    // --- Reduction in shared memory. Only half of the threads contribute to reduction.
    for (unsigned int s=blockDim.x/2; s>0; s>>=1)
    {
        if (tid < s) { temp[tid] = fmaxf(temp[tid], temp[tid   s]); }
        // --- At the end of each iteration loop, we have to make sure that all memory operations have been completed
        __syncthreads();
    }

    // --- Your solution
    //for (unsigned int s = blockDim.x / 2; s > 0; s >>= 1) {
    //    if (tid < s) { if (myId < size) { temp[tid] = fmaxf(d_logLum[myId   s], d_logLum[myId]); } else { temp[tid] = d_logLum[tid]; } }
    //    __syncthreads(); 
    //}

    if (tid == 0) {
        d_out[blockIdx.x] = temp[0];
    }
}

/********/
/* MAIN */
/********/
int main()
{
    const int N = 10;

    thrust::device_vector<float> d_vec(N,3.f); d_vec[4] = 4.f;

    int NumThreads  = (N < BLOCKSIZE) ? nextPow2(N) : BLOCKSIZE;
    int NumBlocks   = (N   NumThreads - 1) / NumThreads;

    // when there is only one warp per block, we need to allocate two warps
    // worth of shared memory so that we don't index shared memory out of bounds
    int smemSize = (NumThreads <= 32) ? 2 * NumThreads * sizeof(int) : NumThreads * sizeof(int);

    // --- reduce2
    thrust::device_vector<float> d_vec_block(NumBlocks);
    reduce_max_kernel<<<NumBlocks, NumThreads, smemSize>>>(thrust::raw_pointer_cast(d_vec_block.data()), thrust::raw_pointer_cast(d_vec.data()), N);

    // --- The last part of the reduction, which would be expensive to perform on the device, is executed on the host
    thrust::host_vector<float> h_vec_block(d_vec_block);
    float result_reduce0 = -FLT_MAX;
    for (int i=0; i<NumBlocks; i  ) result_reduce0 = fmax(h_vec_block[i], result_reduce0);
    printf("Result = %fn",result_reduce0);

}