Метки

, , , , , , , ,

Одной из часто встречающихся операций с массивами является редукция: некая операция со всеми элементами массива. Например, суммирование, выбор максимального или минимального элемента и т. п. В MPI и OpenMP стандартные операции редукции реализованы на уровне системы и программисту заморачиваться с ними не нужно. Но в CUDA её нужно делать самому.

Здесь находится пдф отличной статьи в форме презентации разработчиков NVIDIA, где рассмотрены идеи реализации редукции и их оптимизация.

Следующий код реализует суммирование массива x произвольной длины N.

const unsigned int B_SIZE = 5; // 64
const unsigned int BLOCK_SIZE = 1 << B_SIZE;
typedef float TYPE;

// Host code
TYPE norm(int N, TYPE * vec)
{
	TYPE	ret	= 0;
	size_t	size	= N * sizeof(TYPE);
	dim3	dimBlock(BLOCK_SIZE);
	dim3	dimGrid((N + BLOCK_SIZE - 1) / BLOCK_SIZE);

	TYPE * d_x;
	cudaMalloc((void**)&d_x, size);

	cudaMemcpy(d_x, vec, size, cudaMemcpyHostToDevice);
	int N1 = N;
	while (N1 > 1)
	{
		normKernel<BLOCK_SIZE><<<dimGrid, dimBlock>>>(N1, d_x, d_x);
		N1 = (N1 + BLOCK_SIZE - 1) >> B_SIZE;
	}
	cudaMemcpy(&ret, d_x, sizeof(TYPE), cudaMemcpyDeviceToHost);

	cudaFree(d_x);
	return ret;
}
template<unsigned int blockSize>
// Kernel code
__global__ void reduceKernel(int N, TYPE * d_xin, TYPE * d_xout)
{
	__shared__ TYPE sdata[BLOCK_SIZE];

	size_t tid	= threadIdx.x;
	size_t i		= blockDim.x*2*blockIdx.x + threadIdx.x;
	size_t gridSize	= blockSize*2*gridDim.x;
	sdata[tid]	= 0;
	while (i < N)
	{
		if (i + blockDim.x < N) sdata[tid] += d_xin[i] + d_xin[i + blockDim.x];
		else if (i < N) sdata[tid] += d_xin[i];
		i += gridSize;
	}
	__syncthreads();

	if (blockSize >= 512)
	if (tid < 256) sdata[tid] += sdata[tid + 256];
	__syncthreads();
	if (blockSize >= 256)
	if (tid < 128) sdata[tid] += sdata[tid + 128];
	__syncthreads();
	if (blockSize >= 128)
	if (tid < 64) sdata[tid] += sdata[tid + 64];
	__syncthreads();

	if (tid < 32)
	{
		if (blockSize >= 64 ) sdata[tid] += sdata[tid + 32];
		if (blockSize >= 32 ) sdata[tid] += sdata[tid + 16];
		if (blockSize >= 16 ) sdata[tid] += sdata[tid + 8];
		if (blockSize >=  8 ) sdata[tid] += sdata[tid + 4];
		if (blockSize >=  4 ) sdata[tid] += sdata[tid + 2];
		if (blockSize >=  2 ) sdata[tid] += sdata[tid + 1];
	}

	if (tid == 0) d_xout[blockIdx.x] = sdata[0];
}
Реклама