数值分析实验报告-插值、逼近

  • 格式:doc
  • 大小:284.97 KB
  • 文档页数:10

下载文档原格式

  / 10
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

实验报告:函数逼近&插值多项式补充

问题1:对于给函数2

1()1+25f x x =

,取点21

cos 22k k x n π+=+,k 取0,1,…,n 。n 取10或20。试画出拟合曲线并打印出方程,与第二章计算实习题2的结果进行比较。 问题2:对于给函数2

1

()1+25f x x =

在区间[-1,1]上取x i =-1+0.2i (i=0,1,2,…,10),试

求3次曲线拟合,试画出拟合曲线并打印出方程,与第二章计算实习题2的结果进行比较。

实验目的:通过编程实现牛顿插值方法和函数逼近,加深对多项式插值的理解。应用所编程序解决实际算例。

实验要求:

1. 认真分析问题,深刻理解相关理论知识并能熟练应用; 2. 编写相关程序并进行实验; 3. 调试程序,得到最终结果; 4. 分析解释实验结果; 5. 按照要求完成实验报告。

实验原理:

详见《数值分析 第5版》第二章、第三章相关容。

实验容:

(1)问题1:

这里我们可以沿用实验报告一的代码,对其进行少量修改即可。 当n=10时,代码为:

clear all clc k=0:10; n=length(k);

x1=cos((2*k+1)/2/n*pi); y1=1./(1+25.*x1.^2); f=y1(:); for j=2:n for i=n:-1:j

f(i)=(f(i)-f(i-1))/(x1(i)-x1(i-j+1)); end end syms F x p ; F(1)=1;p(1)=y1(1); for i=2:n

F(i)=F(i-1)*(x-x1(i-1)); p(i)=f(i)*F(i);

end

syms P

P=sum(p);

P10=vpa(expand(P),5);

x0=-1:0.001:1;

y0=subs(P,x,x0);

y2=subs(1/(1+25*x^2),x,x0);

plot(x0,y0,x0,y2)

grid on

xlabel('x')

ylabel('y')

由此我们可以得到P10(x)=-46.633*x^10+3.0962e-14*x^9+130.11*x^8-7.2714e-14*x^7-133.44*x^6+7.1777e-14*x^5+61.443*x^4-1.5805e-14*x^3-12.477*x^2-1.6214e-16*x+1.0

并可以得到牛顿插值多项式在[-1,1]上的图形,并和原函数进行对比,得Fig. 1。

Fig.1 牛顿插值多项式(n=10)函数和原函数图形

当n=20,将上述代码中的“k=0:10;”改为“k=0:20;”即可。

由此我们可以得到P20(x)=6466.6*x^20+8.0207e-13*x^19-34208.0*x^18-3.5038e-12*x^17+77754.0*x^16-99300.0*x^14+3.7253e-9*x^13+78236.0*x^12-

39333.0*x^10+12636.0*x^8-4.6566e-10*x^7-2537.3*x^6+306.63*x^4-21.762*x^2+1.0 并可以得到牛顿插值多项式在[-1,1]上的图形,并和原函数进行对比,得Fig. 2。

Fig.2牛顿插值多项式(

n=20)函数和原函数图形

回顾一下实验一的结果(见Fig. 3),我们不难发现,仅仅是改变了x的取值,结果发生了很大的变化。实验一中,插值多项式与原函数产生了很大的偏差,并且随着分的段数的增加,其误差不断变大,但是在本次实验中,我们不难发现,虽然多项式依旧存在震荡现象,但是误差小了很多,而且随着分的段数的增加,插值多项式曲线与原函数曲线已经十分接近了。

Fig.3实验一结果

这个例子说明:采用切比雪夫节点替代等距节点可以消除龙格现象。

(2)问题2:

分析问题,发现在这个问题中,我们已经知道了原函数,同时它也告诉我们所需取的11个点的值,所以这里可以用两种方法进行函数逼近得到拟合曲线。

首先采用最小二乘法来考虑这个问题,编写代码如下(这里没有直接调用polyfit函数):

clear all

clc

n=3;

x1=-1:0.2:1;

y1=1./(1+25.*x1.^2);

syms S G d a x;

for i=1:n+1;

for j=1:n+1;

G(i,j)=sum(x1.^(i+j-2));

end

end

for i=1:n+1;

d(i)=sum(x1.^(i-1).*y1);

end

a=G^-1*d';

for i=1:n+1

X(i)=x^(i-1);

end

S=vpa(X*a,5)

x0=-1:0.001:1;

y0=subs(S,x,x0);

y2=subs(1/(1+25*x^2),x,x0);

plot(x0,y0,x0,y2)

grid on

xlabel('x')

ylabel('y')

我们可以得到一个三次多项式:S3=1.1665e-16*x^3 - 0.57518*x^2 - 9.4553e-17*x + 0.48412。同时我们也可以得到它与原函数的图形,如图Fig. 4。

Fig.4 最小二乘法n=3的结果

我们发现得到的结果和原函数产生了巨大的误差。首先观察得到的多项式,我们发现它的3次项系数非常小,原因是原函数是一个偶函数,这将导致奇次项系数基本为0。

这里我们调整n,对结果进行观察,取n=4,6,8,10,20。我们可以得到Fig.5-Fig.9。

S4=1.4852*x^4+1.3703e-16*x^3-2.0604*x^2-1.1769e-16*x+0.65522

S6=-4.633*x^6+4.0789e-14*x^5+8.4769*x^4-5.28e-14*x^3-4.5969*x^2+1.3229e-

14*x+0.78461

S8=20.466*x^8-3.8972e-12*x^7-43.601*x^6+6.9014e-12*x^5+30.817*x^4-3.4363e-12*x^3-8.5318*x^2+4.2796e-13*x+0.88802

S10=-220.94*x^10-5.1978e-9*x^9+494.91*x^8+1.0649e-8*x^7-381.43*x^6-6.9693e-9*x^5+123.36*x^4+1.6139e-9*x^3-16.855*x^2-9.6021e-11*x+1.0

S20=-318.82*x^20+74.132*x^19+43.205*x^18-

83.871*x^17+91.867*x^16+68.562*x^15+29.364*x^14-56.393*x^13+260.42*x^12-

32.957*x^11+79.822*x^10+1.8279*x^9-139.85*x^8+49.564*x^7-121.95*x^6-

23.918*x^5+90.922*x^4+3.1437*x^3-15.933*x^2-0.090653*x+1.0

Fig.5 最小二乘法n=4的结果