曾 玲 孫曉光 楊 展 趙 瑩

(1.山西省煤炭地質(zhì)水文勘查研究院，山西 030006；2.中石油煤層氣有限責(zé)任公司忻州分公司，山西 036600)

山西柿莊區(qū)塊煤層氣井產(chǎn)水量數(shù)值法預(yù)測(cè)

曾 玲1孫曉光2楊 展1趙 瑩1

(1.山西省煤炭地質(zhì)水文勘查研究院，山西 030006；2.中石油煤層氣有限責(zé)任公司忻州分公司，山西 036600)

針對(duì)沁水盆地柿莊區(qū)塊煤層氣井普遍高產(chǎn)水的現(xiàn)象，對(duì)柿莊區(qū)塊進(jìn)行了詳細(xì)的水文地質(zhì)研究。建立了柿莊區(qū)塊水文地質(zhì)概念模型：以對(duì)煤層氣井產(chǎn)水影響較大的含水層為目的層，分析出區(qū)塊東南側(cè)為補(bǔ)給區(qū)，西南和東北角為排泄區(qū)。根據(jù)概念模型建立了相應(yīng)的數(shù)學(xué)模型，利用Modflow軟件對(duì)柿莊區(qū)塊進(jìn)行了數(shù)值模擬，并對(duì)煤層氣井產(chǎn)出水量進(jìn)行預(yù)測(cè)，預(yù)測(cè)出三個(gè)高產(chǎn)水區(qū)域。

柿莊區(qū)塊 水文地質(zhì) 模擬 產(chǎn)水量

1 水文地質(zhì)概況

(1)區(qū)塊概況

柿莊區(qū)塊位于沁水盆地南部，區(qū)塊面積711.764km2。本次模擬的目的是用數(shù)值法對(duì)煤層氣井的產(chǎn)水量進(jìn)行預(yù)測(cè)。柿莊區(qū)塊主要抽采的3號(hào)、15號(hào)煤層的煤層氣，對(duì)煤層氣井產(chǎn)出水量影響較大的是煤系地層的含水層。因此，將山西組砂巖裂隙水與石炭系太原組巖溶裂隙含水層合并為一層含水層進(jìn)行分析。含水層上下之間都有較穩(wěn)定的隔水層，在無構(gòu)造發(fā)育地段基本不會(huì)與上下部含水層發(fā)生水力聯(lián)系，是一個(gè)相對(duì)獨(dú)立的含水層。根據(jù)研究區(qū)含水層的特點(diǎn)，確定了本次數(shù)值模擬的三套層系：目的含水層上部弱透水層(第四系松散層+二疊系石盒子組地層)、含水層組(山西組砂巖裂隙含水層+太原組薄層灰?guī)r裂隙含水層)、含水層下部隔水層(本溪組鋁土質(zhì)泥巖)。據(jù)收集到的鉆孔數(shù)據(jù)，山西組頂面標(biāo)高-122～486m，太原組底面標(biāo)高-321～332m，含水層平均厚度141.99m。隔水層本溪組底面標(biāo)高-342～344m,平均厚度17.79m。

圖1 模擬區(qū)平面圖

(2)邊界條件

本次模擬的研究區(qū)域位于長(zhǎng)子縣西南柿莊區(qū)塊，考慮到研究區(qū)的補(bǔ)徑排關(guān)系，將模擬區(qū)范圍在柿莊區(qū)塊原有范圍基礎(chǔ)上向東南部擴(kuò)大，東北部以南漳鎮(zhèn)為界，東南部以高平市米山鎮(zhèn)、河西鎮(zhèn)為界，西南以端氏鎮(zhèn)為界，北及西北部以方山為界，模型區(qū)形狀呈矩形，面積約2565km2，柿莊區(qū)塊是主要的滲流區(qū)(圖1)。

邊界條件表示的是研究區(qū)域的地下水系統(tǒng)與其周圍邊界環(huán)境之間的關(guān)系，周圍環(huán)境的變化對(duì)研究系統(tǒng)的影響。本模型區(qū)的外圍邊界全部為人為邊界(即第二類邊界條件)，即邊界上單位面積或單位長(zhǎng)度上的流量是已知的。

