基于參數(shù)隨機分布的潛流帶地下水?dāng)?shù)值模擬研究

摘 要：含水層水文地質(zhì)參數(shù)的合理性是決定地下水?dāng)?shù)值模型模擬精度的關(guān)鍵因素之一。以淮河流域清水河典型河段為研究區(qū)，通過含水層土體野外監(jiān)測和室內(nèi)試驗，建立了滲透系數(shù)、給水度和貯水率的概率密度分布函數(shù)，并利用參數(shù)中值和置信度為95%的置信區(qū)間，建立了基于多參數(shù)概率密度分布的地下水?dāng)?shù)值模型。結(jié)果表明：相較于采用地質(zhì)勘探建議值的傳統(tǒng)地下水?dāng)?shù)值模型，基于多參數(shù)概率密度分布的地下水?dāng)?shù)值模型的地下水位模擬結(jié)果的精確度得到了顯著提升；河流水位、河床滲透系數(shù)等直接影響河流與地下水之間的水量交換強度，以及地下水的滲流和儲存特征。

關(guān)鍵詞：地下水；數(shù)值模擬；參數(shù)隨機分布；清水河

中圖分類號：TV62；TV882.1 文獻標(biāo)志碼：A doi：10.3969/j.issn.1000-1379.2025.02.016

引用格式：王輝，劉軍生，王荊，等.基于參數(shù)隨機分布的潛流帶地下水?dāng)?shù)值模擬研究[J].人民黃河，2025，47（2）：107-112，155.

基金項目：引江濟淮工程（河南段）工程科研服務(wù)項目（HN-YJJH/JS/FWKY-2021001）；黃河水利科學(xué)研究院基礎(chǔ)研究業(yè)務(wù)費專項（HKY-JBYW-2023-07）；水利部堤防安全與病害防治工程技術(shù)研究中心開放課題（LSDP202402）

NumericalSimulationofGroundwaterinSubsurfaceFlowZoneBasedonRandom DistributionofParameters

WANGHui1，LIUJunsheng2，3，WANGJing2，3，4，ZHAOShougang2，3，SUNRuidong2，3，WANGTianye4，ZHANGYibo2，3，LANYan2，3，YANGHaoming2，3，HOUJiaojian2，3

（1.HenanYangtze？to？HuaiheWaterDiversionCo.，Ltd.，Shangqiu476200，China；2.YellowRiverInstituteofHydraulicResearch，YRCC，Zhengzhou450003，China；3.ResearchCenteronLeveeSafetyandDisasterPrevention，MWR，Zhengzhou450003，China；4.SchoolofWaterConservancyandTransportation，ZhengzhouUniversity，Zhengzhou450001，China）

Abstract：Thereasonablenessofaquiferhydrogeologicalparametersisoneofthekeyfactorstodeterminethesimulationaccuracyofgroundw？ aternumericalmodels.Takingthetypicalcross？sectionoftheQingshuiRiverintheHuaiheRiverBasinasthestudyarea，theprobability densitydistributioncharacteristicfunctionsofthepermeabilitycoefficient，thedegreeofwatersupplyandtherateofwaterstoragewereestab？ lishedthroughthefieldmonitoringofaquifersoilsandindoorexperiments，andthenumericalgroundwatermodelbasedonthemulti？parame？ terprobabilitydensitydistributionwassetupbyusingthemedianvalueoftheparameterandtheconfidenceintervalwithaconfidencelevelof 95%.Theresultsshowthattheaccuracyofthegroundwatersimulationresultsofthegroundwaternumericalmodelbasedonmultiparameter probabilitydensitydistributionissignificantlyimprovedcomparedwiththetraditionalgroundwaternumericalmodelusingthesuggestedvalues ofgeologicalexploration.Thedifferentsimulationscenarios，suchastheriverlevelandthepermeabilitycoefficientoftheriverbed，directly affecttheintensityofwaterexchangebetweentheriverandthegroundwater，aswellastheseepageandstoragecharacteristicsoftheground？ water.

Keywords：groundwater；numericalmodeling；parameterstochasticscheme；QingshuiRiver

