郭旭洋，金 衍，林伯韜，盧運(yùn)虎，訾敬玉

(中國(guó)石油大學(xué)(北京)石油工程學(xué)院，北京 102249)

天然氣水合物儲(chǔ)量巨大，并具有清潔高效的特點(diǎn)，在中國(guó)國(guó)南海已經(jīng)開展2次試采，其中第2次為水平井試采[1]。水平井降壓開采過程中的工程擾動(dòng)會(huì)誘發(fā)水合物分解，導(dǎo)致井周沉積物力學(xué)性質(zhì)劣化，井周應(yīng)力時(shí)空演化特征復(fù)雜，存在出砂和海底滑坡的可能性[2]。因此明確試采過程中的沉積物變形破壞機(jī)制和流-固-熱-化多場(chǎng)耦合機(jī)制對(duì)于南海水合物安全、高效鉆采具有意義。水合物沉積物膠結(jié)性較差，水合物的分解對(duì)其力學(xué)性質(zhì)影響顯著。研究人員采用人工合成樣品和現(xiàn)場(chǎng)保壓取芯等手段開展了大量三軸力學(xué)試驗(yàn)，確定多種水合物沉積物試樣的強(qiáng)度參數(shù)和變形特性，為現(xiàn)場(chǎng)尺度建模分析提供關(guān)鍵力學(xué)參數(shù)。水合物和沉積物基質(zhì)共同決定水合物沉積物的力學(xué)性質(zhì)，但是基質(zhì)占據(jù)主導(dǎo)作用。增加水合物飽和度或減小有效圍壓均能夠使沉積物應(yīng)力-應(yīng)變關(guān)系由應(yīng)變硬化轉(zhuǎn)向應(yīng)變軟化。沉積物內(nèi)聚力受水合物飽和度影響作用明顯，而內(nèi)摩擦角僅在水合物飽和度較高時(shí)受水合物飽和度影響較大[3-9]。針對(duì)中國(guó)南海北部深水淺層天然氣水合物沉積物的力學(xué)試驗(yàn)顯示，有效圍壓增加的情況下，沉積物性質(zhì)會(huì)呈塑性，有效圍壓和水合物飽和度綜合影響沉積物力學(xué)性質(zhì)，導(dǎo)致水合物沉積物鉆采過程中地層失穩(wěn)風(fēng)險(xiǎn)較大[10]。力學(xué)試驗(yàn)獲取的水合物沉積物力學(xué)性質(zhì)為更大尺度的水合物沉積物鉆采過程建模分析提供了基礎(chǔ)數(shù)據(jù)支持，結(jié)合水合物分解動(dòng)力學(xué)、質(zhì)量守恒和能量守恒等，能夠?qū)λ衔锍练e物鉆采過程中的多種物理場(chǎng)演化開展定量建模分析。通過有限元建模分析發(fā)現(xiàn)，較小的井眼半徑更易導(dǎo)致水合物水平井發(fā)生井壁失穩(wěn)，這種影響隨鉆井時(shí)間加劇[11]。水合物相變、傳熱傳質(zhì)和地層變形的耦合建模有助于確定水合物分解前緣和解釋出砂和穩(wěn)產(chǎn)時(shí)間短等現(xiàn)象，為商業(yè)化開發(fā)提供理論參考，而生產(chǎn)壓力、滲透率、初始飽和度和溫度等都與地層變形有一定聯(lián)系[12-13]。基于TOUGH等平臺(tái)，能夠定量表征水合物開采過程中的傳熱傳質(zhì)問題和水合物分解規(guī)律[14]。在此基礎(chǔ)上通過孔隙彈性力學(xué)理論進(jìn)行流固耦合研究，認(rèn)為固體力學(xué)場(chǎng)的耦合方法對(duì)模擬結(jié)果會(huì)產(chǎn)生顯著影響。全耦合對(duì)收斂性和準(zhǔn)確性的提升有較大幫助，而如果采用順序耦合則需要更嚴(yán)密的耦合條件以確保數(shù)值求解的穩(wěn)定性和準(zhǔn)確性[15-17]。筆者根據(jù)南海水合物沉積物力學(xué)特征，建立水平井降壓開采過程的流-固-熱-化多場(chǎng)耦合有限元模型，分析沉積物內(nèi)的變形、破壞和海底沉降規(guī)律。

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

