1、Z,數(shù)學(xué)模型之 插值,插 值,一維插值 1.插值的基本原理； 2.三種插值方法： 分段線性插值; 三次樣條插值; 拉格朗日插值。 3. 用MATLAB作插值計(jì)算 4. 應(yīng)用實(shí)例,二維插值 1.兩類(lèi)插值方法： 網(wǎng)格節(jié)點(diǎn)插值法： 雙線性插值等。 散點(diǎn)數(shù)據(jù)插值法： 修正Shephard法 2. 用MATLAB作二維插值計(jì)算 3. 應(yīng)用實(shí)例,函數(shù)f(x)的產(chǎn)生辦法：插值和擬合。 第一步:適當(dāng)選擇函數(shù)的形式； 第二步:確定函數(shù)的參數(shù)。,設(shè)兩個(gè)變量x與y之間有某種關(guān)系，y=g(x)（很復(fù)雜、或未知或無(wú)封閉形式）。一般可以基于一組實(shí)驗(yàn)觀測(cè)數(shù)據(jù)(xi,yi) ,i=0,1,，構(gòu)造一個(gè)相對(duì)較簡(jiǎn)單的函數(shù)y=f(x

2、)來(lái)近似表示x與y之間的關(guān)系。,一維插 值,引例. 函數(shù)查表問(wèn)題,標(biāo)準(zhǔn)正態(tài)分布函數(shù),標(biāo)準(zhǔn)正態(tài)分布函數(shù)值(1.014)等于多少?,先對(duì)自變量作近似，1.0141.01， 查表得到 (1.01)= 0.8438, 所以 (1.014) (1.01)= 0.8438。 改進(jìn)：函數(shù)值取二者的中點(diǎn)： (1.014) (1.01)+ (1.02)/2 = 0.8438+ 0.8461/2 = 0.84495 還可改進(jìn)嗎？ 利用一個(gè)表格給出的函數(shù)值，近似計(jì)算表格中未給出的函數(shù)值。這實(shí)質(zhì)上就是插值問(wèn)題。,插 值 問(wèn) 題 的 提 法,已知 n+1個(gè)節(jié)點(diǎn),其中,互不相同，不妨設(shè),求任一插值點(diǎn),處的插值,可認(rèn)為節(jié)點(diǎn)

3、滿足函數(shù)關(guān)系y=g(x), 表達(dá)式復(fù)雜, 或無(wú)封閉形式, 或未知。,求 解 插 值 問(wèn) 題 的 基 本 思 路,設(shè)函數(shù)g(x)在n+1個(gè)點(diǎn)x0,x1,xn處的函數(shù)值已知，為 y0,y1,yn 。求一個(gè)分段( 共 n段)線性函數(shù)q(x)，使其滿足： q(xj)=yj, j=0,1,n.,這n+1個(gè)點(diǎn)（ xj, yj) (j=0,1,2,.,n) 為插值節(jié)點(diǎn)，q(x)稱(chēng)為插值函數(shù)。,插值方法,1.分段線性插值,f(x)為被插值函數(shù)。,j=1,2,n。,根據(jù)直線的兩點(diǎn)式方程變形得到q(x)在第j段xj-1, xj 上的表達(dá)式為：,可以證明，,比分段線性插值更光滑。,2. 三次樣條插值,在數(shù)學(xué)上，光滑