0 引言

潛流帶作為地表水和地下水相互作用最強的區(qū)域，水量相對豐富和穩(wěn)定，其水文動態(tài)對水資源評估、開發(fā)和管理具有重要影響。此外，潛流帶孕育的河岸生態(tài)系統(tǒng)具有重要生態(tài)功能和價值，尤其在干旱地區(qū)[1]。水文地質(zhì)參數(shù)是表征含水層或含水層介質(zhì)滲透性、儲水或釋水特性的指標(biāo)，是模擬地下水定量運移的關(guān)鍵參數(shù)。潛流帶地層沉積環(huán)境復(fù)雜[2]，水文地質(zhì)參數(shù)具有很強的空間異質(zhì)性，表現(xiàn)為大尺度的隨機分布[3-4]。為了描述水文地質(zhì)參數(shù)的空間異質(zhì)性，有學(xué)者基于參數(shù)隨機分布理論，探討了參數(shù)統(tǒng)計特征與概率密度分布函數(shù)之間的關(guān)系。隨機場理論最初由Vanmarcke[5]于1977年提出，此后被廣泛應(yīng)用于描述巖土參數(shù)的空間分布特征。近年來，描述土壤空間自相關(guān)特性的隨機場模型得到了廣泛應(yīng)用。Srivastava等[6]采用隨機有限元法模擬了土壤滲透系數(shù)的隨機場，結(jié)合蒙特卡羅模擬方法，研究了滲透系數(shù)的空間變化對邊坡滲流和邊坡安全系數(shù)的影響；Santoso等[7]利用局部平均法對土壤滲透系數(shù)的隨機場進行了離散化處理；Zhu等[8]采用快速傅里葉變換方法模擬了滲透系數(shù)的空間變化；Cho[9]通過蒙特卡羅模擬方法，研究了滲透系數(shù)的空間異質(zhì)性對土石壩滲流和滲流梯度均值、方差等統(tǒng)計特征參數(shù)的影響。

雖然學(xué)者們對巖土參數(shù)的隨機場理論進行了大量研究，但相關(guān)研究成果的實際應(yīng)用主要是在工程地質(zhì)和邊坡穩(wěn)定性領(lǐng)域。在地下水?dāng)?shù)值模擬中，通常給含水層指定一個固定參數(shù)，往往忽略了其空間異質(zhì)性和頻率分布特征。Tang等[10]發(fā)現(xiàn)，河床水文地質(zhì)參數(shù)的隨機分布極大地影響了地下水?dāng)?shù)值模擬的精度，這表明傳統(tǒng)的均質(zhì)含水層假設(shè)在研究小尺度、高精度的透水帶方面存在局限性。本研究以淮河流域清水河典型河段野外監(jiān)測和室內(nèi)試驗為基礎(chǔ)，通過顆粒分析試驗，對土樣進行分類，并將地質(zhì)調(diào)查獲得的土壤參數(shù)與實驗室土壤水力性質(zhì)測試結(jié)果相結(jié)合，建立滲透系數(shù)、給水度和貯水率的概率密度分布模型。同時，基于建立的MODFLOW地下水?dāng)?shù)值模型，分析參數(shù)隨機分布對地下水?dāng)?shù)值模擬精度的影響，探究地下水在含水層非均質(zhì)情況下的運移規(guī)律。

1 數(shù)據(jù)與方法

1.1 研究區(qū)概況

研究區(qū)位于引江濟淮工程清水河河段，引江濟淮工程溝通長江、淮河兩大水系，潤澤安徽、惠及河南、造福淮河、輻射長江，具有保障供水、發(fā)展航運、改善水環(huán)境等綜合效益。研究區(qū)屬黃淮沖積平原，清水河全長137.3km，流域面積901km2，屬暖溫帶半濕潤大陸性季風(fēng)氣候區(qū)，氣溫、降水和風(fēng)向隨季節(jié)變化。年平均降水量為720～820mm，降水集中在6—9月，通常占全年降水量的70%以上。上游年蒸發(fā)能力1200～1400 mm，中下游則達1866mm。數(shù)值模擬區(qū)位于清水河河段重點監(jiān)測斷面16。地下水類型為第四系松散層孔隙潛水，主要賦存于砂壤土、輕粉質(zhì)壤土及粉細(xì)砂層內(nèi)，下部粉細(xì)砂層中地下水具有承壓性，砂壤土、粉細(xì)砂、輕粉質(zhì)壤土具有中等透水性，重粉質(zhì)壤土具有弱透水性（為相對隔水層）。

