CUDA编程04: 线程束基本函数与协作组

YiQi 管理员
  2023-01-14 21:30:50 2023-01-14 21:30 创建   2023-03-11 22:57:53 2023-03-11 22:57 更新      3.1k 字  13 分钟  

一个线程块中连续的32个线程为一个线程束 (warp)。一个SM可以处理一个或多个线程块，一个线程块又可分为若干个线程束。

单指令多线程执行模式

在伏特架构之前，一个线程束中的线程拥有同一个程序计数器。一个线程束同一时刻只能执行一个共同的指令或闲置，即单指令多线程执行模式（single instruction multiple thread, SIMT）。

例如：

1
2
3
4
if (condition) 
    A
else
    B

满足condition的线程会执行语句A，其他的线程闲置，反之亦然，即发生了分支发散（branch divergence）。

从伏特架构开始，引入了独立线程调度机制（independent thread scheduling），每个线程有自己的程序计数器，代价是增加了寄存器负担：单个线程的程序计数器一般需要使用两个寄存器。另外，独立线程调度机制使得假设了线程束同步的代码变得不再安全，需要显式指定同步。

线程束内的线程同步函数

在我们的归约问题中,当所涉及的线程都在一个线程束内时，可以将线程块同步函数__synthreads()换成一个更加廉价的线程束同步函数__syncwarp()，它的原型为：

1
void __syncwarp(unsigned mask = 0xffffffff)

该函数有一个可选的参数，该参数是一个代表掩码的无符号整型数，默认值的全部32个二进制位都为1，代表线程束中的所有线程都参与同步。如果要排除一些线程，可以用一个对应的二进制位为0的掩码参数，例如，掩码Oxfffffffe代表排除第0号线程。基于此，可将上一章的归约核函数修改为

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
void __global__ reduce_syncwarp(const real *d_x, real *d_y, const int N) {
    const int tid = threadIdx.x;
    const int bid = blockIdx.x;
    const int n = bid * blockDim.x + tid;
    extern __shared__ real s_y[];
    s_y[tid] = (n < N) ? d_x[n] : 0.0;
    __syncthreads();

    for (int offset = blockDim.x >> 1; offset >= 32; offset >>= 1) {
        if (tid < offset) {
            s_y[tid] += s_y[tid + offset];
        }
        __syncthreads();
    }

    for (int offset = 16; offset > 0; offset >>= 1) {
        if (tid < offset) {
            s_y[tid] += s_y[tid + offset];
        }
        __syncwarp();
    }

    if (tid == 0) {
        atomicAdd(d_y, s_y[0]);
    }
}

完整代码见 reduceGPU_warp.cu

当offset >= 32时还是使用线程块同步，当offset <= 16是启用线程束同步。使用

1
nvcc -O3 -arch=sm_75 reduce.cu

进行编译后，程序大概耗时3.5ms，比起原函数大概快了10%。注意不能写成

1
2
3
4
for (int offset = 16; offset > 0; offset >>= 1) {
    s_y[tid] += s_y[tid + offset];
    __syncwarp();
}

比如对tid = 0和tid = 16，分别有s_y[0] += s_y[16]和s_y[16] += s_y[32]，既读又写，会有读写竞争（race condition）。

更多线程束内的基本函数

线程束表决函数：

1
2
3
unsigned __ballot_sync(unsigned mask, int predicate);
int __all_sync(unsigned mask, int predicate);
int __any_sync(unsigned mask, int predicate);

线程束洗牌函数：

1
2
3
4
T __shfl_sync(unsigned mask, T v, int srcLane, int w = warpSize);
T __shfl_up_sync(unsigned mask, T v, unsigned d, int w = warpSize);
T __shfl_down_sync(unsigned mask, T v, unsigned d, int w = warpSize);
T __shfl_xor_sync(unsigned mask, T v, int laneMask, int w = warpSize);
  • w只能取2、4、8、16、32
  • 获取束内指标，可使用int lane_id = threadIdx.x % w;或者使用按位与int lane_id = threadIdx.x & (w - 1);
  • 参数mask见__syncwarp
  • 各种函数返回的结果对被掩码排除的线程来说是没有定义的。所以不要尝试在这些被排除的线程中使用函数的返回值

