《CUDA C Programming Guide》

第一章 介绍

GPU在进行浮点数运算时超高的精确度和速度,其原因在于GPU很擅长解决计算密集型和高并行计算。

多核CPU和GPU的一个挑战在于怎么编写并行的程序来利用这些多核。

CUDA并行变成模式就是希望用一个较低的学习曲线来解决这个问题。

三个关键的抽象:线程组等级、共享内存、障碍同步

第二章 编程模型

2.1 Kernels

kernel程序用__global__定义。并且用一种新的执行配置语法<<<...>>>来决定CUDA的线程数量。每一个线程都会在给定的“线程ID”下执行kernel,线程ID可以通过kernel内部的threadIdx变量来获取。

下面的代码显示了向量加法。 总共启动了N个线程,每个都执行一个对应位相加运算。

1
2
3
4
5
6
7
8
9
10
11
12

//Kernel definition
__global__ void VecAdd(float* A, float* B, float* C){
int i = threadIdx.x;
C[i] = A[i] + B[i];
}

int main(){
...
// Kernel invocation with N threads
VecAdd<<<1,N>>>(A,B,C);
}

2.2 线程结构(Thread Hierarchy)

为了方便,threadIdx包含三个组成向量,因此,线程可以使用一维,二维或者三维的thread index,由此,可以表示一维,二维或者三维的线程块(“thread block”)。

线程的索引和它的线程ID是直接关联的。对于一维线程块来说,它们是一样的。对于二维线程块来说,索引为 $(x,y)$ 的线程,其线程ID为:$(x+yD_x)$ 。对于三维线程块来说,索引为 $(x,y,z)$ 的线程,其线程ID为: $(x+ yD_x +zD_xD_y)$ 。

下面的代码显示了矩阵加法。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

//Kernel definition
__global__ void MatAdd(float A[N][N], float B[N][N], float C[C][C]){
int i = threadIdx.x;
int j =threadIndx.y;
C[i][j] = A[i][j] + B[i][j];
}

int main(){
...
//Kernel invocation with one block of N*N*1 threads
int numBlocks = 1;
dim3 threadsPerBlock(N,N); //<<<>>>语法可以接受int类型或者dim3类型的数据
MatAdd<<<numBlocks, threadsPerBlock>>>(A,B,C);
}

thread block的索引可以通过blockIdx变量获得,thread block的维度的size可以通过blockDim变量获得。

扩展上面的MatAdd()代码,使其可以处理多blocks,代码如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

//Kernel definition
__global__ void MatAdd(float A[N][N], float B[N][N], float C[N][N]){
int i = blockIdx.x * blockDim.x + threadIdx.x;
int j = blockIdx.y * blockDim.y + threadInx.y;
if(i<N && j<N){
C[i][j] = A[i][j] + B[i][j];
}
}

int main(){
...
//Kernel invocation
dim3 threadsPerBlock(16, 16);
dim3 numBlocks(N / threadsPerBlock.x, N / threadsPerBlock.y);
MatAdd<<<numBlocks, threadsPerBlock>>>(A, B, C);
...
}

线程块必须独立执行,因此它们必须能够在任意的顺序下并行执行。

CUDA使用__syncthreads()来协调共享内存和各个线程之间的执行。

2.3 内存结构(Memory Hierarchy)

由两种额外的只读内存空间可以被所有线程访问到:常量内存空间(constant memory spaces)和纹理内存空间(texture memory spaces)。

2.4 异构编程

CUDA编程模型假设所有线程都是在物理上单个的设备上运行的,每个设备都当做是主机的从处理程序。比如,kernel在GPU上执行,而剩下的C程序在CPU上执行。

CUDA编程模型还假设主机和设备都DRAM中各自维护自己的内存,对应这“主机内存(host memory)”和“设备内存(divice memory)”。因此,程序通过调用CUDA运行时来管理主机和设备的内存。

2.5 计算能力(Compute Capability)

