今天刷知乎看到個挺有意思的問題：
我想這有啥難的，還能寫出花來不成？結果看到高贊回答，感覺自己的智商有點不夠用了。
随便來看一個高贊回答是怎麼寫的：
這個其實還算比較簡單的，沒啥難度，還有更晦澀的：
這個乍一看根本看不懂在寫啥，當然平時也很少會寫這種晦澀的代碼。
CUDA花式整活！今天我就教大家用cuda來計算一下20的階乘，就當作是CUDA的一個入門例子。
先來看看我的回答：
如何優雅地利用c 編程從1乘到20？13 贊同 · 2 評論回答
我提供了兩種CUDA的實現方法。
方法1
#include <iostream>
typedef unsigned long long int ull;
const int N = 20;
__device__ ull atomicMul(ull* address, ull val) {
ull *address_as_ull = (ull *)address;
ull old = *address_as_ull, assumed;
do {
assumed = old;
old = atomicCAS(address_as_ull, assumed, val * assumed);
} while (assumed != old);
return old;
}
__global__ void mul(ull *x) {
ull i = threadIdx.x 1;
atomicMul(x, i);
}
int main() {
ull *x;
cudaMallocManaged(&x, sizeof(ull));
x[0] = 1;
mul<<<1, N>>>(x);
cudaDeviceSynchronize();
std::cout << x[0] << std::endl;
cudaFree(x);
return 0;
}
這個方法使用原子操作，一共20個線程，每個線程負責一個乘數，然後一起乘到x[0]上。
但是由于并行執行，線程之間沒有先後順序，會導緻同時乘的時候産生沖突，所以需要使用原子操作。在某一個線程将它的乘數乘到x[0]上時，不會被其他線程打斷。也就是會加鎖，同一時刻隻會有一個線程在進行乘法操作。
但是由于CUDA隻提供了加法和減法的原子操作（atomicAdd和atomicSub），所以得自己實現乘法的原子操作atomMul，利用的是atomicCAS操作，也就是compare and swap ，如果目标地址元素和待比較的元素相同，就進行元素的交換，否則不進行任何操作。
可以看出，在atomicMul函數的do while循環中，先用old變量保存x[0]處的當前值，這時候如果有其他線程在x[0]處寫入了新的值，那麼接下來該線程的atomicCAS操作就會判斷元素不相同，不進行任何操作，重新執行下一輪循環。
方法2
#include <iostream>
typedef unsigned long long int ull;
const int N = 20;
const int WARP_SIZE = 32;
__global__ void mul(ull *x) {
int i = threadIdx.x;
ull val = x[i];
for (int mask = WARP_SIZE / 2; mask > 0; mask >>= 1)
val *= __shfl_xor_sync(WARP_SIZE - 1, val, mask, WARP_SIZE);
x[i] = val;
}
int main() {
ull *x;
cudaMallocManaged(&x, WARP_SIZE * sizeof(ull));
for (int i = 0; i < WARP_SIZE; i)
x[i] = i < N ? i 1 : 1;
mul<<<1, WARP_SIZE>>>(x);
cudaDeviceSynchronize();
std::cout << x[0] << std::endl;
cudaFree(x);
return 0;
}
這種方法使用線程束原語__shfl_xor_sync，隻要線程在同一個線程束中（32個線程），就可以獲取其他線程的值，異或運算後寫入指定地址。詳細原理這裡就不解釋了，可以簡單理解為：
所以隻需要在開始時，分配一個大小為32的數組，前20個元素分别保存1-20，後面12個元素是為了滿足線程束大小32的條件，賦值為1就行了。
方法2改進方法2需要額外開辟大小為32的數組，空間存在浪費，并且數組賦值也需要時間。
感謝@NekoDaemon老哥提供的優化建議，隻需要在計算的時候根據線程号計算對應乘積元素就行，但是線程數仍然需要分配32個。
#include <iostream>
typedef unsigned long long int ull;
const int N = 20;
const int WARP_SIZE = 32;
__global__ void mul(ull *x) {
int i = threadIdx.x;
ull val = i < N ? i 1 : 1;
for (int mask = WARP_SIZE / 2; mask > 0; mask >>= 1)
val *= __shfl_xor_sync(WARP_SIZE - 1, val, mask, WARP_SIZE);
if (!i) x[i] = val;
}
int main() {
ull *x;
cudaMallocManaged(&x, sizeof(ull));
mul<<<1, WARP_SIZE>>>(x);
cudaDeviceSynchronize();
std::cout << x[0] << std::endl;
cudaFree(x);
return 0;
}
代碼保存為run.cu，然後執行nvcc run.cu -o run，最後執行./run，就會出來結果2432902008176640000。
今天沒有講解CUDA編程的基礎概念，适合有一定基礎的同學閱讀，如果有對CUDA感興趣的同學，可以在評論區留言，下次專門寫一個CUDA入門系列教程。關注我的公衆号【算法碼上來】，每天跟進最新文章。
,更多精彩资讯请关注tft每日頭條，我们将持续为您更新最新资讯!
zpostcode 
Recruit 
weather 
mreligion 
Yellowpages 
sport 
constellation 
shopping 
name 
game 
directory 
literature 
Word 
tour 
furnish 
Lottery 
tftnews 
lyrics 
News 
digital 
car 
dir 
Edu 
Finance 