1.2 數(shù)據(jù)收集和監(jiān)測

研究區(qū)內(nèi)河道沿線的水文地質(zhì)條件、工程地質(zhì)條件及數(shù)字高程模型（DEM）由《引江濟淮河南省受水區(qū)供水配套工程可行性研究階段工程地質(zhì)勘察報告及圖冊》與現(xiàn)場地質(zhì)鉆孔數(shù)據(jù)綜合分析得出，下墊面條件主要通過衛(wèi)星遙感影像與現(xiàn)場無人機航拍獲取。通過ArcGIS軟件生成研究區(qū)DEM數(shù)據(jù)，為地下水?dāng)?shù)值模擬提供準(zhǔn)確的地表數(shù)據(jù)。

為獲取地層滲透系數(shù)空間分布特征，沿水平縱橫向布設(shè)3×3個地下水監(jiān)測及試驗鉆孔，見圖1。同時，在監(jiān)測孔取樣，垂直向揭露透水砂層，每孔按垂直深度間距0.5～1.0m取樣，分層采集試驗土樣，每層土樣開展密度、顆粒分析、給水度及滲透系數(shù)等試驗。

為了研究河流與地下水之間的相互作用關(guān)系，在距離河道50、100、150m處分別設(shè)置自動化水位監(jiān)測井。這些監(jiān)測井采用VWP-0.35G型地下水位傳感器，每24h監(jiān)測1次。鑒于農(nóng)田是研究區(qū)域的主要土地利用類型，且附近有大量灌溉機井，因此除自動化水位監(jiān)測井外，通過人工巡查的方式定期評估周邊地區(qū)的地下水開采活動。

1.3 數(shù)學(xué)模型

由于目前的勘察技術(shù)無法準(zhǔn)確獲取研究區(qū)的水文地質(zhì)條件，因此使用地下水運動數(shù)學(xué)模型來表達復(fù)雜的水文地質(zhì)體，建立數(shù)學(xué)公式來表達水文地質(zhì)條件與相關(guān)參數(shù)的函數(shù)關(guān)系。可用下列數(shù)學(xué)模型描述潛水水流運動：

此數(shù)學(xué)模型利用VisualMODFLOWFlex中的有限差分法求解。

1.4 隨機統(tǒng)計方法及隨機參數(shù)分析理論

隨機變量可以分為離散型和連續(xù)型兩類。常見的連續(xù)分布密度函數(shù)包括均勻分布、指數(shù)分布、正態(tài)分布、t分布等。在統(tǒng)計和概率論中，均勻分布是一種常見的概率分布，隨機變量均勻分布指在一定區(qū)間內(nèi)，隨機變量每個值出現(xiàn)的概率都是相同的；指數(shù)分布具有無記憶性，是特殊的伽馬分布；正態(tài)分布是一種極為常見的連續(xù)型概率分布，均值μ決定了分布的位置，標(biāo)準(zhǔn)差σ決定了分布的幅度和形狀；t分布用于根據(jù)小樣本來估計呈正態(tài)分布且方差未知的總體的均值，如果總體方差已知或者樣本數(shù)量足夠多，則可以用正態(tài)分布來估計總體均值。

通過單樣本K-S檢驗法檢驗樣本是否來自于特定的理論分布，單樣本K-S檢驗的判斷標(biāo)準(zhǔn)為顯著性。

1.5 模型精度評估

采用納什效率系數(shù)SNSE、均方根誤差SRMSE對地下水運移數(shù)值模型精度進行評估。納什效率系數(shù)計算公式為

2 地下水?dāng)?shù)值模型的構(gòu)建

