頁巖油立體開發誘發的儲層地應力動態響應特征

郭旭洋 金衍 林伯韜

1.中國石油大學(北京)油氣資源與探測國家重點實驗室；2.中國石油大學(北京)石油工程學院

在我國原油對外依存度逐年攀升的背景下，頁巖油大規模建產已經成為保證能源安全的一項重要手段，鄂爾多斯盆地、準噶爾盆地等都已成為我國頁巖油的主要產區［1］。由于頁巖基質滲透率低，自然產能無法滿足工業化生產需求，常通過水平井和水力壓裂等手段進行商業化開采［2-4］。Eagle Ford、Bakken 等北美地區的頁巖油開采已有一些成熟經驗［5-6］，而我國頁巖油具有分布面積和規模相對較小、但累計厚度大的特征。因此，北美頁巖油開采更傾向于在同一埋深內進行橫向上的大規模開發，而我國頁巖油開采更傾向于在水平面上較小的范圍內充分進行立體開發。頁巖油立體開發的關鍵是盡可能增加壓裂縫網與低滲儲層的接觸面積，增加流體的流動能力。立體開發需要在縱向上進行多層水平井布井［7］，使得垂向上水平井井間距較小，投產后會造成儲層壓降，進而通過壓力波誘發儲層地應力的復雜演化。儲層地應力對立體開發中投產井的鄰近井位的鉆完井、壓裂都有直接影響。因此，壓降誘發的儲層地應力響應特征就成為頁巖油高效立體開發需要研究的一個關鍵問題。

金衍等［8］、鄧金根等［9］認為綜合水壓致裂法、Kaiser 效應法、差應變法、多級子測井方法可以較好確定頁巖層的原始地應力。在此基礎上，常可通過建立數學模型的方法對開采過程中的地應力動態演化進行分析。由Biot 等［10］創立并發展的孔隙彈性力學理論是描述儲層壓降誘發的地應力演化這一物理過程的基礎，該理論闡述了一維至三維空間下飽和多孔介質中總應力、有效應力和孔隙壓力的關系。據此，科研人員針對頁巖特征不斷完善相關方法和理論，研究了低滲基質中的流體流動特征、吸附作用、擴散現象、應力敏感性、孔隙結構等相關因素的作用，考慮了基于動量平衡的頁巖骨架變形效應［11-13］。Gupta 等［14］通過水平面內地應力演化計算結果確定了最優井工廠加密井井距，并指導了重復壓裂。郭旭洋等［15］認為，生產壓降誘導的儲層地應力演化可能導致井工廠的水平井壓裂出現復雜裂縫形態，進而造成裂縫串通，干擾鄰近水平井生產，影響Eagle Ford 頁巖油井工廠的開發效率。Roussel 等［16］則提出了一種有限元流固耦合模型，表征儲層二維平面內生產壓差作用下的孔壓、地應力變化規律，并通過水力壓裂模擬證明了儲層開采誘發的地應力演化會顯著改變同一井工廠內鄰近水平井的裂縫密切割形態。

現有成熟的頁巖層地應力表征方法僅適用于未經開發的原始狀態，進入開采狀態后的地應力演化過程一般通過數學模型進行計算。目前已有的數模研究主要針對水平面內橫向尺度情況，較少針對立體情況下縱向延展較大的頁巖油儲層。然而我國頁巖油儲層累計厚度大且采用立體開發的情況較多，因此針對性地結合現場數據建立三維計算模型，開展立體開發過程中的儲層地應力響應特征研究。

1 儲層地應力響應數學模型

生產壓差造成的壓降誘發頁巖油儲層中地應力演化涉及到2 個力學問題：滲流力學和固體力學。其中，滲流力學描述水平井生產過程中地下流體在頁巖中的流動過程；固體力學描述流體流動造成的巖石骨架變形過程。傳統意義上的研究局限于單獨研究壓力場(滲流力學問題)或單獨研究應力場(固體力學問題)，因而無法表征生產誘發地應力演化這一壓力場、應力場協同作用的物理過程。

