微分方程的數(shù)值解法

微分方程的數(shù)值解法第一頁，共六十二頁，2022年，8月28日定義：當求泛函在一個函數(shù)集合K中的極?。ɑ驑O大）問題，則該問題稱為變分問題。

變分問題與微分方程的定解問題有一定的聯(lián)系。（2）初等變分原理①一元二次函數(shù)的變分原理考察J(x)的極值情況。變分原理：設(shè)求，使與求解方程Lx=f

等價。第二頁，共六十二頁，2022年，8月28日對稱正定②多元二次函數(shù)的變分原理求J(x)取極小值的駐點，其中設(shè)設(shè)則J(x)可表示為：第三頁，共六十二頁，2022年，8月28日變分原理：設(shè)矩陣A對稱正定，則下列兩個命題等價：求，使(a)(b)

是方程的解上述兩個例子表明：其中

求二次函數(shù)的極小值問題和求線性代數(shù)方程（組）的解是等價的。第四頁，共六十二頁，2022年，8月28日（1）弦平衡的平衡原理與極小位能原理2兩點邊值問題的變分原理

考察一根長為l的弦，兩端固定在點A(0,0)和B(l,0)。當沒有外力作用時，它的位置沿水平方向與X軸重合。設(shè)有強度為f(x)的外荷載垂直向下作用在弦上，于是弦發(fā)生形變。假定荷載很小，因而發(fā)生的形變也很小。用u(x)表示在荷載f(x)的作用下弦的平衡位置。第五頁，共六十二頁，2022年，8月28日求弦的平衡位置歸結(jié)為求解兩點邊值問題：設(shè)弦處于某一位置u=u(x)，可得到其總位能為

極小位能原理：其中T是弦的張力。平衡原理

弦的平衡位置(記為)將在滿足邊值條件u(0)=0，u(l)=0的一切可能位置中，使位能取極小值。弦的平衡位置是下列變分問題的解第六頁，共六十二頁，2022年，8月28日

在數(shù)學上，要將某個微分方程的定解問題轉(zhuǎn)化為一個變分問題求解，必須針對已給的定解問題構(gòu)造一個相應(yīng)的泛函，并證明定解問題的解與泛函極值問題的解等價。

有限元方法正是利用這種等價性（邊值問題與變分問題的等價性），先將微分方程定解問題轉(zhuǎn)化為變分問題（或變分方程）的求解問題，然后再設(shè)法近似求解變分問題（或變分方程）。第七頁，共六十二頁，2022年，8月28日（2）兩點邊值問題的變分原理①構(gòu)造泛函考察二階常微分方程邊值問題：引入泛函算子則第八頁，共六十二頁，2022年，8月28日②變分問題

與前述二階常微分方程邊值問題相應(yīng)的變分問題是其中求，使第九頁，共六十二頁，2022年，8月28日③變分原理（變分問題與邊值問題的等價性）設(shè)，是邊值問題的解，則使J(u)

達到極小值；

反之，若使J(u)

達到極小值，則是邊值問題的解。其中

是強制邊界條件，是自然邊界條件，區(qū)別這兩類邊界條件在用有限元方法求解邊值問題時很重要。第十頁，共六十二頁，2022年，8月28日（3）虛功原理對兩點邊值問題：其中虛功原理，且滿足變分方程：

設(shè)，以v乘方程兩端，沿[a,b]積分，并利用，得變分方程對任意在力學里，表示虛功設(shè)，則

是邊值問題解的充要條件是：第十一頁，共六十二頁，2022年，8月28日

對于復雜的邊界條件，邊值問題的求解一般是困難的。若將微分方程化為相應(yīng)的變分問題或變分方程，則只需處理強加邊界條件，無需處理自然邊界條件（自然邊界條件已包含于變分問題中泛函的構(gòu)造或已包含于給出的變分方程之中）。這一特點對研究微分方程離散化方法及其數(shù)值解帶來了極大的方便。第十二頁，共六十二頁，2022年，8月28日3二階橢圓邊值問題的變分原理（1）極小位能原理模型方程其中G是平面有界區(qū)域。①構(gòu)造泛函引入泛函算子則第十三頁，共六十二頁，2022年，8月28日②變分問題與前述二階橢圓邊值問題相應(yīng)的變分問題是求，使其中第十四頁，共六十二頁，2022年，8月28日③變分原理（變分問題與邊值問題的等價性）

