ヘキサボナッチ数列の一般式など

よくまとめられてるなおい・・・

フィボナッチ数の類似の数列を参照すると、類推でどこまでもいけるようにできてる。


たとえばL=6ならヘキサボナッチ数列の一般式が出せる。


この6次方程式を解いた6つの複素解がa1～a6に当たるわけだけど
これを使って

こういう風に定義してやれば、たちまち、
nが0～4までのFnがゼロで、n=5になったところからポッと有限の値が出て、
それまでの数列6つの和を取る、まるで規則性のあるようなないようなヘキサボナッチ数列

を、いきなりn=40とかやっても即座に出せる。

ほんともうびっくりだよ。理論上はL=500とかでもできるし
そしたら500番目くらいまで数列がゼロで、501番目とかでいきなり有限の値が現れるんだ
そういうのが無数に作れる
狂気か！

すごーい！たーのしー！

君はトリボナッチ数列の一般式の分母のフレンズなんだね！

方程式メーカーベンチマーク＠scilab

やる気が出ないのでぼーっと、行列2重ループに頼らない方程式メーカーを考えてました。
2378次方程式の時点で、スタックが足りないとか言い出して計算できなくなりました。
ハードの問題なのか、それともコードが要改良なのか、知るか。

m=2577
x=poly(0,'x')
y=zeros(m+1,1);
for n=0:1:m
　y(n+1)=x^n;
end

z=y(1);
for n=1:1:m
　z=z+y(n+1);
end

z=z-2*x^m;

w=roots(z);
re=real(w);im=imag(w);
w2=[re,im];


出てきたw2をExcelのcomplex関数で複素数(文字型)c=a+ibにしてやって、imargument(c)でソートしてやれば
線のグラフにしても編み物にならずシンプルにグラフ化されますよっと

scilabにおけるテンソルのおまけ感

いちおう多次元配列も扱えるみたいだけど、3次元以上の配列にしたらビューアに反映されない件

まあ、2次元配列(行列)にスライスしてやりゃぁ見れるだろうけど
3Dビューアとかほしくね？うん・・・ほしくもないのか・・・。

テンソルの演算ライブラリも特に見つからない感じがする。

いちおう「scilab　テンソル」でググったけどかんばしくない



一般相対論を具体的に理解できるようになるためのツールとしてのプログラミング言語というのはないだろうか
すごく抽象的で、理解したいと思えるまでのハードルがすでに高いんですけど。

なんとかならんのですか

いや別に数値シミュレーションまで期待はしてないんで、
せめてなんかこうイメージをつかみやすくするためだけでいいんで、そういう具体化ツールないかなあ

線形代数の演習中にvitaでウルフラムα使う人

その技教えた人って絶対、「ライバルを蹴散らす」計画立ててますわ。

「行列のべき乗はこうやればいいよー」

結局固有値の意味わからず。なんてひでえやつだ。自分以外のクラスメイトほぼ全員に吹き込むのか。最悪ですわ

身に覚えがないといえばうそになる・・・

自分だけがのし上がったって誰得俺得でしかねーっつの。当時は気づかなかったのかもしれない

なんだこのscilabこええ

a1=ones(1,500)
a2=eye(499,499)
a3=zeros(499,1)
A=[a1;a2,a3]
B=A^10


500次のフィボナッチ系行列Aの10乗が
B=A^10でわりとすぐに出た・・・じゃあ対角化の意味ってなんだったんだってくらいすぐに出た

おっかねえ・・・

ぼなっち

scilabには複素行列などの構造的な変数が最初から入れられるため

for文を使う用事にあまり出くわさない気がする。


でもこれならどうだ

この解をピンポイントで求めて、nを上げてってn次方程式の解を次々表示するなんてのはどうだ！？


いまんとこscilabのグラフ機能はいまいちわかってないから、Excelにデータを移そうと思ってんだけど
その際にfor文を使う用事ができそうな気がする！



たとえば、そうだな。疑似言語で書くと

