1、西南交大数值分析上机报告 数值分析实验报告 班 级:桥梁二班姓 名: 学 号:电 话:2013年12月序言本次数值分析实验采用MATLAB语言。MATLAB是当今最流行的科学计算语言之一,与FORTRAN、C语言一样,在科学计算领域被广泛的应用。本次实验选择MATLAB语言,首先因为MATLAB具有强大的科学计算及数据处理能力,尤其在数值计算与优化方面,具有无可比拟的优越性。MATLAB拥有600多个工程中要用到的数学运算函数,可以方便地实现用户所需的各种计算功能。函数中所使用的算法都是科研和工程计算中的最新研究成果,而且经过了各种优化及容错处理,因此使用起来稳定性和可靠性非常高,在通常情况下
2、,可以用它来代替底层编程语言,如C和C+等。在计算要求相同的情况下,使用MATLAB的编程工作量会大大减少。MATLAB函数所能解决的问题包括矩阵运算、多维数组操作(阵列运算)、复数的各种运算、三角函数和其他初等数学函数运算、非线性方程求根、线性方程组的求解、微分方程及偏微分方程组的求解、符号运算、傅立叶变换和数据的统计分析、工程中的优化问题、稀疏矩阵运算、建模和动态仿真等等。其次MATLAB具有出色的图形处理功能。在科学计算中,往往需要用各种图形把数值计算的结果形象地表现出来,以帮助人们更好地理解、认识和发现其中的科学规律。MATLAB不仅提供数值计算功能和符号运算功能,而且自诞生之日起就具
3、有方便的数据可视化功能,使计算结果的可视化要求得到充分满足。MATLAB在二维曲线和三维曲面的绘制和处理等方面的功能比一般数据可视化软件更加完善。在MATLAB中可以对图形对象的属性进行编辑,以提高输出图形的效果。对一些特殊的可视化要求,例如图形动画等,MATLAB也有相应的功能函数,保证了用户不同层次的要求。MATLAB最突出的优点是程序语言简单易用。早期用于科学计算的计算机语言,由于计算机内存容量和运算速度的限制等原因,常常要定义常量、变量、向量和矩阵等的不同的数据类型,结果导致编程过于复杂化。和这些语言不一样,MATLAB语言对他们进行了高度抽象,实现了数据类型的高度统一,即常量、变量、
4、向量和矩阵等都具有相同的数据类型。因此,用户不需要事先分别定义常量、变量、向量和矩阵等的数据类型就可以直接使用他们。MATLAB语言是一种“数学形式的语言”,它的操作和功能函数指令就是用平时计算机和数学书上的英文单词和符号来表达的,比BASIC、FORTRAN和C等语言更接近于人们书写的数学计算公式、更接近于人们进行科学计算的思维方式,因此,MATLAB语言简单自然,学习和使用更容易。MATLAB程序文件是一个纯文本文件,扩展名为.m,用任何字处理软件都可以对它进行编辑。MATLAB本身就像一个解释系统,对其中的函数程序的执行以一种解释执行的方式进行,程序不必经过编译就可以直接运行,而且能够及
5、时报告出现的错误,进行出错原因分析。因此,程序调试容易、编程效率高。目录一、 解线性方程组 11. 计算与分析 1二、 微分方程初值问题 21. 计算与分析 2三、 非线性方程根的求解 31. 分析与计算 3四、 数值积分 41. 计算与分析 4总结 6附录 7一、 Gauss消元法程序 7二、 Gauss-Seidel迭代法程序 7三、 Runge-Kutta法程序 8四、 Newton迭代法程序 9五、 复合梯形公式程序 9六、 复合Simpson公式程序 10七、 Romberg算法程序 10一、 解线性方程组写出对一般的线性方程组通用的Gauss消元, Gauss-Seidel迭代程序
6、。并以下面的线性方程组为例进行计算,讨论所得到的计算结果是否与理论一致。(1)(2)(3)1. 计算与分析(1)利用选列主元的Gauss消元法编写MATLAB程序,解上述5个方程组。结果如表1.1所示:表1.1 Gauss消元法计算结果x(1)x(2)x(3)x(4)x(5)x1-0.7272727336.363636365.7692307732.69230769-0.63636364x20.80808081-2.070707070.769230777.692307691.54545455x30.25252525114.04040404-4.23076923-42.30769231Gauss消
7、元法为直接解法,从消元的过程可知其结果为精确解,误差来源于机器误差,即计算机内部的存储精度。MATLAB用二进制双精度浮点格式来表达一个数,其计算结果至少可以有12位十进制有效数字的可信度。(2)利用Gauss-Seidel法编写MATLAB程序,用|x(k)-x(k-1)|1,故Gauss-Seidel迭代发散。将方程组的两个方程顺序对换一下,变为,则其系数矩阵为行对角占优阵,Gauss-Seidel迭代收敛。Gauss-Seidel迭代法计算结果如表1.2所示:表1.2 Gauss-Seidel迭代计算结果x1x2x3迭代次数x101 1 1 15x1k-0.72727289 0.8080
8、8085 0.25252512 x2011120x2k36.36363650 -2.07070701 114.04040413 x3011145x3k5.76923210 0.76922890 -4.23076879 x4011151x4k32.69230977 7.69230745 -42.30769378 x501116x5k-0.63636362 1.54545454 二、 微分方程初值问题给定初值问题(精确解 为),用Runge-Kutta 4阶算法按步长求解,分析其中遇到的现象及问题。1. 计算与分析编写4级4阶的Runge-Kutta法MATLAB程序。取步长0.2,计算结果列于表
9、2.1中。表中y(xk)是当x=xk时的准确值,yk是由Runge-Kutta法 得到的当x=xk时的值。表2.1 Runge-Kutta法数值计算结果步长0.20110.251.83156E-020.4253.35463E-040.61256.14421E-060.86251.12535E-071.0 31252.06115E-091.2156253.77513E-111.4781256.91440E-131.63906251.26642E-141.819531252.31952E-162.0 97656254.24835E-18从表中可以看出步长取0.2时Runge-Kutta法是发散的。
10、 4级4阶Runge-Kutta法的绝对稳定区间为(-2.78,0),对于实验方程,h=-20h,当h=0.2时h=-4,故Runge-Kutta法无法收敛。步长取0.1,0.05时,h分别为-2,-1,位于绝对稳定区间,预计Runge-Kutta法收敛。表2.2列出了步长为0.1和0.05时Runge-Kutta法的计算结果和误差表2.2 Runge-Kutta法数值计算结果步长0.1步长0.05011001100.13.33333E-011.35335E-01-1.97998E-010.05 3.75000E-013.67879E-01-7.12056E-030.21.11111E-011
11、.83156E-02-9.27955E-020.10 1.40625E-011.35335E-01-5.28972E-030.33.70370E-022.47875E-03-3.45583E-020.15 5.27344E-024.97871E-02-2.94731E-030.41.23457E-023.35463E-04-1.20102E-020.20 1.97754E-021.83156E-02-1.45975E-030.54.11523E-034.53999E-05-4.06983E-030.25 7.41577E-036.73795E-03-6.77824E-040.61.37174E
12、-036.14421E-06-1.36560E-030.30 2.78091E-032.47875E-03-3.02162E-040.74.57247E-048.31529E-07-4.56416E-040.35 1.04284E-039.11882E-04-1.30961E-040.81.52416E-041.12535E-07-1.52303E-040.40 3.91066E-043.35463E-04-5.56034E-050.95.08053E-051.52300E-08-5.07900E-050.45 1.46650E-041.23410E-04-2.32400E-051.0 1.6
13、9351E-052.06115E-09-1.69330E-050.50 5.49937E-054.53999E-05-9.59374E-061.15.64503E-062.78947E-10-5.64475E-060.55 2.06226E-051.67017E-05-3.92092E-061.21.88168E-063.77513E-11-1.88164E-060.60 7.73348E-066.14421E-06-1.58927E-061.36.27225E-075.10909E-12-6.27220E-070.65 2.90006E-062.26033E-06-6.39727E-071.
14、42.09075E-076.91440E-13-2.09074E-070.70 1.08752E-068.31529E-07-2.55993E-071.56.96917E-089.35762E-14-6.96916E-080.75 4.07820E-073.05902E-07-1.01918E-071.62.32306E-081.26642E-14-2.32306E-080.80 1.52933E-071.12535E-07-4.03975E-081.77.74352E-091.71391E-15-7.74352E-090.85 5.73498E-084.13994E-08-1.59504E-
15、081.82.58117E-092.31952E-16-2.58117E-090.90 2.15062E-081.52300E-08-6.27618E-091.98.60390E-103.13913E-17-8.60390E-100.95 8.06481E-095.60280E-09-2.46201E-092.0 2.86800E-104.24835E-18-2.86800E-101.00 3.02430E-092.06115E-09-9.63150E-101.05 1.13411E-097.58260E-10-3.75854E-101.10 4.25293E-102.78947E-10-1.
16、46346E-101.15 1.59485E-101.02619E-10-5.68660E-111.20 5.98068E-113.77513E-11-2.20554E-111.25 2.24275E-111.38879E-11-8.53960E-121.30 8.41033E-125.10909E-12-3.30124E-121.35 3.15387E-121.87953E-12-1.27434E-121.40 1.18270E-126.91440E-13-4.91262E-131.45 4.43513E-132.54367E-13-1.89147E-131.50 1.66318E-139.
17、35762E-14-7.27413E-141.55 6.23691E-143.44248E-14-2.79443E-141.60 2.33884E-141.26642E-14-1.07242E-141.65 8.77065E-154.65889E-15-4.11176E-151.70 3.28899E-151.71391E-15-1.57509E-151.75 1.23337E-156.30512E-16-6.02861E-161.80 4.62515E-162.31952E-16-2.30563E-161.85 1.73443E-168.53305E-17-8.81126E-171.90 6
18、.50411E-173.13913E-17-3.36498E-171.95 2.43904E-171.15482E-17-1.28422E-172.00 9.14640E-184.24835E-18-4.89805E-18从表中可以看出,当步长取0.1时,Runge-Kutta法是收敛的,但其最终全局误差尽管绝对在很小,相对误差却很大。比真实值要大数个量级。从微分方程的显式解()可以看出,该函数是急遽下降的,其导数,即斜率则是急遽增大的。Runge-Kutta法可由泰勒方法推导而来,其思想即是用xk点处的斜率,xk+1点处的斜率及中点斜率的两个估计值的线性组合来代替,利用公式来预测的值。如果变
19、化很剧烈,而步长又取得较大,则可能在中间过程中会出现“大数吃小数”的现象,使结果精度降低。当步长减半,取0.05时,从表中可以看出其计算精度要好很多,但相对误差仍较大。对于浮点数而言,越靠近0,精度就越高;越远离0,精度就越低。MATLAB中的浮点误差值eps为1到下一个能表示的比1大的浮点数之间的距离。换言之,MATLAB是无法区分1与1+eps之间的数的,它们将被舍入到1或1+eps。在MATLAB命令窗口键入eps可得到其值为2.220446049250313e-16。由显示函数计算的结果从1降到1e-18量级,小于eps,中间过程必然有舍入误差。三、 非线性方程根的求解利用牛顿法求方程
20、于区间的根,考虑不同初值下牛顿法的收敛情况。1. 分析与计算令,转化为求在区间上的零点。通过观察方程易知零点在3附近。取初值x0=3进行牛顿法迭代。编写MATLAB程序,计算结果如表3.1所示,可见迭代收敛很快。如果从区间端点开始迭代,初值分别取2和4,计算结果如表3.2所示,同样,迭代收敛很快。表3.1 Newton迭代法计算结果k03.00000000 -0.09861228913.1479184330.1479184330.00117701423.1461934410.0017249921.50195E-0733.1461932212.20177E-072.66454E-15表3.2 N
21、ewton迭代法计算结果k02.00000000 -0.693147184.00000000 0.61370563913.3862943611.3862943610.1665581463.1817258150.8182741850.02430205623.1499383940.2363559670.0025554993.1462848450.035440976.25023E-0533.1461942570.0037441377.06991E-073.1461932219.16234E-054.24029E-1043.1461932211.03641E-065.41789E-14对于实验方程可用
22、MATLAB画出其在区间2,4上的图形,如图3.1所示。可见其在区间内只有单根,对于单根,Newton法具有2阶收敛速度。图3.1 f(x)在2,4上图形四、 数值积分已知,利用复化梯形公式、复化Simpson公式和Romberg算法求的近似值;并观察实际计算结果,比较它们的收敛速度。1. 计算与分析(1)复合梯形公式与复合Simpson公式分别编写复合梯形公式与复合Simpson公式的MATLAB程序,计算结果如表4.1所示。表4.1 复合梯形公式与复合Simpson公式计算结果子区间数n误差误差13.000000000000000.141592653589793.133*3330.0082
23、593202564623.10000000000000.0415*3.141568627450980.0000240261388143.131*824.010*3.141592502458710.0000001511310983.138*1090.002604159098703.141592651224820.00000000236497163.140941612041390.000651041548403.141592653552840.00000000003696可以看到,相对于复合梯形公式,复合Simpson公式的误差更快地收敛到0。由于梯形公式的误差项为阶,Simpson公式的误差项为
24、,因此Simpson公式比梯形公式具有更高的精度,所需计算的次也更少。(2) Romberg公式编写Romberg算法的MATLAB程序。该方法需要逐次计算梯形公式值,Simpson公式和Cotes公式值,通过其线性组合来提高精度。中间计算值存在矩阵中,如表4.2所示。表4.2 Romberg积分表k03.00000000000000 13.10000000000000 3.13333333333333 23.13117647058824 3.14156862745098 3.14211764705882 33.13898849449109 3.14159250245871 3.1415940
25、9412589 3.14158578376187 43.14094161204139 3.14159265122482 3.14159266114256 3.14159263839680 53.14142989317497 3.14159265355284 3.14159265370804 3.14159265359003 表4.3为Romberg积分的误差表。表4.3 Romberg积分误差表k00.14159265358979 10.04159265358979 0.00825932025646 20.01041618300156 0.00002402613881 0.0005249934
26、690330.00260415909870 0.00000015113109 0.000001440536100.00000686982792 40.00065104154840 0.00000000236497 0.000000007552770.00000001519300 50.00016276041482 0.00000000003696 0.000000000118240.00000000000024可以看到Romberg积分表中第二序列正是复合Simpson公式计算结果。同样容易验证,第三序列是复合Cotes公式计算结果。第四序列(Romberg序列)则是n=8时的Newton-C
27、otes公式。一般而言,右边的列收敛速度要快于左边的列,但并不是说所用公式的阶数越高,结果越精确。因为高阶Newton-Cotes公式的很大,导致数值稳定性较差。从误差表中可以看出,同一行的Cotes序列并不比Simpson序列精度更高。但对于相同子区间数,右边列比左边列精度要高。Romberg序列收敛最快,同时精度也最好。总结本次数值分析实验采用了MATLAB程序设计语言,分别编写了Gauss消元法、Gauss-Seidel迭代法程序求解线性方程组;编写Newton迭代法程序求解非线性方程;编写Runge-Kutta法程序求解线性微分方程;编写复合梯形公式、复合Simpson法程序和Romberg公式程序求解数值积分问题。对Gauss消元法采用了选列主元的方法,增加了程序的通用性;对Gauss-Seidel迭代法讨论了应用条件,要求其迭代矩阵满足一定的条件;对Newton迭代比较了不同初值的收敛速度;对Runge-Kutta法讨论了收敛的稳定区间,要求步长满足一定要求。比较了复合梯形公式、复合Simpson公式和Romberg公式的收敛效率
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1