模型區(qū)的東部邊界為灰?guī)r和含水層地層淺埋區(qū)和出露區(qū)，晉獲斷裂位于此，接受大氣降水和斷層越流補(bǔ)給，地下水由東南部流入研究區(qū)，概化為第二類流量邊界；模型西北邊界為流入邊界，西南、東北邊界均為地下水流出邊界，均概化為第二類流量邊界。研究區(qū)南部的寺頭斷層為隔水?dāng)鄬?，又將研究區(qū)分為東西兩個(gè)區(qū)域(圖2)。 本模型中的第二類邊界條件，采用有流入或者流出邊界的流量，即把這個(gè)流量邊界看做是注水井或抽水井造成的，用水井子程序包(WEL)處理。

圖2 含水層邊界條件概化

(3)地下水水流特征

模型區(qū)山西組砂巖裂隙含水層和太原組灰?guī)r巖溶裂隙含水層均為承壓水，含水層總厚度較大，宏觀上概化為均質(zhì)各向同性含水層，對(duì)于同一參數(shù)分區(qū)的含水層可以看成是均質(zhì)的；地下水水位呈下降趨勢(shì)，地下水流隨時(shí)間和空間變化，故地下水流概化為非穩(wěn)定流；地下水符合達(dá)西定律，從空間上看地下水流近似水平運(yùn)動(dòng)，概化為三維流(縱向運(yùn)動(dòng)微弱)；地下水流態(tài)為均質(zhì)各向同性三維非穩(wěn)定承壓水流。

2 數(shù)學(xué)模型的建立

結(jié)合研究區(qū)的水文地質(zhì)概況，依據(jù)滲流的連續(xù)性方程和達(dá)西定律，建立了與水文地質(zhì)概念模型相對(duì)應(yīng)的二維非穩(wěn)定承壓水流地下水運(yùn)動(dòng)數(shù)學(xué)模型如下：

式中： kxx，kyy——分別為滲透系數(shù)在X和Y方向的分量，單位(LT-1)；

w——源匯項(xiàng)(m3/d)；

Ss——含水層的貯水率(L-1)；

t——時(shí)間(d)；

h——承壓水水頭(m)；

Γ2——給定流量邊界；

q——給定流量邊界上的流量；

cos(n,x),cos(n,y)——邊界上外法線方向與x，y軸的夾角的余弦。

3 數(shù)值解法

數(shù)值法的求解原理是利用分割原理來量化邊界條件以及含水層結(jié)構(gòu)，量化時(shí)要控制各節(jié)點(diǎn)的參數(shù)以及水位水量來構(gòu)造所要模擬的水文地質(zhì)模型。在Modflow中利用格點(diǎn)差分法，將含水層系統(tǒng)劃分為網(wǎng)格系統(tǒng)，將整個(gè)含水層分割為若干層、若干行、若干列，分割出來的每個(gè)小長(zhǎng)方體即為一個(gè)計(jì)算單元。每個(gè)剖分出來的小長(zhǎng)方體的中心稱為節(jié)點(diǎn)，一個(gè)計(jì)算單元格的水頭由水頭在該節(jié)點(diǎn)的值所表示。將連續(xù)函數(shù)離散成為節(jié)點(diǎn)上的函數(shù)值，將偏微分方程轉(zhuǎn)化為代數(shù)方程，通過求解代數(shù)方程來獲得求解函數(shù)。

4 數(shù)值模擬

4.1 模型離散

對(duì)模型區(qū)進(jìn)行網(wǎng)格剖分，將平面分成50行，50列，含水層一層，隔水層兩層?？偣财史至?500個(gè)矩形網(wǎng)格單元。含水層的有效計(jì)算單元2500個(gè)。含水層、隔水層的頂?shù)装甯叱虒?dǎo)入模型中，自動(dòng)生成。在柿莊區(qū)塊以及設(shè)有抽水、注水井的地區(qū)網(wǎng)格加密。地層總體由東北向西南傾斜，西南處地層埋藏最深。

4.2 參數(shù)的確定

(1)滲透系數(shù)

根據(jù)收集到的水文孔抽水試驗(yàn)數(shù)據(jù)及煤層氣井停泵后的恢復(fù)曲線計(jì)算確定含水層參數(shù)(表1)，將含水層滲透系數(shù)劃分為9個(gè)區(qū)域(圖3)。本溪組和含水層上部地層概化為弱透水含水層。

圖3 含水層滲透系數(shù)分區(qū)圖

井名TS-002TS-003TS-006TS-007TS-008SN-004SX-006SX-011SX-012滲透系數(shù)000720022200029000420003900061000250000200013分區(qū)號(hào)986571342

