抽水試驗反演水文地質(zhì)參數(shù)的多線全程加權(quán)法

王旭東，呂 樂，時 俊，譚衛(wèi)佳

（南京工業(yè)大學(xué) 巖土工程研究所，江蘇 南京 210009）

1 研究背景

水文地質(zhì)參數(shù)是合理開展地下水流動分析的基本條件，現(xiàn)場抽水試驗是獲取場地水文地質(zhì)參數(shù)的重要途徑，傳統(tǒng)的標(biāo)準(zhǔn)曲線擬合法（配線法）僅將抽水階段的降深數(shù)據(jù)用于水文地質(zhì)參數(shù)反演，對水位恢復(fù)階段的降深數(shù)據(jù)利用較少。為了充分利用抽水階段和水位恢復(fù)階段的降深數(shù)據(jù)，肖長來等［1］提出了基于Theis 公式的全程曲線擬合法，用于承壓含水層模型的水文地質(zhì)參數(shù)反演。對于工程實踐中廣泛存在的越流承壓含水層模型，楊建民等［2］在Hantush 公式的基礎(chǔ)上，推導(dǎo)了越流承壓含水層剩余降深計算公式，為建立越流承壓含水層模型的全程曲線擬合法奠定了理論基礎(chǔ)。Li 等［3-4］將全程曲線擬合法用于越流承壓含水層水文地質(zhì)參數(shù)反演，研究表明全程曲線擬合法可獲得更合理的含水層水文地質(zhì)參數(shù)。

為了準(zhǔn)確掌握場地含水層系統(tǒng)的地下水流動特性，工程建設(shè)中常采用多觀測井抽水試驗確定場地水文地質(zhì)參數(shù)。王麒等［5］利用多觀測井抽水試驗資料，通過不同觀測井實測數(shù)據(jù)分別求取水文地質(zhì)參數(shù)，再以平均值作為場地水文地質(zhì)參數(shù)，但各觀測井反演求得的水文地質(zhì)參數(shù)差異明顯，平均值對場地含水層水文地質(zhì)特性的反映不盡合理。事實上，受成井工藝、洗井方法、試驗設(shè)備等抽水試驗外界因素的影響，同一抽水試驗中不同觀測井實測數(shù)據(jù)的可靠性存在差異。因此，在多觀測井抽水試驗反演水文地質(zhì)參數(shù)的過程中，有必要考慮觀測井實測數(shù)據(jù)的可靠性，以及對多個觀測井、抽水階段與水位恢復(fù)階段實測數(shù)據(jù)進行全程擬合。

多觀測井抽水試驗采集的數(shù)據(jù)量大且數(shù)據(jù)信息豐富，但存在多觀測井?dāng)M合目標(biāo)函數(shù)求解的困難。仿生算法［6］為復(fù)雜優(yōu)化問題的求解提供了有效途徑，適于解決多觀測井抽水試驗的全程數(shù)據(jù)擬合問題。利用仿生算法求解地下水逆問題方面已取得一定的研究成果［7］，其中粒子群算法因原理簡單、計算效率高，適用于解決大數(shù)據(jù)、多參數(shù)的擬合問題，但在計算過程中易出現(xiàn)局部最優(yōu)現(xiàn)象［8-9］。然而，Mahinthakumar 和Sayeed［10］發(fā)現(xiàn)對于存在多個局部極值的地下水逆問題，遺傳算法相比于其他算法更易獲得最優(yōu)解，但存在計算效率不高的缺點。因此，為保留原粒子群算法計算高效的優(yōu)勢，同時避免出現(xiàn)局部最優(yōu)現(xiàn)象，王瑞峰等［11］在粒子群算法中引入遺傳算法的雜交策略，獲得了較好的尋優(yōu)效果和較高的搜索效率。利用遺傳算法與粒子群算法優(yōu)勢結(jié)合形成的雜交粒子群算法，為多觀測井抽水試驗反演復(fù)雜水文地質(zhì)模型的含水層參數(shù)提供了有效手段。

