徐 亞,胡立堂,劉玉強,王健媛,董 路,王 琪
(1.北京師范大學水科學研究院/地下水污染控制與修復教育部工程研究中心,北京 100875;2.中國環(huán)境科學研究院,北京 100012)
觀測井作為地下水環(huán)境動態(tài)監(jiān)測的主要手段,以及數(shù)值模型反演求取水文地質(zhì)參數(shù)和水質(zhì)運移參數(shù)的主要對象,在水文地質(zhì)學中具有非常重要的地位。地下水觀測井可分兩類:一類是測壓計式觀測孔,這類觀測孔井徑較小(小于1cm),可近似視為只有一個點與含水層接觸,能很好反映接觸點含水層的水質(zhì)和水頭情況;一類是常規(guī)觀測孔,即有限井徑觀測孔。
常規(guī)觀測孔井筒的存在,相當于在含水層中制造出一個具有很大垂向滲透系數(shù)的“管道”。在存在三維流的單層含水層系統(tǒng)中,這個“管道”內(nèi)部會產(chǎn)生垂向的水流,而井筒內(nèi)部水流的垂向流動又會導致井周水頭的重分配乃至對含水層上下部的溶質(zhì)分布產(chǎn)生擾動。對于井筒穿過二層或以上含水層的觀測井,不同含水層之間的水頭差在井管內(nèi)部產(chǎn)生垂向水流將會更加明顯,對含水層天然動態(tài)的擾動也會更加劇烈。因此,常規(guī)觀測井的觀測數(shù)據(jù)很難反映實際的地下水水頭和水質(zhì)情況。然而常規(guī)觀測孔成井工藝相對簡單、成本低廉,在很多發(fā)展中國家仍大量使用。即使在發(fā)達國家,大量歷史的水位和水質(zhì)數(shù)據(jù)都是通過觀測井獲得,因而有必要對常規(guī)觀測井的井筒效應開展“定量”研究,分析觀測數(shù)據(jù)的可靠性及適用范圍。
Bennett和Patten[1]等人首先提出了觀測井中存在的豎直流問題。Giddings[2]分析了觀測井濾管長度對水質(zhì)采樣的影響。Church和Granato[3]通過分析同一場地的測壓計式觀測井和常規(guī)觀測井數(shù)據(jù),從實踐上證實了常規(guī)觀測井可能導致的觀測誤差和對周邊流場的影響。國內(nèi)學者陳崇希等人[4~6]提出了雙重介質(zhì)模擬的等效滲透系數(shù)法,并基于該方法建立了單層抽水的觀測井-含水層耦合模型,定量分析了觀測井井筒存在對含水層天然流場分布的影響,并對井筒效應的各影響因素進行了敏感性分析。
實際上,單層抽水是最簡單和最基本的井筒-含水層系統(tǒng)。大多數(shù)情況下,含水層系統(tǒng)由兩層或兩層以上的含水層組成,含水層之間由滲透系數(shù)相對較小的弱透水層或隔水層分隔,形成多層含水層系統(tǒng)。在實際的開采作業(yè)中,抽水井通常穿過多層含水層以增大井孔涌水量,這即是混合抽水。鑒于混合抽水的普遍存在性,有必要對此種情形下觀測孔中的水位及其井筒效應進行分析。本文將借鑒陳崇希等人耦合模型模擬的基本思路和分析方法,建立有代表性的混合抽水概念模型,分析混合含水層情景下觀測井井筒對地下水天然流場的擾動效果,并對井筒效應的各影響因素進行了敏感性分析。
抽水井附近非達西流模擬思路可分為兩種[7]:第一種是用一個水流控制方程同時描述達西流和非達西流問題。最典型的方法是陳崇希等人的等效滲透系數(shù)法[8~11],基于管流的能量損失方程,將井管中的水流運動方程改寫成達西滲流方程的形式,從而將井孔和含水層統(tǒng)一到一個控制方程下。第二種思路是用不同流態(tài)方程描述不同流態(tài)區(qū)域的水流運動[12~15],需要解決流態(tài)分區(qū)以及不同分區(qū)邊界的耦合問題。
等效滲透系數(shù)方法最初由陳崇希提出,最先應用于廣西北海市龍?zhí)秲蓪踊旌铣樗囼瀳觯?]求取水文地質(zhì)參數(shù),隨后應用于干旱內(nèi)陸區(qū)地下水系統(tǒng)模擬[9~10]和水平井研究[11]中,并取得了較好的效果。其基本思路是觀測井中的管流部分視為滲透系數(shù)很大的圓柱形透鏡體,并對管流引入形式上服從達西定律的“等效滲透系數(shù)”,從而將滲流-管流模型視為含有圓柱形“透鏡體”的新的“滲流”系統(tǒng)。對于這個新“滲流”系統(tǒng),關鍵問題是如何確定井管的“等效滲透系數(shù)”。
流體力學圓柱管中水流的水頭損失:
式中:ΔH——水頭損失/m;
F——摩擦系數(shù);
l——管長/m;
d——管內(nèi)直徑/m;
u——管內(nèi)平均流速/(m·d-1);
g——重力加速度/(m·d-2)。
可以分別得到管流為層流和紊流狀態(tài)時,井管中水流速度的表達式[11]:
井管水流運動為層流時:
井管水流運動為紊流時:
式中:v——流體的運動黏度/(m2·s-1);
μ——流體的動力黏度/(Pa·s);
J——水力坡度;
γ——流體的重度/(N·m-3)。
因為對于井管而言,孔隙度n=1,因此井管內(nèi)的“滲透流速v”等于實際流速u。將式2和式3與達西定律比較分別得到管流中等效滲透系數(shù)的表達式[11]:
觀測井-含水層系統(tǒng)的水流運動可耦合為一個方程,即承壓水非穩(wěn)定運動的基本微分方程[12]:
式中:s——(x,y,z)處水頭降深/m;
W——源匯項/(s-1);
Ss——比彈性儲水系數(shù)/(m-1)。
混合抽水含水層系統(tǒng)由2個及以上的含水層組成,本文僅考慮較簡單的雙層抽水情景。典型的雙層抽水情景可概化如圖1所示。整個系統(tǒng)包含1個有限井徑的完整觀測井、2個承壓含水層與1個弱透水層。同時為模擬自然界中普遍存在的三維流,在含水層中心設置一個非完整抽水井,抽水井以定流量Q抽(注)水,從而在含水層中形成垂向分流速的水流。本文的目的就是研究加入觀測井前后,觀測井內(nèi)及其附近局部區(qū)域內(nèi)水頭和流場的變化。
以抽水井與第一層含水層頂板的交點為坐標原點,抽水井與觀測井連線方向為x坐標軸方向,豎直方向為z軸方向,水平面上x軸的垂向方向為y軸方向,建立坐標系統(tǒng)如圖1(a)所示。此時混合抽水含水層系統(tǒng)中的水頭降深可以用式(4)~(5)描述。
圖1 概念模型及坐標系統(tǒng)示意圖Fig.1 Conceptual model and coordinate system
初始時刻含水層各處降深為0,第一含水層頂板和第二含水層底板均為隔水邊界,水平方向上假設模擬時段內(nèi)水頭降深為0,因此混合抽水含水層系統(tǒng)模型的初始條件和邊界條件表示為:
式中:S0(x,y,z)——表示模擬區(qū)域初始時刻的降深/m;
Sout(t)——表示模擬時段內(nèi)模型外邊界上的降深/m。
含水層系統(tǒng)如圖1所示,由兩個含水層和1個弱透水層組成,各層均為均質(zhì)各向同性介質(zhì),其滲流參數(shù)取值采用參考文獻[14](表1)。抽水井為非完整井,濾管設置在第二含水層頂部,長度為抽水井總長度的1%。
表1 模型基本參數(shù)Table 1 Value of basic model parameters
在抽水初期,含水層水頭降深隨時間變化很明顯,隨著抽水繼續(xù)進行,一段時間后,含水層中的降深變化逐漸放緩。為此,設置初始時間步長為0.001 d,變時間步長因子為1.05,模擬時間步長為127,總模擬時間10.289d。此時含水層中水頭變化很小,降深相對變化率小于5‰/d:
式中:ΔSt——表示某一位置處的降深相對變化率/(d-1);
St——某一位置處t時刻的降深/m;
St+1——某一位置處t+1時刻的降深/m。
平面上采用不等距的結點設置,結點距隨著與抽水井的徑向距離增大而增大,x軸向總計設置49個結點,其中徑距小于40m的結點21個。垂向上,設置了46個不等距結點,其中1~21和26~46結點所在層為含水層,21~26結點所在層為滲透系數(shù)相對較小的弱透水層(隔水層)。單層結點個數(shù)為289個,單層單元個數(shù)為584個;整個模型結點數(shù)為289×46個,單元數(shù)為584×46個。
利用數(shù)值法(平面上基于伽遼金有限單元法,垂向上基于有限差分法)對上述模型進行求解,并基于Fortran平臺編制了程序代碼。圖2為距離抽水井不同位置處,水頭降深的分布。從圖中可知,在模擬時段末,距離抽水井10 000m處水頭降深為零,可推斷模型邊界處降深也應為0,可認為外邊界沒有對內(nèi)部流場產(chǎn)生影響,邊界范圍的設定是合理的。
圖2 模型頂層不同徑距處水頭降深分布Fig.2 Drawdown distribution in different distance at the top layer
2.5.1 井筒內(nèi)降深和流量分析
從圖3可以看出,觀測井井筒中不同深度處降深不同,底層降深小,頂層降深大;頂?shù)讓拥慕瞪畈铍S著井徑增大而減小,如圖3(a)中井徑為0m的觀測井,其含水層頂板和底板的降深差為6.31m,相對應的井徑為0.02、0.05、0.1和0.2m 時,降深差分別為6.19、5.23、2.18和0.15m。頂?shù)讓拥慕瞪畈铍S著徑距增大而減小,如圖3中井徑為0.1m的觀測井,位于2.18、4.47、38.53和79.98m位置時,井筒內(nèi)頂?shù)装彘g的降深差分別為2.18、1.78、0.42和0.13m。
圖3 雙層含水層中,觀測井井筒內(nèi)的水頭降深Fig.3 Drawdown in the wellbore of observation well in double-aquifer system
在含水層中存在三維水流時,同一水平位置不同高度處存在著水頭梯度,含水層底部水頭高,越往上水頭降深越大,水頭越小。若在某一位置處插入一口觀測井,相當于增大了此處的滲透系數(shù),在垂向水頭梯度的作用下,地下水從高水頭水往低水頭處流動,在井管內(nèi)形成明顯的垂向水流。根據(jù)達西定律垂向水流的流量表示為:
式中:J——相鄰井節(jié)點間的水力梯度;
Ke——管流等效滲透系數(shù)/(m·s-1),根據(jù)式(4)求得;
dw——井筒直徑。
由此計算得到的井管內(nèi)部垂直向上的流量分布如圖4所示。
圖4 雙層含水層中,觀測井井筒內(nèi)的垂向水流Fig.4 Vertical flow in the wellbore of observation well in double-aquifer system
若假定流過任意井單元上的垂向流量為Qup_e,流過該單元相鄰的下單元的垂向流量為Qlow_e。則流量經(jīng)過該單元的增加值為Δ Q。令單位長度上垂向流量的增值為,則:
式中:L——單元上下結點間的結點距。
顯然,觀測井內(nèi)單位長度上垂向流量的增值即為該深度上含水層對觀測井的水平補給量(若增值為負則為觀測井向含水層的補給值)。圖5為觀測井徑距相同時不同觀測井井徑條件下,含水層對觀測井的補給流量與深度關系曲線。
分析圖4可以得出以下結論:觀測井中的垂向水流隨著井徑增大而增大,隨著徑距增大而減小。同時分析圖5可知以中間的弱透水層頂?shù)装鍨榻?,弱透水層底板以下含水層對觀測井進行補給,弱透水層頂板以上觀測井向含水層進行排泄。補給和排泄量也與觀測井井徑大小有關,井徑越大補給量(或排泄量)越大。
2.5.2 井筒效應分析
含水層中加入高滲透性的觀測井筒后,井筒內(nèi)會產(chǎn)生明顯的垂向水流,導致井筒周圍含水層的水頭發(fā)生重分布。若定義井徑為零時,含水層中任意位置(x,y,z)處的水頭降深為S0,安裝觀測井井筒后該位置處的水頭降深為SL,則可認為S=SL-S0為井筒對含水層天然水頭分布的影響,將其定義為井筒影響下的降深畸變。圖6為模型頂層不同徑距處的降深畸變,分析模型模擬結果,可以得到降深畸變的以下規(guī)律:距離觀測井距離越遠,受井筒效應影響越小,降深畸變越小。以位于x=2.18 m處,井徑等于0.2 m的觀測井為例,從圖6中可以看出x=2.18 m,即觀測井處降深畸變最大為3.0,距離觀測井越遠井筒的影響越小,以觀測井為中心,觀測井兩側降深的畸變近似呈對稱分布(距離較遠時不對稱)。
圖5 含水層對觀測井的水平補給(排泄)Fig.5 Horizontal recharge(discharge)of aquifer to observation well
圖6 表層降深畸變沿x軸的分布Fig.6 Drawdown difference along the X axis at the surface layer
井筒效應的影響程度隨著觀測井內(nèi)徑增大而增大,井徑為0.2 m的觀測井造成的降深畸變影響最大,井徑0.02 m的觀測井最小。如圖6(a)中,位于x=2.18 m處,井徑為0.02 m的觀測井井筒中的畸變?yōu)?.06 m,而同樣位置的觀測井,井徑為0.05、0.1和0.2 m時,其內(nèi)部的畸變依次為0.522、1.996和2.99 m,明顯呈現(xiàn)隨井徑增大而增大的規(guī)律。
(1)在雙層含水層中利用不完整井抽水時,同一水平位置不同高度處的水頭降深不一樣,由此導致的垂直方向的水頭差會導致觀測井內(nèi)產(chǎn)生垂向水流,水流的流量受觀測井井徑和徑距影響,井徑越大流量越大,徑距越大流量越小。
(2)在三維含水層系統(tǒng)中,利用完整觀測井在雙層含水層中進行水質(zhì)采樣時,采樣結果更能代表下方含水層的水質(zhì),而對于上方含水層而言不具有代表性。
(3)觀測井中的垂向流量會導致井筒周邊含水層水頭的重分布,此即井筒效應。距離觀測井越近,受井筒效應影響越大。
(4)井筒效應的影響程度和影響范圍受觀測井井徑和徑距決定,井徑越大影響程度和范圍越大,徑距越大影響程度和范圍越小。
[1] BENNETT G D,PATTENE P.Constant-head pumpingtest ofa multiaquiferwellto determine characteristics ofindividualaquifers[M]. U S:Geological Survey Water-Supply Paper,1962:181-203.
[2] GIDDINGS T.What is an adequate screen length for monitoringwells:Opinion I[J]. Ground Water Monitoring Review,1987,7(2):96-97.
[3] CHURCH P E,GRANATO G E.Bias in groundwater data caused by well-bore flow in long-screen wells[J].Groundwater ,1996,34(2):262-273.
[4] 陳崇希,蔣健民,林敏,等.地下水不穩(wěn)定混合抽水的滲流-管流耦合模型及其應用[R].武漢:中國地質(zhì)大學(武漢)環(huán)境地質(zhì)研究所,1992.[CHEN C X,JIANG JM,LIN M,etal.Seepageof groundwater pumping unstable hybrid pipe coupling model and its application[R].Wuhan:Environmental Institute of Geology,China University of Geosciences(Wuhan),1992.(in Chinese)]
[5] 陳崇希,胡立堂.滲流-管流耦合模型及其應用綜述[J].水文地質(zhì)工程地質(zhì),2008,35(3):70-75.[CHEN C X,HU L T.A review of the seepage-pipe coupling model and its application[J].Hydrogeology& Engineering Geology,2008,35(3):70-75.(in Chinese)]
[6] CHEN C X,Jiao J J.Numerical simulationof pumping test in multilayer wells with non-darcian flow in the wellbore[J].Ground Water,1999,37(3):465-474.
[7] 徐亞,胡立堂,儀彪奇.井孔-含水層系統(tǒng)數(shù)值模擬方法研究進展[J].水文地質(zhì)工程地質(zhì),2011,38(4):26-29.[ XU Y,HU L T,YI B Q.Advances in numerical simulation methods of a“well-aquifer”system[J].Hydrogeology & Engineering Geology,2011,38(4):26-29.(in Chinese)]
[8] 成建梅,陳崇希.廣西北山巖溶管道-孔隙-裂隙地下水流數(shù)值模擬初探[J].水文地質(zhì)工程地質(zhì),1998,25(4):50-54.[CHENG J M,CHEN C X. GuangxiBeishan karstpipeline - pore :Groundwater Flow Numerical simulation [J].Hydrogeology & Engineering Geology,1998,25(4):50-54.(in Chinese)]
[9] 胡立堂.干旱內(nèi)陸河地區(qū)地表水和地下水集成模型及應用[J].水利學報,2008(4):410-418.[HU L T.Integrated model for surface water and groundwater in arid inland river regions and its application[J].Journal of Hydraulic Engineering,2008(4):410-418.(in Chinese)]
[10] 胡立堂,王忠靜,ROBIN,等.改進的WEAP模型在水資源管理中的應用[J].水利學報,2009(2):171-179.[HU L T,WANG Z J,ROBIN W,et al.Application ofimproved WEAP modelin water resources management[J].Journal of Hydraulic Engineering,2009(2):171-179.(in Chinese)]
[11] Chen C X,WAN J W,ZHAN H B.Theoretical and experimental studies of coupled seepage-pipe flow to a horizontal well[J].Journal of Hydrology,2003,281:163-175.
[12] EUNGYU PARK,ZHAN H B.Hydraulics of a finite-diameter horizontal well with wellbore storage and skin effect[J].Advances in Water Resources,2002,25(4):389-400.
[13] 萬軍偉,沈仲智,潘歡迎.水平井的水力特征及其解析解的適用條件[J].地球科學(中國地質(zhì)大學學報,2003(5):537-543.[WAN J W,SHEN Z Z,PAN H Y.Hydraulic performance of horizontal well and applicability of its analytical solutions[J].Earth Science (JournalofChina University of Geosciences),2003(5):537-543.(in Chinese)]
[14] HU L T,CHEN C C,CHEN X H.Simulation of groundwater flow within observation boreholesfor confined aquifers[J].Journal of Hydrology,2011,398(11/12):101-108.
[15] 徐亞,薛祥山,劉玉強,等.地下水觀測井井筒效應的多場耦合數(shù)值模擬[J].地球科學(中國地質(zhì)大學學報),2014,39(9):1349-1356.[XU Y,XUE X S,LIU Y Q,et al.A coupled Darcy-Brinkman-NS simulation model of wellbore effect of a monitor well[J].Earth Science(Journal of China University of Geosciences,2014,39(9):1349-1356.(in Chinese))]