数值积分和数值微分yj.ppt

《数值积分和数值微分yj.ppt》由会员分享，可在线阅读，更多相关《数值积分和数值微分yj.ppt（54页珍藏版）》请在三一办公上搜索。

1、科大研究生学位课程,计算定积分有微积分基本公式,但很多函数找不到原函数，如,等。而实际上，有很多函数只知一些离散点的函数值，并无表达式，这就需要利用已知条件求出近似值。,第5章 数值积分与数值微分,科大研究生学位课程,数值积分/Numerical Integration/,定义数值积分如下：是离散点上的函数值的线性组合,称为求积系数,与f(x)无关,与积分区间和求积节点有关,称为数值积分公式,数值积分问题可分解为下述三个问题：,1、求积公式的具体构造问题；(包括xi的选取和Ai的构造),3、精确性程度的衡量标准问题。,2、余项估计问题(亦即误差估计问题)；,求积公式的误差 RfI*fIf,科大

2、研究生学位课程,1、解决第一个问题；节点xi 和系数Ai如何选取，即选取原则,两个目标：,1、余项估计问题；求积公式的误差 RfI*fIf尽可能小。,2、求积公式的代数精度尽可能高。,2、解决第二个问题；依赖插值多项式的余项估计公式。,3、对于第三个问题；引进代数精度的概念,科大研究生学位课程,定义5.1 若求积公式,对(x)=xj(j=0,1,2,m)都精确成立,但对(x)=xm+1不精确成立,即,则称此公式具有m次代数精度.,可见,若公式具有m次代数精度,则公式对所有次数不超过m的多项式都精确成立.,注意：,、求积公式的误差是计算精度的度量标志，而代数精度是求积公式优良性能的标志。,2、求

3、积公式的误差小，不代表代数精度高。代数精度高，也不代表求积公式的误差小。它们没有必然联系。,科大研究生学位课程,例 1 确定形如,的求积公式，使其代数精度尽可能高。,数值求积公式为,解 令公式对(x)=1,x,x2都精确成立,则,A0+A1+A2=3,A1+3A2=4.5,A1+9A2=9,解之得:A0=0,A1=9/4,A2=3/4.,例2 试确定参数A0,A1,A2,使求积公式,具有尽可能高的代数精度,并问代数精度是多少?,科大研究生学位课程,解得:A0=A2=1/3,A1=4/3.,求积公式为,当(x)=x3时,左=0,右=0,公式也精确成立.,解 令公式对(x)=1,x,x2 都精确成

4、立,则,A0+A1+A2=2,-A0+A2=0,A0+A2=2/3,当(x)=x4时,左=2/5,右=2/3,公式不精确成立.,所以,此公式的代数精度为3.,科大研究生学位课程,例3 试确定参数A0,A1和x0,x1,使求积公式,具有尽可能高的代数精度,并问代数精度是多少?,解 令公式对(x)=1,x,x2,x3都精确成立,则,A0+A1=2,A0 x0+A1x1=0,A0 x02+A1x12=2/3,A0 x03+A1x13=0,解得:,求积公式为,求积公式的代数精度为3。,科大研究生学位课程,5.1 插值型求积公式,在a,b上取 a x0 x1 xn b，做 f 的 n 次插值多项式，即得