為此，依據(jù)南京江北新區(qū)多觀測井現(xiàn)場抽水試驗，由觀測井水位理論值與實測值的方差確定觀測井權(quán)重，用以定量描述各觀測井實測數(shù)據(jù)的可靠性，并綜合運用多個觀測井抽水階段與水位恢復(fù)階段的全部觀測數(shù)據(jù)，提出考慮觀測井權(quán)重的水文地質(zhì)參數(shù)反演多線全程加權(quán)目標(biāo)函數(shù)。在此基礎(chǔ)上，利用雜交粒子群算法實現(xiàn)多線全程加權(quán)目標(biāo)函數(shù)的最優(yōu)求解，構(gòu)建基于多觀測井抽水試驗反演水文地質(zhì)參數(shù)的多線全程加權(quán)法，實現(xiàn)場地水文地質(zhì)參數(shù)的合理反演。

2 多線全程加權(quán)法

2.1 觀測井權(quán)重受外界試驗因素的影響，多觀測井抽水試驗中不同觀測井的實測數(shù)據(jù)在反映場地地下水流動特性的可靠性上存在差異。傳統(tǒng)方法進行參數(shù)反演前，常常會對不同觀測井的水位降深曲線做出一定篩選，選擇與理論曲線一致性好的降深曲線用于參數(shù)反演的計算分析，其篩選過程在本質(zhì)上就是一種確定不同觀測井?dāng)?shù)據(jù)可靠性的行為。為了定量描述觀測井?dāng)?shù)據(jù)的篩選過程，提出觀測井?dāng)?shù)據(jù)利用的權(quán)重概念。

以單觀測井全程擬合目標(biāo)函數(shù)［1］為基礎(chǔ)，得到單觀測井全程擬合方差：

式中：N 和N ′分別為觀測井抽水階段和水位恢復(fù)階段的觀測時段數(shù)； si,j和s′i,l分別為第i 個觀測井抽水階段和水位恢復(fù)階段的水位降深理論值；和分別為第i 個觀測井抽水階段和水位恢復(fù)階段的水位降深實測值； Δi為第i 個觀測井全程擬合方差。

由單觀測井全程擬合方差，確定其在多線全程加權(quán)法中的權(quán)重大小，單觀測井權(quán)重定義為：

式中：wi為第i 個觀測井權(quán)重；M 為觀測井?dāng)?shù)。

式（1）和式（2）表明，單觀測井全程擬合方差越小，水位降深實測值與理論值的一致性越好，則該觀測井?dāng)?shù)據(jù)在多線全程加權(quán)法中權(quán)重越大。權(quán)重反映了觀測井?dāng)?shù)據(jù)的可靠性，為抽水試驗數(shù)據(jù)的篩選提供了有效方法，體現(xiàn)了單個觀測井?dāng)?shù)據(jù)的可靠性對水文地質(zhì)參數(shù)反演的影響。

2.2 多線全程加權(quán)目標(biāo)函數(shù)多線全程加權(quán)目標(biāo)函數(shù)利用了抽水試驗中多個觀測井在抽水階段與水位恢復(fù)階段的實測數(shù)據(jù)，并通過觀測井權(quán)重考慮了觀測井?dāng)?shù)據(jù)的可靠性。

結(jié)合單觀測井全程擬合目標(biāo)函數(shù)［1］和觀測井權(quán)重wi，則多觀測井抽水試驗的多線全程加權(quán)目標(biāo)函數(shù)為：

式中Z 為抽水試驗多線全程加權(quán)目標(biāo)函數(shù)。

多線全程加權(quán)目標(biāo)函數(shù)利用觀測井權(quán)重考慮了觀測井?dāng)?shù)據(jù)可靠性，實現(xiàn)了多觀測井、不同試驗階段水位降深的多線全程數(shù)據(jù)擬合，奠定了多觀測井抽水試驗合理反演水文地質(zhì)參數(shù)最優(yōu)解的基礎(chǔ)。

3 多線全程加權(quán)目標(biāo)函數(shù)求解

3.1 雜交粒子群算法粒子群算法（particle swarm optimization，PSO）［12］模仿鳥群捕食現(xiàn)象，每個粒子根據(jù)自身飛行經(jīng)驗與其它粒子的飛行經(jīng)驗對自身速度與位置進行動態(tài)調(diào)整，整個粒子種群組成群體尋優(yōu)機制獲取全局最優(yōu)解。

