数值积分与数值微分课件

第5章数值积分与数值微分方法关于定积分计算，已经有较多方法，如公式法、分步积分法等，但实际问题中，经常出现不能用通常这些积分方法计算的定积分问题。

怎样把这些通常方法失效的定积分在一定精度下快速计算出来，特别是通过计算机编程计算出来就是本章研究的内容。

此外，怎样根据函数在若干个点处的函数值去求该函数的导数近似值也是本章介绍的内容。

本章涉及的方法有Newton-Cotes求积公式、Gauss求积公式、复化求积公式、Romberg求积公式和数值微分。

5.1 引例人造地球卫星轨道可视为平面上的椭圆。

我国的第一颗人造地球卫星近地点距离地球表面439km，远地点距地球表面2384km，地球半径为6371km，求该卫星的轨道长度。

本问题可用椭圆参数方程cos ,,0sin x a t a b y b tπ=⎧≤≤>⎨=⎩ (0t 2) 来描述人造地球卫星的轨道，式中a, b 分别为椭圆的长短轴，该轨道的长度L 就是如下参数方程弧长积分但这个积分是椭圆积分，不能用解析方法计算。

5.2问题的描述与基本概念要想用计算机来计算()ba f x dx ⎰，应对其做离散化处理。

注意到定积分是如下和式的极限 01()lim ()nbi i a i f x dx f x λξ→==∆∑⎰ 要离散化，做去掉极限号lim将i ξ取为具体的i x 值为减少离散化带来的误差，将i x ∆用待定系数i A 代替 于是就得到定义 5.1 若存在实数1212,,,;,,,,n n x x x A A A 且任取()[,],f x C a b ∈都有1()()nb i i a i f x dx A f x =≈∑⎰ （5.1） 则称式(5.1)为一个数值求积公式。

式中i A 称为求积系数，i x 称为求积节点，而称 1()()()nb i i a i R f f x dx A f x ==-∑⎰ （5.2） 为求积余项或求积公式（5.1）的截断误差。

从定义可以看到，数值求积公式依赖于求积节点个数n 、求积节点{}i x 和求积系数{}i A ，这三个量有一个发生变化，则产生不同的求积公式。

定义5.2 若求积公式1()()nbi i a i f x dx A f x =≈∑⎰ 对所有不超过m 次的多项式()m p x 有求积余项()0m R p =，而对某一个m+1次多项式1()m p x +有1()0m R p +≠，则称该求积公式的代数精度为m 。

1()()()nb i i a i R f f x dx A f x ==-∑⎰一般，一个求积公式的代数精度越大，则该求积公式越好。

确定代数精度的方法依次取()(0,1,)k f x x k ==代入公式1()()()nb i i a i R f f x dx A f x ==-∑⎰ 并验证()0k R x =是否成立。

若第一个使()0k R x=不成立的k 值为m ，则对应的代数精度为m-1。