5、到,节点,f(x),插值型积分公式,误差,科大研究生学位课程,梯形公式/*trapezoidal rule*/,解：逐次检查公式是否精确成立,代入 P0=1：,=,代入 P1=x：,=,代入 P2=x2：,代数精度=1,定理：形如 的求积公式至少有 n 次代数精度 该公式为插值型（即：）,科大研究生学位课程,为了简化计算,取等距节点xk=a+kh,(k=0,1,2,n,则有,令,则有,称为Newton-Cotes公式.Ck(n)称为Cotes系数.,(5.6),它不仅与函数f(x)无关，而且与积分区间a,b无关。,科大研究生学位课程,设(x)C2a,b,取n=1时的Newton-Cotes公式

6、并估计误差.,计算Cotes系数,于是有,5.2 几个常用的求积公式,从几何上看:用梯形的面积近似曲边梯形的面积。,所以公式,=T,也称为梯形公式,记为T.,5.2.1 梯形公式及其误差,科大研究生学位课程,称之为Simpson公式或抛物线公式，记为S.,构造三次多项式H3(x),使满足 H3(a)=(a),H3(b)=(b),于是有,证明Simpson公式对不高于三次的多项式精确成立,即,这时插值误差为,=S.,设(x)C4a,b,取n=2时的Newton-Cotes公式并估计误差.,解 计算Cotes系数,5.2.2 辛普森公式及其误差,科大研究生学位课程,于是有,科大研究生学位课程,由于

7、构造Newton-Cotes公式需要Cotes系数,将其列表如下:,科大研究生学位课程,牛顿求积公式：,代数精度=3,牛顿公式及其误差,科大研究生学位课程,取n=4的Newton-Cotes公式及误差.,查表可得,于是有,称之为Cotes公式，记为C。其误差为,其中,xk=a+kh,k=0,1,2,3,4,h=(b-a)/4.,代数精度=5,5.2.3 科茨公式及其误差,科大研究生学位课程,一般地,Newton-Cotes公式的截断误差为,例1 用梯形公式、Simpson公式和Cotes公式求积分.,的近似值。,解,IT=1/2*(4+2)=3,IS=1/6*(4+12.8+2)=3.1333

8、3,IC=1/90*(28+14)=3.14212,科大研究生学位课程,5.3 复化求积公式,高次插值有Runge 现象，故采用分段低次插值 分段低次合成的 Newton-Cotes 复化求积公式。,复化梯形公式：,在每个 上用梯形公式：,=Tn,科大研究生学位课程,可见,复化梯形公式是收敛的。而且,要使|RTn|,只要,如果记M2=,复化梯形公式的误差为,，则有,若在每个小区间上的积分采用Simpson公式,则可得到复化Simpson公式:,科大研究生学位课程,复化 Simpson 公式：,=Sn,误差为,如果记M4=,，则有,复化Simpson公式也是收敛的,而且,要使|RSn|,只要,科

9、大研究生学位课程,例 已知函数,分别用复化梯形公式、复化Simpson公式计算积分,解,的数据表,的近似值。,I精确到数点后7位的值是0.9460831。,科大研究生学位课程,例 利用复化梯形公式和复化Simpson公式分别计算上例中定积分,若使精度=10-6,问各需取n为多少?,解 因为(x)=,所以有,于是有,对复化梯形公式,若使|RTn|10-6,只要,故应取n=167.,对复化Simpson公式,若使|RSn|10-6,只要,故只需取n=3.,实际上,S3=0.9460838.,科大研究生学位课程,变步长求积方法,实际的积分计算问题，很难根据误差|Rf|,确定步长 h。,利用两种步长计

10、算积分时，为了减少计算函数f(x)的次数，通常取 h*=h/2.,即先任取步长h 进行计算，然后取较小步长h*进行计算，如果两次计算结果相差较大，则取更小步长进行计算，如此下去，直到相邻两次计算结果相差不大为止，取最小步长算出的结果作为积分值。这种方法称为变步长积分法。,不过由误差表达式可见,只要公式中涉及的高阶导数有界，当h-0时,总有Rf-0.这说明,只需 h 充分小,必可满足误差要求.,因此为计算积分，通常采取逐步缩小步长的办法。,例如应用复化梯形求积公式时，注意当前步长为h时,有,科大研究生学位课程,可见步长减半时,这表明算出T(h)后，为算T(h/2)，只需计算新增节点xi-1/2=

11、a+(i-1/2)h(i=1,n)处的函数值f(xi-1/2)，将它们的和乘新步长h/2，再加上T(h)的一半。,利用T(h)和T(h*)还可近似误差估计,称之事后误差估计.,科大研究生学位课程,对于复化梯形公式,n等分区间h=(b-a)/n,2n等分区间,近似有：,由此引入龙贝格求积方法。,3,事后误差估计公式,科大研究生学位课程,由此得,记T(h)=Tn,T(h/2)=T2n,一方面，若|T2n-Tn|3,则有近似误差|I*-T2n|.,5.4 Romberg求积公式,所以有,另一方面，(4T2n-Tn)/3应比Tn和T2n的近似程度更好.事实上,有,其中,xk=a+kh,k=0,1,2,

12、n,h=(b-a)/n,科大研究生学位课程,而且有,于是有,因此有逐次分半的复化梯形公式的递推公式：,而且，要使,=Sn,只要,科大研究生学位课程,复化Simpson公式能加工成更高精度的公式吗?,由复化Simpson公式的误差估计式有：,科大研究生学位课程,所以有,由此得,一方面，若|S2n-Sn|15,则有近似误差|I*-S2n|.,另一方面,(16S2n-Sn)/15应比Sn和S2n的近似程度更好.,(16S2n-Sn)/15=Cn,类似地,由于,事实上,有,科大研究生学位课程,所以有,由此得,一方面，若|C2n-Cn|63,则有近似误差|I*-C2n|.,另一方面,(64C2n-Cn)

13、/63应比Cn和C2n的近似程度更好.,记(64C2n-Cn)/63=Rn,称为Romberg求积公式.,科大研究生学位课程,用Tm(k)(m=1,2,3,4)分别表示把区间2k等分的复化梯形公式,复化Simpson公式,复化Cotes公式和Romberg求积公式.,而且，要使|I*-Tm(k)|,只要|Tm(k)-Tm(k-1)|(4m-1)(m=1,2,3,4).,则有,若对Romberg求积公式作组合也有,科大研究生学位课程,Romberg算法：,?,?,?,科大研究生学位课程,实际计算可按下表顺序进行,例 利用Romberg积分公式计算积分,科大研究生学位课程,解 按递推公式计算,结果

14、如下,可见,由于|T1(4)-T1(3)|=0.0019531,应有|I*-T1(4)|0.000651033.,由于|T2(3)-T2(2)|=0.0000001,应有|I*-T2(3)|0.00000000667.,由于|T3(2)-T3(1)|=0.0000015,应有|I*-T3(2)|0.00000002381.,由于|T4(1)-T4(0)|=0.0000068,应有|I*-T4(1)|0.00000026667.,科大研究生学位课程,理查德森外推法,利用低阶公式产生高精度的结果。,设对于某一 h 0，有公式 T0(h)近似计算某一未知值 I。由Taylor展开得到：T0(h)I=

15、1 h+2 h2+3 h3+,i 与 h 无关,Q：如何将公式精度由 O(h)提高到 O(h2)？,即：,科大研究生学位课程,前面介绍的 n+1个节点的 Newton-Cotes求积公式，其特征是节点是等距的。这种特点使得求积公式便于构造，复化求积公式易于形成。但同时也限制了公式的精度。n是偶数时，代数精度为n+1，n是奇数时，代数精度为n。,我们知道 n+1个节点的插值型求积公式的代数精确度不低于n。设想：能不能在区间a,b上适当选择n+1个节点 x0,x2,xn,使插值求积公式的代数精度高于n？,答案是肯定的，适当选择节点，可使公式的精度最高达到2n+1，这就是本节所要介绍的高斯求积公式。

16、,5.5 高斯(Gauss)求积公式,科大研究生学位课程,构造具有2n+1次代数精度的求积公式,将节点 x0 xn 以及系数 A0 An 都作为待定系数。令 f(x)=1,x,x2,x2n+1 代入可求解，得到的公式具有2n+1 次代数精度。这样的节点称为Gauss 点，公式称为Gauss 型求积公式。,例：求 的 2 点 Gauss 公式。,代入 f(x)=1,x,x2,x3,不是线性方程组，不易求解。,科大研究生学位课程,这样,高斯求积公式是,由于非线性方程组较复杂,通常对于n=2时就很难求解.,故一般不通过解非线性方程求，,而从分析高斯点的特性来构造高斯求积公式.,此方法称为用待定系数法

17、构造高斯求积公式.,利用正交多项式构造高斯求积公式.,科大研究生学位课程,定理,是高斯点的充分必要条件是以这些节点为零点的多项式,与任何次数不超过 的多项式 带权 正交，,即,求积公式(*)的节点,为了构造高斯求积公式的节点,有下述结论.,Gsuss点应是n+1次正交多项式的零点。,科大研究生学位课程,(2)求出pn+1(x)的n个零点x0,x1,xn 即为Gsuss点.,(1)求出区间a,b上权函数为(x)的正交多项式pn+1(x).,(3)再计算积分系数,Gauss型求积公式的构造方法,一是采用施密特正交化方法.,如何求正交多项式pn+1(x).,另是借用现成的正交多项式函数组.,科大研究

18、生学位课程,解 按 Schemite 正交化过程作出正交多项式:,例：,=x,科大研究生学位课程,故两点Gauss公式为,积分系数为,P2(x)的两个零点为,一些现成的正交多项式组有,1.Legendre多项式,科大研究生学位课程,2.Chebyshev多项式,Tn(x)=cos(narccosx)x-1,1,n=0,1,2,是区间-1,1上权函数(x)=,正交多项式。,3.Laguere多项式,是区间0,+)上权函数(x)=e-x 的正交多项式。,科大研究生学位课程,区间-1,1上权函数(x)=1的Gauss型求积公式,称为Gauss-Legendre求积公式,其Gauss点为Legendr

