为什么 C++ 编译器不做更好的常量折叠？答案

【问题标题】：Why don't C++ compilers do better constant folding?为什么 C++ 编译器不做更好的常量折叠？
【发布时间】：2018-08-31 20:25:56
【问题描述】：

我正在研究加速大部分 C++ 代码的方法，这些代码具有用于计算雅可比的自动导数。这涉及在实际残差中做一些工作，但大部分工作（基于分析的执行时间）是计算雅可比。

这让我很惊讶，因为大多数雅可比是从 0 和 1 向前传播的，所以工作量应该是函数的 2-4 倍，而不是 10-12 倍。为了模拟大量的 jacobian 工作是什么样的，我做了一个超级最小的例子，只有一个点积（而不是真实情况下的 sin、cos、sqrt 等），编译器应该能够优化为单个返回值：

#include <Eigen/Core>
#include <Eigen/Geometry>

using Array12d = Eigen::Matrix<double,12,1>;

double testReturnFirstDot(const Array12d& b)
{
    Array12d a;
    a.array() = 0.;
    a(0) = 1.;
    return a.dot(b);
}

应该和

一样

double testReturnFirst(const Array12d& b)
{
    return b(0);
}

我很失望地发现，在没有启用快速数学的情况下，GCC 8.2、Clang 6 或 MSVC 19 都无法对充满 0 的矩阵的天真点积进行任何优化。即使使用快速数学 (https://godbolt.org/z/GvPXFy)，GCC 和 Clang 中的优化也很差（仍然涉及乘法和加法），并且 MSVC 根本不做任何优化。

我没有编译器的背景，但这有什么原因吗？我相当肯定，在大部分科学计算中，能够进行更好的常数传播/折叠将使更多优化变得明显，即使常数折叠本身并没有导致加速。

虽然我对解释为什么不在编译器方面这样做很感兴趣，但我也对我可以在实际方面做些什么来使我自己的代码在面对这些类型的模式时更快。

【问题讨论】：

浮点数不是实数，它们有严格的正确性要求，明显的优化违反了这些要求。例如。 (1.0 / 3.0) * 3.0 != (1.0 * 3.0)/3.0 因为舍入行为是完全指定的，所以你不能简单地取消 3.
答案取决于dot的实现。可能它不仅仅是一个带有累积的for 循环，而是涉及重新缩放。难怪编译器不能优化它。
-ffast-math 的意思是说“没有必要遵守标准”。 MSVC 等价于 fast-math 是 /fp:fast 如果你指定它，你可能会发现它做了一些优化。
添加-ffast-math 后，剩下的“问题”是显式矢量化，请参阅我的答案。
你可以在godbolt中看到选项。 -O3 用于 gcc/clang，/Ox 用于 MSVC。

标签： c++ compiler-construction eigen automatic-differentiation ceres-solver

【解决方案1】：

这是因为 Eigen 将您的代码显式矢量化为剩余 4 个组件寄存器中的 3 个 vmulpd、2 个 vaddpd 和 1 个水平缩减（假设为 AVX，仅使用 SSE，您将获得 6 个 mulpd 和 5 个 addpd）。使用-ffast-math GCC 和 clang 可以删除最后 2 个 vmulpd 和 vaddpd（这就是他们所做的），但它们不能真正替换 Eigen 显式生成的剩余 vmulpd 和水平缩减。

那么如果你通过定义EIGEN_DONT_VECTORIZE 来禁用Eigen 的显式矢量化呢？然后你会得到你所期望的 (https://godbolt.org/z/UQsoeH) 但其他代码段可能会变得慢得多。

如果你想在本地禁用显式矢量化并且不怕弄乱 Eigen 的内部，你可以在Matrix 中引入一个DontVectorize 选项，并通过专门针对这个traits<> 类型的traits<> 禁用矢量化：

static const int DontVectorize = 0x80000000;

namespace Eigen {
namespace internal {

template<typename _Scalar, int _Rows, int _Cols, int _MaxRows, int _MaxCols>
struct traits<Matrix<_Scalar, _Rows, _Cols, DontVectorize, _MaxRows, _MaxCols> >
: traits<Matrix<_Scalar, _Rows, _Cols> >
{
  typedef traits<Matrix<_Scalar, _Rows, _Cols> > Base;
  enum {
    EvaluatorFlags = Base::EvaluatorFlags & ~PacketAccessBit
  };
};

}
}

using ArrayS12d = Eigen::Matrix<double,12,1,DontVectorize>;

完整示例：https://godbolt.org/z/bOEyzv

【讨论】：

为什么编译器不能优化剩下的向量指令？是 QoI 问题还是技术原因？
@Rakete1111 大概是因为没有人坐下来编写足够详细的规则/模型，编译器将通过这些规则/模型跟踪通过向量指令的常量传播。一些规则（例如乘以或加 0.0）显然已经包含在内，但可能很难使它们像标量规则那样具有包容性。
这在技术上可以通过“取消矢量化”代码来实现，但这会违背用户明确要求的内容，因此这是否合理仍有待商榷。
您对编译器的要求非常高……要让它做您想做的事情，就需要它真正开发一些机器洞察力来了解问题的细节。这并非不可能，但不是编译器作者关注的那种思想。对我们人类来说，很明显，N 维中的点积，其中除了一个向量的第一个元素之外的所有元素都是零，这是一个微不足道的乘法，但这不是编译器的重点。此外，如上所述，为了保持一致性，浮点必须做它所做的事情。出于这个原因，Python 使用了许多已有 30 年历史的 Fortran 库。
您能否就为什么 MSVC 无法优化此代码提供任何见解？也许有解决方法？

【解决方案2】：

我很失望地发现，在没有启用快速数学的情况下，GCC 8.2、Clang 6 或 MSVC 19 都无法对充满 0 的矩阵的幼稚点积进行任何优化。

不幸的是，他们别无选择。由于 IEEE 浮点数有符号零，添加 0.0 不是恒等运算：

-0.0 + 0.0 = 0.0 // Not -0.0!

同样，乘以零并不总是得出零：

0.0 * Infinity = NaN // Not 0.0!

因此，编译器根本无法在点积中执行这些常量折叠，同时保持 IEEE 浮点合规性 - 据他们所知，您的输入可能包含有符号零和/或无穷大。

您必须使用-ffast-math 来获得这些折叠，但这可能会产生不良后果。您可以使用特定标志获得更细粒度的控制（来自http://gcc.gnu.org/wiki/FloatingPointMath）。根据上面的解释，添加以下两个标志应该允许常量折叠：
-ffinite-math-only, -fno-signed-zeros

确实，您可以通过这种方式获得与-ffast-math 相同的程序集：https://godbolt.org/z/vGULLA。您只放弃有符号的零（可能不相关）、NaN 和无穷大。据推测，如果您仍要在代码中生成它们，您会得到未定义的行为，因此请权衡您的选择。

至于为什么即使使用-ffast-math，您的示例也没有得到更好的优化：那是关于 Eigen。大概他们对矩阵运算进行了矢量化，这对于编译器来说更难看穿。使用这些选项正确优化了一个简单的循环：https://godbolt.org/z/OppEhY

【讨论】：

只有 clang 优化了 for 循环，gcc 不这样做。

【解决方案3】：

强制编译器优化 0 和 1 乘法的一种方法是手动展开循环。为简单起见，让我们使用

#include <array>
#include <cstddef>
constexpr std::size_t n = 12;
using Array = std::array<double, n>;

然后我们可以使用折叠表达式（如果不可用则递归）实现一个简单的dot 函数：

<utility>
template<std::size_t... is>
double dot(const Array& x, const Array& y, std::index_sequence<is...>)
{
    return ((x[is] * y[is]) + ...);
}

double dot(const Array& x, const Array& y)
{
    return dot(x, y, std::make_index_sequence<n>{});
}

现在让我们看看你的函数

double test(const Array& b)
{
    const Array a{1};    // = {1, 0, ...}
    return dot(a, b);
}

使用-ffast-math gcc 8.2 produces:

test(std::array<double, 12ul> const&):
  movsd xmm0, QWORD PTR [rdi]
  ret

clang 6.0.0 也是如此：

test(std::array<double, 12ul> const&): # @test(std::array<double, 12ul> const&)
  movsd xmm0, qword ptr [rdi] # xmm0 = mem[0],zero
  ret

例如，对于

double test(const Array& b)
{
    const Array a{1, 1};    // = {1, 1, 0...}
    return dot(a, b);
}

我们得到

test(std::array<double, 12ul> const&):
  movsd xmm0, QWORD PTR [rdi]
  addsd xmm0, QWORD PTR [rdi+8]
  ret

添加。 Clang 展开 for (std::size_t i = 0; i < n; ++i) ... 循环，没有所有这些折叠表达式技巧，gcc 没有，需要一些帮助。

【讨论】：