在水文地質(zhì)參數(shù)反演計算中，每個粒子所在的位置代表一組水文地質(zhì)參數(shù)，參數(shù)個數(shù)即為搜索空間的維度。對于越流承壓含水層模型，涉及導(dǎo)水系數(shù)T、貯水系數(shù)S 和越流因數(shù)B 三個水文地質(zhì)參數(shù)，因此，可以通過在三維搜索空間中設(shè)定一定數(shù)量的粒子，實現(xiàn)搜索計算模型水文地質(zhì)參數(shù)最優(yōu)解的目的。

假設(shè)第i 個粒子自身所在位置為Xi={Ti，Si，Bi}，通過追蹤個體極值Pi={TPi，SPi，BPi}與全局極值Pg={TPg，SPg，BPg}不斷更新粒子的速度與位置，每個粒子所在的位置都是水文地質(zhì)參數(shù)反演問題的潛在解。在進行第t+1 次迭代計算時，粒子根據(jù)如下方式更新自己的速度和位置：

式中：v 為粒子速度；x 為粒子位置；ω為慣性權(quán)重；c1和c2為學(xué)習(xí)因子；rand1和rand2分別為［0，1］區(qū)間內(nèi)變化的隨機數(shù)；pik為第i 個粒子的個體極值Pi在第k 維的分量；pgk為全局極值Pg在第k 維的分量；上標(biāo)t 為迭代次數(shù)；下標(biāo)i 為第i 個粒子；下標(biāo)k 為第k 維。

為避免粒子群算法計算后期易陷入局部極值［13-14］，從而影響水文地質(zhì)參數(shù)反演的準(zhǔn)確性和合理性，引入遺傳算法中的交叉策略，提高粒子多樣性［15］。在粒子群算法的每次迭代過程中，按雜交池大小比率選取一定數(shù)量的粒子進行隨機雜交。子代粒子的位置child（x）和速度child（v）計算方式如下［11］：

式中：rand3為［0，1］區(qū)間內(nèi)變化的隨機數(shù)；parent1（x）和parent2（x）為親代粒子位置；parent1（v）和par?ent2（v）為親代粒子速度。

雜交粒子群算法因其計算高效且不易陷入局部極值的優(yōu)點，適用于多參數(shù)水文地質(zhì)模型、多觀測井抽水試驗的目標(biāo)函數(shù)求解問題，為多觀測井抽水試驗反演水文地質(zhì)參數(shù)提供了有效的求解方法。

3.2 多線全程加權(quán)法求解流程在多觀測井抽水試驗水文地質(zhì)參數(shù)反演計算中，利用雜交粒子群算法實現(xiàn)多線全程加權(quán)目標(biāo)函數(shù)最優(yōu)求解，采用迭代法保證觀測井權(quán)重計算和多線全程加權(quán)目標(biāo)函數(shù)計算中水文地質(zhì)參數(shù)的一致性。在迭代法中，假設(shè)觀測井初始權(quán)重相等，取為w0i =1 M ，迭代法的收斂標(biāo)準(zhǔn)設(shè)定為觀測井權(quán)重、多線全程加權(quán)目標(biāo)函數(shù)和水文地質(zhì)參數(shù)前后兩次迭代值的相對誤差均小于1%，以滿足水文地質(zhì)參數(shù)的精度要求。

4 多觀測井現(xiàn)場抽水試驗

圖1 鉆孔地層柱狀圖與抽水井結(jié)構(gòu)

4.1 場地水文地質(zhì)條件和抽水試驗設(shè)計抽水試驗場地位于南京江北新區(qū)核心區(qū)，屬長江漫灘地貌單元，第四系地層的賦水性、透水性各不相同，主要存在潛水和承壓水兩種類型地下水，鉆孔資料揭示的地層結(jié)構(gòu)與抽水井結(jié)構(gòu)如圖1所示。