對第一邊值問題，無論齊次或非齊次邊界條件，泛函是一樣的，只是邊界條件要作為強加邊值條件加在所取的函數(shù)類上。

設(shè)，是二階橢圓邊值問題的解，則使J(u)達到極小值；

反之，若使J(u)

達到極小值，則是二階橢圓邊值問題的解。其中

對第二、三類邊值問題，無論齊次或非齊次邊界條件，二次泛函形式相對于第一邊值問題有所改變，但函數(shù)類的選取與邊界條件無關(guān)。第十五頁，共六十二頁，2022年，8月28日（2）虛功原理問題其中

設(shè)，以v乘方程兩端后在G上積分，并利用Green公式，得變分方程第十六頁，共六十二頁，2022年，8月28日虛功原理在力學里，表示虛功

設(shè)是邊值問題的解，則對任意，滿足變分方程。

反之，若，且對任意滿足變分方程，則為邊值問題的解。

與極小位能原理類似，第一類邊界條件為強加邊界條件，第二、三類邊界條件為自然邊界條件。

虛功原理比極小位能原理應(yīng)用更廣。第十七頁，共六十二頁，2022年，8月28日目的：求解相應(yīng)的變分問題或相應(yīng)的變分方程。

Ritz方法是近似求解變分問題(即二次泛函極小值)的算法。Galerkin方法是近似求解變分方程的算法，這兩種算法統(tǒng)稱為Ritz-Galerkin方法。Ritz-Galerkin方法的基本思想

以下用V表示等Sobolev空間，L表示微分算子，(u,v)為由L及邊值條件決定的雙線性泛函。4Ritz-Galerkin方法

用有限維空間的函數(shù)代替變分問題(或變分方程)中無限維空間的函數(shù)，從而在有限維函數(shù)空間中求變分問題(或變分方程)的近似解，并要求當有限維空間的維數(shù)不斷增加時，有限維近似解逼近原變分問題(或變分方程)的解。第十八頁，共六十二頁，2022年，8月28日由極小位能原理得出的變分問題為：Ritz方法：求變分問題的近似解。（1）Ritz方法求，使其中，

設(shè)是V

的n維子空間，是的一組基底(稱為基函數(shù))。中任一元素可表示為即選擇適當?shù)?，使取極小值。求，使Ritz方法：第十九頁，共六十二頁，2022年，8月28日展開令則滿足解出代入，則得第二十頁，共六十二頁，2022年，8月28日Ritz方法步驟為：

根據(jù)最小位能原理構(gòu)造相應(yīng)于微分方程或物理問題的變分問題；

取作為的一組基底，即用近似代替無窮維空間V；

根據(jù)二次函數(shù)取極值的必要條件，得到中所滿足的方程組：

求解關(guān)于的線性代數(shù)方程組。第二十一頁，共六十二頁，2022年，8月28日由虛功原理得出的變分方程為：Galerkin方法：求變分方程的近似解。（2）Galerkin方法

設(shè)是V

的n維子空間，是的一組基底(稱為基函數(shù))。中任一元素可表示為即選擇適當?shù)?，使取極小值。Galerkin方法：求，使對，滿足第二十二頁，共六十二頁，2022年，8月28日由的任意性，取作為v

，則得將代入變分方程，則解出代入，則得第二十三頁，共六十二頁，2022年，8月28日Galerkin步驟為：

根據(jù)虛功原理構(gòu)造相應(yīng)于微分方程或物理問題的變分方程；

取作為的一組基底，即用近似代替無窮維空間V；

求解關(guān)于的線性代數(shù)方程組。取作為v

