数值分析与计算方法第三章数值积分与数值微分

1、 本章内容 3.1 3.1 数值积分的基本思想与代数精确度数值积分的基本思想与代数精确度 3.2 3.2 牛顿牛顿- -柯茨公式柯茨公式 3.3 3.3 龙贝格算法龙贝格算法 3.4 3.4 高斯公式高斯公式 3.5 3.5 数值微分数值微分dxxfba)(第三章第三章 数值积分与数值微分数值积分与数值微分实际问题中，往往需要实际问题中，往往需要计算定积分计算定积分dxxfba)(Newton-Leibniz公式公式： ，；，其中)()()()()(xfxFaFbFdxxfba困难：困难： 某些函数的原函数不能用初等函数来表示；某些函数的原函数不能用初等函数来表示；2222sin1sinsin

2、exxxex，如：如： 等；等； 简单，但简单，但 过于复杂，不便计算；过于复杂，不便计算； 以数据形式给出。以数据形式给出。 )(xf)(xf)(xF数值方法由由 由由积分中值定理积分中值定理： )()(abfdxxfba定积分的值可能期望用定积分的值可能期望用被积函数的值被积函数的值来来直接决定直接决定，可见，可见dxxfabfba)()(1)(给出给出平均值平均值：的的某种近似算法某种近似算法，便能相应地获得一种数值积分方法。，便能相应地获得一种数值积分方法。 ；则有：，似平均值：取区间中点处的高度近)2()()()2()()2(bafabdxxfbaffba；则有：，近似平均值：取区间

3、端点高度的平均)()(2)()()(21)() 1 (bfafabdxxfbfaffba如：如：梯形公式梯形公式 中矩形公式中矩形公式 ()( )()22( )ababyfyf axayf x中矩形公式和梯形公式的本质就是在区间a,b上分别用直线和代替曲线,然后积分.( )f x更一般地可以用在区间上次数更高的插值多项式来近似,再对插值多项式积分以取得更理想的近似公式.图4-2图4-1010( )( )( ) ()nnnkkkaxxxbLagrangef xL xlx f x给定一组节点，作插值多项式：，则有：)( fIn称称 为计算为计算 的的插值型求积公式插值型求积公式, ,称称Ak为求积

4、系数求积系数, xk为求积节点求积节点. )( fIn)( fI000( )( )( ) ()()( )( )nnbbkkkkaakknbjkkkajkjj kI ff x dxlx f xdxA f xxxlxAlx dxxx，其中：，易得；(1)1( )( )( )( )( )(1)!nnnnfR xf xL xxn插值余项：，(1)110 ( )( )( )( )( )( )( )(1)!( )bbnnnaabnanbnanniiRfI fIff x dxL x dxR x dxfx dxnxxx 积分余项：；其中，定义：定义：若求积公式若求积公式nkkknxfAfI0)()(对于对于次

5、数次数 m均能准确成立均能准确成立，而，而至少对于一个至少对于一个 的多项式的多项式1m立立，则称求积公式具有，则称求积公式具有 次多项式不能准确成次多项式不能准确成m次次代数精度代数精度。插值型求积公式误差为：插值型求积公式误差为： banndxxnffR)()!1()()(1)1(若若 为次数为次数 的多项式，则的多项式，则 ，)(xfn0)(fR，从而，从而 个节个节 1n点的插值型点的插值型求积公式求积公式至少有至少有 次代数精度次代数精度。n中矩形公式、梯形公式中矩形公式、梯形公式有有1 1次次代数精度代数精度。3.1THn插值型求积公式至少有次代数精度；0( )()( )1Lbkk

6、akmf x dxA f xmf xxx一般地，至少有次代数精度该公式对， ，准确成立。L0L220L11012.1kkkkkmmmkkkAbaA xbaA xbam即，的代数精度；求)43(2)21()41(231)(10fffdxxf例：例：10(1)( )11f xdx ，右边；解：解：3次精度次精度140137(5)( )( )5192f xxf x dx，右边；1301(4)( )( )4f xxf x dx，右边；1201(3)( )( )3f xxf x dx，右边；101(2)( )2f xxx dx，右边；( )0( )()()nbnkkakf x dxbaCf x由此构造的

7、插值型求积公式为：，( )()kkkbana bhxakhnf xxf x， ，步长，节点：；设对应于节点的值为等分， 称为称为Newton-CotesNewton-Cotes公式公式 CotesCotes系数系数001()nnjjx a thkbtjdtnkjanh ；( )011( )nbbjnkkaajkjj kxxClx dxdxbabaxx )() 1() 1() 1()(0nkkkkkkkjknkjj)!() 1(!knkkn多项式积分多项式积分00( 1)()0!()!n knnjj ktj dtknn knk，；，等距节点：ban1)()(21)(21)()(21!0!11)

