原文地址，純翻譯

https://developer.nvidia.com/blog/easy-introduction-cuda-c-and-c/

　　這是cuda平行計算平臺 c和c++介面系列的第一篇文章。學習前要求熟練掌握c，針對cuda fortran程式設計的帖子也會同步更新。這兩個系列將涵蓋cuda平臺上平行計算基本概念。從這裡開始，除非有特別說明，否則我將使用屬於"cuda c"作為cuda c和c++的簡寫。cuda c本質上是帶有一些拓展的c/c++，允許使用多個並行執行緒在gpu上執行的函式。

CUDA程式設計模型基礎

　　在接觸cuda c程式碼之前，那些剛接觸cuda的人最好先了解cuda程式設計模型基本描述和其中的一些術語。

　　cuda程式設計模型是一種同時使用cpu和gpu的異構模型。在cuda中，host指cpu及其記憶體，device指gpu及其記憶體，host上執行的程式碼可以管理host和device上的記憶體、啟動kernel(核函式)，這些kernel是device上執行的函式，他們由gpu上的許多執行緒併發執行。

　　鑑於cuda程式設計模型的異構性質，cuda c程式的普遍操作順序是：

　　1.宣告和分配host端，device端記憶體

　　2.初始化host端資料

　　3.將資料從host端傳送到device端

　　4.執行一個或多個核函式

　　5.將結果從device端傳回host端

　　記住這個操作流程，讓我們來看一個cuda c的例子

第一個cuda c程式

　　上一篇文章中，我介紹了六種SAXPY(Scalar Alpha X Plus Y)的方法,其中就包括了cuda c版本，SAXPY表示單精度A*X+Y，對於平行計算來說是一個很好的hello world程式。在這篇文章中，我將展示一個cuda SAXPY的更完整版本，詳細說明做了什麼以及為什麼這樣做，完整的SAXPY程式碼如下：

#include <stdio.h>

__global__ void saxpy(int n,float a,float *x,float *y)
{
    int i=blockIdx.x*blockDim.x+threadIdx.x;
    
 
if(i<n)
    {
        y[i]=a*x[i]+y[i];
    }
}

int main(void)
{
    int N=1<<20;    // 1左移20
    float *x,*y,*d_x,*d_y;
    x=(float*)malloc(N*sizeof(float));
    y=(float*)malloc(N*sizeof(float));

    cudaMalloc(&d_x,N*sizeof(float));
    cudaMalloc(&d_y,N*sizeof(float));

    for(int i=0;i<N;++i)
    {
        x[i]=1.0f;
        y[i]=2.0f;
    }

    cudaMemcpy(d_x,x,N*sizeof(float),cudaMemcpyHostToDevice);
    cudaMemcpy(d_y,y,N*sizeof(float),cudaMemcpyHostToDevice);

    // perform SAXPY in 1M elements
    int threads=256;
    int blocks=(N+threads-1)/threads;   //  上取整
    saxpy<<<blocks,threads>>>(N,2.0f,d_x,d_y);
    
    cudaMemcpy(y,d_y,N*sizeof(float),cudaMemcpyDeviceToHost);

    float maxError=0.0f;
    for(int i=0;i<N;++i)
    {
        maxError=max(maxError,abs(y[i]-4.0f));
    }
    printf("Max error: %f\n",maxError);

    cudaFree(d_x);
    cudaFree(d_y);
    free(x);
    free(y);

    return 0;
}

　　函式saxpy就是在gpu上執行的核函式，main函式是host端程式碼，讓我們先從host端程式碼開始討論這個程式

Host端程式碼

　　主函式聲明瞭兩對陣列

float *x,*y,*d_x,*d_y;
x=(float*)malloc(N*sizeof(float));
y=(float*)malloc(N*sizeof(float));

cudaMalloc(&d_x,N*sizeof(float));
cudaMalloc(&d_y,N*sizeof(float));

　　指標x和y指向host端陣列，使用malloc函式申請記憶體空間，d_x,d_y陣列指向device端陣列，使用cuda runtime api中的cudaMalloc函式申請空間。cuda程式設計中host端和device端具有獨立的記憶體空間，兩者都可以通過host端程式碼進行管理（cuda c核心還可以在支援他的裝置上分配裝置記憶體）。

　　host端程式碼隨後初始化host端陣列。在此，我們將x設為元素全部為1的陣列，y設為元素全為2的陣列

for(int i=0;i<N;++i)
{
    x[i]=1.0f;
    y[i]=2.0f;
}

　　要初始化device端陣列，我們只需使用cudaMemcpy函式將資料從x，y拷貝到device端對應陣列d_x,d_y上，這個過程就像c語言中的memcpy函式，唯一的區別就是cudaMemcpy需要地四個引數來指定資料拷貝到方向（host端到device端還是device端到host端），在此例中，我們使用cudaMemcpyHostToDevice表示資料從host拷貝到device端。

cudaMemcpy(d_x,x,N*sizeof(float),cudaMemcpyHostToDevice);
cudaMemcpy(d_y,y,N*sizeof(float),cudaMemcpyHostToDevice);