潛水賦存于①-1 層雜填土中，該土層成份復(fù)雜、結(jié)構(gòu)松散、其透水性較好、賦水性較差。潛水水位埋深約1.5 m，主要接受大氣降水及地表水體的滲透補給，排泄以垂直蒸發(fā)和徑流方式為主。②-2 層淤泥質(zhì)粉質(zhì)黏土、②-3 層淤泥質(zhì)粉質(zhì)黏土夾薄層粉土、粉砂以及②-3a 層淤泥質(zhì)粉質(zhì)黏土與粉土、粉砂互層，與下覆含水層相比，這三個地層的透水性和給水性較差，屬弱透水層。

承壓含水層由②-4 中密粉細(xì)砂夾粉質(zhì)黏土及②-5 密實粉細(xì)砂組成，賦水性和透水性好，水量較豐富，分布較均勻，承壓水水頭埋深約1.5 m。天然條件下，承壓水的主要補給來源為地下水徑流的側(cè)向補給，排泄方式以側(cè)向徑流為主。

抽水試驗采用的抽水井與觀測井孔徑均為650 mm，管井外徑均為273 mm，（0.00 ～30.00）m為實管，管外填充黏土；管底部1.00 m 為沉淀管，用黏土球夯實，以減小下覆含水層對抽水試驗的影響。Y4、Y6、Y5 井濾管長度14 m，G1 井濾管長度9 m。

由圖1抽水井結(jié)構(gòu)可知，抽水試驗主要針對②-4 和②-5 含水層進行，以Y4、Y6、Y5 和G1 組成多觀測井抽水系統(tǒng)，其中Y4 為抽水井，Y6、Y5和G1 為觀測井，抽水井與觀測井布置如圖2所示。

圖2 抽水井與觀測井平面布置

4.2 抽水試驗結(jié)果抽水試驗的主要目的是獲得場地淺部含水層的水文地質(zhì)參數(shù)，為降水工程設(shè)計提供可靠的依據(jù)。抽水試驗以定流量80 m3/h 進行抽水，抽水試驗全程時間2880 min，其中（0～1440）min 為抽水階段，（1440～2880）min 為水位恢復(fù)階段。抽水試驗全程對觀測井水位降深進行觀測，三口觀測井獲取的水位降深數(shù)據(jù)如表1所示。

表1 抽水試驗實測水位降深數(shù)據(jù)

圖3 抽水試驗場地水文地質(zhì)模型

4.3 多觀測井抽水試驗反演水文地質(zhì)參數(shù)

4.3.1 水文地質(zhì)模型概化 現(xiàn)場勘察揭示，試驗場地含水層產(chǎn)狀近乎水平，厚度均勻，抽水試驗位置遠(yuǎn)離補給和隔水邊界。相對于②-4 層、②-5 層承壓含水層，②-2 層、②-3 層與②-3a 層屬弱透水層。抽水條件下，上覆①-1 層雜填土中的潛水可以通過弱透水層對承壓含水層進行垂向越流補給?，F(xiàn)場實測潛水含水層初始水位與承壓含水層初始水位基本相等。

假設(shè)含水層為均質(zhì)各向同性，在抽水過程中忽略弱透水層的彈性釋放和潛水水位變化，結(jié)合抽水井設(shè)置，抽水試驗場地的水文地質(zhì)計算模型可概化為越流承壓含水層完整井流問題，如圖3所示。

（1）抽水階段水位降深。Hantush 與Jacob 提出了越流承壓含水層完整井抽水的水位降深公式［16］：

式中：s（r,t）為抽水t 時刻、距離抽水井r 處的承壓含水層水位降深，m；Q 為抽水量，m3/d；T 為導(dǎo)水系數(shù)，m2/d；S 為貯水系數(shù)，無量綱；B 為越流因數(shù)，m；W（u,r/B）為第一類越流系統(tǒng)井函數(shù)；D′為弱透水層厚度，m；K ′為弱透水層滲透系數(shù)，m/d。

（2）水位恢復(fù)階段水位降深。楊建民等［2］在Hantush 公式基礎(chǔ)上，吸收了Theis 水位恢復(fù)公式的建立方法，推導(dǎo)了越流承壓含水層模型水位恢復(fù)階段剩余水位降深公式：

