《并行与分布式程序设计》课程教学参考书：CUDA C PROGRAMMING（CUDA编程指南4.0中文版）

第三章编程接口 19 3.2.11节引入了多种运行时提供的函数，以和两大主要的图形API OpenGL 和Direct3.D互操作。 3.2.1初始化 运行时没有显式的初始化函数：在初次调用运行时函数（更精确地，不在参 考手册中设备和版本管理节中的任何函数)时初始化。在计算运行时函数调用的 时间和解析初次调用运行时产生的错误码时必须牢记这点。 在初始化时，运行时为系统中的每个设备建立一个上下文(3.3节提供了上 下文的更多细节)。这个上下文作为设备的主要上下文，被应用中的主机线程共 享。这些都是隐式发生的，运行时并没有将主要上下文展示给应用。 当主机线程调用cudaDeviceReset()时，这销毁了主机线程操作的设备的主上 下文。主机线程调用的运行时函数将为设备重新建立一个主上下文。 3.2.2设备存储器 正如2.4节所提到的，CUDA编程模型假定系统包含主机和设备，它们各有 自己独立的存储器。内核不能操作设备存储器，所以运行时提供了分配，释放， 拷贝设备存储器和在设备和主机间传输数据的函数。 设备存储器可被分配为线性存储器或CUDA数组。 CUDA数组是不透明的存储器层次，为纹理获取做了优化。它们的细节在 3.2.10节。 计算能力1.x的设备，其线性存储器存在于32位地址空间内，计算能力2.0 的设备，其线性存储器存在于40位地址空间内，所以独立分配的存储器实体能 够通过指针引用，如二叉树。 典型地，线性存储器使用cudaMalloc()分配，通过cudaFree()释放，使用 cudaMemcpy()在设备和主机间传输。在2.1节的向量加法代码中，向量要从主机 存储器复制到设备存储器： ∥Device code global void VecAdd(float*A,float*B,float*C,int N){ int i=blockDim.x blockldx.x threadldx.x; if(i<N) C[门=A[]+B; ∥Host code int main() int N=...; size t size=N*sizeof(float); //Allocate input vectors h A and h B in host memory float*h A=(float*)malloc(size);

第三章 编程接口 19 3.2.11 节引入了多种运行时提供的函数，以和两大主要的图形 API OpenGL 和 Direct3D 互操作。 3.2.1 初始化 运行时没有显式的初始化函数；在初次调用运行时函数（更精确地，不在参 考手册中设备和版本管理节中的任何函数）时初始化。在计算运行时函数调用的 时间和解析初次调用运行时产生的错误码时必须牢记这点。 在初始化时，运行时为系统中的每个设备建立一个上下文（3.3 节提供了上 下文的更多细节）。这个上下文作为设备的主要上下文，被应用中的主机线程共 享。这些都是隐式发生的，运行时并没有将主要上下文展示给应用。 当主机线程调用 cudaDeviceReset()时，这销毁了主机线程操作的设备的主上 下文。主机线程调用的运行时函数将为设备重新建立一个主上下文。 3.2.2 设备存储器 正如 2.4 节所提到的，CUDA 编程模型假定系统包含主机和设备，它们各有 自己独立的存储器。内核不能操作设备存储器，所以运行时提供了分配，释放， 拷贝设备存储器和在设备和主机间传输数据的函数。 设备存储器可被分配为线性存储器或 CUDA 数组。 CUDA 数组是不透明的存储器层次，为纹理获取做了优化。它们的细节在 3.2.10 节。 计算能力 1.x 的设备，其线性存储器存在于 32 位地址空间内，计算能力 2.0 的设备，其线性存储器存在于 40 位地址空间内，所以独立分配的存储器实体能 够通过指针引用，如二叉树。 典型地，线性存储器使用 cudaMalloc()分配，通过 cudaFree()释放，使用 cudaMemcpy()在设备和主机间传输。在 2.1 节的向量加法代码中，向量要从主机 存储器复制到设备存储器： // Device code __global__ void VecAdd(float* A, float* B, float* C, int N) { int i = blockDim.x * blockIdx.x + threadIdx.x; if (i < N) C[i] = A[i] + B[i]; } // Host code int main(){ int N = ...; size_t size = N * sizeof(float); // Allocate input vectors h_A and h_B in host memory float* h_A = (float*)malloc(size);

CUDA编程指南4.0中文版 20 float*h B=(float*)malloc(size); /Initialize input vectors /Allocate vectors in device memory float*d A; cudaMalloc((void**)&d A,size); float*d B: cudaMalloc((void**)&d B,size); float*d C; cudaMalloc((void**)&d C,size); //Copy vectors from host memory to device memory cudaMemcpy(d_A,h_A,size,cudaMemcpyHostToDevice); cudaMemcpy(d_B,h_B,size,cudaMemcpyHostToDevice); ∥Invoke kernel int threadsPerBlock=256; int blocksPerGrid=(N+threadsPerBlock-1)/threadsPerBlock; VecAdd<<<blocksPerGrid,threadsPerBlock>>>(d A,d B,d C,N); /Copy result from device memory to host memory //h C contains the result in host memory cudaMemcpy(h_C,d_C,size,cudaMemcpyDeviceToHost); ∥Free device memory cudaFree(d_A); cudaFree(d_B); cudaFree(d_C); ∥Free host memory … } 线性存储器也可以通过cudaMallocPitch()和cudaMalloc.3DO分配.在分配2D 或3D数组的时候，推荐使用，因为这些分配增加了合适的填充以满足5.3.2.1节 的对齐要求，在按行访问时或者在二维数组和设备存储器的其它区域间复制（用 cudaMemcpy22DO和cudaMemcpy.3DO函数)时，保证了最佳性能。返回的步长 (pitch,.stride)必须用于访问数组元素。下面的代码分配了一个尺寸为 width*height的二维浮点数组，同时演示了怎样在设备代码中遍历数组元素

CUDA 编程指南 4.0 中文版 20 float* h_B = (float*)malloc(size); // Initialize input vectors ... // Allocate vectors in device memory float* d_A; cudaMalloc((void**)&d_A, size); float* d_B; cudaMalloc((void**)&d_B, size); float* d_C; cudaMalloc((void**)&d_C, size); // Copy vectors from host memory to device memory cudaMemcpy(d_A, h_A, size, cudaMemcpyHostToDevice); cudaMemcpy(d_B, h_B, size, cudaMemcpyHostToDevice); // Invoke kernel int threadsPerBlock = 256; int blocksPerGrid = (N + threadsPerBlock – 1) / threadsPerBlock; VecAdd<<<blocksPerGrid, threadsPerBlock>>>(d_A, d_B, d_C, N); // Copy result from device memory to host memory // h_C contains the result in host memory cudaMemcpy(h_C, d_C, size, cudaMemcpyDeviceToHost); // Free device memory cudaFree(d_A); cudaFree(d_B); cudaFree(d_C); // Free host memory ... } 线性存储器也可以通过 cudaMallocPitch()和 cudaMalloc3D()分配。在分配 2D 或 3D 数组的时候，推荐使用，因为这些分配增加了合适的填充以满足 5.3.2.1 节 的对齐要求，在按行访问时或者在二维数组和设备存储器的其它区域间复制（用 cudaMemcpy2D()和 cudaMemcpy3D()函数）时，保证了最佳性能。返回的步长 （ pitch,stride) 必 须 用 于 访 问 数 组 元 素 。 下 面 的 代 码 分 配 了 一 个 尺 寸 为 width*height 的二维浮点数组，同时演示了怎样在设备代码中遍历数组元素

第三章编程接口 21 ∥Host code int width=64,height=64; float*devPtr; int pitch; cudaMallocPitch((void**)&devPtr,&pitch,width sizeof(float),height); MyKernel<<<100,512>>>(devPtr,pitch,width,height); I∥Device code global void MyKernel(float*devPtr,int pitch,int width,int height){ for (int r=0;r<height;++r){ float*row =(float*)((char*)devPtr +r pitch); for (int c=0;c<width;++c){ float element row[c]; } 下面的代码分配了一个尺寸为width*height*depth的三维浮点数组，同时演 示了怎样在设备代码中遍历数组元素。 I∥Host code cudaPitchedPtr devPitchedPtr; cudaExtent extent=make_cudaExtent(64,64,64); cudaMalloc3D(&devPitchedPtr,extent); MyKernel<<<100,512>>>(devPitchedPtr,extent); ∥Device code global void MyKernel(cudaPitchedPtr devPitchedPtr,cudaExtent extent){ char*devPtr=devPitchedPtr.ptr; size_t pitch=devPitchedPtr.pitch; size t slicePitch=pitch extent.height; for (int z=0;z<extent.depth;++z){ char*slice=devPtr+z*slicePitch; for (int y=0;y<extent.height,++y){ float*row=(float*)(slice +y pitch), for (int x=0;x<extent.width;++x){ float element rowx|; 参考手册列出了在cudaMalloc()分配的线性存储器，cudaMallocPitch()或 cudaMalloc.3DO分配的线性存储器，CUDA数组和为声明在全局存储器和常量存 储器空间分配的存储器之间拷贝的所有各种函数。 下面的例子代码复制了一些主机存储器数组到常量存储器中：

第三章 编程接口 21 下面的代码分配了一个尺寸为 width*height*depth 的三维浮点数组，同时演 示了怎样在设备代码中遍历数组元素。 // Host code cudaPitchedPtr devPitchedPtr; cudaExtent extent = make_cudaExtent(64, 64, 64); cudaMalloc3D(&devPitchedPtr, extent); MyKernel<<<100, 512>>>(devPitchedPtr, extent); // Device code __global__ void MyKernel(cudaPitchedPtr devPitchedPtr, cudaExtent extent) { char* devPtr = devPitchedPtr.ptr; size_t pitch = devPitchedPtr.pitch; size_t slicePitch = pitch * extent.height; for (int z = 0; z < extent.depth; ++z) { char* slice = devPtr + z * slicePitch; for (int y = 0; y < extent.height; ++y) { float* row = (float*)(slice + y * pitch); for (int x = 0; x < extent.width; ++x) { float element = row[x]; } } } 参考手册列出了在 cudaMalloc()分配的线性存储器，cudaMallocPitch()或 cudaMalloc3D()分配的线性存储器，CUDA 数组和为声明在全局存储器和常量存 储器空间分配的存储器之间拷贝的所有各种函数。 下面的例子代码复制了一些主机存储器数组到常量存储器中： // Host code int width = 64, height = 64; float* devPtr; int pitch; cudaMallocPitch((void**)&devPtr, &pitch, width * sizeof(float), height); MyKernel<<<100, 512>>>(devPtr, pitch, width, height); // Device code __global__ void MyKernel(float* devPtr, int pitch, int width, int height){ for (int r = 0; r < height; ++r) { float* row = (float*)((char*)devPtr + r * pitch); for (int c = 0; c < width; ++c) { float element = row[c]; } } }

CUDA编程指南4.0中文版 22 constant float constData[256]; float data[256]; cudaMemcpyToSymbol(constData,data,sizeof(data)); cudaMemcpyFromSymbol(data,constData,sizeof(data)); device float devData; float value =3.14f, cudaMemcpyToSymbol(devData,&value,sizeof(float)); _devicefloat*devPointer; float*ptr; cudaMalloc(&ptr,256 sizeof(float)); cudaMemcpyToSymbol(devPointer,&ptr,sizeof(ptr)); 为声明在全局存储器空间的变量分配的存储器的地址，可以使用 cudaGetSymbolAddress()函数检索到。分配的存储器的尺寸可以通过 cudaGetSymbolSize()函数获得。 3.2.3共享存储器 共享存储器使用shared限定词分配，详见B.2节。 正如在2.2节提到的，共享存储器应当比全局存储器更快，详见5.3.2.3节。 任何用访问共享存储器取代访问全局存储器的机会应当被发掘，如下面的矩阵相 乘例子展示的那样。 下面的代码是矩阵相乘的一个直接的实现，没有利用到共享存储器。每个线 程读入A的一行和B的一列，然后计算C中对应的元素，如图31所示。这样， A读了B.width次，B读了A.height次。 /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

CUDA 编程指南 4.0 中文版 22 为声明在全局存储器空间的变量分配的存储器的地址，可以使用 cudaGetSymbolAddress() 函 数 检 索 到 。 分 配 的 存 储 器 的 尺 寸 可 以 通 过 cudaGetSymbolSize()函数获得。 3.2.3 共享存储器 共享存储器使用__shared__限定词分配，详见 B.2 节。 正如在 2.2 节提到的，共享存储器应当比全局存储器更快，详见 5.3.2.3 节。 任何用访问共享存储器取代访问全局存储器的机会应当被发掘，如下面的矩阵相 乘例子展示的那样。 下面的代码是矩阵相乘的一个直接的实现，没有利用到共享存储器。每个线 程读入 A 的一行和 B 的一列，然后计算 C 中对应的元素，如图 3-1 所示。这样， A 读了 B.width 次，B 读了 A.height 次。 // 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 __constant__ float constData[256]; float data[256]; cudaMemcpyToSymbol(constData, data, sizeof(data)); cudaMemcpyFromSymbol(data, constData, sizeof(data)); __device__ float devData; float value = 3.14f; cudaMemcpyToSymbol(devData, &value, sizeof(float)); __device__ float* devPointer; float* ptr; cudaMalloc(&ptr, 256 * sizeof(float)); cudaMemcpyToSymbol(devPointer, &ptr, sizeof(ptr));

第三章编程接口 23 Matrix d_A; d A.width=A.width;d A.height=A.height; size t size =A.width A.height sizeof(float); cudaMalloc((void**)&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((void**)&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((void**)&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 blockldx.y blockDim.y threadldx.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;

第三章 编程接口 23 Matrix d_A; d_A.width = A.width; d_A.height = A.height; size_t size = A.width * A.height * sizeof(float); cudaMalloc((void**)&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((void**)&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((void**)&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; }