壓縮加載條件下含水合物沉積物蠕變特性分析

李輝，張旭輝，陸程，謝鵬飛，魯曉兵

1.中國科學院力學研究所，北京 100190

2.中國科學院大學，北京 100049

3.中國地質調查局油氣資源調查中心，北京 100083

4.中國地質大學（北京）工程技術學院，北京 100083

天然氣水合物（以下簡稱水合物）儲量豐富，是最有可能代替化石燃料的清潔能源之一[1-3]。水合物沉積物層存在于低溫、高壓的自然環境中，陸地與海洋均有分布。中國海域水合物主要分布于南海地區[4]，據估算其儲量約為6.4×1013m3，與中國常規天然氣資源總量基本相當[5]，具有巨大的開發潛力。隨著中國對于海域水合物研究的不斷深入，中國地質調查局分別在2017 年和2020 年進行了海域水合物的試采[6-7]，實現短期穩定產氣。

國內外學者針對含水合物沉積物的力學特性開展了豐富的力學試驗研究[8-11]。Winters 等[12-13]通過核磁共振和三軸壓縮試驗分別對原狀樣品和人工合成樣品進行測量，得到了水合物沉積物的強度和壓縮波速。Masui 等[14]通過實驗證實了原狀樣品和人工合成樣品具有相似的力學強度，但是應力-應變關系隨初始孔隙率和顆粒級配的變化而變化。Yun 等[15]通過室內試驗證明了四氫呋喃水合物（THF）和甲烷水合物具有相似的力學強度，建議使用THF 進行力學性質研究。對于孔隙型水合物的研究較為廣泛，Ren[16]等通過核磁共振發現水合物首先生長于小孔隙中，黏土會延緩整個生長過程。相較于均勻分布的含水合物沉積物樣品，非均勻分布會對沉積物的物理和力學性質產生影響。Xie等[17-18]考慮了水合物沉積物的非均勻性，在人工樣品中加入不同形狀的水合物塊體模擬不同的賦存形式，得到了不同形式下的水合物沉積物的強度和描述模型。含水合物沉積物的滲透率和破裂壓力等均與水合物飽和度及分布形式有關[19-20]。然而大規模商業化開采水合物需要中長期穩定產氣，因此還需要考慮儲層的蠕變特性[21]。

目前，針對粉黏土蠕變特性的研究已經較為成熟[22-24]。余云燕等[25]針對非飽和鹽漬土開展了一系列三軸固結不排水蠕變試驗，在Singh-Mitchell 模型的基礎上進行修正，提出了適合非飽和鹽漬土的經驗模型。羅慶姿等[26]通過三軸不排水壓縮試驗對汕頭的軟土進行了分析。Lai 等[27]針對三峽滑土帶中的非飽和土進行三軸蠕變試驗，引入了考慮基質吸力的相關參數，建立了適用于非飽和土的Singh-Mitchell 蠕變模型。劉業科等[28]通過梯級載荷加載方法進行室內三軸固結不排水的蠕變試驗，在此基礎上對Singh-Mitchell 模型進行修正使其能描述零應力零應變水平。關于常規粉黏土的蠕變特性試驗已經開展了較多研究，積累了較多經驗[29]。區別于常規軟黏土的蠕變特性，含水合物粉質黏土在開采過程中存在水合物相變、滲流等過程，水合物分解前后沉積物的孔隙結構也相應發生改變，目前針對含水合物粉黏土的蠕變特性研究和相關模型的研究較少[30]。

本工作采用固結儀進行蠕變試驗，研究了含水合物粉黏土的蠕變特性，并提出了描述模型。首先，通過試驗獲得了含水合物粉黏土的蠕變特性曲線，分析了含水合物粉黏土分解后的蠕變特性；在Singh-Mitchell 模型中引入水合物飽和度參數，構建并驗證修正后的Singh-Mitchell 模型，得到了可以描述含水合物粉黏土的蠕變特性模型。

1 區域概況

