如何使浮点计算具有确定性？答案

【问题标题】：How can floating point calculations be made deterministic?如何使浮点计算具有确定性？
【发布时间】：2011-09-09 18:16:23
【问题描述】：

浮点计算在处理器上既不关联也不分布。所以，

(a + b) + c 不等于a + (b + c)

而a * (b + c) 不等于a * b + a * c

是否有任何方法可以执行不会给出不同结果的确定性浮点计算。当然，在单处理器上它是确定性的，但在多线程程序中，如果线程相加，例如，它就不是确定性的，因为线程可能存在不同的交错。

所以我的问题是，如何在多线程程序中实现浮点计算的确定性结果？

【问题讨论】：

好问题。虽然答案可能是“你不能”或“使用任意精度算术”，但问它是合理的。
另外，你需要这个什么样的领域？在某些学科中这是一个真正的问题（例如计算几何），而在其他学科中，浮点计算没有问题（实际上是大多数领域，有一些真正重要的变通方法）。
“a + (b*c) 不等于 a*b + a*c”——它们曾经相等吗？
Alexander，我需要它来进行确定性执行，也就是说，如果我多次运行我的程序，它会给出相同的输出。这简化了调试。
奥斯汀，是的，现在更正了：）！

标签： c floating-point deterministic

【解决方案1】：

浮点是确定性的。相同的浮点运算，在相同的硬件上运行，总是产生相同的结果。没有黑魔法、噪音、随机性、模糊或人们通常归因于浮点的任何其他东西。牙仙没有出现，拿走你的结果的低位，在你的枕头下留下四分之一。

也就是说，某些通常用于大规模并行计算的阻塞算法在执行浮点计算的顺序方面是不确定的，这可能导致跨运行的非位精确结果。

你能做些什么呢？

首先，请确保您确实无法忍受这种情况。您可能尝试在并行计算中强制执行排序的许多事情都会损害性能。就是这样。

我还要指出，尽管阻塞算法可能会引入一定程度的不确定性，但它们通常提供的结果比天真的未阻塞串行算法更小的舍入误差（令人惊讶但确实如此！）。如果你能忍受朴素串行算法产生的错误，你可能也能忍受并行阻塞算法的错误。

现在，如果您真的非常需要跨运行的精确再现性，这里有一些建议，它们往往不会对性能产生太大的负面影响：

不要使用可以重新排序浮点计算的多线程算法。问题解决了。这并不意味着您根本不能使用多线程算法，只是您需要确保每个单独的结果仅由同步点之间的单个线程接触。请注意，如果处理得当，这实际上可以提高某些架构的性能，方法是减少内核之间的 D$ 争用。
在归约操作中，您可以让每个线程将其结果存储到数组中的索引位置，等待所有线程完成，然后按顺序累积数组的元素。这会增加少量的内存开销，但通常是可以容忍的，尤其是当线程数“很少”时。
想办法提升并行度。不是计算 24 个矩阵乘法，每个乘法使用并行算法，而是并行计算 24 个矩阵乘积，每个乘法使用串行算法。这也可能对性能有益（有时非常有益）。

还有很多其他方法可以处理这个问题。他们都需要思考和关心。并行编程通常可以。

【讨论】：

它只是机器代码级别的确定性。一个给定的 C 程序编译两次可能会产生不同的结果。
@Antimony 我很困惑，看来你和斯蒂芬在说不同的事情。那么 a+b+c 是不是和 b+c+a 一样呢？

【解决方案2】：

编辑：我已经删除了我的旧答案，因为我似乎误解了 OP 的问题。如果您想查看它，您可以阅读编辑历史记录。

我认为理想的解决方案是切换到每个线程都有一个单独的累加器。这避免了所有锁定，这会对性能产生巨大影响。您可以在整个操作结束时简单地对累加器求和。

或者，如果您坚持使用单个累加器，一种解决方案是使用“定点”而不是浮点。这可以通过在您的累加器中包含一个巨大的“偏差”项以将指数锁定在一个固定值来使用浮点类型来完成。例如，如果您知道累加器永远不会超过 2^32，则可以在 0x1p32 处启动累加器。这会将您锁定在小数点左侧的 32 位精度和 20 位小数精度（假设 double）。如果精度不够，您可以使用较小的偏差（假设累加器不会增长太大）或切换到long double。如果long double 是 80 位扩展格式，则 2^32 的偏差将给出 31 位的小数精度。

然后，每当您想实际“使用”累加器的值时，只需减去偏置项即可。

【讨论】：