對研究區(qū)邊界條件進行概化，西邊界以清水河為界，水頭邊界根據(jù)清水河水位設(shè)定；模型南北跨度為650m，南北邊界垂直于河流方向，考慮地下水流的特點，將其定為零流量邊界；東邊界由序列化監(jiān)測井劃定為給定水頭邊界。根據(jù)鉆孔披露的信息，研究區(qū)主要由淤泥質(zhì)黏土、重淤泥質(zhì)壤土和壤土組成，黏土和砂質(zhì)壤土呈薄層或透鏡狀分布，在垂直方向上模型含水層分為4層。

考慮計算工作量和精度要求，將模型空間分辨率設(shè)定為5m，模擬時間步長設(shè)定為1d。模型被細(xì)分為121行、185列，共22385個矩形網(wǎng)格單元，計算節(jié)點位于網(wǎng)格單元中心。地下水?dāng)?shù)值模型模擬時間從2022年8月26日開始，時間跨度為8個月。

2.1 水文地質(zhì)參數(shù)

研究區(qū)水文地質(zhì)參數(shù)具有空間異質(zhì)性，通過對試驗結(jié)果進行統(tǒng)計分析，可以得到相關(guān)參數(shù)的概率密度分布函數(shù)。第一層含水層為重粉質(zhì)壤土，垂直滲透系數(shù)服從指數(shù)分布，水平滲透系數(shù)和給水度服從t分布，貯水率符合指數(shù)分布。由于這些參數(shù)在數(shù)量級上存在較大差異，因此采用對數(shù)變換進行統(tǒng)計分析，結(jié)果見圖2。

根據(jù)不同土層滲透系數(shù)、給水度、貯水率概率密度分布曲線，得出各土層滲透系數(shù)、給水度和貯水率地質(zhì)勘探建議值、中位數(shù)和置信度為95%的置信區(qū)間，見表1。

2.2 源匯項的確定和處理

模型中地下水運動主要受匯源項輸入的影響，模擬區(qū)源匯項有點、線、面3種要素。點狀要素為農(nóng)業(yè)開采等源匯項；線狀要素主要為河流補給項，通過River模塊處理；面狀要素主要為降雨入滲、灌溉入滲等補給項，由Recharge模塊給出。潛水蒸發(fā)等其他源匯項資料通過Well模塊代入模型計算，均處理成開采井或者補給井。

2.3 初始條件

由于研究區(qū)水文地質(zhì)資料不足，因此無法直接獲得初始流場。根據(jù)已有DEM數(shù)據(jù)、降雨資料、流域流入和流出數(shù)據(jù)，以及研究區(qū)各含水層水文地質(zhì)參數(shù)，模擬2022年地下水流場穩(wěn)定流，獲得穩(wěn)定流初始流場。

2.4 數(shù)值模型

根據(jù)水文地質(zhì)概念模型，將水文地質(zhì)參數(shù)及其他邊界條件輸入模型，基于土體參數(shù)空間隨機分布理論中參數(shù)分區(qū)與取值方法，針對模擬區(qū)水文地質(zhì)條件和源匯項特征，利用初始條件和邊界條件對地下水流動微分方程組進行求解，得到研究區(qū)地下水流場數(shù)值模型，地質(zhì)結(jié)構(gòu)與水位分布見圖3。

3 模型驗證與識別

3.1 模型驗證

為探究多參數(shù)概率密度分布在構(gòu)建地下水?dāng)?shù)值模型方面是否具有優(yōu)勢，建立了兩組地下水?dāng)?shù)值模型。兩組模型保持邊界、源匯項、時間步長等條件不變，模型Ⅰ作為對照組，采用地質(zhì)勘探建議的水文地質(zhì)參數(shù)（平均值）；模型Ⅱ?qū)?shù)中值作為初始模型參數(shù)，同時將95%置信區(qū)間作為參數(shù)反演范圍。使用Visual MODFLOWFlex中內(nèi)嵌的PEST反演模塊對模型參數(shù)進行校準(zhǔn)，最終確定模型參數(shù)。2022年8月24日至2023年2月24日為模擬期，2023年2月24日至2023年4月24日為驗證期，模型Ⅰ和模型Ⅱ模擬的部分鉆孔的地下水位與實測水位見圖4。

