数学の具体的な計算にMaximaを使って、数学もMaximaも同時に学んでしまいましょう。今回はMaximaを使って行列の計算をしてみたいと思います。行列の固有値・固有ベクトルを求め、対角化する方法をみてみます。
Maximaのごく基本的な使い方については以下の記事を参照してください:
また、Maximaでの行列の扱いの基本については
にまとめていますので、参考にしてみてください。
固有値・固有ベクトル
をみたす数 とベクトル のことをいいます。
非自明な固有ベクトル が存在するためには
としたとき、行列
の逆元が存在してはいけません。なぜなら、もし逆元が存在すれば、それを両辺に作用させて が導かれてしまい、 が非自明であることに反してしまうからです。
以上から、固有値に対して以下の固有方程式
が要請されます。この固有方程式を について解くことにより、行列 の固有値が求まります。
一般に、 次正方行列に対しては固有方程式は の 次方程式になるので、複素数の範囲では(縮重度を含めて) 個の固有値があります。
行列の対角化
以下では固有値は互いに全て異なるものとします。
を導入すると、この行列 は正則であることがわかります。このことは、互いに異なる固有値の固有ベクトル が線型独立であることから結論できます。
一番簡単な場合として、2次正方行列の場合を考えてみましょう。このとき一般に2つの固有値と固有ベクトルがあります:
ここで2つの固有値が互いに異なるとき、 と は線型独立です。
これを背理法で示してみましょう。すなわち、2つの固有値が互いに異なるときに と が線型独立でない、と仮定してみます。すると、ある数 が存在して
と書けることになります。仮定から、 と は行列 の固有ベクトルでしたから
です。以上から
となります。いま、固有ベクトル は非自明ですから、この式が成り立つためには
が導かれます。ところがこれは2つの固有値が互いに異なるという前提と矛盾します。よって と は線型独立であることがわかります。
同様の考察によって一般の場合についても、互いに異なる固有値の固有ベクトル が線型独立であることがわかります。したがって線型独立なベクトルから構成された行列 は、その行列式がゼロではあり得ず、その逆行列 が存在すること、つまり正則であることがわかります。
さて、 が行列 の固有ベクトルであることから
より、これに左から行列 を作用させて
という対角行列を作ることができます。このような対角行列を作ることのご利益ですが、行列 のベキ乗が
によって簡単に計算できるということがあげられます。
具体例:2次正方行列の場合
を考えてみます。行列をMaximaに入力するには
A:matrix([1,2],[3,2]);
とします。
固有方程式 は行列式
D:determinant(A-lambda*ident(2));
をゼロとして lambda について方程式を解きます:
L:solve(D=0, lambda);
により、固有値は
と求まります。これらの固有値に対する固有ベクトルを求めていきます。
:
u1:matrix([x1],[y1])$ B1:A.u1$ solve([B1[1][1]=rhs(L[1])*u1[1][1],B1[2][1]=rhs(L[1])*u1[2][1]],[x1,y1]);
より、ベクトルの成分ついて
という条件が課されます。よって、例えば
ととることができます。同様にして固有値 に対する固有ベクトル は
u2:matrix([x2],[y2])$ B2:A.u2$ X:solve([B2[1][1]=rhs(L[2])*u2[1][1],B2[2][1]=rhs(L[2])*u2[2][1]],[x2,y2]);
より
ととれます。
さて、行列 を対角化するために行列 を
によって導入します:
P:matrix([2,1],[3,-1]);
すると
invert(P).A.P;
より
と計算されます。ただし
invert(P);
より
です。
ここで を考えると
と計算できます。実際
Apow[n]:=1/5*matrix([2*4^n+3*(-1)^n, 2*4^n-2*(-1)^n],[3*4^n-3*(-1)^n, 3*4^n+2*(-1)^n])$
から
Apow[0]; Apow[1]; Apow[10]; A^^10;
などによって確認することができます。
行列の対角化(Maximaの組み込み)
固有値、固有ベクトルは原理的には以上のようにして求められるのですが、これらを求める処理がMaximaにはすでに組み込まれています:まず固有値ですが
eigenvalues(A);
により、その値と縮重度を知ることができます。さらに
eigenvectors(A);
では、固有値とその縮重度に加え、それらに対応する固有ベクトルを出力します。
よってこれまでの計算は、例えば、上で考えた とはまた違った行列
に対して
B:matrix([1,2],[2,1])$ E:eigenvectors(B)$ v1:matrix([E[2][1][1][1]],[E[2][1][1][2]])$ v2:matrix([E[2][2][1][1]],[E[2][2][1][2]])$ Q:matrix([v1[1][1],v2[1][1]],[v1[2][1],v2[2][1]]); Lamb:ratsimp(invert(Q).B.Q);
などとして計算できます。
入力する行列を変えることで、それに応じた対角行列が出力されます。
エラーが出る場合は、固有値に縮重度があって線型独立な固有ベクトルがとれない場合になります。
縮重度が2で線型独立な固有ベクトルがない場合は例えば
C1:matrix([1,2],[-2,5])$ eigenvectors(C1);
などです。詳しくは以下の記事で解説しています:
3次正方行列への拡張
これまで 行列についてみてきましたが、 行列に拡張するには
B:matrix([1,2,1],[2,1,1],[0,2,1])$ E:eigenvectors(B)$ v1:matrix([E[2][1][1][1]],[E[2][1][1][2]],[E[2][1][1][3]])$ v2:matrix([E[2][2][1][1]],[E[2][2][1][2]],[E[2][2][1][3]])$ v3:matrix([E[2][3][1][1]],[E[2][3][1][2]],[E[2][3][1][3]])$ Q:expand(ratsimp(matrix( [v1[1][1],v2[1][1],v3[1][1]], [v1[2][1],v2[2][1],v3[2][1]], [v1[3][1],v2[3][1],v3[3][1]]))); Lamb:expand(ratsimp(invert(Q).B.Q));
などとして同様に計算できます。
本記事で学んだ内容の応用として、次の記事もあわせて読んでみてください: