第四章 数值积分与数值微分 第四章 数值积分与数值微分 4.1 Newton-Cotes公式§ 4.2 复合求积法§ 4.3 Romberg算法§ 4.4* Gauss求积法§ 4.5 数值微分§ 本章要点 公式 :近似值的几个基本求积计算定积分 从而导出代替被积函数本章将用插值多项式 ∫ba dxxf xfxP )( ),()( (1) 等距节点下的 :Newton-Cotes公式和 Romberg公式 (2) 数值微分公式 本章应用题 : 为了计算瑞士国土的面积 ,首先对地图作了如下测量 :以西 向东方向为 x轴 ,由南向北方向为 y轴 ,选择方便的原点 ,并将 从最西边界到最东边界在 x轴上的区间适当地划分为若干 段 ,在每个分点的 y方向测出南边界点和北边界点的 y坐标 , 数据如表 (单位 mm): x 7.0 10.5 13.0 17.5 34.0 40.5 44.5 48.0 56.0 y1 44 45 47 50 50 38 30 30 34 y2 44 59 70 72 93 100 110 110 110 x 61.0 68.5 76.5 80.5 91.0 96.0 101.0 104.0 106.5 y1 36 34 41 45 46 43 37 33 28 y2 117 118 116 118 118 121 124 121 121 x 111.5 118.0 123.5 136.5 142.0 146.0 150.0 157.0 158.0 y1 32 65 55 54 52 50 66 66 68 y2 121 122 116 83 81 82 86 85 68 0 20 40 60 80 100 120 140 160 20 40 60 80 100 120 140 瑞士地图的外形如图 (比例尺 18mm:40km) 试由测量数据计算瑞士国土的近似面积 ,并与其精确值 41288平方公里比较 4.1 Newton-Cotes公式§ ∫= ba dxxffI )()(对于积分 公式有则由的原函数如果知道 LeibnizNewtonxFxf ?),()( ∫ba dxxf )( )()()( aFbFxF ba ?== 但是在工程技术和科学研究中 ,常会见到以下现象 : 的一些数值只给出了的解析式根本不存在 )(,)()1( xfxf 不是初等函数如求不出来的原函数 )(,)()()2( xFxFxf 求原函数较困难的表达式结构复杂 ,)()3( xf 以上这些现象 ,Newton-Leibniz很难 发挥作用 只能建立积分的近似计算方法 这类方法很多 ,但为方便起见 ,最常用的一种方法是利用 插值多项式来构造数值求积公式 ,具体步骤如下 : 上取一组节点在积分区间 ],[ ba bxxxa n ≤<<<≤ L10 次插值多项式的作 nxf )( ∑ = = n k kkn xlxfxL 0 )()()( 为插值基函数),,1,0)(( nkxlk L= 不同的 插值方法 有不同的 基函数 有的近似作为被积函数用 ,)()( xfxLn ∫ba dxxf )( ∫≈ ba n dxxL )( ∫ ∑ = = b a n k kk dxxlxf 0 )()( ∑ ∫ = = n k b a kk dxxlxf 0 )()( 则,若计 ∫= b a kk dxxlA )( ∫= ba dxxffI )()( ∑ = ≈ n k kk xfA 0 )( 这就是数值求积公式 称为求积系数其中 kA 为了使一个求积公式能对更多的积分具有较好的实际计 算意义 ,就要求它对尽可能多的被积函数都准确地成立 )( fIn= 因此定义代数精度的概念 : 定义 1. 若求积公式 ∫= ba dxxffI )()( )()( 0 fIxfA n n k kk =≈ ∑ = 即都准确成立次的代数多项式对任意次数不超过 ,))(( mixPm i ≤ 即只要立次多项式却不能准确成但对 ,1+m ∫ba i dxxP )( ∑ = = n k kik xPA 0 )( mi ,,1,0 L= ∫ +ba m dxx 1 ∑ = +≠ n k m kk xA 0 1 则称该求积公式具有 m次的代数精度 代数精度也称 代数精确度 例 1. 试确定下面积分公式中的参数使其代数精确度尽量高 . )()]()0([)]()0([2)()( 12 0 fIhffahhffhdxxffI h =′?′++≈= ∫ ∫= h dxxI 0 0解 : 2 2 1 hI = ]20[2 2 3 1 hah hI ?+= 0)( xxf =对于 hI = 1h= ∫= h dxxI 0 11)( xxf =对于 2 2h = ∫= h dxxI 0 2 2)( xxf =对于 3 3h = 3)221( ha?= 1II =令 12 1=a ]30[2 22 4 1 hah hI ?+=∫= h dxxI 0 3 3)( xxf =对于 4 4h = 4 4h = ]40[2 32 5 1 hah hI ?+=∫= h dxxI 0 4 4)( xxf =对于 5 5h = 6 5h = 3,2,1,0)()( 1 == jxIxI jj )()( 414 xIxI ≠ 因此 所以该积分公式具有 3次代数精确度 一 、 Newton-Cotes数值求积公式 Newton-Cotes公式是指等距节点下使用 Lagrange插值 多项式建立的数值求积公式 ],[)( baCxf ∈设函数 为插值多项式及余项分别的 Lagrangexf )( 等份分割为将积分区间 nba ],[ nkkhaxk ,,1,0, L=+= 为步长其中 n abh ?= 各节点为 ∑ = = n k kkn xlxfxL 0 )()()( )()!1( )()( 1 )1( xnfxR n n n + + += w x ],[ ba∈x ∏ = + ?= n i in xxx 0 1 )()(w其中 ∏ ≠ ≤≤ ? ?= kj nj jk j k xx xxxl 0 )( 而 )()()( xRxLxf nn += 因此对于定积分 ∫= b a dxxffI )()( ∫ += ba nn dxxRxL )]()([有 ∫= ba dxxffI )()( ∫ ∑ = = b a n k kk dxxlxf 0 )()( ∫+ b a n dxxR )( ∑ = = n k kk xfA 0 )( ∫+ b a n dxxR )( 令 ∑ = = n k kkn xfAfI 0 )()( ∫= ba nn dxxRIR )()( ∫= ba dxxffI )()( )()()( nn IRfIfI +=即有 ∫= ba kk dxxlA )(其中 dxxx xxba kj nj jk j∫ ∏ ≠ ≤≤ ? ?= 0 n阶 Newton-Cotes求积公式 Newton-Cotes公式的余项 (误差 ) )()( fIfI n≈ ∫= ba kk dxxlA )( dxxx xxba kj nj jk j∫ ∏ ≠ ≤≤ ? ?= 0 :的计算kA 注意是等距节点 thax +=假设 ],[ bax ∈由 ],0[ nt ∈可知 kA dxxx xxb a kj nj jk j∫ ∏ ≠ ≤≤ ? ?= 0 dthhjk hjtn kj nj ?? ?? ? ? ? ?? ? ? ? ? ?= ∫ ∏ ≠ ≤≤0 0 )( )( dtjtknkh n kj nj kn ∫ ∏ ≠ ≤≤ ? ??? ??= 0 0 )()!(! )1( dtjtknknab n kj nj kn ∫ ∏ ≠ ≤≤ ? ???? ???= 0 0 )()!(! )1()( )()(? n kk CabA ??= ∑ = = n k kkn xfAfI 0 )()( ∑ = ?= n k k n k xfCab 0 )( )()( 所以 Newton-Cotes公式化为 系数称为 CotesC nk )( 思考 使用 n次 Lagrange插值多项式的 Newton-Cotes 公式至少具有 n次代数精度 ,并且 n为偶数时至 少具有 n+1次代数精度 ,试以 n=1,2,4为例说明 该结果 二 、 低阶 Newton-Cotes公式及其余项 在 Newton-Cotes公式中 ,n=1,2,4时的公式是最常用也 最重要三个公式 ,称为低阶公式 1.梯形 (trapezia)公式及其余项 abhbxaxn ?==== ,,,1 10则取 dtt∫ ??= 1 0 )1()1(0CCotes系数为 21= dtt∫= 1 0 )1( 1C 2 1= 求积公式为 )(1 fI ∑ = ?= 1 0 )1( )()( k kk xfCab )]()([2 10 xfxfab +?= )]()([2 bfafab +?=)(1 fI即 上式称为 梯形求积公式 ,也称 两点公式 , 记为 -0.5 0 0.5 1 1.5 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 )]()([2 )( bfafab +?=)(1 fIT = 梯形公式的余项为 )()( 1IRTR = ∫= b a dxxR )(1 dxbxaxfTR b a∫ ??′′= ))((2 )()( x dxbxaxf b a∫ ??′′= ))((2 )(h ],[ ba∈h 第二积分 中值定理 6 )( 2 )( 3abf ?′′?= h )()!1( )()( 1 )1( xnfxR n n n + + += w x )(12 )( 3 hfab ′′??= 2 3 12 )(|)(| MabTR ?≤ |)(|max ],[2 xfM bax ′′= ∈ 梯形 (trapezia)公式具有 1次代数精度 故 2.Simpson公式及其余项 2,,2,,2 210 abhbxabxaxn ?==+=== 则取 Cotes系数为 dtttC ∫ ??= 2 0 )2( 0 )2)(1(4 1 6 1= dtttC ∫ ??= 2 0 )2( 1 )2(2 1 6 4= dtttC ∫ ?= 2 0 )2( 2 )1(4 1 6 1= 求积公式为 )(2 fI ∑ = ?= 2 0 )2( )()( k kk xfCab )](61)(64)(61)[( 210 xfxfxfab ++?= )]()2(4)([6 bfbafafab +++?= )(2 fI -0.5 0 0.5 1 1.5 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 上式称为 Simpson求积公式 , 也称 三点公式或抛物线公式 记为 )(2 fIS = Simpson公式的余项为 )()( 2IRSR = ∫= ba dxxR )(2 )()2(180 )4(4 hfabab ???= Simpson公式具有 3次代数精度 3.Cotes公式及其余项 4,4,,1,0,,4 abhkkhaxn k ?==+== L则取 Cotes系数为 dtttttC )4)(3()2)(1(!44 1 4 0 )4( 0 ?????= ∫ 90 7= dtttttC )4)(3()2(!34 1 4 0 )4( 1 ?????= ∫ 90 32= dtttttC )4)(3()1(!2!24 1 4 0 )4( 2 ?????= ∫ 90 12= dtttttC )4)(2()1(!34 1 4 0 )4( 3 ?????= ∫ 90 32= dtttttC )3)(2()1(!44 1 4 0 )4( 4 ?????= ∫ 90 7= 求积公式为 )(4 fI ∑ = ?= 4 0 )4( )()( k kk xfCab )](907)(9032)(9012)(9032)(907)[( 43210 xfxfxfxfxfab ++++?= )](7)(32)(12)(32)(7[90 43210 xfxfxfxfxfab ++++?= 上式称为 Cotes求积公式 , 也称 五点公式 记为 )(4 fIC = Cotes公式的余项为 )()( 4IRCR = ∫= ba dxxR )(4 )()4(945 )(2 )6(6 hfabab ???= Cotes公式具有 5次代数精度 三 、 Newton-Cotes公式的稳定性 (舍入误差 ) dtjtknknC n kj nj kn n k ∫ ∏ ≠ ≤≤ ? ???? ?= 0 0 )( )( )!(! )1( 考察 Cotes系数 无关与函数的划分有关的节点只与积分区间 )(,],[ xfxba j 因此用 Newton-Cotes公式计算积分的舍入误差主要由 的计算引起函数值 )( kxf 其值可以精确给定 响的舍入误差对公式的影只需讨论 )( kxf )()()(,)( 计算值的近似值作为而以为精确值假设 kkk xfxfxf 为误差)()( kkk xfxf ?=e )( fIn ∑ = ?= n k k n k xfCab 0 )( )()(记 )(计算值的近似值为 nI 而理论值为 )( fIn ∑ = ?= n k k n k xfCab 0 )( )()( 的误差为与 nn II )()( fIfI nn ? ∑ = ??= n k kk n k xfxfCab 0 )( )]()([)( nn II ? ∑ = ?= n k k n kCab 0 )()( e ∑ = ?≤ n k k n kCab 0 )()( e nn II ? ∑ = ?≤ n k n kCab 0 )()( e |}max{| kee = 有若 ,0, )( >≤? nkCnk nn II ? ∑ = ?≤ n k n kCab 0 )()( e ∑ = ??= n k n kCab 0 )( 1)( e ∑ = ??= n k k n k xgCab 0 )( )()( e )1)(( ≡xg ∫= ba dxxg )(e ∫= badxe e)( ab ?= 1 0 )( =∑ = n k n kC性质 : 即 nn II ? e)( ab ?≤ Newton-Cotes公式的舍入误差只是函数值误差的 倍)( ab ? 时 , 公式都是稳定的当事实上 8, <n 公式是稳定的时即 CotesNewtonCnk nk ?>≤? ,0, )( ∑ = ? n k n kCab 0 )()( e 有有正有负若 ,)(nkC ∑ = ?≥ n k n kCab 0 )()( e e)( ab ?= 此时 ,公式的稳定性将无法保证 因此 ,在实际应用中一般不使用高阶 Newton-Cotes公式 而是采用低阶复合求积法 (下节 ) 思考 1. n=0时的 Newton-Cotes公式称为 矩形公式 , 试求出该公式 2. 试编写 trapezia公式 、 Simpson公式 、 Cotes公式的模块程序