微分方程式の数値シミュレーションに 現われる連立1次方程式の係数行列は、 ほとんどの場合、 成分に 0 が多い。 そのような行列をそ疎 (sparse) であるという。 疎行列はその性質をうまく利用することで、効率的な計算が可能になり、 大規模な問題が解けるようになる。 行列が疎であっても、 その逆行列は疎であるとは限らない (大抵の場合、疎ではない) ことに注意する。
MATLAB のコマンド・メモ (2) | ||||||||||||||||||
|
以下の例では、
という三重対角行列 (tridiagonal matrix) について、 連立1次方程式を解いたり、 逆行列を計算したりしている。
>> n=4 n = 4 >> 1:n ans = 1 2 3 4 >> (1:n)'*(1:n) ans = 1 2 3 4 2 4 6 8 3 6 9 12 4 8 12 16 >> eye(n,n) ans = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 >> zeros(n,n) ans = 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 >> ones(n-1,1) ans = 1 1 1 >> diag(ones(n-1,1),1) ans = 0 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 >> diag(ones(n-1,1),-1) ans = 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 0 >> a=2*eye(n,n)-diag(ones(n-1,1),1)-diag(ones(n-1,1),-1) a = 2 -1 0 0 -1 2 -1 0 0 -1 2 -1 0 0 -1 2 >> x=(1:n)' x = 1 2 3 4 >> b=a*x b = 0 0 0 5 >> ia=inv(a) ia = 0.8000 0.6000 0.4000 0.2000 0.6000 1.2000 0.8000 0.4000 0.4000 0.8000 1.2000 0.6000 0.2000 0.4000 0.6000 0.8000 >> ia*a ans = 1.0000 -0.0000 -0.0000 0.0000 0 1.0000 -0.0000 0.0000 0 0 1.0000 0.0000 0 0 0 1.0000 >> n=8;a=2*eye(n,n)-diag(ones(n-1,1),1)-diag(ones(n-1,1),-1) a = 2 -1 0 0 0 0 0 0 -1 2 -1 0 0 0 0 0 0 -1 2 -1 0 0 0 0 0 0 -1 2 -1 0 0 0 0 0 0 -1 2 -1 0 0 0 0 0 0 -1 2 -1 0 0 0 0 0 0 -1 2 -1 0 0 0 0 0 0 -1 2 >> inv(a) ans = 0.8889 0.7778 0.6667 0.5556 0.4444 0.3333 0.2222 0.1111 0.7778 1.5556 1.3333 1.1111 0.8889 0.6667 0.4444 0.2222 0.6667 1.3333 2.0000 1.6667 1.3333 1.0000 0.6667 0.3333 0.5556 1.1111 1.6667 2.2222 1.7778 1.3333 0.8889 0.4444 0.4444 0.8889 1.3333 1.7778 2.2222 1.6667 1.1111 0.5556 0.3333 0.6667 1.0000 1.3333 1.6667 2.0000 1.3333 0.6667 0.2222 0.4444 0.6667 0.8889 1.1111 1.3333 1.5556 0.7778 0.1111 0.2222 0.3333 0.4444 0.5556 0.6667 0.7778 0.8889 >>
が疎であっても、 の成分には つも 0 がないことを確認しよう。
任意の正則行列 に対して、
であれば、 であるから、
となり7、 , の掛け算は非常に簡単なので (付録 D.5.2 参照)、 , , が求まっていれば連立1次方程式は簡単に解ける。 (1), (2) を満たす , , を 求めることを を LU 分解するという。
MATLAB のコマンド・メモ (3) | ||||||||||
|
>> [L u]=lu(a) L = 1.0000 0 0 0 0 0 0 0 -0.5000 1.0000 0 0 0 0 0 0 0 -0.6667 1.0000 0 0 0 0 0 0 0 -0.7500 1.0000 0 0 0 0 0 0 0 -0.8000 1.0000 0 0 0 0 0 0 0 -0.8333 1.0000 0 0 0 0 0 0 0 -0.8571 1.0000 0 0 0 0 0 0 0 -0.8750 1.0000 u = 2.0000 -1.0000 0 0 0 0 0 0 0 1.5000 -1.0000 0 0 0 0 0 0 0 1.3333 -1.0000 0 0 0 0 0 0 0 1.2500 -1.0000 0 0 0 0 0 0 0 1.2000 -1.0000 0 0 0 0 0 0 0 1.1667 -1.0000 0 0 0 0 0 0 0 1.1429 -1.0000 0 0 0 0 0 0 0 1.1250 >> x=(1:n)' x = 1 2 3 4 5 6 7 8 >> b=a*x b = 0 0 0 0 0 0 0 9 >> u\(L\b) ans = 1.0000 2.0000 3.0000 4.0000 5.0000 6.0000 7.0000 8.0000 >>
を LU 分解したときの因子 , , が確かに 置換行列、下三角行列、上三角行列であることを確かめよう (このケースでは、 は単位行列である)。 , が疎行列であることも確かめよう。
今回は使わなかったが、よく使われるコマンドを紹介しておく。
MATLAB のコマンド・メモ (4) | ||||||||||||||||||||||
|
余裕があれば試してみよう | |
|