Exp姜啟源數(shù)學(xué)建模珍藏PPT學(xué)習(xí)教案

1、會(huì)計(jì)學(xué)1Exp姜啟源數(shù)學(xué)建模珍藏姜啟源數(shù)學(xué)建模珍藏計(jì)算機(jī)會(huì)“算”嗎？靠得住嗎？例：把例：把4開開n次方，再平方次方，再平方n次，結(jié)果是次，結(jié)果是4？存在誤差？存在誤差？英國著名數(shù)值分析學(xué)家英國著名數(shù)值分析學(xué)家 Higham (1998): Higham (1998): Can you count on computers?Can you count on computers?精確計(jì)算：精確計(jì)算：解析結(jié)果解析結(jié)果 (Analytical)近似近似計(jì)計(jì)算：算：數(shù)值數(shù)值結(jié)果結(jié)果(Numerical)?422n=55左右：結(jié)果變成左右：結(jié)果變成1計(jì)計(jì)算功效算功效=計(jì)算工具計(jì)算工具*計(jì)算方法計(jì)算方法(算

2、法算法)浮點(diǎn)運(yùn)算：舍入誤差浮點(diǎn)運(yùn)算：舍入誤差第1頁/共42頁實(shí)驗(yàn)3的基本內(nèi)容3.3.數(shù)值積分的數(shù)值積分的梯形公式、辛普森公式和梯形公式、辛普森公式和高斯公式。高斯公式。1.1.插值的基本原理；插值的基本原理； 三種插值方法：拉格朗日插三種插值方法：拉格朗日插 值，分段線值，分段線性性 插值，三次樣條插值。插值，三次樣條插值。2.2.插值的插值的 MATLAB 實(shí)現(xiàn)實(shí)現(xiàn)及插值的應(yīng)用及插值的應(yīng)用。4.4.數(shù)值積分的數(shù)值積分的 MATLAB 實(shí)現(xiàn)實(shí)現(xiàn)及數(shù)值積分的應(yīng)用及數(shù)值積分的應(yīng)用。第2頁/共42頁什么是插值什么是插值(Interpolation)？從查函數(shù)表說？從查函數(shù)表說起起查查 函函 數(shù)數(shù) 表

