微積分模型ppt

1、模型分析和建立模型分析和建立先看如下三個(gè)例子先看如下三個(gè)例子1. 衛(wèi)星軌道的長(zhǎng)度2. 射擊命中概率3. 人口增長(zhǎng)模型1、衛(wèi)星軌道的長(zhǎng)度 人造地球衛(wèi)星的軌道可以視為平面上的橢圓，中國第一顆人造衛(wèi)星近地點(diǎn)距離地球表面439km，遠(yuǎn)地點(diǎn)距離地球表面2384km，地球半徑為6371m，試求該衛(wèi)星的軌道長(zhǎng)度？分析：2、射擊命中概率 炮彈射擊的目標(biāo)為一正橢圓形區(qū)域，當(dāng)瞄準(zhǔn)目標(biāo)的中心發(fā)射時(shí)，在縱多因素的影響下，彈著點(diǎn)與目標(biāo)中心有隨機(jī)偏差。可以合理地假設(shè)彈著點(diǎn)圍繞中心呈二維正態(tài)分布，且偏差在x方向和y方向相互獨(dú)立。若橢圓區(qū)域在x方向半軸長(zhǎng)120m,y方向半軸80m，設(shè)彈著點(diǎn)偏差的均方差在x方向和y方向均為10

2、0m,試求 炮 彈 落 在 橢 圓 形 區(qū) 域 內(nèi) 的 概 率 ？分析：3、人口增長(zhǎng)率20世紀(jì)美國人口統(tǒng)計(jì)數(shù)據(jù)如表6.5所示，計(jì)算表中這些年份的人口增長(zhǎng)率。又已知某地區(qū)20世紀(jì)70年代的人口增長(zhǎng)率如表6.6所示，且1970年人口為210萬，試估計(jì)該地區(qū)1980年的人口？分析:小小結(jié)結(jié)模型求解u 數(shù)值微分相關(guān)知識(shí)回顧u 數(shù)值積分相關(guān)知識(shí)u Matlab 相關(guān)命令u 程序?qū)崿F(xiàn)u 數(shù)值微分相關(guān)知識(shí)回顧hhxfhxfhhxfxfhxfhxfxfhhh2)()(lim)()(lim)()(lim)( 0001. 函數(shù)f(x)以離散點(diǎn)列給出時(shí)，而要求我們給出導(dǎo)數(shù)值，2. 函數(shù)f(x)過于復(fù)雜這兩種情況都要

3、求我們用數(shù)值的方法求函數(shù)的導(dǎo)數(shù)值微積分中，關(guān)于導(dǎo)數(shù)的定義如下：自然，而又簡(jiǎn)單的方法就是，取極限的近似值，即差商 向前差商公式()( )( ),0,1iiif xhf xfxih 向后差商公式 中心差商公式( )()( ),1,2iiif xf xhfxih111()()()2f xhf xhfxh給定點(diǎn)列20)(,(iiixfx且2110,xxxxh)( ),( ),( 012xfxfxf求兩點(diǎn)差商公式：兩點(diǎn)差商公式：三點(diǎn)差商公式三點(diǎn)差商公式0000012201()3 ()4 ()(2 )21 3 ()4 ()()fxf xf xhf xhhf xf xf xxx1002201()()(2 )

4、21 ()()fxf xf xhhf xf xxx2000012201()()4 ()3 (2 )21 =()4 ()3 ()fxf xf xhf xhhf xf xf xxxu 數(shù)值積分相關(guān)知識(shí)回顧)()()(aFbFdxxfba關(guān)于積分，有Newton-Leibniz公式但是，在很多情況下，還是要數(shù)值積分：1、函數(shù)有離散數(shù)據(jù)組成2、F(x)求不出3、F(x)非常復(fù)雜定義數(shù)值積分如下：是離散點(diǎn)上的函數(shù)值的線性組合)()(0iniinxfafI稱為求積系數(shù)，與f(x)無關(guān)，與積分區(qū)間和積分點(diǎn)有關(guān)矩形公式矩形公式左矩形公式左矩形公式( )() ( )baf x dxba f a右矩形公式右矩形公

5、式( )() ( )baf x dxb a f b中矩形公式中矩形公式( )() ()2babaf x dxba f梯形公式梯形公式( )( )( )2babafx dxf af b辛普森辛普森(Simpson)公式公式( ) ( )4 ()( )62babaabf x dxf aff b復(fù)合梯形公式復(fù)合梯形公式11011101( )( ) ( )2()( )22kknbxaxknnkkkkkfx dxfx dxhhfxfxf afxf b復(fù)合辛普森復(fù)合辛普森(Simpson)公式公式111012112( ) ( )4()2()( )6=2nnbkakkkkkkbaf x dxf af xf

6、xf bxxx其中 將將 分為分為 等份，步長(zhǎng)等份，步長(zhǎng) ，分點(diǎn)，分點(diǎn), a bbahnn,0,1,kxakh knMatlab相關(guān)命令 x=linspace(0,pi,50);y=sin(x);R1=trapz(x,y)；% 梯形公式計(jì)算R2=quad(sin,0,pi,1.0e06)；% 辛普森公式計(jì)算結(jié)果R1=1.999314849324062R2=1.9999934965349640sin( )x dx程序?qū)崿F(xiàn)模型實(shí)例講解 森林救火問題 水箱的水流問題 森林救火問題森林救火問題森林失火后，要確定派出消防隊(duì)員的數(shù)量。森林失火后，要確定派出消防隊(duì)員的數(shù)量。隊(duì)員多，森林損失小，救援費(fèi)用大；隊(duì)

7、員多，森林損失小，救援費(fèi)用大；隊(duì)員少，森林損失大，救援費(fèi)用小。隊(duì)員少，森林損失大，救援費(fèi)用小。所以需要綜合考慮森林損失費(fèi)和救援費(fèi)與消防隊(duì)員所以需要綜合考慮森林損失費(fèi)和救援費(fèi)與消防隊(duì)員之間的關(guān)系，以總費(fèi)用最小來決定派出隊(duì)員的數(shù)量。之間的關(guān)系，以總費(fèi)用最小來決定派出隊(duì)員的數(shù)量。問題問題分析分析問題問題記隊(duì)員人數(shù)記隊(duì)員人數(shù)x, 失火時(shí)刻失火時(shí)刻t=0, 開始救火時(shí)刻開始救火時(shí)刻t1, 滅火時(shí)刻滅火時(shí)刻t2, 時(shí)刻時(shí)刻t森林燒毀面積森林燒毀面積B(t). 損失費(fèi)損失費(fèi)f1(x)是是x的減函數(shù)的減函數(shù), 由燒毀面積由燒毀面積B(t2)決定決定. 救援費(fèi)救援費(fèi)f2(x)是是x的增函數(shù)的增函數(shù), 由隊(duì)員人數(shù)

8、和救火時(shí)間決定由隊(duì)員人數(shù)和救火時(shí)間決定.找到找到恰當(dāng)?shù)那‘?dāng)?shù)膞，使，使f1(x), f2(x)之和最小之和最小 關(guān)鍵是對(duì)關(guān)鍵是對(duì)B(t)作出合理的簡(jiǎn)化假設(shè)作出合理的簡(jiǎn)化假設(shè).問題問題分析分析失火時(shí)刻失火時(shí)刻t=0, 開始救火時(shí)刻開始救火時(shí)刻t1, 滅火時(shí)刻滅火時(shí)刻t2, 畫出時(shí)刻畫出時(shí)刻 t 森林燒毀面積森林燒毀面積B(t)的大致圖形的大致圖形t1t20tBB(t2)分析分析B(t)比較困難比較困難,轉(zhuǎn)而討論森林燒毀轉(zhuǎn)而討論森林燒毀速度速度dB/dt. （其中（其中dB/dt表示表示單位時(shí)間燒毀的面單位時(shí)間燒毀的面積）積）模型假設(shè)模型假設(shè) 3）f1(x)與與B(t2)成正比，系數(shù)成正比，系數(shù)c

