数值积分的matlab实现

实验10 数值积分

实验目的：

1．了解数值积分的基本原理；

2．熟练掌握数值积分的MATLAB 实现； 3．会用数值积分方法解决一些实际问题。 实验内容：

积分是数学中的一个基本概念，在实际问题中也有很广泛的应用。同微分一样，在《微积分》中，它也是通过极限定义的，由于实际问题中遇到的函数一般都以列表形式给出，所以常常不能用来直接进行积分。此外有些函数虽然有解析式，但其原函数不是初等函数，所以仍然得不到积分的精确值，如不定积分?1

0 d sin x x x

。这时我们一般考虑用数值方法计算其

近似值，称为数值积分。

10.1 数值微分简介

设函数()y f x =在*

x 可导，则其导数为

h

x f h x f x f h )

()(lim )(**0*

-+='→ （10.1）

如果函数()y f x =以列表形式给出（见表10-1），则其精确值无法求得，但可由下式求得其近似值

h

x f h x f x f )

()()(***

-+≈' （10.2）

一般的，步长h 越小，所得结果越精确。（10.2）式右端项的分子称为函数()y f x =在

*x 的差分，分母称为自变量在*x 的差分，所以右端项又称为差商。数值微分即用差商近似

代替微商。常用的差商公式为：

000()()

()2f x h f x h f x h

+--'≈

（10.3）

h y y y x f 243)(2

100-+-≈

' （10.4）

h

y y y x f n

n n n 234)(12+-≈

'-- （10.5）

其误差均为2

()O h ，称为统称三点公式。

10.2 数值微分的MATLAB 实现

MATLAB 提供了一个指令求解一阶向前差分，其使用格式为： dx=diff(x) 其中x 是n 维数组，dx 为1n -维数组[]21321,,

,n x x x x x x ---,这样基于两点的数值导

数可通过指令diff(x)/h 实现。对于三点公式，读者可参考例1的M 函数文件diff3.m 。

例1 用三点公式计算()y f x =在=x 1.0,1.2,1.4处的导数值，()f x 的值由下表给出。

解：建立三点公式的M 函数文件diff3.m 如下：

function f=diff3(x,y) n=length(x);h=x(2)-x(1); f(1)=(-3*y(1)+4*y(2)-y(3))/(2*h); for j=2:n-1

f(j)=(y(j+1)-y(j-1))/(2*h); end

f(n)=(y(n-2)-4*y(n-1)+3*y(n))/(2*h);

在MATLAB 指令窗中输入指令：

x=[1.0,1.1,1.2,1.3,1.4];y=[0.2500,0.2268,0.2066,0.1890,0.1736];diff3(x,y)

运行得各点的导数值为：-0.2470，-0.2170，-0.1890，-0.1650，-0.0014。所以()y f x =在=x 1.0,1.2,1.4处的导数值分别为-0.2470，-0.1890和-0.0014。

对于高阶导数，MATLAB 提供了几个指令借助于样条函数进行求导，详细使用步骤如下：

step1：对给定数据点（x,y ），利用指令pp=spline(x,y)，获得三次样条函数数据pp ，供后面ppval 等指令使用。其中，pp 是一个分段多项式所对应的行向量，它包含此多项式的阶数、段数、节点的横坐标值和各段多项式的系数。

step2：对于上面所求的数据向量pp ，利用指令[breaks,coefs,m,n]=unmkpp(pp)进行处理，生成几个有序的分段多项式pp 。

step3：对各个分段多项式pp 的系数，利用函数ppval 生成其相应导数分段多项式的系数，再利用指令mkpp 生成相应的导数分段多项式

step4：将待求点xx 代入此导数多项式，即得样条导数值。 上述过程可建立M 函数文件ppd.m 实现如下：

function dy=ppd(pp)

[breaks,coefs,m]=unmkpp(pp); for i=1:m

coefsm(i,:)=polyder(coefs(i,:)); end

dy=mkpp(breaks,coefsm);

于是，如果已知节点处的值x,y ，可用下面指令计算xx 处的导数dyy ：

pp=spline(x,y),dy=ppd(pp);dyy=ppval(dy,xx);

例2 基于正弦函数sin y x 的数据点，利用三点公式和三次样条插值分别求导，并与解析所求得的导数进行比较。

解:编写M 脚本文件bijiao.m 如下：

h=0.1*pi;x=0:h:2*pi;y=sin(x); dy1=diff3(x,y);

pp=spline(x,y);dy=ppd(pp);dy2=ppval(dy,x); z=cos(x);

error1=norm(dy1-z),error2=norm(dy2-z) plot(x,dy1,'k:',x,dy2,'r--',x,z,'b')

运行得结果为：error1 =0.0666，error2 =0.0025，生成图形见图10.1。

图10.1 三点公式、三次样条插值与解析求导比较图

显然利用三次样条插值求导所得误差比三点公式求导小很多，同时由图2.15可知利用三次样条插值求导所得曲线与解析求导曲线基本重合，而三点公式在极值点附近和两个端点附近误差较大，其它点吻合的较好。

10.3 应用示例：湖水温度变化问题

问题：湖水在夏天会出现分层现象，其特点是接近湖面的水的温度较高，越往下水的温度越低。这种现象会影响水的对流和混合过程，使得下层水域缺氧，导致水生鱼类死亡。

对某个湖的水温进行观测得数据见表10-2。

试找出湖水温度变化最大的深度。 1．问题的分析

湖水的温度可视为关于深度的函数，于是湖水温度的变化问题便转化为温度函数的导数问题，显然导函数的最大绝对值所对应的深度即为温度变化最大的深度。对于给定的数据，可以利用数值微分计算各深度的温度变化值，从而得到温度变化最大的深度，但考虑到所给的数据较少，由此计算的深度不够精确，所以采用插值的方法计算加密深度数据的导数值，以得到更准确的结果。

2．模型的建立及求解

记湖水的深度为h （m ），相应的温度为T （℃），且有)(h T T ，并假定函数)(h T 可导。

对给定的数据进行三次样条插值，并对其求导，得到)(h T 的插值导函数；然后将给定的深度数据加密，搜索加密数据的导数值的绝对值，找出其最大值及其相应的深度，相应的MATLAB 指令如下：

h=[0 2.3 4.9 9.1 13.7 18.3 22.9 27.2];T=[22.8 22.8 22.8 20.6 13.9 11.7 11.1 11.1]; hh=0:0.1:27.2;

pp=spline(h,T);dT=ppd(pp);dTT=ppval(dT,hh); [dTTmax,i]=max(abs(dTT)),hh(i) plot(hh,dTT, 'b ',hh(i),dTT(i), 'r. '),grid on

运行得导函数绝对值的最大值点为：h =11.4，

最大值为1.6139,即湖水在深度为11.4m

时温度变化最大，如图10.2所示（黑点为温度变化最大的点）。

图10.2 湖水温度变化曲线图

10.4 数值积分简介

考虑定积分

()d b

a

f x x ?

（10.6）

如果被积函数()f x 是以列表形式给出，则其求解思想同数值微分类似，即用逼近多项式()n P x 近似地代替被积函数()f x ，然后计算积分

()d b

n a

P x x ?

，得（10.6）式的近似值；

如果被积函数的原函数不是初等函数，则将积分区间进行细分，对每个小区间，用一个近似函数代替被积函数()f x ，然后积分得（10.6）式的近似值。这两种类型最终都可归结为函数()f x 在节点k x 上的函数值()k f x 的某种线性组合，即下面数值求积公式： 0()d ()n

b

k k a

k I f x x A f x ==

≈∑?

或 （10.7）

[]0

()d ()n

b

k k a

k f x x A f x R f ==+∑?

（10.8）

其中[]R f 为截断误差。此误差可用代数精度衡量，代数精度越高，误差越小；反之误差越大。

代数精度是用来衡量数值积分公式近似程度的办法，如果)(x f 是一个次数不超过m 的代数多项式，（10.7）式等号成立；而当)(x f 是一个1+m 次多项式时，（10.7）式不能精确成立，则称（10.7）式的代数精度为m 。

选取不同的近似函数，可产生不同的数值求积公式，常见的有：梯形公式、辛普森公式和高斯公式。

10.5 数值积分的MATLAB 实现

MATLAB 提供了下面几个函数计算积分，其使用格式分别为：

（1）trapz(x) 采用梯形公式计算积分（1=h ）,x 为),,1,0(n k f k = （2）quad('fun',a,b,tol) 采用自适应Simpson 法计算积分 （3）quadl('fun',a,b,tol) 采用自适应Gauss-Lobatto 法计算积分 其中fun 为被积函数；tol 是可选项，表示绝对误差，a ，b 为积分的上、下限。

例1分别利用梯形公式、Simpson 公式和Gauss-Lobatto 法计算?

+1

2d 1x x ，并与

其精确值比较。

解：先对积分作符号运算，然后将其计算结果转换为数值型，再将其与这三种方法求得的数值解比较，其MATLAB 指令为：

syms xx

z0=simple(int('sqrt(1+xx^2)',0,1)) z=double(z0);z=vpa(z,8) x=0:0.01:1;y=sqrt(1+x.^2);

z1=trapz(y)*0.01;z1=vpa(z1,8),err1=z-z1;err1=vpa(err1,8)

z2=quad('sqrt(1+x.^2)',0,1);z2=vpa(z2,8),err2=z-z2;err2=vpa(err2,8) z3=quadl('sqrt(1+x.^2)',0,1);z3=vpa(z3,8),err3=z-z3;err3=vpa(err3,8)

运行得精确值为

=--))12ln(2(2

1

1.1477936，三种公式计算得数值积分值分别为1.1477995，1.1477935和1.1477936，其相应误差分别为-.59e-5，.1e-6和0.，由三者误差可见，Gauss-Lobatto 法计算最为精确，Simpson 公式次之，梯形公式最差，但它也

能精确到小数点后5位数。

例2人造地球卫星轨道可视为平面上的椭圆。我国第一颗人造地球卫星近地点距地球表面439km ，远地点距地球表面2384km ，地球半径为6371km ，求该卫星的轨道长度。

解：卫星轨道椭圆的参数方程为cos ,sin (02),,x a t y b t t a b π==≤≤分别是长、短半轴，则根据所给数据知a =6371+2384=8755,b =6371+439=6810。

由对弧长的曲线积分知识知，椭圆的长度为

12222

22

4

(sin cos )d L a t b t t π

=+?

上积分称为椭圆积分，它无法用解析方法计算，可用计算其数值解，编写M 函数文件如下：

function y=y(t) a=8755;b=6810;

y=4*sqrt(a^2*sin(t).^2+b^2*cos(t).^2);

在MATLAB 指令窗中输入以下指令：

l=quad('y',0,pi/2)

运行得结果为：l=4.9090e+004。即卫星轨道长度为49090km 。

对于用列表形式给出的函数，上述方法不再适用，可利用指令spline 构造三次样条插值函数，再计算积分，具体步骤可参考例2。

例3 在桥梁的一端每隔一段时间记录1min 有几辆车过桥，得到数据见表10-3：

试估计一天通过桥梁的车流量。

解：记记录时刻为t 时，相应的车辆数为)(t x ，则所求车流量即为计算积分?

24

d )(t t x ，

则在MATLAB 指令窗中输入下面指令：

x=[0,2,4,5,6,7,8,9,10.3,11.3,12.3,14,16,17,18,19,20,21,22,23,24]; y=[2,2,0,2,5,8,25,12,5,10,12,7,9,28,22,10,9,11,8,9,3];

pp=spline(x,y);s1=quadl('fun',0,24,[],[],pp) %利用三次样条插值计算积分 s2=trapz(x,y) %利用梯形公式计算积分

其中M 函数文件fun.m 为：

function vf=fun(x,pp) vf=ppval(pp,x);

运行得三次样条插值计算所得车流量为212辆，梯形公式计算所得车流量为216辆。

10.6 数值积分的建模示例：煤炭储量计算问题

问题：某煤矿为估计其煤炭的储量，在该矿区内进行勘探，得到数据如表10-4。试估算出该矿区（51,41≤≤≤≤y x ）煤炭的储量。

1．问题的分析

问题给出了很多点对应的煤炭厚度，显然整个煤矿可以看作是一个巨大的曲顶柱体（见图10.3，此图经过插值得到），而煤炭的储量即为此立体的体积。要计算此立体的体积，可

以利用插值得到曲顶柱体的顶面函数，再对其积分；也可以将此曲顶柱体分割成若干个细的曲顶柱体，用数值方法计算这些细曲顶柱体的体积，再对其求和即得原曲顶柱体的体积。

图10.3 煤炭厚度曲面图 2．模型的建立及求解

以煤炭的厚度为三维空间中的z 坐标建立空间坐标系。记煤炭的厚度为h ，则它是坐标

y x ,的二元函数，即),(y x z ?=，则由二重积分的知识可知，此煤矿的煤炭储量W 为

??=D

y x y x W d d ),(? （10.9）

其中}51,41|),{(≤≤≤≤=y x y x D 。

由于函数),(y x ?只给出了一些离散点上的函数值，无法直接计算上述二重积分，所以下面采用数值积分的方法计算其值。

由数值积分的知识知，计算定积分有复合梯形公式为

])(2)()([2d )(1

1

∑?

-=++≈n k k b

a

x f b f a f h

x x f （10.10）

其中h 为步长，),,1,0(n k x k =为节点，且有kh a x k +=。

由（10.9）式得

???=??

????=

b a b

a d c x g(x)x y y x W d d d ),(? （10.11）

其中?

=

d

c

y y x x g d ),()(?，则由（10.10）式可得

])(2)()([2d )(1

1

∑?

-=++≈=

n j j x

b

a

x g b g a g h x x g W （10.12）

而

]),(2),(),([2d ),()(1

1

∑?

-=++≈

=

m k k y d

c

y a d a c a h y y a a g ????

]),(2),(),([2

d ),()(1

1

∑?-=++≈

=m k k y d

c y b

d b c b h y y b b g ????

]),(2),(),([2

d ),()(1

1

∑?-=++≈

=m k k k k k y d

c

k k y x d x c x h y y x x g ????

所以有

[]

[]∑∑∑∑∑∑∑====-=-=-=+++++----=

++++++++≈

n

j m

k k j n

j j j m k k k y x n j m k k j j j m k k k y x y x d x c x y b y a d b c b d a c a h h y x d x c x y b y a d b c b d a c a h h W 00

1

1

11

1

1

)

,(4)],(),([2),(),(2)

,(5),(5),(5),(5[4

]

]),(2),(),([2),(),(2),(),(),(),([4

??????????????????（10.13）

考虑到给定的数据较少，由此产生的误差较大，所以利用插值后的数据计算（10.13）式，相应的MATLAB 计算指令如下：

x=[1 1 1 1 1 2 2 2 2 2 3 3 3 3 3 4 4 4 4 4]*1000;y=[1 2 3 4 5 1 2 3 4 5 1 2 3 4 5 1 2 3 4 5]*1000;

z=-[13.72 25.80 8.47 25.27 22.32 15.47 21.33 14.49 24.83 26.19 23.28 26.48 29.14

12.04 14.58 19.95 23.73 15.35 18.01 16.29]; hx=10;hy=10;cx=1000:hx:4000;cy=1000:hy:5000; [X,Y]= meshgrid(cx,cy);n=length(cx);m=length(cy);

Z=griddata(x,y,z,X,Y,'v4'); %插值 surf(X,Y,Z) %绘制图10.3

W=-hx*hy*(-5*Z(1,1)-5*Z(1,n)-5*Z(m,1)-5*Z(m,n)+2*(sum(Z(1,:)+Z(m,:))+sum(Z(:,1)

+Z(:,n)))+4*sum(sum(Z)))/4

运行得W =2.5242×108，即煤矿的煤炭储量约为2.5242×108m 3. 实验任务：

1．一个物体的运动距离)(t D D =如下表所示：

（1）利用三点公式和三次样条插值方法分别计算此物体在各个时刻的运动速率)(t V ； （2）与)(t D 的解析式10

/70770)(t e

t t D -++-=比较计算结果。

2．已知20世纪美国人口统计数据如表10-5，试计算1900-1990年各年份的人口相对增长率，并以此分析美国在这些年的人口变化情况。

3．计算由表10-6数据给出的积分

?

5

.1 1

.0 d )(x x y 。

4．已知某铸件为曲顶柱体，现要对一个铸件的曲顶面进行涂漆，单位表面的费用为100元，该铸件曲顶面数据如表10-7所示，试估算该项处理的费用。

数值积分算法与MATLAB实现陈悦5133201讲解

东北大学秦皇岛分校 数值计算课程设计报告 数值积分算法及MATLAB实现 学院数学与统计学院 专业信息与计算科学 学号5133201 姓名陈悦 指导教师姜玉山张建波 成绩 教师评语： 指导教师签字： 2015年07月14日

1 绪论 数值分析是计算数学的一个主要部分，计算数学是数学科学的一个分支，它研究用计算机求解各种数学问题的数值检索方其理论与软件的实现.而数值分析主要研究数值计算. 现科学技术的发展与进步提出了越来越多的复杂的数值计算问题，这些问题的圆满解决已远人工手算所能胜任，必须依靠电子计算机快速准确的数据处理能力.这种用计算机处理数值问题的方法，成为科学计算.今天，科学计算的应用范围非常广泛，天气预报、工程设计、流体计算、经济规划和预测以及国防尖端的一些科研项目，如核武器的研制、导弹和火箭的发射等，始终是科学计算最为活跃的领域. 1.1 数值积分介绍 数值积分是数值分析的重要环节，实际问题当中常常需要计算积分，有些数值方法，如微分方程和积分方程的求解，也都和积分计算相联系. 求某函数的定积分时，在多数情况下，被积函数的原函数很难用初等函数表达出来，因此能够借助微积分学的牛顿-莱布尼兹公式计算定积分的机会是不多的.另外，许多实际问题中的被积函数往往是列表函数或其他形式的非连续函数，对这类函数的定积分，也不能用不定积分方法求解.由于以上原因，数值积分的理论与方法一直是计算数学研究的基本课题.对微积分学做出杰出贡献的数学大师，如I.牛顿、L.欧拉、C.F.高斯、拉格朗日等人都在数值积分这个领域作出了各自的贡献，并奠定了这个分支的理论基础. 构造数值积分公式最通常的方法是用积分区间上的n 次插值多项式代替被积函数，由此导出的求积公式称为插值型求积公式.特别在节点分布等距的情形称为牛顿-科特斯公式，例如梯形公式(Trapezoidal Approximations)与抛物线公式(Approximations Using Parabolas)就是最基本的近似公式.但它们的精度较差.龙贝格算法是在区间逐次分半过程中，对梯形公式的近似值进行加权平均获得准确程度较高的积分近似值的一种方法，它具有公式简练、计算结果准确、使用方便、稳定性好等优点，因此在等距情形宜采用龙贝格求积公式(Rhomberg Integration).当用不等距节点进行计算时，常用高斯型求积公式计算，它在节点数目相同情况下，准确程度较高，稳定性好，而且还可以计算无穷积分.数值积分还是微分方程数值解法的重要依据.许多重要公式都可以用数值积分方程导出.现探讨数值积分算法以及运用MATLAB软件的具体实现

数值积分的matlab实现

实验10 数值积分 实验目的： 1．了解数值积分的基本原理； 2．熟练掌握数值积分的MATLAB 实现； 3．会用数值积分方法解决一些实际问题。 实验内容： 积分是数学中的一个基本概念，在实际问题中也有很广泛的应用。同微分一样，在《微积分》中，它也是通过极限定义的，由于实际问题中遇到的函数一般都以列表形式给出，所以常常不能用来直接进行积分。此外有些函数虽然有解析式，但其原函数不是初等函数，所以仍然得不到积分的精确值，如不定积分?1 0 d sin x x x 。这时我们一般考虑用数值方法计算其 近似值，称为数值积分。 10.1 数值微分简介 设函数()y f x =在* x 可导，则其导数为 h x f h x f x f h ) ()(lim )(**0* -+='→ （10.1） 如果函数()y f x =以列表形式给出（见表10-1），则其精确值无法求得，但可由下式求得其近似值 h x f h x f x f ) ()()(*** -+≈' （10.2） 表 10-1 一般的，步长h 越小，所得结果越精确。（10.2）式右端项的分子称为函数()y f x =在 *x 的差分，分母称为自变量在*x 的差分，所以右端项又称为差商。数值微分即用差商近似 代替微商。常用的差商公式为： 000()() ()2f x h f x h f x h +--'≈ （10.3） h y y y x f 243)(2 100-+-≈ ' （10.4）

h y y y x f n n n n 234)(12+-≈ '-- （10.5） 其误差均为2 ()O h ，称为统称三点公式。 10.2 数值微分的MATLAB 实现 MATLAB 提供了一个指令求解一阶向前差分，其使用格式为： dx=diff(x) 其中x 是n 维数组，dx 为1n -维数组[]21321,, ,n x x x x x x ---,这样基于两点的数值导 数可通过指令diff(x)/h 实现。对于三点公式，读者可参考例1的M 函数文件diff3.m 。 例1 用三点公式计算()y f x =在=x 1.0,1.2,1.4处的导数值，()f x 的值由下表给 解：建立三点公式的M 函数文件diff3.m 如下： function f=diff3(x,y) n=length(x);h=x(2)-x(1); f(1)=(-3*y(1)+4*y(2)-y(3))/(2*h); for j=2:n-1 f(j)=(y(j+1)-y(j-1))/(2*h); end f(n)=(y(n-2)-4*y(n-1)+3*y(n))/(2*h); 在MATLAB 指令窗中输入指令： x=[1.0,1.1,1.2,1.3,1.4];y=[0.2500,0.2268,0.2066,0.1890,0.1736];diff3(x,y) 运行得各点的导数值为：-0.2470，-0.2170，-0.1890，-0.1650，-0.0014。所以()y f x =在=x 1.0,1.2,1.4处的导数值分别为-0.2470，-0.1890和-0.0014。 对于高阶导数，MATLAB 提供了几个指令借助于样条函数进行求导，详细使用步骤如下： step1：对给定数据点（x,y ），利用指令pp=spline(x,y)，获得三次样条函数数据pp ，供后面ppval 等指令使用。其中，pp 是一个分段多项式所对应的行向量，它包含此多项式的阶数、段数、节点的横坐标值和各段多项式的系数。 step2：对于上面所求的数据向量pp ，利用指令[breaks,coefs,m,n]=unmkpp(pp)进行处理，生成几个有序的分段多项式pp 。 step3：对各个分段多项式pp 的系数，利用函数ppval 生成其相应导数分段多项式的系数，再利用指令mkpp 生成相应的导数分段多项式 step4：将待求点xx 代入此导数多项式，即得样条导数值。 上述过程可建立M 函数文件ppd.m 实现如下： function dy=ppd(pp) [breaks,coefs,m]=unmkpp(pp);

用递推公式计算定积分(matlab版)

用递推公式计算定积分 实验目的： 1．充分理解不稳定的计算方法会造成误差的积累，在计算过程中会导致误差的迅速增加，从而使结果产生较大的误差。 2．在选择数值计算公式来进行近似计算时，应学会选用那些在计算过程中不会导致误差迅速增长的计算公式。 3．理解不稳定的计算公式造成误差积累的来源及具体过程； 4．掌握简单的matlab语言进行数值计算的方法。 实验题目： 对n=0,1,2,…,20，计算定积分： 实验原理： 由于y(n)= = – 在计算时有两种迭代方法，如下： 方法一： y(n)=– 5*y(n-1),n=1,2,3, (20) 取y(0)= = ln6-ln5 ≈ 0.182322 方法二： 利用递推公式：y(n-1)=-*y(n),n=20,19, (1) 而且，由 = * ≤≤* =

可取：y(20)≈*()≈0.008730. 实验容： 对算法一,程序代码如下： function [y,n]=funa() syms k n t; t=0.182322; n=0; y=zeros(1,20); y(1)=t; for k=2:20 y(k)=1/k-5*y(k-1); n=n+1; end y(1:6) y(7:11) 对算法二，程序代码如下： %计算定积分； %n--表示迭代次数； %y用来存储结果； function [y,n]=f(); syms k y_20;

y=zeros(21,1); n=1; y_20=(1/105+1/126)/2; y(21)=y_20; for k=21:-1:2 y(k-1)=1/(5*(k-1))-y(k)/5; n=n+1; end 实验结果： 由于计算过程中，前11个数字太小，后9个数字比较大，造成前面几个数字只显示0.0000的现象，所以先输出前6个，再输出7—11个，这样就能全部显示出来了。 算法一结果： [y,n]=funa %先显示一y(1)—y(6) ans = 0.1823 -0.4116 2.3914 -11.7069 58.7346

【数学建模学习】Matlab的数值积分

Matlab 的数值积分问题 （１）求和命令sum 调用格式. 如果x 是向量,则sum(x) 给出x 的各个元素的累加和;如果x 是矩阵,则sum(x)是一个元素为x 的每列列和的行向量. 例3.1 调用命令sum 求向量x 的各个元素的累加和。 解:输入 x=[1,2,3,4,5,6,7,8,9,10]; sum(x) 得到 ans=55 例3.2 调用命令sum 求矩阵x 的各列元素的累加和。 解:输入 x=[1,2,3;4,5,6;7,8,9] x= 1 2 3 4 5 6 7 8 9 sum(x) 得到 ans=12 15 18 2．定积分的概念. 定积分是一个积分和的极限. 例如取x e x f =)(,求定积分?10dx e x 的近似值。 积分区间为[0,1],等距划分为20个子区间, x=linspace(0,1,21); 选取每个子区间的端点,并计算端点处的函数值. y=exp(x); 取区间的左端点处的函数值乘以区间长度全部加起来. y1=y(1:20); s1=sum(y1)/20 s1=1.6757 s1可作为定积分?10dx e x 的近似值。 若选取右端点: y2=y(2:21); s2=sum(y2)/20 s2=1.7616 s2也可以作为定积分?10dx e x 的近似值。 下面我们画出图象. plot(x,y);hold on for i=1:20 fill([x(i),x(i+1),x(i+1),x(i),x(i)],[0,0,y(i),y(i),0],'b')

end 如果选取右端点,则可画出图象. for i=1:20 fill([x(i),x(i+1),x(i+1),x(i),x(i)],[0,0,y(i+1),y(i+1),0],'b') hold on end plot(x,y,'r') 在上边的语句中,for … end 是循环语句,执行语句体内的命令20次,fill 命令可以填充多边形,在本例中,用的是兰色(blue)填充. 可试取50个子区间看一看结果怎样.下面按等分区间计算。 syms k n s=symsum(exp(k/n)/n,k,1,n); limit(s,n,inf) 得结果 ans=exp(1)-1 3．计算定积分 例3.6 计算?10dx e x . 解:输入命令: syms x; int(exp(x),0,1) 得结果 ans=exp(1)-1. 这与我们上面的运算结果是一致的. ⒈ 由给定数据进行梯形求积 假设已经建立起向量T N T N y y y y x x x x ],,,[,],,,[2121 ==，则可用以下语句进行梯形求积： sum((2*y(1:end-1,:)+diff(y)).*diff(x))/2 MATLAB 提供的trapz()函数也可直接用梯形法求解积分问题，该函数调用格式为 S=trapz(x,y) [例1-6-17] 试用梯形法求出),0(π∈x 区间内，函数sin(x),cos(x),sin(x/2)的定积分值。 [求解] >> x1=[0:pi/30:pi]'; y=[sin(x1) cos(x1) sin(x1/2)]; x=[x1 x1 x1]; S=sum((2*y(1:end-1,:)+diff(y)).*diff(x))/2 >> S1=trapz(x1,y) [例1-6-18] 用定步长方法求解积分?2 /30)15cos(πdx x 。 [求解] 鉴于求解区域内被积函数有很强的振荡，可先用下述语句绘制被积函数的曲线。 >> x=[0:0.01:3*pi/2,3*pi/2]; y=cos(15*x); plot(x,y) 采用不同的步距，可分别得到积分近似结果。 >> syms x, A=int(cos(15*x),0,3*pi/2) % 求理论值 >> h0=[0.1,0.01,0.001,0.0001,0.00001,0.000001]; v=[]

matlab实现数值分析插值及积分

Matlab实现数值分析插值及积分 摘要： 数值分析(numerical analysis)是研究分析用计算机求解数学计算问题的数值计算方法及其理论的学科，是数学的一个分支，它以数字计算机求解数学问题的理论和方法为研究对象。在实际生产实践中，常常将实际问题转化为数学模型来解决，这个过程就是数学建模。学习数值分析这门课程可以让我们学到很多的数学建模方法。 分别运用matlab数学软件编程来解决插值问题和数值积分问题。题目中的要计算差值和积分，对于问题一，可以分别利用朗格朗日插值公式，牛顿插值公式，埃特金逐次线性插值公式来进行编程求解，具体matlab代码见正文。编程求解出来的结果为：=+。 其中Aitken插值计算的结果图如下： 对于问题二，可以分别利用复化梯形公式，复化的辛卜生公式，复化的柯特斯公式编写程序来进行求解，具体matlab代码见正文。编程求解出来的结果为： 0.6932 其中复化梯形公式计算的结果图如下：

问题重述 问题一：已知列表函数 表格 1 分别用拉格朗日，牛顿，埃特金插值方法计算。 问题二：用复化的梯形公式，复化的辛卜生公式，复化的柯特斯公式计算积分，使精度小于5。 问题解决 问题一：插值方法 对于问题一，用三种差值方法：拉格朗日，牛顿，埃特金差值方法来解决。 一、拉格朗日插值法： 拉格朗日插值多项式如下： 首先构造1+n 个插值节点n x x x ,,,10 上的n 插值基函数，对任一点i x 所对应的插值基函数 )(x l i ，由于在所有),,1,1,,1,0(n i i j x j +-=取零值，因此)(x l i 有因子 )())(()(110n i i x x x x x x x x ----+- 。又因)(x l i 是一个次数不超过n 的多项式，所以只 可能相差一个常数因子，固)(x l i 可表示成： )())(()()(110n i i i x x x x x x x x A x l ----=+- 利用1)(=i i x l 得：

数值积分用matlab实现

数值积分用m a t l a b实 现

东北大学秦皇岛分校 数值计算课程设计报告 数值积分及Matlab实现 学院数学与统计学院 专业信息与计算科学 学号5133117 姓名楚文玉 指导教师张建波姜玉山 成绩 教师评语： 指导教师签字： 2015年07月14日

1 绪论 在科研计算中，经常会碰到一些很难用公式定理直接求出精确解的积分问题，对于这类问题，我们一般转化为数值积分问题，用计算机来实现求解问题． 1.1 课题的背景 对于定积分()b a f x dx ?在求某函数的定积分时，在一定条件下，虽然有牛顿-莱布里 茨公式()()()b a I f x dx F b F a ==-?可以计算定积分的值，但在很多情况下的原函数() f x 不易求出或非常复杂．被积函数的原函数很难用初等函数表达出来，例如 2 sin (),x x f x e x -= 等；有的函数()f x 的原函数()F x 存在，但其表达式太复杂，计算量太大，有的甚至无法有解析表达式．因此能够借助牛顿-莱布尼兹公式计算定积分的情形是不多的．另外，许多实际问题中的被积函数()f x 往往是列表函数或其他形式的非连续函数，对这类函数的定积分，也不能用不定积分方法求解，只能设法求其近似值．因此，探讨近似计算的数值积分方法是有明显的实际意义的，即有必要研究定积分的数值计算方法，以解决定积分的近似计算．而数值积分就是解决此类问题的一种有效的方法，它的特点是利用被积函数在一些节点上的信息求出定积分的近似值．微积分的发明是人类科学史上一项伟大的成就，在科学技术中，积分是经常遇到的一个重要计算环节数值积分是数学上重要的课题之一，是数值分析中重要的内容之一．随着计算机的出现，近几十年来，对于数值积分问题的研究已经成为一个很活跃的研究领域．现在，数值积分在计算机图形学，积分方程，工程计算，金融数学等应用科学领域都有着相当重要的应用，所以研究数值积分问题有着很重要的意义．国内外众多学者在数值积分应用领域也提出了许多新方法．在很多实际应用中，只能知道积分函数在某些特定点的取值，比如天气测量中的气温、湿度、气压等，医学测量中的血压、浓度等等．通过这个课题的研究，我们将会更好地掌握运用数值积分算法求出特殊积分函数的定积分的一些基本方法、理论基础；并且通过Matlab 软件编程的实现，应用于实际生活中． 1.2 课题的主要内容框架

数值积分算法与MATLAB实现

数值积分算法与MATLAB实现 本文从网络收集而来，上传到平台为了帮到更多的人，如果您需要使用本文档，请点击下载按钮下载本文档（有偿下载），另外祝您生活愉快，工作顺利，万事如意！ 摘要：在求一些函数的定积分时，由于原函数十分复杂难以求出或用初等函数表达，导致积分很难精确求出，只能设法求其近似值，因此能够直接借助牛顿-莱布尼兹公式计算定积分的情形是不多的。数值积分就是解决此类问题的一种行之有效的方法。积分的数值计算是数值分析的一个重要分支；因此，探讨近似计算的数值积分方法是有着明显的实际意义的。本文从数值积分问题的产生出发，详细介绍了一些数值积分的重要方法。 本文较详细地介绍了牛顿-科特斯求积公式，以及为了提高积分计算精度的高精度数值积分公式，即龙贝格求积公式和高斯-勒让德求积公式。除了研究这些数值积分算法的理论外，本文还将这些数值积分算法在计算机上通过MATLAB软件编程实现，并通过实例用各种求积公式进行运算，分析比较了各种求积公式的计算误差。 【关键词】数值积分牛顿-科特斯求积公式高精度求积公式MATLAB软件

前言 对于定积分，在求某函数的定积分时，在一定条件下，虽然有牛顿-莱布里茨公式可以计算定积分的值，但在很多情况下的原函数不易求出或非常复杂。被积函数的原函数很难用初等函数表达出来，例如等；有的函数的原函数存在，但其表达式太复杂，计算量太大，有的甚至无法有解析表达式。因此能够借助牛顿-莱布尼兹公式计算定积分的情形是不多的。另外，许多实际问题中的被积函数往往是列表函数或其他形式的非连续函数，对这类函数的定积分，也不能用不定积分方法求解，只能设法求其近似值。因此，探讨近似计算的数值积分方法是有明显的实际意义的，即有必要研究定积分的数值计算方法，以解决定积分的近似计算。而数值积分就是解决此类问题的一种有效的方法，它的特点是利用被积函数在一些节点上的信息求出定积分的近似值。 微积分的发明是人类科学史上一项伟大的成就，在科学技术中，积分是经常遇到的一个重要计算环节。数值积分是数学上重要的课题之一，是数值分析中重要的内容之一，也是应用数学研究的重点。随着计算机的出现，近几十年来，对于数值积分问题的研究已经成为一个很活跃的研究领域。现在，数值积分在计算

用MATLAB算多元函数积分

用MATLAB 计算多元函数的积分 三重积分的计算最终是化成累次积分来完成的，因此只要能正确的得出各累次积分的积分限，便可在MA TLAB 中通过多次使用int 命令来求得计算结果。但三重积分的积分域Ω是一个三维空间区域，当其形状较复杂时，要确定各累次积分的积分限会遇到一定困难，此时，可以借助MATLAB 的三维绘图命令，先在屏幕上绘出Ω的三维立体图，然后执行命令 rotate3d on ↙ 便可拖动鼠标使Ω的图形在屏幕上作任意的三维旋转，并且可用下述命令将Ω的图形向三个坐标平面进行投影： view(0,0),向XOZ 平面投影； view(90,0),向YOZ 平面投影； view(0,90),向XOY 平面投影. 综合运用上述方法，一般应能正确得出各累次积分的积分限。 例11.6.1计算zdv Ω ???，其中Ω是由圆锥曲面222z x y =+与平面z=1围成的闭区域 解 首先用MA TLAB 来绘制Ω的三维图形，画圆锥曲面的命令可以是： syms x y z ↙ z=sqrt(x^2+y^2); ↙ ezsurf(z,[-1.5,1.5]) ↙ 画第二个曲面之前，为保持先画的图形不会被清除，需要执行命令 hold on ↙ 然后用下述命令就可以将平面z=1与圆锥面的图形画在一个图形窗口内： [x1,y1]=meshgrid(-1.5:1/4:1.5); ↙ z1=ones(size(x1)); ↙ surf(x1,y1,z1) ↙ 于是得到Ω的三维图形如图:

由该图很容易将原三重积分化成累次积分： 111zdv dy -Ω=???? 于是可用下述命令求解此三重积分： clear all ↙ syms x y z ↙ f=z; ↙ f1=int(f,z.,sqrt(x^2+ y^2),1); ↙ f2=int(f1,x,-sqrt(1- y^2), sqrt(1- y^2)); ↙ int(f2,y,-1,1) ↙ ans= 1/4*pi 计算结果为4 π 对于第一类曲线积分和第一类曲面积分，其计算都归结为求解特定形式的定积分和二重积分，因此可完全类似的使用int 命令进行计算，并可用diff 命令求解中间所需的各偏导数。 例11.6.2用MATLAB 求解教材例11.3.1 解 求解过程如下 syms a b t ↙ x=a*cos(t); ↙ y=a*sin(t); ↙ z=b*t; ↙ f=x^2 +y^2+z^2; ↙ xt=diff(x,t); ↙ yt=diff(y,t); ↙ zt=diff(z,t); ↙ int(f*sqrt(xt^2 +yt^2+zt^2),t,0,2*pi) ↙ ans= 2/3*( a^2 +b^2)^1/2*a^2*pi+8/3*( a^2 +b^2)^1/2*b^2*pi^3 对此结果可用factor 命令进行合并化简： factor （ans ） ans= 2/3*( a^2 +b^2)^1/2*pi*(3* a^2 +4*b^2*pi^2) 例11.6.3用MATLAB 求解教材例11.4.1 解 求解过程如下 syms x y z1 z2↙ f= x^2 +y^2; ↙ z1=sqrt(x^2 +y^2); ↙ z2=1; ↙ z1x=diff(z1,x); ↙ z1y=diff(z1,y); ↙ z2x=diff(z2,x); ↙ z2y=diff(z2,y); ↙

数值积分的算法比较及其MATLAB实现

编号： 审定成绩： 重庆邮电大学 毕业设计（论文） 设计（论文）题目：数值积分算法与MATLAB实现 学院名称：数理学院 学生姓名： 专业：数学与应用数学 班级： 学号： 指导教师： 答辩组负责人： 填表时间：年月 重庆邮电大学教务处制

摘要 在求一些函数的定积分时，由于原函数十分复杂难以求出或用初等函数表达，导致积分很难精确求出，只能设法求其近似值，因此能够直接借助牛顿-莱布尼兹公式计算定积分的情形是不多的。数值积分就是解决此类问题的一种行之有效的方法。积分的数值计算是数值分析的一个重要分支；因此，探讨近似计算的数值积分方法是有着明显的实际意义的。本文从数值积分问题的产生出发，详细介绍了一些数值积分的重要方法。 本文较详细地介绍了牛顿-科特斯求积公式，以及为了提高积分计算精度的高精度数值积分公式，即龙贝格求积公式和高斯-勒让德求积公式。除了研究这些数值积分算法的理论外，本文还将这些数值积分算法在计算机上通过MATLAB软件编程实现，并通过实例用各种求积公式进行运算，分析比较了各种求积公式的计算误差。 【关键词】数值积分牛顿-科特斯求积公式高精度求积公式MATLAB软件

ABSTRACT When the solution of the definite integral of some function values，because the original function is very complex and difficult to find the elementary function expression, the integral is difficult to accurately calculate, only managed to find the approximate value, and the case is small that allows to direct interface with the Newton - Leibniz formula to calculate the definite integral. Numerical integration is an effective method to solve such problems. The numerical integration is an important branch of numerical analysis; therefore, exploring the approximate calculation of the numerical integration method has obvious practical significance. This article departure from the numerical integration problem, described in detail some important numerical integration methods. This paper has introduced detail the Newton - Coates quadrature formula, and in order to improve the calculation accuracy of numerical integration formulas, More precise formulas have Romberg quadrature formulas and the Gauss - Legendre quadrature formula. In addition to the study of these numerical integration algorithm theory, the article also involve what these numerical integration algorithm be programmed by matlab software on the computer, and an example is calculated with a variety of quadrature formulas, finally analysis and comparison to various quadrature formulas calculation error. 【Key words】Numerical integration Newton-Cotes quadrature formula High-precision quadrature formula Matlab software

利用Matlab实现Romberg数值积分算法----系统建模与仿真结课作业

利用Matlab 实现Romberg 数值积分算法 一、内容摘要 针对于某些多项式积分，利用Newton —Leibniz 积分公式求解时有困难，可以采用数值积分的方法，求解指定精度的近似解，本文利用Matlab 中的.m 文件编写了复化梯形公式与Romberg 的数值积分算法的程序，求解多项式的数值积分，比较两者的收敛速度。 二、数值积分公式 1.复化梯形公式求解数值积分的基础是将区间一等分时的Newton —Cotes 求积公式: I =(x)[f(a)f(b)]2 b a b a f dx -≈ +? 其几何意义是，利用区间端点的函数值、与端点构成的梯形面积来近似(x)f 在区间[a,b]上的积分值，截断误差为： 3" (b a)()12 f η-- (a,b)η∈ 具有一次的代数精度，很明显，这样的近似求解精度很难满足计算的要求，因而，可以采用将积分区间不停地对分，当区间足够小的时候，利用梯形公式求解每一个小区间的积分近似值，然后将所有的区间加起来，作为被求函数的积分，可以根据计算精度的要求，划分对分的区间个数，得到复化梯形公式： I =1 1 (b a)(b a) (x)dx [f(a)f(b)2(a )]2n b a k k f f n n -=--≈+++∑? 其截断误差为：

2" (b a)h ()12 R f η--= (a,b)η∈ 2.Romberg 数值积分算法 使用复化的梯形公式计算的数值积分，其收敛速度比减慢，为此，采用Romberg 数值积分。其思想主要是，根据I 的近似值2n T 加上I 与2n T 的近似误差，作为新的I 的近视，反复迭代，求出满足计算精度的近似解。 用2n T 近似I 所产生的误差可用下式进行估算： 12221 ()3 n n n I T T T -?=-=- 新的I 的近似值： 122 n n j T T -=?+ j =（0 1 2 ….） Romberg 数值积分算法计算顺序 i=0 (1) 002T i=1 (2) 102T (3) 012T i=2 (4) 202T (5) 112T (6) 022T i=3 (7) 302T (8) 212T (9) 122T (10) 032T i=4 (11) 402T (12) 312T (13) 222T (14) 132T … … … … 其中，第一列是二阶收敛的，第二列是四阶收敛的，第三列是六阶收敛的，第四列是八阶收敛的，即Romberg 序列。

[整理]Matlab积分.

一.数值积分的实现方法 1．变步长辛普生法 基于变步长辛普生法，MA TLAB给出了quad函数来求定积分。该函数的调用格式为：[I,n]=quad('fname',a,b,tol,trace) 其中fname是被积函数名。a和b分别是定积分的下限和上限。tol用来控制积分精度，缺省时取tol=0.001。trace控制是否展现积分过程，若取非0则展现积分过程，取0则不展现，缺省时取trace=0。返回参数I即定积分值，n为被积函数的调用次数。 例8-1 求定积分。 (1) 建立被积函数文件fesin.m。 function f=fesin(x) f=exp(-0.5*x).*sin(x+pi/6); (2) 调用数值积分函数quad求定积分。 [S,n]=quad('fesin',0,3*pi) S = 0.9008 n = 77 2．牛顿－柯特斯法 基于牛顿－柯特斯法，MA TLAB给出了quad8函数来求定积分。该函数的调用格式为：[I,n]=quad8('fname',a,b,tol,trace) 其中参数的含义和quad函数相似，只是tol的缺省值取10-6。?该函数可以更精确地求出定积分的值，且一般情况下函数调用的步数明显小于quad函数，从而保证能以更高的效率求出所需的定积分值。 (1) 被积函数文件fx.m。 function f=fx(x) f=x.*sin(x)./(1+cos(x).*cos(x)); (2) 调用函数quad8求定积分。 I=quad8('fx',0,pi) I = 2.4674 分别用quad函数和quad8函数求定积分的近似值，并在相同的积分精度下，比较函数的调用次数。 调用函数quad求定积分： format long; fx=inline('exp(-x)'); [I,n]=quad(fx,1,2.5,1e-10) I = 0.28579444254766 n = 65 调用函数quad8求定积分： format long; fx=inline('exp(-x)'); [I,n]=quad8(fx,1,2.5,1e-10) I = 0.28579444254754 n = 33

matlab求定积分之实例说明

一、符号积分 符号积分由函数int来实现。该函数的一般调用格式为： int(s)：没有指定积分变量和积分阶数时，系统按findsym函数指示的默认变量对被积函数或符号表达式s求不定积分； int(s,v)：以v为自变量，对被积函数或符号表达式s求不定积分； int(s,v,a,b)：求定积分运算。a,b分别表示定积分的下限和上限。该函数求被积函数在区间[a,b]上的定积分。a和b可以是两个具体的数，也可以是一个符号表达式，还可以是无穷(inf)。当函数f关于变量x在闭区间[a,b]上可积时，函数返回一个定积分结果。当a,b中有一个是inf时，函数返回一个广义积分。当a,b中有一个符号表达式时，函数返回一个符号函数。 例： 求函数x^2+y^2+z^2的三重积分。内积分上下限都是函数，对z积分下限是sqrt(x*y)，积分上限是x^2*y；对y积分下限是sqrt(x)，积分上限是x^2；对x的积分下限1，上限是2，求解如下： >>syms x y z %定义符号变量 >>F2=int(int(int(x^2+y^2+z^2,z,sqrt(x*y),x^2*y),y,sqrt(x),x^2),x,1,2) %注意定积分的书写格式 F2 = 1610027357/6563700-6072064/348075*2^(1/2)+14912/4641*2^(1/4)+64/225*2 ^(3/4) %给出有理数解 >>VF2=vpa(F2) %给出默认精度的数值解 VF2 = 224.92153573331143159790710032805 二、数值积分 1.数值积分基本原理 求解定积分的数值方法多种多样，如简单的梯形法、辛普生(Simpson)?法、牛顿－柯特斯(Newton-Cotes)法等都是经常采用的方法。它们的基本思想都是将整个积分区间[a,b]分成n个子区间[xi,xi+1]，i=1,2,…,n，其中x1=a，xn+1=b。这样求定积分问题就分解为求和问题。 2.数值积分的实现方法 基于变步长辛普生法，MATLAB给出了quad函数来求定积分。该函数的调用格式为： [I,n]=quad('fname',a,b,tol,trace) 基于变步长、牛顿－柯特斯(Newton-Cotes)法，MATLAB给出了quadl函数来求定积分。该函数的调用格式为： [I,n]=quadl('fname',a,b,tol,trace) 其中fname是被积函数名。a和b分别是定积分的下限和上限。tol用来控制积分精度，缺省时取tol=0.001。trace控制是否展现积分过程，若取非0则展现积分过程，取0则不展现，缺省时取trace=0。返回参数I即定积分值，n为被积函数的调用次数。 例： 求函数'exp(-x*x)的定积分，积分下限为0，积分上限为1。 >>fun=inline('exp(-x.*x)','x'); %用内联函数定义被积函数fname

详解Matlab求积分的各种方法

详解Matlab求积分的各种方法 一、符号积分由函数int来实现。 该函数的一般调用格式为： int(s)： 没有指定积分变量和积分阶数时，系统按findsym函数指示的默认变量对被积函数或符号表达式s求不定积分；int(s,v)： 以v为自变量，对被积函数或符号表达式s求不定积分；int(s,v,a,b)： 求定积分运算。 a,b分别表示定积分的下限和上限。 该函数求被积函数在区间[a,b]上的定积分。 a和b可以是两个具体的数，也可以是一个符号表达式，还可以是无穷(inf)。 当函数f关于变量x在闭区间[a,b]上可积时，函数返回一个定积分结果。 当a,b中有一个是inf时，函数返回一个广义积分。 当a,b中有一个符号表达式时，函数返回一个符号函数。 例： 求函数x^2+y^2+z^2的三重积分。 内积分上下限都是函数，对z积分下限是sqrt(x*y)，积分上限是x^2*y；对y积分下限是sqrt(x)，积分上限是x^2；对x的积分下限1，上限是2，求解如下： >>syms x y z %定义符号变 量>>F2=int(int(int(x^2+y^2+z^2,z,sqrt(x*y),x^2*y),y,sqrt(x),x^2),x,1,2) %注意定积分的书写格式F2 =57/-

/348075*2^(1/2)+14912/4641*2^(1/4)+64/225*2^(3/4) %给出有理数 解>>VF2=vpa(F2) %给出默认精度的数值解VF2 = 224.9 232805二、数值积分 1.数值积分基本原理求解定积分的数值方法多种多样，如简单的梯形法、辛普生(Simpson)?法、牛顿－柯特斯(Newton-Cotes)法等都是经常采用的方法。 它们的基本思想都是将整个积分区间[a,b]分成n个子区间[xi,xi+1]， i=1,2,…,n，其中x1=a，xn+1=b。 这样求定积分问题就分解为求和问题。 2.数值积分的实现方法基于变步长辛普生法，MATLAB给出了quad函数来求定积分。 该函数的调用格式为： [I,n]=quad('fname',a,b,tol,trace)基于变步长、牛顿－柯特斯(Newton-Cotes)法，MATLAB给出了quadl函数来求定积分。 该函数的调用格式为： [I,n]=quadl('fname',a,b,tol,trace)其中fname是被积函数名。 a和b分别是定积分的下限和上限。 tol用来控制积分精度，缺省时取tol= 0.0 01。 trace控制是否展现积分过程，若取非0则展现积分过程，取0则不展现，缺省时取trace= 0。

MATLAB计算积分

函数的积分和椭圆的周长 1．正弦函数的积分 [问题]求正弦函数从0到π的积分 y = sin x 当x = 0时，积分为0，画出积分的函数曲线。 [数学模型] 定积分的结果为 ππ00 sin d cos 2S x x x ==-=? 不定积分的结果为 sin d cos I x x x C ==-+? 其中C 是积分常量，由初始条件决定。当x = 0时，积分为I = 0，必有C = 1。结果为 I = -cos x + 1 [算法]根据积分的基本概念，将积分区域分为多份，用矩形法求曲线下的近似面积表示积分的近似值 1()n i i S f x x ==?∑ 矩形法的函数是sum(f)。 用梯形法求曲线下的近似面积表示积分的近似值 1 101[()()]2 n i i i S f x f x x -+==+?∑ 梯形法的函数是trapz(f)。 用数值积分的函数是quad 和quadl ，常用使用格式是 S = quad(f,a,b) 其中，f 表示被积函数，a 表示积分的下限，b 表示积分的下限。 用符号的函数是int ，常用使用格式是 S = int(f,a,b) [程序]zqy4_1.m 如下。 %正弦函数的积分 clear %清除变量 x=linspace(0,pi); %自变量向量 dx=x(2); %间隔 y=sin(x); %被积函数 s1=sum(y)*dx %矩形法积分 s2=trapz(y)*dx %梯形法积分 f=inline('sin(x)'); %被积的内线函数 s3=quad(f,0,pi) %数值定积分

s4=int('sin(x)',0,pi) %符号积分 sc1=cumsum(y)*dx; %矩形法累积积分(精度稍差) sc2=cumtrapz(y)*dx; %梯形法累积积分 figure %创建图形窗口 plot(x,-cos(x)+1,x,sc1,'.',x,sc2,'o') %画解析式和矩阵法以及梯形法积分曲线 s=int('sin(x)') %符号积分 sc3=subs(s,'x',x); %替换数值求符号积分的值 C=-sc3(1) %求积分常数 hold on %保持图像 plot(x,sc3+C,'c*') %画符号法积分曲线 grid on %加网格 fs=16; %字体大小 xlabel('\itx','FontSize',fs) %横坐标 ylabel('\intsin\itx\rmd\itx','FontSize',fs)%纵坐标 title('正弦函数的积分','FontSize',fs) %标题 legend('解析解','矩形法','梯形法','符号法')%图例 zqy4.1图 zqy4.2图 2.三角函数和指数的积分 [问题]求如下函数的积分 y = e ax sin bx 其中a = 0.5，b = 2。积分下限为0。画出积分的函数曲线。 [数学模型] 设 11e sin d sin de {e sin e cos d }ax ax ax ax I bx x bx bx b bx x a a == =-??? 11{e sin cos de }{e sin [e cos e sin d ]}ax ax ax ax ax b b bx bx bx bx b bx x a a a a =-=-+?? 因此不定积分为 221e (sin cos )ax I a bx b bx C a b =-++ 当x = 0时，I 应该为零，所以 22b C a b =+

详解Matlab求积分的各种方法

详解Matlab求积分地各种方法 一、符号积分 符号积分由函数int来实现.该函数地一般调用格式为： int(s)：没有指定积分变量和积分阶数时，系统按findsym函数指示地默认变量对被积函数或符号表达式s求不定积分； int(s,v)：以v为自变量，对被积函数或符号表达式s求不定积分； int(s,v,a,b)：求定积分运算.a,b分别表示定积分地下限和上限.该函数求被积函数在区间[a,b]上地定积分.a和b可以是两个具体地数，也可以是一个符号表达式，还可以是无穷(inf).当函数f关于变量x在闭区间[a,b]上可积时，函数返回一个定积分结果.当a,b中有一个是inf时，函数返回一个广义积分.当a,b中有一个符号表达式时，函数返回一个符号函数. 例： 式 F2 = 1. 2. 其中 取tol=0.001.trace控制是否展现积分过程，若取非0则展现积分过程，取0则不展现，缺省时取trace=0.返回参数I即定积分值，n为被积函数地调用次数. 例： 求函数'exp(-x*x)地定积分，积分下限为0，积分上限为1. >>fun=inline('exp(-x.*x)','x'); %用内联函数定义被积函数fname >>Isim=quad(fun,0,1) %辛普生法 Isim = 0.746824180726425 IL=quadl(fun,0,1) %牛顿－柯特斯法

IL = 0.746824133988447 三、梯形法求向量积分 trapz(x,y)—梯形法沿列方向求函数Y关于自变量X地积分(向量形式，数值方法). >>d=0.001; >>x=0:d:1; >>S=d*trapz(exp(-x.^2)) S= 0.7468 或： S =