数值分析上机实习报告文档格式.docx
- 文档编号:20612134
- 上传时间:2023-01-24
- 格式:DOCX
- 页数:19
- 大小:177.65KB
数值分析上机实习报告文档格式.docx
《数值分析上机实习报告文档格式.docx》由会员分享,可在线阅读,更多相关《数值分析上机实习报告文档格式.docx(19页珍藏版)》请在冰豆网上搜索。
2松弛因子对SOR法收敛速度的影响。
12
2.1SOR法12
2.2小结13
3Runge-Kutta4阶算法求解常微分方程14
3.1问题分析14
3.2实验结果及分析14
3.3结论17
4总结17
附录18
1Jacobi迭代法与Gauss-seidel迭代法的比较
2.用雅格比法与高斯-赛德尔迭代法解下列方程组Ax=b,研究其收敛性,上机验证理论分析是否正确,比较它们的收敛速度,观察右端项对迭代收敛有无影响。
(1)A行分别为A1=[6,2,-1],A2=[1,4,-2],A3=[-3,1,4];
b1=[-3,2,4]T,b2=[100,-200,345]T,
(2)A行分别为A1=[1,0,8,0.8],A2=[0.8,1,0.8],A3=[0.8,0.8,1];
b1=[3,2,1]T,b2=[5,0,-10]T,
(3)A行分别为A1=[1,3],A2=[-7,1];
b=[4,6]T,
1.1雅格比迭代法原理
雅可比迭代法的分量形式:
即:
考虑非奇异线性方程组:
AX=b
将A分解:
A=D-L-U
其中:
D=diag(diag(A))为对角矩阵;
L=-tril(A)为下三角矩阵;
U=-triu(A)为上三角矩阵;
B=D-1(L+U)为迭代矩阵。
那么方程组AX=b可写成:
X=D-1(L+U)X+D-1b
写成迭代形式为:
X(K)=D-1(L+U)X(K-1)+D-1b
式中上标K=1,2,3…n。
系数矩阵A=
62-1
14-2
-314
可以看出矩阵A是一个严格对角占优矩阵所以用雅格比法肯定收敛。
1.00000.80000.8000
0.80001.00000.8000
0.80000.80001.0000
迭代矩阵:
B=
0-0.8000-0.8000
-0.80000-0.8000
-0.8000-0.80000
由B可知迭代矩阵的谱半径为:
1.6>
1所以该迭代不收敛。
13
-71
0-3
70
4.5826>
1所以该迭代不收敛
1.2雅格比迭代法的运行结果
以下为Matlab输出结果:
(1)系数矩阵A=
右端列向量b=[-324]T
x=
3.0000
4.0000
-5.0000
n=
69
右端列向量b=[100-200345]T
119.3750
-125.8333
54.7917
77
37.0175
-3.5965
114.9123
29
(2)系数矩阵A=
右端列向量b=[321]T
waring:
迭代次数太多,可能不收敛!
1.0e+040*
1.5388
200
右端列向量b=[50-10]T
1.0e+041*
1.0942
(3)系数矩阵A=
右端列向量b=[46]T
1.0e+132*
2.7278
-0.9093
1.3高斯-赛德尔迭代法原理
高斯—赛德尔迭代公式的矩阵形式为:
x(k+1)=BGxk+fG
即:
AX=b
B=(D-L)-1U为迭代矩阵。
(D-L)X=UX+b
X=(D-L)-1UX+(D-L)-1b
写成迭代格式为:
X(K)=(D-L)-1UX(k-1)+(D-L)-1b
可以看出矩阵A是一个严格对角占优矩阵所以用高斯-赛德尔迭代法肯定收敛。
b1=[3,2,1]T,b2=[5,0,-10]T
00.6400-0.1600
00.12800.7680
0.7155<
1所以该迭代收敛。
0-21
21>
1.4高斯-赛德尔迭代法的运行结果
-0.9825
1.4035
-0.0877
22
5.7692
0.7692
-4.2308
45
32.6923
7.6923
-42.3077
52
warning:
迭代次数太多,可能不收敛!
1.0e+264*
-0.6135
-4.2945
1.5小结
1.不管是雅格比法还是高斯-赛德尔迭代法迭代是否收敛都只和B有关或者说只与系数矩阵A有关而与b无关。
只要矩阵的谱半径小于1迭代就收敛。
2.雅克比法不收敛的情况下高斯-赛德尔迭代法不一定发散。
另外,高斯-赛德尔迭代法收敛性也与方程组右端项无关。
3.对同一个方程组,高斯-赛德尔迭代法比雅可比法迭代次数少,收敛更快。
3.松弛因子对SOR法收敛速度的影响。
用SOR法求解方程组Ax=b,其中
要求程序中不存系数矩阵A,分别对不同的阶数取w=1.1,1.2,...,1.9进行迭代,记录近似解x(k)达到||x(k)-x(k-1)||<
10-6时所用的迭代次数k,观察松弛因子对收敛速度的影响,并观察当w0或w2会有什么影响?
2.1SOR法
多数情况下,Gauss-Seidel迭代法比Jacobi迭代法收敛速度快,分析Gauss-Seidel迭代公式和Jacobi迭代公式,可以发现:
Gauss-Seidel迭代法实际上,是每一个分量的第k次迭代值,由k-1次迭代的值再加上一个残量的倍数构成。
此时,得到新的迭代公式
Xik=xi(k-1)+wΔxik
当w=1时,为Gauss-Seidel迭代,w>1时,称作超松弛迭代法,而w<1,则为低松弛法。
调用自定义功能函数function[x,n]=SOR(A,b,x0,w,eps,M):
可以得到如下表所示在矩阵的阶数m不同和松弛因子w的情况下的结果。
W值
m值
4
6
8
1.1
9
11
1.2
12
13
1.3
14
15
16
1.4
18
19
1.5
23
24
1.6
31
33
1.7
43
44
1.8
68
69
71
1.9
144
143
147
2.2小结
从表格中的数据可以看出:
随着松弛因子w的不断增大,迭代次数在不断增加,收敛速度越来越慢。
由Kahan必要条件可知0<
w<
2,因此可以判断当w0或w2时,迭代不收敛。
当在程序中输入一个小于等于零或大于等于二的w值是结果显示错误。
因此可以得到松弛因子w在(12)范围内,w越大SOR法的收敛速度越慢;
而当w>
=2及=<
0时,SOR法对线性方程组的求解不收敛。
3Runge-Kutta4阶算法求解常微分方程
3.1问题分析
用Runge-Kutta4阶算法对初值问题y/=-20*y,y(0)=1按不同步长求解,用于观察稳定区间的作用,推荐两种步长h=0.1,0.2。
注:
此方程的精确解为:
y=e-20x
Runge-Kutta法的基本思想:
通过f(x,y)某些点函数值的适当线性组合替换Euler法中的
,可能使得方法的精确度更高。
在这种思想的指导下,标准的四阶Runge-Kutta的算法如下(h为求解步长):
算法稳定性的定义:
如果
是某方法第k步的近似值,
是其准确值,其绝对误差为
,即有
。
假定第k步之后的计算不再有舍入误差,只是由
引起的扰动
,都有
,则称此方法是绝对稳定的。
数值解法的稳定性一般与步长h以及f(x,y)有关。
3.2实验结果及分析
在MATLAB中按照标准四阶Runge-Kutta的算法编写程序代码,具体代码见附录六。
并运行程序。
解:
首先探讨不同步长下算法稳定的问题,取不同步长得到求解的微分方程的解同精确解的误差分别如表1、表2、表3、表4、表5。
表1h=0.0001下的Runge-Kutta四阶算法
表2h=0.001下的Runge-Kutta四阶算法
表3h=0.01下的Runge-Kutta四阶算法
表4h=0.1下的Runge-Kutta四阶算法
表5h=0.2下的Runge-Kutta四阶算法
从上面得表中可以看到,在h=0.0001,0.001,0.01,0.1时,绝对误差在逐步变小,因此在这种情况下认为此方法是稳定的。
当h=0.2时,绝对误差逐渐扩大,因此在h=0.2时,算法不稳定。
现在来探讨不同h下面的精度问题,取不同的步长求y(0.1)的值得表6。
表6不同步长h下的y(0.1)值
从表6中可以看到步长越小,标准四阶Runge-Kutta算法得到的值精度越高。
计算所需的步骤也越长。
3.3结论
(1)用标准的四阶Runge-Kutta算法计算常微分方程时,不同的步长对算法的稳定性有影响。
不够稳定的步长下面的计算,误差会越来越大,结果失真严重。
(2)一般情况下,步长越小,标准的四阶Runge-Kutta算法的稳定性越高,精度也越高。
4总结
本次实验使用MATLAB软件编程解决线性方程组的Jacobi迭代法和Gauss-Seidel迭代法以及线性方程组求解的SOR迭代法程序,并对得到的结果进行分析。
通过对MATLAB软件的使用,进一步熟悉了MATLAB程序编写过程中的一些注意事项及规则,对以后MATLAB软件的应用打下好的基础。
通过实验使得我对Jacobi迭代法、Gauss-Seidel迭代法及SOR法解线性方程组有了更加深刻的认识,同时还了解到这些方法在使用过程中可能出现的一些现象和问题,可以更加熟练的应用这些方法来解决工程中的遇到数值问题。
1.不管是雅格比法还是高斯-赛德尔迭代法迭代是否收敛都只和B由管或者说只与系数矩阵A有关而与b无关。
2.雅可比法收敛条件情况和高斯-赛德尔迭代法的收敛情况没有必然联系。
4.松弛法迭代过程中,松弛因子w越大SOR法的迭代次数越多,收敛速度越慢;
附录
Jacobi迭代法的迭代程序:
function[x,n]=jacobi(A,b,x0,esp,varargin)
ifnargin==3
esp=1.0e-6;
M=200
elseifnargin<
3
error
return
elseifnargin==5
M=varargin{1};
end
D=diag(diag(A));
L=-tril(A,-1);
U=-triu(A,1);
B=D\(L+U);
f=D\b;
x=B*x0+f;
n=1;
whilenorm(x-x0)>
=esp
x0=x;
x=B*x0+f;
n=n+1;
if(n>
=M)
disp('
'
);
return;
end
Gauss-seidel迭代法的迭代程序:
function[x,n]=gaussseidel(A,b,x0,esp,M)
ifnargin==3
esp=1.0e-6;
M=200;
elseifnargin==4
elseifnargin<
error
return;
G=(D-L)\U;
f=(D-L)\b;
x=G*x0+f;
x0=x;
x=G*x0+f;
if(n>
disp('
SOR法的迭代程序:
function[x,n]=SOR(A,b,x0,w,eps,M)
eps=1.0e-6;
B=inv(D-L*w)*((1-w)*D+w*U);
f=w*inv((D-L*w))*b;
=eps
矩阵A和b的输入程序:
clearall
symsm
m=input('
pleaseinputm'
)
A=-4*eye(m);
fori=1:
m-1
A(i+1,i)=1;
A(i,i+1)=1;
A
b(1,1)=-3;
b(m,1)=-3;
fori=2:
m-1;
b(i,1)=-2;
b
m;
X(i,1)=0;
X;
x0=X;
[x,n]=SOR(A,b,x0,1.8)
Runge-Kutta4阶算法
h=[0.0001,0.001,0.01,0.1]
k=1
forxk=0.1:
0.1:
1%取不同的x进行计算
1:
4%取不同的步长计算每个x
hj=h(i)
yk=1
forxkp=hj:
hj:
xk%四阶标准Runge-Kutta迭代算法
k1=hj*(-20*yk);
formatlong
k1;
k2=hj*(-20*(yk+k1*0.5));
k3=hj*(-20*(yk+k2*0.5));
k4=hj*(-20*(yk+k3));
ykm=yk+(k1+2*k2+2*k3+k4)*(0.16666666666667);
yk=ykm;
yeal=exp(-20*xkp);
%求出精确值
ck=abs(yk-yeal);
%精确值和计算值之间的绝对误差
ykeal(k)=yeal;
ckk(k)=ck;
yxk(k)=yk;
xkk(k)=xk;
k=k+1;
D=[yxk'
ykeal'
ckk'
];
%将迭代得到的结果及误差存到矩阵D中
10%将同一个步长下不同X得到的结果集中到一个矩阵
forj=1:
Ah1(i,j)=D(4*i-3,j);
Ah2(i,j)=D(4*i-2,j);
Ah3(i,j)=D(4*i-1,j);
Ah4(i,j)=D(4*i,j);
Axk(i)=xkk(4*i-3);
Ch1=[Axk'
Ah1]
Ch2=[Axk'
Ah2]
Ch3=[Axk'
Ah3]
Ch4=[Axk'
Ah4]
Cheng=[Ch1(1,[234]);
Ch2(1,[234]);
Ch3(1,[234]);
Ch4(1,[234])]
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数值 分析 上机 实习 报告