南大数值分析ppt课件第六章曲线拟合与函数逼近.ppt

《南大数值分析ppt课件第六章曲线拟合与函数逼近.ppt》由会员分享，可在线阅读，更多相关《南大数值分析ppt课件第六章曲线拟合与函数逼近.ppt（34页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第六章第六章曲线拟合与函数逼近曲线拟合与函数逼近/*Approximation Theory*/仍然是已知仍然是已知x1xm;y1ym,求一个简单易求一个简单易算的近似函数算的近似函数P(x)f(x)。但是但是m很大；很大；yi 本身是测量值，不准确，即本身是测量值，不准确，即yi f(xi)这时没必要取这时没必要取P(xi)=yi,而要使而要使P(xi)yi 总体上总体上尽可能小。尽可能小。常见做法：常见做法：使使最小最小/*minimax problem*/太太复杂复杂 使使最小最小不可导，求解困难不可导，求解困难 使使最小最小/*Least-Squares method*/1最小二乘拟合

2、最小二乘拟合多项式多项式/*L-S approximating polynomials*/确定多项式确定多项式，对于一组数据，对于一组数据(xi,yi)(i=1,2,n)使得使得达到达到极小极小，这里，这里nm。naaa10 实际上是实际上是a0,a1,an的多元函数，即的多元函数，即 =+=miinininyxaxaaaaa121010.),.,(在在 的极值点应有的极值点应有kiminjijijxyxa =102=+njmikiimikjijxyxa0112记记 =mikiikmikikxycxb11,法方程组法方程组(或或正规方程组正规方程组)/*normal equations*/回归

3、系数回归系数/*regression coefficients*/1L-SApproximatingPolynomials定理定理L-S拟合多项式拟合多项式存在唯一存在唯一(n 0,b0)线性化：线性化：由由可做变换可做变换xbay lnlnbBaAxXyY=,ln,1,lnBXAY+就是个就是个线性问题线性问题将将化为化为后易解后易解A和和B),(iiYX),(iiyxHW:p.233#7，#9,#10,#11例例用用来拟合来拟合。2正交多项式与最小二乘拟合正交多项式与最小二乘拟合/*Orthogonal Polynomials&Least-Squares Approximation */已

4、知已知x1xm;y1ym,求一个简单易算的近求一个简单易算的近似函数似函数P(x)f(x)使得使得最小。最小。已知已知a,b上定义的上定义的f(x)，求一个简单易算的求一个简单易算的近似函数近似函数P(x)使得使得最小。最小。线线性性无无关关/*linearly independent*/函函数数族族 0(x),1(x),n(x),满足条件：其中任意函数的线性组合满足条件：其中任意函数的线性组合a0 0(x)+a1 1(x)+an n(x)=0对任意对任意x a,b成立成立当且仅当当且仅当a0=a1=an=0。2OrthogonalPolynomials&L-SApproximation考考虑

5、虑一一般般的的线线性性无无关关函函数数族族=0(x),1(x),n(x),，其其有有限限项项的的线线性性组组合合称称为为广广义义多多项项式式/*generalized polynomial*/.常见多项式：常见多项式：j(x)=x j对应对应代数代数多项式多项式/*algebraic polynomial*/j(x)=cosjx、j(x)=sinjx j(x),j(x)对对应应三角三角多项式多项式/*trigonometric polynomial*/j(x)=e kj x,ki kj对应对应指数指数多项式多项式/*exponential polynomial*/2OrthogonalPoly

6、nomials&L-SApproximation权函数：权函数：离散型离散型/*discrete type*/根据一系列离散点根据一系列离散点拟合时，在每一误差前拟合时，在每一误差前乘一正数乘一正数wi，即即误差函数误差函数，这个，这个wi就称就称作作权权/*weight*/，反映该点的重要程度。反映该点的重要程度。=niiiiyxPw12)(连续型连续型/*continuous type*/在在a,b上用广义多项式上用广义多项式P(x)拟合连续函数拟合连续函数f(x)时，定义时，定义权函数权函数(x)Ca,b，即误差函数即误差函数=。权函数必须。权函数必须(x)满足：非负、可积，且在满足：非

7、负、可积，且在a,b的任何子的任何子区间上区间上(x)0。2OrthogonalPolynomials&L-SApproximation广义广义L-S拟合：拟合：离散型离散型/*discrete type*/在点集在点集x1xm 上测得上测得y1ym，在一组权系数在一组权系数w1wm 下求广义多项式下求广义多项式P(x)使得使得误差函数误差函数 最小。最小。=niiiiyxPw12)(连续型连续型/*continuous type*/已知已知y(x)Ca,b以及权函数以及权函数(x)，求广义多项式求广义多项式P(x)使使得误差函数得误差函数=最小最小。dxxyxPxba2)()()(内积内积与

8、与范数范数离散型离散型连续型连续型则易证则易证(f,g)是是内积内积，而而是是范数范数。(f,g)=0表示表示f与与g带权正交带权正交。广义广义L-S问题可叙述为：求广义多项式问题可叙述为：求广义多项式P(x)使得使得最小。最小。2OrthogonalPolynomials&L-SApproximationnkyaknjjjk,.,0,),(),(0=设设则完全类似地有：则完全类似地有：)(.)()()(1100 xaxaxaxPnn +=法法方程组方程组/*normal equations*/定理定理Ba=c 存在唯一解存在唯一解 0(x),1(x),n(x)线性无线性无关。关。即：即：),

