数值分析实验二讲解Word下载.docx
《数值分析实验二讲解Word下载.docx》由会员分享，可在线阅读，更多相关《数值分析实验二讲解Word下载.docx（20页珍藏版）》请在冰点文库上搜索。
3.理解插值方法与数据拟合的区别，掌握数据拟合方法解决实际问题的基本步骤和求解理论，并能通过数值实验进行验证。
实验准备：
1.在开始本实验之前，请回顾教科书的相关内容；
2.需要一台准备安装WindowsXPProfessional操作系统和装有数学软件的计算机。
实验内容及要求
A题在某海域测得一些点
处的水深
由表1给出，在矩形区域
内画出海底曲面的图形；
若船的吃水深度为5英尺，请问在矩形区域
里的哪些地方船要避免进入。
表1海域数据表
129
140
103.5
88
185.5
195
105
7.5
141.5
23
147
22.5
137.5
85.5
4
8
6
157.5
107.5
77
81
162
117.5
-6.5
-81
3
56.5
-66.5
84
-33.5
9
B题设从某一实验中测得2个变量
和
的一组数据表2所示。
表2变量关系数据表
1
2
5
7
10
求一代数多项式曲线，使其最好地拟合这组给定数据。
要求：
请用多项式和指数模型进行曲线拟合，观察其结果的精度，再直接调用polyfit命令进行比较。
C题考察函数
（1）设计插值方案，构造龙格（Runge）现象。
（2）设计插值方法消除龙格（Runge）现象，要求至少用分段线性插值和分段三次插值两种算法进行设计，并比较两种算法的精度。
D题已知数据（表3）
表3数据表
0.56160
0.56280
0.56401
0.56521
0.82741
0.82659
0.82577
0.81495
试用Lagrange插值多项式求
时的函数近似值．
说明：
实验过程应包括对问题的简要分析、求解方法、求解步骤、程序及其必要的图表等内容。
实验过程：
A题的实验过程
1、实验中问题的重述
在某海域测得一些点
处的水深由表1给出，在矩形区域
2、对实验的分析
2.1对软件的选择
对于此题来说，利用Matlab10.0软件编程计算。
2.2对“如何画海底曲面图形”的分析
要画出海底曲面图形，首先就应该有该海域海底是平滑的假设；
然后，由于测量点是散乱分布的，所以在平面上先作出测量点的分布图；
其次，利用二维插值方法griddadt函数补充一些点的水深；
最后，再利用surf函数作出海底曲面图。
2.3对“哪些地方需避免进入”的分析
由于船只的吃水深度为5，所以在水深为5以下的区域都是危险区域，船只都需避免。
于是，先作出等高线图，再通过等高线图得到水深小于5的区域。
2.4对结果呈现的分析
由于题中明确要求画出海底曲面，所以利用插值，以曲面图显示结果。
题目中还要求给出哪些地方需要避免进入，所以可以利用图形显示出那一部分区域，也可以给出该危险区域中点的坐标范围。
3、实验求解过程（程序见附录）
首先，利用plot函数，画出所给14个点的位置分布图，分布图见下图1；
然后，利用meshgrid函数，生成一些采样点，便于接下来根据这些采样点绘制出整个海底曲面图形；
其次，利用griddadt函数进行插值；
再次，利用surf函数，画出经过插值后的海底曲面图形，见下图2。
图1测量点的分布图
图2海底曲面图
最后，利用contour函数，画出海水的等高线图，见下图3。
并且，利用clabel函数标记出各个等高线的水深高度。
图3等高线图
于是得到该危险区域为：
图3中水深为-5的等高线区域内部所有的位置。
B题的实验过程
一、用多项式拟合曲线
观察这一实验中测得的变量x、y的数据并进行分析，当i=[1,5]时y随x的增大而减小，当i=[6,9]时，y随x的增大呈线性增大。
那么我们就可以采取n>
=2次多项式进行曲线拟合,然后观察图像的变化及其精度，并分析。
用matlab多项式拟合曲线图像如下：
图4二次多项式拟合图像
图5三次多项式拟合图像
图6四次多项式拟合图像
图7五次多项式拟合图像
图8六次多项式拟合图像
图9七次多项式拟合图像
通过上面六个拟合图像我们可以观察到当拟合多项式的次数越高，我们就猜想：
拟合多项式的次数越高那么拟合曲线的精度就越高，这一结论是否正确呢？
如果一直这样下去会产生什么样的情况？
图10八次多项式拟合图像
从上一图像我们可以看出，在x=[1,3]这个区间时，图像幅度变化的非常大。
但是我们可以观察到下图是当n=12时的图像，可以看出当x=[910]范围内时变化幅度大。
图11次多项式拟合图像
尽管曲线会因n的取值的增大，会不同程度的发生变化。
但是通过上面的图像我们可以总结出当多项式的次数越高时，曲线穿过点的个数是越多的，当n=8时，基本上曲线能穿过每一个点，当n=12时也一样。
二、用polyfit命令拟合曲线
我们直接调用命令进行作图：
p=polyfit（x,y,2）和p=polyfit（x,y,8），分别得到以下结果：
图12polyfit（x,y,2）拟合图像
图13polyfit（x,y,8）拟合图像
我们发现直接用polyfit命令非常方便，只需要改变维数就行了，结果和多项式拟合的图像一样，因此我们可以用这个命令来判断上述的结论，当n=100时拟合图像如下：
图14polyfit（x,y,100）拟合图像
结论：
并非多项式拟合的次数越高越好，比如当n=13时图像就如下图：
图14polyfit（x,y,13）拟合图像
三、用指数模型拟合曲线
同样的道理，指数模型拟合也是一样，得到以下拟合图像：
C题的实验过程
一、实验分析
（1）、要求构造龙格现象，则采用多项式插值的方法呈现。
（2）、针对于
（1）中的龙格现象，采用分段线性插值和分段三次插值的两种算法进行改进。
二、实验求解过程（代码见附录）
（1）、利用多项式插值的方法对函数
进行插值，分别取5、10个点进行插值进而构造出Runge现象，并呈现图像。
（2）、分别采用分段线性插值和分段三次插值的两种算法对
（1）中的Runge现象进行改进计算。
并呈现图像作比较。
三、实验结果
（1）MATLAB编程实现Runge现象的呈现如下图：
（2）分别采用分段线性插值和分段三次插值的两种算法对
（1）中的Runge现象进行改进计算。
并呈现图像如下图：
由图像中的关系可以看出，分段三次插值方法的误差相对低些。
D题的实验过程
针对问题试用Lagrange插值多项式求
时的函数近似值。
首先求出Lagrange插值基函数的值，然后输入表格的数据及插值点，最后计算插值点所对应的函数值。
二、实验求解过程
用Matlab编程可以求出插值点的函数值。
Step1：
在M文件中编写程序，求出
的字符长度；
Step2：
判断m与n是否相等，如果不等，程序结束，反之继续编写；
Step3：
用for循环求插值基函数的值；
Step4：
在命令窗口中调用函数求出插值点对应的函数值；
具体求解的流程图如下：
三、实验求解结果
时的函数近似值分别为：
0.8265,0.8268,0.8231。
附录
A题实验程序
clc
clear
x0=[129140103.588185.5195105157.5107.57781162162117.5];
y0=[7.5141.52314722.5137.585.5-6.5-81356.5-66.584-33.5];
z0=-[48686889988949];
%原始数据x0、y0、z0
[x,y]=meshgrid（75:
5:
200,-50:
150）;
z=griddata（x0,y0,z0,x,y,'
v4'
）;
%散点数据的插值函数griddata（）
figure
（1）
plot（x0,y0,'
r*'
）
figure
（2）
plot3（x0,y0,z0,'
%已知节点的三维散点图
holdon
surf（x,y,z）%画出经插值的图形
figure（3）
[c,h]=contour（x,y,z）;
%画等高线
clabel（c,h）%标明等高线的高程，高程为5的区域避免进入
B题实验程序
%%数据文件（data）的读取
[vd]=xlsread（'
data.xls'
x=v（2,:
y=v（3,:
n=length（x）;
%读取文件中的点的个数
%%最小二乘法
symsab
y1=a*n+b*sum（x）-sum（y）;
y2=a*sum（x）+b*sum（x.^2）-sum（x.*y）;
[ab]=solve（y1,y2,'
a'
'
b'
%误差的计算
y1=a+b.*x;
formatlong
w1=sum（abs（y1-y））/sum（y）
%%指数模型的拟合
%y=A*exp（B/x）
%%matlab内部函数的命令
a1=polyfit（x,y,1）;
y2=a1
（1）.*x+a1
（2）;
w2=sum（abs（y2-y））/sum（y）
C题实验程序
function[]=Runge（a1,a2,n）
%%
（1）构造龙格现象
%n:
构造n次的多项式插值
x1=a1:
（a2-a1）/n:
a2;
y1=1./（1+x1.^2）;
symsx
y=0;
fori=1:
n+1
ss=1;
forj=1:
n+1
ifj~=i
c=（x-x1（j））/（x1（i）-x1（j））;
ss=ss.*c;
else
continue
end
y=y+ss.*y1（i）;
end
%插值多项式
x=-5:
0.01:
5;
yy=eval（y）;
plot（x,yy）
%龙格现象图
End
%%c题问题主程序
%Runge现象的呈现的程序
a=[510];
a1=-5;
a2=5;
length（a）
Runge（a1,a2,a（i））;
holdon
xlabel（'
Runge现象图'
%%分段线性插值
（a2-a1）/a
（2）:
symsx
s1=[];
a
（2）
b0=（1-（x-x1（i））/（x1（i+1）-x1（i）））*y1（i）;
b1=（（x-x1（i））/（x1（i+1）-x1（i）））*y1（i+1）;
s1{1,i}=b0+b1;
%%分段三次插值
y2=-2.*（x1）./（（1+x1.^2）.^2）;
s3=[];
h=x1（i+1）-x1（i）;
c0=（（（x-x1（i））/（x1（i+1）-x1（i））-1）^2）*（2*（x-x1（i））/（x1（i+1）-x1（i））+1）*y1（i）;
c1=（（x-x1（i））/（x1（i+1）-x1（i）））^2*（-2*（x-x1（i））/（x1（i+1）-x1（i））+3）*y1（i+1）;
c2=（（x-x1（i））/（x1（i+1）-x1（i）））*（（x-x1（i））/（x1（i+1）-x1（i））-1）^2*（y2（i））;
c3=（（x-x1（i））/（x1（i+1）-x1（i）））^2*（（x-x1（i））/（x1（i+1）-x1（i））-1）*y2（i+1）;
s3{1,i}=c0+c1+c2+c3;
%%画图
figure
x11=a1:
（a2-a1）/（a
（2）*100）:
y11=1./（1+x11.^2）;
%分段线性插值图
x=x1（i）:
（a2-a1）/（10*a
（2））:
x1（i+1）;
yy=eval（s1{1,i}）;
plot（x,yy,'
r'
holdon
plot（x11,y11）
%axis（[a1,a2,-0.5,2]）;
分段线性图红色：
插值蓝色：
原值'
%分段三次插值
yy=eval（s3{1,i}）;
'
D题实验程序
function[y0,N]=lagrange_eval（X,Y,x0）
%%lagrange插值
%X，Y是一直的插值点的坐标点
%x0是插值点
%y0是lagrange多项式在x0处的值
%N是lagrange插值函数的权系数
m=length（X）;
N=zeros（m,1）;
y0=0;
m
N（i）=1;
ifj~=i
N（i）=N（i）*（x0-X（j））/（X（i）-X（j））;
y0=y0+Y（i）*N（i）;
%%lagrange插值（D题）
x=v（1,:
y=v（2,:
%数据的读取
x1=[0.56260.56350.5645];
%需插值的点
y1=[];
length（x1）
[y0,N]=lagrange_eval（x,y,x1（i））;
y1=[y1y0];
fprintf（'
lagrange插值函数的结果'
y1%lagrange插值函数的结果
实验总结（由学生填写）：
升级会员 