8、1(21) 1(!1!01) 1(100)1(1101)1(0fTbfafabdxxfdttCdttCba；，dtjtknknCnnkjjknnk 00)()(!)( !) 1(梯形公式梯形公式nk,1,0- - 线性插值线性插值直线近似替代曲线直线近似替代曲线 xyoba)(af)(bf)(xf)0)(时xf，梯形曲边梯形SS；，等距节点：bbaan2222(2)0012(2)1002(2)20( 1)1(1)(2)2 0! 2!6( 1)4(2)2 1! 1!6( 1)1(1)2 2! 0!6CttdtCt tdtCt tdt；SimpsonSimpson公式公式()( )( )4 ()(

9、 )62( )babaabf x dxf afffbSxyoba)(af)(bf)(xf)2(baf2ba616164)0)(时xf形抛物线为曲边的曲边梯为曲边的曲边梯形SSxf)(- - 抛物线插值抛物线插值抛物线近似替代曲线抛物线近似替代曲线 ；，等距节点：4321044kkabaxnk；)()(7)(32)(12)(32)(790)()(43210fCxfxfxfxfxfabdxxfba注：注：( )( )( )12837nnnkkn kCCCnCotesNewtonCotesnRunge不难看出具有某种对称性：；当时，系数有正有负，这时稳定性没有保证， 因而实际计算不用高阶公式；高次可

10、能产生现象；CotesCotes公式公式907903290129032907818383816164612121)4()3()2()1(kkkkCCCC一般，一般， 梯形公式梯形公式 SimpsonSimpson公式公式 NewtonNewton公式公式 CotesCotes公式公式 )( fT)( fS)( fN)( fC次代数精度；具有公式证明3)()2(4)(6)(bfbafafabfSSimpson例：例：(1)( )1( )( )(14 1)6bsaf xIfdxbabaS fba，；证明：证明：22221(2)( )( )()21( )2()()62bsaf xxIfx dxbab

11、aS fababba，；2233222331(3)( )( )()31( )()()63bsaf xxIfx dxbabaS faabbba，；44554441(5)( )( )()5()( )( )64bsasf xxIfx dxbabaabS fabIf，；33441(4)( )( )()4( )( )bsasf xxIfx dxbaS fIf，；.3次代数精度具有公式即SimpsonTH3.2 当阶当阶 n 为为偶数偶数时，牛顿时，牛顿-柯特斯公式柯特斯公式至少有至少有 n+1 次代数精度次代数精度 . nkknknxfCabI0)()()( 2210200220( )()()/ 2()

12、0nnnbnniiannninnninf xxRfxx dxxath hti dttunhui du证明：只要证 为偶数时，对的余项为零即可。事实上，因为被积函数是奇函数。(1)111( )( )( )( )(1)!1( ) ()()2!bbnTnnaanbaRfR x dxfx dxnfxa xb dx （不变号）0)()(2bxaxx广义积分中值定理广义积分中值定理 3()( )()12bafa b ， 。1( )()()2bafxa xb dx；，)()(290)()()4(55bafabfRs方法一：方法一：事实上，由于事实上，由于SimpsonSimpson公式具有公式具有三次代数精

13、度三次代数精度，因而，因而对对三次多项式准确成立三次多项式准确成立，作满足条件的，作满足条件的3 3次次HermiteHermite插值插值 ，)(3xp件：处有插值条件及导数条，：在节点bbaa2；，)()()2()2()2()2()()(3333bfbpbafbapbafbapafap(4)231( )( )( ) () () ()4!2abf xp xfxaxxb则有：；3 3次代数精度次代数精度插值条件插值条件333( )4()( )62baabp app b ( )4 ()( )62baabf aff b( )S f33( )( ) )bap x dxS p x；，)()(290)(

14、)2() 1() 1(!4)()()2()()(!41)()2()()(!41)()()()()()4(551152)4(2)4(2)4(3bafabdtabtttfdxbxbaxaxfdxbxbaxaxfdxxpdxxffSfIfRbabababastabbax22令方法二：方法二：22( ), ,()()()22ababababR xfx abxa xxb以节点 ， 的二次插值多项式的余项为 , ,()()()22bsaababRffx abxa xxb dx221, ,() () )42baabfx ab dxaxb221() (), , ,42badabxaxbfx x ab dxdx

