用MATLAB求数值积分的方法

形公式 ①, ②计算积分的近似值.
3. 2 梯形公式命令 t rapz
命令格式 :
t rapz ( x ,y) % 输入向量 x = [ x0 ,x1 , …,xn ] ,
输出同维数的向量 y = [f0 ,f1 , …,f n ] ,按梯形公式计
算积分近似值 ;
t rapz (y) % 按梯形公式计算积分 ,但取步长
值 ,非常接近精确值 1. (3) 用辛普森公式命令 输入命令 : > > I4 = quad (′sin ( x ) ′,0 , pi/ 2) 结果为 I4 = 1. 000 0. 该值就是精确值 1. 但是 quad 和 quad8 需给出
函数表达式 ,不能像 sum 和 trapz 那样 ,能对用数值 给出的 x , y 进行积分. 4 数值积分法的应用实例
y1 = [ 44 45 47 50 50 38 30 30 34 36 34 41 45 46 43 37 33 28 32 65 55 54 52 50 66 66 68 ] 3 (40/ 18) ;
y2 = [ 44 59 70 72 93 100 110 110 110 117 118 116 118 118 121 124 121 121 121 122 116 83 81 82 86 85 68 ] 3 (40/ 18) ;
© 1994-2011 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net
© 1994-2011 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net
第6期
陈佩宁等 :用 MA TL A B 求数值积分的方法
59
(2) 梯形公式 如果将二者求平均值 ,则每个小区间上的小矩 形变为小梯形 ,整个区间上的值变为 :
积分 Sn ,即
m- 1
∑ Sn
=
h 3
[f
( x0)
+ f (x2m)
+ 4 f ( x2k+1)
k=0
+
m- 1
∑ 2 f (x2k) ] .
④
k=1
④式称为辛普森公式.
3 求数值积分的 MATLAB 积分命令
3. 1 矩形公式命令 su m
命令格式 :
su m (y) % 输出一个向量 y 的分量的和 ,按矩
取左端点和右端点的函数值为小矩形的高 ,则分别
得到两个曲边梯形面积的近似计算公式 :
n- 1
∑ L n
=
h f (xk)
k=0
,h
=
b
n
a
①
n
∑ Rn
=
h f (xk)
k=1
,h
=
bn
a
②
称公式 ①, ②分别为左 、右矩形公式 ,两个矩形
面积分别小于和大于所求曲边梯形的面积.
Ξ 收稿日期 :2007212224 作者简介 :陈佩宁 (19712 ) ,女 ,河北望都人 ,石家庄职业技术学院讲师.
∫ ∫ 2
函数表示 ,如 ex dx ,
si n x
x
d
x
等
.
(2) 函数 f (x) 结构复杂 ,求其原函数非常困难.
(3) 函数 f (x) 的结构虽然简单且其原函数存
在 ,但其原函数的结构相对复杂.
(4) 函数 f (x) 没有具体的表达式 ,只有一些由
试验测试数据形成的表格或图形.
而在这些情况下 ,可采用“数值积分”的方法求
出定积分 (近似值) .
2 数值积分法
用数值积分的方法求一个函数在区间[ a ,b ] 上
的定积分 ,可利用定积分的定义来求解 :
∫ ∑ I =
b
f (x)
a
dx
=
n
l
n
im
→∞wenku.baidu.com
k
=
1
f
(ξk)
bn
a,
∑ 设 In
=
n
f (ξk)
k=1
b
n
a
,则
I
=
l
n
im
→∞
In
.
此时称 In 为数值积分. 显然 ,数值积分 In 就是 I
60
石家庄职业技术学院学报
第 20 卷
表 1 瑞士地图测量数据
x
7 10. 5 13 17. 5 34
40. 5 44. 5 48
56
61
68. 5 76. 5 80. 5 91
y1
44 45
47 50
50
38
30
30
34
36
34
41
45
46
y2
44 59
70 72
93 100 110 110 110 117 118 116 118
Abstract :This paper discusses ,wit h illust rations ,t he reasons and met hods of using matlab to solve t he nu2 merical integration.
Key words :MA TLAB ; numerical integration ; rectangular formula ; t rapezoid formula ; Simpson formula
1 引言
在一元微积分学中 , 若已知函数 f ( x ) 在闭区
间 [ a , b ] 上连续且其原函数为 F ( x ) ,求 f ( x ) 在该
区间上的定积分可用牛顿 - 莱布尼兹公式求解 , 即
∫b
f ( x) d x =
F( x)
|
b a
=
F ( b)
-
F( a) . 而 在
a
∫b
MA TLAB 中可以用符号积分命令 int 求 f ( x ) d x ,
n- 1
∑ Tn
=
h f (xk)
k=1
+
h 2
[f
( x0)
+ f (xn) ]
③
其中 ,h
=
b
n
a.
③式称为梯形公式.
(3) 辛普森公式
为了提高计算结果的精度 ,用分段二次插值函
数代替 f (x) . 由于每段要用到相邻两个小区间端点
的三个函数值 ,所以 ,小区间的数目必须是偶数 ,记
n = 2m ,k = 0 ,1 , …,m - 1 ,在第 k 个小区间上用三
的近似值 ,并且当 n 越大 , In 就越接近于精确值 I. 由
于ξk 取值不同 ,数值积分 In 的结果会有所不同. [2 ] 数值积分的计算公式也有多种 :
(1) 矩形公式
将积分区间[ a ,b ]n 等分 ,每个小区间宽度均为
h = (b - a) / n ,h 称为积分步长.
记 a = x0 < x1 < … < xk … < xn = b ,在小区 间上用小矩形面积近似小曲边梯形的面积 ,若分别
118
x
96 101
104 106. 5 111. 5 118 123. 5 136. 5 142 146 150 157 158
y1
43 37
33 28
32
65
55
54
52
50
66
66
68
y2 121 124
121 121 121 122 116
83
81
82
86
85
68
x = [ 7 10. 5 13 17. 5 34 40. 5 44. 5 48 56 61 68. 5 76. 5 80. 5 91 96 101 104 106. 5 111. 5 118 123. 5 136. 5 142 146 150 157 158 ] 3 (40/ 18) ;
(1) 用左 、右矩形公式命令 输入命令 : > > x = linspace (0 , pi/ 2 ,50) ; y = sin ( x ) > > I1 = sum ( y (1 :49) ) 3 ( pi/ 2) / (50 - 1) > > I2 = sum ( y (2 :50) ) 3 ( pi/ 2) / (50 - 1) 结果为 I1 = 0. 983 9 , I2 = 1. 015 9. (2) 用梯形公式命令 输入命令 : > > I3 = t rapz ( x , y) 结果为 I3 = 0. 999 9. 该值就是 I1 和 I2 的平均
参考文献 : [ 1 ] 柯善军. 高等数学与应用实验 [ M ] . 北京 : 北京航空航天大学
出版社 ,2007.
[ 2 ] 萧树铁. 大学数学数学实验 [ M ] . 北京 : 高等教育 出 版 社 , 1997.
责任编辑 :金 欣
Methods to solve numerical integration by using MATLAB
陈佩宁a , 刘 竞b
(石家庄职业技术学院 a. 信息工程系 ; b. 机电工程系 ,河北 石家庄 050081)
摘 要 :介绍了数值积分法的几种计算公式及相应的 MA TLAB 命令 ,并给出了用 MA TLAB 编程求数值积分 的实例.
关键词 :MA TLAB ;数值积分 ;矩形公式 ;梯形公式 ;辛普森公式 中图分类号 :O172 文献标识码 :A
h = 1. 3. 3 辛普森公式命令 quad
命令格式 : quad (′f u n′, a , b) %计算函数 f u n 从 a 到 b 的定积分 ,自动选择步长 ,误差为 10 - 3 ; quad8 (′f un′, a , b) % 用高精度计算 ,效果比 quad 更好. 就例 1 中的定积分 ,分别介绍三种命令的用法 , 并对计算结果进行比较.
根据地图的比例 ,18 mm 相当于 40 km ,由测量 数据根据下列程序计算瑞士国土的近似面积 :
图 1 瑞士地图
© 1994-2011 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net
t rapz ( x , y2) - t rapz ( x , y1) 运行该程序得结果为 42 414 km2 . 5 结束语 使用 MA TLAB 命令求积分不仅可对已知解析 表达式的函数进行积分 ,还可以对一些离散点的数 据进行积分 (近似值) ,通过输入简单的命令就可得 到结果 ,轻松完成手工计算难以处理的含有大量数 据的工作.
CHEN Pei2ninga , L IU J ingb
(a. Depart ment of Information Technology ; b. Depart ment of Electric and Electronic Engineering , Shijiazhuang Vocational Technology Institute ,Shijiazhuang , Hebei 050081 ,China)
2008 年 12 月 第 20 卷第 6 期
石家庄职业技术学院学报 Journal of Shijiazhuang Vocational Technology Institute
文章编号 :100924873 (2008) 0620058203
Dec. 2008 Vol. 20 No. 6
用 MATL AB 求数值积分的方法Ξ
a
该命令格式为 :
int ( f , x , a , b) % 求函数 f 在区间[ a , b ] 上的
定积分.
π
∫ 例 1 求 2 si nx dx. 0
解 输入命令 : > > sy ms x
> > I = i nt ( si n (x) ,x ,0 ,pi/ 2) , 结果显示为 : I = 1. 用牛顿 - 莱布尼兹公式计算定积分的方法在 理论上和解决实际问题中起到了很大的作用 ,但它 并不能解决定积分计算的所有问题. 在工程技术领 域常遇到十分复杂的情况而无法用牛顿 - 莱布尼 兹公式求解. 其可能出现的情况[1 ] 有 : (1) 某些被积函数 f (x) , 其原函数无法用初等
个节点 ( x2k ,f2k) , ( x2k+1 ,f2k+1) , ( x2k+2 ,f 2k+2) 作二次
插值函数 Sk (x) ,然后积分可得 :
∫x2k+1 Sk ( x) d x = x2k
h 3
[f
( x2k)
+ 4f ( x2k+1)
+
f (x2k+2) ] ,求出 m 段的和就得到整个区间上的近似
对于不知道函数具体表达式 , 只知道 x , y 的 一些实测数据的函数 , 可以用数值积分法的梯形公 式求定积分. 以求瑞士国土面积为例说明其应用. 图 1 是瑞士地图 , 为了计算它的国土面积 , 首先对地 图作如下测量 : 以由西向东方向为 x 轴 , 由南向 北方向为 y 轴 , 选择方便的原点 , 并将从最西边 界点到最东边界点在 x 轴上的区间适当地划分为 若干段 , 在每个分点的 y 方向测出南边界点和北 边界点的 y 坐标 y1 和 y2 , 这样就得到表 1 中的数 据 (单位 mm) .