3、表xtdtex2221)(x0121.0 0.8413 0.8438 0.84611.1 0.8643 0.8665 0.86861.2 0.8849 0.8869 0.8888標(biāo)準(zhǔn)正態(tài)分布函數(shù)表標(biāo)準(zhǔn)正態(tài)分布函數(shù)表求求 (1.114) (1.114)=0.8665 (0.8686 0.8665) 0.4=0.8673插值插值插值在圖像處理插值在圖像處理/數(shù)控加工數(shù)控加工/外觀設(shè)計(jì)等領(lǐng)域有重要應(yīng)用外觀設(shè)計(jì)等領(lǐng)域有重要應(yīng)用第3頁/共42頁插值的基本原理插值的基本原理插值問題的提法插值問題的提法已知已知 n+1n+1個(gè)節(jié)點(diǎn)個(gè)節(jié)點(diǎn), 1 , 0(),(njyxjj其中其中jx互不相同，不妨設(shè)互不相同，

4、不妨設(shè)),10bxxxan求任一插值點(diǎn)求任一插值點(diǎn))(*jxx 處的插值處的插值.*y0 x1xnx0y1y節(jié)點(diǎn)可視為由節(jié)點(diǎn)可視為由)(xgy 產(chǎn)生產(chǎn)生,g表達(dá)式復(fù)雜表達(dá)式復(fù)雜,甚至無表達(dá)甚至無表達(dá)式式*x*y第4頁/共42頁0 x1xnx0y1y求解插值問題的基本思路求解插值問題的基本思路構(gòu)造一個(gè)構(gòu)造一個(gè)( (相對(duì)簡(jiǎn)單的相對(duì)簡(jiǎn)單的) )函數(shù)函數(shù)),(xfy 通過全部節(jié)點(diǎn)通過全部節(jié)點(diǎn), ,即即), 1 ,0()(njyxfjj再用再用)(xf計(jì)算插值，計(jì)算插值，即即).(*xfy *x*y插值的插值的基本原理基本原理第5頁/共42頁1.1.拉格朗日拉格朗日(Lagrange)(Lagrange

5、)多項(xiàng)式插值多項(xiàng)式插值1.0 1.0 插值多項(xiàng)式插值多項(xiàng)式) 1 ()(0111axaxaxaxLnnnnnnnnnnnnnyyYaaAxxxxX001100,11在什么條件下）(0)det(X), 1 , 0()(njyxLjjn)2(YXA 求ia三種插值三種插值方法方法有唯一解)2(第6頁/共42頁1.1 1.1 拉格朗日插值多項(xiàng)拉格朗日插值多項(xiàng)式式nixxxxxxxxxxxxxxxxxlniiiiiiniii1 , 0,)()()()()()()(110110)3()()(0 xlyxLiniinjjnjiyxLjijixl)(,0, 1)(則若又又（2）有唯一解，故有唯一解，故（3）

6、與與（1）相同相同。 基函數(shù)基函數(shù)( )ilx) 1 ()(0111axaxaxaxLnnnnn) 2(YXA三種插值三種插值方法方法第7頁/共42頁),(),()!1()()()()(0)1(baxxngxLxgxRnjjnnn1)1()(nnMg減小（粗略地看）如何使誤差)(xRn平緩gjxx 接近njjnnxxnMxR01)!1()(三種插值三種插值方法方法1.2 1.2 誤差估計(jì)誤差估計(jì)增加n第8頁/共42頁1.3 1.3 拉格朗日插值多項(xiàng)式的振蕩拉格朗日插值多項(xiàng)式的振蕩?)(?)(xRxLnnn55,11)(2xxxg63. 363. 3),()(limxxgxLnnRunge現(xiàn)現(xiàn)象

7、象取n=2,4,6,8,10,計(jì)算Ln(x), 畫出圖形-505-1.5-1-0.500.511.52y=1/(1+x2)n=2n=4n=6n=8n=10三種插值三種插值方法方法Runge.m第9頁/共42頁2.2.分段線性插值分段線性插值xjxj-1xj+1x0 xn其它,0,)()()(1111110jjjjjjjjjjjnjjjnxxxxxxxxxxxxxxxlxlyxI計(jì)算量與計(jì)算量與n n無關(guān)無關(guān); ;n n越大，誤差越小越大，誤差越小. .nnnxxxxgxI0),()(lim三種插值方法三種插值方法第10頁/共42頁機(jī)翼下輪廓線3. 3. 三次樣條插值三次樣條插值樣條函數(shù)的由來樣

8、條函數(shù)的由來飛機(jī)、船體、汽車外形等的放樣（設(shè)計(jì)飛機(jī)、船體、汽車外形等的放樣（設(shè)計(jì)）細(xì)木條：樣條細(xì)木條：樣條第11頁/共42頁3. 3. 三次樣條插值三次樣條插值, 1,),()(1nixxxxsxSiii,)()3), 1 ,0()()2), 1()()10223niiiiiiixxCxSniyxSnidxcxbxaxs數(shù)學(xué)樣條（數(shù)學(xué)樣條（spline）iiiidcban,4 個(gè)待定系數(shù)3）) 1, 1()()()()(),()(111 nixsxsxsxsxsxsiiiiiiiiiiii3）2），3）共 4n-2個(gè)方程三種插值方法三種插值方法第12頁/共42頁自然邊界條件）(0)()()40

9、 nxSxS)(,)4)3)2xSdcbaiiii三次樣條插值確定三次樣條插值確定4 4n n個(gè)系數(shù)需增加個(gè)系數(shù)需增加 2 2個(gè)條件個(gè)條件思考1）自然邊界條件的幾何意義是什么？2）樣條插值為什么普遍用3次多項(xiàng)式，而不是2或4次？三次樣條插值三次樣條插值).()(limxgxSn第13頁/共42頁三種插值方法小結(jié)三種插值方法小結(jié) 拉格朗日插值（高次多項(xiàng)式插值）：拉格朗日插值（高次多項(xiàng)式插值）：曲線光滑；誤差估計(jì)有表達(dá)式；收斂性不能保證。曲線光滑；誤差估計(jì)有表達(dá)式；收斂性不能保證。用于理論分析，實(shí)際意義不大用于理論分析，實(shí)際意義不大。 分段線性和三次樣條插值（低次多項(xiàng)式插值）：分段線性和三次樣條插

10、值（低次多項(xiàng)式插值）：曲線不光滑（三次樣條插值已大有改進(jìn)）；誤差估曲線不光滑（三次樣條插值已大有改進(jìn)）；誤差估計(jì)較難（對(duì)三次樣條插值）；收斂性有保證。計(jì)較難（對(duì)三次樣條插值）；收斂性有保證。簡(jiǎn)單實(shí)用，應(yīng)用廣泛簡(jiǎn)單實(shí)用，應(yīng)用廣泛。 其他：其他：Hermite插值、分段三次插值、二維插值等插值、分段三次插值、二維插值等根據(jù)需要，根據(jù)需要，各取所需各取所需。第14頁/共42頁1. 1. 拉格朗日插值拉格朗日插值: :自編程序自編程序, ,如名為如名為 lagr.m 的的M文件，文件， 第一行為第一行為 function y=lagr(x0,y0,x) 輸入輸入: :節(jié)點(diǎn)節(jié)點(diǎn)x0,y0, 插值點(diǎn)插值點(diǎn)

11、x ( (均為均為數(shù)組，長度自定義數(shù)組，長度自定義) )）；）； 輸出輸出: :插值插值y ( (與與x同長度數(shù)組同長度數(shù)組) )）。）。 應(yīng)用時(shí)輸入應(yīng)用時(shí)輸入x0,y0,x后后, ,運(yùn)行運(yùn)行 y=lagr(x0,y0,x)2. 2. 分段線性插值分段線性插值: :已有程序已有程序 y=interp1(x0,y0,x) y=interp1(x0,y0,x,linear)3. 3. 三次樣條插值三次樣條插值: :已有程序已有程序 y=interp1(x0,y0,x,spline) 或或 y=spline(x0,y0,x)用MATLAB作插值計(jì)算注：注：MATLAB有樣條工具箱（有樣條工具箱（Sp

12、line Toolbox）第15頁/共42頁用MATLAB作插值計(jì)算55,11)(2xxxg為例，作三種插值的比較為例，作三種插值的比較以以 0 1.0000 1.0000 1.0000 1.0000 0.5000 0.8000 0.8434 0.7500 0.8205 1.0000 0.5000 0.5000 0.5000 0.5000 1.5000 0.3077 0.2353 0.3500 0.2973 2.0000 0.2000 0.2000 0.2000 0.2000 2.5000 0.1379 0.2538 0.1500 0.1401 3.0000 0.1000 0.1000 0.1

13、000 0.1000 3.5000 0.0755 -0.2262 0.0794 0.0745 4.0000 0.0588 0.0588 0.0588 0.0588 4.5000 0.0471 1.5787 0.0486 0.0484 5.0000 0.0385 0.0385 0.0385 0.0385 x y y1 y2 y3 用用n=11個(gè)節(jié)點(diǎn)，個(gè)節(jié)點(diǎn)，m=21個(gè)插值點(diǎn)，三種方法作插值，畫圖。個(gè)插值點(diǎn)，三種方法作插值，畫圖。chazhi1第16頁/共42頁插值的應(yīng)用加工時(shí)需要加工時(shí)需要x每改變每改變0.05時(shí)的時(shí)的y值值chazhi2圖1 零件的輪廓線 （x間隔0.2）表1 x間隔0.2的加

14、工坐標(biāo)x,y（圖1右半部的數(shù)據(jù)）數(shù)控機(jī)床加工零件數(shù)控機(jī)床加工零件 0.0,5.00 0.2,4.710.4,4.31 0.6,3.68 0.8,3.051.0,2.50 1.2,2.051.4,1.69 1.6,1.40 1.8,1.182.0,1.00 2.2,0.862.4,0.74 2.6,0.64 模型模型 將圖1逆時(shí)針方向轉(zhuǎn)90度，輪廓線上下對(duì)稱，只需對(duì)上半部計(jì)算一個(gè)函數(shù)在插值點(diǎn)的值。 圖2 逆時(shí)針方向轉(zhuǎn)90度的結(jié)果-5-4-3-2-101234500.511.522.533.544.55uv令v=x, u= -y 第17頁/共42頁為什么要作數(shù)值積分為什么要作數(shù)值積分 許多函數(shù)許多

15、函數(shù)“積不出來積不出來”, ,只能用數(shù)值方法，如只能用數(shù)值方法，如dxxxdxebabaxsin,22 積分是重要的數(shù)學(xué)工具，是微分方程、概率積分是重要的數(shù)學(xué)工具，是微分方程、概率論等的基礎(chǔ)；在實(shí)際問題中有直接應(yīng)用。論等的基礎(chǔ)；在實(shí)際問題中有直接應(yīng)用。 對(duì)于用離散數(shù)據(jù)或者圖形表示的函數(shù)對(duì)于用離散數(shù)據(jù)或者圖形表示的函數(shù)，(可以先做插值然后積分；或者直接利用點(diǎn)做數(shù)值積分)計(jì)算積分只有求助于數(shù)值方法。計(jì)算積分只有求助于數(shù)值方法。數(shù)值積分?jǐn)?shù)值積分第18頁/共42頁nabfIIdxxfInkknnnba)(,lim)(1數(shù)數(shù) 值值 積積 分分 的的 基基 本本 思思 路路回回 憶憶 定定 積積 分分 的

16、的 定定 義義各種數(shù)值積分方法研究的是各種數(shù)值積分方法研究的是k),(ba如何取值，區(qū)間如何取值，區(qū)間如何劃分，如何劃分，使得既能保證一定精度，計(jì)算量又小。使得既能保證一定精度，計(jì)算量又小。n n充分大時(shí)充分大時(shí)I In n就是就是I I的數(shù)值積分的數(shù)值積分（計(jì)算功效：算得準(zhǔn)，算得快）（計(jì)算功效：算得準(zhǔn)，算得快）第19頁/共42頁1.1.從矩形公式到梯形公式數(shù)值積分?jǐn)?shù)值積分yy=f(x)xbao)1(10nkknfhL)(,10kknkxffnabhbxxxxa)2(1nkknfhRnnRL,平均，得到梯形公式) 3()(2011nnkknffhfhTxk+1xkxk-1fk第20頁/共42頁

17、2.2.辛普森辛普森(Simpson)(Simpson)公式公式（拋物線公式）（拋物線公式） 梯形公式相當(dāng)于用分段線性插值函數(shù)分段線性插值函數(shù)代替)(xf每段要用相鄰兩小區(qū)間兩小區(qū)間端點(diǎn)的三個(gè)函數(shù)值端點(diǎn)的三個(gè)函數(shù)值拋物線拋物線公式公式提高精度提高精度分段二次插值函數(shù)分段二次插值函數(shù)2221212222(,),(,),(,)0,1, ,1kkkkkkxfxfxfkm數(shù)值積分?jǐn)?shù)值積分yy=f(x)xbaox2kf2kx2k+1x2k+2f2k+1f2k+2區(qū)間數(shù)必須為偶數(shù)區(qū)間數(shù)必須為偶數(shù)mn2第21頁/共42頁) 4(2),24(3112101220mabhffffhSmkkmkkmm 對(duì)對(duì)k求和

18、求和(共共m段段) )，得（復(fù)合），得（復(fù)合）辛普森公式辛普森公式：)4(3)(22122222kkkxxkfffhdxxskk二次插值函數(shù)sk(x)構(gòu)造用),(),(),(2222121222kkkkkkfxfxfx2.2.辛普森辛普森(Simpson)(Simpson)公式（拋物線公式）公式（拋物線公式）第22頁/共42頁bannnTdxxfTITfR)(),(梯形公式在每小段上是用梯形公式在每小段上是用線性插值函數(shù)線性插值函數(shù)T T( (x) )代替代替 f( (x) )11()( )( )()(),(,)2kkkkkkff xT xxxxxxx x （拉格朗日插值余項(xiàng)）梯形公式梯形公式

19、的誤差估計(jì)的誤差估計(jì))(2011nnkknffhfhTbadxxf)()(12)(2)()()(3111kxxkkkxxfhdxxxxxfdxxTxfkkkk 因?yàn)椋阂驗(yàn)椋?x-xk)(x-xk+1)在在(xk,xk+1)不變號(hào)，所以：不變號(hào)，所以：第23頁/共42頁)5()(12|),(|22abMhTfRn梯形公式梯形公式Tn的的誤差是誤差是h2階的階的),(, )(max2baxxfM 估計(jì)估計(jì)habn因?yàn)?103)(12|),(|nkknfhTfR 梯形公梯形公式的誤差式的誤差)( )( (121)(121 2afbfdxxfhTIban 103)(12)(nkkbannfhTdxxf

20、TI) 5()()(122afbfhTIn第24頁/共42頁同理可得：同理可得：) 6()(180| ),(|44abMhSfRn其中其中),(,)(max)4(4baxxfM辛普森公式辛普森公式Sn的誤差是的誤差是h4階的階的。辛普森公式的誤差估辛普森公式的誤差估計(jì)計(jì)第25頁/共42頁梯形公式和辛普森公式的收斂性若對(duì)若對(duì)I某個(gè)數(shù)值積分某個(gè)數(shù)值積分In有有chIIpnnlim（非零常數(shù)）（非零常數(shù)）則稱則稱 In是是 p 階收斂的階收斂的。梯形公式梯形公式 2 2 階收斂，辛普森公式階收斂，辛普森公式 4 4 階收斂。階收斂。c=0: 至少至少p階收斂（超階收斂（超p階收斂）階收斂）第26頁/

21、共42頁積分步長的自動(dòng)選取積分步長的自動(dòng)選取選定數(shù)值積分公式后，如何確定步長選定數(shù)值積分公式后，如何確定步長h以滿足給定的誤差以滿足給定的誤差 )()(122afbfhTIn梯形公式)(412nnTITInnTT2用二分法只要用二分法只要其中其中fk+1/2是原是原分點(diǎn)分點(diǎn)xk,xk+1的中點(diǎn)的中點(diǎn)（記記xk+1/2）的函數(shù)值的函數(shù)值1021222nkknnfhTT且且T T2n2n可在可在T Tn n基礎(chǔ)上計(jì)算基礎(chǔ)上計(jì)算)(3122nnnTTTInTI2)2(2nnhh第27頁/共42頁高斯高斯(Gauss)(Gauss)求積公式求積公式矩形公式矩形公式(1)、(2)梯形公式梯形公式(3)辛

22、普森公式辛普森公式(4)A Ak k是與是與f f無關(guān)的常數(shù)無關(guān)的常數(shù)代數(shù)代數(shù)精度精度設(shè)設(shè),)(kxxf用用(7)計(jì)算計(jì)算,)(badxxfI若對(duì)于若對(duì)于mk,1 ,0都有都有, IIn而當(dāng)而當(dāng), 1IImkn則則稱稱In的代數(shù)精度為的代數(shù)精度為m.)7()(1nkkknxfAINewton-Cotes方法方法第28頁/共42頁梯形公式的代數(shù)精度（考察梯形公式的代數(shù)精度（考察T1）k=1f(x)=x222abxdxIba2)()()(21baabbfafhT3332abdxxIba2)(221baabTk=2f(x)=x2IT1IT 1梯形公式的代數(shù)精度為梯形公式的代數(shù)精度為1辛普森公式的代數(shù)

23、精度為辛普森公式的代數(shù)精度為3第29頁/共42頁高斯公式的思路高斯公式的思路取消對(duì)節(jié)點(diǎn)的限制，按照代數(shù)精度最大取消對(duì)節(jié)點(diǎn)的限制，按照代數(shù)精度最大的原則，同時(shí)確定節(jié)點(diǎn)的原則，同時(shí)確定節(jié)點(diǎn)xk和系數(shù)和系數(shù)Ak構(gòu)造求積公式構(gòu)造求積公式)()(22112xfAxfAG對(duì)于對(duì)于11)(dxxfI使使G G2 2的代數(shù)精度為的代數(shù)精度為3 332, 1)(xxxxf)()()(221111xfAxfAdxxf確定確定2121,AAxx第30頁/共42頁03/202322311222211221121xAxAxAxAxAxAAA將將f(x)f(x)代入計(jì)算得代入計(jì)算得1,3/1,3/12121AAxx)3/

24、1 ()3/1(2ffG用用n個(gè)節(jié)點(diǎn)，個(gè)節(jié)點(diǎn)，Gn的代數(shù)精度可達(dá)的代數(shù)精度可達(dá)2n-1, 但是需解但是需解復(fù)雜的非線性方程組，實(shí)用價(jià)值不大。復(fù)雜的非線性方程組，實(shí)用價(jià)值不大。第31頁/共42頁常 用 的 高 斯 公 式將將( (a,ba,b) )分小，把小區(qū)間變換為分小，把小區(qū)間變換為(-1(-1，1), 1), 再用再用G G2 2mkkkbazfzfhdxxf1)2()1 ()()(2)(1)(2)11,222 32 3kkkkkkxxxxhhzz伸縮變換公式mkkhaxmabhk, 1, 0,/ )(代數(shù)精度為代數(shù)精度為3節(jié)點(diǎn)加密時(shí)，原計(jì)算信息無法利用節(jié)點(diǎn)加密時(shí)，原計(jì)算信息無法利用第32

25、頁/共42頁思路思路：將積分區(qū)間分小，在小區(qū)間上用：將積分區(qū)間分小，在小區(qū)間上用n不太不太 大大的的 。而在節(jié)點(diǎn)加密一倍時(shí)能夠利用原節(jié)點(diǎn)的。而在節(jié)點(diǎn)加密一倍時(shí)能夠利用原節(jié)點(diǎn)的函數(shù)值，可以把區(qū)間的端點(diǎn)作為固定節(jié)點(diǎn)。函數(shù)值，可以把區(qū)間的端點(diǎn)作為固定節(jié)點(diǎn)。改進(jìn)的高斯公式nG)()()(121bfAxfAafAGnnkkknGauss-Lobatto求積公式求積公式 其中a, b為小區(qū)間的端點(diǎn)，nnAAxx,112為2n-2個(gè)參數(shù)，代數(shù)精度可達(dá)到代數(shù)精度可達(dá)到2n-3注意：實(shí)際計(jì)算中一般采用自適應(yīng)方法確定步長注意：實(shí)際計(jì)算中一般采用自適應(yīng)方法確定步長第33頁/共42頁用用MATLAB 作數(shù)值積分作數(shù)值

26、積分10nkknfhLnkknfhR1矩形矩形公式公式Sum(x)輸入數(shù)組x(即fk),輸出x的和(數(shù))cumsum(x)輸入數(shù)組x,輸出x的依次累加和(數(shù)組)梯形梯形公式公式)(2011nnkknffhfhTTrapezoidal (梯形)trapz(x)輸入數(shù)組x,輸出按梯形公式x的積分(單位步長)trapz(x,y)輸入同長度數(shù)組 x,y,輸出按梯形公式y(tǒng)對(duì)x的積分(步長不一定相等)第34頁/共42頁用用MATLAB 作數(shù)值積分作數(shù)值積分mabhffffhSmkkmkkmn2),24(3112101220辛普森公式辛普森公式(quadrature: 求積) quad(fun,a,b,to

27、l,trace)I,fn=quad()用自適應(yīng)辛普森公式計(jì)算tol為絕對(duì)誤差，缺省時(shí)為10-6Gauss-Lobatto公公式式)()()(121bfAxfAafAGnnkkknquadl(fun,a,b,tol,trace)I,fn=quadl()用自適應(yīng)Gauss-Lobatto公式計(jì)算 tol為絕對(duì)誤差，缺省時(shí)為10-6注意：注意：fun.m中應(yīng)以自變量為矩陣的形式輸入中應(yīng)以自變量為矩陣的形式輸入(點(diǎn)運(yùn)算點(diǎn)運(yùn)算)第35頁/共42頁矩形域上計(jì)算二重積分的命令：矩形域上計(jì)算二重積分的命令：dblquad(fun,xmin,xmax,ymin,ymax,tol)廣義積分、二重和三重積分長方體上計(jì)算三重積分的命令：長方體上計(jì)算三重積分的命令：triplequad(fun,xmin,xmax,ymin,ymax, zmin,zmax,tol)注：注：fun是被積函數(shù)，本身可以有自己的參數(shù)是被積函數(shù)，本身可以有自己的參數(shù)廣義積分：廣義積分：通過分析和控制誤差，轉(zhuǎn)換成普通積分通過分析和控制誤差，轉(zhuǎn)換成普通積分quadv(fun,a,b,tol,trace)向量值積分：向量值積分：第36頁/共42頁用用MATLAB 作數(shù)值積分作數(shù)值積分例例. 計(jì)算計(jì)算4011s i nd xx1 1）矩形公式和梯形公