，將代入變分方程，得到滿足的方程組：第二十四頁，共六十二頁，2022年，8月28日有限元法廣泛應(yīng)用的原因Ritz-Galerkin方法應(yīng)用的困難①基函數(shù)選取必須滿足強加邊界條件，因此選取困難；②計算量、存儲量巨大；③方程組求解病態(tài)嚴重。

充分發(fā)揮了變分形式和Ritz-Galerkin方法的優(yōu)點；②擺脫了傳統(tǒng)的基函數(shù)取法；③各種問題的結(jié)構(gòu)程序格式統(tǒng)一。第二十五頁，共六十二頁，2022年，8月28日

有限元方法基于變分原理，又具有差分方法的一些特點，并且適于較復雜的區(qū)域和不同粗細的網(wǎng)格。二橢圓型方程的有限元方法

差分法解偏微分方程，解得的結(jié)果就是準確解u在節(jié)點上的近似值；Ritz-Galerkin方法得到近似的解析解，但對一般區(qū)域，卻往往難以實現(xiàn)。

有限元方法與傳統(tǒng)Ritz-Galerkin方法的差別在于有限維函數(shù)空間的構(gòu)造方法。Ritz-Galerkin方法選用的基函數(shù)在整個定解區(qū)域上整體光滑，有限元則取分段或分片連續(xù)且局部非零的基函數(shù)。第二十六頁，共六十二頁，2022年，8月28日考慮兩點邊值問題：1一維問題的線性元

將區(qū)間[a,b]分割為n個子區(qū)間。第i個單元記為，其長度。（1）試探函數(shù)與試探函數(shù)空間設(shè)則稱為試探函數(shù)空間，稱為試探函數(shù)。第二十七頁，共六十二頁，2022年，8月28日（2）用單元形狀函數(shù)表示試探函數(shù)設(shè)在節(jié)點上試探函數(shù)在節(jié)點上的一組值為最簡單的試探函數(shù)空間由分段線性函數(shù)組成。在第i個單元上的線性插值函數(shù)為即

當時，的（線性）插值公式稱為（線性）單元形狀函數(shù)。第二十八頁，共六十二頁，2022年，8月28日

把每個單元形狀函數(shù)合并起來，就得到整個區(qū)間[a,b]上都有定義的函數(shù)：第二十九頁，共六十二頁，2022年，8月28日為使分段插值標準化，通常用仿射變換顯然把變到，令則變?yōu)榛虻谌?，共六十二頁?022年，8月28日定義基函數(shù)系（3）用節(jié)點基函數(shù)表示試探函數(shù)第三十一頁，共六十二頁，2022年，8月28日

線性無關(guān)，它們可組成試探函數(shù)空間的基，常稱為節(jié)點基函數(shù)。幾何形狀如圖ab任一試探函數(shù)可表示為

用這類插值型基函數(shù)，可以構(gòu)造出適合各種邊界條件的試探函數(shù)。第三十二頁，共六十二頁，2022年，8月28日若借助前述放射變換節(jié)點基函數(shù)可用變量表示為第三十三頁，共六十二頁，2022年，8月28日①直接形成有限元方程(a)把表達式代入泛函；（4）從Ritz方法出發(fā)形成有限元方程(b)將泛函表達式中積分區(qū)間[a,b]變到[0,1]；(c)由達到極小值的條件得到含的有限元方程這兒(d)解出有限元方程的數(shù)值解，就得到使二次泛函取極小的近似函數(shù)(有限元解)第三十四頁，共六十二頁，2022年，8月28日有限元方程可用矩陣表示為其中稱為總剛矩陣。第三十五頁，共六十二頁，2022年，8月28日

工程中形成有限元方程時，通常先在每個單元上形成單元矩陣（稱為單元剛度矩陣），然后由單元剛度矩陣形成總剛度矩陣（稱為總體合成）。②用單元剛度分析形成有限元方程(a)把按單元組織，則在第i個單元上，令其中稱為單元剛度矩陣。各元素可計算得到。第三十六頁，共六十二頁，2022年，8月28日

