吳寶楊,李全生,曹志國,郭 強,孔軍峰,武書泉,張宗鮮
(1.中國礦業大學(北京) 力學與建筑工程學院,北京 100083; 2.國家能源集團 煤炭開采水資源保護與利用國家重點實驗室,北京 102209; 3.北京低碳清潔能源研究院,北京 102211; 4.國家能源集團 寧夏煤業公司靈新煤礦,寧夏 靈武 751410)
煤炭開采過程中的水資源保護與利用一直是采礦行業面臨的一大難題,尤其是在西部干旱半干旱的生態脆弱礦區,由于采動引起含水層的破壞,導致地下水的流失和污染?;诖?,煤礦地下水庫儲水技術應運而生[1]。建立煤礦地下水庫,并采取相應的礦井水處理措施,利用采空區及采動裂隙區實現礦井水的存儲和分質利用,已在神東礦區得到了大規模應用。然而通過在復雜地質條件下建設煤礦地下水庫來封存高鹽礦井水的工程實踐還尚無先例[2-4]。
在煤礦地下水庫建設和運行過程中,煤層工作面重復開采引起的覆巖裂隙導通含水層,可能會對礦區地下水資源造成初次影響[5-9],其次在高濃度鹽水注入地下水庫的過程中,由于濃度梯度的作用,礦物離子將主要沿著裂隙通道進行擴散,進而可能對礦區水資源造成二次影響[10-12],這2種影響方式都與采空區覆巖中的采動裂隙是否形成貫通的滲流通道密切相關。當前,對煤層覆巖裂隙滲流通道的研究主要針對煤礦水害防治方面,且已取得了較豐富的成果。煤層在開采后,頂板巖層先后經歷“彎曲下沉-裂隙發育-破斷垮落”的動態演化過程,在應力場的不斷重分布過程中,覆巖裂隙將經歷重復的張開與閉合循環過程,導致覆巖產生持續變化的破斷和離層裂隙[13-16],巖層斷裂和離層裂隙在復雜的演化過程中不斷形成透水通道并突破含水層,進而引發覆巖裂隙網絡的涌水災害[17-18],覆巖中的垂直斷裂裂隙和離層裂隙的發育高度是決定含水層水資源流失的主要因素[19-21],而對于煤礦地下水庫,在礦井水分質利用后產生的高礦化度水注入地下水庫的過程中,圍巖中的斷裂-離層裂隙仍然是滲流主要通道[22-24]。因此,通過分析采空區覆巖裂隙分布形態及高鹽礦井水擴散規律,進而研判煤礦地下水庫建設和運行期間對臨近含水層的影響,就成為研究煤礦地下水庫高鹽礦井水封存的核心問題之一。
靈新礦位于寧夏回族自治區靈武市寧東鎮境內,北距銀川市45 km,西距靈武市40 km。靈新煤礦為生產礦井,目前礦井涌水量為450 m3/h,最大涌水量550 m3/h。靈新礦生產過程中產生的礦井水屬于高礦化度礦井水,不能直接外排地表,針對這一問題該煤礦擬計劃在井下建設一座礦井水處理廠,處理后產生的高濃鹽水,封存于煤礦地下水庫中,因此需要利用現有的采空區建設一座庫容不小于300×104m3的煤礦地下水庫。
靈新礦井田位于磁窯堡向斜的西翼,呈一東傾的單斜構造,發育有規模不大的斷層,本井田構造復雜程度屬簡單。全礦井共劃分為6個采區,其中一采區主要開采14,15,16煤,煤層采厚分別為2.78,3.18,4.28 m,煤層傾角10°左右,14煤~15煤,15煤~16煤層間距分別為20和18 m。煤層群頂板充水水源主要為延安組K2,K3,K4含水層,3個含水層距14煤的距離分別為220,120,40 m。工作面采用傾斜長壁采煤方法,后退式回采,全部采用冒落法管理頂板。
一采區北翼下組14,15,16煤現已全部開采完畢,該區域地質構造簡單,沒有大斷層,地質水文條件相對簡單。因此,靈新礦煤礦地下水庫建設擬選取一采區北翼的若干采空區(分別為14號煤層中的L1614采空區、L1814采空區,15號煤層中的L1615采空區、L1815采空區,16號煤層中的L1616采空區、L1816采空區,如圖1所示)作為地下儲水空間建設煤礦地下水庫。上述各采空區的工作面采寬180 m左右,兩工作面間隔煤柱約為25 m。