研究區域位于中國南海北部的西北大陸邊緣，呈NE 向延展，面積達4.5×104km2，具備穩定的高壓低溫環境，可以形成天然氣水合物，也是我國重要的常規石油富集區域[31-32](圖1)。試驗采用研究區海域水合物儲層的粉黏土沉積物，水合物賦存于海底7～158 m 范圍內的第四系沉積物中，該區域平均水深約為1 000 m。

圖1 研究區域位置[31]Fig.1 Location of the research area[31]

2 試驗樣品來源及方法

2.1 試驗材料

土體的基本物理和力學性質如下：土體顆粒比重Gs=2.65，孔隙度n=0.509，平均粒徑d50=0.0067，其中黏質組分（d≤4 μm）約占50%，砂質組分（4 μm＜d≤63 μm）約占45%，顆粒級配曲線如圖2 所示。根據《土的工程分類標準》，試驗用土粒徑小于0.075 mm的顆粒質量分數大于50%，屬于細粒土；不均勻系數Cu= 13.56，屬于不均勻土，曲率系數Cc= 0.57，級配不連續。綜上，試驗土為級配不良的細粒土。

圖2 沉積物顆粒級配曲線圖Fig.2 Curves of sediment particle size distribution

試驗所用的四氫呋喃溶液采購于中國國藥試劑網，濃度為99.9%，通過加入蒸餾水配置成不同濃度的四氫呋喃溶液。質量分數為19%的四氫呋喃溶液在低溫常壓條件下可以完全合成四氫呋喃水合物，考慮到四氫呋喃的揮發性，試驗中使用四氫呋喃溶液的質量分數為21%。根據試驗樣品所需的水合物飽和度，計算并配置一定質量分數的四氫呋喃溶液。

2.2 試驗樣品制備與試驗設置

試驗介質為粉黏土，基于試驗區域水合物儲層實際工況和土工試驗結果，粉黏土比重Gs= 2.65，制備水合物沉積物樣品時選取干密度 ρd=1.3 g/cm3，研究在恒定的上覆層荷載下水合物儲層分解后對水合物儲層的蠕變力學特性的影響。具體制備過程如下：

土骨架制備：將洗鹽后的粉黏土沉積物置于烘箱中在105℃條件下烘干12 h 后取出冷卻，研磨土體至顆粒完全粉碎（土體顆粒粒徑小于0.2 mm）。根據干密度稱取一定質量的沉積物平均二等分，將每份土體加入適量水后填入環刀內砸實至所要求的高度，在每層樣品制備完成后將試樣拉毛，以降低其分層的影響，制備成面積為30 cm2、高2 cm樣品，置于烘箱中在105℃條件下烘干12 h 后取出冷卻。

水合物合成：根據試驗中水合物飽和度的設置配置不同濃度的四氫呋喃溶液，將制備好的土骨架置于密封罐內，密封罐與真空裝置和儲液裝置相連，開啟抽真空裝置對樣品抽取真空兩小時，完成后通過儲液裝置自下而上澆注配置好的四氫呋喃溶液，保持密封狀態一天至土樣含液量穩定。最后取出浸泡完成的沉積物樣品置于0～2 ℃的冷柜中1～2 d 合成水合物，如圖3 所示。

圖3 試驗前后的試樣Fig.3 Sample before and after test

在試驗前取出含水合物沉積物樣品置于密封罐內分解，從而得到不同水合物飽和度分解后的沉積物樣品。采用分別加載的方式進行，在荷載的作用下，連續觀測。加載應力分別為：①號0.5 MPa、②號2 MPa、③號3 MPa 條件下進行試驗；分解前水合物飽和度分別為 0%、20%、60%、80%。試驗儀器使用固結儀如圖4，通過采用單杠桿固結儀試驗獲得水合物沉積物的蠕變力學響應參數。

圖4 高壓固結儀實物圖Fig.4 High-pressure consolidation apparatus

