求解系数矩阵由16阶Hilbert方程组构成的线性方程组,右端项为
[2877/851,3491/1431,816/409,2035/1187,2155/1423,538/395,1587/1279,573/502,…
947/895,1669/1691,1589/1717,414/475,337/409,905/1158,1272/1711,173/244]
要求:
1.Gauss_Sedel迭代法
2.最速下降法
3.共轭梯度法
一、 Gauss_Sedel迭代法
1.Gauss_Sedel迭代法函数
function [ x,k] =gauss_seidel_method( a,b,x0 )
n=length(x0);%x0的维数
k=1;%迭代次数
N=1000;%限制最大迭代次数
x=x0;
while k<N
for i=1:n%迭代过程
if i==1
x(1,1)=(b(1,1)-a(1,2:n)*x0(2:n,1))/a(1,1);
else if i<n
x(i,1)=(b(i,1)-a(i,1:i-1) *x(1:i-1,1)-a(i,i+1:n)*x0(i+1:n,1))/a(i,i);
else
x(n,1)=(b(n,1)-a(n,1:n-1)*x(1:n-1,1))/a(n,n);
end
end
end
if norm(x-x0)<10^(-4)%精度
break;
else
k=k+1;
x0=x;
end
end
end
2.Gauss_Sedel迭代法执行文件
a=hilb(16);
b=[2877/851,3491/1431,816/409,2035/1187,2155/1423,538/395,1587/1279,573/502,...
947/895,1669/1691,1589/1717,414/475,337/409,905/1158,1272/1711,173/244]';
x0=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]';%初始值
[x,k]=gauss_seidel_method(a,b,x0);
3.Gauss_Sedel迭代法执行结果
x =
1.0021
0.9630
1.1095
0.9648
0.9377
0.9530
0.9838
1.0146
1.0295
1.0517
1.0403
1.0246
1.0351
0.9759
0.9719
0.9388
k =
859
二、 最速下降法
1.最速下降法函数
function [ x,k ] = max_speed_method( a,b,x0 )%最速下降法
k=1;%迭代次数
N=1000;%限制最大迭代次数
x=x0;
while k<N%迭代过程
r=b-a*x;%步长
l=dot(r,r)/dot(a*r,r);%方向
x=x+l*r;
if norm(x-x0)<10^(-4)%精度
break;
else
k=k+1;
end
x0=x;
end
end
2.最速下降法执行文件
a=hilb(16);
b=[2877/851,3491/1431,816/409,2035/1187,2155/1423,538/395,1587/1279,573/502,... 947/895,1669/1691,1589/1717,414/475,337/409,905/1158,1272/1711,173/244]';
x0=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]';初始值
[x,k]=max_speed_method(a,b,x0);
3.最速下降法结果
x =
0.9986
1.0159
0.9831
0.9800
0.9907
1.0030
1.0128
1.0187
1.0204
1.0191
1.0144
1.0072
0.9988
0.9880
0.9766
0.9642
k =
423
三、 共轭梯度法
1.共轭梯度法函数
function [ x,k ] = CG_method( a,b,x0 )%共轭梯度法函数
x=x0;%初始赋值
r0=b-a*x0;
p0=r0;
k=1;%迭代次数
N=1000;%限制最大迭代次数
while norm(p0)>10^(-4)
m=dot(r0,r0)/dot(a*p0,p0);
x=x0+m*p0;
r1=r0-m*a*p0;
n=dot(r1,r1)/dot(r0,r0);
p1=r1+n*p0;
p0=p1;
r0=r1;
x0=x;
if r0==0
break;
else if dot(p0,a*p0)==0
break;
else if k>=N
break;
else
k=k+1;
end
end
end
end
end
2.共轭梯度法执行文件
a=hilb(16);
b=[2877/851,3491/1431,816/409,2035/1187,2155/1423,538/395,1587/1279,573/502,...
947/895,1669/1691,1589/1717,414/475,337/409,905/1158,1272/1711,173/244]';
x0=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]';%初始值
[x,k]=CG_method(a,b,x0);
3.共轭梯度法结果
x =
0.9965
1.0271
0.9776
0.9732
0.9861
1.0013
1.0132
1.0205
1.0230
1.0216
1.0168
1.0092
0.9998
0.9884
0.9761
0.9628
k =
5