SIMD minmag和maxmag

Question

我想实现SIMD minmag和maxmag函数.据我所知,这些功能是

minmag(a,b) = |a|<|b| ? a : b
maxmag(a,b) = |a|>|b| ? a : b

我希望这些浮动和双重,我的目标硬件是Haswell.我真正需要的是计算两者的代码.以下是我对SSE4.1 for double的看法(AVX代码几乎完全相同)

static inline void maxminmag(__m128d & a, __m128d & b) {
    __m128d mask    = _mm_castsi128_pd(_mm_setr_epi32(-1,0x7FFFFFFF,-1,0x7FFFFFFF));
    __m128d aa      = _mm_and_pd(a,mask);
    __m128d ab      = _mm_and_pd(b,mask);
    __m128d cmp     = _mm_cmple_pd(ab,aa);
    __m128d cmpi    = _mm_xor_pd(cmp, _mm_castsi128_pd(_mm_set1_epi32(-1)));
    __m128d minmag  = _mm_blendv_pd(a, b, cmp);
    __m128d maxmag  = _mm_blendv_pd(a, b, cmpi);
    a = maxmag, b = minmag;
}

但是,这并不像我想的那样高效.是否有更好的方法或至少值得考虑的替代方案？我想尝试避免端口1,因为我已经使用该端口进行了许多添加/减少.该_mm_cmple_pd禀进入端口1.

我感兴趣的主要功能是:

//given |a| > |b|
static inline doubledouble4 quick_two_sum(const double4 & a, const double4 & b)  {
    double4 s = a + b;
    double4 e = b - (s - a);
    return (doubledouble4){s, e};
}

所以我真正追求的就是这个

static inline doubledouble4 two_sum_MinMax(const double4 & a, const double4 & b) {
    maxminmag(a,b);       
    return quick_to_sum(a,b);
}

编辑:我的目标是two_sum_MinMax比two_sum下面更快:

static inline doubledouble4 two_sum(const double4 &a, const double4 &b) {
        double4 s = a + b;
        double4 v = s - a;
        double4 e = (a - (s - v)) + (b - v);
        return (doubledouble4){s, e};
}

编辑:这是我追求的终极功能.它有20个add/subs,所有这些都将转到Haswell的端口1.使用我two_sum_MinMax在这个问题中的实现将它降低到端口1上的16个add/subs,但是它的延迟更差并且仍然更慢.你可以看到这个函数的程序集,并在优化快速乘法 - 但是 - 慢 - 加法 - fma - 和 - doubledouble中阅读更多关于我为什么关心它的原因

static inline doublefloat4 adddd(const doubledouble4 &a, const doubledouble4 &b) {
        doubledouble4 s, t;
        s = two_sum(a.hi, b.hi);
        t = two_sum(a.lo, b.lo);
        s.lo += t.hi;
        s = quick_two_sum(s.hi, s.lo);
        s.lo += t.lo;
        s = quick_two_sum(s.hi, s.lo);
        return s;
        // 2*two_sum, 2 add, 2*quick_two_sum = 2*6 + 2 + 2*3 = 20 add
}

Answer 1

这是一个使用较少指令的替代实现:

static inline void maxminmag_test(__m128d & a, __m128d & b) {
    __m128d cmp     = _mm_add_pd(a, b); // test for mean(a, b) >= 0
    __m128d amin    = _mm_min_pd(a, b);
    __m128d amax    = _mm_max_pd(a, b);
    __m128d minmag  = _mm_blendv_pd(amin, amax, cmp);
    __m128d maxmag  = _mm_blendv_pd(amax, amin, cmp);
    a = maxmag, b = minmag;
}

它使用了一个有点微妙的算法(见下文),结合我们可以使用符号位作为选择掩码的事实.

它还使用@ EOF建议只使用一个掩码并切换操作数顺序,这样可以保存指令.

我用少量案例测试了它,它似乎与你原来的实现相匹配.

算法:

 if (mean(a, b) >= 0)       // this can just be reduced to (a + b) >= 0
 {
     minmag = min(a, b);
     maxmag = max(a, b);
 }
 else
 {
     minmag = max(a, b);
     maxmag = min(a, b);
 }