【问题标题】:Fast dot product of a bit vector and a floating point vector位向量和浮点向量的快速点积
【发布时间】:2013-04-09 16:52:49
【问题描述】:

我试图在 i7 上以最有效的方式计算浮点向量和位向量之间的点积。 实际上,我在 128 或 256 维向量上执行此操作,但为了说明,让我编写 64 维的代码来说明问题:

// a has 64 elements. b is a bitvector of 64 dimensions.
float dot(float *restrict a, uint64_t b) {
    float sum = 0;
    for(int i=0; b && i<64; i++, b>>=1) {
        if (b & 1) sum += a[i];
    }
    return sum;
}

这当然可行,但问题是,这是整个程序的时间关键点(占用 50 分钟运行的 95% 的 CPU 时间),所以我迫切需要让它更快。

我的猜测是上面的分支是游戏杀手(防止乱序执行,导致错误的分支预测)。 我不确定向量指令是否可以在这里使用并有所帮助。 使用带有 -std=c99 -march=native -mtune=native -Ofast -funroll-loops 的 gcc 4.8,我目前正在得到这个输出

    movl    $4660, %edx
    movl    $5, %ecx
    xorps   %xmm0, %xmm0
    .p2align 4,,10
    .p2align 3
.L4:
    testb   $1, %cl
    je  .L2
    addss   (%rdx), %xmm0
.L2:
    leaq    4(%rdx), %rax
    shrq    %rcx
    testb   $1, %cl
    je  .L8
    addss   4(%rdx), %xmm0
.L8:
    shrq    %rcx
    testb   $1, %cl
    je  .L9
    addss   4(%rax), %xmm0
.L9:
    shrq    %rcx
    testb   $1, %cl
    je  .L10
    addss   8(%rax), %xmm0
.L10:
    shrq    %rcx
    testb   $1, %cl
    je  .L11
    addss   12(%rax), %xmm0
.L11:
    shrq    %rcx
    testb   $1, %cl
    je  .L12
    addss   16(%rax), %xmm0
.L12:
    shrq    %rcx
    testb   $1, %cl
    je  .L13
    addss   20(%rax), %xmm0
.L13:
    shrq    %rcx
    testb   $1, %cl
    je  .L14
    addss   24(%rax), %xmm0
.L14:
    leaq    28(%rax), %rdx
    shrq    %rcx
    cmpq    $4916, %rdx
    jne .L4
    ret

编辑可以对数据进行置换(只要所有参数的置换相同),顺序无关紧要。

我想知道是否有什么东西可以以 Chris Dodd 的 SSE2 代码的 3 倍以上的速度运行。

新说明:也欢迎使用 AVX/AVX2 代码!

编辑 2 给定一个位向量,我必须将它与 128 个(或 256,如果它是 256 位)不同的浮点向量相乘(因此一次也可以涉及多个浮点向量)。这是整个过程。任何可以加快整个过程的东西也是受欢迎的!

【问题讨论】:

  • 这似乎是在 GPU 上重写和执行的好代码,或者至少是多线程
  • GPU 解决方案也很受欢迎 :) 但是 dot() 需要多次执行并反馈给大 CPU 端算法(无法想象它可以写在 GPU 上运行) ,意味着 GPU 和 CPU 之间的数据传输过多,这可能会导致巨大的损失。该程序已经是多线程的,并行运行 dot() 顺便说一句。
  • 我在考虑一个流程:xmm0=0; xmm1 = *next_4_floats++; xmm0 = conditional_copy(xmm1,mask); xmm2+=xmm0; mask&lt;&lt;=1; 最后,必须对 xmm2 元素进行一次水平求和。基于掩码的复制可能比形成全 1/0 的 32 位掩码并使用并行和短一条指令。 OTOH 需要对面具进行预处理才能订购 vec_128 = |048c.....|159d....|26ae....|37bf....|
  • 对不起,我对SSE不是很熟悉,但是这里的conditional_copy和mask是什么?而且mask和b有关系吗?
  • 掩码将基于b;并且条件复制指令很可能是来自 SSE4.1 的blendvps

标签: c optimization assembly x86-64 sse


【解决方案1】:

这里有一些可以尝试的方法。

尝试让编译器使用CMOV 而不是分支。 (Note that using a union this way is well-defined in C11 but undefined in C++11.)

union {
    int i;
    float f;
} u;
u.i = 0;
if (b & 1) {
    u.f = a[i];
}
sum += u.f;

使用乘法而不是分支。

sum += (b & 1) * a[i];

保留几个和并在末尾添加它们以减少数据流依赖性。 (您可以将上述任一建议与此建议结合使用。)

float sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
for (int i = 0; i < 64; i += 4; b >>= 4) {
    if (b & 1) sum0 += a[i];
    if (b & 2) sum1 += a[i+1];
    if (b & 4) sum2 += a[i+2];
    if (b & 8) sum3 += a[i+3];
}
return sum0 + sum1 + sum2 + sum3;

通过一次处理几个位来减少分支的数量:

for (int i = 0; i < 64; i += 4, b >>= 4) {
    switch (b & 0xf) {
        case 0:
            break;
        case 1:
            sum += a[i];
            break;
        case 2:
            sum += a[i + 1];
            break;
        case 3:
            sum += a[i] + a[i+1];
            break;
        case 4:
            sum += a[i+2];
            break;
        // etc. for cases up to and including 15
    }
}

您可以为每个求和过程一次保留多个求和。在这种情况下,您可能希望使用宏或内联函数并调用它四次。

【讨论】:

  • !!(b &amp; 1)(b &amp; 1) 一样,不是吗?
  • @CarlNorum 哈哈,是的。好决定。我已经编辑了我的答案。
  • 没问题,只是打字少了点。 :-) 很棒的建议。
  • 如果您想保留多个运行总和,您仍然需要!!(b &amp; 2)(b &amp; 2) &gt;&gt; 1
  • 感谢您的意见。看起来不错,但实际上它们的表现更差(使用与以前相同的数据集)。乘法:花费 201% 的时间。在写这篇文章之前,我实际上用不同的总和进行了基准测试,大约慢了 158%。使用开关:一次移动两位得到原始的 %154 时间,一次移动四位得到 130%。
【解决方案2】:

最好的办法是使用一次对 4 个浮点数进行操作的 SSE ps 指令。您可以利用 float 0.0 全为 0 位这一事实,使用 andps 指令屏蔽不需要的元素:

#include <stdint.h>
#include <xmmintrin.h>

union {
    uint32_t i[4];
    __m128   xmm;
} mask[16] = {
 {  0,  0,  0,  0 },
 { ~0,  0,  0,  0 },
 {  0, ~0,  0,  0 },
 { ~0, ~0,  0,  0 },
 {  0,  0, ~0,  0 },
 { ~0,  0, ~0,  0 },
 {  0, ~0, ~0,  0 },
 { ~0, ~0, ~0,  0 },
 {  0,  0,  0, ~0 },
 { ~0,  0,  0, ~0 },
 {  0, ~0,  0, ~0 },
 { ~0, ~0,  0, ~0 },
 {  0,  0, ~0, ~0 },
 { ~0,  0, ~0, ~0 },
 {  0, ~0, ~0, ~0 },
 { ~0, ~0, ~0, ~0 },
};

float dot(__m128 *a, uint64_t b) {
    __m128 sum = { 0.0 };
    for (int i = 0; i < 16; i++, b>>=4)
        sum += _mm_and_ps(a[i], mask[b&0xf].xmm);
    return sum[0] + sum[1] + sum[2] + sum[3];
}

如果您预计掩码中有很多 0,那么短循环 0 可能会更快:

for (int i = 0; b; i++, b >>= 4)
    if (b & 0xf)
        sum += _mm_and_ps(a[i], mask[b&0xf].xmm);

但如果 b 是随机的,这会更慢。

【讨论】:

  • 谢谢!基准测试结果并不引人注目,但改进肯定是存在的:运行时间约为原来的 86%。但是,我仍然想知道是否有更好的答案。
  • GCC 报告它无法向量化循环。我浏览了gcc.gnu.org/projects/tree-ssa/vectorization.html,但到目前为止找不到原因。
  • @JohnSmith:循环已经矢量化(使用显式矢量类型和指令),因此无需进一步矢量化。使用 -O3,gcc 也完全展开上述循环。
  • 我明白了。我真的不熟悉 SSE :)
  • @ChrisDodd:你的代码是一个很好的例子,说明不同编译器的汇编器输出会有多么不同 :-) 我能看到的唯一变化是将它变成 map-reduce,因为这样你就赢了'不依赖于后续的addps。我怀疑效果是否巨大。
【解决方案3】:

我发现为 Chris Dodd 的代码生成的程序集非常依赖于编译器; clang 将其变成一个循环,而 gcc(4.6 和 4.7)和 Intel icc(12.x 和 13.x)改为展开循环。不过,可以通过将其转换为 map-reduce 来减少依赖关系(需要等待前一个 sum +=),

float dot(__m128 *a, uint64_t b) {
    __m128 sum[8];
    int i;

    for (i = 0; i < 8; i++) {
            sum[i] = _mm_add_ps(
                    _mm_and_ps(a[2*i], mask[b & 0xf].xmm),
                    _mm_and_ps(a[2*i+1], mask[(b & 0xf0) >> 4].xmm));
            b >>= 8;
    }
    for (i = 0; i < 4; i++) {
            sum[i] = _mm_add_ps(sum[2*i], sum[2*i+1]);
    }
    sum[0] = _mm_add_ps(sum[0], sum[1]);
    sum[2] = _mm_add_ps(sum[2], sum[3]);

    sum[0] = _mm_add_ps(sum[0], sum[2]);

    sum[0] = _mm_hadd_ps(sum[0], sum[0]);
    sum[0] = _mm_hadd_ps(sum[0], sum[0]);

    i = _mm_extract_ps(sum[0], 0);
    return *((float*)(&i));
}

这会创建明显 inferior 带有 clang 的程序集(它将 sum[] 存储在堆栈上)但 更好 代码(不依赖于后续的 addpsgccicc。有趣的是,只有gcc 最后才知道sum[0] 中较低的float 可以原地返回...

关于如何针对特定编译器进行调整的很好的练习 ...

【讨论】:

  • 对于这么短的向量来说,累加器可能有点矫枉过正,但至少有 2 个累加器会是一个更好的主意。特别是如果使用标量整数内容和 LUT 处理掩码,展开它来处理整个字节是很好的。如果编译器很聪明,他们甚至会注意到mask[b&amp;0xf0)&gt;&gt;4] 可以优化右移,因为mask 的每个元素都是一个16 字节的向量,所以b&amp;0xf0 已经是一个字节索引。
  • 顺便说一句,你的标量提取是邪恶的。永远不要使用_mm_extract_ps,使用_mm_cvtss_f32(sum[0]);。请参阅Fastest way to do horizontal float vector sum on x86 以获得比 2x haddps 更好的 hsum。
【解决方案4】:

如果允许在数据float data[128]或者进行稍微不同的排列,或者在__m128 mask; 的位掩码中进行相应的排列,则可以稍微改进上面 Chris Dodd 建议的算法. (不计算置换掩码所需的时间,这个实现(+开销)大约快 25%)。这当然只是我在 cmets 中提供的想法的草稿。

union {
    unsigned int i[4];
    float f[4];
    __m128   xmm;
} mask = { 0xFF00FF00, 0xF0F0F0F0, 0xCCCCCCCC, 0xAAAAAAAA };

float dot2(__m128 *a, __m128 mask);
// 20M times 1.161s

float dotref(__m128 *a, unsigned int *mask) // 20M times 8.174s
{
    float z=0.0f;
    int i;
    for (i=0;i<32;i++) {
       if (mask[0] & (0x80000000U >> i)) z+= a[i][0];
       if (mask[1] & (0x80000000U >> i)) z+= a[i][1];
       if (mask[2] & (0x80000000U >> i)) z+= a[i][2];
       if (mask[3] & (0x80000000U >> i)) z+= a[i][3];       
    }
   return z;
}

相应的汇编器实现将是:

dot2:
    // warm up stage: fill in initial data and
    // set up registers
    pxor %xmm1, %xmm1      ;; // clear partial sum1
    pxor %xmm2, %xmm2      ;; // clear partial sum2
    movaps (%rdi), %xmm3   ;; // register warm up stage1
    movaps 16(%rdi), %xmm4 ;; // next 4 values
    pxor %xmm5, %xmm5
    pxor %xmm6, %xmm6
    lea 32(%rdi), %rdi
    movl $16, %ecx            ;; // process 2x4 items per iteration (total=128)
a:  ;; // inner loop -- 2 independent data paths
    blendvps %xmm3, %xmm5
    pslld $1, %xmm0
    movaps (%rdi), %xmm3   
    blendvps %xmm4, %xmm6
    pslld $1, %xmm0
    movaps 16(%rdi), %xmm4
    addps %xmm5, %xmm1
    pxor  %xmm5, %xmm5
    addps %xmm6, %xmm2
    pxor  %xmm6, %xmm6
    lea 32(%rdi), %rdi
    loop a
 ;; // cool down stage: gather results (xmm0 = xmm1+xmm2)
 ;; // in beautiful world this stage is interleaved
 ;; // with the warm up stage of the next block
    addps %xmm2, %xmm1
    movaps  %xmm1, %xmm2
    movaps  %xmm1, %xmm0
    shufps  $85, %xmm1, %xmm2
    addss   %xmm2, %xmm0
    movaps  %xmm1, %xmm2
    unpckhps %xmm1, %xmm2
    shufps  $255, %xmm1, %xmm1
    addss   %xmm2, %xmm0
    addss   %xmm1, %xmm0
    ret

【讨论】:

  • 还应该复制掩码xmm1 = xmm0 &lt;&lt; 1;并将两者都移动 2。这将删除内部循环中的最后一个依赖项。也可以尝试找到更好的交错偏移。这使用 16 的偏移量,但 64 或 256 的偏移量可能会更好。此外,任一选项都可能导致更容易计算位掩码排列(在 xmm0 和 xmm1 中)。
【解决方案5】:

如果你有 i7,那么你有 SSE4.1,你可以使用 _mm_dp_ps 内在函数,它做一个点积。使用您的示例代码,它看起来像

#include <stdint.h>
#include <immintrin.h>

const float fltmask[][4] =
  {{0, 0, 0, 0}, {0, 0, 0, 1}, {0, 0, 1, 0}, {0, 0, 1, 1},
   {0, 1, 0, 0}, {0, 1, 0, 1}, {0, 1, 1, 0}, {0, 1, 1, 1},
   {1, 0, 0, 0}, {1, 0, 0, 1}, {1, 0, 1, 0}, {1, 0, 1, 1},
   {1, 1, 0, 0}, {1, 1, 0, 1}, {1, 1, 1, 0}, {1, 1, 1, 1}};

// a has 64 elements. b is a bitvector of 64 dimensions.
float dot(float * restrict a, uint64_t b) {
     int i;
    float sum = 0;
    for(i=0; b && i<64; i+=4,b>>=4) {
        __m128 t0 = _mm_load_ps (a);
        a += 4;
        __m128 t1 = _mm_load_ps (fltmask[b & 15]);
        sum += _mm_cvtss_f32 (_mm_dp_ps (t0, t1, 15));
    }
    return sum;
}

PS。数组 afltmask 最好是 16 字节对齐!

PPS。当使用gcc -std=c99 -msse4 -O2 编译时,循环看起来像:

.L3:
    movq    %rdx, %rax
    movaps  (%rcx), %xmm1
    shrq    $4, %rdx
    andl    $15, %eax
    addq    $16, %rcx
    addl    $4, %r8d
    salq    $4, %rax
    testq   %rdx, %rdx
    dpps    $15, (%r9,%rax), %xmm1
    addss   %xmm1, %xmm0
    jne .L13

当然,-O3 是展开的。

【讨论】:

  • 感谢您的意见!抱歉,我这里只有 i3 (SSE3),所以我现在不能做基准测试。不过这个周末我可以访问 i7。我想知道基准数字是多少(对于 ~1M uint64_t 数字)?
  • i5 也有 SSE4.1。我认为自 Core 2 (Penryn) 以来的所有台式机 Intel CPU 都具有 SSE4.1。不知道基准数字,它还取决于许多其他因素(内存、缓存)。
  • 对不起,我检查了 /proc/cpu,结果是 i3 + SSE3。不管怎样,我明天试试。我有一个广泛的数据集,因此这些因素将被消除。
  • 已知点积和其他水平指令不如(速度方面)更长和更复杂的序列。
  • 确实,这比我使用 SSE4.2 的 i5 机器上的 SSE2 代码慢了约 2 倍 :(
【解决方案6】:

为了扩展 Aki Suihkonen 的答案,重塑位串有助于有条件地移动浮点数。在下面的解决方案中,使用 SSE 指令 PMOVMASKB 和 PSHUFB 的两级位置换,加上指令 BLENDVPS 已用于在 Core 2 Duo 2.26GHz 上实现 1.25 个元素处理/周期,这是我参考速度的 20 倍C 代码。

[编辑:添加了 AVX2 实现。性能未知,因为我自己无法测试,但预计速度会提高一倍。 ]

这是我的实现和测试平台,解释如下。

SSE4.1(旧,AVX2 见下文)

代码

/* Includes */
#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <smmintrin.h>  /* SSE 4.1 */
#include <time.h>



/* Defines */
#define ALIGNTO(n) __attribute__((aligned(n)))
#define USE_PINSRW 1
#define NUM_ITERS  2260000



/**
 * Bit mask shuffle.
 * 
 * This version uses a loop to store eight u16 and reloads them as one __m128i.
 */

__m128 bitMaskShuffleStoreAndReload(__m128i mask){
    const __m128i perm    ALIGNTO(16) = _mm_set_epi8(15, 7, 14, 6, 13, 5, 12, 4,
                                                     11, 3, 10, 2, 9,  1, 8,  0);
    int i;
    uint16_t interMask[8] ALIGNTO(16);

    /* Shuffle bitmask */
        /* Stage 1 */
    for(i=7;i>=0;i--){
        interMask[i] = _mm_movemask_epi8(mask);
        mask = _mm_slli_epi32(mask, 1);
    }

        /* Stage 2 */
    return _mm_castsi128_ps(
                    _mm_shuffle_epi8(
                        _mm_load_si128((const __m128i*)interMask),
                        perm)
                );
}

/**
 * Bit mask shuffle.
 * 
 * This version uses the PINSTRW instruction.
 */

__m128 bitMaskShufflePINSRW(__m128i mask){
    const __m128i perm    ALIGNTO(16) = _mm_set_epi8(15, 7, 14, 6, 13, 5, 12, 4,
                                                     11, 3, 10, 2, 9,  1, 8,  0);
    __m128i  imask ALIGNTO(16);

    /* Shuffle bitmask */
        /* Stage 1 */
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 7);
    mask = _mm_slli_epi16(mask, 1);
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 6);
    mask = _mm_slli_epi16(mask, 1);
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 5);
    mask = _mm_slli_epi16(mask, 1);
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 4);
    mask = _mm_slli_epi16(mask, 1);
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 3);
    mask = _mm_slli_epi16(mask, 1);
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 2);
    mask = _mm_slli_epi16(mask, 1);
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 1);
    mask = _mm_slli_epi16(mask, 1);
    imask = _mm_insert_epi16(imask, _mm_movemask_epi8(mask), 0);

        /* Stage 2 */
    return _mm_castsi128_ps(
                    _mm_shuffle_epi8(
                        imask,
                        perm)
                );
}

