SIMD优化应用

SIMD优化应用

SIMD 介绍

费林分类法(Flynn's Taxonomy)

Ÿ   单一指令流单一数据流（SISD）

Ÿ   单一指令流多数据流（SIMD）

Ÿ   多指令流单一数据流（MISD）

Ÿ   多指令流多数据流（MIMD）

 

 

 SIMD(Single Instruction Multiple Data) 是一种处理器指令类型，即单个指令可以同时处理多个数据。这里的多个，一般指的就是就是多个 CPU 寄存器。

例如，以下为加法的标量(Scalar)与SIMD(Vector)两种执行方式:

 

 

 现代大多数处理器指令为 SISD / SIMD, 而多核计算方式被认为是 MIMD。

处理器如何实现 SIMD

 

 

 

 

 有一段经典SIMD应用代码 ClickHouse: Columns/ColumnsCommon.cpp, filterArraysImplGeneric()，这段代码的注释是这样形容的:

/** A slightly more optimized version.

* Based on the assumption that often pieces of consecutive values

*  completely pass or do not pass the filter.

* Therefore, we will optimistically check the parts of `SIMD_BYTES` values.

*/

在这段代码中，filterArraysImplGeneric 函数需要完成下面一段 SQL Where 子句的部分逻辑:

Ÿ   有一个数组 data, 里面存放了数据

Ÿ   有一个 int8 类型数组 filter，size等于 data，存放了一些布尔值，用于标识 data 里面的数据是否应该被筛选出来

Ÿ   生成一个新的数组 res，根据 filter 数组，对 data 数组进行筛选

这个逻辑正常实现很简单，伪码如:

let res = []

for (let i = 0; i < data.size(); i ++) {

    if (filter[i] != 0) {

        res.append(data[i])

    }

}

有什么优化点呢 ？如果考虑到 data 的长度很长的话 (在 ClickHouse 的应用场景，这个 Array 的长度，在一次查询中，所有执行累加在一起，通常在亿级以上)

ClickHouse 的实现思路如下:

Ÿ   上述代码耗时因素在于循环次数非常多，等于 data 数组的长度

Ÿ   如果可以降低循环次数，同时保证单次循环耗时变化不大，总体执行效率更高

Ÿ   要满足上述条件，需要在同一个指令周期做更多运算，SIMD 指令可以满足诉求

Ÿ   考虑到目前 SIMD 指令最大支持 512 比特数据，filter 本身是 8 比特，单循环处理数据量可为 512 / 8 = 64 个 (考虑到当前大多数计算机处理器为 64 位)

Ÿ   每次取 64 个 filter 数组项，组合为一个 64 位无符号整型值 mask，这样每个 filter 数组项占用一个比特位

Ÿ   考虑两种特殊情况: 如果 mask 64 位比特均为 1，那么这个循环中，64 个 data 项都应该拷贝到 res，反之如果 mask 64 位比特均为 0，可以直接跳过循环。这两种特殊情况在业务场景中都经常发生

Ÿ   如果没有命中特殊情况，需要循环 64 次嘛 ？有没有进一步优化方法 ？

Ÿ   有多少项需要复制到新的数组，取决于 mask 中比特位为 1 的个数

Ÿ   __builtin_ctzll 编译器内置函数，可以获取最低比特位为 1 的 index，以此可知道哪项数据应该复制

Ÿ   最低置1比特位的项已经复制后，需要将其比特位置 0，这里有两个都比较高效的方法: _blsr_u64 指令 及 mask = mask & (mask-1)，_blsr_u64 可以查指令集，但是 mask & (mask-1) 具体怎么才可以想出来，需要深厚的编程功底了。

Ÿ   循环设置最低置1比特位，直到 mask 中所有比特位均为 0

Ÿ   完成

分析一下每个 64项 Loop 性能可以提升多少 ？

这里为了方便，每个指令的指令周期相等，现实情况中并非如此
下面的分析，更多是定性而非定量

Ÿ   如果按照最基础实现，需要经历 64 个循环，每个循环需要一次条件判断，以及可能的一次数据读取与写入，大约 3 个指令周期，大约 194 指令周期

Ÿ   下面考虑 ClickHouse 的 SIMD 优化实现

Ÿ   Bytes64MaskToBits64Mask 需要 3 个指令周期

Ÿ   最佳情况, filter 不命中: 无额外指令周期，性能显著高于标量版本

Ÿ   次佳情况, filter 全部命中: 需要大约 10 个指令周期复制数据，性能显著高于标量版本

Ÿ   其余情况，指令周期与 mask 中 1 比特位个数相关

Ÿ   在最差情况下，SIMD 优化版本耗时可能略高于普通实现

代码如下:

    static constexpr size_t SIMD_BYTES = 64;

   

    while (filt_pos < filt_end_aligned)

    {

        // 将 64 字节的 filter 数组，转变为 64 比特的无符号整数

        uint64_t mask = Bytes64MaskToBits64Mask(filt_pos);

 

        // 是否为边界情况: filter 全部命中

        if (0xffffffffffffffff == mask)

        {

            // copy all

            ...

        }

        else

        {

            while (mask)

            {

                // 找出最低位置1比特的 index

                size_t index = __builtin_ctzll(mask);

                copy_array(offsets_pos + index);

 

            #ifdef __BMI__

                // _blsr_u64 将最低置一比特位置 0

                mask = _blsr_u64(mask);

            #else

                // 一种不依赖原生指令支持的方案

                mask = mask & (mask-1);

            #endif

            }

        }

 

        filt_pos += SIMD_BYTES;

        offsets_pos += SIMD_BYTES;

    }

其中，Bytes64MaskToBits64Mask 的实现同样很讲究:

/// Transform 64-byte mask to 64-bit mask

inline UInt64 Bytes64MaskToBits64Mask(const UInt8 * bytes64)

{

#if defined(__AVX512F__) && defined(__AVX512BW__)

    static const __m512i zero64 = _mm512_setzero_epi32();

    UInt64 res = _mm512_cmp_epi8_mask(_mm512_loadu_si512(reinterpret_cast<const __m512i *>(bytes64)), zero64, _MM_CMPINT_EQ);

#elif defined(__AVX__) && defined(__AVX2__)

    static const __m256i zero32 = _mm256_setzero_si256();

    UInt64 res =

        (static_cast<UInt64>(_mm256_movemask_epi8(_mm256_cmpeq_epi8(

        _mm256_loadu_si256(reinterpret_cast<const __m256i *>(bytes64)), zero32))) & 0xffffffff)

        | (static_cast<UInt64>(_mm256_movemask_epi8(_mm256_cmpeq_epi8(

        _mm256_loadu_si256(reinterpret_cast<const __m256i *>(bytes64+32)), zero32))) << 32);

#elif defined(__SSE2__) && defined(__POPCNT__)

    static const __m128i zero16 = _mm_setzero_si128();

    UInt64 res =

        (static_cast<UInt64>(_mm_movemask_epi8(_mm_cmpeq_epi8(

        _mm_loadu_si128(reinterpret_cast<const __m128i *>(bytes64)), zero16))) & 0xffff)

        | ((static_cast<UInt64>(_mm_movemask_epi8(_mm_cmpeq_epi8(

        _mm_loadu_si128(reinterpret_cast<const __m128i *>(bytes64 + 16)), zero16))) << 16) & 0xffff0000)

        | ((static_cast<UInt64>(_mm_movemask_epi8(_mm_cmpeq_epi8(

        _mm_loadu_si128(reinterpret_cast<const __m128i *>(bytes64 + 32)), zero16))) << 32) & 0xffff00000000)

        | ((static_cast<UInt64>(_mm_movemask_epi8(_mm_cmpeq_epi8(

        _mm_loadu_si128(reinterpret_cast<const __m128i *>(bytes64 + 48)), zero16))) << 48) & 0xffff000000000000);

#else

    UInt64 res = 0;

    const UInt8 * pos = bytes64;

    const UInt8 * end = pos + 64;

    for (; pos < end; ++pos)

        res |= ((*pos == 0)<<(pos-bytes64));

#endif

    return ~res;

}

按照优先级尽可能利用位数多的指令集进行计算，同时在所有指令集都不可用及未64字节对齐(align)的部分进行了保底处理。利用了以下指令集:

Ÿ   AVX512F / AVX512BW

Ÿ   AVX/AVX2

Ÿ   SSE2

 

参考链接：

https://zhuanlan.zhihu.com/p/449154820