Я хочу сделать умножение разреженной матрицы на плотный вектор.Предположим, что единственным форматом хранения для сжатия записей в матрице является сжатое хранилище строк CRS.
Мое ядро выглядит следующим образом:
__global__ void
krnlSpMVmul1(
float *data_mat,
int num_nonzeroes,
unsigned int *row_ptr,
float *data_vec,
float *data_result)
{
extern __shared__ float local_result[];
local_result[threadIdx.x] = 0;
float vector_elem = data_vec[blockIdx.x];
unsigned int start_index = row_ptr[blockIdx.x];
unsigned int end_index = row_ptr[blockIdx.x + 1];
for (int index = (start_index + threadIdx.x); (index < end_index) && (index < num_nonzeroes); index += blockDim.x)
local_result[threadIdx.x] += (data_mat[index] * vector_elem);
__syncthreads();
// Reduction
// Writing accumulated sum into result vector
}
Как вы можете видеть, ядро должнобыть настолько наивным, насколько это возможно, и даже делает несколько вещей неправильно (например, vector_elem не всегда является правильным значением).Я знаю об этих вещах.
Теперь к моей проблеме: предположим, я использую размер блока из 32 или 64 потоков.Как только строка в моей матрице содержит более 16 ненулевых элементов (например, 17), только первые 16 умножений выполняются и сохраняются в общей памяти.Я знаю, что значение local_result[16], которое является результатом 17-го умножения, просто равно нулю.Использование размера блока из 16 или 128 потоков решает описанную проблему.
Поскольку я довольно новичок в CUDA, я, возможно, упустил простейшую вещь, но я не могу придумать больше ситуаций, на которые можно посмотреть.
Помощь очень ценится!
Редактировать в комментарии к комментариям:
Я напечатал значения, которые были в local_result[16] непосредственно после вычисления.Это было 0.Тем не менее, вот отсутствующий код:
Часть сокращения:
int k = blockDim.x / 2;
while (k != 0)
{
if (threadIdx.x < k)
local_result[threadIdx.x] += local_result[threadIdx.x + k];
else
return;
__syncthreads();
k /= 2;
}
и как я записываю результаты обратно в глобальную память:
data_result[blockIdx.x] = local_result[0];
Вот и всеГот.
Сейчас я тестирую сценарий с матрицей, состоящей из одной строки с 17 элементами, которые все не равны нулю.Буферы выглядят так в псевдокоде:
float data_mat[17] = { val0, .., val16 }
unsigned int row_ptr[2] = { 0, 17 }
float data_vec[17] = { val0 } // all values are the same
float data_result[1] = { 0 }
И это отрывок из моей функции-обертки:
float *dev_data_mat;
unsigned int *dev_row_ptr;
float *dev_data_vec;
float *dev_data_result;
// Allocate memory on the device
HANDLE_ERROR(cudaMalloc((void**) &dev_data_mat, num_nonzeroes * sizeof(float)));
HANDLE_ERROR(cudaMalloc((void**) &dev_row_ptr, num_row_ptr * sizeof(unsigned int)));
HANDLE_ERROR(cudaMalloc((void**) &dev_data_vec, dim_x * sizeof(float)));
HANDLE_ERROR(cudaMalloc((void**) &dev_data_result, dim_y * sizeof(float)));
// Copy each buffer into the allocated memory
HANDLE_ERROR(cudaMemcpy(
dev_data_mat,
data_mat,
num_nonzeroes * sizeof(float),
cudaMemcpyHostToDevice));
HANDLE_ERROR(cudaMemcpy(
dev_row_ptr,
row_ptr,
num_row_ptr * sizeof(unsigned int),
cudaMemcpyHostToDevice));
HANDLE_ERROR(cudaMemcpy(
dev_data_vec,
data_vec,
dim_x * sizeof(float),
cudaMemcpyHostToDevice));
HANDLE_ERROR(cudaMemcpy(
dev_data_result,
data_result,
dim_y * sizeof(float),
cudaMemcpyHostToDevice));
// Calc grid dimension and block dimension
dim3 grid_dim(dim_y);
dim3 block_dim(BLOCK_SIZE);
// Start kernel
krnlSpMVmul1<<<grid_dim, block_dim, BLOCK_SIZE>>>(
dev_data_mat,
num_nonzeroes,
dev_row_ptr,
dev_data_vec,
dev_data_result);
Надеюсь, это просто, но объясню, если это будет интересно..
Еще одна вещь: я только что понял, что использование BLOCK_SIZE из 128 и 33 ненулевых значений также приводит к сбою ядра.Опять же, только последнее значение не вычисляется.