模型Ⅰ對QSH5-1、QSH5-2、QSH5-3的水位計算值相對于實測值的均方根誤差分別為0.47、0.65、1.13 m，確定系數(shù)R2分別為0.47、0.59、0.22，納什效率系數(shù)分別為0.73、0.02、-2.96。模型Ⅱ?qū)SH5-1、QSH5-2、QSH5-3水位計算值的均方根誤差明顯較小，分別為0.13、0.16、0.14m；R2顯著提高，分別為0.81、0.76、0.77；納什效率系數(shù)顯著增大，分別為0.98、0.94、0.94。比較模擬水位與實測水位可知，模型Ⅰ計算水位與實測水位最大差值出現(xiàn)在2022年10月24日QSH5-3處，差值為127cm，最小差值出現(xiàn)在2023年1月24日QSH5-3處，差值為69cm；模型Ⅱ計算水位與實測水位最大差值出現(xiàn)在2022年10月24日QSH5-3處，差值為27cm，較模型Ⅰ降低了78.7%，最小差值出現(xiàn)在2022年12月24日QSH5-2處，差值為0.9cm，較模型Ⅰ降低了98.7%。因此，考慮參數(shù)隨機分布的地下水?dāng)?shù)值模型計算結(jié)果與實測值更接近，模型擬合效果更好。

3.2 農(nóng)田輸水過程中地下水流場的動態(tài)變化

清水河2022年11月15日首次通水后，研究區(qū)地下水位受河道補給影響顯著，其中距離河道約50m處的觀測孔QSH5-1地下水位由34.389m（2022年12月2日）上升到最高水位36.160m（2022年12月18日），16d抬升1.771m；距離河道約100m處的觀測孔QSH5-2地下水位由34.554m（2022年12月5日）上升到最高水位35.667m（2023年1月5日），31d抬升1.113m；距離河道約150m處的監(jiān)測井QSH5-3地下水位由34.651m（2022年12月9日）上升到最高水位35.708m（2023年2月19日），72d抬升1.057m。

根據(jù)監(jiān)測結(jié)果可知，通水后地下水流場由周邊農(nóng)田補給河道轉(zhuǎn)換為河道補給周邊農(nóng)田；由于各監(jiān)測井與河道距離不同，因此其地下水位變化也存在較大差異。河道補給地下水存在一定滯后性，以監(jiān)測井QSH5-1、QSH5-3為例，垂向水位變化速率分別為0.295、0.015m/d，距離河道較近處地下水受河道補給影響遠(yuǎn)大于距離河道較遠(yuǎn)處的。

3.3 不同輸水條件下地下水流場變化

基于上述反演校正后的地下水?dāng)?shù)值模型，以模擬期第80d為例，分別把河流水位34.5、35.0、35.5、36.0、36.5、37.0m輸入模型，對比分析河流水位變化對地下水流場的影響，不同河流水位浸潤線見圖5。當(dāng)河流水位≤35.5m時，周邊農(nóng)田地下水補給河流；當(dāng)河流水位>35.5m時，河流補給地下水，地下水位在距離河道150～300m處逐漸趨于平穩(wěn)。其中，河流水位為35.5～36.0m時，河流對地下水位的影響較小，地下水位在距離河道約100m處已趨于平穩(wěn)；河流水位為36.5～37.0m時，河流補給地下水，且河流對地下水位影響距離較遠(yuǎn)，地下水位在距離河道250～400m處趨于平穩(wěn)。

3.4 河床滲透系數(shù)對河流滲流量的影響

分別選取各地層概率密度分布函數(shù)中滲透系數(shù)置信區(qū)間95%的小值、中值及大值組合（見表2），利用VisualMODFLOWFlex建立地下水?dāng)?shù)值模型，模擬計算河道正常調(diào)度水位條件下，河道運行期單側(cè)年滲流量。計算結(jié)果表明，最大、最小滲流量分別為57.2萬、25.2萬m3，滲流量變幅是最小滲流量的126.9%，由此可見河床滲透系數(shù)隨機性對河流滲流量影響十分顯著。