19、e多项式的零点.,几种Gauss型求积公式,(1)Gauss-Legendre求积公式,公式的Gauss点和求积系数可在数学用表中查到.,科大研究生学位课程,例 用3点Gauss公式计算积分,解 查表得x0=-0.7745966692,x1=0,x2=0.7745966692,A0=A2=0.5555555556,A1=0.8888888889,所以有,Gauss-Legendre求积公式的余项为,误差为,实际上,I*=2sin1=1.68294197,误差为|R|=6.15810-5.,用Simpson公式,则有I*1.69353487,误差为|R|=1.0610-2.,科大研究生学位课程,

20、由于,因此,a,b上权函数(x)=1的Gauss型求积公式为,例 用3点Gauss公式计算积分,结果远比Simpson公式的结果精确.,解 这里Gauss点和积分系数与上例相同,所以,求积误差可表示为,科大研究生学位课程,区间0,)上权函数(x)=e-x的Gauss型求积公式,称为Gauss-Laguerre求积公式,其Gauss点为Laguerre多项式的零点.,(2)Gauss-Laguerre求积公式,公式的Gauss点和求积系数可在数学用表中查到.,科大研究生学位课程,Gauss-Laguerre求积公式为,求积公式的误差为,由于,所以,对0,+)上权函数(x)=1的积分,也可以构造类

