第5章兩點(diǎn)邊值問題求解方法

1、航空航天中的計(jì)算方法航空航天中的計(jì)算方法授課教師：陳琪鋒授課教師：陳琪鋒中南大學(xué)航空航天學(xué)院中南大學(xué)航空航天學(xué)院第二部分第二部分 邊值問題求解方法邊值問題求解方法第第5章章 兩點(diǎn)邊值問題求解方法兩點(diǎn)邊值問題求解方法航空航天中的計(jì)算方法Page 42022-3-17內(nèi)容提要內(nèi)容提要5.1常微分方程邊值問題的概念常微分方程邊值問題的概念5.2打靶法打靶法5.3有限差分法有限差分法5.4有限元法有限元法1 Part 3: Two-Point Boundary Value Problems.2 David L. Darmofal, Computational Methods in Aerospace

2、Engineering (Lecture Notes), MIT, 2005. Chap11,12.3 清華大學(xué)數(shù)學(xué)系編，現(xiàn)代應(yīng)用數(shù)學(xué)手冊清華大學(xué)數(shù)學(xué)系編，現(xiàn)代應(yīng)用數(shù)學(xué)手冊計(jì)算方法分冊（計(jì)算方法分冊（第十一章，常微分方程邊值問題的數(shù)值方法），北京出版第十一章，常微分方程邊值問題的數(shù)值方法），北京出版社，社，1990.航空航天中的計(jì)算方法Page 52022-3-175.1 常微分方程邊值問題的概念常微分方程邊值問題的概念 對于常微分方程：對于常微分方程：其中其中 ，x為標(biāo)量，為標(biāo)量， y和和 f 為為m維向量維向量。在。在上求解之需要上求解之需要m個(gè)定解條件，若定解條件的形式為：個(gè)定解條件，若

