GPGPU

2个矩阵乘法的例子：

1. 分块：每个block负责目标矩阵中的一块，好处：（读显存的数据量：计算次数=1：小块的边长）；如果直接每个thread负责1个目标值，读显存数据量：计算次数=1：1，而且读column显存时可能无法连续读显存；（？？？好像不对）

2. 分块矩阵把数据读到了shared-memory，充分利用高速缓存进行矩阵乘法；减少了访问显存的量，使得计算：访存是O(N^3) : O(N^3/K) ; 只要shared-memory能放的下，K越大越好；

 

1. A*X+Y的例子：https://devblogs.nvidia.com/easy-introduction-cuda-c-and-c/

注意：应对向量长度不是block大小的整数倍，1. block个数向上取整 2. kernel里用if判断是否有活儿干

2. 卡时间：https://devblogs.nvidia.com/how-implement-performance-metrics-cuda-cc/

同步操作(cudaMemcpy)和异步操作(kernel)解释的很好；（cudaMemcpy是同步的（阻塞的），即上面的kernel执行完毕，cudaMemcpy才能开始，cudaMemcpy执行完毕，后面的kernel才能开始；kernel是异步的，即发起调用后，控制立即返回CPU端继续往下执行）

CPU端计时：调cudaDeviceSynchronize()等kernel执行完再打点。缺点：等待期间CPU啥也干不了，浪费了；

GPU端计时：event, 会记下这个stream执行到这里的时间点；

global-mem带宽计算：1. 理论带宽（同时读写）是查硬件手册计算得到的，MHZ(每秒多少次传输，memoryClockRate）* 内存接口宽度（384bit，memoryBusWidth) * 2(DDR RAM可以读写同时）2. 实际带宽，跑程序卡时间得到（也要把读的和写的都算进去）；（2012年的GPU，理论双向带宽148GB/s)

GFLOPS是针对float数据计算的，double数据的性能一般是GFLOPS的一半。A*X+Y看成是2个Float计算；这个例子读写12个byte才做2个float计算，显然是memory-bound的，不是compute-bound。复杂程序性能优化建议用工具来profiling，看看bottleneck在哪儿。

3. 查GPU属性，错误处理：https://mp.csdn.net/postedit/88734169

通过调API来获取GPU参数来计算memory带宽的例子；

nvcc编译选项可以指定编译成在X.X版本上跑；

每个block多少threads，很重要：太少，则每个SM上的block数（硬件规定有上限的）填不满SM；太多，不能超过block上限thread总数（报错）；

同步操作，会返回错误码；kernel启动等异步操作，在出错的时候系统会往系统的变量里写入错误码，host可以同步一下（cudaDeviceSynchronize），然后调GetLastError得到错误码（于此同时错误码被清除了）；kernel之后立马调的GetLastError只能返回launch相关的错误，不同步完的话无法返回执行期间的错误。cudaDeviceSynchronize拖慢速度所以尽量只在debug的时候用，release以后别用；

 

 

1. 优化数据传输：https://devblogs.nvidia.com/how-optimize-data-transfers-cuda-cc/

host<-->device的PCI带宽8GB/s; device memory<-->GPU的带宽144GB/s; 相差一个数量级；

优化原则：1. 尽量少在host和device之间传输数据；2.使用page-locked(pinned)主存（双刃剑）；3.多次小传输打包成一次大传输，可减少每次传输的额外开销；4.该数据传输可以和kernel执行或者其他数据传输并行起来（stream)

卡时间的方法：1. cudaEventElapsedTime()； 2. nvprof工具查看数据传输耗时

对某个任务来说，GPU和CPU谁快，应该把GPU和device的数据传输时间也计算进去，否则不准。

host<-->device之间copy数据，CUDA无论如何会使用host的page-locked(pinned)主存，如果用户开辟的不是这个则CUDA自动开辟一份先做一次copy，如果用户开辟的是这个则省了这一步，会有百分之小几十的速度提升。

开辟pinned内存，把小传输打包到这个大pinned内存里，再一次性传输；2D/3D的数据调用cudaMemcpy2D()/cudaMemcpy3D()会更快

