如何高效地并行设置位向量中的位？

Question

考虑一个N位的比特向量（N很大）和一个数字数组M（M适中，通常远小于N），每个数字在范围0..N-1内，指示向量的哪个比特必须设置为1。后面的数组未排序。比特向量只是整数数组，具体来说是__m256i，其中256个比特被打包进每个__m256i结构中。

如何有效地将这项工作分割到多个线程中？

首选语言为C++（MSVC++2017工具集v141），汇编也可以。首选CPU为x86_64（使用内部函数也可以）。如果有任何益处，则希望使用AVX2。

- Serge Rogatch

---

嗯...看起来主要是内存带宽的问题。我不确定是否真的有比显而易见的方法更好的方式。一种方法可能是先对数组进行排序，这样你就可以按顺序设置位，从而使缓存更加高效。 - fuz

M已经排序了吗？如果没有，你几乎肯定想要为单个线程进行优化。 - zzxyz

使用一些算法和典型数据来衡量性能...请展示你的代码。通过位向量，你是指std::bitset还是std::vector<bool>或其他什么东西。另请参阅：如何证明std::bitset比std::vector<bool>更快？。如果你的数据尚未排序且非常大，则很难进行优化。同时避免过早优化。只有在你可以证明显而易见的方法不足以解决问题时才需要优化。对于小数据量，线程或复杂算法的开销会使代码变慢。 - Phil1970

如果数组未排序，请考虑使用二叉搜索树。这样，您就不必进行任何内存地址算术或位移操作；只需直接使用位数即可。 - prl

我同意@prl的观点，如果M远小于N，则争用会很少。但是，如果更新被排序并组成一个掩码“m”，该掩码被“ored”到位置“w = w | m;”中，并且（硬核心）如果T个线程索引为t [0 ... T-1]，则每个线程都被赋予缓存行，使得addr％T == t。因此，每个线程都在自己拥有的行中查找更新。我并不是说这两种方法都会提高性能，因为存在开销且M << N，但我只是指出，在规模上，您希望切分目标空间（N），而不是域空间（M）。对齐或过度对齐也可能有所帮助。 - Persixty

Answer 1

@IraBaxter发布了一个有趣但存在缺陷的想法，可以通过付出巨大代价使其实现。我怀疑@BeeOnRope关于对M数组进行部分排序/分区的想法会表现得更好（特别是对于具有大型私有高速缓存的CPU来说，它们可以保持N的某些部分热度）。我将总结我在他删除的答案评论中描述的修改版本。（该答案提供了一些关于何时值得多线程处理的建议，需要N有多大。）

---

每个写线程都会获得一块未排序/分区的M。

这个想法是，冲突非常罕见，因为N相对于可以同时进行的存储数量来说很大。由于设置位是幂等的，所以我们可以通过检查内存中的值来处理冲突（当两个线程想要在同一个字节中设置不同的位时），以确保我们想要的位实际上已经被设置了，如RMW操作：or [N + rdi], al（没有lock前缀）。

例如，线程1尝试存储0x1并覆盖了线程2的存储0x2。线程2必须注意并重试读取-修改-写入操作（可能使用lock or使其简单化并使多次重试变得不可能），最终在冲突字节中得到0x3。

我们需要在读取之前加上一个“mfence”指令。否则，存储转发将会给我们刚刚写入的值，在其他线程看到我们的存储之前。换句话说，一个线程可以比它们在全局顺序中出现更早地观察到自己的存储。x86确实有存储的总顺序，但没有加载的总顺序。因此，我们需要“mfence”来防止存储-加载重排序。（英特尔的“加载不会与较旧的存储重新排序到相同位置”保证并不像听起来那么有用：存储/重新加载不是内存屏障；他们只是谈论乱序执行保留程序顺序语义。）

mfence是昂贵的，但比仅使用lock or [N+rdi], al更好的技巧在于我们可以批量操作。例如，执行32个or指令，然后进行32次读取操作。这是mfence开销与每个操作之间的平衡，还是增加了误共享（读取已被另一个CPU声明的缓存行）的机会。

