很久没有写最近学习的一些内容了，有些小忙，也因为业余时间活动安排地太满了，时间不足。其实写了很多笔记，但是规划得不是很工整，零零散散，只有自己看得懂，就不发出来了hhh。
最近因为一个project的需要，连续肝了好几周的cuda代码，把CPU的代码转到GPU上去实现。目前结果也挺好，即使是在我笔记本辣鸡的GTX1050Ti里也提高了约120倍的速度。当然也得益于我优秀的"设计"hhhh不自夸了。
一维的数据
在我目前写的Cuda代码中，我把所有数据都一维化，因为嫌弃在Cuda里处理高维数据时对齐指针很麻烦。如下一个简单的例子(具体内容需要参考于[1])，如果要在Cuda里实现一个二维数组的相加看起来像下面
C[idy][idx] = A[idy][idx] + B[idy][idx];
那么在配置时，host(CPU)端你需要做的是设定一个二维指针并分配空间
int **A = (int **)malloc(sizeof(int*) * Row);
int **B = (int **)malloc(sizeof(int*) * Row);
int **C = (int **)malloc(sizeof(int*) * Row);
int *dataA = (int *)malloc(sizeof(int) * Row * Col);
int *dataB = (int *)malloc(sizeof(int) * Row * Col);
int *dataC = (int *)malloc(sizeof(int) * Row * Col);
device(GPU)端设定二维指针并分配空间
cudaMalloc((void**)&d_A, sizeof(int **) * Row);
cudaMalloc((void**)&d_B, sizeof(int **) * Row);
cudaMalloc((void**)&d_C, sizeof(int **) * Row);
cudaMalloc((void**)&d_dataA, sizeof(int) *Row*Col);
cudaMalloc((void**)&d_dataB, sizeof(int) *Row*Col);
cudaMalloc((void**)&d_dataC, sizeof(int) *Row*Col);
注意dataA/d_dataA是CPU/GPU实际储存数据的变量。而A,B,C/d_A,d_B,d_C是储存一个二维矩阵的每一行第一个元素的变量。在原博客中，提取了d_data的每一行的首地址，赋值给了A,B,C
for (int i = 0; i < Row; i++) {
A[i] = d_dataA + Col * i;
B[i] = d_dataB + Col * i;
C[i] = d_dataC + Col * i;
}
最后再把数据从host拷贝到device。
cudaMemcpy(d_A, A, sizeof(int*) * Row, cudaMemcpyHostToDevice);
cudaMemcpy(d_B, B, sizeof(int*) * Row, cudaMemcpyHostToDevice);
cudaMemcpy(d_C, C, sizeof(int*) * Row, cudaMemcpyHostToDevice);
cudaMemcpy(d_dataA, dataA, sizeof(int) * Row * Col, cudaMemcpyHostToDevice);
cudaMemcpy(d_dataB, dataB, sizeof(int) * Row * Col, cudaMemcpyHostToDevice);
上面比较重要的一步就是把从机端的data地址赋值给A,B,C，这样A[0]就代表了第一行的地址A[0][0]就代表了一行第一列的数据了。在核函数里就可以比较直观地相加了。
C[idy][idx] = A[idy][idx] + B[idy][idx];
对于高维数据，为了在核函数直观地相加，在Cuda里我们通常需要人为地先对齐数据。二维可能比较好对齐，但更高维可能就比较烧脑袋。
如果真的使用过一些计算机视觉库(当然不止)的底层的数据的话，你会发现其实底层的数据都是一维储存的。比如Opencv 的cv mat。
cv::Mat m = cv::Mat::zeros(7, 7, CV_8UC1);;
char* ptr = m.data;
成员变量data返回数据指针，上面的ptr[0],ptr[1]代表第一行第一列，第一行第二列的数据，ptr[7*k+n]代表第k行第n列的数据。
广为使用的C++线性代数库Eigen也是如此。
Eigen::Matrix3f m;
m << 1, 2, 3,
4, 5, 6,
7, 8, 9;
float* ptr = m.data();
成员函数data()返回数据指针，同样ptr[0],ptr[1]代表第一行第一列，第一行第二列的数据，ptr[7*k+n]代表第k行第n列的数据。
c++的标准库std vector同样如此，不管几维的向量，vector.data()提供储存的一维的数据。
可见我们一直都在处理一维的数据，当然这也不奇怪，毕竟数据在RAM里储存的方式就是按着地址从小到大排列，肯定是一维的。只是当上面的库为我们提供了很方便的高维接口时我们忽略了这点。这也让我对从最基本的数据来操作计算产生了兴趣，又加上cudas里不能直接使用Eigen,CV Mat等，我只能获取对应变量的一维数据在核函数中使用，如果自定义结构体/类，需要手动对齐数据。在我的project里，有无数的Eigen, Cv Mat, vector变量，所以我最终决定，直接在他们的底层一维数据上操作吧。比如下面是我一个Cuda代码核函数的参数
(bool calculate_der, double* im0, double* im1, double* points3d, double* bs_value_ref, int*bs_index_ref, double* pose, double* in, int thread_work, int bin_num, int bs_degree, int rows, int cols, int cell, double* d_d_sum_bs_pose, double* d_d_sum_joint_bs_pose, double* d_pro_target, double* d_pro_joint)
当然你肯定不感兴趣他们是什么，不过你大概可以看到，只有基本的数据类型，他们全是从eigen, CV Mat, vector等中提取出来的。我也没有定义一些复杂的类甚至结构体了，一切从最原始的数据指针的开始操作。当然我上面的操作仅限于参数比较少的，几个十几个，如果有几十个参数可能就麻烦了，还是得把类似的参数归纳到自定义的结构体或者类里去。数据都一维化的操作使得我处理数据更直接，但是代价是代码的可读性变差了不少。比如一个三维的变量，
//三维变量var，有m行n列k层，如今我们获取了它的一维数据，尝试获取第五行第六列第七层
var[6*m*n+4*n+5];
每一层有m行n列，所以第７层开始的数据地址是6mn，对于该层来讲，4*n是第５行数据的其实位置，加５得到第六列的数据。可能这个还比较直观，但一旦你要面对的不单单是m行，可能是m+k+s行的数据，上面代码的获取方式还是挺麻烦的。
虽然以后对于高维数据我可能还是会设计出一个strcut来获取特定行/列的元素，但是这次的project针对一维数据的操作确实让我获益匪浅，就算最基础的c风格的代码方法也能获得可观的效果。
Cuda函数耗时
cuda函数的调用是比较耗时的，我曾经做过实验，在我的机子上调用10000次cudaMalloc花费了大概30ms的时间，也就是１次大概3us。cudaFree()同样会消耗几微秒的时间。之所以做一个实验，是因为我当初有部分代码(我以为)需要大量调用相关函数，当时有点懵，因为希望代码能100ms内出个结果的，结果没想到还没拷贝数据做任何计算就花费了10几毫秒。如果数据量大，从host到device的数据拷贝本来就很耗时了，这就很麻烦。所以在设计Cuda代码时，将数据归类，尽量使用少的次数把数据分配/复制完成。
递归
递归其实是我不太喜欢的编程方式，栈溢出先不说，一旦函数内容繁杂了之后，出了问题很不好debug。所以其实我没有在自己项目里写过递归代码，虽然刷题(也没刷几道)感觉大家很喜欢用的样子= =。但好巧不巧，这次代码有一小部分需要使用别人写好的递归函数，并且是在GPU里，第一次移植进去之后没什么大问题，程序正常运行，但是后来cuda程序时不时出现"illegal memory access"。后来经过查找，发现是在Cuda中，能分配给一个线程的栈空间是很有限的，至少在我的1050Ti中，运行下列代码
size_t limit = 0;
cudaDeviceGetLimit(&limit, cudaLimitStackSize);
printf("cudaLimitStackSize: %u\n", (unsigned)limit);
print出来的结果是1024 byte。这样如果栈稍微深层一点，就无法继续了。于是我手动把栈限制提高到2048byte，这样我的程序基本都能运行了。
size_t limit = 2048;
cudaDeviceSetLimit(cudaLimitStackSize, limit);
我并不太清楚能手动设置的上限是多少，但肯定不大，不然不会默认给吝啬的1024byte了。但由此可见，由于栈尺寸的限制，在Cuda里入栈需要更加的小心谨慎了。在CUDA论坛里逛时很多人不建议在CUDA中使用递归，或者只使用很简单的递归。
原子操作
这一次也使用了一定次数的GPU内的原子操作，为了多个线程同时修改一个变量时不起冲突。比如一千个线程同时执行下面操作
__global__ void Add(double* A){
A[0] += 1;
}
如果原本A[0]是0，那么结果不会是１０００，因为很多线程会同时基于原来的０加１得到１，又有不少线程会基于原来的数加１得到２．最终结果会小于１０００。比较值得一提的是，不像CPU里的多线程，同时读写一个数，如果没有mutex会报错，由于GPU本身的结构设计，它是不会报错的，计算正常进行，只是结果不是１０００而已。
这时如果我们想得到１０００，我们需要使用atomicAdd
__global__ void Add(double* A){
atomicAdd(&A[0], 1);
}
可能很多同学也熟悉了c++里的原子操作。只是比较令我惊讶的是，这个操作几乎没有让我的程序增加时间消耗。
最后atomicAdd要能对double变量进行操作的话，需要GPU结构sm_60以上。CMake里可以输入下列语句来帮助编译
set(CUDA_NVCC_FLAGS -arch=compute_60)
另外附上一张表格方便查看自己GPU的版本[2]，基本上1000+系列的显卡都是支持sm_60的。
本次笔记写地匆忙，如有不足之处请指出并交流
[1] https://www.cnblogs.com/skyfsm/p/9673960.html
[2]https://arnon.dk/matching-sm-architectures-arch-and-gencode-for-various-nvidia-cards/
