數(shù)學(xué)建模數(shù)值建模

1、1擬合與插值擬合與插值第第5 5章章 數(shù)值分析法建模數(shù)值分析法建模作業(yè)：習(xí)題5：1，2上機實驗3：習(xí)題5：32擬合函數(shù)v已知函數(shù)f在離散樣本點集離散樣本點集M上的函數(shù)值y1,y2,yn，擬合是指選取含參數(shù)的函數(shù)形式f(x;a1,a2 ,am) ，通過調(diào)整待定參數(shù)a1,a2 ,am, 使得在樣本點集上的函樣本點集上的函數(shù)計算值與實際值誤差最小數(shù)計算值與實際值誤差最小。v如果選取線性函數(shù)f，則稱線性擬合或者線性回歸，否則稱為非線性擬合或者非線性回歸。v擬合特點：函數(shù)不要求經(jīng)過樣本數(shù)據(jù)點。因此當(dāng)測量值與真實值有誤差時，一般用數(shù)據(jù)擬合。3插值函數(shù)v插值是指已知某函數(shù)的在離散點集離散點集M上的函數(shù)值或者

2、導(dǎo)數(shù)信息，通過求解函數(shù)形式，使得函數(shù)在給定離散點上滿足一定約束條件。v如果約束條件中只有函數(shù)值的約束，叫作Lagrange插值，否則叫作Hermite插值 v插值函數(shù)特點：必須經(jīng)過所有樣本數(shù)據(jù)點。因此當(dāng)測量值是準(zhǔn)確的，沒有誤差時，一般用插值。4曲線擬合步驟1.確定經(jīng)驗公式形式確定經(jīng)驗公式形式2.計算經(jīng)驗公式中的參數(shù)計算經(jīng)驗公式中的參數(shù)3.檢驗經(jīng)驗公式有效性檢驗經(jīng)驗公式有效性5 利用已知的結(jié)論曲線改直技術(shù)作散點圖法多項式近似6v最小二乘法原理就是找一組參數(shù) 使最小二乘誤差 最小。12,naa aa2( )Syf x數(shù)學(xué)模型TTA AaA y曲線擬合問題最常用的解法曲線擬合問題最常用的解法線性最小

