优化从子位掩码生成父位掩码

Question

给定一个 64 位子掩码输入，例如：

10000000 01000000 00100000 00010000 00001000 00000100 00000010 00000000

8 位父掩码为：

11111110

父掩码中的单个位映射到子掩码字符串中的 8 位，并且当 8 个子位中的一个设置为 1 时，父掩码中的位设置为 1。计算此值的简单算法将如下：

unsigned __int64 childMask = 0x8040201008040200; // The number above in hex
unsigned __int8 parentMask = 0;
for (int i = 0; i < 8; i++)
{
    const unsigned __int8 child = childMask >> (8 * i);
    parentMask |= (child > 0) << i;
}

我想知道上面的代码中是否还有任何优化要做。代码将在 CUDA 上运行，我希望尽可能避免使用分支。对于答案，C++/C 中的代码就可以了。 for 循环可以展开，但我宁愿将其留给编译器进行优化，在必要时使用#pragma unroll 提供提示。

Answer 1

一种可能的方法是使用__vcmpgtu4 进行逐字节比较，将结果作为打包掩码返回，可以与 0x08040201（高半部分为 0x80402010）进行与运算以将它们转换为最终结果，但随后需要对它们进行水平求和，这似乎没有得到很好的支持，但可以使用普通的旧 C 样式代码来完成。

例如，

unsigned int low = childMask;
unsigned int high = childMask >> 32;
unsigned int lowmask = __vcmpgtu4(low, 0) & 0x08040201;
unsigned int highmask = __vcmpgtu4(high, 0) & 0x80402010;
unsigned int mask = lowmask | highmask;
mask |= mask >> 16;
mask |= mask >> 8;
parentMask = mask & 0xff;

Answer 2

在 CUDA 支持的至少某些 GPU 架构上，这种基于经典位旋转技术的解决方案可能比 accepted answer 更快，因为 __vcmp* 内在函数并非在所有架构上都快。

由于 GPU 基本上是 32 位架构，因此 64 位 childMask 被处理为两半，hi 和 lo。

处理包括三个步骤。在第一步中，我们将每个非空字节设置为0x80，否则该字节保持不变。换句话说，如果字节非零，我们设置每个字节的最高有效位。一种方法是使用 Alan Mycroft 在 1980 年代设计的空字节检测算法的修改版本，该算法通常用于C-string 处理。或者，我们可以使用hadd (~0, x) 仅在x != 0 设置最高有效位这一事实，其中hadd 是减半添加：hadd (a, b) = (a + b) / 2，没有在中间计算中溢出。 Peter L. Montgomery 在 2000 年发表了一个有效的实现。

在第二步中，我们将每个字节的最高有效位收集到最高半字节。为此，我们需要将第 7 位移动到第 28 位，将第 15 位移动到第 29 位，将第 23 位移动到第 30 位，将第 31 位移动到第 31 位，对应的移位因子为 21、14、7 和 0。为了避免单独的班次，我们将班次因子组合成一个“魔术”乘数，然后乘以它，从而并行执行所有班次。

在第三步中，我们组合包含结果的半字节并将它们移动到正确的位位置。对于 hi 字，这意味着将位 中的半字节移动到位 中，对于 lo 字，这意味着将位 中的半字节移动到位 。这种组合可以通过按位或或加法来执行。哪种变体更快可能取决于目标架构。

#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>

#define USE_HAROLDS_SOLUTION  (0)

#define USE_MYCROFT_ZEROBYTE  (0)
#define USE_TWO_MASKS         (1)
#define USE_ADD_COMBINATION   (1)

