在没有硬件SIMD的情况下，通过SWAR并行将64位整数中打包的8位整数减1。

Question

在没有硬件SIMD的情况下，通过SWAR并行将64位整数中打包的8位整数减1。

80

如果我有一个64位整数，我把它解释为一个包含8个元素的紧密打包的8位整数数组。我需要从每个打包的整数中减去常量1，同时处理溢出，而不是让一个元素的结果影响另一个元素的结果。

目前我有这个代码，它可以工作，但我需要一种并行地执行每个打包的8位整数减法且不进行内存访问的解决方案。在x86上，我可以使用类似psubb的SIMD指令，以并行方式减去打包的8位整数，但我编写代码的平台不支持SIMD指令（在这种情况下是RISC-V）。

因此，我正在尝试手动实现SWAR（寄存器内的SIMD），取消掉uint64_t字节之间的进位传递，做类似于以下内容的操作：

uint64_t sub(uint64_t arg) {
    uint8_t* packed = (uint8_t*) &arg;

    for (size_t i = 0; i < sizeof(uint64_t); ++i) {
        packed[i] -= 1;
    }

    return arg;
}

我认为您可以使用位运算符完成这个问题，但我不确定。我正在寻找一种不使用SIMD指令的解决方案。我正在寻找一种在C或C++中非常可移植的解决方案，或者只是其背后的理论，以便我可以实现自己的解决方案。

- cam-white

5

它们需要是8位的还是可以改成7位呢？ - tadman

