付 鵬，范榮亮，黃 旭，王光磊

（1.松遼水利委員會水文局，吉林 長春 130021；2.徐州市水務局，江蘇 徐州 221018）

洮兒河流域洪水流量較大，但缺少地表水控制性工程，致使洪水資源難以利用。下游扇形地潛水含水層埋藏淺、厚度大，具有較大的地下調蓄空間，是洮兒河流域最具開發(fā)利用價值的地下水庫。為研究地下水庫的調蓄功能，實現(xiàn)地表水與地下水的聯(lián)合配置，構建、開發(fā)了洮兒河流域扇形地下水庫調控模型。

1 洮兒河流域概況

洮兒河是嫩江下游右側的一大支流，流經(jīng)內蒙古興安盟的科右前旗、突泉縣和吉林省的白城、洮南、鎮(zhèn)賚、大安等市縣，流域面積約為3萬km2。干流及其主要支流歸流河和蛟流河均發(fā)源于大興安嶺東南麓，自西北向東南流入嫩江，干流總長395 km。流域形狀狹長，東西寬，南北窄。東為嫩江，西為大興安嶺，南為霍林河，北與綽爾河相鄰。流域整體上位于松嫩平原西部向大興安嶺山地的生態(tài)過渡帶，屬于典型的大陸性干旱-半干旱季風氣候區(qū)。多年平均降雨量為458.4 mm，6～9月降雨量占全年降水量的85.9%。多年平均蒸發(fā)量為1780～1910 mm。

洮兒河扇形地位于大興安嶺東麓山前傾斜平原，處于松嫩平原西緣，占地面積約2924 km2。洮兒河與蛟流河自大興安嶺東坡進入平原，形成了沖積洪積扇，構成了山前傾斜平原。地面由扇頂向前緣傾斜，坡降0.1%～0.3%。扇頂海拔200～220 m，前緣海拔145～150 m，部分隱伏于松嫩平原之下，形成埋藏扇形地。洮兒河扇形地是一個巨大的天然地下水庫，具有良好的調蓄條件。

2 模型簡介

2.1 基本方程

洮兒河扇形地主要地下水類型為第四系孔隙潛水，含水層厚10～40 m，是模擬研究的目標含水層。該層孔隙潛水在天然狀態(tài)下的運動規(guī)律基本符合達西定律。水流運動形式可概化為非均質各向同性的平面二維流，用數(shù)學方程表示為：

式中：K——潛水含水層滲透系數(shù)；μ——潛水含水層給水度；W，P——分別為單位面積、單位時間上的垂向補給量和排泄量；h，H——分別為地下水位和潛水水位；B——含水層底板高程；h0，h1——分別為初始水和一類邊界點水位；q——二類邊界單寬流量；D——模擬區(qū)范圍。

2.2 基礎軟件

洮兒河流域地下水庫調控模型利用了美國楊百翰大學環(huán)境模型研究實驗室的GMS地下水模擬系統(tǒng)。該系統(tǒng)包含了MODFLOW、FEMWATER、MT3DMS、RT3D等計算模塊。其中的MODFLOW模塊是專門用于孔隙介質地下水的三維有限差分數(shù)值模擬軟件。因其程序結構的模塊化、離散方法的簡單化和求解方法的多樣化等優(yōu)點，已被廣泛用于模擬井流、河流、排泄、蒸發(fā)和補給等對非均質和復雜邊界條件的水流系統(tǒng)的影響。洮兒河流域地下水庫調控模型主要利用了GMS系統(tǒng)中的MODFLOW模塊及其輔助模塊。

3 扇形地GMS模型的構建

3.1 含水層的水力特征

考慮到地下水庫調控的含水層應具有良好的開啟性，同時綜合考慮含水層的供水意義和水資源開發(fā)利用現(xiàn)狀，將模擬目標含水層定為第四系潛水含水層。該層孔隙潛水在天然狀態(tài)下的水力坡度不大，一般小于0.1%，其運動規(guī)律基本符合達西定律，水流運動形式可概化為平面二維流。

3.2 邊界條件

通過研究，近20年的調查統(tǒng)計資料及地下水等水位線圖可以看出，扇形地西部邊界和北部邊界與等水位線近于垂直，可視為隔水邊界；由扇形地頂部流入的洮兒河和蛟流河河谷較寬，滲透性較好，含水層厚度為10 m左右，兩河入口可視為流量邊界；扇形地前緣地帶含水層逐漸由單層過渡為多層，為亞粘土和砂礫石夾層，地形開闊平坦，地下水的水力坡度很小，可視為弱透水邊界。

3.3 源匯項處理

模型的主要源項包括降水入滲補給量、河道與渠道滲漏補給量、側向徑流補給量和灌溉回歸入滲補給量；主要匯項包括潛水蒸發(fā)量、人工開采量和側向流出量。其中，降水入滲補給量、潛水蒸發(fā)量和灌溉回歸水量根據(jù)歷年的降水、徑流、蒸發(fā)等實測資料采用SWAT模型模擬計算。其它各源匯項的處理方法如下：

1）側向徑流補給量。地下水側向徑流補給量按下式估算：

式中：Qlr——山前側向補給量（104m3/a）；K——含水層滲透系數(shù)（m/d）；I——垂直于剖面方向上的水力坡度；M——含水層厚度（m）；L——剖面長度（km）。