以水合物沉積物為對(duì)象構(gòu)建多物理場(chǎng)耦合模型，表征沉積物內(nèi)水合物分解、熱傳遞、氣水滲流和變形破壞過程。根據(jù)質(zhì)量守恒，多孔介質(zhì)內(nèi)的流動(dòng)過程可以表示為

(1)

式中，φ為孔隙度；Si為流體飽和度；ρi為流體密度，kg/m3；vi為流體速度，m/s；si為匯/源質(zhì)量速度，kg/(m3·s)；下角i代表相態(tài)類型。

沉積物中流體滲流速度低，根據(jù)達(dá)西定律[18]，多孔介質(zhì)內(nèi)流體流速表示為

(2)

式中，k為絕對(duì)滲透率，m2；Kri為相對(duì)滲透率；μi為黏度，Pa·s；pi為孔隙流體壓力，Pa；g為重力加速度，m/s2。

水合物降壓開采過程中的傳熱問題需要考慮沉積物固相骨架內(nèi)的熱傳導(dǎo)、多孔介質(zhì)內(nèi)流體的熱對(duì)流和水合物分解產(chǎn)生的吸熱過程[19]。該過程中的傳熱問題可以表示為

(hC+hE)=Q.

(3)

式中，E為比內(nèi)能，J/kg；S為飽和度；h為對(duì)流傳熱項(xiàng)；Q為熱匯/源，J/(m3·s)；下角L、G、H分別代表液相、氣相和水合物相3種相態(tài)。

降壓開采過程中，水平井井底壓力降低，導(dǎo)致天然氣水合物穩(wěn)定區(qū)域平衡被破壞，水合物開始分解[20-21]。根據(jù)水合物分解動(dòng)力學(xué)模型，降壓過程中水合物分解速率表示為

RMH=-kdMMH(pe-pG)As.

(4)

式中，RMH為分解速率；kd為反應(yīng)速率，Pa·s；MMH為水合物摩爾質(zhì)量，kg/mol；pe為相平衡壓力，Pa；pG為氣相壓力，Pa；As為反應(yīng)比表面積，1/m2。

由于溫度和壓力的變化會(huì)導(dǎo)致多孔介質(zhì)及多孔介質(zhì)內(nèi)流體發(fā)生體積變化，具體可以表示為

(5)

(6)

式中，χ為壓縮系數(shù)，1/Pa；β為膨脹系數(shù)，1/K；T為溫度，K。

由于降壓，沉積物內(nèi)水合物分解，水合物飽和度下降，導(dǎo)致沉積物力學(xué)性質(zhì)出現(xiàn)劣化[22]。

水合物分解會(huì)改變內(nèi)聚力，c可以表示為水合物飽和度的函數(shù)

(7)

式中，c0為不考慮水合物部分的沉積物內(nèi)聚力，Pa；θ內(nèi)摩擦角，(°)；SMH為水合物飽和度；α和β為系數(shù)。

為表征降壓開采過程中的沉積物力學(xué)響應(yīng)，需要根據(jù)動(dòng)量平衡對(duì)應(yīng)力張量建立表達(dá)式

(8)

式中，σ為總應(yīng)力張量；ρb為體積密度，kg/m3。

根據(jù)無窮小變換，應(yīng)變張量ε和位移矢量u的關(guān)系可以表示為

(9)

在多孔介質(zhì)中，本構(gòu)關(guān)系表示為

δσ=C:δ(ε-εp)-αδpI.

(10)

式中，C為彈性張量；εp為塑性應(yīng)變；α為耦合系數(shù)；I為二階單位張量。

