許增光, 王　雪, 柴軍瑞, 覃　源

(西安理工大學 省部共建西北旱區(qū)生態(tài)水利國家重點實驗室, 陜西 西安 710048)

目前,潔凈水資源短缺已成為一個世界性問題,如何準確掌握地下水的運動規(guī)律,利用地下水的天然可再生性,并通過人工調(diào)控方法對地下水資源進行有效管理和科學利用已成為當前水文地質(zhì)學研究的前沿及今后的發(fā)展方向[1-3]。因滲流產(chǎn)生的滲透變形、滑坡、潰壩以及地下水開采引起地面下沉、庫區(qū)滲漏,庫岸滑坡等情形,一旦發(fā)生,往往后果嚴重,因此掌握滲流場分布規(guī)律尤為重要。滲流分析通常都是基于傳統(tǒng)的確定性模型,但實際滲流場常常受到諸多隨機性因素的影響,具有很強的不確定性。造成不確定性的因素主要包括:滲流參數(shù)、初始條件、邊界條件、源匯項等[4-5]。而滲透系數(shù)是影響滲流場的重要因素,因此獲取較為精確的滲透系數(shù)的分布規(guī)律至關(guān)重要。

滲透系數(shù)是表征含水層特性的一個關(guān)鍵要素。通過壓水試驗計算滲透系數(shù)是實際工程的一貫做法,可以獲取直接的滲透系數(shù),精度高,但同時也存在成本高及操作繁瑣的問題。加入間接測量數(shù)據(jù)電阻率來獲取滲透系數(shù),與壓水試驗得到直接的滲透系數(shù)結(jié)合,可以在一定程度上減少打井鉆孔的數(shù)量,在操作和成本上都是有益的。特別是在大范圍區(qū)域情況下,效果會更加明顯。電阻率與滲透系數(shù)之間的關(guān)系一直在被研究,有很多經(jīng)驗、半經(jīng)驗公式。電阻率可由電測深、電剖面、高密度電阻率儀測量,操作簡單且數(shù)據(jù)量大。文獻[6]利用ERT數(shù)據(jù)推求出了非均質(zhì)多孔介質(zhì)滲透系數(shù)的分布。諸多國外學者對電阻率與滲透系數(shù)之間的關(guān)系進行了深入的研究[7-18],最后皆得到了二者之間的聯(lián)系。

本文采用兩種方法,利用電阻率測量數(shù)據(jù)求滲透系數(shù)。一種是經(jīng)驗公式法,Archie公式利用電阻率計算孔隙率,Kozeny-Carman公式利用孔隙率計算滲透系數(shù);另一種是人工神經(jīng)網(wǎng)絡(luò)法,利用MATLAB訓練電阻率與滲透系數(shù)的關(guān)系。然后借助Surfer軟件畫出滲透系數(shù)分布的等值線圖。最后對這兩種方法進行對比分析。

1　方法原理

經(jīng)驗公式法是將Archie公式與Kozeny-Carman公式結(jié)合,充分利用了電阻率、孔隙率和滲透系數(shù)之間的關(guān)系。人工神經(jīng)網(wǎng)絡(luò)法是利用數(shù)學方法直接建立電阻率與滲透系數(shù)的數(shù)學關(guān)系式。

1.1　經(jīng)驗公式法

Archie公式:

ρ0=aρwφ-m

(1)

式中,ρ0和ρw分別是地層電阻率和地層流體電阻率(Ω·m),φ是介質(zhì)的孔隙率,a和m是與介質(zhì)有關(guān)的常數(shù)值。Archie公式只有在純砂巖地層或者泥層含量小于5%的地層中適用[19]。

Kozeny-Carman公式:

(2)

式中,k是滲透系數(shù)(m/s),δw是水的密度(1 000 kg/m3),g是重力加速度(9.81 m/s2),μ是運動粘滯系數(shù)(kg/m·s),d是地層巖性顆粒粒徑值(m)。

通過公式(1)求出孔隙率,代入公式(2)就可得到滲透系數(shù)。

1.2　人工神經(jīng)網(wǎng)絡(luò)

傳輸函數(shù):

a=f(wp+b)

(3)

隱含層節(jié)點的輸出:

(4)

輸出節(jié)點的計算輸出:

(5)

誤差函數(shù):

(6)

在上述公式中,p是輸入數(shù)據(jù),w和v是權(quán)重,b是閾值,t是目標數(shù)據(jù)。誤差函數(shù)分別對隱含層節(jié)點、輸出層節(jié)點和閾值求導(dǎo)。求導(dǎo)的目的是修正變量Δvlj和Δwji。根據(jù)方程(4),沿著負梯度方向修正w和b。訓練直到滿足精度要求為止。

