基于滲流―固體相互作用的三維水文地質參數估計研究

［摘要］為準確反映地下水系統的真實情況，并提高水文地質參數的估計精度，研究提出基于三維流固耦合數學模型，對三維水文地質參數進行反演估計。對研究構建的三維水文地質參數反演算法進行有效性研究，發現其導水系數反演誤差為5.54%，貯水系數計算誤差為2.74%，水位降深反演誤差為0.2 m，運算時間為34.2 min，其運算速度與精度性能較其他對比算法更高。綜上結果可知，研究構建的三維水文地質參數反演算法可應用于地下水流動特性和水文地質參數的研究中，為未來水資源管理和土地利用規劃提供理論參考。

［關鍵詞］滲流—固體；相互作用；水文地質；參數估計；水資源管理

地下水資源是人類生活和經濟發展的重要組成部分。為了更好地管理和保護地下水資源，需對地下水的流動特性和水文地質參數進行準確估計。但地下水流動是一個復雜得多的物理過程，其涉及水力、滲透性、彈性和化學等多個方面的相互作用[1，2]。傳統的水文地質參數估計方法往往忽略了地下水系統中滲流和固體之間的相互作用，導致估計結果的不準確性[3，4]。三維數學模型在描述滲流―固體相互作用中具有良好的適用性，其可基于物理方程和實驗數據，模擬地下水的流動和固體介質的響應[5]。同時通過對水文地質參數的估計，可以更好地理解地下水流動的機理和特性[6]。故研究提出基于固體―流體相互作用的理論，通過構建固流耦合數學模型，并利用紐曼函數對該數學模型中的貯水率及滲透系數進行反演，實現水文地質參數的精確估計，提高數學模型的擬合能力。

1 試驗工程概況

研究針對貴州某工程進行現場抽水試驗，該試驗區域的面積約7000 m2，周長為377 m。試驗地區的底層及含水層劃分概況如圖1所示。

如圖1所示，試驗區的地質主要為松散巖類，地下水類型主要為松散巖類孔隙水。根據地下水的成因及年代，可將其分為淺部土層的潛水，及粉細砂中細砂層中的承壓水。實驗場景內共布置三口試驗觀察井，首先，研究在淺部土層的潛水區域布置一口試驗觀察井。這口試驗井將用于觀測和采集淺部土層的潛水水位和水質數據。隨后，研究在粉細砂中細沙層中的承壓水區域布置兩口試驗井。這兩口試驗井將用于觀測和采集承壓水的水位和水質數據。承壓水是深部地下水，在巖石或土層中受到一定壓力的水。通過分析承壓水的特性和壓力變化，研究可分析地下水的流動性和水文特征，進而分析基于滲流―固體相互作用的三維水文地質參數。

2 基于滲流―固體相互作用的三維水文地質參數估計

2.1 滲流―固體相互作用的三維數學模型構建

滲流―固體相互作用是指在多孔介質中，流體（水或氣體）與固體顆粒之間的相互作用過程[7]。在這個過程中，流體通過多孔介質的孔隙和通道，與固體顆粒表面發生接觸和相互作用[8]。在傳統的陸地工程數字仿真模型中，為保證計算的復雜度及計算時間，一般將固體假定為絕對剛體，對固流相互作用不作考慮[9]。但在地下水開采等水資源管理工程中，固體―流體的相互作用則不可被忽視。在流體―固體的相互作用中，流體和固體彼此影響，并通過相互作用來實現能量和質量的傳遞，該過程又稱流固耦合[10]。此外，研究還對多孔介質含水層中的流固耦合進行考慮，構建三維流固耦合數學模型，以更加準確地反映地下水流動實際運動情況。三維流固耦合數學模型的基本內容如圖2所示。

如圖2，三維流固耦合數學模型包括水的壓縮性方程、多孔介質的壓縮性方程、固體力學方程、多孔介質的壓縮方程及多孔介質的含水層參數。其中地下水控制方程是描述地下水彈性變形的基本方程，研究利用虎克定律和Taylor級數展開來表示其水流彈性變形，其計算公式如（1）所示。

