matlab 插值与拟合

实验 2
插值与拟合
实验内容：1. 三种插值方法 2．用 Matlab 计算插值 3．拟合的基本原理 4．用 Matlab 拟合曲线 实验目的： 掌握插值与拟合方法
一、
概念的引入
1. 插值与拟合在现实生活中的应用 机械制造：汽车外观设计 采样数据的重新建构：电脑游戏中场景的显示，地质勘探，医学领域（CT）
二、 插值理论
设函数 y=f(x)在区间[a,b]上连续，在[a,b]上有互异点 x0,x1,…,xn 处取值 y0,y1,…,yn 。如果函数φ(x)在点 xi 上满足φ(xi)=yi (i=0,1,2,…,n)，则称φ(x)是函数 y=f(x)的插值函数，x0,x1,…,xn 是插值节点。若此时φ(x)是代数多项式 P(x)，则称 P(x) 为插值多项式。显然 f(x)≈φ(x)，x∈[a,b]
y ( x) hi[( xi x)(2aiyi yi ' ) yi ]
i 1
n
其中： hi
(
j 1ai xi xj i 1 xi xj i j
2) Matlab 实现 Hermite.m function y=hermite(x0,y0,y1,x) n=length(x0);m=length(x);
3. 分段插值------子区间内，避免函数在某些区间失真 1） 线性插值 已知 n+1 个不同节点 x0,x1,…,xn ,构造分段一次线性多项式 P(x)，使之满足 P(x)在[a,b]上连续 P(xk)=yk P(x)在[xi,xi+1]上是线性函数，P(x)=
y l ( x)
2) Matlab 实现函数
Lagrange.m
function y=lagrange(x0,y0,x) n=length(x0);m=length(x); for i=1:m z=x(i); s=0.0;
for k=1:n p=1.0; for j=1:n if j~=k p=p*(z-x0(j))/(x0(k)-x0(j)); end end s=p*y0(k)+s; end y(i)=s; end 3) 例题的 Matlab 实现 x=[1 2 3 4]; y=[0 –5 –6 3]; lagrange(x,y,2.5)
很显然，当 n=1、插值节点只有两个 xk,xk+1 时 P1(x)=yklk(x)+yk+1lk+1(x) 其中基函数 lk(x)= 2. 牛顿插值 构造 n 次多项式 Nn(x)=f(x0)+f(x0,x1)(x-x0)+f(x0,x1,x2)(x-x0)(x-x1)+… +f(x0,x1,x2,…,xn)(x-x0)(x-x1)…(x-xn) 称为牛顿插值多项式，其中
此时 i ( x ) = ui ( x ) + vi ( x ) ，其中 hi xi 1 xi
4.
三次样条插值------二阶可导 对于给定 n+1 个不同节点 x0,x1,…,xn 及函数值 y0,y1,…,yn，其中 a=x0<x1<…<xn=b。构
造三次样条插值函数 S(x)。S(x)称为三次样条函数时需满足： S(x)在[a,b]上二阶导数连续 S(xk)=yk (k=0,1,…,n) 每个子区间[xk,xk+1]上 S(x)是三次多项式(k=0,1,…,n)
显然 lk (x)满足 lk (xi)=
1(i k ) 0(i k ) f ( n 1) ( ) n 1(x) (n 1)!
此时 Pn(x)≈f(x),误差 Rn(x)=f(x)-Pn(x)=
其中 ∈(a,b)且依赖于 x， n 1(x) =（x-x0）(x-x1)…(x-xn)
i i i 0
n
2） 两点带导数插值---避免尖点、一阶连续
区间[a,b]上两个互异节点 xi,xi+1，已知实数 y i,y i+1,m i,m i+1，为了构造次数不大于 3 的 多项式 i ( x ) 满足条件
i ( xi ) yi ' i ( xi ) mi
i ( xi 1) yi 1 ' i ( xi 1) mi 1 vi ( xi ) 0 vi ( xi ) 0 vi ( xi 1) yi 1 ' vi ( xi 1) mi 1
N3(x)=f(x0)+f(x0,x1)(x-x0)+f(x0,x1,x2)(x-x0)(x-x1)+f(x0,x1,x2,x3)(x-x0)(x-x1)(x-x2) =0+(-5)(x-1)+2*(x-1)(x-2)+1*(x-1)(x-2)(x-3) =x3-4x2+3
三、 Matlab 在插值中的应用
5. 例题 已知函数 y=f(x)的观测值 X Y 1 0 2 -5 3 -6 4 3
求三次插值多项式φ(x)及φ(2.5). 解： (1)拉格朗日插值
l 0( x ) l1( x) l 2( x )
( x 2)( x 3)( x 4) ( x 2)( x 3)( x 4) (1 2)(1 3)(1 4) 6 ( x 1)( x 3)( x 4) ( x 1)( x 3)( x 4) (2 1)(2 3)(2 4) 2 ( x 1)( x 2)( x 4) ( x 1)( x 2)( x 4) (3 1)(3 2)(3 4) 2
3
+ (-6)
( x 1)( x 2)( x 4) ( x 1)( x 2)( x 3) +3 2 6
2
=x -4x +3 φ(x)≈P3(x)= x -4x +3 φ(2.5)=2.5 -4*2.5 +3=-6.375
3
(2)牛顿插值 xi 1 2 3 4 f(xi) 0 –5 –6 3 -5 -1 9 2 5 1 一阶差商 二阶差商 三阶差商
2. Runge 现象及分段线性插值 １）Runge 现象 Runge 在本世纪初发现：在[-1,1]上用 n+1 个等距结点作插值多项式 Pn(x)，使其在个结 点的值与函数 y(x)=1/(1+25x2)在结点的值相等。但在 n→∞时，插值多项式 Pn(x)在区间中 部趋于 y(x)。但对于 0.726≤∣x∣≤1 的 x，Pn(x)严重发散。 通过下面的例子，以图形的方式体会 Runge 现象（ f(x)=1/(1+x ) ） x=[-5:1:5]; y=1./(1+x.^2); x0=[-5:0.1:5]; y0=lagrange(x,y,x0); y1=1./(1+x0.^2); %绘制图形 plot(x0,y0,'--r') hold on plot(x0,y1,'-b')
ui ( xi ) yi ui ' ( xi ) mi 引入 ui ( x ) , vi ( x ) 使之满足 ui ( xi 1) 0 ' ui ( xi 1) 0
2 x xi 1 2 ui ( x) [ yi (mi yi )( x xi )]( ) hi hi 可以求出 vi ( x) [ yi 1 (mi 1 2 yi 1)( x xi 1)]( x xi ) 2 hi hi
for k=1:m yy=0.0; for i=1:n h=1.0; a=0.0; for j=1:n if j~=i h=h*((x(k)-x0(j))/(x0(i)-x0(j)))^2; a=1/(x0(i)-x0(j))+a; end end yy=yy+h*((x0(i)-x(k))*(2*a*y0(i)-y1(i))+y0(i)); end y(k)=yy; end 3） 例题 如下数据表，构造 Hermite 多项式，并求出 Sin0.34 的近似值 x Sin(x) Cos(x) x0=[0.3 0.32 0.35]; y0=[0.29552 0.31457 0.34290]; y1=[0.95534 0.94924 0.93937]; x=[0.3:0.005:0.35]; y=hermite(x0,y0,y1,x); plot(x,y) y=hermite(x0,y0,y1,0.34); y sin(0.34) y2=sin(x); hold on plot(x,y2,'--r') 4. 三次样条插值 y3=interp1(x,y,x0,'*spline'); y3=spline(x,y,x0); plot(x0,y3,'-g') 0.30 0.29552 0.95534 0.32 0.31457 0.94924 0.35 0.34290 0.93937
-0.5 -5 -4 -3 -2 -1 0 1 2 3 4 5 0 1 2
2
1.5
0.5
2）Matlab 实现分段插值------一维插值 interp1 yi=interp1(x,y,xi) yi=interp1(y,xi) yi=interp1(x,y,xi,’method’) linear spline cubic 线性插值 三次样条插值 三次插值 对(x,y)进行插值，计算插值点 xi 的函数值 默认 x=1:n，n 是向量 y 的元素个数 指定特定算法插值，method 可以是如下字符串
x xk 1 xk xk 1
lk+1(x)=
x xk xk 1 xk
f ( x 0, x1) f ( x 0, x1, x 2)

f ( x 0) f ( x1) x 0 x1 f ( x 0, x1) f ( x1, x 2) x0 x2 f ( x 0, x1,..., xn 1) f ( x1, x 2,..., xn ) x 0 xn