它们必须是8位的，很抱歉:( - cam-white

12

这种技术通常称为SWAR。 - harold

1

相关：通过分段比较64位整数 - phuclv

1

你是否期望一个包含零的字节会被转换为0xff？ - Alnitak

显示剩余5条评论

8个回答

17

对于 RISC-V，你可能会使用 GCC/clang。

有趣的事实是：GCC 知道一些 SWAR 位运算技巧（在其他答案中展示），并且可以在为没有硬件 SIMD 指令的目标编译代码时，使用 GNU C 本地向量为您执行这些技巧。（但是针对 RISC-V 的 clang 将只是天真地展开它到标量操作，因此如果要在不同的编译器上获得良好的性能，您必须自行完成）。

本地向量语法的一个优点是，在针对具有硬件 SIMD 的机器进行目标定位时，它将使用该硬件而不是自动向量化您的位运算或类似于此的可怕东西。

它使编写 vector -= scalar 操作变得容易；语法工作正常，隐式广播（也称为 splatting）标量。

此外，请注意，从 uint8_t 数组[] 加载 uint64_t* 是严格别名 UB，因此请小心处理。（另请参见为什么 glibc 的 strlen 需要如此复杂才能快速运行？关于在纯 C 中使 SWAR 位运算技巧严格别名安全的内容）。您可能需要像这样声明一个 uint64_t，以便将其指针转换为访问任何其他对象，就像 ISO C / C++ 中的 char* 一样。

使用以下内容将 uint8_t 数据转换为 uint64_t，以便与其他答案一起使用：

// GNU C: gcc/clang/ICC but not MSVC
typedef uint64_t  aliasing_u64 __attribute__((may_alias));  // still requires alignment
typedef uint64_t  aliasing_unaligned_u64 __attribute__((may_alias, aligned(1)));

另一种安全地进行别名加载的方法是使用 memcpy 到一个 uint64_t 中，这也消除了 alignof(uint64_t) 对齐要求。但对于没有有效未对齐加载的 ISA，gcc/clang 在无法证明指针对齐时不会内联和优化掉 memcpy，这将对性能产生灾难性的影响。

简而言之，最好的方法是将您的数据声明为 uint64_t array[...] 或动态分配为 uint64_t，或者更好的是 alignas(16) uint64_t array[];。这确保至少对齐到8个字节，如果您指定了 alignas，则对齐到16个字节。

由于 uint8_t 几乎肯定是指向 unsigned char*，所以可以通过 uint8_t* 安全访问 uint64_t 的字节（但对于 uint8_t 数组来说不能反过来访问）。因此，在这种特殊情况下，窄元素类型为 unsigned char，您可以绕开严格别名问题，因为 char 是特殊的。

GNU C 本地向量的语法示例：

GNU C 本地向量始终允许与其基础类型别名（例如，int __attribute__((vector_size(16))) 可以安全地别名为 int，但不能是 float 或 uint8_t 或其他任何类型）。

#include <stdint.h>
#include <stddef.h>

// assumes array is 16-byte aligned
void dec_mem_gnu(uint8_t *array) {
    typedef uint8_t v16u8 __attribute__ ((vector_size (16), may_alias));
    v16u8 *vecs = (v16u8*) array;
    vecs[0] -= 1;
    vecs[1] -= 1;   // can be done in a loop.
}

对于没有任何硬件SIMD的RISC-V，您可以使用vector_size（8）来表达您可以有效使用的粒度，并做两倍数量的较小向量。

但是，在x86上，vector_size（8）无论使用GCC还是clang都编译得很愚蠢：GCC在GP-integer寄存器中使用SWAR位操作，clang则拆分为2字节元素以填充16字节的XMM寄存器，然后重新打包。（MMX已经过时，至少在x86-64上GCC / clang甚至不会使用它。）

但是，对于vector_size（16）（Godbolt），我们得到了预期的movdqa / paddb。（由pcmpeqd same，same生成的全1矢量）。使用-march = skylake 仍然会得到两个单独的XMM操作，而不是一个YMM，所以目前的编译器也不会将矢量操作“自动矢量化”为更宽的矢量：/

对于AArch64，使用vector_size（8）并不太糟糕（Godbolt）； ARM / AArch64可以使用d或q寄存器以8或16字节块的方式进行本地处理。

因此，如果您希望在x86，RISC-V，ARM / AArch64和POWER之间具有可移植性的性能，则可能需要使用vector_size（16）进行编译。 但是，其他一些ISA在64位整数寄存器内执行SIMD，例如MIPS MSA。

vector_size(8)使得查看汇编更容易（只有一个寄存器的数据）：

# GCC8.2 -O3 for RISC-V for vector_size(8) and only one vector

dec_mem_gnu(unsigned char*):
        lui     a4,%hi(.LC1)           # generate address for static constants.
        ld      a5,0(a0)                 # a5 = load from function arg
        ld      a3,%lo(.LC1)(a4)       # a3 = 0x7F7F7F7F7F7F7F7F
        lui     a2,%hi(.LC0)
        ld      a2,%lo(.LC0)(a2)       # a2 = 0x8080808080808080
                             # above here can be hoisted out of loops
        not     a4,a5                  # nx = ~x
        and     a5,a5,a3               # x &= 0x7f... clear high bit
        and     a4,a4,a2               # nx = (~x) & 0x80... inverse high bit isolated
        add     a5,a5,a3               # x += 0x7f...   (128-1)
        xor     a5,a4,a5               # x ^= nx  restore high bit or something.

        sd      a5,0(a0)               # store the result
        ret

我认为这与其他非循环答案的基本思路相同；避免进位然后修正结果。

我认为这是5个ALU指令，比顶部答案更差。但是它看起来关键路径延迟只有3个周期，有两个由2个指令组成的链条导致XOR。@Reinstate Monica - ζ--的答案编译为4个周期的依赖链（对于x86）。 5个周期的循环吞吐量也受到瓶颈的影响，包括一个天真的sub在关键路径上，并且循环会受到延迟的限制。

然而，使用clang这是无用的。它甚至不按照加载的顺序添加和存储，因此它甚至没有进行良好的软件流水线处理！

# RISC-V clang (trunk) -O3
dec_mem_gnu(unsigned char*):
        lb      a6, 7(a0)
        lb      a7, 6(a0)
        lb      t0, 5(a0)
...
        addi    t1, a5, -1
        addi    t2, a1, -1
        addi    t3, a2, -1
...
        sb      a2, 7(a0)
        sb      a1, 6(a0)
        sb      a5, 5(a0)
...
        ret

- Peter Cordes

13

我要指出的是，一旦您开始处理超过单个uint64_t的内容，您编写的代码实际上会进行向量化。 https://godbolt.org/z/J9DRzd

- robthebloke

1

你能解释一下或者提供一个参考，说明那里正在发生什么吗？这看起来很有趣。 - n314159

2

我一开始尝试不使用SIMD指令来完成这个，但我还是觉得这很有趣 :) - cam-white

