1、樣條插值函數(shù)及應(yīng)用摘要樣條函數(shù)具有廣泛的應(yīng)用,是現(xiàn)代函數(shù)論的一個十分活潑的分支,是計算方法的主要根底和工具之一,由于生產(chǎn)和科學(xué)技術(shù)向前開展的推動以及電子計算機廣泛應(yīng)用的需要,人們便更多地應(yīng)用這個工具,也更深刻的熟悉了它的本質(zhì).在實際問題中所遇到許多函數(shù)往往很復(fù)雜,有些甚至是很難找到解析表達(dá)式的.根據(jù)函數(shù)已有的數(shù)據(jù)來計算函數(shù)在一些新的點處的函數(shù)值,就是插值法所需要解決的問題.插值法是數(shù)值逼近的重要方法之一,它是根據(jù)給定的自變量值和函數(shù)值,求取未知函數(shù)的近似值.早在一千多年前,我國科學(xué)家就在研究歷法時就用到了線性插值和二次插值.而在實際問題中,有許多插信函數(shù)的曲線要求具有較高的光滑性,在整個曲線中

2、,曲線不但不能有拐點,而且曲率也不能有突變.因此,對于插值函數(shù)必須二次連續(xù)可微且不變號,這就需要用到三次樣條插值.關(guān)鍵詞三次樣條函數(shù)；插值法引言1第一章三次樣條插值21.1 樣條插值函數(shù)簡介21.2 三次樣條函數(shù)應(yīng)用3第二章AMCM91粘計水塔水流量52.1 理論分析及計算62.2 運用MATLA歆件計算9參考文獻(xiàn)14引言樣條函數(shù)具有廣泛的應(yīng)用,是現(xiàn)代函數(shù)論的一個十分活潑的分支,是計算方法的主要根底和工具之一,由于生產(chǎn)和科學(xué)技術(shù)向前開展的推動以及電子計算機廣泛應(yīng)用的需要,人們便更多地應(yīng)用這個工具,也更深刻的熟悉了它的本質(zhì).上世紀(jì)四十年代,在研究數(shù)據(jù)處理的問題中引出了樣條函數(shù),例如,在1946年

3、Schoenberg將樣條引入數(shù)學(xué),即所謂的樣條函數(shù),直到五十年代,還多應(yīng)用于統(tǒng)計數(shù)據(jù)的處理方面,從六十年代起,在航空、造船、汽車等行業(yè)中,開始大量采用樣條函數(shù).在我國,從六十年代末開始,從船體數(shù)學(xué)放樣到飛機外形設(shè)計,逐漸出現(xiàn)了一個使用樣,逐漸出現(xiàn)了一個使用樣條函數(shù)的熱潮,并推廣到數(shù)據(jù)處理的許多問題中.在實際生活中有許多計算問題對插值函數(shù)的光滑性有較高的要求,例如飛機機翼外形、發(fā)動機進(jìn)、排氣口都要求有連續(xù)的二階導(dǎo)數(shù),用三次樣條繪制的曲線不僅有很好的光滑度,而且當(dāng)節(jié)點逐漸加密時其函數(shù)值整體上能很好地逼近被插函數(shù),相應(yīng)的導(dǎo)數(shù)值也收斂于被插函數(shù)的導(dǎo)數(shù)值,不會發(fā)生“龍格現(xiàn)象.現(xiàn)在國內(nèi)外學(xué)者對這方面的研

4、究也越來越重視,根據(jù)我們的需要來解決不同的問題,而且函數(shù)的形式也在不斷地改良,長期以來很多學(xué)者致力于樣條插值的研究,對三次樣條的研究已相當(dāng)成熟.第一章三次樣條插值1.1樣條插值函數(shù)簡介在實際問題中所遇到許多函數(shù)f(x)往往很復(fù)雜,有些甚至是很難找到解析表達(dá)式的.有時通過實驗或者數(shù)值計算所得到的也只是一些離散的點X(i=0,1,2.n)上的函數(shù)值,即yi=f(xi),i=0,1,2.n.根據(jù)函數(shù)f(x)已有的數(shù)據(jù)來計算函數(shù)f(x)在一些新的點x處的函數(shù)值,就是插值法所需要解決的問題.插值法的根本思想就是,首先根據(jù)已有的函數(shù)值來構(gòu)造一個簡單的函數(shù)y(x)作為f(x)的近似表達(dá)式,然后用y(x)來計

