李 飛，李松達，單圣蘇，唐石振，王霞雨，劉佩貴

(合肥工業(yè)大學土木與水利工程學院，安徽 合肥 230009)

0 引言

地下水在維持沙漠綠洲區(qū)生態(tài)環(huán)境健康方面發(fā)揮著至關重要的作用，但受地理位置和參數(shù)空間變異性等因素的制約，地下水數(shù)值模型建立過程中存在眾多不確定性因素，嚴重影響到地下水資源評價結果及數(shù)值模型結果的可靠度。因此，需要開展數(shù)值模型的參數(shù)靈敏度分析，目前局部靈敏度分析方法較成熟，也廣泛應用于地下水數(shù)值模型中[1]。全局靈敏度分析方法的研究成果主要有多元回歸法[2]、Morris法[3]、擴展傅里葉幅度靈敏度檢驗法(EFAST)[4]以及Sobo’l法[5]等。為進一步對比分析局部與全局靈敏度分析結果的異同，本文以塔克拉瑪干沙漠南緣的綠洲為例，在參數(shù)局部靈敏度分析的基礎上結合全局靈敏度分析方法(Morris法)，以地下水水位變幅作為分析對象，對該綠洲地區(qū)的地下水數(shù)值模型進行局部和全局靈敏度分析，以期為指導沙漠綠洲區(qū)地下水資源的開發(fā)利用提供技術支撐。

1 研究區(qū)地下水數(shù)值模型

1.1 研究區(qū)概況

研究區(qū)位于塔克拉瑪干大沙漠的西南緣的皮山縣，區(qū)內主要河流為皮山河。流域內賦存有第四系松散堆積層孔隙潛水；第三系碎屑巖裂隙孔隙水及層間承壓水；以泥盆系、志留系、二迭系、侏羅系和白堊系為主的碎屑巖裂隙潛水；地下水的賦存、分布具有明顯的分帶性，水文地質結構如圖1所示。

圖1 山前水文地質剖面示意圖

皮山河發(fā)源于喀喇昆侖山北麓，上游由阿克肖河及康艾孜河兩大支流組成，兩支流在瑙阿巴提塔吉克民族鄉(xiāng)匯合，由南向北流經克里陽注入雅普泉水庫，后流經皮山縣城消失于塔克拉瑪干大沙漠，全長160 km，皮山河流域形似一羽毛狀。皮山河水文站以上集水面積為1 894.7 km2，徑流系數(shù)為17.52×104m3/km2，海拔高程在2 300～5 500 m之間，河流長度為72 km，河源最高點海拔6 000 m以上。在5 500 m以上的高山區(qū)分布有大量的冰川和積雪。多年平均降水量僅171.2 mm左右，雪線以上冰川和永久積雪區(qū)年降水量可達500 mm左右，出山口后平原區(qū)多年平均降水量僅為48.2 mm左右，高山區(qū)冰雪消融水和中淺山區(qū)降雨是皮山河徑流的主要補給源。

1.2 地下水數(shù)值模型

因研究區(qū)所在的水文地質單元范圍較大，結合數(shù)據(jù)及資料情況，模擬區(qū)范圍南界和北界分別以1 340 m和1 500 m等高線為界，東界為兵團農場水庫和干渠，西界垂直于地下水水位埋深等值線，面積約412 km2。模擬區(qū)的東側、北側和南側均概化為水頭邊界，西側為零流量邊界；勝利干渠和兵團農場干渠為一混合邊界，在灌溉期間，接受干渠的滲漏補給，非灌溉期為零流量邊界。垂向上接受大氣降水的入滲補給及水面蒸發(fā)，模擬區(qū)范圍及邊界條件如圖2所示。

圖2 模擬區(qū)范圍及邊界條件

根據(jù)水文地質概念模型，將模擬區(qū)地下水流概化成非均質各向同性非穩(wěn)定二維地下水流系統(tǒng)，并建立相應的數(shù)學模型。

(1)

式中：K為含水層的滲透系數(shù)(m/d)；h為地下水水位(m)；z為潛水含水層底板(m)；W為單位體積流量，用以代表流進源或流出匯的水量；μ為給水度；H0為初始水位(m)；t為時間(d)；D為模擬區(qū)范圍；Γ1為一類邊界；Γ2為二類邊界；q為含水層側向單寬補排量(m2/ d)，隔水邊界為零。

2 靈敏度分析

靈敏度分析可以對模型中參數(shù)的不確定性對模型的輸出結果的影響程度進行評價與量化。故而，在地下水數(shù)值模型中，參數(shù)靈敏度分析是必不可少的[6]。

參數(shù)靈敏度分析主要包括局部靈敏度分析與全局靈敏度分析[7-8]。局部靈敏度分析是用來分析單個參數(shù)變化對數(shù)值模型輸出結果的影響程度，不考慮參數(shù)之間的耦合作用。從數(shù)學角度來說 ，在某設計點Xk 處某個設計函數(shù) Fj(X)對設計變量xi的靈敏度可由式(2)計算得到[9]：

(2)

其中|Sji|為函數(shù)Fj(X)對Xi的靈敏程度，其值越大，則說明靈敏度越大，對數(shù)值模型的影響程度也就越大。

