有限差分法求解偏微分方程

1、范文范例指導(dǎo)參考有限差分法求解偏微分方程摘要：本文主要使用有限差分法求解計(jì)算力學(xué)中的系統(tǒng)數(shù)學(xué)模型，推導(dǎo)了有限差分法的理論基礎(chǔ)，并在此基礎(chǔ)上給出了部分有限差分法求解偏微分方程的算例驗(yàn)證了推導(dǎo)的正確性及 操作可行性。關(guān)鍵詞：計(jì)算力學(xué)，偏微分方程，有限差分法Abstract : This dissertation mainly focuses on solving the mathematic model ofcomputation mechanics with finite-difference method. The theoretical basis of finite-difference

2、is derived in the second part of the dissertation, and then I use MATLAB to program the algorithms to solve some partial differential equations to confirm the correctness of the derivation and the feasibility of the method.Key words : Computation Mechanics, Partial Differential Equations, Finite-Dif

3、ference Methodword版整理1引言機(jī)械系統(tǒng)設(shè)計(jì)常常需要從力學(xué)觀點(diǎn)進(jìn)行結(jié)構(gòu)設(shè)計(jì)以及結(jié)構(gòu)分析，而這些分析的前提就是建立工程問(wèn)題的數(shù)學(xué)模型。通過(guò)對(duì)機(jī)械系統(tǒng)應(yīng)用自然的基本定律和原 理得到帶有相關(guān)邊界條件和初始條件的微分積分方程，這些微分積分方程構(gòu)成了 系統(tǒng)的數(shù)學(xué)模型。求解這些數(shù)學(xué)模型的方法大致分為解析法和數(shù)值法兩種，而解析法的局限性眾所周知，當(dāng)系統(tǒng)的邊界條件和受載情況復(fù)雜一點(diǎn)， 往往求不出問(wèn)題的解析解或 近似解。另一方面，計(jì)算機(jī)技術(shù)的發(fā)展使得計(jì)算更精確、更迅速。因此，對(duì)于絕 大多數(shù)工程問(wèn)題，研究其數(shù)值解法更具有實(shí)用價(jià)值。 對(duì)于微分方程而言，主要分 為差分法和積分法兩種，本論文主要討論差分

4、法。2有限差分法理論基礎(chǔ)2.1 有限差分法的基本思想當(dāng)系統(tǒng)的數(shù)學(xué)模型建立后，我們面對(duì)的主要問(wèn)題就是微分積分方程的求解。 基本思想是用離散的只含有限個(gè)未知量的差分方程組去近似地代替連續(xù)變量的 微分方程和定解條件，并把差分方程組的解作為微分方程定解問(wèn)題的近似解。將 原方程及邊界條件中的微分用差分來(lái)近似， 對(duì)于方程中的積分用求和或及機(jī)械求 積公式來(lái)近似代替，從而把原微分積分方程和邊界條件轉(zhuǎn)化成差分方程組。有限 差分法求解偏微分方程的步驟主要有以下幾步：區(qū)域離散，即把所給偏微分方程的求解區(qū)域細(xì)分成由有限個(gè)格點(diǎn)組成的網(wǎng)格，這些離散點(diǎn)稱作網(wǎng)格的節(jié)點(diǎn)；近似替代，即采用有限差分公式替代每一個(gè)格點(diǎn)的導(dǎo)數(shù)；逼近求

5、解，換而言之，這一過(guò)程可以看作是用一個(gè)插值多項(xiàng)式及其微分來(lái)代替偏微分方程的解的過(guò)程。從原則上說(shuō)，這種方法仍然可以達(dá)到任意滿意的計(jì)算精度。因?yàn)榉匠痰倪B續(xù) 數(shù)值解可以通過(guò)減小獨(dú)立變量離散取值的問(wèn)格， 或者通過(guò)離散點(diǎn)上的函數(shù)值進(jìn)行 插值計(jì)算來(lái)近似得到。理論上，當(dāng)網(wǎng)格步長(zhǎng)趨近于零時(shí)，差分方程組的解應(yīng)該收 斂于精確解，但由于機(jī)器字節(jié)的限制，網(wǎng)格步長(zhǎng)不可能也沒(méi)有必要取得無(wú)限小，那么差分法的收斂性或者說(shuō)算法的穩(wěn)定性就顯得至關(guān)重要。因此，在運(yùn)用有限差分法時(shí)，除了要保證精度外，還必須要保證其收斂性。2.2 系統(tǒng)微分方程的一般形式由于大多數(shù)工程問(wèn)題都是二維問(wèn)題， 所以得到的微分方程一般都是偏微分方 程，對(duì)于一維問(wèn)

