ガウス消去法によって、行列のLU分解ができます。以下、$4\times 4$の行列$A$を例にして説明します。
Aの第i列の対角成分以下の要素を消去するために、以下の作業行列$E_i$が用意されます。
$$ E_1= \left( \begin{array}{cccc} 1 & 0 & 0 & 0\\ f_{21} & 1 & 0 & 0 \\ f_{31} & 0 & 1 & 0\\ f_{41} & 0 & 0 & 1 \end{array} \right) $$
$$ E_2= \left( \begin{array}{cccc} 1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 \\ 0 & f_{32} & 1 & 0\\ 0 & f_{42} & 0 & 1 \end{array} \right) $$
$$ E_3= \left( \begin{array}{cccc} 1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0\\ 0 & 0& f_{43} & 1 \end{array} \right) $$ ガウス消去法によって、以下の分解ができます。 $$ E_3E_2E_1 A = U $$ よって、 $$ A = E_1^{-1} E_2^{-1} E_3^{-1} U = LU \quad (L := E_1^{-1} E_2^{-1} E_3^{-1} ) $$
Lの各成分を計算すると、以下の形が分かります。
$$ L= \left( \begin{array}{cccc} 1 & 0 & 0 & 0\\ -f_{21} & 1 & 0 & 0 \\ -f_{31} & -f_{32} & 1 & 0\\ -f_{41} & -f_{42} & -f_{43} & 1 \end{array} \right) $$
上記の計算過程によって、ガウス消去法の各列の作業とともに、Lの各列の成分を計算でき、Aの対角成分以下に格納することができます。 即ち、ガウス消去法の計算後、以下の行列が得られます。
$$ \left( \begin{array}{cccc} u_{11} & u_{12} & u_{13} & u_{14}\\ -f_{21} & u_{22} & u_{23} & u_{24} \\ -f_{31} & -f_{32} & u_{33} & u_{34}\\ -f_{41} & -f_{42} & -f_{43} & u_{44} \end{array} \right) $$
行列のLU分解を行い、上記の説明のようにLとUを求めてください。
$$ A= \left( \begin{array}{cccc} 2& 1 & -2 & 3\\ 2 & 0 & 4 & 1 \\ 3 & 0 & 5& 2\\ 1 & -1 & 2 & 1 \end{array} \right) $$
Aの行列の第１行を使って、行列の一番目の列の対角以下の成分を０に変形します。
まず、第１行のalpha倍を第2,3,4行に加えます。 作業行列$E_1$の成分 $f_{j1}$ (j=2,3,4)$をAの第１列格納する。
A = [2, 1, -2, 3; 2, 0, 4, 1; 3, 0, 5, 2; 1, -1, 2, 1];
n = 4;
%対角の第１成分からの作業
i = 1;
%----------- 第１行によるAの第１列の対角以下の成分を消去する -----------
for j=2:n
%作業の行の番号
alpha = -A(j,1)/A(1,1);
for k=2:n %ｋ：各列の番号
A(j,k) = A(j,k) + alpha*A(1,k);
end
%作業行列のf_ｊ1をAの第１列格納する。
A(j,1) = - alpha ;
end
%変形後のA
A
Aの行列の第2行を使って、2番目の列の対角以下の成分を０に変形します。
第2行のalpha倍を第3,4行に加えます。 作業行列$E_2$の成分 $f_{j2}$ (j=3,4)$をAの第2列の対角以下に格納する。
i=2;
%----------- 第２行によるAの第２列の対角以下の成分を消去する -----------
for j=2:n
%作業の行の番号
alpha = ??;
for k=??:n %ｋ：各列の番号
A(j,k) = A(j,k) + ??;
end
%作業行列のf_ｊ2をAの第2列に格納する。
A(??,??) = ?? ;
end
%変形後のA
A
Aの行列の第3行を使って、3番目の列の対角以下の成分を０に変形します。
第2行のalpha倍を第4行に加えます。 作業行列$E_3$の成分 $f_{43}$をAの第3列の対角以下に格納する。
i=3;
%----------- 第3行によるAの第3列の対角以下の成分を消去する -----------
%ここにコードを書いてください。
上記の計算をここにまとめて書いてください。行の番号iに関するループを使ってください。
%ここにコードを書いてください。
AのLU分解が分かると、以下の２つの三角行列に関する方程式を解いてから、Ax=bの解は容易に計算できます。
$$ Ly=b, Ux=y $$
Lは以下の形を持っています。 $$ L= \left( \begin{array}{cccc} 1 & 0 & 0 & 0\\ l_{21} & 1 & 0 & 0 \\ l_{31} & l_{32} & 1 & 0\\ l_{41} & l_{42} & l_{43} & 1 \end{array} \right) $$
よって、$Ly=b$の解$y=(y_1, y_2, y_3, y_4)$は以下のように$y_i$を順番に算出することができます。
for loopでまとめます。
for i=1:n
y(i) = b(i)
for j = 1:(i-1)
y(i) = y(i) - y(j)*L(i,j)
end
end
講義番号：１８１S1520 水準 ０４
上記のAのLU分解に利用して、以下の$b_1$, $b_2$に関する連立一次方程式$Ax=b_1, Ax=b_2$を解いてください。 $$ b_1= \left( \begin{array}{cccc} 11\\ -5\\ -5\\ -4 \end{array} \right),\quad b_2= \left( \begin{array}{cccc} 10\\ 18\\ 26\\ 9 \end{array} \right) $$
以下の行列$A_1$,$A_2$のLU分解を計算してください。
$A_1$は以下のように生成されます。
A1=diag(2*ones(10,1))+diag(ones(9,1),1)+diag(ones(9,1),-1)
A1=diag(2*ones(10,1))+diag(ones(9,1),1)+diag(ones(9,1),-1)
$A_2 x=b$を解いてください。ただし、$b$はすべて成分が１である$6\times 1$ベクトルです。