基于Schwarz交替法求解分析淺埋近距離煤層 開采覆巖結構演化規律

李 磊，李壽君

( 1. 煤炭科學研究總院，北京 100013；2. 天地科技股份有限公司 開采設計事業部，北京 100013；3. 中煤科工開采研究院有限公司，北京 100013；4. 扎賚諾爾煤業有限責任公司，內蒙古 滿洲里 021412 )

隨著我國中東部地區煤炭資源的逐漸枯竭，國家煤炭資源開采重心逐漸向西部轉移。同時，西部煤田一些礦井的上部煤層逐漸開采結束，礦井進入下部煤層開采階段，比如隆德煤礦、活雞兔礦井和布爾臺煤礦等，該類型礦井煤層埋藏較淺，上下煤層間距較小，覆蓋層為含水砂層。該類型礦井在開采過程中，不僅面臨著第一主煤層采空區積水的威脅，而且可能面臨下部煤層重復開采引起原導水裂縫帶二次發育造成的礦井災害問題，例如波及含水層造成工作面涌水潰砂等。

基于板或梁的力學分析模型，國內外專家學者已經對工作面開采覆巖結構演化規律做了大量的科學研究。1867年，德國學者舒里茲認為采空區上方巖層類似于梁結構，提出采場上覆巖層結構的“懸梁假說”[1]；蘇聯學者里特捷爾將采空區上方巖層類比于梁或板，提出了“壓力拱假說”[2-4]；錢鳴高院士于1981年提出了“砌體梁”結構模型，研究了覆巖的“O-X”型破斷規律[5-6]；1988年，宋振騏院士提出了“傳遞巖梁”的理論模型[7]；20世紀90年代，錢鳴高院士[8-9]等提出了巖層控制的“關鍵層理論”，認為巖層中存在控制巖層移動變化規律的關鍵巖層，并系統分析了采動損害現象和防止采動災害；趙德深[10]提出了上覆巖層運動的“拱 板式”和“板式”兩種力學平衡模型，研究了開采空間在覆巖中的傳播規律；姜福興[11]等將基本頂結構劃分為“類拱式”、“拱梁式”和“梁式”；黃慶 享[12-13]提出了“拱梁”結構數學模型，并研究了 厚砂土層產生破裂的判據；竇林名[14]等分析了工作面巖層的邊界差異，將覆巖運動結構劃分為OX，F和T型結構；吳立新[15]等在研究大面積開采巖層 變形破壞時提出了覆巖托板控制理論；閆少宏等提出了頂板短懸臂梁-鉸接巖梁結構，并研究了煤 壁前方破壞面積和支承壓力規律[16]；尹希文[17]提出了淺埋煤層工作面覆巖“切落體”結構模型，并 研究了切落體結構的穩定條件及失穩類型；任艷 芳[18]提出了淺埋深工作面覆巖中“懸臂梁-鉸接 巖梁”組合巖層結構，認為“鉸接巖梁”結構在 煤壁處出現滑落失穩是造成淺埋長壁工作面頂 板沿煤壁切落、形成壓架事故的根本原因。目前，大部分淺埋近距離煤層重復開采的研究是基于 板或梁的力學分析模型，但該類力學模型更適合分析單煤層開采情況。淺埋近距離煤層重復開采的大部分研究成果主要集中在采場礦壓規律方面，而對于重復開采垮落帶和導水裂縫帶二次發育高度的研究并不多[18-24]。因此，本文以淺埋近距離煤層開采工作面為研究對象，提出一種更適合于近距 離煤層開采擾動的“無限平面雙矩形孔口”理論分析模型，通過理論分析和相似材料模擬試驗研究近距離煤層采煤工作面的覆巖結構演化規律，為 淺埋近距離煤層安全開采提供理論基礎和技術保障，從而避免工作面充水潰砂等重大災害事故的發生。

1 近距離煤層重復采動應力場理論解析