我们可以将一组中的最后一个or作为lock or而不是实际的mfence指令。这对于AMD和Intel的吞吐量都更好。例如，根据Agner Fog's tables，在Haswell/Skylake上，mfence的吞吐量为33c中的一个，而lock add（与or具有相同的性能）的吞吐量为18c或19c。或者对于Ryzen，约为70c（mfence）vs.约17c（lock add）。

如果我们将每个围栏的操作数量保持得非常低，那么数组索引(m[i]/8) + 掩码 (1<<(m[i] & 7)) 可以在所有操作中都保持在寄存器中。这可能不值得；因为围栏的代价太高，不能像每6个or操作一样经常进行。使用bts和bt位串指令意味着我们可以将更多的索引保持在寄存器中(因为不需要移位结果)，但可能不值得，因为它们很慢。

使用向量寄存器来保存索引可能是一个好主意，以避免在屏障后重新从内存加载它们。我们希望加载地址在读回加载uops可以执行时就准备好(因为它们正在等待屏障前的最后一个存储提交到L1D并变为全局可见)。

使用单字节读取-修改-写入使实际冲突尽可能不太可能发生。每个字节的写入只对相邻的7个字节进行非原子RMW操作。当两个线程在同一个64B缓存行中修改字节时，性能仍会受到false-sharing的影响，但至少我们避免了必须重新执行许多or操作。32位元素大小会使某些事情更有效率（例如使用xor eax, eax / bts eax, reg生成只有2个uops或1个BMI2 shlx eax，r10d，reg（其中r10d = 1）的1<<（m [i] & 31））。

避免使用类似于bts [N]，eax的位串指令：它的吞吐量比为or [N + rax]，dl进行索引和掩码计算要差。这是它的“完美”用例（除了我们不关心内存中位的旧值，我们只想设置它），但它的CISC负担太大。

在C语言中，函数可能如下所示：

/// UGLY HACKS AHEAD, for testing only.

//    #include <immintrin.h>
#include <stddef.h>
#include <stdint.h>
void set_bits( volatile uint8_t * restrict N, const unsigned *restrict M, size_t len)
{
    const int batchsize = 32;

    // FIXME: loop bounds should be len-batchsize or something.
    for (int i = 0 ; i < len ; i+=batchsize ) {
        for (int j = 0 ; j<batchsize-1 ; j++ ) {
           unsigned idx = M[i+j];
           unsigned mask = 1U << (idx&7);
           idx >>= 3;
           N[idx] |= mask;
        }

        // do the last operation of the batch with a lock prefix as a memory barrier.
        // seq_cst RMW is probably a full barrier on non-x86 architectures, too.
        unsigned idx = M[i+batchsize-1];
        unsigned mask = 1U << (idx&7);
        idx >>= 3;
        __atomic_fetch_or(&N[idx], mask, __ATOMIC_SEQ_CST);
        // _mm_mfence();

        // TODO: cache `M[]` in vector registers
        for (int j = 0 ; j<batchsize ; j++ ) {
           unsigned idx = M[i+j];
           unsigned mask = 1U << (idx&7);
           idx >>= 3;
           if (! (N[idx] & mask)) {
               __atomic_fetch_or(&N[idx], mask, __ATOMIC_RELAXED);
           }
        }
    }
}

这在gcc和clang中大致符合我们的要求。汇编代码（Godbolt）可以通过多种方式更加高效，但尝试这个可能会有趣。这不是安全的：我只是在C语言中拼凑出了这个独立函数所需的汇编代码，而没有将其内联到调用者或其他任何地方。__atomic_fetch_or并不是非原子变量的适当编译器屏障，就像asm("":::"memory")一样。（至少C11的stdatomic版本不是。）我应该使用传统的__sync_fetch_and_or，它对所有内存操作都是一个完整的屏障。