各函数功能：

  • unsigned __ballot_sync(mask, predicate)：
    如果线程束内第n个线程参与计算且predicate值非零，则将所返回无符号整数的第n个二进制位取为1，否则取为0
  • int __all_sync(mask, predicate)
    线程束内所有参与线程的predicate值都不为零才返回1，否则返回0
  • int __any_sync(mask, predicate)
    线程束内所有参与线程的predicate值有一个不为零则返回1，否则返回0
  • T __shfl_sync(mask, v, srcLane, w)
    参与线程返回标号为srcLane的线程中变量v的值。这是一种广播式数据交换，即将一个线程中的数据广播到所有（包括自己）线程
  • T __shfl_up_sync(mask, v, d, w)
    标号为t的参与线程返回标号为t - d的线程中变量v的值。标号满足t - d < O的线程返回原来的v。例如，当w＝8，d = 2时，该函数将第0 ~ 5号线程中变量v的值传送到第2 ~ 7号线程，而第0 ~ 1号线程返回它们原来的v。形象地说，这是—种将数据向上平移的操作
  • T __shfl_down_sync(mask, v, d, w)
    标号为t的参与线程返回标号为t + d的线程中变量v的值。标号满足t + d >= O的线程返回原来的v。例如，当w＝8，d = 2时，该函数将第2 ~ 7号线程中变量v的值传送到第0 ~ 5号线程，而第6 ~ 7号线程返回它们原来的v。形象地说，这是—种将数据向下平移的操作
  • T __shfl_xor_sync(mask, v, laneMask, w)
    标号为t的参与线程返回标号为t ^ 1aneMask的线程中变量v的值。这里，t ^ 1aneMask表示两个整数按位异或运算的结果。例如，当w = 8，1aneMask = 2时，第0 ~ 7号线程的按位异或运算t ^ 1aneMask分别如下： 
    1
    2
    3
    4
    5
    6
    7
    8
    0 ^ 2 = 0000 ^ 0010 = 0010 = 2
    1 ^ 2 = 0001 ^ 0010 = 0011 = 3
    2 ^ 2 = 0010 ^ 0010 = 0000 = 0
    3 ^ 2 = 0011 ^ 0010 = 0001 = 1
    4 ^ 2 = 0100 ^ 0010 = 0110 = 6
    5 ^ 2 = 0101 ^ 0010 = 0111 = 7
    6 ^ 2 = 0110 ^ 0010 = 0100 = 4
    7 ^ 2 = 0111 ^ 0010 = 0101 = 5
    有一测试程序，见warp_func.cu 
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    threadIdx.x:  0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 
    lane_id:      0  1  2  3  4  5  6  7  0  1  2  3  4  5  6  7 
    FULL_MASK = ffffffff
    mask1     = fffe
    mask2     = 1
    all_sync (FULL_MASK): 0
    all_sync     (mask1): 1
    any_sync (FULL_MASK): 1
    any_sync     (mask2): 0
    shfl:       2  2  2  2  2  2  2  2 10 10 10 10 10 10 10 10 
    shfl_up:    0  0  1  2  3  4  5  6  8  8  9 10 11 12 13 14 
    shfl_down:  1  2  3  4  5  6  7  7  9 10 11 12 13 14 15 15 
    shfl_xor:   1  0  3  2  5  4  7  6  9  8 11 10 13 12 15 14

利用线程束洗牌函数进行归约计算

几个线程束洗牌函数中，T __shfl_down_sync()比较适合进行我们的归约计算，函数如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
void __global__ reduce_shfl(const real *d_x, real *d_y, const int N) {
    const int tid = threadIdx.x;
    const int bid = blockIdx.x;
    const int n = bid * blockDim.x + tid;
    extern __shared__ real s_y[];
    s_y[tid] = (n < N) ? d_x[n] : 0.0;
    __syncthreads();

    for (int offset = blockDim.x >> 1; offset >= 32; offset >>= 1) {
        if (tid < offset) {
            s_y[tid] += s_y[tid + offset];
        }
        __syncthreads();
    }

    real y = s_y[tid];

    for (int offset = 16; offset > 0; offset >>= 1) {
        y += __shfl_down_sync(FULL_MASK, y, offset);
    }

    if (tid == 0) {
        atomicAdd(d_y, y);
    }
}