Biot 孔隙彈性理論［11］指出，多孔介質的孔隙壓力和應力狀態可以通過Biot 系數和有效應力概念聯系起來，建立了應力作用下巖石變形與流體壓力的關系為

式中，σij′為有效應力，Pa；σij為總應力，Pa；α為Biot 系數；δij為狄拉克德爾塔函數；p為孔隙壓力，Pa；Kd為巖石排水模量，Pa；Ks為巖石固相體積模量，Pa。式（1）和式（2）說明了孔隙壓力和應力的求解是對模型展開研究的關鍵參數。

為了表征儲層開采造成的壓力場演化，基于質量守恒建立的流動方程為［17］

式中，mi為i相(油或水)在體積元內的運移質量，kg/m3；fi為i相在面積元上的流動質量，kg/(m2· s)；n為單位法向向量；t為時間，s；qi為水平井設計產量，kg/(m3· s)；Ω為儲層內方程的域；Γ為方程求解區域的邊界。式(3)表征了質量守恒過程中質量累計項、質量流動項、匯/源項的平衡。流動項可進一步表示為

式中，ρi為i相密度，kg/m3；k為滲透率，μm2；μi為i相黏度，Pa· s ；p為孔隙壓力，Pa。

為了表征頁巖巖石變形造成的應力場演化，采用準靜態假設，根據動量平衡建立應力平衡方程為

式中，σ為應力張量；ρb為固體密度，kg/m3。

式(5)描述了三維空間內各方向是主應力和切應力的分量與外界作用的應力平衡狀態。根據無窮小變換假設，位移和應變的關系表示為

式中，ε為應變；?s為對稱梯度；u為位移，m。

得到孔隙介質流動和巖石骨架應力平衡方程后，需要通過Biot 孔隙彈性理論將2 個物理場進行流固耦合，具體通過孔隙度變化實現

式中，φ為孔隙度；εv為體積應變。

綜合式(3)至式(7)則可以求解壓降誘發的三維地應力演化。對時間采用有限差分方法進行隱式歐拉法求解，以確保數值解的準確性和收斂性，對空間離散采用伽遼金方法(Galerkin Method)進行有限元分析，并采用Newton-Raphson 法處理系統的非線性特征。頁巖油儲層水力裂縫在三維模型中采用局部網格加密(Local Grid Refinement)和儲層改造體積(Stimulated Reservoir Volume)的方法進行表征。

2 三維建模

F 油田頁巖油儲層厚度較大，適合立體開發手段開采。研究區域的流體和儲層物性見表1。前期已投產分段分簇水力壓裂水平井的某段裂縫微地震監測數據如圖1 所示，綠色線段代表水平井段，紫色代表共計13 個微地震事件及其位置，灰色輪廓代表該段壓裂的儲層改造監測區域(縫網)，作為對單段壓裂規模的測算，對應的裂縫幾何形態見表2。重點研究儲層三維空間內的孔壓和地應力演化。

表1 模型流體、儲層參數Table 1 Model fluid and reservoir parameters

圖1 某段壓裂裂縫微地震數據Fig.1 Microseismic data of one certain hydraulic fracture

表2 單段裂縫幾何形態Table 2 Geometry of single fracture

根據上述儲層、裂縫數據，建立如圖2 所示三維數值計算模型，研究水平井一段裂縫投產后的儲層壓降和地應力演化，重點研究壓力場、應力場在縱向隨時間發生立體演化的過程。三維建模中，通過2 條離散裂縫及SRV(儲層改造區域15 m×180 m×15 m)，對圖1 中單段水力壓裂裂縫進行表征。在縱向上，研究的儲層厚度為縫網高度的3 倍左右，以表征地應力響應在縱向尺度上的變化特征。SRV 位于模型正中心，裂縫高度在模型中z方向位置是15 m至30 m。模擬縫網生產時，使用固定井底流壓5 MPa，生產時長為2 年。通過模擬，可以獲得三維空間內孔壓、地應力的演化特征，進而進行分析。

