连分数

What · 是什么

简单连分数是把实数 x 表示成 [a_0,a_1,a_2,...]=a_0+1/(a_1+1/(a_2+...)) 的形式，其中 a_0 为整数，a_1,a_2,... 为正整数。任意有理数对应有限简单连分数，任意实数对应（可能无限）简单连分数，无理数对应无限连分数。截断前 k+1 项得到第 k 个渐近分数 x_k=p_k/q_k，是实数 x 的最佳有理逼近。二次无理数（形如 (P+√D)/Q）的连分数是循环连分数。核心用途：求最大公约数（欧几里得算法的连分数解释）、丢番图逼近、佩尔方程、Stern-Brocot 树、线性丢番图方程、以及『找与直线最近整点』类问题。

Why · 为什么

How · 实现

伪代码

# 实数 -> 连分数系数
def to_cf(x):
  a = []; 
  while x is not rational-exact:
    a.append(floor(x))
    x = 1 / (x - floor(x))

# 渐近分数递推
p_prev2, p_prev1 = 0, 1
q_prev2, q_prev1 = 1, 0
for a in coeffs:
  p = a*p_prev1 + p_prev2
  q = a*q_prev1 + q_prev2
  yield (p, q)
  p_prev2, p_prev1 = p_prev1, p
  q_prev2, q_prev1 = q_prev1, q

Python 实现

from fractions import Fraction
from math import isqrt


def frac_to_cf(x: Fraction):
    """有理数 -> 有限简单连分数系数列表."""
    a = []
    while x != 0:
        n = x.numerator // x.denominator
        a.append(n)
        x = x - n
        if x != 0:
            x = 1 / x
    return a


def cf_to_frac(a):
    """连分数系数 -> 有理数(精确)."""
    f = Fraction(a[-1])
    for k in range(len(a) - 2, -1, -1):
        f = Fraction(a[k]) + 1 / f
    return f


def convergents(coeffs):
    """返回所有渐近分数 (p_k, q_k)."""
    p_prev2, p_prev1 = 0, 1
    q_prev2, q_prev1 = 1, 0
    for a in coeffs:
        p = a * p_prev1 + p_prev2
        q = a * q_prev1 + q_prev2
        yield (p, q)
        p_prev2, p_prev1 = p_prev1, p
        q_prev2, q_prev1 = q_prev1, q


def sqrt_cf(D: int, n: int):
    """sqrt(D) 的前 n 个连分数系数 (D 非完全平方)."""
    r = isqrt(D)
    assert r * r != D
    P, Q = 0, 1
    a0 = r
    res = [a0]
    for _ in range(n - 1):
        P2 = a0 * Q - P
        Q2 = (D - P2 * P2) // Q
        a0 = (P2 + r) // Q2
        P, Q = P2, Q2
        res.append(a0)
    return res


if __name__ == "__main__":
    # 415/93 = 4 + 1/(2+1/(6+1/7)) = [4;2,6,7]
    cf = frac_to_cf(Fraction(415, 93))
    assert cf == [4, 2, 6, 7], cf
    assert cf_to_frac(cf) == Fraction(415, 93)

    # 黄金比例 phi = [1;1,1,1,...] 渐近分数 = 斐波那契比
    gold = [1] * 8
    cg = list(convergents(gold))
    assert cg[0] == (1, 1)
    assert cg[1] == (2, 1)
    assert cg[2] == (3, 2)
    assert cg[3] == (5, 3)
    assert cg[4] == (8, 5)

    # sqrt(2) = [1;2,2,2,...], 渐近分数 3/2,7/5,17/12 (佩尔解)
    c2 = sqrt_cf(2, 5)
    assert c2 == [1, 2, 2, 2, 2]
    cv = list(convergents(c2))
    assert cv[1] == (3, 2) and cv[1][0] ** 2 - 2 * cv[1][1] ** 2 == 1
    assert cv[3] == (17, 12)

    # 交叉验证: 随机有理数往返一致
    import random
    random.seed(11)
    for _ in range(200):
        x = Fraction(random.randint(1, 10**6), random.randint(1, 10**6))
        c = frac_to_cf(x)
        assert cf_to_frac(c) == x
        # 渐近分数逼近性质: |x - p/q| < 1/q^2
        for (p, q) in convergents(c):
            assert abs(x - Fraction(p, q)) < Fraction(1, q * q) or x == Fraction(p, q)
    print("All tests passed.")

易错点

• 浮点展开丢精度：必须用有理数 Fraction 或整数 PQa
• 渐近分数递推初值 p_{-2}=0,p_{-1}=1,q_{-2}=1,q_{-1}=0 记错会导致全体错位
• 误差估计 |x-p/q|<1/q^2 的常数在不同表述下不同（渐近分数用 1/(q_k q_{k+1}) 更紧）
• 奇数/偶数下标渐近分数逼近方向相反，用于『找一侧最近分数』时需选对方向
• √D 展开中 P,Q 为负值要正确归一（Q 整除性由二次无理数性质保证，勿断言整除而用整除计算）

例题

信源

• OI Wiki✓ https://oi-wiki.org/math/number-theory/continued-fraction/
• CP-Algorithms✓ https://cp-algorithms.com/algebra/continued-fractions.html