定积分的近似计算

《定积分的近似计算》由会员分享，可在线阅读，更多相关《定积分的近似计算（9页珍藏版）》请在装配图网上搜索。

1、实验二 定积分的近似计算一、问题背景与实验目的利用牛顿莱布尼兹公式虽然可以精确地计算定积分的值，但它仅适用于被 积函数的原函数能用初等函数表达出来的情形如果这点办不到或者不容易办到 这就有必要考虑近似计算的方法在定积分的很多应用问题中，被积函数甚至没 有解析表达式，可能只是一条实验记录曲线，或者是一组离散的采样值，这时只 能应用近似方法去计算相应的定积分本实验将主要研究定积分的三种近似计算算法：矩形法、梯形法、抛物线法对于定积分的近似数值计算，Matlab有专门函数可用.厢巴二、相关函数(命令)及简介1. sum(a)：求数组a的和.2. format long:长格式，即屏幕显示15位有效数

2、字.(注：由于本实验要比较近似解法和精确求解间的误差，需要更高的精度).3. double。：若输入的是字符则转化为相应的ASCII码；若输入的是整型数值 则转化为相应的实型数值.4. quad() ：抛物线法求数值积分.格式： quad(fun,a,b) ，注意此处的 fun 是函数，并且为数值形式的，所以 使用*、/、八等运算时要在其前加上小数点，即*、/、等.例：Q = quad(1./(x3-2*x-5),0,2);5. trapz():梯形法求数值积分.格式： trapz(x,y)其中 x 为带有步长的积分区间； y 为数值形式的运算(相当于上面介绍的函数 fun)例：计算曲x=0:

3、pi/100:pi;y=sin(x);trapz(x,y)6. dblquad()：抛物线法求二重数值积分.格式：dblquad(fun,xmin,xmax,ymin,ymax), fun 可以用 inline 定义，也可 以通过某个函数文件的句柄传递例1： Q1 = dblquad(inline(y*sin(x), pi, 2*pi, 0, pi)顺便计算下面的Q2,通过计算，比较Q1与Q2结果(或加上手工验算)， 找出积分变量x、y的上下限的函数代入方法.Q2 = dblquad(inline(y*sin(x), 0, pi, pi, 2*pi)例 2： Q3 = dblquad(inte

4、grnd, pi, 2*pi, 0, pi) 这时必须存在一个函数文件 integrnd.m： function z = integrnd(x, y)z = y*sin(x);7. fprintf (文件地址，格式，写入的变量)：把数据写入指定文件.例： x = 0:.1:1;y = x; exp(x);fid = fopen(exp.txt,w);%打开文件fprintf(fid,%6.2f %12.8fn,y);%写入fclose(fid)%关闭文件8. syms变量1变量2：定义变量为符号.9. sym(表达式)：将表达式定义为符号.解释： Matlab 中的符号运算事实上是借用了 Ma

5、ple 的软件包，所以当在 Matlab中要对符号进行运算时，必须先把要用到的变量定义为符号.10. int(f,v,a,b):求f关于v积分，积分区间由a到b.11. subs(f，x，a):将a的值赋给符号表达式f中的x，并计算出值若 简单地使用subs(f)，则将f的所有符号变量用可能的数值代入，并计算出 值。三、实验内容一、问题的提出计算定积分的方法：(1) 求原函数；(2) 利用牛顿莱布尼茨公式计算结果。问题：(1) 被积函数的原函数不能用初等函数表示；(2) 被积函数难于用公式表示，而是用图形或表格给出的(3) 被积函数虽然能用公式表示，但计算其原函数很困难 解决办法：建立定积分的

6、近似计算方法。思路：在数值上表示曲边梯形的面积，只要近似地算出相应的曲边梯形的面积，就可得 到所给定积分的近似值。常用方法：矩形法、梯形法、抛物线法二、矩形法用分点a=x ,x，x =b将区间a,bn等分，取小区间左端点的函数值01 ny (i=0,l,2,n-1)作为窄矩形的高，如图：i则有：取小区间右端点的函数值yj(i=l,2,n)作为窄矩形的高，如图:则有：以上两公式称为矩形法公式。例：用矩形法求，并与用牛顿莱布尼茨公式计算的结果进行比较程序如下： #include #include void main()double result,a=0,b=1,i,n=1000000,h;prin

7、tf( 按牛顿公式计算得到的结果： %fn,sin(b)-sin(a);result=0;h=(b-a)/n;/计算区间高度for(i=l;i二n;i+)/求和result=result+cos(a+i*h);resul t二h *resul t; /乘以区间高度printf（“用近似公式计算得到的结果：fn,result）;三、梯形法y珀梯形法就是在每个小区间上，以窄梯形的面积近似代替窄曲边梯形的面积，如图：则有：皿宅+仇+乃）+旳）心+几-i + xJAx=儿）十1十必十十儿例：用梯形法求，并与用牛顿莱布尼茨公式计算的结果进行比较。#include #include void main()

8、double result,a=0,b=1,i,n=1000000,h;printf( 按牛顿公式计算得到的结果： %fn,sin(b)-sin(a);result=0;h=(b-a)/n;/ 计算区间高度for(i=1;i=n-1;i+)/ 求和result=result+cos(a+i*h);result+=(cos(a)+cos(b)/2;result=h*result;/ 乘以区间高度printf( 用近似公式计算得到的结果： %fn,result);四、抛物线法此法就是将曲线分成许多小段，用对称轴平行于 y 轴的二次抛物级上的一段弧来 近似替代原来的曲线弧，从而得到定积分的近似值。用

9、分点a二x,x，x =b将区间a,bn等分(偶数)，这些分点对应曲线上的点为 01 nM(x,y)(其中 y=f(x),i=0,l,2,，n)，如图：i i i i i因为经过三个不同的点可以唯一确定一条抛物线，可将这些曲线上的点M互相i衔接地分成n/2组M,M,M,M,M,M,，M ,M ,M，即每相邻两个区间为012234n-2n-1 n一组。在每组 M , M , M (k=l,2,n/2)所对应的子曲间x ,x 上，用2k-2 2k-1 2k 2k-2 2k经过点M , M , M的二次抛物线近似代替曲线弧。2k-22k-1 2k下面讨论如何计算积分二。设h为区间高度，即h=x -x

10、=x -x。根据积分性质(积分在数值上表示曲2k 2k-12k-12k-2边梯形的面积)有如下等式成立：j/ _ (/T-= j ： (/r- +毀 + f)dL即将区间x ,x 平移到区间-h,h上，计算所得的定积分的值与原区间上的相 2k-2 2k同。计算在-h,h上过三点 ; ;的抛物线 n；m 为曲边的面积。抛物线二-1中的可由下列方程组确定：员=斫+吐+F.由此得：-： = ； 于是所求面积为：显然，曲边梯形的面积只与二；的纵坐标 -及底边所在的区间的长度2h 有关。由此可知 n/2 组梯形的面积为：启冷心1+牧+曲4冷咫+匹+从十枪】十矗十十人例：用抛物法求JcosMdx.，并与用

11、牛顿莱布尼茨公式计算的结果进行比较。#include #include void main()double result,a=0,b=1,i,n=1000000,h;printf( 按牛顿公式计算得到的结果： %fn,sin(b)-sin(a);result=0;h=(b-a)/n;/ 计算区间高度for(i=1;i=n/2;i+)/ 求和result=result+2*cos(a+2*i*h)+4*cos(a+(2*i-1)*h);result+=cos(a)+cos(b);result=h*result/3;/ 乘以区间高度printf( 用近似公式计算得到的结果： %fn,result)

12、;注意：对于以上三种方法当 n 取得越大时近似程度就越好。练习题： 用三种积分法近似计算如下定积分的值：4. 直接应用 Matlab 命令计算结果(1)数值计算方法 1： int(l/(l+x辽)，x,O,l)(符号求积分)方法2： quad(l./(l+x.2),0,l)(抛物线法求数值积分)方法 3： x=0:0.001:1;y=l./(l+x.2);trapz(x,y)(梯形法求数值积分)数值计算EL*皿方法 1： int(int(x+y2,y,T,l),x,0,2)(符号求积分)方法 2： dblquad(inline(x+y2),0,2,-l,l)(抛物线法二重数值积分)四、自己动手

13、r 1 d?r1. 实现实验内容中的例子，即分别采用矩形法、梯形法、抛物线法计算产, 取科=2方，并比较三种方法的精确程度.r 2 dx2. 分别用梯形法与抛物线法，计算，取归辽并尝试直接使用函数trapz() 、 quad() 进行计算求解，比较结果的差异.3. 试计算定积分(注意：可以运用trapz()、quad ()或附录程序求解吗？为什么？)r 14. 将的近似计算结果与Mat lab中各命令的计算结果相比较，试猜测Matlab 中的数值积分命令最可能采用了哪一种近似计算方法？并找出其他例 子支持你的观点.5 通过整个实验内容及练习，你能否作出一些理论上的小结，即针对什么类 型的函数（具有某种单调特性或凹凸特性），用某种近似计算方法所得结果 更接近于实际值？6 学习 fulu2sum.m 的程序设计方法，尝试用函数 sum 改写附录 1 和附录 3 的程序，避免 for 循环。