圖1 煤礦地下水庫建設區域剖面Fig.1 Sectional view of construction area of underground coal mine reservoir
研究高鹽礦井水在煤礦地下水庫中的滲流規律時,覆巖裂隙在重復開采后形成的裂隙網絡和滲流通道是關鍵。為了獲取覆巖裂隙網絡和滲流動通道特征,采用二維相似模型試驗,考慮含水層對覆巖裂隙的影響,相似模擬試驗中考慮了固-液耦合過程。同時,基于相似模擬試驗分析滲流通道特征,為分析高鹽礦井水的滲流和擴散過程提供基礎。
二維相似模擬試驗模型現場如圖2(a)所示。試驗模型架子尺寸為3.2 m×0.25 m×1.6 m(長×寬×高)。根據相似模擬試驗臺尺寸及礦井地質資料,確定模型的相似比見表1。

表1 物理模擬試驗相似比Table 1 Similarity ratio of the physical experiment

圖2 固-液耦合相似模擬Fig.2 Solid-liquid coupling physical simulation
試驗中采用了注水系統模擬含水層水壓,如圖2(b)所示,該系統包括一個裝N2的高壓罐和一個裝堿水的圓柱型容器罐。采用N2加壓作為水壓控制的動力源,通過解壓閥控制壓力。最大壓力:80 kPa,即8 m高水柱。
由于實驗架高度限制,未能把巖層從開挖煤層水平一直模擬至地表,設計模型累高為1.12 m(原型280 m),對模型未能模擬的上覆巖層厚度為0.48 m(原型120 m),應采用等效荷載方式實現,模型中未模擬巖層平均容重取1 666.7 kg/m3(原型2 500 kg/m3)。因此,模型需要施加的重力補償荷載為2.67 kPa。工作面巖層主要力學性能參數及分層厚度見表2。
溪蓀鳶尾株型整齊,花大艷麗,可用作花壇、花境布置,彌補早春開花植物種類少,以及夏季花壇用花時冷色系花不足的遺憾。作為花境種植,既可以鑲邊,也可以與其他灌木和草本花卉搭配,進行分區域、重復種植[5]。

表2 巖層主要力學性能參數及分層厚度Table 2 Mechanical parameters and thickness of strata
根據靈新煤礦巖層主要力學性能參數,結合以往配比經驗,選擇河砂為骨料,石膏作為膠結物,云母片作為分層材料以模擬巖層層理[25-27]。地層材料分為隔水層、含水層和煤層等,根據室內測定滲透系數和強度(包括抗壓和抗拉強度)等數據以及結合相似比,確定了隔水層和含水層中碎石和凡士林配比,模型分層與配比情況見表3。