9、1 (燒毀單位面積損失費(fèi)）燒毀單位面積損失費(fèi)） 1）0 t t1, dB/dt 與與 t成正比，系數(shù)成正比，系數(shù) (火勢(shì)蔓延速度）火勢(shì)蔓延速度） 2）t1 t t2, 降為降為 - x ( 為隊(duì)員的平均滅火為隊(duì)員的平均滅火速度）速度） 4）每個(gè)）每個(gè)隊(duì)員的單位時(shí)間滅火費(fèi)用隊(duì)員的單位時(shí)間滅火費(fèi)用c2, 一次性費(fèi)用一次性費(fèi)用c3 (運(yùn)運(yùn)送消防隊(duì)員和器材等一次性支出送消防隊(duì)員和器材等一次性支出)假設(shè)假設(shè)1）的解釋的解釋 rB火勢(shì)以失火點(diǎn)為中心，火勢(shì)以失火點(diǎn)為中心，均勻向四周呈圓形蔓延，均勻向四周呈圓形蔓延，半徑半徑 r與與 t 成正比成正比面積面積 B與與 t2成正比，成正比， dB/dt與與 t成

10、正比成正比.xbtt12202)()(tdttBtB模型建立模型建立dtdBb0t1tt2x假設(shè)假設(shè)1）,1tbxcttxcxftBcxf31222211)()(),()(目標(biāo)函數(shù)目標(biāo)函數(shù)總費(fèi)用總費(fèi)用)()()(21xfxfxC假設(shè)假設(shè)3）4）xttt112假設(shè)假設(shè)2）)(222212212xttbt0dxdCxcxxtcxtctcxC3122121211)(22)(模型建立模型建立目標(biāo)函數(shù)目標(biāo)函數(shù)總費(fèi)用總費(fèi)用模型求解模型求解求求 x使使 C(x)最小最小231221122ctctcx結(jié)果解釋結(jié)果解釋 / 是火勢(shì)不繼續(xù)蔓延的最少隊(duì)員數(shù)是火勢(shì)不繼續(xù)蔓延的最少隊(duì)員數(shù)dtdBb0t1t2tx其中其中

11、 c1,c2,c3, t1, , 為已知參數(shù)為已知參數(shù)模型模型應(yīng)用應(yīng)用c1,c2,c3已知已知, t1可估計(jì)可估計(jì), c2 x c1, t1, x c3 , x 結(jié)果結(jié)果解釋解釋231221122ctctcxc1燒毀單位面積損失費(fèi)燒毀單位面積損失費(fèi), c2每個(gè)每個(gè)隊(duì)員單位時(shí)間滅火費(fèi)隊(duì)員單位時(shí)間滅火費(fèi), c3每個(gè)每個(gè)隊(duì)員一次性費(fèi)用隊(duì)員一次性費(fèi)用, t1開始救火時(shí)刻開始救火時(shí)刻, 火火勢(shì)蔓延速度勢(shì)蔓延速度, 每個(gè)每個(gè)隊(duì)員平均滅火隊(duì)員平均滅火速度速度.為什么為什么? ? , 可可設(shè)置一系列數(shù)值設(shè)置一系列數(shù)值由模型決定隊(duì)員數(shù)量由模型決定隊(duì)員數(shù)量x練習(xí)題 在森林救火模型中，如果考慮消防隊(duì)員的滅火速度 與

12、開始救火時(shí)的火勢(shì)b有關(guān)，試假設(shè)一個(gè)合理的函數(shù)關(guān)系，重新求解模型。程序?qū)崿F(xiàn) 美國某州的用水管理機(jī)構(gòu)要求各社區(qū)提供以每小時(shí)多少加美國某州的用水管理機(jī)構(gòu)要求各社區(qū)提供以每小時(shí)多少加侖計(jì)的用水量以及每天所用的總水量。許多社區(qū)沒有測(cè)量流入侖計(jì)的用水量以及每天所用的總水量。許多社區(qū)沒有測(cè)量流入或流出水塔水量的裝置，只能代之以每小時(shí)測(cè)量水塔中的水位，或流出水塔水量的裝置，只能代之以每小時(shí)測(cè)量水塔中的水位，其誤差不超過其誤差不超過5%。需要注意的是，當(dāng)水塔中的水位下降到最低。需要注意的是，當(dāng)水塔中的水位下降到最低水位水位L時(shí)，水泵就自動(dòng)向水塔輸水直到最高水位時(shí)，水泵就自動(dòng)向水塔輸水直到最高水位H，此期間不能，

13、此期間不能測(cè)量水泵的供水量，因此，當(dāng)水泵正在輸水時(shí)不容易建立水塔測(cè)量水泵的供水量，因此，當(dāng)水泵正在輸水時(shí)不容易建立水塔中水位和用水量之間的關(guān)系。水泵每天輸水一次或兩次，每次中水位和用水量之間的關(guān)系。水泵每天輸水一次或兩次，每次約約2小時(shí)。小時(shí)。水箱的水流問題水箱的水流問題 試估計(jì)任何時(shí)刻（包括水泵正在輸水時(shí)間）從水塔試估計(jì)任何時(shí)刻（包括水泵正在輸水時(shí)間）從水塔流出的水流量流出的水流量f(t)，并估計(jì)一天的總用水量。已知該水塔，并估計(jì)一天的總用水量。已知該水塔是一個(gè)高為是一個(gè)高為40ft（英尺），直徑為（英尺），直徑為57ft（英尺）的正圓（英尺）的正圓柱，表柱，表6-1給出了某個(gè)小鎮(zhèn)一天水塔水

14、位的真實(shí)數(shù)據(jù)，給出了某個(gè)小鎮(zhèn)一天水塔水位的真實(shí)數(shù)據(jù)，水位降至約水位降至約27.00ft水泵開始工作，水位升到水泵開始工作，水位升到35.50ft停止停止工作。（注：工作。（注：1ft（英尺）（英尺）=0.3024m（米）（米）表表6.1 某小鎮(zhèn)某天水塔水位某小鎮(zhèn)某天水塔水位時(shí)間時(shí)間/s水位水位/0.01ft時(shí)間時(shí)間/s水位水位/0.01ft031754663633503316311049953326066353054539363167106192994572543087139372947605743012179212892645542927212402850685352842252232795

15、71854276728543275275021269732284269779254水泵開動(dòng)水泵開動(dòng)35932水泵開動(dòng)水泵開動(dòng)82649水泵開動(dòng)水泵開動(dòng)39332水泵開動(dòng)水泵開動(dòng)8596834753943535508995333974331834459327033402.2 問題分析問題分析 模型中所指流量可視為單位時(shí)間內(nèi)流出水的模型中所指流量可視為單位時(shí)間內(nèi)流出水的體積體積。由于水塔是正圓柱形，橫截面積是常數(shù)，所以在水泵由于水塔是正圓柱形，橫截面積是常數(shù)，所以在水泵不工作時(shí)段，流量很容易根據(jù)水位相對(duì)時(shí)間的變化算不工作時(shí)段，流量很容易根據(jù)水位相對(duì)時(shí)間的變化算出。問題的難點(diǎn)在于如何估計(jì)水泵供水時(shí)段

16、的流量。出。問題的難點(diǎn)在于如何估計(jì)水泵供水時(shí)段的流量。 水泵供水時(shí)段的流量只能靠供水時(shí)段前后的流量經(jīng)水泵供水時(shí)段的流量只能靠供水時(shí)段前后的流量經(jīng)插值或擬合得到。作為用于插值或擬合的原始數(shù)據(jù)，我插值或擬合得到。作為用于插值或擬合的原始數(shù)據(jù)，我們希望水泵不工作時(shí)段的流量越準(zhǔn)確越好。這此流量大們希望水泵不工作時(shí)段的流量越準(zhǔn)確越好。這此流量大體上可由兩種方法計(jì)算，一是直接對(duì)表體上可由兩種方法計(jì)算，一是直接對(duì)表6-1中的水量用中的水量用數(shù)值微分算出各時(shí)段的流量，用它們擬合其它時(shí)刻或連數(shù)值微分算出各時(shí)段的流量，用它們擬合其它時(shí)刻或連續(xù)時(shí)間的流量；二是先用表中數(shù)據(jù)擬合水位一時(shí)間函數(shù)，續(xù)時(shí)間的流量；二是先用表

17、中數(shù)據(jù)擬合水位一時(shí)間函數(shù)，求導(dǎo)數(shù)即可得到連續(xù)時(shí)間的流量。求導(dǎo)數(shù)即可得到連續(xù)時(shí)間的流量。 有了任何時(shí)刻的流量，就不難計(jì)算一天的總用水有了任何時(shí)刻的流量，就不難計(jì)算一天的總用水量。其實(shí)，水泵不工作時(shí)段的用水量可以由測(cè)量記錄量。其實(shí)，水泵不工作時(shí)段的用水量可以由測(cè)量記錄直接得到，由表直接得到，由表6-1中下降水位乘以水搭的截面積就是中下降水位乘以水搭的截面積就是這一時(shí)段的用水量。這個(gè)數(shù)值可以用來檢驗(yàn)數(shù)據(jù)插值這一時(shí)段的用水量。這個(gè)數(shù)值可以用來檢驗(yàn)數(shù)據(jù)插值或擬合的結(jié)果?；驍M合的結(jié)果。 在具體給出本問題的解答之前，先介紹一個(gè)簡(jiǎn)單在具體給出本問題的解答之前，先介紹一個(gè)簡(jiǎn)單的數(shù)據(jù)插值方法。的數(shù)據(jù)插值方法。2.

18、3 拉格朗日插值拉格朗日插值1、線性插值、線性插值)(xf假設(shè)已知假設(shè)已知在區(qū)間在區(qū)間,ba上的兩個(gè)結(jié)點(diǎn)上的兩個(gè)結(jié)點(diǎn)10, xx和它們的函數(shù)值和它們的函數(shù)值 10 xxx )()()(10 xfxfxf求一個(gè)一次多項(xiàng)式求一個(gè)一次多項(xiàng)式xaaxP101)(，使得多項(xiàng)式，使得多項(xiàng)式)(1xP在結(jié)點(diǎn)上滿足條件在結(jié)點(diǎn)上滿足條件 1, 0),()(1ixfxPii 這種插值方法稱為線性插值方法（也稱兩點(diǎn)插值）。這種插值方法稱為線性插值方法（也稱兩點(diǎn)插值）。 可以求出：可以求出： )()()(101001011xfxxxxxfxxxxxP2、拋物插值、拋物插值)(xf已知已知在區(qū)間在區(qū)間,ba上的三個(gè)結(jié)點(diǎn)

19、上的三個(gè)結(jié)點(diǎn)210,xxx和它們的函數(shù)值和它們的函數(shù)值 210 xxxx )()()()(210 xfxfxfxf求一個(gè)次數(shù)不超過求一個(gè)次數(shù)不超過2的多項(xiàng)式的多項(xiàng)式22102)(xaxaaxP，使得它在結(jié)點(diǎn)上滿足條件，使得它在結(jié)點(diǎn)上滿足條件2, 1, 0),()(2ixfxPii這種插值方法稱為拋物線插值法，這種插值方法稱為拋物線插值法， 可求出：可求出：)()( )()( )()()( )()( )()(121012002010212xfxxxxxxxxxfxxxxxxxxxP)()( )()( )(2021210 xfxxxxxxxx3、n次拉格朗日插值次拉格朗日插值,ba假設(shè)取區(qū)間假設(shè)取

20、區(qū)間 上的上的1n個(gè)結(jié)點(diǎn)個(gè)結(jié)點(diǎn)nxxx,10，并且已知函數(shù)，并且已知函數(shù))(xf在這此點(diǎn)的函數(shù)值在這此點(diǎn)的函數(shù)值 nxxxxx210 )()()()()(210nxfxfxfxfxf現(xiàn)在求一個(gè)次數(shù)不超過現(xiàn)在求一個(gè)次數(shù)不超過n的多項(xiàng)式的多項(xiàng)式nnnxaxaxaaxP2210)(，使得滿足條件，使得滿足條件 nixfxPiin, 1, 0),()( 這種插值方法稱為這種插值方法稱為n次多項(xiàng)式插值（或稱代數(shù)插值），次多項(xiàng)式插值（或稱代數(shù)插值）， 利用拉格朗日插值插值方法可得利用拉格朗日插值插值方法可得)()()()()()()()()(,2, 21, 10, 0nnnnnnnxfxlxfxlxfxl

21、xfxlxP0( )0,1,njknjkjj kxxlxknxx其中上述多項(xiàng)式稱為上述多項(xiàng)式稱為 n次拉格朗日（次拉格朗日（Lagrange）插值多項(xiàng)式，）插值多項(xiàng)式， 函數(shù)函數(shù)nkxlnk, 1, 0),(,稱為拉格朗日插值基函數(shù)。稱為拉格朗日插值基函數(shù)。當(dāng)當(dāng)n=1,2時(shí)，時(shí)，n次拉格朗日（次拉格朗日（Lagrange）插值多項(xiàng)式即）插值多項(xiàng)式即為線性插值多項(xiàng)式和拋物插值多項(xiàng)式。為線性插值多項(xiàng)式和拋物插值多項(xiàng)式。4321x例例6.1 已知函數(shù)發(fā)已知函數(shù)發(fā)f(x)的函數(shù)表如下：的函數(shù)表如下： 31 . 28 . 15 . 2)(xf求其拉格朗日插值多項(xiàng)式，并求求其拉格朗日插值多項(xiàng)式，并求)5

22、. 2(f的近似值。的近似值。 解解 由于給出了由于給出了4個(gè)插值結(jié)點(diǎn)，所以可做出次數(shù)不超個(gè)插值結(jié)點(diǎn)，所以可做出次數(shù)不超過過3的拉格朗日插值多項(xiàng)式。的拉格朗日插值多項(xiàng)式。43132361)41 ( )31 ( )21 ()4( )3( )2()(233 , 0 xxxxxxxl6219421)42( )32( )22()4( )3( )2()(233 , 1xxxxxxxl472721)43( )33( )23()4( )3( )2()(233 , 2xxxxxxxl161161)44( )34( )24()4( )3( )2()(233 , 3xxxxxxxl 將上列將上列4式代入式代入n=

23、3 的拉格朗日插值公式，可得所的拉格朗日插值公式，可得所要求的插值多項(xiàng)式為要求的插值多項(xiàng)式為6 . 4933. 29 . 00667. 0)(233xxxxP將將x=2.5代入可得代入可得f(2.5)的近似值為的近似值為1.8496。 拉格朗日插值法適合節(jié)點(diǎn)較少的情況，當(dāng)節(jié)點(diǎn)較多拉格朗日插值法適合節(jié)點(diǎn)較少的情況，當(dāng)節(jié)點(diǎn)較多的大范圍高次插值的逼近效果往往并不理想且當(dāng)插值的大范圍高次插值的逼近效果往往并不理想且當(dāng)插值結(jié)點(diǎn)增加時(shí)，計(jì)算越來越繁。為了提高精度和減少計(jì)結(jié)點(diǎn)增加時(shí)，計(jì)算越來越繁。為了提高精度和減少計(jì)算，還有牛頓插值法下、三次樣條插值等，具體可參算，還有牛頓插值法下、三次樣條插值等，具體可參

24、閱有關(guān)書籍。閱有關(guān)書籍。1一維插值命令一維插值命令interp1的具體使用格式的具體使用格式y(tǒng)y=interp1(x,y,xx,method) 其中其中x,y是插值結(jié)點(diǎn)的數(shù)據(jù)向量，如果是插值結(jié)點(diǎn)的數(shù)據(jù)向量，如果y是矩陣，是矩陣，則對(duì)矩陣則對(duì)矩陣y的每一列相對(duì)的每一列相對(duì)x進(jìn)行插值，進(jìn)行插值，xx是待求函數(shù)值是待求函數(shù)值的插值結(jié)點(diǎn)向量，可以缺省。的插值結(jié)點(diǎn)向量，可以缺省。method是可選項(xiàng)，是可選項(xiàng)，說明插值使用的方法。對(duì)于一維插值，說明插值使用的方法。對(duì)于一維插值，MATLAB提供提供可選的方法有：可選的方法有：nearest，linear，spline，cubic，它們分別表示最近插值、線

25、性插值、三次樣條插值和它們分別表示最近插值、線性插值、三次樣條插值和三次插值。三次插值。2.4 MATLAB軟件實(shí)現(xiàn)插值法軟件實(shí)現(xiàn)插值法 MATLAB軟件提供了專門做各種插值的命令：軟件提供了專門做各種插值的命令：interp1（一維插值），（一維插值），interp2（二維插值），（二維插值），interp3（三維插值），（三維插值），interpn（n維插值），維插值），spline（樣條插值）等。（樣條插值）等。 2二維插值命令二維插值命令interp2的具體使用格式的具體使用格式zz=interp2(x,y,z,xx,yy,method) 該指令的意思是根據(jù)數(shù)據(jù)向量該指令的意思是根據(jù)數(shù)

26、據(jù)向量x，y，z按按method指定的方法來做插值，然后將指定的方法來做插值，然后將xx，yy處插值處插值函數(shù)的插值結(jié)點(diǎn)向量，如果函數(shù)的插值結(jié)點(diǎn)向量，如果xx，yy在插值范圍之內(nèi)，在插值范圍之內(nèi)，則返回值在則返回值在zz中，否則返回值為空中，否則返回值為空NaN。method是插值方法可選項(xiàng)，具體要求同一是插值方法可選項(xiàng)，具體要求同一維插值的情況。維插值的情況。 該命令還有以下幾種省略格式：該命令還有以下幾種省略格式：zz=interp2(z,xx,yy)zz=interp2(x,y,z,xx,yy)zz=interp2(z,ntimes)3三維插值命令三維插值命令interp3的具體使用格式

27、的具體使用格式vi=interp3(x,y,z,v,xi,yi,zi,method) 它的具體含義跟前面的一、二維插值是相似的，在它的具體含義跟前面的一、二維插值是相似的，在此不作解釋，讀者可在此不作解釋，讀者可在MATLAB工作空間中用工作空間中用help interp3命令獲得。命令獲得。4樣條插值命令樣條插值命令spline的具體使用格式的具體使用格式y(tǒng)y=spline(x,y,xx)它的意思等同于命令它的意思等同于命令yy=interp1(x,y,xx,spline)例例6.2 在用外接電源給電容器充電時(shí)，電容器兩端在用外接電源給電容器充電時(shí)，電容器兩端的電壓的電壓V將會(huì)隨著充電時(shí)間將

28、會(huì)隨著充電時(shí)間t發(fā)生變化，已知在某一次發(fā)生變化，已知在某一次實(shí)驗(yàn)時(shí)，通過測(cè)量得到下列觀測(cè)值，分別用拉格朗日實(shí)驗(yàn)時(shí)，通過測(cè)量得到下列觀測(cè)值，分別用拉格朗日插值法、分段線性插值法、三次樣條插值法畫出插值法、分段線性插值法、三次樣條插值法畫出V隨隨著時(shí)間著時(shí)間t變化的曲線圖，分別計(jì)算當(dāng)時(shí)間變化的曲線圖，分別計(jì)算當(dāng)時(shí)間t=7s時(shí)，三時(shí)，三種插值法各自算得電容器兩端電壓的近似值。種插值法各自算得電容器兩端電壓的近似值。1295 . 64321st 4 .101 .106 . 90 . 92 . 83 . 72 . 6vV解解:由于由于MATLAB沒有提供現(xiàn)成的拉格朗日插值命令，沒有提供現(xiàn)成的拉格朗日插值

29、命令，我們可以編寫一個(gè)函數(shù)我們可以編寫一個(gè)函數(shù)lglrcz.m來完成，其他兩種插值來完成，其他兩種插值法可用現(xiàn)成的命令。法可用現(xiàn)成的命令。 用用MATLAB軟件進(jìn)行三種插值計(jì)算的程序?yàn)檐浖M(jìn)行三種插值計(jì)算的程序?yàn)門hree_interpmethods.m。程序程序lagrange.m：function y=lagrange(x0,y0,x) n=length(x0)； m=length(x) ； for i=1:m z=x(i) ； s=0.0； for k=1:n p=1.0； for j=1:n if j=k p=p*(z-x0(j)/(x0(k)-x0(j) ； end end s=p*

30、y0(k)+s； end y(i)=s； end程序?yàn)槌绦驗(yàn)門hree_interpmethods.mt=1,2,3,4,6.5,9,12 ；v=6.2,7.3,8.2,9.0,9.6,10.1,10.4 ；t0=0.2:0.1:12.5；lglr=lagrange(t,v,t0) ；laglr=lagrange(t,v,7) ；fdxx=interp1(t,v,t0) ；fendxx=interp1(t,v,7) ；scyt=interp1(t,v,t0,spline) ；sancyt=interp1(t,v,7,spline)plot(t,v,*,t0,lglr,r,t0,fdxx,g,t

31、0,scyt,b)gtext(lglr)gtext(fdxx)gtext(scyt)執(zhí)行結(jié)果是執(zhí)行結(jié)果是laglr= 9.52988980716254fendxx= 9.70000000000000sancyt= 9.67118039327444圖形如圖圖形如圖6.1所示。所示。圖中曲線圖中曲線lglr表示拉格朗日插值曲線，表示拉格朗日插值曲線，scyt表示三次樣表示三次樣條插值曲線，條插值曲線，fdxx表示分段線性插值曲線。表示分段線性插值曲線。 2.5 問題求解問題求解 為了表示方便，我們將為了表示方便，我們將2.1節(jié)問題中所給表節(jié)問題中所給表6-1中的中的數(shù)據(jù)全部化為國際標(biāo)準(zhǔn)單位（表數(shù)據(jù)

32、全部化為國際標(biāo)準(zhǔn)單位（表6-2），時(shí)間用小時(shí)），時(shí)間用小時(shí)（h），高度用米（），高度用米（m）：）：表表6-2 一天內(nèi)水塔水位記錄一天內(nèi)水塔水位記錄(單位：米單位：米)時(shí)間時(shí)間/h水位水位/m時(shí)間時(shí)間/h水位水位/m09.6812.9510.210.929.4513.889.941.849.3114.989.652.959.1315.909.413.878.9816.839.184.988.8117.938.925.908.6919.048.667.008.5219.968.437.938.3920.848.228.978.2222.02水泵開動(dòng)水泵開動(dòng)9.98水泵開動(dòng)水泵開動(dòng)22.96水泵開

33、動(dòng)水泵開動(dòng)10.93水泵開動(dòng)水泵開動(dòng)23.8810.5910.9510.8224.9910.3512.0310.5025.9110.181模型假設(shè)模型假設(shè)（1）流量只取決于水位差，與水位本身無關(guān)，故由物）流量只取決于水位差，與水位本身無關(guān)，故由物理學(xué)中理學(xué)中Torriceli定律：小孔流出的液體的流速正比于定律：小孔流出的液體的流速正比于水面高度的平方根。問題中給出水塔的最低和最高水水面高度的平方根。問題中給出水塔的最低和最高水位分別是位分別是8.1648m)3024. 027(和和10.7352m)3024. 050.35(（設(shè)出口的水位為零），因?yàn)椋ㄔO(shè)出口的水位為零），因?yàn)閟qrt1467