8

另一方面，SIMD 代码糟糕透了。编译器完全误解了这里发生的情况。即：这是一个“显然是由编译器完成的示例，因为没有人会如此愚蠢”。 - harold

1

@PeterCordes：我更倾向于使用__vector_loop(index, start, past, pad)结构，实现可以将其视为for(index=start; index<past; index++)[这意味着任何实现都可以通过定义宏来处理使用它的代码]，但是它的语义会更加宽松，以便编译器可以按照任何2的幂次块大小处理事物，最多扩展到pad，如果起始点和结束点不是块大小的倍数，则向下扩展起始点并向上扩展结束点。每个块内的副作用将不被排序，如果循环中发生break，则其他重复... - supercat

1

@PeterCordes：虽然restrict很有帮助（如果标准认可“至少可能基于”的概念，并直接明确定义“基于”和“至少可能基于”，而不是使用愚蠢且无法工作的边角情况，那么它将更有帮助），但我的建议还允许编译器执行比请求更多的循环次数——这将极大地简化向量化，但标准没有做出任何规定。 - supercat

显示剩余5条评论

11

您可以确保减法不会溢出，然后修复高位：

uint64_t sub(uint64_t arg) {
    uint64_t x1 = arg | 0x80808080808080;
    uint64_t x2 = ~arg & 0x80808080808080;
    // or uint64_t x2 = arg ^ x1; to save one instruction if you don't have an andnot instruction
    return (x1 - 0x101010101010101) ^ x2;
}

- Falk Hüffner

我认为它适用于一个字节的所有256个可能值; 我将其放在Godbolt上（使用RISC-V clang）https://godbolt.org/z/DGL9aq，以查看各种输入（如0x0、0x7f、0x80和0xff（移位到数字中间））的常量传播结果。看起来很好。我认为顶部答案归结为同样的东西，但是它用更复杂的方式解释了它。 - Peter Cordes

编译器在这里构建常量寄存器时可以做得更好。clang花费了很多指令来构建splat(0x01)和splat(0x80)，而不是通过移位从另一个中获取一个。即使在源代码https://godbolt.org/z/6y9v-u中以这种方式编写它也不能手把手地引导编译器生成更好的代码；它只是进行常量传播。 - Peter Cordes

我想知道为什么它不直接从内存中加载常量；这是 Alpha 编译器（一种类似的架构）所做的。 - Falk Hüffner

RISC-V的GCC确实会从内存中加载常量。看起来clang需要进行一些调整，除非数据缓存未命中是预期的，并且与指令吞吐量相比代价高昂。（自Alpha以来，这种平衡肯定已经发生了变化，而且RISC-V的不同实现可能也不同。如果编译器意识到它是一个重复模式，他们可以在开始时进行移位/OR扩展，以便使用一个LUI/add获取20 + 12 = 32位立即数据。AArch64的位模式立即数甚至可以将这些作为AND/OR/XOR的立即数，智能解码与密度选择） - Peter Cordes

添加了一个答案，展示了GCC的本地向量SWAR用于RISC-V。 - Peter Cordes

7

不确定这是否符合您的需求，但它可以并行地执行8个减法：

#include <cstdint>

constexpr uint64_t mask = 0x0101010101010101;

uint64_t sub(uint64_t arg) {
    uint64_t mask_cp = mask;
    for(auto i = 0; i < 8 && mask_cp; ++i) {
        uint64_t new_mask = (arg & mask_cp) ^ mask_cp;
        arg = arg ^ mask_cp;
        mask_cp = new_mask << 1;
    }
    return arg;
}