再把擴展成nn矩陣，使其第i1行、第i行和第i1列、第i列交叉位置的元素就是單元剛度矩陣的四個元素，其余全為零（只是第一行，第一列元素非零）。即記則其中稱為總剛矩陣。第三十七頁，共六十二頁，2022年，8月28日(b)由達到極小值的條件(c)解出有限元方程的數(shù)值解，就得到使二次泛函取極小的近似函數(shù)(有限元解)得到有限元方程。第三十八頁，共六十二頁，2022年，8月28日（5）從Galerkin方法出發(fā)形成有限元方程把表達式代入變分方程對前面的兩點邊值問題，變分方程變?yōu)槠渲?/p>

與Ritz方法相比，Galerkin方法形成的有限元方程其系數(shù)矩陣就是總剛矩陣。該方程即為Galerkin法形成的有限元方程。

由Galerkin方法推導有限元方程更加方便直接，且適用面廣。第三十九頁，共六十二頁，2022年，8月28日

若希望在每個單元上提高逼近的精確度，則可通過提高插值多項式次數(shù)來實現(xiàn)，

在單元上可構(gòu)造一、二、三及高次插值多項式，其方法有兩種：2一維問題的高次元

整個問題計算的全過程除分析單元插值外，均與前面框架類似。①Lagrange型：在單元內(nèi)部增加一些插值節(jié)點。②Hermite型：在節(jié)點引進一階、二階乃至更高階導數(shù)。第四十頁，共六十二頁，2022年，8月28日①線性元(Lagrange型)要求：在每一個單元上是一次多項式，在單元節(jié)點處連續(xù)。插值條件：在單元的兩個端點取指定值。②二次元(Lagrange型)要求：在每一個單元上是二次多項式，在單元節(jié)點處連續(xù)。插值條件：在單元的兩個端點及單元中點取指定值。③三次元（Hermite型）要求：在每一個單元上是三次多項式，在單元節(jié)點處連續(xù)。插值條件：在兩個端點取指定的函數(shù)值和一階導數(shù)值。第四十一頁，共六十二頁，2022年，8月28日

采用高次元，有限元方程形成的方法和線性元類似，但工作量增加。一是計算積分的復雜性增加，二是矩陣的帶寬增加。

高次元的主要優(yōu)點是收斂階高，且提高了函數(shù)逼近的光滑性。第四十二頁，共六十二頁，2022年，8月28日

假定區(qū)域G可以分割成有限個矩形的和，且每個小矩形（單元）的邊和坐標軸平行。3二維問題的矩形元通過仿射變換采用矩形剖分后，任一個矩形總可變成單位正方形

如果在上造出單元形狀函數(shù)，就可得到試探函數(shù)。而上的形狀函數(shù)可通過先在上造出形狀函數(shù)，再通過仿射變化而得到。第四十三頁，共六十二頁，2022年，8月28日

在上構(gòu)造形狀函數(shù)，也采用Lagrange型和Hermite型插值。Lagrange型：根據(jù)若干插值節(jié)點處的函數(shù)值決定插值函數(shù)。Hermite型：根據(jù)若干插值節(jié)點處的函數(shù)值、一階偏導數(shù)乃至更高階偏導數(shù)決定插值函數(shù)。第四十四頁，共六十二頁，2022年，8月28日（1）Lagrange型公式①雙一次插值插值條件：給定頂點上的函數(shù)值求：雙線性函數(shù)滿足設(shè)令由為雙線性函數(shù)，可求得第四十五頁，共六十二頁，2022年，8月28日令則

通過仿射變換消去、

，就得到上的形狀函數(shù)。把這些函數(shù)按單元疊加，即對所有單元求和，就得到G上的試探函數(shù)。