它使用GNU C的原子内建函数（atomic builtins）来在非atomic_uint8_t变量上执行原子RMW操作。同时从多个线程运行此函数将会产生C11 UB，但我们只需要它在x86上工作。我使用volatile来获取atomic中允许异步修改的部分，而不强制N[idx] |= mask;成为原子操作。这样做的目的是确保读取检查不被优化掉。

我使用__atomic_fetch_or作为内存屏障，因为我知道它会在x86上工作。虽然使用seq_cst，它在其他ISA上也可能会工作，但这都是一种大型的hack。

- Peter Cordes

Answer 2

假设您想将此工作分配给T个线程。这是一个非常有趣的问题，因为它不能通过分区轻松并行化，不同大小的N和M可能需要不同的解决方案。

完全并发基准

您可以简单地将数组M分成T个部分，并让每个线程在其自己的M分区上使用共享的N。主要问题是，由于M未排序，所有线程都可以访问N的任何元素，因此会干扰彼此的工作。为了避免这种情况，您必须对共享的N数组的每次修改使用原子操作，例如std::atomic::fetch_or，否则就得想出一些锁定方案。这两种方法都可能降低性能（即使用原子操作设置位可能比等效的单线程代码慢一个数量级）。

让我们看看可能更快的思路。

私有N

避免“共享N”问题的一个相对明显的想法是，简单地为每个T提供一个私有副本的N，并通过or在最后合并它们。

不幸的是，这个解决方案的复杂度为O(N) + O(M/T)，而原始的单线程解决方案的复杂度为O(M)，上面的“原子”解决方案大约为O(M/T)⁴。由于我们知道N >> M，因此在这种情况下，这可能是一个糟糕的权衡。但是，值得注意的是，每个项中的隐藏常数非常不同：来自合并步骤⁰的O(N)项可以使用256位宽的vpor指令，这意味着吞吐量接近200-500位/周期（如果缓存），而设置位的步骤，即O(M/T)，我估计更接近1位/周期。因此，即使N的大小是M的10或100倍，对于适度的T，这种方法肯定可以是最佳的。

M的分区

基本思路是将M中的索引分成几部分，以便每个工作线程可以在N数组的不同部分上工作。如果M排序了，那就很简单了，但它没有排序，所以...

一个简单的算法，如果M分布得很好，将会很有效，首先将M的值划分为T个桶，桶中的值在范围[0, N/T), [N/T, 2N/T], ..., [(T-1)N/T, N)内。也就是说，将N分为T个不相交的区域，然后找到落入每个区域的M的值。您可以通过将每个线程分配一个相等大小的M块，并让它们每个创建T个分区，然后在最后逻辑合并¹它们，以便您拥有M的T个分区来分散这项工作。

第二步是实际设置所有位：您将一个分区分配给每个线程T，它可以以“单线程”方式设置位，即不必担心并发更新，因为每个线程都在处理N的不相交分区²。

两个步骤O(M)，第二步与单线程情况相同，因此并行化的开销是第一步。我怀疑第一步的速度范围将从与第二步相同到可能慢2-4倍，具体取决于实现和硬件，因此您可以期望在具有许多核心的机器上加速，但只有2或4个核心时可能不会更好。

如果M的分布不是平滑的，例如第一步创建的分区大小差异很大，它将表现不佳，因为某些线程将获得更多的工作。一个简单的策略是创建大约10 * T个分区，而不仅仅是T，并且在第二遍中，所有线程都从相同的分区队列中消耗，直到完成。通过这种方式，您可以更均匀地分散工作，除非数组M非常紧密。在这种情况下，您可能需要考虑第一步的改进，该步骤首先基本上创建元素的分桶直方图，然后进行减少阶段，该阶段查看组合的直方图以创建良好的分区。

本质上，我们只是逐渐将第一阶段精细化为一种并行排序/分区算法，对此已经有很多文献。您甚至可能发现完整的（并行）排序最快，因为它将极大地帮助位设置阶段，因为访问将有序并具有最佳的空间局部性（分别帮助预取和缓存）。