近距離煤層開采過程中，沿工作面推進方向垂直剖開形成無限個平面，每個平面的邊界處受到相同的分布荷載，因此，此類問題可以看作平面應變問題。在這些平面中，上下兩個采煤工作面類似于在一個半無限面內開2個矩形孔，這屬于彈性力學的雙連域問題—“無限平面雙矩形孔口問題”。

本節采用彈性力學和復變函數理論構建近距離煤層開采擾動的“無限平面雙矩形孔口”理論模型，并對其有效性進行驗證。通過該理論模型分析近距離煤層開采覆巖應力場分布情況。

1.1 多連域問題的復變函數解法

對于彈性力學的雙連域問題，常規彈性力學解法難以求解，因此需要復變函數法進行變換求解。復變函數解法是求出雙調和方程的通解，并在不同的應力邊界條件下尋求特解的過程。在復變函數解法中，Airy應力函數U可用復數z的2個解析函數φ1( z )和θ1( z )來表示[25-27]。但是，直接求解函數φ1( z )和ψ1( z )仍然比較困難。因此，本節利用復變函數保角變換法[26]，把z平面上具有復雜形狀邊界( 如矩形、橢圓形 )的一個區域變換為另一復平面ξ 上具有簡單形狀的邊界( 如單位圓內或單位圓外 )，如圖1所示，然后把在簡單邊界下求得的解再變換回去，從而求得原問題的解。

圖1 復平面內邊界形狀變換示意 Fig. 1 Schematic diagram of boundary shape transformation in complex plane

令ψ1( z )=θ1′( z )，以z為自變量的2個解析函數φ1( z )和ψ1( z )通過映射函數 z = ω （ξ ）變成以ξ 為自變量的函數。變換關系為

利用上述變換關系，應力分量的復變函數可改寫為

位移分量的復變函數可表示為

則應力邊界條件為

式中，E為彈性模量；μ為泊松比；Xij和 Yij分別為邊界上面力的2個分量。

對于多孔問題，多連通域不可能單值地保角映射到單連通域，無法直接獲得其映射函數。因此，本文將多連通域轉化成一系列單連通域后，采用Schwarz交替法進行逐一求解。

1.2 近距離煤層開采擾動覆巖應力場理論解析方法

將近距離煤層開采簡化為平面應變問題，并抽象為無限平面雙矩形孔口問題，如圖2所示。為便于分析，設上下煤層工作面尺寸相同，均為長2a，寬2b，中間巖層高為h，水平和垂直原巖應力分別為p1和p2，在上下工作面分別建立直角坐標系x1O1y1和x2O2y2。z1和z2分別表示x1O1y1和x2O2y2坐標系下的坐標。

圖2 近距離煤層開采擾動平面模型 Fig. 2 Plane model of near-range coal mining disturbance

要將多連通域問題轉化成一系列單連通域問題，需要用到Schwarz交替法，然后針對每個單連通域問題采用柯西積分解法求解。具體的求解流程為：① 只開采上煤層，利用柯西積分法求得只有上工作面1的應力場；② 在上工作面1的應力場作用下，下工作面2的邊界必然產生多余的面力( 洞邊不為零的面力 )，而洞邊是要滿足零面力條件的，因此，在下工作面2邊界處施加與多余面力大小相等且方向相反的反面力，模擬下工作面2，此時只有下工作面2沒有上工作面1，可求得只有下工作面2的應力場；③ 在上工作面1邊界位置也會產生多余的面力；④ 反復進行上述模擬開挖過程，直到兩洞洞邊均滿足零面力條件，將每次的迭代求解結果疊加，即為最終的求解結果。

1.2.1 上煤層回采后工作面復應力函數

把z平面上矩形孔洞外部區域變換成ξ平面上單位圓外部區域，根據Schwarz-Christoffel積分公式，其映射函數為

式中，n為正奇數；R為反映洞形大小的正實數，

式中，k的值與工作面的寬長比b/a有關，需進行試算，控制方程為式( 8 )。