6、題得到的是常微分方程， 解法與偏微分方程類似，故為了不是一 般性，這里只討論偏微分方程。由于工程中高階偏微分較少出現(xiàn)， 所以本文僅僅 給出二階偏微分方程的一般形式， 對(duì)于高階的偏微分，可進(jìn)行類似地推廣。二階 偏微分方程的一般形式如下：（1）其中，為彈性體上的某一特征物理量（連續(xù)函數(shù)）。當(dāng)A、B、C都是常數(shù)時(shí)，（1）式稱為準(zhǔn)線性，有三種準(zhǔn)線性方程形式：如果，則稱為橢圓型方程；如果，則稱為拋物型方程；如果，則稱為雙曲型方程。橢圓型方程主要用來(lái)處理穩(wěn)態(tài)或靜態(tài)問(wèn)題，如熱傳導(dǎo)等問(wèn)題；拋物線方程主 要用來(lái)處理瞬態(tài)問(wèn)題，如滲透、擴(kuò)散等問(wèn)題；雙曲型方程主要用來(lái)處理振動(dòng)問(wèn)題， 如玄震動(dòng)、薄膜震動(dòng)等問(wèn)題。除了上述

7、微分方程外，必須給出定解條件，通常有如下三類：第一類邊界條件（Dirichlet 條件）：；第二類邊界條件（Neuman除件）：一；第三類邊界條件（Robin條件）：一；其中，為求解域 的邊界，為的單位外法矢，。第二類和第三類邊界條件統(tǒng)稱為導(dǎo)數(shù)邊界條件。2.3有限差分方程的數(shù)學(xué)基礎(chǔ)2.3.1 一元函數(shù)導(dǎo)數(shù)的差分公式一個(gè)函數(shù)在x點(diǎn)上的導(dǎo)數(shù)，可以近似地用它所臨近的兩點(diǎn)上的函數(shù)值的差分范文范例指導(dǎo)參考來(lái)表示。函數(shù)在 處的泰勒展式如下:(2)（10）word版整理將區(qū)間等距劃分，我們得到一系列節(jié)點(diǎn)：, , （一），與節(jié)點(diǎn)相鄰的節(jié)點(diǎn)有和 ，因此在點(diǎn)處對(duì)一個(gè)單變量函數(shù)，以步長(zhǎng)， ，然后求出在這些節(jié)點(diǎn)上的近

8、似值可以構(gòu)造如下形式的展開(kāi)式：（3）一（4）由式（3）和式（4）可得到：一階向前差分： （5）一階向后差分： （6）一階中心差分： 不妨，記 ,則式（5）、（6）、（7）分別簡(jiǎn)寫(xiě)為：一階向前差分：（8）一階向后差分：（9）一階中心差分：根據(jù)式（8）、式（9）和式（10）,可得二階差分：范文范例指導(dǎo)參考二階向前差分： (11)二階向后差分： (12)二階中心差分：(13)差分公式(13)是以相隔2h的兩結(jié)點(diǎn)處的函數(shù)值來(lái)表示中間結(jié)點(diǎn)處的一階 導(dǎo)數(shù)值，可稱為中點(diǎn)導(dǎo)數(shù)公式。式(11)和式(12)是以相鄰三結(jié)點(diǎn)處的函數(shù)值 來(lái)表示一個(gè)端點(diǎn)處的一階導(dǎo)數(shù)值， 可稱為端點(diǎn)導(dǎo)數(shù)公式。應(yīng)當(dāng)指出：中點(diǎn)導(dǎo)數(shù)公 式與端點(diǎn)