34、. 1)1648. 87352.10(，約為，約為1，所以可忽略水位對(duì)流速的影響。，所以可忽略水位對(duì)流速的影響。（2）將流量看作是時(shí)間的連續(xù)函數(shù)，為計(jì)算簡(jiǎn)單，）將流量看作是時(shí)間的連續(xù)函數(shù)，為計(jì)算簡(jiǎn)單，不妨將不妨將流量定義成單位時(shí)間流出水的高度流量定義成單位時(shí)間流出水的高度，即水位對(duì)，即水位對(duì)時(shí)間變化率的絕對(duì)值（水位是下降的），水塔截面積時(shí)間變化率的絕對(duì)值（水位是下降的），水塔截面積為為3475.2334*)3024. 0*57(2S（m2），得到結(jié)果），得到結(jié)果后乘以后乘以s即可得真實(shí)的流量。即可得真實(shí)的流量。 2流量估計(jì)方法流量估計(jì)方法 首先依照表首先依照表6-2所給數(shù)據(jù)，用所給數(shù)據(jù)，用MA

35、TLAB作出時(shí)間作出時(shí)間水水位散點(diǎn)圖（圖位散點(diǎn)圖（圖6.2）。）。下面來計(jì)算水箱下面來計(jì)算水箱流量流量與與時(shí)間時(shí)間的關(guān)系。的關(guān)系。 根據(jù)圖根據(jù)圖6.2，一種簡(jiǎn)單的處理方法為，將表，一種簡(jiǎn)單的處理方法為，將表6-2中的中的數(shù)據(jù)分為三段，然后對(duì)每一段的數(shù)據(jù)做如下處理：數(shù)據(jù)分為三段，然后對(duì)每一段的數(shù)據(jù)做如下處理： 設(shè)某段數(shù)據(jù)設(shè)某段數(shù)據(jù)),(,),(),(1100nnyxyxyx，相鄰數(shù)據(jù)中點(diǎn)的平均流速用下面的公式（流速，相鄰數(shù)據(jù)中點(diǎn)的平均流速用下面的公式（流速=（左（左端點(diǎn)的水位右端點(diǎn)的水位）端點(diǎn)的水位右端點(diǎn)的水位）/區(qū)間長(zhǎng)度）：區(qū)間長(zhǎng)度）：iiiiiixxyyxxv111)2(每段數(shù)據(jù)首尾點(diǎn)的流速

36、用下面的公式計(jì)算：每段數(shù)據(jù)首尾點(diǎn)的流速用下面的公式計(jì)算：)()43()(022100 xxyyyxv)()43()(221nnnnnnxxyyyxv用以上公式求得時(shí)間與流速之間的數(shù)據(jù)如表用以上公式求得時(shí)間與流速之間的數(shù)據(jù)如表6-3。時(shí)間/h流速/cmh-1時(shí)間/h流速/ cmh-1029.8912.4931.520.4621.7413.4229.031.3818.4814.4326.362.39516.2215.4426.093.4116.3016.3724.734.42515.3217.3823.645.4413.0418.4923.426.4515.4519.5025.007.46513.

37、9820.4023.868.4516.3520.8422.178.9719.2922.02水泵開動(dòng)9.98水泵開動(dòng)22.96水泵開動(dòng)10.93水泵開動(dòng)23.8827.0910.9533.5024.4321.6211.4929.6325.4518.4825.9113.30由表由表6-3作出時(shí)間作出時(shí)間流速散點(diǎn)圖如圖流速散點(diǎn)圖如圖6.3。1）插值法）插值法 由表由表6-3，對(duì)水泵不工作時(shí)段，對(duì)水泵不工作時(shí)段1，2采取插值方法，可以得采取插值方法，可以得到任意時(shí)刻的流速，從而可以知道任意時(shí)刻的流量。到任意時(shí)刻的流速，從而可以知道任意時(shí)刻的流量。 我們分別采取拉格朗日插值法，分段線性插值法及我們分別采

38、取拉格朗日插值法，分段線性插值法及三次樣條插值法；對(duì)于水泵工作時(shí)段三次樣條插值法；對(duì)于水泵工作時(shí)段1應(yīng)用前后時(shí)期應(yīng)用前后時(shí)期的流速進(jìn)行插值，由于最后一段水泵不工作時(shí)段數(shù)據(jù)的流速進(jìn)行插值，由于最后一段水泵不工作時(shí)段數(shù)據(jù)太少，我們將它與水泵工作時(shí)段太少，我們將它與水泵工作時(shí)段2合并一同進(jìn)行插值合并一同進(jìn)行插值處理（該段簡(jiǎn)稱混合時(shí)段）。處理（該段簡(jiǎn)稱混合時(shí)段）。 我們總共需要對(duì)四段數(shù)據(jù)（第我們總共需要對(duì)四段數(shù)據(jù)（第1，2未供水時(shí)段，第未供水時(shí)段，第1供水時(shí)段，混合時(shí)段）進(jìn)行插值處理，下面以第供水時(shí)段，混合時(shí)段）進(jìn)行插值處理，下面以第1未供未供水時(shí)段數(shù)據(jù)為例分別用三種方法算出水時(shí)段數(shù)據(jù)為例分別用三種方

39、法算出流量函數(shù)流量函數(shù)和和用水用水量量（用水高度）。（用水高度）。 下面是用下面是用MATLAB實(shí)現(xiàn)該過程的程序?qū)崿F(xiàn)該過程的程序(第第1未供水時(shí)段未供水時(shí)段)。 t=0,0.46,1.38,2.395,3.41,4.425,12.44,6.45,7.465,8.45,8.97；v=29.89,21.74,18.48,16.22,16.30,15.32,13.04,15.45,13.98,16.35,19.29； t0=0:0.1:8.97；lglr=lagrange(t,v,t0); /*注：注：lagrange為一函數(shù)，程序同為一函數(shù)，程序同 lagrange.m*/lglrjf=trapz

