数值逼近实验报告2
- 格式:docx
- 大小:285.82 KB
- 文档页数:11
实验报告
实验项目名称 函数逼近与快速傅里叶变换
实 验 室 数学实验室
所属课程名称 数值逼近
实 验 类 型 算法设计
实 验 日 期 2013年10月16日
班 级 11信息与计算科学
学 号 2011119404
姓 名 冯学宁
成 绩
实验概述:
【实验目的及要求】
本次实验的目的是熟练《数值分析》第三章“函数逼近与快速傅里叶变换”的相关内容,掌握切比雪夫多项式、勒让德多项式、n次曲线拟合以及快速傅里叶变换。
本次试验要求编写牛顿多项式插值,三次样条插值,拉格朗日插值的程序编码以及画图,并在MATLAB软件中去实现。
【实验原理】
《数值分析》第三章“函数逼近与快速傅里叶变换”的相关内容,包括:切比雪夫多项式、勒让德多项式、n次曲线拟合以及快速傅里叶变换的相应算法和相关性质。
【实验环境】(使用的软硬件)
软件:
MATLAB 2012a
硬件:
电脑型号:联想 Lenovo 昭阳E46A笔记本电脑
操作系统:Windows 8 专业版
处理器:Intel(R)Core(TM)i3 CPU M 350 @2.27GHz 2.27GHz
实验内容:
【实验方案设计】
第一步,将书上关于切比雪夫多项式、勒让德多项式、n次曲线拟合以及快速傅里叶变换的内容转化成程序语言,用MATLAB实现;第二步,分别用切比雪夫多项式、勒让德多项式、n次曲线拟合以及快速傅里叶变换求解不同的问题。
【实验过程】(实验步骤、记录、数据、分析)
实验的主要步骤是:首先分析问题,根据分析设计MATLAB程序,利用程序算出问题答案,分析所得答案结果,再得出最后结论。
实验一:编写程序实现[-1,1]上n阶切比雪夫多项式,并作画(n=0,1,…,10
在一个figure中)。要求:输入Chebyshev(-1,1,n),输出如anxn+an-1xn-1+…多项式。
在MATLAB的Editor中建立一个M-文件,输入程序代码,实现切比雪夫多项式的程序代码如下:
function Pn=Chebyshev(n,x)
syms x;
if n==0
Pn=1;
else if n==1
Pn=x;
else Pn=expand(2*x*Chebyshev(n-1)-Chebyshev(n-2));
end
end
x=[-1:0.01:1];
A=sym2poly(Pn);
yn=polyval(A,x);
plot (x,yn);
hold on
end
end
在command Windows中输入命令:Chebyshev(10),得出的结果为:
Chebyshev(10) ans =
512*x^10 - 1280*x^8 + 1120*x^6 - 400*x^4 + 50*x^2 - 1
并得到Figure,图像如下:
实验二:编写程序实现[-1,1]上n阶勒让德多项式,并作画(n=0,1,…,10 在一个figure中)。要求:输入Legendre(-1,1,n),输出如anxn+an-1xn-1+…多项式。
在MATLAB的Editor中建立一个M-文件,输入程序代码,实现勒让德多项式的程序代码如下:
function Pn=Legendre(n,x)
syms x;
if n==0
Pn=1;
else if n==1
Pn=x;
else Pn=expand((2*n-1)*x*Legendre(n-1)-(n-1)*Legendre(n-2))/(n);
end
x=[-1:0.1:1];
A=sym2poly(Pn);
yn=polyval(A,x);
plot (x,yn,'-o');
hold on end
end
在command Windows中输入命令:Legendre(10),得出的结果为:
Legendre(10)
ans =
(46189*x^10)/256 - (109395*x^8)/256 + (45045*x^6)/128 - (15015*x^4)/128 +
(3465*x^2)/256 - 63/256
并得到Figure,图像如下:
实验三:利用切比雪夫零点做拉格朗日插值,并与以前拉格朗日插值结果比较。
在MATLAB的Editor中建立一个M-文件,输入程序代码,实现拉格朗日插值多项式的程序代码如下:
function [C,D]=lagr1(X,Y)
n=length(X);
D=zeros(n,n);
D(:,1)=Y';
for j=2:n
for k=j:n
D(k,j)=(D(k,j-1)- D(k-1,j-1))/(X(k)-X(k-j+1));
end
end
C=D(n,n);
for k=(n-1):-1:1
C=conv(C,poly(X(k)));
m=length(C);
C(m)= C(m)+D(k,k);
end
在command Windows中输入如下命令:
clear,clf,hold on;
k=0:10;
X=cos(((21-2*k)*pi)./22); %这是切比雪夫的零点
Y=1./(1+25*X.^2);
[C,D]=lagr1(X,Y);
x=-1:0.01:1;
y=polyval(C,x);
plot(x,y,X,Y,'.');
grid on;
xp=-1:0.01:1;
z=1./(1+25*xp.^2);
plot(xp,z,'r')
得到Figure,图像如下所示:
比较后发现,使用切比雪夫零点做拉格朗日插值不会发生龙格现象。
实验四:对于给定函数22511)(xxf在区间1][-1,上取)10,,1,0(2.01iixi,试求3次曲线拟合,试画出拟合曲线并打印出方程,与第2章计算实习题2的结果比较。 在MATLAB的Editor中输入程序代码,实现3次曲线拟合的程序代码如下:
x=[-1:0.2:1];
y=1./(1+25*x.^2);
p=polyfit(x,y,3);
yp=polyval(p,x);
plot(x,y,'o-',x,yp,'*-')
fp=poly2sym(p)
command Windows中输出:
fp =
(8299760523663595*x)/649037107316853453566312041152512 - (3305*x^2)/5746
- (572862092712169*x^3)/81129638414606681695789005144064 +
4360609662300613/9007199254740992
并得到Figure,图像如下:
第2章计算实习题2的结果如下:
牛顿插值多项式: 三次样条多项式:
比较下来可知:三次样条的拟合效果最好,其次是牛顿插值多项式,最次是3次曲线拟合。
实验五:由实验给出数据表
x 0.0 0.1 0.2 0.3 0.5 0.8 1.0 y
1.0 0.41 0.50 0.61 0.91 2.02 2.46
试求3次、4次多项式的曲线拟合,再根据数据曲线形状,求一个另外函数的拟合曲线,用图示数据曲线以及相应的三种拟合曲线。
在MATLAB的Editor中建立一个M-文件,输入程序代码,实现3次曲线拟合的程序代码如下:
x=[0 0.1 0.2 0.3 0.5 0.8 1];
y=[1 0.41 0.5 0.61 0.91 2.02 2.46];
p1=polyfit(x,y,3)
p2=polyfit(x,y,4)
y1=polyval(p1,x);
y2=polyval(p2,x);
plot(x,y,'g-',x,y1,'r-',x,y2,'b-')
hold on
p3=polyfit(x,y,2)%观察图形与抛物线接近,故采用2次曲线拟合
y3=polyval(p3,x);
plot(x,y3,'c-')
得到Figure,图像如下:
实验六:使用快速傅里叶变换确定函数f(x)=x^2*cosx在[-π,π]上的16次三角插值多项式。 在MATLAB的Editor中建立一个M-文件,输入程序代码,实现3次曲线拟合的程序代码如下:
%三角插值多项式程序
function [an,bn,f]=fseries(fx,x,n,a,b)
if nargin==3
a=-pi;
b=pi;
end
l=(b-a)/2;
if a+b
fx=subs(fx,x,x+l+a);
end
an=int(fx,x,-l,l)/l;
bn=[];
f=an/2;
for ii=1:n
ann=int(fx*cos(ii*pi*x/l),x,-l,l)/l;
bnn=int(fx*sin(ii*pi*x/l),x,-l,l)/l;
an=[an,ann];
bn=[bn,bnn];
f=f+ann*cos(ii*pi*x/l)+bnn*sin(ii*pi*x/l);
end
if a+b
f=subs(f,x,x-l-a);
end
在command Windows中输入如下命令:
>> syms x;
fx=x^2*cos(x)
[an,bn,f]=fseries(fx,x,16,-pi,pi)
得到结果如下:
an =
[ -4, pi^2/3 + 1/2, -20/9, 5/8, -68/225, 13/72, -148/1225, 25/288, -260/3969,
41/800, -404/9801, 61/1800, -580/20449, 85/3528, -788/38025, 113/6272,
-1028/65025]
bn =
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]