第七讲定积分的近似计算

第七讲定积分的近似计算第1页，共32页，2023年，2月20日，星期三

问题背景和实验目的定积分的近似计算

定积分计算的基本公式是牛顿－莱布尼兹公式。但当被积函数的原函数不知道时，如何计算？这时就需要利用近似计算。特别是在许多实际应用中，被积函数甚至没有解析表达式，而是一条实验记录曲线，或一组离散的采样值，此时只能用近似方法计算定积分。本实验主要研究定积分的三种近似计算算法：矩形法、梯形法和抛物线法。同时介绍Matlab计算定积分的相关函数。第2页，共32页，2023年，2月20日，星期三

矩形法

定积分的定义：定积分的近似计算第3页，共32页，2023年，2月20日，星期三矩形法n

充分大，xi

充分小

定积分的近似：

通常我们取左点法右点法中点法

点可以任意选取，常见的取法有：

左端点，右端点和中点。第4页，共32页，2023年，2月20日，星期三步长节点

右点法：

中点法：

左点法：左点法、右点法和中点法第5页，共32页，2023年，2月20日，星期三解：矩形法举例==>h=1/100=0.01,xi=i*h,a=0,b=1,n=100例：用不同的矩形法计算下面的定积分(取n=100)，

并比较这三种方法的相对误差。

左点法：

右点法：

中点法：(i=0,1,2,...,100)第6页，共32页，2023年，2月20日，星期三

理论值：左点法相对误差：误差分析矩形法举例右点法相对误差：中点法相对误差：不同的方法有不同的计算精度有没有更好的近似计算定积分的方法

?第7页，共32页，2023年，2月20日，星期三定积分几何意义第8页，共32页，2023年，2月20日，星期三

曲边小梯形的面积可以由直边小梯形的面积来近似整个曲边梯形的面积：梯形法第9页，共32页，2023年，2月20日，星期三

如果我们n

等分区间[a,b]，即令：则==>梯形公式梯形法梯形公式与中点公式有什么区别

?第10页，共32页，2023年，2月20日，星期三解：==>例：用梯形法计算下面定积分(取n=100)，

并计算相对误差梯形法举例a=0,b=1,n=100,f(x)=1/(1+x2)==>h=1/100=0.01,xi=i*h,yi=f(xi)

相对误差：第11页，共32页，2023年，2月20日，星期三

2n

等分区间[a,b]，得该直线用抛物线代替，计算精度是否会更好？

计算每个节点上的函数值：抛物线法

在区间[x0,x2]上，用过以下三点的抛物线来近似原函数f(x)。第12页，共32页，2023年，2月20日，星期三

设过以上三点的抛物线方程为：则在区间[x0,x2]上，有y=

x2+x

+

=p1(x)

抛物线法第13页，共32页，2023年，2月20日，星期三

同理可得：

相加即得：抛物线法第14页，共32页，2023年，2月20日，星期三

整理后可得：或辛卜生(Simpson)公式抛物线法公式抛物线法第15页，共32页，2023年，2月20日，星期三==>例：用抛物线法计算下面定积分(取n=100)，

并计算相对误差解：a=0,b=1,n=100,yi

=f(xi)=1/(1+xi2)

相对误差：抛物线法第16页，共32页，2023年，2月20日，星期三矩形法：cumsumcumsum(x)

x

若为向量，cumsum(x)返回一个向量，此向量的第n个元素为向量x的前n个元素的和；x

若为矩阵，cumsum(x)返回一个矩阵，此矩阵的列为对应的x列的累积和

Matlab近似计算定积分的相关函数Matlab计算定积分函数介绍cumsum(x，dim)

dim决定求和的方向，dim=1，按列的方向求，dim=2，按行的方向求。第17页，共32页，2023年，2月20日，星期三>>x=[012;345]；>>cumsum(x)ans=012357>>cumsum(x,1)ans=012357>>cumsum(x,2)ans=0133712a=[123]；>>cumsum(a,1)ans=123>>cumsum(a,2)ans=136>>cumsum(a)ans=136第18页，共32页，2023年，2月20日，星期三梯形法：trapztrapz(x,y)

x

为分割点（节点）组成的向量，