40、(t0,lglr)fdxx=interpl(t,v,t0)；fdxxjf=trapz(t0,fdxx)scyt=interpl(t,v,t0,spline)；sancytjf=trapz(t0,scyt)plot(t,v,o,t0,lglr,r,t0,fdxx,g,t0,scyt,b)gtext(lglr)gtext(fdxx)gtext(spline)運(yùn)行結(jié)果lglrjf=145.6231 ； fdxxjf=147.1430； sancytjf=145.6870 圖圖6.4是對(duì)第是對(duì)第1未供水段數(shù)據(jù)用三種不同方法得到的未供水段數(shù)據(jù)用三種不同方法得到的插值函數(shù)圖，插值函數(shù)圖， 圖中曲線圖中曲線

41、lglr、fdxx和和spline分別表示用拉格朗日插分別表示用拉格朗日插值法，分段線性插值法及三次樣條插值法得到的曲線。值法，分段線性插值法及三次樣條插值法得到的曲線。 由表由表6-2知，第知，第1未供水時(shí)段的總用水高度為未供水時(shí)段的總用水高度為146(=968822)，可見上述三種插值方法計(jì)算的結(jié)果，可見上述三種插值方法計(jì)算的結(jié)果與實(shí)際值（與實(shí)際值（146）相比都比較接近?？紤]到三次樣條插）相比都比較接近?？紤]到三次樣條插值方法具有更加良好的性質(zhì)，建議采取該方法。值方法具有更加良好的性質(zhì)，建議采取該方法。 其他三段的處理方法與第其他三段的處理方法與第1未供水時(shí)段的處理方法未供水時(shí)段的處理方