上述分區(qū)方法在分區(qū)范圍和數(shù)值受人為因素影響比較大。為了能夠得到相對(duì)比較準(zhǔn)確的值，把一定范圍的水文地質(zhì)參數(shù)輸入計(jì)算機(jī)模型中進(jìn)行反演，經(jīng)過不斷調(diào)試，最終得出符合實(shí)際的分區(qū)和滲透系數(shù)值。

(2)初始水頭

根據(jù)收集到的水文孔和煤層氣井，分析模型區(qū)水文地質(zhì)特征，確定含水層地下水流場(chǎng)，將含水層流場(chǎng)矢量化后導(dǎo)入模型中生成含水層初始水頭(圖4)。

(3)儲(chǔ)水系數(shù)

根據(jù)辛安泉域區(qū)域資料進(jìn)行賦值(表2)。

表2 模型中各層滲透系數(shù)及儲(chǔ)水系數(shù)

圖4 含水層初始流場(chǎng)

4.3 源匯項(xiàng)概化

模型計(jì)算區(qū)內(nèi)的源匯項(xiàng)主要是大氣降水、河流入滲補(bǔ)給和側(cè)向補(bǔ)給；礦井排水和工農(nóng)業(yè)城市用水。

本報(bào)告模型范圍內(nèi)補(bǔ)給源有兩個(gè)來源：一是主要以模型區(qū)東南部奧灰地層出露，接受大氣降水補(bǔ)給；二是接受含水層的西北側(cè)向補(bǔ)給。

(1)大氣降水補(bǔ)給項(xiàng)

地下水由模型范圍東南部接受補(bǔ)給后向西徑流，在模擬區(qū)中部形成向西南和東北不同方向進(jìn)行徑流。由于季節(jié)性因素的影響，且各月分配不均勻，補(bǔ)給量也不同，降水入滲補(bǔ)給量采用以下公式計(jì)算：

式中：Q——降水入滲補(bǔ)給量，萬m3/a；

α——降水入滲補(bǔ)給系數(shù)；

Ri——第i個(gè)月的降水量，mm/月；

F——入滲區(qū)面積，km2。

模型區(qū)東南部的降水入滲區(qū)的面積約為3km2，入滲區(qū)屬巖溶化較強(qiáng)的灰?guī)r，入滲系數(shù)取0.35。根據(jù)長(zhǎng)治縣和壺關(guān)兩縣降水量近10年的降水量，結(jié)合公式計(jì)算得到模型區(qū)平均每年降水入滲補(bǔ)給量為57.5萬m3/a，豐水期(2011年)降水入滲補(bǔ)給量為72.2萬m3/a，枯水期(2012年)降水入滲補(bǔ)給量為44.3萬m3/a，將計(jì)算結(jié)果分別加到相應(yīng)的計(jì)算單元上。

(2)側(cè)向徑流補(bǔ)給

補(bǔ)給量由以下公式得出：

Q=KJML

式中：Q——地下水側(cè)向補(bǔ)給、排泄量，m3/d；

K——計(jì)算斷面含水層平均滲透系數(shù)，m/d；

J——計(jì)算斷面平均水力坡度；

M——計(jì)算斷面含水層平均厚度，m；

L——計(jì)算斷面的長(zhǎng)度，m。

根據(jù)西北部含水層平均厚度10.11m，水力坡度1.5‰，斷面長(zhǎng)度5km，平均滲透系數(shù)131m/d，最后計(jì)算求得側(cè)向徑流量Q=9933m3/d，計(jì)算值加到模型中相應(yīng)的計(jì)算單元上。

4.4 模型的識(shí)別和校驗(yàn)

為了確保模型求解的唯一性，在模型調(diào)試過程中充分利用各種定解條件，也就是用那些靠得住的實(shí)測(cè)資料，如邊界流量、生產(chǎn)井開采量等來約束模型對(duì)原型的擬合。在模型調(diào)試過程中，還充分利用水文地質(zhì)勘探中獲得的有關(guān)信息及計(jì)算者對(duì)水文地質(zhì)條件的認(rèn)識(shí)，來約束模型的調(diào)試和識(shí)別。這樣就能使識(shí)別后的含水層模型參數(shù)、地下水流場(chǎng)及水位動(dòng)態(tài)、地下水資源量三者達(dá)到唯一的最佳匹配，使模型識(shí)別結(jié)果唯一、正確、可靠(圖5)。

圖5 水位歷史擬合

圖6 產(chǎn)水量預(yù)測(cè)分布圖