利用Mohr-Coulomb和Drucker-Prager模型表征彈塑性變形下的剪切破壞。屈服函數(shù)f和Drucker-Prager模型的塑性勢(shì)函數(shù)g可分別表示為

(11)

(12)

式中，I1為有效應(yīng)力第一不變量；J2為有效偏應(yīng)力第二不變量；βf、κf、βg和κg為破壞包絡(luò)線系數(shù)[23-25]。

Mohr-Coulomb模型可表示為

(13)

(14)

(15)

(16)

2 水平井降壓開采模型

根據(jù)上述控制方程，建立水合物沉積物降壓開采模型。由于水平井鉆采是南海試采工作的一項(xiàng)重要內(nèi)容，也是提高產(chǎn)氣效率的關(guān)鍵手段，因此建立單個(gè)水平井模型和雙水平井模型，其中雙水平井的水平段深度不同。圖1為單個(gè)水平井筒降壓開采模型和雙水平井筒降壓開采模型。由于水平井筒沿x方向貫穿模型，井筒長(zhǎng)度較長(zhǎng)，可以采用平面應(yīng)變假設(shè)簡(jiǎn)化數(shù)值運(yùn)算，僅分析y-z截面內(nèi)的二維流-固-熱-化耦合問題。對(duì)于單井模型，井筒位于y-z截面正中(25 m，25 m)；對(duì)于雙井模型，井筒分別位于y-z截面的(25 m，15 m)和(25 m，35 m)。對(duì)于二維y-z平面模型，y方向?yàn)镾Hmax，z方向?yàn)镾v，x方向施加Shmin。在厚度為50 m的沉積物層內(nèi)，假設(shè)地層壓力梯度和地應(yīng)力梯度的影響較小，地層壓力為14 MPa，垂向主應(yīng)力為16.6 MPa，水平最大主應(yīng)力為14.76 MPa，不隨深度而增大。沉積物密度為2 300 kg/m3，初始飽和度為40%，孔隙度為0.15，滲透率為1.8×10-2μm2，地層溫度為284.15 K，干濕導(dǎo)熱系數(shù)分別為1和3 W/(m·K-1)，降壓生產(chǎn)壓力為3 MPa。

圖1 單個(gè)水平井筒和雙水平井筒降壓開采模型Fig.1 Sketch of models for depressurization in a single horizontal well and two horizontal wells

3 結(jié)果分析

根據(jù)圖1所示模型，開展水平井筒降壓開采過程中的數(shù)值模擬，分別針對(duì)單個(gè)水平井筒和雙水平井筒開展分析，模擬降壓開采時(shí)長(zhǎng)為15 d。降壓開采模型能夠表征該過程中井周沉積物內(nèi)的水合物分解范圍、分解吸熱引發(fā)溫度變化、孔隙壓力降低和彈塑性演化，并能夠預(yù)測(cè)水平井產(chǎn)能。在二維y-z平面內(nèi)，由于數(shù)值解關(guān)于z軸對(duì)稱，僅展示y軸正半軸數(shù)值結(jié)果。

圖2為采用單個(gè)水平井和2口水平井降壓開采5和15 d時(shí)的井周孔隙流體壓力分布。近井區(qū)域壓力降最明顯，這是由井筒降壓開采使用的較低生產(chǎn)壓力造成的。隨著距離向遠(yuǎn)場(chǎng)移動(dòng)，壓力越接近初始?jí)毫?。隨著降壓開采時(shí)間的增加，壓降前緣的波及范圍更廣。對(duì)于單井模型，降壓開采5 d時(shí)遠(yuǎn)場(chǎng)壓力仍等于初始?jí)毫?4 MPa，但開采15 d后遠(yuǎn)場(chǎng)壓力也開始下降。相對(duì)應(yīng)地，采用雙井模型能夠增加壓降前緣的波及范圍。對(duì)比開采5 d時(shí)的單井模型和雙井模型的孔隙流體壓力結(jié)果，發(fā)現(xiàn)雙井模型井周壓降范圍更大，壓降更加劇烈。開采15 d后雙井模型的遠(yuǎn)場(chǎng)區(qū)域壓力普遍降至低于11 MPa，單井模型原廠區(qū)域壓力仍可達(dá)13 MPa，說明采用兩口水平井降壓開采能夠更有效地建立地層內(nèi)的壓降，為甲烷流動(dòng)提供促進(jìn)作用。