例 5.1确定求积公式11()((33f x dx f f -≈-+⎰ 的代数精度。

解 取()kf x x =代入求积公式有11()[((]331[(](1(1))13k k k k k k R x x dx k -=--+=-+-+⎰易验证 0123()()()()0R x R x R x R x ====，但48()045R x =≠，故本题求积公式代数精度为3。

例 5.2确定下面求积公式22()()(0)(2)h h f x dx Af h Bf Cf h -≈-++⎰ 的参数A ，B ，C ，使它具有尽可能高的代数精度，并指出相应的代数精度。

解 本题要先求出具体的求积公式，然后再判断所求公式的代数精度。

公式有3个待定参数，h不是求积公式的参数，故利用3个条件得到的3个等式关系就可以解决求出具体求积公式的问题。

依次取2()1,,f x x x =代入求积公式并取等号，有解之得故所求的求积公式为 221648()()(0)(2)939hh h h h f x dx f h f f h -≈-++⎰ 为确定其代数精度，再取3()f x x =代入求出的公式继续计算，有3416()03R x h =-≠，故所求的求积公式具有二阶代数精度。

5.3 插值型求积公式借助多项式插值函数来构造的求积公式称为插值型求积公式。

一般选用不同的插值公式就可以得到不同的插值型求积公式。

基本思想利用被积函数()f x的插值函数代替()f x做定积分的近似计算来构造求积公式。

1.构造原理考虑()f x 在n 个节点12,,,n x x x 上的n-1次Lagrange 插值多项式1()n L x -与()f x 的余项，有()()11()()()[,]!n ni in n i f f x f x l x x x a b n ξω-==∑+∈（） 这里()()()111,nnk in n k kk i k k ix x l x x x x x x ω-==≠⎛⎫-=∏=- ⎪-⎝⎭∏。

两边取积分,有()()11()()()!n nbbbiin n aa ai f f x dx f x lx dx x dx n ξω-==∑+⎰⎰⎰（）记1()b i in aA l x dx -=⎰ （5.3）则有（5.4） 若舍去()()()!n bn af x dx n ξω⎰，得求积公式 1()()nbi i ai f x dx A f x =≈∑⎰求积系数1()bi in aA l x dx -=⎰的求积公式就是插值型求积公式。

插值型求积公式的求积余项()()()()!n bn af R f x dx n ξω=⎰当()f x 为次数小于n 次的多项式时，有()()0n f x ≡，对应的()0R f =。

因此插值型求积公式的代数精度至少为n-1。

若取()1f x =，代入式(5.4)，可得插值型求积公式的求积系数之和为1nkk Ab a ==-∑下面具体介绍常用的几个插值型求积公式。

2. Newton-Cotes 求积公式1) n 点的Newton-Cotes 公式的构造将求积节点i x 取为[a,b]上的等距节点2,,n做积分变量变换:x a th =+则当[],xa b ∈时,有[]0,1t n ∈-,于是有插值型求积公式的求积系数为()(),1,2,,1(5.5)n i i A b a C i n =-=-()n i C 常称为Cotes 系数，易验证()11nn i i C ==∑通常称()()1()()n bn i i ai f x dx b a C f x =≈-∑⎰（5.6）为n 点的Newton-Cotes 公式。

由于求积节点i x 是等距的，因此也称式（5.6）为等距节点求积公式。

Newton-Cotes公式A) 2 点的Newton-Cotes公式（5.7）这正是我们熟悉的梯形公式。

B) 3点的Newton-Cotes公式为（5.8）称它为Simpson公式或抛物线公式。

表5.1n ()niC 23456789例5.3 试分别用梯形公式和Simpson公式计算1 0sin xI dxx=⎰解用梯形公式计算，有用Simpson公式计算，有2）n 点Newton-Cotes 公式的代数精度定理5.1 当求积节点个数n 为奇数时，对应的Newton-Cotes 求积公式的代数精度至少为n 。

证明 由于是插值型求积公式，故有()0,0,1,,1kR xk n ==-对nx有()()()()()()11112121111!1n nnnx a thb n x nn kak k n l l s t ll n lk x R xx x dxht k dtn h s l k dsξ=+-=+===++=-+-==-=-+=+-+∏∏⎰⎰∏⎰记()()()21111l lm l k k m ls s l k s m ϕ+=-+==-=+-+=+∏∏，易知()()s s ϕϕ-=-，故是奇函数，得()()10lnn lR x h s ds ϕ+-==⎰，得证。

3） 梯形公式与Simpson 公式的余项引理 5.1 (积分中值定理)设()()[],,f x g x C a b ∈，且()g x 在[],a b 上不变号，则有()()()()[],bb aaf xg x dx f g x dxa b ξξ=∈⎰⎰梯形公式余项为()()()()1[,]2!baf R f x a x b dx a b ξξ''=--∈⎰()()x a x b --在[a,b]不变号，()f x 有2阶连续导数，由引理 5.1，有()()()()()()31[,]2!12baf b a R f x a x b dx f a b ηηη''-''=--=-∈⎰梯形公式余项（5.9）抛物线公式的余项（5.10）4）Newton-Cotes 公式的数值稳定性设计算函数()k f x 时产生舍入误差为k ε实际在计算机中参加计算的是()k f x 的近似值()()1,2,k k kf x f x k n ε=-=故Newton-Cotes 公式在计算机中产生的误差为()()()()()()111()()nnn n n i i i i i i nn i ii b a C f x b a C f x b a C ηε====---=-∑∑∑若记1max ii nεε≤≤=，则有()()1nn n i i b a C ηε=≤-∑由Cotes 表5.1,当时, ()0niC >，从而有 ()()()1nn n i i b a C b a ηεε=≤-=-∑说明此时计算舍入误差可以控制，从而Newton-Cotes 公式是数值稳定的。

但当n>8时，()n iC 有正有负，()1nn ii C =∑ 随n 增大而增大，从而导致舍入误差增加。

故n>8时，Newton-Cotes 公式是数值不稳定的。

因而一般不用n>8的Newton-Cotes 公式来做定积分计算。

3. Gauss 求积公式1）Gauss 求积公式的构造与概念n 点的插值型求积公式的代数精度至少是n-1 ,那么是否还能提高其代数精度呢？若能，其代数精度最大能是多少？为回答这个问题，观察一下插值求积公式的构造方法，发现其至少具有n-1次代数精度的结论是在限定求积节点i x 的前提下得出的，若让求积节点i x 也可以自由取值，则就给提高代数精度创造了条件。

为使问题讨论更具一般性，这里考虑带权的定积分求积公式()()()1nbkkak x f x dx A f x ρ=≈∑⎰ （5.11）式中()x ρ是已知的非负函数，为区间[a,b]上的权函数，[a,b]可以取为∞。