2.把默认stream和"stream们各个操作提交顺序造成的性能差异"讲的很透彻的文章：https://devblogs.nvidia.com/how-overlap-data-transfers-cuda-cc/   

默认stream(NULL)就是他要执行必须之前提交的所有stream的操作都完成，他执行的时候后面所有stream提交的操作都不能开始；

对于host而言，cudaMemcpy是阻塞的，kernel提交是非阻塞的（可以和接下来的CPU计算并行），对device而言，kernel执行也是阻塞的，他执行完，后面的操作才能开始（同stream内的）；

想copy和计算并行，必须用cudaMemcpyAsync版本，且pinned host内存；

提交完一个stream的3个操作再提交下一个stream的3个操作那个例子，讲的很清楚。两种情况：copy和kernel各有1个engine对应1个队列，H2D和D2H和kernel各有1个engine对应1个队列；按提交顺序queue到不同的engine队列里（必须按先入先出来执行），且满足同一个stream内的顺序；Hyper-Q之后这些trick没用了。

3. global memory:  https://devblogs.nvidia.com/how-access-global-memory-efficiently-cuda-c-kernels/

用__device__定义变量，或者用cudaMalloc()开辟动态数组；

每个硬件warp上是32个threads；warp是SIMD执行（齐步走）；早期1代显卡要求访问global memory要对齐且按thread顺序访问，不这样的话就有1个数量级的带宽损失（分成多次访存了）；2代显卡之后有cache，凑齐了才去访存，所以没有什么带宽损失

row-major的矩阵，如果这批threads要访问矩阵某一列，涉及到跨stride访问，造成带宽严重下降（GPU没办法将这些访问请求合并到一起）；data-locality很重要！解决方法：1. 按column-major存放矩阵；2. 使用shared-memory(没有stride降带宽这说）；

4. shared memory: https://devblogs.nvidia.com/using-shared-memory-cuda-cc/

bank处理好了，shared memory延迟比global memory好100倍；

同一份shared memory可以被同一个block里所有线程访问；

block里的线程们不一定并行执行（不同warp时），先写后读的例子写的快的线程读写的慢的线程的目标数据会race, 要用__syncthreads()来barrier; 注意：有分支的时候要确保所有threads都调了，否则会死锁！

分配语法：1. 编译之前知道大小，直接在kernel里定义数组即可；2. 不知道大小，则在kernel<<<>>>里告诉大小，在kernel里声明extern __shared int s[]即可，当多个数组用的话手动设置那些数组的指针即可；

划分成32个bank（同warp大小）; 每个bank带宽是4B/cycle; 3.0之后可以配置bank为8B; bank之间的地址是连续的，threads访问连续shared-memory可以跨越多个bank，从而达到最佳带宽。多个threads同时访问同一个bank，会被拆成多次请求（访问同一地址例外，广播和多播确保其一次完成）

64KB/SM的shared memory，可以被user配置成多种L1-cache和shared-memory比例；如果系统需要更多shared-memory则会强制调整。

5. 矩阵转置例子：https://devblogs.nvidia.com/efficient-matrix-transpose-cuda-cc/

用矩阵copy来做对比，来看出来各种不同的转置实现有多少潜力。每个sub-matrix由每个线程搬4轮，可以将kernel里计算index的开销"均摊薄"

全用global-memory来转置，比copy慢了一两倍。原因：写入的时候threads没有“连续”访存（典型的跨stride访问）

使用32*32的shared-memory来中转，读的时候正常读入（读一行写一行），写的时候读（读shared-mem)一列写(写global-mem)一行, 使得对global-mem的访问都是“连续”的了；速度还是比copy慢了百分之大几十

为了调查是不是__syncthreads()造成的拖累，对矩阵copy也用shared-mem来实现并用上__syncthreads(), 发现并没有慢，说明不是这个的原因；

真正原因：threads读shared-memory（一列）的时候，访问的是同一个bank; 解决：shared-mem开辟[32][32+1]，从而让同一列的相邻元素都分布在不同的bank; 完美达到copy例子的带宽！