k值應使得ε趨近于0。R的值亦與工作面尺寸有關，用式( 9 )表示。

單位圓外包括無窮遠點在內的2個解析函數φ11( ξ1)和ψ11( ξ1)通式為

式 中，B =（ p1+ p2）4，B′ = （p2- p1）2，φ11( ξ1)和ψ11( ξ1)為ξ1平面單位圓外整個區域的解析函數，可表示為

在地應力p1和p2作用下只存在工作面1時，φ11( ξ1)和ψ11( ξ1)為

1.2.2 上煤層回采引起的下煤層工作面邊界面力

x1O1y1坐標面簡稱為z1平面，x2O2y2坐標面簡稱為z2平面。假定φ11( ξ1)和ψ11( ξ1)在ξ2平面下的形式為φ12( ξ2)和ψ12( ξ2)，φ22( ξ2)和ψ22( ξ2)在ξ1平面下的形式為φ21( ξ1)和ψ21( ξ1)。工作面2邊界點在ξ2平面下坐標為σ2，在ξ1平面下坐標為γ1。則φ11( ξ1)和ψ11( ξ1)在工作面2邊界引起的多余面力f12( σ2)的表達式為

由式( 13 )求得的多余面力 f12（σ2）是下工作面2邊界面力值，其分布的顯式表達可利用復數級數式( 14 )進行逼近，Dk是復數逼近式中 σ2的各冪次項系數。

如圖3所示，由于上煤層回采引起的下煤層工作面邊界多余面力的求解步驟為：① σ2通過映射變換t2=ω( σ2)得到對應z2平面下的點t2；② t2通過坐標平移T1=t2+c得到z1平面下的坐標T1；③ T1通過映射變換γ1=ω-1( T1)得到ξ1平面下的坐標γ1；④ 將γ1代入式( 13 )即可求得洞2周邊的多余面力值f12( σ2)；⑤ 對多余面力f12( σ2)進行復數級數逼近。

圖3 坐標點變換示意 Fig. 3 Diagram of coordinate point transformation

1.2.3 反面力作用下只存在下煤層時的復應力函數

1.2.4 應力分量迭代求解

ξ1按照式( 10 )中的坐標變換得到ξ2與ξ1的關系，ξ2=f( ξ1)。將該式代入φ22( ξ2)和ψ22( ξ2)得到φ21( ξ1)和ψ21( ξ1)。因此，Schwarz交替法的第1次迭代結果為

第1次迭代求解后，下工作面2邊界滿足零面力條件，而上工作面1在φ21( ξ1)和ψ21( ξ1)作用下仍可能存在多余面力。再次按照1.2.2節和1.2.3節，所述，進行反復迭代，直到2工作面的多余面力小到一定數值，從而最終將所有迭代結果疊加得到φ1( ξ1)和ψ1( ξ1)。

1.2.5 應力分量

z平面上的應力分量組合式，經變換以ξ為自變量的復變函數表達式，利用式( 3 )可得：

1.3 理論模型有效性分析

通過對比分析雙矩形孔口問題的理論計算和數值模擬結果，驗證該理論方法的可靠性。模型的相關參數見表1。

表1 參數設置 Table 1 Parameters setting

1.3.1 映射函數確定的合理性

利用式( 8 )，當k=0.069時，ε ≈0。此時的映射函數表達式為

根據該映射函數表達式分析可得，映射后的矩形除了在4個邊角處微翹外，其余各處與原矩形的形狀尺寸基本一致，說明所選映射函數是合理的。

1.3.2 零面力條件驗證

無支護情況下，工作面邊界應滿足零面力條件，Schwarz交替法求解本身就是一個逐步消除多余面力的過程。經反復迭代后，在反面力作用下工作面2周邊正應力和剪應力分布如圖4所示。從圖4可以看出，上煤層1回采后下煤層2邊界多余面力分量雖然未全部為0，但這些面力分量均在零值附近小幅度波動，這說明Schwarz交替求解法基本滿足零面力條件，該理論方法滿足工程計算的要求。