计算能力用设备的版本号标识,有时也称为“SM version”。有major revision number“X”和minor revision number“Y”组成,记为X.Y。

第三章 编程接口(Programming Interface)

CUDA提供了诸多扩展语法来并行变成(如kernel),任何包含这些扩展语法的源文件都需要需要nvcc编译器来编译。

3.1 Compilation with NVCC

Kernels can be written using the CUDA instruction set architecture, called PTX, which
is described in the PTX reference manual.

3.1.1 编译工作流

3.1.1.1 离线编译 Offline Compilation

nvcc编译器可以编译包含host code和device code的混合代码。nvcc首先会将device code从host code中分离出来,然后,会进行以下两步:

  • 编译device code到assembly form(PTX code)或者binary form(cubin objec)
  • 修改host code,将<<<...>>>语法替换成必要的CUDA C运行时函数。

被修改的host code要么输出成C code,要么把直接输出成object code。

3.1.1.2 即时编译 Just-in-Time Compilation

二进制兼容性 Binary Compatibility

Binary code is architecture-specific。通过-code=sm_35的形式来指定特定的architecture。注意,编译好的二进制代码是向上兼容的,也就是如果二进制代码设定的计算能力版本号为 $X.y$ ,那么该代码就只能运行在 $X.z$上,其中 $z\ge y$ 。

3.1.3 PTX 兼容性

一些PTX指令仅仅支持在高计算机能力版本的GPU中使用。

指针某些计算能力版本的PTX code总是可以转换成binary code,进而都更高计算能力的版本中使用。

其他有些版本,则不能直接使用。

3.1.4 应用兼容性

为了满足bianry兼容性和PTX兼容性,推荐使用即时编译。

通过-gencode-arch-code编译选项来指定计算能力版本号及其兼容性。

3.1.5 C/C++ 兼容性

CUDA源文件的前端遵循C++语法规则。host code可以完美支持C++语言发。但是,device code支持C++语法中的一个子集,详细可以参见C/C++ Language Support。

3.1.6 64-Bit Compatibility

64位的nvcc会将device code编译成64位模式(即,指针占64位)。此时host code必须在32位模式下执行。

32位device、host同理

-m64 选项可以切换nvcc的位数。

3.2 CUDA C Runtime

运行时使用cudart库实现。 或者使用cudart.liblibcudart.a进行静态链接,或者使用cudart.dllcudart.so进行动态链接。

正如前面异构编程提到的。CUDA编程模型将系统看走是由device和host组成的一个整体,二者管理各自的内存。

3.2.1 初始化

没有专门的初始化函数
在初始化阶段,runtime为系统中的每一个device创建上下文。

3.2.2 设备内存 Device Memory

Kernels操控device memory之外的代码。所以runtime会提供allocate,deallocate和copy等函数来复制device memory,同时在host和device内存之间传输数据。

Divice memory可以申请linear memory或者CUDA arrays。

CUDA arrays是不透明的内存形式,专门用于优化texture fetching。

Linear memory存在于40位的地址空间中,可以使用cudaMalloc()函数申请,用cudaFree()函数释放,同时,可以用cudaMemcpy()函数来将数据在host和device之间交换。

cudaMallocPitch()cudaMalloc3D()推荐用来申请2D和3D数组。cudaMemcpy2D

cudaGetSymbolAddress()用来检索指向内存的地址。

cudaGetSymbolSize()可以得到申请内存的大小。

3.2.3 共享内存 Shared Memory

Shared Memory在Thread Hierarchy上要比global memory更快。