式（1）中，ρ 表示流體密度；β 表示水體積的可變形系數；p 表示壓力。多孔介質的含水層可能存在非彈性變形及彈性變形的情況，因此，研究利用虎克定律對其進行計算，計算公式如（2）所示。

式（2）中，α 表示多孔介質體積的壓縮系數；δ 表示多孔介質骨架上的有效應力；Vb表示多孔介質的體積；Vs表示固體部分的體積；Vv表示空隙體積；n 表示多孔介質的孔隙度。水和多孔介質的壓縮性計算公式如（3）所示。

σ = λσs + （1 - λ） p （3）

式（3）中，σ 表示一維平衡狀態的垂直總壓應力；σs表示粒間應力的垂直分量；λ 表示單位面積中的粒間接觸面積水平投影。含水層參數包括貯水率及滲透系數，其是用于描述多孔介質中水的分布和流動特性的重要參數。其中貯水率是指多孔介質中孔隙的體積與總體積之比。表示多孔介質中水可以存儲的容量。貯水率越大，說明多孔介質中可以存儲的水越多。貯水率的計算公式如（4）所示。

式（3）中，σ 表示一維平衡狀態的垂直總壓應力；σs表示粒間應力的垂直分量；λ 表示單位面積中的粒間接觸面積水平投影。含水層參數包括貯水率及滲透系數，其是用于描述多孔介質中水的分布和流動特性的重要參數。其中貯水率是指多孔介質中孔隙的體積與總體積之比。表示多孔介質中水可以存儲的容量。貯水率越大，說明多孔介質中可以存儲的水越多。貯水率的計算公式如（4）所示。

式（4）中，μs表示貯水率；g 表示重力系數；H 表示地下水單位面積的剖面高度；dVb表示含水層體積降低排放的水量；dV 表示遇水體積增加排放的水量。滲透系數作為描述多孔介質中水流動性質的參數，其可表示為單位時間內單位面積的水通過多孔介質的能力。滲透系數越大，說明多孔介質中水的滲透性越好，流動性越強。而導水系數是多孔介質對流體的滲透能力的衡量指標，其表示單位時間內單位面積的介質能夠通過的水量。導水系數越大，說明介質滲透性能越好，水分能夠更快地通過。導水系數的大小受到介質孔隙度、孔隙連通性、孔隙大小和土壤顆粒組成等因素的影響。研究利用滲透試驗對導水系數的值進行確定。滲透系數和導水系數的計算公式如式（5）所示。

T = KM （5）

式（5）中，T 表示導水系數；M 表示滲透系數。隨后，研究基于以上參數進行固流耦合的數學模型的構建。

2.2 基于流固耦合模型的三維水文地質參數反演

完成固流耦合數學模型的構建后，研究利用紐曼函數對數學模型中的貯水率及滲透系數進行反演，獲取三維水文地質參數，以更加精確地描述多孔介質中水的流動變化。紐曼函數標準曲線分為A組曲線及B組曲線，分別代表標準曲線的上凹部分及下凸部分[11，12]。其反演的基本思想為通過對A組和B組的曲線與觀察井的實際時間―降深曲線（s-t曲線）進行最優配線，以獲取最優水文地質變量[13]。但由于固流耦合模型的復雜性，傳統的優化算法往往面臨著維度高、非線性強等問題，導致求解難度較大。而隨機權重粒子群優化算法（Random Weighted Particle Swarm Optimization，RandWPSO）采用了粒子群優化的思想，通過模擬鳥群的行為來搜索最優解，能夠有效地應對流固耦合模型的高維度和非線性問題，提高優化的效率[14，15]。故研究選擇RandWPSO算法作為流固耦合模型的優化配線方法，其配線流程如圖3所示。

