-175-
第九章 插值与拟合
插值：求过已知有限个数据点的近似函数。
拟合：已知有限个数据点，求近似函数，不要求过已知数据点，只要求在某种意义下它在这些点上的总偏差最小。
插值和拟合都是要根据一组数据构造一个函数作为近似，由于近似的要求不同，二者的数学方法上是完全不同的。而面对一个实际问题，究竟应该用插值还是拟合，有时容易确定，有时则并不明显。
§1 插值方法
下面介绍几种基本的、常用的插值：拉格朗日多项式插值、牛顿插值、分段线性插值、Hermite 插值和三次样条插值。
1.1 拉格朗日多项式插值 1.1.1 插值多项式
用多项式作为研究插值的工具，称为代数插值。其基本问题是：已知函数)(x f 在区间],[b a 上1+n 个不同点n x x x ,,,10L 处的函数值)(i i x f y =),,1,0(n i L =，求一个至多n 次多项式
n n n x a x a a x +++=L 10)(? （1）
使其在给定点处与)(x f 同值，即满足插值条件
),,1,0()()(n i y x f x i
i i n L ===? （2） )(x n ?称为插值多项式，),,1,0(n i x i L =称为插值节点，简称节点，],[b a 称为插值区间。从几何上看，n 次多项式插值就是过1+n 个点))(,(i i x f x ),,1,0(n i L =，作一条多项式曲线)(x y n ?=近似曲线)(x f y =。
n 次多项式（1）有1+n 个待定系数，由插值条件（2）恰好给出1+n 个方程
??????
?=++++=++++=++++n n n n n n n n n
n y
x a x a x a a y x a x a x a a y x a x a x a a L L L L L L L L L L L L L L L 2210
1
12
1211000202010 （3） 记此方程组的系数矩阵为A ，则
n
n
n n n
n
x x x x x x x x x A L L
L L L L L L L L 212110200111)det(=
是范德蒙特(Vandermonde)行列式。当n x x x ,,,10L 互不相同时，此行列式值不为零。因此方程组（3）有唯一解。这表明，只要1+n 个节点互不相同，满足插值要求（2）的
插值多项式（1）是唯一的。
插值多项式与被插函数之间的差
)()()(x x f x R n n ??=
-176-
称为截断误差，又称为插值余项。当)(x f 充分光滑时，
),(),()!
1()
()()()(1)1(b a x n f x L x f x R n n n n ∈+=
?=++ξωξ 其中∏=+?=
n
j j
n x
x x 0
1)()(ω。
1.1.2 拉格朗日插值多项式
实际上比较方便的作法不是解方程（3）求待定系数，而是先构造一组基函数
)())(()()
())(()()(110110n i i i i i i n i i i x x x x x x x x x x x x x x x x x l ????????=
+?+?L L L L
)10( , 0,n ,,i x x x
x n i
j j j
i j
L =??=∏≠= )(x l i 是n 次多项式，满足
??
?=≠=i
j i
j x l j i 1
0)( 令
∑∑∏==≠=?
?
??????
????==n
i n
i n
i
j j j i j i i i n x x x x y x l y x L 000)()( （4） 上式称为n 次Lagrange 插值多项式，由方程（3）解的唯一性，1+n 个节点的n 次
Lagrange 插值多项式存在唯一。
1.1.3 用Matlab 作Lagrange 插值
Matlab 中没有现成的Lagrange 插值函数，必须编写一个M 文件实现Lagrange 插值。 设n 个节点数据以数组0,0y x 输入（注意Matlat 的数组下标从1开始），m 个插值点以数组x 输入，输出数组y 为m 个插值。编写一个名为lagrange.m 的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
-177-
1.2 牛顿（Newton ）插值
在导出Newton 公式前，先介绍公式表示中所需要用到的差商、差分的概念及性质。 1.2.1 差商
定义 设有函数L ,,,),(210x x x x f 为一系列互不相等的点，称
j
i j i x x x f x f ??)()()(j i ≠为)(x f 关于点j i x x ,一阶差商（也称均差）记为],[j i x x f ，即
j
i j i j i x x x f x f x x f ??=
)()(],[
称一阶差商的差商
k
i k j j i x x x x f x x f ??]
,[],[
为)(x f 关于点k j i x x x ,,的二阶差商，记为],,[k j i x x x f 。一般地，称
k
k k x x x x x f x x x f ???021110]
,,,[],,,[L L
为)(x f 关于点k x x x ,,,10L 的k 阶差商，记为
k
k k k x x x x x f x x x f x x x f ??=
?02111010]
,,,[],,,[],,,[L L L 容易证明，差商具有下述性质： ],[],[i j j i x x f x x f =
],,[],,[],,[k i j j k i k j i x x x f x x x f x x x f == 1.2.2 Newton 插值公式 线性插值公式可表成
],[)()()(10001x x f x x x f x ?+=?
称为一次Newton 插值多项式。一般地，由各阶差商的定义，依次可得
],[)()()(000x x f x x x f x f ?+=
],,[)(],[],[101100x x x f x x x x f x x f ?+=
],,,[)(],,[],,[210221010x x x x f x x x x x f x x x f ?+=
L L
],,,[)(],,,[],,,[01010n n n n x x x f x x x x x f x x x f L L L ?+=?
将以上各式分别乘以,),)((),(,1100L x x x x x x ???)())((110????n x x x x x x L ，然
后相加并消去两边相等的部分，即得
]
,,,,[)())(( ],,,[)())(( ],[)()()(1010101101000n n n n x x x x f x x x x x x x x x f x x x x x x x x f x x x f x f L L L L L
???+???++?+=?
记
-178-
]
,,,,[)( ]
,,,,[)())(()(],,,[)())(( ],[)()()(1011010101101000n n n n n n n n x x x x f x x x x x f x x x x x x x R x x x f x x x x x x x x f x x x f x N L L L L L L
+?=???=???++?+=ω
显然，)(x N n 是至多n 次的多项式，且满足插值条件，因而它是)(x f 的n 次插值多项式。这种形式的插值多项式称为Newton 插值多项式。)(x R n 称为Newton 插值余
项。
Newton 插值的优点是：每增加一个节点，插值多项式只增加一项，即
],,,[)()()()(11001++??+=n n n n x x x f x x x x x N x N L L
因而便于递推运算。而且Newton 插值的计算量小于Lagrange 插值。
由插值多项式的唯一性可知，Newton 插值余项与Lagrange 余项也是相等的，即
)
,()()!
1()(],,,,[)()(1)1(101b a x n f x x x x f x x R n n n n n ∈+==+++ξωξωL
由此可得差商与导数的关系
!
)
(],,,[)(10n f x x x f n n ξ=L
其中}{max },{min ),,(00i n
i i n
i x x ≤≤≤≤==∈βαβαξ。
1.2.3 差分
当节点等距时，即相邻两个节点之差（称为步长）为常数，Newton 插值公式的形式会更简单。此时关于节点间函数的平均变化率（差商）可用函数值之差（差分）来表示。
定义 设有等距节点),,1,0(0n k kh x x k L =+=，步长h 为常数，)(k k x f f =。称相邻两个节点1,+k k x x 处的函数值的增量)1,,1,0(1?=?+n k f f k k L 为函数)(x f 在点k x 处以h 为步长的一阶差分，记为k f Δ，即
),,1,0(1n k f f f k k k L =?=Δ+
类似地，定义差分的差分为高阶差分。如二阶差分为
)2,,1,0(12
?=Δ?Δ=Δ+n k f f f k
k k L
一般地，m 阶差分为 ),3,2(111
L =Δ?Δ
=Δ?+?k f f f k
m k m k m
，
上面定义的各阶差分又称为向前差分。常用的差分还有两种： 1??=?k k k f f f
称为)(x f 在k x 处以h 为步长的向后差分； ?????
?
???????
?+
=22h x f h x f f k k k δ 称为)(x f 在k x 处以h 为步长的中心差分。一般地，m 阶向后差分与m 阶中心差分公式为
-179-
111
??????=?k m k m k m f f f
2
112
11?
?+??=k m k m k m f
f
f δδδ
差分具有以下性质：
（i ）各阶差分均可表成函数值的线性组合，例如
j m k m
j j k m
f j m f ?+=∑?????????=
Δ )1(0
j k m
j j k m f j m f ?=∑???
??????=? )1(0
（ii ）各种差分之间可以互化。向后差分与中心差分化成向前差分的公式如下： m k m
k m
f f ?Δ=?
2
m
k m k m f
f ?
Δ=δ
1.2.4 等距节点插值公式
如果插值节点是等距的，则插值公式可用差分表示。设已知节点kh x x k +=0),,2,1,0(n k L =，则有
)
())((!)( )())(](,,,[)
](,[)()(1100
000110100100?????Δ++?Δ+=???++?+=n n n n n n x x x x x x h n f x x h
f f x x x x x x x x x f x x x x f x f x N L L L L L
若令th x x +=0，则上式又可变形为
0000!
)1()1()(f n n t t t f t f th x N n
n Δ+??+
+Δ+=+L L 上式称为Newton 向前插值公式。
1.3 分段线性插值
1.3.1 插值多项式的振荡
用Lagrange 插值多项式)(x L n 近似))((b x a x f ≤≤，虽然随着节点个数的增加，
)(x L n 的次数n 变大，多数情况下误差|)(|x R n 会变小。但是n 增大时，)(x L n 的光滑性变坏，有时会出现很大的振荡。理论上，当∞→n ，在],[b a 内并不能保证)(x L n 处处收敛于)(x f 。Runge 给出了一个有名的例子：
]5,5[,11
)(2
?∈+=x x
x f 对于较大的||x ，随着n 的增大，)(x L n 振荡越来越大，事实上可以证明，仅当63.3||≤x 时，才有)()(lim x f x L n n =∞
→，而在此区间外，)(x L n 是发散的。
高次插值多项式的这些缺陷，促使人们转而寻求简单的低次多项式插值。 1.3.2 分段线性插值 简单地说，将每两个相邻的节点用直线连起来，如此形成的一条折线就是分段线性插值函数，记作)(x I n ，它满足i i n y x I =)(，且)(x I n 在每个小区间],[1+i i x x 上是线性
-180-
函数),,1,0(n i L =。
)(x I n 可以表示为 ∑==n
i i i n x l y x I 0)()(
????
????
?=∈??=∈??=+++???其它时舍去时舍去，）（，）
0 ],[ 0( ],[ ,)(111
11
1
n i x x x x x x x i x x x x x x x x l i i i i i i i i i
i i )(x I n 有良好的收敛性，即对于],[b a x ∈有，
)()(lim x f x I n n =∞
→。
用)(x I n 计算x 点的插值时，只用到x 左右的两个节点，计算量与节点个数n 无关。但n 越大，分段越多，插值误差越小。实际上用函数表作插值计算时，分段线性插值就足够了，如数学、物理中用的特殊函数表，数理统计中用的概率分布表等。
1.3.3 用Matlab 实现分段线性插值
用Matlab 实现分段线性插值不需要编制函数程序，Matlab 中有现成的一维插值函数interp1。
y=interp1(x0,y0,x,'method')
method 指定插值的方法，默认为线性插值。其值可为： 'nearest' 最近项插值 'linear' 线性插值
'spline' 逐段3次样条插值 'cubic' 保凹凸性3次插值。 所有的插值方法要求x0是单调的。
当x0为等距时可以用快速插值法，使用快速插值法的格式为'*nearest'、'*linear'、'*spline'、'*cubic'。
1.4 埃尔米特(Hermite)插值 1.4.1 Hermite 插值多项式 如果对插值函数，不仅要求它在节点处与函数同值，而且要求它与函数有相同的一阶、二阶甚至更高阶的导数值，这就是Hermite 插值问题。本节主要讨论在节点处插值函数与函数的值及一阶导数值均相等的Hermite 插值。
设已知函数)(x f y =在1+n 个互异节点n x x x ,,,10L 上的函数值)(i i x f y =
),,1,0(n i L =和导数值)(''i i x f y = ),,1,0(n i L =，要求一个至多12+n 次的多项式)(x H ，使得
i i y x H =)( i i y x H ')('= ),,1,0(n i L = 满足上述条件的多项式)(x H 称为Hermite 插值多项式。
Hermite 插值多项式为
-181-
∑=+??=n
i i i i i i i y y y a x x h x H 0
])'2)([()(
其中∏≠=????
??
????=n i
j j j i j
i x x x x h 02
，∑≠=?=n i
j j j i i x x a 01。 1.4.2 用Matlab 实现Hermite 插值
Matlab 中没有现成的Hermite 插值函数，必须编写一个M 文件实现插值。
设n 个节点的数据以数组0x （已知点的横坐标）
，0y （函数值），1y （导数值）输入（注意Matlat 的数组下标从1开始），m 个插值点以数组x 输入，输出数组y 为m 个插值。编写一个名为hermite.m 的M 文件： function y=hermite(x0,y0,y1,x); n=length(x0);m=length(x); 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
1.5 样条插值
许多工程技术中提出的计算问题对插值函数的光滑性有较高要求，如飞机的机翼外形，内燃机的进、排气门的凸轮曲线，都要求曲线具有较高的光滑程度，不仅要连续，而且要有连续的曲率，这就导致了样条插值的产生。
1.5.1 样条函数的概念
所谓样条（Spline ）本来是工程设计中使用的一种绘图工具，它是富有弹性的细木条或细金属条。绘图员利用它把一些已知点连接成一条光滑曲线（称为样条曲线），并使连接点处有连续的曲率。
数学上将具有一定光滑性的分段多项式称为样条函数。具体地说，给定区间],[b a 的一个分划
b x x x x a n n =<<<<=Δ?110:L 如果函数)(x s 满足：
（i ）在每个小区间)1,,1,0](,[1?=?n i x x i i L 上)(x s 是k 次多项式；
（ii ）)(x s 在],[b a 上具有1?k 阶连续导数。
则称)(x s 为关于分划Δ的k 次样条函数，其图形称为k 次样条曲线。n x x x ,,,10L 称为样条节点，121,,,?n x x x L 称为内节点，n x x ,0称为边界点，这类样条函数的全体记做
-182-
),(k S P Δ，称为k 次样条函数空间。
显然，折线是一次样条曲线。 若),()(k S x s P Δ∈，则)(x s 是关于分划Δ的k 次多项式样条函数。k 次多项式样条函数的一般形式为
∑
∑
?=+=?+=1
1
)(!
!
)(n j k
j j
k
i i
i k x x k i x x s βα，
其中),,1,0(k i i L =α和)1,,2,1(?=n j j L β均为任意常数，而
?????<≥?=?+
j
j
k
j k j x x x x x x x x ,0,)()(，（1,,2,1?=n j L ）
在实际中最常用的是2=k 和3的情况，即为二次样条函数和三次样条函数。 二次样条函数：对于],[b a 上的分划b x x x a n =<<<=ΔL 10:，则
∑
?=+Δ∈?++
+=1
1
22
2
102)2,()(!
2!
2)(n j P j j
S x x x x x s βααα，
（5） 其中??
??
?<≥?=?+j j j j x x x x x x x x ,0,)()(2
2
，（1,,2,1?=n j L ）。 三次样条函数：对于],[b a 上的分划b x x x a n =<<<=ΔL 10:，则
∑
?=+Δ∈?++
+
+=1
1
33
3
2
2
103)3,()(!
3!
3!
2)(n j P j j
S x x x x x x s βαααα，
（6） 其中??
??
?<≥?=?+j j j j x x x x x x x x ,0,)()(3
3
，（1,,2,1?=n j L ）。 利用样条函数进行插值，即取插值函数为样条函数，称为样条插值。例如分段线性插值
是一次样条插值。下面我们介绍二次、三次样条插值。 1.5.2 二次样条函数插值 首先，我们注意到)2,()(2Δ∈P S x s 中含有2+n 个特定常数，故应需要2+n 个插值条件，因此，二次样条插值问题可分为两类： 问题（1）： 已知插值节点i x 和相应的函数值),,1,0(n i y i L =以及端点0x （或n x ）处的导数值0'y （或n y '），求)2,()(2Δ∈P S x s 使得
??
?====)
)('()(')
,,2,1,0()(0022n n n i i y x s y x s n i y x s 或L （7）
问题（2）：
已知插值节点i x 和相应的导数值),,2,1,0('n i y i L =以及端点0x （或n x ）处的函
数值0y （或n y ），求)2,()(2Δ∈P S x s 使得
??
?====))(()()
,,2,1,0(')('200
22n n i i y x s y x s n i y x s 或L （8）
-183-
事实上，可以证明这两类插值问题都是唯一可解的。 对于问题（1），由条件（7）
????????
??
?=+===?+++==++==++=∑?=0
0210211
2
2210212121101
202
020100
2')(')
,,3,2()(2121)(21)(2
1)(y x x s n j y x x x x x s y
x x x s y x x x s j i j
i j i j j j ααβαααααααααL
引入记号T
n X ),,,,,(11210?=ββαααL 为未知向量，)',,,,(010y y y y C n L =为已
知向量。
??
??????????
?
????
??
??
??
???
???=?001
0)(21)(2
1
2110
)(2
1
2110
2110
02110
2
1212212222
2112
00L
L M
M M M M L L
L x x x x x x x x x x x x x x x A n n n n n
于是，问题转化为求方程组C AX =的解T n X
),,,,,(11210?=ββαααL 的问题，
即可得到二次样条函数)(2x s 的表达式。 对于问题（2）的情况类似。 1.5.3 三次样条函数插值 由于)3,()(3Δ∈P S x s 中含有3+n 个待定系数，故应需要3+n 个插值条件，已知插值节点i x 和相应的函数值),,2,1,0()(n i y x f i i L ==，这里提供了1+n 个条件，还需要2个边界条件。
常用的三次样条函数的边界条件有3种类型：
（i ）n y b s y a s ')(',')('303==。由这种边界条件建立的样条插值函数称为)(x f 的完备三次样条插值函数。
特别地，0''0==n y y 时，样条曲线在端点处呈水平状态。
如果)('x f 不知道，我们可以要求)('3x s 与)('x f 在端点处近似相等。这时以
3210,,,x x x x 为节点作一个三次Newton 插值多项式)(x N a ，以321,,,???n n n n x x x x 作一个三次Newton 插值多项式)(x N b ，要求
)(')('),(')('b N b s a N a s b a ==
由这种边界条件建立的三次样条称为)(x f 的Lagrange 三次样条插值函数。
（ii ）3303")(",")("y b s y a s ==。特别地0""==n n y y 时，称为自然边界条件。
-184-
（iii ）)0(")0("),0(')0('3333?=+?=+b s a s b s a s ，(这里要求=+)0(3a s
)0(3?b s )此条件称为周期条件。
1.5.4 三次样条插值在Matlab 中的实现
在Matlab 中数据点称之为断点。如果三次样条插值没有边界条件，最常用的方法，就是采用非扭结（not-a-knot ）条件。这个条件强迫第1个和第2个三次多项式的三阶导数相等。对最后一个和倒数第2个三次多项式也做同样地处理。
Matlab 中三次样条插值也有现成的函数： y=interp1(x0,y0,x,'spline')； y=spline(x0,y0,x)；
pp=csape(x0,y0,conds)，y=ppval(pp,x)。
其中x0,y0是已知数据点，x 是插值点，y 是插值点的函数值。
对于三次样条插值，我们提倡使用函数csape ，csape 的返回值是pp 形式，要求插值点的函数值，必须调用函数ppval 。
pp=csape(x0,y0)：使用默认的边界条件，即Lagrange 边界条件。 pp=csape(x0,y0,conds)中的conds 指定插值的边界条件，其值可为： 'complete' 边界为一阶导数，即默认的边界条件 'not-a-knot' 非扭结条件 'periodic' 周期条件
'second' 边界为二阶导数，二阶导数的值[0, 0]。 'variational' 设置边界的二阶导数值为[0,0]。
对于一些特殊的边界条件，可以通过conds 的一个21×矩阵来表示，conds 元素的取值为1，2。此时，使用命令
pp=csape(x0,y0_ext,conds)
其中y0_ext=[left, y0, right]，这里left 表示左边界的取值，right 表示右边界的取值。
conds(i)=j 的含义是给定端点i 的j 阶导数，
即conds 的第一个元素表示左边界的条件，第二个元素表示右边界的条件，conds=[2,1]表示左边界是二阶导数，右边界是一阶导数，对应的值由left 和right 给出。
详细情况请使用帮助help csape 。 例1 机床加工
待加工零件的外形根据工艺要求由一组数据),(y x 给出（在平面情况下），用程控铣床加工时每一刀只能沿x 方向和y 方向走非常小的一步，这就需要从已知数据得到加工所要求的步长很小的),(y x 坐标。
表1中给出的y x ,数据位于机翼断面的下轮廓线上，假设需要得到x 坐标每改变0.1时的y 坐标。试完成加工所需数据，画出曲线，并求出0=x 处的曲线斜率和
1513≤≤x 范围内y 的最小值。
表1
x 0 3 5 7 9 11 12 13 14 15 y 0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6
要求用Lagrange 、分段线性和三次样条三种插值方法计算。
解 编写以下程序： clc,clear
x0=[0 3 5 7 9 11 12 13 14 15];
-185-
y0=[0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6]; x=0:0.1:15;
y1=lagrange(x0,y0,x); %调用前面编写的Lagrange 插值函数 y2=interp1(x0,y0,x);
y3=interp1(x0,y0,x,'spline');
pp1=csape(x0,y0); y4=ppval(pp1,x);
pp2=csape(x0,y0,'second'); y5=ppval(pp2,x); fprintf('比较一下不同插值方法和边界条件的结果:\n') fprintf('x y1 y2 y3 y4 y5\n') xianshi=[x',y1',y2',y3',y4',y5'];
fprintf('%f\t%f\t%f\t%f\t%f\t%f\n',xianshi')
subplot(2,2,1), plot(x0,y0,'+',x,y1), title('Lagrange')
subplot(2,2,2), plot(x0,y0,'+',x,y2), title('Piecewise linear') subplot(2,2,3), plot(x0,y0,'+',x,y3), title('Spline1') subplot(2,2,4), plot(x0,y0,'+',x,y4), title('Spline2') dyx0=ppval(fnder(pp1),x0(1)) %求x=0处的导数 ytemp=y3(131:151);
index=find(ytemp==min(ytemp)); xymin=[x(130+index),ytemp(index)]
计算结果略。
可以看出，拉格朗日插值的结果根本不能应用，分段线性插值的光滑性较差（特别是在14=x 附近弯曲处），建议选用三次样条插值的结果。
1.6 B 样条函数插值方法 1.6.1 磨光函数
实际中的许多问题，往往是既要求近似函数（曲线或曲面）有足够的光滑性，又要求与实际函数有相同的凹凸性，一般插值函数和样条函数都不具有这种性质。如果对于一个特殊函数进行磨光处理生成磨光函数（多项式），则用磨光函数构造出样条函数作为插值函数，既有足够的光滑性，而且也具有较好的保凹凸性，因此磨光函数在一维插值（曲线）和二维插值（曲面）问题中有着广泛的应用。
由积分理论可知，对于可积函数通过积分会提高函数的光滑度，因此，我们可以利用积分方法对函数进行磨光处理。
定义 若)(x f 为可积函数，对于0>h ，则称积分
∫+
?=22
,1)(1)(h
x h x h dt t f h x f 为)(x f 的一次磨光函数，h 称为磨光宽度。
同样的，可以定义)(x f 的k 次磨光函数为
∫+??=2
2
,1,)(1)(h
x h x h k h k dt t f h x f （1>k ）
事实上，磨光函数)(,x f h k 比)(x f 的光滑程度要高，且当磨光宽度h 很少时)
(,x f h k 很接近于)(x f 。
1.6.2 等距B 样条函数
对于任意的函数)(x f ，定义其步长为1的中心差分算子δ如下：
-186-
)2
1()21()(??+=x f x f x f δ
在此取0
)(+=x x f ，则
2121+++?????????????
?+=x x x δ 是一个单位方波函数（如图1），记0
0)(+=Ωx x δ。并取1=h ，对)(0x Ω进行一次磨光
得
∫
∫
+
?+++
????
???????????????????+=Ω=Ω212
10021
2
1012121)()(x x x x dt t t dt t x +++?+++
?+?+=?=
∫∫
)1(2)1(1
1
0x x x dt t dt t x x x x
显然)(1x Ω是连续的（如图1）。
图1 )(0x Ω和)(1x Ω的图形
类似地可得到k 次磨光函数为
∑+=++???????++?=Ω1
121!)1()(k j k
j
k j k j k x k C x 实际上，可以证明：)(x k Ω是分段k 次多项式，且具有1?k 阶连续导数，其k 阶
导数有2+k 个间断点，记为2
1
+?=k j x j （1,,2,1,0+=k j L ）。从而可知)(x k Ω是对应于分划+∞<<<<
条函数，简称为k 次B 样条。由于样条节点为2
1
+?=k j x j （1,,2,1,0+=k j L ）是等距的，故)(x k Ω又称为k 次等距B 样条函数。 对于任意函数)(x f 的k 次磨光函数，由归纳法可以得到：
dt t f h t x h x f k h k ∫∞+∞???
?
?
????Ω=)(1)(1,（22h x t h x +≤≤?） 特别地，当1)(=x f 时，有111=??
?
????Ω∫∞+∞??dt h t x h k ，从而1)(=Ω∫+∞∞?dx x k ，且当1≥k 时有递推关系
-187-
????????????
?Ω???????????????+Ω?????
?++=
Ω??212121211)(1
1x x k x k x k x k k k 1.6.3 一维等距B 样条函数插值
等距B 样条函数与通常的样条有如下的关系：
定理 设有区间],[b a 的均匀分划jh x x j +=Δ0:（n j ,,2,1,0L =），n
a
b h ?=，则对任意k 次样条函数),()(k S x s P k Δ∈都可以表示为B 样条函数族
1
021?=?=???
????????
?+???Ωn j k j k
k j h x x 的线性组合。
根据定理，如果已知曲线上一组点),(j j y x ，其中jh x x j +=0（n j h ,,2,1,0,0L =>），则可以构造出一条样条磨光曲线（即为B 样条函数族的线性组合）
∑??=????
????Ω=
1
0)(n k
j k j k j h x x c x s 其中j c （1,,1,?+??=n k k j L ）为待定常数。用它来逼近曲线，既有较好的精度，
又有良好的保凸性。 实际中，最常用的是3=k 的情况，即一般形式为
∑+?=???
?????Ω=
1
1
03
3)(n j j j h x x c x s 其中3+n 个待定系数j c （1,,0,1+?=n j L ）可以由插值条件确定。
对于插值条件
????
?====)
,0(')(')
,2,1,0()(33n j y x s n j y x s j j j j L 有
????
?
?
??
???=?Ω===?Ω==?Ω=∑∑∑+?=+?=+?=n
n j j n i n j j i n j j y j n c h x s n i y j i c x s y j c h x s ')('1)(',,2,1,0,)()(')('1)('31
1311330
1
1
303L （9）
注意到)(3x Ω的局部非零性及其函数值：32)0(3=
Ω，6
1
)1(3=±Ω，当2≥x 时0)(3=Ωx ；且由)21()2
1
()('223?Ω?+Ω=Ωx x x 知，0)0('3=Ω，21)1('3m =±Ω，
当2≥x 时0)('3=Ωx 。则（9）式的每一个方程只有三个非零系数，具体为
-188-
?????=+?==++=+?+?+??n n n i i i i hy c c n i y c c c hy c c '
2,,2,1,0,64'21
111011L （10）
由方程组（10）容易求解出j c （1,,0,1+?=n j L ），即可得到三次样条函数)
(3x s 表达式。 1.6.4 二维等距B 样条函数插值 设有空间曲面),(y x f z =（未知），如果已知二维等距节点=),(j i y x
),(00τj y ih x ++（0,>τh ）上的值为ij z （m j n i ,,1,0;,,1,0L L ==），则相应的B
样条磨光曲面的一般形式为
???
?????Ω????????Ω=
∑∑??=??=j y y i h x x c y x s l
n k i k m l
j ij τ0101
),( 其中ij c （m j n i ,,1,0;,,1,0L L ==）为待定常数，l k ,可以取不同值，常用的也是
2,=l k 和3的情形。这是一种具有良好保凸性的光滑曲面（函数），在工程设计中是常
用的，但只能使用于均匀划分或近似均匀划分的情况。
1.7 二维插值
前面讲述的都是一维插值，即节点为一维变量，插值函数是一元函数（曲线）。若节点是二维的，插值函数就是二元函数，即曲面。如在某区域测量了若干点（节点）的高程（节点值），为了画出较精确的等高线图，就要先插入更多的点（插值点），计算这些点的高程（插值）。
1.7.1 插值节点为网格节点
已知n m ×个节点：),,(ij j i z y x （n j m i ,,2,1;,,2,1L L ==），且m x x < n y y < Matlab 中有一些计算二维插值的程序。如 z=interp2(x0,y0,z0,x,y,'method') 其中x0,y0分别为m 维和n 维向量，表示节点，z0为m n ×维矩阵，表示节点值，x,y 为一维数组，表示插值点，x 与y 应是方向不同的向量，即一个是行向量，另一个是列向量，z 为矩阵，它的行数为y 的维数，列数为x 的维数，表示得到的插值，'method'的用法同上面的一维插值。 如果是三次样条插值，可以使用命令 pp=csape({x0,y0},z0,conds,valconds)，z=fnval(pp,{x,y}) 其中x0,y0分别为m 维和n 维向量，z0为n m ×维矩阵，z 为矩阵，它的行数为x 的维数，列数为y 的维数，表示得到的插值，具体使用方法同一维插值。 例2 在一丘陵地带测量高程，x 和y 方向每隔100米测一个点，得高程如2表，试插值一曲面，确定合适的模型，并由此找出最高点和该点的高程。 解 编写程序如下： clear,clc x=100:100:500; y=100:100:400; z=[636 697 624 478 450 698 712 630 478 420 -189- 680 674 598 412 400 662 626 552 334 310]; pp=csape({x,y},z') xi=100:10:500;yi=100:10:400 cz1=fnval(pp,{xi,yi}) cz2=interp2(x,y,z,xi,yi','spline') [i,j]=find(cz1==max(max(cz1))) x=xi(i),y=yi(j),zmax=cz1(i,j) 1.7.2 插值节点为散乱节点 已知n 个节点：),,2,1)(,,(n i z y x i i i L =，求点),(y x 处的插值z 。 对上述问题，Matlab 中提供了插值函数griddata ，其格式为： ZI = GRIDDATA(X,Y ,Z,XI,YI) 其中X 、Y 、Z 均为n 维向量，指明所给数据点的横坐标、纵坐标和竖坐标。向量XI 、YI 是给定的网格点的横坐标和纵坐标，返回值ZI 为网格（XI ，YI ）处的函数值。XI 与YI 应是方向不同的向量，即一个是行向量，另一个是列向量。 例3 在某海域测得一些点（x,y ）处的水深z 由下表给出，在矩形区域（75,200）×(-50,150) 内画出海底曲面的图形。 表3 x 129 140 103.5 88 185.5195 105157.5107.57781 162 162 117.5y 7.5 141.5 23 147 22.5137.585.5-6.5-81 3 56.5-66.5 84 -33.5z 4 8 6 8 6 8 8 9 9 8 8 9 4 9 解 编写程序如下： x=[129 140 103.5 88 185.5 195 105 157.5 107.5 77 81 162 162 117.5]; y=[7.5 141.5 23 147 22.5 137.5 85.5 -6.5 -81 3 56.5 -66.5 84 -33.5]; z=-[4 8 6 8 6 8 8 9 9 8 8 9 4 9]; xi=75:1:200; yi=-50:1:150; zi=griddata(x,y,z,xi,yi','cubic') subplot(1,2,1), plot(x,y,'*') subplot(1,2,2), mesh(xi,yi,zi) §2 曲线拟合的线性最小二乘法 2.1 线性最小二乘法 曲线拟合问题的提法是，已知一组（二维）数据，即平面上的n 个点),(i i y x ， n i ,,2,1L =，i x 互不相同，寻求一个函数（曲线）)(x f y =，使)(x f 在某种准则下 与所有数据点最为接近，即曲线拟合得最好。 -190- 线性最小二乘法是解决曲线拟合最常用的方法，基本思路是，令 )()()()(2211x r a x r a x r a x f m m +++=L （11） 其中)(x r k 是事先选定的一组线性无关的函数，k a 是待定系数),,,2,1(n m m k <=L 。拟合准则是使i y ，n i ,,2,1L =， 与)(i x f 的距离i δ的平方和最小，称为最小二乘准则。 2.1.1 系数k a 的确定 记 21 1 2 1])([),,(i i n i n i i m y x f a a J ?==∑∑==δL （12） 为求m a a ,,1L 使J 达到最小，只需利用极值的必要条件0=??k a J ),,1(m k L =，得到关于m a a ,,1L 的线性方程组 ),,1(, 0])()[(1 1 m j y x r a x r n i m k i i k k i j L ==?∑∑== 即 ),,1(, )()]()([1 1 1 m j y x r x r x r a m k n i i n i i j i k i j k L ==∑∑∑=== （13） 记 m n n m n m x r x r x r x r R ×? ???? ?? ???=)()()()(1111L M M M L ， T m a a A ],,[1L =，T n y y Y ),,(1L = 方程组（13）可表为 Y R RA R T T = （14） 当)}(,),({1x r x r m L 线性无关时，R 列满秩，R R T 可逆，于是方程组（14）有唯一解 Y R R R A T T 1 )(?= 2.1.2 函数)(x r k 的选取 面对一组数据n i y x i i ,,2,1),,(L =，用线性最小二乘法作曲线拟合时，首要的、也是关键的一步是恰当地选取)(,),(1x r x r m L 。如果通过机理分析，能够知道y 与x 之间应该有什么样的函数关系，则)(,),(1x r x r m L 容易确定。若无法知道y 与x 之间的关系，通常可以将数据n i y x i i ,,2,1),,(L =作图，直观地判断应该用什么样的曲线去作拟合。人们常用的曲线有 （i ）直线 21a x a y += （ii ）多项式 11++++=m m m a x a x a y L （一般3,2=m ，不宜太高） （iii ）双曲线（一支） 21 a x a y += -191- （iv ）指数曲线 x a e a y 21= 对于指数曲线，拟合前需作变量代换，化为对21,a a 的线性函数。 已知一组数据，用什么样的曲线拟合最好，可以在直观判断的基础上，选几种曲线分别拟合，然后比较，看哪条曲线的最小二乘指标J 最小。 2.2 最小二乘法的Matlab 实现 2.2.1 解方程组方法 在上面的记号下， 21||||),,(Y RA a a J m ?=L Matlab 中的线性最小二乘的标准型为 ,min 22 Y RA A ? 命令为 Y R A \=。 例4 用最小二乘法求一个形如2 bx a y +=的经验公式，使它与表4所示的数据拟合。 表4 x 19 25 31 38 44 y 19.0 32.3 49.0 73.3 97.8 解 编写程序如下 x=[19 25 31 38 44]'; y=[19.0 32.3 49.0 73.3 97.8]'; r=[ones(5,1),x.^2]; ab=r\y x0=19:0.1:44; y0=ab(1)+ab(2)*x0.^2; plot(x,y,'o',x0,y0,'r') 2.2.2 多项式拟合方法 如果取},,,1{)}(,),({11m m x x x r x r L L =+，即用m 次多项式拟合给定数据，Matlab 中有现成的函数 a=polyfit(x0,y0,m) 其中输入参数x0,y0为要拟合的数据，m 为拟合多项式的次数，输出参数a 为拟合多项式y=a m x m +…+a 1x+a 0系数a=[ a m , …, a 1, a 0]。 多项式在x 处的值y 可用下面的函数计算 y=polyval(a,x)。 例5 某乡镇企业1990-1996年的生产利润如表5。 表5 年份 1990 1991 1992 1993 1994 1995 1996 利润（万元） 70 122 144 152 174 196 202 试预测1997年和1998年的利润。 解 作已知数据的的散点图， x0=[1990 1991 1992 1993 1994 1995 1996]; y0=[70 122 144 152 174 196 202]; plot(x0,y0,'*') -192- 发现该乡镇企业的年生产利润几乎直线上升。因此，我们可以用01a x a y +=作为拟合函数来预测该乡镇企业未来的年利润。编写程序如下： x0=[1990 1991 1992 1993 1994 1995 1996]; y0=[70 122 144 152 174 196 202]; a=polyfit(x0,y0,1) y97=polyval(a,1997) y98=polyval(a,1998) 求得201=a ，4 01007054×?=.a ，1997年的生产利润y97=233.4286，1998年的生产利润y98=253.9286。 §3 最小二乘优化 在无约束最优化问题中，有些重要的特殊情形，比如目标函数由若干个函数的平方和构成。这类函数一般可以写成： ∑== m i i x f x F 1 2)()(,n R x ∈ 其中T n x x x ),,(1L =，一般假设n m ≥。我们把极小化这类函数的问题： ∑==m i i x f x F 1 2)()(min 称为最小二乘优化问题。 最小二乘优化是一类比较特殊的优化问题，在处理这类问题时，Matlab 也提供了一些强大的函数。在Matlab 优化工具箱中，用于求解最小二乘优化问题的函数有：lsqlin 、lsqcurvefit 、lsqnonlin 、lsqnonneg ，用法介绍如下。 3.1 lsqlin 函数 求解 2 22 1min d Cx x ? s.t. ?? ? ??≤≤=≤ub x lb beq x Aeq b x A ** 其中Aeq A C ,,为矩阵，x ub lb beq b d ,,,,,为向量。 Matlab 中的函数为： x=lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0) 例6 用lsqlin 命令求解例4。 解 编写程序如下： x=[19 25 31 38 44]'; y=[19.0 32.3 49.0 73.3 97.8]'; r=[ones(5,1),x.^2]; ab=lsqlin(r,y) x0=19:0.1:44; y0=ab(1)+ab(2)*x0.^2; plot(x,y,'o',x0,y0,'r') -193- 3.2 lsqcurvefit 函数 给定输入输出数列ydata xdata ,,求参量x ,使得 ()2 22 ),(21),(2 1 min ∑?= ?i i i x ydata xdata x F ydata xdata x F Matlab 中的函数为 X=LSQCURVEFIT(FUN,X0,XDATA,YDATA,LB,UB,OPTIONS) 其中FUN 是定义函数),(xdata x F 的M 文件。 例7 用下面表6中的数据拟合函数kt be a t c 02.0)(?+=中的参数k b a ,,。 表6 j t 100 200 300 400 500 600 700 800 900 1000 j c 4.54 4.99 5.35 5.65 5.90 6.10 6.26 6.39 6.50 6.59 解 该问题即解最优化问题： ∑=??+= 10 1 202.0)(),,(min i j kt c be a k b a F j （1）编写M 文件fun1.m 定义函数),(tdata x F ： function f=fun1(x,tdata); f=x(1)+x(2)*exp(-0.02*x(3)*tdata); %其中x(1)=a ，x(2)=b ，x(3)=k （2）调用函数lsqcurvefit，编写程序如下： td=100:100:1000; cd=[4.54 4.99 5.35 5.65 5.90 6.10 6.26 6.39 6.50 6.59]; x0=[0.2 0.05 0.05]; x=lsqcurvefit(@fun1,x0,td,cd) 3.3 lsqnonlin 函数 已知函数向量[]T k x f x f x F )(,),()(1L =，求x 使得 2 2)(2 1min x F x Matlab 中的函数为 X=LSQNONLIN(FUN,X0,LB,UB,OPTIONS) 其中FUN 是定义向量函数)(x F 的M 文件。 例8 用lsqnonlin 函数求解例7。 解 这里 [ ] T kt kt c be a c be a t x F x F 1002.0102.0101 ,,),()(?+?+==??L []k b a x ,,= （1）编写M 文件fun2.m 如下： function f=fun2(x); td=100:100:1000; cd=[4.54 4.99 5.35 5.65 5.90 6.10 6.26 6.39 6.50 6.59]; f=x(1)+x(2)*exp(-0.02*x(3)*td)-cd; -194- （2）调用函数lsqnonlin，编写程序如下： x0=[0.2 0.05 0.05]; %初始值是任意取的 x=lsqnonlin(@fun2,x0) 3.4 lsqnonneg 函数 求解非负的x ，使得满足 2 22 1min d Cx x ?， Matlab 中的函数为 X = LSQNONNEG(C,d,X0,OPTIONS) 例9 已知?? ????? ?? ???=6170.06344.06245.06233.07071.06861 .02869.00372.0C ，?????? ????? ?=8405.00747.01781.08587.0d ，求)0(≥x x 满足 2 22 1min d Cx x ? 最小。 解 编写程序如下： c=[0.0372 0.2869;0.6861 0.7071;0.6233 0.6245;0.6344 0.6170]; d=[0.8587;0.1781;0.0747;0.8405]; x=lsqnonneg(c,d) 3.5 曲线拟合的用户图形界面求法 Matlab 工具箱提供了命令cftool ，该命令给出了一维数据拟合的交互式环境。具体执行步骤如下： （1）把数据导入到工作空间； （2）运行cftool ，打开用户图形界面窗口； （3）对数据进行预处理； （4）选择适当的模型进行拟合； （5）生成一些相关的统计量，并进行预测。 可以通过帮助（运行doc cftool ）熟悉该命令的使用细节。 §4 曲线拟合与函数逼近 前面讲的曲线拟合是已知一组离散数据},,1),,{(n i y x i i L =，选择一个较简单的函数)(x f ，如多项式，在一定准则如最小二乘准则下，最接近这些数据。 如果已知一个较为复杂的连续函数],[),(b a x x y ∈，要求选择一个较简单的函数 )(x f ，在一定准则下最接近)(x f ，就是所谓函数逼近。 与曲线拟合的最小二乘准则相对应，函数逼近常用的一种准则是最小平方逼近，即 ∫?=b a dx x y x f J 2)]()([ (15) 达到最小。与曲线拟合一样，选一组函数},,1),({m k x r k L =构造)(x f ，即令 )()()(11x r a x r a x f m m ++=L 插值和拟合的定义 1.定义：若x为自变量，y为因变量，则x与y之间有一个确定的函数表达式f（x），现实中，这个函数关系式很难确定，运用逼近的方法处理：取得一组数据点（xi,yi,i=1,2,3...n）,构造一个简单函数p(x)作为f(x)的近似表达式，且p(x)满足： p(xi)=f(xi)=yi i=1,2,3,4...n 这就是插值问题。 若不要求p(x)通过所有数据点，而是要求曲线在某种准则下整体与数据点接近，例如运用最小二乘法得到p(x),这种问题称为拟合。 2插值类型：一维插值是指对一维函数y=f(x)进行插值，二维插值就是对二维函数y=f(x,y)进行插值. 3.插值的matlab函数及其应用 （1）.一维插值：yi=interp1(x,Y,xi，’method’)，对一组节点（x,Y）进行插值，计算插值点xi 的函数值。method包括了一下几种类型： Nearest：线性最邻近插值（速度最快，平滑性最差） Linear：线性插值（默认项）（生成效果连续，但是顶点处有坡度变化） Spline：三次样条插值（运行时间较长，插值数据和导数均连续） Pchip：分段三次艾米尔特（hermite）插值 Cubic：双三次插值（较高版本的matlab不能运用，v5cubic能够运行） 运行的代码及插值效果： clear; clc; x=0:0.2:2; y=(x.^2-3*x+5).*exp(-3*x).*sin(x); xi=0:0.03:2; yi_nearest=interp1(x,y,xi,'nearest'); yi_linear=interp1(x,y,xi,'linear'); yi_spline=interp1(x,y,xi,'spline'); yi_pchip=interp1(x,y,xi,'pchip'); yi_cubic=interp1(x,y,xi,'v5cubic'); figure; hold on; subplot(2 ,3, 1); plot(x,y,'r*'); title('已知数据值'); subplot(2, 3 ,2); plot(x,y,'r*',xi,yi_nearest,'b-'); title('最邻近插值'); 插值法和拟合实验报告 一、 实验目的 1.通过进行不同类型的插值，比较各种插值的效果，明确各种插值的优越性； 2.通过比较不同次数的多项式拟合效果，了解多项式拟合的原理； 3.利用matlab 编程，学会matlab 命令； 4.掌握拉格朗日插值法； 5.掌握多项式拟合的特点和方法。 二、 实验题目 1.、插值法实验 将区间[-5,5]10等分，对下列函数分别计算插值节点 k x 的值，进行不同类型 的插值，作出插值函数的图形并与)(x f y =的图形进行比较： ;11)(2x x f += ;a r c t a n )(x x f = .1)(42 x x x f += (1) 做拉格朗日插值； (2) 做分段线性插值； (3) 做三次样条插值. 2、拟合实验 给定数据点如下表所示： 分别对上述数据作三次多项式和五次多项式拟合，并求平方误差，作出离散函数 ),(i i y x 和拟合函数的图形。 三、 实验原理 1.、插值法实验 ∏∑∏∏∏∑∑≠==≠=≠=≠=+-==--= =-= ==-=-=----==++==j i j j i i i i i n i i n n j i j j n j i j j i i n j i j j n i i i n i i n n n o i n i i n x x x x x y x l x L x x c n i x x c x x x c x x x x x x x x c y x l x L y x l y x l y x l x L ,00 ,0,0,01100 00 )(l )()() (1 ,1,0, 1)()(l ) ()())(()()()()()()()(， 故， 得 再由，设 2、拟合实验 习题5.1： Matlab程序如下： clc,clear x=1:0.5:10; y=x.^3-6*x.^2+5*x-3; y0=y+rand; f1=polyfit(x,y0,1) y1=polyval(f1,x); plot(x,y,'+',x,y1); grid on title('一次拟合曲线'); figure(2); f2=polyfit(x,y0,2) y2=polyval(f2,x); plot(x,y,'+',x,y2); grid on title('二次拟合曲线'); figure(3); f4=polyfit(x,y0,4) y3=polyval(f4,x); plot(x,y,'+',x,y3); grid on title('四次拟合曲线'); figure(4); f6=polyfit(x,y0,6) y4=polyval(f6,x); plot(x,y,'+',x,y4); grid on title('六次拟合曲线'); 计算结果及图如下 f1 = 43.2000 -148.8307 f2 = 10.5000 -72.3000 90.0443 f4 = 0.0000 1.0000 -6.0000 5.0000 -2.3557 f6 = -0.0000 0.0000 -0.0000 1.0000 -6.0000 5.0000 -2.3557 5.2高程数据问题解答如下：matlab程序: clc,clear x0=0:400:5600 y0=0:400:4800 z0=[1350 1370 1390 1400 1410 960 940 880 800 690 570 430 290 210 150 1370 1390 1410 1430 1440 1140 1110 1050 950 820 690 540 380 300 210 1380 1410 1430 1450 1470 1320 1280 1200 1080 940 780 620 460 370 350 1420 1430 1450 1480 1500 1550 1510 1430 1300 1200 980 850 750 550 500 利用matlab实现插值与拟合实验 张体强1026222 张影 晁亚敏 [摘要]：在测绘学中，无论是图形处理，还是地形图处理等，大多离不开插值与拟合的应用，根据插值与拟合原理，构造出插值和拟合函数，理解其原理，并在matlab平台下，实现一维插值，二维插值运算，实现多项式拟合，非线性拟合等，并在此基础上，联系自己所学专业，分析其生活中特殊例子，提出问题，建立模型，编写程序，以至于深刻理解插值与拟合的作用。 [关键字]: 测绘学插值多项式拟合非线性拟合 [ Abstract]: in surveying and mapping, whether the graphics processing, or topographic map processing and so on, are inseparable from the interpolation and fitting application, according to the interpolation and fitting theory, construct the fitting and interpolation function, understanding its principle, and MATLAB platform, achieve one-dimensional interpolation, two-dimensional interpolation, polynomial fitting, non-linear fitting, and on this basis, to contact their studies, analysis of their living in a special example, put forward the question, modeling, programming, so that a deep understanding of interpolation and fitting function. [ Key words]: Surveying and mapping interpolation polynomial fitting nonlinear 实验三 插值法与拟合实验 一、实验目的 1. 通过本实验学会利用程序画出插值函数，并和原图形相比较 2. 通过本实验学会拟合函数图形的画法，并会求平方误差 二、实验题目 1. 插值效果的比较 实验题目：区间[]5,5-10等分，对下列函数分别计算插值节点k x 的值，进行不同类型的插值，作出插值函数的图形并与)(x f y =的图形进行比较： 2 11)(x x f +=； x x f arctan )(=； 4 41)(x x x f += （1） 做拉格朗日插值； （2） 做三次样条插值. 2. 拟合多项式实验 实验题目：给定数据点如下表所示： 分别对上述数据作三次多项式和五次多项式拟合，并求平方误差，作出离散函数),(i i y x 和拟合函数的图形. 三、实验原理 本实验应用了拉格朗日插值程序、三次样条插值程序、多项式拟合程序等实验原理. 四、实验内容 1（1） figure x=-5:0.2:5; y=1./(1+x.^2); plot(x,y,'r'); hold on %拉格朗日插值 x1=-5:1:5; y1=1./(1+x1.^2); xx=-4.5:0.5:4.5; yy=malagr(x1,y1,xx); plot(xx,yy,'+') %三次样条插值 dy0=1./(1+25); dyn=1./(1+25); m=maspline(x1,y1,dy0,dyn,xx); plot(xx,m,'ok') 1（2） x=-5:0.2:5; y=atan(x); plot(x,y,'r'); hold on %拉格朗日插值 x1=-5:1:5; y1=atan(x1); xx=-4.5:0.5:4.5; yy=malagr(x1,y1,xx); plot(xx,yy,'+') %三次样条插值 dy0=1./(1+25); dyn=1./(1+25); m=maspline(x1,y1,dy0,dyn,xx); plot(xx,m,'ok') 1（3） x=-5:0.2:5; y=x.^2./(1+x.^4); plot(x,y,'r'); hold on %拉格朗日插值 x1=-5:1:5; y1=x1.^2./(1+x1.^4); xx=-4.5:0.5:4.5; yy=malagr(x1,y1,xx); plot(xx,yy,'+') %三次样条插值 dy0=1./(1+25); dyn=1./(1+25); m=maspline(x1,y1,dy0,dyn,xx); plot(xx,m,'ok') 2. x=[-1.5 -1.0 -0.5 0.0 0.5 1.0 1.5]'; y=[-4.45 -0.45 0.55 0.05 -0.44 0.54 4.55]'; plot(x,y,'or'); hold on %三次多项式拟合 p1=mafit(x,y,3); MATLAB中的插值、拟合与查表 插值法是实用的数值方法，是函数逼近的重要方法。在生产和科学实验中，自变量x与因变量y的函数y = f(x)的关系式有时不能直接写出表达式，而只能得到函数在若干个点的函数值或导数值。当要求知道观测点之外的函数值时，需要估计函数值在该点的值。 如何根据观测点的值，构造一个比较简单的函数y=φ(x)，使函数在观测点的值等于已知的数值或导数值。用简单函数y=φ(x)在点x处的值来估计未知函数y=f(x)在x点的值。寻找这样的函数φ(x)，办法是很多的。φ(x)可以是一个代数多项式，或是三角多项式，也可以是有理分式；φ(x)可以是任意光滑（任意阶导数连续）的函数或是分段函数。函数类的不同，自然地有不同的逼近效果。在许多应用中，通常要用一个解析函数（一、二元函数）来描述观测数据。 根据测量数据的类型： 1．测量值是准确的，没有误差。 2．测量值与真实值有误差。 这时对应地有两种处理观测数据方法： 1．插值或曲线拟合。 2．回归分析（假定数据测量是精确时，一般用插值法，否则用曲线拟合）。 MATLAB中提供了众多的数据处理命令。有插值命令，有拟合命令，有查表命令。 2.2.1 插值命令 命令1 interp1 功能一维数据插值（表格查找）。该命令对数据点之间计算内插值。它找出一元函数f(x)在中间点的数值。其中函数f(x)由所给数据决定。各个参量之间的关系示意图为图2-14。 格式 yi = interp1(x,Y,xi) %返回插值向量yi，每一元素对应于参量xi，同时由向量x 与Y的内插值决定。参量x指定数据Y的点。若Y为一矩阵，则按Y的每列计算。yi是阶数为length(xi)*size(Y,2)的输出矩阵。 yi = interp1(Y,xi) %假定x=1:N，其中N为向量Y的长度，或者为矩阵Y的行数。 yi = interp1(x,Y,xi,method) %用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； 数值计算之美SHU ZHI JI SUAN ZHI MEI 胡家威 http://hujiaweibujidao.github.io/ 清华大学逸夫图书馆·北京 内容简介 本书是我对数值计算中的若干常见的重要算法及其应用的总结，内容还算比较完整。 本人才疏学浅，再加上时间和精力有限，所以本书不会详细介绍很多的概念，需要读者有一定的基础或者有其他的参考书籍，这里推荐参考文献中的几本关于数值计算的教材。 本书只会简单介绍下算法的原理，对于每个算法都会附上我阅读过的较好的参考资料以及算法的实现(Matlab或者其他语言)，大部分代码是来源于参考文献[1]或者是经过我改编而成的，肯定都是可以直接使用的，需要注意的是由于Latex对代码的排版问题，导致中文注释中的英文字符经常出现错位，对于这种情况请读者自行分析，不便之处还望谅解。写下这些内容的目的是让自己理解地更加深刻些，顺便能够作为自己的HandBook，如有错误之处，还请您指正，本人邮箱地址是：hujiawei090807@https://www.doczj.com/doc/5212708916.html,。 目录 第四章曲线拟合和多项式插值1 4.1曲线拟合 (1) 4.1.1使用线性方程进行曲线拟合 (1) 4.1.2非线性方程进行曲线拟合 (2) 4.1.3使用二次或者高次多项式进行曲线拟合[最小二乘问题].3 4.2多项式插值 (4) 4.2.1拉格朗日插值多项式 (4) 4.2.2牛顿插值多项式 (5) 4.2.3分段线性插值 (7) 4.2.4保形分段三次插值 (8) 4.2.5三次样条插值 (10) 4.3Matlab函数解析 (13) 参考文献14 第十章 插值与拟合方法建模 在生产实际中，常常要处理由实验或测量所得到的一批离散数据，插值与拟合方法就是要通过这些数据去确定某一类已经函数的参数，或寻求某个近似函数使之与已知数据有较高的拟合精度。插值与拟合的方法很多，这里主要介绍线性插值方法、多项式插值方法和样条插值方法，以及最小二乘拟合方法在实际问题中的应用。相应的理论和算法是数值分析的内容，这里不作详细介绍，请参阅有关的书籍。 §1 数据插值方法及应用 在生产实践和科学研究中，常常有这样的问题：由实验或测量得到变量间的一批离散样点，要求由此建立变量之间的函数关系或得到样点之外的数据。与此有关的一类问题是当原始数据 ),(,),,(),,(1100n n y x y x y x 精度较高，要求确定一个初等函数)(x P y =（一般用多项式或分段 多项式函数）通过已知各数据点（节点），即n i x P y i i ,,1,0,)( ==，或要求得函数在另外一些点（插值点）处的数值，这便是插值问题。 1、分段线性插值 这是最通俗的一种方法，直观上就是将各数据点用折线连接起来。如果 b x x x a n =<<<= 10 那么分段线性插值公式为 n i x x x y x x x x y x x x x x P i i i i i i i i i i ,,2,1,,)(11 1 11 =≤<--+--= ----- 可以证明，当分点足够细时，分段线性插值是收敛的。其缺点是不能形成一条光滑曲线。 例1、已知欧洲一个国家的地图，为了算出它的国土面积，对地图作了如下测量：以由西向东方向为x 轴，由南向北方向为y 轴，选择方便的原点，并将从最西边界点到最东边界点在x 轴上的区间适当的分为若干段，在每个分点的y 方向测出南边界点和北边界点的y 坐标y1和y2，这样就得到下表的数据（单位：mm ）。 插值和拟合区别 运输1203黎文皓通过这个学期的《科学计算与数学建模》课程的学习，使我掌握了不少数学模型解决实际问题的方法，其中我对于插值与拟合算法这一章，谈一谈自己的看法可能不是很到位，讲得不好的地方也请老师见谅。 首先，举一个简单的例子说明一下这个问题。 如果有100个平面点，要求一条曲线近似经过这些点，可有两种方法：插值和拟合。 我们可能倾向于用一条（或者分段的多条）2次、3次或者说“低次”的多项式曲线而不是99次的曲线去做插值。也就是说这条插值曲线只经过其中的3个、4个（或者一组稀疏的数据点）点，这就涉及到“滤波”或者其他数学方法，也就是把不需要90多个点筛选掉。如果用拟合，以最小二乘拟合为例，可以求出一条（或者分段的多条）低次的曲线（次数自己规定），逼近这些数据点。具体方法参见《数值分析》中的“线性方程组的解法”中的“超定方程的求解法”。经过上面例子的分析，我们可以大致的得到这样一个结论。插值就是精确经过，拟合就是逼近。 插值和拟合都是函数逼近或者数值逼近的重要组成部分。他们的共同点都是通过已知一些离散点集M上的约束，求取一个定义在连续集合S(M包含于S)的未知连续函数，从而达到获取整体规律目的，即通过"窥几斑"来达到"知全豹"。 所谓拟合是指已知某函数的若干离散函数值{f1,f2,…,fn}，通过调 整该函数中若干待定系数f(λ1, λ2,…,λ3), 使得该函数与已知点集的差别(最小二乘意义)最小。如果待定函数是线性，就叫线性拟合或者线性回归(主要在统计中)，否则叫作非线性拟合或者非线性回归。表达式也可以是分段函数，这种情况下叫作样条拟合。 而插值是指已知某函数的在若干离散点上的函数值或者导数信息，通过求解该函数中待定形式的插值函数以及待定系数，使得该函数在给定离散点上满足约束。 从几何意义上将，拟合是给定了空间中的一些点，找到一个已知形式未知参数的连续曲面来最大限度地逼近这些点；而插值是找到一个(或几个分片光滑的)连续曲面来穿过这些点。 不过是插值还是拟合都是建立在一定的数学模型的基础上进行的。多项式插值虽然在一定程度上解决了由函数表求函数的近似表达式的问题，但是在逼近曲线上有明显的缺陷，很可能不能很好的表示函数的走向，存在偏差，在实际问题中我们往往通过函数近似表达式的拟合法来得到一个较为准却的表达式。 计算方法实验报告 专业班级：医学信息工程一班姓名：陈小芳学号：201612203501002 实验成绩： 1．【实验题目】 插值法与曲线拟合 2．【实验目的】 3．【实验内容】 4. 【实验要求】 5. 【源程序（带注释）】 （1）拉格朗日插值 #include 第四、五讲作业题参考答案 一、填空题 1、拉格朗日插值基函数在节点上的取值是（ 0或1 ）。 2、当1,1,2x =-，时()034f x =-， ，，则()f x 的二次插值多项式为 （ 2527 633 x x +- ）。 3、由下列数据 所确定的唯一插值多项式的次数为（ 2次 ）。 4、根据插值的定义，函数()x f x e -=在[0,1]上的近似一次多项式1()P x = （ 1(1)1e x --+ ），误差估计为（ 18 ）。 5、在做曲线拟合时，对于拟合函数x y ax b = +，引入变量变换y =（ 1 y ），x =（ 1 x ）来线性化数据点后，做线性拟合y a bx =+。 6、在做曲线拟合时，对于拟合函数Ax y Ce =，引入变量变换（ ln()Y y = ）、 X x =和B C e =来线性化数据点后，做线性拟合Y AX B =+。 7、设3()1f x x x =+-，则差商[0,1,2,3]f =（ 1 ）。 8、在做曲线拟合时，对于拟合函数()A f x Cx =，可使用变量变换（ln Y y =）（ln X x = ）和B C e =来线性化数据点后，做线性拟合Y AX B =+。 9、设(1)1,(0)0,(1)1,(2)5,()f f f f f x -====则的三次牛顿插值多项式为 （ 3211 66x x x +-），其误差估计式为（ 4()(1)(1)(2),(1,2)24f x x x x ξξ+--∈-） 10、三次样条插值函数()S x 满足：()S x 在区间[,]a b 内二阶连续可导， (),,0,1,2,,,k k k k S x y x y k n ==（已知）且满足()S x 在每一个子区间1[,] k k x x +上是（ 三次多项式 ）。 实验报告三 一、实验目的 通过本实验的学习，各种插值法的效果，如多项式插值法，牛顿插值法，样条插值法，最小二乘法拟合（即拟合插值），了解它们各自的优缺点及插值。 二、实验题目 1、 插值效果比较 实验题目：将区间[]5,5-10等份，对下列函数分别计算插值节点k x 的值，进行不同类型的插值，作出插值函数的图形并与)(x f y =的图形进行比较： 211)(x x f +=；x x f arctan )(=；4 2 1)(x x x f +=。 （1） 做拉格朗日插值； （2） 做三次样条插值。 2、 拟合多项式实验 实验题目：给定数据点如下表所示： i x -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 i y -4.45 -0.45 0.55 0.05 -0.44 0.54 4.55 分别对上述数据做三次多项式和五次多项式拟合，并求平方误差，作出离散函数()i i y x ,和拟合函数的图形。 三、实验原理 n 阶拉格朗日插值 设已知n x x x ,,,10 及()()()x L n i x f y n i i ,,,1,0 ==为不超过n 次的多项式，且满足 插值条件()().,,1,0n i y x L i i n ==由对()x L 2的构造经验，可设 ()()()()(),11000 n n n i i i n y x l y x l y x l y x l x L +++==∑= 其中，()()n i x L i ,,1,0 =均为n 次多项式且满足() .,,1,0,, ,0, ,1n j i j i j i x l j i =?? ?≠==不难验 证，这样构造出的()x L n 满足插值条件。因此问题归结为求()()n i x l i ,,1,0 =的表达式。因 ()i j x i ≠是n 次多项式()x l i 的n 个根，故可设 第七讲插值方法与数据拟合 § 7.1 引言 在工程和科学实验中，常常需要从一组实验观测数据(x i , y i ) (i= 1, 2, …, n) 揭示自变量x与因变量y 之间的关系，一般可以用一个近似的函数关系式y = f (x) 来表示。函数f (x) 的产生办法因观测数据与要求的不同而异，通常可采用两种方法：插值与数据拟合。 § 7.1.1 插值方法 1．引例1 已经测得在北纬32.3?海洋不同深度处的温度如下表： 根据这些数据，我们希望能合理地估计出其它深度（如500米、600米、1000米…）处的水温。 解决这个问题，可以通过构造一个与给定数据相适应的函数来解决，这是一个被称为插值的问题。 2．插值问题的基本提法 对于给定的函数表 其中f (x) 在区间[a, b] 上连续，x0，x1，…，x n为[a, b] 上n + 1个互不相同的点，要求在一个性质优良、便于计算的函数类{P(x)} 中，选出一个使 P(x i ) = y i，i= 0, 1, …, n(7.1.1) 成立的函数P(x) 作为 f (x) 的近似，这就是最基本的插值问题（见图7.1.1）。 为便于叙述，通常称区间[a, b] 为插值区间，称点x0，x1，…，x n为插值节点，称函数类{P(x)} 为插值函数类，称式(7.1.1) 为插值条件，称函数P(x) 为插值函数，称f (x) 为被插函数。求插值函数P(x) 的方法称为插值法。 § 7.1.2 数据拟合 1．引例2 在某化学反应中，已知生成物的浓度与时间有关。今测得一组数据如下： 根据这些数据，我们希望寻找一个y = f (t) 的近似表达式（如建立浓度y与时间t之间的经验公式等）。从几何上看，就是希望根据给定的一组点（1, 4.00），…,（16, 10.60），求函数y = f (t) 的图象的一条拟合曲 实验 2 插值与拟合 系班姓名学号 【实验目的】 1、掌握用MATLAB计算拉格朗日、分段线性、三次样条三种插值的方法，改变节点的数目， 对三种插值结果进行初步分析。 2、掌握用MATLAB作线性最小二乘的方法。 3、通过实例学习如何用插值方法与拟合方法解决实际问题，注意二者的联系和区别。 【实验内容】 预备：编制计算拉格朗日插值的M文件： 以下是拉格朗日插值的名为y_lagrl的M文件： function y=y_lagr1(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 第1题（d） 选择函数y=exp(-x2) (-2≤x≤2),在n个节点上（n不要太大，如5~11）用拉格朗日、分段线性、三次样条三种插值方法，计算m个插值点的函数值（m要适中，如50~100）。通过数值和图形输出，将三种插值结果与精确值进行比较。适当增加n，在作比较，由此作初步分析。 运行如下程序： n=7;m=61;x=-2:4/(m-1):2; y=exp(-x.^2); z=0*x; x0=-2:4/(n-1):2; y0=exp(-x0.^2); y1=y_lagr1(x0,y0,x); y2=interp1(x0,y0,x); y3=interp1(x0,y0,x,'spline'); [x'y'y1'y2'y3'] plot(x,z,'w',x,y,'r--',x,y1,'b:',x,y2,'m',x,y3,'b') gtext('y=exp(-x^2)'),gtext('Lagr.'),gtext('Piece.-linear.'),gtext ('Spline'), 将三种插值结果y1,y2,y3与精确值y 项比较，显然y1在节点处不光滑,拉格朗日插值出现较大的振荡,样条插值得结果是最好的.增加n 值（使n=11），再运行以上程序，得到的图形如右图所示，比较这两个图可发现，节点增加后，三种插值方法结果的准确度均有所提高，因此可近似地认为：增加节点个数可以提高插值结果的准确程度。 第3题 用给定的多项式，如y=x 3-6x 2+5x-3，产生一组数据（x i ,y i ,i=1,2,…,n ），再在yi 上添加随机干扰（可用rand 产生（0，1）均匀分布随机数，或用randn 产生N （0，1）分布随机数），然后用x i 和添加了随机干扰的y i 作3次多项式拟合，与原系数比较。如果作2或4次多项式拟合，结果如何？ 解：2 编制y_2_3.m 文件 n=15; x=0:8/(n-1):8; y=x.^3-6*x.^2+5*x-3; z=0*x; y0=y+rand(1,15); f=polyfit(x,y0,m); r=polyval(f,x) pl2ot(x,z,'k',x,y,'r:'r,'b') 程序及运行结果如下：m=2 ,y_2_3 f = 5.9888 -31.9916 17.6679 m=3 ,y_2_3 1 山区地貌：在某山区测得一些地点的高程如下表：(平面区域1200<=x<=4000,1200<=y<=3600)，试作出该山区的地貌图和等高线图，并对几种插值方法进行比较。 2 用给定的多项式，如y=x3-6x2+5x-3，产生一组数据(xi,yi，i=1,2,…,n),再在yi上添加随机干扰(可用rand产生(0,1)均匀分布随机数,或用rands产生N(0,1)分布随机数)，然后用xi 和添加了随机干扰的yi作的3次多项式拟合，与原系数比较。如果作2或4次多项式拟合，结果如何？ 3 用电压V=10伏的电池给电容器充电，电容器上t时刻的电压为，其中V0是电容器的初始电压，是充电常数。试由下面一组t，V数据确定V0，。 2用给定的多项式，如y=x3-6x2+5x-3，产生一组数据(xi,yi，i=1,2,…,n),再在yi上添加随机干扰(可用rand产生(0,1)均匀分布随机数,或用rands产生N(0,1)分布随机数)，然后用xi和添加了随机干扰的yi作的3次多项式拟合，与原系数比较。 分别作1、2、4、6次多项式拟合，比较结果，体会欠拟合、过拟合现象。 解：程序如下： x=1:0.5:10; y=x.^3-6*x.^2+5*x-3; y0=y+rand; f1=polyfit(x,y0,1)%输出多项式系数 y1=polyval(f1,x);%计算各x点的拟合值 plot(x,y,'+',x,y1) grid on title('一次拟合曲线'); figure(2); f2=polyfit(x,y0,2)%2次多项式拟合 y2=polyval(f2,x); plot(x,y,'+',x,y2); 《数值分析》课程实验一：插值与拟合 一、实验目的 1. 理解插值的基本原理，掌握多项式插值的概念、存在唯一性； 2. 编写MA TLAB 程序实现Lagrange 插值和Newton 插值，验证Runge 现象； 3. 通过比较不同次数的多项式拟合效果，理解多项式拟合的基本原理； 4. 编写MA TLAB 程序实现最小二乘多项式曲线拟合。 二、实验内容 1. 用Lagrange 插值和Newton 插值找经过点(-3, -1), (0, 2), (3, -2), (6, 10)的三次插值公式，并编写MATLAB 程序绘制出三次插值公式的图形。 2. 设 ]5,5[,11 )(2 -∈+= x x x f 如果用等距节点x i = -5 + 10i /n (i = 0, 1, 2, …, n )上的Lagrange 插值多项式L n (x )去逼近它。不妨取n = 5和n = 10，编写MATLAB 程序绘制出L 5(x )和L 10(x )的图像。 (2) 编写MA TLAB 程序绘制出曲线拟合图。 三、实验步骤 1. (1) Lagrange 插值法：在线性空间P n 中找到满足条件： ?? ?≠===j i j i x l ij j i , 0, , 1)(δ 的一组基函数{}n i i x l 0)(=，l i (x )的表达式为 ∏ ≠==--= n i j j j i j i n i x x x x x l ,0),,1,0()( 有了基函数{}n i i x l 0)(=，n 次插值多项式就可表示为 ∑==n i i i n x l y x L 0 )()( (2) Newton 插值法：设x 0, x 1, …, x n 是一组互异的节点，y i = f (x i ) (i = 0, 1, 2, …, n )，f (x )在处的n 阶差商定义为 计算方法实验 题目： 班级： 学号： 姓名： 目录 计算方法实验 (1) 1 实验目的 (3) 2 实验步骤 (3) 2.1环境配置： (3) 2.2添加头文件 (3) 2.3主要模块 (3) 3 代码 (4) 3.1主程序部分 (4) 3.2多项式方程部分 (4) 3.3核心算法部分 (8) 3.4数据结构部分 (13) 4运行结果 (19) 4.1拉格朗日插值法运行结果 (19) 4.2牛顿插值法运行结果 (20) 4.3多项式拟合运行结果 (20) 5总结 (21) 拉格朗日插值法 (21) 牛顿插值法 (21) 多项式拟合 (21) 6参考资料 (22) 1 实验目的 1.通过编程对拉格朗日插值法、牛顿插值法以及多项式拟合数据的理解 2.观察上述方法的计算稳定性和求解精度并比较各种方法利弊 2 实验步骤 2.1环境配置： VS2013，C++控制台程序 2.2添加头文件 #include "stdio.h" #include "stdlib.h" #include "stdafx.h" 2.3主要模块 程序一共分成三层，最底层是数据结构部分，负责存储数据，第二层是交互部分，即多项式方程部分，负责输入输出获得数据，最上层是核心的算法部分，负责处理已获得的数据。具体功能如下： ●数据结构部分 数据结构部分是整个程序的最底层，负责存储部分。因方程系数作为数据元素插入和删除操作较少，而顺序表空间利用率大且查看方便，故此程序选用顺序表保存系数。数据结构文件中写的是有关顺序表的所有基本操作以供其他文件调用。本次实验使用列主元高斯消元法作为求解方程组的方法，所以也用了二维顺序表存储数组。综上，数据结构部分文件是前两个试验的文件内容和，稍作修改。 ●常系数微分方程部分 多项式方程部分是程序的第二层，内容主要是常系数微分方程导数的计算和显示菜单部分。 ●算法部分 算法部分分为两个文件，一个是插值部分，一个是拟合部分。 插值部分文件负责有关插值的核心算法，处于整个程序最上层部分，负责拉格朗日插值法和牛顿插值法的具体实现过程。调用方程文件的函数，将获得的数据进行处理运算，将结果返回给方程主函数和输出的第二层。每种方法有两个函数，一个为仅仅实现一次插值的算法，另一个是和方程部分联系的 Matlab中插值拟合函数汇总和使用说明 命令1 interp1 功能一维数据插值（表格查找）。该命令对数据点之间计算内插值。它找出一元函数f(x)在中间点的数值。其中函数f(x)由所给数据决定。 x：原始数据点 Y：原始数据点 xi：插值点 Yi：插值点 格式 (1)yi = interp1(x,Y,xi) 返回插值向量yi，每一元素对应于参量xi，同时由向量x 与Y 的内插值决定。参量x 指定数据Y 的点。 若Y 为一矩阵，则按Y 的每列计算。yi 是阶数为 length(xi)*size(Y,2)的输出矩阵。 (2)yi = interp1(Y,xi) 假定x=1:N，其中N 为向量Y 的长度，或者为矩阵Y 的行数。(3)yi = interp1(x,Y,xi,method) 用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； ’linear’：线性插值（缺省方式），直接完成计算； ’spline’：三次样条函数插值。对于该方法，命令interp1 调用 函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函数。命令spline 用它们执行三次样条函数插值； ’pchip’：分段三次Hermite 插值。对于该方法，命令interp1 调用函数pchip，用于对向量x 与y 执行分段三次内插值。该方法保留单调性与数据的外形； ’cubic’：与’pchip’操作相同； ’v5cubic’：在MATLAB 5.0 中的三次插值。 对于超出x 范围的xi 的分量，使用方 法’nearest’、’linear’、’v5cubic’的插值算法，相应地将返回NaN。对其他的方法，interp1 将对超出的分量执行外插值算法。 (4)yi = interp1(x,Y,xi,method,'extrap') 对于超出x 范围的xi 中的分量将执行特殊的外插值法extrap。(5)yi = interp1(x,Y,xi,method,extrapval) 确定超出x 范围的xi 中的分量的外插值extrapval，其值通常取NaN 或0。 例1 1. 2.>>x = 0:10; y = x.*sin(x); 3.>>xx = 0:.25:10; yy = interp1(x,y,xx); 4.>>plot(x,y,'kd',xx,yy) 例2 1. 关于几种曲线拟合基本方法的比较 学院：材料科学与工程学院专业：材料学（博）姓名：郑文静学号：1014208040 在实际工作中，变量之间的关系未必都是线性关系，更多时候，它们之间呈现出了曲线关系，在科学实验或社会活动中，通过实验或观测得到一些x和y数据，为了对位置点进行研究，很多时候，我们通过曲线拟合的方式，将这些离散点近似为一条连续的曲线，从而来预测或者得到所需结果。曲线拟合的方法很多，本文中，主要讨论了曲线拟合的三种基础方法--插值法、磨光法、最小二乘法的特点，并对其在科学实验和生产实践中的应用性进行了比较。 插值法是函数逼近的一种基本方法，插值法就是通过函数在有限个点处的取值情况，估算出函数在其他点处的近似值。插值法中，选取不同的插值公式，来满足实际或运算需求，得到拟合的函数。其中，最基础的插值方法是三弯矩法，该方法是利用拉格朗日插值为基础，已知平面中的n+1个不同点，寻找一条n次多项式曲线通过这些点。该曲线具有唯一性。另外，还有三转角法，该方法是利用Henmiter插值为基础，其思路与三弯矩法相同，已知条件有所差别，在Henmiter插值中，不仅已知函数在一些点的函数值，而且，还知道它在这些点的导数值，甚至知道其高阶导数值，要求所求函数不仅满足过这些点，同时也要求其导函数，甚至高阶导函数满足条件。采用Henmiter插值法求得的多项式比拉格朗日法求得的多项式有较高的光滑逼近要求。此外，还有以分段和B-样条函数为基础的δ-基函数法，其中，样条函数是：对于[a,b]上的划分，称函数S(x)为[a,b]上关于划分△的k次样条函数，记做S k，△[a,b]。该方法避免了高次插值可能引起的大幅度波动现象，在实际中通常采用分段低次插值来提高近似程度。插值法常用于填充图像变换时像素之间的空隙。 磨光法是适应保凸性要求的数据拟合方法。积分可以改变函数的光滑度，而微商是积分的逆运算，对函数进行积分，然后在微商，可以将函数还原。而差商近似为微商，对函数积分后差商，可以将函数近似还原，同时可以更光滑，这种变换就是磨光。可以采用其他方法拟合得到函数，对于不光滑的点采用一次或多次磨光，得到更加光滑连续的函数。这种方法常用于外形设计。 最小二乘法也是函数逼近的一种基本方法。该方法不要求拟合曲线通过已知点，而是通过最小化误差的平方和寻找数据的最佳函数匹配。其解题步骤是：首先通过数据点，确定其可能所属的函数类型；然后，设出函数，并求出误差平方和的表达式；之后，由表达式对函 插值法与曲线拟合上机报告 江帆2010042020009 电子科技大学，物理电子学院,四川,成都 摘要：本文对插值法和曲线拟合进行了研究和分析，对拉格朗日插值法、分段线性插值法，就实例计算了拉格朗日插值法的五次插值函数和分段线性插值插值函数，画出了两者图像，比较了他们的插值余项。对于牛顿前插公式的插值余项进行了研究。就实例对给定的数据进行了二次拟合和指数函数拟合，比较了他们的拟合优劣。 关键字：雅各比迭代法，高斯迭代法，LU分解法，追赶法，线性方程组 引言 在生产实际以及科学实验中经常要研究变量之间的函数关系，但是在很多情况下很难找到具体的函数表达式，往往只能用测量或者观察获得一张数据表，即根据这种用表格形式给出的函数无法得到不在表中的点的函数值，也不能进一步研究函数的分析性质，。有的虽然能给出一个函数的分析表达式，但是式子很复杂,不适合使用。为了解决这些问题，我们设法通过这张表格求出一个简单函数P(x)，使这种球P(x)的方法称为插值法。 对于插值法由于节点上的函数值一般测量或实验得到的数据，如果个别点的误差较大插值函数保留了这些误差，影响逼近的精度。为此我们希望用另外的方法构造逼近函数，使得从总体上趋势上更能反映被逼近函数的特性，希望求得的逼近函数与已给的函数总体上来说其偏差按某种方法度量能达到最小。这就是曲线拟合。 1.拉格朗如插值法 1.1 理论原理 通常我们采用构造法，即直接构造一个满足差值条件式的n次插值多项式。若令则有，和被称为一次差值的基本插值多项式。利用基本插值多项式容易得出满足条件的n次插值多项式。常n次拉格朗日插值多项式，常记为，即 1.2 算法描述 设方程组的系数矩阵的对角线元虚为迭代次数容许的最大值，为容许误差。 （1）取Xi初值0.25到0.75，步长为0.001，令。插值和拟合
插值法和拟合实验报告(数值计算)
第五章插值与拟合答案—牟善军
插值与拟合(使用插值还是拟合)
计算方法--插值法与拟合实验
matlab中插值拟合与查表
第四章曲线拟合和多项式插值 - Hujiawei Bujidao
数学建模案例分析插值与拟合方法建模1数据插值方法及应用
插值和拟合区别
实验四 插值法与曲线拟合
第4、5讲 插值与拟合 作业参考答案
数值计算插值法与拟合实验
插值法与数据拟合法
清华大学_计算方法(数学实验)实验2插值与拟合
插值与拟合例题
数值分析实验插值与拟合
计算方法上机作业插值与拟合实验报告
Matlab中插值拟合函数汇总和使用说明
关于几种曲线拟合基本方法的比较
数值方法-插值与拟合
相关主题
文本预览