圖2 單個(gè)水平井筒和雙水平井筒降壓開采5和15 d孔隙流體壓力分布Fig.2 Pore pressure distribution after 5 and 15 d of depressurization for single- and two-well scenarios

圖3為單井模型和雙井模型內(nèi)降壓開采5和15 d時(shí)的水合物飽和度分布特征，水合物飽和度為0的區(qū)域表示水合物完全分解。結(jié)果顯示水合物分解區(qū)域近似呈圓形，這是因?yàn)閥-z平面內(nèi)的井周壓降前緣也近似呈圓形，壓降導(dǎo)致水合物穩(wěn)定區(qū)的相平衡遭到破壞而出現(xiàn)水合物分解，分解區(qū)域受到壓降前緣影響較為直接。根據(jù)圖3數(shù)據(jù)，5 d時(shí)單井水合物分解前緣距水平井筒擴(kuò)展1.01 m，雙井為0.96 m；15 d時(shí)水合物分解前緣擴(kuò)展的距離在單井模型中為1.74 m，在雙井模型中則為1.58 m。說明增加水平井?dāng)?shù)量和增加降壓開采時(shí)間均會(huì)對(duì)水合物分解范圍產(chǎn)生影響。隨著降壓開采時(shí)間的增加，水合物分解前緣向遠(yuǎn)場(chǎng)擴(kuò)展，15 d內(nèi)能夠擴(kuò)展1.58～1.74 m。

圖3 單個(gè)水平井筒和雙水平井筒降壓開采5 d和15 d水合物飽和度分布Fig.3 Hydrate saturation distribution after 5 and 15 d of depressurization for single- and two-well scenarios

圖4為單井模型和雙井模型降壓開采5和15 d后的溫度分布。由于水合物分解是吸熱過程，導(dǎo)致分解區(qū)域及鄰近沉積物內(nèi)溫度降低。降溫前緣也近似呈圓形，與壓降前緣和水合物分解前緣形狀類似。造成這一現(xiàn)象的原因除了壓降與水合物分解直接相關(guān)外，也與式(1)、(3)有關(guān)，兩者均為拋物型微分方程，呈現(xiàn)的傳質(zhì)傳熱特征具有類似性。圖4表明，單井模型中的降溫前緣在5 d時(shí)為1.25 m，在15 d時(shí)為2.71 m；雙井模型中的降溫前緣在5 d時(shí)為1.21 m，在15 d時(shí)為2.67 m。這說明開采時(shí)間仍是影響降溫范圍的主要因素，水平井?dāng)?shù)量對(duì)降溫范圍的影響相對(duì)較小。

圖4 單個(gè)水平井筒和雙水平井筒降壓開采5和15 d溫度分布Fig.4 Temperature distribution after 5 and 15 d of depressurization for single- and two-well scenarios

圖5 水平井累積產(chǎn)氣對(duì)比Fig.5 Comparison of cumulative gas production between studied scenarios

