3. CUDA简介

　　nVidia GPU使用了CUDA编程模型，对硬件进行抽象，包括三个基本概念：线程组的层次，共享内存，同步。提供了细粒度数据并行和线程并行，循环的粗粒度数据并行和任务并行。CUDA有可扩展的多线程流处理器(SMs)阵列构成。当执行CUDA的kernel时，grid的块就被分布到多处理器上。一个多处理器由8个标量处理器(SP)。

　　CUDA线程模型：GPU上的thread有三个层次。最高层是grid，grid内的thread执行相同的kernel，每个grid可以有2^16-1个blocks， blocks组织成一维或者二维形式。中间层是block，每个block内部的thread可以组织成三维，最多可以有512个threads。同一block内部的thread可以通过一个低延迟的片上共享内存共享数据，对数据进行原子操作，并通过__syncthreads原语进行同步。最后，硬件将thread分为warp执行以提高访存性能。

　　CUDA kernel函数是C函数，在被调用时按照指定的grid维数和线程块维数并行的执行。CUDA可以通过内置的threadIdx变量获得线程ID。kernel的限制：不能包含循环，静态变量，参数的个数必须固定。threadIdx是一个三元矢量，所以可以将线程按照一维，二维或者三维来划分。在kernel内可以通过内置的blockIdx变量和blockDim变量访问某个线程块。同一网格中不同线程块中的线程不能互相通信和同步。Block内的线程通过shared memory，atomic operations和barrier synchronization协同工作，不同block内的线程不能通信。Kernel启动时内置变量和函数参数存放在shared memory中。

　　每个线程有私有本地内存，每个线程块有共享内存，所有线程可以访问全局内存。所有线程还可以访问二个只读内存空间，常量和纹理内存空间。Kernel内的数组存储在local memory中。Shared memory，G80上为16个bank，4bytes为单位进行寻址，bank ID=4-byte address%16，相邻的4-byte地址映射为相邻的bank，每一个bank的带宽为4bytes per clock cycle。对同一个bank的同时访问导致bank conflict，只能顺序处理。解决方法，padding，transpose。Coalesced global memory accesses，在half-warp层对访问global memory进行协调，

　　CUDA将32个标量线程称为warp，以warp为单位来进行创建，管理，调度，执行。一个warp共享和执行相同的指令，由于每个SM有8个core，执行一个warp指令需要4个cycles，类似于一组矢量指令的流，所以标量线程可以看成矢量处理单元。一个cycle一条指令从L1指令cache加载到指令buffer中，当warp的数据都可用时被选择执行。Block内的thread按顺序分配到不同的warp中，但是warp的顺序可能随着GPU的更新而变化。GPU调度线程零开销。执行结构为单指令多线程SIMT，SM将线程块中的每一个线程映射到一个标量处理器，每一个标量线程独立的在自己的指令空间和寄存器状态运行。SIMT与SIMD的不同指出在于后者指定了数据宽度，然而SIMT中的每一个线程可以执行不同的代码路径。SIMT可以使得程序员编写指令级并行代码，

　　主要的编译制导语句：

　　l __device__均为inline函数，__noinline__指示编译器不要inline。

　　l #pragma unroll n告知编译器展开循环5次，#pragma unroll 1禁止循环展开，如果不加数字，#pragma unroll，如果循环计数为常数，则全部展开，否则不展开。

　　l __restrict__告诉编译器指针没有重名，即写入一个指针不会影响其他指针指向的数据，所有的指针参数都要指定__restrict__，因此编译器可以使用指令重排和公共子表达式消除优化方法。这样可能增加了register需求数量，因此需要平衡。

　　l 编译时可以配置global内存为L1和L2都 (-Xptxas -dlcm=ca)，或者只有L2进行缓存(-Xptxas -dlcm=cg)。

　　l __align__(n)，在n字节位置对齐。例如struct __align__(16){float x,y,z;}，主机和设备端内存分配函数返回的首地址均在256字节对齐。

　　主要的编译选项：

　　l 加入--ptxas-options=-v编译选项可以看各存储空间分配大小。

　　l 加入-cubin生成.cubin

　　l 加入-ptx生成. Ptx文件。

　　l 加入--keep，保留中间结果，中间文件。

　　4. CUDA优化方法

　　4.1 从硬件角度观察

　　共享内存作用：

　　l 同一线程块内线程共享数据，减少内存带宽需求，增加计算访存比。例如，不利用共享内存的矩阵乘代码，假设每个线程计算一个矩阵C元素，内存带宽需求为2n2+2n3，分块矩阵乘代码[PPOPP08，CUBLAS1.1]在共享内存中缓存矩阵A和B大小为16*16子块，内存带宽需求降为2n2+2n3/B。

　　l 可以缓存内存全局数据，优化warp内线程内存访问，使得coalesced加载内存数据。例如考虑矩阵转置，计算访存比为1，不存在数据重用，似乎利用不到共享内存，如果要求源矩阵A和目标矩阵B均为行主或列主存储，每个线程处理一个元素，则无论怎样组织一个warp内的线程，对其中一个矩阵中元素一定是非coalesced访问，如果利用共享内存，一个线程块处理矩阵一个子块，则可以进行优化使得所有内存访问均为coalesced，且共享内存没有bank冲突。

　　共享内存优化方法，防止bank冲突：

　　l 数组填充。

　　l 矩阵转置。