42、法類似，這里不再詳細(xì)敘述，只給出數(shù)值結(jié)果和函數(shù)圖類似，這里不再詳細(xì)敘述，只給出數(shù)值結(jié)果和函數(shù)圖像（圖像（圖6.5圖圖6.7），圖中曲線標(biāo)記同圖），圖中曲線標(biāo)記同圖6.4。 圖圖6.5 第一供水段時(shí)間第一供水段時(shí)間流速示意圖流速示意圖圖圖12.6 第第2未供水段時(shí)間未供水段時(shí)間流速示意圖流速示意圖圖圖6.7 混合時(shí)段時(shí)間混合時(shí)段時(shí)間流速示意圖流速示意圖圖圖6.8是用分段線性及三次樣條插值方法得到的整個(gè)過是用分段線性及三次樣條插值方法得到的整個(gè)過程的時(shí)間程的時(shí)間流速函數(shù)示意圖。流速函數(shù)示意圖。表表6-4 各時(shí)段及一天的總用水量（用水高度）各時(shí)段及一天的總用水量（用水高度）第第1未供水段未供水段第第

43、2未供水段未供水段第第3供水段供水段混合時(shí)段混合時(shí)段全天全天拉格朗日拉格朗日插值法插值法145.6231258.866454.068992.1337550.6921分段線性分段線性插值法插值法147.1430258.969749.605176.4688532.1866三次樣條三次樣條插值法插值法145.6870258.654753.333481.7699539.4450 表表6-5是對(duì)一天中任取的是對(duì)一天中任取的4個(gè)時(shí)刻分別用個(gè)時(shí)刻分別用3種方法得到種方法得到的水塔水流量近似值。的水塔水流量近似值。時(shí)間時(shí)間6.8810.8815.8822.8815.9826671234851433.74260

