我正在尝试将浮点双精度值x转换为具有12位（正确四舍五入）有效数字的小数。我假设x在10^110和10^111之间，这样它的十进制表示形式将是x.xxxxxxxxxxxE110。而且，为了好玩，我尝试只使用浮点运算。
我得到了下面的伪代码，其中所有的运算都是双精度运算，符号1e98表示最接近数学10^98的二重，1e98_2表示最接近算术减法10^98-1e98结果的二重。符号fmadd(X * Y + Z)用于具有操作数X、Y、Z的融合乘加运算。
y = x * 2^-1074; // exact
q = y / 1e98; // q is denormal and the significand of q interpreted
// as an integer is our candidate for the 12 decimal
// digits of x
r = fmadd(q * 1e98 - y); // close to 1e98 * (error made during the division)
// If 1e98_2 >= 0, we divided by a number that was smaller than we wished
// The correct answer may be q or q+1.
if (r and 1e98_2 have opposite signs)
{
return the significand of q;
}
s = copysign(2^-1074, r);
r1 = abs(r);
r2 = abs(1e98_2);
h = 1e98 * 0.5 * 2^-1074;
Set rounding mode to downwards
r3 = fmadd(r2 * q + r1);
if (r3 < h)
{
return the significand of q;
}
else
{
return significand of (q + s)
}
我对上述伪代码中的混乱表示歉意，但我还不太清楚，因此提出了以下问题：
第一个fmadd是否按预期工作（计算1e98*（除法过程中出现的错误））？
标志。我无法说服自己他们是对的。但我也不能说服自己他们错了。
关于这种算法可能产生错误结果的频率，有什么想法，也许是争论吗？
如果它真的有效，如果"q=y/1e98"改为"q=y*1e-98"（使所有其他指令保持不变），算法是否有可能继续工作？
我还没有测试过这个算法。我没有任何一台带有fmadd指令的计算机，尽管我希望能找到一台这样我就可以执行上面的指令。
设y/d为精确运算，q=rnd(y/d)为四舍五入到最接近浮点值的结果
那么乘以d的真实误差是rt=(rnd(y/d)-y/d)*d=q*d-y，我们用fmadd执行的运算是r=rnd(q*d-y)
为什么q*d-y是精确的（fmadd不进行最后舍入）还不太清楚，但假设q*d的位数有限（<nbits(q)+nbits(d)），y的指数是q*d的指数（+/-1），并且由于误差是|rt|<0.5*ulp(q)*d，这意味着第一个nbits(q)正在消失。。。这回答了问题1。
所以q*1e98 - y = r，其中|r|*2^1074 <= 0.5e98 < 5*10^98（第二不等式是幸运的）
q*(10^98) - y = r + (10^98-1e98)*q，其中|10^98-1e98|*q*2^1074 <= 0.5e95（假设精度至少为15位，log(2^53)/log(10) > 15）
所以你问|q*(10^98)-y|*2^1074>5*10^97
你有一个|q*(10^98)-y|的近似值，即r+1e98_2*q
由于|r| < 5*10^98和|r+(10^98-1e98)*q|<|r|（如果符号相反），我认为这肯定回答了问题2。但如果1e98_2是<0.
如果r和1e98_2具有相同的符号，则可能超过5*10^97，因此您将进一步讨论r3 = 1e98_2*q + r与h=0.5e98*2^-1074
对于问题3，乍一看，我认为有两件事可能会使算法失败：
-
1e98_2不准确（约为10^98-1e98-1e98_2 = -3.6e63） -
并且CCD_ 39不是CCD_。
真正的误差r3t大约是(1e98_2-3e63)*q + r < r3（并且只有当>0时我们才感兴趣，因为1e98_2>0）。
因此，当真实误差r3t低于真实平局ht时，误差r3的近似值落在近似平局h之上可能导致不正确的舍入。有可能吗？如果有，你的问题3的频率有多高？
为了减轻上述不平等风险，您尝试截断r3的大小，从而截断r3 <= 1e98_2*q + r。我觉得有点累了，无法对误差范围进行真正的分析。。。
因此，我扫描了一个错误，我发现的第一个失败的例子是1.000000001835e110（我假设正确地四舍五入到最接近的双，但实际上是1000000000183.4999998415379982112042942630528225695526491963291846957919215885146546696544423465444844268032e98）。
在这种情况下，r和1e98_2具有相同的符号，并且
-
(x/1e98) > 1000000000183.50000215 -
q有效位因此四舍五入为1000000000184 -
r3>h（r3*2^1074约为5000001584620017e97），我们错误地增加了q+s，而它本应是q-s，肯定是一个错误。
我的答案是：
-
是的，
r=fmadd(q * 1e98 - y)正好是1e98*（除法时出错），但我们不在乎除法，它只是提供了一个猜测，重要的是减法是精确的。 -
是的，符号是正确的，因为
|r| < 5*10^98，如果符号相反，则为|r+(10^98-1e98)*q|<|r|。但如果1e98_2是<0. -
以第一个失败的例子
(1.0000000001835e110 - 1.0e110)/1.0e110 ulp -> 1.099632e6为例，一个非常天真的猜想是，一百万分之一的情况下，r3正在下降到h……因此，一旦q+s校正为q-s，r3>h而r3t<ht的出现在任何情况下都远小于1/100000。。。在感兴趣的范围内有超过10^15的双打，所以考虑到这不是一个严肃的答案。。。 -
是的，上面的讨论只是关于猜测q，与它的产生方式无关，以及1中的减法。仍然是准确的。。。