浅说佩尔方程(1)：计算连分数的渐近分数

连分数(Continued Fraction, 繁分数)是一种十分好玩的东西。一个连分数可以写成下面的表达式:

$c=a_0+\cfrac{1}{a_1+\cfrac{1}{a_2+\cfrac{1}{a_3+\cfrac{1}{a_4+\cfrac{1}{a_5+\cfrac{1}{a_6+\cfrac{1}{a_7+\ldots}}}}}}}\tag{1}$

也可以简记为:

$c=a_0+\cfrac{1}{a_1+}\cfrac{1}{a_2+}\cfrac{1}{a_3+}\cfrac{1}{a_4+}\cfrac{1}{a_5+}\cfrac{1}{a_6+}\cfrac{1}{a_7+}\ldots\tag{2}$

或者:

$c=[a_0;a_1,a_2,a_3,a_4,a_5,a_6,a_7\ldots]\tag{3}$

$(3)$ 这种记法还可以表示循环的情况：

$c=[a_0;\dot{a_1},a_2,a_3,a_4,a_5,a_6,\dot{a_7}]\tag{4}$

许多无理数都可以写成连分数的形式，比较著名的有:

$\sqrt{2} = 1+\cfrac{1}{2+\cfrac{1}{2+\cfrac{1}{2+\cfrac{1}{2+\cfrac{1}{2+\cfrac{1}{2+\cfrac{1}{2+\ldots}}}}}}}\tag{5}$

$\frac{\sqrt{5}-1}{2} = \cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\ldots}}}}}}}\tag{6}$

$e = 2+\cfrac{1}{1+\cfrac{1}{2+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{4+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{6+\ldots}}}}}}}}\tag{7}$

连分数的资料可以参见()()()()()()

为了玩一下连分数的性质，首先来实现一个计算连分数的渐近分数的算法。

考虑如下广义连分数：

$c=a_0+\cfrac{b_0}{a_1+}\cfrac{b_1}{a_2+}\cfrac{b_2}{a_3+}\cfrac{b_3}{a_4+}\cfrac{b_4}{a_5+}\cfrac{b_5}{a_6+}\cfrac{b_6}{a_7+}\ldots\cfrac{b_{n-1}}{a_n}\tag{8}$

设其第 $n$ 项的渐进分数为:

$c_n = \frac{p_n}{q_n}\tag{9}$

作一个新的数据结构来表述形如 $(8)$ 式的连分数：

data Cf = Cf [Integer] [Integer]

于是 $(8)$ 式表示为:

$Cf \quad [a_0,a_1,a_2,a_3\ldots,a_{n}] \quad [b_0,b_1,b_2,b_3\ldots,b_{n-1}]\tag{12}$

注意到

$a_0+\cfrac{b_0}{a_1+\cfrac{b_1}{a_2+\cfrac{b_2}{a_3+\cfrac{b_3}{a_4}}}} = a_0+\cfrac{b_0}{a_1+\cfrac{b_1}{a_2+\cfrac{b_2a_4}{a_3a_4+b_3}}}\tag{10}$

也即

$a_0+\cfrac{b_0}{a_1+}\cfrac{b_1}{a_2+}\cfrac{b_2}{a_3+}\ldots\cfrac{b_{n-1}}{a_n}=a_0+\cfrac{b_0}{a_1+}\cfrac{b_1}{a_2+}\cfrac{b_2}{a_3+}\ldots\cfrac{b_{n-3}}{a_{n-2}}\cfrac{a_nb_{n-2}}{a_{n}a_{n-2}+b_{n-2}}\tag{11}$

Haskell 实现:

data Cf = Cf [Integer] [Integer]
    
pq :: Cf -> (Integer, Integer)
pq (Cf xs ys) = simple $ pq' (reverse xs) (reverse ys)
    
pq' :: [Integer] -> [Integer] -> (Integer, Integer)
pq' (x0:[]) []            = (x0, 1)
pq' (x0:x1:[]) (y0:[])    = (x0*x1+y0, x0)
pq' (x0:x1:xs) (y0:y1:ys) = pq' xs' ys' where
    xs' = x0*x1+y0 : xs
    ys' = y1*x0    : ys

simple (x, y) = let d = gcd x y in (quot x d, quot y d)

其中， simple 的作用是约分一个分数。为了递归的方便，我们用 pq' 这个函数来处理倒序的输入。这里偷懒了，没有检验输入的合法性。

