数值计算方法I试验手册

《数值计算方法I试验手册》由会员分享，可在线阅读，更多相关《数值计算方法I试验手册（41页珍藏版）》请在装配图网上搜索。

1、计算方法实验手册江西财经大学信息管理学院数学与决策科学系2003.1- 1 -目录MATLAB268Gauss10Runge1315171820Monte-Carlo23Euler2628293133343536- 2 -预备实验MATLAB使用练习一、实验目的MATLAB是本实验课的主要软件，尽快熟悉其基本命令和语法，以适应后续实验对MATLAB的使用要求掌握 MATLAB中常用的基本命令和函数；掌握 MATLAB的简单作图方法；掌握 MATLAB的简单编程技巧二、实验内容 . 矩阵和向量的输入和运算在 MATLAB的 Commond窗口中输入以下命令，并仔细观察执行后的结果（ 1）矩阵的输

2、入：A=1 2 3;4 5 6（分号表示矩阵换行，表示回车）A=1,2,3;4,5,6（空格或逗号表示同一行元素的分隔）A=1 2;3 4;（在一行的末尾加分号，运行结果将不显示）Z=zeros(2,3)（ 23 阶零矩阵）E=ones(2)（ 22阶全 1 矩阵）I=eye(2,4)（ 24阶对角线元素为1 的矩阵）R=rand(2,3)（ 23 阶 (0,1) 均匀分布随机矩阵）H=hilb(4)（4 阶 hilbert矩阵）另外还有标准正态分布矩阵randn(m,n),n阶幻方矩阵 magic(n)等 .（ 2）矩阵元素的调用和赋值：a=A(1,3)（矩阵 A 的第行第列元素赋值给a）A

3、(1,3)（查看矩阵 A 的第行第列元素）A(1,3)=0（将 A 的第行第列元素赋值为）A(3,4)=1（注意矩阵 A 原本是 23 阶矩阵）B(2,3)=2（注意前面没有给出矩阵B）（ 3）矩阵的裁剪和拼接：A（显示矩阵 A）A(2,:)（ A 的第行）A(:,3)（A的第 3列）C1=A(1:2,:)（将 A 的第 , 行赋值给矩阵 C1）- 3 -C =A(:,: )（将 A 的第列赋值给矩阵C2）C=A(: , : )（将 A 的第 , 行 , 第列赋值给 C3）C4= A(1: ,4:-1:2)（将 A 的第 1、行 , 第、列赋值给C3）A（显示矩阵 A 为 34 ）B1=A,C

4、2（左右拼接矩阵A和 C2为 B1）B2=A； C1（上下拼接矩阵A和 C1为 B2）（）矩阵的运算MATLAB为矩阵的运算提供了下列矩阵运算符：（加法）；（减法）； （转置）；（乘法） ；（乘幂）；（左除）； / （右除）. （点乘法） ；. （点乘幂） ； . （点左除） ； ./ （点右除）矩阵的“” 、“”、“”、“”、“”必须遵循矩阵的运算规则，点运算“. ”、“.”、“ . ”、“ ./ ”是指矩阵对应元素之间进行相应的运算在 MATLAB的 Commond窗口中继续输入以下命令，并仔细观察结果：C=zeros(3)D=ones(3)比较 C+1和 C+DA=1 2;4 5 6;7

5、 8 9B=1 2 3（ B 为一列矩阵即列向量）X=AB（左除表示左边乘以A 的逆矩阵 ,X 为方程组 AX=B的解）X=B/A（右除表示左边乘以A 的逆矩阵 ,X 为矩阵方程 XA=B的解）C=AA*CA .*C（ A,C 的对应分量相乘）A2A .2（ A 的每个分量平方）A/A（相当与 AA 1 ）A./A（结果为全矩阵）在 MATLAB中，数组、向量与矩阵的形式完全一致，并且运算与矩阵也基本类似，但要注意行列向量之分除了与矩阵一样的输入方式外，数组（向量）还有以下的输入方式：a=1:5b=1:2:8- 4 -c=10:-2:0b=0:2:8,zeros(1,3)(相当与矩阵的拼接)c