9、導(dǎo)數(shù)公式相比，精度較高。因?yàn)榍罢叻从沉私Y(jié)點(diǎn)兩邊的函數(shù)變化，而后 者卻只反映了結(jié)點(diǎn)一邊的函數(shù)變化。 因此，我們總是盡可能應(yīng)用前者，而只有在 無(wú)法應(yīng)用前者時(shí)才不得不應(yīng)用后者。但是，由于式(11)中的各階導(dǎo)數(shù)均使用的是向前差分，導(dǎo)致用到的節(jié)點(diǎn)不 相鄰，同時(shí)為了均衡誤差，將節(jié)點(diǎn) 處用到的一階差分換成向后差分，則式(11) 修正為： (14)同理，根據(jù)上述推導(dǎo)過(guò)程，可得到任意階的差分公式：n階向前差分： (15)n階向后差分： (16)n階中心差分： (17)說(shuō)明，上述公式中各節(jié)點(diǎn)處前一階導(dǎo)數(shù)的代入可能存在不一致，可能是向前 差分、向后差分或者中心差分，從而使最終的公式在系數(shù)上存在差別。當(dāng)然，也 可以對(duì)

10、各相鄰節(jié)點(diǎn)進(jìn)行需要階數(shù)的泰勒展開(kāi)，從而建立方程組直接求各階導(dǎo)數(shù)。word版整理范文范例指導(dǎo)參考2.3.2微分方程轉(zhuǎn)化為線性方程由于三種類型的微分方程解法類似，故這里僅以橢圓型微分方程為例，將微 分方程轉(zhuǎn)化為代數(shù)方程，對(duì)于雙曲型和拋物型方程依次類推即可。不妨記：（稱為拉普拉斯算子），和是求解域上的連續(xù)函數(shù)。假設(shè)求解區(qū)域?yàn)椋?，將求解區(qū)域劃分成個(gè)網(wǎng)格,其中： ，如圖1所示。記 ，則根據(jù)式（14）可得到：(18)444卜圖1五點(diǎn)差分公式式（18）也稱為五點(diǎn)差分公式，同理根據(jù)式（12）和式（13）可分別得到向 前差分公式（19）和向后差分公式（20）,如圖（2所示）。向前差分(19)向后差分'1

11、1F-七圖2向前差分（左）圖3中心差分、向前差分和向后差分的拉普拉斯算子表示利用中心差分公式（18）,由于式（18）在點(diǎn)處具有二階精度（）,所以式（18）可近似改寫(xiě)成下式：（21）根據(jù)橢圓方程的具體形式可以將其分為以下三種形式：拉普拉斯（Laplace）方程：泊松（Poison）方程：赫耳墨次（Helmholtz ）方程：根據(jù)式（21）,可建立三種不同形式橢圓方程的代數(shù)方程如下：拉普拉斯方程：化簡(jiǎn)后得到拉普拉斯方程的計(jì)算公式:（22）泊松方程:(23)word版整理范文范例指導(dǎo)參考赫耳墨次方程：（24）2.3.3建立有限差分方程組根據(jù)式（22） （24）建立方程組，但是需要知道對(duì)應(yīng)的邊界條件才

12、能使方 程組存在定解，根據(jù)2.2中可知，邊界條件一般分為狄利克雷邊界條件和導(dǎo)數(shù)邊 界條件兩種，下面分別給出這兩種邊界條件的有限差分方程組的建立過(guò)程：狄利克雷邊界條件：對(duì)于狄氏條件而言，給出了邊界上各節(jié)點(diǎn)出的函數(shù)計(jì)算公式，直接代入節(jié)點(diǎn)值 計(jì)算即可，如下所示為矩形區(qū)域的邊界點(diǎn)計(jì)算：（左邊界）（右邊界）（25）（下邊界）（上邊界）導(dǎo)數(shù)邊界條件：一以右邊界點(diǎn)為例，對(duì)于右邊界點(diǎn)，根據(jù)Neuman僚件可得下式：一 一（26）對(duì)于拉普拉斯方程，根據(jù)計(jì)算公式（22）,對(duì)于邊界上的點(diǎn) 可得：（27）顯然，上式中的在求解域外，是未知量。根據(jù)中心差分公式（10）可得到：（28）根據(jù)式（28）可得到逼近表示：，并且具