表3 模型地層材料配比Table 3 Model formation material ratio
按照各分層尺寸自下而上鋪設模型,對各分成間撒云母粉,并夯實各層相似模擬材料。隔水層相似材料制作選用河砂作為骨料,碳酸鈣、石膏作為膠料,石蠟、凡士林作為隔水添加劑,該配料制作的隔水層相似材料可以實現對原型地層的低強度、大變形和抗水性的相似模擬。含水層相似材料制作選用塊石(粒徑>10 mm)、卵石(粒徑2 mm)、粗砂(粒徑<1 mm)為骨料,碳酸鈣、石膏作為膠料,通過3種骨料材料的級配可以實現對不同滲透性能含水層的相似模擬。為了防止含水層前后兩側出現滲水的情況,在鋪設含水層時,在含水層前后兩側采用了蠟封的措施,以保證含水層中的水不從含水層兩側滲出。將模型自然養護20~30 d。待模型完全干燥、定型后,拆除槽鋼,并對含水層加載水壓,即可對模型進行開挖。模型中煤層的開采,根據現場工作面開采進度進行模擬。工作面開采順序為:L1614—L1615—L1616—L1814—L1815—L1816綜采工作面。
開采14煤層L1614工作面,如圖3(a)所示,覆巖自下而上相繼垮落、斷裂、離層,覆巖破壞整體形成一梯形破壞帶。覆巖破壞帶高度發育至K4含水層底部,裂縫導通K4含水層,分別向工作面上山邊緣(圖3(a)中Ⅰ處)、中央(圖3(a)中Ⅱ處)、下山邊緣(圖3(a)中Ⅲ處)涌入3股水流。
開采15煤層L1615工作面,如圖3(b)所示,由于覆巖重復采動,破壞范圍顯著向上擴展發育。覆巖破壞帶高度發育至K4含水層與K3含水層之間巖層內,在L1615工作面上山邊緣新增2處涌水點,涌水點達到5處,分布于工作面下山邊緣(圖3(b)Ⅳ處)和上山邊緣(圖3(b)Ⅴ處),其中以下山邊緣的涌水現象最為明顯。開采16煤層L1616工作面,如圖3(c)所示,覆巖破壞范圍進一步向上擴展發育。覆巖破壞帶高度發育至K3含水層下邊界,L1616工作面上山邊緣新增2處涌水點,如圖3(c)所示。接續開采14煤層L1814工作面,如圖3(d)所示,L1814采空區上方覆巖形成獨立梯形破壞帶。破壞帶高度發育至K4含水層底部,但因上山原采空區覆巖破壞,水源徑流路徑被破壞,L1814涌水主要集中在采空區上山方向(圖3(d)中Ⅷ處),采空區下山邊緣(圖3(d)中Ⅸ處)只出現少量涌水。開采15煤層L1815工作面,如圖3(e)所示,當L1815傾向推進至L1614,L1814工作面間隔離煤柱正下方時,發生強烈的頂板來壓現象,L1814與L1614工作面隔離煤柱塑化;覆巖垮落異常發育,L1814新破壞帶與上山方向破壞帶貫通,連成一更大范圍的破壞帶。L1815工作面開采上方覆巖離層發育至K3含水層底部,L1814采空區新增涌水點及涌水量有所減少,涌水點集中在采空區上山邊緣(圖3(e)中Ⅹ處)。開采16煤層L1816工作面,如圖3(f)所示,覆巖垮落異常發育,發育高度達到150 m左右,覆巖破壞帶高度發育至K3含水層內,覆巖呈內凹型,下山裂縫連線開裂程度比上山方向大;下山采空區無涌水,上山采空區多處涌水(圖3(f)中Ⅺ處)。

圖3 不同工作面采動后裂隙及涌水模型試驗Fig.3 Modelling results of fractures and water bursting during the repeated mining of different mining faces
在所有工作面開采完成后得到的裂隙分布及涌水結果如圖4所示。其中,圖4(a)為模型試驗現場結果,圖4(b)為根據現場照片進行像素二值化處理后得到的裂隙網絡分布圖。總體上,在工作面鄰近巖層有較密集的離層裂隙,采空區上山覆巖、下山覆巖和中部有明顯的斷裂裂隙,L1614和L1616采空區上山和下山邊緣有較為明顯的涌水現象,且涌水來源與K4含水層。因此,滲流通道主要集中于上述區域,且連通K4含水層。離層裂隙和巖層斷裂裂隙擴展到了K3和K2含水層,但K3和K2含水層未出現明顯的滲流過程,表明滲流通道不通暢。圖4(b)提供了裂隙分布特征,但由于像素識別困難,不能提供更為精確的裂隙結構,尤其巖層斷裂和內部裂隙通道。因此,需要采用數值模擬進一步分析采動裂隙通道,并基于數值模擬結果進行高鹽礦井水滲流規律的分析。

圖4 開采完成后試驗現場圖及裂隙網絡Fig.4 Overall drawing of modelling experiment and the characterized fracture network
根據靈新煤礦地質資料,對于L1614~L1615和L1814~L1815工作面寬度,以及結合相似模型試驗尺寸,建立靈新煤礦采動數值模型??紤]到開采后應力重分布及變形特征,將采動地質體簡化為平面問題,基于3DEC數值模擬軟件,建立的數值計算模型如圖5(a)所示。數值模型厚度為10 m。考慮采動對工作面臨近巖層的強擾動,對工作面臨近巖層的模擬節理劃分進行了局部加密。如圖5(b)所示,為了使巖層垮落更加合理,不同巖層之間的模擬節理面進行了錯落式層疊,節理間距從L1616和L1816工作面向上,逐漸放大,其中,底層節理間距為5 m,地表最大節理間距28 m。根據實際開采步驟,從L1614~L1616,L1814~L1816依次進行開采計算,工作面開采區域如圖5(c)所示。