完整代码见 reduceGPU_warp.cu

  • 在进行线程束内的循环之前，这里将共享内存中数据复制到了寄存器，寄存器更高效
  • 因为洗牌函数能够自动处理同步与读写竞争问题，所以去掉了同步函数

在本人机器上测试，大概耗时2.8ms.

协作组

使用协作组功能要包含源文件

1
#include <cooperative_groups.h>

并使用cooperative_groups命名空间

线程块级别的协作组

（这里东西有点多，这几个类型间的关系有点迷惑）

协作组编程模型中最基本的类型（基类？）是线程组thread_group，有如下成员：

  1. void sync(): 同步组内所有线程
  2. unsigned size(): 组的大小
  3. unsigned thread_rank(): 当前调用该函数的线程在组内的标号
  4. bool is_valid(): 如果定义的组违反了任何CUDA限制，则为假，否则为真

线程组类型有一个称为thread_block的导出类型，有额外的两个函数：

  1. dim3 group_index(): 等价于blockIdx
  2. dim3 thread_index(): 等价于threadIdx

可以使用

1
thread_block g = this_thread_block();

定义并初始化一个线程块对象。g.sync()完全等价于__syncthreads()，g.group_index()完全等价于blockIdx，g.thread_index()完全等价于threadIdx。

可以使用函数tiled_partition()将一个线程块划分为若干片（tile），每一片构成一个新的线程组。目前仅仅可以将片的大小设置为2的正整数次方且不大于32。例如，下面语句将一个线程块分割为我们熟知的线程束：

1
thread_group g32 = tiled_partition(this_thread_block(), 32);

还可以进一步细分，比如把一个线程束再分割为包含4个线程的线程组：

1
thread_group g4 = tiled_partition(g32, 4);

如果分割在编译期就已知，可使用模板：

1
2
thread_block_tile<32> g32 = tiled_partition<32>(this_thread_block());
thread_block_tile<4>  g4  = tiled_partition<4>(this_thread_block());

这样定义的线程组称为线程块片（thread block tile）。线程块片还有一系列类似线程束内函数的方法：

1
2
3
4
5
6
7
unsigned __ballot_sync(int predicate);
int __all_sync(int predicate);
int __any_sync(int predicate);
T __shfl_sync(T v, int srcLane);
T __shfl_up_sync(T v, unsigned d);
T __shfl_down_sync(T v, unsigned d);
T __shfl_xor_sync(T v, int laneMask);

线程块片的函数有两点不同:

  1. 线程块片的函数少了第亿个代表掩码的参数，因为线程组内的所有线程都必须参与相关函数的运算
  2. 线程块片的洗牌函数（上述函数中的后4个）少了最后一个代表宽度的参数，因为该宽度就是线程块片的大小，即定义线程块片的模板参数。

利用协作组进行归约计算

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
void __global__ reduce_cp(const real *d_x, real *d_y, const int N) {
    const int tid = threadIdx.x;
    const int bid = blockIdx.x;
    const int n = bid * blockDim.x + tid;
    extern __shared__ real s_y[];
    s_y[tid] = (n < N) ? d_x[n] : 0.0;
    __syncthreads();

    for (int offset = blockDim.x >> 1; offset >= 32; offset >>= 1) {
        if (tid < offset)
        {
            s_y[tid] += s_y[tid + offset];
        }
        __syncthreads();
    }

    real y = s_y[tid];

    thread_block_tile<32> g = tiled_partition<32>(this_thread_block());
    for (int i = g.size() >> 1; i > 0; i >>= 1) {
        y += g.shfl_down(y, i);
    }

    if (tid == 0) {
        atomicAdd(d_y, y);
    }
}

完整代码见 reduceGPU_warp.cu
计算耗时与使用线程束洗牌函数的一样。

数组归约程序的进一步优化