　　l 共享数据大小要考虑同步作用，由于在共享内存的数据通常被线程块内所有部分线程访问，所以要在加载数据之后进行同步，共享数据太小会导致同步次数增多，降低性能。

　　l 利用共享内存数据通常要在kernel内循环利用不用数据，共享数据太大会导致循环次数增多，可以利用循环展开降低循环等待时间。

　　全局内存优化方法：合并访问，利用共享内存降低带宽需求。解决非coalesced内存访问方法：

　　l 使用structure of array(SOA)而不是AOS，

　　l 使用__align(x)对齐数据，

　　l 或者使用shared memory来协作访问。

　　l Coalescing float3 access说明如何使用shared memory进行Coalesced access。

　　寄存器作用：

　　l GPU的设计理念是大量并行多线程隐藏延迟，因此在一个SM上配置了大量寄存器。

　　l 开发ILP可以在寄存器重用数据，而TLP只能在共享内存重用数据。[SC08]GEMM对矩阵A和C元素在寄存器重用以开发ILP，对矩阵B子块在共享内存重用以开发TLP。ILP-TLP trade-off，寄存器比共享内存空间大，

　　4.2 从软件角度观察

　　由于GPU的设计理念是大量并行多线程隐藏延迟，因此在一个SM上配置了大量寄存器，由于SM还限制了thread block，warp和threads数目，所以通常利用循环展开增加寄存器利用率。Tesla C1060上每个SM最多容纳768个threads，而每个SM有64KB大小寄存器，所以每个线程可以使用85B，[SC08]中的GEMM算法每个线程使用30个寄存器，而CUBLAS 1.1使用15个寄存器。

　　GPU的目标是高带宽而不是低延迟。GPU的线程切换低开销，较高访存延迟和较小共享内存和寄存器空间都限制了GPU上应用以细粒度并行为主。GPU带宽比CPU高，可以通过线程间调度隐藏延迟，因此GPU优化主要在访存对齐以及利用共享内存，优化：最大化独立并行性，最大化算术计算密度，重复计算优于访问存储器

　　对于CUDA线程模型三个层次，[PPOPP10]进行了分析并提出了相应指标计算方法评估算法在每一层中特性：在线程级提出ILP，表示在寄存器一级共享数据以及隐藏延迟;在warp级开发DLP，保证全局内存coalesced访问，减少共享内存bank冲突;在线程块一级开发线程间数据共享。由于一个warp内线程同步运行，我们将GPU上运行的一个warp视为CPU中一个线程，则本文所指GPU线程ILP(Instruction-Level parallelism，指令级并行)为一个GPU上单线程内部指令级并行性，而GPU线程TLP(Thread-Level Parallelism，线程级并行性)指一个SM上活动warp间并行性。优化计算：循环展开减少动态指令数。[SC08]过多的线程并不能提高性能，增加线程内计算任务，提高ILP。

　　ILP vs.TLP：访存延迟为400cycles，一个warp在G200上需要4个cycles，如果每个n个指令有一个访存操作，则容忍延迟需要400/4n个warp。如果增加一个register能够增加独立指令个数，但是降低了block数量或者block内线程数量，即降低TLP，需要比较两者孰优孰劣。

　　4.3 Fermi基本优化方法

　　Global内存可以配置成L1和L2缓存，或者只有L2缓存，而local内存必须进行L1和L2缓存。如果kernel使用较多local内存，则设置为48KB L1 cache。计算能力1.x可以利用texture内存加速数据访问，然而由于2.0加入了带宽更高的L1 cache，使用texture cache并不一定能改进性能。

　　计算能力1.x的global内存访问以half-warp为单位，而2.0则以一个warp为单位。

　　计算能力1.x的shared内存访问以half-warp为单位，有16个bank，而2.0则以一个warp为单位，有32个bank，依然以32bit为单位存储，即连续32bit存储在连续bank，每个bank的带宽为32bit/(2 cycles)。只有在两个或更多线程同时访存一个bank内不同32bit字时产生bank conflict，多个线程访问同一数据不会产生bank conflict，而且2.0增加多播使得多个字可以在一个transaction内完成，多个线程可以访问同一个32bit字内的不同字节。对于64bit double寻址，只要一个half-warp内线程不访问同一bank内的不同32bit数据，就不存在bank conflict，因此double data = shared[BaseIndex + tid]没有bank conflict。

　　l 计算能力1.x中整数乘以mul24为基本实现，而2.0则32bit实现整数乘。

　　l 浮点运算配置成-ftz=true (denormalized numbers are flushed to zero)，-prec-div=false (less precise division), and -prec-sqrt=false (less precise square root)会生成性能更高的代码。

　　l nvcc设置–m64或者–m32将生成64位和32位代码，如果GPU不用更多内存，可以分别编译设备代码和主机代码。

　　l 标准函数可以在主机或者设备端调用，intrinsic函数只运行在设备端调用，并且具有更高性能，准确率降低。intrinsic函数以__为前缀。

　　l 2.0提供了一些高级intrinsic函数。

　　为了隐藏延迟为L的指令，1.x需要L/4个warp，而2.0需要L/2个warp，因此2.0上需要更多warp隐藏延迟。

　　由于2.0支持并发kernel，可以将并行度较低或者串行部分代码交由GPU执行，以减少主机和设备端数据交换，对于未来主机和设备融合共享内存，这种优化则不必要。