44、09085346325.5662241818047734.7099621055169414.8272413793103432.9976262626262725.4465591397849525.471578947368415.0527858196582033.7089553595325925.5490892055742329.41733175863551注：注：拉格朗日插值法拉格朗日插值法分段線性插值法分段線性插值法三次樣條插值法三次樣條插值法2）擬合法）擬合法（1）擬合水位）擬合水位時(shí)間函數(shù)時(shí)間函數(shù) 從表從表6-2中的測(cè)量記錄看，一天有兩次供水時(shí)段和三中的測(cè)量記錄看，一天有兩次供水時(shí)段和三次

45、未供水時(shí)段，分別對(duì)第次未供水時(shí)段，分別對(duì)第1，2未供水時(shí)段的測(cè)量數(shù)據(jù)未供水時(shí)段的測(cè)量數(shù)據(jù)直接作多項(xiàng)式擬合，可得到水位函數(shù)（注意，根據(jù)多直接作多項(xiàng)式擬合，可得到水位函數(shù)（注意，根據(jù)多項(xiàng)式擬合的特點(diǎn)，此處擬合多項(xiàng)式的次數(shù)不宜過高，項(xiàng)式擬合的特點(diǎn)，此處擬合多項(xiàng)式的次數(shù)不宜過高，一般以一般以36次為宜）。對(duì)第次為宜）。對(duì)第3未供水時(shí)段來說，數(shù)據(jù)未供水時(shí)段來說，數(shù)據(jù)過少不能得到很好的擬合。過少不能得到很好的擬合。 設(shè)設(shè)t，h分別為已輸入的時(shí)刻和水位測(cè)量記錄（由表分別為已輸入的時(shí)刻和水位測(cè)量記錄（由表6.2提供，水泵啟動(dòng)的提供，水泵啟動(dòng)的4個(gè)時(shí)刻不輸入），這樣第個(gè)時(shí)刻不輸入），這樣第1未未供水時(shí)段各時(shí)刻的

