《数值分析讲义第四章数值积分》由会员分享,可在线阅读,更多相关《数值分析讲义第四章数值积分(77页珍藏版)》请在金锄头文库上搜索。
1、第四章第四章 数值积分与数值微分数值积分与数值微分 /* Numerical Integration and differentiation*/近似计算近似计算1 引言引言 对对f()采用不同的近似计算方法,从而得到各采用不同的近似计算方法,从而得到各种不同的求积公式。种不同的求积公式。 以上三种方法都是用被积函数值的以上三种方法都是用被积函数值的线性组合线性组合来表示积来表示积分值分值。推广,一般地有推广,一般地有求积节求积节点点求积系数,与被求积系数,与被积函数无关积函数无关像这样,将积分用若干节点上被积函数值的像这样,将积分用若干节点上被积函数值的线性组合线性组合来表示来表示的数值积分公
2、式称为的数值积分公式称为机械求积公式机械求积公式。 求积误差求积误差 机械型求积公式的构造归结为,确定求积节点机械型求积公式的构造归结为,确定求积节点xk和求积系和求积系数数Ak,使在某种意义下精确度较高。总之,要解决三个问,使在某种意义下精确度较高。总之,要解决三个问题:题:1.精确度的度量标准;精确度的度量标准;2.如何构造具体的求积公式;如何构造具体的求积公式;3.具体求积公式构造出来后,误差如何估计?具体求积公式构造出来后,误差如何估计? 以上三种方法都是用被积函数值的以上三种方法都是用被积函数值的线性组合线性组合来表示积来表示积分值分值。推广,一般地有推广,一般地有求积节求积节点点求
3、积系数,与被求积系数,与被积函数无关积函数无关像这样,将积分用若干节点上被积函数值的像这样,将积分用若干节点上被积函数值的线性组合线性组合来表示来表示的数值积分公式称为的数值积分公式称为机械求积公式机械求积公式。 求积误差求积误差 机械型求积公式的构造归结为,确定求积节点机械型求积公式的构造归结为,确定求积节点xk和求积系和求积系数数Ak,使在某种意义下精确度较高。总之,要解决三个问,使在某种意义下精确度较高。总之,要解决三个问题:题:1.精确度的度量标准;精确度的度量标准;2.如何构造具体的求积公式;如何构造具体的求积公式;3.具体求积公式构造出来后,误差如何估计?具体求积公式构造出来后,误
4、差如何估计? 以上三种方法都是用被积函数值的以上三种方法都是用被积函数值的线性组合线性组合来表示积来表示积分值分值。推广,一般地有推广,一般地有求积节求积节点点求积系数,与被求积系数,与被积函数无关积函数无关像这样,将积分用若干节点上被积函数值的像这样,将积分用若干节点上被积函数值的线性组合线性组合来表示来表示的数值积分公式称为的数值积分公式称为机械求积公式机械求积公式。 求积误差求积误差 机械型求积公式的构造归结为,确定求积节点机械型求积公式的构造归结为,确定求积节点xk和求积系和求积系数数Ak,使在某种意义下精确度较高。总之,要解决三个问,使在某种意义下精确度较高。总之,要解决三个问题:题
5、:1.精确度的度量标准;精确度的度量标准;2.如何构造具体的求积公式;如何构造具体的求积公式;3.具体求积公式构造出来后,误差如何估计?具体求积公式构造出来后,误差如何估计?若若某某个个求求积积公公式式对对次次数数 m 阶阶的的多多项项式式准准确确成成立立,而而对对m+1 阶阶的的多多项项式式不不一一定定准准确确成成立立。即即对对应应的的误误差差满满足足:R Pk =0 对对任任意意 k m 阶阶的的多多项项式式成成立立,且且 R Pm+1 0 对对某某个个 m+1 阶多项式成立,则称此求积公式的阶多项式成立,则称此求积公式的代数精度代数精度为为 m 。代数精度与误差的关系代数精度与误差的关系
6、:代数精度越高,求积误差越小。:代数精度越高,求积误差越小。结论:结论:问题问题1由上面代数精度条件确定求积公式可分两种情形:由上面代数精度条件确定求积公式可分两种情形:1.若事先给定求积节点若事先给定求积节点xk(k=0,n),例如被积函数以表的形式例如被积函数以表的形式给出时给出时xk确定,可令确定,可令m=n,由上式确定由上式确定n+1个系数个系数Ak即可即可-待定系数法和插值法待定系数法和插值法。2.若若xk和和Ak都可选择,令都可选择,令m=2n +1,确定,确定xk和法和法Ak -Gauss法法要使求积公式具有要使求积公式具有m阶代数精度,则它对阶代数精度,则它对1,x,xm均准确
7、成立,均准确成立,即即m+1个方程,个方程,2n+2个未知数个未知数问题问题2Case 1-方法方法11 插值型求积插值型求积 公式公式思思路路利用利用插值多项式插值多项式 则积分易算。则积分易算。 在在a, b上取上取 a x0 x1 xn b,做,做 f 的的 n 次插值次插值多项式多项式 ,即得到,即得到Ak由由 决定,决定,与与 无关。无关。节点节点 f (x)插值型积分公式插值型积分公式/*interpolatory quadrature*/误差误差Case 1-方法方法21 Newton-Cotes Formulae例:例:对于对于a, b上上1次插值,有次插值,有考察其代数精度。
8、考察其代数精度。f(x)abf(a)f(b)梯形公式梯形公式/* trapezoidal rule*/解:解:逐次检查公式是否精确成立逐次检查公式是否精确成立代入代入 P0 = 1:=代入代入 P1 = x :=代入代入 P2 = x2 : 代数精度代数精度 = 11 Newton-Cotes FormulaeTh1.形如形如 的求积公式至少有的求积公式至少有 n 次代数精度次代数精度 该该公式为公式为插值型插值型(即:(即: ) 当节点当节点等距分布等距分布时:时:令令Cotes系数系数注:注:Cotes 系数仅取决于系数仅取决于 n 和和 i,可查表得到。与,可查表得到。与 f (x) 及
9、及区间区间a, b均无关。均无关。 2 Newton-Cotes 公式公式NewtonCotes formula1 Newton-Cotes Formulaen = 1:Trapezoidal Rule/* 令令 x = a+th, h = b a, 用中用中值定理值定理 */代数精度代数精度 = 1n = 2:Simpsons Rule代数精度代数精度 = 3n = 4: Cotes Rule, 代数精度代数精度 = 5,偶数阶偶数阶N-C公式具公式具有有n+1阶代数精度阶代数精度n = 3: Simpsons 3/8-Rule, 代数精度代数精度 = 3,对称节点的系数相同对称节点的系数相
10、同Cotes公式是公式是用不同节点用不同节点的函数值的函数值(高度)的(高度)的加权平均来加权平均来近似区间的近似区间的平均高度平均高度注:当注:当n8时,时,Cotes系数有负,造成公式不稳定,因此常用系数有负,造成公式不稳定,因此常用低阶低阶Cotes公式。公式。证明:只需证明n为偶数时,为偶数时, N-C公式对公式对f(x)=xn+1的余项的余项R(f)=0即可。即可。因因 f(n+1)(x)=(n+1)!, 由余项公式得由余项公式得Th2. n为偶数时为偶数时, N-C公式至少具有公式至少具有n+1阶代数精度。阶代数精度。注:当注:当n 为偶数时,为偶数时,Cotes公式具有公式具有n
11、+1阶精度,与阶精度,与n+1阶阶Cotes公式精度相同,但少计算一个节点上的函数值,公式精度相同,但少计算一个节点上的函数值,因此一般常用偶数阶因此一般常用偶数阶Cotes公式。公式。偶数阶偶数阶N-C公式具公式具有有n+1阶代数精度阶代数精度N-C公式具有公式具有n阶代数精度阶代数精度余项余项R=o(h n+2)Hint:construct a interpolation polynomial of order 5, H(x), satisfying H(a)=f(a), H(b)=f(b), H (k)(a+b)/2) = f (k)(a+b)/2).HW: p.151-152 #1-6
12、数值稳定性的一般概念数值稳定性的一般概念N-C的稳定性的稳定性3 复合求积复合求积 /* Composite Quadrature */Havent we had enough formulae? Whats up now?Oh come on, you dont seriously consider h=(b a)/2 acceptable, do you?Why cant you simply refine the partition if you have to be so picky?Dont you forget the oscillatory nature of high-degr
13、ee polynomials!Uh-oh高次插值有高次插值有Runge 现象现象,故采用分段低次插值,故采用分段低次插值 分段低次合成的分段低次合成的 Newton-Cotes 复合复合求积公式。求积公式。 复合梯形公式:复合梯形公式:在每个在每个 上用梯形公式:上用梯形公式:= Tn/*中值定理中值定理*/2 Composite Quadrature 复化复化 Simpson 公式:公式:44444= Sn注:注:为方便编程,可采用另一记法:令为方便编程,可采用另一记法:令 n = 2n 为偶数,为偶数, 这时这时 ,有,有2 Composite Quadrature 收敛速度与误差估计:收
14、敛速度与误差估计: 若若一一个个复复化化积积分分公公式式的的误误差差满满足足 且且C 0,则,则称该公式是称该公式是 p 阶收敛阶收敛的。的。/*中值定理中值定理*/类似的,可得类似的,可得2阶收敛阶收敛4阶收敛阶收敛6阶收敛阶收敛例例1:计算计算解:解:其中其中= 3.138988494其中其中= 3.141592502运算量基运算量基本相同本相同Q: 给定精度给定精度 ,如何取,如何取 n ?例如:要求例如:要求 ,如何判断,如何判断 n = ?上例中若要求上例中若要求 ,则,则即:取即:取 n = 4092 Composite Quadrature事后误差估事后误差估计式,可用计式,可用
15、来判断迭代来判断迭代是否停止。是否停止。始步长始步长h0.510-24 龙贝格龙贝格积分积分 /* Romberg Integration */复化梯形公式算法简单,但精度较差,收敛速度(复化梯形公式算法简单,但精度较差,收敛速度(2阶收敛)阶收敛)较慢,如何提高收敛速度?较慢,如何提高收敛速度?注:按上面规律,可以构造线性组合系数为注:按上面规律,可以构造线性组合系数为的新的积分公式,但当的新的积分公式,但当m4时,前一个系数接近于时,前一个系数接近于1,后一,后一个系数接近于个系数接近于0,这样构造出的新公式与前一个公式结果差别,这样构造出的新公式与前一个公式结果差别不大,反而增加计算量,
16、因此实际上常做到不大,反而增加计算量,因此实际上常做到Romberg公式为公式为止。止。例:例:计算计算已知对于已知对于 = 10 6 须将区间对分须将区间对分 9 次,得到次,得到 T512 = 3.14159202由由 来计算来计算 I 效果是否好些效果是否好些?考察考察= 3.141592502 = S4一般有:一般有:Romberg 序列序列 Romberg 算法:算法: ? ? 0, k=0,1,n定理定理求积系数的求积系数的另一种计算另一种计算方法方法结论:结论:Gauss型积分公式是数值稳定型积分公式是数值稳定的。(稳定性的。(稳定性分析类似于分析类似于n 7的的N-C公式)公式
17、) 区间-1,1上权函数(x)=1的Gauss型求积公式,称为Gauss-Legendre求积公式求积公式,其Gauss点为Legendre多项式的零点. 几种几种Gauss型求积公式型求积公式 (1) Gauss-Legendre求求积公式公式 公式的Gauss点和求积系数可在数学用表中查到 .nxkAknxkAk10260.93246951420.66120938650.23861918610.17132449240.36076157300.467913934620.5773502692130.774596669200.55555555560.888888888970.9491079123
18、0.74153118560.405845151400.12948496620.27970539150.38183005050.417959183740.86113631160.33998104360.34785484510.652145154980.96028985650.79666647740.52553240990.18343464250.10122853630.22238103450.31370664590.362683783450.90617984590.538469310100.23692688510.47862867050.5688888889 用3点Gauss公式计算积分 解解
19、查表得x1=-0.7745966692,x2=0,x3=0.7745966692, A1=A3=0.5555555556,A2=0.8888888889, 所以有 Gauss-Legendre求积公式的余项为 例例13误差为 实际上,I=2sin1=1.68294197, 误差为|R|=6.15810-5 . 用Simpson公式,则有I1.69353487, 误差为|R|=1.0610-2 . 由于因此,a,b上权函数(x)=1的Gauss型求积公式为 用3点Gauss公式计算积分结果远比Simpson公式的结果精确.例例14 解解 这里Gauss点和积分系数与上例相同,所以 求积误差可表示
20、为 区间0,)上权函数(x)=e-x的Gauss型求积公式,称为Gauss-Laguerre求积公式求积公式,其Gauss点为Laguerre多项式的零点. (2) Gauss-Laguerre求求积公式公式 公式的Gauss点和求积系数可在数学用表中查到 .nxkAknxkAk20.58588643763.41421356230.85355339050.146446609450.26356031971.41340305913.59642577107.085810005812.64080084420.52175561050.39866681100.07594244970.00361175870
21、.000023370030.41577455672.29428036026028994508290.71109300990.27851773350.010389256560.22284660411.18893210162.99273632605.77514356919.837467418315.98287398060.45896467930.41700083070.11337338200.01039919750.00026101720.000000898540.32254768961.74576110114.53662029699.39507091230.60315410430.3574186
22、9240.03888790850.0005392947Gauss-Laguerre求积公式为 求积公式的误差为 由于 所以,对0, +)上权函数(x)=1的积分,也可以构造类似的Gauss-Laguerre求积公式: 区间(-,)上权函数(x)= (3) Gauss-Hermite求求积公式公式 公式的Gauss点和求积系数可在数学用表中查到 .nxkAknxkAk20.70710678110.886226925460.43607741191.33584907042.35060497360.72462959520.15706732030.004530009931.224744871300.29
23、540897511.816359000640.52464762321.65068012380.80491409000.081312835470.81628788281.67355162872.651961356300.42560725260.05451558280.00097178120.810264617550.95857246462.020182870400.39361932310.01995324210.9453087204 的Gauss型求积公式,称为Gauss-Hermite求积公式求积公式, 其Gauss点为Hermite多项式的零点. Gauss-Hermite求积公式为6 Numerical Differentiationh=0.8,精,精度度像这样,将微分的计算归结为在若干节点上函数值的像这样,将微分的计算归结为在若干节点上函数值的线性组线性组合合的数值微分方法称为的数值微分方法称为机械求导方法机械求导方法。二次插值二次插值f(xn)7.Gauss型积分公式型积分公式8.HW: p.152-153 #7-12