式中：s′( r,t′ )為水位恢復(fù)階段剩余降深，m；t′為水位恢復(fù)階段的持續(xù)時間，d；sconst為穩(wěn)定降深，m；Qn為穩(wěn)定前最后一階段流量，即穩(wěn)定流量，m3/d；u′n為穩(wěn)定流量停抽后水位恢復(fù)階段持續(xù)時間代入所得的無量綱變量。

式（8）和式（11）分別用以描述越流承壓含水層抽水階段和水位恢復(fù)階段的水位降深，為越流承壓含水層模型的多線全程加權(quán)法建立奠定了理論基礎(chǔ)。

4.3.2 雜交粒子群算法參數(shù)設(shè)定 通過matlab 軟件實現(xiàn)雜交粒子群算法對多線全程加權(quán)目標(biāo)函數(shù)的求解，綜合考慮計算精度和收斂速度設(shè)定雜交粒子群算法參數(shù)。粒子群算法參數(shù)取學(xué)習(xí)因子c1=c2=2.0［17］，慣性權(quán)重ω=0.5；雜交策略參數(shù)取雜交概率Pc=0.9，雜交池大小比率Sp=0.2［11］。針對越流承壓含水層計算模型試算，雜交粒子群算法計算參數(shù)取種群數(shù)量N=50，最大迭代次數(shù)M=100。

4.3.3 水文地質(zhì)參數(shù)反演 對于三口觀測井的抽水試驗，取觀測井權(quán)重初始值為根據(jù)多線全程加權(quán)目標(biāo)函數(shù)，經(jīng)4 次迭代后，觀測井權(quán)重、目標(biāo)函數(shù)值和水文地質(zhì)參數(shù)均達(dá)到收斂標(biāo)準(zhǔn)，計算結(jié)果如表2所示，抽水試驗實測水位降深和水位降深理論值如圖4所示。

表2 多線全程加權(quán)法迭代計算結(jié)果

圖4 多觀測井抽水試驗水位降深擬合結(jié)果

從表2給出的計算結(jié)果看，隨著迭代次數(shù)增加，觀測井權(quán)重趨于收斂，目標(biāo)函數(shù)減小，反演的水文地質(zhì)參數(shù)趨于穩(wěn)定。觀測井權(quán)重大小反映了試驗數(shù)據(jù)與計算模型理論值的一致性，為觀測井實測數(shù)據(jù)篩選評估和合理利用提供了有效方法，計算結(jié)果表明，G1 觀測井的權(quán)重最大，故其實測數(shù)據(jù)更具可靠性和可利用性。

為了對比分析，表3給出了配線法確定的水文地質(zhì)參數(shù)，可觀察到不同觀測井、同一觀測井不同抽水試驗階段的反演結(jié)果都存在明顯的差異。以導(dǎo)水系數(shù)為例，抽水階段反演結(jié)果的平均值為631.69 m2/d，極差為405.61 m2/d，達(dá)到平均值的64.21%，表明采用平均值作為反演參數(shù)存在不合理性。從單觀測井不同試驗階段的反演結(jié)果看，Y6 觀測井抽水階段的導(dǎo)水系數(shù)是水位恢復(fù)階段的3.99 倍，而G1 觀測井僅為0.83 倍，顯然G1 觀測井的數(shù)據(jù)更為合理可靠。由此可知，多觀測井抽水試驗中不同觀測井、同一觀測井不同抽水試驗階段實測數(shù)據(jù)的可靠性存在差異，采用單一觀測井?dāng)?shù)據(jù)確定場地水文地質(zhì)參數(shù)并用于評價場地地下水流動特性存在一定的局限性。因此，在多觀測井抽水試驗反演水文地質(zhì)參數(shù)時，有必要引入觀測井權(quán)重以考慮觀測井實測數(shù)據(jù)的可靠性。

從圖4中可以看出，無論是在抽水階段還是水位恢復(fù)階段Y5 和G1 觀測井的計算值與實測值均有較好的一致性，但Y6 觀測井在（0 ～500）min 抽水階段的水位降深理論值與實測值存在明顯差異，實測值隨時間的變化規(guī)律偏離了越流承壓含水層模型的地下水流動特性，表明了外界試驗因素對觀測井?dāng)?shù)據(jù)的影響，因此在多線全程加權(quán)分析中該觀測井?dāng)?shù)據(jù)的權(quán)重最小，體現(xiàn)了對抽水試驗觀測井?dāng)?shù)據(jù)的有效篩選。

