搜索
东南大学计算方法与实习实验报告学院：电子科学与工程学院学号：06A12528姓名：陈毓锋同组人员：罗关生，肖洋指导老师：李元庆1、Lagrange插值多项式的表达式为：𝑃(𝑥)=∑𝑓(𝑥𝑗)∏𝑥−𝑥𝑗𝑥𝑗−𝑥𝑖𝑖≠𝑗𝑛𝑗=1=∑𝑓(𝑥𝑗)𝑙𝑗(𝑥)𝑛𝑗=1但上述公式不适合实际计算：-对x的求值，每次需要的工作量是O（n2）的；-增加（或删减）和修改节点，都需要重新计算表达式。但这只是通常的论断，只要稍作修改，它一样可以适合计算。解析：1）对Lagrange插值多项式稍作修改：运用多项式𝑙(𝑥)=∏(𝑥−𝑥𝑖)𝑛𝑖=1可以将拉格朗日基本多项式重新写为：𝑃(𝑥)=𝑙(𝑥)𝑥−𝑥𝑗1∏（𝑥𝑗−𝑥𝑖）𝑛𝑖=0,𝑖≠𝑗令ω𝑗=1∏（𝑥𝑗−𝑥𝑖）𝑛𝑖=0,𝑖≠𝑗，则𝑃(𝑥)=𝑙(𝑥)∑𝜔𝑖(𝑥−𝑥𝑖)𝑛𝑖=1则它的优点是当插值点的个数增加一个时，每个ωj都除以（xj-xn+1），就可以得到新的ωn+1，则计算的工作量O（n），比重新计算每个多项式所需要的复杂度O（n2）减少了一个量级。算法：①输入xi,yi(i=0,1,2,…,n)，令Ln(x)=0;②对i=1,2,…,n计算𝑙(𝑥)=∏(𝑥−𝑥𝑖)𝑛𝑖=1ω𝑗=1∏（𝑥𝑗−𝑥𝑖）𝑛𝑖=0,𝑖≠𝑗𝐿𝑛(𝑥)=𝐿𝑛(𝑥)+𝑙(𝑥)𝜔𝑖(𝑥−𝑥𝑖)算法程序代码：#includeiostream#includecmath#includectime#includeiomanip#includecstdlibusingnamespacestd;floatBarycentricLagrange(floatx[],floaty[],floatxx,intn){inti,j;floatlx=1,*w,*L,W=0;w=newfloat[n];L=newfloat[n];floatP=0;for(i=0;in;i++){//多项式l(x)=(x-x0)(x-x1)…(x-xn)lx*=xx-x[i];}for(i=0;in;i++){//将1/(x-xj)看做一个整体L()L[i]=1/(xx-x[i]);}for(i=0;in;i++){w[i]=L[i]*y[i];for(j=0;jn;j++){if(j!=i){w[i]*=1/(x[i]-x[j]);}}W+=w[i];}P=lx*W;returnP;}对比拉格朗日基本表达式算法情况：#includeiostream#includecmath#includectime#includeiomanip#includecstdlibusingnamespacestd;floatLagrange(floatx[],floaty[],floatxx,intn){inti,j;float*l,L=0;l=newfloat[n];for(i=0;in;i++){l[i]=y[i];for(j=0;jn;j++){if(j!=i){l[i]*=(xx-x[j])/(x[i]-x[j]);}}L+=l[i];}deletel;returnL;}2）当选择均匀节点时，取xi=x0+ih,可以得到ω𝑗=(−1)𝑛−1ℎ𝑛−1(n−1)!(𝑛−1𝑖)经过验证，即ω𝑗=(−1)𝑛−1ℎ𝑛−1(n−1)!𝐶𝑛−1𝑗此时该均匀节点重新定义的ωj与之前定义的ωj值相同，验证为正确定义。算法：③输入xi,yi(i=0,1,2,…,n)，令Ln(x)=0;④对i=1,2,…,n计算𝑙(𝑥)=∏(𝑥−𝑥𝑖)𝑛𝑖=1ω𝑗=(−1)𝑛−1ℎ𝑛−1(n−1)!𝐶𝑛−1𝑗𝐿𝑛(𝑥)=𝐿𝑛(𝑥)+𝑙(𝑥)𝜔𝑖(𝑥−𝑥𝑖)算法程序代码：#includeiostream#includecmath#includectime#includeiomanip#includecstdlibusingnamespacestd;floatCorrctAverageCodeBarycentricLagrange(floatx0,floaty[],floatxx,intn,floath){inti,k=0;floatlx=1,jc=1,xs=1,*x,*w,*L,*XS,W=0;x=newfloat[n];w=newfloat[n];L=newfloat[n];XS=newfloat[n];floatP=0;for(i=0;in;i++){//给x[i]赋值x[i]=x0+i*h;}for(i=0;in;i++){//多项式l(x)=(x-x0)(x-x1)…(x-xn)lx*=xx-x[i];}for(i=0;in-1;i++){//计算(n-1)!jc*=n-(i+1);}for(i=0;in;i++){//将1/(x-xj)看做一个整体L()L[i]=1/(xx-x[i]);}for(i=1;in;i++){//(n-1)(n-2)…(n-i)/ixs*=(float)(n-i)/i;XS[i]=xs;}h=pow(h,n-1);for(i=0;in;i++){w[i]=L[i]*y[i];if((n-i)%2==0){if(i==0){w[i]*=(-1)/(h*jc);}else{w[i]*=(-1)*XS[i]/(h*jc);}}elseif((n-i)%2!=0){if(i==0){w[i]*=1/(h*jc);}else{w[i]*=XS[i]/(h*jc);}}W+=w[i];}P=lx*W;returnP;}运行结果代码：voidmain(){inta=0,b=0;clock_tstart1,start2,end1,end2;doubleduration1,duration2;floatx1[5]={0.4,0.55,0.65,0.8,0.9};//非均匀节点测试floaty1[5]={0.41075,0.57815,0.69675,0.88811,1.02652};floatxx1=0.596,L,BL1;floatx2[4]={0.3,0.4,0.5,0.6};//均匀节点测试floaty2[5]={1.2222,1.2681,1.3033,1.3293};floatxx2=0.45,BL2,ACBL;start1=clock();do{a++;L=Lagrange(x1,y1,xx1,4);}while(a=100000);end1=clock();duration1=(double)(end1-start1)/CLOCKS_PER_SEC;cout所用时间为：duration1endl;start2=clock();do{b++;BL1=BarycentricLagrange(x1,y1,xx1,4);}while(b=100000);end2=clock();duration2=(double)(end2-start2)/CLOCKS_PER_SEC;cout所用时间为：duration2endl;BL2=BarycentricLagrange(x2,y2,xx2,4);ACBL=CorrctAverageCodeBarycentricLagrange(0.3,y2,xx2,4,0.1);coutx=setprecision(6)xx1,L[4]=setprecision(7)Lendl;coutx=setprecision(6)xx1,BL1[4]=setprecision(7)BL1endl;coutx=setprecision(6)xx2,BL2[4]=setprecision(7)BL2endl;coutx=setprecision(6)xx2,ACBL[4]=setprecision(7)ACBLendl;}运行结果：3）实验分析：利用for循环将拉格朗日基本插值多项式和修正后的拉格朗日插值多项式循环了10万次，通过程序运行如图可以看到分别所用时间为0.22ms和0.248ms，则每个程序所用时间为0.22*10-8s和0.248*10-8s。利用多次计算，得到多次所用时间数据，取平均值则可以发现修正后的插值多项式比基本多项式运行速度更快一点。另外由于代码段事实上调用进行运算时修正所要经过的循环比基本多，仍能跟基本持平，甚至更快，所以该方法所需要的运算量更快。在关于稳定性的取值时在初值+0.0001/rand()取微小量，可以得出：则得出该方法的稳定性良好。通过该程序对该算法的验证，可以看出修正后的拉格朗日插值多项式比拉格朗日基本插值多项式计算量较小，稳定性良好，并且均匀节点时算法也具有这样的优点。