提高线程利用率

之前折半归约只有1/2，1/4，1/8, … 的线程在工作，其他闲置。线程利用率较低。

提升的中心思想是，一个线程以跨度为grid_size * block_size去处理整个数组的数据，再用洗牌函数归约到长度为grid_size的数组，再调用一次函数，此时grid_size = 1，这样就可以将剩下的数归约到d_y[0]了。完整代码见 reduceGPU_stride.cu

归约函数为：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
void __global__ reduce_cp(const real *d_x, real *d_y, const int N) {
    const int tid = threadIdx.x;
    const int bid = blockIdx.x;
    extern __shared__ real s_y[];

    real y = 0.0;
    const int stride = blockDim.x * gridDim.x;
    for (int n = bid * blockDim.x + tid; n < N; n += stride) {
        y += d_x[n];
    }
    s_y[tid] = y;
    __syncthreads();

    for (int offset = blockDim.x >> 1; offset >= 32; offset >>= 1) {
        if (tid < offset)
        {
            s_y[tid] += s_y[tid + offset];
        }
        __syncthreads();
    }

    y = s_y[tid];

    thread_block_tile<32> g = tiled_partition<32>(this_thread_block());
    for (int i = g.size() >> 1; i > 0; i >>= 1) {
        y += g.shfl_down(y, i);
    }

    if (tid == 0) {
        d_y[bid] = y;
    }
}

调用函数为：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
real reduce(const real *d_x) {
    const int ymem = sizeof(real) * GRID_SIZE;
    const int smem = sizeof(real) * BLOCK_SIZE;

    real h_y[1] = {0};
    real *d_y;
    CHECK(cudaMalloc(&d_y, ymem));

    reduce_cp<<<GRID_SIZE, BLOCK_SIZE, smem>>>(d_x, d_y, N);
    reduce_cp<<<1, 1024, sizeof(real) * 1024>>>(d_y, d_y, GRID_SIZE);

    CHECK(cudaMemcpy(h_y, d_y, sizeof(real), cudaMemcpyDeviceToHost));
    CHECK(cudaFree(d_y));

    return h_y[0];
}

本机执行时间约为1.3ms左右，有大幅度提升。更为重要的是，该计算结果为123000064.0，相比精确结果123000000.0有七位有效数字。之前使用原子函数所得到的结果为123633392.0，仅有3位准确的有效数字。这是因为，在使用两个核函数时，将数组d_y归约到最终结果的计算也使用了折半求和，比直接累加要稳健。

避免反复分配与释放设备内存

设备内存的分配与释放是比较耗时的。可以使用静态全局内存为d_y提前分配好空间（编译期），在使用运行时API函数cudaGetSymbolAddress()将d_y与静态全局内存地址绑定。完整代码见 reduceGPU_static.cu 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
__device__ real static_y[GRID_SIZE]; 

real reduce(const real *d_x) {
    real *d_y;
    CHECK(cudaGetSymbolAddress((void**)&d_y, static_y));
    
    const int smem = sizeof(real) * BLOCK_SIZE;

    reduce_cp<<<GRID_SIZE, BLOCK_SIZE, smem>>>(d_x, d_y, N);
    reduce_cp<<<1, 1024, sizeof(real) * 1024>>>(d_y, d_y, GRID_SIZE);

    real h_y[1] = {0};
    CHECK(cudaMemcpy(h_y, d_y, sizeof(real), cudaMemcpyDeviceToHost));
    // CHECK(cudaMemcpyFromSymbol(h_y, static_y, sizeof(real)); // also ok

    return h_y[0];
}

本机执行时间约为1.1ms，有进一步的提升。

作者的测试结果如为下表：

此页目录
CUDA编程04: 线程束基本函数与协作组
  1. 单指令多线程执行模式
  2. 线程束内的线程同步函数
  3. 更多线程束内的基本函数
    1. 利用线程束洗牌函数进行归约计算
  4. 协作组
    1. 线程块级别的协作组
    2. 利用协作组进行归约计算
  5. 数组归约程序的进一步优化
    1. 提高线程利用率
    2. 避免反复分配与释放设备内存