我們已經知道了threads/blocks在CUDA端的組織方式，接下來我們學學多維度空間下的多執行緒模型，下面以矩陣乘法為例。

1. 行優先

儲存方式

二維矩陣在記憶體中的儲存方式受到程式語言的影響，主要可以分為兩種：行優先和列優先。對於程式語言諸如C/C++/CUDA而言，資料在記憶體中的組織方式是行優先。舉例說明行優先的儲存方式，如下，

給定一個3 $\times$ 3大小的矩陣 $A$

$A_{3\times 3}=\begin{pmatrix} a_{1,1} & a_{1,2} & a_{1,3} & a_{2,1} & a_{2,2} & a_{2,3} & a_{3,1} & a_{3,2} & a_{3,3} \end{pmatrix}$

位置索引

基於給定的 $n \times m$大小的矩陣 $A$，則 $a_{i,j}$在記憶體中的位置索引為 $i\times m + j$，其中， $0 \leq i \leq n-1$， $0 \leq j \leq m-1$.

2. CPU實現和分析

CPU核心演算法

根據行優先的位置索引，我們可以容易地計算矩陣的乘法如下，

for(int row = 0; row < matrix_width; row++)
{
    for(int col = 0; col < matrix_width; col++)
    {
        int single_element = 0;
        for(int k = 0; k < matrix_width; k++)
        {
            single_element += matrix_a_host[row * matrix_width + k] * matrix_b_host[matrix_width  * k + col];
        }
        matrix_c_host[row * matrix_width + col] = single_element;
    }
}

注：矩陣乘法需要保證第一個矩陣的列數等於第二個矩陣的行數。

演算法複雜度

基於CPU實現的矩陣乘法的演算法複雜度為 $O\left ( n^{3} \right )$。

3. GPU實現和分析

Threads/Blocks組織方式

出於程式設計上的方便，我們採用一維的blocks組織方式，二維的threads組織方式，如下，

dim3 dimGrid(1, 1, 1);
dim3 dimBlock(3, 3, 1);

注：預設情況下threads/blocks均採用一維的組織方式。

GPU核心演算法

__global__ void matrixs_1D_multiplication(int *matrix_a_dev, int *matrix_b_dev, int *matrix_c_dev, int matrix_width)
{
    int row = threadIdx.y;
    int col = threadIdx.x;

    if(row < matrix_width && col < matrix_width)
    {
        for(int k = 0; k < matrix_width; k++)
        {
            matrix_c_dev[row * matrix_width + col] += matrix_a_dev[row * matrix_width + k] * matrix_b_dev[k * matrix_width + col]; 
        }
    }
}

演算法複雜度

基於GPU實現的矩陣乘法的演算法複雜度為 $O\left ( n \right )$。

4. 編譯除錯

CPU版本

原始碼：matrix_mul_cpu.cu

編譯

nvcc matrix_mul_cpu.cu -o main

執行

./main

GPU版本

原始碼：matrix_mul_gpu.cu

編譯

nvcc matrix_mul_gpu.cu -o main

執行

./main