前の記事(Boost LambdaとCUDA/GPUについての雑感)が結構読まれているので、続編を書くことにした。
ラムダ式は無名関数を提供する機能。
→ 個人的には、「いかにforループ式を書かないか」というコンセプトの１つ。
思えば、CUDA/GPU/Thrustも、いかに「C/C++的なforループを書かないか」という設計になっているとも言える。
であるが、ラムダ式（Boost.Lambda)もCUDA/GPUも、割と同じ感覚（イメージ）で書けるということでもある。
ラムダ式のコード
1#include <boost/lambda/lambda.hpp>
2#include <boost/lambda/if.hpp>
3#include <vector>
4#include <algorithm>
5#include <iostream>
6
7int main()
8{
9 std::vector<int> v{1, 3, 2, 0, 4, 5, 2};
10 std::for_each(v.begin(), v.end(),
11 boost::lambda::if_then(boost::lambda::_1 > 2,
12 std::cout << boost::lambda::_1 << "\n"));
13}
これは、リストv{1, 3, 2, 0, 4, 5, 2}に対し、2よりも大きい値を表示する。
実行してみる。。
# ./a.out
3
4
5
個人的には、リストの各要素にif_thenという関数を張り付ける（割り当てる）イメージである。
実は、CUDA/GPUも書いているときの感覚は同じである。
少し長いがコードを見てみる。。
1#include
2#include
3
4void validateResult(float *hostRef, float *gpuRef, const int N)
5{
6 double epsilon = 1.0E-8;
7 bool match = 1;
8
9 for (int i = 0; i < N; i++)
10 {
11 if (abs(hostRef[i] - gpuRef[i]) > epsilon)
12 {
13 match = 0;
14 printf("Arrays do not match!\n");
15 printf("host %5.2f gpu %5.2f at current %d\n", hostRef[i], gpuRef[i], i);
16 break;
17 }
18 }
19
20 if (match) printf("Arrays match.\n\n");
21
22 return;
23}
24
25void genData(float *ip, int size)
26{
27 time_t t;
28 srand((unsigned) time(&t));
29
30 for (int i = 0; i < size; i++)
31 {
32 ip[i] = (float)( rand() & 0xFF ) / 10.0f;
33 }
34
35 return;
36}
37
38void sumArraysOnHost(float *A, float *B, float *C, const int N)
39{
40 for (int idx = 0; idx < N; idx++)
41 {
42 C[idx] = A[idx] + B[idx];
43 }
44}
45__global__ void sumArraysOnGPU(float *A, float *B, float *C, const int N)
46{
47 int i = blockIdx.x * blockDim.x + threadIdx.x;
48
49 if (i < N) C[i] = A[i] + B[i];
50}
51
52int main(int argc, char **argv)
53{
54 struct timespec startTime, endTime, sleepTime;
55
56 printf("%s Starting...\n", argv[0]);
57
58 int dev = 0;
59 cudaDeviceProp deviceProp;
60 cudaGetDeviceProperties(&deviceProp, dev);
61 printf("Using Device %d: %s\n", dev, deviceProp.name);
62 cudaSetDevice(dev);
63
64 int nElem = 1 << 28;
65 printf("Vector size %d\n", nElem);
66
67 size_t nBytes = nElem * sizeof(float);
68
69 float *host_A, *host_B, *hostRef, *gpuRef;
70 host_A = (float *)malloc(nBytes);
71 host_B = (float *)malloc(nBytes);
72 hostRef = (float *)malloc(nBytes);
73 gpuRef = (float *)malloc(nBytes);
74
75 genData(host_A, nElem);
76 genData(host_B, nElem);
77
78 memset(hostRef, 0, nBytes);
79 memset(gpuRef, 0, nBytes);
80
81 clock_gettime(CLOCK_REALTIME, &startTime);
82 sleepTime.tv_sec = 0;
83 sleepTime.tv_nsec = 123; 84
85 sumArraysOnCPU(host_A, host_B, hostRef, nElem);
86
87 clock_gettime(CLOCK_REALTIME, &endTime);
88 if (endTime.tv_nsec < startTime.tv_nsec) {
89 printf("%ld.%09ld", endTime.tv_sec - startTime.tv_sec - 1
90 ,endTime.tv_nsec + 1000000000 - startTime.tv_nsec);
91 } else {
92 printf("%ld.%09ld", endTime.tv_sec - startTime.tv_sec
93 ,endTime.tv_nsec - startTime.tv_nsec);
94 }
95 printf(" sec\n");
96
97 float *device_A, *device_B, *device_C;
98 cudaMalloc((float**)&device_A, nBytes);
99 cudaMalloc((float**)&device_B, nBytes);
100 cudaMalloc((float**)&device_C, nBytes);
101
102 cudaMemcpy(device_A, host_A, nBytes, cudaMemcpyHostToDevice);
103 cudaMemcpy(device_B, host_B, nBytes, cudaMemcpyHostToDevice);
104 cudaMemcpy(device_C, gpuRef, nBytes, cudaMemcpyHostToDevice); 105
106 int iLen = 1024;
107 dim3 block (iLen);
108 dim3 grid ((nElem + block.x - 1) / block.x);
109
110 clock_gettime(CLOCK_REALTIME, &startTime);
111 sleepTime.tv_sec = 0;
112 sleepTime.tv_nsec = 123;
113
114 sumArraysOnGPU<<>>(device_A, device_B, device_C, nElem);
115 cudaDeviceSynchronize();
116 cudaGetLastError() ;
117
118 clock_gettime(CLOCK_REALTIME, &endTime);
119 if (endTime.tv_nsec < startTime.tv_nsec) {
120 printf("%ld.%09ld", endTime.tv_sec - startTime.tv_sec - 1
121 ,endTime.tv_nsec + 1000000000 - startTime.tv_nsec);
122 } else {
123 printf("%ld.%09ld", endTime.tv_sec - startTime.tv_sec
124 ,endTime.tv_nsec - startTime.tv_nsec);
125 }
126 printf(" sec\n");
127
128 cudaMemcpy(gpuRef, device_C, nBytes, cudaMemcpyDeviceToHost);
129
130 validateResult(hostRef, gpuRef, nElem);
131
132 cudaFree(device_A);
133 cudaFree(device_B);
134 cudaFree(device_C);
135
136 free(host_A);
137 free(host_B);
138 free(hostRef);
139 free(gpuRef);
140
141 return(0);
142}
ここで、この２つの関数であるが、なぜGPUの方の関数はforループがないかというと、上のラムダ式と似たような仕組みになっているからである。
38void sumArraysOnCPU(float *A, float *B, float *C, const int N)
39{
40 for (int idx = 0; idx < N; idx++)
41 {
42 C[idx] = A[idx] + B[idx];
43 }
44}
45__global__ void sumArraysOnGPU(float *A, float *B, float *C, const int N)
46{
47 int i = blockIdx.x * blockDim.x + threadIdx.x;
48
49 if (i < N) C[i] = A[i] + B[i];
50}
図にするとこんな感じであらう。
つまり、GPUデバイス内の各スレッドに対して、sumArrayOnGPUを張り付ける（割り当てる）イメージである。
実行してみる。。
# ./a.out
./a.out Starting...
Using Device 0: Quadro GV100
Vector size 268435456
1.719520959 sec
0.005400906 sec
（｀ー´）b
実は、Google発のMapRuduceも同じ感覚で書ける。。
アレは、各データの要素にMap関数を張り付けていくイメージである。
このあたりのことは、LISPをやりまくると開眼するらしい。。
(ﾟДﾟ;)