21、似的Gauss-Laguerre求积公式:,科大研究生学位课程,数值微分就是用离散方法近似地求出函数在某点的导数值.,5.6 数值微分,后一种数值微分方法称为中点方法，它其实是前两种方法的算术平均.,但它的误差阶却由 提高到,按照Taylor展开原理可得,其中h步长。,5.6.1.差商法,科大研究生学位课程,设Ln(x)是(x)以a=x0 x1xn=b为节点的n次Lagrange插值多项式,则取,当(x)Cn+1+ka,b时,有,5.6.2 插值型求导公式,特别,当k=1时有,如果仅限定在节点xi处求导,则有,科大研究生学位课程,如取n=1的线性插值L1(x)=(x-x0)(x1)-(x-x1

22、)(x0)/h,(其中h=x1-x0)可得数值微分的二点公式:,如取n=2的等距节点(x2-x1=x1-x0=h)抛物线插值:,L2(x)=(x-x1)(x-x2)(x0)-2(x-x0)(x-x2)(x1)+(x-x0)(x-x1)(x2)/2h2,则有,L2(x)=(2x-x1-x2)(x0)-2(2x-x0-x2)(x1)+(2x-x0-x1)(x2)/2h2,L2(x)=2(x0)-4(x1)+2(x2)/2h2,L2(x)=(x0)-2(x1)+(x2)/h2,科大研究生学位课程,可得数值微分的三点公式:,科大研究生学位课程,由向后差商公式有,由三点公式有,精确值。,例 已知函数f(x)=ex的数据位：f(2.6)=13.4637,f(2.7)=14.8797f(2.8)=16.4446.用二点、三点公式计算f(x)在x=2.7处的一阶、二阶导数的近似值。,解：由向前差商公式得,由中心差商公式有,科大研究生学位课程,练 习 题,第116页 习题55.1-5.4,5.13-5.14,