在第四篇提到費氏數列的公式解會遇到浮點數問題,只能精準算到 122 位。因為浮點數是利用二進位的小數來做儲存,對於非二進位的數字會有誤差。而且其儲存位數有上限,對於無理數的運算,超過一定的大小就會出現精確度問題。
所以這篇就要來簡單使用 Python 內建的 Decimal module 來拉高浮點數運算的精確度。
以費氏數列公式解為例
首先回顧一下公式:
$$
F_n = \frac{1}{\sqrt{5}}\left( (\frac{1 + \sqrt{5}}{2})^n - (\frac{1 - \sqrt{5}}{2})^n \right)
$$
可以注意到裡面有無理數,也就是黃金比例的部分。對無理數做 n 次方運算,當 n 放大後很明顯誤差也會樂乘越大。為了提高精確度,我們可以提高儲存的位元。
Decimal module
Python 提供了一個內建 module:Decimal
,是一個十進位運算模型,可以避免一些浮點數誤差問題:
In [1]: 0.1 * 3 == 0.3
Out[1]: False
In [2]: from decimal import Decimal
...:
In [3]: Decimal('0.1') * 3 == Decimal('0.3')
Out[3]: True
注意我們是使用 Decimal('0.1')
而不是 Decimal(0.1)
,因為後者會把不精確的 float 0.1 丟到 Decimal,還是變得不精確:
In [4]: Decimal(0.1)
Out[4]: Decimal('0.1000000000000000055511151231257827021181583404541015625')
套用到費氏數列公式解
原本的 code:
def fib_binet(n):
golden_ratio = (1 + 5 ** 0.5) / 2
golden_ratio_alt = (1 - 5 ** 0.5) / 2
return round((golden_ratio ** n - golden_ratio_alt ** n) / 5 ** 0.5)
我們把數字用 Decimal 包起來:
def accurate_fib_binet(n):
sqrt5 = Decimal(5) ** Decimal('0.5')
golden_ratio = (1 + sqrt5) / Decimal(2)
golden_ratio_alt = (1 - sqrt5) / Decimal(2)
return round((golden_ratio ** Decimal(n) - golden_ratio_alt ** Decimal(n)) / sqrt5)
並拿矩陣解法的 fib_fast_double_iter
來比較:
x = 0
while True:
x += 1
diff = decimal_fib_binet(x) - fib_fast_double_iter(x)
if diff != 0:
print(x, diff)
break
$ python binet_de.py
123 -1
嗯?跟原本一樣?
主要是因為我們的問題最主要是儲存位數不夠,畢竟我們原本的數值就是無理數,用十進位也無法精確表示。
然而 Decimal module 提供了我們提升精確度的方法:
以空間換取精確度
from decimal import getcontext
getcontext().prec = 10000 # default is 28
把上面這段 code 加進去之後,原本的 while loop 就會跑很久…,也就代表數字到很大的時候,還是可以算準。我自己測試約可以到四萬附近。
不過精確度越大,就會算越久,可以從 MAX_PREC
這個常數看最大可以設到多少:
In [1]: from decimal import MAX_PREC
...: MAX_PREC
Out[1]: 999999999999999999
不要真的設到這麼大,不然會卡住 …
但實際上還是要看應用情境吧,對我來說如果不知道進來的數字會多大,不如直接使用公式解以外的算法比較保險。