6、=linspace(1,10,10)(从到 10 的共 10 个数值的等差数组)x=0:pi/4:pi(pi是圆周率 ) . 变量和函数MATLAB中的变量区分大小写字母，第一个必须是字母，不需要说明变量类型或维数，系统自动分配内存空间 , 字符串是用单引号括起来的字符集合，可以像向量或矩阵一样进行拼接与剪裁MATLAB有几个特殊的变量：pi( 圆周率 ) ；eps（最小浮点数） ；inf(正无穷大， 特指 1 ) ；0NaN(不定值，特指0 ) ； i,j( 虚数单位 ) 观察下列结果：0a=0 1 0,b=1 0 0,c=a./bMATLAB提供了大量的函数，按照分类分为标量函数、向量函数

7、和矩阵函数三种类型请用 help 命令查看下列各种函数的调用格式：标量函数： sin cos tan cot sec csc asin acos atan acot asec acsc sqrt exp log log10 abs round floor ceil fix sign real imag向量函数： max min sum length mean median prod sort矩阵函数： zeros ones eye rand randn diag det rank inv eig trace expm poly norm cond lu qr svd . 基本图象功能二维图象最

8、常用的命令有：基本作图命令： plot fplot polar图象注释命令： title xlabel ylabel textgtext legend坐标管理命令： axis窗口管理命令： hold onhold offgrid subplot figure二维图象最常用的线型和颜色有：线型（线方式） ：实线：点线虚点线波折线线型（点方式） ：圆点加号星号xX 形o 小圆颜色 ：y 黄r 红g 绿b 蓝w 白k 黑m 紫c 青 . 程序设计关系运算符：大于 =大于或等于=等于=不等于- 5 -关系运算比较两个数值， 当指出的关系成立时结果为， 否则为， 关系运算可以作用于两个同型的矩阵或向量，

9、 比较时逐个元素或分量比较， 结果是一由和组成的矩阵或向量逻辑运算符：& 与运算（表示且）|或运算（表示或）非运算（表示不）流控制语句： if （条件语句）；while （循环语句）；for（循环语句）；switch（条件语句），它们都用end 结束，具体形式用help 命令查看文本 M 文件：由语句组成的MATLAB命令集，运行时不必输入参数函数 M 文件：可以方便自己调用的函数文件，第一行有特殊要求，其形式必须为function函数名（自变量）该文件的存储名必须与函数名相同三、实验要求本实验主要是为了熟悉MA TLAB 而设计的 MATLAB是本实验课的主要软件，它将计算、可视化和编程功能

10、集成在非常便于使用的环境中，是一个交互式的、以矩阵计算为基础的科学和工程软件，具有编程效率高、计算功能强、使用简便和易于扩充等特点要求通过本实验，尽快熟悉其基本命令和语法，以适应后续实验对MA TLAB 的使用要求- 6 -实验一数值计算中误差的传播规律一、相关原理1设 x * 是 x 的一个近似值，则有如下概念绝对误差E( x*)x * x ，简称误差，简记为E误 差 限( x*)E( x*) ，简记为E r ( x*)E( x*)x *x相对误差xx，简记为 ErE r* ( x*)E(x*)x*x或x*x *，简记为 E r *相对误差限r(x*)Er (x*) ，简记为r或r * (

11、x*)E r * ( x*) ，简记为r *2设 yf ( x1, x2 ) ， x1* 是 x1的一个近似值，x2* 是 x2 的一个近似值， y*f ( x1* , x2* ) ，则 y * 为相应的 y f ( x1 , x2 ) 的近似值设 E1 , E2 分别为 x1* , x*2 的绝对误差；1 , 2 分别为 x1* , x2* 的绝对误差限；Er*1 , Er*2 分别为 x1*, x2* 的相对误差；r*1, r*2 分别为 x1* , x2* 的相对误差限；f*f*则y * 的绝对误差为E( y*)E1E2x1x2f*f*y * 的绝对误差限为12x1x2x1*x2*y *

12、 的相对误差为Er ( y*)fEr*1fEr*2y *x1y *x2x1*x2*y * 的相对误差限为*f*f*ry *x1r 1y *x2r 23.设 x*0.a1a2am10k （ a10 ）是 x 的一个近似值，如果x * 有 n （ nm ）位有效数字，则Er ( x*)1101n2a1- 7 -如果 x * 的相对误差满足Er ( x*)1101 n2(a11)则有 x * 至少有 n （ nm ）位有效数字二、实验目的观察并初步分析数值计算中误差的传播；观察有效数字与误差传播的关系三、实验内容使用 MATLAB的 help命令学习 MATLAB命令 digtis和 vpa 的用途

