数组求和的高速方法（利用cuda的共享内存）-第二部分之程序完善

上一篇提到，那份源码的使用是有限制的。

这次来完善一下。其实就是迭代多次，使得最后一次刚好在一个线程块可以求和。

完善部分：

template<class DType>DType array_sum_gpu(DType *dev_array,const int array_size,DType *dev_result){	//const size_t max_block_size  = 512;//目前有些gpu的线程块最大为512，有些为1024.	const size_t block_size = 512;//线程块的大小。	size_t num_elements = array_size;	size_t num_blocks = (num_elements/block_size) + ((num_elements%block_size) ? 1 : 0);	double *dev_input_array = 0;	double *dev_block_sums = 0;//一个线程块一个和。	while(num_elements > block_size)	{		if(dev_block_sums == 0)//第一次		{			dev_input_array = dev_array;		}		else //除了第一次		{			if(dev_input_array != dev_array)				cudaFree(dev_input_array);			dev_input_array = dev_block_sums;		}		num_blocks = (num_elements/block_size) + ((num_elements%block_size) ? 1 : 0);		//给输出结果分配内存		cudaMalloc((void**)&dev_block_sums, sizeof(double) * (num_blocks ));		// launch one kernel to compute, per-block, a partial sum//把每个线程块的和求出来		block_sum<<<num_blocks,block_size,block_size * sizeof(double)>>>(dev_input_array, dev_block_sums, num_elements);		num_elements = num_blocks;	}	block_sum<<<1,num_elements,num_elements * sizeof(double)>>>(dev_block_sums, dev_result, num_elements);	double result = 0;	cudaMemcpy(&result, dev_result, sizeof(double), cudaMemcpyDeviceToHost);	cudaFree(dev_block_sums);	return result;}

核函数block_sum还是原来的代码。

下面是测试我的代码;

void test_sum2(){	// create array of 256k elements	//const int num_elements = 1<<18;//=512*512=262144	const int num_elements = 1<<20;	// generate random input on the host	std::vector<double> h_input(num_elements);	for(int i = 0; i < h_input.size(); ++i)	{		h_input[i] = 1;//random_num<double>();	}	const double host_result = std::accumulate(h_input.begin(), h_input.end(), 0.0f);	std::cerr << "Host sum: " << host_result << std::endl;	// move input to device memory//分配内存	double *d_input = 0;	cudaMalloc((void**)&d_input, sizeof(double) * num_elements);	cudaMemcpy(d_input, &h_input[0], sizeof(double) * num_elements, cudaMemcpyHostToDevice);	double *dev_result=0;	cudaMalloc((void**)&dev_result,sizeof(double));	double sum = array_sum_gpu(d_input,num_elements,dev_result);	std::cout << "Device sum: " << sum << std::endl;}

其实这个程序还是有点限制的。

请注意第一次求num_blocks.

size_t num_blocks = (num_elements/block_size) + ((num_elements%block_size) ? 1 : 0);

万一第一次求出的num_blocks大于线程块的最大数量，一般是65535，那就不行了。

所以如果数组的元素数量大于1024*65535，那就无法计算了。

解决这中问题的通常方法，是让一个线程串行执行多个相同任务。

由于求解的问题暂时没有超过这个数量级（6-7千万），所以先这样。