2　工程案例

本文案例取自于文獻[8],研究區(qū)域位于德國波鴻魯爾大學南部的魯爾河谷,在凱姆納德湖西南部河岸,覆蓋面積大約占199 016 m2。這個區(qū)域有6口已經(jīng)存在的鉆井,6個滲透系數(shù)值作為實測值。用電測深法測量20個點的電阻率。已知的地質(zhì)和水文地質(zhì)資料詳見文獻[8],測量資料見表1[8]。

表1　測量資料Tab.1　Survey data

2.1　經(jīng)驗公式法

根據(jù)已有6口井的滲透系數(shù)反推常數(shù)a和m,最終求得a=0.8,m=1.4。結(jié)合式(1)(2)求出20個滲透系數(shù)。計算的滲透系數(shù)與電阻率的關(guān)系見圖1。

圖1　經(jīng)驗公式法求解的滲透系數(shù)與電阻率的關(guān)系Fig.1　Relationship between hydraulic conductivity and resistivity calculated by empirical formula method

地層因子F=ρ/ρw。圖2是計算滲透系數(shù)與地層因子的關(guān)系。

圖2　經(jīng)驗公式法計算的滲透系數(shù)與地層因子的關(guān)系Fig.2　Relationship between hydraulic conductivity and formation factor calculated by empirical formula method

由圖1和圖2可以看出用地層因子F表達滲透系數(shù)更為合適。經(jīng)驗公式法滲透系數(shù)的計算值與滲透系數(shù)的觀測值之間的關(guān)系見圖3。因為w6與ves1相差很大,可能是由于測量上的誤差或者是這一點上變異大非均質(zhì)性強,所以在圖中未畫出點。

圖3　經(jīng)驗公式法滲透系數(shù)的計算值與滲透系數(shù)的觀測值比較Fig.3　Comparison between the calculated values and the observed values of hydraulic conductivity from the empirical formula

最后滲透系數(shù)與地層因子的關(guān)系式為:

k=129.44F-2.445+1.8522

(7)

平均誤差百分比為:

(8)

式中,N是滲透系數(shù)觀測數(shù)據(jù)的個數(shù),obsi是滲透系數(shù)觀測值;compi是滲透系數(shù)計算值。

最后用6個點的測量數(shù)據(jù)和計算數(shù)據(jù)計算經(jīng)驗公式法的平均誤差百分比。PE=10.82%。除去ves1的值,計算的PE=5.9%。誤差較小,證明預(yù)測的效果是較為理想的。

2.2　人工神經(jīng)網(wǎng)絡(luò)法

利用BP神經(jīng)網(wǎng)絡(luò)對電阻率進行訓練。公式(3)是訓練的基礎(chǔ)公式,變量由輸入數(shù)據(jù)、權(quán)重和閾值組成。輸入變量為ρ和ρw,輸出變量是6口井觀測的滲透系數(shù)。訓練出合適的函數(shù),然后再求20個電阻率對應(yīng)的滲透系數(shù)。

人工神經(jīng)網(wǎng)絡(luò)法滲透系數(shù)的計算值與滲透系數(shù)的觀測值之間的關(guān)系見圖4。兩種方法得到的20個滲透系數(shù)數(shù)據(jù)見表2。

圖4　人工神經(jīng)網(wǎng)絡(luò)法滲透系數(shù)的計算值與滲透系數(shù)的觀測值比較Fig.4　Comparison between the calculated values and the observed values of hydraulic conductivity from artificial neural network

點數(shù)經(jīng)驗公式法k/(10-2m/s)人工神經(jīng)網(wǎng)絡(luò)法k/(10-2m/s)1292283220522932352544334268536525762052247196215825727194553551024326911199220122262301331824214259238152132441621018317226276182341841920322020438242

按公式(8)計算人工神經(jīng)網(wǎng)絡(luò)法的平均誤差百分比,PE=9.23%。除去ves1的值,計算的PE=3.2%。這個誤差百分比也比較小,預(yù)測效果比較理想。

2.3　地質(zhì)統(tǒng)計法

利用地質(zhì)統(tǒng)計法中的克里格插值分別對這兩種方法計算出的滲透系數(shù)分布。Surfer軟件有很強大的插值功能,并且可以繪制等值線圖。圖5和圖6分別是經(jīng)驗公式法和人工神經(jīng)網(wǎng)絡(luò)得到的滲透系數(shù)分布圖。