圖4 上煤層1回采后下煤層2邊界處多余的面力 Fig. 4 Excess surface force at boundary of coal seam 2 after mining in coal seam 1

1.3.3 應力理論解與數值模擬結果比較

如圖5所示，本文采用ABAQUS軟件建立平面雙矩形孔口計算模型，模型尺寸為400 m×400 m。上下工作面開挖尺寸均為100 m×4 m，中間巖層厚度為20 m。單元屬性為彈性，單元類型為CPE4R。模型的應力邊界條件為：上下邊界施加3 MPa均布載荷，左右兩側施加2 MPa均布載荷。然后在上工作面的右側和上側選取2條應力分量監測線，從而對比分析數值模擬和理論計算2種方法獲得的應力分量分布規律的差異性。

圖5 有限元模型及監測點 Fig. 5 Finite element model and monitoring method

圖6是2條監測線上的數值模擬和理論計算的應力分量分布規律。橫坐標是距洞邊的距離。對于測線1，σy的理論值和模擬值變化趨勢均逐漸減小，但在工作面及其附近位置，理論值小于模擬值。σx的變化趨勢是先增大后減小，模擬值的變化幅度明顯。對于測線2，σx的變化趨勢仍然是先增大后減小。σy的趨勢則為逐漸增大，理論值普遍略大于模擬值。綜上所述，理論計算和數值模擬的σx和σy的變化趨勢和數值大致相同，相似程度高，說明在彈性條件下該理論方法求解雙矩形孔口的應力分布問題的可靠性較好。

圖6 2條測線上應力分布規律對比 Fig. 6 Comparison of stress distribution laws on two measuring lines

2 單一煤層開采擾動覆巖應力場及破壞區分析

通過理論解析方法計算近距離煤層開采覆巖應力場分布后，根據不同位置覆巖所受的應力狀態分析覆巖破壞范圍及其演化過程。判斷標準為：垮落帶為采空區上方拉張破壞區及局部拉張區，裂縫帶為垮落帶上部的拉裂隙區。本節將垮落帶和裂縫帶作為覆巖破壞范圍進行分析。

由于彈性力學理論求解近距離煤層開采擾動覆巖應力場問題的復雜性，求解不同力學性質和垂直各向異性的層狀煤巖層力學模型變得非常困難，甚至無法得到理論解析解。因此，將實際的層狀煤巖層地質模型簡化為均勻各向同性的彈性體模型。單煤層開采理論計算的相關力學參數見表2，分析其應力場分布及其破壞范圍。

基于Mohr-Coulomb準則的推導式( 19 )，建立覆巖塑性極限平衡區計算方法，并以此作為判斷覆巖的破壞區。當 0f ＜ 時，圍巖應力狀態處于抗剪強度包線下方，此時圍巖未破壞；當 =0f 時，圍巖應力狀態在包線上，處于極限平衡狀態；當 0f ＞ 時，圍巖應力狀態超出包線，圍巖發生破壞。

表2 理論模型的參數設置 Table 2 Parameter setting of mechanical model

2.1 單一煤層開采覆巖應力場分布規律

圖7為單一煤層開采覆巖應力分量σx和σy的分布云圖。在圖7( a )中，覆巖σx整體分布大致呈“X”狀，頂底板處于拉應力狀態，呈拱形分布，隨著遠離頂底板，σx由拉應力變為壓應力。左右兩側應力呈橢圓狀分布，存在壓應力集中區。在圖7( b )中，工作面附近σy應力場分為3個區域：頂底板鐘型分布區、左右兩側耳狀分布區、外部原巖應力分布區。頂底板、左右兩側均為壓應力，左右兩側耳狀分布區壓應力值要普遍大于頂底板鐘型分布區，越靠近工作面的左右兩側，壓應力值越大。

2.2 單一煤層開采覆巖破壞區分布

