有效移位或大比特量矢量

Question

我有一个大内存数组作为指针uint64_t * arr (加号)，它表示普通位.我需要非常有效地(大多数表演性/快速)将这些位元从0移到63。

通过移动整个数组，我的意思不是移动每个元素(如a[i] <<= Shift)，而是将其转换为单个大位向量。换句话说，对于每个中间位置i (除了第一个和最后一个元素)，我可以在一个循环中执行以下操作：

dst[i] = w | (src[i] << Shift);
w = src[i] >> (64 - Shift);

其中w是一些临时变量，保存前一个数组元素的右移值。

上述解决方案简单明了。但我需要更高效的东西，因为我有千兆字节的数据。

理想的情况是使用一些SIMD指令，所以我正在寻找来自专家的SIMD建议。我需要实现所有四种流行指令集的转换代码- SSE-SSE4.2 / AVX / AVX-2 / AVX-512。

但据我所知，例如，对于SSE2，只有si128()内禀/指令，它的移位量仅为8倍(换句话说，字节移动)。我需要任意大小的移位，而不仅仅是字节移位。

没有SIMD，我也可以通过使用shld reg, reg, reg指令同时移位128位，这允许进行128位移位。它是在MSVC中作为内部shiftleft128()实现的，并产生可以是在这里看到的的汇编程序代码。

顺便说一下，我需要所有MSVC/GCC/CLang的解决方案。

此外，在单循环迭代中，我可以在顺序操作中移动4或8个字，这将使用CPU流水线来加速几条指令的并行无序执行。

如果需要的话，我的位向量可以对齐到内存中的任意字节，如果这将有助于提高SIMD的速度，例如通过执行对齐读/写。此外，源和目标位矢量存储器是不同的(不重叠)。

换句话说，我正在寻找关于如何在不同的Intel CPU上最有效(最高效)地解决我的任务的所有建议。

Note，为了澄清，我实际上需要做几次换档，而不仅仅是单班。我有大位矢量X，以及几百个移位大小s0, s1, ..., sN，其中每个移位大小不同，也可以是大的(例如，100 K比特的移位)，然后我想要计算产生的大位矢量Y = (X << s0) | (X << s1) | ... | (X << sN)。我只是将StackOverflow的问题简化为移动单个向量。但是，这个关于原始任务的细节可能是非常重要的。

根据@Jake‘’Alquimista‘’LEE的要求，我决定实现一个现成的玩具最小可复制示例，计算输入位向量src的移位，以生成或编辑最终的dst位向量。这个例子根本没有优化，只是一个简单的简单变体如何解决我的任务。为了简单起见，这个例子的输入向量的大小很小，而不是像我这样的千兆字节。这是一个玩具示例，我没有检查它是否正确地解决了任务，它可能包含一些小错误：

在网上试试！

#include <cstdint>
#include <vector>
#include <random>

#define bit_sizeof(x) (sizeof(x) * 8)

using u64 = uint64_t;
using T = u64;