圖5　經(jīng)驗公式法得到的滲透系數(shù)分布Fig.5　Hydraulic conductivity distribution obtained by empirical formula method

圖6　神經(jīng)網(wǎng)絡(luò)法得到的滲透系數(shù)分布Fig.6　Hydraulic conductivity distribution obtained by artificial neural network method

圖5與圖6分布圖,證明兩種預(yù)測方法是可行的,基于實測的電阻率數(shù)據(jù)可以得到滲透系數(shù)的分布。不僅操作簡單還可以獲取大量的數(shù)據(jù),效果很理想。從中也可以看出經(jīng)驗公式法得到的滲透系數(shù)值比人工神經(jīng)網(wǎng)絡(luò)法大,而人工神經(jīng)網(wǎng)絡(luò)法分布較為均勻。在本工程案例上,人工神經(jīng)網(wǎng)絡(luò)預(yù)測的效果最佳。計算值與觀測值相差較小,其平均誤差小于經(jīng)驗公式法的平均誤差。

通過對兩種預(yù)測方法的對比分析,經(jīng)驗公式法是依據(jù)電阻率、孔隙率和滲透系數(shù)之間的關(guān)系而定,預(yù)測值比較大,這與a和m的確定和地質(zhì)條件都有關(guān)系。神經(jīng)網(wǎng)絡(luò)的預(yù)測是依賴目標數(shù)據(jù),由于觀測數(shù)據(jù)分布均勻(1.83～2.8 cm/s)而且數(shù)量較少(6個),這導(dǎo)致神經(jīng)網(wǎng)絡(luò)預(yù)測的整體效果也比較均勻。但是總體的趨勢是一致的。從平均誤差來看,本案例預(yù)測效果最佳的是人工神經(jīng)網(wǎng)絡(luò)法。

3　結(jié)　論

本文主要是采取兩種方法利用電阻率數(shù)據(jù)得到滲透系數(shù)數(shù)據(jù)。從而減少壓水試驗帶來的高成本和復(fù)雜操作問題?？梢越?jīng)濟高效地獲取滲透系數(shù)的分布。

經(jīng)驗公式法是結(jié)合Archie公式和Kozeny-Carman公式,結(jié)合前人的經(jīng)驗,理論方法簡單。但是在a和m的確定上存在著誤差,a和m的值影響著最后的計算結(jié)果。Archie公式不適合在黏土中使用,所以經(jīng)驗公式法是不適合黏土地層地質(zhì)情況。同時也不適合太復(fù)雜地層,影響滲透系數(shù)的因素不僅僅是孔隙率,影響的因素越多,預(yù)測精度就會越低。

人工神經(jīng)網(wǎng)絡(luò)法是把工程問題數(shù)學化,采用人工神經(jīng)網(wǎng)絡(luò),利用計算機龐大的計算能力訓練出電阻率與滲透系數(shù)之間的關(guān)系,但人工神經(jīng)網(wǎng)絡(luò)法不適合在沒有已知的滲透系數(shù)情況下使用。

兩種方法得到的ves1的滲透系數(shù)與測量值相差都很大,可能這一點上滲透系數(shù)的觀測值存在著誤差或者是此處的非均質(zhì)變異性較其他地方強。根據(jù)平均誤差百分比的計算,兩種方法的平均誤差都較小。根據(jù)最后的插值等值線圖,預(yù)測的結(jié)果很理想。

在這個工程案例下兩種預(yù)測方法都是可行的,而且預(yù)測的效果都是很理想的。雖然兩者之間有一些偏差,可能是由于研究區(qū)域的面積大,已知數(shù)據(jù)太少造成的,再加上測量出現(xiàn)一些誤差等情況都會導(dǎo)致結(jié)果的不準確。

參考文獻:

[1] JORDI B A. Importance of stream infiltration data for modelling surface water-groundwater interactions [J]. Journal of Hydrology, 2015, 528: 683-693.

[2] 張博, 李國秀, 程品, 等. 基于隨機理論的地下水環(huán)境風險評價 [J]. 水科學進展, 2016, 27(1): 100-106.

ZHANG Bo, LI Guoxiu, CHENG Pin, et al. Groundwater environment risk assessment based on stochastic theory [J]. Advances in Water Science, 2016, 27(1): 100-106.

[3] 王林, 徐青. 基于蒙特卡羅隨機有限元法的三維隨機滲流場研究 [J]. 巖土力學, 2014, 35(1): 287-292.