矩阵乘法

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
// Matrices are stored in row-major order:
// M(row, col) = *(M.elements + row * M.width + col)
typedef struct {
int width;
int height;
float* elements;
} Matrix;
// Thread block size
#define BLOCK_SIZE 16
// Forward declaration of the matrix multiplication kernel
__global__ void MatMulKernel(const Matrix, const Matrix, Matrix);
// Matrix multiplication - Host code
// Matrix dimensions are assumed to be multiples of BLOCK_SIZE
void MatMul(const Matrix A, const Matrix B, Matrix C)
{
// Load A and B to device memory
Matrix d_A;
d_A.width = A.width; d_A.height = A.height;
size_t size = A.width * A.height * sizeof(float);
cudaMalloc(&d_A.elements, size);
cudaMemcpy(d_A.elements, A.elements, size,
cudaMemcpyHostToDevice);
Matrix d_B;
d_B.width = B.width; d_B.height = B.height;
size = B.width * B.height * sizeof(float);
cudaMalloc(&d_B.elements, size);
cudaMemcpy(d_B.elements, B.elements, size,
cudaMemcpyHostToDevice);
// Allocate C in device memory
Matrix d_C;
d_C.width = C.width; d_C.height = C.height;
size = C.width * C.height * sizeof(float);
cudaMalloc(&d_C.elements, size);
// Invoke kernel
dim3 dimBlock(BLOCK_SIZE, BLOCK_SIZE);
dim3 dimGrid(B.width / dimBlock.x, A.height / dimBlock.y);
MatMulKernel<<<dimGrid, dimBlock>>>(d_A, d_B, d_C);
// Read C from device memory
cudaMemcpy(C.elements, Cd.elements, size,
cudaMemcpyDeviceToHost);
}
// Free device memory
cudaFree(d_A.elements);
cudaFree(d_B.elements);
cudaFree(d_C.elements);
// Matrix multiplication kernel called by MatMul()
__global__ void MatMulKernel(Matrix A, Matrix B, Matrix C)
{
// Each thread computes one element of C
// by accumulating results into Cvalue
float Cvalue = 0;
int row = blockIdx.y * blockDim.y + threadIdx.y;
int col = blockIdx.x * blockDim.x + threadIdx.x;
for (int e = 0; e < A.width; ++e)
Cvalue += A.elements[row * A.width + e]
* B.elements[e * B.width + col];
C.elements[row * C.width + col] = Cvalue;
}

3.2.4 Page-Locked Host Memory

plhm有诸多好处,如:

  • 可以使copy更快
  • 可以直接映射到地址空间,消除copy到device memory的需求
  • 在fron-side bus上, plhm的bandwidth更高。

但是,plhm是稀缺资源,因此,在申请时一定要注意不能过度使用。

3.2.4.1 移动内存 Portable Memory

一块page-locked内存可以备用在与任意一个设备的连接上。但是默认情况下,使用page-locked内存的好处仅仅会在当前block所申请的设备(或者那些共享了内存的设备)上显现。

为了使这些好处对于所有设备来说都是可行的,block需要将flag cudaHostAllocPortable传递到cudaHostAlloc中去,或者将flag cudaHostRegisterPortable传递到cudaHostRegister中去。

3.2.4.2 Write-Combining Memory

该内存通过释放L1和L2的缓存资源,使得更多缓存资源可以给后面的应用使用。

从host中读取write-combining 内存比较慢,所以常常只用于写。

3.2.4.3 Mapped Memory

使用cduaHostAllocMapped或者cudaHostRegisterMapped可以使一块page-locked host内存映射到设备的地址空间中。

因此,这样的内存块往往有两个地址,分别是host地址和device地址。

直接kernel中访问host内存有以下几点优势:

  • 无需在device中申请内存块,也无需在device和host中来还传输数据
  • 在执行kernel时,无需使用“流”来overlap数据的传输

3.2.5 异步并发执行 Asynchronous Concurrent Execution

CUDA将下列操作当作一个独立的任务,它们可以互相并发执行:

  • 在host上计算
  • 在device上计算
  • 将host上的内存传输到device上
  • 将device上的内存传输到host上

3.2.5.1 在Host和Device之间并发执行