¹ 合并的最简单的概念形式是将每个线程对M的分区复制，以便您拥有所有M的连续分区，但实际上，如果分区很大，您可以将分区保留在原处并将它们链接在一起，这会增加消费代码的复杂性，但避免了紧凑步骤。

² 为了使它从线程角度真正不相交，您需要确保N的分区落在“字节边界”上，甚至可能是缓存行边界，以避免错误共享（尽管后者可能不是一个大问题，因为它仅发生在每个分区的边缘，并且处理顺序意味着您不太可能出现争用）。

⁴ 实际上，使用共享N的基线并发解决方案的确切“顺序”很难定义，因为会有争用，因此O（M/T）扩展将在足够大的T的情况下崩溃。如果我们假设N相当大，而T限于最多十几个核心的典型硬件并发性，那么这可能是一个可以接受的近似值。

- BeeOnRope

---

如果在循环外部有一个初始化为1的寄存器，则可以使用shlx代替xor和bts。 - BeeOnRope

这可以解释为存储转发。如果读/写现在是8字节，则下一次迭代的读取会命中上一次迭代的存储。虽然在我的心理模式中实际上没有任何存储转发，因为锁定操作的隐含栅栏不应允许后续加载在SB为空之前进行，但是谁知道它在实践中如何运作。一堆连续的原子操作并不是很常见。 - BeeOnRope

我尝试了使用times 10 imul ecx，ecx并注释掉（或不注释）lock or块。 差异（如果有的话）低于测量噪声水平，在25M个迭代中为约750.4Mc。 - Peter Cordes

是的，英特尔明确表示HT 静态地 分区存储缓冲区，因此每个逻辑线程都有自己的缓冲区。(https://dev59.com/lIbca4cB1Zd3GeqPYaRx#27902942) - Peter Cordes

尝试在相同地址情况下使用“lock or”的好建议。它不会变慢，仍然是18c。这使得“lock bts”看起来更奇怪。 - Peter Cordes

Answer 3

集合中涉及到几种操作（A，B=集合，X=集合中的元素）：

Set operation           Instruction
---------------------------------------------
Intersection of A,B     A and B
Union of A,B            A or B
Difference of A,B       A xor B
A is subset of B        A and B = B     
A is superset of B      A and B = A       
A <> B                  A xor B <> 0
A = B                   A xor B = 0
X in A                  BT [A],X
Add X to A              BTS [A],X
Subtract X from A       BTC [A],X

鉴于您可以使用布尔运算符替换集合操作，因此可以使用VPXOR、VPAND等。
要设置、重置或测试单个位，只需使用

mov eax,BitPosition
BT [rcx],rax

您可以使用以下代码设置一个集合是否为空（或其他内容）

vpxor      ymm0,ymm0,ymm0       //ymm0 = 0
//replace the previous instruction with something else if you don't want
//to compare to zero.
vpcmpeqqq  ymm1,ymm0,[mem]      //compare mem qwords to 0 per qword
vpslldq    ymm2,ymm1,8          //line up qw0 and 1 + qw2 + 3
vpand      ymm2,ymm1,ymm2       //combine qw0/1 and qw2/3
vpsrldq    ymm1,ymm2,16         //line up qw0/1 and qw2/3
vpand      ymm1,ymm1,ymm2       //combine qw0123, all in the lower 64 bits.
//if the set is empty, all bits in ymm1 will be 1.
//if its not, all bits in ymm1 will be 0.     

（我相信使用混合/聚集等指令可以改进此代码） 从这里，您只需扩展到更大的集合或其他操作即可。

请注意，带有内存操作数的bt、btc、bts不限于64位。
以下内容将正常工作。

mov eax,1023
bts [rcx],rax   //set 1024st element (first element is 0).

- Johan

---

问题在于如何高效地并行（多线程）地将位设置为“1”，给定一个要设置为“1”的位索引数组（并保留其他位不变）。 - Serge Rogatch

“and's”和“or's”是你的好朋友，就像上面所详细说明的那样。 - Johan