c++ - 执行小插入/移位的并行算法

Question

假设我有一个包含 8 个数字的数组 A，我有另一个包含数字的数组 B 来确定 A 中的数字应该向右移动多少位

A 3, 6, 7, 8, 1, 2, 3, 5

B 0, 1, 0, 0, 0, 0, 0, 0

0 表示有效，1 表示这个数字应该在第 1 位之后，输出数组应该在 3 之后插入 0，输出数组 C 应该是:

C: 3,0,6,7,8,1,2,3

是插入0还是其他什么并不重要，关键是3之后的所有数字都移动了一位。出站号码将不再在数组中。

另一个例子:

A 3, 6, 7, 8, 1, 2, 3, 5

B 0, 1, 0, 0, 2, 0, 0, 0

C 3, 0, 6, 7, 8, 0, 1, 2

........................................

A 3, 6, 7, 8, 1, 2, 3, 5

B 0, 1, 0, 0, 1, 0, 0, 0

C 3, 0, 6, 7, 8, 1, 2, 3

我正在考虑使用 scan/prefix-sum 或类似的东西来解决这个问题。这个数组也很小，我应该能够将数组放入一个经线(<32 个数字)并使用随机播放指令。有人有想法吗？

Answer 1

一种可能的方法。

由于您的移位 (0, 1, 0, 1, 0, 1, 1, 1 和 0, 1, 0 ， 0 都产生相同的数据偏移模式，例如)不可能只创建移位模式的前缀和来在每个位置产生相对偏移。然而，我们可以观察到，如果移位模式中的每个零都被其左侧的第一个非零移位值替换，则将创建一个有效的偏移模式:

0, 1, 0, 0   (shift pattern)
0, 1, 1, 1   (offset pattern)

或

0, 2, 0, 2   (shift pattern)
0, 2, 2, 2   (offset pattern)

那么如何做到这一点呢？假设我们有第二个测试用例转换模式:

      0, 1, 0, 0, 2, 0, 0, 0

我们想要的偏移模式是:

      0, 1, 1, 1, 2, 2, 2, 2

1. 对于给定的移位模式，创建一个二进制值，如果移位模式中相应索引处的值为零，则每个位为一，否则为零。我们可以使用 warp vote指令，为此调用了 __ballot()。每条车道将从选票中获得相同的值:

     1  0  1  1  0  1  1  1  (this is a single binary 8-bit value in this case)
   

2. 现在每个 warp lane 都将采用该值，并在 warp lane 位置添加一个值为 1 的值。在示例的其余部分使用泳道 1:

   + 0  0  0  0  0  0  1  0  (the only 1 bit in this value will be at the lane index)
   = 1  0  1  1  1  0  0  1
   

3. 我们现在将第 2 步的结果与第 1 步的结果按位异或:

   = 0  0  0  0  1  1  1  0
   

4. 我们现在计算这个值中 1 的位数(有一个 __popc() intrinsic 为此)，并从结果中减去 1。因此对于上面的 channel 1 示例，此步骤的结果将是 2，因为设置了 3 个位。这给出了到我们左边第一个值的距离，它在原始移位模式中是非零的。因此对于车道 1 示例，车道 1 左侧的第一个非零值比车道高 2 个车道，即车道 3。

5. 对于每条车道，我们使用第 4 步的结果为该车道获取适当的偏移值。我们可以使用 __shfl_down() 一次处理所有车道 warp shuffle说明。

     0, 1, 1, 1, 2, 2, 2, 2
   

   从而产生我们想要的“偏移模式”。

一旦我们有了所需的偏移模式，让每个扭曲 channel 使用其偏移值来适本地移动其数据项的过程就很简单了。

这是一个完整的示例，使用了您的 3 个测试用例。上面的步骤 1-4 包含在 __device__ 函数 mydelta 中。内核的其余部分正在执行步骤 5 洗牌，适本地索引数据并复制数据。由于使用了 warp shuffle 指令，我们必须为 cc3.0 或更高版本的 GPU 编译它。 (但是，用其他允许在 cc2.0 或更高版本设备上运行的索引代码替换 warp shuffle 指令并不难。)此外，由于使用了各种内在函数，此函数不能用于超过 32 个数据项，但这是您问题中所述的先决条件。

$ cat t475.cu
#include <stdio.h>
#define DSIZE 8

#define cudaCheckErrors(msg) \
    do { \
        cudaError_t __err = cudaGetLastError(); \
        if (__err != cudaSuccess) { \
            fprintf(stderr, "Fatal error: %s (%s at %s:%d)\n", \
                msg, cudaGetErrorString(__err), \
                __FILE__, __LINE__); \
            fprintf(stderr, "*** FAILED - ABORTING\n"); \
            exit(1); \
        } \
    } while (0)


