I want to convert an openMP program to cuda c.
I try to find my way on the web and the sdk. But the material is beyond my level.
My c program loop over n=2^30 index and add the weight of each index.
1) What is the correct grid_size and block_size?
My guess is to replicate openMP and do
grid_size=n/max_number_of_cuda_threads;
block_size=1;
2) How can I implement openMP reduction in cuda?
I try a cudaMemcpy and then reduce the array in standard c, but it seems slow.
I look at the thrust library and its reduce operator. But I don't see how to integrate it with my current code.
program.c
#include <math.h>
#include <omp.h>
float get_weigth_of_index(long index,float* data){
int i;
float v=0;
for(i=0;i<4;i++)
v+=index*data[i];
return v;
}
int main(){
long i;
float r=0;
long n=pow(2,30);
float data[4]={0,1,2,3};
#pragma omp parallel for reduction (+:r)
for(i=0;i<n;i++)
r+=get_weigth_of_index(i,data);
return 0;
}
program.cu
#include <stdlib.h>
#include <stdio.h>
#include <omp.h>
#include <math.h>
__device__ float get_weigth_of_index(long index,float* data){
int i;
float v=0;
for(i=0;i<4;i++)
v+=index*data[i];
return v;
}
__global__ void looper(long max_number_of_cuda_threads, float* data,float* result){
long bid=blockIdx.x;
long start=bid*max_number_of_cuda_threads;
long end=start+max_number_of_cuda_threads;
long i;
float r=0;
for(i=start;i<end;i++)
r+=get_weigth_of_index(i,data);
result[bid]=r;
}
int main(){
long n=pow(2,30);
int max_number_of_cuda_threads=1024; //I'm not sure it's correct
long grid_size=n/max_number_of_cuda_threads;
long block_size=1;
float data_host[4]={0,1,2,3};
float* data_device=0;
float* result_device=0;
cudaMalloc((void**)&data_device, sizeof(int)*4);
cudaMemcpy(data_device, data_host, sizeof(int)*4, cudaMemcpyHostToDevice);
cudaMalloc((void**)&result_device, sizeof(float)*grid_size);
looper<<<grid_size,block_size>>>(max_number_of_cuda_threads,data_device,result_device);
//reduction with standard c: cudaMemcpy seems slow
float* result_host=(float*)malloc(sizeof(float)*grid_size);
cudaMemcpy(result_host, result_device, sizeof(float)*grid_size, cudaMemcpyDeviceToHost);
long i;
float v=0;
#pragma omp parallel for reduction(+:v)
for(i=0;i<grid_size;i++)
v+=result_host[i];
printf("result:%f",v);
return 0;
}
my gpu card
Device 0: "Tesla M2050"
Number of multiprocessors: 14
Number of cores: 448
Total amount of constant memory: 65536 bytes
Total amount of shared memory per block: 49152 bytes
Total number of registers available per block: 32768
Warp size: 32
Maximum number of threads per block: 1024
Maximum sizes of each dimension of a block: 1024 x 1024 x 64
Maximum sizes of each dimension of a grid: 65535 x 65535 x 1
Maximum memory pitch: 2147483647 bytes
Texture alignment: 512 bytes
I think that thrust::transform_reduce can solve your problem. This code shows how you can use it:
#include <thrust/transform_reduce.h>
#include <thrust/functional.h>
#include <thrust/device_vector.h>
#include <thrust/host_vector.h>
#include <cmath>
struct get_weigth_of_index
{
get_weigth_of_index(float* data, size_t n)
{
cudaMalloc((void**)&_data,n * sizeof(float));
cudaMemcpy(_data, data, n * sizeof(float), cudaMemcpyHostToDevice);
_n = n;
}
float* _data;
size_t _n;
__host__ __device__
float operator()(const int& index) const
{
float v=0;
for(size_t i=0; i<_n; i++)
v += index * _data[i];
return v;
}
};
int main(void)
{
float x[4] = {1.0, 2.0, 3.0, 4.0};
size_t len = 1024; // init your value
float * index //init and fill you array here
// transfer to device
thrust::device_vector<float> d_index(index, index + len);
get_weigth_of_index unary_op(x, 4);
thrust::plus<float> binary_op;
float init = 0;
float sum = thrust::transform_reduce(d_x.begin(), d_x.end(), unary_op, init, binary_op);
std::cout << sum<< std::endl;
return 0;
}