解释：位掩码在每个8位数字中都以1开头。我们将其与参数进行异或运算。如果此处有1，则减去1并停止。这是通过在new_mask中将相应的位设置为0来完成的。如果此处为0，则将其设置为1并进行进位，因此该位保持为1并且我们将掩码向左移动。最好自己检查新掩码的生成是否按预期工作，我认为是这样，但第二个意见也不会有坏处。

附注：实际上，我不确定循环中对mask_cp的非空检查是否会减慢程序速度。如果没有它，代码仍然正确（因为0掩码什么也不做），并且编译器更容易执行循环展开。

- n314159

for不会并行运行，你是不是混淆了for_each？ - LTPCGO

3

@LTPCGO 不，我的意图不是将这个for循环并行化，这实际上会破坏算法。但是这段代码是在64位整数中并行处理不同的8位整数，也就是所有8个减法操作同时进行，但它们需要最多8个步骤。 - n314159

我意识到我的要求可能有点不合理，但这已经非常接近我所需要的了，谢谢 :) - cam-white

4

int subtractone(int x) 
{
    int f = 1; 

    // Flip all the set bits until we find a 1 at position y
    while (!(x & f)) { 
        x = x^f; 
        f <<= 1; 
    } 

    return x^f; // return answer but remember to flip the 1 at y
}

您可以使用上述位运算来完成此操作，您只需要将整数分成8位一组，8次发送到此函数中即可。以下内容摘自如何将一个64位数字拆分为八个8位的值？，我在其中添加了上述函数。

uint64_t v= _64bitVariable;
uint8_t i=0,parts[8]={0};
do parts[i++] = subtractone(v&0xFF); while (v>>=8);

不管是怎样接触到这个技术，它都是有效的C或C++代码。

- LTPCGO

5

这并没有实现工作的并行化，这是原帖作者的问题。 - nickelpro

是的，@nickelpro说得对，这会一个接一个地执行每个减法，我想同时减去所有8位整数。不过还是感谢你的回答，谢谢兄弟。 - cam-white

2

