实验报告.docx
《实验报告.docx》由会员分享,可在线阅读,更多相关《实验报告.docx(13页珍藏版)》请在冰豆网上搜索。
实验报告
2011——2012学年第一学期
合肥学院数理系
实验报告
课程名称:
数值分析
实验项目:
插值和数据拟合实验
实验类别:
综合性□设计性□验证性□
专业班级:
09级数学与应用数学
(2)班
姓名:
徐国华学号:
06
实验地点:
7#604
实验时间:
2011年11月17日
指导教师:
孙梅兰成绩:
1.实验目的:
1、熟悉利用拉格朗日多项式求函数的近似值
2、学会编写用newton前差值公式计算函数值的程序
3、学会利用最小二乘法对数据进行拟合
2、实验内容:
实验一:
1、
1
用四次拉格朗日多项式求
的函数近似值
2、函数y=f(x)由下列数据给出:
k
0
1
2
3
4
5
6
7
8
试用k=3、4、5、6这四组构造三次多项式,求出y=f(x)在x=处的近似值,并用matlab中的polyfit命令验证结果
3、编写一个用newton前差值公式计算函数值的程序,要求先输出差分表,再计算x的函数值,并应用下面的问题:
20
21
22
23
24
求x=时的三次差值多项式的值。
实验二:
1、由化学实验得到某物质浓度与时间的关系如下:
时间t
1
2
3
4
5
6
7
8
浓度y
时间t
9
10
11
12
13
14
15
16
浓度y
求浓度与时间的二次拟合曲线。
2、试求下表中y与t的拟合曲线
t
0
5
10
15
20
25
30
35
40
45
50
55
Y/(*10^(-4))
0
要求:
(1)用plot画出原始数据分布趋势图
(2)用最小二乘法拟合,近似表达式为y=a2x^3+a1x^2+a0x
(3)打印出y
(4)另取方程y=ax^b拟合,并比较
3、用形如a*e^x+bsin(x)+cln(x)+dcos(x)的函数在最小二乘法的意义下拟合数据表
x
y
4、
使用次数x
2
3
5
6
7
9
10
容积y
使用次数x
11
12
14
16
17
19
20
容积y
选用双曲线1/y=a+b/x对数据进行拟合,使用最小二乘法求出拟合函数,做出曲线图
三、实验方案(程序设计说明)
实验一:
1、建立的lagrange1文件:
functionyy=lagrange1(x,y,xi)
m=length(x);n=length(y);
ifm~=n,error('向量x与y的长度必须一致');end
s=0;
fori=1:
n
z=ones(1,length(xi));
forj=1:
n
ifj~=i
z=z.*(xi-x(j))/(x(i)-x(j));
end
end
s=s+z*y(i);
end
yy=s;
在命令窗口调用Lagrange函数如下
>>x=[,,,,1];
>>y=[,,,,];
>>xi=[];
>>yi=lagrange1(x,y,xi)
yi=
>>plot(x,y,'o',xi,yi,'g^')
2、建立的lagrange1文件:
functionyy=lagrange1(x,y,xi)
m=length(x);n=length(y);
ifm~=n,error('向量x与y的长度必须一致');end
s=0;
fori=1:
n
z=ones(1,length(xi));
forj=1:
n
ifj~=i
z=z.*(xi-x(j))/(x(i)-x(j));
end
end
s=s+z*y(i);
end
yy=s;
在命令窗口调用Lagrange函数如下
x=[,,,];
>>y=[,,,];
>>xi=[];
>>yi=lagrange1(x,y,xi)
yi=
>>plot(x,y,'o',xi,yi,'g^')
3、
functionf=Newtonforward(x,y,x0)
symst;
if(length(x)==length(y))
n=length(x);
c(1:
n)=;
else
disp('x和y的维数不相等!
');
return;
end
f=y
(1);
y1=0;
xx=linspace(x
(1),x(n),(x
(2)-x
(1)));
if(xx~=x)
disp('节点之间不是等距的!
');
return;
end
for(i=1:
n-1)
for(j=1:
n-i)
y1(j)=y(j+1)-y(j);
end
c(i)=y1
(1);
l=t;
for(k=1:
i-1)
l=l*(t-k);
end;
f=f+c(i)*l/factorial(i);
simplify(f);
y=y1;
if(i==n-1)
if(nargin==3)
f=subs(f,'t',(x0-x
(1))/(x
(2)-x
(1)));
else
f=collect(f);
f=vpa(f,6);
end
end
end
l=t;
for(k=1:
i-1)
l=l*(t+k);
end;
f=f+c(i)*l/factorial(i);
simplify(f);
y=y1;
if(i==n-1)
if(nargin==3)
f=subs(f,'t',(x0-x(n))/(x
(2)-x
(1)));
else
f=collect(f);
f=vpa(f,6);
end
end
end
在命令窗口调用Lagrange函数如下
>>x=[20,21,22,23,24];
>>y=[,,,,];
>>x0=;
>>f=Newtonforward(x,y,x0)
实验二:
1、>>xi=[12345678910111213141516];
>>yi=[];
plot(xi,yi,'o')
>>A=[ones(size(xi));xi;xi.^2]'
A=
111
124
139
1416
1525
1636
1749
1864
1981
110100
111121
112144
113169
114196
115225
116256
>>a=A\yi'
a=
>>b=[]
b=
>>y=poly2str(b,'x')
y=
x^2+x+
>>f2=polyval(flipud(a),xi);
>>plot(xi,yi,'bo',xi,f2,'r-')
2、
(1)画出数据分布趋势图
>>ti=[0510152025303540455055];
>>yi=[0].*(10^(-4));
>>plot(ti,yi,'o')
(2)建立数学模型y=a3x^3+a2x^2+a1x+a0
建立超定方程组系数矩阵
>>A=[ones(size(ti));ti;ti.^2;ti.^3]'
A=
1525125
1101001000
1152253375
1204008000
12562515625
13090027000
135122542875
140160064000
145202591125
1502500125000
1553025166375
(3)求超定方程组的最小二乘解
>>a=A\yi'
a=
*
(4)求拟合曲线方程
>>b=[a
(1)a
(2)a(3)a(4)]
b=
*
y=poly2str(b,'x')
y=
x^3+x^2-x+
f3=polyval(flipud(a),ti);
(5)用方程y=ax^b拟合
t=[ones(size(ti));log(ti)];
>>aa=t'\log(yi)'
aa=
>>yy=exp(*xi.^(;
>>plot(ti,yi,'bo',ti,yy,'r--',ti,f3,'b--')
legend('原始数据','指数拟合曲线','3阶拟合曲线')
title('拟合曲线比较')
3、编写命令M文件t1如下
data=[;
;
;
;
;
;
;
;
;
];
x=data(:
1);y=data(:
2);
A(:
1)=exp(x);A(:
2)=sin(x);
A(:
3)=log(x);A(:
4)=cos(x);
c=A\y
xx=:
:
;
g=c
(1)*exp(xx)+c
(2)*sin(xx)+c(3)*log(xx)+c(4)*cos(xx);
plot(x,y,'bo',xx,g,'r--');
c=
g=*exp(xx)+*sin(xx)+*log(xx)+*cos(xx);
4、
>>xi=[2356791011121416171920];
>>yi=[];
>>xj=1./xi;
yj=1./yi;
>>A=[ones(size(xj));xj]'
A=
>>a=A\yj'
a=
>>b=[a
(2)a
(1)]
b=
>>y=poly2str(b,'x')
y=
x+
>>f2=polyval(flipud(a),xj);
plot(xj,yj,'bo',xj,f2,'r-')
四、程序运行结果
实验一:
第一题图第二题图
第三题运行结果:
f=
实验二:
第一题图第二题图
2、
3、
第三题图第四题图
五、实验总结:
通过matlab的运用,我基本熟悉了插值多项式构造,并通过多次调试和实验解决了实验问题;掌握了运用matlab数据插值的思想方法;明白了小二乘法的基本原理。
学生签名:
徐国华
2011年11月17日
六.教师评语及成绩
教师签名:
年月日