解析OneFlow Element-Wise算子实现方法

0x0. 前言
由于cuda水平太菜,所以一直没写过这方面的笔记。现在日常的工作中已经不能离开写cuda代码,所以准备学习zzk随缘做一做cuda的笔记记录一下学习到的知识和技巧。这篇文章记录的是阅读oneflow的element-wise系列cuda算子实现方案学习到的技巧,希望可以帮助到一起入门cuda的小伙伴们。elemet-wise算子指的是针对输入tensor进行逐元素操作,比如relu就是针对输入tensor的每个值进行判断是否大于0,大于0的话输出就是输入否则就是0。用cuda来表达最简单的写法就是:
__global__ void relu_kernel(float* input, float* output){  int32_t idx = blockidx.x * blockdim.x + threadidx.x;  output[idx] = input[idx] < 0 ? 0 : input[idx];}int main(){  float* input;  float* output;  int32_t elem_cnt = 3*224*224;    cudamalloc(&input, sizeof(float)*elem_cnt);  cudamalloc(&output, sizeof(float)*elem_cnt);  int32_t thread_num = 256;  int32_t grid_size = (elem_cnt + thread_num -1) / thread_num;  relu_kernel<>(src, dst);     cudadevicesynchronize();  cudafree(src);  cudafree(dst);  return 0;} 虽然这种写法非常简单明了,但却存在明显的性能问题。所以这篇文章将基于oneflow开源的element-wise cuda算子方案来解释如何写一个高性能的element-wise cuda算子。
0x1. 性能
以gelu激活函数为例子,分别测试 dtype = float32,不同shape下的前向耗时以及带宽利用率(nvidia a100-pcie-40gb)。性能情况如下图所示:
在这里插入图片描述
在这里插入图片描述
可以看到对于 gelu 来说,无论是性能还是带宽 oneflow 的实现都是更优的,接下来我们就来了解一下为什么 oneflow 的 element-wise 算子性能可以做到更优。
0x2. 用法
oneflow在 elementwise.cuh 文件中分别针对一元,二元,三元运算的 element-wise 操作实现了模板函数。在包含这个头文件之后我们可以使用  cuda::unary/binary/ternary  这几个模板函数来针对我们自己定义的 element-wise 操作进行计算。注意,这里说的一元,二元,三元代表的是这个 element-wise 操作有几个输入 tensor。
我们举个例子,假设我们要做的 element-wise 操作是逐点乘法,也即有 2 个输入tensor x 和 y,然后 x 和 y的形状和数据类型都是一致的。那么我们可以定义一个模板类:
templatestruct multiplyfunctor {  of_device_func t operator()(t x, t y) const {    return x*y;  }}; 这里 of_device_func 表示我们定义的这个函数既可以运行在 cpu 又可以运行在 gpu 上,它的定义是:
#if defined(__cudacc__)#define of_device_function __device__ __host__ __forceinline__#else#define of_device_function inline#endif 然后我们就可以使用 cuda::binary 这个模板函数来完成这个二元的 element-wise 算子了。示例代码如下:
    const user_op::tensor* x = ctx->tensor4argnameandindex(x, 0);    const user_op::tensor* y = ctx->tensor4argnameandindex(y, 0);    user_op::tensor* out = ctx->tensor4argnameandindex(out, 0);    const int64_t elem_cnt = x->shape().elem_cnt();    of_cuda_check(cuda::binary(multiplyfunctor(), elem_cnt, out->mut_dptr(),                                            x->dptr(),                                             y->dptr(),                                            ctx->device_ctx()->cuda_stream())); 这里的 x, y, out 分别代表这个 element-wise 操作的输入输出 tensor,然后 element_cnt 表示 tensor 的元素个数,输出张量的数据首地址 out->mut_dptr(), 输入张量的数据首地址 x->dptr() && y->dptr() ,最后一个参数则是当前 kernel 运行的 cuda stream对象。
0x3. 原理&&代码实现解析
我个人认为这里有几个要点,分别是一个线程处理多个数据,向量化数据访问提升带宽,设置合理的block数量(gridsize)和线程数量(blocksize)以及在合适的地方进行循环展开(unrool)以及一些编程上的技巧。
0x3.1 给 element-wise 操作设置合理的 gridsize 和 blocksize
下面这段代码展示了 oneflow 针对 element-wise 算子是如何设置 gridsize 和 blocksize 的。对应的源码地址为:https://github.com/oneflow-inc/oneflow/blob/master/oneflow/core/cuda/elementwise.cuh#l30-l52 。
constexpr int kblocksize = 256;constexpr int knumwaves = 32;inline cudaerror_t getnumblocks(int64_t n, int* num_blocks) {  int dev;  {    cudaerror_t err = cudagetdevice(&dev);    if (err != cudasuccess) { return err; }  }  int sm_count;  {    cudaerror_t err = cudadevicegetattribute(&sm_count, cudadevattrmultiprocessorcount, dev);    if (err != cudasuccess) { return err; }  }  int tpm;  {    cudaerror_t err = cudadevicegetattribute(&tpm, cudadevattrmaxthreadspermultiprocessor, dev);    if (err != cudasuccess) { return err; }  }  *num_blocks = std::max(1, std::min((n + kblocksize - 1) / kblocksize,                                                   sm_count * tpm / kblocksize * knumwaves));  return cudasuccess;} 这个地方 blocksize 直接被设置为了 256 ,对应 constexpr int kblocksize = 256; 这行代码,也就是说每个 block 有 256 个线程。为什么是 256 ?大家不妨读一下俊丞大佬这篇经典的 给cuda kernel设置合适的 gridsize 和 block size 的文章 。文章中通过对 sm 的资源分析确定在主流的gpu上将 blocksize 设置为 128 或者 256 是比较合适,在这里直接设置为了 256 。
确定了 blocksize 之后需要确定 kernel 启动线程块的数量,我一直觉得上述文章中对这一段的分析是尤其精彩的,这里再截图展示一下:
选自oneflow cuda kernel 中 grid_size 和 block_size 应该怎么设置 一文
根据这里的分析,对于 element-wise 操作要设置合适的 gridsize 不仅需要考虑元素的数量还要考虑由于 sm 硬件本身带来的限制。如下公式所述:
*num_blocks = std::max(1, std::min((n + kblocksize - 1) / kblocksize,                                                   sm_count * tpm / kblocksize * knumwaves)); 这里的 (n + kblocksize - 1) / kblocksize 就是根据 element-wise 操作的元素个数来计算需要启动多少个线程块,比如在文章开头的例子中有 = 个元素,那么就一共需要 个线程块。然后这里以gtx 3080ti为例,它的sm个数也就是sm_count=80,每个sm最多调度的线程数tpm=1536,那么sm_count * tpm / kblocksize * knumwaves = 80 * 1536 / 256 * 32 = 15360,所以在这个例子中我们最终设置的线程块个数为 588 个。
通过上述讲解和分析我们已经确定了启动 element-wise cuda kernel 的 gridsize 和 blocksize。
0x3.2 向量化数据访问提升带宽
对于大多数 element-wise 算子来说,一般它们的计算量不会太大,所以它们的瓶颈一般在gpu的带宽上。在 nvidia 的性能优化博客 https://developer.nvidia.com/blog/cuda-pro-tip-increase-performance-with-vectorized-memory-access/ 中提到,对于很多 cuda 核函数我们都可以通过向量化数据访问的方式来提升带宽受限的 kernel 的性能,特别是对于架构比较新的 gpu 向量化数据访问的效果会更加明显。
在 oneflow 的 element-wise 系列算子中,为了更好的进行向量化的数据访问,俊丞设计了如下的 pack 数据结构(代码位置:https://github.com/oneflow-inc/oneflow/blob/master/oneflow/core/cuda/elementwise.cuh#l54-l70):
templatestruct getpacktype {  using type = typename std::aligned_storage::type;};templateusing packtype = typename getpacktype::type;templateunion pack {  static_assert(sizeof(packtype) == sizeof(t) * pack_size, );  __device__ pack() {    // do nothing  }  packtype storage;  t elem[pack_size];}; 对getpacktype理解有误请看知乎的修改后正确版本用了 std::aligned_storage 先声明了一个内存对齐的数据类型 type ,注意这个 type 的内存长度为 pack_size * sizeof(t) 。然后这里的 t 是我们需要进行 pack 的数据类型,而 pack_size 则表示我们需要 pack 的元素个数。接下来我们看到 pack 联合体中声明了 storage 和 elem 两个数组,它们公用同一段对齐的内存。然后 pack 联合体的入口有一个检查: static_assert(sizeof(packtype) == sizeof(t) * pack_size, ); 这是用来判断我们之前声明的 type 的内存长度是否符合预期。
接下来我们从 https://github.com/oneflow-inc/oneflow/blob/master/oneflow/core/cuda/elementwise.cuh#l155-l194 这里可以看到这个 pack 联合体主要是用在 kernel 启动之前判断 element-wise 操作的输入输出 tensor 对应的数据指针地址是否满足内存对齐的条件,如果不满足则这个 element-wise 操作无法执行数据 pack 。对应下图2个画红色框的地方。
接下来,oneflow 定义了真正要执行数据 pack 的数据结构 packed 并且定义了计算 packsize 的工具函数。代码位置为:https://github.com/oneflow-inc/oneflow/blob/master/oneflow/core/cuda/elementwise.cuh#l72-l95 。
templatestruct alignas(sizeof(t) * pack_size) packed {  __device__ packed() {    // do nothing  }  union {    t elem[pack_size];  };};constexpr int kmaxpackbytes = 128 / 8;constexpr int kmaxpacksize = 8;constexpr int min(int a, int b) { return a < b ? a : b; }templateconstexpr int packsize() {  return min(kmaxpackbytes / sizeof(t), kmaxpacksize);}templateconstexpr int packsize() {  return min(packsize(), packsize());} 这里需要注意的是对于 cuda 来说,最多支持 128 个 bit 的访问粒度,也就是说 packsize 的大小不能超过 128 个bit。然后对于各种数据类型来说,half 数据类型的 bit 数是最少的即 16,所以一次性可以支持 pack 8个half类型的数据,4个float32的数据,以此类推。所以这里的定义的 kmaxpacksize 表示 128/16=8 ,然后 kmaxpackbytes 则表示最大可以 pack 的 byte 数 。
请注意区分 bit 和 byte 。
接下来 https://github.com/oneflow-inc/oneflow/blob/master/oneflow/core/cuda/elementwise.cuh#l97-l144 则是真正的为 element-wise 操作完成数据 pack 并执行计算。
首先来看这段充满技巧的代码:
在这里插入图片描述
首先这里定义了一个 hasapply2 类用来判断是否可以支持一次性pack 2个 char/int8/half2 类型的元素,这个地方是一个针对 int8/half2/char 数据类型的特殊处理,某些 element-wise 算子 kernel 确实需要支持这种数据类型的计算。也就是说对于 half2 的话,在一个内存访问粒度里我们其实是可以 pack 128 / 8 = 16个的。然后用了c++模板元编程的 std::enable_if 来控制针对 half2 类型的特殊 pack 处理,也就是上图代码中的两个 applypack 函数。可以看到对于 half2 类型的 element-wise 操作我们需要给对应的 functor 定义一个 apply2 函数,比如对于 cast 操作的 functor 定义如下:
templatestruct castfunctor {  __device__ to operator()(from from) const { return static_cast(from); }};templatestruct castfunctor {  __device__ to operator()(half from) const { return static_cast(static_cast(from)); }  __device__ void apply2(to* to, const half* from) const {    const float2 f2 = __half22float2(*reinterpret_cast(from));    to[0] = static_cast(f2.x);    to[1] = static_cast(f2.y);  }}; 0x3.3 启动 kernel
我们接下来看一下 element-wise 的 kernel 实现:https://github.com/oneflow-inc/oneflow/blob/master/oneflow/core/cuda/elementwise.cuh#l133-l144 。
在这里插入图片描述
在 kernel 中我们发现每一个线程实际上处理了多个 pack 后的数据,也即:for (int64_t i = global_tid; i  1026-1024 = 2  int num_blocks;  {    cudaerror_t err = getnumblocks(n_pack, &num_blocks); // 计算线程块数目    if (err != cudasuccess) { return err; }  }  applygeneric<>(      factory, n_pack, reinterpret_cast(r),      (reinterpret_cast(in))..., n_tail, r + tail_offset,      (in + tail_offset)...);   return cudapeekatlasterror();} 0x4. 总结
以上就是我对 oneflow element-wise 系列 cuda 算子实现的解析,后续有空会持续更新学习到的新知识。


基于FPGA 的FMC 接口应用实例
IC生产流程步骤详解 从上游到下游解说
LabVIEW Tab选项卡控件XTab的使用方法
P-HJT新纪录浮出水面,N型难说替代P型?
GC1101射频前端单芯片概述、应用优势及主要特性
解析OneFlow Element-Wise算子实现方法
全球半导体材料市场呈现增长趋势,中国台湾市场位居全球第一
一个静电保护ESD的工艺电路解析
网约车新政内容解读:哪些人可以加入打车司机行列?
成为“全球第一”需要三年吗?小米全球化战略再升级,越南工厂交付第一批手机
需求迫切的医疗器械产品上市 加快了高端医疗器械国产化替代的步伐
国内新能源汽车增长迅速 但短板也非常明显
云原生计算对这三个热门市场的影响
Agilent安捷伦E4438C射频发生器6GHz
rt-thread SDIO驱动框架分析(贴片SD卡flash驱动\SD Nand flash驱动)
使用RX单片机实现数字电源控制的示例
国风设计 北通新品SWITCH游戏机收纳包发布
亚马逊旗下的Ring推出了其第一款视频门铃的更新版本
rfid在智慧医疗上可以怎样融入进去
用SwitchLaboVR玩游戏是什么样的体验