通过异步调用,许多device操作都可以排成队列,进而利用CUDA driver来执行。这减轻了host线程管理device的负担,使得它可以执行其他任务。下列device操作相对于host来说是异步的。

  • 启动Kernel
  • 将内存拷贝到一个单一的deivce内存中
  • 将内存从host拷贝到device中小于64将内存从host拷贝到device中小于64KB的内存块中
  • 用带有Async后缀的函数来进行内存拷贝
  • 内存设定函数的调用(Memory set function calls)

3.2.5.2 Concurrent Kernel Execution

计算能力告诉2.x的版本可以并发启动多个kernels。

可以通过concurrentKernel设备属性来查看是否支持并发启动。

3.2.5.3 Overlap of Data Tranfer and Kernel Execution

有一些设备可以从GPU利用kernel执行来进行异步的内存拷贝,应用需要想设备询问是否具有这种功能,可以通过检查asyncEngineCoune设备属性来查看,如果大于0就是支持的。

3.2.5.4 并发数据传输 Concurrent Data Transfers

一些计算能力在2.x以上的设备可以在进行数据传输是进行overlap,应用可以通过asyncEngineCount来检查设备是否支持该功能

3.2.5.5 流 Streams

应用通过“流”来管理并发的操作。一个“流”代表这一串由“命令”组成的序列

3.2.5.5.1 创建和销毁:Creation and Destruction

下面的代码创建了两个流,同时在page-locked内存中申请了flaot类型的hosrtPtr

1
2
3
4
5
6

cudaStream_t stream[2];
for (int i = 0; i < 2; ++i)
cudaStreamCreate(&stream[i]);
float* hostPtr;
cudaMallocHost(&hostPtr, 2 * size);

这两条流都经过下面的代码定义,执行一个从host到device的内存拷贝操作,一个kernel启动操作,以及一个从device到host的的内存拷贝操作。

1
2
3
4
5
6
7
8
9

for (int i = 0; i < 2; ++i) {
cudaMemcpyAsync(inputDevPtr + i * size, hostPtr + i * size,
size, cudaMemcpyHostToDevice, stream[i]);
MyKernel <<<100, 512, 0, stream[i]>>>
(outputDevPtr + i * size, inputDevPtr + i * size, size);
cudaMemcpyAsync(hostPtr + i * size, outputDevPtr + i * size,
size, cudaMemcpyDeviceToHost, stream[i]);
}

通过cudaStreamDestroy()可以释放流:

1
2
3
4

for (int i =0;i<2;i++){
cudaStreamDestroy(stream[i]);
}
3.2.5.5.2 默认流 Default Stream

对于使用--default-stream per-thread编译选项的代码来说,默认流是一条常规的流,并且每个主线程都有它自己的默认流

对于使用--default-stream legacy编译选项的代码来说,默认流是一条特殊的流,被称为NULL stream ,并且每一个设备都具有一条单一的NULL stream供所有的host线程使用。它的特殊性源自于它会隐式的进行同步操作。

3.2.5.5.3 显式同步

以下几种方式可以对流进行显式同步操作:

  • cudaDeviceSynchronize()
  • cudaStreamSynchronize()
  • cudaStreamWaitEvent()
  • cudaStreamQuery()

为了避免不必要的速率降低,以上所有的同步函数通常用于timing purpoese或者用于孤立一个拷贝或启动的失败。

3.2.5.5.4 隐式同步

从两个不同流传过来的指令中的其中一条出现了问题,那么就无法并发运行。

同步操作越晚执行越好。

3.2.5.5.5 Overlapping 行为

对于两个流里面的操作,如一个是从host拷贝内存到device,另一个是从device拷贝内存到host,则这两条指令直接存在Overlap。

3.2.5.5.6 回调函数Callbacks

运行时提供了可以在流中的任何位置插入回调函数的指令:cudaStreamAddCallback()

3.2.5.5.7 流优先级 Stream Priorities

相对优先级使用cudaStreamCreateWithPriority()