/**
 * SSE 4.1 implementation.
 */

float dotSSE41(__m128 f[32], unsigned char maskArg[16]){
    int i, j, k;
    __m128i  mask      ALIGNTO(16) = _mm_load_si128((const __m128i*)maskArg);
    __m128   shufdMask ALIGNTO(16);
    __m128   zblended  ALIGNTO(16);
    __m128   sums      ALIGNTO(16) = _mm_setzero_ps();
    float    sumsf[4]  ALIGNTO(16);

    /* Shuffle bitmask */
#if USE_PINSRW
    shufdMask = bitMaskShufflePINSRW(mask);
#else
    shufdMask = bitMaskShuffleStoreAndReload(mask);
#endif

    /* Dot product */
    for(i=1;i>=0;i--){
        for(j=1;j>=0;j--){
            for(k=7;k>=0;k--){
                zblended  = _mm_setzero_ps();
                zblended  = _mm_blendv_ps(zblended, f[i*16+j+k*2], shufdMask);
                sums      = _mm_add_ps(sums, zblended);
                shufdMask = _mm_castsi128_ps(_mm_slli_epi32(_mm_castps_si128(shufdMask), 1));
            }
        }
    }

    /* Final Summation */
    _mm_store_ps(sumsf, sums);
    return sumsf[0] + sumsf[1] + sumsf[2] + sumsf[3];
}

/**
 * Reference C implementation
 */


float dotRefC(float f[128], unsigned char mask[16]){
    float sum = 0.0;
    int i;

    for(i=0;i<128;i++){
        sum += ((mask[i>>3]>>(i&7))&1) ? f[i] : 0.0;
    }

    return sum;
}

/**
 * Main
 */

int main(void){
    /* Variables */
        /* Loop Counter */
    int           i;

        /* Data to process */
    float         data[128] ALIGNTO(16);
    unsigned char mask[16]  ALIGNTO(16);
    float         refCSum, sseSum;

        /* Time tracking */
    clock_t       t1, t2, t3;
    double        refCTime, sseTime;



    /* Initialize mask and float arrays with some random data. */
    for(i=0;i<128;i++){
        if(i<16)
            mask[i]=rand();

        data[i] = random();
    }


    /* RUN TESTS */
    t1 = clock();
    for(i=0;i<NUM_ITERS;i++){
        refCSum = dotRefC(data, mask);
    }
    t2 = clock();
    for(i=0;i<NUM_ITERS;i++){
        sseSum  = dotSSE41((__m128*)data, mask);
    }
    t3 = clock();


    /* Compute time taken */
    refCTime = (double)(t2-t1)/CLOCKS_PER_SEC;
    sseTime  = (double)(t3-t2)/CLOCKS_PER_SEC;


    /* Print out results */
    printf("Results:\n"
           "RefC: Time: %f    Value: %f\n"
           "SSE:  Time: %f    Value: %f\n",
           refCTime, refCSum,
           sseTime,  sseSum);

    return 0;
}

说明

BLENDVPS 使用 128 位寄存器 XMM0 的所有四个 32 位通道中的最高位来确定是否将其源操作数的相应通道中的值移动到其目标操作数中。使用 MOVAPS 加载数据时,会获得 4 个连续的浮点数:例如,第 8、第 9、第 10 和第 11 个浮点数。当然,它们的选择或取消选择必须由相应的一组位来控制:例如,位串中的第 8、9、10 和 11 位。

问题是,当第一次加载掩码时,这些集合的位彼此相邻(在第 8、第 9、第 10 和第 11 位),而实际上它们应该相隔 32 位;请记住,在某些时候,它们必须占据每个通道的第 31 位位置(XMM0 寄存器中的第 31、63、95 和 127 位)。

理想情况下会发生位转置,将位 0、4、8、12,... 带入通道 0,位 1, 5, 9, 13, ... 在通道一中,位 2, 6, 10, 14, ... 在通道二和位 3, 7, 11, 15, ... 在通道三。因此,以前连续的所有 4 位集合现在都间隔 32 位,四个 32 位通道中的每一个通道中都有一个。那么它所需要的只是一个循环,迭代 32 次,每次将一组新的 4 位转移到每个通道的最高位位置。

不幸的是,x86 不具备良好的位操作指令,因此由于缺乏一种干净的方式来进行完美的转置,因此这里是一种合理的折衷方案。

在掩码中,128 位

  0   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15
 16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31
 32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47
 48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63
 64  65  66  67  68  69  70  71  72  73  74  75  76  77  78  79
 80  81  82  83  84  85  86  87  88  89  90  91  92  93  94  95
 96  97  98  99 100 101 102 103 104 105 106 107 108 109 110 111
112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127

由 8 个 PMOVMASKB 和 8 个 PSLLW 指令排列,首先到

  0   8  16  24  32  40  48  56  64  72  80  88  96 104 112 120
  1   9  17  25  33  41  49  57  65  73  81  89  97 105 113 121
  2  10  18  26  34  42  50  58  66  74  82  90  98 106 114 122
  3  11  19  27  35  43  51  59  67  75  83  91  99 107 115 123
  4  12  20  28  36  44  52  60  68  76  84  92 100 108 116 124
  5  13  21  29  37  45  53  61  69  77  85  93 101 109 117 125
  6  14  22  30  38  46  54  62  70  78  86  94 102 110 118 126
  7  15  23  31  39  47  55  63  71  79  87  95 103 111 119 127

并通过一条 PSHUFB 指令来

  0   8  16  24  32  40  48  56   4  12  20  28  36  44  52  60
 64  72  80  88  96 104 112 120  68  76  84  92 100 108 116 124
  1   9  17  25  33  41  49  57   5  13  21  29  37  45  53  61
 65  73  81  89  97 105 113 121  69  77  85  93 101 109 117 125
  2  10  18  26  34  42  50  58   6  14  22  30  38  46  54  62
 66  74  82  90  98 106 114 122  70  78  86  94 102 110 118 126
  3  11  19  27  35  43  51  59   7  15  23  31  39  47  55  63
 67  75  83  91  99 107 115 123  71  79  87  95 103 111 119 127

。我们现在迭代四个“运行”,每个“运行”包含八组四位,以 32 位的间隔分布(如我们所愿),使用这些组作为 BLENDVPS 的掩码控制。 bit shuffle 固有的笨拙导致dotSSE41() 中看起来笨拙的三重嵌套循环,但是

clang -Ofast -ftree-vectorize -finline-functions -funroll-loops -msse4.1 -mssse3 dot.c -o dottest

无论如何,循环都已展开。内循环迭代由 16 次重复组成

blendvps    0x90(%rsi),%xmm1
addps       %xmm4,%xmm1
pslld       $0x1,%xmm2
movdqa      %xmm2,%xmm0
xorps       %xmm4,%xmm4

.

顺便说一句,我无法确定我的两个 shuffle 版本中哪个最快,所以我在回答中给出了这两种实现。

AVX2(新,但未经测试)

代码

/* Includes */
#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <immintrin.h>  /* AVX2 */
#include <time.h>



/* Defines */
#define ALIGNTO(n) __attribute__((aligned(n)))
#define NUM_ITERS  2260000


/**
 * Bit mask shuffle.
 * 
 * This version uses the PINSTRW instruction.
 */

__m256 bitMaskShufflePINSRW(__m256i mask){
    __m256i  imask ALIGNTO(32);

    /* Shuffle bitmask */
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 7);
    mask = _mm256_slli_epi32(mask, 1);
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 6);
    mask = _mm256_slli_epi32(mask, 1);
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 5);
    mask = _mm256_slli_epi32(mask, 1);
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 4);
    mask = _mm256_slli_epi32(mask, 1);
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 3);
    mask = _mm256_slli_epi32(mask, 1);
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 2);
    mask = _mm256_slli_epi32(mask, 1);
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 1);
    mask = _mm256_slli_epi32(mask, 1);
    imask = _mm256_insert_epi32(imask, _mm256_movemask_epi8(mask), 0);

    /* Return bitmask */
    return _mm256_castsi256_ps(imask);
}

/**
 * AVX2 implementation.
 */

float dotAVX2(__m256 f[16], unsigned char maskArg[16]){
    int i, j, k;
    /* Use _mm_loadu_si128 */
    __m256i  mask      ALIGNTO(32) = _mm256_castsi128_si256(_mm_load_si128((const __m128i*)maskArg));
    __m256   shufdMask ALIGNTO(32);
    __m256   zblended  ALIGNTO(32);
    __m256   sums      ALIGNTO(32) = _mm256_setzero_ps();
    float    sumsf[8]  ALIGNTO(32);

    /* Shuffle bitmask */
    shufdMask = bitMaskShufflePINSRW(mask);
    shufdMask = _mm256_castsi256_ps(_mm256_slli_epi32(_mm256_castps_si256(shufdMask), 16));

    /* Dot product */
    for(i=15;i>=0;i--){
        zblended  = _mm256_setzero_ps();
        /* Replace f[i] with _mm256_loadu_ps((float*)&f[i]) */
        zblended  = _mm256_blendv_ps(zblended, f[i], shufdMask);
        sums      = _mm256_add_ps(sums, zblended);
        shufdMask = _mm256_castsi256_ps(_mm256_slli_epi32(_mm256_castps_si256(shufdMask), 1));
    }

    /* Final Summation */
    _mm256_store_ps(sumsf, sums);
    return sumsf[0] + sumsf[1] + sumsf[2] + sumsf[3] + sumsf[4] + sumsf[5] + sumsf[6] + sumsf[7];
}

/**
 * Reference C implementation
 */

float dotRefC(float f[128], unsigned char mask[16]){
    float sum = 0.0;
    int i;

    for(i=0;i<128;i++){
        sum += ((mask[i>>3]>>(i&7))&1) ? f[i] : 0.0;
    }

    return sum;
}

/**
 * Main
 */

int main(void){
    /* Variables */
        /* Loop Counter */
    int           i;

        /* Data to process */
    float         data[128] ALIGNTO(32);
    unsigned char mask[16]  ALIGNTO(32);
    float         refCSum, sseSum;

        /* Time tracking */
    clock_t       t1, t2, t3;
    double        refCTime, sseTime;



    /* Initialize mask and float arrays with some random data. */
    for(i=0;i<128;i++){
        if(i<16)
            mask[i]=rand();

        data[i] = random();
    }


    /* RUN TESTS */
    t1 = clock();
    for(i=0;i<NUM_ITERS;i++){
        refCSum = dotRefC(data, mask);
    }
    t2 = clock();
    for(i=0;i<NUM_ITERS;i++){
        sseSum  = dotAVX2((__m256*)data, mask);
    }
    t3 = clock();


    /* Compute time taken */
    refCTime = (double)(t2-t1)/CLOCKS_PER_SEC;
    sseTime  = (double)(t3-t2)/CLOCKS_PER_SEC;


    /* Print out results */
    printf("Results:\n"
           "RefC: Time: %f    Value: %f\n"
           "SSE:  Time: %f    Value: %f\n",
           refCTime, refCSum,
           sseTime,  sseSum);

    return 0;
}

说明

