У меня следующая проблема. Я реализовал несколько разных алгоритмов параллельного сокращения, и все они работают правильно, если я уменьшаю только одно значение для каждого ядра. Но теперь мне нужно уменьшить несколько (21), и я просто понятия не имею, почему иногда это работает, а иногда нет.
Выполняемые шаги:
- вычисление соответствующих значений для каждого потока (в примере я просто установил их равными 1, поскольку он показывает такое же поведение)
- загрузить их в общую память
- синхронизировать мои темы внутри блока
- уменьшить значения в разделяемой памяти
Вот полный код, который вы можете просто скопировать и запустить.
#include <stdio.h>
#include <cuda_runtime.h>
// switch the compiler flag if you don't have the sdk's helper_cuda.h file
#if 1
#include "helper_cuda.h"
#else
#define checkCudaErrors(val) (val)
#define getLastCudaError(msg)
#endif
#ifdef __CDT_PARSER__
#define __global__
#define __device__
#define __shared__
#define __host__
#endif
// compute sum of val over num threads
__device__ float localSum(const float& val, volatile float* reductionSpace, const uint& localId)
{
reductionSpace[localId] = val; // load data into shared mem
__syncthreads();
// complete loop unroll
if (localId < 128) reductionSpace[localId] += reductionSpace[localId + 128];
__syncthreads();
if (localId < 64) reductionSpace[localId] += reductionSpace[localId + 64];
__syncthreads();
// within one warp (=32 threads) instructions are SIMD synchronous
// -> __syncthreads() not needed
if (localId < 32)
{
reductionSpace[localId] += reductionSpace[localId + 32];
reductionSpace[localId] += reductionSpace[localId + 16];
reductionSpace[localId] += reductionSpace[localId + 8];
reductionSpace[localId] += reductionSpace[localId + 4];
reductionSpace[localId] += reductionSpace[localId + 2];
reductionSpace[localId] += reductionSpace[localId + 1];
}
## Edit: Here we need to sync in order to guarantee that the thread with ID 0 is also done... ##
__syncthreads();
return reductionSpace[0];
}
__global__ void d_kernel(float* od, int n)
{
extern __shared__ float reductionSpace[];
int g_idx = blockIdx.x * blockDim.x + threadIdx.x;
const unsigned int linId = threadIdx.x;
__shared__ float partialSums[21];
float tmp[6] =
{ 0, 0, 0, 0, 0, 0 };
// for simplification all computations are remove - this version still shows the same behaviour
if (g_idx < n)
{
tmp[0] = 1.0f;
tmp[1] = 1.0f;
tmp[2] = 1.0f;
tmp[3] = 1.0f;
tmp[4] = 1.0f;
tmp[5] = 1.0f;
}
float res = 0.0f;
int c = 0;
for (int i = 0; i < 6; ++i)
{
for (int j = i; j < 6; ++j, ++c)
{
res = tmp[i] * tmp[j];
// compute the sum of the values res for blockDim.x threads. This uses
// the shared memory reductionSpace for calculations
partialSums[c] = localSum(res, reductionSpace, linId);
}
}
__syncthreads();
// write back the sum values for this block
if (linId < 21)
{
atomicAdd(&od[linId], partialSums[linId]);
}
}
int main()
{
int w = 320;
int h = 240;
int n = w * h;
// ------------------------------------------------------------------------------------
float *d_out;
checkCudaErrors(cudaMalloc(&d_out, 21 * sizeof(float)));
float* h_out = new float[21];
int dimBlock = 256;
int dimGrid = (n - 1) / dimBlock + 1;
int sharedMemSize = dimBlock * sizeof(float);
printf("w: %d\n", w);
printf("h: %d\n", h);
printf("dimBlock: %d\n", dimBlock);
printf("dimGrid: %d\n", dimGrid);
printf("sharedMemSize: %d\n", sharedMemSize);
int failcounter = 0;
float target = (float) n;
int c = 0;
// ------------------------------------------------------------------------------------
// run the kernel for 200 times
for (int run = 0; run < 200; ++run)
{
cudaMemset(d_out, 0, 21 * sizeof(float));
d_kernel<<<dimGrid, dimBlock, sharedMemSize>>>(d_out, n);;
getLastCudaError("d_kernel");
checkCudaErrors(cudaMemcpy(h_out, d_out, 21 * sizeof(float), cudaMemcpyDeviceToHost));
// check if the output has target value
// since all threads get value 1 the kernel output corresponds to counting the elements which is w*h=n
bool failed = false;
for (int i = 0; i < 21; ++i)
{
if (abs(h_out[i] - target) > 0.01f)
{
++failcounter;
failed = true;
}
}
// if failed, print the elements to show which one failed
if (failed)
{
c = 0;
for (int i = 0; i < 6; ++i)
{
for (int j = i; j < 6; ++j, ++c)
{
printf("%10.7f ", h_out[c]);
}
printf("\n");
}
}
}
printf("failcounter: %d\n", failcounter);
// ------------------------------------------------------------------------------------
delete[] h_out;
checkCudaErrors(cudaFree(d_out));
// ------------------------------------------------------------------------------------
return 0;
}
Некоторые комментарии:
BlockSize всегда 256, поэтому развернутый цикл в localSum() проверяет правильные идентификаторы потоков. Как упоминалось в начале, из 200 прогонов иногда это полностью правильно, иногда только 2 значения неверны, а иногда 150 или около того неверны.
И ему не нужно ничего делать с точностью с плавающей запятой, поскольку только 1x1 умножается и сохраняется в переменной res в d_kernel(). Я ясно вижу, что иногда просто некоторые потоки или блоки не запускаются, но я не знаю, почему. :/
Просто глядя на результаты, должно быть очевидно, что есть какое-то состояние гонки, но я просто не вижу проблемы.
Кто-нибудь знает, в чем проблема?
Редактировать:
Я протестировал сейчас много вещей и увидел, что это должно что-то делать с BlockSize. Если я уменьшу его до чего-то ‹=64 и соответствующим образом изменю localSum(), тогда все всегда будет работать так, как ожидалось.
Но это просто не имеет смысла для меня?! Я по-прежнему не делаю здесь ничего другого, кроме обычного параллельного редукции с разделяемой памятью с той лишь разницей, что делаю это 21 раз за поток.
Редактировать 2:
Теперь я в полном замешательстве. Проблема заключается в развертывании цикла!! Или, лучше сказать, синхронизации деформации. Работает следующий код localSum():
// compute sum of val over num threads
__device__ float localSum(const float& val, volatile float* reductionSpace, const uint& localId)
{
reductionSpace[localId] = val; // load data into shared mem
__syncthreads();
for (unsigned int s = blockDim.x / 2; s > 0; s >>= 1)
{
if (localId < s)
{
reductionSpace[localId] += reductionSpace[localId + s];
}
__syncthreads();
}
return reductionSpace[0];
}
Но если я разворачиваю последнюю деформацию и не синхронизирую потоки, я снова иногда получаю 2 или 3 неверных результата из 2000 прогонов. Таким образом, следующий код НЕ работает:
// compute sum of val over num threads
__device__ float localSum(const float& val, volatile float* reductionSpace, const uint& localId)
{
reductionSpace[localId] = val; // load data into shared mem
__syncthreads();
for (unsigned int s = blockDim.x / 2; s > 32; s >>= 1)
{
if (localId < s)
{
reductionSpace[localId] += reductionSpace[localId + s];
}
__syncthreads();
}
if (localId < 32)
{
reductionSpace[localId] += reductionSpace[localId + 32];
reductionSpace[localId] += reductionSpace[localId + 16];
reductionSpace[localId] += reductionSpace[localId + 8];
reductionSpace[localId] += reductionSpace[localId + 4];
reductionSpace[localId] += reductionSpace[localId + 2];
reductionSpace[localId] += reductionSpace[localId + 1];
}
return reductionSpace[0];
}
Но какой в этом смысл, если CUDA одновременно выполняет один варп (32 потока) и __syncthreads() не требуется?!
Мне не нужно, чтобы кто-то публиковал мой рабочий код здесь, но я действительно прошу кого-то с большим опытом и глубокими знаниями в программировании CUDA описать мне здесь основную проблему. Или хотя бы подсказать.