由圖7可知，隨著采煤工作面開挖，工作面前方的覆巖破壞范圍也隨著擴大。采空區上方巖層處于拉應力狀態，極值拉應力為0.56 MPa，上覆巖層破壞范圍達到32 m。

3 近距離煤層開采擾動覆巖應力場及破壞區分析

基于上述理論力學模型，分析在不同中間巖層厚度( h=12，20，36 m )條件下覆巖應力場、破壞區，以 及其與中間巖層厚度和力學參數的關系。除中間巖層厚度不同外，近距離煤層開采擾動分析的相關參數見表2。

圖7 單一煤層開采覆巖應力場云圖 Fig. 7 Overburden stress field nephogram in single coal seam mining

3.1 覆巖應力場隨中間巖層厚度的演化規律

中間巖層h=20 m時，重復開采覆巖內水平和垂直應力場如圖8( c )，( d )所示。σx整體呈現“X”型分布，左右兩側呈橢圓狀分布，上下頂板附近呈橄欖狀分布。與單一煤層不同的是，左右兩側存在的4個集中區分布于遠離中間巖層一側角點位置。σy呈中心對稱分布，頂底板附近呈現鐘型分布，左右兩側呈耳狀分布。工作面兩側存在應力集中區，兩煤層集中區相互連接。

圖8( a )，( c )，( e )為不同中間巖層厚度下σx云圖。隨著中間巖層厚度的增大，“X”型分布區域有所變化，上煤層區域的水平應力σx的拉應力范圍縮??；當h＞20 m時，X型區域出現斷連，橢圓狀分布區域和橄欖狀分布區域略有增大。中間巖層左右兩側的壓應力逐漸減小，逐漸趨于中間巖層內部應力值。圖8( b )，( d ) ，( f )為不同中間巖層厚度下σy云圖。隨著中間巖層厚度的增加，頂底板附近壓應力區沒有明顯變化。

圖8 不同中間巖層厚度下應力演化規律 Fig. 8 Evolution law of stress field with different intermediate rock thickness

3.2 覆巖破壞區隨中間巖層厚度的演化規律

由圖8可知不同中間巖層厚度( 12，20，36 m )下破壞區演化規律，當中間巖層為12～20 m時，上下2工作面破壞區向外發育并交匯在一起，形成大面積貫通性破壞區，此時中間巖層垮落，并引起上煤層工作面原有破壞區的增大。當中間巖層為36 m時，破壞區有一定發育，但并未貫穿中間巖層，說明隨著中間巖層厚度增大，2工作面相互作用導致破壞效果受到限制，中間巖層的支撐作用增強。

中間巖層h=20 m時，隨著采煤工作面開挖，工作面前方的覆巖破壞區也隨著擴大。采空區上方巖層處于拉應力狀態，極值拉應力為0.82 MPa，上覆巖層破壞范圍達到58 m。與上煤層開采引起的32 m裂縫帶破壞高度相比，破壞區高度明顯增大。說明下工作面采動影響范圍到達上工作面的煤層采空區，對已破壞覆巖產生二次影響，上工作面原破壞區發生二次發育。

4 近距離煤層重復開采覆巖裂隙演化物理模擬

為了驗證理論解析方法對分析近距離煤層上覆巖層破壞區域范圍問題的可靠性，對靈露煤礦 Ⅱ2-1和Ⅱ3近距離煤層開采進行相似材料模擬試驗，并對比分析重復采動覆巖應力場解析方法的破壞范圍與相似材料模擬試驗結果，從而證明理論解析法的有效性。

4.1 模型建立

以綜合鉆孔柱狀圖作為參考依據，模擬分析 Ⅱ2-1，Ⅱ3兩層煤開采條件下上覆巖層的變形破壞情況( 包括覆巖垮落帶和導水裂縫帶高度 )。其中，上煤層Ⅱ2-1的模擬開采厚度為4 m，下煤層Ⅱ3模擬開采厚度為4 m，煤層間距20 m。鋪設模型高度為139 m，用等效荷載代替74.5 m的巖層( 第四紀和泥巖 )厚度。Ⅱ2-1和Ⅱ3近距離煤層的地質條件和煤層厚度與中間巖層h=20 m的理論分析模型相近。

