樣條插值函數(shù)及應(yīng)用

1、樣條插值函數(shù)及應(yīng)用摘要樣條函數(shù)具有廣泛的應(yīng)用，是現(xiàn)代函數(shù)論的一個十分活躍的分支，是計算方法的主要基礎(chǔ)和工具之一，由于生產(chǎn)和科學(xué)技術(shù)向前發(fā)展的推動以及電子計算機廣泛應(yīng)用的需要，人們便更多地應(yīng)用這個工具，也更深刻的認(rèn)識了它的本質(zhì)。在實際問題中所遇到許多函數(shù)往往很復(fù)雜，有些甚至是很難找到解析表達式的。根據(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第二章AMCM91A 估計水塔水流量52.1 理論分析及計算62.2運用MATLAB軟件計算9參考文獻14引 言樣條函數(shù)具有廣泛的應(yīng)用，是現(xiàn)代函數(shù)論的一個十分活躍的分支，是計算方法的主要基礎(chǔ)和工具之一，由于生產(chǎn)和科學(xué)技術(shù)向前發(fā)展的推動以及電子計算機廣泛應(yīng)用的需要，人們便更多地應(yīng)用這個工具，也更深刻的認(rèn)識了它的本質(zhì)。上世紀(jì)四十年代，在研究數(shù)據(jù)處理的問題中引出了樣條函數(shù)，

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

4、內(nèi)外學(xué)者對這方面的研究也越來越重視，根據(jù)我們的需要來解決不同的問題，而且函數(shù)的形式也在不斷地改進，長期以來很多學(xué)者致力于樣條插值的研究，對三次樣條的研究已相當(dāng)成熟。第一章三次樣條插值1.1 樣條插值函數(shù)簡介 在實際問題中所遇到許多函數(shù)往往很復(fù)雜，有些甚至是很難找到解析表達式的。有時通過實驗或者數(shù)值計算所得到的也只是一些離散的點上的函數(shù)值，即。根據(jù)函數(shù)已有的數(shù)據(jù)來計算函數(shù)在一些新的點處的函數(shù)值，就是插值法所需要解決的問題。插值法的基本思想就是，首先根據(jù)已有的函數(shù)值來構(gòu)造一個簡單的函數(shù)作為的近似表達式，然后用來計算新的點上的函數(shù)值作為的近似值。通?？梢赃x多項式函數(shù)作為近似函數(shù)，因為多項式具有各階導(dǎo)

5、數(shù)，求值也比見方便。常用的有Lagrange插值、Newton插值、Hermite插值和樣條插值。線性插值在分段點上僅連續(xù)而不可導(dǎo)，三次埃爾米特插值有連續(xù)的一階導(dǎo)數(shù)，這樣的光滑程度常不能滿足物理問題的需要，樣條函數(shù)可以同時解決這兩個問題,使插值函數(shù)既是低階分段函數(shù),又是光滑的函數(shù)，并且只需在區(qū)間端點提供某些導(dǎo)數(shù)信息。三次樣條函數(shù)定義：設(shè)在區(qū)間上取個節(jié)點，函數(shù)在各個節(jié)點處的函數(shù)值為若滿足： (1) (2) 在區(qū)間上，具有連續(xù)的二階導(dǎo)數(shù)； (3) 在區(qū)間上,是三次的多項式；則稱是函數(shù)在區(qū)間上的三次樣條插值函數(shù)。 由以上定義可以看出，雖然每個子區(qū)間上的多項式可以各不相同，但在相鄰子區(qū)間的連接處卻是光

6、滑的。因此，樣條插值也稱為分段光滑插值。從定義知要求出，在每一個小區(qū)間上確定4個待定系數(shù)，共有n個小區(qū)間，故應(yīng)有4n個參數(shù)。根據(jù)在上二階導(dǎo)數(shù)連續(xù)，在節(jié)點出滿足連續(xù)性條件共有3n-3個條件，再加上滿足插值條件共有4n-2個條件，因此還需要2個條件才能確定。通?？稍趨^(qū)間端點上各加一個條件（稱為邊界條件），可根據(jù)實際問題的要求給定，通常有以下三種：（1）已知端點的一階導(dǎo)數(shù)值，即（2）倆端點的二階導(dǎo)數(shù)已知，即其特殊情況稱為自然邊界條件。（3）當(dāng)是以為周期的函數(shù)時，則要求也是周期函數(shù)。這時邊界條件應(yīng)滿足而此時。這樣確定的樣條函數(shù)，稱為周期函數(shù)。1.2 三次樣條函數(shù)應(yīng)用作函數(shù)在取間隔為0.1的點圖，用插值

7、進行實驗。使用MATLAB軟件程序代碼如下：%產(chǎn)生原始數(shù)據(jù)x=0:0.1:1;y=(x.2-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);%作圖subplot(1,2,2)plot(x,y,ro,xx,yy,b) 運行截圖圖1.1 運行結(jié)果第二章AMCM91A 估計水塔水流量美國某洲的各用水管理機構(gòu)要求各社區(qū)提供以每小時多少加侖計的用水率以及每天總的用水量,但許多社區(qū)并沒有測量水流入或流出當(dāng)?shù)厮乃康脑O(shè)備,他們只能代之以每

