第3章微积分实验-4节数值积分

3.4数值积分3.4.1数值积分计算方法 3.4.2误差估计和收敛性

3.4.3数值积分的MATLAB实现3.3.4应用实例3.4.1数值积分计算方法给定函数

，对，有Newton-Leibniz

公式

但是，在下列情况下：函数仅在离散点处给出；被积函数的原函数无法用初等函数表示；被积函数的原函数虽有初等函数表达式，但过于复杂；就必须借助数值方法来求函数的积分3.4.1数值积分计算方法用数值方法近似地求一个函数在区间[a,b]上的定积分的基本思路，可归结到定积分定义：(3-8)3.4.1数值积分计算方法取等距步长，当n充分大时，就是的数值积分，

—第k小区间中x的取值。显然，取值不同，数值积分的结果就不同。这种做法相当于用相对简单阶梯函数

代替作积分。实际上各种不同的数值积分方法就在于，研究用什么样的简单函数代替，使得既能保证结果有一定的精度，计算量又小。3.4.1数值积分计算方法牛顿-柯特斯(Newton-Cote’s)公式设为给定的求积结点，将其作为插值结点，作的拉格朗日插值多项式，然后，利用拉格朗日插值多项式替代作积分，选取不同的多项式，就得到不同的求积公式。设，将区间n等分，记，，以这个等距结点为求积结点的插值型求积公式通常被称为牛顿-柯特斯(Newton-Cote’s)公式3.4.1数值积分计算方法常用的Newton-Cote’s公式：

n=0分别，

和近似可得

几何解释：用以点（左矩形公式）

为顶点的矩形面积近似所要求的积分（矩形面积）右矩形公式左矩形公式中矩形公式3.4.1数值积分计算方法n=1此时，，相当于用

和的算术平均值近似，则得

几何解释：用以点为顶点的梯形面积近似所要求的积分（曲边梯形的面积）梯形公式

图3.10矩形公式、梯形公式的几何意义3.4.1数值积分计算方法3.4.1数值积分计算方法

n=2此时，，，，求积公式为通常称此公式为辛普森（Simpson）公式（抛物线公式）。3.4.1数值积分计算方法图3.11辛普森公式的几何意义3.4.1数值积分计算方法复化的牛顿-柯特斯(Newton-Cote’s)公式为构造高精度的数值积分公式，可采用分段低次多项式替代整体高次多项式，这就导出了复化的牛顿-柯特斯公式，其基本思想是：先把积分区间分成一些长度较小的子区间，在每个子区间上使用低阶的牛顿-柯特斯公式，最常用的是下面的复化梯形公式和复化辛普森公式。设在上，定积分表示曲线下的面积，我们先从图形上看看如何近似计算这块面积。3.4.1数值积分计算方法图3.12定积分的复化矩形公式和复化梯形式3.4.1数值积分计算方法将

n等分，—积分步长。记在每个小区间上用矩形面积近似下面曲边梯形的面积，在整个区间

内构成台阶形。3.4.1数值积分计算方法

容易看出，两个台阶形面积分别为在图3.7中，两个台阶形分别小于和大于所求面积故(3-10)、(3-11)就是计算定积分的复化矩形公式。(3-10)

(3-11)

或将两者平均，则每个小区间上的小矩形变为小梯形，

整个区间上的结果为

视为结点，3-7式相当于用分段线性插值函数作为的近似，称为复化梯形求积公式。(3-12)

3.4.1数值积分计算方法3.4.1数值积分计算方法用分段二次插值函数代替，记，在第k段的两个小区间上，用三个结点作二次插值函数，然后积分，求m段之和可得整个区间上的近似积分

(3-13)3-13式称为复化Simpson求积公式（抛物线公式）。3.4.2误差估计和收敛性误差估计有了求积公式，如何度量它对原积分的近似程度呢？一种方式是考察另一种方式是用使得的函数类的大小来度量。人们称它为求积余项引出了代数精确度的概念3.4.2误差估计和收敛性