WANG Lin, XU Qing. Analysis of three dimensional random seepage field based on Monte Carlo stochastic finite element method [J]. Rock and Soil Mechanics, 2014, 35(1): 287-292.

[4] 覃榮高, 曹廣祝, 仵彥卿. 非均質(zhì)含水層中滲流與溶質(zhì)運移研究進展 [J]. 地球科學進展, 2014, 29(1): 30-41.

QIN Ronggao, CAO Guangzhu, WU Yanqing. Review of the study of groundwater flow and solute transport in heterogeneous aquifer [J]. Advances in Earth Science, 2014, 29(1): 30-41.

[5] 梁婕, 曾光明, 郭生練, 等. 滲透系數(shù)的非均質(zhì)性對地下水溶質(zhì)運移的影響 [J]. 水利學報, 2008, 39(8): 900-906.

LIANG Jie, ZENG Guangming, GUO Shenglian, et al. Effect of hydraulic conductivity heterogeneity on solute transport in groundwater [J]. Journal of Hydraulic Engineering, 2008, 39(8): 900-906.

[6] 喻永祥, 吳吉春. 利用ERT數(shù)據(jù)推求非均質(zhì)多孔介質(zhì)滲透系數(shù)初探 [J]. 水文地質(zhì)工程地質(zhì), 2006, 33(2): 41-44.

YU Yongxiang, WU Jichun. Application of electrical resistivity tomography data to estimate hydraulic conductivity in porous medium [J]. Hydrogeology and Engineering Geology, 2006, 33(2): 41-44.

[7] SUN K, GOLTZ M N. Direct estimation of hydraulic parameters relating to steady state groundwater flow [J]. Environmental Modelling &Software, 2016, 86: 50-55.

[8] NIWAS S, CELIK M. Equation estimation of porosity and hydraulic conductivity of Ruhrtal aquifer in Germany using near surface geophysics [J]. Journal of Applied Geophysics, 2012, 84(9): 77-85.

[9] OLATUNJI S, MUSA A. Estimation of aquifer hydraulic characteristics from surface geoelectrical methods: case study of the Rima Basin, north western Nigeria [J]. Arabian Journal for Science and Engineering, 2014, 39(7): 5475-5487.

[10] OKIONGBO K, MEBINE P. Estimation of aquifer hydraulic parameters from geoelectrical method-a case study of Yenagoa and environs, Southern Nigeria [J]. Arabian Journal of Geosciences, 2015, 8(8): 6085-6093.

[11] MONDAL N C, DEVI A B, RAJ P A, et al. Estimation of aquifer parameters from surfacial resistivity measurement in granitic area in Tamil Nadu [J]. Current Science, 2016, 111(3): 524-534.

[12] SINGH S, SINGH V S. Estimation of hydraulic characteristics from electrical resistivity data in coastal aquifers of southern India [J]. Journal Geological Society of India, 2016, 88(1): 77-86.

[13] TIZRO A T, VOUDOURIS K, BASAMI Y. Estimation of porosity and specific yield by application of geoelectrical method-A case study in western Iran [J]. Journal of Hydrology, 2012, 454-455(4): 160 -172.

[14] KAZAKIS N, VARGEMEZIS G, VOUDOURIS K. Estimation of hydraulic parameters in a complex porous aquifer system using geoelectrical methods [J]. Science of the Total Environment,2016,550:742 -750.

[15] PERDOMO S, AINCHIL J, KRUS E. Hydraulic parameters estimation from well logging resistivity and geoelectrical measurements [J]. Journal of Applied Geophysics, 2014, 105(2): 50-58.

[16] BOUDERBALA A, REMINI B, HAMOUDI S. Geoelectrical investigation of saline water intrusion into freshwater aquifers: A case study of Nador coastal aquifer, Tipaza, Algeria [J]. Geofísica Internacional, 2016, 55(4): 239-253.

[17] TIZRO A T, VOUDOURIS K, SALEHZADE M, et al. Hydrogeological framework and estimation of aquifer hydraulic parameters using geoelectrical data: a case study from West Iran [J]. Hydrogeology Journal, 2010, 18(4): 917-929.

[18] ASFAHANI J. Hydraulic parameters estimation by using an approach based on vertical electrical soundings (VES) in the semi-arid Khanasser valley region, Syria [J]. Journal of African Earth Sciences, 2016, 117: 196-206.

[19] 楊飛. 地球物理資料綜合解釋 [M]. 北京: 石油工業(yè)出版社, 1981: 158-159.