トリボナッチ数列～複素無理数が正の整数に化ける不思議～

トリボナッチ数列というのは以下のように

n番目の数とn+1番目の数を足してn+2番目の数列の値にするというもので
式にするとこのような漸化式となるのですが



これを行列で表現すると、このようになります。


この漸化式を、一般式に直すのは、行列では簡単です。

F1=1,F0=0とすると

こうすればいいだけです。

しかし、この式の中の行列のべき乗がややこしくなってくるため、
固有値・固有ベクトル・対角化の技術が必要になってくるのです。

この行列はフィボナッチ数列とは違って、対称ではないため、固有値は複素数になりえます。
ただし、非負なので、ペロンフロベニウスの定理には従います。


このような3次方程式を解くことになります。


フリーのプログラミング言語、scilabで解いてみましょう。

scilab上で、
A=[1,1,1;1,0,0;0,1,0]
p=spec(A)

と書くと、解いてくれます。

ただし、scilabはあくまでフリーの言語なので
Tを対角化された行列
Dを対角化するための行列として
[D,T]=bdiag(A)

というコマンドを指定しても、非対称(非エルミート)行列では望み通りの対角化をしてくれないかもしれません
(ただ、T*D*inv(T)を計算させたら矛盾なくAになったので、道理は通っているみたいです)


ここでは、spec(A)で得られた固有値を用いることにします。

p=spec(A)
というコマンドで、
λ1=1.8392868、λ2,3=− 0.4196434 ± 0.6062907i

という固有値を得ました。
p=
このように、列ベクトルの形で入ってます。

マセマティカによると、この3つの固有値で対角化された行列
J=
を算出する、対角化のための行列は
P=
だそうなので、

scilabに

P=[p(1)^2,p(2)^2,p(3)^2;p(1),p(2),p(3);1,1,1]
Q=inv(P)
J=clean(Q*A*P)

というコマンドを打ちこむと、対角化された行列Jを算出できます。
クリーンってコマンドで、誤差範囲の純虚数成分をカットしてくれます。

そして、行列Aのべき乗は

P*J^n*Qで算出可能なので

forループを使って、以下のように算出してみました。
Cという列ベクトルに、逐次トリボナッチ数列を吐き出してます。
forでnを増やすたびに、これまでのCの値が逐一出てしまうので、コマンドの後ろにセミコロンを入力して、計算結果の表示を省略してます。
最後の行でCとだけ書いているのは、「列ベクトルCの値を表示せよ」という意味です。

roundは、実数の範囲で出てきてしまった誤差を取り除くためのものです。
床関数floorじゃダメです。0.99999とかいうのが0になってしまいます。
型変換int32もダメです。Cという入れ物が一度整数型になってしまうと、動作テストのために「やっぱりint32なしで計算させてみたい」ってなっても「型の不一致」というエラーが出ます。融通が利くように、roundにしたほうがいいようです。

for n=1:31
B=round(clean(P*J^n*Q));
C(n)=B(1,1);
end
C

算出してみると

1.
2.
4.
7.
13.
24.
44.
81.
149.
274.
504.
927.
1705.
3136.
5768.
10609.
19513.
35890.
66012.
121415.
223317.
410744.
755476.
1389537.
2555757.
4700770.
8646064.
15902591.
29249425.
53798080.
98950096.


0,1に続く3番目からの数列ですが、ちゃんと出ているのがわかるかと思います。
1e+10とかって指数表示になるような大きな数になる手前までにしておきました。

少々無理強いしましたが、しっかりと複素無理数が打ち消されて、整数に化けています。


そう考えると、デカボナッチ数列を実装する際には、最初9つの初期値が必要で
もはや数列の体をなしていないともいえそうです。
実際、フィボナッチの初期値を変えたものとして、リュカ数が有名？ですよね


じゃあ501次方程式なんか500個の初期値を持つわけだから、いったい数列として何の意味を成すんだと思っちゃいます。

が、これは、数列と行列と高次方程式がつながっているからこそ
意味があるのではないかと思いました。


scilabは最初見たときはコマンドプロンプトのように見えました。
しかし、scipad(今のscinotes)があることを知って、ちゃんと閉じた格好して部品として保存できることがわかりました。
また、scinotesではないscilab本体にfor文を入れても、ちゃんとendと打つまで演算を待ってくれます。

ただ、新たな疑問がわいてきました。
こいつは実行ファイルを作ることができるのかと。

たとえばC言語だったら実行ファイルを作ることができて
コード次第ではいろんな拡張子のファイルを作ってそこに有用なデータを入れることができます。
bmp仕様のbmpファイルを作る、実行ファイルも作れることでしょう。

しかし、このscilabはそんなことをするように設計されてはいないのかもしれません