y为被积函数在节点上的函数值组成的向量。第19页，共32页，2023年，2月20日，星期三前面的做法例：用梯形法计算下面定积分(取n=100解：a=0,b=1,n=100,yi

=f(xi)=1/(1+xi2)>>

x=0:1/100:1;>>

y=1./(1+x.^2);>>

trapz(x,y)trapz函数trapz(x,1./(1+x.^2))trapz举例第20页，共32页，2023年，2月20日，星期三quad(f,a,b,tol)f=f(x)为被积函数，[a,b]为积分区间，tol

为计算精度将自变量看成是向量抛物线法：quad不用自己分割积分区间可以指定计算精度，若不指定，缺省精度是10-6精度越高，函数运行的时间越长此处的函数

f是数值形式，应该使用数组运算，即

点运算：.*，./，.\，.^

注：抛物线法第21页，共32页，2023年，2月20日，星期三解：>>

quad('1./(1+x.^2)',0,1)>>

quad('1./(1+x.^2)',0,1,10e-10)>>

quad('1./(1+x.^2)',0,1,10e-16)函数表达式一定要用单引号括起来！涉及的运算一定要用数组运算！例：用quad

计算定积分：quad举例第22页，共32页，2023年，2月20日，星期三fun.mfunctionf=fun(x)f=x./(4+x.^2);quad(‘fun’，0，1）quad('x./(4+x.^2);',0,1)第23页，共32页，2023年，2月20日，星期三抛物线法计算二重积分：dblquaddblquad(f,a,b,c,d,tol)

tol

为计算精度，若不指定，则缺省精度为10-6

f(x,y)

可以由inline

定义，或通过一个函数句柄传递

[a,b]

是第一积分变量的积分区间，[c,d]

是第二积分变量

的积分区间按字母顺序，大写字母排在小写字母的前面二重积分的计算第24页，共32页，2023年，2月20日，星期三>>

f=inline('4*x*y+3*y^2');>>

I=dblquad(f,-1,1,0,2)

f(x,y)

中关于第一自变量的运算是数组运算，

即把x

看成是向量，y

看成是标量。也可以全部采用数组运算例2：计算二重积分>>

dblquad(inline('4*x*y+3*x^2'),-1,1,0,2)>>

dblquad(inline('4*x*y+3*x.^2'),-1,1,0,2)X例1：计算二重积分dblquad举例第25页，共32页，2023年，2月20日，星期三例：计算二重积分>>

dblquad(@(x,y)4*x*y+3*x.^2,-1,1,0,2)指定x、y

分别是第一和第二积分变量>>

dblquad(inline('4*x*y+3*x.^2'),-1,1,0,2)被积函数f(x,y)

的另一种定义方法：匿名函数>>

dblquad(@(y,x)4*x*y+3*x.^2,-1,1,0,2)下面的命令运行结果和上面的一样吗？dblquad举例第26页，共32页，2023年，2月20日，星期三int(f,a,b)

计算

f

关于默认自变量

的定积分，积分区间为[a,b]。int(f)

计算

f

关于默认自变量

的不定积分。int(f,v,a,b)

计算函数f

关于自变量v

的定积分，积分区间为[a,b]int(f,v)

计算函数

f

关于自变量

v

的不定积分findsym(f,1)符号积分：

intint符号积分第27页，共32页，2023年，2月20日，星期三>>

symsxy;>>

f=y*sin(x);>>

int(f,x)>>

int(f,y)>>

int(f)>>

int('a+b')ans=-y*cos(x)ans=1/2*y^2*sin(x)ans=-y*cos(x)ans=a*b+1/2*b^2例：指出下面各条语句的输出结果int举例第28页，共32页，2023年，2月20日，星期三例：用int

函数计算定积分：解：>>

symsx;>>

f=1/(1+x^2);>>

int(f,x,0,1)>>

f=sym('1/(1+x^2)');>>

int(f,x,0,1)>>

int('1/(1+x^2)',x,0,1)或>>

int('1/(1+x^2)',0,1)或或int举例第29页，共32页，2023年，2月20日，星期三double(a)将a

转化为双精度型，若a

是字符，则取对应的ASCII码>>

a=3;>>

double(a)>>

double('a')例：ans=3ans=97其它相关函数第30页，共32页，2023年，2月20日，星期三>>

x=1:0.001:2;>>

y=exp(x.^(-2));>>

trapz(x,y)

梯形法：

抛物线法：>>

quad('exp(x.^(-2))',1,2,10e-10)

符号积分法：>>