圖2 考慮裂縫的三維計算模型Fig.2 3D calculation model considering fractures

3 計算結果與分析

縫網生產2 年后三維模型不同深度層的孔隙壓力分布情況如圖3 所示，5 個深度層裂縫均未穿透。第1 層(a)深度為0～5 m，孔隙壓力分布為12.3～13.0 MPa；第2 層(b)深度為5～10 m，孔隙壓力分布為11.8～13.0 MPa；第3 層(c)深度為10～15 m，孔隙壓力分布為9.8～13.0 MPa；第4 層(d)深度為15～20 m，孔隙壓力分布為5.1～13.0 MPa；第5 層(e)深度為20～25 m，孔隙壓力分布為5.0～13.0 MPa。(a)～(e)層對應的深度為三維模型中的z方向對應的高度，而不是埋深。

每個深度層平面內初始最大主應力方向與裂縫方向相同，生產誘發儲層壓降后引起的最大主應力方向發生變化，與裂縫方向產生的夾角用 θp表示為

式中，θp為主應力方向變化角度，rad；τxy為剪應力，Pa；σx為x方向主應力，Pa；σy為y方向主應力，Pa。σx、σy和 τxy均因為縫網生產而發生變化，各個深度層內主應力角度也會因而發生變化。

圖3 不同深度層孔隙壓力分布Fig.3 Pore pressure distribution in different layers

圖3 中，遠離SRV 的區域，孔壓變化幾乎可以忽略，這是由基質極低的滲透率阻礙了生產井壓降向儲層遠處傳播而造成的。這些區域的最大主應力轉向也幾乎可以忽略，均為指向初始最大主應力方向，這是由于壓降的壓力波幾乎未傳播到這些區域。各層的SRV 均可觀測到孔壓和最大主應力的角度變化。(c)～(e)層由于被水力裂縫穿過，造成的壓降和主應力轉向較為明顯，而且這種趨勢越接近水力裂縫中間層越明顯，并在(e)層達到最明顯。(a)、(b)層由于未被水力裂縫穿過，雖然可以觀測到一定的壓降和最大主應力轉向，但是與(c)～(e)相比較為輕微。需要注意的是，根據理想均質假設，主應力轉向數值計算結果均以裂縫為中心呈左右對稱，但由于繪圖采集的樣本點不以裂縫左右對稱，使得圖中主應力轉向沒有嚴格按裂縫左右對稱。

圖3 中各個深度層最大主應力 σH的分布情況如圖4 所示。受立體開發影響，σH的變化主要集中在SRV 區域以及SRV 沿y方向延伸數10 m 的區域。σH演化最明顯的區域是SRV 中的(e)層，計算得到(e)層的 σH最大值為23 MPa。類似地，如圖5 所示是各個深度層最小主應力 σh的分布情況。σh的變化也是主要集中在SRV 區域內，在(e)層 σh最低低至20 MPa。在SRV 區域外，在x方向和y方向均有σh的演化。(a)～(b)層由于未被水力裂縫穿過，壓降不明顯，但(c)～(e)層壓降造成的應力擾動較明顯地傳遞到(a)～(b)層，因此造成少部分區域最小主應力在生產2 年后上升的現象。

圖4 生產2 年后各個深度層最大主應力分布Fig.4 Distribution of the maximum principal stress in each formation at different depth after 2 years’ production

圖5 生產2 年后各個深度層最小主應力分布Fig.5 Distribution of the minimum principal stress in each formation at different depth after 2 years’ production

