3

我想构建一个表示多个（比如说N）算术类型的数据类型，并使用运算符重载提供与算术类型相同的接口，这样我就得到了像 Agner Fog 的vectorclass这样的数据类型。

请看这个例子：Godbolt

#include <array>

using std::size_t;

template<class T, size_t S>
class LoopSIMD : std::array<T,S>
{
public:
    friend LoopSIMD operator*(const T a, const LoopSIMD& x){
        LoopSIMD result;
        for(size_t i=0;i<S;++i)
            result[i] = a*x[i];
        return result;
    }

    LoopSIMD& operator +=(const LoopSIMD& x){
        for(size_t i=0;i<S;++i){
            (*this)[i] += x[i];
        }
        return *this;
    }
};

constexpr size_t N = 7;
typedef LoopSIMD<double,N> SIMD;

SIMD foo(double a, SIMD x, SIMD y){
    x += a*y;
    return x;
}

对于一定数量的元素，这似乎工作得很好，gcc-10 为 6，clang-11 为 27。对于大量元素，编译器不再使用 FMA（例如vfmadd213pd）操作。相反，它们分别进行乘法（例如vmulpd）和加法（例如vaddpd）。

问题：

  • 这种行为有充分的理由吗？
  • 是否有任何编译器标志可以让我可以为 gcc 增加 6 和为 clang 增加 27 的上述值？

谢谢！

4

3 回答 3

0

对于 gcc 10.2，我做了以下操作，并获得了一些相当不错的结果，与-Ofast -march=skylake -ffast-math您的 Godbolt 链接相同。

friend LoopSIMD operator*(const T a, const LoopSIMD& x) {
    LoopSIMD result;
    std::transform(x.cbegin(), x.cend(), result.begin(),
                   [a](auto const& i) { return a * i; });
    return result;
}

LoopSIMD& operator+=(const LoopSIMD& x) {
    std::transform(this->cbegin(), this->cend(), x.cbegin(), this->begin(),
                   [](auto const& a, auto const& b) { return a + b; });
    return *this;
}

std::transform有一些疯狂的超载，所以我想我需要解释一下。

第一个重载捕获a，将每个值相乘，并将其存储回结果的开头。

第二个重载的作用是zip将两个值从x和相加this并将结果存储回this。

如果你没有结婚operator+=，你可以像这样operator*创建自己的fma

    LoopSIMD& fma(const LoopSIMD& x, double a ){
        std::transform_inclusive_scan(
            x.cbegin(),
            x.cend(),
            this->begin(),
            std::plus{},
            [a](auto const& i){return i * a;},
            0.0);
        return *this;
    }

这需要 c++17，但会循环将 SIMD 指令保留在

foo(double, LoopSIMD<double, 40ul>&, LoopSIMD<double, 40ul> const&):
        xor     eax, eax
        vxorpd  xmm1, xmm1, xmm1
.L2:
        vfmadd231sd     xmm1, xmm0, QWORD PTR [rsi+rax]
        vmovsd  QWORD PTR [rdi+rax], xmm1
        add     rax, 8
        cmp     rax, 320
        jne     .L2
        ret
于 2020-11-04T15:07:48.153 回答
0

我发现给定示例的改进。

在循环之前添加#pragma omp simdGCC 设法使 FMA 优化达到N=71.

https://godbolt.org/z/Y3T1rs37W

如果使用 AVX512，尺寸可能会更大：

https://godbolt.org/z/jWWPP7W5G

于 2021-09-20T10:29:11.580 回答
0

您也可以简单地制作自己的 fma 函数：

template<class T, size_t S>
class LoopSIMD : std::array<T,S>
{
public:
    friend LoopSIMD fma(const LoopSIMD& x, const T y, const LoopSIMD& z) {
        LoopSIMD result;
        for (size_t i = 0; i < S; ++i) {
            result[i] = std::fma(x[i], y, z[i]);
        }
        return result;
    }
    friend LoopSIMD fma(const T y, const LoopSIMD& x, const LoopSIMD& z) {
        LoopSIMD result;
        for (size_t i = 0; i < S; ++i) {
            result[i] = std::fma(y, x[i], z[i]);
        }
        return result;
    }
    // And more variants, taking `const LoopSIMD&, const LoopSIMD&, const T`, `const LoopSIMD&, const T, const T`, etc
};

SIMD foo(double a, SIMD x, SIMD y){
    return fma(a, y, x);
}

但是为了首先进行更好的优化，您应该对齐阵列。如果您这样做，您的原始代码会很好地优化：

constexpr size_t next_power_of_2_not_less_than(size_t n) {
    size_t pow = 1;
    while (pow < n) pow *= 2;
    return pow;
}

template<class T, size_t S>
class LoopSIMD : std::array<T,S>
{
public:
    // operators
} __attribute__((aligned(next_power_of_2_not_less_than(sizeof(T[S])))));

// Or with a c++11 attribute
/*
template<class T, size_t S>
class [[gnu::aligned(next_power_of_2_not_less_than(sizeof(T[S])))]] LoopSIMD : std::array<T,S>
{
public:
    // operators
};
*/

SIMD foo(double a, SIMD x, SIMD y){
    x += a * y;
    return x;
}
于 2020-11-04T17:39:54.400 回答