西南交通大学研究生数值分析上机实习Word文件下载.doc
《西南交通大学研究生数值分析上机实习Word文件下载.doc》由会员分享，可在线阅读，更多相关《西南交通大学研究生数值分析上机实习Word文件下载.doc（10页珍藏版）》请在冰点文库上搜索。
（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
迭代次数n
6
收敛失败
Steffensen加速法
9
7
求sinx=0的根。
1.4
1.6
1.8
3
3.1416
31.4159
6.2832
5
8
-3.14159
25.1327
6.28319
3.14159
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矩阵带入
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
36.36363656
36.3636365
-2.070706833
-2.070707008
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
迭代失败
5.769232169
0.769228821
-4.230768792
44
32.69230969
7.692307507
-42.30769376
51
2.2结果分析
Jacobi迭代迭代矩阵
特征值
-0.333333333
0.166666667
-0.542663188
-0.25
0.5
0.271331593969688+0.370849252830665i
0.75
0.271331593969688-0.370849252830665i
-0.8
-1.6
0.8
-3
0+4.58257569495584i
0-4.58257569495584i
Causs-Seidel迭代矩阵
0.083333333
0.458333333
0.046875+0.350432210812591i
-0.270833333
0.010416667
0.046875-0.350432210812591i
0.64
-0.16
0.704+0.128i
0.128
0.768
0.704-0.128i
-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.2
0.3
0.4
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
0.6
0.7
0.9
6.14E-06
8.32E-07
1.13E-07
1.52E-08
2.06E-09
0.00137174
0.0004572
0.00015242
5.08E-05
1.69E-05
0.0013656
0.0004564
0.0001523
当h=0.2时
0.01831564
0.0003355
25
125
625
3125
4.98168436
24.999665
124.999994
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这一软件的基本功能与应用。
通过上机中我深刻明白，在所面临的实际问题，大多烦琐复杂的计算必须基于一定的数值算法，用计算机编程来实现。
而人的手工计算是远远不可行的。
这就要求我们寻找一种稳定性好，收敛速度快，又满足要求精度的数值算法。
所以上机实习和这门课程使我了解了一些普遍应用的基本算法，使我在以后工程计算中能够进行广泛的应用。
通过这次上机练习，让我对数值分析所介绍的迭代求解方法及其理论有了更深层次的理解，了解了各种方法之间的优缺点，并且认识到了自己在以前的学习中所存在的问题，及时的修补了自己知识上的漏洞。
同时也提高了我在编写程序上的熟练程度，所以，我认为这次上机实习是非常有收获的，给予我学习数值分析的帮助也是非常大的。
附录
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
function[x,n]=Newton1（xi）%牛顿法求解
sin（x）'
error
function[x,n]=Steffensen（xi）%Steffensen加速法
0.000001）
y=x-f（x）/f2（x）
z=y-f（y）/f2（y）
x=x-（y-x）*（y-x）/（z-2*y+x）
function[x,n]=Steffensen1（xi）%Steffensen加速法
x-sin（x）/cos（x）'
y=f（x）
z=f（y）
题目二
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
x0=x
B=I-inv（D）*A
x=B*x+inv（D）*b
n=n+1
2、Causs-Seidel迭代法
function[x,n]=gaussseidel（A,b）
M=triu（A）
N=tril（A）
L=A-M
U=A-N
200
C=D+L
G=-inv（C）*U
x=G*x+inv（C）*b
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
升级会员 