@nickelpro 当我开始回答时，编辑还没有被做出来 其中包括问题的并行部分，所以我在提交后才注意到它，将其保留以防对其他人有用，因为它至少回答了位运算的部分，并且可以通过利用 for_each(std::execution::par_unseq,... 来并行工作，而不是使用 while 循环。 - LTPCGO

2

我的错，我提交了问题后才意识到我没有说需要并行处理，所以进行了编辑。 - cam-white

2

不打算编写代码，但对于减1操作，您可以通过减少一组8个1并检查结果的LSB是否已经“翻转”来实现。任何未切换的LSB都表示从相邻的8位发生了进位。应该可以计算出一系列AND / OR / XOR序列来处理此问题，而无需使用任何分支。

- Hot Licks

这可能有效，但要考虑一种情况，即进位传播穿过一个8位组并进入另一个8位组。在好的答案中（首先设置MSB或其他内容），确保不会传播进位的策略可能至少与此相同效率。目前要超越的目标（即良好的非循环无分支答案）是5个RISC-V asm ALU指令，指令级并行性使关键路径仅为3个周期，并使用两个64位常量。 - Peter Cordes

0

专注于每个字节的工作，然后将其放回原处。

uint64_t sub(uint64_t arg) {
   uint64_t res = 0;

   for (int i = 0; i < 64; i+=8) 
     res += ((arg >> i) - 1 & 0xFFU) << i;

    return res;
   }

- user12450543

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- nanofarad · Accepted Answer

如果您拥有具有高效SIMD指令的CPU，SSE/MMX paddb (_mm_add_epi8)也是可行的。Peter Cordes' answer还描述了GNU C (gcc/clang)向量语法，并且对于严格别名 UB 的安全性进行了说明。我强烈建议您也查看该答案。

使用uint64_t自己完成完全可移植，但仍需要小心以避免访问uint64_t*中的uint8_t数组时出现对齐问题和严格别名 UB。您通过在uint64_t中开始使用数据而将该部分省略，但对于GNU C，may_alias typedef可以解决该问题（请参见Peter的答案或memcpy）。

否则，您可以将数据分配/声明为uint64_t，并在需要单个字节时通过uint8_t*访问它。 unsigned char*允许别名任何内容，因此针对8位元素的特定情况可以避开这个问题。（如果根本不存在uint8_t，那么很可能可以安全地假设它是一个unsigned char。）

请注意，这是从之前错误的算法更改而来（请参见修订历史记录）。

对于任意减法，不需要循环即可实现，并且在每个字节中使用已知常数（如1）会更有效率。主要技巧是通过设置高位来防止每个字节的进位，然后纠正减法结果。

我们将稍微优化给出的减法技巧here。它们定义：

SWAR sub z = x - y
    z = ((x | H) - (y &~H)) ^ ((x ^~y) & H)

将H定义为0x8080808080808080U（即每个打包整数的MSB）。对于递减，y是0x0101010101010101U。

我们知道y的所有MSB都已清除，因此我们可以跳过一个掩码步骤（即y&~H在我们的情况下与y相同）。计算如下：

我们将x的每个组件的MSB设置为1，这样借位就无法传播到下一个组件。称之为调整后的输入。
我们从矫正后的输入中减去每个分量的1，通过减去0x01010101010101来实现。由于步骤1，这不会导致分量间的借位。称之为调整后的输出。
我们现在需要更正结果的MSB。我们将调整后的输出与原始输入的反转MSB异或以完成修复结果。

该操作可写为：

#define U64MASK 0x0101010101010101U
#define MSBON 0x8080808080808080U
uint64_t decEach(uint64_t i){
      return ((i | MSBON) - U64MASK) ^ ((i ^ MSBON) & MSBON);
}

最好由编译器内联（使用编译器指令强制执行），或将表达式作为另一个函数的一部分内联编写。

测试用例：

in:  0000000000000000
out: ffffffffffffffff

in:  f200000015000013
out: f1ffffff14ffff12

in:  0000000000000100
out: ffffffffffff00ff

in:  808080807f7f7f7f
out: 7f7f7f7f7e7e7e7e

in:  0101010101010101
out: 0000000000000000

性能细节

这里是单次调用函数的x86_64汇编代码。为了获得更好的性能，应该将其内联，希望常量尽可能长时间地存储在寄存器中。在常量存储在寄存器中的紧密循环中，实际的减少需要五个指令：优化后的或+非+与+加+xor。我没有看到比编译器优化更好的替代方案。

uint64t[rax] decEach(rcx):
    movabs  rcx, -9187201950435737472
    mov     rdx, rdi
    or      rdx, rcx
    movabs  rax, -72340172838076673
    add     rax, rdx
    and     rdi, rcx
    xor     rdi, rcx
    xor     rax, rdi
    ret

通过对以下代码片段进行IACA测试：

// Repeat the SWAR dec in a loop as a microbenchmark
uint64_t perftest(uint64_t dummyArg){
    uint64_t dummyCounter = 0;
    uint64_t i = 0x74656a6d27080100U; // another dummy value.
    while(i ^ dummyArg) {
        IACA_START
        uint64_t naive = i - U64MASK;
        i = naive + ((i ^ naive ^ U64MASK) & U64MASK);
        dummyCounter++;
    }
    IACA_END
    return dummyCounter;
}

我们可以展示，在Skylake机器上，执行减法、异或和比较跳转操作只需要不到5个时钟周期：

Throughput Analysis Report
--------------------------
Block Throughput: 4.96 Cycles       Throughput Bottleneck: Backend
Loop Count:  26
Port Binding In Cycles Per Iteration:
--------------------------------------------------------------------------------------------------
|  Port  |   0   -  DV   |   1   |   2   -  D    |   3   -  D    |   4   |   5   |   6   |   7   |
--------------------------------------------------------------------------------------------------
| Cycles |  1.5     0.0  |  1.5  |  0.0     0.0  |  0.0     0.0  |  0.0  |  1.5  |  1.5  |  0.0  |
--------------------------------------------------------------------------------------------------

（当然，在x86-64上，您只需加载或使用movq将其加载到XMM寄存器中进行paddb，因此更有趣的是看看它在像RISC-V这样的ISA上如何编译。）