二、幾種主要的差分格式

（一）顯式差分格式

為了便于說(shuō)明，以下列一維河間地塊均質(zhì)各向同性承壓含水層中的地下水流問(wèn)題

為例來(lái)加以說(shuō)明。

首先將研究區(qū)域[0，L]用直線等分為l份，步長(zhǎng)Δx=L/l，把時(shí)間段[0，Tsum]用直線等分成m份，時(shí)間步長(zhǎng)Δt=Tsum/m，構(gòu)成如圖2-1所示的網(wǎng)格，結(jié)點(diǎn)坐標(biāo)xi=iΔx，tκ=κΔt（i=0，1，…，l；κ=0，1，…，m），簡(jiǎn)記為（i，κ），并以表示H（iΔx，κΔt），以表示原方程的差分方程解（即H的近似值）。

式（2-5）中的導(dǎo)數(shù)，用差商代替，在典型結(jié)點(diǎn)（i，κ）處表示為

略去O（Δt）和O（[Δx]2），可得和式（2-5）以及圖2-1中x，t平面上的網(wǎng)格對(duì)應(yīng)的差分方程為

截?cái)嗾`差為O（Δt）+O（[Δx]2）。即當(dāng)Δt和Δx很小時(shí)，這一誤差是Δt階的一個(gè)量與（Δx）2階的一個(gè)量之和。若定義

則式（2-10）可改寫(xiě)為

由此可知，只要知道某一時(shí)段κ開(kāi)始時(shí)刻tκ各結(jié)點(diǎn)的值，利用上式便能算出tκ+1時(shí)刻，即κ時(shí)段終了時(shí)刻的值（1≤i≤l-1，1≤κ≤m）。所以稱(chēng)這一方法為顯式方法。邊界結(jié)點(diǎn)的水頭值則由邊界條件給出，即

這樣利用t=0時(shí)刻時(shí)各結(jié)點(diǎn)的值已由初始條件式（2-6）

給出的情況，直接計(jì)算t1時(shí)刻各內(nèi)部結(jié)點(diǎn)的水頭值，并利用邊界條件補(bǔ)充邊界結(jié)點(diǎn)上的水頭值。再把求得的t1時(shí)刻各結(jié)點(diǎn)的水頭作為初值，重復(fù)上述過(guò)程求t2時(shí)刻各結(jié)點(diǎn)的水頭值。如此一個(gè)時(shí)間水平、一個(gè)時(shí)間水平地做下去，就能求得計(jì)算區(qū)Ω上所有時(shí)刻的水頭值。

前面我們給出了求的方法，但必須回答一個(gè)問(wèn)題，即差分方程的解是不是很逼近原微分方程的解在相應(yīng)結(jié)點(diǎn)上的值？為此，需要從兩個(gè)方面，即差分方程的收斂性和穩(wěn)定性來(lái)回答上述問(wèn)題。

如果差分方程的解在步長(zhǎng)Δx、Δt取得充分小時(shí)，和微分方程的解析解在某種意義上很接近的話(huà)，便說(shuō)這種差分格式是收斂的。研究收斂性就是討論當(dāng)Δx→0、Δt→0時(shí)，差分方程的解和微分方程解的差（在一維條件下為）的絕對(duì)值在什么條件下趨近于零。其次，實(shí)際計(jì)算中由于只能用有限位計(jì)算，每一步都會(huì)有舍入誤差，而且它還影響以后的計(jì)算結(jié)果。于是要考慮一個(gè)問(wèn)題，當(dāng)某一步結(jié)果本身有誤差時(shí)，利用它去計(jì)算，若Δx和Δt固定，隨著計(jì)算時(shí)間或計(jì)算次數(shù)的增加，誤差是逐漸消除？還是逐步積累，愈變愈大？如是后者，則當(dāng)t→∞時(shí)（或計(jì)算次數(shù)無(wú)限增多時(shí)），盡管某一步的誤差很小，但其影響最終有可能達(dá)到十分可觀的程度，使所得解面目全非。這時(shí)所考慮的差分格式便是不穩(wěn)定的。顯然，不收斂和不穩(wěn)定的差分格式是沒(méi)有實(shí)用價(jià)值的。

顯式格式經(jīng)證明，只有當(dāng)0<λ≤1/2時(shí)才收斂和穩(wěn)定。因此，Δt不能取的太大。

（二）隱式差分格式

如式（2-5）左端二階導(dǎo)數(shù)取在tκ+1時(shí)間水平上[即用κ時(shí)段末t=（κ+1）Δt時(shí)刻的水頭值]，便得隱式差分方程為

截?cái)嗾`差為O（Δt）+O（[Δx]2）。仍用式（2-11）定義的λ，則式（2-12）化為

式（2-13）左端包含3個(gè)未知數(shù)，不能直接解出，所以稱(chēng)為隱式方法。必須對(duì)所有內(nèi)結(jié)點(diǎn)（本例共有l(wèi)-1個(gè)）都列出與式（2-13）相應(yīng)的方程，并把邊界條件

代入第一個(gè)和最后一個(gè)方程，形成由l-1個(gè)方程組成的方程組。聯(lián)立求解，可得l-1個(gè)內(nèi)結(jié)點(diǎn)上tκ+1時(shí)刻的水頭值。所得代數(shù)方程組的系數(shù)矩陣只在三條對(duì)角線上有值，其余均為零，故稱(chēng)三對(duì)角線方陣。可用追趕法求解。

可以證明，隱式方法對(duì)任何λ都是收斂、穩(wěn)定的，也就是說(shuō)它的收斂、穩(wěn)定是無(wú)條件的，Δt的取值不受Δx的嚴(yán)格限制。

（三）中心式差分格式

如果對(duì)在tκ和tκ+1的中點(diǎn)t=tκ+Δt/2處取中心差

把沿t的正向在tκ+1/2處用Talyor級(jí)數(shù)展開(kāi)，沿t的負(fù)向在tκ+1/2處用Talyor級(jí)數(shù)展開(kāi)，各取兩項(xiàng)，兩式相加，得

于是，式（2-5）可寫(xiě)成下列差分格式

截?cái)嗾`差為O（[Δt]2）+O（[Δx]2）稱(chēng)為Crank-Nicolson中心式差分格式或六點(diǎn)格式，形成的代數(shù)方程組的系數(shù)矩陣也是三對(duì)角線方陣。它和隱式方法一樣也是無(wú)條件收斂、穩(wěn)定的。

（四）加權(quán)顯式-隱式格式

前面介紹的幾種差分格式可以統(tǒng)一到下列一般形式中

θ為權(quán)因子。上述格式稱(chēng)為加權(quán)顯式-隱式格式。當(dāng)θ=0時(shí)即為顯式方法；θ=1時(shí)即為隱式方法；θ=1/2時(shí)即為中心式差分方法。不難證明，如取則式（2-15）是無(wú)條件穩(wěn)定的，截?cái)嗾`差為O[（Δt）2+（Δx）4]（Lapidus和Pinder，1982）。