通过正确舍入获得精确的双精度答案

【问题标题】：Exact double precision by correct rounding通过正确舍入获得精确的双精度
【发布时间】：2017-06-10 08:15:27
【问题描述】：

虽然我的问题听起来微不足道，但事实并非如此。希望你能帮助我。

我想在我的 .NET (C#) 项目中实现区间运算。这意味着每个数字都由下限和上限定义。这有助于解决诸如

之类的问题

1 / 3 = 0.333333333333333 (15 significant digits)

因为你会拥有

1 / 3 = [ 0.33333333333333 , 0.333333333333334 ] (14 significant digits each)

，所以我现在可以肯定，正确的答案在这两个数字之间。如果没有区间表示，我已经有了舍入误差（即 0.0000000000000003）。

为了实现这一点，我编写了自己的 Interval 类型，它重载了所有标准运算符，如 +-*/ 等。为了使这种类型正常工作，我需要能够在两个方向上舍入 1 / 3 的结果。将结果向下舍入将为我的区间提供下限，将结果向上舍入将为我提供区间的上限。

.NET 具有 Math.Round(double,int) 方法，可将 double 舍入到 int 小数位。看起来不错，但不能强制向上/向下舍入。 Math.Round(1.0/3.0,14) 会向下舍入，但也无法像这样实现向上舍入到 0.33...34。

但是你可能会说Math.Ceil和Math.Floor！好的，这些方法四舍五入到下一个较低或较高的整数。所以如果我想四舍五入到小数点后 14 位，我首先需要修改我的结果：

1 / 3 = 0.333333333333333 -> *E14 -> 33333333333333.3

所以现在我可以调用 Math.Ceil 和 Math.Floor 并在重整后得到两个舍入结果

33333333333333 & 33333333333334 -> /E14 -> 0.33333333333333 & 0.33333333333334

看起来不错，但是：假设我的号码接近double.MaxValue。我不能只给*E14 一个接近double.MaxValue 的值，因为这会给我一个OverflowException。所以这也不是解决办法。

而且，最重要的是：在尝试舍入 0.9999999999999999999999999 (more than 15 digits) 时，所有这些都更加失败，因为在我开始尝试向下舍入之前，内部表示已经四舍五入为 1。

我可以尝试以某种方式解析包含双精度的字符串，但这无济于事，因为 (1/3 * 3).ToString() 已经打印出 1 而不是 0.99...9。

Decimal 也不起作用，因为我不想要那么深的精度，14 位就足够了；但我仍然想要那个双倍范围！

在存在多个区间算术实现的 C++ 中，可以通过动态通知处理器将其循环模式切换为例如“始终向下”或“始终向上”来解决此问题。我在 .NET 中找不到任何方法。

那么，你有什么想法吗？提前致谢！

【问题讨论】：

我可能误解了你的问题，但你不能简单地将1/(Math.Pow(10, nDecimalPlaces)) 添加到Math.Round 结果中吗？
当结果为 0.999999999999999 时失败，因为 C# 将其四舍五入为 1。然后我会得到 [ 0.99999999999999 , 1.00000000000001 ] 这没关系，但不是那么好，因为我失去了上限的精度（1作为上限就足够了）。
我不确定interval arithmetic 是这样的。但是当你想向上或向下舍入时，为什么要四舍五入到 14 个十进制数字呢？对于像double 这样的二进制类型没有多大意义。相反，您可以四舍五入说 46 个有效位。这似乎更“统一”。 编辑：不要忘记像0.333333333333334（14 位小数）这样的数字不能用像double 这样的二进制类型精确表示。相反，如果您使用二进制“截止”，例如 46 位或其他任何值，则间隔的端点将是二进制意义上的“圆形”。
14 来自这样一个事实，即双精度最高可达 15 位有效数字。因为第 15 位可能会四舍五入，所以我想确保这个数字不算数，所以只取 14。
不确定这是否有效，但这个问题显示了如何从 C# 调用 _control87（可用于更改舍入模式）：stackoverflow.com/questions/18811466/…

标签： c# floating-point rounding intervals

【解决方案1】：

这可能更像是一个长评论而不是真正的答案。

这段代码根据截断最低七位返回一个“间隔”（我只使用Tuple<,>，你可以使用自己的Interval类型）：

static Tuple<double, double> GetMinMaxIntervalBasedOnBinaryNumbersThatAreRoundOnLastSevenBits(double number)
{
  if (double.IsInfinity(number) || double.IsNaN(number))
    return Tuple.Create(number, number); // maybe treat this case differently

  var i = BitConverter.DoubleToInt64Bits(number);

  const int numberOfBitsToClear = 7; // your seven, can change this value, must be below 52
  const long precision = 1L << numberOfBitsToClear;
  const long bitMask = ~(precision - 1L);

  //truncate i
  i &= bitMask;

  return Tuple.Create(BitConverter.Int64BitsToDouble(i), BitConverter.Int64BitsToDouble(i + precision));
}

免责声明：我不确定这是否对任何目的有用。特别是不确定它对区间算术有用。

使用此代码，GetMinMaxIntervalBasedOnBinaryNumbersThatAreRoundOnLastSevenBits(1.0 / 3.0) 返回元组 (0.333333333333329, 0.333333333333336)。

此代码，就像您在问题中要求的代码一样，具有明显的“问题”，即如果原始值接近（甚至等于）我们使用的“整数”之一，则返回间隔是“倾斜的”，原始数字接近间隔的末端之一。例如，输入42.0（已经舍入），你得到元组(42, 42.0000000000009)。

这段代码的一个好处是我希望它非常快。

【讨论】：

感谢您的努力！也玩了一段时间的 BitConverter。您的方法因大量输入而失败。例如，999999999999 返回 999999999999.016。错误 0.016 已经相当大了...
@selmaohneh 这就是我想要的！不，0.016 在这种情况下是 small 的！你在你的问题中说你不想乘以1E14，因为它可能会溢出到无穷大。这意味着您要处理大量数字。例如，您想包含号码3.14E298。我的方法可以很好地解决这个问题。您认为3.14E298 的错误是什么？

【解决方案2】：

假设nextDown(x) 是一个返回小于x 的最大双精度的函数，而nextUp(x) 是一个返回大于x 的最小双精度的函数。请参阅Get next smallest Double number 了解实施思路。

如果您将下限结果向下舍入，则使用舍入到最近结果的nextDown。如果您需要向上取整，请使用取整结果的 nextUp。

此方法确保区间继续包含精确的实数结果。它引入了额外的舍入误差 - 在某些情况下，下限将比应有的值小一 ULP，和/或上限将大一 ULP。但是，它是间隔的最小扩展，比您使用十进制或通过抑制低有效位得到的扩展要小得多。

【讨论】：