西南交通大学研究生数值分析上机实习.doc
- 文档编号:655355
- 上传时间:2023-04-29
- 格式:DOC
- 页数:10
- 大小:173KB
西南交通大学研究生数值分析上机实习.doc
《西南交通大学研究生数值分析上机实习.doc》由会员分享,可在线阅读,更多相关《西南交通大学研究生数值分析上机实习.doc(10页珍藏版)》请在冰点文库上搜索。
目录
解题:
1
题目一:
1
1.1计算结果 1
1.2结果分析 1
题目二:
2
2.1计算结果 2
2.2结果分析 3
题目三:
4
3.1计算结果 4
3.2结果分析 5
总结 5
附录 6
Matlab程序:
6
题目一:
6
第一问Newton法:
6
第二问Newton法:
6
第一问Steffensen加速法:
7
第二问Steffensen加速法:
7
题目二 8
1、Jacobi迭代法 8
2、Causs-Seidel迭代法 8
题目三:
9
解题:
题目一:
分别用牛顿法,及基于牛顿算法下的Steffensen加速法
(1)求ln(x+sinx)=0的根。
初值x0分别取0.1,1,1.5,2,4进行计算。
(2)求sinx=0的根。
初值x0分别取1,1.4,1.6,1.8,3进行计算。
分析其中遇到的现象与问题。
1.1计算结果
求ln(x+sinx)=0的根,可变行为求解x-sinx-1=0的根。
结果 :
牛顿法
初值x0
0.1
1
1.5
2
4
收敛解x
0.51097
0.51097
0.51097
迭代次数n
6
6
收敛失败
6
收敛失败
Steffensen加速法
收敛解x
0.51097
0.51097
0.51097
0.51097
0.51097
迭代次数n
9
4
6
9
7
求sinx=0的根。
初值x0分别取1,1.4,1.6,1.8,3进行计算。
牛顿法
初值x0
1
1.4
1.6
1.8
3
收敛解x
0
3.1416
31.4159
6.2832
3.1416
迭代次数n
5
7
8
4
3
Steffensen加速法
收敛解x
0
-3.14159
25.1327
6.28319
3.14159
迭代次数n
3
4
6
3
3
1.2结果分析
从结果对比我们可发现牛顿—Steffensen加速法比牛顿法要收敛的快,牛顿法对于初值的选取特别重要,比如第
(1)问中的初值为4的情况,100次内没有迭代出来收敛解,而用Steffensen加速法,7次迭代可得;在第
(2)问中的初值为1.6的情况,收敛解得31.4159,分析其原因应该是,x0=1.6,;迭代式在迭代过程中会出现分母趋近于0,程序自动停止迭代的情况,此时得到的x往往非常大,而在第一问中我们如果转化为用x+sinx=1,则可以收敛到结果。
题目二:
用雅格比法与高斯-赛德尔迭代法解下列方程组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
2.1计算结果
初值均为0矩阵带入
(1)A行分别为A1=[6,2,-1],A2=[1,4,-2],A3=[-3,1,4];b1=[-3,2,4]T,b2=[100,-200,345]T
Jacobi迭代
Causs-Seidel迭代
b1
x1
-0.727272981
-0.727272584
x2
0.808080524
0.808080904
x3
0.25252547
0.252525336
迭代次数
23
15
b2
x1
36.36363656
36.3636365
x2
-2.070706833
-2.070707008
x3
114.0404039
114.0404041
迭代次数
30
20
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
Jacobi迭代
Causs-Seidel迭代
b1
x1
迭代失败
5.769232169
x2
0.769228821
x3
-4.230768792
迭代次数
44
b2
x1
迭代失败
32.69230969
x2
7.692307507
x3
-42.30769376
迭代次数
51
(3)A行分别为A1=[1,3],A2=[-7,1];b=[4,6]T
Jacobi迭代
Causs-Seidel迭代
b1
x1
迭代失败
迭代失败
x2
x3
迭代次数
2.2结果分析
Jacobi迭代迭代矩阵
特征值
1
0
-0.333333333
0.166666667
-0.542663188
-0.25
0
0.5
0.271331593969688+0.370849252830665i
0.75
-0.25
0
0.271331593969688-0.370849252830665i
2
0
-0.8
-0.8
-1.6
-0.8
0
-0.8
0.8
-0.8
-0.8
0
0.8
3
0
-3
0+4.58257569495584i
7
0
0-4.58257569495584i
Causs-Seidel迭代矩阵
特征值
1
0
-0.333333333
0.166666667
0
0
0.083333333
0.458333333
0.046875+0.350432210812591i
0
-0.270833333
0.010416667
0.046875-0.350432210812591i
2
0
-0.8
-0.8
0
0
0.64
-0.16
0.704+0.128i
0
0.128
0.768
0.704-0.128i
3
0
-3
0
0
-21
-21
第一小题的经计算谱半径为小于1,故方程组雅可比迭代收敛。
而经计算高斯-赛德尔迭代的谱半径为0.3535小于1,故原方程组高斯-赛德尔迭代矩阵收敛。
第二小题的经计算谱半径为1.6大于1,故方程组雅可比迭代发散。
而经计算高斯-赛德尔迭代的谱半径为0.7155小于1,故原方程组高斯-赛德尔迭代矩阵收敛。
第三小题谱半径为4.5826大于1,故方程组雅可比迭代发散。
而高斯-赛德尔迭代的谱半径为21大于1,故原方程组高斯-赛德尔迭代矩阵发散。
从2.1中的结果列表中可以看到,Causs-Seidel迭代法比Jacobi迭代法收敛速度要快很多。
题目三:
用Runge-Kutta4阶算法对初值问题y/=-20*y,y(0)=1按不同步长求解,用于观察稳定区间的作用,推荐两种步长h=0.1,0.2。
注:
此方程的精确解为:
y=e-20x
3.1计算结果
当h=0.1时
Xk
0.1
0.2
0.3
0.4
0.5
y(xk)
0.13533528
0.0183156
0.00247875
0.00033546
4.54E-05
yk
0.33333333
0.1111111
0.03703704
0.01234568
0.00411523
0.19799805
0.0927955
0.03455828
0.01201022
0.00406983
Xk
0.6
0.7
0.8
0.9
1
y(xk)
6.14E-06
8.32E-07
1.13E-07
1.52E-08
2.06E-09
yk
0.00137174
0.0004572
0.00015242
5.08E-05
1.69E-05
0.0013656
0.0004564
0.0001523
5.08E-05
1.69E-05
当h=0.2时
Xk
0.2
0.4
0.6
0.8
1
y(xk)
0.01831564
0.0003355
6.14E-06
1.13E-07
2.06E-09
yk
5
25
125
625
3125
4.98168436
24.999665
124.999994
625
3125
3.2结果分析
h=0.2时,h=-4,而Runge-Kutta4阶算法的绝对稳定区间是[-2.78,0],故h=0.2时计算不稳定;而h=0.1时,h=-2,在绝对稳定区间内,计算稳定,结果可靠。
总结
本次上机实习使用matlab编写了牛顿法、牛顿-Steffensen法方程求解的程序和雅格比法、高斯-赛德尔迭代法求解方程组的程序、Runge-Kutta4阶算法。
深刻了解了matlab这一软件的基本功能与应用。
通过上机中我深刻明白,在所面临的实际问题,大多烦琐复杂的计算必须基于一定的数值算法,用计算机编程来实现。
而人的手工计算是远远不可行的。
这就要求我们寻找一种稳定性好,收敛速度快,又满足要求精度的数值算法。
所以上机实习和这门课程使我了解了一些普遍应用的基本算法,使我在以后工程计算中能够进行广泛的应用。
通过这次上机练习,让我对数值分析所介绍的迭代求解方法及其理论有了更深层次的理解,了解了各种方法之间的优缺点,并且认识到了自己在以前的学习中所存在的问题,及时的修补了自己知识上的漏洞。
同时也提高了我在编写程序上的熟练程度,所以,我认为这次上机实习是非常有收获的,给予我学习数值分析的帮助也是非常大的。
附录
Matlab程序:
题目一:
第一问Newton法:
function[x,n]=Newton(xi)%牛顿法求解
x=xi
xi=xi-0.1
n=0
y=sym('log(x+sin(x))')
f=inline(y)
f1=diff(y)
f2=inline(f1)
while(abs(x-xi)>=0.000001)
xi=x
x=x-f(x)/f2(x)
n=n+1
ifn>100
error
break
end
end
end
第二问Newton法:
function[x,n]=Newton1(xi)%牛顿法求解
x=xi
xi=xi-0.1
n=0
y=sym('sin(x)')
f=inline(y)
f1=diff(y)
f2=inline(f1)
while(abs(x-xi)>=0.000001)
xi=x
x=x-f(x)/f2(x)
n=n+1
ifn>100
error
break
end
end
end
第一问Steffensen加速法:
function[x,n]=Steffensen(xi)%Steffensen加速法
x=xi
xi=xi-0.1
n=0
y=sym('log(x+sin(x))')
f=inline(y)
f1=diff(y)
f2=inline(f1)
while(abs(x-xi)>0.000001)
xi=x
y=x-f(x)/f2(x)
z=y-f(y)/f2(y)
x=x-(y-x)*(y-x)/(z-2*y+x)
n=n+1
ifn>100
error
break
end
end
end
第二问Steffensen加速法:
function[x,n]=Steffensen1(xi)%Steffensen加速法
x=xi
xi=xi-0.1
n=0
y=sym('x-sin(x)/cos(x)')
f=inline(y)
while(abs(x-xi)>0.000001)
xi=x
y=f(x)
z=f(y)
x=x-(y-x)*(y-x)/(z-2*y+x)
n=n+1
ifn>100
error
break
end
end
end
题目二
1、Jacobi迭代法
function[x,n]=jacobi(A,b)
[m,n]=size(A)
x=zeros(n,1)
e=ones(n,1)
x0=x+e
I=eye(n)
a=diag(A)
D=diag(a)
n=0
whilenorm(x0-x,inf)>0.000001
ifn>200;
error
break
end
x0=x
B=I-inv(D)*A
x=B*x+inv(D)*b
n=n+1
end
2、Causs-Seidel迭代法
function[x,n]=gaussseidel(A,b)
[m,n]=size(A)
x=zeros(n,1)
e=ones(n,1)
x0=x+e
a=diag(A)
D=diag(a)
M=triu(A)
N=tril(A)
L=A-M
U=A-N
n=0
whilenorm(x0-x,inf)>0.000001
ifn>200
error
break
end
x0=x
C=D+L
G=-inv(C)*U
x=G*x+inv(C)*b
n=n+1
end
题目三:
functionA=Runge(x,h)
a=x/h
y=1
x=0
fori=1:
a
k1=h*-20*y
k2=h*-20*(y+1/2*k1)
k3=h*-20*(y+1/2*k2)
k4=h*-20*(y+k3)
y=y+(k1+2*k2+2*k3+k4)/6
A(i)=y
x=x+h
f=exp(-20*x)
B(i)=f
d=abs(f-y)
C(i)=d
end
9
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 西南交通大学 研究生 数值 分析 上机 实习