Numpy矩阵乘方得到负值

Question

我想在一个斐波那契问题中使用NumPy，因为它在矩阵乘法方面很有效率。你知道有一种方法可以使用矩阵[[1, 1], [1, 0]]来找到斐波那契数。

我写了一些非常简单的代码，但是在增加n之后，矩阵开始出现负数。

import numpy
def fib(n):
    return (numpy.matrix("1 1; 1 0")**n).item(1)

print fib(90)
# Gives -1581614984

这可能的原因是什么？

注意：`linalg.matrix_power` 也会给出负值。

注意2：我尝试了从0到100的数字。在47之后，它开始给出负值。这是否是一个大整数问题，因为NumPy是用C编码的？如果是这样，我该如何解决？

编辑：使用带有 `linalg.matrix_power` 的常规Python `list` 矩阵也给出了负结果。还要补充一点，不是所有的结果在47之后都是负数，它们是随机出现的。

编辑2：我尝试使用 @AlbertoGarcia-Raboso 建议的方法。它解决了负数问题，但出现了另外一个问题。它将答案作为 `-5.168070885485832e+19` 给出，而我需要的是 `-51680708854858323072L`。所以我尝试使用 `int()`，它将其转换为了 `L`，但现在似乎由于精度损失，答案是不正确的。

- Rockybilly

---

@mgilson 我使用了 numpy.__version__，它返回了 1.11.1。 - Rockybilly

听起来你使用的是32位版本的NumPy。fib(47)应该是2971215073，这个数字略大于有符号32位整数的最大值。 - Alex Riley

@ajcr 我使用了这个链接 https://dev59.com/t1wX5IYBdhLWcg3wpw1h 来检查我的 NumPy。它显示是 64 位的。 - Rockybilly

@Rockybilly：你可以使用NumPy来完成这个任务，但是如果你要计算大量的值，你需要使用Python内置的整数类型作为矩阵（这种类型比机器宽度整数慢）。只需使用(np.matrix("1 1; 1 0", dtype=np.object)**n).item(1)（即声明矩阵应使用“object”类型）。 - Alex Riley

@ajcr，您能否将此作为答案添加，并解释np.object的含义？ - Rockybilly

Answer 1

你看到负值出现的原因是因为NumPy默认使用np.int32数据类型来表示矩阵。

这种数据类型能够表示的最大正整数是2的31次方减1，即2147483647。不幸的是，这个数字小于第47个斐波那契数2971215073。导致溢出产生负数:

>>> np.int32(2971215073)
-1323752223

使用更大的整数类型（比如np.int64）可以解决这个问题，但只是暂时的：如果你继续请求更多、更大的斐波那契数列的话仍然会遇到问题。

唯一确保解决的方法是使用无限大小的整数类型，例如Python的int类型。要做到这一点，请将矩阵修改为np.object类型：

def fib_2(n):
    return (np.matrix("1 1; 1 0", dtype=np.object)**n).item(1)

np.object 类型允许矩阵或数组保存任何混合本地 Python 类型。本质上，矩阵不再保存机器类型，而是像 Python 列表一样，仅由指向内存中整数对象的指针组成。Python 整数将在计算斐波那契数时使用，而溢出不是一个问题。

>>> fib_2(300)
222232244629420445529739893461909967206666939096499764990979600

这种灵活性是以性能降低为代价的：NumPy的速度源于整数/浮点类型的直接存储，这些类型可以由您的硬件进行操作。

- Alex Riley