4.4 反演方法性能比較雜交粒子群算法為基于多觀測井抽水試驗的多線全程加權(quán)目標(biāo)函數(shù)求解提供了有效求解方法，實現(xiàn)了水文地質(zhì)模型的多參數(shù)反演。圖5給出了水文地質(zhì)參數(shù)和目標(biāo)函數(shù)值隨雜交粒子群算法迭代次數(shù)的收斂過程，計算結(jié)果表明，50 個粒子迭代至60 次時水文地質(zhì)參數(shù)和目標(biāo)函數(shù)均趨于收斂。

為探討雜交粒子群算法（BreedPSO）的適用性和性能，開展了與粒子群算法（PSO）［12］、線性遞減慣性權(quán)重的粒子群算法（LinWPSO）［18］、基于模擬退火的粒子群算法（SimuAPSO）［19］的對比分析。4 種計算方法均取種群數(shù)量N=50，最大迭代次數(shù)M=100，學(xué)習(xí)因子c1=c2=2.0。其中PSO 算法的慣性權(quán)重ω=0.5，LinWPSO 算法的慣性權(quán)重最大值ωmax=0.9，最小值ωmin=0.2；SimuAPSO 算法的退火常數(shù)λ=0.5。計算結(jié)果如表4所示。

表4中PSO 算法、BreedPSO 算法和LinWPSO 算法的目標(biāo)函數(shù)值表明，目標(biāo)函數(shù)值的大小不僅反映擬合精度，而且影響水文地質(zhì)參數(shù)的取值，在參數(shù)反演中合理選擇擬合精度，有助于提高反演參數(shù)的可靠性。

圖6給出的目標(biāo)函數(shù)收斂過程顯示，PSO 算法、LinWPSO 算法和BreedPSO 算法的目標(biāo)函數(shù)值較小，且均趨于收斂，但BreedPSO 算法的收斂速度明顯優(yōu)于PSO 算法和LinWPSO 算法。SimuAPSO 算法收斂速度快，但目標(biāo)函數(shù)值明顯偏大，反演得到的水文地質(zhì)參數(shù)不合理，疑似陷入局部最優(yōu)解。

圖5 水文地質(zhì)參數(shù)與目標(biāo)函數(shù)隨雜交粒子群算法迭代次數(shù)的收斂過程

表4 不同算法的水文地質(zhì)參數(shù)反演結(jié)果

圖6 多線全程加權(quán)目標(biāo)函數(shù)收斂過程

對比分析表明，雜交粒子群算法具有計算效率高、全局搜索能力強、收斂速度快等優(yōu)點，適用于求解多觀測井抽水試驗反演水文地質(zhì)參數(shù)問題。

5 結(jié)論

（1）針對多觀測井抽水試驗受外界因素影響導(dǎo)致的實測數(shù)據(jù)可靠性差異，提出了觀測井權(quán)重概念，用以定量描述各觀測井實測數(shù)據(jù)的可靠性。通過單觀測井的全程擬合方差確定了觀測井權(quán)重，實現(xiàn)了對抽水試驗數(shù)據(jù)的有效篩選。

（2）基于多個觀測井的抽水階段與水位恢復(fù)階段的實測數(shù)據(jù)，提出了考慮觀測井權(quán)重的多線全程加權(quán)目標(biāo)函數(shù)，構(gòu)建了多觀測井抽水試驗水文地質(zhì)參數(shù)反演的多線全程加權(quán)法。結(jié)合實際工程案例，驗證了多線全程加權(quán)法的可行性和水文地質(zhì)參數(shù)的合理性。

（3）運用遺傳算法雜交策略改進的雜交粒子群算法，實現(xiàn)了對多觀測井抽水試驗反演水文地質(zhì)參數(shù)問題的求解。對于越流承壓含水層計算模型，雜交粒子群算法的粒子數(shù)取50、迭代次數(shù)取60 次即可滿足計算收斂要求，雜交粒子群算法為多線全程加權(quán)目標(biāo)函數(shù)最優(yōu)求解提供了有效方法。