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

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

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

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

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

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

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

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

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

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

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

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

達(dá)到極小值；

反之，若使J(u)

達(dá)到極小值，則是邊值問(wèn)題的解。其中

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

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

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

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

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

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

反之，若使J(u)

達(dá)到極小值，則是二階橢圓邊值問(wèn)題的解。其中

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

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

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

反之，若，且對(duì)任意滿足變分方程，則為邊值問(wèn)題的解。

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

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

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

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

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

設(shè)是V

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

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

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

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

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

設(shè)是V

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

工程中形成有限元方程時(shí)，通常先在每個(gè)單元上形成單元矩陣（稱為單元?jiǎng)偠染仃嚕?，然后由單元?jiǎng)偠染仃囆纬煽倓偠染仃嚕ǚQ為總體合成）。②用單元?jiǎng)偠确治鲂纬捎邢拊匠?a)把按單元組織，則在第i個(gè)單元上，令其中稱為單元?jiǎng)偠染仃?。各元素可?jì)算得到。第三十六頁(yè)，共六十二頁(yè)，2022年，8月28日

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

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

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

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

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

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

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

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

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

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

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

通過(guò)仿射變換消去、

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

實(shí)際計(jì)算時(shí)，并不消去中間變量、，因?yàn)橛?jì)算剛度矩陣元素（定積分）用、作自變量更為方便。第四十六頁(yè)，共六十二頁(yè)，2022年，8月28日插值條件：給定II上九個(gè)插值節(jié)點(diǎn)(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ù)滿足②雙二次插值第四十七頁(yè)，共六十二頁(yè)，2022年，8月28日故

通過(guò)仿射變換消去、

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

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

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

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

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

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

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

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

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