2）河道滲漏補給量。通過模型校準，確定河床的水力傳導系數(shù)。MODFLOW可根據(jù)河道斷面數(shù)據(jù)、河道水位及河床的水力傳導系數(shù)計算河道滲漏補給量。

3）渠道滲漏補給量。渠道滲漏補給量計算采用滲漏系數(shù)法。渠道滲漏系數(shù)根據(jù)引洮工程渠段實測流量資料分析確定。渠道在引水初期滲漏系數(shù)較大，然后快速減小，5～6日后基本趨于穩(wěn)定，穩(wěn)定滲漏系數(shù)為0.15～0.20。

4）側向流出量。計算公式同側向徑流補給量公式。

5）人工開采量。根據(jù)各鄉(xiāng)鎮(zhèn)地下水實際開采量進行統(tǒng)計，分農業(yè)用水和其它用水兩部分。

按照上述方法，可得到近20年扇形地多年平均地下水的補給量與排泄量（見表1）。由表1可知，地下水多年平均總補給量為5.0823×108m3/a，扣除地下水灌溉回歸水量，地下水資源量為4.3182×108m3/a。降水入滲補給量與河道滲漏補給量占地下水資源量的89.6%。

表1 扇形地多年平均地下水補給、排泄量 108m3/a

3.4 計算區(qū)剖分及參數(shù)分區(qū)

1）計算區(qū)剖分。扇形地地下水雖為潛水，但在扇形地前緣夾有第三系粘性土層，導致下層含水層具有一定的承壓性。因此綜合考慮水文地質及地形地貌等條件，將扇形地含水層概化為三層，每一層剖分為4640個單元，單元平均面積為0.632 km2。共有結點4871個，其中邊界結點117個。邊界結點中一類邊界結點75個，有水位觀測資料的邊界結點18個；二類邊界結點42個。共利用觀測井137個，其中邊界結點井18個，內部結點井119個。扇形地網(wǎng)格剖分情況如圖1所示。

圖1 研究區(qū)網(wǎng)格剖分圖

2）參數(shù)分區(qū)。對每個含水層的參數(shù)進行分區(qū)。其中第一層和第二層分為3個大區(qū)、6個小區(qū)；第三層為扇形地前緣承壓含水層，劃為6個分區(qū)。在MODFLOW模型中，降水入滲補給量與潛水蒸發(fā)量可通過SWAT模擬結果確定，另外還需要確定的參數(shù)主要是滲透系數(shù)及含水層給水度。這兩個參數(shù)初始值（見表2），根據(jù)前述參數(shù)確定方法給定，再根據(jù)實測資料進行率定。

表2 扇形地各含水層分區(qū)初始參數(shù)表

4 模型的識別與檢驗

在模型的識別和檢驗過程中，為了如實反映1998年洪水后洮兒河和蛟流河持續(xù)斷流的情況，將洮兒河和蛟流河的滲漏轉化為注水井的形式。沿河道設置一系列注水井，注水量與河流流量的變化一致。

1）模型的識別。選擇2000年9月1日～2001年6月1日，作為模型的識別期，共274 d。由于模型識別期較短，且處于枯水期，源匯項較少，數(shù)據(jù)擬合相對比較容易。根據(jù)所掌握的識別期實測資料，以一個月為一個應力期，將識別期劃分為9個應力期。每個應力期平均劃分為3個時段，平均10天為一個時段。

由于模擬區(qū)域條件比較簡單，且識別期實測資料翔實，加之輸入源匯項和參數(shù)能較真實地反映實際的水文地質條件，因而模型識別取得了較好的結果（見圖3）。地下水位的模型計算值與實測值擬合誤差小于0.5 m的結點觀測井數(shù)和結點數(shù)均達到總數(shù)的75%以上。由圖3可見，實測與計算水位的等值線在整體上達到了較好的擬合，說明所構建的數(shù)學模型真實地刻畫了扇形地的水文地質特征。

圖2 模型識別期水位擬合圖

圖3 模型檢驗期水位擬合圖

2）模型的檢驗。為進一步驗證所建立的數(shù)值模型和模型參數(shù)的可靠性，利用2001年與2002年地下水位動態(tài)觀測資料對模型進行檢驗。模型檢驗期為1年，檢驗期時段的劃分同識別期一致。

將檢驗期的源匯項輸入已識別好的數(shù)學模型，以2001年6月1日實際觀測水位作為初始水位運行模型。檢驗期末（2002年6月1日）計算結果顯示：地下水位的模型計算值與實測值擬合誤差小于0.5 m的結點觀測井數(shù)和結點數(shù)均占總數(shù)的75%以上。計算水位與實測水位等值線的整體擬合程度很好，說明模型的構建是成功的。

5 結論

1）洮兒河扇形地下水庫調控模型基于成熟的MODFLOW軟件，經(jīng)綜合分析當?shù)厮牡刭|條件并合理確定源匯項及分區(qū)參數(shù)后得以構建。檢驗結果說明所建立的數(shù)學模型、邊界條件、水文地質參數(shù)選取和源匯項處理是合理和可靠的，

2）所建立的數(shù)值計算模型能夠真實反映洮兒河扇形地地下水動態(tài)變化特征，可以將該模型應用于洮兒河扇形地下水庫調控，并進一步應用于洮兒河流域地表水與地下水資源的聯(lián)合配置與利用。