for(n++){

a1=ones(1,n)
a2=eye(n-1,n-1)
a3=zeros(n,1)
A=[a1;a2,a3]
x=poly(0,'x')
q=det(x*eye(A)-A)
p=roots(q)
re(n,:)=real(p)
im(n,:)=imag(p)
}


みたいなことをやりたい

reとimはベクトルじゃなく行列になるけど、領域の半分ぐらいが空になる感じ。

もう複素二分法のためのスツルム列とか考える意味なくなっちゃうよね・・・そこなんだよな、プログラミングのちょっと虚しいところ

せっかくこの図のソースになるExcelファイルを改良しようと思ってたのに、意欲が削がれる削がれる

scilabで計算するフィボナッチ系501次方程式

a1=ones(1,500)
a2=eye(500,500)←2重ループなので多少時間かかります
a3=zeros(500,1)
A=[a1,1;a2,a3]←2重
x=poly(0,'x')
q=det(x*eye(A)-A)←ここで一番時間食いそうです
p=roots(q)
re=real(p)
im=imag(p)

これを慣れたExcelに張り付けて、
re
im
complex(re,im)
imargument←これでソートして、re対imでグラフ化


さすがにちょっと時間かかりましたね。パソコンが

ベンチマークとしていかがでしょうｗｗｗｗｗｗ


もうね、解が、ほぼほぼ線です。吹き出しみたくなってます

1001次は終わる気しなかったので止めました。


実数だけども、非対称行列だとspecとbdiagの結果が食い違いますねえ。
specのほうがもっともらしく見えるけど、
もしかしてどっちも合ってるなんてことあるかな？？
今度、実際に対角化して確かめてみましょうか。

二分法の複素版

こないだ、高次方程式の複素解を求めてる最中にふと思ったんだけど
ニュートン法や二分法にも複素版がほしいよなぁって。


どうも、重解とかガウス平面をてさぐる範囲とかの問題で、一般的な方法論はないみたいなんだけど
スツルム列とかいうのを使って、暫定的に使える方法論とかはあるみたい。

スツルム・・・どっかで聞いたよな。スツルム・リウビル問題とか。なんだっけあれ



やっぱ効率よく精度ほしいもんな。まあ手動でできなくもないけどね


10次方程式
x^10-x^9-x^8-x^7-x^6-x^5-x^4-x^3-x^2-x-1=0

の解を求めた際に、全部足しても0.5にしかならないんだもんな。ホントは1になるはずなのに。
10個の解の積のほうがマイナス1に近かった気がする

そりゃあメッシュの刻み幅0.05固定だったら限度あるって。

あのExcelファイル、要改良だよなぁ色々と

ノナッチ行列とか

といきなり方程式の形で出されると、この高次方程式の解の法則性はよく見えないように思えますが、
この方程式が

こういう規則性のうえになり立った行列の固有値を求める特性方程式であるとわかると
得られる知識は一気に広がります。


まず、この行列のトレースは1なので、高次方程式の解の和は1になることがわかりますし

中身が全部正の実数であることから、ペロンフロベニウスの定理が活きて、ガウス平面における解の分布が右寄り横長であることもわかります。

det(デターミナントつまり行列式)も±1ということから、解を全部掛け算したものも±1であることがわかります。
行列式は、行列の大きさが偶数次だとマイナス1、奇数次だとプラス1だということもわかります。
行列式は、たとえば5次行列のペンタボナッチ数列なら、このような再帰的アルゴリズムで求めることができます。



また、ペロンフロベニウスの定理から、絶対値の一番大きな固有値は正の実数であることがわかるため、


この式のλに1や2を入れると、
解の1つが1以上2未満であることがわかるかと思います。

λ=1を代入した時

λ=2を代入した時

λが1と2のときで符号が変わっていますね。

このことから、固有値λをてさぐる範囲は、最大でも実部および虚部が-2から+2の間までであることがわかります。


それでは実際にこの多項式の零点を探してみましょう。

固有値の和も積も1あるいは±1ということは
おおよそ、ガウス平面で上下左右対称気味の配置になりそうな気がしますね。