如圖3所示，研究首先對RandWPSO算法的粒子維度、粒子種群大小、隨機權重平均值的范圍及最大迭代次數等參數進行設置。隨后，研究對粒子群的速度及位置進行初始化操作，并利用適應度函數對粒子群進行評價并對比，獲取最小適應度的粒子作為群體最優粒子，以此更新粒子速度及位置，直至到達最大迭代次數。隨后，研究對B組進行配線，獲取最優曲線的配線參變量。紐曼函數的配線方式包括兩種方式，研究分為方法1及方法2。以B組曲線為例，方法1的具體過程表示為，首先對s-t曲線進行繪制，并將該曲線的后半段對B組曲線進行擬合，并確定其最優配線參變量。隨后選取任意擬合點，計算B組導水系數和導水度。配線參變量通過計算曲線的導數和二階導數，找到曲線上的極值點和拐點來實現。根據最優配線參變量，即可獲得與B組曲線相反的曲線，即A組曲線。將A組曲線與s-t曲線的前半段進行最優擬合，計算A組導水系數和貯水系數。方法2與方法1除先后順序外，其計算方法均一致。方法2首先對A組曲線進行配線，獲取其配線參變量，隨后以該參變量進行B組的配線。A組和B組的導水系數應大體相近。研究將對方法1和方法2獲取的AB兩組導水系數差值進行計算，選取性能更優的配線方法作為參數的計算方法。最后，研究通過A組和B組配線的曲線及結果圖，輸出水文地質參數值。

3 基于三維流固耦合數學模型的水文地質參數反演算法實證分析

3.1 三維流固耦合數學模型有效性驗證

為驗證研究構建的三維固耦合模型的有效性，研究利用抽水試驗中觀察井的降深―時間觀測值作為數據集，對三維流固耦合模型進行有效性驗證。對比模型為一維流固耦合數學模型和二維流固耦合數學模型。性能評價指標為與實際水位降深值的偏差及水位降深值的計算精度。實驗環境為多物理場仿真軟件COMSOL。該三維流固耦合數學模型的水位降深模擬數據及實測水位降升對比結果如圖4所示。

圖4（a）為各三維固流耦合模型的擬合準確度對比結果。如圖4（a）所示，研究提出的三維固流耦合數學模型計算出的水位降深值與實際值更為接近，其與實際值的最大偏差為0.06 m，較一維流固耦合模型的最大偏差低0.09 m，較二維固流模型的最大偏差低0.07 m，故研究提出的三維固流模型擬合準確度更優。圖4（b）為各三維固流耦合模型的擬合精確度對比結果。如圖4（b）所示，研究提出的基于固流耦合模型計算出的水位降深值計算精確度為90.1%，較其他模型更高。綜上結果可知，研究提出的固流耦合模型擬合準確度及精度更高。

3.2 水文地質參數反演算法有效性驗證

為選擇參數反演算法中，性能更優的配線方法。研究分別對配線方法1及配線方法2進行計算結果分析。實驗數據集為抽水試驗的觀測值，實驗環境為MATLABE。兩種方法的計算結果如表1所示。

如表1所示，方法1中，A組的導水系數為176.73，B組的導水系數為130.82，A組與B組的導水系數差值為45.91。在方法2中，A組的導水系數為176.28，B組的導水系數為137.55，A組與B組的導水系數差值為39.73。綜上結果可知，使用配線方法2的導水系數差值較配線方法1的導水系數低，故配線方法2更符合實際情況，可應用于水文地質參數的反演中。最后，為驗證研究提出的基于固流耦合的參數反演算法性能，研究基于COMSOL及MATLABE對抽水試驗觀察井的觀測值進行計算。對比方法為基于一維固流耦合的參數反演算法、基于二維流固耦合的參數反演算法及基于三維非流固耦合的參數反演法。對比指標為導水系數、貯水系數、觀察井的水位降深及運行時間。各反演算法的計算結果對比如圖5所示。

圖5（a）為各反演算法的擬合準確度對比結果。如圖5（a）所示，研究提出的基于三維耦合模型的參數反演算法導水系數反演誤差最小，為5.54%，較其他反演算法的誤差更低。此外，該算法的貯水系數計算誤差為2.74%，其性能較其他模型更優。圖5（b）為各反演算法的水位降深與運行效率對比結果。如圖5（b）所示，研究提出的基于三維耦合模型的參數反演算法水位降深反演誤差最小，為0.2 m，較其他模型的誤差更小，此外，該算法的運算時間為34.2 min，較一維和二維的反演算法較高，但較三維非固流耦合低52 min。綜上結果可知，研究構建的基于三維固流耦合的參數反演算法可平衡水文地質參數的求解性能及求解速度，可對實際水位進行預測和求解。

