清华大学数值分析报告某实验报告材料Word文件下载.docx
《清华大学数值分析报告某实验报告材料Word文件下载.docx》由会员分享,可在线阅读,更多相关《清华大学数值分析报告某实验报告材料Word文件下载.docx(48页珍藏版)》请在冰点文库上搜索。
![清华大学数值分析报告某实验报告材料Word文件下载.docx](https://file1.bingdoc.com/fileroot1/2023-4/29/bd8021bb-c7dc-4ded-b1fb-a6ca1d3da654/bd8021bb-c7dc-4ded-b1fb-a6ca1d3da6541.gif)
最终解为x=(1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000001,0.999999999999998,1.000000000000004,0.999999999999993,1.000000000000012,0.999999999999979,1.000000000000028)T
使用无穷范数衡量误差,得到=2.842170943040401e-14,可以发现,采用顺序高斯消元法求得的解与精确解之间误差较小。
通过进一步观察,可以发现,按照顺序高斯消去法计算时,其选取的主元值和矩阵中其他元素大小相近,因此顺序高斯消去法方式并没有对结果造成特别大的影响。
若采用列主元高斯消元法,则结果为:
最终解为x=(1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000)T
同样使用无穷范数衡量误差,有=0;
若使用完全主元高斯消元法,则结果为
最终解x=(1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000,1.000000000000000)T
(2)
若每步都选取模最小或尽可能小的元素为主元,则计算结果为
最终解x=(1.0000000000000001.0000000000000001.0000000000000001.0000000000000010.9999999999999981.0000000000000040.9999999999999931.0000000000000120.9999999999999791.000000000000028)T
使用无穷范数衡量误差,有为2.842170943040401e-14;
而完全主元消去法的误差为=0。
从
(1)和
(2)的实验结果可以发现,列主元消去法和完全主元消去法都得到了精确解,而顺序高斯消去法和以模尽量小的元素为主元的消去法没有得到精确解。
在后两种消去法中,由于程序计算时的舍入误差,对最终结果产生了一定的影响,但由于方程组的维度较低,并且元素之间相差不大,所以误差仍比较小。
为进一步分析,计算上述4种方法每步选取的主元数值,并列表进行比较,结果如下:
第n次消元
顺序
列主元
完全主元
模最小
1
6.000000000000000
8
2
4.666666666666667
3
4.285714285714286
4
4.133333333333333
5
4.064516129032258
6
4.031746031746032
7
4.015748031496063
4.007843137254902
9
4.003913894324853
10
4.001955034213099
0.015617370605469
从上表可以发现,对这个方程组而言,顺序高斯消去选取的主元恰好事模尽量小的元素,而由于列主元和完全主元选取的元素为8,与4在数量级上差别小,所以计算过程中的累积误差也较小,最终4种方法的输出结果均较为精确。
在这里,具体解释一下顺序法与模最小法的计算结果完全一致的原因。
该矩阵在消元过程中,每次选取主元的一列只有两个非零元素,对角线上的元素为4左右,而其正下方的元素为8,该列其余位置的元素均为0。
在这样的情况下,默认的主元也就是该列最小的主元,因此两种方法所得到的计算结果是一致的。
理论上说,完全高斯消去法的误差最小,其次是列主元高斯消去法,而选取模最小的元素作为主元时的误差最大,但是由于方程组的特殊性(元素相差不大并且维度不高),这个理论现象在这里并没有充分体现出来。
(3)
时,重复上述实验过程,各种方法的计算结果如下所示,在这里,仍采用无穷范数衡量绝对误差。
顺序高斯消去法
列主元高斯消去
完全主元高斯消去
选取模最小或尽可能小元素作为主元消去
X
1.0000000000000001.0000000000000001.0000000000000001.0000000000000010.9999999999999981.0000000000000040.9999999999999931.0000000000000140.9999999999999721.0000000000000570.999999999999886
1.0000000000002270.9999999999995471.0000000000009020.9999999999982091.0000000000035240.9999999999931791.0000000000127320.9999999999781731.000000000029102
2.910205409989430e-11
可以看出,此时列主元和完全主元的计算结果仍为精确值,而顺序高斯消去和模尽可能小方法仍然产生了一定的误差,并且两者的误差一致。
与n=10时候的误差比相比,n=20时的误差增长了大约1000倍,这是由于计算过程中舍入误差的不断累积所致。
所以,如果进一步增加矩阵的维数,应该可以看出更明显的现象。
(4)
不同矩阵维度下的误差如下,在这里,为方便起见,选取2-条件数对不同维度的系数矩阵进行比较。
维度
条件数
顺序消去
模尽量小
1.7e+3
2.84e-14
1.8e+6
2.91e-11
5.7e+7
9.31e-10
1.8e+9
2.98e-08
1.9e+12
3.05e-05
3.8e+16
3.28e+04
3.88e-12
8.5e+16
3.52e+13
4.2e-3
从上表可以看出,随着维度的增加,不同方法对计算误差的影响逐渐体现,并且增长较快,这是由于舍入误差逐步累计而造成的。
不过,方法二与方法三在维度小于40的情况下都得到了精确解,这两种方法的累计误差远比方法一和方法四慢;
同样地,出于与前面相同的原因,方法一与方法四的计算结果保持一致,方法二与方法三的计算结果保持一致。
4.结论
本文矩阵中的元素差别不大,模最大和模最小的元素并没有数量级上的差异,因此,不同的主元选取方式对计算结果的影响在维度较低的情况下并不明显,四种方法都足够精确。
对比四种方法,可以发现采用列主元高斯消去或者完全主元高斯消去法,可以尽量抑制误差,算法最为精确。
不过,对于低阶的矩阵来说,四种方法求解出来的结果误差均较小。
另外,由于完全选主元方法在选主元的过程中计算量较大,而且可以发现列主元法已经可以达到很高的精确程度,因而在实际计算中可以选用列主元法进行计算。
附录:
程序代码
clear
clc;
formatlong;
%方法选择
n=input('
矩阵A阶数:
n='
);
disp('
选取求解方式'
1顺序Gauss消元法,2列主元Gauss消元法,3完全选主元Gauss消元法,4模最小或近可能小的元素作为主元'
a=input('
求解方式序号:
'
%赋值A和b
A=zeros(n,n);
b=zeros(n,1);
fori=1:
n
A(i,i)=6;
ifi>
A(i,i-1)=8;
end
ifi<
n
A(i,i+1)=1;
end
n
forj=1:
b(i)=b(i)+A(i,j);
end
给定系数矩阵为:
A
右端向量为:
b
%求条件数及理论解
线性方程组的精确解:
X=(A\b)'
fprintf('
矩阵A的1-条件数:
%f\n'
cond(A,1));
矩阵A的2-条件数:
cond(A));
矩阵A的无穷-条件数:
cond(A,inf));
%顺序Gauss消元法
ifa==1
A1=A;
b1=b;
fork=1:
ifA1(k,k)==0
disp('
主元为零,顺序Gauss消元法无法进行'
break
end
fprintf('
第%d次消元所选取的主元:
%g\n'
k,A1(k,k))
%disp('
此次消元后系数矩阵为:
%A1
forp=k+1:
l=A1(p,k)/A1(k,k);
A1(p,k:
n)=A1(p,k:
n)-l*A1(k,k:
n);
b1(p)=b1(p)-l*b1(k);
x1(n)=b1(n)/A1(n,n);
fork=n-1:
-1:
forw=k+1:
n
b1(k)=b1(k)-A1(k,w)*x1(w);
end
x1(k)=b1(k)/A1(k,k);
disp('
顺序Gauss消元法解为:
disp(x1);
所求解与精确解之差的无穷-范数为'
norm(x1-X,inf)
%列主元Gauss消元法
ifa==2
A2=A;
b2=b;
[max_i,max_j]=find(A2(:
k)==max(abs(A2(k:
n,k))));
ifmax_i~=k
A2_change=A2(k,:
A2(k,:
)=A2(max_i,:
A2(max_i,:
)=A2_change;
b2_change=b2(k);
b2(k)=b2(max_i);
b2(max_i)=b2_change;
ifA2(k,k)==0
主元为零,列主元Gauss消元法无法进行'
k,A2(k,k))
%A2
l=A2(p,k)/A2(k,k);
A2(p,k:
n)=A2(p,k:
n)-l*A2(k,k:
b2(p)=b2(p)-l*b2(k);
x2(n)=b2(n)/A2(n,n);
b2(k)=b2(k)-A2(k,w)*x2(w);
x2(k)=b2(k)/A2(k,k);
列主元Gauss消元法解为:
disp(x2);
norm(x2-X,inf)
%完全选主元Gauss消元法
ifa==3
A3=A;
b3=b;
VV=eye(n);
[max_i,max_j]=find(A3(k:
n,k:
n)==max(max(abs(A3(k:
n)))));
ifnumel(max_i)==0
[max_i,max_j]=find(A3(k:
n)==-max(max(abs(A3(k:
W=eye(n);
W(max_i
(1)+k-1,max_i
(1)+k-1)=0;
W(k,k)=0;
W(max_i
(1)+k-1,k)=1;
W(k,max_i
(1)+k-1)=1;
V=eye(n);
V(k,k)=0;
V(max_j
(1)+k-1,max_j
(1)+k-1)=0;
V(k,max_j
(1)+k-1)=1;
V(max_j
(1)+k-1,k)=1;
A3=W*A3*V;
b3=W*b3;
VV=VV*V;
ifA3(k,k)==0
主元为零,完全选主元Gauss消元法无法进行'
k,A3(k,k))
%A3
l=A3(p,k)/A3(k,k);
A3(p,k:
n)=A3(p,k:
n)-l*A3(k,k:
b3(p)=b3(p)-l*b3(k);
x3(n)=b3(n)/A3(n,n);
b3(k)=b3(k)-A3(k,w)*x3(w);
x3(k)=b3(k)/A3(k,k);
完全选主元Gauss消元法解为:
disp(x3);
norm(x3-X,inf)
%模最小或近可能小的元素作为主元
ifa==4
A4=A;
b4=b;
AA=A4;
AA(AA==0)=NaN;
[min_i,j]=find(AA(k:
n,k)==min(abs(AA(k:
ifnumel(min_i)==0
[min_i,j]=find(AA(k:
n,k)==-min(abs(AA(k:
n))));
W(min_i
(1)+k-1,min_i
(1)+k-1)=0;
W(min_i
(1)+k-1,k)=1;
W(k,min_i
(1)+k-1)=1;
A4=W*A4;
b4=W*b4;
ifA4(k,k)==0
主元为零,模最小Gauss消元法无法进行'
k,A4(k,k))
%A4
l=A4(p,k)/A4(k,k);
A4(p,k:
n)=A4(p,k:
n)-l*A4(k,k:
b4(p)=b4(p)-l*b4(k);
x4(n)=b4(n)/A4(n,n);
b4(k)=b4(k)-A4(k,w)*x4(w);
x4(k)=b4(k)/A4(k,k);
模最小Gauss消元法解为:
disp(x4);
norm(x4-X,inf)
二、实验3.3
考虑方程组的解,其中系数矩阵H为Hilbert矩阵:
这是一个著名的病态问题。
通过首先给定解(例如取为各个分量均为1)再计算出右端的办法给出确定的问题。
(1)选择问题的维数为6,分别用Gauss消去法(即LU分解)、J迭代法、GS迭代法和SOR迭代法求解方程组,其各自的结果如何?
将计算结果与问题的解比较,结论如何。
(2)逐步增大问题的维数,仍用上述的方法来解它们,计算的结果如何?
计算的结果说明的什么?
(3)讨论病态问题求解的算法。
1.算法设计
对任意线性方程组,分析各种方法的计算公式如下,
(1)Gauss消去法:
首先对系数矩阵进行LU分解,有,则原方程转化为,令,则原方程可以分为两步回代求解:
具体方法这里不再赘述。
(2)J迭代法:
首先分解,再构造迭代矩阵,其中,进行迭代计算,直到误差满足要求。
(3)GS迭代法:
首先分解,再构造迭代矩阵,其中,进行迭代计算,直到误差满足要求。
(4)SOR迭代法:
首先分解,再构造迭代矩阵,其中,进行迭代计算,直到误差满足要求。
2.实验过程
一、根据维度n确定矩阵H的