実際、以上のように、2付近の正の実数の固有値が1つあり
1付近にはなく
あとは大方、複素平面上で半径1の円より少し小さい半径のところを等しく分割したような感じの配置になります。
半径1より少し小さいのは、2に漸近していく実数の固有値との兼ね合いだと思います。

多項式の次数が大きくなるにつれて、真円に近付いていくでしょうし、配置もほぼ円周上びっしりに漸近していくはずです。
が、ちょっと直感からズレるかもしれませんが、どれだけ次数が高くなっても、この固有値の和は1なんです。

どういうことかというと、正の実数1のところには解がなく、その代わりに正の実数2のところにあるからなのです。
単位円の円周すべての点を足し合わせたら打ち消しあってゼロになるはずです。
そこから、正の実数1を取り除くと、マイナス1になります。
どれだけ円周上に固有値を敷き詰めた円でも、正の実数1を取り除くとポッカリと空いてすべての和が-1になるのです。
そこに、正の実数1の代わりに正の実数2を付け足すと、-1+2で1になるのです。これがトレース＝1の由来になります。


DL用Excelファイル を載せておきます。


モノボナッチ　フィボナッチ　トリボナッチ　テトラボナッチ　ペンタボナッチ　ヘキサボナッチ　ヘプタボナッチ　オクタボナッチ　ノナボナッチ　デカボナッチ

モノボッチ　ジボナッチ　テトラナッチ　ノナッチ　ノナボッチ　ノナッチ





にほんブログ村

(1対)^2応

あーよかった。scilabちゃんはちゃんと複素行列の逆行列は普通に出してくれた。
ようやく、Excelより手間が省ける実感がわいてきたというか

縦ベクトルを横に並べる入れ子もうまくいった。

気を付けなきゃいけないのは、行列AのA'ってやると自動的にエルミート共役になってしまうことだね。ただの転置にしたければA.'とするらしい。(ダッシュとの間にピリオドをいれた)


で、なんとかナッチ数列の行列モデルは
有向グラフのごく狭い部分集合なんだな。
いや、こじつけ以上の意味があるかどうかはわかんないけど。

非負バイナリ行列というか非負ブーリアン行列というか
有向グラフの他のバージョンには、それぞれなにかしらの数列が1対1対応してるんだろうか


あ、でもあれだな、等比数列だったらもはやバイナリとかブーリアンじゃなくなるな。
等差数列だったらなんか非線形っぽいっつーかアフィン変換みたいなことになりそう

きのこは新しい「物理的な」高次方程式のネタを得た！

トリボナッチ以上のテトラナッチとかの数列の一般式を求めるための行列って
非負だけど対称じゃないんですよね。


そうすると、エルミート行列の性質「固有値がすべて実数」が当てはまらなくて
でもペロンフロベニウスの定理は生きるんで、大きなLボナッチ数列の概算はほぼ実数の指数関数に近似するんですけど

途中計算に複素の無理数を含みながら、数列としての結果は整数っていう奇妙な状態になりましてね
フィボナッチの無理数に飽き足らず複素数まで拡張するんですよｗｗｗ


それで、ちょっと計算させてみようと思ってscilabを立ち上げたら
さすが無償！どうもエルミート行列以外の行列に対する固有値の求め方が、数値解析でもscilabちゃんはできないみたいでして＾＾；
結局ウルフラムアルファとExcelで対処しましたｗｗｗｗｗあほかｗｗｗｗｗ
せめてガウス平面を固有値てさぐったりしてなんとかならないんですかねえ・・・まあてさぐる範囲があらかじめわかんなかったら即手詰まりな感じはありますけどね・・・？



今回はほぼ確認のためだったので、サラスの公式が使える3次の行列にとどめまして
対角化のための行列Pの逆行列も、有効数字3桁でウルフラムアルファに求めてもらいました＞＜
だって複素数なんだもんｗｗｗｗｗ

あ、でもそうですね、Excelをscilabちゃんに取りこむ方法が見つかれば
逆行列ぐらいは出してくれますよね！そのくらいできてもらわないと困る！だって信号処理のために複素行列扱うプログラミング言語ですもんね！


ウルフラムアルファの弱点は、無償バージョンではデータがコピペで移動できない点です。まーた狡猾なお手前で。＾＾


あれ？でも結局サラス使わなかったじゃん。今さっきウルフラムアルファに頼ったって言ったよね・・・？
じゃあテトラナッチでもいけたやん。手間めっちゃかかるけど。


あ、それでな
ユニタリ化できたら楽なのになーって思ってサラス使ったんでした。
でもユニタリ化失敗ｗｗｗｗｗウルフラム先生がユニタリ化して出してくれなかったら不可能なんだよたぶん俺たちの中ではなｗｗｗｗｗｗ
よく見たら固有ベクトルがなんか、線形従属っぽいんですよ。だめだこりゃー



あーアホみたいなことしたｗたのしかったー







でもこの場合、固有ベクトルには規則性があって
たとえば3次方程式だったらa1～a3までの複素根があるとすると
n番目の固有ベクトルが、以下みたいにかけるってわかったのはすごいね。

して、それを即座に見抜くウルフラムさんも化けもんですわ。

ネタ切れから逃亡！Lボナッチ数列の一般式の一般化

うおおおおお！！！ネタが切れているときに友達からネタが降ってきた！親方ぁ！

一般化したLボナッチ数列の一般式を行列で求める際には
このような行列の固有値と固有ベクトルを解けばいいらしい！

特性方程式はこうなるっぽい！

ちなみにEは単位行列、nは整数、Lは自然数

時間がない！もうこれでいいや！

せめて今日何か残しといてよかった。



4次の特殊ユニタリSU(4)生成子の最大4つの固有値は、生成子がエルミート行列だから全部実数なんです

固有値をλとして


λ^4-λ^2+p4×λ+q4=0

の、p4とq4は適した値を取りますが、常に実数解が4つ得られるんです。

そういう、最初(定義)と最後(性質)がわかってるのに、途中がよくわからない状況なんです！

フェラーリの方法に出てくる

√(u)が常に実数ってのはわかるんですよ。

中身はu=A+Bcosδって感じで、A>Bだから、u>0の実数なんですね。


ただ、λ=±√(u)±√(D)

みたいな恰好をしていて
このDがよくわからんのですよね。D>0なのは確かなんですけどねえ

うまく可視化したイメージができない。
このDの中に、±√(u)が入ってるからこそ、λは非対称になりえるってのに

わっかんねえんだよなぁ

ってか手間かかりすぎて翌日になると意欲もごっそり忘れてるから困る


そんなわけで片・両平方根方眼紙の中に円を入れてみたわけです。


悩みの種はほかにもあって
平方根方眼紙は楽なんですよ。縦横のx乗のxを0.5～1.0でなかば連続的に動かせばちゃんと機能するから。


ほかにぜひやりたいのは対数・リニア変換のニューンって動くバージョン！
できるのか！？

テイラー展開使おうって目論見までは定まってるんだけど
収束半径的なあたりが主な理由で、なんかこう漠然と不安がよぎる。

動いたらプランクの法則とかカルノーサイクル3Dとかをぬるぬる動かしてみたいよね。
もちろん回転もさせたいし、トリプル対数・リニア変換をぬるぬる動かすってのはぜひしたい

プランクの法則はなあ
4乗法則と反比例法則があるからなー、両対数方眼紙はもってこいなのよー
3Dにしがいもあるし。

っていうか前にプロトを一度作ったこともあるんだけど、いまいち伝えられてない感があったからリベンジしたい

もちろん対数・リニアのぬるぬるもそうなんだけど
てっぺん同士を比較したり、積分したものをプロットしたりするところ
あれをもうちょっと視覚に訴える何かがほしい


もしこれからやるんなら、Excelグラフはもっぱら2Dを自作3D化だね。ワイヤーフレームになるのは仕方ない。

今後プログラミングの腕があがったら、ポリゴンとかでも表現したいな。
そこはやっぱ、プログラミング言語なんだから、かゆいところに手が届いてなんぼだよね