使用与 SSE4.1 相同的概念。不同之处在于,现在我们一次处理 8 个浮点数,并使用 AVX2 的 256 位寄存器和来自 ymm 寄存器的 PMOVMASKB(收集 256/8 = 32 位)。正因为如此,我们现在有了一个更简单的位掩码 shuffle 和一个更简单的循环。

在掩码中,256 位

  0   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15
 16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31
 32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47
 48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63
 64  65  66  67  68  69  70  71  72  73  74  75  76  77  78  79
 80  81  82  83  84  85  86  87  88  89  90  91  92  93  94  95
 96  97  98  99 100 101 102 103 104 105 106 107 108 109 110 111
112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127
128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143
144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159
160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175
176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191
192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207
208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223
224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239
240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255

使用 8 条 PMOVMASKB 和 8 条 PSLLW 指令置换为

  0   8  16  24  32  40  48  56  64  72  80  88  96 104 112 120
128 136 144 152 160 168 176 184 192 200 208 216 224 232 240 248
  1   9  17  25  33  41  49  57  65  73  81  89  97 105 113 121
129 137 145 153 161 169 177 185 193 201 209 217 225 233 241 249
  2  10  18  26  34  42  50  58  66  74  82  90  98 106 114 122
130 138 146 154 162 170 178 186 194 202 210 218 226 234 242 250
  3  11  19  27  35  43  51  59  67  75  83  91  99 107 115 123
131 139 147 155 163 171 179 187 195 203 211 219 227 235 243 251
  4  12  20  28  36  44  52  60  68  76  84  92 100 108 116 124
132 140 148 156 164 172 180 188 196 204 212 220 228 236 244 252
  5  13  21  29  37  45  53  61  69  77  85  93 101 109 117 125
133 141 149 157 165 173 181 189 197 205 213 221 229 237 245 253
  6  14  22  30  38  46  54  62  70  78  86  94 102 110 118 126
134 142 150 158 166 174 182 190 198 206 214 222 230 238 246 254
  7  15  23  31  39  47  55  63  71  79  87  95 103 111 119 127
135 143 151 159 167 175 183 191 199 207 215 223 231 239 247 255

。对于 128 个元素的位与浮点点积,我们然后在八组 16 个元素上并行迭代。通过迭代 32 个元素的集合,可以轻松地将这个实现扩展到 256 个元素的 DP。现在只需要一个循环。

具体来说,要将其更改为适用于 256 元素点积,您需要

  • 函数参数的大小加倍。 __m256 f[32], unsigned char maskArg[32]
  • 将掩码负载 (= _mm256_castsi128_si256(_mm_load_si128((const __m128i*)maskArg));) 替换为 = _mm256_load_si256((const __m256i*)maskArg);
  • 在对bitMaskShufflePINSRW 的调用正下方删除补偿左移16。
  • 从第 31 组而不是第 15 组向下运行循环:for(i=31;i&gt;=0;i--)

由于我的 CPU 是 SSE4.1,我既不能测试也不能运行代码,但是 Clang with

clang -Ofast -ftree-vectorize -finline-functions -funroll-loops -mavx2 -msse4.1 -mssse3 dotavx2.c -o dottest

编译干净(无需展开即可获得更快的代码),生成以下代码:

(gdb) disas dotAVX2
Dump of assembler code for function dotAVX2:
0x0000000100001730 <+0>:    push   %rbp
0x0000000100001731 <+1>:    mov    %rsp,%rbp
0x0000000100001734 <+4>:    vmovdqa (%rsi),%xmm0
0x0000000100001738 <+8>:    vpslld $0x1,%ymm0,%ymm1
0x000000010000173d <+13>:   vpslld $0x1,%ymm1,%ymm2
0x0000000100001742 <+18>:   vpmovmskb %ymm2,%eax
0x0000000100001746 <+22>:   vpslld $0x1,%ymm2,%ymm2
0x000000010000174b <+27>:   vpmovmskb %ymm2,%ecx
0x000000010000174f <+31>:   vxorps %ymm3,%ymm3,%ymm3
0x0000000100001753 <+35>:   vmovd  %ecx,%xmm4
0x0000000100001757 <+39>:   vpinsrd $0x1,%eax,%xmm4,%xmm4
0x000000010000175d <+45>:   vpmovmskb %ymm1,%eax
0x0000000100001761 <+49>:   vpinsrd $0x2,%eax,%xmm4,%xmm1
0x0000000100001767 <+55>:   vpslld $0x1,%ymm2,%ymm2
0x000000010000176c <+60>:   vpslld $0x1,%ymm2,%ymm4
0x0000000100001771 <+65>:   vpslld $0x1,%ymm4,%ymm5
0x0000000100001776 <+70>:   vpmovmskb %ymm0,%eax
0x000000010000177a <+74>:   vpinsrd $0x3,%eax,%xmm1,%xmm0
0x0000000100001780 <+80>:   vpmovmskb %ymm5,%eax
0x0000000100001784 <+84>:   vpslld $0x1,%ymm5,%ymm1
0x0000000100001789 <+89>:   vpmovmskb %ymm1,%ecx
0x000000010000178d <+93>:   vmovd  %ecx,%xmm1
0x0000000100001791 <+97>:   vpinsrd $0x1,%eax,%xmm1,%xmm1
0x0000000100001797 <+103>:  vpmovmskb %ymm4,%eax
0x000000010000179b <+107>:  vpinsrd $0x2,%eax,%xmm1,%xmm1
0x00000001000017a1 <+113>:  vpmovmskb %ymm2,%eax
0x00000001000017a5 <+117>:  vpinsrd $0x3,%eax,%xmm1,%xmm1
0x00000001000017ab <+123>:  vinserti128 $0x1,%xmm0,%ymm1,%ymm0
0x00000001000017b1 <+129>:  vpslld $0x10,%ymm0,%ymm0
0x00000001000017b6 <+134>:  vblendvps %ymm0,0x1e0(%rdi),%ymm3,%ymm1
0x00000001000017c0 <+144>:  vpslld $0x1,%ymm0,%ymm0
0x00000001000017c5 <+149>:  vblendvps %ymm0,0x1c0(%rdi),%ymm3,%ymm2
0x00000001000017cf <+159>:  vaddps %ymm2,%ymm1,%ymm1
0x00000001000017d3 <+163>:  vpslld $0x1,%ymm0,%ymm0
0x00000001000017d8 <+168>:  vblendvps %ymm0,0x1a0(%rdi),%ymm3,%ymm2
0x00000001000017e2 <+178>:  vaddps %ymm2,%ymm1,%ymm1
0x00000001000017e6 <+182>:  vpslld $0x1,%ymm0,%ymm0
0x00000001000017eb <+187>:  vblendvps %ymm0,0x180(%rdi),%ymm3,%ymm2
0x00000001000017f5 <+197>:  vaddps %ymm2,%ymm1,%ymm1
0x00000001000017f9 <+201>:  vpslld $0x1,%ymm0,%ymm0
0x00000001000017fe <+206>:  vblendvps %ymm0,0x160(%rdi),%ymm3,%ymm2
0x0000000100001808 <+216>:  vaddps %ymm2,%ymm1,%ymm1
0x000000010000180c <+220>:  vpslld $0x1,%ymm0,%ymm0
0x0000000100001811 <+225>:  vblendvps %ymm0,0x140(%rdi),%ymm3,%ymm2
0x000000010000181b <+235>:  vaddps %ymm2,%ymm1,%ymm1
0x000000010000181f <+239>:  vpslld $0x1,%ymm0,%ymm0
0x0000000100001824 <+244>:  vblendvps %ymm0,0x120(%rdi),%ymm3,%ymm2
0x000000010000182e <+254>:  vaddps %ymm2,%ymm1,%ymm1
0x0000000100001832 <+258>:  vpslld $0x1,%ymm0,%ymm0
0x0000000100001837 <+263>:  vblendvps %ymm0,0x100(%rdi),%ymm3,%ymm2
0x0000000100001841 <+273>:  vaddps %ymm2,%ymm1,%ymm1
0x0000000100001845 <+277>:  vpslld $0x1,%ymm0,%ymm0
0x000000010000184a <+282>:  vblendvps %ymm0,0xe0(%rdi),%ymm3,%ymm2
0x0000000100001854 <+292>:  vaddps %ymm2,%ymm1,%ymm1
0x0000000100001858 <+296>:  vpslld $0x1,%ymm0,%ymm0
0x000000010000185d <+301>:  vblendvps %ymm0,0xc0(%rdi),%ymm3,%ymm2
0x0000000100001867 <+311>:  vaddps %ymm2,%ymm1,%ymm1
0x000000010000186b <+315>:  vpslld $0x1,%ymm0,%ymm0
0x0000000100001870 <+320>:  vblendvps %ymm0,0xa0(%rdi),%ymm3,%ymm2
0x000000010000187a <+330>:  vaddps %ymm2,%ymm1,%ymm1
0x000000010000187e <+334>:  vpslld $0x1,%ymm0,%ymm0
0x0000000100001883 <+339>:  vblendvps %ymm0,0x80(%rdi),%ymm3,%ymm2
0x000000010000188d <+349>:  vaddps %ymm2,%ymm1,%ymm1
0x0000000100001891 <+353>:  vpslld $0x1,%ymm0,%ymm0
0x0000000100001896 <+358>:  vblendvps %ymm0,0x60(%rdi),%ymm3,%ymm2
0x000000010000189d <+365>:  vaddps %ymm2,%ymm1,%ymm1
0x00000001000018a1 <+369>:  vpslld $0x1,%ymm0,%ymm0
0x00000001000018a6 <+374>:  vblendvps %ymm0,0x40(%rdi),%ymm3,%ymm2
0x00000001000018ad <+381>:  vaddps %ymm2,%ymm1,%ymm1
0x00000001000018b1 <+385>:  vpslld $0x1,%ymm0,%ymm0
0x00000001000018b6 <+390>:  vblendvps %ymm0,0x20(%rdi),%ymm3,%ymm2
0x00000001000018bd <+397>:  vaddps %ymm2,%ymm1,%ymm1
0x00000001000018c1 <+401>:  vpslld $0x1,%ymm0,%ymm0
0x00000001000018c6 <+406>:  vblendvps %ymm0,(%rdi),%ymm3,%ymm0
0x00000001000018cc <+412>:  vaddps %ymm0,%ymm1,%ymm0
0x00000001000018d0 <+416>:  vpshufd $0x1,%xmm0,%xmm1
0x00000001000018d5 <+421>:  vaddss %xmm1,%xmm0,%xmm1
0x00000001000018d9 <+425>:  vmovhlps %xmm0,%xmm0,%xmm2
0x00000001000018dd <+429>:  vaddss %xmm1,%xmm2,%xmm1
0x00000001000018e1 <+433>:  vpshufd $0x3,%xmm0,%xmm2
0x00000001000018e6 <+438>:  vaddss %xmm1,%xmm2,%xmm1
0x00000001000018ea <+442>:  vextracti128 $0x1,%ymm0,%xmm0
0x00000001000018f0 <+448>:  vaddss %xmm1,%xmm0,%xmm1
0x00000001000018f4 <+452>:  vpshufd $0x1,%xmm0,%xmm2
0x00000001000018f9 <+457>:  vaddss %xmm1,%xmm2,%xmm1
0x00000001000018fd <+461>:  vpshufd $0x3,%xmm0,%xmm2
0x0000000100001902 <+466>:  vmovhlps %xmm0,%xmm0,%xmm0
0x0000000100001906 <+470>:  vaddss %xmm1,%xmm0,%xmm0
0x000000010000190a <+474>:  vaddss %xmm0,%xmm2,%xmm0
0x000000010000190e <+478>:  pop    %rbp
0x000000010000190f <+479>:  vzeroupper 
0x0000000100001912 <+482>:  retq   
End of assembler dump.

正如我们所见,内核现在是 3 条指令(vblendvps、vaddps、vpslld)。