9、(),(),(00yyaabnnjiij =c证明：证明：若存在一组系数若存在一组系数 i使得使得0.1100=+nn 则等式两边分别与则等式两边分别与 0,1,n作内积，得到：作内积，得到：即：即：B =02OrthogonalPolynomials&L-SApproximation例：例：用用来拟合来拟合，w 1解：解：0(x)=1，1(x)=x，2(x)=x2Itissoooosimple!Whatcanpossiblygowrong?7623)(463|484,|1=BcondBB2OrthogonalPolynomials&L-SApproximation例：例：连续型拟合中，取连续

10、型拟合中，取则则Hilbert阵！阵！改进：改进：若能取函数族若能取函数族=0(x),1(x),n(x),，使得任意一对使得任意一对 i(x)和和 j(x)两两两两（带权）正（带权）正交交，则，则B 就化为就化为对角阵对角阵！这时直接可算出这时直接可算出ak=Well,nofreelunchanyway 正交正交多项式多项式的构造：的构造：将正交函数族中的将正交函数族中的 k 取为取为k阶阶多项式多项式，为简单起见，可取，为简单起见，可取 k 的的首项系数为首项系数为1。有递推有递推关系式：关系式：其中其中证明略证明略2OrthogonalPolynomials&L-SApproximatio

11、n例：例：用用来拟合来拟合，w 1解：解：通过正交多项式通过正交多项式 0(x),1(x),2(x)求解求解设设)()()(221100 xaxaxay +=1)(0=x 229),(),(0000=ya25),(),(00001=x25)()()(011=xxxx 537),(),(1111=ya25),(),(11112=x45),(),(00111=55)(45)()25()(2012+=xxxxxx 21),(),(2222=ya与前例结果一致。与前例结果一致。注：注：手算时也可手算时也可用待定系数法确用待定系数法确定函数族。定函数族。2OrthogonalPolynomials&L-

12、SApproximation Algorithm:Orthogonal Polynomials Approximation Toapproximateagivenfunctionbyapolynomialwitherrorboundedbyagiventolerance.Input:numberofdatam;xm;ym;weightwm;toleranceTOL;maximumdegreeofpolynomialMax_n.Output:coefficientsoftheapproximatingpolynomial.Step 1Set 0(x)1;a0=(0,y)/(0,0);P(x)=a

13、0 0(x);err=(y,y)a0(0,y);Step 2Set 1=(x 0,0)/(0,0);1(x)=(x 1)0(x);a1=(1,y)/(1,1);P(x)+=a1 1(x);err=a1(1,y);Step 3Setk=1;Step 4While(k 0distinctpoints.YouaresupposedtowriteafunctionvoidOPA(double(*f)(),doublex,doublew,intm,doubletol,FILE*outfile)toapproximatethefunctionfbyanorthogonalpolynomialusingth

14、eexactfunctionvaluesatthegivenmpointsx.Thearraywmcontainsthevaluesofaweightfunctionatthegivenpointsx.Thetotalerrormustbenolargerthantol.2OrthogonalPolynomials&L-SApproximationInputThereisnoinputfile.Instead,youmusthandinyourfunctionina*.hfile.Theruleofnamingthe*.hfileisthesameasthatofnamingthe*.cor*

15、.cppfiles.Output(representsaspace)Foreachtestcase,youaresupposedtooutputthefollowinginformation:The 1st linecontainsthe integer6 n 0 which isthe degree of thepolynomialintheformat:fprintf(outfile,%dn,n);The 2nd line contains the n+1 coefficients of the approximationpolynomialwhere.Eachofthecoefficie

16、ntistobeprintedasinCprintf:fprintf(outfile,%8.4e,coefficient);The3rdlinecontainsthetotalerrorintheformat:fprintf(outfile,error=%12.8en,err);Note:Ifthetotalerrorisstillnotsmallenoughwhenn=6,simplyoutputtheresultobtainedwhenn=6.Theoutputsoftwotestcasesmustbeseperatedbyablankline.2OrthogonalPolynomials

17、&L-SApproximationSample Judge ProgramSample Judge Program#include#include#defineMAX_m200#defineMAX_n6#include98115001_12.hdoublef1(doublex)returnsin(x);doublef2(doublex)returnexp(x);voidmain()FILE*outfile=fopen(out.txt,w);intm,i;doublexMAX_m,wMAX_m,tol;m=90;for(i=0;im;i+)xi=3.1415926535897932;xi=xi*

18、(double)(i+1)/180.0;wi=1.0;tol=0.001;OPA(f1,x,w,m,tol,outfile);m=200;for(i=0;im;i+)xi=0.01*(double)i;wi=1.0;tol=0.001;OPA(f2,x,w,m,tol,outfile);fclose(outfile);2OrthogonalPolynomials&L-SApproximationSample OutputSample Output(represents a space)represents a space)3 2.5301e 0031.0287e+000 7.2279e 002

19、 1.1287e 001error=6.33097847e 00541.0025e+0009.6180e 0016.2900e 0017.0907e 0031.1792e 001error=1.61711536e 0042函数的最佳逼近函数的最佳逼近/*Optimal Approximation*/最佳最佳平方平方逼近：即连续型逼近：即连续型L-S逼近，在逼近，在意义意义下，下，使得使得最小。最小。最佳最佳一致一致逼近逼近/*uniform approximation*/在在意义下，使得意义下，使得最小。也称为最小。也称为minimaxproblem。偏差偏差/*deviation*/若若，则

20、称，则称x0为为 偏差点偏差点。Didntyousayitsaverydifficultproblem?Takeiteasy.Itsnotsodifficultifweconsiderpolynomialsonly.3OptimalApproximationv1.0最佳一致逼近多项式最佳一致逼近多项式/*optimal uniform approximating polynomial*/的构造：求的构造：求n 阶多项式阶多项式Pn(x)使得使得|Pn y|最小。最小。直接构造直接构造OUAP 的确比较困难，不妨换个角度，先的确比较困难，不妨换个角度，先考察它应该具备的考察它应该具备的性质性质。

21、有如下结论：。有如下结论：OUAP 存在，且必同时有存在，且必同时有 偏差点。偏差点。证明：证明：存在性证明略。后者用反证法，设只有正偏差点。存在性证明略。后者用反证法，设只有正偏差点。设设而而对于所有的对于所有的x a,b都有都有是是n阶多项式阶多项式是是误差更小的多项式误差更小的多项式3OptimalApproximation（Chebyshev定理）定理）Pn是是y的的OUAP Pn关于关于y在在定义域上定义域上至少有至少有n+2个个交错交错的的 偏差点。偏差点。即存在点集即存在点集a t1tn+2 b 使得使得tk称为称为切比雪夫交错组切比雪夫交错组/*Chebyshev altern

22、ating sequence*/若若且且y不是不是n次多项式，则次多项式，则n 次次OUAP唯一唯一。证明：证明：反证，设有反证，设有2个个OUAPs，分别是分别是Pn 和和Qn 。则则它们的平均函数它们的平均函数也是一个也是一个OUAP。2)()()(xQxPxRnnn+=对于对于Rn有有Chebyshev交错组交错组t1,tn+2使得使得nkknkknkknnEtytQtytPtytRE +=|)()(|21|)()(|21|)()(|nkknkknEtytQtytP=|)()(|)()(|则至少在一个点上必须有则至少在一个点上必须有)()()()(knkkkntQtytytP=0)()(

23、=kkntytR0=nE3OptimalApproximation由由Chebyshev定理可推出：定理可推出：Pn(x)y(x)在定义域上至少变号在定义域上至少变号 次，故至少有次，故至少有个个根根。xy0yy x=()yy xEn=+()yy xEn=()yP xn=()n+1n+1可见可见Pn(x)是是y(x)的的某一个某一个插值多项式插值多项式如何确定插如何确定插值节点值节点x0,xn 的位置，使得的位置，使得Pn(x)刚好是刚好是y 的的OUAP？即，使插值余项即，使插值余项v2.0达到极小？达到极小？3OptimalApproximationv2.1在在 1,1上求上求x1,xn使

24、得使得的的|wn|最小。最小。=niinxxxw1)()(注意到注意到，要使，要使|wn|最小就意味着最小就意味着)()(1xPxxwnnn =v3.0在在 1,1上求函数上求函数xn的的n 1阶阶 OUAP。由由Chebyshev定理可推出：定理可推出：Pn 1(x)关于关于xn有有n+1个个偏差点，即偏差点，即wn(x)在在n+1个点上个点上交错取极大、极小值。交错取极大、极小值。v3.1在在 1,1上求上求切比雪夫交错组切比雪夫交错组t1,tn+1。切比雪夫多项式切比雪夫多项式/*Chebyshev polynomials*/3OptimalApproximation考虑三角函数考虑三角

25、函数cos(n )在在0,上的上的个极值点。个极值点。n+1当当时，时，cos(n )交错达到极大值交错达到极大值1和极小值和极小值 1，且存在系数，且存在系数a0,an使得使得令令x=cos()，则则x 1,1。)cosarccos()cos()(xnnxTn=称为称为Chebyshev多多项式项式Tn的重要性质：的重要性质：当当时，时，交错取到极大值交错取到极大值1和极小和极小值值 1，即，即1当当时时，即，即x1,xn为为Tn(x)的的n个零点。个零点。3OptimalApproximationTn(x)满足递推关系：满足递推关系：T0(x)=1，T1(x)=x，Tn(x)为为n次多项式

26、，首项系数为次多项式，首项系数为。且。且T2n(x)只含只含x的的次幂，次幂，T2n+1(x)只含只含x的的次幂。次幂。2n 1偶偶奇奇T0(x),T1(x),是是 1,1 上上关于权关于权正交的正交的函数族。即在内积函数族。即在内积的意义下有的意义下有OKOK,IthinkitsenoughforusWhatsourtargetagain?v3.1在在 1,1上求上求切比雪夫交错组切比雪夫交错组t1,tn+1。v3.0在在 1,1上求函数上求函数xn的的n 1阶阶 OUAP。Tn(x)的的n个个零点零点。3OptimalApproximation可见：若取可见：若取，则，则wn在在 1,1

27、上有上有n+1个个极值极值点点tk，也即，也即Pn 1(x)=xn wn(x)关于关于xn在在 1,1 上有上有n+1个交错个交错偏差点偏差点tk。v3.0OKv2.1在在 1,1上求上求x1,xn使得使得的的|wn|最小。最小。=niinxxxw1)()(取最小值取最小值 n=首项系数为首项系数为1的的n阶多项式阶多项式/*monic polynomials of degree n*/x1,xn即为即为如何确定插值节点如何确定插值节点x0,xn 的位置，使得的位置，使得Pn(x)刚好是刚好是y 的的OUAP？即，使插值余项即，使插值余项达到极小？达到极小？v2.0取取x0,xn为为Tn+1(

28、x)的的n+1个个零点零点，做，做y 的的插值多项式插值多项式Pn(x)，则，则插值余项的上界可达极小插值余项的上界可达极小。3OptimalApproximation注：注：上界上界最小不表示最小不表示|Rn(x)|最小，故最小，故Pn(x)严格意义上只是严格意义上只是y(x)的的近似近似最佳逼近多项式；最佳逼近多项式；对于一般区间对于一般区间x a,b，可作变量替换可作变量替换，则，则t 1,1，这时这时即以即以为插值节点为插值节点(k=0,n)，得，得Pn(x)，余项余项有最小上界。有最小上界。3OptimalApproximation例：例：求求f(x)=ex在在0,1上的近似最佳逼近

29、多项式，使其误差上的近似最佳逼近多项式，使其误差不超过不超过0.5 10 4。解：解：根据误差上界确定根据误差上界确定n:n=4计算计算T5(t)的根：的根：以以x0,x4为节点作为节点作L4(x)3OptimalApproximation Chebyshev多项式的其它应用多项式的其它应用多项式降次多项式降次/*reduce the degree of polynomial with a minimal loss of accuracy*/设设f(x)Pn(x)。在降低在降低Pn(x)次数的同时次数的同时,使因此增加的误差尽可能小使因此增加的误差尽可能小,也叫也叫economiza-tion

30、 of power series。从从Pn中去掉一个含有其最高次项的中去掉一个含有其最高次项的,结果降结果降次为次为,则：则：PnPn 1|)(|max|)()(|max|)()(|max1,11,111,1xPxPxfxPxfnnn +因降次而增的误差因降次而增的误差设设Pn 的的首项系数为首项系数为an，则取则取可使精度可使精度尽可能少损失。尽可能少损失。12)()(=nnnnxTaxP3OptimalApproximation例：例：f(x)=ex在在 1,1上的上的4阶阶Taylor展开为展开为，此时误差此时误差请将其请将其降为降为2阶阶多项式多项式。解：解：取取（查表知（查表知）取取（查表知（查表知）若简单取若简单取，则误差，则误差另类另类解法可查解法可查p.163表表7-2，将将x3和和x4中的中的T3和和T4删除。删除。注：注：对一般区间对一般区间a,b，先将先将x换为换为 t，考虑考虑f(t)在在 1,1上的逼近上的逼近Pn(t)，再将再将t换回换回x，最后得到最后得到Pn(x)。HW:p.164#3