根據靈露煤礦Ⅱ2-1、Ⅱ3煤層采煤工作面地質條件與開采技術條件，相似模型的幾何相似比設置為CL=180∶1，容重比為Cγ=1.6∶1，則其應力比為Cσ=288∶1；煤層工作面開采時間與相似模型開挖時間的時間相似常數為：。靈露煤礦真實的煤巖層巖性和厚度及力學性質見表3。

表3 靈露煤礦煤巖層的力學性質 Table 3 Mechanical properties of coal and rock strata

4.2 覆巖垮落帶和裂縫帶特征

上煤層開采破壞影響區域：當上煤層開挖空間達到42 m時，上覆巖層基本頂巖層發生大范圍垮落，垮落高度達到7 m，同時，工作面上方的巖層斷裂角很明顯，認定42 m為工作面的初次來壓。當開挖空間達到63 m時，工作面上覆頂板再次垮落，判定為第2次頂板的周期來壓步距( 13 m )。如圖9( a )所示，垮落高度達到13 m，裂隙富集區主要集中在圖中的橢圓形內，裂隙發育高度32 m。當上煤層開挖空間達到99 m時，工作面上方的巖層垮落角為50°，如圖9( b )所示，垮落高度達到12 m，裂隙發育高度35 m。

下煤層開采破壞影響區域：當下煤層工作面推進到45 m時，如圖10( a )所示，上覆巖層離層頂板垮落，工作面上方的巖層垮落角為55°，垮落高度發育至16 m，中間巖層未完全破壞。當下煤層工作面推進至終采線位置，如圖10( b )所示，上覆巖層垮落高度發育至20 m，即整個中間巖層發生破壞。下煤層開采引起的裂隙富集區主要集中在圖中的橢圓形內，與上煤層開采引起的35 m裂縫帶高度相比，導水裂縫帶發育高度明顯增大，達到61 m。Ⅱ3煤層頂板平均周期來壓步距為12 m。說明下煤層采動影響范圍達到Ⅱ2-1煤層采空區，對已破壞覆巖產生二次影響，上覆巖層裂縫帶發生二次發育。重復采動覆巖應力場解析方法的破壞區求解結果與相似材料模擬試驗結果相差不超過10%，滿足工程應用要求。

圖9 上部煤層開采覆巖破壞情況 Fig. 9 Fractured zone of overlying strata of upper coal seam

圖10 下部煤層開采覆巖破壞情況 Fig. 10 Fractured zone of overlying strata of lower coal seam

5 結 論

( 1 ) 針對我國西部淺埋近距離煤層重復開采問題，提出一種新的“無限平面雙矩形孔口”理論解析模型，并采用Schwarz交替求解法推導了近距離煤層開采擾動覆巖應力場理論計算公式。該理論模型的建立和求解為求解分析近距離煤?層開采擾動的覆巖結構演化規律提供了理論基礎。

( 2 ) 近距離煤層重復開采覆巖應力場分布特征：① 近距離煤層重復開采后，煤巖層內的水平應力σx呈現“X”型分布，頂底板均處于拉應力狀態，采空區兩側巖體應力呈橢圓狀分布；② 隨著厚度的增大，中間巖層內應力集中范圍和程度降低，下工作面開采引起的覆巖破壞區域二次發育情況也不同。

( 3 ) 近距離煤層重復開采覆巖破壞范圍分布規律：當中間巖層厚度達到一定程度，下工作面開采并不會引起原覆巖破壞高度的明顯增大。當中間巖層厚度h=20 m時，下工作面開采后上覆巖層破壞區貫穿整個中間巖層，采空區上方巖層處于拉應力狀態，下工作面開采引起的原裂縫帶發育高度明顯增大，下煤層工作面采動引起了上工作面原破壞區域 的二次發育。