實際計算時，并不消去中間變量、，因為計算剛度矩陣元素（定積分）用、作自變量更為方便。第四十六頁，共六十二頁，2022年，8月28日插值條件：給定II上九個插值節(jié)點(0,0)、(1/2,0)、(1,0)、(0,1/2)、(1/2,1/2)、(1,1/2)、(0,1)、(1/2,1)、(1,1)的函數(shù)值。求：雙二次函數(shù)滿足②雙二次插值第四十七頁，共六十二頁，2022年，8月28日故

通過仿射變換消去、

，就得到上的形狀函數(shù)。令由為二次函數(shù)，可求得設(shè)第四十八頁，共六十二頁，2022年，8月28日插值條件：給定II上十六個插值節(jié)點(見圖)。求：雙三次函數(shù)滿足設(shè)③雙三次插值第四十九頁，共六十二頁，2022年，8月28日故令由為三次函數(shù)，可求得第五十頁，共六十二頁，2022年，8月28日

可以在四個頂點分別給定函數(shù)值、兩個一階偏導數(shù)的值和二階混合偏導數(shù)的值(共十六條件)，確定一個雙三次多項式的十六個系數(shù)。（2）Hermite型公式Lagrange型公式中不出現(xiàn)導數(shù)，這樣的試探函數(shù)只屬于。為了得到屬于的試探函數(shù)，需要Hermite型插值公式。雙三次多項式含有十六項：

簡單且常用的是不完全的雙三次多項式插值。它去掉雙三次多項式中的項。第五十一頁，共六十二頁，2022年，8月28日插值條件：給定II上四個插值節(jié)點。求：不完全雙三次函數(shù)滿足四個頂點處的函數(shù)值等于在該點的函數(shù)值；四個頂點處的值等于在該點的值；四個頂點處的值等于在該點的值。根據(jù)仿射變換則可將原插值問題轉(zhuǎn)化為II上的插值問題。第五十二頁，共六十二頁，2022年，8月28日滿足四個頂點處的函數(shù)值等于在該點的函數(shù)值；四個頂點處的值等于在該點的值乘以x；四個頂點處的值等于在該點的值乘以y。插值條件：給定II上四個插值節(jié)點(0,0)、(1,0)、(0,1)、(1,1)。求：不完全雙三次函數(shù)

類似于Lagrange型公式的構(gòu)造，可以求得上的形狀函數(shù)。第五十三頁，共六十二頁，2022年，8月28日

在三角形元的有限元方法中，先將定解區(qū)域G化分為若干個小三角形（稱作單元）。然后在每個單元上構(gòu)造插值型函數(shù)，并用分片函數(shù)（但整體連續(xù)的函數(shù)）代替變分問題或變分方程中所需求解的函數(shù)。4二維問題的三角形元

用有限元求解二維橢圓邊值問題時，應(yīng)用最廣的是三角形元。第五十四頁，共六十二頁，2022年，8月28日（1）三角剖分

將定解區(qū)域化分成若干個小三角形單元時應(yīng)注意：③為了保證有限元解的精確度和收斂性，并避免其離散后代數(shù)方程組系數(shù)矩陣的病態(tài)性，網(wǎng)格剖分中疏密的過渡不要太陡。錯誤為了保證有限元解有較好的精度，每個單元中應(yīng)盡量避免出現(xiàn)大的鈍角。應(yīng)避免④單元頂點的編號順序可以任意，但節(jié)點編號順序?qū)⒂绊懹邢拊匠探M系數(shù)矩陣的結(jié)構(gòu)（帶寬）。①為了方便構(gòu)造插值型函數(shù)，要求每個單元的頂點是相鄰單元的頂點。第五十五頁，共六十二頁，2022年，8月28日（2）面積坐標及有關(guān)公式

在三角形單元上構(gòu)造插值型函數(shù)，并不簡單類同于矩形單元。①面積坐標

考慮一個面積為S的三角形單元，其頂點按反時針順序記為i,j,k。在此單元內(nèi)部任取一點p(x,y)，連接p和三個頂點，此單元則被分成三個小三角形它們的面積記為和。ijkP(x,y)記單元內(nèi)任一點p(x,y)的位置與三維數(shù)一一對應(yīng)，稱