将单精度浮点数转换为双精度以进行除法

Question

作为高性能计算人员，我们倾向于尽可能默认使用单精度浮点数（float 或 real）。这是因为如果每个操作单独执行得更快，您每秒可以执行更多操作。

然而，与我一起工作的一位资深人士总是坚持（当需要准确性时）您应该暂时将单精度数据转换为双精度数据以执行除法。即：

float a, b;
float ans = ((double)a)/((double)b);

或

real :: a, b, ans
ans = real(dble(a)/dble(b))

取决于您使用的语言。在我看来，这看起来真的很难看，老实说，我什至不知道 ans 中的答案是否会比如果更准确您只是简单地以单点精度编写了 ans = a/b。

谁能告诉我在算术之前转换你的数字，专门用于执行除法，是否真的会得到更准确的答案？这是一个 language/compiler 的具体问题，还是由 IEEE 决定？这种准确性提高在什么数值下最明显？

任何有启发性的 comments/answers 将不胜感激。

Answer 1

是的，转换为双精度会给你更好的除法精度（或者，我应该说，精度）。可以说这取决于 IEEE，但这只是因为 IEEE 定义了格式和标准。 doubles 本质上比 floats 更精确，存储数字和除法。

要回答你的最后一个问题，这在大 a 和小（小于 1）b 时最为明显，因为这样你最终会得到一个非常大的商，在范围内所有浮点数的粒度都较小。

Answer 2

运行在 x86 (GCC 4.9.3) 上：

#include "stdio.h"
int main(int arc, char **argv)
{
    float a=73;
    float b=19;

    float ans1 = (a*a*a/b/b/b);
    float ans2 = ((double)a*(double)a*(double)a/(double)b/(double)b/(double)b);
    printf("plain: %f\n", ans1);
    printf("cast:  %f\n", ans2);
    return 0;
}

输出：

plain: 56.716282
cast:  56.716286

Windows计算器中的相同操作return:

56.716285172765709287068085726782

显然，第二个结果更准确。

Answer 3

这在很大程度上取决于所使用的平台。

使用非 SSE 指令的 80x86（或 1980 年代的 8087）使用 80 位精度（long double 或 real*10）执行所有运算。正是 "store" 指令将结果从数字处理器移动到内存中，从而失去了精度。

除非它是一个非常愚蠢的编译器，否则最大精度应该出现在

float a = something, b = something_else;
float ans = a/b;

由于要进行除法运算，单精度操作数加载后会被扩展精度，结果也会被扩展精度。

如果您正在做一些更复杂的事情并希望保持最大精度，请不要将中间结果存储在较小的变量中：

float a, b, c, d;

float prod_ad = a * d;
float prod_bc = b * c;
float sum_both = prod_ad + prod_bc;   // less accurate

与一次完成所有操作相比，这给出的结果不太精确，因为大多数编译器会生成将所有中间值保持在扩展精度的代码：

float a, b, c, d;

float sum_both = a * d + b * c;   // more accurate

基于 Eugeniu Rosca 的示例程序：

#include "stdio.h"
void main(void)
{
    float a=73;
    float b=19;

    long double a1 = a;
    long double b1 = b;

    float ans1 = (a*a*a/b/b/b);
    float ans2 = ((double)a*(double)a*(double)a/(double)b/(double)b/(double)b);
    float ans3 = a1*a1*a1/b1/b1/b1;
    long double ans4 = a1*a1*a1/b1/b1/b1;

    printf ("plain:  %.20g\n", ans1);
    printf ("cast:   %.20g\n", ans2);
    printf ("native: %.20g\n", ans3);
    printf ("full:   %.20Lg\n", ans4);
}

提供，无论优化级别

plain:  56.716281890869140625
cast:   56.71628570556640625
native: 56.71628570556640625
full:   56.716285172765709289

这表明对于微不足道的操作，差别不大。但是，将常量更改为更精确的挑战：

float a=0.333333333333333333333333;
float b=0.1;

提供

plain:  37.03704071044921875
cast:   37.037036895751953125
native: 37.037036895751953125
full:   37.037038692721614131

精度差异显示出更明显的效果。

Answer 4

float ans = ((double)a)/((double)b);

此 article 证明 ans 始终与 IEEE 754 算术的单精度除法和 FLT_EVAL_METHOD=0.

的计算结果相同

当FLT_EVAL_METHOD=1时，同样的属性也是平凡的。

当FLT_EVAL_METHOD=2时，我不确定。人们可能会将规则解释为 a/b 的 long double 计算必须首先四舍五入到 double，然后到 float。在这种情况下，它可能不如直接从 long double 舍入到 float 准确（后者产生正确舍入的结果，而前者在极少数情况下可能无法这样做，除非另一个定理像Figueroa 的申请表明这从未发生过）。

长话短说，对于任何现代合理的浮点计算平台（*），float ans = ((double)a)/((double)b);有任何好处都是迷信。你应该让你在问题中提到的资深人士展示一对 a, b 结果不同的值，更不用说更准确了。当然，如果他们坚持认为这样更好，那么提供一对值对他们来说应该没有问题。

(*) 记得将 -fexcess-precision=standard 与 GCC 一起使用以保持理智

将单精度浮点数转换为双精度以进行除法

Converting single-precision floating point numbers to double-precision for division

c

floating-point

fortran

hpc