15、 221() (), , ,42baabxaxbfx x ab dx 22(4)1() ()( )44!baxaxbfdx 广义积分中值定理广义积分中值定理 (4)22( )() ()96bafxaxb dx 5(4)()( ) , 2880bafa b (4)221( )() ()44!bafxaxbdx 同理可得，同理可得， )( fC具有具有5 5次代数精度次代数精度； 。，)()(2945)()()6(117bafabfRc 由误差公式可知由误差公式可知区间过大区间过大，误差亦大误差亦大；为避免可；为避免可选取适当多的节点，即选取相对选取适当多的节点，即选取相对高阶高阶的的Newton

16、-cotesNewton-cotes公式；但由稳定性分析又知：当公式；但由稳定性分析又知：当 稳定稳定的现象；的现象； 8n时，会出现时，会出现不不 111100101( )( )( )()2()2( )()2)iinnbxiiaxiinninihf x dxf x dxf xf xhf xf xf xTf；0101nia b naxxxbxaibhinahn将，等分，记，则， ， ，21032( )1()()( )( )( )( )( )1212niibnaf xCabTha bffnnhbaf x dxTfffh 若， 则由连续函数的介值， ，；。3103110( )( )( )12( )

17、12nbniainiiiiihf x dxTffhfxx ，；3()( )( )( )12TbaT fRffa b 由的误差：，T Tn n的积分余项为的积分余项为0111121()2niiiiiibana baxxxbxaihhnxxxxxSimpson等分 ， ，；记， 对各子区间，用公式得：111110021101112( )( ) ( )4 ()()6 ()2( )()4()( )6iinnbxiiaxiiinninniiihf x dxf x dxf xf xf xhf xf xf xf xSf；5(4)515(4)1015(4)0()( )( )( )90 21( )( )( )(

18、 )9021( )( )902snbniiiiainiibaS fRffa bhf x dxSffxxhf 由的误差：，；41(4)(4)05(4)(4)4( )1()( )( )( )( )( )( ) ( )9021802niibnaf xCa bTha bffnnhbahf x dxSfff 若， ，则由连续函数介值， ，；10413124( )7 ( )32 ()9012 ()32 ()7 ()( )nbiaiiiniihf x dxf xf xf xf xf xCf则124iiiiiiibaa b nhxaihxxhnxxhxxhin将，等分，记， ，；，)

