【问题标题】:Open MP: SIMD compatible function in SIMD loops?打开 MP:SIMD 循环中的 SIMD 兼容功能?
【发布时间】:2019-01-14 16:33:44
【问题描述】:

通常,我可能会编写一个 simd 循环,例如:

float * x = (float *) malloc(10 * sizeof(float));
float * y = (float *) malloc(10 * sizeof(float));

for(int i = 0; i < 10; i++)
    y[i] = 10;

#pragma omp simd
for(int i = 0; i < 10; i++)
    x[i] = y[i]*y[i];

假设我有两个任务:

float square(float x) {
    return x * x;
}
float halve(float x) {
    return x / 2.;
}

还有一个 omp 循环原语:

void apply_simd(float * x, float * y, int length, float (*simd_func)(float c)){
    #pragma omp simd
    for(int i = 0; i < length; i++)
         x[i] = simd_func(y[i])
}

这在 SIMD 的参数范围内是否合法?或者编译器生成的代码是否会比我显式内联所有代码效率低?

是否写作:

float inline square(float x){ ... } 

改变什么?或者只有当我明确地写下本机函数/运算符的操作时,我才能期望从 SIMD 中受益?

【问题讨论】:

    标签: c++ gcc g++ openmp simd


    【解决方案1】:

    是的,启用优化 (-O3 -march=native),现代编译器可以通过函数指针可靠地内联如果满足这些条件:

    • 函数指针有一个编译时常量值
    • 它指向一个函数,其定义编译器可以看到

    这听起来很容易确保,但是如果此代码用于 Unix/Linux 上的共享库(使用 -fPIC 编译),那么符号插入规则意味着 float halve(float x) { return x * 0.5f; } 1 即使在同一个翻译单元中也不能内联。见Sorry state of dynamic libraries on Linux

    即使在构建共享库时也可以使用inline 关键字允许内联;就像static 一样,它允许编译器根本不发出函数的独立定义,如果它决定在每个调用点内联。

    halvesquareapply_simd 上使用inline。 (因为apply_simd 需要内联到将halve 作为函数arg 传递的调用者。apply_simd 的独立定义是无用的,因为它不能内联未知函数。)如果他们在@ 987654340@ 而不是 .h,您也可以将它们设为 static,或者改为将它们设为 inline


    一次完成尽可能多的工作

    我怀疑你想写这样一些非常低效的东西:

    apply_simd(x, y, length, halve);   // copy y to x
    apply_simd(x, x, length, square);  // then update x in-place
    // NEVER DO THIS, make one function that does both things
    // with gcc and clang, compiles as written to two separate loops.
    

    只进行0.5f 的复制和乘法的循环通常会成为内存带宽的瓶颈。像 Haswell(或 Skylake)这样的现代 CPU 的 FMA/mul(或添加)吞吐量(每个时钟 2 个 256 位向量)是存储带宽(每个时钟到 L1d 的 1 个 256 位向量)的两倍。 计算强度很重要。不要通过编写多个执行单独的琐碎操作的循环来削弱您的代码

    如果展开任何循环,或者如果数据不适合 L1d,则 SIMD x[i] = 0.25f * y[i]*y[i] 的吞吐量将与单独使用其中任何一个操作的吞吐量相同。

    我检查了 g++ 8.2 和 clang++ 6.0 on the Godbolt compiler explorer 的 asm 输出。即使__restrict 告诉它 x 和 y 不重叠,编译器仍然创建了 2 个单独的循环。


    将 lambda 作为函数指针传递

    我们可以使用 lambda 轻松地将任意操作组合成单个函数,并将其作为函数指针传递。这解决了上述创建两个单独循环的问题,同时仍然为您提供了将循环包装在函数中的所需语法。

    如果您的 halve(float) 函数是不平凡事物的占位符,您可以在 lambda 中使用它来与其他事物组合。例如square(halve(a))

    在早期的 C++ 标准中,您需要将 lambda 分配给函数指针。 (Lambda as function parameter)

    // your original function mostly unchanged, but with size_t and inline
    inline  // allows inlining even with -fPIC
    void apply_simd(float * x, const float *y, size_t length, float (*simd_func)(float c)){
        #pragma omp simd
        for(size_t i = 0; i < length; i++)
             x[i] = simd_func(y[i]);
    }
    

    C++11 调用者:

    // __restrict isn't needed with OpenMP, but you might want to assert non-overlapping for better auto-vectorization with non-OpenMP compilers.
    void test_lambda(float *__restrict x, const float *__restrict y, size_t length)
    {
        float (*funcptr)(float) = [](float a) -> float {
             float h=0.5f*a; // halve first allows vmulps with a memory source operand
             return h*h;    // 0.25 * a * a doesn't optimize to that with clang :/
        };
    
        apply_simd(x, y, length, funcptr);
    }
    

    在 C++17 中它更容易,并且可以使用字面匿名 lambda:

    void test_lambda17(float *__restrict x, const float *__restrict y, size_t length)
    {
        apply_simd(x, y, length, [](float a) {
            float h = 0.5f*a;
            return h * h;
          }
        );
    }
    

    它们都可以使用 gcc 和 clang 高效编译到像这样的内部循环 (Godbolt compiler explorer)

    .L4:
        vmulps  ymm0, ymm1, YMMWORD PTR [rsi+rax]
        vmulps  ymm0, ymm0, ymm0
        vmovups YMMWORD PTR [rdi+rax], ymm0
        add     rax, 32
        cmp     rax, rcx
        jne     .L4
    

    clang 展开一些,并且可能接近每个时钟加载+存储一个 256 位向量,并进行 2 次乘法运算。 (非索引寻址模式可以通过展开隐藏两个指针增量来实现。愚蠢的编译器。:/)


    Lambda 或函数指针作为模板参数

    使用本地 lambda 作为模板参数(在函数内部定义),编译器绝对可以始终内联。但是(由于 gcc 错误)这目前不可用。

    但只有一个函数指针,它实际上并不能帮助捕捉您忘记使用inline 关键字或破坏编译器内联能力的情况。这仅意味着函数地址必须是动态链接时间常量(即直到动态库的运行时绑定才知道),因此它不会使您免于符号插入。在使用-fPIC编译时,编译器仍然不知道它可以看到的全局函数的版本是否是在链接时实际解析的版本,或者LD_PRELOAD 或主可执行文件中的符号将覆盖它。所以它只是发出从 GOT 加载函数指针的代码,并在循环中调用它。 SIMD 当然是不可能的。

    它确实可以阻止你通过传递函数指针而不总是内联的方式来射击自己。不过,可能使用 constexpr 在模板中使用之前,您仍然可以将它们作为 args 传递。 因此,如果不是因为 gcc 错误阻止您将其与 lambda 一起使用,您可能想要使用它。

    C++17 允许将没有捕获的自动存储 lambda 作为函数对象传递。 (以前的标准需要外部或内部 (static) 链接作为模板参数传递的函数。)

    template <float simd_func(float c)>
    void apply_template(float *x, const float *y, size_t length){
        #pragma omp simd
        for(size_t i = 0; i < length; i++)
             x[i] = simd_func(y[i]);
    }
    
    
    void test_lambda(float *__restrict x, const float *__restrict y, size_t length)
    {
        // static // even static doesn't help work around the gcc bug
        constexpr auto my_op = [](float a) -> float {
             float h=0.5f*a; // halve first allows vmulps with a memory source operand
             return h*h;    // 0.25 * a * a doesn't optimize to that with clang :/
        };
    
        // I don't know what the unary + operator is doing here, but some examples use it
        apply_lambda<+my_op>(x, y, length); // clang accepts this, gcc doesn't
    }
    

    clang 编译得很好,但 g++ 错误地拒绝了它,即使使用 -std=gnu++17

    不幸的是,gcc 在以这种方式使用 lambda 时存在错误 (83258)。 详情请参阅 Can I use the result of a C++17 captureless lambda constexpr conversion operator as a function pointer template non-type argument?

    不过,我们可以在模板中使用常规函数。

    // `inline` is still necessary for it to actually inline with -fPIC (in a shared lib)
    inline float my_func(float a) { return 0.25f * a*a;}
    
    void test_template(float *__restrict x, const float *__restrict y, size_t length)
    {
        apply_lambda<my_func>(x, y, length);   // not actually a lambda, just a function
    }
    

    然后我们从 g++8.2 -O3 -fopenmp -march=haswell 得到一个像这样的内部循环。请注意,我使用0.25f * a * a; 而不是首先使用halve,以查看我们得到了什么样的错误代码。这就是 g++8.2 所做的。

    .L25:
        vmulps  ymm0, ymm1, YMMWORD PTR [rsi+rax]   # ymm0 = 0.25f * y[i+0..7]
        vmulps  ymm0, ymm0, YMMWORD PTR [rsi+rax]   # reload the same vector again
        vmovups YMMWORD PTR [rdi+rax], ymm0        # store to x[i+0..7]
        add     rax, 32
        cmp     rax, rcx
        jne     .L25
    

    如果 gcc 不使用索引寻址模式(在 Haswell/Skylake 上 stops it from micro-fusing),则重新加载相同的向量两次以保存指令可能是一个好主意。所以这个循环实际上是 7 个微指令,每次迭代最多运行 7/4 个循环。

    根据英特尔的优化手册,随着展开,对于宽向量而言,接近每个时钟 2 次读取 + 1 次写入的限制显然是持续运行的问题。 (他们说 Skylake 可能会维持每个时钟 82 个字节,而峰值是 96 个加载 + 存储在一个时钟中。)如果数据不知道是对齐的,这是特别不明智的,并且 gcc8 已经切换到未知的乐观策略 -对齐数据:使用未对齐的加载/存储,并让硬件处理没有 32 字节对齐的情况。 gcc7 及更早版本在主循环之前对齐指针,并且只加载一次向量。


    脚注 1:幸运的是 gcc 和 clang 可以将 x / 2. 优化为 x * 0.5f,避免升级为 double

    在没有-ffast-math 的情况下可以使用乘法而不是除法,因为0.5f 可以精确地表示为float,这与分母不是2 的幂的分数不同。

    但请注意,0.5 * x 确实没有优化为 0.5f * x; gcc 和 clang 实际上会扩展到 double 并返回。我不确定这是否是与x / 2. 相比错过的优化,或者当它可以完全表示为float 时,是否存在真正的语义差异阻止它将双常量优化为浮点数。

    【讨论】:

    • 这是我遇到过的任何主题的最全面的知识部署之一。惊人的。谢谢你。关于所有这些,目标是重用语法和代码标记——但你的反应是深刻的。我会检查这些这个神螺栓。
    猜你喜欢
    • 2015-05-23
    • 1970-01-01
    • 2013-11-03
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2017-12-16
    • 2017-10-13
    • 1970-01-01
    相关资源
    最近更新 更多