4、程度的定量描述是：函數(shù)(曲線)的k階導(dǎo)數(shù)存在且連續(xù)，則稱(chēng)該曲線具有k階光滑性。 光滑性的階次越高，則越光滑。是否存在較低次的分段多項(xiàng)式達(dá)到較高階光滑性的方法？三次樣條插值就是一個(gè)很好的例子。,如果三次樣條函數(shù)S(x)在a,b這兩個(gè)端點(diǎn)處的二階導(dǎo)數(shù)(當(dāng)然只能是單邊導(dǎo)數(shù))等于零，則稱(chēng)其為自然三次樣條插值函數(shù)。,三次樣條函數(shù) 記為S(x), 它是定義在區(qū)間a,b 上的函數(shù)， 滿足以下兩個(gè)條件： 1) S(x) 在每一個(gè)小區(qū)間xi-1,xi上是一個(gè)三次多項(xiàng)式函數(shù) ； 2) 在整個(gè)區(qū)間a,b上，其二階導(dǎo)數(shù)存在且連續(xù)。 即在每個(gè)節(jié)點(diǎn)處的二階導(dǎo)數(shù)連續(xù)。,問(wèn)題：給定n+1個(gè)節(jié)點(diǎn)（x0, y0 )，(x1, y

5、1 ) , (xn, yn), 求一個(gè)三次樣條函數(shù)S(x)，使其滿足： S(xi)=yi,i=0,1,n.,三次樣條插值函數(shù),應(yīng)滿足的條件：,如何確定三次樣條函數(shù)在每一個(gè)小區(qū)間上的三次多項(xiàng)式函數(shù)的系數(shù)呢？,參數(shù)：每個(gè)小段上4個(gè)，n個(gè)小段共計(jì)4n個(gè)。 方程：1) 每個(gè)小段上由給定函數(shù)值得到2個(gè)，n個(gè)小段共計(jì)2n個(gè)； 2) 光滑性要求每一個(gè)內(nèi)部節(jié)點(diǎn)的一階二階導(dǎo)數(shù)連續(xù)，得出其左右導(dǎo)數(shù)相等，因此，每個(gè)節(jié)點(diǎn)產(chǎn)生2個(gè)方程，共計(jì)2(n1) 個(gè) 。現(xiàn)在得到了4n2個(gè)方程，還差兩個(gè)。為此，常用的方法是對(duì)邊界節(jié)點(diǎn)除函數(shù)值外附加要求，這就是所謂的邊界條件。需要兩個(gè)，正好左右兩個(gè)端點(diǎn)各一個(gè)。,3) 周期性邊界條件 在

6、兩個(gè)邊界點(diǎn)的一階導(dǎo)數(shù)值以及二階導(dǎo)數(shù)值均相等：即 S(x0)=S(xn)； S(x0)= S(xn) 。,常用如下三類(lèi)邊界條件： 1) m邊界 條件 給定兩個(gè)邊界節(jié)點(diǎn)的一階導(dǎo)數(shù)值：m0,mn, 即： S(x0)=m0, S(xn)=mn。 2) M邊界條件 給定兩個(gè)邊界節(jié)點(diǎn)的二階導(dǎo)數(shù)值：M0,Mn, 即： S (x0)=M0, S(xn)=Mn。,特別地，當(dāng) M0和Mn都為零時(shí)，稱(chēng)為自然邊界條件。,Method可取： nearest ：最鄰近插值； spline ： 三次樣條插值； cubic ： 立方插值。 缺省時(shí)： 分段線性插值。 注意：所有的插值方法都要求x是單調(diào)的，并且xi不能夠超過(guò)x的

7、范圍。,用MATLAB作插值計(jì)算,一維插值函數(shù)：,yi=interp1(x，y，xi，method),用MATLAB作分段線性插值計(jì)算,分段線性插值: y=interp1(x0,y0,x),例：在1-12的11小時(shí)內(nèi)，每隔1小時(shí)測(cè)量一次溫度，測(cè)得的溫度依次為：5，8，9，15，25，29，31，30，22，25，27，24。試估計(jì)每隔1/10小時(shí)的溫度值。,To MATLAB (xiancha1.m),To MATLAB (yangCha),用MATLAB作三次樣條插值計(jì)算,稱(chēng)為拉格朗日插值基函數(shù)。,已知函數(shù)f(x)在n+1個(gè)點(diǎn)x0,x1,xn處的函數(shù)值為 y0,y1,yn 。求一n次多項(xiàng)式函