获取优先级范围[highest priority,lowest priority]使用cudaDeviceGetStramPriorityRange

下面的代码包含了获取当前设备的优先级范围:

1
2
3
4
5
6
7
// get the range of stream priorities for this device
int priority_high, priority_low;
cudaDeviceGetStreamPriorityRange(&priority_low, &priority_high);
// create streams with highest and lowest available priorities
cudaStream_t st_high, st_low;
cudaStreamCreateWithPriority(&st_high, cudaStreamNonBlocking, priority_high);
cudaStreamCreateWithPriority(&st_low, cudaStreamNonBlocking, priority_low);

3.2.5.6 事件 Events

runtime 同时还提供了密切监视device进程的方法,主要是通过appliocation来异步记录事件完成时的事件点。

3.2.5.6.1 创建和销毁

下面的代码创建了2个事件:

1
2
3
cudaEvent_t event1, event2;
cudaEventCreate(&event1);
cudaEventCreate(&event2);

下面的代码销毁了2个事件:

1
2
cudaEventDestroy(event1);
cudaEventDestroy(event2);

3.2.5.6.2 运行时间

下面的代码可以用来监视事件的运行时间(从start到end)

1
2
3
4
5
6
7
8
9
10
11
12
13
14

cudaEventRecord(start, 0);
for (int i = 0; i < 2; ++i) {
cudaMemcpyAsync(inputDev + i * size, inputHost + i * size,
size, cudaMemcpyHostToDevice, stream[i]);
MyKernel<<<100, 512, 0, stream[i]>>>
(outputDev + i * size, inputDev + i * size, size);
cudaMemcpyAsync(outputHost + i * size, outputDev + i * size,
size, cudaMemcpyDeviceToHost, stream[i]);
}
cudaEventRecord(stop, 0);
cudaEventSynchronize(stop);
float elapsedTime;
cudaEventElapsedTime(&elapsedTime, start, stop);

3.2.5.7 同步调用

可以使用cudaSetDeviceFlags()配合一些特定的flags??

3.2.6 Multi-Device System

3.2.6.1 Device Enumeration

一个host系统可以有多个devices,下面的代码展示了如何枚举这些设备,并查询它们的属性,决定可启用CUDA的设备数量。

1
2
3
4
5
6
7
8
9
int deviceCount;
cudaGetDeviceCount(& deviceCount);
int device;
for( device = 0; device<deviceCount; device++){
cudaDeviceProp deviceProp;
cudaGetDeviceProperties(&deviceProp, device);
printf("Device %d has compute capability %d.%d. \n"),
device, deviceProp.major, deviceProp.minor);
}

3.2.6.2 Device Selection

host线程可以在任何时候使用cudaSetDevice()来设置device。device的内存申请和kernel启动都会在当前设置的device上进行,如果没有调用该函数,则当前的设备默认为0.

下面的代码展示了如何设置当前的device,以及申请内存和执行kernel

1
2
3
4
5
6
7
8
9
size_t size = 1024 * sizeof(float);
cudaSetDevice(0);
float* p0;
cudaMalloc(&p0, size);// Allocate memory on device 0
MyKernel<<<1000, 128>>>(p0);// Launch kernel on device 0
cudaSetDevice(1);// Set device 1 as current
float* p1;
cudaMalloc(&p1,size);// Allocate memory on device 1
MyKernel<<<1000,128>>>(p1);// Launch kernel on device 1

3.2.6.3 Stream and Event Behavior

如果stream没有绑定到当前的device,那么kernel launch就会失败

即使stream没有绑定到当前的device,memory copy也会成功

如果input event和input stream 绑定到了不同的diveces上面,那么cudaEventRecord() will fail

如果两个input events绑定到了不同的devices上面,那么cudaEventElapsedTime() will fail

即使input event绑定到了不同于当前device的device上面,cudaEventSynchronize()cudaEventQuery()仍然will succeed