46、水位可由如下供水時(shí)段各時(shí)刻的水位可由如下MATLAB程序完成：程序完成：t=0,0.92,1.84,2.95,3.87,4.98,5.90,7.00,7.93,8.97,10.95,12.03,12.95,13.88,14.98,15.90,16.83,17.93,19.04,19.96,20.84,23.88, 24.99,25.91h=9.68,9.45,9.31,9.13,8.98,8.81,8.69,8.52,8.39,8.22,10.82,10.50,10.21,9.94,9.65,9.41,9.18,8.92,8.66,8.43,8.22,10.59,10.35,10.18；c1=

47、polyfit(t(1:10),h(1:10),3)；% 用用3次多項(xiàng)式擬合次多項(xiàng)式擬合tp1=0:0.1:8.9；x1=polyval(c1,tp1)；plot(tp1,x1) 圖圖6.9給出的是第給出的是第1未供水時(shí)段的時(shí)間未供水時(shí)段的時(shí)間水位擬合水位擬合函數(shù)圖形。函數(shù)圖形。變量變量x1中存放了以中存放了以0.1為步長(zhǎng)算出的各個(gè)時(shí)刻的水位高為步長(zhǎng)算出的各個(gè)時(shí)刻的水位高度。同樣地，第度。同樣地，第2未供水時(shí)段時(shí)間未供水時(shí)段時(shí)間水位圖可由如下水位圖可由如下MATLAB程序完成，讀者可自己上機(jī)運(yùn)行查看。程序完成，讀者可自己上機(jī)運(yùn)行查看。c2=polyfit(t(11:21),h(11:21),3

48、)；tp2=10.9:0.1:20.9；x2=-polyval(c2,tp2)；plot(tp2,x2)（2）確定流量）確定流量時(shí)間函數(shù)時(shí)間函數(shù)對(duì)于第對(duì)于第1，2未供水時(shí)段的流量可直接對(duì)水位函數(shù)求導(dǎo)未供水時(shí)段的流量可直接對(duì)水位函數(shù)求導(dǎo)，程序如下：，程序如下：c1=polyfit(t(1:10),h(1:10),3)；c2=polyfit(t(11:21),h(11:21),3)；a1=polyder(c1)；a2=polyder(c2)； % 對(duì)擬合得到的水位函數(shù)求導(dǎo)，得流量對(duì)擬合得到的水位函數(shù)求導(dǎo)，得流量tp1=0:0.01:8.97；tp2=10.95:0.01:20.84；x13=-po

