n 的第 N 个斐波那契数等于 10^19?
Nth Fibonacci number for n as big as 10^19?
我正在尝试编写一个程序来查找 1 < n < 10^19 的第 n 个斐波那契数。
这是我使用动态规划的代码。
memo = {}
def fib(n):
if n in memo:
return memo[n]
if n <= 2:
f = 1
else:
f = fib(n-1) + fib(n-2)
memo[n]=f
return f
print fib(input()) % 1000000007
我的代码似乎不适用于大数。我收到无效响应错误。
有什么建议吗?
Python 的默认值 recursion limit 为 1000(通常)。要了解您系统的确切限制:
>>> import sys
>>> sys.getrecursionlimit()
首先,如果你想递归地写这个并且你正在使用 Python 3.2 及更高版本(它看起来不像你,从 print
语句判断)那么你可以像这样使用 @functools.lru_cache(maxsize=128, typed=False)
:
import functools
@functools.lru_cache()
def fib(n):
if n <= 2:
return 1
else:
return fib(n-1) + fib(n-2)
话虽如此,这对于大量数据来说仍然不会很快。更好的方法是编写一个迭代解决方案,在任何给定时间,您只需要"memoize"就是最后两个数字。
您当然可以使用 matrix form 以获得更好的性能。
最终,因为 n
和 10**19
一样大,你将很难编写任何在 Python 中运行的东西而不给你一个 OverflowError
.
以 O(n) 的效率,您将永远无法到达那里。与代码无关,但 Dijkstra's note "In honor of Fibonacci" 描述了一种在 O(log(n)) 效率中找到 F(n) 的方法。
F(2n-1) = F(n-1)^2 + F(n)^2
F(2n) = (2*F(n-1)+F(n))*F(n)
你不仅可以这样做,而且还可以递归地做。
当 N 为 10^19 时获得第 N 个斐波那契数如果你以天真的方式来做是行不通的(至少我猜它不会行得通)。
有一种多更好的方法来做到这一点。这种技术适用于很多像这样的系列。它被称为 Fibonacci Q Matrix.
哪里
可以这样想:
你有一些矩阵将向量 A 转换为 B:
填写这些条目很容易。特殊之处在于它现在是一个矩阵运算符,所以如果我们想要第 1000 个斐波那契数,我们只需要进行矩阵乘法即可。
你可以用一个循环来做到这一点,但它会花费你相当长的时间才能达到 10^19,并且做 10^19 矩阵乘法(即使它们很小)正在进行也要花点时间。
相反,我们走另一条捷径。 x^N 可以重写为幂的乘积,它们总和为 N,即
x**100 == x**90 * x**10
所以目标是在不进行大量计算的情况下在索引中获得大量数字:
x**2
与 x*x
一样困难 - 他们花费的时间相同。但是 x*x*x*x
给出与 (x**2)**2
相同的答案,同时需要额外的乘法。当您获得更高的权力时,收益会增加。因此,如果您将指数分解为 2 的幂(任何幂都有效,但这是最简单的情况),
X**100 == X**64 * X**32 * X**4
即
X**100 == (((((X**2)**2)**2)**2)**2)**2 + ...
所以你要做的是计算你想要达到的总功率的两个幂,然后取 Q
矩阵的两个幂的乘积。
这似乎对我有用:
fib_matrix = [[1,1],
[1,0]]
def matrix_square(A, mod):
return mat_mult(A,A,mod)
def mat_mult(A,B, mod):
if mod is not None:
return [[(A[0][0]*B[0][0] + A[0][1]*B[1][0])%mod, (A[0][0]*B[0][1] + A[0][1]*B[1][1])%mod],
[(A[1][0]*B[0][0] + A[1][1]*B[1][0])%mod, (A[1][0]*B[0][1] + A[1][1]*B[1][1])%mod]]
def matrix_pow(M, power, mod):
#Special definition for power=0:
if power <= 0:
return M
powers = list(reversed([True if i=="1" else False for i in bin(power)[2:]])) #Order is 1,2,4,8,16,...
matrices = [None for _ in powers]
matrices[0] = M
for i in range(1,len(powers)):
matrices[i] = matrix_square(matrices[i-1], mod)
result = None
for matrix, power in zip(matrices, powers):
if power:
if result is None:
result = matrix
else:
result = mat_mult(result, matrix, mod)
return result
print matrix_pow(fib_matrix, 10**19, 1000000007)[0][1]
然后,你可以更进一步 - 它只是一个 2x2 矩阵,所以我们可以对角化它,然后得到第 n 个斐波那契数的公式,就像 n 的函数 - 没有递归.像这样:
如上所述,我们计算使我们从一步到下一步的矩阵:
然后是从一组数字到下一组数字的关系:
我们可以在哪里链接这些矩阵乘法:
没有什么能阻止我们一路回到第一个斐波那契数列:
现在游戏变成了 "how do we raise that matrix to the power n" - 这正是上面代码中所做的。但是有比我上面提出的解决方案更好的方法。我们可以将 Q-matrix 分解为特征值和向量,这样写:
其中 U 是包含 Q[= 的特征值的酉矩阵106=],而Λ是对应的特征值矩阵。这些特征值和向量是:
然后你使用这种分解方式的标准优势之一,当你将它提升到一个幂时,相邻的 U 矩阵和它的逆结合给出酉矩阵,留下一个 U 和它在末端是逆矩阵,中间有一串对角矩阵,将它们提升到一个幂是微不足道的:
所以现在我们拥有了根据单个公式编写第 n 个斐波那契数所需的一切,无需递归。我会在本周晚些时候 tomorrow/some 完成它...
我不认为你可以用这个达到 1E19,但这是避免双重溢出和递归深度限制的方法:
import decimal
import operator
def decimal_range(start, stop, step=1):
"""Provides an alternative to `xrange` for very high numbers."""
proceed = operator.lt
while proceed(start, stop):
yield start
start += step
def fib(n):
"""
Computes Fibonacci numbers using decimal.Decimal for high
precision and without recursion
"""
a, b = decimal.Decimal(0), decimal.Decimal(1)
for i in decimal_range(0, n):
a, b = b, a + b
return a
在我的机器上,计算 1E6 需要 26.5 秒,但我不能保证结果的正确性:
In [26]: %time f2(n)
CPU times: user 26.4 s, sys: 130 ms, total: 26.5 s
Wall time: 26.5 s
Out[26]: Decimal('1.953282128707757731632014830E+208987')
迭代器取自this SO thread with minimal alterations, while the fib
function can be found in this other thread。
我正在尝试编写一个程序来查找 1 < n < 10^19 的第 n 个斐波那契数。
这是我使用动态规划的代码。
memo = {}
def fib(n):
if n in memo:
return memo[n]
if n <= 2:
f = 1
else:
f = fib(n-1) + fib(n-2)
memo[n]=f
return f
print fib(input()) % 1000000007
我的代码似乎不适用于大数。我收到无效响应错误。 有什么建议吗?
Python 的默认值 recursion limit 为 1000(通常)。要了解您系统的确切限制:
>>> import sys
>>> sys.getrecursionlimit()
首先,如果你想递归地写这个并且你正在使用 Python 3.2 及更高版本(它看起来不像你,从 print
语句判断)那么你可以像这样使用 @functools.lru_cache(maxsize=128, typed=False)
:
import functools
@functools.lru_cache()
def fib(n):
if n <= 2:
return 1
else:
return fib(n-1) + fib(n-2)
话虽如此,这对于大量数据来说仍然不会很快。更好的方法是编写一个迭代解决方案,在任何给定时间,您只需要"memoize"就是最后两个数字。
您当然可以使用 matrix form 以获得更好的性能。
最终,因为 n
和 10**19
一样大,你将很难编写任何在 Python 中运行的东西而不给你一个 OverflowError
.
以 O(n) 的效率,您将永远无法到达那里。与代码无关,但 Dijkstra's note "In honor of Fibonacci" 描述了一种在 O(log(n)) 效率中找到 F(n) 的方法。
F(2n-1) = F(n-1)^2 + F(n)^2
F(2n) = (2*F(n-1)+F(n))*F(n)
你不仅可以这样做,而且还可以递归地做。
当 N 为 10^19 时获得第 N 个斐波那契数如果你以天真的方式来做是行不通的(至少我猜它不会行得通)。
有一种多更好的方法来做到这一点。这种技术适用于很多像这样的系列。它被称为 Fibonacci Q Matrix.
哪里
可以这样想:
你有一些矩阵将向量 A 转换为 B:
填写这些条目很容易。特殊之处在于它现在是一个矩阵运算符,所以如果我们想要第 1000 个斐波那契数,我们只需要进行矩阵乘法即可。
你可以用一个循环来做到这一点,但它会花费你相当长的时间才能达到 10^19,并且做 10^19 矩阵乘法(即使它们很小)正在进行也要花点时间。
相反,我们走另一条捷径。 x^N 可以重写为幂的乘积,它们总和为 N,即
x**100 == x**90 * x**10
所以目标是在不进行大量计算的情况下在索引中获得大量数字:
x**2
与 x*x
一样困难 - 他们花费的时间相同。但是 x*x*x*x
给出与 (x**2)**2
相同的答案,同时需要额外的乘法。当您获得更高的权力时,收益会增加。因此,如果您将指数分解为 2 的幂(任何幂都有效,但这是最简单的情况),
X**100 == X**64 * X**32 * X**4
即
X**100 == (((((X**2)**2)**2)**2)**2)**2 + ...
所以你要做的是计算你想要达到的总功率的两个幂,然后取 Q
矩阵的两个幂的乘积。
这似乎对我有用:
fib_matrix = [[1,1],
[1,0]]
def matrix_square(A, mod):
return mat_mult(A,A,mod)
def mat_mult(A,B, mod):
if mod is not None:
return [[(A[0][0]*B[0][0] + A[0][1]*B[1][0])%mod, (A[0][0]*B[0][1] + A[0][1]*B[1][1])%mod],
[(A[1][0]*B[0][0] + A[1][1]*B[1][0])%mod, (A[1][0]*B[0][1] + A[1][1]*B[1][1])%mod]]
def matrix_pow(M, power, mod):
#Special definition for power=0:
if power <= 0:
return M
powers = list(reversed([True if i=="1" else False for i in bin(power)[2:]])) #Order is 1,2,4,8,16,...
matrices = [None for _ in powers]
matrices[0] = M
for i in range(1,len(powers)):
matrices[i] = matrix_square(matrices[i-1], mod)
result = None
for matrix, power in zip(matrices, powers):
if power:
if result is None:
result = matrix
else:
result = mat_mult(result, matrix, mod)
return result
print matrix_pow(fib_matrix, 10**19, 1000000007)[0][1]
然后,你可以更进一步 - 它只是一个 2x2 矩阵,所以我们可以对角化它,然后得到第 n 个斐波那契数的公式,就像 n 的函数 - 没有递归.像这样:
如上所述,我们计算使我们从一步到下一步的矩阵:
然后是从一组数字到下一组数字的关系:
我们可以在哪里链接这些矩阵乘法:
没有什么能阻止我们一路回到第一个斐波那契数列:
现在游戏变成了 "how do we raise that matrix to the power n" - 这正是上面代码中所做的。但是有比我上面提出的解决方案更好的方法。我们可以将 Q-matrix 分解为特征值和向量,这样写:
其中 U 是包含 Q[= 的特征值的酉矩阵106=],而Λ是对应的特征值矩阵。这些特征值和向量是:
然后你使用这种分解方式的标准优势之一,当你将它提升到一个幂时,相邻的 U 矩阵和它的逆结合给出酉矩阵,留下一个 U 和它在末端是逆矩阵,中间有一串对角矩阵,将它们提升到一个幂是微不足道的:
所以现在我们拥有了根据单个公式编写第 n 个斐波那契数所需的一切,无需递归。我会在本周晚些时候 tomorrow/some 完成它...
我不认为你可以用这个达到 1E19,但这是避免双重溢出和递归深度限制的方法:
import decimal
import operator
def decimal_range(start, stop, step=1):
"""Provides an alternative to `xrange` for very high numbers."""
proceed = operator.lt
while proceed(start, stop):
yield start
start += step
def fib(n):
"""
Computes Fibonacci numbers using decimal.Decimal for high
precision and without recursion
"""
a, b = decimal.Decimal(0), decimal.Decimal(1)
for i in decimal_range(0, n):
a, b = b, a + b
return a
在我的机器上,计算 1E6 需要 26.5 秒,但我不能保证结果的正确性:
In [26]: %time f2(n)
CPU times: user 26.4 s, sys: 130 ms, total: 26.5 s
Wall time: 26.5 s
Out[26]: Decimal('1.953282128707757731632014830E+208987')
迭代器取自this SO thread with minimal alterations, while the fib
function can be found in this other thread。