六、MATLAB在计算方法中的应用(插值拟合)讲解

样条（11节点）
x=linspace(-6,6,11); y=1./(1+x.^2); xi=linspace(-6,6,50); yi = interp1(x, y, xi,'spline'); yii=1./(1+xi.^2); plot(xi,yi,'r') hold on plot(xi,yii)
分析下边的程序
x = -4:4; y = [0 .15 1.12 2.36 2.36 1.46 .49 .06 0]; cs = spline(x,y); xx = linspace(-4,4,101); plot(x,y,'o',xx,ppval(cs,xx),'-')
ppval(cs,xx)表示由样条函数cs上对应于横坐标xx的 纵坐标的值
Matlab在计算方法中的应用(1)
一 插值
1,插值问题的提法
已知n 1个节点(x j , y j ), j 0,1,, n，其中 x j互不 相同，不妨设a x0 x1 xn b，求任一插 值点 x( x j ) 处的插值 y，(x j , y j ) 可以看成由某 个函数 y g(x)产生的，g 的解析表达式可能十分
样条
x=linspace(-6,6,11); y=1./(1+x.^2); xi=linspace(-6,6,50); yi = interp1(x, y, xi,'spline'); yii=1./(1+xi.^2); plot(xi,yi,'r') hold on plot(xi,yii)
7，其他的插值方法命令
一维插值: yi = interp1(x, y, xi, ‘method ’)
method
nearest linear spline cubic
最近点插值
线性插值 样条插值 立方插值
二维插值 zi=interp2(x, y, z, xi, yi, ‘method’) 三维插值vi = interp3(x,y,z,v,xi,yi,zi, ‘method’)
yi = interp1(x, y, xi,‘linear’); ％由节点计算插值函数在xi处的函数值
yii=1./(1+xi.^2); ％计算出xi对应的精确函数值
plot(xi,yi,‘r’) ％绘出插值函数的图形
hold on ％保持图形窗口
plot(xi,yii) ％绘出原函数的精确图形 x=linspace(-6,6,11);
立方插值
x=linspace(-6,6,5); y=1./(1+x.^2); xi=linspace(-6,6,50); yi = interp1(x, y, xi,'cubic'); yii=1./(1+xi.^2); plot(xi,yi,'r') hold on plot(xi,yii)
8，各种一维插值的比较
最邻近 x=linspace(-6,6,5); y=1./(1+x.^2); xi=linspace(-6,6,50); yi = interp1(x, y, xi,'nearst'); yii=1./(1+xi.^2); plot(xi,yi,'r') hold on plot(xi,yii)
例 对下面一组数据作二次多项式拟合
xi 0.1 0.2 0.4 0.5 0.6 0.7 0.8 0.9 1 yi 1.978 3.28 6.16 7.34 7.66 9.58 9.48 9.30 11.2
12
10
8
6
4
2
0
-2
0
0.2
0.4
0.6
0.8
1
1）输入以下命令： x=0:0.1:1; y=[-0.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2]; A=polyfit(x,y,2) z=polyval(A,x); plot(x,y,'k+',x,z,'r') %作出数据点和拟合曲线的图形 2）计算结果： Ａ = -9.8108 20.1293 -0.0317
yi = interp1(x, y, xi, ‘linear ’)
命令解释： interp1（一维插值函数） x, y,（用于存放已知节点的坐标） ‘linear ’（可选参数，linear表示用分段线性插值的方法做一维值） xi ，yi（得到插值函数后并不是返回具体的函数而是对于给定的xi返 回他所对应的插值函数的函数值yi）
Ln (x) anxn an1xn1 a1x a0 对于节点 (x j , y j )
应有 Ln (x j ) y j j 0,1,, n
4,拉格朗日多项式插值的matlab实现
拉格朗日多项式插值没有现成的matlab命令 可用，如果要求拉格朗日插值多项式则必须自己 编制matlab函数。
复杂，或不存在封闭形式，也可以未知。
2,求解的基本思路
构造一个相对简单的函数y f (x) ，使 f 通过全部节点， 即 f (x j ) y j ( j 0,1,, n)，再用 f (x) 计算插值，即 y f (x)
f
3,拉格朗日多项式插值
从理论和计算角度看，多项式是最简单的函数，设f(x) 是n次多项式，记作
线性插值 较快，有足够精度
linear
立方插值 较慢，精度高，平滑性好。 cubic
样条插值 最慢，精度高，平滑性最好 spline
最邻近插值 最快，精度低，不平滑
arst
9，样条插值（spline）
x = 0:10; y = sin(x); xx = 0:0.25:10; yy = spline(x,y,xx); ％直接用样条插值命令求插值 plot(x,y,'o',xx,yy)
x = lsqcurvefit(@(x,xdata) x(1)*sin(xdata)+x(2),[2 7],xdata,ydata);
xdd=4:0.1:6;
ydd= x(1)*sin(xdd)+x(2);
plot(xdd,ydd)
hold on
plot(xdata,ydata,'ro')
xdd=4:0.1:6; ydd= x(1)*sin(xdd)+x(2); plot(xdd,ydd) hold on plot(xdata,ydata,'ro')
或写为
xdata = [5;4;6];
% example xdata
ydata = 3*sin([5;4;6])+6; % example ydata
例
g (x)

1 1 x2
,
5 x 5
采用拉格朗日多项式插值：选取不同插值 节点个数n+1，其中n为插值多项式的次数，当n 分别取2,4,6,8,10时，绘出插值结果图形.
5，分段线性插值
插值函数为分段的线性函数，即用线段把相 邻的两个节点连接起来而得的函数
6，分段线性插值的matlab命令
二 拟合
1，离散数据的多项式拟合
已知数据表 x f(x)
x1 x2 ·········· xm y1 y2 ·········· ym
求拟合函数: (x) = a0 + a1x + ······+ anxn
使得
m
[(a0 a1xk an xkn ) yk ]2
k 1
例题
建立myfun.m内容为
function F = myfun(x,xdata) F = x(1)*sin(xdata)+x(2)
运行下列命令
xdata = [5;4;6]; xdata
% example
ydata = 3*sin([5;4;6])+6; % example ydata
x = lsqcurvefit(@myfun,[1 6],xdata, ydata);
达到最小
2，多项式拟合命令
1. 作多项式f(x)=a1xm+ …+amx+am+1拟合,可利用已有程序: a=polyfit(x,y,m)
输出拟合多项式系数
输入同长度
a=[a1, …am , am+1] (数组)） 的数组X，Y
拟合多项 式次数
2.多项式在x处的值y可用以下命令计算： y=polyval（a，x）
f (x) 9.8108x2 20.1293x 0.0317
4，可转换为多项式拟合
例题：数据为x=[-70 -50 -30 -10 20] y=[0.182 0.219 0.077 0.332 1]
拟合方程为： y=exp（ax+b）
将拟合方程两边取对数则得 lny＝ax＋b 令lny＝z 得 Z＝ax＋b 由y值算出z的值对z与x用一次多项式拟合。
x1=[-70 -50 -30 -10 20]; y1=[-1.704 -1.519 -2.564 -1.103 0]; subplot(1,2,1); plot(x1,y1,'o'); hold on;
％将y中的数据取对数y1=lny
p1=polyfit(x1,y1,1); yy1=polyval(p1,-70:20); plot(-70:20,yy1,'r');
y=1./(1+x.^2);
xi=linspace(-6,6,50);
yi = interp1(x, y, xi,'linear');