3、定解條件的形式為：其中其中 g為為m維向量。則該問題稱為兩點(diǎn)邊值問題（維向量。則該問題稱為兩點(diǎn)邊值問題（TPBVP，Two Point Boundary Value Problem）。）。 如果邊值條件形式可寫為：如果邊值條件形式可寫為：其中其中g(shù)L和和gR的維數(shù)之和等于的維數(shù)之和等于m，則邊界條件為分離的。，則邊界條件為分離的。 xxx( )( , ( )yfy xydy d xa b, ab( ( ), ( )0g yy 5.1 常微分方程邊值問題的概念常微分方程邊值問題的概念aa( ( )0,( ( )0LRgygy 航空航天中的計(jì)算方法Page 62022-3-175.2 打靶法打靶法

4、 以二階系統(tǒng)為例，考慮邊值問題：以二階系統(tǒng)為例，考慮邊值問題：變換：變換：考慮初值問題：考慮初值問題：初值問題的解為：初值問題的解為：找到找到滿足：滿足： xxxxxa bab( )( , ( ),( ),( ),( )yfyyyAyB 5.2 打靶法打靶法12yyyy xxxxxa bab1221211( )( )( )( ,),( ),( )yyyfyyyAyB xxa baa1221212,( ,),( ),( )yyyfyyyAy xx12( ; ),( ; )yyb1( ; )yB 如何求如何求？航空航天中的計(jì)算方法Page 72022-3-17打靶法的幾何解釋：打靶法的幾何解釋：5

5、.2 打靶法打靶法打靶：求解初值問題打靶：求解初值問題航空航天中的計(jì)算方法Page 82022-3-175.1.1 割線法割線法 以兩個(gè)不同的以兩個(gè)不同的值求解初值問題，得到兩個(gè)解：值求解初值問題，得到兩個(gè)解：根據(jù)初值條件知：根據(jù)初值條件知：假設(shè)假設(shè) 是是的線性函數(shù)，可取的線性函數(shù)，可取 為：為：迭代求解公式：迭代求解公式：結(jié)束條件：結(jié)束條件：5.2 打靶法打靶法 bbb1020101110( ;)( ;)( ;)Byyy b1( ; )yxx1011( ;),( ;)yyaa1011( ;)( ;)yyA bbb11111111( ;)( ;)( ;)mmmmmmmByyy b111( ;)

6、myB 航空航天中的計(jì)算方法Page 92022-3-17割線法的幾何解釋：割線法的幾何解釋：5.2 打靶法打靶法線性近似：按割線求根線性近似：按割線求根航空航天中的計(jì)算方法Page 102022-3-175.1.2 牛頓法牛頓法 求解非線性方程（組）：求解非線性方程（組）：在已知初值在已知初值0的處的處Taylor展開：展開：線性近似：線性近似：迭代求解公式：迭代求解公式：結(jié)束條件：結(jié)束條件：5.2 打靶法打靶法 bb110100( ;)( ;)yBy b1( ; )yB bbb21111001010( ;)( ;)( ;)yyyOB bb111( ;)( ;)mmmmyBy b111( ;

7、)myB b1( ; ) ?y 航空航天中的計(jì)算方法Page 112022-3-17差分法求偏導(dǎo)數(shù)差分法求偏導(dǎo)數(shù)或采用其它數(shù)值微分方法?；虿捎闷渌鼣?shù)值微分方法。f 可微時(shí)解偏導(dǎo)數(shù)微分方程可微時(shí)解偏導(dǎo)數(shù)微分方程微分方程對微分方程對求偏導(dǎo)：求偏導(dǎo)：5.2 打靶法打靶法bbb11101010( ;)( ;)( ;)yyy xxa baa1221212,( ,),( ),( )yyyfyyyAy xx12( ; ),( ; )yy xa baa122121212,( ; )0,( ; )1yyyyyffyyyy 初值問題，可解！初值問題，可解?。ㄅc割線法等價(jià)）（與割線法等價(jià)）割線代替切線割線代替切線航

8、空航天中的計(jì)算方法Page 122022-3-17每一步迭代求解初值問題每一步迭代求解初值問題其中：其中：解得：解得：得到的終端值和對得到的終端值和對的偏導(dǎo)數(shù)：的偏導(dǎo)數(shù)：5.2 打靶法打靶法axaaa1212122121212212,( )( ,),( ),( )0,( )1yyyAyfyyyzzzffzzzzyy xa b, 1212,yyzz xxxx1212( ; ),( ; ),( ; ),( ; )yyzzbb11( ; ),( ; )yy 航空航天中的計(jì)算方法Page 132022-3-17作業(yè)題作業(yè)題5：用牛頓打靶法求解兩點(diǎn)邊值問題用牛頓打靶法求解兩點(diǎn)邊值問題迭代初始條件取迭代初

9、始條件取 。5.2 打靶法打靶法xxxxx2222sin(ln ),12(1)1,(2)2yyyyy (1)0y 航空航天中的計(jì)算方法Page 142022-3-175.3 有限差分法有限差分法 以二階系統(tǒng)為例，邊值問題：以二階系統(tǒng)為例，邊值問題：有限差分近似有限差分近似將區(qū)間將區(qū)間 等分為等分為N個(gè)子區(qū)間個(gè)子區(qū)間將將 在在xi處處Taylor展開：展開： xxxxxa bab( )( , ( ),( ),( ),( )yfyyyAyB 5.3 有限差分法有限差分法 baa,0,1,2,ihxihiNN x( )y xa b, 2311( )()()()()23!iiiiiiiy xy xy

10、xxxyxxxyxxx用差分近似代替微分，將微用差分近似代替微分，將微分方程化為代數(shù)方程求解分方程化為代數(shù)方程求解航空航天中的計(jì)算方法Page 152022-3-17若取若取x=xi+1=x+ih：忽略二階以上部分，得一階導(dǎo)數(shù)的忽略二階以上部分，得一階導(dǎo)數(shù)的前向差分前向差分近似：近似：若取若取x=xi-1=x-ih：忽略二階以上部分，得一階導(dǎo)數(shù)的忽略二階以上部分，得一階導(dǎo)數(shù)的后向差分后向差分近似：近似：5.3 有限差分法有限差分法1()()()iiiy xy xy xh 23111()()()()()26iiiiiy xy xy x hh yxh yx 23111()()()()()26iii

11、iiy xy xy x hh yxh yx 1()()()iiiy xy xy xh 一階精度一階精度一階精度一階精度航空航天中的計(jì)算方法Page 162022-3-17xi+1和和xi-1在在xi處的處的Taylor展開相減，忽略三階以上部分，得展開相減，忽略三階以上部分，得一階導(dǎo)數(shù)一階導(dǎo)數(shù)的中心差分的中心差分近似：近似：xi+1和和xi-1在在xi處的處的Taylor展開相加，忽略四階以上部分，得展開相加，忽略四階以上部分，得二階導(dǎo)數(shù)二階導(dǎo)數(shù)中心差分中心差分近似：近似：三階導(dǎo)數(shù)的中心差分近似？三階導(dǎo)數(shù)的中心差分近似？5.3 有限差分法有限差分法11()()()2iiiy xy xy xh

12、二階精度二階精度112()()2 ()()iiiiy xy xy xyxh 二階精度二階精度航空航天中的計(jì)算方法Page 172022-3-17xi+1和和xi-1在在xi處的處的Taylor展開相減，忽略五階以上部分：展開相減，忽略五階以上部分：xi+2和和xi-2在在xi處的處的Taylor展開相減，忽略五階以上部分：展開相減，忽略五階以上部分：三階導(dǎo)數(shù)的中心差分近似：三階導(dǎo)數(shù)的中心差分近似：四階導(dǎo)數(shù)的中心差分近似：四階導(dǎo)數(shù)的中心差分近似：5.3 有限差分法有限差分法21123()2 ()2 ()()()iiiiiy xy xy xy xyxh 二階精度二階精度35111()()2()()

13、()3iiiiy xy xy x hh yxO h 35228()()4()()()3iiiiy xy xhy xh yxO h (4)21124()4 ()6 ()4 ()()()iiiiiiy xy xy xy xy xyxh 二階精度二階精度航空航天中的計(jì)算方法Page 182022-3-17有限差分法解微分方程兩點(diǎn)邊值問題有限差分法解微分方程兩點(diǎn)邊值問題微分方程微分方程離散化，將區(qū)間離散化，將區(qū)間 等分為等分為N個(gè)子區(qū)間：個(gè)子區(qū)間：在節(jié)點(diǎn)上應(yīng)用中心差分公式，得到代數(shù)方程組：在節(jié)點(diǎn)上應(yīng)用中心差分公式，得到代數(shù)方程組： xxxxxa bab( )( , ( ),( ),( ),( )yfy

14、yyAyB 5.3 有限差分法有限差分法 x1111202(,),1,2,12,iiiiiiiNyyyyyfyiNhhyAyB baa,0,1,2,ihxihiNN xa b, 航空航天中的計(jì)算方法Page 192022-3-17有限差分法解微分方程兩點(diǎn)邊值問題的幾何解釋有限差分法解微分方程兩點(diǎn)邊值問題的幾何解釋5.3 有限差分法有限差分法離散點(diǎn)：微分用有限差分近似離散點(diǎn)：微分用有限差分近似航空航天中的計(jì)算方法Page 202022-3-17例例5.1：用有限差分法求解兩點(diǎn)邊值問題：用有限差分法求解兩點(diǎn)邊值問題取離散化區(qū)間取離散化區(qū)間h=0.1，N=10。xxxxx2222sin(ln ),1

15、2(1)1,(2)2yyyyy 5.3 有限差分法有限差分法 x1111202(,),1,2,12,iiiiiiiNyyyyyfyiNhhyAyB xxxxx111122202sin(ln)2221,21,1,2,1iiiiiiiiiiNiyyyyyyhhyyihiN 航空航天中的計(jì)算方法Page 212022-3-17線性方程組：線性方程組：即：即：5.3 有限差分法有限差分法 xxxxxx2211220sin(ln)12 111,21,1,2,1iiiiiiiiNihhhhyyyyyihiNxxxxxxxx10112222111()1()sin(ln ),( )()1NNNNHYGhgyy

16、yghYGgyhgy 航空航天中的計(jì)算方法Page 222022-3-175.4 有限元法有限元法 以二階系統(tǒng)為例，考慮邊值問題：以二階系統(tǒng)為例，考慮邊值問題：5.4.1 投影類方法的基本思想投影類方法的基本思想 以一簡單函數(shù)以一簡單函數(shù) 近似近似y(x)，給出連續(xù)近似解，例如：，給出連續(xù)近似解，例如：一般形式：一般形式： ， 已知，已知， 待定。待定。殘差：殘差：某種意義上使殘差最小，則得到某種準(zhǔn)則下最佳的近似解。某種意義上使殘差最小，則得到某種準(zhǔn)則下最佳的近似解。 xxxxxa bab( )( , ( ),( ),( ),( )yfyyyAyB 5.4 有限元法有限元法x( ) yxxxx

17、x( )( )( , ( ),( )Ryfyyx10( )Nkkkyc x x0( )kkkyc x x 10( )( )Nkkkyx k( )kx航空航天中的計(jì)算方法Page 232022-3-17區(qū)間殘差平方和最小：最小二乘法區(qū)間殘差平方和最?。鹤钚《朔ㄈ舾商囟c(diǎn)處殘差為零：配點(diǎn)法若干特定點(diǎn)處殘差為零：配點(diǎn)法加權(quán)殘差為零：加權(quán)殘差法加權(quán)殘差為零：加權(quán)殘差法Galerkin法：法： 。5.4 有限元法有限元法x012,min( )NbaRdx 12()0,1,2,iNRiNabxx( ) ( )0,1,2,biawRdxiN ( )iw x計(jì)算復(fù)雜，不常用計(jì)算復(fù)雜，不常用為權(quán)函數(shù)為權(quán)函數(shù)(

18、 )( )iiw xx 航空航天中的計(jì)算方法Page 242022-3-17例例5.2：考慮兩點(diǎn)邊值問題：考慮兩點(diǎn)邊值問題解析解為：解析解為：試用配點(diǎn)法和加強(qiáng)殘差法求解該問題近似解。試用配點(diǎn)法和加強(qiáng)殘差法求解該問題近似解。x50,11( 1)100,(1)100 xyeyy 5050sinh(1),10050cosh(1)11sinh(),cosh()22xxxxxyeaxbabxeexee 5.4 有限元法有限元法航空航天中的計(jì)算方法Page 252022-3-175.4 有限元法有限元法解析解解析解航空航天中的計(jì)算方法Page 262022-3-17設(shè)近似解的形式：設(shè)近似解的形式：基函數(shù)的

19、選擇示例：基函數(shù)的選擇示例：為滿足邊值條件要求為滿足邊值條件要求取二次函數(shù)取二次函數(shù)以及三次項(xiàng)以及三次項(xiàng)N=2（1）配點(diǎn)法）配點(diǎn)法近似解的殘差近似解的殘差令令N個(gè)點(diǎn)處殘差為零求解系數(shù)，如個(gè)點(diǎn)處殘差為零求解系數(shù)，如x 1( )100( )Niiiyx 5.4 有限元法有限元法( 1)0i 1( )(1)(1)xxx 2( )(1)(1)xxxx 線性函數(shù)不滿足線性函數(shù)不滿足xx12( )502650 xxRyee x1 3 121125cosh,25sinh33 配點(diǎn)？配點(diǎn)？航空航天中的計(jì)算方法Page 272022-3-17（2）加權(quán)殘差法）加權(quán)殘差法要求：要求：Galerkin法，取法，取即

20、：即：5.4 有限元法有限元法111212587575,22eee 1122( )( ),( )( )wxxwxx xx( ) ( )0,1,2biawRdxi x1121(1)(1)( 2650)0 xxxedx x1121(1)(1)( 2650)0 xxxxedx 航空航天中的計(jì)算方法Page 282022-3-175.4 有限元法有限元法-1-0.8-0.6-0.4-0.200.20.40.60.8195100105110115120125130 exactCollocationGalerkin配點(diǎn)法、配點(diǎn)法、Galerkin加權(quán)殘差法與精確解的比較加權(quán)殘差法與精確解的比較航空航天中的

21、計(jì)算方法Page 292022-3-175.4.2 有限元法的基本思想有限元法的基本思想 將區(qū)域（區(qū)間）劃分為小的單元，在單元上表示近似解將區(qū)域（區(qū)間）劃分為小的單元，在單元上表示近似解以及求殘差加權(quán)積分。以及求殘差加權(quán)積分。 第第i個(gè)單元，個(gè)單元，Ei， ， 在每個(gè)單元上解用多項(xiàng)式近似；在每個(gè)單元上解用多項(xiàng)式近似； 在每個(gè)單元上計(jì)算加權(quán)殘差；在每個(gè)單元上計(jì)算加權(quán)殘差； 根據(jù)各單元滿足的方程確定多項(xiàng)式近似解的系數(shù)。根據(jù)各單元滿足的方程確定多項(xiàng)式近似解的系數(shù)。5.4 有限元法有限元法xxx1ii xxx1iii 局部近似，分段光滑局部近似，分段光滑可以用簡單的低階近似可以用簡單的低階近似航空航天

22、中的計(jì)算方法Page 302022-3-175.4.3 有限元法：線性元為例有限元法：線性元為例 解在每個(gè)單元上采用解在每個(gè)單元上采用x的線性函數(shù)近似表示。的線性函數(shù)近似表示。5.4 有限元法有限元法航空航天中的計(jì)算方法Page 312022-3-17 線性函數(shù)具有線性函數(shù)具有2個(gè)自由度：由兩個(gè)端點(diǎn)的函數(shù)值確定；個(gè)自由度：由兩個(gè)端點(diǎn)的函數(shù)值確定； N個(gè)線性單元，近似連續(xù)函數(shù)，個(gè)線性單元，近似連續(xù)函數(shù)，N+1個(gè)自由度：由個(gè)自由度：由N+1個(gè)節(jié)點(diǎn)的函數(shù)值唯一確定。個(gè)節(jié)點(diǎn)的函數(shù)值唯一確定。 設(shè)近似解表達(dá)為：設(shè)近似解表達(dá)為： 由由 可知可知5.4 有限元法有限元法x 11( )( )Niiiyx xx

23、 11()(),1,2,1NjiijiyjN x121,1,2,1(),1,2,1iNjjL yyyiNyyjN 線性函數(shù)線性函數(shù)航空航天中的計(jì)算方法Page 322022-3-17 近似解可以用節(jié)點(diǎn)基函數(shù)表示為：近似解可以用節(jié)點(diǎn)基函數(shù)表示為：節(jié)點(diǎn)基函數(shù)節(jié)點(diǎn)基函數(shù) 在節(jié)點(diǎn)在節(jié)點(diǎn) j 處：處： 由于由于 近似表達(dá)式中近似表達(dá)式中 取值的任意性，可知：取值的任意性，可知：5.4 有限元法有限元法xx11( )( )Niiiyy xxxx111111()()()()NjiijijNiijjjjiijiijyyyyy 節(jié)點(diǎn)基函數(shù)的特性節(jié)點(diǎn)基函數(shù)的特性x( ) yx1,()0,ijijij 航空航天中的

24、計(jì)算方法Page 332022-3-17 對于線性元，節(jié)點(diǎn)基函數(shù)在每個(gè)單元內(nèi)是線性函數(shù)：對于線性元，節(jié)點(diǎn)基函數(shù)在每個(gè)單元內(nèi)是線性函數(shù)： 5.4 有限元法有限元法xxxxxxxxxxxxxxxxx11111110()( )()0(2,3,)iiiiiiiiiiiiN 第第i 個(gè)節(jié)點(diǎn)基函個(gè)節(jié)點(diǎn)基函數(shù)的幾何表示數(shù)的幾何表示xxxxxx111iiiiii xxxxxxxxx211212()( )0 xxxxxxxxx110( )NNNNNN 航空航天中的計(jì)算方法Page 342022-3-17節(jié)點(diǎn)函數(shù)值的求解：加權(quán)殘差節(jié)點(diǎn)函數(shù)值的求解：加權(quán)殘差Galerkin法法 近似解的節(jié)點(diǎn)基函數(shù)表示：近似解的節(jié)點(diǎn)

25、基函數(shù)表示： Galerkin法求法求 ，即解方程組：，即解方程組：其中：其中： （對于示例二階系統(tǒng)）（對于示例二階系統(tǒng)）即：即：分部積分：分部積分：5.4 有限元法有限元法xxx11( )() ( )Niiiyy x()iyxx( ) ( )0,1,2,1biaRdxiN xxxxx( )( )( , ( ),( )Ryfyy()0,1,2,1bbbiiiaaayydxfdxiN 0,1,2,1bbiiaadyfdxiN 航空航天中的計(jì)算方法Page 352022-3-17在第在第i個(gè)單元內(nèi)，個(gè)單元內(nèi)，Ei： 5.4 有限元法有限元法單元內(nèi)單元內(nèi) 為連續(xù)函數(shù)為連續(xù)函數(shù)x( ) yxxxxxx

26、xxxxxxx1111101( )10(2,3,)iiiiiiiiiiN xxxxxxx112121( )0 xxxxxxx110( )1NNNNN xxxxxx11111( )() ( )( )( )Niiiiiiiiiiyyyyyy xxx1ii 航空航天中的計(jì)算方法Page 362022-3-17方程組：方程組： 當(dāng)當(dāng) 時(shí)方程為：時(shí)方程為： 當(dāng)當(dāng) 時(shí)方程為：時(shí)方程為： 當(dāng)當(dāng) 時(shí)方程為：時(shí)方程為：5.4 有限元法有限元法1i xxx211122()0 xaxyyyfdx ()0,1,2,1bbbiiiaaayydxfdxiN 2,3,iN 11110iiiixxiixxydxfdx 221

27、1110 xxaxxyydxfdx xxxxxxxxxx111111111()()0iiiiiiiiiiixxiiiixxyyyfdxfdx 1iN11110NNNNxxbNNxxyydxfdx xxx11()0NNxNbNNNxyyyfdx 航空航天中的計(jì)算方法Page 372022-3-17 給定邊值條件：給定邊值條件： 上述上述N+1個(gè)方程可解出個(gè)方程可解出N+1個(gè)未知量：個(gè)未知量： 最終得到問題的解：最終得到問題的解：5.4 有限元法有限元法xx11(),()NyAyB xaxb11,N xxx11( )() ( )Niiiyy xxx23,(),(),(),aNbyyyyy 航空航天

28、中的計(jì)算方法Page 382022-3-17例例5.2（續(xù)）：兩點(diǎn)邊值問題（續(xù)）：兩點(diǎn)邊值問題用線性元和加權(quán)殘差用線性元和加權(quán)殘差Galerkin法得到法得到N+1個(gè)方程：個(gè)方程：給定：給定： ，可求解剩余，可求解剩余N+1個(gè)未知數(shù)。個(gè)未知數(shù)。 x50,11( 1)100,(1)100 xyeyy 5.4 有限元法有限元法 xxxxxxxxxxxxxxxx2111111221111111150()50()50(),2,3,50()0iiiiNNxxaxiiiiiiixxxxiiiixxxxNbNNNxyyye dxyyye dxe dxiNyyye dx xxx( , ( ),( )50 xfyye 11100,100Nyy 航空航天中的計(jì)算方法Page 392022-3-17高斯求積高斯求積 將積分表示為被積函數(shù)在若干點(diǎn)處的函數(shù)值加權(quán)和：將積分表示為被積函數(shù)在若干點(diǎn)處的函數(shù)值加權(quán)和：若適當(dāng)選取若適當(dāng)選取 和和 ，可使公式對次數(shù)，可使公式對次數(shù) 2N+1的多項(xiàng)式被積的多項(xiàng)式被積函數(shù)均精確成立（具有函數(shù)均精確成立（具有2N+1次代數(shù)精度）