3、二乘法線性最小二乘法的基本思路的基本思路選定一組基函數(shù)選定一組基函數(shù) r1(x), r2(x), rn(x), nm, 選取擬合函數(shù)形式為選取擬合函數(shù)形式為 f(x)=a1r1(x)+a2r2(x)+ +anrn(x) ，其中其中 a1,a2, an 為待定系數(shù)。為待定系數(shù)。正規(guī)方程正規(guī)方程特例：多項式擬合模型求解1TTaA AA yMatlab求解：求解：aAy7v差分與差商概念差分與差商概念v一階前向差分 v二階前向差分 vvm階前向差分 v一階差商 v二階差商 vvm階差商 )()(11kkkkkxfxfyyy21kkkyyy , 2 , 1,111myyykmkmkmkkkkkxxyx

4、xf11,kkkkkkkkkxxxxfxxfxxxf212121,kmkmkkkmkkkmkkkxxxxxfxxxfxxxf,112118差商與導(dǎo)數(shù)的聯(lián)系v若y=f(x)在a,b上m階可微，且 則bxxxxxamkkkkk321!)(,)(21mfxxxxfmmkkkk),(mkkxxv若節(jié)點為等距分布時，有 ，h為節(jié)點距，則khxxk0)(,!.)(1mmmkkkmkmfhxxxfhmyv定擬合多項式的最佳階數(shù)： 差分表/插商表波動最小原則v等距節(jié)點用差分表，不等距節(jié)點用差商表9差分表差分表0123nxxxxx0123nf xf xf xf xf x0121nyyyy202122nyyy30

5、33nyy0ny23nxfffff10.xk f(xk) 一階差商一階差商 二階差商二階差商 三階差商三階差商 n 階差商階差商差商表差商表0123nxxxxx0123nf xf xf xf xf x0112231,nnf xxf x xf xxf xx01212321,nnnf xx xf x xxf xxx0123321 , , ,nnnnf x x xxf xxxx01 , , , nf x xx11與擬合有關(guān)的與擬合有關(guān)的MATLAB 函數(shù)函數(shù)vpolyfit:多項式擬合多項式擬合vpoly2sym:由多項式系數(shù)向量得多項式符號表達式由多項式系數(shù)向量得多項式符號表達式vpolyval:

6、計算多項式函數(shù)在指定處的函數(shù)值計算多項式函數(shù)在指定處的函數(shù)值vpoly:計算過固定點的多項式計算過固定點的多項式vlsqcurvefit, lsqnonlin非線性最小二乘擬合非線性最小二乘擬合vfit,fittype一般非線性擬合一般非線性擬合vcftool 擬合工具箱擬合工具箱 作多項式作多項式f(x)=anxn+ +a1x+a0擬合擬合, ,可利用可利用MatlabMatlab命令命令:a=polyfit(x,y,n)yy=polyval(a,xx)yy=polyval(a,xx). .計算多項式計算多項式a a在在xxxx處的值處的值12 當(dāng)精確函數(shù) y = f(x) 非常復(fù)雜或未知時

7、，在區(qū)間a,b上一系列節(jié)點 x0 xm 處測得函數(shù)值 y0 = f(x0), , ym = f(xm)，由此構(gòu)造一個簡單易算的近似函數(shù) g(x) f(x)，滿足條件 g(xj) = f(xj) (j = 0, m) (*)這個問題稱為“插值問題”插值問題g(x) 稱為稱為f(x) 的的插值函數(shù)插值函數(shù),一般取一般取多項式函數(shù)多項式函數(shù)。x0 xm稱為插值節(jié)點稱為插值節(jié)點, 插值節(jié)點互不相同插值節(jié)點互不相同條件條件(*)稱為稱為插值條件插值條件，區(qū)間，區(qū)間a,b稱為稱為插值區(qū)間插值區(qū)間13x0 x0 x2xm-1xm xf(x)g(x)精確函數(shù) y = f(x) 非常復(fù)雜或未知g(x) 稱為稱為

8、f(x) 的的插值函數(shù)插值函數(shù) 節(jié)點節(jié)點 x0 xm稱為插值節(jié)點稱為插值節(jié)點14基本思想：基本思想：在n次多項式空間Pn中找一組合適的基函數(shù)中找一組合適的基函數(shù) 0(x), 1(x), n(x),使使pn(x)=a0 0(x) +a1 1(x) +an n(x)不同的基函數(shù)的選取導(dǎo)致不同的不同的基函數(shù)的選取導(dǎo)致不同的插值方法插值方法 無論是從理論和計算的角度,還是從應(yīng)用的角度看,多項式都是最簡單的函數(shù),因此,多項式插值是最基本的插值方法15根據(jù)插值條件，應(yīng)有0111)(axaxaxaxLnnnnnniyxLiin, 2 , 1 , 0,)(1001112111011,1nnnnnnnnnnnx

9、xayayxxAaYyaxx A aY0det()0ijj i nAxx 拉格朗日拉格朗日(Lagrange)插值插值16 稱為拉格朗日插值基函數(shù)拉格朗日插值基函數(shù)。0( )( )nniiiL xl xy拉格朗日插值多項式公式也可以通過構(gòu)造插值基函數(shù)方法直接得到:其中l(wèi)i(x) 為n次多項式 0,njjjixxx記拉格朗日拉格朗日(Lagrange)插值插值 0,1, 2,( )njjjjiiiixlxxxinxxx17拉格朗日拉格朗日(Lagrange)插值插值特別地特別地:兩點一次兩點一次(線性線性)插值多項式插值多項式: 011010110 xxxxLxyyxxxx三點二次三點二次(拋物

10、線拋物線)插值多項式插值多項式: 1202012012010210122021x xx xx xx xx xx xL xyyyxxxxx xx xxxxx 18 拉格朗日多項式插值的這種振蕩現(xiàn)象叫 Runge現(xiàn)象現(xiàn)象21( ),551f xxx 采用拉格朗日多項式插值：選取不同插值節(jié)點個數(shù)n+1，其中n為插值多項式的次數(shù)，當(dāng)n分別取2,4,6,8,10時，繪出插值結(jié)果圖形.例例1ch.m1ch.m-505-1.5-1-0.500.511.52y=1/(1+x2)n=2n=4n=6n=8n=10( )( )( )?nnnRxLxf x 19 Newtons Interpolation基本原理基本

11、原理Lagrange Lagrange 插值雖然易算，但若要增加一個節(jié)點時，全部插值雖然易算，但若要增加一個節(jié)點時，全部基函數(shù)基函數(shù) li( (x x) ) 都需要重新計算都需要重新計算能否重新在能否重新在Pn中尋找新的基中尋找新的基函數(shù)？希望每加一個節(jié)點時，函數(shù)？希望每加一個節(jié)點時，只附加一項只附加一項上去即可上去即可2001020101( )()()() .().()nnnN xAA x xA x xx xA x xx x 選取選取1,1,x - x0 , (x - x0)(x - x1),(x-x0)(x-x1) (x-xn-1)構(gòu)成構(gòu)成Pn的一組基函數(shù)的一組基函數(shù) Newtons In

12、terpolation基本原理基本原理利用插值條件利用插值條件Nn(xj)=f(xj), j=0,1，n代入上式，代入上式，得關(guān)于得關(guān)于Ak (k=0,1，n)的線性代數(shù)方程組的線性代數(shù)方程組0010111000100()10()1()()nninniAf xxxAf xxxxxAf x當(dāng)當(dāng)xj 互異時，系數(shù)矩陣互異時，系數(shù)矩陣非奇異，方程有唯一解非奇異，方程有唯一解0100,0mi mkkkkmmiif xAf x xxknxx21Lagrange插值與Newton插值的異同點v兩者都是通過給定n+1個互異的插值節(jié)點，求一條n次多項式曲線近似地表示待插值的函數(shù)曲線vLagrange插值和Ne

13、wton法插值的結(jié)果和余項都是一致的，因為都是利用n次多項式插值v區(qū)別： Lagrange插值法在求每個函數(shù)的時候要用到所有結(jié)點，因此如果需要再多加一個結(jié)點的話，需要重新求出函數(shù)才可，而這需要大工作量，于是數(shù)學(xué)家們就發(fā)明了Newton法。vLagrange插值法是通過構(gòu)造n+1個n次基函數(shù)，作線性組合而得到 Newton法插值是通過求各階差商，遞推得到公式 f(x)=f(x0)+(x-x0)fx0,x1+(x-x0)(x-x1)fx0,x1,x2.+(x-x0).(x-xn-1)fx0,x1.xn22 在數(shù)學(xué)上，光滑程度的定量描述是：函數(shù)(曲線)的k階導(dǎo)數(shù)存在且連續(xù)，則稱該曲線具有k階光滑性。

14、 光滑性的階次越高，則越光滑。是否存在較低次的 分 段 多 項 式 達 到 較 高 階 光 滑 性 的 方 法 ？三次樣條插值就是一個很好的例子。三次樣條插值三次樣條插值23Cubic Spline LagrangeInterpolation Interpolation24 , 1,),()(1nixxxxsxSiii,)()3), 1 ,0()()2), 1()()10223niiiiiiixxCxSniyxSnidxcxbxaxs(0)(0)(0)(0) ,(1,1)(0)(0)iiiiiiS xS xSxSxinSxSx04)()()0 (nSxSx自 然 邊 界 條 件 ）)(,)4)

