计算方法实验报告

Matelab 实习三

一．题目

1.对于给函数f(x)=1/(1+25x^2)在区间[-1,1]上取Xi=-1+0.2i （i=0,1,…10),试求三次曲线拟合，试画出拟合曲线并打印方程，与第2章计算实习题2的结果比较。

试求3次、4次多项式的曲线拟合，再根据数据曲线形状，求一个另外函数的拟合曲线，用图示数据曲线及相应的三种拟合曲线。

3.使用快速傅里叶变换确定函数)cos()(2x x x f =在[-pi,pi]上的16次三角插值多项式。

二．实验理论 1.，由此求出格拉姆矩阵G ，可写成矩阵形式：

Ga=d 。实验中，只要求出矩阵G 和向量d ，即可由Ga=d 得出ai(i=0,1,2 …n)，进而可求出S(x)=∑ai*x^i.

2.，可写成矩阵形式： Ga=d 。实验中，只要求出矩阵G 和向量d ，即可由Ga=d 得出ai(i=0,1,2 …n)

进而可求出S(x)=∑ai*x^i.

3.FFT 计算公式如下：

其中q=1,2,…p,k=0,1,…2^(p-q)-1,j=0,1,…2^(q-1)-1。经过不断的二分迭代，可得到Ap(j).

即Ck ，而ak=real(Ck),bk=image(Ck),

∑-=+++=11)sin cos (2

1n k n k k a ky b ky a S

三．实验程序

1.

function [ output_args ] = triRunge(m)

n=m+1;

x=zeros(1,11);

y=zeros(1,11);

x(1)=-1;

for i =1:10

x(i+1)=-1+0.2*i;

end ```·,1,0,]),(==∑k d a k j j k ?? ```·,1,0,]),(==∑k d a k j j k ?? ```·,1,0,]),(==∑k d a k j j k ?? ```·,1,0,]),(==∑k d a k j j k ?? []12111111-q 11111)22(-)2()22(A )

22()2()2(A -----------+++=++++++=+q k p q q q q q q p q q q q q q w j k A j k A j k j k A j k A j k

for i=1:11

y(i)=1/(1+25*x(i)^2);

end

d=zeros(n,n); f=zeros(1,n); for i=1:n

for j=1:n

for k=1:11

d(i,j)=d(i,j)+x(k)^(i+j-2); end

end

end

for i=1:n

for k=1:11

f(i)=f(i)+x(k)^(i-1)*y(k);

end

end

a=f/d; syms x;

g=a(n)*x+a(n-1);

for i=1:n-2

g=g*x+a(n-2-i+1);

end

vpa(g);

g=collect(g);

g=vpa(g,3); ezplot(g,[-1,1]);

gtext('拟合曲线');

k=char(g);

hold on

ezplot('1/(1+25*x^2)',[-1,1]);

gtext('龙格函数1/(1+25x^2)');

gtext(k)

title('三次曲线拟合');

end

注：蓝色为3次拟合曲线，红色为龙格函数曲线

可以看出：三次拟合曲线与龙格函数基本不重合，拟合效果差。

多次实验可以得到结果：

10次及10次以上的曲线拟合，拟合曲线与龙哥函数曲线拟合度较高：

注：红色为龙格函数曲线，蓝色为10次曲线拟合（拟合曲线表达式过长，未能完全显示）。与第2章实习题2的结果进行比较：

3次曲线拟合效果较差。

2.

function [ output_args ] = tri(m) n=m+1;

x=[0 0.1 0.2 0.3 0.5 0.8 1.0];

y=[0.1 0.41 0.5 0.61 0.91 2.02 2.46]; plot(x,y,'g');

d=zeros(n,n);

f=zeros(1,n);

for i=1:n

for j=1:n

for k=1:7

d(i,j)=d(i,j)+x(k)^(i+j-2); end

end

end

for i=1:n

for k=1:7

f(i)=f(i)+x(k)^(i-1)*y(k);

end

end

a=f/d; syms x;

g=a(n)*x+a(n-1);

for i=1:n-2

g=g*x+a(n-2-i+1);

end

vpa(g);

g=collect(g);

g=vpa(g,3);

hold on;

ezplot(g,[0,1]);

gtext('4次拟合曲线');

k=char(g);

gtext(k)

end

注：红色为3次拟合曲线，蓝色为4次拟合曲线，绿色为原函数曲线。分析结果图，可以看出，4次曲线拟合与原函数重合较好。

3.function [ ] = t(m)

n=2*m;

a=-pi;

b=pi;

f=x^2*cos(x);

h=(b-a)/n;

a1=zeros(1,n);

a2=zeros(1,n);

w=zeros(1,m);

syms x;

for d=0:n-1

x=-pi+h*d;

a1(d+1)=eval(f);

end

for k=1:m

w(k)=(exp(-i*2*pi/n))^(k-1);

end

p=log2(n)

for q=1:p

if rem(q,2)

for k=0:2^(p-q)-1

for j=1:2^(q-1)

a2(k*2^q+j)=a1(k*2^(q-1)+j)...

+a1(k*2^(q-1)+j+2^(p-1));

a2(k*2^q+j+2^(q-1))=(a1(k*2^(q-1)+j)...

-a1(k*2^(q-1)+j+2^(p-1)))*w(k*2^(q-1)+1);

end

end

else

for k=0:2^(p-q)-1

for j=1:2^(q-1)

a1(k*2^q+j)=a2(k*2^(q-1)+j)...

+a2(k*2^(q-1)+j+2^(p-1));

a1(k*2^q+j+2^(q-1))=(a2(k*2^(q-1)+j)...

-a2(k*2^(q-1)+j+2^(p-1)))*w(k*2^(q-1)+1);

end

end

end

end

if rem(p,2)

c=a2;

else

c=a1;

end

syms x;

ak=zeros(n); bk=zeros(n);

for k=1:n

ak(k)=real(c(k)*exp(-i*(k-1)*(-pi)))/m;

bk(k)=imag(c(k)*exp(-i*(k-1)*(-pi)))/m;

end

s=ak(1)/2;

for k=1:m

s=s+ak(k+1)*cos(k*x)+bk(k+1)*sin(k*x);

s=vpa(s,3);

end

ezplot(s,[-pi,pi]);

gtext('16次三角插值多项式');

hold on

ezplot(f,[-pi,pi]);

gtext('f(x)=x^2*cos(x)');

end