根據(jù)降壓開采多物理場(chǎng)耦合模擬可以進(jìn)行水平井產(chǎn)氣預(yù)測(cè)(圖5)。對(duì)比了單井模型水平井產(chǎn)量、雙井模型井1(z=35 m)產(chǎn)量、雙井模型井2(z=15 m)產(chǎn)量和雙井模型兩口水平井平均產(chǎn)量?？傮w累產(chǎn)曲線顯示單井累產(chǎn)氣差異不明顯，這是因?yàn)閮煽谒骄噍^遠(yuǎn)，盡管井周壓降區(qū)域互相波及，但水合物分解區(qū)域和溫度降低區(qū)域仍相隔一段距離，而水平井產(chǎn)氣主要是由水合物分解區(qū)釋放的甲烷提供，因此在井周水合物分解范圍相似的情況下，獲得的水平井累產(chǎn)氣結(jié)果也較為近似。比較單井模型和雙井模型內(nèi)兩口井的累積產(chǎn)氣曲線，發(fā)現(xiàn)15 d內(nèi)累產(chǎn)氣由高到低排序?yàn)殡p井模型井1>單井模型>雙井模型井2。分析發(fā)現(xiàn)這一關(guān)系能夠和三口水平井井位深度建立起聯(lián)系，累產(chǎn)氣最高得雙井模型井1深度最小，單井模型其次，雙井模型井2深度最大。由于在y-z平面內(nèi)多場(chǎng)耦合模型考慮了重力的作用，導(dǎo)致降壓開采后深度更大的區(qū)域壓力相對(duì)更高、壓降相對(duì)更弱。據(jù)此，在深度更大的區(qū)域進(jìn)行水平井降壓開采形成的壓降差異會(huì)略小于在深度較小的區(qū)域進(jìn)行降壓開采，一定程度上減緩水合物分解和甲烷釋放，減小了對(duì)應(yīng)水平井累產(chǎn)氣。此外，還計(jì)算了雙井模型兩口井累產(chǎn)的平均值，發(fā)現(xiàn)其數(shù)值略高于單井模型。這說明在本研究設(shè)定的條件下，增加水平井?dāng)?shù)量能夠適當(dāng)提升單井累產(chǎn)效果。4條累產(chǎn)曲線均呈現(xiàn)出開采初期產(chǎn)氣效率很高、開采一段時(shí)間后產(chǎn)量急劇下降的特點(diǎn)，與前期國(guó)內(nèi)外部分試采工作結(jié)果具有一定的一致性。

降壓開采誘發(fā)的沉積物變形與破壞也是需要研究的重點(diǎn)之一。圖6為采用單井模型和雙井模型降壓開采5和15 d后的沿垂直方向的總應(yīng)力分布情況。由于近井區(qū)域受到降壓開采影響，孔隙流體壓力下降明顯，導(dǎo)致垂直方向總應(yīng)力在井筒和近井區(qū)域處于較低值。在同一時(shí)間，增加水平井?dāng)?shù)量進(jìn)一步加劇了垂直方向總應(yīng)力的降低。隨著開采時(shí)間增加，垂直方向總應(yīng)力降低的區(qū)域逐漸增大，而且遠(yuǎn)場(chǎng)區(qū)域也出現(xiàn)一定的垂向總應(yīng)力降低。通過與圖2對(duì)比，發(fā)現(xiàn)垂直方向總應(yīng)力的演化規(guī)律與孔隙流體壓力額演化規(guī)律具有一定關(guān)聯(lián)性，只是因?yàn)樵诳紫稄椥粤W(xué)理論中，孔隙流體壓力是總應(yīng)力中各正應(yīng)力分量的重要組成部分。在15 d時(shí)，雙井模型內(nèi)遠(yuǎn)場(chǎng)垂向總應(yīng)力比單井模型更低，也是由雙井模型遠(yuǎn)場(chǎng)孔隙流體壓力更低造成的。

圖6 單個(gè)水平井筒和雙水平井筒降壓開采5 d和15 d垂向應(yīng)力分布Fig.6 Vertical stress distribution after 5 and 15 d of depressurization for single- and two-well scenarios

