使用双精度浮点运算

Question

我有一个函数，它需要两个字符串（浮点），运算和浮点位宽：

EvaluateFloat(const string &str1, const string &str2, enum operation/*add,subtract, multiply,div*/, unsigned int bit-width, string &output)

输入 str1 和 str2 可以是 float（32 位）或 double（64 位）。

如果将输入存储为双精度并执行双精度操作而不考虑位宽并且取决于位宽类型转换它是否为 32 位，是否可以。 例如

double num1 = atof(str1);
double num2 = atof(str2);
double result = num1 operation num2; //! operation will resolved using switch
if(32 == bit-width)
{
 float f_result = result;
 output = std::to_string(f_result);
}
else
{
 output = std::to_string(result);
}

如果我使用浮点类型执行浮点操作，我可以安全地假设 f_result 将完全相同，即

float f_num1 = num1;
float f_num2 = num2;
float f_result = f_num1 operation f_num2

PS：

1. 我们假设不会有任何级联操作，即 out = a + b + c 相反，它将转换为： temp = a +b out = temp + c
2. 我不关心 inf 和 nan 值。
3. 我正在尝试编码冗余，否则我有两个执行相同的操作 两次浮动一次，其他两次两次

Answer 1

C++ 没有指定 float 或 double 使用哪些格式。如果使用 IEEE-754 binary32 和 binary64，则 +、-、*、/ 或 sqrt 不会发生双舍入错误。给定float x 和float y，以下成立（左边是float 算术，右边是double）：

• x+y = (float) ((double) x + (double) y)。
• x-y = (float) ((double) x - (double) y)。
• x*y = (float) ((double) x * (double) y)。
• x/y = (float) ((double) x / (double) y)。
• sqrt(x) = (float) sqrt((double) x)。

这是根据 Samuel A. Figueroa del Cid 的论文A Rigorous Framework for Fully Supporting the IEEE Standard for Floating-Point Arithmetic in High-Level Programming Languages，2000 年 1 月，纽约大学。从本质上讲，double 在float 之外有很多数字（位），因此舍入到double 永远不会隐藏正确舍入到float 所需的信息以获取这些操作的结果。 （这不适用于一般操作；它取决于这些操作的属性。）在第 57 页，Figueroa del Cid 给出了一个表格，表明如果float 格式具有 p 位，那么，为避免双舍入错误，double 必须有 2p+1 位用于加法或减法，2p 用于乘法和除法，2p i>+2 表示sqrt。由于 binary32 的有效位有 24 位，double 有 53 位，因此这些都满足。 （详见论文。有一些注意事项，例如对于各种操作，p 必须至少为 2 或 4。）

Answer 2

根据标准对double 的浮点运算相当于以无限精度进行运算。如果我们将其转换为float，我们现在已经将其四舍五入了两次。一般来说，这并不等同于首先舍入到float。例如。 0.47 舍入到 0.5 舍入到 1，但 0.47 直接舍入到 0。正如 chtz 所提到的，两个浮点数的乘法应该总是恰好是两倍（使用 IEEE 数学，其中 double 的精度是 float 的两倍以上） ，所以当我们转换为 float 时，我们仍然只损失了一次精度，所以结果应该是一样的。同样，加法和减法应该不是问题。

除法不能精确地用双精度表示（甚至不是 1/3），所以我们可能认为除法有问题。但是，我已经在夜间运行了示例代码，尝试了超过 3 万亿个案例，但没有发现任何将原始除法作为 double 运行的案例会给出不同的答案。

#include <iostream>

int main() {
        long i=0;
        while (1) {
                float x = static_cast <float> (rand()) / static_cast <float> (RAND_MAX);
                float y = static_cast <float> (rand()) / static_cast <float> (RAND_MAX);
                float  f = x / y;
                double d = (double)x / (double)y;
                if(++i % 10000000 == 0) { std::cout << i << "\t" << x << "," << y << std::endl; }
                if ((float(d) !=  f)) {
                        std::cout << std::endl;
                        std::cout << x << "," << y << std::endl;
                        std::cout << std::hex << *(int*)&x << "," << std::hex << *(int*)&y << std::endl;
                        std::cout << float(d) - f << std::endl;
                        return 1;
                }
        }
}