13、和使用格式；在 4 位浮点数下解二次方程x 262 x 1 0 ；计算下列个函数在点x2处的近似值（） y0(x1)6 ，（）y116 ，( x1)（） y2(32x)3 ，（）y31，(32 x)3（）y49970x 四、实验要求本次实验包含三个相对独立的内容在内容中，请解释两个命令的格式和作用；求解方程 x262 x 10 时，分别使用求根公式和韦达定理两种方法，并比较其有效数字和相对误差；实验内容中的个函数在x2 处的精确值都是相等的，若取2 1.4 进行计算，计算各函数的结果，作图观察并比较它们的绝对误差（作图区间可取1.4,1.42 甚至更小） , 并从算法设计原则上说明原因- 8

14、-实验二数值计算中的算法稳定性一、实验目的误差扩张的算法是不稳定的，是我们所不期望的；误差衰竭的算法是稳定的，是我们努力寻求的， 这也是贯穿本课程的目标本实验的目的是通过实验，体会稳定性在选择算法中的地位二、实验内容考虑一个简单的由积分定义的序列I n1n 0,1,2,xnex 1 dx ,（ E）0显然 I n 0, n 当 n0 时，显然有1e 1I 0ex 1dx 10而对 n1，由分部积分法不难得到I n1nex 1dxxn x 1 11n 1 x 10xenxedx00即I n 1nI n 1,n 1,2,3,（E）同样由（ E）式，我们还可以得到I n1xn ex 1dx1100x

15、 n dx（E）n 1由递推关系（ E）式，可以得到计算积分序列（E）的两种算法：第一种是（E）的直接应用，即E01e 1（E）En1nEn 1 ,n1,2,3,另一种是利用（E）和（ E）的变形得到I N0I n1 I n（ E）1n, n N 1, N 2, ,3,2,1三、实验要求分别利用算法（E）和（ E）式计算，并在计算中分别采用位、位和位有效数字，请判断哪种算法能给出更精确的结果；- 9 -两种算法的优劣，与你的第一感觉是否吻合？请从理论上证明你实验得出的结果，解释实验的结果设算法（E）中 I 0 的计算误差为 e0 ，由 I 0 递推计算 I n 的误差为en ；设算法（ E）中

16、 I N 的计算误差为N ，由 I N 向前递推计算 I n （ nN ）的误差为n ；如果在上面两种算法中都假定后面的计算不会再引入其他误差，试给出en 与 e0 、 n 与 N的关系；算法（ E）中 e0通常会很小，当 n 增大时，误差 en 的变化趋势如何？算法（E）中 N 通常相对比较大，当 n 减小时，误差 n 又是如何传播的？即比较上述两个算法，当某一步产生误差后，该误差对后面的影响是衰减还是扩张？通过理论分析和数值实验，针对算法（E）和（ E）的稳定性，给出你的结论- 10-实验三Gauss消去法主元素的选取与算法的稳定性一、相关原理解线性方程组Axb , ARn n ,bRn其

