列规约

加速计算 CUDA
1 
__global__ void SAD(int *A,int *B,int *C)
{
	int stride;
	int tx=blockDim.x*blockIdx.x+threadIdx.x;
	int ty=blockDim.y*blockIdx.y+threadIdx.y;
	
		C[ty*COL+tx]=fabs((float)(A[ty*COL+tx]-B[ty*COL+tx]));
	
	__syncthreads();
	
		for(stride=ROW/2;stride>0;stride=stride>>1){
			__syncthreads();
			if (ty < stride)
   C[ty*COL+tx] += C[(ty+stride)*COL+tx];
		}
	
	
}

代码如下，其实我就想实现的就是列规约的，但是结果有时对，有时对不的。我猜想原因是不是多个线程同时读写了，可是自己看了半天就是没看出来。求万能的版主大人帮帮我。谢谢。

2

LZ您好，非万能的斑竹来回答一下您的问题：

由于您没有给出部分参数的确切定义，先大致推断如下：
COL表示列数，而不是列内元素数；ROW表示行数，而不是行内元素数。比如一个3行5列的数组/矩阵，这里ROW=3，COL=5。并且COL和ROW是宏或者__device__全局变量等，无需参数传递就能在kernel函数中正常使用。
假设您*A,*B,*C对应的分配空间大小是合适的，以及正常初始化过。
假设您的数据是行优先存储的。
以及您没有给出您kernel调用时的形状，这个在下面将继续讨论。

再次基础上，通过目测您的代码，我发现问题如下：

1：您定义tx，ty之后，直接使用了tx，ty进行寻址，这里可能有隐藏的问题。如果您的线程发布形状并非和您的数组一一对应，那么需要if判断将超出部分的线程提前return，（即将tx>=COL,ty>=ROW的线程提前return;）否则寻址会出现问题。

2：您在生成C数组内容之后，使用了__syncthreads()同步，这个同步只能对block内部进行同步的，而您的kernel很可能是多block的，这样各个block之间并不保证同步，而您在后面的规约操作中并非是一个block只使用自己生成的数据，所以这里是有问题的。
简而言之，这里并无法保证C数组完整地生成了。

同理，您下面循环中的__syncthreads()也是同样的问题。

3：您的for循环，stride=ROW/2。这里ROW如果是奇数，那么最后规约的次数可能和您预期的次数是不相符的。我们比较ROW=2和ROW=3两种情况，stride均为1，那么均只运行下面的循环一次。但是对于ROW=3的情况，这个行为和预期的行为并不一致。以及，实际上ROW只要含有奇数因子，都会有类似的问题，而只有在ROW为2的幂次的时候，这样做才是完全无问题的。因为您在1#并未给出如此保证，所以这里提醒您注意这一点。

您代码逻辑上的问题，目测结果如上，供您参考。

祝您编码顺利~

3

多谢ice的解答。有些地方我是没写清楚的。看了您的回答，我想应该是第二个问题。多谢您的解答！

4

LZ您好，不客气的。

其实规约的实现有很多方法的，并不局限于全部铺上线程。
您也可以用中等数量的线程数（即并不和原数据点数一对一对应，但数量足够GPU的各SM跑满），每个线程多干一些活。

假如您有M行N列的数据需要进行列规约，以及这里M,N都较大。
那么您可以使用一维的block，每个block的线程数为N/2或者N的几分之一这样。
每个线程在列方向上完成M个元素的几分之一。
也就是说利用一个一维的block完成一个二维的块的规约。

此时每个线程只需要循环叠加自己这一列的数据即可，只需要从global里面读取一次，并累加到自己的寄存器变量里面，而无需反复读写global memory。

每个block计算完毕的时候，可以向global memory的结果矩阵（1行N列）的对应位置使用atomicAdd累加本block的结果即可。（也可以先直接保存于global memory某处，最后再对block结果进行一次小规模的规约操作。）

或者如果按照您之前大规模铺线程的思路，那么您每个2维的block在规约的时候，只对自己这个block的数据规约，得到一个局部数据，最后再合并，这样也是可以的。

因为一般来说除了结束kernel没有其他的全局同步的手段，所以您之前的写法是不可以的。

大致如此，祝您好运~