第四章 地下水流、熱流及污染物傳輸模式建置
4.2 數值離散說明
在說明各應用模組如進行數值離散之前，先說明本研究數值離散之 方法。以第三章之圖3.7(a)為例，該圖為應用凡諾依圖所建立之模擬 格網，這是一組不規則網格，圖上環繞計算結點i之相鄰結點分別為 A、B 與 C，而圖上虛線為中央結點與相鄰結點連結iA、iB與iC，實 線則為虛線之中垂線，中垂線所形成之區域即為計算結點i之控制體 積，以下離散方法均以此示意圖進行說明。
1. 微分運算子：方程式中如有對空間之微分運算子，在此以最簡單 的差分方法進行離散。式(4-29)是針對iB方向之空間偏微分，其 為節點B之數值( fB)減去節點i之數值( fi)，再除以兩點距離(iB)。
114
115
116
(4-4)、(4-5)、 (4-6)、(4-9)、(4-10)與(4-11)本身為顯示形式之方程式 (explicit form)，不存有微分或積分運算子，亦即將自變數代入即可直
117
接計算得出其應變數，因此在計算上無需藉由數值技巧特殊處理即可 計算。式(4-1)為連續方程式，為水質量守恆之雷諾傳輸定律，依據式 (3-19)可將其離散，其離散方程式如下所示：
M ∆, − M , + ω −q ∆, + M ,∆
∈
× ∆t
+ (1 − ω) × −q , + M ,
∈
× ∆t = 0
(4-37)
M , = n ρ ,S , V , (4-38)
M , = m , A (4-39)
另外，式(4-7)為計算穿越水質量流率之計算方程式，搭配達西定律 進行計算，其中包含微分運算子，因此以簡單差分(simple difference) 來離散微分運算子(如式(4-29)所示)，如式(4-40)所示。ij 代表格點 i 與 格點 j之直線距離，正值代表水流由格點 i 流向格點 j 。此外，K , 代表線段ij之等效水力傳導係數，其必須利用K , 與K ,來估算，估 算方式包含算數平均(如式(3-20)所示)、幾何平均(如式(3-21)所示)、
調和平均(如式(3-22)所示)與上風法(如式(3-23)所示)四種，可適性計 算架構可在背景直接擇一進行推估。傳統上水力傳導係數 (K , ) 建
118
議以調和平均進行計算。另外，溶液密度 (ρ , ) 則可應用算數平 均進行推估。
m , = −ρ , K , h − h
ij (4-40)
式(4-8)內具有積分運算子，積分運算子可應用梯形積分法進行離散 (如式(4-30)所示)，其將定積分區間切割成 N 段，故各段長度應符合 z = z + NΔz ，梯形積分法將各段端點之函數值乘上梯形積分權重 (M , )再乘上各段長度(Δz)，如 i 等於 0 或 N 時，其權重值為 0.5 ， 其他權重值則為1 。因此利用前述辛普森積分法，可以將式(4-8)離 散成式(4-41)。其中 N 代表節點集合，涵蓋 z 至 z 之間的所有節 點。
h = p + ζ = p + ρ gdz
ρ° g ≅ p + ∑ ∈ C ρ , g × ∆z
ρ° g (4-41)
4.2.2 熱流傳輸應用模組
H ∆ − H + ω −q ,∆ + H ∆
∈
× ∆t + (1 − ω)
× −q , + H
∈
× ∆t = 0
(4-42)
119
H = n S , ρ , ∅ + (1 − n°)ρ , ∅ , V , (4-43)
H = h , + h , A (4-44)
h , = −K , , T − T ij
(4-45a)
h , = −K , , T − T
ij − n D , T − T
ij − n D , ∂T
∂n
(4-45b)
在熱流應用模組中，式(4-16)為藉由雷諾傳輸定律所建立的守恆方程 式，其中具有積分運算子。(4-17)為擴散項之熱傳輸計算方程式，其 中具有微分運算子，因此這兩個方程式必須要進行數值離散處理。式 (4-42)至(4-44)為守恆方程式的離散式，其處理流程同式(4-37)一般。
式(4-45a)為擴散項之離散方程式。ij代表格點 i 與格點 j 之直線距離，
若h , 為正值代表熱流由格點 i 流向格點 j 。式(4-45b)為延散項之 離散方程式。在第一項中之說明已於式(4-45a)說明中敘述，在第二項 中，若為該項為正值則代表熱流由格點 i 流向格點 j 。在第三項中，
對於橫向溫度梯度之處理，請詳見附錄E。
h , = m , ϕ, = m , ϕ, for m , > 0
m , ϕ, for m , < 0 (4-46) 式(4-46)為對流項之離散方程式，其中 m , 為水流質量流率，可由 地下水流應用模組估算； ϕ, 為水流之單位質量熱容量，傳統數值
120
離散上，由於熱能對流是透過水流帶動，因此在線段 ij 之單位質量 熱容量之推估，建議以上風法進行推估，意即直接取水流流動之上風 方位之單位質量熱容量數值，因此式(4-46)以 m , 之正負號進行上 風方向之判斷。
4.2.3 污染傳輸應用模組
M ,∆ − M , = ω −q ,∆ + M ,∆
∈
× ∆t
+(1 − ω) × −q , + M ,
∈
× ∆t
(4-47)
M , = n S , ρ C, V , (4-48)
M , = m , , + m , , A (4-49)
在溶質傳輸應用模組中，式(4-22)為藉由雷諾傳輸定律所建立的守恆 方程式，其中具有積分運算子。式(4-23)為延散項之溶質傳輸計算方 程式，其中具有微分運算子，因此這兩個方程式必須要進行數值離散 處理。式(4-47)至(4-49)為守恆方程式的離散式，其處理流程同式(4-37) 一般。式(4-50)為延散項之離散方程式。在第一項中， ij 代表格點 i 與格點 j 之直線距離，第一項為正值代表污染物由格點 i 流向格點
j 。在第二項中，對於橫向濃度梯度之處理，請詳見附錄 E。
121
m , , = −n D , C, − C,
ij − n D , ∂C,
∂n
(4-50)
式(4-51)為對流項之離散方程式，其中 m , 為水流質量流率，可由
地下水流應用模組估算；傳統數值離散上，由於溶質對流是透過水流 帶動，因此在線段 ij 之濃度(C , )與密度(ρ , )推估上，建議以上風 法進行推估，意即直接取水流流動之上風方位之濃度數值，因此式 (4-51)以 m 之正負號進行上風方向之判斷。
m , , = m , C, ρ , for m , > 0
m , C, ρ , for m , < 0 (4-51)
122