数学实验 5线性代数方程组的数值解法概要.docx
- 文档编号:2037646
- 上传时间:2022-10-26
- 格式:DOCX
- 页数:16
- 大小:98.61KB
数学实验 5线性代数方程组的数值解法概要.docx
《数学实验 5线性代数方程组的数值解法概要.docx》由会员分享,可在线阅读,更多相关《数学实验 5线性代数方程组的数值解法概要.docx(16页珍藏版)》请在冰豆网上搜索。
数学实验5线性代数方程组的数值解法概要
实验5:
线性代数方程组的数值解法
习题3:
已知方程组,其中,定义为:
试通过迭代法求解此方程组,认识迭代法收敛的含义以及迭代初值和方程组系数矩阵性质对收敛速度的影响。
实验要求:
(1)选取不同的初始向量x0和不同的方程组右端向量b,给定迭代误差要求,用雅可比迭代法和高斯-赛德尔迭代法计算,观测得到的迭代向量序列是否均收敛?
若收敛,记录迭代次数,分析计算结果并得出结论;
(2)取定右端向量b和初始向量x0,将A的主对角线元素成倍的增长若干次,非主对角元素不变,每次用雅可比迭代法计算,要求迭代误差满足,比较收敛速度,分析现象并得出结论。
1、程序设计(可直接粘贴运行)
1)Jacobi迭代法
functiony=jacobi(a,b,x0,e,m)
%定义jacobi函数,其中:
a,b为线性方程组中的矩阵和右端向量;x0为初始值;
%e和m分别为人为设定的精度和预计迭代次数;运行结果y为迭代的结果和所有中间值组成的
%矩阵
y=0;%对y初始化
d=diag(diag(a));%按雅可比迭代标准形形式取主对角元素作为矩阵D
u=-triu(a,1);%取上三角矩阵u
l=-tril(a,-1);%取下三角矩阵l
bj=d^-1*(l+u);
fj=d^-1*b;
x=[x0,zeros(20,m-1)];%初始化x,其中x1=x0,即初始值
fork=1:
m%人为规定迭代次数,防止不收敛迭代导致死循环
x(:
k+1)=bj*x(:
k)+fj;%jacobi迭代
ifnorm(x(:
k+1)-x(:
k),inf) %判断迭代后是否满足迭代中止条件: y=x(: 1: k+1);%赋给y所有中间值和迭代结果 sizej=k;%若去掉;号,则输出迭代次数 break%并结束迭代 end%若不成立,继续迭代 end %以下部分为验证迭代公式收敛的方法,仅需运行一次即可,因为收敛性完全由A矩阵决定,而A %在本题是固定不变的;通过判断中B的谱半径或范数大小(B在jacobi迭代法中%为矩阵bj),即可得知收敛性: e=eig(bj) %输出全部特征值,若,则收敛 n1=norm(bj,1);%计算1-范数 n2=norm(bj);%计算2-范数 nn=norm(bj,inf);%计算-范数 q=min([n1n2nn]) %由于谱半径不超过人以一种范数,所以只要范数的最小值q<1,也可判断迭代法收敛 2)Gauss迭代法: 与Jacobi程序结构相同,不再注释 functiony=gauss(a,b,x0,e,m) y=0; d=diag(diag(a)); u=-triu(a,1); l=-tril(a,-1); x=[x0,zeros(20,m-1)]; bgs=(d-l)^-1*u; fgs=(d-l)^-1*b; fork=1: m x(: k+1)=bgs*x(: k)+fgs; ifnorm(x(: k+1)-x(: k),inf) y=x(: 1: k+1); sizeg=k; break end end e=eig(bgs) n1=norm(bgs,1); n2=norm(bgs); nn=norm(bgs,inf); min([n1n2nn]) 3)操作函数: %构造矩阵A n=20; a1=sparse(1: n,1: n,3,n,n);%按稀疏矩阵的输入法构造,比较方便 a2=sparse(1: n-1,2: n,-0.5,n,n); a3=sparse(1: n-2,3: n,-0.25,n,n); a=a1+a2+a3+a2'+a3'; a=full(a);%还原为满矩阵 %通过给定不同的初始向量x0或者右端项b,以及规定不同的误差要求,进行jacobi和gauss %迭代,得到的结果y1、y2位两种迭代的次数,同时输出迭代结果,便于分析 b= x0= e= m= y1=jacobi(a,b,x0,e,m); y2=gauss(a,b,x0,e,m); 4)改变矩阵A: 先对jacobi函数作一定修正,方便分析,命名为jacobi2,如下: functiony=jacobi2(a,b,x0,e,m) d=diag(diag(a)); u=-triu(a,1); l=-tril(a,-1); bj=d^-1*(l+u); fj=d^-1*b; x=[x0,zeros(20,m-1)]; n1=norm(bj,1);%计算范数 n2=norm(bj); nn=norm(bj,inf); q=min([n1n2nn]); y (1)=q;%输出结果1: 范数的最小值,判断收敛速度的方法 fork=1: m x(: k+1)=bj*x(: k)+fj; ifnorm(x(: k+1)-x(: k),inf) y (2)=k;%输出结果2: 迭代次数 break end end %改变A矩阵主对角元素的值,比较jacobi迭代的收敛速度,即迭代误差满足%时的迭代次数 b=(1: 20)';%取定b x0=20*ones(20,1);%取定x0 e=10^-5;%取定精度 r=20;%设置主对角元素扩大的最大倍数 y=0; m=50; fork=1: r; a1=k*sparse(1: n,1: n,3,n,n);%只需更改A的主对角元素 a=a1+a2+a3+a2'+a3';%重新构造A a=full(a); y(k,1: 3)=[k,jacobi2(a,b,x0,e,m)]; end y%输出对角线扩大倍数\最小范数\迭代次数 2、运行结果分析 1)不同初值(x0)和右端项(b)条件下的解的情况 表一: 收敛性判断 1-范数 2-范数 -范数 Jacobi 0.0163 0.0167 0.0163 0.0167 Gauss 0.0008 0.0084 0.0084 0.0084 可以看到,矩阵A无论是谱半径或是任意范数的值都小于1,可知在A不变的情况下,Jacobi和Gauss法必然收敛。 2)b取不同的值,x0=20*ones(20,1),e=10^-5,m=50条件下的情况对比 迭代次数 B=[1: 20]’ B=[10: 10: 200]’ B=[20: -1: 1]’ B=20*ones(20,1) B=2000*ones(20,1) Jacobi 24 26 24 23 30 Gauss 16 17 16 15 20 根据1)分析的结果,可以证明无论b取任何值,采用两种方法迭代均收敛,但b的值的变化会影响迭代的次数;且Gauss迭代法总是比Jacobi迭代法收敛速度更快。 下表列出的是B=[1: 20]’情况下部分结算结果,可以很明显的看到两种迭代法的收敛速度不同: Jacobi K=1 K=2 K=3 K=4 K=5 K=6 … K=22 K=23 K=24 标准值 X1 5.3333 2.7500 1.5394 1.0874 0.8858 0.7982 … 0.7247 0.7247 0.7247 0.7247 X2 9.0000 4.3333 2.6644 1.9248 1.6082 1.4652 … 1.3444 1.3444 1.3444 1.3444 X3 11.0000 5.8056 3.7199 2.7801 2.3626 2.1717 … 2.0072 2.0072 2.0072 2.0072 Gauss K=1 K=2 K=3 K=4 K=5 K=6 … K=14 K=15 K=16 标准值 X1 5.3333 2.0540 1.1354 0.8539 0.7657 0.7377 … 0.7247 0.7247 0.7247 0.7247 X2 6.5556 2.9432 1.8459 1.5033 1.3949 1.3605 … 1.3444 1.3444 1.3444 1.3444 X3 7.5370 3.7386 2.5553 2.1815 2.0627 2.0249 … 2.0072 2.0072 2.0072 2.0072 由于精度问题,在迭代的最后几次中从显示的数位已经不能看出标准值与计算值得差别,但是若采用long显示设定,就可以看到更多位小数的显示,其结果符合最初设定的精度e,数据繁琐,略。 3)x0取不同的值,b=20*ones(20,1),e=10^-5,m=50条件下的情况对比 迭代次数 X0=[1: 20]’ X0=[10: 10: 200]’ X0=[20: -1: 1]’ X0=20*ones(20,1) X0=2000*ones(20,1) Jacobi 22 27 22 23 31 Gauss 15 18 15 15 20 根据1)分析的结果,可以证明无论X0取任何值,采用两种方法迭代均收敛,但x0的值的变化会影响迭代的次数;且Gauss迭代法总是比Jacobi迭代法收敛速度更快。 下表列出的是x0=[1: 20]’情况下部分结算结果,可以很明显的看到两种迭代法的收敛速 度不同: Jacobi K=1 K=2 K=3 K=4 K=5 K=6 … K=20 K=21 K=22 标准值 X1 7.2500 8.6250 9.2228 9.4536 9.5541 9.5974 … 9.6327 9.6327 9.6327 9.6327 X2 7.6667 9.9583 10.814 11.183 11.341 11.411 … 11.4683 11.4683 11.4683 11.4683 X3 8.1667 10.757 11.816 12.282 12.487 12.579 … 12.6560 12.6560 12.6560 12.6560 Gauss K=1 K=2 K=3 K=4 K=5 K=6 … K=13 K=14 K=15 标准值 X1 7.2500 8.9352 9.4267 9.5720 9.6150 9.6276 … 9.6327 9.6327 9.6327 9.6327 X2 8.7083 10.654 11.228 11.398 11.448 11.463 … 11.4683 11.4683 11.4683 11.4683 X3 9.8056 11.813 12.408 12.584 12.636 12.650 … 12.6560 12.6560 12.6560 12.6560
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数学实验 5线性代数方程组的数值解法概要 数学 实验 线性代数 方程组 数值 解法 概要