圖5 采動裂隙數值模擬Fig.5 Numerical simulation of mining-induced fractures
圖6(a)為分布開挖后計算得到的位移云圖,由圖6(a)可知,最大位移發生在采空區中部,對比圖6(a)和圖4(a),數值模擬得到的巖層垮落外輪廓線與相似模型試驗結果近似。圖6(b)為覆巖裂隙分布圖,密集裂隙主要分布于采空區鄰近巖層,表現為離層裂隙和部分斷裂裂隙,近地面裂隙主要為離層裂隙且擴展至K3含水層。基于數值模擬結果,后處理時設置為只顯示裂隙輪廓線,得到裂隙網絡如圖6(c)所示。對比圖6(c)和圖4(b)可知,數值模擬的裂隙網絡圖與相似模型試驗結果盡管在裂隙分布上具有一定的差異,但總體上,裂隙網絡分布區域較為接近??傮w裂隙主要分布于采空區鄰近地層、上山位置和下山位置并延伸至K3含水層。由于試驗場地高度限制,相似模擬中模型高度地表軟土層及以下小部分區域未進行模擬,模型試驗最大高度大約位于圖6(c)中的虛線位置。數值與試驗結果顯示,模型右側最大裂隙分布高度與數值模擬近似,且均出現了較大的離層裂隙。圖6(c)中的Ⅰ~Ⅺ分別對應于圖3中的涌水通道特征,可見數值模擬得到的裂隙網絡能夠全部涵蓋相似模擬試驗中的涌水通道特征點,此外,數值模擬的裂隙網絡能夠提供更多的斷裂裂隙細節。因此,數值模擬的裂隙網絡較為可靠,可以作為后續裂隙滲流過程數值分析的裂隙模型。

圖6 重復開采后的位移云圖、裂隙圖和裂隙網絡提取圖Fig.6 After repeated mining,displacement cloud diagram,fracture diagram and fracture network extraction diagram
基于圖6(c),得到裂隙網絡中單個裂隙節點坐標數據后,在GMSH軟件中,首先建立裂隙網絡線單元,隨后建立地層模型,最后進行布爾運算,建立含裂隙網絡的二維平面模型。采用三角形單元并在裂隙附近進行加密,其中裂隙為一維線單元。最終建立了采動裂隙網絡的數值計算模型,如圖7所示。

圖7 滲流數值模擬模型Fig.7 Seepage numerical simulation model
整個模型為原型的1∶1比例,考慮細砂巖與粉砂巖滲透率相對裂隙較小,且數值模擬重點關注裂隙中滲流情況,因此,細砂巖與粉砂巖巖層在數值模型中統一為一種材料組。模型中材料組分為底部泥巖、粉細砂巖、含水層K2、含水層K3、含水層K4和3種裂隙。單元尺寸最小為0.5 m,最大為30 m。裂隙單元尺寸控制在1 m以內。
根據相似模擬試驗的裂隙結果(圖4(a))和數值模擬裂隙(圖6(b))將裂隙開度大體上分成了3類,如圖8所示,紅色裂隙為開度較大的裂隙,灰色裂隙開度其次,綠色裂隙開度最小。

圖8 裂隙分類示意Fig.8 Fracture classification diagram
采用OpenGeoSys數值模擬平臺進行計算,該軟件在涉及溫度-滲流-應力-化學等多物理場耦合問題的多個領域中均有應用[28]。本次計算所采用的數值模擬參數見表4,5。采用濃鹽水模擬高礦化度礦井水,表4中鹽離子在水中擴散系數參考了文獻[29]。表5中泥巖、砂巖和含水層的孔隙率、曲度和滲透率參考文獻[30]中的參數,對于彌散系數,主要參考文獻[31]中的參數。

表4 流體物理參數Table 4 Fluid physical parameter

