三次样条插值算法详解概要

1三次样条插值鉴于高次插值不收敛又不稳定的特点，低次插值既具有收敛性又具有稳定性，因此低次值更具有实用价值，但是低次插值的光滑性较差，比如分段线性插值多项式在插值区间中仅具有连续性，在插值节点处有棱角，一阶导数不存在；分段三次Hermite插值多项式在插值区间中仅具有一阶导数即一阶光滑性但不具备二阶光滑性，不能满足某些实际应用比如汽车、轮船、飞机等的外形中流线形设计。样条是在二十世纪初期经常用于图样设计的一种富有弹性的细长条，多个样条互相弯曲连接后沿其边缘画出的曲线就是三次样条曲线。后来数学上对其进行了抽象，定义了m次样条函数，并成为数值逼近的重要研究分枝，进一步扩大了样条函数的应用范围。2样条函数的定义定义4.1设区间[a,b]上给定一个节点划分a=x0x1……xn-1xn=b如果存在正整数k使得[a,b]上的分段函数s(x)满足如下两条：(1)在[a,b]上有直到k-1阶连续导数。(2)在每个小区间[xi,xi+1]上是次数不大于k的多项式。则称分段函数s(x)是以(2.6)为节点集的k次样条函数。3处的函数值为在节点如果函数nxxxxf,,,)(10njyxfjj,,1,0,)(njyxSjj,,1,0,)(三次样条插值函数的定义并且关于这个节点集的三次样条函数s(x)满足插值条件：则称这个三次样条函数s(x)为三次样条插值函数。4则其必满足的三次样条插值函数是如果,)()(xfxSnjyxSjj,,1,0,)(1,,1,)()(limnjmxSxSjjxxj1,,1),()(limnjxSxSjxxj1,,1,)()(limnjyxSxSjjxxj三次样条插值函数的边界条件插值条件：连续性条件：一阶导数连续条件：二阶导数连续条件：5（1）因为s(x)在每个小区间上是一个次小于三次的多项式，故有四个未知系数；（2）因为s(x)有n分段，从而共有4n个未知系数！（3）但插值条件与样条条件仅给出4n-2个条件，无法定出4n个未知系数，还差2个条件！这2个条件我们用边界条件给出！6通常我们对插值多项式在两端点的状态加以要求也就是所谓的边界条件:第一边界条件：由区间端点处的一阶导数给出即30003()(),()(),nnnsxmfxsxmfx7第二边界条件：由区间端点处的二阶导数给出即30003()(),()(),nnnsxMfxsxMfx特殊情况为自然边界条件：由区间端点处的二阶导数恒为0给出即3003()0()0nnsxMsxM8这样三次样条插值问题就分成三类！其实不止这三类!第三类又称周期边界条件：由区间端点处的函数值或导数值满足周期条件给出303303303(0)(0)(0)(0)(0)(0)nnnsxsxsxsxsxsx9样条函数的例子容易验证：3215)18218(2115)1427163(1015)152611()(232323xxxxxxxxxxxxxS是满足如下数据的第一类边界样条插值问题解：x0123y0000y’1010样条函数的例子3215)18218(2115)1427163(1015)152611()(232323xxxxxxxxxxxxxS11通常有三转角法、三弯矩法、B样条基函数法。三次样条插值函数的求法这三种方法的基本思想是类似的，都是通过待定某些参数来确定插值函数，但肯定不是待定4n个参数。而是利用已知条件将待定参数减小到最少。比如：待定一阶导数、待定二阶导数、采用基函数方法来确定插值函数。12三转角法:待定一阶数为了确定三次样条插值函数的表达式S(x)，我们采用待定系数法来求解，我们待定什么系数呢？考虑到带一阶导数的分段三次Hermite插值多项式)()()()()(11111010)(3iiiiiiiihxxiihxxiihxxihxxiihmhmyyxH1,,1,0,1nixxhiiinnnxxxxHxxxxHxxxxHxH1)1(321)1(310)0(33)()()()(13njmxSjj,,1,0,)(我们采用待定一阶导数的方法即设因为分段三次Hermite插值多项式已经至少是一阶连续可导了，为了让它成为三次样条函数只需确定节点处的一阶导数使这些节点处的二阶导数连续即可！1,,1),0()0(nixSxSii)()()()()(11111010iiiiiiiihxxiihxxiihxxihxxihmhmyyxS1,,1,0,],,[11nixxhxxxiiiii212)1()(,)1)(12()(xxxxxx14并整理后得求二阶导数对,)(xSi)()2(6)(131iiiiiiyyhxxxxS11122624642iiiiiiiixxxxxxmmhh112642lim()(0)()iiiiiixxiiiSxSxyymmhhh112111624lim()(0)()iiiiiixxiiiSxSxyymmhhh15(0)(0),1,2,,1iiSxSxin由于在内部节点处二阶导数连续条件：11111)11(21iiiiiiimhmhhmh)(321121iiiiiihyyhyy1,,1ni111112iiiiiiiiimhhhmmhhh)(3111111iiiiiiiiiiiihyyhhhhyyhhh整理化简后得：16112,1,,1iiiiiimmmgin1iiiihhh111iiiiihhh1113()iiiiiiiiiyyyyghh个未知量个方程共个1,1nn称为三转角法基本方程组以上推导还没有考虑边界条件！针对不同类型的三次样条问题，就可以导出不同的方程组！17第一类三次样条插值问题方程组00)(mxSnnmxS)(基本方程组化为n-1阶方程组0112112fgmmkkkkkkgmmm1122,,3,2nknnnnnnfgmm111212由于已知：化为矩阵形式181111022223333422221111222222nnnnnnnnnmgmmgmgmgmgm这是一个严格对角占优的三对角方程组，用追赶法可以求解！19第二类三次样条插值问题的方程组00)(MxSnnMxS)(时，称为自然边界条件00nMM由于已知：)(6)(01200yyhxS004mh102mh0M)(6)(121nnnnyyhxS112nnmhnnmh14nM故得：200000110232MhhyymmnnnnnnnMhhyymm23211110gng稍加整理得212222121132211nnnnmmmmm1210nnggggg1210联合基本方程组得一个n+1阶三对角方程组，化成矩阵形式为：仍然是严格对角占优21第三类样条插值问题的方程组立即可得下式：nnmmxSxS00)0()0(由于：)0()0(0nxSxS)2(2)(6)0(01001200mmhyyhxS)2(2)(6)0(11121nnnnnnnmmhyyhxS22nnnnngmmm211110013nnnnnnhyyhyygnnnnnnhhhhhh1,101100其中：联合基本方程得一个广义三对角或周期三对角方程组：1111112222311222222nnnnnnmgmgmg这个方程组的系数矩阵仍然是严格对角占优阵！23求解这些不同类型的样条插值问题的方程组，我们可得所要待定的一阶导数：称为三次样条插值问题三转角公式！nmmmm,,,,210再代入S(x)的每一段表达式，就求得三次样条函数的表达式!)()()()()(11111010iiiiiiiihxxiihxxiihxxihxxihmhmyyxS1,......,2,1,0ni24例1.对于给定的节点及函数值1245()1342kkxfx的近似值并求插值函数的三次样条求满足自然边界条件)3(),(0)()(0fxSxSxSn解:这是自然边界条件下的样条问题。3213123113221kkkkhhh11kkkkhhhk125291g272g60g63g我们可以将上述计算列于表中：)(3111kkkkkkkkkhyyhyyg00001023MhhyygnnnnnnMhhyyg23111k0123xk1245yk1342mk????Mk0??0hk121*λk*2/31/3*μk*1/32/3*gk69/2-7/2-626213/223/13/123/2123210mmmm3210gggg由些得如下方程组：819,45,47,8173210mmmm利用三转角公式：)()()()()(11111010iiiiiiiihxxiihxxiihxxihxxihmhmyyxS27)()()()()(0100010101100010000hxxhxxhxxhxxhmhmyyxS)2()1()2(3)1(147181700xxxx21147838123xxxx)()()()()(1211121111211102011hxxhxxhxxhxxhmhmyyxS)()()(4)(32412522127240220xxxx42147838123xxxx同样可以求得第三段表达式！28323232137112884137()1248843451033345884xxxxSxxxxxxxxx417)3()3(Sf2930三弯矩法:待定二阶导数()(0,1,2,,)kkMSxkn选择二阶导数作为待定参数:(1)1111()()(),[,]()kkkkkkkkkkkkkkSxMMMSxxxMxxxSxMh由于三次样条S(x)是三次多项式,故它的二阶导数是一次多项式,从而1(0,1,,1)kkkhxxkn思考:(1)的原因?31(2)21()()()2kkkkkkkkMMSxxxMxxph(3)321()()()()62kkkkkkkkkkMMMSxxxxxpxxqh