　　在執行完核函式後，為了將資料d_y從device端拷貝回host端的y，我們使用cudaMemcpy函式，第四個引數指定為cudaMemcpyDeviceToHost

cudaMemcpy(y,d_y,N*sizeof(float),cudaMemcpyDeviceToHost);

啟動核函式

　　saxpy核函式由以下程式碼啟動

saxpy<<<blocks,threads>>>(N,2.0f,d_x,d_y);

　　<<<>>>中間的資訊是核函式的執行配置，具體指有多少裝置執行緒並行執行核心。在cuda中，軟體有一個執行緒層次結構，它模仿執行緒處理器在gpu上的分組方式。在cuda程式設計模型中，我們通常說通過grid來啟動執行緒block。執行配置中的第一個引數指定grid中的block數目，第二個引數指定一個block中的執行緒數。

　　block和grid可以通過給dim3傳值構造成1維，2維或3維度資料結構（由cuda定義的具有x,y,z策劃那個圓的簡單結構），但是對於這個簡單的例子，我們只需要一個維度，所以我們傳遞一個整數。在這個例子中，我們設定一個block中包含256個執行緒，使用取整演算法來計算得到處理陣列N個元素所需要的block數目（(N+256-1)/256）。

　　對於陣列中元素的數量不能被執行緒塊大小整除的情況，核心程式碼必須檢查越界記憶體訪問。（我認為就是檢測取整多出來的那部分）

釋放記憶體

　　程式結束時，我們需要釋放掉所有申請的記憶體空間。對於用cudaMalloc()申請的device端記憶體，使用cudaFree()來釋放。對於host端申請單記憶體，使用free()來釋放。

cudaFree(d_x);
cudaFree(d_y);
free(x);
free(y);

device端程式碼

　　現在我們來看核函式程式碼

__global__ void saxpy(int n,float a,float *x,float *y)
{
    int i=blockIdx.x*blockDim.x+threadIdx.x;
    if(i<n)
    {
        y[i]=a*x[i]+y[i];
    }
}

　　在cuda中，我們使用__global__宣告符來宣告核函式。device程式碼中定義的變數不需要指定為device變數，因為他們被指定駐留在device上。（我的理解是核函式裡宣告的變數，他的生命週期就是該執行緒上執行的核函式生命週期）。在本例中，n，a和i變數將由每個執行緒儲存在一個暫存器中，並且指標x和y必須是指向device記憶體地址空間的指標。因為當我們從host端啟動核函式時，我們要將d_x,d_y傳遞給kernel。前兩個變數n和a，我們沒有明確的程式碼將它們從host端傳到device端，由於函式引數在c/c++中預設按值傳遞，因此cuda執行時可以自動處理這些值到device的傳輸。cuda runtime api這個特性使得gpu上啟動核函式非常自然和容易，幾乎與呼叫c函式一樣。

　　核函式saxpy的程式碼只有兩行，前面我們提到過，核函式是由多執行緒併發執行的。如果我們希望每個執行緒處理陣列的一個元素，那麼我們需要一種區分和識別每個執行緒的方法。cuda定義了變數blockDim，blockIdx和threadIdx。這些預定義變數的型別為dim3，類似於主機程式碼中的執行配置引數。預定義變數blockDim包含在核函式啟動的第二個執行配置引數中，指定block中的執行緒個數（維度）。預定義變數threadIdx和blockIdx分別表示block內的執行緒索引和grid內地block索引。表達方式：

int i=blockIdx.x*blockDim.x+threadIdx.x;

　　生成用於訪問陣列元素的全域性索引。我們在這個例子中沒有使用gridDim，gridDim表示核函式啟動的第一個執行配置引數中指定的grid尺寸（一個grid裡面有多少個block）

　　在使用索引訪問陣列元素之前，會根據元素個數n檢查索引是否合法，以確保沒有記憶體越界訪問。如果陣列中的元素不能被block中的執行緒數整除，（由於向上取整）核函式啟動的執行緒數會大於陣列大小，此時需要進行檢查。核函式的第二行逐個元素執行SAXPY操作，除了邊界檢查外，他的結果等於host端通過迴圈實現的SAXPY

if(i<n) y[i]=a*x[i]+y[i];

編譯和執行程式碼

　　cuda c編譯器，nvcc是NVIDIA CUDA Toolkit（http://www.nvidia.com/content/cuda/cuda-toolkit.html）的一部分。為了編譯我們的SAXPY例子，我們將程式碼儲存為.cu格式，檔名是saxpy.cu。我們可以使用nvcc編譯它

nvcc -o saxpy saxpy.cu

　　我們可以執行程式碼

./saxpy

總結與結論

　　通過使用cuda c對saxpy簡單實現，我們現在瞭解了cuda程式設計的基礎知識。將c程式碼移植到cuda c上只需要幾個c擴充套件：核函式的__global__宣告符；核函式啟動時的執行配置（<<<blocks,threads>>>）；用於識別和區分並行執行核函式的gpu執行緒的內建裝置變數blockDim,blockIdx和threadIdx。

　　異構cuda程式設計模型的一個優點是，可以增量的將現有程式碼從c移植到cuda c，一次一個核心。

　　在本系列的下一篇文章中，我麼將瞭解一些效能測量和指標。