（k阶代数精度）

一个数值积分有k阶代数精度：如果当是次数小于或等于k的多项式时，，而对于

k+1次多项式，。即对任意次数不高于k次的多项式，数值积分没有误差。定义：由拉格朗日插值多项式的性质可知，如上构造的求积公式的代数精确度至少是n.梯形公式：求积余项为

Simpson公式：求积余项为

1阶代数精度3阶代数精度3.4.2误差估计和收敛性收敛性

若对的某个数值积分，有（非零常数），则称是p阶收敛的。按此定义可以判断复化梯形公式是2阶收敛，类似，

复化辛普森公式是4阶收敛。定义：3.4.2误差估计和收敛性（p阶收敛）3.4.3数值积分的MATLAB实现对于向量x，cumsum(x)返回一结果向量，此向量x的第n个元素为向量的前个元素之和，如>>x=[1,2,4,-1];>>I=cumsum(x)I=1376矩形求积指令cumsum(x)1.一元函数定积分的数值计算x是由每个小区间左端点的函数值构成的向量

x是由每个小区间右端点的函数值构的向量

x是由每个小区间中间点的函数值构成的向量cumsum(x)*h其中h为子区间步长左矩形公式右矩形公式中矩形公式计算矩形积分公式3.4.3数值积分的MATLAB实现>>x=linspace(0,pi,50);h=1:49；t(h)=x(h);y=sqrt(sin(t).^3-sin(t).^5);T=cumsum(y)*((pi-0)/49);I=T(49)>>I=0.8003（指令x=linspace(a,b,n):在[a,b]区间中的n个等分点（包括端点）构成的向量）3.4.3数值积分的MATLAB实现利用复化左矩形公式计算

用梯形方法计算Y的积分近似值。对于向量Y，Y,应为

此时步长相同且固定h=1。使用梯形法计算Y对X的积分。其中输入量：梯形求积指令trapz；Z=trapz(Y)Z=trapz(X,Y)3.4.3数值积分的MATLAB实现>>x=linspace(0,450);>>formatlong>>y=(x+2)./sqrt(2*x+1)>>T=trapz(x,y)T=

7.3339490782752983.4.3数值积分的MATLAB实现利用trapz计算辛普森公式求积指令quad用simpson公式求函数fun(x)在[a,b]的积分近似值,自动选取步长,绝对误差为1.e-6，输出积分值。同上,但绝对误差为tol。quad(‘fun’,a,b)quad(‘fun’,a,b,tol)3.4.3数值积分的MATLAB实现

用辛普森公式计算积分的值，并估计误差，再用MATLAB验证所求的积分值。解：辛普森公式为此时，。从而有误差为3.4.3数值积分的MATLAB实现MATLAB程序：>>quad('exp(-x)',0,1)ans=0.6321默认误差是，可见与上面所得结果一致。3.4.3数值积分的MATLAB实现几种求积指令的比较

计算积分（很显然这个积分的值是2）编写M文件如下：

x=linspace(0,pi,50);formatlongy=sin(x);T1=cumsum(y(1:49))*(pi/49); %左矩形法T1=T1(49)

T2=trapz(x,y); %梯形法

T3=quad('sin',0,pi); %辛普森法3.4.3数值积分的MATLAB实现%r1r2r3是几种方法的绝对误差r1=abs(T1-2);r2=abs(T2-2);r3=abs(T3-2);

计算结果：r1=6.851506759375514e-004r2=6.851506759371073e-004r3=3.601569042999131e-009

3.4.3数值积分的MATLAB实现计算3.4.3数值积分的MATLAB实现计算近似计算二重积分2.广义积分的数值计算3.重积分的数值计算

3.4.4应用实例人造地球卫星轨道计算问题：人造地球卫星轨道可视为平面上的椭圆。我国第一颗人造地球卫星近地点距地球表面439km，远地点距地球