uint8_t parentMask (uint64_t childMask)
{
#if USE_TWO_MASKS
    const uint32_t LSB_MASK = 0x01010101;
#endif // USE_TWO_MASKS
    const uint32_t MSB_MASK = 0x80808080;
    const uint32_t MAGICMUL = (1 << 21) | (1 << 14) | (1 << 7) | (1 << 0);
    uint32_t lo, hi;

    /* split 64-bit argument into two halves for 32-bit GPU architecture */
    lo = (uint32_t)(childMask >>  0);
    hi = (uint32_t)(childMask >> 32);

#if USE_MYCROFT_ZEROBYTE
    /* Set most significant bit in each byte that is not zero. Adapted from Alan 
       Mycroft's null-byte detection algorithm (newsgroup comp.lang.c, 1987/04/08,
       https://groups.google.com/forum/#!original/comp.lang.c/2HtQXvg7iKc/xOJeipH6KLMJ):
       null_byte(x) = ((x - 0x01010101) & (~x & 0x80808080))
    */
#if USE_TWO_MASKS
    lo = (((lo | MSB_MASK) - LSB_MASK) | lo) & MSB_MASK;
    hi = (((hi | MSB_MASK) - LSB_MASK) | hi) & MSB_MASK;
#else // USE_TWO_MASKS
    lo = (((lo & ~MSB_MASK) + ~MSB_MASK) | lo) & MSB_MASK;
    hi = (((hi & ~MSB_MASK) + ~MSB_MASK) | hi) & MSB_MASK;
#endif // USE_TWO_MASKS

#else // USE_MYCROFT_ZEROBYTE

    /* Set most significant bit in each byte that is not zero. Use hadd(~0,x).
       Peter L. Montgomery's observation (newsgroup comp.arch, 2000/02/11,
       https://groups.google.com/d/msg/comp.arch/gXFuGZtZKag/_5yrz2zDbe4J):
       (A+B)/2 = (A AND B) + (A XOR B)/2.
    */
#if USE_TWO_MASKS
    lo = (((~lo & ~LSB_MASK) >> 1) + lo) & MSB_MASK;
    hi = (((~hi & ~LSB_MASK) >> 1) + hi) & MSB_MASK;
#else // USE_TWO_MASKS
    lo = (((~lo >> 1) & ~MSB_MASK) + lo) & MSB_MASK;
    hi = (((~hi >> 1) & ~MSB_MASK) + hi) & MSB_MASK;
#endif // USE_TWO_MASKS

#endif // USE_MYCROFT_ZEROBYTE

    /* collect most significant bit of each byte in most significant nibble */
    lo = lo * MAGICMUL;
    hi = hi * MAGICMUL;

    /* combine nibbles with results for high and low half into final result */
#if USE_ADD_COMBINATION    
    return (uint8_t)((hi >> 24) + (lo >> 28));
#else // USE_ADD_COMBINATION
    return (uint8_t)((hi >> 24) | (lo >> 28));
#endif // USE_ADD_COMBINATION
}

uint8_t parentMask_ref (uint64_t childMask)
{
    uint8_t parentMask = 0;
    for (uint32_t i = 0; i < 8; i++) {
        uint8_t child = childMask >> (8 * i);
        parentMask |= (child > 0) << i;
    }
    return parentMask;
}

uint32_t build_mask (uint32_t a)
{
    return ((a & 0x80808080) >> 7) * 0xff;
}

uint32_t vcmpgtu4 (uint32_t a, uint32_t b)
{
    uint32_t r;
    r = ((a & ~b) + (((a ^ ~b) >> 1) & 0x7f7f7f7f));
    r = build_mask (r);
    return r;
}

uint8_t parentMask_harold (uint64_t childMask)
{
    uint32_t low = childMask;
    uint32_t high = childMask >> 32;
    uint32_t lowmask = vcmpgtu4 (low, 0) & 0x08040201;
    uint32_t highmask = vcmpgtu4 (high, 0) & 0x80402010;
    uint32_t mask = lowmask | highmask;
    mask |= mask >> 16;
    mask |= mask >> 8;
    return (uint8_t)mask;
}

/*
  From: geo <gmars...@gmail.com>
  Newsgroups: sci.math,comp.lang.c,comp.lang.fortran
  Subject: 64-bit KISS RNGs
  Date: Sat, 28 Feb 2009 04:30:48 -0800 (PST)

  This 64-bit KISS RNG has three components, each nearly
  good enough to serve alone.    The components are:
  Multiply-With-Carry (MWC), period (2^121+2^63-1)
  Xorshift (XSH), period 2^64-1
  Congruential (CNG), period 2^64
*/
static uint64_t kiss64_x = 1234567890987654321ULL;
static uint64_t kiss64_c = 123456123456123456ULL;
static uint64_t kiss64_y = 362436362436362436ULL;
static uint64_t kiss64_z = 1066149217761810ULL;
static uint64_t kiss64_t;
#define MWC64  (kiss64_t = (kiss64_x << 58) + kiss64_c, \
                kiss64_c = (kiss64_x >> 6), kiss64_x += kiss64_t, \
                kiss64_c += (kiss64_x < kiss64_t), kiss64_x)
#define XSH64  (kiss64_y ^= (kiss64_y << 13), kiss64_y ^= (kiss64_y >> 17), \
                kiss64_y ^= (kiss64_y << 43))
#define CNG64  (kiss64_z = 6906969069ULL * kiss64_z + 1234567ULL)
#define KISS64 (MWC64 + XSH64 + CNG64)

int main (void)
{
    uint64_t childMask, count = 0;
    uint8_t res, ref;

    do {
        childMask = KISS64;
        ref = parentMask_ref (childMask);
#if USE_HAROLDS_SOLUTION
        res = parentMask_harold (childMask);
#else // USE_HAROLDS_SOLUTION
        res = parentMask (childMask);
#endif // USE_HAROLDS_SOLUTION
        if (res != ref) {
            printf ("\narg=%016llx  res=%02x  ref=%02x\n", childMask, res, ref);
            return EXIT_FAILURE;
        }
        if (!(count & 0xffffff)) printf ("\r%llu", count);
        count++;
    } while (1);
    return EXIT_SUCCESS;
}