前面介紹的計算平方和的程式，似乎沒有什麼實用價值。所以我們的第二個 CUDA 程式，要做一個確實有（某些）實用價值的程式，也就是進行矩陣乘法。而且，這次我們會使用浮點數。
雖然矩陣乘法有點老套，不過因為它相當簡單，而且也可以用來介紹一些有關 CUDA 的有趣性質。
矩陣乘法
為了單純起見，我們這裡以方形的矩陣為例子。基本上，假設有兩個矩陣 A 和 B，則計算 AB = C 的方法如下：
for(i = 0; i < n; i++) { for(j = 0; j < n; j++) { C[i][j] = 0; for(k = 0; k < n; k++) { C[i][j] += A[i][k] * B[k][j]; } } }
一開始，我們先準備好產生資料、設定 CUDA 等等的工作：
int main() { float *a, *b, *c, *d; int n = 1000; if(!InitCUDA()) return 0; a = (float*) malloc(sizeof(float) * n * n); b = (float*) malloc(sizeof(float) * n * n); c = (float*) malloc(sizeof(float) * n * n); d = (float*) malloc(sizeof(float) * n * n); srand(0); matgen(a, n, n); matgen(b, n, n); clock_t time = matmultCUDA(a, n, b, n, c, n, n); matmult(a, n, b, n, d, n, n); compare_mat(c, n, d, n, n); double sec = (double) time / CLOCKS_PER_SEC; printf("Time used: %.2f (%.2lf GFLOPS)\n", sec, 2.0 * n * n * n / (sec * 1E9)); return 0; }
InitCUDA 函式和第一個 CUDA 程式一樣，可以直接參考前面的文章。以下是上面用到的一些其它的函式：
產生矩陣：
void matgen(float* a, int lda, int n) { int i, j; for(i = 0; i < n; i++) { for(j = 0; j < n; j++) { a[i * lda + j] = (float) rand() / RAND_MAX + (float) rand() / (RAND_MAX * RAND_MAX); } } }
這個函式只是利用亂數產生器把矩陣填滿 0 ~ 1 之間的數字。特別注意到因為 C 語言中無法宣告變動大小的二維矩陣，所以我們使用 i * lda + j 的方式。
進行矩陣乘法：
void matmult(const float* a, int lda, const float* b, int ldb, float* c, int ldc, int n) { int i, j, k; for(i = 0; i < n; i++) { for(j = 0; j < n; j++) { double t = 0; for(k = 0; k < n; k++) { t += a[i * lda + k] * b[k * ldb + j]; } c[i * ldc + j] = t; } } }
這是以 CPU 進行矩陣乘法、用來進行驗證答案正確與否的程式。特別注意到它用 double 來儲存暫時的計算結果，以提高精確度。
驗證結果：
void compare_mat(const float* a, int lda, const float* b, int ldb, int n) { float max_err = 0; float average_err = 0; int i, j; for(i = 0; i < n; i++) { for(j = 0; j < n; j++) { if(b[i * ldb + j] != 0) { float err = fabs((a[i * lda + j] - b[i * ldb + j]) / b[i * ldb + j]); if(max_err < err) max_err = err; average_err += err; } } } printf("Max error: %g Average error: %g\n", max_err, average_err / (n * n)); }
這個函式計算兩個矩陣的最大相對誤差和平均相對誤差，並把結果印出來。
最後是 CUDA 的矩陣乘法的部份：
#define NUM_THREADS 256 clock_t matmultCUDA(const float* a, int lda, const float* b, int ldb, float* c, int ldc, int n) { float *ac, *bc, *cc; clock_t start, end; start = clock(); cudaMalloc((void**) &ac, sizeof(float) * n * n); cudaMalloc((void**) &bc, sizeof(float) * n * n); cudaMalloc((void**) &cc, sizeof(float) * n * n); cudaMemcpy2D(ac, sizeof(float) * n, a, sizeof(float) * lda, sizeof(float) * n, n, cudaMemcpyHostToDevice); cudaMemcpy2D(bc, sizeof(float) * n, b, sizeof(float) * ldb, sizeof(float) * n, n, cudaMemcpyHostToDevice); int blocks = (n + NUM_THREADS - 1) / NUM_THREADS; matMultCUDA<<<blocks * n, NUM_THREADS>>>(ac, n, bc, n, cc, n, n); cudaMemcpy2D(c, sizeof(float) * ldc, cc, sizeof(float) * n, sizeof(float) * n, n, cudaMemcpyDeviceToHost); cudaFree(ac); cudaFree(bc); cudaFree(cc); end = clock(); return end - start; }
這個函式相當單純，就是在顯示記憶體中配置存放矩陣的記憶體，然後把主記憶體中的矩陣資料複製到顯示記憶體上。不過，因為我們的矩陣乘法函式可以指定 pitch（即 lda、ldb、和 ldc），所以如果用一般的 cudaMemcpy 函式來複製記憶體的話，會需要每個 row 都分開複製，那會需要呼叫很多次 cudaMemcpy 函式，會使效率變得很差。因此，在這裡我們用了一個新的 cudaMemcpy2D 函式，它是用來複製二維陣列，可以指定陣列的 pitch。這樣就可以透過一次函式呼叫就可以了。
進行計算的 kernel 如下：
__global__ static void matMultCUDA(const float* a, size_t lda, const float* b, size_t ldb, float* c, size_t ldc, int n) { const int tid = threadIdx.x; const int bid = blockIdx.x; const int idx = bid * blockDim.x + tid; const int row = idx / n; const int column = idx % n; int i; if(row < n && column < n) { float t = 0; for(i = 0; i < n; i++) { t += a[row * lda + i] * b[i * ldb + column]; } c[row * ldc + column] = t; } }
這個函式一開始先從 bid 和 tid 計算出這個 thread 應該計算的 row 和 column，在判斷 row 和 column 在範圍內之後，就直接進行計算，並把結果寫到 c 矩陣中，是非常單純的函式。
在 GeForce 8800GT 上實際執行的結果如下：
Max error: 2.01484e-006 Average error: 3.36637e-007 Time used: 1.1560 (1.73 GFLOPS)
可以看到兩個問題：
- 很明顯的，執行效率相當低落。
- 最大相對誤差偏高。理想上應該要低於 1e-6。
計算結果的誤差偏高的原因是，在 CPU 上進行計算時，我們使用 double（即 64 bits 浮點數）來累進計算過程，而在 GPU 上則只能用 float（32 bits 浮點數）。在累加大量數字的時候，由於累加結果很快會變大，因此後面的數字很容易被捨去過多的位數。
由於 CUDA 的浮點數運算，在進行加、減、乘法時是符合 IEEE 754 規定的精確度的，因此，我們可以利用 Kahan’s Summation Formula 來提高精確度。把程式改成：
if(row < n && column < n) { float t = 0; float y = 0; for(i = 0; i < n; i++) { float r; y -= a[row * lda + i] * b[i * ldb + column]; r = t - y; y = (r - t) + y; t = r; } }
修改後的程式的執行結果是：
Max error: 1.19209e-007 Average error: 4.22751e-008 Time used: 1.1560 (1.73 GFLOPS)
可以看到相對誤差有很大的改善，效率則沒什麼變化。
由於 Kahan’s Summation Formula 需要的運算量提高，但是效率卻沒有什麼改變，可以看出這個 kernel 主要的瓶頸應該是在記憶體的存取動作上。這是因為有大量的記憶體讀取是重覆的。例如，矩陣 a 的一個 row 在每次進行計算時都被重覆讀入，但這是相當浪費的。這樣的計算方式，總共需要讀取 2*n3 次記憶體。如果讓一個 row 只需要讀入一次的話，就可以減到為 n3+n2 次。
第一個改良
和我們的第一個 CUDA 程式一樣，我們可以利用 shared memory 來儲存每個 row 的資料。不過，因為只有同一個 block 的 threads 可以共用 shared memory，因此現在一個 row 只能由同一個 block 的 threads 來進行計算。另外我們也需要能存放一整個 row 的 shared memory。因此，把先把呼叫 kernel 的部份改成：
matMultCUDA<<<n, NUM_THREADS, sizeof(float) * n>>>(ac, n, bc, n, cc, n, n);
kernel 的部份則改成：
__global__ static void matMultCUDA(const float* a, size_t lda, const float* b, size_t ldb, float* c, size_t ldc, int n) { extern __shared__ float data[]; const int tid = threadIdx.x; const int row = blockIdx.x; int i, j; for(i = tid; i < n; i += blockDim.x) { data[i] = a[row * lda + i]; } __syncthreads(); for(j = tid; j < n; j += blockDim.x) { float t = 0; float y = 0; for(i = 0; i < n; i++) { float r; y -= data[i] * b[i * ldb + j]; r = t - y; y = (r - t) + y; t = r; } c[row * ldc + j] = t; } }
第一個部份先把整個 row 讀到 shared memory 中，而第二個部份則進行計算，並沒有太大的變化。主要的差別是現在一個 row 只由一個 block 進行計算。
在 GeForce 8800GT 上，執行的結果是：
Max error: 1.19209e-007 Average error: 4.22751e-008 Time used: 0.4220 (4.74 GFLOPS)
很明顯的，計算的結果並沒有改變，不過速度則提高了超過一倍。雖然如此，但是這樣的效率仍不盡理想，因為理論上 GeForce 8800GT 有超過 300GFLOPS 的運算性能。即使是把 Kahan’s Summation Formula 所需要的額外運算考慮進去，這樣的效率仍然連理論最大值的十分之一都不到。
會有這樣的結果，原因其實還是同樣的：對記憶體的存取次數太多了。雖然現在 A 矩陣的 row 的資料已經不再需要重覆讀取，但是 B 矩陣的 column 的資料仍然一直被重覆讀取。
另一個問題比較不是那麼明顯：對 B 矩陣的讀取，雖然看起來不連續，但實際上它是連續的。這是因為不同的 thread 會讀取不同的 column，因此同時間每個 thread 讀取的各個 column 加起來，就是一個連續的記憶體區塊。那麼，為什麼效率還是不佳呢？這是因為，GPU 上的記憶體控制器，從某個固定的倍數位址開始讀取，才會有最高的效率（例如 16 bytes 的倍數）。由於矩陣大小並不是 16 的倍數（這裡使用的是 1000×1000 的矩陣），所以造成效率不佳的情形。
要解決這個問題，我們可以在 cudaMalloc 的時候稍微修改一下，讓寬度變成 適當的倍數就可以了。但是，適當的倍數是多少呢？幸運的是，我們並不需要知道這些細節。CUDA 提供了一個 cudaMallocPitch 的函式，可以自動以最佳的倍數來配置記憶體。因此，我們可以把 cudaMalloc 的部份改成：
size_t pitch_a, pitch_b, pitch_c; cudaMallocPitch((void**) &ac, &pitch_a, sizeof(float) * n, n); cudaMallocPitch((void**) &bc, &pitch_b, sizeof(float) * n, n); cudaMallocPitch((void**) &cc, &pitch_c, sizeof(float) * n, n);
cudaMallocPitch 函式會以適當的倍數配置記憶體，並把配置的寬度傳回。因此，在把矩陣複製到顯示記憶體上時，要使用它傳回的寬度：
cudaMemcpy2D(ac, pitch_a, a, sizeof(float) * lda, sizeof(float) * n, n, cudaMemcpyHostToDevice); cudaMemcpy2D(bc, pitch_b, b, sizeof(float) * ldb, sizeof(float) * n, n, cudaMemcpyHostToDevice);
呼叫 kernel 的部份也需要修改：
matMultCUDA<<<n, NUM_THREADS, sizeof(float) * n>>>( ac, pitch_a / sizeof(float), bc, pitch_b / sizeof(float), cc, pitch_c / sizeof(float), n);
同樣的，把計算結果複製回到主記憶體時，也要使用傳回的寬度值：
cudaMemcpy2D(c, sizeof(float) * ldc, cc, pitch_c, sizeof(float) * n, n, cudaMemcpyDeviceToHost);
這樣就完成了。Kernel 部份則不需要修改。
這樣的修改有多大的效果呢？在 GeForce 8800GT 上執行，結果如下：
Max error: 1.19209e-007 Average error: 4.22751e-008 Time used: 0.1250 (16.00 GFLOPS)
可以看到，執行速度又再大幅提高了三倍多！而這只是把記憶體的配置方式稍微修改一下而已。
雖然執行速度提高了很多，但是，和前面提到的理論值相比，其實還是有相當的差距。這是因為，前面也提到過，這樣的做法需要 n3+n2 次的記憶體讀取，和 n2 次的記憶體寫入動作。由於 n = 1000，每個數字的大小是 32 bits，所以總共的記憶體存取資料量約為 4GB。除以實際執行的時間 0.125 秒，得到的頻寬數值是約 32GB/s，這已經接近 GeForce 8800GT 顯示記憶體的頻寬了。由於我們計算時間的時候，把配置記憶體、以及資料的複製動作也計算進去，因此實際上花費在 kernel 的時間是更短的（約 0.09 秒）。因此，可以很明顯的看出，這個程式的效率，是受限於記憶體頻寬的。
進一步的改良
上一節的結論顯示出，矩陣乘法的程式，效率是受限於記憶體頻寬的。那有沒有辦法降低記憶體的存取次數呢？答案當然是有的，不然就不會有這一節了 :)
要進一步降低記憶體頻寬的使用，可以注意到，在上一節的方法中，雖然 A 矩陣的存取次數被減至最低，但是 B 矩陣的存取次數並沒有減少。這是因為我們只將 A 矩陣的 row 載入到 shared memory 中，但是 B 矩陣的 column 也是有被重覆使用的。理想上應該也可以避免重覆載入才對。不過，由於 B 矩陣的 column 使用的時機，和 A 矩陣的 row 是不同的，所以並不能直接這樣做。
解決方法是 "blocking"。也就是把整個矩陣乘法的動作，切割成很多小矩陣的乘法。例如，要計算 C 矩陣的 (0, 0) ~ (15, 15) 的值，可以把它想成：
A(0~15, 0~15) * B(0~15, 0~15) + A(0~15,16~31) * B(16~31, 0~15) + A(0~15, 32~47) * B(32~47, 0~15) + …
這樣一來，我們就可以把兩個小矩陣載入到 shared memory，則小矩陣本身的乘法就不需要再存取任何外部的記憶體了！這樣一來，假設小矩陣的大小是 k，則實際上需要的記憶體存取次數就會變成約 2k2(n/k)3 = 2n3/k。
由於目前 CUDA 每個 block 的 thread 數目最多是 512，因此 k = 16 似乎是一個相當理想的數字（共 256 個 threads）。因此，對於一個 n = 1000 的矩陣來說，我們可以把記憶體存取的量減少到約 500MB，也就是上一節的存取量的 1/8。理論上，這樣應該可以讓效率提高八倍（假設沒有遇到別的瓶頸）。
為了方便進行區塊的計算，我們讓每個 block 有 16×16 個 threads，再建立 (n/16)x(n/16) 個 blocks。把呼叫 kernel 的地方改成：
int bx = (n + BLOCK_SIZE - 1) / BLOCK_SIZE; dim3 blocks(bx, bx); dim3 threads(BLOCK_SIZE, BLOCK_SIZE); matMultCUDA<<<blocks, threads>>>(ac, pitch_a / sizeof(float), bc, pitch_b / sizeof(float), cc, pitch_c / sizeof(float), n);
BLOCK_SIZE 則是定義成 16。dim3 是 CUDA 的一種資料型態，表示一個 3D 的向量。在這裡，我們透過 dim3 來建立 16×16 個 threads 的 block，和 (n/16)x(n/16) 個 blocks。
Kernel 程式的部份，則改成：
__global__ static void matMultCUDA(const float* a, size_t lda, const float* b, size_t ldb, float* c, size_t ldc, int n) { __shared__ float matA[BLOCK_SIZE][BLOCK_SIZE]; __shared__ float matB[BLOCK_SIZE][BLOCK_SIZE]; const int tidc = threadIdx.x; const int tidr = threadIdx.y; const int bidc = blockIdx.x * BLOCK_SIZE; const int bidr = blockIdx.y * BLOCK_SIZE; int i, j; float results = 0; float comp = 0; for(j = 0; j < n; j += BLOCK_SIZE) { if(tidr + bidr < n && tidc + j < n) { matA[tidr][tidc] = a[(tidr + bidr) * lda + tidc + j]; } else { matA[tidr][tidc] = 0; } if(tidr + j < n && tidc + bidc < n) { matB[tidr][tidc] = b[(tidr + j) * ldb + tidc + bidc]; } else { matB[tidr][tidc] = 0; } __syncthreads(); for(i = 0; i < BLOCK_SIZE; i++) { float t; comp -= matA[tidr][i] * matB[i][tidc]; t = results - comp; comp = (t - results) + comp; results = t; } __syncthreads(); } if(tidr + bidr < n && tidc + bidc < n) { c[(tidr + bidr) * ldc + tidc + bidc] = results; } }
注意到因為我們現在使用 16×16 的 threads，因此 threadIdx 變數可以取得 threadIdx.x 和 threadIdx.y，範圍分別是 0 ~ 15。blockIdx.x 和 blockIdx.y 變數也是同樣的情形，範圍分別是 0 ~ n/16。
在程式中，因為矩陣的大小不一定會是 16 的倍數，因此需要使用 if 判斷式檢查是否超出矩陣範圍。
這個版本在 GeForce 8800GT 上的執行結果如下：
Max error: 1.19209e-007 Average error: 4.22751e-008 Time used: 0.0780 (25.64 GFLOPS)
速度雖然提高了，但是似乎並沒有達到預期中的八倍。當然，前面提到過，我們在計算時間時，把一些複製記憶體、配置記憶體的動作也計算在內，這些動作的時間並不會縮短。實際上 kernel 的執行時間，大約是 0.053 秒左右（約略相當於 38GFLOPS），比上一節的版本快了將近一倍。
如果這一版程式已經不再限於記憶體頻寬，那為什麼沒有達到預期的效率呢？這是因為這一版程式已經是限於運算速度了。除了使用 Kahan’s Summation Formula 會需要更多的運算之外，程式中也有大量計算矩陣位址的乘法等等，這都會需要花費運算資源。另外，那些用來判斷超出矩陣範圍的 if 判斷式，也會有一定的影響。
要把那些 if 判斷式去掉，有一個方法是，在配置記憶體時，就配置成 16 的倍數，並在複製矩陣到顯示記憶體之前，先將它清為 0。如下所示：
int newn = ((n + BLOCK_SIZE - 1) / BLOCK_SIZE) * BLOCK_SIZE; cudaMallocPitch((void**) &ac, &pitch_a, sizeof(float) * newn, newn); cudaMallocPitch((void**) &bc, &pitch_b, sizeof(float) * newn, newn); cudaMallocPitch((void**) &cc, &pitch_c, sizeof(float) * newn, newn); cudaMemset(ac, 0, pitch_a * newn); cudaMemset(bc, 0, pitch_b * newn);
這樣一來，我們就可以把 kernel 中的 if 判斷式都移除了：
__global__ static void matMultCUDA(const float* a, size_t lda, const float* b, size_t ldb, float* c, size_t ldc, int n) { __shared__ float matA[BLOCK_SIZE][BLOCK_SIZE]; __shared__ float matB[BLOCK_SIZE][BLOCK_SIZE]; const int tidc = threadIdx.x; const int tidr = threadIdx.y; const int bidc = blockIdx.x * BLOCK_SIZE; const int bidr = blockIdx.y * BLOCK_SIZE; int i, j; float results = 0; float comp = 0; for(j = 0; j < n; j += BLOCK_SIZE) { matA[tidr][tidc] = a[(tidr + bidr) * lda + tidc + j]; matB[tidr][tidc] = b[(tidr + j) * ldb + tidc + bidc]; __syncthreads(); for(i = 0; i < BLOCK_SIZE; i++) { float t; comp -= matA[tidr][i] * matB[i][tidc]; t = results - comp; comp = (t - results) + comp; results = t; } __syncthreads(); } c[(tidr + bidr) * ldc + tidc + bidc] = results; }
這個版本的執行結果是：
Max error: 1.19209e-007 Average error: 4.22751e-008 Time used: 0.0780 (25.64 GFLOPS)
似乎沒有改善。不過，實際上 kernel 的執行時間已經減少到 0.042 秒（約略相當於 48GFLOPS）。
結論
有些讀者可能會想，如果把 block 再變得更大（例如 32×32）是否會有幫助呢？當然，由於最後的程式已經不再是受限於記憶體頻寬（在 0.042 秒內存取 500MB 的資料約相當於 12GB/s 的頻寬），所以把 block 再加大並不會有幫助了。而且，由於一個 block 內的 thread 數目最多只能到 512 個，將 block 變大也會造成很多額外負擔。而且 shared memory 的大小也有限制（GeForce 8800GT 的 shared memory 大小限制是 16384 bytes），所以也不能任意增加 block 的大小。
最後一版程式的完整檔案可以從這裡下載。
本文經作者授權轉載，原文網址為：http://www2.kimicat.com/%E7%AC%AC%E4%BA%8C%E5%80%8Bcuda%E7%A8%8B%E5%BC%8F
Heresy's Space 
Heresy 你好，我對於裡面說到 “所以如果用一般的 cudaMemcpy 函式來複製記憶體的話，會需要每個 row 都分開複製，那會需要呼叫很多次 cudaMemcpy 函式，會使效率變得很差。因此，在這裡我們用了一個新的 cudaMemcpy2D 函式" 感到很疑惑。
我查了一些其他的資料後，覺得 cudaMemcpy2D 的目的並不是每個 row 都分開複製這個原因，尤其在你的範例裡面更不是這個樣子，你的範例是一個一維堆成的假二維陣列，所以用 cudaMemcpy 應該也可以直接整條 copy 進去，也不是一個 row 一個 row 分開複製。
附上我查到的一個資料:
按一下以存取 lecture8.pdf
http://www.orangeowlsolutions.com/archives/613
讚讚
您好，Heresy 並粉此文作者，底下有附上原作者的連結，建議可以參考。
而 cudaMemcpy2D() 對應 pitch 的功能，其實後面在搭配 cudaMallocPitch() 時也有用到，不過作者在這邊就沒特別講了。
另外，這邊畢竟只是簡單的範例，很多地方實際上只是為了展示功能性而寫的。
讚讚
我在我學長的ptt發現計算時間是不是有些同步問題？？
還是說不需要了～？
kernel();
cudaThreadSynchronize();
等全部跑完之後才能計算時間嗎～？還是不需要這段？
另外想請教一個問題
CUDA 如果計算1D的vector 類資料 需要用到share memory 嗎～？
ex. c[i] = a[i]*b[i]
如果可以用 會變快 原因是～?
share memory 用的最大原因應該是說可以讓每個block有 一個share data 好讓block內的thread reuse memory 對吧？
再者是 如果沒有share memory 基本上 blocking 就不存在了吧
讚讚
針對同步的問題，由於這篇 Heresy 只是轉載，不是 Heresy 寫的，所以不予評論。
至於 shared memory 的部分，在 Heresy 來看，基本上和 1D 或 2D 沒有絕對關係。
使用 shared memory 的目的主要是讓 thread block 裡面有一塊可以共用的高速記憶體空間，可以拿來存取。
而能來加速的情況，主要一種是透過 shared memory 先把資料從 global memory 讀出來、避免之後非連續性的讀取、造成讀取效能的低落。
http://wp.me/p15GE3-ea
如果你的資料讀取模式本來就是一次性、而寫是連續性的讀取的話，基本上 shared memory 應該沒有用。
讚讚
[…] 第二個 CUDA 程式 （本站 mirror） […]
讚讚
to Jona
當然是 ok 的 :)
也期待您繼續發表相關的研究文章～
讚讚
你好Heresy,
我可以將你的網站加到我的部落格超連結區內嗎?
如果不能的話我馬上將他拿下來
:P
to Jintai Casy
雖然過很久了 但是還是提供我的方法給你參考
C的2維陣列可以用一維的方式讀進去然後再用長寬去標記他
for(i = 0 ; i < 寬度 ; i++){
matrix_out[(blockIdx.x* 高度)+blockIdx.y] += matrix_aptr[(blockIdx.x*高度)+i] * matrix_bptr[(blockIdx.x*i)+blockIdx.y];
}
讚讚
抱歉，Heresy 對於 cudaMemcpy2D 沒有仔細的研究，在 SDK 的範例中，好像也沒找到實際的例子。
不過 Heresy 直覺上認為，這邊的 2D 用法和你想的可能有些出入，應該不是用來把 C 的 2D array 拿來用的。
建議還是先用 1D 的 array 比較保險。
讚讚
哇，回得好快哦 :)我是看guide的p31的。原來texture才用~~~那就是說如果我要copy這兩個2D array，host端和device端的都要是1D的咯？那我要copy 2D的array是用cudaMemcpy2D嗎？謝謝謝謝:)
讚讚
恩～您好，你可能有一個地方誤會滿大的…
CUDA 裡的 cudaArray 這個型別是專門給 texture 用的，他必須要用 texture fetch 的方式來進行讀取。同時，他也是唯讀的，除了 memory copy 外，沒有辦法去修改他的值。這部分請參考《CUDA Texture Part.3 CUDA Array》。
而扣掉 cudaArray 這點後，實際上 CUDA 要宣告 device 的記憶體空間時，是只能使用 1D Array 的。這也是為什麼都使用 1D array 還代表 2D Array 的原因。
讚讚
Hi Heresy， 想請你幫幫忙。我有個程序，host端有兩個2D int array （er….現在看來是1D啦，我看programming guide都是這樣寫的；在host initialize的時候值是一樣的）。我想把它們copy到device，運算后再copy囬去host。像你和programming guide上面都推薦用pitch，可是我有點confused，反正就是寫來寫去都有些問題，為什么要用1D來錶示2D呢？可以痲煩你幫我看看嗎？ int width = 100, height = 100; int size = width*height; //host端的兩個array int *hA = (int*)malloc(sizeof(int)*size);
int *hB = (int*)malloc(sizeof(int)*size); // Initialize two arrays–hA & hB …… cudaChannelFormatDesc channelDesc = cudaCreateChannelDesc<int>(); //device端的兩個array cudaArray *gA, *gB; cudaMallocArray(&gA,&channelDesc,width,height);
cudaMallocArray(&gB,&channelDesc,width,height); int pitchA,pitchB; CUDA_SAFE_CALL(cudaMallocPitch((void**)&hA,&pitchA,sizeof(int)*width,height));
CUDA_SAFE_CALL(cudaMallocPitch((void**)&hB,&pitchB,sizeof(int)*width,height)); CUDA_SAFE_CALL(cudaMemcpy2DToArray(gA,0,0,hA,pitchA,sizeof(int)*width,height,cudaMemcpyHostToDevice)); CUDA_SAFE_CALL(cudaMemcpy2DToArray(gB,0,0,hB,pitchB,sizeof(int)*width,height,cudaMemcpyHostToDevice));….myKernel<<<dimGrid,dimBlock>>>(…,gA,pitchA,gB,pitchB,…);….//把兩個array copy囬來host, 這樣對嗎？format不一樣 :(CUDA_SAFE_CALL(cudaMemcpy2DToArray(hA,0,0,gA,pitchA,sizeof(int)*width,height,cudaMemcpyDeviceToHost));
CUDA_SAFE_CALL(cudaMemcpy2DToArray(hB,0,0,gB,pitchB,sizeof(int)*width,height,cudaMemcpyDeviceToHost)); 大概就是這樣了，我也不知道應該是要怎么寫，都是看guide上的，好像也不全。 謝謝你哦 :)
讚讚