3、 地層下陷模式之研選與測試
3.1 地層下陷模式研選
地下水超抽導致地層下陷主要的原因乃是由於地下水之抽取導致孔隙水壓改變，破壞 地層原有的平衡狀態而發生地層下陷。對於此種因超抽地下水而引致地層下陷之研究最早 始於 Fuller(1908)，至今雖然已有許多研究方法提出，但大致可將其歸納分為經驗法 (empirical method)、半經驗法(semi-empirical method)、半理論法(semi-theoretical method)以 及理論法(theoretical method)等四類。
所謂經驗法即是直接推求地層下陷與時間之迴歸關係，並利用外插方式預測未來可能 之下陷量，如曹以松(民國 58 年)與 Vega(1984)。一般常用的迴歸關係式有二次函數、指數 函數及對數函數等。此法雖然簡單，但是由於不考慮土體與流體之物理特性，且通常需要 較大量之現場實測資料以求得迴歸關係，故較為不經濟且不準確。
半經驗法則是考慮地層下陷量及其相關因子之迴歸函數關係，如地下水位、水位洩降 量或抽水量等。Wadachi(1940)曾探討地層下陷量與地下水位之相關性，Gabrysch(1969)則 提出地層下陷量與水頭降低量之比值與含水量黏土含量之百分比成一簡單正比關係。此 外，廖日昇等(民國 79 年)與柳志錫等(民國 81 年)分析北港沿海地區與彰雲地區，結果發現 該地區平均地下水位下降 1 公尺，會發生 0.1 至 0.8 公分之地層下陷量。一般而言，雖然半 經驗法比經驗法更為合理，但是仍缺乏嚴謹之理論架構，故往往無法說明各物理量間之定 率關係，因此也較難達到準確之要求。
半理論法則是以應力應變等概念模式，配合現場之實測值，求得模式所需參數值後，
即可應用於地層下陷量之推估，如深度孔隙率模式(Schatz(1978))，對數孔隙比有效應力模 式(Gibbs(1950))以及簡單估計模式(盛若磐(民國 78 年))等。半理論法通常應用於小範圍之 地層下陷，較不適用於複雜土層之地層下陷預測。
由上述可知，不論經驗法、半經驗法以及半理論法，一般而言，因未能完整地考慮地 下水超抽導致地層下陷之發生機制，所以預測結果可能較不精確，故應採用較完整之理論 法。而地下水超抽導致地層下陷研究之理論法，大致上可分為一維非耦合與三維耦合計算 兩種。
嚴謹的地層下陷理論始於土壤力學大師 Terzaghi(1925)所提出之有效應力觀念。
Terzaghi 於飽和土壤中，在土體之總應力不變且變形為微小垂向一維，以及流體為不可壓 縮且在土壤中流動滿足達西定律(Darcy’s Law)等假設下，導出著名的 Terzaghi 單向度一維 垂直壓密方程式。因為 Terzaghi 僅考慮流體為一維流動，所以 Jacob(1940,1950)在假設流體 具壓縮性下，將傳統二維或三維之地下水流方程式(Theis(1938))與 Terzaghi 一維壓密理論 結合，推得 Terzaghi-Jacob 單向度壓密二維、三維地下水流方程式，並清楚地交代地下水 流方程式中貯蓄係數(storage coefficient)之物理意義。
上述兩種方法均將孔隙水壓與土體變形分開求解，即所謂兩步驟法(two-step method) 或非耦合法(uncouple method)。也就是，若採用一維非耦合地層下陷模式，僅需先求解地 下水流方程式，再將所得之孔隙水壓變化代入土體位移代數式，即可求得下陷量。一維非 耦合地層下陷模式已被廣泛的使用於地下超抽導致地層下陷之計算，尤其是區域性
(regional)問題，如 Helm(1975,1976)、Gambolati and Freeze(1973)，Neuman et al.(1982)、
Narasimhan and Witherspoon (1977)、Bear and Corapcioglu(1981a)、Mishra et al. (1993)、
Gambolati et al. (1991)、李德河等（民國 86 年）以及歐陽湘等（民國 86 年）。因為地層下 陷量只是孔隙水壓變化在已知土體壓縮係數下之簡單代數計算結果而已，所以對於一維非 耦合地層下陷計算，首要的工作乃是如何準確且有效地模擬地下水之孔隙水壓變化。
由於大自然之沉積作用，使得區域性地下水流域通常由一連串導水性質良好之含水層 夾雜導水性質差之阻水層所形成，此種層狀系統稱為多含水層系統(multi-aquifers system)。
一般而言，多含水層系統之地下水為三維流動狀況，且現今已有許多三維地下水流模式 (three-dimensional groundwater model)可供模擬使用，如 FE3DGW(Gupta(1984))、FLOW3D、
MODFLOW 與 3DFEMWATER 等。但是由於區域性地下水流域中，水平方向之延展距離 遠大於垂直方向之厚度，且相鄰阻水層與含水層間之導水性質差異甚大，以及考量使用直 接三維模擬計算量太大等原因，乃有擬三維水流概念(quasi three-dimensional concept)之提 出，也就是將阻水層與含水層中之水流分別假設為垂直流與水平流，而阻水層與含水層再 以兩者間之滲流量加以耦合（Hantush(1960)）。如此三維水流控制方程式，於含水層中被簡 化為水平二維，而阻水層中則被簡化為垂向一維，將大大的減少計算時間。此外，對於阻 水層與含水層間之滲流量推估法共有兩種，第一種方法為藉由求解阻水層一維垂直水流方 程式來估計滲流量，即考慮所謂阻水層之貯蓄效應(aquitard storage effect) (見
Hantush(1960))，另一種方法則為假設滲流量與阻水層相鄰上、下含水層之水頭差(hydraulic head)成正比，也就是不考慮阻水層貯蓄效應(Polubarinov-Kochina(1962))。
利用擬三維水流概念模擬多含水層系統之地下水流相當多，如 Berdehoft and
Pinder(1970)、Chorley and Frind(1978)、Neumann et al.(1982)與已廣泛被使用之地下水模擬 軟體 MODFLOW。Neuman and Witherspon(1969)曾以雙含水層單井抽水之案例，對擬三維 與三維水流概念之差異加以探討，並指出當相鄰阻水層與含水層間之滲透率(permeability) 相差兩個量階以上，其誤差小於 5﹪。雖然擬三維水流概念對區域性之地下水模擬有不錯 之結果，但是對於如抽水井附近以及土層不連續之局部三維水流流況則無法模擬，而需要 採用三維水流模式。由上述可知，要如何才能建立一套高效率且準確，並可同時適用於區 域與局部地下水模擬之地下水流計算方法，仍值得進一步探討。
雖然一維非耦合地層下陷模式已被廣泛應用於地下水超抽導致區域性地層下陷之計 算，但是由於土體之變形僅假設為垂向一維而不考慮水平位移之影響，仍有不足之處，故 有三維耦合地層下陷理論之提出。Rendulic(1935)、Terzaghi(1943)以及 Barron(1948)等企圖 發展三維 Terzaghi 壓密理論，但因缺乏條件說明應力應變關係，以滿足應力應變之相容性 與唯一解之特性而失敗，此問題直到 Biot(1941,1955)提出孔隙流與彈性體應力應變關係才 獲得解決。Biot 假設土體變形符合彈性理論，且在土體飽和與流體微可壓縮以及達西相對 流速定律適用情況下，將三維應力應變關係結合三維地下水流方程式，建立著名的 Biot 三 維孔隙彈性理論(three-dimensional poroelastic theory)。Biot 孔隙彈性理論之理論特色在於可 以同時求解孔隙水壓與三維土體變形量(水平與垂直變形)，即所謂耦合法(couple method)。
此外，Verruijt(1969)繼續發展 Biot 之模式，在引用多相流之觀念下，分別考慮固相土體與 液相水以及汽相空氣等之質量守恆關係，以及固相土體之力平衡，有系統地建構出三維水 流及三維彈性壓密方程式。Verruijt 模式之特點為可適用於非飽和土體，且在假設土體總應 力不變及一維垂向飽和壓密下，可簡化為 Terzaghi-Jacob 單向度壓密三維水流方程式。利 用 Biot 三維孔隙彈性理論來分析地層下陷者，如 Safai and Pinder(1979,1980)、Lewis and
Schrefler(1978,1991)、Yeh et al.(1995)與 Tarn and Lu(1991)等。至於一維非耦合模式與三維 耦合模式之差異，Lewis and Schrefler(1991)有詳細討論。此外，Ghaboussi and Wilson(1973) 指出，對於三維耦合地層下陷之計算，因為土體與流體之量階與參數差異甚大，所以時間 間距之選用需特別注意，才不至於產生計算不穩定現象。
三維耦合模式在理論上雖較完整，但是實際應用時需要大量計算時間。針對此問題，
Bear and Corapcioglu(1981b)與 Corapcioglu and Bear(1983)提出深度平均二維區域模式 (depth-averaged two-dimensional regional model)概念，即是先將土體依導水性質之不同加以 分層，且在平面應力作用及各層之孔隙水壓與水平位移變形量於垂直深度上均相同的假設 下，將三維耦合地層下陷方程式沿垂項深度作積分。Bear 所提出之二維區域模式，雖然能 在減少計算量下有效地耦合求解孔隙水壓及水平與垂直方向之變形，但是由於不考慮阻水 層之貯蓄現象，且孔隙水壓與水平位移變形的過度簡化假設，可能無法精確地估計地層下 陷量。但是 Bear and Corapcioglu 應用垂向積分技巧將三維問題簡化為二維形式以減少計算 量之概念，卻值得加以進一步探討。
u 為超額孔隙水壓， e C 為壓密係數(coefficient of consolidation)， v k
為滲透係數，M 為 v
體積壓縮係數，γ w
為流體之密度重， HΔ 為下陷量， H 為土層厚度，C 為壓縮指數 c
(Compression index)，C r
為再壓指數(recompression index)，e 為初始孔隙比(void ratio)， 0 σ 0
為初始土體有效應力，σ 為土體之預壓密應力，c σ f
為土體有效應力。模擬時所需輸入之 隙比、壓密係數、壓縮指數(compression index,C )、再壓指數(recompression index， c C r
) 等。所輸入之土壤參數，壓縮比介於 0.1 至 0.2 間、壓密係數介於 1.638 至 93.05(m2/month) 間、壓縮指數介於 0.1855 至 0.3249 間、再壓指數介於 0.0107 至 0.1176 間、初始孔隙比介 於 0.604 至 1.199 等。4.工研院能資所研究團隊
依據 Terzaghi 一維壓密理論，利用 COMPAC 模式(Helm, (1984))模擬雲林縣王功、台 西與宜梧國中等地區之地層下陷。模擬時所需輸入之資料有土層分布情況、地下水位變化 值、水力傳導係數、非彈性儲水係數、彈性儲水係數等。所輸入之土壤參數值，水力傳導
係數介於 0.0001 至 0.004 間、砂土之彈性儲水係數為 0.00001、黏土之非彈性儲水係數介於
應用商用地下水模擬軟體 Modflow 以及一維壓密模式 Interbed(Leake,(1990)；Leake 與 Prudic(1991))，以非耦合方式求解，模擬雲林地區抽水對地層下陷的影響。其地下水流方
h
為地下水水位，K xx
、K yy
與K zz
為水力傳導係數，W
為抽水量或補注量，S 為儲 s
水係數，Δ h
為地下水位變化，S skv
為非彈性儲水係數，S skv
為彈性儲水係數。模擬時所需 輸入之資料分別有土體分層資料、分層地下水位觀測值、水力傳導係數、儲水係數、抽水 量、砂土之彈性模數、非彈性儲水係數、彈性儲水係數，論文中未詳述土壤參數如何給定。7.林美聆教授研究團隊
依據 Biot 軸對稱三維孔隙彈性理論，利用有限元素軟體 ABAQUS(Hibbit, (1998))模擬 台西與麥寮地區之地層下陷。以蚊港為模擬中心，其所使用之水流控制方程式為 數與楊氏模數分別介於 1.6×10-7 至 8.0×10-4(m/sec)間以及 0.825 至 2.489(MPa)間，而多井 分析所使用之滲透係數與楊氏模數分別介於 2.8×10-9 至 1.4×10-5(m/sec)間以及 3.30 至 9.9(MPa)間。模擬時，為了使數值計算穩定，時間間距需滿足