4 結(jié)論

基于清水河典型河段含水層土體野外監(jiān)測和室內(nèi)試驗，分析了滲透系數(shù)、給水度和貯水率等水文地質(zhì)參數(shù)變化，并建立了滲透系數(shù)、給水度和貯水率的概率密度分布函數(shù)，利用參數(shù)中值和置信度為95%的置信區(qū)間，建立了基于多參數(shù)概率密度分布的地下水?dāng)?shù)值模型，并將該模型與采用地質(zhì)勘探建議值的傳統(tǒng)地下水?dāng)?shù)值模型模擬結(jié)果進行了對比。結(jié)果表明，相較于采用地質(zhì)勘探建議值的傳統(tǒng)地下水?dāng)?shù)值模型，基于多參數(shù)概率密度分布的地下水?dāng)?shù)值模型的模擬水位與實測水位的最大差值比傳統(tǒng)模型降低了78.7%，最小差值比傳統(tǒng)模型降低了98.7%，模型模擬結(jié)果精確度明顯提升。

將基于多參數(shù)概率密度分布的地下水?dāng)?shù)值模型計算結(jié)果與實測水位進行比較，在相同外部條件下，當(dāng)河流水位增大至35.5m以上時，研究區(qū)潛水帶由地下水補給河流轉(zhuǎn)換為河流補給地下水；河床滲透系數(shù)由最小值增大到最大值時，河流單側(cè)年滲流量從25.2萬m3增加到57.2萬m3，滲流量變幅是最小滲流量的126.9%。這說明河流水位、河床滲透系數(shù)等直接影響河流與地下水之間的水量交換強度，以及地下水的滲流和儲存特征。

參考文獻：

[1] WANGTianye，WUZening，WANGPing，etal.Plant？Ground？ waterInteractionsinDrylands：AReviewofCurrentResearch andFuturePerspectives[J].AgriculturalandForestMeteorol？ ogy，2023，341：109636.

[2] 劉鑫，王宇，李典慶.考慮土體參數(shù)空間變異性的邊坡大變形破壞模式研究[J].工程地質(zhì)學(xué)報，2019，27（5）：1078-1084.

[3] 紀(jì)文貴，羅躍，劉金輝，等.考慮滲透系數(shù)不確定性的地浸過程溶浸范圍隨機模擬[J].原子能科學(xué)技術(shù)，2023，57（6）：1099-1110.

[4] 張俊，劉天罡，董佳秋，等.含水層層狀非均質(zhì)對地下水流系統(tǒng)的影響[J].中國地質(zhì)，2020，47（6）：1715-1725.

[5] VANMARCKEEH.ProbabilisticModelingofSoilProfiles[J].JournaloftheGeotechnicalEngineeringDivision，1977，103（11）：1227-1246.

[6] SRIVASTAVAA，BABULSG，HALDARS.Influenceof SpatialVariabilityofPermeabilityPropertyonSteadyState SeepageFlowandSlopeStabilityAnalysis[J].Engineering Geology，2010，110（3-4）：93-101.

[7] SANTOSOAM，PHOONKK，QUEKST.EffectsofSoilSpa？ tialVariabilityonRainfall？InducedLandslides[J].ComputersStructures，2011，89（11）：893-900.

[8] ZHUH，ZHANGLM，ZHANGLL，etal.Two？Dimensional ProbabilisticInfiltrationAnalysiswithaSpatiallyVarying PermeabilityFunction[J].ComputersandGeotechnics，2013，48：249-259.

[9] CHOSE.ProbabilisticAnalysisofSeepageThatConsiders theSpatialVariabilityofPermeabilityforanEmbankmenton SoilFoundation[J].EngineeringGeology，2012，133：30-39.

[10] TANGQ，KURTZW，SCHILLINGOS，etal.TheInfluence ofRiverbedHeterogeneityPatternsonRiver？AquiferEx？ changeFluxesUnderDifferentConnectionRegimes[J]. JournalofHydrology，2017，554：383-396.

【責(zé)任編輯 呂艷梅】