49、lyval(a1,tp1)； % Water_total1=100*trapz(tp1,x13); % 計(jì)算第計(jì)算第1未供水時(shí)段的總用水量未供水時(shí)段的總用水量x14=-polyval(a1,7.93,8.97)； % 求某些點(diǎn)的流量求某些點(diǎn)的流量 x23=-polyval(a2,tp2)；Water_total2=100*trapz(tp2,x114); %計(jì)算第計(jì)算第2未供水時(shí)段的總用水量未供水時(shí)段的總用水量x24=-polyval(a2,10.95,12.03)； x25=-polyval(a2,19.96,20.84)；subplot(1,2,1)plot(tp1,x13*100) %與

50、圖與圖6.10單位保持一致單位保持一致 % 倍數(shù)倍數(shù)100是使單位由米變成厘米是使單位由米變成厘米subplot(1,2,2)plot(tp2,x23*100) %與圖與圖6.10單位保持一致單位保持一致 程序運(yùn)行得到第程序運(yùn)行得到第1，2未供水時(shí)段的時(shí)間未供水時(shí)段的時(shí)間流量圖如圖流量圖如圖6.10，可以看到與圖，可以看到與圖6.8中用插值給出的曲線比較吻合。中用插值給出的曲線比較吻合。 如果用如果用5次多項(xiàng)式擬合則得圖次多項(xiàng)式擬合則得圖6.11的圖形，顯然較的圖形，顯然較三次擬合的效果好。三次擬合的效果好。 而第而第1供水時(shí)段的流量則用前后時(shí)期的流量進(jìn)行擬供水時(shí)段的流量則用前后時(shí)期的流量進(jìn)行

51、擬合得到。為使流量函數(shù)在合得到。為使流量函數(shù)在t=9和和t=11連續(xù)，我們只取連續(xù)，我們只取4個(gè)點(diǎn)，用三次多項(xiàng)式擬合得到第個(gè)點(diǎn)，用三次多項(xiàng)式擬合得到第1供水時(shí)段的時(shí)間供水時(shí)段的時(shí)間流量圖形如圖流量圖形如圖6.12，可以看到與圖，可以看到與圖6.8中的相應(yīng)部分中的相應(yīng)部分比較吻合。比較吻合。圖圖6.12圖圖6.8程序如下：程序如下： 在第在第2供水時(shí)段之前取供水時(shí)段之前取t=19.96,20.84兩點(diǎn)的流量，兩點(diǎn)的流量，用第用第3未供水時(shí)段的未供水時(shí)段的3個(gè)記錄做差分得到兩個(gè)流量數(shù)個(gè)記錄做差分得到兩個(gè)流量數(shù)據(jù)據(jù)21.62，18.48，然后用這，然后用這4個(gè)數(shù)據(jù)做三次多項(xiàng)式擬個(gè)數(shù)據(jù)做三次多項(xiàng)式擬合

52、得到第合得到第2供水時(shí)段與第供水時(shí)段與第3未供水時(shí)段的時(shí)間未供水時(shí)段的時(shí)間流量流量圖如圖圖如圖6.13，可以看到與圖，可以看到與圖6.8中的相應(yīng)部分也比較中的相應(yīng)部分也比較吻合吻合。圖圖6.13，圖圖6.8程序如下：程序如下：t3=19.96,20.84,t(22),t(23)；x3=x25*100,21.62,18.48；x_flow3=polyfit(t3,x3,3)；tp3=19.96:0.01:25.91；xx3=polyval(x_flow3,tp3)；tp4=20.84:0.01:23.88；gs2=polyval(x_flow3,tp4)；gsysl2=trapz(tp4,gs2); %該語句計(jì)算第該語句計(jì)算第2供水時(shí)段的總用水量供水時(shí)段的總用水量figure;plot(tp3,xx3)（3）一天總用水量的估計(jì)）一天總用水量的估計(jì) 分別對(duì)供水的兩個(gè)時(shí)段和不供水的兩個(gè)時(shí)段積分分別對(duì)供水的兩個(gè)時(shí)段和不供水的兩個(gè)時(shí)段積分（流量對(duì)時(shí)間）并求和得到一天的總用水量約為（流量對(duì)時(shí)間）并求和得到一天的總用水量約為526.8935（此數(shù)據(jù)是總用水高度，單位為（此數(shù)據(jù)是總用水高度，單位為cm）。表）。表6-6列出各段用水量，與插值法算得的表列出各段用水量，與插值法算得的表6-4相比，二相比，二者較為吻合。者較為吻合。表表6-6時(shí)段時(shí)段第第