int main() {
    std::mt19937_64 rng{123};

    // Random generate source bit vector
    std::vector<T> src(100'000);
    for (size_t i = 0; i < src.size(); ++i)
        src[i] = rng();

    size_t const src_bitsize = src.size() * bit_sizeof(T);

    // Destination bit vector, for example twice bigger in size
    std::vector<T> dst(src.size() * 2);

    // Random generate shifts
    std::vector<u64> shifts(200);
    for (size_t i = 0; i < shifts.size(); ++i)
        shifts[i] = rng() % src_bitsize;

    // Right-shift that handles overflow
    auto Shr = [](auto x, size_t s) {
        return s >= bit_sizeof(x) ? 0 : (x >> s);
    };

    // Do actual Shift-Ors
    for (auto orig_shift: shifts) {
        size_t const
            word_off = orig_shift / bit_sizeof(T),
            bit_off = orig_shift % bit_sizeof(T);

        if (word_off >= dst.size())
            continue;
        
        size_t const
            lim = std::min(src.size(), dst.size() - word_off);

        T w = 0;
        
        for (size_t i = 0; i < lim; ++i) {
            dst[word_off + i] |= w | (src[i] << bit_off);
            w = Shr(src[i], bit_sizeof(T) - bit_off);
        }

        // Special case of handling for last word
        if (word_off + lim < dst.size())
            dst[word_off + lim] |= w;
    }
}

我的实际项目的当前代码与上面的玩具示例不同。这个项目已经正确地解决了一个现实世界的任务.我只需要做额外的优化。我已经做了一些优化，比如使用OpenMP并行化shift，或者在所有核上进行操作。另外，正如注释中所说的，我为每个移位大小创建了专门的模板函数，总共创建了64个函数，并选择了64个函数中的一个来执行实际的移位--或者。每个C++函数都有移位大小的编译时间值，因此编译器会在考虑编译时间值的情况下进行额外的优化。

Answer 1

您可以，而且可能您甚至不需要显式地使用SIMD指令。目标编译器GCC、CLANG和MSVC以及ICC等编译器都支持自动矢量化.虽然手工优化的程序集可以优于编译器生成的矢量化指令，但通常很难实现，您可能需要针对不同架构的几个版本。导致高效自动矢量化指令的通用代码是一种可以跨许多平台移植的解决方案。

例如，一个简单的shiftvec版本

void shiftvec(uint64_t* dst, uint64_t* src, int size, int shift)
{
    for (int i = 0; i < size; ++i,++src,++dst)
    {
        *dst = ((*src)<<shift) | (*(src+1)>>(64-shift));
    }
}

与最近的GCC编译(或CLANG工作，以及)和-O3 -std=c++11 -mavx2导致SIMD指令在核心循环的组装

.L5:
  vmovdqu ymm4, YMMWORD PTR [rsi+rax]
  vmovdqu ymm5, YMMWORD PTR [rsi+8+rax]
  vpsllq ymm0, ymm4, xmm2
  vpsrlq ymm1, ymm5, xmm3
  vpor ymm0, ymm0, ymm1
  vmovdqu YMMWORD PTR [rdi+rax], ymm0
  add rax, 32
  cmp rax, rdx
  jne .L5

见godbolt.org：https://godbolt.org/z/5TxhqMhnK

如果您想在dst中合并多个移位，这也是一个概括：

void shiftvec2(uint64_t* dst, uint64_t* src1, uint64_t* src2, int size1, int size2, int shift1, int shift2)
{
    int size = size1<size2 ? size1 : size2;
    for (int i = 0; i < size; ++i,++src1,++src2,++dst)
    {
        *dst = ((*src1)<<shift1) | (*(src1+1)>>(64-shift1));
        *dst |= ((*src2)<<shift2) | (*(src2+1)>>(64-shift2)); 
    }
    for (int i = size; i < size1; ++i,++src1,++dst)
    {
        *dst = ((*src1)<<shift1) | (*(src1+1)>>(64-shift1));        
    }
    for (int i = size; i < size2; ++i,++src2,++dst)
    {
        *dst = ((*src2)<<shift2) | (*(src2+1)>>(64-shift2));
    }
}

编译成核心循环：

.L38:
  vmovdqu ymm7, YMMWORD PTR [rsi+rcx]
  vpsllq ymm1, ymm7, xmm4
  vmovdqu ymm7, YMMWORD PTR [rsi+8+rcx]
  vpsrlq ymm0, ymm7, xmm6
  vpor ymm1, ymm1, ymm0
  vmovdqu YMMWORD PTR [rax+rcx], ymm1
  vmovdqu ymm7, YMMWORD PTR [rdx+rcx]
  vpsllq ymm0, ymm7, xmm3
  vmovdqu ymm7, YMMWORD PTR [rdx+8+rcx]
  vpsrlq ymm2, ymm7, xmm5
  vpor ymm0, ymm0, ymm2
  vpor ymm0, ymm0, ymm1
  vmovdqu YMMWORD PTR [rax+rcx], ymm0
  add rcx, 32
  cmp r10, rcx
  jne .L38

将多个源组合在一个循环中将减少用于加载/写入目标的内存带宽的总量。您可以组合的数量的限制当然受到可用寄存器的限制。请注意，xmm2和xmm3 for shiftvec包含shift值，因此使用不同版本的编译时已知的shift值可以释放这些寄存器。

此外，对每个指针使用__restrict (由GCC、CLANG、MSVC支持)将告诉编译器范围不重叠。

最初，我在MSVC提供适当的自动矢量化代码时遇到了问题，但似乎添加更多类似SIMD的结构将使其适用于所有三个所需的编译器GCC、CLANG和MSVC：

void shiftvec(uint64_t* __restrict dst, const uint64_t* __restrict src, int size, int shift)
{
    int i = 0;
    // MSVC: use steps of 2 for SSE, 4 for AVX2, 8 for AVX512
    for (; i+4 < size; i+=4,dst+=4,src+=4)
    {
        for (int j = 0; j < 4; ++j)
            *(dst+j) = (*(src+j))<<shift;
        for (int j = 0; j < 4; ++j)
            *(dst+j) |= (*(src+1)>>(64-shift));
    }
    for (; i < size; ++i,++src,++dst)
    {
        *dst = ((*src)<<shift) | (*(src+1)>>(64-shift));
    }    
}

Answer 2

我会试图依靠x64的能力，从未对齐的地址，并这样做，几乎没有可见的损失时，恒星是正确地(联合国)对齐。一个人只需要处理几种情况(shift % 8)或(shift % 16) --所有这些都可以用SSE2指令集完成，用零修正余数，并有一个与数据矢量不对齐的偏移量，并通过memcpy寻址UB。

也就是说，内部循环看起来是这样的：

uint16_t const *ptr;
auto a = _mm_loadu_si128((__m128i*)ptr);
auto b = _mm_loadu_si128((__m128i*)(ptr - 1);
a = _mm_srl_epi16(a, c);
b = _mm_sll_epi16(b, 16 - c);
_mm_storeu_si128((__m128i*)ptr, mm_or_si128(a,b));
ptr += 8;

将这个循环展开几次，您可能可以使用_mm_alignr_epi8 on SSE3+来放松内存带宽(以及那些需要结合来自非对齐内存访问的结果的管道阶段)：

auto a0 = w; 
auto a1 = _mm_load_si128(m128ptr + 1);
auto a2 = _mm_load_si128(m128ptr + 2);
auto a3 = _mm_load_si128(m128ptr + 3);
auto a4 = _mm_load_si128(m128ptr + 4);
auto b0 = _mm_alignr_epi8(a1, a0, 2);
auto b1 = _mm_alignr_epi8(a2, a1, 2);
auto b2 = _mm_alignr_epi8(a3, a2, 2);
auto b3 = _mm_alignr_epi8(a4, a3, 2);
// ... do the computation as above ...
w = a4;   // rotate the context

Answer 3

换句话说，我正在寻找关于如何在不同的Intel CPU上最有效(最高效)地解决我的任务的所有建议。

效率的关键是懒惰。懒惰的关键是撒谎--假装你没有做任何改变。

对于最初的例子(仅说明这一概念)，请考虑：

struct Thingy {
    int ignored_bits;
    uint64_t data[];
}

void shift_right(struct Thingy * thing, int count) {
    thing->ignored_bits += count;
}

void shift_left(struct Thingy * thing, int count) {
    thing->ignored_bits -= count;
}

int get_bit(struct Thingy * thing, int bit_number) {
    bit_number += thing->ignored_bits;
    return !!(thing->data[bit_number / 64] & (1 << bit_number % 64));
}

对于实际的代码，您需要关注各种细节--您可能希望从数组开头的备用位开始(和非零ignored_bits)，这样您就可以假装正确地移动；对于每一个小的移位，您可能都希望清除“移入”位(否则它的行为就像浮点--例如(5.0 << 8) >> 8) == 5.0)；如果/当ignored_bits超出某个范围时，您可能想要一个大的memcpy()；等等。

为了更多的乐趣，滥用低级内存管理--使用VirtualAlloc() (Windows)或mmap() (Linux)来预留一个巨大的空间，然后将数组放在空间的中间，然后根据需要在数组的开始/结束处分配/释放页面；这样您只需要在原始位被“移动”到左边/右边之后才能使用memcpy()。

当然，结果是这会使代码的其他部分变得复杂--例如，对于OR 2位字段，您将不得不做一个棘手的“获取A；移位A以匹配B；结果=A或B”的调整。这并不是破坏业绩的关键。