高壓固結儀放置于20℃的恒溫室中，在試驗過程中保持樣品處于飽和狀態，當沉降值穩定后樣品停止變形，即認為試樣的變形不超過0.001 mm/d 后達到穩定狀態。

3 結果和分析

3.1 試驗結果

通過對原始數據的整理，得到了分別加載方式下，不同水合物飽和度（SH）粉質黏土樣品分解后在軸向壓力σ下試驗時間t和應變ε的蠕變試驗曲線（圖5）。

圖5 不同飽和度下蠕變全過程曲線圖Fig.5 The whole-process curve of creeping under different saturation degrees

由圖5 中可以看出，在分別加載的條件下，含水合物沉積物分解后樣品在各級荷載下的應變-時間關系曲線的變化規律基本相同，只有初始蠕變階段和穩定蠕變階段，并未出現加速蠕變階段。施加荷載的初期，均產生了較大的瞬間變形，曲線呈飛躍式增長。隨著時間的推移，曲線的斜率逐漸減小，表示試樣的變形速率逐漸減小。最終，試樣的變形速率趨于穩定，曲線呈斜率較小的直線，試樣的變形隨時間不斷增大，但該階段變形值較小。

含水合物粉質黏土樣品的蠕變特性曲線可以分為瞬時變形、固結變形和蠕變變形3 個階段[33]。水合物的賦存改變了粉黏土的初始孔隙結構，在施加荷載的瞬間會產生較大的變形進入瞬時變形階段；當樣品孔隙被壓縮，孔隙水被排出，而由于黏土沉積物骨架的滲透率較低因此會產生超靜孔隙水壓力并需要較長時間消散，處于固結變形階段；最后固結完成后，在恒定荷載的作用下，由于恒定軸向壓力的作用，土體骨架持續變形，應變不斷增加但總體應變率較小，最后趨于穩定，處于蠕變變形階段（圖6）。無水合物沉積物應變率和時間的曲線見圖7。瞬時變形認為當軸向荷載加載完成后完成，取1 min 對應的應變作為瞬時變形。將圖7中的轉折點作為固結變形和蠕變變形分界點，即將應變率1×10-5/min 作為固結變形和蠕變變形分界點。

圖6 變形過程三階段示意圖Fig.6 Diagram of three-stage deformation process

圖7 無水合物沉積物應變率和時間的曲線Fig.7 Curves of strain rate and time for hydrate-free sediments

3.2 Singh-Mitchell 蠕變模型

Singh 和Mitchell 等[34]在總結大量土體蠕變的基礎上，使用指數函數及冪函數的形式描述黏土的應力-應變-時間關系，并在此基礎上提出了Singh-Mitchell 蠕變模型。Singh-Mitchell 蠕變模型具有參數少、適用范圍廣等特點，因此本文先擇此模型來描述粉黏土的蠕變特性。Singh-Mitchell 蠕變模型可以表示為：

式(2)中 ε為軸向應變， ε0為初始軸向應變；在不 考 慮 初 始 的 軸 向 應 變 時， ε0=0 ，令 λ=1-m，，得到式(3)：

式(3)即蠕變方程，其中需要確定的參數有3 個分別為B、β、λ。

當t=t1時，有 ε=ε1，可以得到方程(4)，兩邊取對數得到方程(5)：

β和B的值可以直接通過單位參考時間t1的lnε-Dr直線關系圖得到。土體的非線性蠕變是一簇曲線，從上式中可以看出Singh-Mitchell 蠕變模型也是一簇指數形式的曲線，形狀與現有的粉黏土蠕變試驗結果吻合，因此選擇Singh-Mitchell 蠕變模型來描述粉黏土的蠕變特性是合適的。

3.3 試驗用粉黏土蠕變特性分析

將蠕變特性曲線繪制在雙對數坐標軸中，可以得到 lnε-lnt關系圖（圖8），圖中可以看出，在不同的剪切應力水平下樣品的蠕變試驗曲線為一組平行近似的直線，直線的斜率為模型中 λ 值。將t1=1min時的應變帶入式(4)和式(5)中得到 lnεr-Dr關系曲線，擬合后具有良好的線性關系，通過式(5)可以計算出不含水合物沉積物的β和B值（圖9）。