17、中a11a12a1nb1Aa21a22a2n， bb2，an1an 2annbn记 A(1)A ， b(1)b Gauss 列主元消去法通常包含下面四个步骤：选取主元素： 在矩阵 A(k ) 的第 k 列元素 akk(k ) 及其以下的元素中按某种规则（如按模最大或最小） ，寻找到列主元，记录其行号，记为i k ，并且当 ai(kk, k)0 时进行换行，否则退出消去法换行：交换矩阵( A( k) ,b(k ) ) 的第 k 行与 i k 行，即取1101 k行1I ik ,k(ik k)110 ik 行11作乘积 I i k ,k ( A( k ) , b( k) ) ，即为交换矩阵( A(

18、 k ) , b( k ) ) 的第 k 行与 ik 行，当 ikk 时，表示不需换行换行后，仍将 I i k ,k ( A( k ) , b( k) ) 的元素记为 aij( k ) 消元过程：定义行乘数 mikaik(k ),n ），令( k ) （ i k 1,akk- 11-11Lkmk 1,k1mk 2,k1mn ,k1作乘积 Lk I i ,k ( A( k) ,b(k) )( A(k 1) , b( k1) k重复以上过程， k1,2, , n1，当 k n 1时， Gauss 选主元的消元过程结束回代：通过上述过程，增广矩阵( A(1) ,b(1) ( A(n ) ,b(n )

19、 ) ，并且 A( n) 是上三角矩阵a11(1)a12(1)a1(1n)b1( 1)A( n )a22(2)a2( n2)， b(n )b2( 2)ann(n )bn( n)当 ann(n )0 时，利用下面的回代公式即可得到方程组的解：xnbn(n )ann(n)nbi(i )aij(i )xij i1,(i n 1, n 2, ,2,1)aii( i )Gauss 列主元消去法结束二、实验目的简单了解条件数的概念（条件数的理论将在教材的第六章介绍）；通过实验了解算法的稳定性与主元素及条件数的关系三、实验内容Gauss 消去法是我们在线性代数中已经熟悉的但由于计算机的数值运算是在一个有限的

20、浮点数集合上进行的，如何才能确保Gauss 消去法作为数值算法的稳定性呢？Gauss 消去法从理论算法到数值算法，其关键是主元素的选择主元素的选择从数学理论上看起来平凡，它却是数值分析中十分典型的问题考虑线性方程组Axb , ARn n ,bRn编制一个能自动选取主元，又能手动选取主元的求解线性代数方程组的Gauss 消去过程 三、实验要求- 12-6171861151取矩阵 A， b，则方程组有解精确解 x，取86115186141n 10 ，用 MATLAB的 cond 命令计算矩阵 A 的条件数让程序自动选取主元，结果如何？现选择程序中手动选取主元素的功能 每步消去过程总选取按模最小或按

21、模尽可能小的元素作为主元素，观察并记录结果 若每步消去过程总选取按模最大的元素作为主元素，结果又如何？分析实验的结果取矩阵阶数n20 或者更大，重复上述实验过程，观察并记录分析不同的问题及消去过程中选择不同的主元素时计算结果的差异，说明主元素的选取在消去过程中的作用选取其他你感兴趣的问题或者随机生成矩阵，计算其条件数重复上述实验，观察记录并分析实验的结果选择著名的Hilbert矩阵 H 作方程组 Hx b 的系数矩阵，其中H(hij ) nn ， hij1， i , j1,2,ni j1,bn ) Tnb( b1 , b2 ,， bihij ， i1,2, nj 1该方程组的精确解为x(1,1

22、,1) T 分别选择方程的维数为 3，6， 10，计算矩阵 H 的条件数，并用Gauss 列主元消去法求解，记录结果，并给出你的结论- 13-实验四观察 Runge 现象和对非光滑函数进行插值的可能性一、相关原理f ( x) 在节点 ax0x1xn b 处的函数值为 y0 , y1 , yn ，构造其 Lagrange插值多项式 Ln ( x) 的插值基函数为l j ( x)nxxix j， j 0,1,2, , ni0xiijLagrange 插值多项式为nLn ( x)y j l j ( x)j 0其截断误差为f ( n 1) ()1 ( x)Rn (x)n(n1)!n其中(a,b) ，

23、n 1 ( x)( xxi )i 0二、实验目的观察高次 Lagrange插值多项式 Ln (x) 的 Runge 现象；观察非光滑函数进行多项式插值的可能性三、实验内容考虑在一个固定的区间上用Lagrange 插值逼近一个函数显然 Lagrange 插值中使用的节点越多， 插值多项式似的次数就越高我们自然关心插值多项式的次数增加时，Ln (x)是否也更加靠近被逼近的函数设区间 1,1 上函数1f ( x)125x2考虑区间 1,1 上的一个等距分割，节点为xi12i , i 0,1,2, , nn作 f ( x) 在 1,1 上的 Lagrange插值多项式n1Ln ( x)l i ( x)

24、i 0 1 25xi2- 14-其中 li (x) 为 Lagrange 插值基函数连续非光滑函数的几何特性非常差，在几何图象上一般会出现大量的尖点在构造非光滑函数的多项式插值时，由于多项式具有高阶光滑度，两者之间会产生怎样的现象？选择区间 0,1 上的连续非光滑函数g (x)sin k x作 g( x) 区间 0,1 上的 Lagrange 插值多项式四、实验要求选择不断增大的节点数目 n 2,3, ，画出原函数 f ( x) 及插值多项式 Ln ( x) 在区间的 1,1 上的图象，比较并分析实验结果选择其他的函数，例如定义在区间 5,5 上的函数u( x)x,v( x)arctan(x)

25、x41重复上述的实验过程，观察其结果又将如何如果不取等距节点，而改为取如下节点xib ab a cos2i1， i 0,1,2, n222(n1)以 x0 , x1, , xn 为插值节点构造上述函数f (x), u( x), v( x) 的 Lagrange 插值多项式，比较其结果选择不同的 k 和 n ，用等距节点作g( x) 的 n 次 Lagrange插值多项式，观察其误差大小及收敛情况- 15-实验五观察及改善最小二乘拟合的数值不稳定现象一、相关原理最小二乘拟合：设f ( x) 是定义在区间a, b 上的函数， x0 , x1 , xm 是区间 a,b上 的 一 组 节 点 ， yi

26、f ( xi ) （ i0,1,2,m ） 是 函 数 f ( x) 对 应 于 节 点 的 函 数值0 ( x),1 ( x),n ( x) 是定义在区间 a, b 上的线性无关的连续函数，是由拟合基函数0 ( x),1( x),n (x) 生成的函数类，即span0 ( x),1 ( x),n ( x)n在中必存在一个函数S* ( x)a*jj( x) 是函数 f ( x) 在 a,b 上的最小二乘拟合j 0函数，其中拟合系数a*j满足法方程组(0,0)(0,1)( 0 , n )a0( 0 , f )(1,0)(1,1)(1,2)a1( 1 , f )，( n , 0 ) ( n , 1

27、 )( n , n )an( n , f )记为 Aab ，其中(0,0)(0,1)( 0, n )a0( 0, f )A(1,0)(1,1)(1,2)， aa1， b( 1 , f )( n , 0 )( n , 1 )( n , n )an( n , f )mm( i ,j )i (xk ) j ( xk ) ， ( i , f )i ( xk ) yik 0k0最小二乘解的平方误差为2m(S * ( xi ) yi ) 2( f , f )( S*, f )2i 0( f ,f )( a*, b)mnyi2a*j ( j , f )i 0j 0 Legendre （勒让德）正交多项式在区

28、间 1,1 上正交多项式- 16-P0 ( x)1Pn ( x)1d n2n , n 1,2,3,4,nn! dxn ( x1)2称为 Legendre 正交多项式 Legendre 正交多项式具有递推公式：P0 (x)1 , P1 ( x)x(n1)Pn 1 (x)(2n1) xPn (x) nPn 1 ( x) ， n 1,2,二、实验目的1观察最小二乘拟合多项式的数值不稳定现象；2探索改善最小二乘拟合多项式的不稳定现象的可能性三、实验内容在区间 1,1 上取 m20 个等距节点， 计算出相应节点上函数ex的值做为样本数据，以 1, x, x2 , xn 为基函数作出 n3,5,7,9 次

29、的最小二乘多项式；在区间 1,1 上取 m20 个等距节点， 计算出相应节点上函数ex 的值做为样本数据，以 Legendre多项式P0(),P1( ),( ) 为基函数作出 n 3,5,7,9 次的最小二乘多项xxPn x式四、实验要求在区间 1,1 上取 m20 个等距节点，计算出相应节点上函数ex 的值做为样本数据，以2n1, x, x , x为基函数作出n3,5,7,9次的最小二乘多项式用的condMATLAB命令，求出相应的法方程组的系数矩阵A 的条件数，画出 ln( cond ( A) n 之间的曲线计2算出不同次数的最小二乘拟合多项式的平方误差(n) 2 仍在区间 1,1 上取

30、m20个等距节点，计算出相应节点上函数ex 的值做为样本数据，以 Legendre多项式P0(),P1(x), ,( ) 为基函数作出 n3,5,7,9 次的最小二乘多xPn x项式用MATLAB的 cond 命令求出相应的法方程组的系数矩阵A 的条件数，画出2ln( cond ( A) n 之间的曲线 计算出不同次数的最小二乘拟合多项式的平方误差( n) 2 比较使用两类不同基函数得到的拟合结果，给出你的结论- 17-实验六曲线逼近方法的比较一、实验目的掌握各种函数逼近方法的适用范围及实际应用中方法选择时应注意的问题二、实验内容曲线的拟合和插值，是逼近函数的基本方法，每种方法具有各自的特点和

31、特定的适用范围，实际工作中合理选择方法是重要的现在仍然考虑实验四中的著名问题，即用不同的插值和拟合方法作区间 1,1 上函数f ( x)125x 21的逼近函数，并加以比较下面的MATLAB程序给出了该函数的二次和三次拟合多项式x=-1:0.2:1;y=1./(1+25*x.*x);xx=-1:0.02:1;p2=polyfit(x,y,2);yy2=polyval(p2,xx);plot(x,y,0,xx,yy2)xlabel(x);ylabel(y);hold onp3=polyfit(x,y,3);yy3=polyval(p3,xx);plot(xx,yy3)hold off适当修改上述

32、MATLAB程序，也可以拟合其他你感兴趣的函数三、实验要求将拟合的结果与Lagrange 插值及样条插值的结果比较；归纳总结数值实验的结果，试定性地说明函数逼近各种方法的适用范围，及实际应用中选择方法应注意的问题- 18-实验七利用数值微分方法求解偏微分方程一、相关原理常见的数值微分公式及其误差阶在数轴上给定一组点x0 , x1, xn 及在相应点的函数值f ( x0 ), f ( x1), f ( xn ) 假定xk 1xkh ，则f ( xk ) f ( xk ) f ( xk ) f ( xk )f (xk 1)f ( xk )O (h) （两点公式）hf (xk 2 )4 f ( xk

33、 1 ) 3 f ( xk )O(h2 ) （三点公式）2hf ( xk 1 )f (xk 1 )O(h2 )（三点中点公式）2hf ( xk 1 )2 f ( xk )f (xk 1 )O (h2 ) （二阶导数的三点公式）h2常见的数值偏导数公式及其误差阶（补充内容）设 u(x, y) 是定义在矩形区域axb, cy d 上的二元函数，将该区域等距划分为M N 个小格子， h和分别为 x 轴方向和 y 方向的步长， 并把格子线的交点记为 ( xi , y j ) ，记 u( xi , y j ) uij ， i0,1,2, M ;j0,1,2, N 函数 u( x, y) 在点 ( xi

34、, y j ) 处的偏导数和二阶偏导数可表示为差商形式u( xi , y j )xu( xi 1, y j ) u( xi , y j )O (h)hu( xi , y j )y2u(xi , y j )2u(xi , y j )u( xi , y j 1 )u( xi , y j )O( )u(xi 1, yj )2u(xi , y j )u( xi 1 , y j )h2O(h2 )u(xi , y j 1 )2u(xi , y j )u( xi , y j 1 )2 )2O(简单偏微分方程数值解的原理设 偏 微分 方 程 的 未 知函 数 是 定 义 在矩 形 区 域 axb, cyd

35、上 的 二 元 函 数u( x, y) ，将该区域等距划分为MN 个小格子，h和分别为 x 轴方向和y 方向的步长，- 19-并把格子线的交点记为(xi , y j ) ，记 u(xi , y j )uij ，i0,1,2, M ;j0,1,2, N 在每个格点上用上述数值导数表示相应的偏导数，从而可以建立关于（未知函数u( x, y) 在格点 (xi , yj ) 上函数值） uij 的一个方程组，解出这个方程组即可得到偏微分方程的数值近似解多数情形下方程组为稀疏的线性方程组，求解时可用相关的稀疏方程组的知识求解二、实验目的通过实验， 掌握使用数值微分的方法，将简单的偏微分方程化为方程组，进

36、而求出偏微分方程在格点处近似解的方法三、实验内容在单位正方形上有偏微分方程2 u2uxy( x22)x 2y2y其边界条件为u(0, y)u( x,0) 0 ；u(1, y)y3x36;u( x,1).6使用数值微分的方法求解这个偏微分方程四、实验要求分别取 N 5,10,15 把单位正方形划分成NN 个小格子，并把格子线的交点记为( xi , yj ) ，显然 xiih ， y jjh ， h1， i , j0,1,2, , N ，记 u(xi , y j ) uij 写出N用 数 值 微 分 公 式 建 立 的 关 于 uij 的 线 性 代 数 方 程 组 ， 并 求 解 该 方 程 的

37、 精 确 解 是( xy)3u( x, y)将精确解和近似解画在同一张三维图上6画三维图请参考MATLAB的 plot3 ，mesh， surf 等命令- 20-实验八比较不同数值积分公式对积分方程解的影响一、相关原理复合梯形公式原始积分I ( f )bf ( x) dxahn积分公式Tn ( f ) f (a)f (b)2f (xk )2k 1余项E(Tn )b a h2 f ( )12复合 Simpson 公式原始积分I ( f )bf (x)dxaba nn1积分公式Sn ( f )f(a)f()2f(xk) 4f(x1) 6nbkk01k2ban1n 1或 ( )( )2() 4()Sn ( f )6nf af bf x2 kf x2k 1k1k 0b ah4余项E(Sn )f (4 ) ()1802 Gauss Legendre 积分公式（补充）原始积分I ( f )1f ( x