如圖6 所示為SRV 內和SRV 外 σx和 σy隨時間的變化，總體上呈現隨時間遞減的規律，這主要是由于孔壓下降造成的：由于總應力由有效應力和孔壓組成，孔壓下降會造成總應力顯著下降。但SRV 外σx在(a)～(c)層初始階段呈現上升趨勢，這是由于初始階段這些層位的壓降誘發的有效應力上升無法被有限的壓降抵消。此外，越接近SRV 中間層 σx和σy隨時間的演化就越明顯，這也可以說明SRV 對立體開發造成的儲層應力動態演化具有增強作用。

圖6 σx 和σy 時間變化曲線Fig.6 Time dependent curves of σx and σy

圖3 至圖6 的結果表明，立體開發造成的儲層地應力大小和方向的演化與縱向深度有顯著的關系。SRV 穿過的層位中地應力變化的響應更為明顯，SRV 未穿過的層位中雖然地應力大小和方向的變化有限，但是依然造成了一定的地應力場改變，仍可能影響到立體開發井網中縱向上鄰近井位的鉆完井和壓裂效果。

4 層間應力干擾分析

立體開發過程中水平井分段分簇裂縫造成的不同層位中的應力干擾是頁巖油現場開發關注的一個重點，可以通過研究縱向尺度上的水平最小主應力σh和地層破裂壓力pf來表征。破裂壓力表達式為

式中，pf為破裂壓力，Pa；σh為最小主應力，Pa；σH為最大主應力，Pa；St為巖石抗拉強度，Pa。在層間應力干擾分析中，重點研究經過2 年立體開發后相關變量與初始值對比的變化量Δσh和 Δpf。

如圖7 所示為水平井所處x-z垂直平面內的Δσh分布。立體空間內，該垂直平面與水平井井筒在同一平面內。結果顯示分段分簇裂縫生產導致的層間最小主應力干擾可以延伸至SRV 的厚度以外，直至儲層的頂端和底端，并在這些區域造成最高可達1 MPa 的應力值增加。在SRV 內，由于壓降較為明顯，導致最小主應力值可降低最多3 MPa。立體維度下，這種層間干擾可以擴散至SRV 以外15 m 處。這說明如果在圖示水平井位置的深淺±20 m (總厚度40 m)以內布水平井，均可能受到圖示水平井造成的應力干擾。

圖7 x-z垂面內不同層位間最小主應力擾動情況Fig.7 Disturbance of minimum principal stress between different horizons in the x-z vertical plane

如圖8 展示，水平井所處x-z垂直平面內的Δpf分布。圖示水平井對立體層位造成的地層破裂壓力干擾主要集中在與SRV 同寬度的范圍內，在厚度上則波及至儲層頂端和底端，造成了地層破裂壓力增加。這意味著在圖示水平井深淺±20 m 且與SRV 同寬度的區域內進行立體開發時，造成地層拉伸破壞所需的井內流體壓力升高，且最高可升高4 MPa，使地層更難發生破裂。此外，在x方向上，靠近SRV 左右邊緣較近(＜30 m)的區域內，則由于干擾作用導致地層破裂壓力下降，最多可下降1 MPa，說明此區域地層拉伸破壞所需流壓降低，更易發生地層破裂。

圖8 x-z垂面內不同層位間地層破裂壓力擾動情況Fig.8 Disturbance of formation fracturing pressure between different horizons in the x-z vertical plane

5 結論

(1)投產后，在縱向上會影響3 倍SRV 厚度內儲層的地應力大小和方向。該影響在SRV 中間層最大、在SRV 遠端最小。SRV 對立體開發造成的儲層地應力演化具有增強作用。

(2)總應力在演化過程中總體上遞減，但是在壓降誘發的有效應力增加大于壓降的區域，可能出現暫時的總應力增加。

(3)水平井生產造成的應力干擾可波及至深淺±20 m (總厚度40 m)區域內，對最小主應力和地層破裂壓力造成干擾，增加或降低地層破裂難度。