5、算新的點上的函數(shù)值作為f(x)的近似值.通?？梢赃x多項式函數(shù)作為近似函數(shù)y(x),由于多項式具有各階導(dǎo)數(shù),求值也比見方便.常用的有Lagrange插值、Newton插值、Hermite插值和樣條插值.線性插值在分段點上僅連續(xù)而不可導(dǎo),三次埃爾米特插值有連續(xù)的一階導(dǎo)數(shù),這樣的光滑程度常不能滿足物理問題的需要,樣條函數(shù)可以同時解決這兩個問題,使插值函數(shù)既是低階分段函數(shù),又是光滑的函數(shù),并且只需在區(qū)間端點提供某些導(dǎo)數(shù)信息.三次樣條函數(shù)定義：設(shè)在區(qū)間a,b】上取n+1個節(jié)點a=x0<x1<x2<xn=b,函數(shù)y=f(x)在各個節(jié)點處的函數(shù)值為yi=f(xj,i=0,1,.,n,假設(shè)S

6、(x)滿足：(1)S(x)=yi,i=Q1,.n;(2)在區(qū)間Ja,b】上,S(x)具有連續(xù)的二階導(dǎo)數(shù)；(3)在區(qū)間k/i=0,1,.,n-1)上,S(x)是x三次的多項式；那么稱S(x)是函數(shù)y=f(x)在區(qū)間a,b上的三次樣條插值函數(shù).由以上定義可以看出,雖然每個子區(qū)間上的多項式可以各不相同,但在相鄰子區(qū)間的連接處卻是光滑的.因此,樣條插值也稱為分段光滑插值.從定義知要求出S(x),在每一個小區(qū)間k,xi_Ji=0,1,.,n-1)上確定4個待定系數(shù),共有n個小區(qū)間,故應(yīng)有4n個參數(shù).根據(jù)S(x)在口力】上二階導(dǎo)數(shù)連續(xù),在節(jié)點x(i=1,2.n-1)出滿足連續(xù)性條件S(Xi-0)=S(Xi

7、0)S(Xi-0)=S(Xi0).一、,二iS(Xi-0)=S(xi0)共有3n-3個條件,再加上S(x)滿足插值條件S(x)=丫4=0,1,.n;共有4n-2個條件,因此還需要2個條件才能確定S(x)o通?？稍趨^(qū)間端點上各加一個條件(稱為邊界條件),可根據(jù)實際問題的要求給定,通常有以下三種：(1)端點的一階導(dǎo)數(shù)值,即S(x0)-f0,S(xn)=fn(2)倆端點的二階導(dǎo)數(shù),即S(x0)=f0,S(xn)=fn其特殊情況S"(%)=S"(Xn)=0,稱為自然邊界條件.(3)當(dāng)f(x)是以X.為周期的函數(shù)時,那么要求S(X)也是周期函數(shù).這時邊界條件應(yīng)滿足S(Xi-0)=S(

8、Xi0)S(Xi-0)=S(Xi0)S(xi-0)-S(xi0)而此時y°=yn.這樣確定的樣條函數(shù)S(x),稱為周期函數(shù).1.2三次樣條函數(shù)應(yīng)用作函數(shù)y=(x2-3x+7)eXsin(2x)ft10,1取間FB為0.1的點圖,用插值進(jìn)行實驗.使用MATLA歆件程序代碼如下：%產(chǎn)生原始數(shù)據(jù)x=0:0.1:1;y=(x.A2-3*x+7).*exp(-4*x).*sin(2*x);%乍圖subplot(1,2,1);plot(x,y,x,y,'ro')%寺求插值點xx=0:0,02:1;yy=interp1(x,y,xx,'spline');%乍圖sub

9、plot(1,2,2)plot(x,y,'ro',xx,yy,'b')運行截圖運行結(jié)果第二章AMCM91AT計水塔水流量美國某洲的各用水治理機構(gòu)要求各社區(qū)提供以每小時多少加侖計的用水率以及每天總的用水量,但許多社區(qū)并沒有測量水流入或流出當(dāng)?shù)厮乃康脑O(shè)備,他們只能代之以每小時測量水塔中的水位,精度在0.5%以內(nèi),更為重要的是,無論什么時候,只要水塔中的水位下降到某一最低水位L時,水泵就啟動向水塔重新充水至某一最高水位H,但也無法得到水泵的供水量的測量數(shù)據(jù).因此,在水泵工作時,人們?nèi)菀捉⑺械乃慌c水泵工作時的用水量之間的關(guān)系.水泵每天向水塔充水一次或兩次,

10、每次約兩小時.表1白某小鎮(zhèn)某天的水塔水位時間033166635106191393717921212402522328543水位317531103054299429472892285027972752時間322843593539332394354331846636499535393657254水位2697水泵工作水泵工作355034453350326031673087時間605746455468535718547502179154826498596889953水位30122927284227672697水泵工作水泵工作34753397時間93270水位3340單位：時間/秒；水位/0.01英尺試

11、估計在任何時刻,甚至包括水泵正在工作期間內(nèi),水從水塔流出的流量f(t),并估計一天的總用水量,表1中給出了某個真實小鎮(zhèn)某一天的真實數(shù)據(jù).表1中給出了從第一次測量開始的以秒為單位的時刻,以及該時刻的高度單位為百分之一英尺的水塔中水位的測量值,例如,3316秒后,水塔中的水位到達(dá)31.10英尺.水塔是一個垂直圓形柱體,高為40英尺,直徑57英尺,通常當(dāng)水塔的水位降至27.00英尺時水泵開始向水塔充水,而當(dāng)水塔的水位升至35.50英尺時水泵停止工作.2.1 理論分析及計算1 .水塔充水時間確實定(1)第一次充水時間確實定當(dāng)時間t=32284秒時,水位26.97英尺,約低于最低水位27英尺,因此可作為

12、第一次開始充水時刻.當(dāng)t=39435秒時,水塔水位35.5英尺,恰為最高水位,因此可作為第一次充水的結(jié)束時刻.充水時間為dt=(39435-32284)/3600=1.9864小時,也接近充水時間2小時.(2)第二次充水時間確實定當(dāng)時間t=75021秒時,水位26.97英尺,約低于最低水位27英尺,因此可作為第二次開始充水時刻.當(dāng)t=82649秒時,水泵在工作,但充水時間到達(dá)dt=(8264-7502)/3600=2.1189小時；但下一時刻t=85968時,水塔水位34.75英尺,低于最高水位35.50英尺.因此可將t=82649秒作為第二次充水的結(jié)束時刻,且該時刻水位為最大充水高度35.5

13、0英尺.2 .計算各時刻塔內(nèi)水的體積單位轉(zhuǎn)換為1英尺=0.3048米,1升=1/3.785411加侖體積計算公式為v=：.d2h/4表2不同時刻水體積表時間水體積時間水體積時間水體積0(1)6061258.967851487219.03755425540.921159371610.954267771519.95945282361.843158302612.032865767020.83925148722.949757157112.954463953422.9581(3)6777153.871456259913.875862235223.88006633974.978155209914.9822

14、60459824.98696485065.900054408115.903958932525.90836376257.006453396316.82615750087.928652537217.9317558781其中(1)表示第一段開始,(2)表示第二段開始,(3)表示第三段開始;單位：時間(小時),水流量：加侖/小時；3 .計算各時刻點的水流量加侖/小時水流量公式為：|dv(t)I出以上25個時刻處的水流量采用差分的方法得到,共分三段分別處理差分公式為：1對每段前兩點采用向前差分公式f(ti)=-3Vi4Vi-Vi22(tii-1)2對每段最后兩點采用向后差分公式f(ti)=W-4VVc2

15、(t-t)3對每段中間點采用中央差分公式f(ti)=V七+8Vt8/-+V/12(tkti)得到各點水流量表表3不同時刻水流量表時間水流量時間水流量時間水流量0(1)P144048.96781102319.03751166520.92111118210.954(2)1946919.9594144951.84311006312.03282021520.8392146482.9497:1101212.95441894122.958(3)152203.8714879813.87581590323.88001152634.9781999114.98221805524.9869137115.9000:8

16、12415.90391564625.908396347.00641016116.8261137427.9286848717.931714962其中1表示第一段開始,2表示第二段開始,3表示第三段開始：單位：時間小時,水流量：加侖/小時；4.用三次樣條擬合流量數(shù)據(jù)對表3中25個時刻點的流量數(shù)據(jù)采用三次樣條插值得到一條光滑曲線,作為任意時刻的流量曲線,見圖2.1水塔流量圖0.80510152025時間小時28.B421111苣47中縣煤1圖2.1水塔流量圖其中*表示數(shù)據(jù)點,實線為樣條曲線5.一天總用水量計算一大流水總量計算：方法1:直接積分法：24§=ftdt=332986加侖方法2:分

17、段計算法第一次充水前用水Vi=606125514872=91253力口侖第一次充水后第二次充水前用水V2=677715-514872=162843加侖22.9581,23.88期間用水V3=677715663397=14318力口侖10.9542第一次充水期間用水：V4=ftdt=30326加侖8.967822.9581第二次充水期間用水：V5=口392ftdt=31605加侖,2423.88,24期間用水:V6=ftdt=1524加侖23.886總共用水S2=£Vi=331869加侖i4兩種方法結(jié)果相差err=0.34%6,水泵水流量計算第一次充水期間水塔體積增加V1=677715

18、-514872=162843力口侖充水時間：純=10.9542-8.9678=1.9864小時第一次充水期間水泵平均流量10.9542MftdtP1=678=97246力1侖/小時:t1第二次充水期間水塔體積增加忸2=677715-514872=162843力口侖充水時間：t2=22.9581-20.8392=2.1189小時第二次充水期間水泵平均流量22.9581V2ftdtP2=208392-91769力1侖/小時上那么整個充水期間水泵平均流量P=P+P2/2=94507加侖/小時2.2運用MATLA歆件計算1,依據(jù)理論分析,編寫程序代碼.MATLAB程序代碼c=0.3048;%1英尺等于

19、0.3048米p=1.0/3.785;%1升=1/3,785411力口侖d=57*c;h=31,75*c;v=pi*d*d*h/4*1000*p;data=0,3175;3316,3110;6635,3054;10619,2994;13937,2947;17921,2892;21240,2850;25223,2797;28543,2752;32284,2697;39435,3550;43318,3445;46636,3350;49953,3260;53936,3167;57254,3087;60574,3012;64554,2927;68535,2842;71854,2767;75021,26

20、97;82649,3550;85968,3475;89953,3397;93270,3340;%原始數(shù)據(jù)t=data(:,1)/3600;%計算時間(小時為單位)v=pi*d*d*data(:,2)/100*c/4*1000*p;%計算體積%計算差分n=length(v);f=zeros(n,1);%存儲差分值n1=10;%算第一段fori=1:n1ifi<=2%前兩點采用向前差分f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2)/(2*(t(i+1)-t(i);elseifi<=n1-2%采用三點中央差分公式f(i)=abs(-v(i+2)+8*v(i+1)-8*v

21、(i-1)+v(i-2)/(12*(t(i+1)-t(i);elseifi>=n1-1f(i)=abs(3*v(i)-4*v(i-1)+v(i-2)/(2*(t(i)-t(i-1);endendn2=21;%算第二段fori=n1+1:n2ifi<=n1+2%前兩點采用向前差分f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2)/(2*(t(i+1)-t(i);elseifi<=n2-2f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2)/(12*(t(i+1)-t(i);elseifi>=n2-1f(i)=abs(3*v(i)

22、-4*v(i-1)+v(i-2)/(2*(t(i)-t(i-1);endendn3=25;%計算第三段fori=n2+1:n3ifi<=n2+2%前兩點采用向前差分f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2)/(2*(t(i+1)-t(i);elseifi<=n3-2f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2)/(12*(t(i+1)-t(i);elseifi>=n3-1f(i)=abs(3*v(i)-4*v(i-1)+v(i-2)/(2*(t(i)-t(i-1);endendplot(t,f,'r*'

23、);刎原始點圖tmin=min(t);tmax=max(t);tt=tmin:0.1:tmax;%獲得離散的時間點,用于作樣條曲線ff=spline(t,f,tt);%計算三次樣條插值holdonplot(tt,ff,'b');%畫樣條曲線xlabel(,時間(小時),)；ylabel('流量(加侖/小時),)；title(,水塔流量圖,);holdoffdt=0.05;t2=0.5:dt:24.5;%獲得離散白時間點,用于積分nn=length(t2);f2=spline(t,f,t2);s=(f2(1)+f2(nn)+2*sum(f2(2:nn-1)*dt/2;%計

24、算24小時用水量,采用復(fù)化梯形公式fprintf('(全部積分法)1天總水流量s=%8.2fn',s);v10=v(11)-v(10);%第一次水塔增加的水dt1=t(11)-t(10);tp=t(10):dt:t(11);%第一次充水其間流出的水nn=length(tp);yp=spline(t,f,tp);%計算三次樣條插值v11=(yp(1)+yp(nn)+2*sum(yp(2:nn-1)*dt/2;v1=v10+v11;%第一次充水總量p1=v1/dt1;%第一次充水的平均水流量v20=v(22)-v(21);%第二次水塔增加的水dt2=t(22)-t(21);tp1=t(21):dt:t(22);%第二次充水其間流出的水nn=length(tp1);yp1=spline(t,f,tp1);%+算三次樣條插值v21=(yp1(1)+yp1(nn)+2*sum(yp1(2:nn-1)*dt/2;v2=v20+v21;%!二次充水總量p2=v2/dt2;%第二次充水的平均水流量p=(p1+p2)/2;%兩次充水平均水流量fprintf('兩次充水平均水流量p=%8.2fn',p