數(shù)值計(jì)算的資料

數(shù)值計(jì)算的資料第1頁/共51頁WenjianYu2第七章數(shù)值積分與數(shù)值微分

第2頁/共51頁3數(shù)值積分的基本概念第3頁/共51頁4數(shù)值積分目的與用途經(jīng)典問題:算幾何形體的面積、

體積，力學(xué)中

物體的重心位置例:鋁制波紋瓦的長(zhǎng)度問題由一塊平整的鋁板壓制而成.

若每個(gè)波紋的高度(自中心線)

為1英寸,周期為2英寸,做4英尺長(zhǎng)波紋瓦需多長(zhǎng)鋁板?

第二類橢圓積分,無法解析求出!

第4頁/共51頁5數(shù)值積分基本思想

...

積分系數(shù)

積分節(jié)點(diǎn)希望用較少的計(jì)算量得到較準(zhǔn)確的結(jié)果第5頁/共51頁6插值型求積公式

中矩形公式梯形公式第6頁/共51頁7積分余項(xiàng)與代數(shù)精度

反映了計(jì)算的截?cái)嗾`差插值余項(xiàng)的積分衡量求積公式準(zhǔn)確度的另一個(gè)指標(biāo)注意:對(duì)某些情況，代數(shù)精度并不是越高越好第7頁/共51頁8積分余項(xiàng)與代數(shù)精度

(至少0次代數(shù)精度)它至少有n次代數(shù)精度即插值型

求積公式第8頁/共51頁9求積公式的收斂性與穩(wěn)定性

(一系列求積公式的性質(zhì))

積分問題一般不太敏感

第9頁/共51頁10求積公式的收斂性與穩(wěn)定性

這是控制數(shù)值計(jì)算誤差能達(dá)到的最佳情況

要盡量尋求穩(wěn)定的求積公式第10頁/共51頁11牛頓-柯特斯公式第11頁/共51頁12Newton-Cotes公式

就是n階牛頓-柯特斯公式

n=1,1/2,1/2n=2,1/6,2/3,1/6

n=4,7/90,16/45,2/15,16/45,7/90Cotes系數(shù)表

一系列求積公式便于使用第12頁/共51頁13Newton-Cotes公式

n=1,1/2,1/2n=2,1/6,2/3,1/6

n=4,7/90,16/45,2/15,16/45,7/90

n=8,Cotes系數(shù)表思考題

梯形公式Simpson公式Cotes公式中矩形公式可看成是n=0時(shí)的特例.例題(板書)第13頁/共51頁14Newton-Cotes公式

關(guān)鍵看積分:

(n階公式至少有n次代數(shù)精度)一般不用n=3對(duì)應(yīng)的N-C公式第14頁/共51頁15低階N-C公式的積分余項(xiàng)

不保號(hào),無法用積分中值定理

2

詳細(xì)過程自己看書第15頁/共51頁16穩(wěn)定性、收斂性

n=8,實(shí)際只使用n<8的偶數(shù)階N-C公式,也看出代數(shù)精度不是越高越好第16頁/共51頁17復(fù)合求積公式第17頁/共51頁18復(fù)合求積公式

(compositequadrature)

積分誤差:n增大,誤差減小

仍是機(jī)械求積公式第18頁/共51頁19復(fù)合求積公式

2階準(zhǔn)確度

第19頁/共51頁20復(fù)合求積公式

與復(fù)合梯形公式對(duì)比,例7.4第20頁/共51頁21復(fù)合求積公式步長(zhǎng)折半的復(fù)合求積公式計(jì)算積分余項(xiàng)公式包含被積函數(shù)的高階導(dǎo)數(shù),很難應(yīng)用.

常常動(dòng)態(tài)地確定步長(zhǎng)h常用的動(dòng)態(tài)減小步長(zhǎng)策略是:步長(zhǎng)折半,利用已算出結(jié)果復(fù)合梯形公式的情況遞推化的復(fù)合梯形公式:

(逐漸減小,直到滿足精度要求)

只需再計(jì)算新增節(jié)點(diǎn)的函數(shù)值第21頁/共51頁22復(fù)合求積公式步長(zhǎng)折半的復(fù)合求積公式計(jì)算復(fù)合Simpson公式的情況很少使用中矩形公式的原因

與梯形公式有相同的代數(shù)精度/準(zhǔn)確度,計(jì)算量更小

可類似構(gòu)造復(fù)合中矩形公式,但在步長(zhǎng)折半時(shí),無法重用以前的結(jié)果第22頁/共51頁23Remberg積分算法第23頁/共51頁24復(fù)合梯形公式的余項(xiàng)展開式

(可與Richardson外推結(jié)合)Th7.5

所有小區(qū)間的積分求和:第24頁/共51頁25復(fù)合梯形公式的余項(xiàng)展開式

Th7.5

所有小區(qū)間乘h/2求和:

第25頁/共51頁26RichardsonExtrapolation

(“0”代表未經(jīng)外推的原始公式)

第26頁/共51頁27RichardsonExtrapolation

第27頁/共51頁28Romberg算法龍貝格算法列三角形表格,按行依次計(jì)算計(jì)算公式

可證明:

具有2k+1次代數(shù)精度

(類似高階差商的計(jì)算)第28頁/共51頁29Romberg算法

做等距的節(jié)點(diǎn)分布,充分利用其上函數(shù)值得到最高準(zhǔn)確度的結(jié)果要求被積函數(shù)充分光滑!第29頁/共51頁30Romberg算法

011/80.99739781/40.98961583/80.97672671/20.95885115/80.93615563/40.90885177/80.877192610.841471010.9207355

1/20.93979330.9461459

1/220.94451350.94608690.946083001/230.9460833第30頁/共51頁31Romberg算法

10.500000

1/20.4267770.402369

1/220.4070180.4004320.400302

1/230.4018120.4000770.4000540.400050

1/240.4004630.4000140.4000090.4000090.400009

1/250.4001180.4000020.4000020.4000020.4000020.400002

第31頁/共51頁32自適應(yīng)積分算法第32頁/共51頁33自適應(yīng)積分算法(adaptivequadrature)基本思想Romberg算法效果不好的情況積分節(jié)點(diǎn)沒必要均勻分布怎樣自動(dòng)地非均勻取點(diǎn),使計(jì)算

結(jié)果達(dá)到精度要求?1.評(píng)估當(dāng)前區(qū)間積分結(jié)果的準(zhǔn)

確度,若不準(zhǔn)確就將它一分為二,直至小區(qū)間的結(jié)果準(zhǔn)確2.

用兩個(gè)低階求積公式算同一個(gè)積分,它們之差可近似判斷結(jié)果的準(zhǔn)確度

(遞歸計(jì)算過程)第33頁/共51頁34自適應(yīng)積分算法一個(gè)自適應(yīng)求積算法對(duì)每個(gè)區(qū)間,用Simpson公式、復(fù)合Simpson公式計(jì)算無論區(qū)間大小,用相同的誤差閾值;函數(shù)的遞歸調(diào)用原理算法:

實(shí)際的算法需保證函數(shù)值不重復(fù)計(jì)算見課本的quadtx程序

演示模塊7.3,quadgui第34頁/共51頁35自適應(yīng)積分算法一個(gè)自適應(yīng)求積程序quadtx例子:更多討論通過閾值設(shè)置控制相對(duì)誤

差;不連續(xù)函數(shù)的特殊處理還有其他估計(jì)積分誤差的

方法,比如利用中矩形公式,梯形公式的差注意與Romberg算法的不同

[Q,fcnt]=quadtx(@humps,0,1,1e-3)fcnt=69,用了69個(gè)積分點(diǎn)第35頁/共51頁36高斯求積公式第36頁/共51頁37高斯求積公式

解得:

第37頁/共51頁38高斯求積公式

高斯積分有2n+1次代數(shù)精度

插值型求積公式,代數(shù)精度的概念也可擴(kuò)展

第38頁/共51頁39高斯求積公式

Th7.7比自適應(yīng)積分算法使用方便第39頁/共51頁40高斯-勒讓德公式

012345高斯-勒讓德積分表第40頁/共51頁41高斯-勒讓德公式

第41頁/共51頁42數(shù)值微分第42頁/共51頁43數(shù)值微分

(向前差分)(向后差分)(中心差分)利用Taylor展開推出:

第43頁/共51頁44數(shù)值微分

第44頁/共51頁45數(shù)值微分

h0.10.4516049081

0.050.45407616940.4548999231

00250.45469262880.45489811520.4548979947p359

準(zhǔn)確的有效數(shù)字位數(shù)第45頁/共51頁46數(shù)值微分的應(yīng)用

二階中心差分

例如n=2

第46頁/共51頁47Matlab中的積分計(jì)算第47頁/共51頁Matlab中的數(shù)值積分指定被積函數(shù)inline命令Matlabv7推薦使用匿名函數(shù)@輸入?yún)?shù)可以是向量，因此需采用逐項(xiàng)運(yùn)算符號(hào)M文件:可處理含奇異點(diǎn)的積分帶多個(gè)參數(shù):>>f1=inline('1./sqrt(1+x.^4)'

)>>f2=@(x)1./sqrt(1+x.^4)Functionf=sinc(r)ifx==0f=1;elsef=sin(x)./x;end>>f_beta=inline('t.^(z-1).*(1-t).^(w-1)','t','z','w')函數(shù)句柄第48頁/共51頁Matlab中的數(shù)值積分一維積分的命令quad(第7.5節(jié)的自適應(yīng)積分)quadl,quadgk(擴(kuò)展的Gauss自適應(yīng)積分,pp.231)>>Q=quad(f1,0,1)

%inlinefunction>>quad(@sinc,0,pi)%functiondefinedbyM-fileans=0.58949>>[beta2_5,fcnt]=quad(f_beta,0,1,1e-5,0,2,5)beta2_5=0.033333

fcnt=17[q,fcnt]=quad(fun,a,b,tol,trace,p1,p2,...)積分參數(shù)是否輸出函數(shù)計(jì)算次數(shù)等.0/非0準(zhǔn)確度控制：缺省值10-6第49頁/共51頁Matlab中的數(shù)值積分二重、三重積分dblquad,quad2d,triplequad符號(hào)積分定義符號(hào)變量:

sym(),syms,符號(hào)積分:intsimple(表達(dá)式化簡(jiǎn))double離散數(shù)據(jù)點(diǎn)積分復(fù)合梯形法trapz(x,y)>>symsx>>h=1/((x-.3)^2+.01)+1/((x-.9)^2+.0