本次研究采用煤層氣孔抽水試驗(yàn)期間的水文孔水位動(dòng)態(tài)資料來進(jìn)行模型識(shí)別與驗(yàn)證，使計(jì)算和地質(zhì)條件的分析相結(jié)合，及時(shí)指導(dǎo)調(diào)參，推斷出水文孔的水位觀測(cè)值與模擬計(jì)算值基本一致。經(jīng)反復(fù)多次調(diào)參計(jì)算，四個(gè)水位觀測(cè)孔的擬合結(jié)果較為滿意，說明數(shù)學(xué)模型是正確的，數(shù)值模型是基本上是可靠的，所以此數(shù)值模型是可以用于研究區(qū)的煤層氣孔出水量的預(yù)測(cè)。

4.5 產(chǎn)水量預(yù)測(cè)及分析

在進(jìn)行模型內(nèi)涌水量預(yù)測(cè)時(shí)，由于含水層發(fā)育情況、富水性、滲透性極為復(fù)雜，在進(jìn)行預(yù)測(cè)的時(shí)候煤層氣孔涌水量與孔內(nèi)水位變化情況這兩個(gè)參數(shù)都是待求的。

在這種情況下，最好的處理方法就是根據(jù)歷史的涌水量和地下水位觀測(cè)數(shù)據(jù)建立涌水量與地下水位(或者地下水位降深)之間的數(shù)學(xué)關(guān)系。

本次模型預(yù)測(cè)的目的是為了預(yù)測(cè)柿莊區(qū)塊的產(chǎn)水量，由北向南隨機(jī)設(shè)置8個(gè)區(qū)塊來預(yù)測(cè)產(chǎn)水量。通過對(duì)8個(gè)區(qū)塊的預(yù)測(cè)(表3)，找出區(qū)塊9、區(qū)塊8、區(qū)塊5三個(gè)高產(chǎn)水風(fēng)險(xiǎn)區(qū)(圖6)，后續(xù)勘探開發(fā)過程中，應(yīng)盡量避免高產(chǎn)水風(fēng)險(xiǎn)區(qū)。

表3 區(qū)域1與其他8個(gè)區(qū)域水量交換預(yù)測(cè)結(jié)果表

[1] 丁繼紅,周德亮,馬生忠.國(guó)外地下水模擬軟件的發(fā)展現(xiàn)狀與趨勢(shì)[J].勘察科術(shù),2002,(1):37-42.

[2] 吳劍鋒,朱學(xué)愚.由MODFLOW淺談地下水流數(shù)值模擬軟件的發(fā)展趨勢(shì)[J].工程勘察.2000,(2):12-15.

[3] 賈金生,田冰,劉昌明. Visual MODFLOW在地下水模擬中的應(yīng)用[J].河北農(nóng)業(yè)大學(xué)學(xué)報(bào),2003,26(2):71-78.

[4] 武強(qiáng),董東林,武鋼,等.水資源評(píng)價(jià)的可視化專業(yè)軟件(Visual Modflow )與應(yīng)用潛力[J].水文地質(zhì)工程地質(zhì),1999,(5):21-23.

[5] 何杉. Processing MODFLOW 軟件在地下水污染防治中的應(yīng)用[J].水資源保護(hù),1999,57(3):16-18.

[6] 李平,盧文喜,馬洪云,等. Visual MODFLOW在地下水?dāng)?shù)值模擬中的應(yīng)用-以公主嶺市黃龍工業(yè)園水源地為例[J].工程勘察,2006,(3):24-27.

(責(zé)任編輯 桑逢云)

Numerical Method to Predict Water Production of CBM Well in Shizhuang Block of Shanxi Province

ZENG Ling1,SUN Xiaoguang2, YANG Chen1, ZHAO Ying1

(1.Shanxi Coal Geological Prospecting Institute of Hydology, Shanxi 030006;2. Xinzhou Branch, PetroChina Coalbed Methane Company Limited, Shanxi 036600)

In view of generally high water production of CBM wells in the Shizhuang Block of Qinshui Basin, the paper carries out the detailed study of hydrogeology in the block, and the hydrogeological concept model of Shizhuang Block is established. With greater influence on the CBM Wells produce water aquifer strata as the target layer, it is analyzed that the southeast corner of the block is the recharge area, and the southwest and northeast corner the discharge area. According to the conceptual model, the corresponding mathematical model is established. Using the Modflow software to conduct simulation in Shizhuang Block, it can forecast the water production of CBM wells and to predict the three areas with high water production.

Shizhuang Block; hydrogeology; simulation; water production

曾玲，女，碩士研究生，主要從事煤炭水文地質(zhì)研究。