将浮点数转换为有理数，保证转换回原始浮点数

Question

我正在寻找一种将浮点数转换为有理数的算法，以便保证有理数计算回原始浮点数，并且分母被最小化。

一个简单的算法只能将浮点数的实际值返回为 X / 2 ^N，但对于不是有限二进制分数的任何东西，2 ^N往往非常大。例如，数字 0.1，当存储在双精度浮点数中时，实际上近似为 ³??²???????¹???³????????????????? ？？？（分母为 2 ⁵⁵）。但是，将 0.1 转换为 ¹??? 显然更好，而且¹？？？将评估为 ³??²?????¹???³?????????????????????? 在浮点运算下。

一个相关的问题是用最少的数字打印十进制浮点数（本文描述了一些技术），这可以被认为是这个问题的一个特殊版本，附加约束是分母必须是 10 的幂。

有一个现有的问题，并且可能更多，但它们没有转换后的有理数必须评估为原始浮点数的约束。

Answer 1

让我们从一个定义开始，该定义精确地确定我们在任何特定情况下寻找的分数：

定义。 说，一个馏分a/b是更简单的比另一种级分c/d（用写在最低方面，具有正分母两者的级分），如果b <= d，abs(a) <= abs(c)这两个不等式的，并且至少一个是严格。

因此，例如5/7比 更简单6/7，并且比5/7更简单5/8，但两者2/5都不3/4比另一个更简单。（我们这里没有总订单。）

然后根据这个定义，有一个不立即显而易见的定理，它保证我们正在寻找的分数总是存在：

定理。给定一个J包含至少一个分数的实数子区间，J包含一个唯一的最简单 分数。换句话说，有一个独特的分数f使得

• f在J和，
• 任何其他部分g中J，f比简单g。

特别是，根据问题的要求，区间中最简单的分数总是具有最小的分母。定理中的“包含至少一个分数”条件是为了排除像闭区间这样的退化情况[?2, ?2]，它根本不包含任何分数。

我们的工作就是编写一个函数，采用有限浮点输入x，并返回最简单的分数n/d为这x是最接近浮子n/d，目标浮点格式。假设合理合理的浮点格式和舍入模式，舍入到的实数集x将形成实线的非空子区间，具有有理端点。所以我们的问题自然分解为两个子问题：

• 问题 1.给定x目标浮点格式的浮点数，描述x根据该浮点格式规则舍入到的所有值的区间。这涉及识别该区间的端点并确定该区间是开放的、封闭的还是半开放的。

• 问题 2.给定J具有有理端点的实线的非空子区间，计算该子区间中的最简单分数。

第二个问题更有趣，更少依赖平台和语言细节；让我们先解决这个问题。

寻找区间中最简单的分数

假设 IEEE 754 浮点格式和默认的舍入到偶数舍入模式，舍入到给定浮点数的间隔将是开放的或封闭的；对于其他舍入模式或格式，它可能是半开放的（一端开放，另一端封闭）。所以本节我们只看开区间和闭区间，但是适应半开区间并不难。

假设这J是具有有理端点的实线的非空子区间。为简单起见，我们假设它J是正实数线的一个子区间。如果不是，那么它要么包含0——在这种情况下0/1是最简单的分数J——或者它是负实数线的一个子区间，我们可以取反，找到最简单的分数，然后取反。

那么下面给出了一个简单的递归算法，用于在 中找到最简单的分数J：

• 如果 J 包含1，则1/1是 中最简单的分数J
• 否则，如果J是 的子区间(0, 1)，则 中的最简单分数J是1/f，其中f是 中的最简单分数1/J。（这是从“最简单”的定义中直接得出的。）
• 否则，J必须是 的子区间(1, ?)，并且 中的最简单分数J是q + f，其中是仍然在正实数内的q最大整数J - q，并且f是 中的最简单分数J - q。

对于最后一个陈述的证明草图：如果a / b是 中的最简单分数J并且c / d是 中的最简单分数J - q，则a / b比 更简单或等于(c + qd) / d，并且c / d更简单于或等于(a - qb) / b。所以b <= d, a <= c + qd, d <= band c <= a - qb, 并且它跟随着b = dand a = c + qd, so c / d = a / b - q。

在类似 Python 的伪代码中：

def simplest_in_interval(J):
    # Simplest fraction in a subinterval J of the positive reals
    if J < 1:
        return 1 / simplest_in_interval(1/J)
    else if 1 < J:
        q = largest_integer_to_the_left_of(J)
        return q + simplest_in_interval(J - q)
    else:
        # If we get here then J contains 1.
        return 1

要看到算法必须始终终止，不能进入无限循环，请注意，每个求逆步骤都必须后跟一个J - q步骤，并且每个J - q步骤都将间隔的左右端点的分子减少。具体来说，如果区间的端点是a/b和c/d，则总和abs(a) + abs(c) + b + d是一个正整数，随着算法的进行稳步减少。

要将上述内容转换为真正的 Python 代码，我们必须处理一些细节。首先，让我们暂时假设这J是一个闭区间；我们将适应下面的开放间隔。

我们将用它的端点left和 来表示我们的区间right，这两个都是正fraction.Fraction实例。然后下面的Python代码实现了上面的算法。

from fractions import Fraction
from math import ceil

