《矩阵与数值分析》课程数值实验报告.docx
- 文档编号:7732235
- 上传时间:2023-01-26
- 格式:DOCX
- 页数:45
- 大小:100.62KB
《矩阵与数值分析》课程数值实验报告.docx
《《矩阵与数值分析》课程数值实验报告.docx》由会员分享,可在线阅读,更多相关《《矩阵与数值分析》课程数值实验报告.docx(45页珍藏版)》请在冰豆网上搜索。
《矩阵与数值分析》课程数值实验报告
《矩阵与数值分析》课程数值实验报告
一、Gauss消去法与Gauss列主元消去法求解方程组
1.1.问题
给定n阶方程组
,其中
,
则方程组有解
。
对
和
,分别用Gauss消去法和列主元消去法解方程组,并比较计算结果。
1.2.算法描述
1.2.1算法基本思路
1)Gauss消去法基本思路:
设有方程组
,设
是可逆矩阵。
高斯消去法的基本思想就是将矩阵的初等行变换作用于方程组的增广矩阵(A|b),将其化成一个上三角形矩阵,然后通过回代求这个上三角形方程组的解,就是原方程组的解。
2)Gauss列主元消去法基本思路:
设有方程组
,设
是可逆矩阵。
列主元高斯消去法的基本思想就是在Gauss消去法中增加选主元的过程,即在第k步消元时,首先在第k列主对角线以下(含主对角元)元素中挑选绝对值最大的数,并通过初等行变换,使得该数位于主对角线上,然后再继续消元。
1.2.2算法实现步骤:
1.初始化矩阵A和b;
2.选择消元方法:
1.Gauss消去法2.Gauss列主元消去法
3.当选择1时,进行Gauss消去:
进行n-1次消元,并实时输出每次消元结果,以便观测并检查每次消元结果;消元后得到上三角形矩阵,回代求解;
当选择2时,进行Gauss列主元消去:
进行n-1次消元,在每次消元之前查找列主元,并使列主元位于对角线上,然后再继续消元并实时输出每次消元结果,以便观测并检查每次消元结果;
4.输出最后的求解结果。
1.3程序代码
本程序是用C++语言在MicrosoftvisualC++编译环境下编译和运行的。
由于代码太长,程序代码见附录1
1.4计算结果:
1)Gauss消去法求解结果:
1.当矩阵阶数n=10,即附录1程序代码中的N设为10时,运行程序得到如下结果:
消去之前的a和b:
6100000000|7
8610000000|15
0861000000|15
0086100000|15
0008610000|15
0000861000|15
0000086100|15
0000008610|15
0000000861|15
0000000086|14
请选择消元方法:
1.Gauss消去法2.Gauss列主元消去法
1↘(由键盘输入数字1,并回车)
第1次消元结果:
6100000000|7
04.666710000000|5.6667
0861000000|15
0086100000|15
0008610000|15
0000861000|15
0000086100|15
0000008610|15
0000000861|15
0000000086|14
……
第9次消元结果:
6100000000|7
04.666710000000|5.6667
004.28571000000|5.2857
0004.1333100000|5.1333
00004.064510000|5.0645
000004.03171000|5.0317
0000004.0157100|5.0157
00000004.007810|5.0078
000000004.00391|5.0039
0000000004.002|4.002
求解结果:
X[10]={1111111111}
2.当矩阵阶数n=84,即程序代码中的N为84时,将附录1程序中N的值改为84,运行程序得到如下结果:
(由于此时矩阵太大,消元过程和1类似,所以省略写出。
)
求解结果:
X[84]={11111111111111111111111111111111110.99999810.9999921.000020.9999691.000060.9998781.000240.9995121.000980.9980471.003910.9921881.015620.9687521.06250.8750081.249980.5000311.99994-0.9998784.99976-6.9995116.999-30.99864.9961-126.992256.984-510.9691024.94-2046.874096.74-8190.4716383.9-32764.565531-131055262097-5241271.048e+006-2.09498e+0064.18586e+006-8.35533e+0061.66451e+007-3.30281e+0076.50077e+007-1.25821e+0082.34867e+008-4.02629e+0085.36838e+008}
2)Gauss列主元消去法求解结果:
1.当矩阵阶数n=10,即附录1程序代码中的N设为10时,运行程序得到如下结果:
消去之前的A和b:
6100000000|7
8610000000|15
0861000000|15
0086100000|15
0008610000|15
0000861000|15
0000086100|15
0000008610|15
0000000861|15
0000000086|14
请选择消元方法:
1.Gauss消去法2.Gauss列主元消去法
2↘(由键盘输入数字2,并回车)
第1次消元结果:
8610000000|15
0-3.5-0.750000000|-4.25
0861000000|15
0086100000|15
0008610000|15
0000861000|15
0000086100|15
0000008610|15
0000000861|15
0000000086|14
……
第9次消元结果:
8610000000|15
0861000000|15
0086100000|15
0008610000|15
0000861000|15
0000086100|15
0000008610|15
0000000861|15
0000000086|14
000000000-0.015617|-0.015617
求解结果:
X[10]={1111111111}
2.当矩阵阶数n=84,即程序代码中的N为84时,将附录1程序中N的值改为84,运行程序得到如下结果:
(由于此时矩阵太大,消元过程和1类似,所以省略写出。
)
求解结果:
X[84]={1111111111111111111111111111111111111111111111111111111111111111111111111111111110.99999910.999997}
1.5结果分析:
1.当n=10时,用Gauss消去法和Gauss列主元消去法得到的值解和精确解
全吻合,此时无论Gauss消法法还是列主元消去法求解,可以得到理想结果。
2.当n=84,即矩阵阶数过大时,可以看出此时用Gauss消元法得到的数值解结果与精确解结果
相去甚远,误差很大。
这就是因为当矩阵太大时,在消元过程中出现了小数做除数的情形,在计算过程中的舍入误差使解的面目全非了。
而用Gauss列主元消元法得到的数值解与精确解仍然可以很好的吻合,所以比较Gauss消去法与Gauss列主元消去法求解方程组可知,后者能有效避免消元过程中的小主元做除数,从而使求解精度更高。
附录1:
#include
#include
usingnamespacestd;
#defineN10//N表示矩阵大小
voidSet_value(doublea[][N],double*b);
voidGauss_Array(double(*a)[N],double*b,double*X);//Gauss消元求解函数
voidSelect_main(doublea[][N],double*b,inti);
voidSelectmain_Array(double(*a)[N],double*b,double*X);//Gauss列主元消去法
voiddisplay(doublea[][N],double*b);//输出当前(A|b)的函数
intmain()
{
doubleA[N][N],b[N];
Set_value(A,b);//初始化系数矩阵A和右端向量b
doubleX[N];//向量X用于存求解结果
intchoice;
cout<<"请选择消元方法:
1.Gauss消元2.Gauss列主元消元"< cin>>choice; switch(choice) { case1: Gauss_Array(A,b,X);break; case2: Selectmain_Array(A,b,X);break; default: cout<<"输入有误! ";exit(0); } cout<<"求解结果: "< cout<<"X["< for(inti=0;i cout< cout<<"}"< return0; } voidSet_value(doublea[][N],double*b) { for(inti=0;i { for(intj=0;j { a[i][j]=0; if(j==i-1) a[i][j]=8; if(j==i) a[i][j]=6; if(j==i+1) a[i][j]=1; } } for(inti=0;i { b[i]=15; if(i==0) b[i]=7; if(i==N-1) b[i]=14; } cout<<"消去之前的A和b: "< display(a,b); } voidGauss_Array(double(*a)[N],double*b,double*X) { for(inti=0;i { for(intI=i+1;I { if(a[i][i]==0)//如果对角线元素出现为0的情况,则矩阵奇异,退出程序 { cout<<"矩阵奇异! "< exit(0); } doubleL=a[I][i]/a[i][i]; for(intj=i;j a[I][j]=a[I][j]-L*a[i][j]; b[I]=b[I]-L*b[i]; } cout<<"第"< "< display(a,b);//输出每次消元的结果 } for(intk=N-1;k>=0;k--)//倒序求解方程 { if(a[k][k]==0)//如果消元后对角线元素出现为0的情况,则矩阵奇异,退出程序 { cout<<"矩阵奇异! "< exit(0); } doublesum=0; for(intm=N-1;m>k;m--) sum+=a[k][m]*X[m]; X[k]=(b[k]-sum)/a[k][k]; } } voidSelect_main(doublea[][N],double*b,intj) { intmaxi=j; doublemax=a[j][j]; for(inti=j;i if(fabs(a[i][j])>fabs(max)) { max=a[i][j]; maxi=i; } if(maxi! =j)//将列主元所在行放至列首 { doubleTemp=b[j]; b[j]=b[maxi]; b[maxi]=Temp; for(intm=j;m { doubletemp=a[j][m]; a[j][m]=a[maxi][m]; a[maxi][m]=temp; } } } voidSelectmain_Array(double(*a)[N],double*b,double*X) { for(inti=0;i { Select_main(a,b,i); for(intI=i+1;I { if(a[i][i]==0)//如果对角线元素出现为0的情况,则矩阵奇异,退出程序 { cout<<"矩阵奇异! "< } doubleL=a[I][i]/a[i][i]; for(intj=i;j a[I][j]=a[I][j]-L*a[i][j]; b[I]=b[I]-L*b[i]; } cout<<"第"< "< display(a,b);//输出每次消元结果 } for(intk=N-1;k>=0;k--)//倒序求解方程 { if(a[k][k]==0)//如果消元后对角线元素出现为0的情况,则矩阵奇异,退出程序 { cout<<"矩阵奇异! "< } doublesum=0; for(intm=N-1;m>k;m--) sum+=a[k][m]*X[m]; X[k]=(b[k]-sum)/a[k][k]; } } voiddisplay(double(*a)[N],double*b)//输出系数矩阵A和b { for(inti=0;i { for(intj=0;j cout< cout<<"|"< } } 二、用Jacobi和Gauss-Seidel迭代法求解方程组 2.1.问题 设有方程组 ,其中 , (a)选取不同的初始向量 和不同的右端向量 ,给定迭代误差要求,用Jacobi和Gauss-Seidel迭代法计算,观测得出的迭代向量序列是否收敛。 若收敛,记录迭代次数,分析计算结果并得出你的结论。 (b)选定初始向量 和不同的右端向量 ,如取 。 将 的主对角线元素成倍增长若干次,非主对角元素不变,每次用Jacobi法计算,要求迭代误差满足 ,比较收敛速度,分析现象并得出你的结论。 2.2.算法描述 1)算法基本思路 1.Jacobi迭代法基本思路: 设有方程组 ,设A=(aij)nxn是可逆矩阵,对其移移项和变形后可得如下迭代格式: 然后,取定一个初始向量如x(0)=(0,0,…,0)T,带入上式开始迭代,然后如果达到迭代精度则停止迭代,即得到方程的数值解。 2.Guass-Seidel迭代法基本思路: Guass-Seidel迭代法是对Jacobi迭代法 (1)做出如下改进: 其余同Jacobi迭代法。 2)算法实现步骤: 1.用题目中已知的系数矩阵初始化矩阵A,并输入不同的右端向量b和不同的初始向量x(0); 2.选择迭代方式: 1.Jacobi迭代法2.Gauss—Seidel迭代法; 当选择1时,按照迭代格式 (1)开始迭代,并给定迭代误差,同时在每次迭代结束后输出此时的迭代结果,以便观测和检查; 当选择2时,按照迭代格式 (2)开始迭代,并给定迭代误差,同时在每次迭代结束后输出此时的迭代结果,以便观测和检查; 3.如果达到给定迭代误差则迭代停止,得到方程数值解,并输出。 2.3.程序代码 本程序是用C++语言在MicrosoftvisualC++编译环境下编译和运行的,由于代码太长,程序代码见附录2 2.4.计算结果 a)对问题a的求解结果如下: 1.当输入的右端向量b=[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]T;初始向量x(0)=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]T时,用Jacobi迭代法求解,得到如下结果: 请选择迭代方式: 1.Jacobi迭代法2.Gauss—Seidel迭代法 1↘(由键盘输入数字1,并回车) 第1次迭代结果: X[1]={0.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.3333330.333333} …… 第15次迭代结果: X[15]={0.4816320.5734090.6327920.652090.6609320.6642940.6657020.666260.6664830.666560.666560.6664830.666260.6657020.6642940.6609320.652090.6327920.5734090.481632} 第16次迭代结果: X[16]={0.4816340.5734120.6327970.6520960.6609390.6643020.665710.6662690.6664920.6665690.6665690.6664920.6662690.665710.6643020.6609390.6520960.6327970.5734120.481634} 迭代精度达到要求,迭代停止! 迭代结果为X[16] 2.同样输入右端向量b=[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]T;初始向量x(0)=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]T时,用G-S迭代法求解,得到如下结果: 请选择迭代方式: 1.Jacobi迭代法2.Gauss—Seidel迭代法 2↘(由键盘输入数字2,并回车) 第1次迭代结果: X[1]={0.3333330.3888890.4259260.4367280.4416150.443330.4440230.4442810.4443820.444420.4444350.4444410.4444430.4444440.4444440.4444440.4444440.4444440.4444440.444444} …… 第10次迭代结果: X[10]={0.4816290.5734060.6327910.6520910.6609340.6642980.6657080.6662680.6664920.6665710.6665730.6664970.6662740.6657160.6643080.6609440.6521010.6328010.5734150.481636} 第11次迭代结果: X[11]={0.4816340.5734120.6327980.6520980.6609410.6643050.6657150.6662740.6664980.6665760.6665760.6664990.6662760.6657180.6643080.6609450.6521010.6328010.5734150.481636} 迭代精度达到要求,迭代停止! 迭代结果为X[11] b)对问题b)的求解结果如下: 选取初始向量和右端向量分别为 ,将 的主对角线元素成倍增长若干次,非主对角元素不变,每次选用Jacobi迭代法计算,要求迭代误差满足, 将程序中表示倍数的M取1时,迭代次数为23次,能达到迭代误差精度要求; 当M取2时,迭代次数为13次,能达到迭代误差精度要求; 当M取3时,迭代次数为10次,能达到迭代误差精度要求; 当M取4时,迭代次数为9次,能达到迭代误差精度要求; 当M取5时,迭代次数为8次,能达到迭代误差精度要求; 当M取7时,迭代次数为8次,能达到迭代误差精度要求; 当M取8时,迭代次数为7次,能达到
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 矩阵与数值分析 矩阵 数值 分析 课程 实验 报告