Numpy matrix exponentiation gives negative value
<p>I wanted to use <code>NumPy</code> in a Fibonacci question because of its efficiency in matrix multiplication. You know that there is a method for finding Fibonacci numbers with the matrix <code>[[1, 1], [1, 0]]</code>.</p>

<p><a href="http://i.stack.imgur.com/JuBny.gif" rel="nofollow"><img src="http://i.stack.imgur.com/JuBny.gif" alt="Fibo"></a></p>

<p>I wrote some very simple code but after increasing <code>n</code>, the matrix is starting to give negative numbers.</p>

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

print fib(90)
# Gives -1581614984
</code></pre>

<p>What could be the reason for this?</p>

<p><strong>Note:</strong> <code>linalg.matrix_power</code> also gives negative values.</p>

<p><strong>Note2:</strong> I tried numbers from 0 to 100. It starts to give negative values after 47. Is it a large integer issue because NumPy is coded in C ? If so, how could I solve this ?</p>

<p><strong>Edit:</strong> Using regular python <code>list</code> matrix with <code>linalg.matrix_power</code> also gave negative results. Also let me add that not all results are negative after 47, it occurs randomly.</p>

<p><strong>Edit2:</strong> I tried using the method @AlbertoGarcia-Raboso suggested. It resolved the negative number problem, however another issues occured. It gives the answer as <code>-5.168070885485832e+19</code> where I need <code>-51680708854858323072L</code>. So I tried using <code>int()</code>, it converted it to <code>L</code>, but now it seems the answer is incorrect because of a loss in precision.</p>

<p>The reason you see negative values appearing is because NumPy has defaulted to using the <code>np.int32</code> dtype for your matrix.</p>

<p>The maximum positive integer this dtype can represent is 2<sup>31</sup>-1 which is 2147483647. Unfortunately, this is less the 47th Fibonacci number, 2971215073. The resulting overflow is causing the negative number to appear:</p>

<pre><code>&gt;&gt;&gt; np.int32(2971215073)
-1323752223
</code></pre>

<p>Using a bigger integer type (like <code>np.int64</code>) would fix this, but only temporarily: you'd still run into problems if you kept on asking for larger and larger Fibonacci numbers.</p>

<p>The only sure fix is to use an unlimited-size integer type, such as Python's <code>int</code> type. To do this, modify your matrix to be of <code>np.object</code> type:</p>

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

<p>The <code>np.object</code> type allows a matrix or array to hold any mix of native Python types. Essentially, instead of holding machine types, the matrix is now behaving like a Python list and simply consists of pointers to integer objects in memory. Python integers will be used in the calculation of the Fibonacci numbers now and overflow is not an issue.</p>

<pre><code>&gt;&gt;&gt; fib_2(300)
222232244629420445529739893461909967206666939096499764990979600
</code></pre>

<p>This flexibility comes at the cost of decreased performance: NumPy's speed originates from direct storage of integer/float types which can be manipulated by your hardware.</p>

python numpy matrix