4 結論

為更加精確地模擬地下水流動和貯水特性，研究利用滲流―固體相互作用構建三維流固耦合數學模型，并將該模型與RandWPSO算法進行融合，構建基于滲流―固體相互作用的三維水文地質參數反演算法，獲取水文地質參數。對研究構建的三維流固耦合數學模型進行有效性驗證，發現水位降深誤差為0.06 m，精確度為90.1%，模擬準確度及精確度性能較其他模型更優。此外，研究還對基于三維耦合模型的參數反演算法進行有效性驗證，發現其導水系數反演誤差為5.54%，貯水系數計算誤差為2.74%，水位降深反演誤差為0.2 m，顯著低于其他對比算法。此外，該算法的運算時間為34.2 min，較三維非固流耦合低52 min。綜合以上結果可知，研究提出的基于三維耦合模型的參數反演算法可對三維水文地質參數進行準確估計，并為地下水流動特性和水文地質參數的研究提供理論參考，進而對地下水資源管理和環境保護等領域提供有益的支持。

[參考文獻]

［1］高成路，李術才，周宗青，等. 巖體水力壓裂應力―滲流耦合近場動力學模擬［J］. 同濟大學學報：自然科學版，2022，50（4）：470-480.

［2］張培森，趙成業，侯季群，等. 溫度―應力―滲流耦合條件下紅砂巖滲流特性試驗研究［J］. 巖石力學與工程學報，2020，39（10）：1957-1974.

［3］樊輝，呂義清，鄭劭鐸. 基于滲流―應力耦合的采動堆積層滑坡變形分析［J］. 煤礦安全，2020，51（11）：281-289.

［4］劉慶生，肖浩，談成亮，等. 離子型稀土礦浸出過程的滲流―反應―應力耦合模擬研究［J］. 中國稀土學報，2022，40（5）：880-892.

［5］慎國強，劉浩杰，王振濤，等. 致密儲層復雜縫網流固耦合滲流規律及產能模擬［J］. 大慶石油地質與開發，2023，42（4）：157-165.

［6］周亞醒，陳浩，張飛，等. 群孔抽水試驗方法求解含水層水文地質參數［J］. 水文，2022，42（6）：78-82.

［7］劉淑惠，肖長來，梁秀娟. 基于線性、非線性規劃的雙遼市水文地質參數計算方法［J］. 水電能源科學，2020，38（4）：72-75.

［8］Brunetti C，Bianchi M，Pirot G，Linde N. Hydrogeological Model Selection Among Complex Spatial Priors［J］. Water Resources Research，2019，55（8）：6729-6753.

［9］樊冬冬，劉天任，譚勇，等. 富水砂層水文地質參數反演及基坑降水圍護優化［J］. 建筑技術，2023，54（13）：1587-1591.

［10］王旭東，呂樂，時俊，等. 抽水試驗反演水文地質參數的多線全程加權法［J］. 水利學報，2020，51（3）：276-285.

［11］張淵博，鄒德旋，張春韻，等. 自適應慣性權重的粒子群優化算法［J］. 計算機仿真，2023，40（4）：350-357.

［12］Marshall S K，Cook P G，Simmons C T，Konikow L F，Dogramaci S. An Approach to Include Hydrogeologic Barriers With Unknown Geometric Properties in Groundwater Model Inversions［J］. Water Resources Research，2022，58（7）：1-22.

［13］石鴻蕾，郝奇琛，邵景力，等. 基于多源數據的弱透水層水文地質參數反演研究——以呼和浩特盆地某淤泥層為例［J］. 水文地質工程地質，2021，48（2）：1-7.

［14］李磊，陳干，唐沛，等. 基于解析法和數值法反演哈爾濱漫灘區水文地質參數［J］. 城市軌道交通研究，2021，24（10）：54-59.

［15］邱淑偉，吳亞敏，柯昱琪，等. 基于遍歷搜索算法的水文地質參數優化求解［J］. 吉林大學學報：地球科學版，2020，50（6）：1854-1861.