生成斐波那契数列时出现了OverflowError 'Numerical result out of range'错误

Question

可能是重复问题：
如何在Python中处理非常大的数字

我有一个用于生成斐波那契数列的Python函数：

def fib(n):                                                                                                            
        return ((1+math.sqrt(5))**n - (1-math.sqrt(5))**n)/(2**n*math.sqrt(5))

我可以输入由 fib 函数计算的数字，直到 700，此时函数的性能开始下降。

OverflowError: (34, 'Numerical result out of range')

我需要使用一种特殊类型，例如long来解决这个问题吗？

- fergusdawson

Answer 1

问题在于您使用double类型计算值时会溢出。 double类型只能准确地计算到大约第85个斐波那契数。 如果您想要快速且准确地计算，最好使用基于更好递推关系的算法，并使用Python bignum整数。 特别地，您可以使用以下方法:

 fib(2*n) = fib(n)^2 + fib(n-1)^2
 fib(2*n-1) = fib(n)*(2*fib(n-1)+fib(n))

或者等价的矩阵指数公式（请原谅丑陋的格式）

 [ F_n     F_{n-1} ]      [ 1   1 ] ^N 
 [                 ]  =   [       ]
 [ F_{n-1} F_{n-2} ]      [ 1   0 ]

这两种方法都可以得到算法，其计算次数为O(log(N))而非O(N)。

以下是伪代码完整的解决方案

---

如果您想使用双精度和显式公式进行计算，则可以调整公式以获得更快的结果，直到大约第1500个斐波那契数之前不会溢出，并保持与您的版本相同的准确性。我IRC 是：

def fib(n):                                                                                                            
    return round( ((1+math.sqrt(5))/2)**n / math.sqrt(5) )

- Michael Anderson

---

应该使用“round”而不是“floor”。被压制的术语“((1-sqrt(5))/2)**n”具有交替符号。 - Daniel Fischer

谢谢你，Daniel。你是正确的。现在已经修复了。 - Michael Anderson

Answer 2

很容易分离错误

>>> (1+math.sqrt(5))**700
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
OverflowError: (34, 'Numerical result out of range')

这种方法不太可行，因为浮点数没有足够的精度。 例如，在这里：

>>> (1+math.sqrt(5))**600
1.024664165563927e+306

你正在使用的是前15位数字，其余的291位在进行任何算术运算时都将被视为零。

有关浮点数精度问题的更多信息，请参见wikipedia 链接。

- John La Rooy

Answer 3

您可以尝试以下方法：

您始终可以尝试这种方法：

def fib(n, memo={0:0, 1:1}):
    if n not in memo:
        memo[n] = fib(n-1) + fib(n-2)
    return memo[n]

print fib(800)

输出：

69283081864224717136290077681328518273399124385204820718966040597691435587278383112277161967532530675374170857404743017623467220361778016172106855838975759985190398725

- arshajii

---

更好的做法是将 fib(n, memo={0:0, 1:1}) 作为可变默认参数传入。在函数定义时，memo 将被创建，但仍保留在本地作用域中。 - wim

+1. 但是... print fib(1000) 将会引发一个 RuntimeError: maximum recursion depth exceeded（至少在64位Mac CPython 2.7.2上会这样，不同的情况可能有所不同）。因此，这并没有比他的原始代码更进一步。当然，除非你只是通过这种方式来避免 RuntimeError: fib(500); print fib(1000). - abarnert

@wim：我不确定我喜欢那个。你知道，“显式比隐式更好”。 - georg

这个用例有什么不明确的地方吗？当所讨论的函数不太适合作为类对象时，使用它比使用全局变量更加优美。而且标准库中也有这种模式的例子，如果我没记错的话。 - wim

@thg435：Guido 承认 3.x 中的一个缺陷是，“默认变量黑科技”仍然是 Pythonic 的处理方式，尽管它本不应如此。语言中没有更好的解决方案，试图修改语言以添加更好的解决方案似乎都会引起比解决的问题更糟糕的问题。 - abarnert

Answer 4

如果您确实想使用该算法，并且希望超越内置的float限制，那么是的，您需要不同的类型。

如果您只想获得近似答案而不是异常，那很简单；您可以直接使用无穷范围。但是，如果您还想消除舍入误差，则无法具有无限精度（这将花费无限时间/空间），因此您必须知道如何计算所需精度以用于输入范围。（我会把这留给读者作为练习。）

标准库类型decimal.Decimal可能是您所需的全部内容。它根据IEEE-854标准提供任意精度的定点或浮点十进制算术。虽然它不提供足够的数学函数而对许多情况不可用，但您只需要基本算术和sqrt，它就能够胜任。 对于大型数字，它也可能很慢，但是如果您只想在几个三位数上计算fib，则足够了。

当 Decimal 不足时，通常会使用一些第三方模块来包装行业标准的 C 库，例如 gmp/mpfr，比如 bigfloat。

以下是如何获得无限范围，但舍入误差大致与内置浮点数相同的方法：

>>> s5 = decimal.Decimal(5).sqrt()
>>> def fib(n):
...     return ((1+s5)**n - (1-s5)**n)/(2**n*s5)
>>> fib(800)
Decimal('6.928308186422471713629008226E+166')
>>> int(fib(800))
69283081864224717136290082260000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000L
>>> s5 = bigfloat.sqrt(5)
>>> def fib(n):
...     return ((1+s5)**n - (1-s5)**n)/(2**n*s5)
>>> fib(800)
BigFloat.exact('6.9283081864226567e+166', precision=53)
>>> int(fib(800))
69283081864226566841137772774650010139572747244991592044952506898599601083170460360533811597710072779197410943266632999194601974766803264653830633103719677469311107072L

但请注意，这两个答案都不是您在完美计算时得到的答案；您已经因为舍入误差而丢失了24位数字。（这些值不同的原因是bigfloat在二进制中进行舍入，而decimal则在十进制中进行舍入。）

为了解决这个问题，您需要更多的精度。所有库都提供一些改变精度的方法；bigfloat具有比大多数库更方便的选项，但没有太多负担：

>>> decimal.getcontext().prec = 300
>>> s5 = decimal.Decimal(5).sqrt()
>>> def fib(n):
...     return ((1+s5)**n - (1-s5)**n)/(2**n*s5)
>>> fib(800)
69283081864224717136290077681328518273399124385204820718966040597691435587278383112277161967532530675374170857404743017623467220361778016172106855838975759985190398725.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000048
>>> def fibp(n, p):
...     with bigfloat.precision(p):
...         s5 = bigfloat.sqrt(5)
...         return ((1+s5)**n - (1-s5)**n)/(2**n*s5)
>>> fibp(800, 125)
BigFloat.exact('6.92830818642247171362900776814484912138e+166', precision=125)
>>> int(fibp(800, 125))
69283081864224717136290077681448491213794574774712670552070914552025662674717073354503451578576268674564384721027806323979200718479461097490537109958812524476157132800L

- abarnert

---

没有内置的方法：Decimal呢？ - DSM

Decimal 不像 bigfloat 一样传播错误（这对计算机工程师更有意义，对科学家来说则不太有意义），而且它要慢得多，并且提供的数学函数远远不如后者。然而，它确实提供了 sqrt，这就是 OP 所需的全部内容，如果他还没有理解误差传播，那么他不会对 IEEE 模型感到失望，所以唯一的真正问题是速度。我会更新答案。 - abarnert