总结：两个有效提速：global-memory的"连续"访问，shared-memory的bank散开访问，有效

 

6. CUDA7 Thread default stream: https://devblogs.nvidia.com/gpu-pro-tip-cuda-7-streams-simplify-concurrency/

两个例子：1. Thread default stream上执行kernel，不会等之前提交的所有stream上操作完成，也不会阻塞后续提交的所有stream上操作开始；2. 每个thread在自己的default stream上提交任务，可以实现多线程的stream并行执行的效果；

 

7. Unified Memory: https://devblogs.nvidia.com/unified-memory-cuda-beginners/

老卡：cudaMallocManaged()调用的时候就在GPU上分配存储了；CPU这边访问内存时，会触发page fault，从GPU上一页一页搬到host上，CPU改了内存后，GPU上启动kernel时会先自动把全部这部分内存搬到GPU里，再运行kernel； nvprof显示的kernel耗时不包括之前的CPU->GPU内存搬家耗时；

新卡：lazy策略，cudaMallocManaged()调用时也许并不在GPU上分配存储；用的时候才分配；GPU kernel启动后用到的时候触发page fault一页一页吧数据从CPU搬到GPU; nvprof显示的kernel耗时就高了；

该怎么做：CPU的活儿让GPU来干，从而减少CPU干完的东西搬到GPU的开销；尽早prefetch即将访问的数据；

老卡上CPU和GPU不能同时访问同一地址（会报错）；新卡支持page fault，因此可以同时访问同一地址（用户自己确保数据一致性）; 新卡支持Unified Memory上的原子操作；
另一篇好文章：https://devblogs.nvidia.com/beyond-gpu-memory-limits-unified-memory-pascal/

 

8. Cooperative Groups: https://devblogs.nvidia.com/cooperative-groups/

A.折半reduce的例子（每一步都要sync)（很像MPI_reduce)

B.使用上面的例子(reduce_sum)来完成一个数组的加和：每个thread加和一部分数组(stride-style)得到1个数，block内部用折半reduce来把threads个数加成1个数，用原子加操作把每个block的数加到一起；

C.thread group：可以是所在的block，也可使是自定义拆分的；

D. 把B例子的block的thread-group替换成手工切分的32threads的thread-group, 改动很小；group内部的sync必须该group的所有threads都参与！所以有分支的时候，尽量拆分group，避免后续代码的sync造成死锁；

E. 把线程个数放到模板参数里，可以在编译时告诉编译器thread-group大小，从而让编译器去有机会优化代码（比如循环展开，32的时候利用warp来去掉同步）

F. 使用shuffle命令来加速warp（直接使用同warp内其他thread的寄存器数据，不用shared memory了）；B例子使用shfl_down和模板参数分thread-group来实现；https://blog.csdn.net/Bruce_0712/article/details/64926471

G. coalesced_threads可以拿到本warp内当前活跃的线程，编成一个group; 可以用于同步，可以用于选举rank0来做事；(里面的例子用到了这个：https://devblogs.nvidia.com/cuda-pro-tip-optimized-filtering-warp-aggregated-atomics/ ）（shfl可以用来做warp内广播）

 

9. copy_if()的实现：https://devblogs.nvidia.com/cuda-pro-tip-optimized-filtering-warp-aggregated-atomics/

用原子操作递增offset；很慢；

开一个shared memory变量（__shared__)，用原子操作递增该变量，每个线程能有自己的offset，每个block派一个thread去原子更新这段的offset；（代码很巧妙）；稍快一点儿；

用以下代码实现，比thrust还快半个数量级：（如果用shuffle指令集也可实现，只是代码多些）

__device__ int atomicAggInc(int *ctr) {
  auto g = coalesced_threads();
  int warp_res;
  if(g.thread_rank() == 0)
    warp_res = atomicAdd(ctr, g.size());
  return g.shfl(warp_res, 0) + g.thread_rank();
}

调用端：

  if(src[i] > 0)
    dst[atomicAggInc(nres)] = src[i];