圖8 無水合物試樣 lgε vs lgt蠕變特性曲線Fig.8 The lgε vs lgt creep characteristic curve of hydrate-free samples

圖9 lnεr vs Dr關系擬合曲線Fig.9 Fitting curve lnεr vs Dr relationship

考慮含水合物沉積物的剪切應力，偏應力的計算公式為q=(1-K0)*σ ， 其中K0為靜止土壓力系數，σ為軸向應力。根據Simpson 等[35]正常固結黏土K0經驗公式：

將三軸剪切試驗內摩擦角 φ′=5°帶入式(6)后得到K0=0.84，即q=0.16σ，同時注意到不同圍壓下qf=0.15～0.2 MPa[36]，，因此后續取值Dr=σ。

帶入上述參數計算Singh-Mitchell 蠕變模型中的相關參數值，結果如表1 所示。

表1 Singh-Mitchell 蠕變模型中參數計算值Table 1 Parameter values of the Singh-Mitchell creeping model

將表1 中計算的模型參數帶入式(4)中，以軸向壓力σ為0.5 MPa 為例得到Singh-Mitchell 蠕變模型的表達式：

將模型計算曲線和試驗值進行對比（圖10）?？梢钥闯?，Singh-Mitchell 蠕變模型計算曲線和蠕變試驗曲線的發展趨勢具有較好的一致性，在加載瞬間應變快速增長，隨著時間的增加軸向應變迅速衰減。在較高的剪切應力水平下Singh-Mitchell 蠕變模型和試驗曲線具有較好的一致性，在較低的剪切應力水平下模型計算值略大于試驗值，有一定的誤差，總體擬合關系較好。

圖10 無水合物模型計算曲線和試驗值對比Fig.10 Comparison of model calculation curve to the test result of hydrate-free sediments

3.4 含水合物試樣蠕變特性分析

上一節中討論了不含水合物沉積物樣品的蠕變特性曲線，通過結果可以看出Singh-Mitchell 蠕變模型可以描述其蠕變特性，為了使Singh-Mitchell 蠕變模型可以描述含水合物沉積物樣品，需要對模型進行修正。引入水合物飽和度SH來描述水合物含量，SH是指水合物與孔隙體積之比[37]：

水合物在合成時會膨脹從而引起土體孔隙的變化，在經典模型中，模型參數A反映出土的組成、結構和應力歷史等特性，因此考慮通過對參數A進行修正來表征水合物飽和度的影響，在Singh-Mitchell 蠕變模型的基礎上引入模型如式(9)：

式中，α為模型參數，根據試驗數據，本文中取0.37。

將修正后的參數A和不同水合物飽和度經典模型中的參數A進行對比（圖11），可以看出其具有良好的對應關系。

圖11 經典模型和修正模型參數 A對比圖Fig.11 Comparison of classical model and modified model parameter A

由此建立了一種修正的Singh-Mitchell 蠕變模型：

當SH=0 時，模型又退化成經典的Singh-Mitchell蠕變模型。

通過計算，得到不同工況下修正的Singh-Mitchell蠕變模型參數，得到參數（表2）。

表2 修正的 Singh-Mitchell 蠕變模型參數值Table 2 Parameter values of modified Singh-Mitchell creeping model

將表2 中的參數帶入公式(10)中，得到模型計算值，同時和試驗數據對比如圖12 所示。

圖12 不同水合物飽和度模型計算曲線和試驗值對比Fig.12 Comparison of model calculated curves and the test values at different hydrate saturation degrees

