我可以/应该在 GPU 上运行此统计应用程序的代码吗？
Posted
技术标签:
【中文标题】我可以/应该在 GPU 上运行此统计应用程序的代码吗？【英文标题】：Can/Should I run this code of a statistical application on a GPU? 【发布时间】：2012-10-29 08:58:04 【问题描述】：我正在开发一个在数组中包含大约 10 到 3000 万个浮点值的统计应用程序。
几种方法在嵌套循环中对数组执行不同但独立的计算，例如：
Dictionary<float, int> noOfNumbers = new Dictionary<float, int>();
for (float x = 0f; x < 100f; x += 0.0001f)
int noOfOccurrences = 0;
foreach (float y in largeFloatingPointArray)
if (x == y)
noOfOccurrences++;
noOfNumbers.Add(x, noOfOccurrences);
当前应用程序是用 C# 编写的，在 Intel CPU 上运行，需要几个小时才能完成。我不了解 GPU 编程概念和 API，所以我的问题是：
是否可以（并且有意义）利用 GPU 来加速此类计算？ 如果是：有谁知道任何教程或有任何示例代码（编程语言无关紧要）？【问题讨论】：
【参考方案1】：更新GPU版本
__global__ void hash (float *largeFloatingPointArray,int largeFloatingPointArraySize, int *dictionary, int size, int num_blocks)
int x = (threadIdx.x + blockIdx.x * blockDim.x); // Each thread of each block will
float y; // compute one (or more) floats
int noOfOccurrences = 0;
int a;
while( x < size ) // While there is work to do each thread will:
dictionary[x] = 0; // Initialize the position in each it will work
noOfOccurrences = 0;
for(int j = 0 ;j < largeFloatingPointArraySize; j ++) // Search for floats
// that are equal
// to it assign float
y = largeFloatingPointArray[j]; // Take a candidate from the floats array
y *= 10000; // e.g if y = 0.0001f;
a = y + 0.5; // a = 1 + 0.5 = 1;
if (a == x) noOfOccurrences++;
dictionary[x] += noOfOccurrences; // Update in the dictionary
// the number of times that the float appears
x += blockDim.x * gridDim.x; // Update the position here the thread will work
我只是针对较小的输入进行了测试，因为我正在笔记本电脑上进行测试。尽管如此，它仍然有效，但需要进行更多测试。
更新顺序版本
我刚刚做了这个简单的版本，它在不到 20 秒（包括生成数据的函数所用的时间）内为一个包含 30,000,000 个元素的数组执行算法。
这个简单的版本首先对你的浮点数组进行排序。然后，将遍历排序后的数组并检查给定的value 在数组中出现的次数，然后将该值连同它出现的次数一起放入字典中。
您可以使用sorted 映射，而不是我使用的unordered_map。
代码如下：
#include <stdio.h>
#include <stdlib.h>
#include "cuda.h"
#include <algorithm>
#include <string>
#include <iostream>
#include <tr1/unordered_map>
typedef std::tr1::unordered_map<float, int> Mymap;
void generator(float *data, long int size)
float LO = 0.0;
float HI = 100.0;
for(long int i = 0; i < size; i++)
data[i] = LO + (float)rand()/((float)RAND_MAX/(HI-LO));
void print_array(float *data, long int size)
for(long int i = 2; i < size; i++)
printf("%f\n",data[i]);
std::tr1::unordered_map<float, int> fill_dict(float *data, int size)
float previous = data[0];
int count = 1;
std::tr1::unordered_map<float, int> dict;
for(long int i = 1; i < size; i++)
if(previous == data[i])
count++;
else
dict.insert(Mymap::value_type(previous,count));
previous = data[i];
count = 1;
dict.insert(Mymap::value_type(previous,count)); // add the last member
return dict;
void printMAP(std::tr1::unordered_map<float, int> dict)
for(std::tr1::unordered_map<float, int>::iterator i = dict.begin(); i != dict.end(); i++)
std::cout << "key(string): " << i->first << ", value(int): " << i->second << std::endl;
int main(int argc, char** argv)
int size = 1000000;
if(argc > 1) size = atoi(argv[1]);
printf("Size = %d",size);
float data[size];
using namespace __gnu_cxx;
std::tr1::unordered_map<float, int> dict;
generator(data,size);
sort(data, data + size);
dict = fill_dict(data,size);
return 0;
如果你的机器上安装了库推力，你应该使用这个：
#include <thrust/sort.h>
thrust::sort(data, data + size);
而不是这个
sort(data, data + size);
肯定会更快。
原帖
我正在开发一个具有大型数组的统计应用程序 包含 10 - 30 百万个浮点值。
是否有可能（并且有意义）利用 GPU 来加速 这样的计算？
是的，是的。一个月前，我在 GPU 上运行了一个完整的分子动力学模拟。其中一个计算粒子对之间的力的内核作为参数6 数组接收，每个数组带有500,000 双精度数，总计3 百万双精度数(22 MB)。
所以如果你打算放30百万个浮点数，大约是全局内存的114 MB，那是没有问题的。
在您的情况下，计算次数会成为问题吗？根据我在分子动力学 (MD) 方面的经验，我会说不。顺序 MD 版本大约需要 25 小时才能完成，而 GPU 版本需要 45 分钟。您说您的应用程序花费了几个小时，同样基于您的代码示例，它看起来比 MD 更软。
以下是力计算示例：
__global__ void add(double *fx, double *fy, double *fz,
double *x, double *y, double *z,...)
int pos = (threadIdx.x + blockIdx.x * blockDim.x);
...
while(pos < particles)
for (i = 0; i < particles; i++)
if(//inside of the same radius)
// calculate force
pos += blockDim.x * gridDim.x;
一个简单的 CUDA 代码示例可以是两个二维数组的总和：
在 C 中：
for(int i = 0; i < N; i++)
c[i] = a[i] + b[i];
在 CUDA 中：
__global__ add(int *c, int *a, int*b, int N)
int pos = (threadIdx.x + blockIdx.x)
for(; i < N; pos +=blockDim.x)
c[pos] = a[pos] + b[pos];
在 CUDA 中，您基本上将每次 for 迭代并分配给每个线程，
1) threadIdx.x + blockIdx.x*blockDim.x;
每个块都有一个ID，从0 到N-1（N 是块的最大数量），每个块有一个'X' 线程数，从ID 到X-1 .
-
为您提供每个线程将根据其
ID 和线程所在的块ID 计算的for 循环迭代； blockDim.x 是块拥有的线程数。
因此，如果您有 2 个块，每个块带有 10 线程和 N=40，则：
Thread 0 Block 0 will execute pos 0
Thread 1 Block 0 will execute pos 1
...
Thread 9 Block 0 will execute pos 9
Thread 0 Block 1 will execute pos 10
....
Thread 9 Block 1 will execute pos 19
Thread 0 Block 0 will execute pos 20
...
Thread 0 Block 1 will execute pos 30
Thread 9 Block 1 will execute pos 39
查看您当前的代码，我已经制作了您的代码在 CUDA 中的样子的草稿：
__global__ hash (float *largeFloatingPointArray, int *dictionary)
// You can turn the dictionary in one array of int
// here each position will represent the float
// Since x = 0f; x < 100f; x += 0.0001f
// you can associate each x to different position
// in the dictionary:
// pos 0 have the same meaning as 0f;
// pos 1 means float 0.0001f
// pos 2 means float 0.0002f ect.
// Then you use the int of each position
// to count how many times that "float" had appeared
int x = blockIdx.x; // Each block will take a different x to work
float y;
while( x < 1000000) // x < 100f (for incremental step of 0.0001f)
int noOfOccurrences = 0;
float z = converting_int_to_float(x); // This function will convert the x to the
// float like you use (x / 0.0001)
// each thread of each block
// will takes the y from the array of largeFloatingPointArray
for(j = threadIdx.x; j < largeFloatingPointArraySize; j += blockDim.x)
y = largeFloatingPointArray[j];
if (z == y)
noOfOccurrences++;
if(threadIdx.x == 0) // Thread master will update the values
atomicAdd(&dictionary[x], noOfOccurrences);
__syncthreads();
您必须使用atomicAdd，因为来自不同块的不同线程可能同时写入/读取noOfOccurrences，因此您必须确保mutual exclusion。
这只是一种方法；您甚至可以将外部循环的迭代分配给线程而不是块。
教程
Rob Farmer 撰写的 Dobbs 博士期刊系列 CUDA: Supercomputing for the masses 非常出色，在其十四期中几乎涵盖了所有内容。它的启动也相当温和，因此对初学者非常友好。
还有其他人：
Volume I: Introduction to CUDA Programming Getting started with CUDA CUDA Resources List看看最后一条，你会发现很多学习CUDA的链接。
OpenCL：OpenCL Tutorials | MacResearch
【讨论】：
【参考方案2】：我对并行处理或 GPGPU 知之甚少，但对于这个特定示例，您可以通过单次遍历输入数组而不是循环遍历一百万次来节省大量时间。对于大型数据集，如果可能，您通常希望一次性完成。即使您正在执行多个独立的计算，如果它在相同的数据集上，您可能会在同一遍中获得更快的速度，因为您将获得更好的参考位置。但由于代码复杂性的增加，这可能不值得。
此外，您真的不想像那样重复地将少量添加到浮点数，舍入误差会加起来，您将无法得到您想要的结果。我在下面的示例中添加了一个 if 语句来检查输入是否与您的迭代模式匹配，但如果您实际上不需要它，请忽略它。
我不懂任何 C#，但您的示例的单遍实现看起来像这样：
Dictionary<float, int> noOfNumbers = new Dictionary<float, int>();
foreach (float x in largeFloatingPointArray)
if (math.Truncate(x/0.0001f)*0.0001f == x)
if (noOfNumbers.ContainsKey(x))
noOfNumbers.Add(x, noOfNumbers[x]+1);
else
noOfNumbers.Add(x, 1);
希望这会有所帮助。
【讨论】：
您可以通过使用 TryGet 而不是 ContainsKey 和 noOfNumbers[x] 来改进您的代码。使用 TryGet 为您节省了一次字典查找，它是 O(1) 摊销的（即并不总是 O(1)）并且是一个昂贵的 O(1)，因为字典是一种相当复杂的数据类型。无论如何+1 感谢两位的帮助。非常感谢，您的建议将很快添加到我的应用程序中。不幸的是，我有近 100 种其他方法，我认为这些方法无法进一步优化。即使我使用代码优化将此类计算速度提高了 90%，在快速 CPU 上仍可能需要几个小时才能完成。 请将有限数据集（以及您自己的基准）的完整方法发送给我们。这将使我们有能力为您提供更多帮助。根据我目前在您的代码中看到的情况，我很确定在我们开始使用 GPU 之前，我可以将代码速度提高一倍。【参考方案3】：绝对YES，这种算法通常是大规模数据并行处理的理想候选者，GPU 非常擅长这一点。是否有可能（并且有意义）利用 GPU 来加速 这样的计算？
如果是：有谁知道任何教程或有任何示例代码 （编程语言无所谓）？
如果您想采用 GPGPU 方式，您有两种选择：CUDA 或 OpenCL。
CUDA 已经成熟，拥有很多工具，但以 NVidia GPU 为中心。
OpenCL 是在 NVidia 和 AMD GPU 以及 CPU 上运行的标准。所以你真的应该喜欢它。
对于教程，您有 Rob Farber 编写的关于 CodeProject 的精彩系列：http://www.codeproject.com/Articles/Rob-Farber#Articles
对于您的特定用例，有很多使用 OpenCL 构建的直方图样本（请注意，许多是图像直方图，但原理相同）。
在使用 C# 时，您可以使用 OpenCL.Net 或 Cloo 等绑定。
如果您的数组太大而无法存储在 GPU 内存中，您可以对其进行块分区并轻松地为每个部分重新运行 OpenCL 内核。
【讨论】：
关于高效直方图算法的附加资源...users.cecs.anu.edu.au/~ramtin/cuda.htm 感谢您的帮助！非常感激。您对 DirectX 有何看法？似乎有一个很好的 C# SDK www.sharpdx.org 做了一些额外的研究。 OpenCL 非常有趣，因为它还支持 Xeon Phi 和现代 Intel CPU 的集成 GPU，请参阅此处bit.ly/Ta29ab @Mike：至于 DirectCompute，您可以忘记它并使用更高级的 API，例如 C++ AMP，但前提是您的计算必须并且总是在 Windows 上运行。否则，请使用像 OpenCL 这样的标准 API，因为如果您需要在 Linux 农场上运行它，它将使您的代码面向未来。【参考方案4】：除了上述海报的建议之外，在适当的时候使用 TPL（任务并行库）在多个内核上并行运行。
上面的示例可以使用 Parallel.Foreach 和 ConcurrentDictionary，但更复杂的 map-reduce 设置将数组拆分为块，每个块生成一个字典，然后将其缩减为单个字典，从而获得更好的结果。
我不知道您的所有计算是否都正确映射到 GPU 功能，但无论如何您都必须使用 map-reduce 算法将计算映射到 GPU 核心，然后将部分结果缩减为单个结果，因此您最好先在 CPU 上执行此操作，然后再转到不太熟悉的平台。
【讨论】：
感谢您的建议。我已经在使用 TPL，但在更高级别上。这意味着我的应用程序并行调用了几个似乎运行良好的方法（CPU 使用率超过 90%）。【参考方案5】：鉴于此，我不确定使用 GPU 是否合适 'largerFloatingPointArray' 值需要从内存中检索。我的理解是 GPU 更适合自包含计算。
我认为将这个单进程应用程序转变为在许多系统上运行的分布式应用程序并调整算法应该会大大加快速度，具体取决于可用系统的数量。
您可以使用经典的“分而治之”的方法。我会采取的一般方法如下。
使用一个系统将“largeFloatingPointArray”预处理为哈希表或数据库。这将一次性完成。它将使用浮点值作为键，并将数组中出现的次数作为值。最坏的情况是每个值只出现一次，但这不太可能。如果每次运行应用程序时 largeFloatingPointArray 都在不断变化，那么内存中的哈希表是有意义的。如果它是静态的，那么该表可以保存在诸如 Berkeley DB 之类的键值数据库中。我们称之为“查找”系统。
在另一个系统上，我们称其为“主”，创建工作块并将工作项“分散”到 N 个系统中，并在结果可用时“收集”结果。例如，一个工作项可以像两个数字一样简单，表示系统应该工作的范围。当系统完成工作时，它会发回事件数组，并准备好处理另一块工作。
性能得到了提高，因为我们不再对 largeFloatingPointArray 进行迭代。如果查找系统成为瓶颈，那么它可以根据需要在尽可能多的系统上进行复制。
如果有足够多的系统并行工作，应该可以将处理时间减少到几分钟。
我正在为基于多核的系统（通常称为微服务器）使用 C 语言并行编程开发编译器，这些系统是/或将使用系统内的多个“片上系统”模块构建的. ARM 模块供应商包括 Calxeda、AMD、AMCC 等。英特尔可能也会提供类似的产品。
我有一个可以工作的编译器版本，它可以用于这样的应用程序。该编译器基于 C 函数原型，生成 C 网络代码，实现跨系统的进程间通信代码 (IPC)。可用的 IPC 机制之一是 socket/tcp/ip。
如果您在实施分布式解决方案方面需要帮助，我很乐意与您讨论。
于 2012 年 11 月 16 日添加。
我对算法进行了更多思考，我认为这应该一次性完成。它是用 C 编写的，与现有的相比应该非常快。
/*
* Convert the X range from 0f to 100f in steps of 0.0001f
* into a range of integers 0 to 1 + (100 * 10000) to use as an
* index into an array.
*/
#define X_MAX (1 + (100 * 10000))
/*
* Number of floats in largeFloatingPointArray needs to be defined
* below to be whatever your value is.
*/
#define LARGE_ARRAY_MAX (1000)
main()
int j, y, *noOfOccurances;
float *largeFloatingPointArray;
/*
* Allocate memory for largeFloatingPointArray and populate it.
*/
largeFloatingPointArray = (float *)malloc(LARGE_ARRAY_MAX * sizeof(float));
if (largeFloatingPointArray == 0)
printf("out of memory\n");
exit(1);
/*
* Allocate memory to hold noOfOccurances. The index/10000 is the
* the floating point number. The contents is the count.
*
* E.g. noOfOccurances[12345] = 20, means 1.2345f occurs 20 times
* in largeFloatingPointArray.
*/
noOfOccurances = (int *)calloc(X_MAX, sizeof(int));
if (noOfOccurances == 0)
printf("out of memory\n");
exit(1);
for (j = 0; j < LARGE_ARRAY_MAX; j++)
y = (int)(largeFloatingPointArray[j] * 10000);
if (y >= 0 && y <= X_MAX)
noOfOccurances[y]++;
【讨论】：
工作可以在第二次机器网络中拆分；但是恕我直言，使用 GPU 功能进行廉价（通常是巨大的）改进要好得多。至于您的框架，它与 MPI 相比如何？ :) 感谢所有信息和 c 代码。也许我为我的问题找到了一个很好的解决方案：bit.ly/Ta4aSL [PDF] 听起来很有希望......你觉得呢？ Mike，这是一种利用 DirectX 而又不受特定 GPU 限制的有趣方式。我在考虑副作用，如果有的话。在积极使用 DirectX 的同时，对将图形渲染到显示器的其他应用程序是否有任何影响？尝试在您的应用程序运行和不运行的情况下播放 youtube 或 windows 媒体播放器视频，看看您是否注意到正在播放的视频质量有任何恶化。此外，您是否知道将来您是否需要扩展超出工作站的能力？由于它都是 Windows 环境的一部分，我认为值得一试。以上是关于我可以/应该在 GPU 上运行此统计应用程序的代码吗？的主要内容，如果未能解决你的问题，请参考以下文章
GPU cuda代码可以在多个GPU卡上运行而无需任何实现吗？