表5 巖層滲流數值模擬參數Table 5 Numerical simulation parameters of rock seepage
在數值模擬過程中,考慮到煤層開采后,由于采空區排水,上覆巖層裂隙通道導致含水層中水壓接近0,因而設定初始滲流壓力場為0??紤]注水壓力較大(設定最大注水壓力水頭為100 m),因此,可忽略K4含水層中水的補給對滲流場的影響。由于裂隙網絡未到達K2含水層且粉砂巖滲透率較低,也忽略K2含水層中水的補給水壓。而裂隙網絡到達了K3含水層,因此,設定K3含水層接近地表處有補給水源。
在OpenGeoSys數值模擬軟件中考慮質量守恒和壓力平衡,常規滲流無法計算注水后的流體分散區域。因此,在計算滲流場時,采用了物質追蹤的方法分析濃鹽水滲流范圍。設初始滲流場壓力為0,追蹤物質為0。數值計算模型的邊界條件如圖9所示。

圖9 滲流場與濃度場計算數值模型邊界Fig.9 Numerical model boundary of seepage field and concentration field is calculated
在數值模型中,考慮K4含水層在煤層開采時降水釋放全部水壓,因此,在K4含水層的2個邊界處設定水壓為0,但考慮地表水和深部孔隙水緩慢補給影響濃度場擴散,因此設定邊界處濃度為0。首先計算在采空區注水至100 m水頭過程中的滲流場演化過程,隨后計算50 a的濃度場變化過程。
注水數值模擬過程中,計算時間步數采用Neumann自動時間步長法,設定注水壓力達到100 m水頭需要0.5 a,最終計算了6步達到最終平衡。采動裂隙網絡中的注水滲流過程如圖10所示,其中,紅色線條代表濃鹽水沿著采空區裂隙的滲流路徑。

圖10 注水條件下滲流過程Fig.10 Seepage process under water injection condition

圖11 注水后濃度場演化Fig.11 Evolution diagram of concentration field after injection
由圖10(a)可知,在注水后的很短時間內,水流在注水點附近,而注水開始后,鹽水沿著采空區裂隙流動(圖10(b)),隨后,如圖10(c)所示,濃鹽水沿著垮落帶裂隙流動并進入上部采空區及附近裂隙,且濃鹽水進入K4含水層。由圖10(d)可知,濃鹽水在進入K4含水層之后,主要沿著含水層流動,并隨后進入部分不連通裂隙區(圖10(e))以及在采空區上方的覆巖中發生緩慢流動(圖10(f))。
綜上所述,在注水后,采動裂隙中濃鹽水先沿著采空區流進,隨后沿著上方覆巖破裂區進入K4含水層,最終在K4含水層及其鄰近覆巖斷裂裂隙中流動。
在注水1 a后,如圖11(a)所示,濃鹽水主要在裂隙中、斷裂帶中破碎巖層中有部分分布。在注水5~10 a內,如圖11(b),(c)所示,裂隙中的濃鹽水逐步向巖體中擴散,L1614上山斷裂裂隙中濃鹽水擴散較快,主要原因在于濃鹽水深入破裂裂隙中,在濃度梯度的作用下,沿裂隙水擴散。此外,L1814下山方向的K4含水層中有濃鹽水擴散,而采空區中部裂隙中濃鹽水擴散不明顯。在注水20 a后,如圖11(d)~(g)所示,采空區中部破碎巖層中濃鹽水擴散較為明顯,兩側的斷裂裂隙中濃鹽水擴撒速度較慢,而K4含水層中濃鹽水闊山較為顯著??傮w上,在注水50 a時,會有較小相對濃度的濃鹽水進入到K3含水層,主要采空區中部斷裂裂隙,濃鹽水相對濃度0.1~0.3。
(1)數值模擬的采動裂隙網絡分布與相似模擬試驗結果相似,在煤層群開采后,采空區圍巖的離層裂隙和巖層斷裂裂隙相互導通并擴展至K4含水層,采空區涌水來源于距煤層較近的K4含水層。K3含水層與采空區之間未出現明顯的滲流路勁,表明滲流通道不通暢,認為煤層群開采對K3含水層的影響不大。
(2)注水條件下,高濃度鹽水首先在裂隙采空區流動,隨后沿著覆巖斷裂裂隙進入K4含水層,并最終在K4含水層及其鄰近覆巖斷裂裂隙中流動,滲流場離K3含水層較遠。濃度場在自身重力所用下,在10 a后擴散速度減緩,在50 a內,濃度場主要在K4含水層鄰近覆巖、K4含水層和底部泥巖中擴散,K3含水層中有少量鹽水分布。