【讨论】:

  • @John Smith,我已经编辑了我的解决方案以包含(未经测试,因为我没有 Haswell 硬件)AVX2 代码。请尝试一下。
  • @Aki Suihkonen:事实证明,在 AVX2 中可以杀死内核中的 两条 指令!我的内核现在只有 3 条指令。
  • 我也不知道在 Haswell 上展开循环的智慧;循环内核非常小,以至于不展开或限制展开因子可能是有利可图的,这样它就完全适合 uop 缓存和循环缓存。克隆循环体并在两个副本中向左移动两个也可能有所帮助,但我不知道有什么方法可以说服编译器这样做。
  • 嗯...我想你是想测试-fno-unroll-loops? AVX2 在某种程度上是不可避免的,因为 VPMOVMSKB 和 VPSLLD 是 256 位整数指令,而英特尔出于某种原因决定在 AVX 中只有浮点运算才能获得 256 位扩展。整数运算只有通过 AVX2 才能获得该特权。
  • Chris Dodd 的代码,因为它只使用 AND 位运算(在浮点版本中可用),实际上可以在 AVX 中实现。但是,您面临的选择是将他的 LUT 扩展到 256 个元素,每个元素有 8 个浮点数,每个 4 字节 (8KB),或者将 LUT 中的掩码加载到两个 4 元素半部分中,将上半部分洗牌到YMM 寄存器的前 128 位,将其与您的数据相加并累加。
【解决方案7】:

你可以像这样消除分支:

for(int i=0; b && i<64; i++, b>>=1)
    sum += a[i] * (b & 1);

虽然这会增加额外的 mult,但至少不会破坏您的管道。

控制分支的另一种方法是按您的方式使用它,但也使用编译器宏。我猜在 gcc 中它是 likely(if ...) 宏。您将使用分支,但这样您会告诉编译器该分支将被更频繁地执行,而 gcc 将进行更多优化。

另一个可以做的优化是“缓存”点积。因此,您将拥有一个变量,而不是拥有一个计算点积的函数,该变量在开始时将乘积初始化为 0。并且每次插入/删除/更新向量的元素时,也会更新保存结果的变量。

【讨论】:

  • 如果你想鼓励无分支,写sum += b&amp;1 ? a[i] : 0;。实际的 FP 乘以 0.0 或 1.0 通常是个坏主意,而严格的 FP 就是你会得到的。尽管理想情况下编译器会将其压缩为 FMA(如果可用),但仍然不好。
【解决方案8】:

描述这个的另一种方式是掩蔽缩减(sum)。

AVX512 内置了掩码,例如与_mm256_mask_add_ps。 (对于如此短的数组,您可能只想在 Skylake-AVX512 上使用 256 位向量,以避免限制最大 turbo。除非您的代码花费大量时间执行此操作或其他可以从 512 位向量中受益的循环。

它实际上是从内存源到vaddps 的一条asm 指令,在未设置掩码寄存器中的位的位置保持元素不变。例如你的循环体实际上可以像

那样简单
vaddps    zmm0{k1}, zmm0, [rdi]    ; merge-masking into ZMM0
kshiftrq  k1, k1, 16               ; bring the next 16 mask bits down to the bottom

(从初始化 zmm0 的零掩码加载开始。实际上,您需要多个累加器来隐藏 FP-add 延迟,因此使用掩码中的低对位向量进行多个零掩码加载是理想的。)

所以你真的希望编译器能够输出类似这样的 asm 内容。 (用内在函数写起来应该很简单,但更冗长,所以我只写了 asm 而不是花时间查找内在名称。)

   ;float dot(float *restrict a, uint64_t b)
dot:
    kmovq     k1, rsi                  ; __mmask64 m = b
    vmovups   zmm0{k1}{z}, [rdi]       ; zero-masking
    kshiftrq  k2, k1, 16
    vmovups   zmm1{k2}{z}, [rdi+ 16*4]   ; next vector of 16x 4-byte floats
    kshiftrq  k3, k1, 32
    kshiftrq  k4, k1, 48

    vaddps    zmm0{k3}, zmm0, [rdi + 16*4 * 2]   ; merge-masking into ZMM0
    vaddps    zmm1{k4}, zmm1, [rdi + 16*4 * 3]

    ;; if you have more mask data, use it and do some more adds, or maybe run 2 more dep chains of FP add before combining.

    vaddps    zmm0, zmm0, zmm1               ; reduce down to 1 accumulator

    ;then horizontal sum down to scalar
    VEXTRACTF64x4   ymm1, zmm0, 1
    vaddps          ymm0, ymm0, ymm1    ; narrow to 256
    vextractf128    xmm1, ymm0, 1
    vaddps          xmm0, xmm0, xmm1    ; narrow to 128
    vmovshdup / vaddps / vmovhlps / vaddss

    ret

请参阅Fastest way to do horizontal float vector sum on x86,了解如何在 C 中使用内部函数有效地编写水平和,因此它们可以编译为不糟糕的 asm。 (不是haddps,也不是一个简单的存储阵列和循环。)

如果掩码数据来自内存,您可以将其直接加载到掩码寄存器中。不幸的是,在当前的 Skylake-AVX512 上,这需要 3 微秒,因此您想要进行 64 位加载并使用掩码移位。或者让编译器做它的事情。


使用 AVX2,@IWill 的回答看起来是在正确的轨道上,但我怀疑最好将掩码处理与掩码使用交错。

例如从浮点数组的末尾开始,并使用掩码的前 8 位。首先将掩码的前 32 位广播到 YMM 的所有元素。然后vpsllvd左移[0,1,2,3,4,5,6,7],将掩码的最高位留在高元素的顶部,将掩码的第8高位留在低元素中。

然后使用vblendvps 或更好的vmaskmovps 从向量的高端进行屏蔽加载(Skylake 上的微指令更少,与 Haswell 上的 blendv 相同)。 vblendvps 始终至少为 2 微秒。有趣的事实:具有隐式 XMM0 源的非 VEX 版本是 Skylake 上的单微指令。

然后将vpslld 加 8 以将接下来的 8 个掩码位放置到位,并将指针减 32 个字节。

您必须在每 4 个向量中停止并广播 + 可变移位新掩码位,但使用 FP 添加的交错工作可能非常适合隐藏所涉及的所有内容的延迟。

请参阅 is there an inverse instruction to the movemask instruction in intel avx2? 了解有关将位掩码转换为矢量掩码的更多详细信息,包括我上面描述的广播 + 可变移位策略。 (Chris Dodd 的 LUT 适用于 4 位 -> 4 个双字,因为表很小。除此之外,您还需要 ALU。)

【讨论】:

    猜你喜欢
    • 1970-01-01
    • 2014-06-23
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多