他说的是这样的场景：你有两个线程，每个线程都向同一个累加器添加一个浮点值。如果不控制加法发生的顺序，可能会得到(accumulator + thread A value) + thread B value 或(accumulator + thread B value) + thread A value，由于四舍五入可能不相等。
好的，当时问题不是很清楚。无论如何，在我看来，每个线程都应该有自己的累加器，最后它们都可以加在一起。
这也将避免累加器上的锁定，这几乎肯定会使整个过程比单线程慢一个数量级......
@Stephen：我已经用一些新想法替换了我的答案。

【解决方案3】：

即使使用高精度定点数据类型也无法解决使所述方程的结果具有确定性的问题（某些情况除外）。正如 Keith Thompson 在评论中指出的那样，1/3 是一个无法以标准 base-10 或 base-2 浮点表示形式正确存储的值的反例（无论使用的精度或内存如何）。

根据特定需求，可以解决此问题的一种解决方案（它仍然有限制）是使用Rational number 数据类型（一种同时存储分子和分母的数据类型）。 Keith 建议将 GMP 作为此类库之一：

GMP 是一个免费的任意精度算术库，对有符号整数、有理数 和浮点数进行运算。精度没有实际限制...

它是否适合（或足够）这项任务是另一回事......

编码愉快。

【讨论】：

1.0/3.0 的不精确性与非确定性无关。
它与它有关 - 如果不是因为它，那么如果你有 1/3+1/3+1/3 或 (1+1+1)/ 也没关系3，因此这两种方法之间的非确定性选择不会使结果不确定。
@R.. 抓住了我！ :) 我在那里用它来表明并非所有数字都可以使用标准编码存储在有限空间中。但是，由于问题中没有划分，因此它的用处不大，因此可以作为示例。假设有界精度，当 FP 操作的顺序超过 + 和 * 时，仍然很容易提出示例。
我将 1.0/3.0 视为原子量而不是除法。
@R.. 但它不是标准尾数指数浮点或固定精度浮点中的精确原子量。如果存在除法，则很容易争论 Random8321 提供的情况：中间值只是 不能以任何内存总和的形式充分表示。仅通过加法和乘法，我认为这个论点必须被删除为“为了一个有界的精度”。

【解决方案4】：

使用十进制类型或支持这种类型的库。

【讨论】：

它仍然不是“确定性的”。对于任何有限的存储数量，操作的顺序仍然很重要。（但是，使用固定的 高精度 类型会有所帮助。）
pst，固定高精度类型？你这是什么意思？
@Metallic 有一个有限的内存集无法存储的精度水平。
如果您需要计算 1.0/3.0，十进制类型将无济于事。
Decimal 无法解决此问题。非关联性并不是二进制浮点格式所独有的。

【解决方案5】：

尝试将每个中间结果存储在 volatile 对象中：

volatile double a_plus_b = a + b;
volatile double a_plus_b_plus_c = a_plus_b + c;

这可能会对性能产生不利影响。我建议测量两个版本。

编辑：volatile 的目的是抑制即使在单线程环境中也可能影响结果的优化，例如~~更改操作顺序或~~将中间结果存储在更宽的寄存器中。它没有解决多线程问题。

EDIT2：还有一点需要考虑的是

一个浮动表达式可能是收缩，也就是说，被评估为好像这是一个原子操作，因此省略了隐含的舍入错误通过源代码和表达式求值方法。

这可以通过使用来抑制

#include <math.h>
...
#pragma STDC FP_CONTRACT off

参考：C99 standard（大 PDF），第 7.12.2 节和第 6.5 节第 8 段。这是 C99 特定的；有些编译器可能不支持它。

【讨论】：

如果我有多个线程添加到全局总和（关联属性）怎么办？这会造成不确定性，因为线程可能每次都以不同的顺序添加到总和中。
@MetallicPriest：是的，这确实会产生不确定性。考虑让线程将它们的结果附加到某种列表中（使用适当的代码以避免弄乱数据结构）。完成后，对列表进行排序并求和。
我看不出volatile 在这里有什么帮助。他已经在做 lock() 和 unlock() 希望同步内存。问题是不同的线程以不同的顺序完成工作，而 volatile 不会改变这一点。
编译器不能重新排序可能会改变结果的浮点运算，除非用户许可它这样做，不管有没有volatile关键字。
@Stephen：好点子。操作数的计算顺序没有定义，但这只是意味着x + y 可以先计算 x 然后 y，或者 y 然后 x；它不影响结果。操作数与运算符的关联定义明确。

【解决方案6】：

使用packed decimal。

【讨论】：

Decimal 无法解决此问题。非关联性并不是二进制浮点格式所独有的。
每种浮点格式（甚至十进制格式）都会导致舍入（唯一的选择是使用无限位数）。