SIMD минмаг и макмаг

Question

Я хочу реализовать SIMD-функции minmag и maxmag. Насколько я понимаю, эти функции

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

Я хочу использовать их для чисел с плавающей запятой и двойных чисел, а моим целевым оборудованием является Haswell. Что мне действительно нужно, так это код, который рассчитывает оба. Вот что у меня есть для SSE4.1 для двойного (код 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 Instinsic идет в порт 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, все из которых идут в порт 1 на Haswell. Используя мою реализацию two_sum_MinMax в этом вопросе он уменьшается до 16 add/subs на порту 1, но имеет худшую задержку и все еще медленнее. Вы можете увидеть сборку для этой функции и узнать больше о том, почему я забочусь об этом, в разделе оптимизация для быстрого умножения, но медленного сложения fma-and-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);
 }