求解拋物線(xiàn)方程的顯式差分法

1、蘇州大學(xué)本科生畢業(yè)設(shè)計(jì)（論文）目錄摘要2Abstract31. 前言42. 一維熱傳導(dǎo)方程第一邊值問(wèn)題的定義及差分格式的建立52.1 第一邊值問(wèn)題的定義52.2差分格式的建立52.2.1 連續(xù)問(wèn)題離散化：網(wǎng)格剖分52.2.2 顯式差分格式的定義及截?cái)嗾`差72.3 熱傳導(dǎo)方程混合問(wèn)題的差分方法82.3.1顯式差分方程問(wèn)題及向量表示82.3.2 顯式差分格式的穩(wěn)定性92.3.3 顯式差分格式的收斂性133. 數(shù)值算例及圖像153.1 數(shù)值算例154 總結(jié)185 參考文獻(xiàn)196 致謝20摘要在物理學(xué)中，通常采用二階拋物型偏微分方程來(lái)對(duì)熱傳導(dǎo)和擴(kuò)散現(xiàn)象進(jìn)行描述，統(tǒng)稱(chēng)為熱傳導(dǎo)方程。目前常用的求解熱傳導(dǎo)方

2、程的差分格式有六點(diǎn)格式（Crank-Nicolson）、隱式差分格式、顯式差分格式和Richardson格式。本文將使用顯式差分法分析一維熱傳導(dǎo)方程問(wèn)題，給出其截?cái)嗾`差，并證明穩(wěn)定性及收斂性，最后通過(guò)數(shù)值算例和圖像對(duì)誤差進(jìn)行簡(jiǎn)單的分析，從而驗(yàn)證我們的理論分析。關(guān)鍵字：一維熱傳導(dǎo)方程，初邊值條件，顯式差分法，穩(wěn)定性，收斂性。AbstractThe second order parabolic partial differential equation is usually used to describe the phenomena of thermal transmission and dif

3、fusion in physics.At present, there are explicit difference schemes, implicit difference schemes, Crank-Jolson schemes and Richardson schemes for solving thermal transmission equations.In this paper, the explicit difference method is used to analyze the one-dimensional thermal transmission question,

4、 and its steadiness and astringency are proved. The truncation error is given, and the error is analyzed simply by numerical examples and images, so as to verify our theoretical analysis.Keywords:one-dimensional thermal transmission equation, explicit difference method, boundary condition, steadines

5、s, astringency.1. 前言熱傳導(dǎo)是由溫差引起的熱量傳運(yùn)過(guò)程，當(dāng)物體內(nèi)部各點(diǎn)溫度不相同時(shí)，熱量就會(huì)從高溫處向低溫處傳播。熱傳導(dǎo)方程是拋物型偏微分方程最簡(jiǎn)單的例子，用于探究固體內(nèi)部的熱量傳輸，在許多現(xiàn)象的數(shù)學(xué)模型中出現(xiàn)；熱傳導(dǎo)的探究對(duì)工業(yè)和經(jīng)濟(jì)長(zhǎng)遠(yuǎn)發(fā)展也是十分重要的，比如：材料結(jié)構(gòu)熱應(yīng)力計(jì)算、金屬材料在鑄造、鍛壓等過(guò)程中內(nèi)部溫度分布，以及食品的冷凍過(guò)程等等都跟導(dǎo)熱理論關(guān)系密切，同時(shí)它也是目前研究各種傳熱現(xiàn)象及工程熱物理學(xué)科必不可少的工具。所以，探究熱傳導(dǎo)方程具有十分重要的實(shí)踐意義。下面，先考察某個(gè)物體在三維情況下的熱傳導(dǎo)問(wèn)題，物體G在時(shí)刻時(shí)在處的溫度我們用函數(shù)來(lái)表示，假設(shè)關(guān)于具有連續(xù)

6、一階偏導(dǎo)數(shù)，關(guān)于具有連續(xù)二階偏導(dǎo)數(shù)。計(jì)算熱傳導(dǎo)問(wèn)題，還需要一些參數(shù)，那么用表示物體在處的熱傳導(dǎo)系數(shù)，設(shè)物體G的密度為，比熱容為。那么根據(jù)熱量守恒定律、傅里葉定律以及高斯公式得：，假設(shè)此時(shí)物體是均勻的，那么均是常數(shù)，令，化簡(jiǎn)得到：拋物型方程的求解是很有研究意義的，通常情況下，對(duì)于簡(jiǎn)單的方程或定解問(wèn)題，精確解比較容易求得，但是過(guò)于復(fù)雜的方程，精確解很難求甚至求不到，這個(gè)時(shí)候我們就只能想辦法去找拋物型方程定解問(wèn)題的近似解，所以就有必要去研究拋物型方程數(shù)值解。而在目前，已經(jīng)總結(jié)出了求解近似解的方法，即有限差分法和有限元法。本文主要是來(lái)討論有限差分法當(dāng)中的顯式差分。本文主要探究一維熱傳導(dǎo)方程問(wèn)題。全文分

7、為兩部分，第一部分從一維熱傳導(dǎo)方程的定義出發(fā)，介紹了求解熱方程的顯式差分法,所謂顯式差分法，第一步就是把問(wèn)題的定義域劃分成網(wǎng)格，第二步使用數(shù)值微分公式，用差商取代替原來(lái)熱傳導(dǎo)問(wèn)題的微商，就可以把熱傳導(dǎo)問(wèn)題轉(zhuǎn)換成顯式差分格式，再討論其在不同網(wǎng)格上的誤差，第三步證明其穩(wěn)定性及收斂性，第二部分是通過(guò)具體的數(shù)值算例對(duì)誤差進(jìn)行簡(jiǎn)單的分析，驗(yàn)證我們的理論分析。2. 一維熱傳導(dǎo)方程第一邊值問(wèn)題的定義及差分格式的建立2.1 第一邊值問(wèn)題的定義本文討論的一維熱傳導(dǎo)方程第一邊值問(wèn)題如下：在區(qū)域內(nèi)求函數(shù)滿(mǎn)足方程及初始條件邊界條件 并且，其中是給定的函數(shù).2.2差分格式的建立2.2.1 連續(xù)問(wèn)題離散化：網(wǎng)格剖分用差分

8、方法求解偏微分方程的第一步是將連續(xù)問(wèn)題離散化，所以取為自變量的增量，稱(chēng)為空間步長(zhǎng)；取為自變量的增量，稱(chēng)為時(shí)間步長(zhǎng)。接下來(lái)我們用兩組平行線(xiàn)構(gòu)成的長(zhǎng)方形網(wǎng)格覆蓋整個(gè)平面,網(wǎng)格點(diǎn)簡(jiǎn)記為.在上的結(jié)點(diǎn)稱(chēng)為邊界結(jié)點(diǎn)（結(jié)點(diǎn)：網(wǎng)格線(xiàn)的交點(diǎn)），其它的在內(nèi)的結(jié)點(diǎn)稱(chēng)為內(nèi)部結(jié)點(diǎn).下面介紹建立差分格式時(shí)幾個(gè)常見(jiàn)公式：，其中 （1），其中. （2），其中. （3），其中. （4）現(xiàn)就公式（1）（2）做如下證明：將在點(diǎn)處作拉格朗日型泰勒展開(kāi)得到 所以有，公式(2)得證.將在點(diǎn)處作拉格朗日型泰勒展開(kāi)得到 得到 （5）將在處作泰勒展開(kāi)得到 得到 （6）由(5)-(6)得：,公式(1)得證。稱(chēng)為關(guān)于自變量的二階中心差商，是關(guān)于自

