というわけで、Pythonで検算してみた

フィボナッチ数列の場合、Aは正則です。

N+1番目のフィボナッチ数列を調べるばあい、末尾再帰


F :- F( A X)
が N 回呼び出されますが、

行列 A の N 乗 と ベクトル X の積と、Fの結果が一致します。
\begin{equation}F ( \begin{bmatrix}N&0&0&1\end{bmatrix} ) ==  A^{N}  X \end{equation}

今回は素直にnumpyで行列のN乗を使ってますけれども、 スペクトル分解などを使って簡単に行列のべき乗を計算できるようにすれば高速化できるのかも。

書いてて自分でも良くわかってないのだけれども、ここから もうちょい歩いていくと、何か面白いんじゃないかという気がします。


# -*- coding: cp932 -*-
import numpy as np

A=np.array([[1,0,0,-1],
[0,0,1,0],
[0,1,1,0],
[0,0,0,1]])

assert np.linalg.det(A)==-1

def fib(x):
n,x1,x2,one=x
if n==0:
return x
elif n==1:
return x
else :
#n=n-1 x1=x2 x2=x1+x2
return fib(np.dot(A,x))

testdata=[0, 1, 1, 2, 3, 5, 8, 13,
21, 34, 55, 89, 144, 233,
377, 610, 987, 1597, 2584, 4181, 6765, 10946]

assert np.all(fib(np.array([0,0,1,1]))==np.array([0,0,1,1]))

A_=np.matrix(A)
for i in range(2,len(testdata)):
f=fib(np.array([i,0,1,1]))
assert np.all(f==np.array([1,testdata[i-1],testdata[i],1]))
x_=np.matrix([i,0,1,1]).T
#
# A の i-1 乗 と、xの行列積が、上記と同じになる
#
f_=(A_**(i-1))*x_
assert np.all(f.flatten()==f_.flatten())
print f


print "done"