8、數(shù)Pn(x)，使其滿足： Pn(xi)=yi,i=0,1,n.,解決此問(wèn)題的拉格朗日插值多項(xiàng)式公式如下,其中Li(x) 為n次多項(xiàng)式：,3. 拉格朗日(Lagrange)多項(xiàng)式插值,例 將0,/2 n等分，在g(x)=cos(x)上取n+1個(gè)節(jié)點(diǎn)，作Pn(x)（取n=1,2) ,計(jì)算Pn(/6),與 cos(/6)比較, 觀察誤差。,解: n=1, (x0,y0)=(0,1), (x1,y1)=(/2,0),/6,P1(x)=y0L0+y1L1=1-2x/, P1(/6)=0.6667 精確值：cos (/6)=0.8660,n=2時(shí): (x0,y0)=(0,1), (x1,y1)=(/4,0

9、.7071), (x2,y2)=(/2,0),P2(x)=y0L0+y1L1+y2L2 =8(x-/4)(x-/2)/2-16x(x-/2)0.7071/2,P2(/6)=0.8508 精確值：cos (/6)=0.8660,/6,思考:是否n越大,插值的誤差就越小?,拉格朗日多項(xiàng)式插值的振蕩現(xiàn)象,To MATLAB (runge),Runge現(xiàn)象:,采用拉格朗日多項(xiàng)式插值：選取不同插值節(jié)點(diǎn)個(gè)數(shù)n+1，其中n為插值多項(xiàng)式的次數(shù)，當(dāng)n分別為2,4,6,8,10時(shí)的插值計(jì)算結(jié)果如下頁(yè)圖.,由圖中可以看出：當(dāng)插值節(jié)點(diǎn)個(gè)數(shù)n+1增大時(shí)，其插值效果變差；當(dāng)n進(jìn)一步增大時(shí)，插值函數(shù)還會(huì)出現(xiàn)劇烈的振蕩現(xiàn)象,

10、從而產(chǎn)生很大的誤差，例如n=10。這種現(xiàn)象稱(chēng)為龍格(Runge) 振蕩現(xiàn)象。發(fā)生這種現(xiàn)象，不是由于這個(gè)被插值函數(shù)的問(wèn)題，因?yàn)檫@樣的函數(shù)太多，俯首皆拾；而是拉格朗日多項(xiàng)式插值本身所固有的問(wèn)題。而對(duì)于三次樣條插值，情況則與此正好相反，隨著n的增大，插值效果變得越來(lái)越好。 所以，在實(shí)際應(yīng)用中，應(yīng)該根據(jù)具體情況，并結(jié)合各種實(shí)際問(wèn)題的特點(diǎn)靈活選用。,分段線性和三次樣條插值（低次多項(xiàng)式插值）：曲線不光滑（三次樣條插值已大有改進(jìn)）；誤差估計(jì)較難（對(duì)三次樣條插值）；收斂性有保證。簡(jiǎn)單實(shí)用，應(yīng)用廣泛。,一維插值方法小結(jié),拉格朗日插值（高次多項(xiàng)式插值）：其插值函數(shù)在整個(gè)區(qū)間上是一個(gè)解析表達(dá)式，便于再次開(kāi)發(fā)利用；曲

11、線光滑；誤差估計(jì)有表達(dá)式；收斂性不能保證（振蕩現(xiàn)象）。用于理論分析，實(shí)際意義不大。,二維插 值,網(wǎng)格節(jié)點(diǎn)插值法： 最鄰近插值； 分片線性插值; 雙線性插值； 雙三次插值。 散點(diǎn)數(shù)據(jù)插值法： 修正Shephard法,引例1：如何繪制山區(qū)地貌圖,要在某山區(qū)方圓大約27平方公里范圍內(nèi)修建一條公路，從山腳出發(fā)經(jīng)過(guò)一個(gè)居民區(qū)，再到達(dá)一個(gè)礦區(qū)。橫向縱向分別每隔400米測(cè)量一次，得到一些地點(diǎn)的高程：(平面區(qū)域0=x=5600,0=y=4800)，首先需作出該山區(qū)的地貌圖和等高線圖。,引例2：船在該海域會(huì)擱淺嗎？,在某海域測(cè)得一些點(diǎn)(x,y)處的水深z由下表給出，船的吃水深度為5英尺，在矩形區(qū)域（75，200

12、）*（-50，150）里的哪些地方船要避免進(jìn)入。,二維插值的提法,已知 mn個(gè)節(jié)點(diǎn),第一種（網(wǎng)格節(jié)點(diǎn)）：,二維插值的提法,第二種（散亂節(jié)點(diǎn)）：,求 解 二 維 插 值 問(wèn) 題 的 基 本 思 路,或,最鄰近插值； 分片線性插值； 雙線性插值； 雙三次插值。,網(wǎng)格節(jié)點(diǎn)插值方法,二維或高維情形的最鄰近插值，與被插值點(diǎn)最鄰近的節(jié)點(diǎn)的函數(shù)值即為所求。 注意：最鄰近插值一般不連續(xù)。具有連續(xù)性的最簡(jiǎn)單的插值是分片線性插值。,1最鄰近插值,將四個(gè)插值點(diǎn)（矩形的四個(gè)頂點(diǎn)）處的函數(shù)值依次簡(jiǎn)記為： f (xi, yj)=f1，f (xi+1, yj)=f2，f (xi+1, yj+1)=f3，f (xi, yj+

13、1)=f4,2分片線性插值,分兩片的函數(shù)表達(dá)式如下： 第一片(下三角形區(qū)域)： (x, y)滿足,插值函數(shù)為：,第二片(上三角形區(qū)域)：(x, y)滿足,插值函數(shù)為：,注意：(x, y)當(dāng)然應(yīng)該是在插值節(jié)點(diǎn)所形成的矩形區(qū)域內(nèi)。顯然，分片線性插值函數(shù)是連續(xù)的；,雙線性插值是一片一片的空間二次曲面構(gòu)成。 雙線性插值函數(shù)的形式如下：,其中有四個(gè)待定系數(shù)，利用該函數(shù)在矩形的四個(gè)頂點(diǎn)（插值節(jié)點(diǎn)）的函數(shù)值，得到四個(gè)代數(shù)方程，正好確定四個(gè)系數(shù)。,3雙線性插值,二維插值:已有程序,Method可?。?nearest 最鄰近插值 linear 雙線性插值 cubic 雙三次插值 缺省時(shí), 雙線性插值 要求x0,

14、y0單調(diào)；x，y可取為矩陣，或x取行向量，y取為列向量，x,y的值分別不能超出x0,y0的范圍。,z=interp2(x0,y0,z0,x,y,method),用MATLAB作網(wǎng)格節(jié)點(diǎn)數(shù)據(jù)的插值,例：測(cè)得平板表面3*5網(wǎng)格點(diǎn)處的溫度分別為： 82 81 80 82 84 79 63 61 65 81 84 84 82 85 86 試作出平板表面的溫度分布曲面z=f(x,y)的圖形。,用MATLAB作網(wǎng)格節(jié)點(diǎn)數(shù)據(jù)的插值,To MATLAB (wenduqm),width=1:5; depth=1:3; temps=82 81 80 82 84;79 63 61 65 81;84 84 82 85

15、 86; mesh(width,depth,temps);pause di=1:.1:3;di=di; wi=1:.1:5; zlin=interp2(width,depth,temps, wi,di,linear); figure(2); mesh(wi,di,zlin);,M文件wenduqm.m,加密數(shù)據(jù)點(diǎn),xlabel(Width of Plate), ylabel(Depth of Plate) zlabel(Degrees Celsius), axis(ij),grid, pause; zlin=interp2(width,depth,temps,wi,di, cubic); fi

16、gure(3); mesh(wi,di,zlin) xlabel(Width of Plate), ylabel(Depth of Plate) zlabel(Degrees Celsius), axis(ij),grid,范例：繪制山區(qū)地貌圖,要在某山區(qū)方圓大約27平方公里范圍內(nèi)修建一條公路，從山腳出發(fā)經(jīng)過(guò)一個(gè)居民區(qū)，再到達(dá)一個(gè)礦區(qū)。橫向縱向分別每隔400米測(cè)量一次，得到一些地點(diǎn)的高程（網(wǎng)格節(jié)點(diǎn)）：(平面區(qū)域0=x=5600,0=y=4800)，首先需作出該山區(qū)的地貌圖和等高線圖。,To MATLAB (shanqu),通過(guò)此例對(duì)雙線性插值方法和雙三次插值方法的插值效果進(jìn)行比較。,返 回,散點(diǎn)數(shù)據(jù)的插值問(wèn)題,二維散點(diǎn)數(shù)據(jù)插值問(wèn)題的提法,用MATLAB作散點(diǎn)數(shù)據(jù)的插值計(jì)算,注意：x0,y0,z0均為向量，長(zhǎng)度相等。 Method可取nearest,linear,cubic,v4; linear是缺省值。,z=griddata(x0,y0,z0,x,y,method),在某海域測(cè)得一些點(diǎn)(x,y)處的水深z由下表給出，船的吃水深度為5英尺，在矩形區(qū)域（75，200）*（-50，15