9、變量的一階向前差商、是關(guān)于自變量的一階向后差商、是關(guān)于自變量的一階中心差商.上述公式過(guò)于冗長(zhǎng)，為方便起見(jiàn)，我們規(guī)定用下列符號(hào)來(lái)替代關(guān)于自變量的差商： 2.2.2 顯式差分格式的定義及截?cái)嗾`差差分方法的思想：(1)用差商代替微商；(2)求網(wǎng)格點(diǎn)上的近似值.利用關(guān)于的二階中心差商以及關(guān)于的向前差商，也就是公式（1）（2），在點(diǎn)處就有： （7）設(shè) （8） （9）則（7）式可以改寫(xiě)為，其中表示差分算子，表示在點(diǎn)以逼近的截?cái)嗾`差。從（9）中我們能得到，若在所考慮的區(qū)域內(nèi)是保持有界的，那么當(dāng)時(shí)有。用表示的近似值，忽略截?cái)嗾`差，則得到差分方程令，則有.顯式差分格式的圖示為： (k,j+1) j+1 （k-1

10、,j） (k,j) (k+1,j) j2.3 熱傳導(dǎo)方程混合問(wèn)題的差分方法2.3.1顯式差分方程問(wèn)題及向量表示與第一邊值問(wèn)題對(duì)應(yīng)的顯式差分問(wèn)題是：當(dāng)固定時(shí)，令，則方程組（10）可以表示成 （13）方程組（13）右端關(guān)于的系數(shù)是組成了一個(gè)N-1階的三對(duì)角矩陣 （14）再將方程組（10）（11）（12）也寫(xiě)成向量形式，令，則方程組（10）（11）（12）可以簡(jiǎn)單的記作 （15）所以，根據(jù)初始條件（11）以及邊界條件（12），再通過(guò)方程組（10）就可以逐層算出.例如已知層的值，則有，.初始條件和邊界條件已經(jīng)給出，由這兩個(gè)條件即可算出層的值，以此類(lèi)推。雖然顯式差分格式已經(jīng)構(gòu)造出來(lái)，但一個(gè)數(shù)值格式可以被

11、真正用于計(jì)算還需考慮以下兩點(diǎn)：1. 數(shù)值格式的穩(wěn)定性，即顯式差分格式的誤差在計(jì)算中能不能保持有界？2. 數(shù)值格式的收斂性，即當(dāng)無(wú)限小時(shí)，差分方程的精確解是否逼近微分方程解？2.3.2 顯式差分格式的穩(wěn)定性建立差分格式的時(shí)候，我們?cè)诰W(wǎng)格區(qū)域內(nèi)計(jì)算通常是分層進(jìn)行的。所謂的穩(wěn)定性問(wèn)題就是如果在某一層引入了誤差,那么一定會(huì)影響到下一層和之后各層的計(jì)算，但是我們?cè)诮獠罘址匠探M時(shí)必定會(huì)因?yàn)榭陀^因素的影響而引入誤差，比如初始誤差和舍入誤差,那么研究諸如這些誤差在計(jì)算過(guò)程中的影響。如果說(shuō)誤差是在剛開(kāi)始計(jì)算時(shí)引入的，要使差分方程穩(wěn)定，就要滿(mǎn)足在之后逐步計(jì)算過(guò)程中，誤差的影響有界或逐漸消失；不然就表明該差分格式是

12、不穩(wěn)定的。不穩(wěn)定的差分格式即使收斂也沒(méi)有意義,因?yàn)榇藭r(shí)誤差會(huì)隨步數(shù)的增多慢慢積累逐漸擴(kuò)大,最終差分方程的真實(shí)解將被完全淹沒(méi)。下面我們用圖來(lái)證明探究差分格式的穩(wěn)定性是十分有必要的。令時(shí)，差分格式為假定邊界條件的計(jì)算是精確的，只在初始層某個(gè)點(diǎn)上產(chǎn)生了誤差，初始層其他點(diǎn)以及以后各層上都沒(méi)有誤差。誤差應(yīng)該滿(mǎn)足方程組由此可以得到的如下分布表：000000000000000表1.時(shí)的誤差分布表從圖中能明顯看到，誤差是在逐層減少，所以當(dāng)步長(zhǎng)比時(shí)，顯式差分格式穩(wěn)定。再來(lái)看一下當(dāng)時(shí)的情形，誤差方程為此時(shí)誤差分布表如下：表2.時(shí)的誤差分布表從表中可以明顯看出誤差是在逐層增大的，所以當(dāng)時(shí)顯式差分格式不穩(wěn)定。因此研究

13、要使顯式差分格式穩(wěn)定，步長(zhǎng)比應(yīng)該在什么范圍內(nèi)就至關(guān)重要。下面開(kāi)始研究顯式差分格式穩(wěn)定的充要條件。設(shè)方程組（15）的近似解為，則滿(mǎn)足方程組因此誤差應(yīng)滿(mǎn)足方程由此就能得到，故有.由這個(gè)不等式可知，如果，則有，也就是在各層結(jié)點(diǎn)上的范數(shù)都不超過(guò)初始誤差的范數(shù)。在此時(shí)，顯式差分格式是穩(wěn)定的。所以問(wèn)題就轉(zhuǎn)化成了在什么條件下才能使成立。矩陣A也就是公式（14）是一個(gè)實(shí)對(duì)稱(chēng)矩陣，取范數(shù)，就有，其中是A的絕對(duì)值最大的特征值，所以說(shuō)顯式格式穩(wěn)定的充分條件是，,充分條件已證明，接下來(lái)用反證法說(shuō)明這一條件也是必要的。假設(shè)顯式格式穩(wěn)定，但是。如果是A的與對(duì)應(yīng)的特征向量，也就是，那就有，所以，。如果不是A的與對(duì)應(yīng)的特征向

14、量，那么就可以唯一的表示為A的N-1個(gè)線(xiàn)性無(wú)關(guān)的標(biāo)準(zhǔn)直交特征向量的線(xiàn)性組合，即，其中，這些特征向量是一個(gè)標(biāo)準(zhǔn)直交系，即所以，又有，因?yàn)锳是實(shí)對(duì)稱(chēng)矩陣，所以A的特征值都是實(shí)數(shù)，只要，當(dāng)時(shí)就總有。以上這兩種情況都說(shuō)明顯式差分格式是不穩(wěn)定的，與假設(shè)矛盾，所以。因此，顯式差分格式穩(wěn)定。定理1：顯式差分格式穩(wěn)定,.證明：求A的特征值.A的特征多項(xiàng)式為.按第一行展開(kāi)，就有,這個(gè)式子其實(shí)是的一個(gè)遞推公式 其中， .由此就可以依次求出令，那么特征方程為特征方程的根為，因此.設(shè)，那么， 再由， 能解出， 所以就有.接下來(lái)就要求方程的根，也就是求矩陣A的特征值.令那么，因此，那么要求的根，也就轉(zhuǎn)化成求的根。即，再

15、由，我們可以知道，這就是矩陣A的特征值。根據(jù)前面的討論我們已經(jīng)知道顯式差分格式穩(wěn)定的充要條件是，即，這個(gè)公式等價(jià)于 ，不等式右端很明顯對(duì)于一切都成立，要使左端不等式也成立，必須滿(mǎn)足。所以當(dāng)滿(mǎn)足時(shí)，顯式差分格式穩(wěn)定。同樣的再用反證法證明這一條件也是必要的。假設(shè)差分格式穩(wěn)定，但，這個(gè)時(shí)候可以把表示為，。對(duì)于充分大的N來(lái)說(shuō)，一定可以在中找到一個(gè)，使得.所以當(dāng)時(shí)，差分格式不穩(wěn)定。綜上，我們得到顯式差分格式穩(wěn)定的充要條件是，也就是步長(zhǎng)比一定要限制在這一范圍內(nèi)，所以通常也稱(chēng)顯式差分格式是條件穩(wěn)定的差分格式。2.3.3 顯式差分格式的收斂性接下來(lái)我們討論當(dāng)無(wú)限小時(shí)，差分方程精確解是否無(wú)限趨近于微分方程解.定

16、理2：若滿(mǎn)足第一邊值問(wèn)題的解在區(qū)域內(nèi)有連續(xù)的偏導(dǎo)數(shù)，并且，則顯式差分格式一致收斂.證明：由前面的討論我們知道微分方程的解應(yīng)滿(mǎn)足其中，=常數(shù).差分方程精確解要滿(mǎn)足因此，將微分方程解和差分方程精確解作差就能得到上式應(yīng)滿(mǎn)足其中.設(shè) ， .那我們就可以把方程組寫(xiě)成向量的形式所以令 ,.那么.當(dāng)時(shí)，.那我們就能得到 ，所以 ，（對(duì)所有k,j都成立）。因此，在我們考慮的區(qū)域內(nèi)的任意一點(diǎn)（k,j）上都有，當(dāng)時(shí).這就足以證明，要使顯式差分格式一致收斂，必須要滿(mǎn)足第一邊值問(wèn)題的解在區(qū)域內(nèi)有連續(xù)的偏導(dǎo)數(shù)，并且.3. 數(shù)值算例及圖像3.1 數(shù)值算例為了便于觀察和表示微分方程的解，我們將最終解用圖像表示出來(lái)。下面給出

17、了一個(gè)熱傳導(dǎo)方程問(wèn)題具體算例：真解為下面利用顯式格式，取不同的網(wǎng)格比和步長(zhǎng)h來(lái)計(jì)算，得到數(shù)值解與準(zhǔn)確解的三維圖像。 圖1. 時(shí)的數(shù)值解圖 圖2. 時(shí)的真解圖下面列出部分結(jié)點(diǎn)數(shù)值解、精確解和誤差絕對(duì)值的表格，便于我們更直觀的觀察。 數(shù)值解精確解 40（0.5,0.1）0.375690.372712.980e-380（0.5,0.2）0.141140.138912.230e-3120（0.5,0.3）0.0530260.0517731.253e-3160（0.5,0.4）0.0199210.0192966.250e-4200（0.5,0.5）0.00748430.00719192.924e-424

18、0（0.5,0.6）0.00281180.00268051.313e-4280（0.5,0.7）0.00105640.000999035.737e-5320（0.5,0.8）0.000396860.000372352.451e-5360（0.5,0.9）0.000149100.000138781.032e-5表2. 部分結(jié)點(diǎn)數(shù)值解、精確解和誤差的絕對(duì)值（） 取不同空間步長(zhǎng)時(shí)所得的數(shù)值解最大誤差我們記為固定網(wǎng)格比，分別取，利用loglog函數(shù)，探究與的關(guān)系，結(jié)果如下圖： 圖3. ，時(shí)的誤差圖從上圖中我們能很清楚的看到，真解和數(shù)值解的誤差以平方階收斂。4 總結(jié)本文介紹了求解熱傳導(dǎo)方程的顯式差分格式

19、，使用顯式差分法分析一維熱傳導(dǎo)問(wèn)題，把問(wèn)題的定義域劃分成網(wǎng)格，使用適當(dāng)?shù)臄?shù)值微分公式，就可以用差商代替原來(lái)熱傳導(dǎo)問(wèn)題的微商，從而把原問(wèn)題轉(zhuǎn)換成顯式差分格式，證明穩(wěn)定性及收斂性，討論其在不同網(wǎng)格上的誤差，并用實(shí)例和圖像簡(jiǎn)單的驗(yàn)證并說(shuō)明關(guān)于誤差的結(jié)論。熱傳導(dǎo)方程的求解方法多種多樣，比如Crank-Nicolson格式、隱式差分格式、顯式差分格式和Richardson格式，古典顯格式不用求解一個(gè)迭代方程組，計(jì)算速度快，需要的內(nèi)存也比隱式格式要少，使用起來(lái)方便；但是這種格式是條件穩(wěn)定的，也就是它要求計(jì)算的時(shí)間步長(zhǎng)要受到嚴(yán)格的控制，而隱式差分格式則不會(huì)受到任何限制，它是無(wú)條件穩(wěn)定，所以在實(shí)際應(yīng)用問(wèn)題當(dāng)中通常都會(huì)采用隱式差分格式而非顯式，這是古典顯格式的不足之處。5 參考文獻(xiàn)1 蘇煜城，吳啟光：偏微分方程數(shù)值解法 M南京大學(xué)數(shù)學(xué)系計(jì)算數(shù)學(xué)用書(shū) 科學(xué)出版社 2 李榮華：偏微分方程數(shù)值解法 M高等教育出版社 3 馮立偉：熱傳導(dǎo)方程幾種差分格式的MATLAB數(shù)值解法比較 M沈陽(yáng)化工大學(xué)，遼寧沈陽(yáng).2011-64 梁昆淼：數(shù)學(xué)物理方法