实验三 迭代法解线性方程组

实验三 迭代法解线性方程组
实验目的
学会用Jacobi 迭代法、Gauss-Seidel 迭代法和超松驰迭代法求线性方程组解。

学会对迭代法做收敛性分析，研究求方程组解的最优迭代方法。

学会用共轭梯度法求线性方程组的解，研究共轭梯度法的计算效率。

实验要求
按照题目要求完成实验内容。

写出相应的Matlab 程序。

给出实验结果。

对实验结果进行分析讨论。

写出相应的实验报告。

实验步骤
1、研究Jacobi 迭代法求解线性方程组的方法和相应的收敛性。

（1）用Jacobi 迭代法（Jacobi.m ）求解线性方程组
⎪⎪⎩⎪⎪⎨⎧-=+--=-+-=--71912263532311321321321x x x x x x x x x （4.31） 取初始点()()T
x 0,0,00=，精度要求为105
-=ε。

请给出满足精度要求的迭代次数和相应的
计算结果。

function [x,k]=jc(a,b,x0,ep,max)
n=length(a);
k=0;
if nargin<5
max=500;
end
if nargin<4
ep=1e-5;
end
if nargin<3
x0=zeros(n,1);
y=zeros(n,1);
end
x=x0;x0=x+2*ep;
while norm(x0-x,inf)>ep&&k<max
k=k+1;x0=x;
for i=1:n
y(i)=b(i);
for j=1:n
if j~=i
y(i)=y(i)-a(i,j)*x0(j);
end
end
if abs(a(i,i))<1e-10||k==max
warning('a(i,i) ̫С');
return ;
end
y(i)=y(i)/a(i,i);
end
x=y;
end
end
>> a=[11 -3 -2;-1 5 -3;-2 -12 19];
>> b=[3 6 -7]';
>> [x,k]=jc(a,b)
x =
0.999985953146656
1.999977859006902
0.999978180649185
k = 33
研究相应迭代矩阵的谱半径和Jacobi 迭代的渐近收敛速度。

2、研究Gauss-Seidel 迭代法求解线性方程的方法和相应的收敛性
（1）用Gauss-Seidel 迭代法求解线性方程组，取初始点T x )0,0,0()0(=，精度要求为510-=ε，请给出满足精度要求的迭代次数和相应的计算结果。

function [x,k]=gsd(a,b,x0,ep,max)
n=length(a);
k=0;
if nargin<5
max=500;
end
if nargin<4
ep=1e-5;
end
if nargin<3
x0=zeros(n,1);
x=x0;x0=x+2*ep;
while norm(x0-x,inf)>ep&&k<max
k=k+1;
x0=x;
z=x;
for i=1:n
z(i)=b(i);
for j=1:n
if j~=i
z(i)=z(i)-a(i,j)*x(j);
end
end
if abs(a(i,i))<1e-10||k==max
warning('a(i,i) 太小');
return
end
z(i)=z(i)/a(i,i);
x(i)=z(i);
end
end
if k==max
warning('已达上限');
end
end
（2）研究相应迭代矩阵的谱半径和Gauss-Seidel 迭代的渐近收敛速度。

3、研究松弛因子对超松弛迭代法的影响
（1）选取,6.1,4.1,2.1,1,8.0=ω用超松弛（SOR ）迭代法（sor.m)求解线性方程组（4.31），取初始点T x )0,0,0()0(=，精度要求为510-=ε，列出不同松弛因子下的迭代次数。

function [x,k]=sor(a,b,om,ep,Nmax)
n=length(b);
k=0;
if nargin<6
Nmax=500;
end
if nargin<5
ep=1e-5;
if nargin<4
x0=zeros(n,1);
end
if nargin<3
om=1.5;
end
x=x0;
x0=x+2*ep;k=0;
l=tril(a,-1);
u=tril(a,1);
while norm(x0-x,inf)>ep&k<Nmax
k=k+1;
x0=x;
for i=1:n
x1(i)=(b(i)-l(i,1:i-1)*x(1:i-1,1))-u(i,i+1:n)*x0(i+1:n,1))/a(i,i);
x(i)=(1-om)*x0(i)+om*x1(i);
end
end
if k==max
warning('已达上限');
end
end
（2）用算法4.4（求最优松弛因子的0.618法）求相应的最优松弛因子，并与式（4.20）的计算结果相比较
function [wopt,ropt]=sor618(A)
n=length(A);
D=zeros(n);
L=zeros(n);
U=zeros(n);
for i=1:n
D(i,i)=A(i,i);
for j=1:i-1
L(i,j)=-A(i,j);
end
for j=i+1:n
U(i,j)=-A(i,j);
end
end
t=2/(sqrt(5)+1);
ep=1e-5;
a=0;
b=2;
w1=a+(1-t)*(b-a);
B=(D-w1*L)\((1-w1)*D+w1*U);
r1=max(abs(eig(B)));
wr=a+t*(b-a);
B=(D-wr*L)\((1-wr)*D+wr*U);
rr=max(abs(eig(B)));
while abs(b-a)>ep
if r1<rr
b=wr;
wr=w1;
rr=r1;
w1=a+(1-t)*(b-a);
B=(D-wr*L)\((1-wr)*D+w1*U);
r1=max(abs(eig(B)));
else
a=w1;
w1=wr;
r1=rr;
wr=a+t*(b-a);
B=(D-wr*L)\((1-wr)*D+w1*U);
rr=max(abs(eig(B)));
end
if r1<r1
wopt=w1;
ropt=r1;
else
wopt=wr;
ropt=rr;
end
end
（3）取算法 4.4计算出的最优松弛因子opt ω，用SOR 迭代法再次求解线性方程组（4.31）其迭代次数与2.1=ω的迭代次数相比较
4、讨论共轭梯度法求解线性方程组的方法
（1）用共轭梯度法求解Hilbert 系数矩阵的线性方程组，分别取n=60，80，100，b=z H n ，取
T z ）
，1,1,1(⋯=和T n z ),2,1(，⋯=
记录各自的迭代次数、计算结果与精确解的绝对误差和相对误差。

（2）考虑对共轭梯度法的推广，用共轭梯度法求解一般线性方程组
b Ax =
式中矩阵A 是非正定对称阵。

现考虑极小化问题
b)-T(A x b)-min(Ax (4.32)
其极小点满足
b A Ax A T T = （4.33）
注意到，当矩阵A 是非奇异矩阵时，A A T
为正定对称矩阵。

因此，可以对于方程组应用共轭梯度方法求解。

试用上述方法求解线性方程组 ⎪⎩
⎪⎨⎧-=-+=-+=++53.048.104.368.203.130.193.048.105.053.448.151.2321321321x x x x x x x x x
计算结果保留8位小数。

function[x,k]=getd(a,b,x0,ep,nmax)
n=length(a);
if nargin<5
nmax=500;
end
if nargin<4
ep=1e-10;
end
if nargin<3
x0=zeros(n,1);
end
x=x0;
x0=x+2*ep;
r=b'-a*x;
d=r;
k=0;
while norm(x0-x,inf)>ep&&k<nmax
k=k+1;
x0=x;
alpha=(r'*r)/(d'*a*d);
r1=r;
s=alpha*d;
x=x+s;
r=r-a*s;
beta=(r'*r/((r1)'*r1));
d=r+beta*d;
end
if k==nmax
warning('已达上限'); end。