15、3)2xSdcbaiiiif( (x) )為被插值函數(shù)為被插值函數(shù)插值條件連接條件25yi=interp1(x，y，xi，method)插值方法插值方法被插值點被插值點插值節(jié)點插值節(jié)點xixi處的插處的插值結(jié)果值結(jié)果nearest ：最鄰近插值：最鄰近插值linear ： 線性插值；線性插值；spline ： 三次樣條插值；三次樣條插值；cubic ： 立方插值。立方插值。缺省時：缺省時： 線性插值。線性插值。 插值方法要求插值方法要求x x是單調(diào)的，且是單調(diào)的，且xi不能夠超過不能夠超過x的范圍。的范圍。用用MATLABMATLAB作一維插值計算作一維插值計算擴展材料261. 1. 拉格朗日

16、插值拉格朗日插值: :自編程序自編程序2. 2. 分段線性插值分段線性插值: :已有程序已有程序 y=interp1(x0,y0,x)3. 3. 三次樣條插值三次樣條插值: :已有程序已有程序 y=interp1(x0,y0,x,spline) 或或 y=spline(x0,y0,x) 注意：幾乎所有的插值方法都注意：幾乎所有的插值方法都要求要求x0單調(diào)單調(diào)，并且，并且x不不能夠超過能夠超過x0的范圍的范圍。interp2,interp3,interpn多元函數(shù)插值多元函數(shù)插值27pp=spline (x,y) ：樣條函數(shù)的表示（樣條函數(shù)的表示（PP結(jié)構(gòu)）結(jié)構(gòu)）pp=cscvn(x,y): 自