19、()4()(945)(2)(2945)()(6)6(10)6(117bahfabfhfCdxxfniinba？数及步长分别为多少形公式计算时所需节点，用复化梯，要求误差小于给定63110sindxxex例：例：误差事先估计误差事先估计 积分余项为积分余项为31222sin( )( )12122( )( )( )63xnIex dxbaITfh fffnn 记，则；解：解：( )2sin()cos()442sin()2cos2xxxfxexxexex；( )sin( )(sincos )2sin()4xxxf xexfxexxex，；31313362max( )max 2cos22( )2104

20、926.043xxxnfxexeITfenn；24927492849272313231nSimpsonnh取，即取个节点， 步长为；类似可对复化公式加以讨论得：；即取个节点；步长为； 由复化求积公式的截断误差可知，加密节点可以提高由复化求积公式的截断误差可知，加密节点可以提高求积公式的精度，但求积公式的精度，但困难困难在于：使用公式之前需给出在于：使用公式之前需给出合适合适的的步长步长，h h过大，满足不了精度；过大，满足不了精度；h h过小，计算量过大，而事先过小，计算量过大，而事先给出等分数很困难，因此将研究一种自动变步长的算法。给出等分数很困难，因此将研究一种自动变步长的算法。 , a

21、bbahnkxakh0,1,2,.,kn1,kkxx( )1()()2knkkhTf xf x将积分区间n等分， 在子区间上采用梯形公式，则 ，11( )( )01( )2()( )2nnkknnkkkhTTf af xf b从而可得n 等分时的复化梯形公式 121102111100211022( ) ( )2 ()()21 ( )()()2221( )()22nniiiinniiiiinniihTff xf xf xhhf xf xf xhTff x；2110( )( )( )( )()2nnnniiiTfTfhbaTff xf xhn计算的简化利用已得的；，；将子区间1,kkxx再二等分一

22、次并分别采用梯形公式， )()(2)(22121)(2,iiiinxfxfxfhT按这种方式不断计算所得结果将越来越精确，但对于我们预先给定的精度要求，何时停止计算呢？32()112()b annITf 32()2212(2 )()b annITf 12()()24ffnnITIT221( )( )( )( )3nnnI fTfTfTf2( )( )nnTfTf可用作为不太严格的终止准则。 , a b1,2,4,8,.,2 ,.,k通常积分区间的等分数取为递推公式可改成 1121221, 2 , 1 ),2) 12(22),()(2kkkikkkabiafabTTbfafabT用复化梯形公式计

23、算 精确至3位有效数字。12021dxx ，例：110() (0)(1)1.52ibaTfffhi ；令令；解：1211()( )1.55222；TTff24113()( )( )1.56562444；TTfff 4811357() ( )( )( )( )1.5695288888；TTfffff 824()()0.003110129.57TfTfI ；221()()()()3因为；nnnI fTfTfTf 用右端误差来修正用右端误差来修正 得：得：2221( )( )( )( )34( )( )4 1nnnnnI fTfTfT fTfT f 2( )nTf ，()nSf 通过适当组合Tn,T

24、2n能使代数精度提高2次.同理,对Sn,S2n组合 54142880( )()()nbaISffn 542242880 2( )()()nbaISffn (4)(4)12()()224ffnnISIS2241()()()()21nnnI fSfSfSf ；222221()()()()14()(5)41nnnnnI fSfSfSSfSff ()nCf 同样地，对复化Cotes公式进行类似处理，可得代数精度为7次的Romberg公式。2261()()()()21nnnI fCfCfCf ；232231()()()()4()()4()163nnnnnnCfCfRfI fCfCfCf 进一步，对进一步

25、，对Romberg公式进行组合公式进行组合42244()R ()2561()41255255nnnnnRffEfRR 但因系数接近但因系数接近1及误及误差的积累，使组合失差的积累，使组合失去组合意义去组合意义 这种将粗糙的复化梯形公式逐步加工成精度较高的求积这种将粗糙的复化梯形公式逐步加工成精度较高的求积公式的方法称为公式的方法称为 Romberg Romberg方法方法; ;又称又称逐次分半加速收敛法逐次分半加速收敛法. . 1(0)0(1)2( )000(1)( )( )11( )( )2(21)222441kkkkkkmkkkmmmmbaTf af bTbabaTf aiTTT写成便于程

26、序设计的计算格式，即写成便于程序设计的计算格式，即；要求误差不超过的近似值，积分法求：用5102102114dxxIRomberg例：例：248161248124121RCSTRCSTCSTSTT解：解：121415925. 31389885. 31421179. 31415687. 31311765. 33133333. 31 . 33RC；1311765. 3)43()41(41213133333. 331341 . 3)21(21213)()(224121121ffTTTTSfTTbfafabT。；14159. 31021101587. 01511415925. 331341389885

27、. 3)87()85()83()81(81211421179. 315115161415687. 331344552448448121242SISSTTSffffTTSSCTTS注：注：计算计算S1用到用到T2,用到用到3个点的函数值，代数精度为个点的函数值，代数精度为3次；计算次；计算C1，需用到，需用到T4,用到用到5个点的函数值，代数精度为个点的函数值，代数精度为5次；计算次；计算R1时，时，用到用到T8,用到用到9个点函数值，代数精度为个点函数值，代数精度为7次，所以次，所以Romberg公式公式不是插值型求积公式。不是插值型求积公式。 0 21nbkkkkakkxAf x dxA f

28、 xGaussxGa ssnu适当选择 和 ，使求积公式：具有次代数精度， 称为公式。节点 称点。 1-11-1 1,122=222baf x dxabbaxttbaabbaf x dxft dt 不失一般性，仅讨论积分形式，否则令则 100-110-100100-11=2=2,0=0=fx dxA fxdxAAxxdxA x具有 次代数精度即1Gauss点公式故一点Gauss公式为 11( )2 (0)f x dxf即中矩形公式就是一点Gauss公式 10011-1101-1101001 1-1122201001 1-113320013-131=2=1=0=32=33=0=f x dxA f

29、 xA f xdxAAAAxdxA xAxxxx dxA xAxx dxA xAx 具有 次代数精度即2Gauss点公式所以二点Gauss公式为133331( )()()f x dxff 多点多点Gauss求积公式类似可求得，但求解非线性方程组求积公式类似可求得，但求解非线性方程组较复杂，通常较复杂，通常n2就很难求解故一般不通过解方程求就很难求解故一般不通过解方程求 xk 及及 Ak 而从研究而从研究高斯点的基本特性高斯点的基本特性来着手解决来着手解决Gauss 求积求积公式的构造问题公式的构造问题 对于插值型求积公式，对于插值型求积公式，x0 xn 为为 Gauss 点点 与任意次数不大于

30、与任意次数不大于n 的多项式的多项式 P(x) （带权）正交，即（带权）正交，即 nkkxxxw0)()(定理定理11( ) ( )0 x p x dx求求 Gauss 点点 求求w(x)的的零点零点0)(d)()()(1 banxxxxxP 证明：证明： “” x0 xn 为为 Gauss 点点, 则公式则公式 至少有至少有 2n+1 次代数精度。次代数精度。 bankkkxfAdxxfx0)()()( 对任意次数对任意次数不大于不大于n 的多项式的多项式 Pm(x)， Pm(x) w(x)的次数的次数不大于不大于2n+1，则代入公式应则代入公式应精确成立精确成立： nkkkmkbamxwx

31、PAdxxwxPx0)()()()()( 0= 0 “” 要证明要证明 x0 xn 为为 Gauss 点，即要证公式对任意次点，即要证公式对任意次数数不大于不大于2n+1 的多项式的多项式 Pm(x) 精确成立，即证明：精确成立，即证明： nkkmkbamxPAdxxPx0)()()( 设设)()()()(xrxqxwxPm bababamdxxrxdxxqxwxdxxPx)()()()()()()( 0 nkkkxrA0)( nkkmkxPA0)( Legendre 多项式族：多项式族：1)( x 定义在定义在 1, 1上，上，kkkkkxdxdkxP)1(!21)(2 满足：满足： lkl

32、kPPklk1220),(xPP 10, 1由由 有递推有递推11)12()1( kkkkPxPkPk以以 Pn+1 的根为节点的求积公式称为的根为节点的求积公式称为Gauss-Legendre 公式公式。注：注：一般一般a,b上的积分可化为上的积分可化为-1,1上特殊高斯公式进行计算。上特殊高斯公式进行计算。 (22)11( )2(22)!011 Gauss-legendre( )()( )nnfkknkkR ff x dxA f xx dxxn+1Legendre对求积公式，其余项为其中 为次的多定理项式零点。n n0 00 02 21 11 12 20 0(0.8888889)(0.88

33、88889)(0.5555556)(0.5555556)3 30.86113630.86113630.33998100.33998100.34785480.34785480.65214520.65214524 40 00.53846930.53846930.90617980.90617980.56888890.56888890.47862870.47862870.23692690.2369269kxkA)5773503. 0(33)7745967. 0(5159895部分部分Gauss-Legendre公式的节点和系数公式的节点和系数稳定性：稳定性：1 1）因为积分式对函数）因为积分式对函数f

34、(x)=1f(x)=1准确成立，则准确成立，则abAknk02 2）因求积公式对）因求积公式对2n2n次多项式次多项式 也准确成，则也准确成，则)()(2xlxfk220( )( )0bnkk kikiaAA lxlx dx若若 有误差有误差 ,则由此引起最终结果的误差为则由此引起最终结果的误差为 )(kxfk001 10.nnnkkAAAAba即即Gauss求积公式有好的求积公式有好的数值稳定数值稳定性性收敛性：收敛性：0( ),()( )bnnkkkaf xC a bA f xf x dx 可以证明，若有高斯型求积公式Gauss型求积公式除有好的收敛性和稳定性之外，而且是具型求积公式除有好

35、的收敛性和稳定性之外，而且是具有最高代数精度的数值求积公式，即对有最高代数精度的数值求积公式，即对n+1个求积节点构成个求积节点构成的求积公式不可能具有的求积公式不可能具有2n+2次代数精度。次代数精度。 20()niif xxx对对2n+2次多项式次多项式20( )()0bbniiaaf x dxxxdx 2000()0nnnkkkikkiA f xAxx对于确定的函数表达式对于确定的函数表达式 f(x)，函数的导数值总是可以按各种求，函数的导数值总是可以按各种求导的公式和规则或是直接通过导数定义的极限积分求得，但当导的公式和规则或是直接通过导数定义的极限积分求得，但当函数是以数据表的形式给

36、出时，微分学中的方法就不能使用。函数是以数据表的形式给出时，微分学中的方法就不能使用。 x x0 x1 x2 xn y y0 y1 y2 yn 对于列表函数对于列表函数 y = f (x)： 插值多项式插值多项式 y = Pn(x)作为它的近似作为它的近似, 我们取我们取统称统称插值型的求导公式插值型的求导公式)()(xPxfn )(xPn )(xf 的近似值，的近似值，作为作为建立的数值公式建立的数值公式 依据插值余项定理，求导公式的依据插值余项定理，求导公式的余项余项为为)(dd)!1()()()!1()()()()1(11)1( nnnnnfxnxxnfxPxf式中式中. )()(01

37、niinxxx 我们我们限定限定:求某个节点求某个节点 xk 上的导数值，上面的第二项变上的导数值，上面的第二项变为零，这时有余项公式为零，这时有余项公式)()!1()()()(1)1(knnknkxnfxPxf 下面我们仅仅考察下面我们仅仅考察节点处的导数值节点处的导数值为简化讨论为简化讨论, 假定假定所给的节点是等距的所给的节点是等距的 两点公式两点公式 已给两节点已给两节点 x0, x1 上的函数值上的函数值 f (x0), f (x1)，做线性插值，做线性插值)()()(101001011xfxxxxxfxxxxxP 记记 x1 x0 = h，对上式两端求导，有，对上式两端求导，有)(

38、)(1)(101xfxfhxP )()(1)( )()(1)(01110101xfxfhxPxfxfhxP 于是有下列求导公式：于是有下列求导公式：)(2)()(1)()(2)()(1)(011010 fhxfxfhxffhxfxfhxf 而利用余项公式而利用余项公式(3.43)知，带知，带余项的两点公式余项的两点公式是：是： 设已给出三节点设已给出三节点x0, xl=x0+h, x2=x0+2h上的函数值上的函数值,做二次做二次插值插值)()()()()()()()()()(2120210121012002010212xfxxxxxxxxxfxxxxxxxxxfxxxxxxxxxP 令令 x = x0+th, 则则 )()12()()44()()32(21)()()1(21)()2()()2)(1(21)(2100221002xftxftxfththxPxfttxfttxfttthxP 三点公式三点公式 这里撇号这里