根據(jù)式(2)并結合當?shù)刭Y料，本次數(shù)值模型中選用的參數(shù)包括滲透系數(shù)(K)和給水度(S)，以及刻畫邊界條件的水力傳導系數(shù)(C)。首先對滲透系數(shù)K進行靈敏度分析。保持其他兩個參數(shù)數(shù)值不變，將K值變化依次 ±10%、±20%、±30%，每改變一次K值，記錄一次對應的地下水位變幅，其他參數(shù)的靈敏度分析同理，得到參數(shù)局部靈敏度分析結果如圖3所示。

圖3 局部靈敏度分析結果

由圖3可知，不論參數(shù)增大還是減小，滲透系數(shù) K值對模型地下水位影響最大，即K值的靈敏度最高。例如當參數(shù)變化均為-10%時，滲透系數(shù)K、水力傳導系數(shù)C、給水度S變化引起的地下水位變幅分別為-1.083 m、-0.155 m、-0.037 m。此外當參數(shù)變幅為負值，即參數(shù)值減小時，C較S靈敏度更高，變幅為正時，S較C的靈敏度偏大，但因對數(shù)值模型得到的水位值的影響程度均較小。由于該研究區(qū)地下水數(shù)學模型是非線性模型，因此參數(shù)增大和減小相同百分比時，圖像并不是關于y軸完全對稱，即參數(shù)增大和減小相同百分比時，對模型的影響并不完全相同，這也間接反映出參數(shù)的變化幅度對水位輸出結果的影響并不是線性變化的。

局部靈敏度分析方法操作簡便，但未考慮參數(shù)之間的相關性及相互影響的作用，以及參數(shù)之間會相互影響繼而影響數(shù)值模型的輸出結果。而全局靈敏度分析方法可以考慮多個參數(shù)之間的耦合影響，因此可以彌補局部靈敏度分析方法的不足。

根據(jù)局部靈敏度分析選取的三個參數(shù)進行組合，得到KS、KC、SC、KSC四種組合形式，并依次同時變化±10%、±20%、±30%，分別運行相應的數(shù)值模型得到變化后的地下水水位變幅，結果如圖4所示。

圖4 參數(shù)全局靈敏度分析

由圖4可知，三個參數(shù)KSC同時變化時，對模型輸出結果影響最大，與數(shù)值模型的理論分析結果相一致，KC、KS、SC組合對模型水位輸出結果影響幅度相應的有所降低。例如，當參數(shù)組合均減小10%，KSC引起的地下水位變幅為1.315 m，而KC、KS、SC引起的地下水位變幅分別為1.301 m、1.146 m、0.156 m，均小于1.315 m。

將局部靈敏度分析結果與全局靈敏度分析結果對比，K單獨變化+10%時，地下水變幅為0.925 m，K分別與S、C的參數(shù)組合KS、KC引起的地下水位變幅為0.958 m、0.939 m，由K單獨引起的變化量分別占到96.56%和98.55%。K單獨變化-30%時，地下水位變幅為-3.913 m，KS、KC引起的地下水位變幅為-4.339 m、-4.593 m，此時，由K單獨引起的變化量分別占到90.18%和85.19%。由此可以看出，該研究區(qū)地下水數(shù)值模型靈敏度最大的參數(shù)是K，但同時也應注意到，在變幅較小時，K的變化量在參數(shù)組合中占比極大(96%以上)，變幅增大時，K的變化量雖然仍舊占據(jù)絕大部分(85%以上)，但是其比例有所下降。因此，在實際應用中，如果對數(shù)值模型的輸出精度要求較高，則需重點考慮滲透系數(shù)K和C的變化及資料豐富程度，反之，如果對模型輸出結果精度要求不高，考慮到人力、物力、財力等因素的影響，其他靈敏度較低的參數(shù)對模型的影響可暫不考慮。

由上分析可知，全局靈敏度分析與局部靈敏度分析結果類似，滲透系數(shù)K值的影響最大，由于全局靈敏度分析考慮了各個參數(shù)之間的耦合作用，因此更加貼近實際。對參數(shù)進行靈敏度分析，可以使地下水數(shù)值模型更加完善，更好服務生產生活。

3 結語

本文通過對新疆皮山縣塔克拉瑪干沙漠南緣的綠洲建立數(shù)值模擬模型，圍繞數(shù)值模型中的參數(shù)進行局部和全局靈敏度分析，得出在所選的三個參數(shù)含水層滲透系數(shù)K、給水度S、水力傳導系數(shù)C中，含水層滲透系數(shù)K的靈敏度最高，給水度S次之，C影響不明顯。其排序為|K|>|S|>|C|。全局靈敏度分析可以考慮各個參數(shù)之間的相互關系對地下水數(shù)值模型輸出結果的影響，較局部靈敏度分析方法更貼合實際。對參數(shù)進行靈敏度分析，可以得出對所構建的地下水數(shù)值模型影響最大的參數(shù)，在構建模型以及率定相關參數(shù)等方面可以重點關注，為決策者做出更加合理的決策提供依據(jù)。