8、小時測量水塔中的水位,精度在0.5%以內(nèi),更為重要的是,無論什么時候,只要水塔中的水位下降到某一最低水位L時,水泵就啟動向水塔重新充水至某一最高水位H,但也無法得到水泵的供水量的測量數(shù)據(jù)。因此，在水泵工作時，人們?nèi)菀捉⑺械乃慌c水泵工作時的用水量之間的關(guān)系。水泵每天向水塔充水一次或兩次，每次約兩小時。表1 白某小鎮(zhèn)某天的水塔水位時間033166635106191393717921212402522328543水位317531103054299429472892285027972752時間322843593539332394354331846636499535393657254水位2697

9、水泵工作水泵工作355034453350326031673087時間605746455468535718547502179154826498596889953水位30122927284227672697水泵工作水泵工作34753397時間93270水位3340單位：時間/秒；水位/0.01英尺試估計在任何時刻，甚至包括水泵正在工作期間內(nèi)，水從水塔流出的流量，并估計一天的總用水量，表1中給出了某個真實小鎮(zhèn)某一天的真實數(shù)據(jù)。表1中給出了從第一次測量開始的以秒為單位的時刻，以及該時刻的高度單位為百分之一英尺的水塔中水位的測量值，例如，3316秒后，水塔中的水位達到31.10英尺。水塔是一個垂直圓形柱

10、體，高為40英尺，直徑57英尺，通常當(dāng)水塔的水位降至27.00英尺時水泵開始向水塔充水，而當(dāng)水塔的水位升至35.50英尺時水泵停止工作。2.1 理論分析及計算1. 水塔充水時間的確定(1) 第一次充水時間的確定當(dāng)時間t=32284秒時，水位26.97英尺,約低于最低水位27英尺，因此可作為第一次開始充水時刻。當(dāng)t=39435秒時，水塔水位35.5英尺，恰為最高水位，因此可作為第一次充水的結(jié)束時刻。充水時間為=(39435-32284)/3600=1.9864小時，也接近充水時間2小時。(2) 第二次充水時間的確定當(dāng)時間t=75021秒時，水位26.97英尺,約低于最低水位27英尺，因此可作為第

11、二次開始充水時刻。當(dāng)t=82649秒時，水泵在工作，但充水時間達到=(8264-7502)/3600=2.1189小時；但下一時刻t=85968時，水塔水位34.75英尺，低于最高水位35.50 英尺。因此可將t=82649秒作為第二次充水的結(jié)束時刻,且該時刻水位為最大充水高度35.50 英尺。2.計算各時刻塔內(nèi)水的體積單位轉(zhuǎn)換為1英尺=0.3048米, 1升=1/3.加侖體積計算公式為 表2 不同時刻水體積表時間水體積時間水體積時間水體積0 (1)8.967819.03750.921110.9542 (2)19.95941.843112.032820.83922.949712.954422.

12、9581 (3)3.871413.875823.88004.978114.982224.98695.900015.903925.90837.006416.82617.928617.9317其中(1)表示第一段開始，(2) 表示第二段開始，(3) 表示第三段開始；單位：時間(小時)，水流量：加侖/小時；3.計算各時刻點的水流量(加侖/小時)水流量公式為： 以上25個時刻處的水流量采用差分的方法得到，共分三段分別處理。差分公式為：(1) 對每段前兩點采用向前差分公式(2) 對每段最后兩點采用向后差分公式(3) 對每段中間點采用中心差分公式得到各點水流量表表3 不同時刻水流量表時間水流量時間水流量時

13、間水流量0 (1)144048.96781102319.0375166520.92111118210.954(2)1946919.9594144951.84311006312.03282019520.8392146482.94971101212.95441894122.958(3)152203.8714879813.87581590323.8800152634.9781999114.98221805524.9869137115.9000812415.90391564625.908396347.00641016116.8261137427.9286848717.931714962其中(1)表示第

14、一段開始，(2) 表示第二段開始，(3) 表示第三段開始：單位：時間(小時)，水流量：加侖/小時；4用三次樣條擬合流量數(shù)據(jù)對表3中25個時刻點的流量數(shù)據(jù)采用三次樣條插值得到一條光滑曲線，作為任意時刻的流量曲線,見圖2.1。圖2.1 水塔流量圖其中*表示數(shù)據(jù)點，實線為樣條曲線5.一天總用水量計算一天流水總量計算:方法1:直接積分法: 方法2:分段計算法第一次充水前用水(加侖) 第一次充水后第二次充水前用水(加侖)22.9581，23.88期間用水(加侖)第一次充水期間用水: 第二次充水期間用水: 23.88，24期間用水: 總共用水兩種方法結(jié)果相差6.水泵水流量計算第一次充水期間水塔體積增加充水

15、時間:第一次充水期間水泵平均流量第二次充水期間水塔體積增加充水時間:第二次充水期間水泵平均流量則整個充水期間水泵平均流量2.2運用MATLAB軟件計算1. 依據(jù)理論分析，編寫程序代碼。MATLAB程序代碼c=0.3048; %1英尺等于0.3048米p=1.0/3.785; %1升=1/3.加侖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

16、,3550; 43318,3445;46636,3350;49953,3260;53936,3167;57254,3087; 60574,3012;64554,2927;68535,2842;71854,2767;75021,2697; 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; %計算第一段for i=1:

17、n1 if i=2 %前兩點采用向前差分 f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2)/(2*(t(i+1)-t(i); elseif i=n1-1 f(i)=abs(3*v(i)-4*v(i-1)+v(i-2)/(2*(t(i)-t(i-1); end endn2=21; %計算第二段for i=n1+1:n2 if i=n1+2 %前兩點采用向前差分 f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2)/(2*(t(i+1)-t(i); elseif i=n2-1 f(i)=abs(3*v(i)-4*v(i-1)+v(i-2)/(2*(t(i)-t(i-1)

18、; end endn3=25; %計算第三段for i=n2+1:n3 if i=n2+2 %前兩點采用向前差分 f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2)/(2*(t(i+1)-t(i); elseif i=n3-1 f(i)=abs(3*v(i)-4*v(i-1)+v(i-2)/(2*(t(i)-t(i-1); end endplot(t,f,r*); %畫原始點圖tmin=min(t); tmax=max(t);tt=tmin:0.1:tmax; %獲得離散的時間點,用于作樣條曲線ff=spline(t,f,tt); %計算三次樣條插值hold onplot(tt,

19、ff,b); %畫樣條曲線xlabel(時間(小時);ylabel(流量(加侖/小時);title(水塔流量圖);hold offdt=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小時用水量,采用復(fù)化梯形公式fprintf(全部積分法)1天總水流量s= %8.2fn,s); v10=v(11)-v(10); %第一次水塔增加的水 dt1=t(11)-t(10); tp=t(10):dt:t(11); %第一次充水其間流

20、出的水 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); vv1=v(1)-v(10); %第一次充水前總流量 vv2