如何找到为每个双精度值产生较小输出的乘数？

Question

您如何找到最大和正的 IEEE-754 二进制 64 值 C，使得正的标准化二进制 64 值 A 与 C 的每个 IEEE-754 乘积都小于 A ?

我知道它一定接近 0.999999...但我想找到最大的那个。

假设从四舍五入到最近，与偶数相关。

Answer 1

有几种实验方法；这是C = 1 - ε的证明，其中ε是机器ε（即1与大于1的最小可表示数之间的距离。）

我们当然知道C < 1，所以尝试C = 1 - ε/2 是有意义的，因为它是下一个小于1 的可表示数字。 （ε/2 是因为C 在可表示数字的[0.5, 1) 桶中。）让我们看看它是否适用于所有A。

我将在本段中假设1 <= A < 2。如果A 和AC 都在“正常”区域，那么指数是多少并不重要，情况与指数2^0 相同。现在，C 的选择显然适用于A=1，所以我们只剩下区域1 < A < 2。查看A = 1 + ε，我们看到AC（精确值，不是四舍五入的结果）已经大于1；而对于A = 2 - ε，我们看到它小于2。这很重要，因为如果AC 介于1 和2 之间，我们知道AC 和round(AC) 之间的距离（即四舍五入到最接近的可表示值）最多为ε/2。现在，如果A - AC < ε/2，那么我们不想要的round(AC) = A。 （如果 A - AC = ε/2 则 可能 舍入到 A 给定正常 FP 舍入规则的“与偶数”部分，但让我们看看我们是否可以做得更好。）因为我们选择了C = 1 - ε/2，我们可以看到A - AC = A - A(1 - ε/2) = A * ε/2。因为它大于ε/2（记住，A>1），所以它距离A 足够远，可以舍入它。

但是！我们必须检查的A 的另一个值是可表示的最小正常值，因为AC 不在正常范围内，因此我们的“到最近的相对距离”规则不申请。我们发现，在这种情况下，A-AC 正好是该地区机器 epsilon 的一半。 “四舍五入到最近，平局”开始，产品四舍五入回到等于A。德拉特。

对C = 1 - ε 进行同样的处理，我们看到round(AC) < A，并且没有其他任何东西接近于向A 进行四舍五入（我们最终询问是否A * ε > ε/2，当然是哪个）。所以妙语是C = 1-ε/2几乎有效，但法线和非法线之间的界限把我们搞砸了，C = 1-ε 让我们进入了终点区。

Answer 2

由于浮点类型的性质，C 会根据 A 的值大小而有所不同。您可以使用nextafter 获得小于1 的最大值，这将是C 的粗略值

但是，如果 A 太大或太小，A*C 将与 A 相同。我无法从数学上证明 nextafter(1.0, 0) 将适用于所有可能的 A，因此我'我建议这样的解决方案

double largestCfor(double A)
{
    double C = nextafter(1.0, 0);
    while (C*A >= A)
        C = nextafter(C, 0);
    return C;
}

如果您想要一个适用于任何A 的C 值，即使C*A 可能不是最大的可能值，您也需要检查该类型可以表示的每个指数

double C = 1;
for (double A = 0x1p-1022; isfinite(A); A *= 2) // loop through all possible exponents
{
    double c = largestCfor(A);
    if (c < C) C = c;
}

我尝试在Ideone 上运行并得到了结果

C                 = 0.999999999999999777955395074969
nextafter(1.0, 0) = 0.999999999999999888977697537484

编辑：

0.999999999999999777955395074969 是0x1.ffffffffffffep-1，也是1 - DBL_EPSILON。这与上面 Sneftel 的证明一致