测试一下 $\sqrt{2}$ 的展开：

λ> pq (Cf [1] [])
(1,1)
λ> pq (Cf [1,2] [1])
(3,2)
λ> pq (Cf [1,2,2] [1,1])
(7,5)
λ> pq (Cf [1,2,2,2] [1,1,1])
(17,12)
λ> pq (Cf [1,2,2,2,2] [1,1,1,1])
(41,29)
λ> pq (Cf [1,2,2,2,2,2] [1,1,1,1,1])
(99,70)
λ> pq (Cf [1,2,2,2,2,2,2] [1,1,1,1,1,1])
(239,169)
λ> pq (Cf [1,2,2,2,2,2,2,2] [1,1,1,1,1,1,1])
(577,408)
λ> pq (Cf [1,2,2,2,2,2,2,2,2] [1,1,1,1,1,1,1,1])
(1393,985)

这样我们就可以很方便地把 $\sqrt{2}$ 的前 $N$ 项渐近分数算出来：

$c_{\sqrt{2}, 0} = \cfrac{1}{1} , c_{\sqrt{2}, 1} = \cfrac{3}{2} , c_{\sqrt{2}, 2} = \cfrac{7}{5} , c_{\sqrt{2}, 3} = \cfrac{17}{12} , c_{\sqrt{2}, 4} = \cfrac{41}{29} , c_{\sqrt{2}, 5} = \cfrac{99}{70} , c_{\sqrt{2}, 6} = \cfrac{239}{169} , c_{\sqrt{2}, 7} = \cfrac{577}{408} ,$

$\begin{align*} c_{\sqrt{2}, 8} &= \frac{1393}{985} \\\ c_{\sqrt{2}, 9} &= \frac{3363}{2378} \\\ c_{\sqrt{2}, 10} &= \cfrac{8119}{5741} \\\ c_{\sqrt{2}, 11} &= \cfrac{19601}{13860} \\\ c_{\sqrt{2}, 12} &= \cfrac{47321}{33461} \\\ c_{\sqrt{2}, 13} &= \cfrac{114243}{80782} \\\ c_{\sqrt{2}, 14} &= \cfrac{275807}{195025} \\\ c_{\sqrt{2}, 15} &= \cfrac{665857}{470832} \\\ c_{\sqrt{2}, 16} &= \cfrac{1607521}{1136689} \\\ c_{\sqrt{2}, 17} &= \cfrac{3880899}{2744210} \\\ c_{\sqrt{2}, 18} &= \cfrac{9369319}{6625109} \\\ c_{\sqrt{2}, 19} &= \cfrac{22619537}{15994428} \\\ c_{\sqrt{2}, 20} & = \cfrac{54608393}{38613965} \\\ c_{\sqrt{2}, 30} & = \cfrac{367296043199}{259717522849} \\\ c_{\sqrt{2}, 40} & = \cfrac{2470433131948081}{1746860020068409} \\\ c_{\sqrt{2}, 50} & = \cfrac{16616132878186749607}{11749380235262596085} \\\ c_{\sqrt{2}, 60} & = \cfrac{111760107268250945908601}{79026329715516201199301} \\\ c_{\sqrt{2}, 70} & = \cfrac{751698464870122983994500719}{531531081917181734003902441} \\\ c_{\sqrt{2}, 80} & = \cfrac{5055923762956339922096065927393}{3575077977948634627394046618865} \\\ c_{\sqrt{2}, 90} & = \cfrac{34006142477945877445895155433144599}{24045973948151434586670623554583549} \\\ c_{\sqrt{2}, 100} & = \cfrac{228725309250740208744750893347264645481}{161733217200188571081311986634082331709} \\\ c_{\sqrt{2}, 200} & = \cfrac{43339386297227576661095959458572614328030405537398930692163383984677691985393}{30645573943232956180057972969833245887630954508753693529117371074705767728665} \\\ c_{\sqrt{2}, 300} & = \cfrac{8212044442188195711251397481104956818556043349289906025493014758803935392763929364315400117877457284225002419343001}{5806792312476572226540901922582599644888642760862041653307913237758537888853160679392465448940731843733274415796501} \\\ \end{align*}$

简直不费吹灰之力。

Continued fraction on Wikipedia

维基百科 连分数 条目`

Continued Fraction on Wolfram Mathworld

An Introduction to Continued Fractions

Continued Fractions...an Introduction

認識連分數(臺灣大學數學系)