表征降壓開采誘發(fā)的彈塑性變形能夠幫助預(yù)測(cè)氣井出砂和海底沉降。圖7為單井模型水平井、雙井模型中井1和井2同深度的塑性變形時(shí)間演化規(guī)律。由圖7(a)看出，在降壓開采起始時(shí)，僅在井壁處存在塑性應(yīng)變，儲(chǔ)層沉積物均無塑性應(yīng)變發(fā)生。降壓開采0.5～1 d后，由井壁至沉積物層內(nèi)出現(xiàn)隨距離逐漸減小的塑性應(yīng)變，這是由于近井區(qū)域壓降劇烈，沉積物固體骨架變形程度更大。0.5～1 d時(shí)，由于水合物分解尚不明顯，水合物飽和度對(duì)塑性變形的影響不顯著。降壓開采5 d后，由于水合物分解前緣的擴(kuò)展，分解區(qū)向未分解區(qū)過渡的塑性變形呈現(xiàn)出突變特征，分解區(qū)展現(xiàn)的塑性應(yīng)變高于未分解區(qū)。從0.5～15 d的結(jié)果顯示，水解物分解前緣和塑性應(yīng)變突變點(diǎn)具有較強(qiáng)相關(guān)性，分解前緣外的未分解儲(chǔ)層內(nèi)的塑性應(yīng)變則由近及遠(yuǎn)平滑降低。該結(jié)果說明水合物分解區(qū)塑性變形特點(diǎn)明顯，更易出砂和誘發(fā)地層失穩(wěn)。

圖7 不同模型井塑性變形時(shí)間演化Fig.7 Temporal evolution of plastic deformation in different model wells

由圖7(b)、(c)看出雙井模型中井1和井2各時(shí)間點(diǎn)的塑性變形分布特征與單井模型結(jié)果比較接近，均呈現(xiàn)出水合物分解前緣與塑性應(yīng)變突變點(diǎn)吻合的情況。5、10和15 d的塑性應(yīng)變?cè)?0 m以外出現(xiàn)小幅度升高的現(xiàn)象，這是受到了模型的邊界條件的影響。在y-z平面的右側(cè)邊界設(shè)置為固定位移邊界，雙井模型的壓降和變形程度和范圍更大，導(dǎo)致近井區(qū)域壓縮變形程度較高，由于右側(cè)邊界為固定邊界，在靠近右側(cè)邊界的部分區(qū)域出現(xiàn)拉伸變形。說明在本研究中模型尺寸和邊界條件對(duì)塑性變形求解過程有影響。

4 結(jié) 論

(1)由于降壓開采中的傳熱和傳質(zhì)過程控制方程均具有拋物型微分方程的特征，沉積物層內(nèi)的孔隙壓力、溫度和水合物分解前緣的時(shí)空演化規(guī)律具有較強(qiáng)的相關(guān)性。

(2)沉積物層內(nèi)的總應(yīng)力時(shí)空演化規(guī)律受孔隙壓力演化的影響明顯，近井區(qū)域的總應(yīng)力數(shù)值小于遠(yuǎn)場(chǎng)受壓降影響較小區(qū)域的總應(yīng)力數(shù)值。這一現(xiàn)象能夠通過孔隙彈性力學(xué)理論進(jìn)行解釋。

(3)在水平井降壓開采的多場(chǎng)耦合研究時(shí)，重力作用對(duì)水平井產(chǎn)能、力學(xué)響應(yīng)和孔隙壓力分布的影響較為明顯。降壓開采后，重力作用導(dǎo)致深層水平井壓降程度低于淺層水平井，一定程度上抑制水合物分解前緣的擴(kuò)展。本研究數(shù)模顯示，淺層水平井累積產(chǎn)氣量更高。

(4)水合物分解會(huì)明顯地劣化沉積物力學(xué)性質(zhì)，誘發(fā)較強(qiáng)的塑性應(yīng)變。在降壓開采初期水合物分解前緣未形成時(shí)，沉積物層塑性應(yīng)變不明顯；隨著降壓開采時(shí)間的增加，水合物分解區(qū)力學(xué)性質(zhì)劣化顯著，塑性應(yīng)變明顯高于未分解區(qū),內(nèi)聚力等強(qiáng)度參數(shù)變低,在水合物分解前緣附近的塑性應(yīng)變會(huì)出現(xiàn)明顯的階梯式下降特征。此外，建模分析時(shí)的模型尺寸和邊界條件設(shè)置均在一定程度上可能影響塑性變形求解過程。