13、有2階逼近精度，代入式（27）可得下式：同理，對(duì)于其它邊界可獲得如下邊界方程：word版整理（29）范文范例指導(dǎo)參考word版整理(30)（下邊界）（上邊界）對(duì)于泊松方程和赫耳墨次方程同樣根據(jù)上述方法，獲得邊界條件的線性方程，然后將這些方程添加到式（22） （24）所建立的方程組中，從而建立起個(gè)元的線性方程組，解該方程組即可獲得各節(jié)點(diǎn)的函數(shù)值。對(duì)于上述過(guò)程建立的線性方程組的求解，可采用多種方法，比如Jacob迭代法、Gauss-Seidel迭代法、超松弛迭代法（SOFRS）、高斯消元法等方法求解。2.4有限差分法的收斂性和穩(wěn)定性由于迭代法必須保證收斂性，所以在解有限差分方程組時(shí)還應(yīng)保證其收斂性

14、， 也就是通常所說(shuō)的算法穩(wěn)定性。有限差分法的算法穩(wěn)定性可以通過(guò)特征值方法、 傅里葉變換（馮諾依曼條件）以及能量估計(jì)等方法來(lái)判斷，下面給出常用的馮諾 依曼條件：向前差分：，絕對(duì)收斂；向后差分：，絕對(duì)收斂；中心差分：對(duì)任何的 對(duì)不收斂；假設(shè)求解域內(nèi)方向網(wǎng)格劃分的步長(zhǎng)為，方向網(wǎng)格劃分的步長(zhǎng)為，將偏微分方 程化為標(biāo)準(zhǔn)形式，具體來(lái)說(shuō)標(biāo)準(zhǔn)形式如下：雙曲方程：_（31）對(duì)于式（31）所示的雙曲方程，馮諾依曼條件為：拋物方程：-（32）對(duì)于式（32）所示的拋物方程，馮諾依曼條件為：橢圓方程：_ _（33）對(duì)于式（33）所示的橢圓方程，馮諾依曼條件為：為了使算法在任何情況下都能保持穩(wěn)定性， 去掉對(duì)網(wǎng)格劃分的馮諾

15、依曼條件， 通常采用隱式方案，對(duì)五點(diǎn)差分公式中的節(jié)點(diǎn)所在的行做差分， 然后把這些差分 的加權(quán)作為中心點(diǎn)的差分值，則拉普拉斯算子可修正為：）（34）利用式（34）進(jìn)行計(jì)算時(shí)，穩(wěn)定性沒(méi)有任何限制。取不同的值得到不同的差分公式，通常取為了提高計(jì)算精度，很明顯的一個(gè)措施就是網(wǎng)格細(xì)分，但是由于隨著網(wǎng)格步 長(zhǎng)的減小，未知量的數(shù)目將會(huì)呈指數(shù)增長(zhǎng)，網(wǎng)格劃分太細(xì)會(huì)導(dǎo)致計(jì)算量過(guò)于龐大 而無(wú)法計(jì)算。通常，我們可以通過(guò)提高逼近的精度，采用更高精度的差分公式， 例如對(duì)公式（21）進(jìn)行修改，可得到九點(diǎn)差分公式：(35)nLM卜4卜 .L、卜卜圖5九點(diǎn)差分公式3有限差分法求解實(shí)例根據(jù)上述推導(dǎo)有限差分法理論，對(duì)于不同類型的偏

16、微分方程建立有限差分方 程組，采用mat lab編程給出一些計(jì)算實(shí)例如下：1.橢圓型方程拉普拉斯方程：;求解域：下面分別給出拉普拉斯方程在不同的邊界條件下的解。a）狄利克雷邊界條件：下邊界：上邊界：左邊界：右邊界：圖6狄利克雷邊界條件下拉普拉斯方程的解b) Neumanrii 界條件:下邊界:上邊界:左邊界:右邊界:圖7 Neumann邊界條件下拉普拉斯方程的解范文范例指導(dǎo)參考泊松方程：;求解域:狄利克雷邊界條件：下邊界：上邊界：左邊界：右邊界：|網(wǎng)格步長(zhǎng)：h=0.05-40 4-圖8狄利克雷邊界條件下泊松方程的解赫耳墨茲方程：；求解域：狄利克雷邊界條件：下邊界：上邊界：左邊界：右邊界：通過(guò)圖9與圖6的對(duì)比發(fā)現(xiàn)，微分方程 的解與微分方程的具體形式關(guān)系不 大，主要由求解域和邊界條件所決定。word版整理2.雙曲方程: