有限差分法

有限差分法簡單的說便是在離散的向度上，將控制方程式及邊界狀態 中的微分項以差分式替代

(Boore, 1972)

，使得原先的微分方程式變成簡單的 四則運算式

(Mocoz, Kristek and Halada, 2004)

。當我們要將有限差分法應用 至特定的控制方程式或微分方程式中時，總括來看有三個主要步驟：先

(1)

建構離散的有限差分離散模型，亦即建構離散的網格，並對微分方程式導 入有限差分式。然後對導出的有限差分運算式進行

(2)

有限差分離散模型的 分析，確定離散模型的收斂性、穩定性及餘尾誤差級數。最後使用計算系 統協助

(3)

疊代計算

(iteration)

。

網格網格

網格網格系統系統系統系統(Grid)：使用有限差分法的第一步便是要先決定所採用的離散網格 系統。所謂離散網格系統的意思是說，在微分方程式中，其所定義的時間 及空間皆是連續，而為了將有限差分是導入，必須將連續的空間以及時間 向度離散，而以離散的點代替原本連續的物理量。這樣離散後的時間及空 間向度，便稱做為離散網格系統，或簡稱網格。

在 傳 統 的 有 限 差 分 法 使 用 上 ， 皆 是 採 用 均 勻 方 型 網 格

(Uniform

Rectangular Grid, Alterman and Karal, 1968

；

Boore, 1972

；

Kelly et al., 1976)

， 亦即其將空間離散為均勻分布的格點

(

如圖

2.22a

所示

)

。而後為因應不同的 使 用 需 求 而 有 均 勻 三 角 形 網 格

(Uniform Trigular Grid)

、 非 均 勻 網 格

(non-uniform grid)

等網格系統的使用

(

如圖

2.22b,c,d

所示

)(Moczo, 1989

；

Zhang, 1997

；

Pitarka, 1999)

。

圖圖

圖圖 2222.22.22.22.22 有限差分法有限差分法有限差分法網格系統示意圖有限差分法網格系統示意圖網格系統示意圖網格系統示意圖

以均勻方型網格舉例，若有一個具有三變數的物理量其變數為

(x,y,t)

，

x

表示空間上水平向，

y

表示空間上垂直向，

t

表示時間，則在這三個量度上，

其離散後的變數

x

ⁱ、

y

^j、

t

k定義為：

x i x

x

_i

=

₀

+ ∆

，

for i = 0, ± 1, ± 2, … y

j y

y

_j

=

₀

+ ∆ , for j = 0, ± 1, ± 2, … t

k t

t

_k

=

₀

+ ∆ , for k = 0, 1, 2, …

(2. 105)

其中，∆x,∆y,∆t分別為在

x

方向、

y

方向以及時間軸上，離散後相鄰兩點的 間距；

x

₀

, y

₀

, t

₀為網格中各向度上任一點的實際值。假設格點上的物理量為

U(x

i

,y

j

,t

k

)

，則可簡單表示為

U

_i^k_,j。

標準有限差分式 標準有限差分式 標準有限差分式

標準有限差分式(Standard Finite Difference Approximation)：當離散網格系 統選定後，便可開始選用適當的有限差分式導入微分方程式中。習慣上，

對於一開始發展有限差方法時應用在均勻方型網格系統的有限差分式，我 們便稱其為標準有限差分式。假設有一函數

F(x)

存在有連續的一次微分式，

則此微分式可以採用有限差分式

((2. 106)~((2. 108)

式進行導入：

前向差分式 ^{'( ) ( ) ( ) O(h)}

後向差分式 ^{( )}

difference approximation)

。對於更高級餘尾誤差的一次微分差分式或是更高

階的微分方程的有限差分式，都可使用相同的方式獲得，在這不詳加介紹，

可逕行參考

Mitchell and Griffiths (1980)

、

Morton and Mayers (1994)

或

Durran

(1999)

等人的著作，有詳細的推導及應用。

一致性一致性

一致性一致性、、、、穩定性與收斂性分析穩定性與收斂性分析穩定性與收斂性分析穩定性與收斂性分析：在採用有限差分法的時候應要特別注意幾 點有限差分運算式的特性。首先，對於有限差分式與其所欲計算的微分式

的一致性

(consistent)

。一致性的意思是說，當離散的間距趨近於零的時候，

有限差分式將會與微分式幾乎相等。大致而言，這一點是可確認的。經由 泰勒展開式得到的差分式餘尾誤差級數可很容易的證明這點，但在某些情 況下，必須要使空間間距與時間間距符合某種關係時才可成立。若一差分 式具有這樣的特性，我們稱此差分式是有條件的具有一致性

(Conditionally

Consistent)

。惟有當一致性能成立時，才可使用此有限差分式。

再者需注意的是導出的差分運算式

(Finite Differece Equation)

所得到的 解是否具有穩定性。當微分方程的解具有結界

(bounded

)時，差分運算式所 得到的解亦具有結界，則稱此差分運算式是穩定的；但若差分運算是所得 到的解並不具有結界，則稱其是不穩定。大致而言，我們所探求的物理現 象都是具有結界的，這點可由微分方程的解析解得知，因此我們期待使用 差分運算式所獲得的解亦是具有結界。若任意的空間與時間離散間距都可 使 差 分 運 算 式 滿 足 穩 定 性 的 要 求 ， 則 稱 此 運 算 方 程 式 是 無 條 件 穩 定

(unconditionally stable)

；反之，若需使空間與時間的離散間距滿足某種關係

時才可達到穩定性的要求，便稱其是有條件式的穩定

(conditionally stable)

。 若無法達到穩定性的要求，則此差分運算式所獲得的解將無法表現實際解 的行為。穩定性的分析只能對具有線性行為的微分方程進行，若針對對非 線性問題，則通常需要先使問題成為區域性線性

(linerarized locally)

。最常 被使用於穩定性分析的方法是

von Neumann Method

，對於此方法的說明可 參考

Moczo, Kristek, and Hakada (2004)

。

最 後 一 個 需 要 對 採 用 的 差 分 運 算 式 進 行 解 析 的 是 其 收 斂 性

(Convergent)

。收斂性所表示的是，當離散的間距趨近於零的時候，微分方

程的解與差分運算式所得到的解會幾乎相等。在這需要做一個釐清，一致 性與收斂性所關心的是不同的兩個問題：一致性所關心的是差分式本身與 微分這動作的關係，而收斂性所關心的是由差分運算式所得到的解。對於 收斂性的分析是困難的，然而，幸運地，收斂性與一致性以及穩定性有絕 對的關係。根據

Lax equivalent theorem

，在數值架構中，當一致性與穩定性 同時滿足時，其必然亦滿足收斂性

(Moczo, Kristek and Halada, 2004)

。 因此，在採用有限差分式時應當要注意一致性與穩定性的滿足條件。

除此之外，採用不同的有限差分式以及網格系統，會因為不同的問題而有

其需考慮的特定問題，若僅是盲目的使用任一有限差分法於不同的問題 上，而不對其可能造成的誤差進行了解，過度信任所得到的結果，很容易 便有錯誤的使用發生。以震波模擬問題而言，在採用有限差分法時，除了 要注意滿足穩定性的條件，還應當注意因為不同網格系統而造成不同程度 的假頻散現象

(Grid Dispersion

，即在地質材料本身性質外，因為空間離散而 造成波速隨頻率變化的現象

)

，亦應該要注意其整個有限差分系統

(

包含控制 方程式，網格系統，以及使用的有限差分式

)

所能應用的地質材料狀況

(Moczo, Kristek and Halada, 2004)

。

有限差分法震波模擬 有限差分法震波模擬 有限差分法震波模擬

有限差分法震波模擬：在震波數值模擬領域中，有限差分法因為在應用上 及程式撰寫上的便利性成為最為廣為使用的技術。其在震測資料模擬的研 究上，最先是由學者

Z. S. Alterman

在

1968

年代導入使用。她與她的研究團 隊利用標準的有限差分法，施以不同的起始值與邊界值進行一系列震波在 均質均向的半無限域土層中傳遞的模擬（

Alterman and Kornfeld, 1968

；

Alterman and Karal, 1968

；

Alterman and Rotenberg, 1969

；

Alterman and Loewenthal, 1970

；

Ottaviani, 1971

；

Ilan et al., 1975

等）。然而，如前所述，

我們常須了解的是具有不平整介面且有側向變化的土層，因此

Alterman

研 究團隊的數值計算模式發展到

1972

年後，便有針對在非均質土壤中的二階 微分方程的有限差分計算模式出現（

Boore, 1972

），並加以改進擴充（

Kelly

et al., 1976

），使得運算式能自動符合模型內部的地層材料變化，而無須設

定內部的邊界狀態。

然而，在標準有限差分法的架構下所進行的震波模擬，容易因為高波

松比

(Poison’s Ratio)

地層材料的存在而發生即使已滿足穩定性條件，亦會有

不穩定以及嚴重的假頻散現象產生的問題。為解決這問題，目前在有限差 分 法 網 格 系 統 中 受 到 大 量 採 用 的 錯 置 網 格

(Staggered Grid)

系 統 是 由

Madariaga

在

1976

年提出，應用在地震時地層裂隙的動態模擬上。錯置網 格與均勻方型網格的主要差異是在於，在錯置網格中，各物理量的定義位 置不在同一格點上，且相距離散間距的一半

(

如圖

2.23, 2.24

所示

)

。此一格

網系統由

Virieux(1984,1986)

引進至震波模擬上使用，採用二級有限差分式

導入控制方程式。其在控制方程的處理上，與過去以位移做為導入物理量 處理二階微分方程的方式不同，以速度以及應力為其導入物理量，成功的 將控制方程降為一階微分式，因而增加了其導入有限差分式的簡單性。因 其採用速度及應力作為導入物理量，故常以速度應力法

(Velocity-Stress Method)

稱之

(

詳細內容請見

3.1.1

節

)

。經由

Levander (1988)

將四級有限差分 式導入速度應力法，其發展已大致完成。由於其計算上的效率以及準確性，

雖仍有許多新概念的方法出現，但依舊是目前有限差分法在震波模擬上的 主流架構

(Olsen and Schuster, 1992

；

Graves, 1993

；

Rodrigues, 1993

；

Yomogida and Etgen, 1993

；

Graves, 1996

；

Robertsson and Holliger, 1997

；

Hayashi and Suzuki, 2004)

。

圖圖圖圖 2222.23.23.23 .23 速度錯速度錯速度錯速度錯置網格架構置網格架構置網格架構，置網格架構，，時間定，時間定時間定義在時間定義在義在 k+1/2義在k+1/2k+1/2k+1/2

圖圖圖圖 2.242.242.24 2.24

應力錯置網格架構應力錯置網格架構應力錯置網格架構應力錯置網格架構，，，，時間定義在時間定義在時間定義在時間定義在 k+1k+1k+1 k+1

由於速度應力法所得到的是速度值而非位移值，須經由對速度進行時 間積分才可得到位移。然而在這存在著一問題，積分的過程中會有一未知 的常數存在，而使得所得位移會與實際位移有一常數差距，造成震波的位

移歷時

(time history)

有平移現象的可能；又積分過程需要時間，因而降低了

些 許 效 率 。 於 是 ， 在

1990

年

Luo and Schuster

提 出 位 移 應 力 法

(Displacement-Velocity Method)

，依舊採用錯置網格，但以位移及應力為其 導入物理量。在這雖需處理位移的二階微分，但其不須經由積分便可得到 位移資料，進而使其效率增加。

近來對於地層的震波模擬除討論純彈性的非均質土層外，還對具有黏 滯性及非均向的土層感到興趣，更複雜的定義關係使得錯置網格不敷使 用，因此

Magnier et al.(1994)

引進了

Andrews(1973)

提出的部份錯置網格

（

partly-stggered grid

）系統。部份錯置網格與錯置網格的不同是在於，錯

置網格的物理量是定義在不同的位置上，而部分錯置網格則是有些共用格 點而有些不共用。由於其能使得錯置網格無法應用的非均向地層模擬問題

得到解決，而引起不少的注意與應用

(Zhang, 1997

；

Saenger et al., 2000

；

得到解決，而引起不少的注意與應用

(Zhang, 1997

；

Saenger et al., 2000

；

在文檔中 Pseudo-section概念於表面波震測應用之數值模擬探討 (頁 63-71)