今回の授業では、Hilbert行列の例を利用して反復法の収束性を検討してみます。
Hilbert行列を作成します。MATLABではhilb(n)という命令で$n\times n$のHilbert行列を作成することができます。
$$ H=( h_{ij})_{n \times n}, \quad h_{ij}=\frac{1}{j+j-1}\:. $$
n=2
A=hilb(n)
n=3
A=hilb(n)
n=3
A=hilb(n)
printf("行列式= %15.10f \n\n", det(A))
myeig = eig(A);
printf("固有値=\n\n")
disp(myeig)
cond_num = max(myeig)/min(myeig);
printf("\n条件数= %15.10f \n", cond_num)
SOR法を使って、Ax=bの解を求めてください。ここで、$b$のすべての成分を１とする。
以下のコードでは、SOR法が関数として定義されています。
特に、trilとtriuを使って、行列の下三角行列と上三角行列を取ることができます。
function [sol, errlist, iteration_num] = sor(A,b,w,max_N)
E=tril(A,-1);
F=triu(A,1);
D=diag(diag(A));
exact_sol = A\b;
M=D/w+E; N=(w-1)/w*D+F;
x=zeros(length(b),1);
errlist = [];
for k=1:max_N;
x = M\(b-N*x);
errlist(k) = norm(x - exact_sol);
if errlist(k)<1E-5
break
end
end
sol=x;
iteration_num =k;
end
A=hilb(n);
b = ones(n,1); % Vector b can also be inputed by: b=[1,1,1,1]';
[sol, errlist,num] = sor(A,b,1.5,1000);
sol
num
plot(log10(errlist))
grid on
$N=3$の場合、$1.1\le w \le 1.9$の中の$w$について、近似解のL2ノルムの誤差が$10^{-6}$以下になるまでに必要な反復計算の回数を調べて、一番良い$w$を選んでください。
$N=3,4,5$の場合、上記の計算で定めた$w$の値を利用して、SOR法で近似解を求めてください。近似解のL2ノルムの誤差が$10^{-6}$以下ことが必要です。特に、反復計算回数が20000になる可能性があるので、大きい反復計算回数で試してみてください。
行列の条件と反復計算の回数の関係を調べてください。
ガウス消去法を使って、Ax=bの解を求めてください。ここで、$b$のすべての成分を１とする。
n=3
A=hilb(n)
b = ones(n,1); % Vector b can also be inputed by: b=[1,1,1,1]';
sol_x = A\b
%ここでガウス消去法によって近似解を求めてください。（ピボット選択の使用は自由です。）
$N=3,4,5,6$の場合、Hilbert行列の条件数を調べて、近似解の誤差との関係を検討してみてください。