def simplest_in_closed_interval(left, right):
    """
    Simplest fraction in [left, right], assuming 0 < left <= right < ?.
    """
    if right < 1:  # J ? (0, 1)
        return 1 / simplest_in_closed_interval(1 / right, 1 / left)
    elif 1 < left:  # J ? (1, ?):
        q = ceil(left) - 1  # largest q satisfying q < left
        return q + simplest_in_closed_interval(left - q, right - q)
    else:  #  left <= 1 <= right, so 1 ? J
        return Fraction(1)

这是一个示例运行：

>>> simplest_in_closed_interval(Fraction("3.14"), Fraction("3.15"))
Fraction(22, 7)

原则上，开区间的代码同样简单，但实际上有一个复杂的问题：我们可能需要处理无限区间。例如，如果我们的原始间隔是J = (2, 5/2)，那么第一步将该间隔移动2获得(0, 1/2)，然后将该间隔反转为(2, ?)。

所以对于开放区间，我们将继续用它的一对(left, right)端点来表示我们的区间，但现在right要么是一个fractions.Fraction实例，要么是一个特殊的常量INFINITY。而不是简单地能够用来1 / left取左端点的倒数，我们需要一个辅助函数来计算分数或 的倒数INFINITY，以及另一个用于减法的辅助函数，确保INFINITY - q给出INFINITY。以下是这些辅助功能：

#: Constant used to represent an unbounded interval
INFINITY = "infinity"

def reciprocal(f):
    """ 1 / f, for f either a nonnegative fraction or ? """
    if f == INFINITY:
        return 0
    elif f == 0:
        return INFINITY
    else:
        return 1 / f


def shift(f, q):
    """ f - q, for f either a nonnegative fraction or ? """
    if f == INFINITY:
        return INFINITY
    else:
        return f - q

这是主要功能。注意if和elif条件中不等式的变化，以及我们现在想要使用floor(left)而不是ceil(left) - 1找到q位于区间左侧的最大整数的事实：

from fractions import Fraction
from math import floor

def simplest_in_open_interval(left, right):
    """
    Simplest fraction in (left, right), assuming 0 <= left < right <= ?.
    """
    if 1 <= left:  # J ? (1, ?):
        q = floor(left)
        return q + simplest_in_open_interval(shift(left, q), shift(right, q))
    elif right != INFINITY and right <= 1:  # J ? (0, 1)
        return 1 / simplest_in_open_interval(reciprocal(right), reciprocal(left))
    else:  #  left < 1 < right, so 1 ? J
        return Fraction(1)

上面的代码是为了清晰而不是效率而优化的：它在 big-Oh 复杂性方面相当有效，但在实现细节方面却没有。我把它留给读者来转换成更有效的东西。第一步是使用整数分子和分母而不是fractions.Fraction实例。如果您对它的外观感兴趣，请查看我在 PyPI 上的simplefractions包中的实现。

查找舍入到给定浮点数的间隔

现在我们可以找到给定区间中的最简单分数，我们需要解决问题的另一半：找到四舍五入到给定浮点数的区间。这样做的细节将在很大程度上取决于语言、使用的浮点格式，甚至是使用的舍入模式等。

在这里，我们概述了在 Python 中执行此操作的一种方法，假设 IEEE 754 binary64 浮点格式具有默认的舍入到偶数舍入模式。

为简单起见，假设我们的输入浮点数x是正的（并且是有限的）。

Python >= 3.9 提供了一个函数math.nextafter，允许我们从x. 这是对最近的浮点数 ? 执行此操作的示例：

>>> import math
>>> x = 3.141592653589793
>>> x_plus = math.nextafter(x, math.inf)
>>> x_minus = math.nextafter(x, 0.0)
>>> x_plus, x_minus
(3.1415926535897936, 3.1415926535897927)

（请注意，一般来说，要做到这一点，我们还需要处理x最大可表示浮点数并math.nextafter(x, math.inf)给出无穷大的特殊情况。）

舍入到的间隔的边界x是x和相邻浮动之间的中间。Python 允许我们将浮点数转换为相应的精确值作为分数：

>>> from fractions import Fraction
>>> left = (Fraction(x) + Fraction(x_minus)) / 2
>>> right = (Fraction(x) + Fraction(x_plus)) / 2
>>> print(left, right)
14148475504056879/4503599627370496 14148475504056881/4503599627370496

我们还需要知道我们有一个闭区间还是一个开区间。我们可以查看位表示来弄清楚（这取决于浮点数的最低有效位是0还是1），或者我们可以测试一下我们的区间端点是否四舍五入x：

>>> float(left) == x
True
>>> float(right) == x
True

他们这样做，所以我们有一个闭区间。通过查看浮点数的十六进制表示可以证实这一点：

>>> x.hex()
'0x1.921fb54442d18p+1'

所以我们可以找到四舍五入为xusing的最简单分数simplest_in_closed_interval：

>>> simplest_in_closed_interval(left, right)
Fraction(245850922, 78256779)
>>> 245850922 / 78256779 == x
True

把这一切放在一起

虽然核心算法很简单，但有足够多的极端情况需要处理（负值、开区间与闭区间sys.float_info.max等），以至于一个完整的解决方案最终有点过于混乱，无法在此答案中完整发布。前段时间我整理了一个 Python 包，称为simplefractions处理所有这些极端情况；它在 PyPI上可用。这是在行动：

>>> from simplefractions import simplest_from_float
>>> simplest_from_float(0.1)
Fraction(1, 10)
>>> simplest_from_float(-3.3333333333333333)
Fraction(-10, 3)
>>> simplest_from_float(22/7)
Fraction(22, 7)
>>> import math
>>> simplest_from_float(math.pi)
Fraction(245850922, 78256779)