modflow計算地下水流程式原理.docVIP

1.1 Modflow計算地下水流程式原理 　　運用PMWIN(processing modflow)程式中之Modflow來地下水模擬，模式採用有線差分法(FDM)來求得近似數值解。模擬之前需將模擬區域劃分成網格(Cell)，每一網格假設為均勻物質，以網格節點或網格中心之水頭代表網格水頭，含水層厚度在水平方向上可以變化，但在數值上會有較大的誤差。 　　在Modflow中使用的固定密度(Constant Density)地下水流偏微分方程式(McDonald and Harbaugh，1988，p.2-1)如下所示：　　　　　　式1-1 式中、、分別為沿著x、y、z座標軸的水力傳導係數。 　　h=壓力水頭(Potentiometric Head) 　　WW0.0表示流出地下水系統(Out Flow of Ground-Water System)，反之，W0.0表示流進地下水系統。 　　Ss=孔隙介質之儲水率。　　t = 時間(T)。 　　當式1-1結合邊界及初始條件時，描述在非均質性及非等向性介質(提供主軸方向的水力傳導係數)中的暫態三維地下水流。左側代表含水層在不同方向之入滲與水頭變化及地下水補給與抽水造成之三向度地下水流動模式；在均向性穩定流含水層系統中，因為流率為一常數，多孔介質之儲水量為定值，此式即為Laplace方程式，可由給定之含水層流況、初始水位狀況與邊界條件獲得正確解。惟在大多數非穩定流系統中，含水層之水位變化在整個模式中是個變數，因此需將水層系統分隔成小區塊，以求得滿足各小區塊間之解與整個系統之近似解；在此假設環境下每一小區塊在單位時間內體積流率即可設定是一常數(McDonald and Harbaugh，1988)。 　　而由連續方程式(Continuity Equation)來發展有限差分形式之地下水流方程式，亦表示網格流入及流出的總和要等於網格中儲存量的變化率。在地下水密度固定的假設下，連續方程式表示網格中地下水流平衡方程式為式1-2所示： 　　　　　　　　式1-2 式中：流入及流出網格之地下水流量之總和 　　　：在有限差分法方程式中的比儲水係數。比儲水係數的定義為每單位含水層結構體積及每單位壓力水頭變化之抽出或補助地下水的變化。 　　 ：網格體積大小。 　　 ：在單位時間間隔中，壓力水頭的變化。 　　在式1-2等號右邊的項表示在已知的水頭變化及的時間間隔中，網格中儲水體積的變化率。因此，式1-2等號的左邊表示流入及流出網格的項及體積總和。 　　而外流力的作用如井水抽水、補注、蒸發、排水、河流入滲等，對網格之流量率影響亦可由下式1-3表示之(McDonald and Harbaugh，1988)： 　　　　　　式1-3 式中：由外界流入網格(i,j,k)之總體積流量率。 　　 ：外流力介質之導水係數。　　 ：網格中之水位。　　 ：直接注入網格之流量率。 　　 ｎ　：外流力數。 　　因此整個網格(i,j,k)之三向度水流有限差分連續方程式，即可由網格間之流量率變化與外流力作用比表示成式1-4： 　　　　　　　　　　　　　　　　　　　　　　　　　　　式1-4 　　其中右側之為不同時間之水位變化，亦可以有限差分法表示成式1-5： 　　　　　　　　　　　　　　　　　　　　　式1-5 　　 　　而將單位網格整體地下水流有限差分方程統整如下式1-6： 　　　　　　　　　　　　　　　　　　　　　　　　　式1-6 　　式1-6即為McDonald與Harbaugh(1988)發展之三向度地下水流有限差分連續方程式，可由起始水頭與邊界條件以反覆迭代(Interative Method)的方式，來得到每個時間階段(Time Step)之有限差分方程式系統的解。在每個時間階段的計算期間，初始壓力水頭分佈會在迭代法的迭代過程中，連續地產生許多中間的壓力水頭分佈。 1.2 地層下陷分析方法 　　因地下水位下降造成之地層下陷，其下陷量大小與土壤種類有關。砂質土層的沉陷速度較快但沉陷量較小；而黏性土層之壓縮則較具延滯性且沉陷量較大。 　　因地下水位下降所引致的地表壓密沉陷之分析方法大致可分為兩類，一為基於Terzaghi 壓密理論之＂分離式方法＂(decoupled approach)與另一為基於Biot理論之＂耦合式方法＂(coupled approach)。前者是先求出孔隙水壓力之分佈，再應用有效應力觀念計算土層之應變及地表壓密沉陷量。後者是基於孔隙水與介質之間存在某種互制關係，以介質位移和孔隙水壓力為基本變數之耦合壓密理論。 Terzaghi 壓密理論 　　本模式以Terzaghi(1943)之單向壓密理論分析地層下陷問題。該理論係利用飽和土壤內孔隙水在穩定層流情況下之連續條件，及有效應力原理和土壤本身之應力－應變