17、然樣條函數(shù)的表示（不要求x單調(diào)?。﹜i=ppval(pp,xx)或 yi=spline(x,y,xx)：樣條函數(shù)求值樣條函數(shù)求值fprime = fnder(pp), fprime = fnder(pp,dorder) 樣條函數(shù)求導(dǎo)樣條函數(shù)求導(dǎo)inters = fnint(pp) , intgrf = fnint(pp,value)樣條函數(shù)積分樣條函數(shù)積分fnplt(pp)樣條函數(shù)繪圖樣條函數(shù)繪圖樣條函數(shù)的相關(guān)樣條函數(shù)的相關(guān)MATLABMATLAB命令命令28v用程控銑床加工機翼斷面的下輪廓線時用程控銑床加工機翼斷面的下輪廓線時v每一刀只能沿每一刀只能沿x x方向和方向和y y方向走非常小的一

18、步。方向走非常小的一步。v表表3-13-1給出了下輪廓線上的部分數(shù)據(jù)給出了下輪廓線上的部分數(shù)據(jù)v但工藝要求銑床沿但工藝要求銑床沿x x方向每次只能移動方向每次只能移動0.10.1單位單位. . v這時需求出當(dāng)這時需求出當(dāng)x x坐標(biāo)每改變坐標(biāo)每改變0.10.1單位時的單位時的y y坐標(biāo)。坐標(biāo)。v試完成加工所需的數(shù)據(jù)，畫出曲線試完成加工所需的數(shù)據(jù)，畫出曲線. .29x0=0 3 5 7 9 11 12 13 14 15 ;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;y=interp1(x0,y0,x,spline);plot(x0,y0

19、,k+,x,y,r)grid on05101500.511.522.530 例例1 1：每隔：每隔1 1小時測量一次溫度，測得的小時測量一次溫度，測得的1212個溫度個溫度依次為：依次為：5 5，8 8，9 9，1515，2525，2929，3131，3030，2222，2525，2727，2424。試估計每隔。試估計每隔6 6分鐘的溫度值。分鐘的溫度值。Temp.mhours=1:12;temps=5 8 9 15 25 29 31 30 22 25 27 24;h=1:0.1:12;t=interp1(hours,temps,h,spline); (直接輸出數(shù)據(jù)將是很多的)plot(hou

20、rs,temps,+,h,t,hours,temps,r:) %作圖xlabel(Hour),ylabel(Degrees Celsius)31 要求要求x0,y0 x0,y0單調(diào)；單調(diào)；x x，y y可取為矩陣，或可取為矩陣，或x x取行向量，取行向量，y y取為列向量，取為列向量，x,yx,y的值分別的值分別不能超出不能超出x0,y0 x0,y0的范圍。的范圍。z=interp2(x0,y0,z0,x,y,method)被插值點插值方法插值節(jié)點被插值點的函數(shù)值nearestnearest 最鄰近插值最鄰近插值linearlinear 雙線性插值雙線性插值cubiccubic 雙三次插值雙三

21、次插值缺省時缺省時, , 雙線性插值雙線性插值用用MATLABMATLAB作網(wǎng)格節(jié)點數(shù)據(jù)的插值作網(wǎng)格節(jié)點數(shù)據(jù)的插值( (二維二維) )32 cz =griddata（x，y，z，cx，cy，method） 要求要求cxcx取行向量，取行向量，cycy取為列向量。取為列向量。被插值點插值方法插值節(jié)點被插值點的函數(shù)值nearestnearest 最鄰近插值最鄰近插值linearlinear 雙線性插值雙線性插值cubiccubic 雙三次插值雙三次插值v4 Matlab提供的插值方法提供的插值方法缺省時缺省時, , 雙線性插值雙線性插值用用MATLABMATLAB作散點數(shù)據(jù)的插值計算作散點數(shù)據(jù)的插值計算33v某海域上頻繁地有各種噸位的船只經(jīng)過。v為保證船只的航行安全，有關(guān)機構(gòu)在低潮時對水深進行了測量,下表是測量數(shù)據(jù)：表3 水道水深的測量數(shù)據(jù)x 129.0 140.0 103.5 88.0 185.5 195.0 105.5y 7.5 141.5 23.0 147.0 22.5 137.5 85.5z 4 8 6 8 6 8 8x 157.5 107.5 77.0 81.0 162.0 162.0 117.5y -6.5 -81.0 3.0 56.5 -66.5 84.0 -33.5z 9 9 8 8 9 4 934v其中（其中