從圖12 可以看出，修正后的Singh-Mitchell 蠕變模型在不同水合物飽和度和剪切應力下都能夠較好地預測均勻分布的含水合物粉質黏土樣品的蠕變特性，因此可以使用此修正的Singh-Mitchell 蠕變模型來預測含水合物沉積物樣品的蠕變特性曲線。修正的Singh-Mitchell 蠕變模型參數計算簡便、參數較少并且容易獲取，使用較為方便，因此本文建立的修正Singh-Mitchell 蠕變模型對描述含水合物沉積物樣品的蠕變特性曲線更為合適，對我國南海水合物儲層長期開采后的沉降預測具有一定的指導意義。

4 討論

4.1 蠕變試驗結果分析

南海海域水合物主要賦存于海底以下到200 m范圍內的粉黏土等未固結的沉積物內，具有含水高、孔隙度高的特點，水合物飽和度最高可達90%以上[38]。因此本文根據儲層的實際賦存條件，以水合物儲層的粉黏土為介質，在0.5、2 和3 MPa 的壓縮加載條件下進行了0、20%、60%、80%水合物飽和度分解后的蠕變特性試驗。試驗表明：加載過程中，含水合物沉積物經歷瞬時變形、固結變形和蠕變變形3 個階段；在加載的短時間內，瞬時變形階段完成，引起的軸向變形約占總變形量的70%～80%，最大值為85.8%，在固結變形階段，孔隙水在穩定的應力條件下從孔隙中排出，一般數小時內完成，引起的軸向應變占應變量的15%～25%，最后的蠕變變形階段一般數天內完成，引起的軸向應變小于總應變量的5%。軸向壓力和水合物飽和度均對軸向變形有影響。軸向壓力越大，蠕變變形越大，軸向壓力在較低水平時變化對軸向變形的影響大，隨著軸向壓力的增加，這種影響程度減小。水合物飽和度主要通過改變樣品的微觀孔隙結構進而影響蠕變特性，分解前的水合物飽和度越高，蠕變變形越大。水合物飽和度在低軸壓條件下對變形影響小，在軸壓為0.5 MPa 時，從無水合物到水合物飽和度增加到80%時，軸向變形從17.25%增加到18.8%，而在軸壓為3 MPa 時，軸向變形從26.83%增加到了34.73%。

4.2 基于Singh-Mitchell 的修正模型

數學模型是定量描述各個試驗變量之間相互關系的常用方法。之前的海洋地質調查海域水合物儲層粉質黏土進行蠕變特性的研究不足，本文利用蠕變試驗結果，通過建立修正的Singh-Mitchell 蠕變模型，對研究區域內水合物儲層土的蠕變特性進行初步分析。從蠕變模型中可以看出，通過軸向壓力和水合物飽和度兩個物理指標能夠建立出較好的蠕變模型，反映出這兩個物理指標是蠕變特性的主要影響因素。模型中參數α和土體初始條件有關，主要根據試驗獲取，反映出土體的蠕變特性與土體孔隙微觀結構、顆粒級配、應力歷史等多種因素相關。

5 結論

（1）含水合物沉積物樣品蠕變特性曲線可以分為瞬時變形、固結變形和蠕變變形3 個階段。瞬時變形數分鐘內完成，引起較大的變形，固結變形數小時內完成，引起較小的變形，蠕變變形數天內完成，引起極小的變形。

（2）蠕變變形隨著加載應力水平、水合物飽和度的增大而增加，水合物飽和度主要通過改變沉積物樣品的孔隙結構從而影響蠕變變形，水合物飽和度越高，孔隙結構受影響程度也越大，試樣的物理和力學性質越差。在低應力水平時，水合物飽和度對蠕變變形的影響不明顯。

（3）基于應力水平和水合物飽和度建立的修正Singh-Mitchell 蠕變模型，可以有效地預測含水合物沉積物的蠕變特性，對水合物儲層長期開采后的沉降變形預測具有指導意義，便于實際工程中推廣。但由于試驗組數少等原因本文未對水合物沉積物蠕變變形機理進行深入研究，針對非均勻分布的含水合物沉積物蠕變特性也需要進一步研究。