__device__ int mydelta(const int shift){
  unsigned nz = __ballot(shift == 0);
  unsigned mylane = (threadIdx.x & 31);
  unsigned lanebit = 1<<mylane;
  unsigned temp = nz + lanebit;
  temp = nz ^ temp;
  unsigned delta = __popc(temp);
  return delta-1;
}
__global__ void mykernel(const int *data, const unsigned *shift, int *result, const int limit){ // limit <= 32
  if (threadIdx.x < limit){
    unsigned lshift = shift[(limit - 1) - threadIdx.x];
    unsigned delta = mydelta(lshift);
    unsigned myshift = __shfl_down(lshift, delta);
    myshift = __shfl(myshift, ((limit -1) - threadIdx.x)); // reverse offset pattern
    result[threadIdx.x] = 0;
    if ((myshift + threadIdx.x) < limit)
    result[threadIdx.x + myshift] = data[threadIdx.x];
  }
}

int main(){
  int A[DSIZE]         = {3, 6, 7, 8, 1, 2, 3, 5};
  unsigned tc1B[DSIZE] = {0, 1, 0, 0, 0, 0, 0, 0};
  unsigned tc2B[DSIZE] = {0, 1, 0, 0, 2, 0, 0, 0};
  unsigned tc3B[DSIZE] = {0, 1, 0, 0, 1, 0, 0, 0};

  int *d_data, *d_result, *h_result;
  unsigned *d_shift;
  h_result = (int *)malloc(DSIZE*sizeof(int));
  if (h_result == NULL) { printf("malloc fail\n"); return 1;}
  cudaMalloc(&d_data, DSIZE*sizeof(int));
  cudaMalloc(&d_shift, DSIZE*sizeof(unsigned));
  cudaMalloc(&d_result, DSIZE*sizeof(int));
  cudaCheckErrors("cudaMalloc fail");
  cudaMemcpy(d_data, A, DSIZE*sizeof(int), cudaMemcpyHostToDevice);
  cudaMemcpy(d_shift, tc1B, DSIZE*sizeof(unsigned), cudaMemcpyHostToDevice);
  cudaCheckErrors("cudaMempcyH2D fail");
  mykernel<<<1,32>>>(d_data, d_shift, d_result, DSIZE);
  cudaDeviceSynchronize();
  cudaCheckErrors("kernel fail");
  cudaMemcpy(h_result, d_result, DSIZE*sizeof(int), cudaMemcpyDeviceToHost);
  cudaCheckErrors("cudaMempcyD2H fail");
  printf("index: ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", i);
  printf("\nA:     ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", A[i]);
  printf("\ntc1 B: ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", tc1B[i]);
  printf("\ntc1 C: ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", h_result[i]);
  cudaMemcpy(d_shift, tc2B, DSIZE*sizeof(unsigned), cudaMemcpyHostToDevice);
  cudaCheckErrors("cudaMempcyH2D fail");
  mykernel<<<1,32>>>(d_data, d_shift, d_result, DSIZE);
  cudaDeviceSynchronize();
  cudaCheckErrors("kernel fail");
  cudaMemcpy(h_result, d_result, DSIZE*sizeof(int), cudaMemcpyDeviceToHost);
  cudaCheckErrors("cudaMempcyD2H fail");
  printf("\ntc2 B: ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", tc2B[i]);
  printf("\ntc2 C: ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", h_result[i]);
  cudaMemcpy(d_shift, tc3B, DSIZE*sizeof(unsigned), cudaMemcpyHostToDevice);
  cudaCheckErrors("cudaMempcyH2D fail");
  mykernel<<<1,32>>>(d_data, d_shift, d_result, DSIZE);
  cudaDeviceSynchronize();
  cudaCheckErrors("kernel fail");
  cudaMemcpy(h_result, d_result, DSIZE*sizeof(int), cudaMemcpyDeviceToHost);
  cudaCheckErrors("cudaMempcyD2H fail");
  printf("\ntc3 B: ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", tc3B[i]);
  printf("\ntc2 C: ");
  for (int i = 0; i < DSIZE; i++)
    printf("%d, ", h_result[i]);
  printf("\n");
  return 0;
}
$ nvcc -arch=sm_35 -o t475 t475.cu
$ ./t475
index: 0, 1, 2, 3, 4, 5, 6, 7,
A:     3, 6, 7, 8, 1, 2, 3, 5,
tc1 B: 0, 1, 0, 0, 0, 0, 0, 0,
tc1 C: 3, 0, 6, 7, 8, 1, 2, 3,
tc2 B: 0, 1, 0, 0, 2, 0, 0, 0,
tc2 C: 3, 0, 6, 7, 8, 0, 1, 2,
tc3 B: 0, 1, 0, 0, 1, 0, 0, 0,
tc2 C: 3, 0, 6, 7, 8, 1, 2, 3,
$