水合物開發(fā)儲層溫度場擾動研究

董京楠，班凡生，車陽，曲帥，王子健，劉智強

1 中國石油集團工程技術研究院有限公司，北京 100081

2 中國石油大學(北京)石油工程學院，北京 102249

0 引言

天然氣水合物是以甲烷為主要成分，在高壓力、低溫度的條件下形成的一種籠型晶體，在全球范圍內(nèi)有著巨大的儲量，也是公認的21世紀繼天然氣和石油能源之后最具開發(fā)潛力的清潔能源。在海底深水區(qū)和永久凍土深處廣泛存在，有著豐富的儲量[1-4]。在我國，除了東海、南海等海域均有大量的天然氣水合物礦藏，青藏髙原等凍土區(qū)域也已證實了天然氣水合物的存在，并成功獲取了水合物天然巖心樣品[5]。同時，據(jù)估計青藏高原凍土區(qū)天然氣水合物資源具有1.2×1012～2.4×1014m3的巨大的儲量[6]。

海洋天然氣水合物的開采不僅受溫度和壓力的影響，同時開采過程的地質(zhì)構(gòu)造作用、海平面的升降以及海底工程都會使得開發(fā)十分復雜?；谔烊粴馑衔镌诔爻合聵O其不穩(wěn)定性的特性，目前所提出的天然氣水合物開釆方法基本都是利用降壓、升溫等手段打破沉積層水合物的相平衡條件，促使其分解，然后再收集溢出的天然氣。具體的開發(fā)手段主要包括降壓法、熱激發(fā)法、化學試劑法和CO2置換開采法等[7-13]。水合物開采過程，溫度壓力的控制決定著開采的效率，而開采對儲層的影響也不可忽視。由于海底表層沉積物大多屬于黏土礦物，其強度較低，儲層內(nèi)的水合物晶體維持著基本骨架的穩(wěn)定性。開采時溫度和壓力的變化會導致水合物分解成氣體和水，水合物對骨架的膠結(jié)作用降低，從而降低凍土骨架的強度，進而導致井壁失穩(wěn)。而無論是海上還是陸地，水合物的大規(guī)模開釆必然會導致水合物沉積層的強度降低，繼而發(fā)生沉降，進而可能造成蓋層的失穩(wěn)，甚至會引發(fā)地震等危害[14-17]。

對于水合物分解引起的一系列問題，其中包含相變的傳熱問題是其中首先要面對的。在經(jīng)典的傳熱問題中，主要是考慮材料的本身熱導率、熱容對傳熱過程的影響。同時，大多數(shù)研究都是在常溫和高溫環(huán)境下的傳熱過程。而對于水合物的開發(fā)，第一其環(huán)境溫度較低，第二在開發(fā)過程中水合物會發(fā)生相變，同時吸收大量的熱，會嚴重影響儲層的溫度場。這些問題相對于傳統(tǒng)傳熱問題更加復雜，研究也相對較少，甚至是阻礙了水合物儲層開發(fā)的關鍵理論基礎之一。

本文主要研究在海洋區(qū)域天然氣水合物開采時，井筒周圍的地層溫度場變化。通過理論分析建立了低溫條件下含相變的巖石傳熱控制方程。在控制方程中，將比熱和潛熱合并成一個等效熱容參數(shù)?；贕ibbs-Thomson方程[18-19]，提出了一種利用壓汞數(shù)據(jù)得到孔隙內(nèi)含冰量隨溫度變化的方法。采用所提出的模型，通過類比水合物開發(fā)過程水合物相變與冰水相變的本質(zhì)，以水為介質(zhì)來模擬水合物開采過程中井筒溫度升高對周圍地層的溫度擾動，為進一步研究開采過程中儲層的穩(wěn)定性提供儲層內(nèi)相變傳熱的理論基礎。

1 理論模型的建立

通常情況下巖石內(nèi)的熱量的傳遞由熱導率決定，可以寫出一般形式的熱傳導控制方程(即Poisson方程[20])：

其 中T為 溫 度(℃)，t為 時 間(s)，ρ為 巖 石 密 度(kg·m-3)，c為 巖 石 比 熱(J·kg-1·K-1)，k是 巖 石 熱 導 率(W·m-1·K-1)，?是梯度算子。而當孔隙中發(fā)生相變時，熱量不僅會以比熱的形式傳遞(體系溫度發(fā)生變化，但組分相態(tài)不變)，還會以潛熱形式發(fā)生轉(zhuǎn)移(體系組分發(fā)生相變，但溫度不變)。因此，含相變的控制方程可以寫為：

其中L是水的潛熱(3.33×105J·kg-1)，ρi是冰的密度(kg·m-3)，θi孔隙中冰的體積分數(shù)。當假設凍融過程都是瞬時的，將dθi/dt= (dθi/dT)(dT/dt)帶入上式，則有

其中ρc-Lρidθi/dT這項定義為等效熱容，從等式可以看出，孔隙含冰量隨溫度的變化使得傳熱過程變得復雜。

為了求解這個非線性方程，首先需要確定含冰量與溫度間的關系(即dθi/dT)。冰在孔隙中發(fā)生相變，其相變溫度不是一個常數(shù)，而是與孔隙的尺寸有關。巖石孔隙中的相變可以由Gibbs-Thomson方程表示[18]

其中ri是被冰晶所侵入的最小入孔半徑，S是物質(zhì)的摩爾熵(J·mol-1·K-1)，Vi表示冰的摩爾體積(m3·mol-1)，γiw是冰水間的界面應力(0.04 J·m-2)，α是冰水間的接觸角(在這里取α=0)。Gibbs-Thomson方程揭示出在一定溫度下，水只有在大于某一臨界尺寸的巖石孔隙中才能夠發(fā)生相變，孔隙越小，其間發(fā)生相變所需的溫度越低。

由Gibbs-Thomson方程可以得到冰晶侵入孔徑尺寸與溫度的關系。但是需要注意的是，由于分離壓的存在，在孔隙中，冰晶與壁面會存在一層較薄的液膜。這層液膜的厚度同樣與溫度有關[18]

其中ξ是表示分子間力的作用的特征長度(在這里取ξ～2.3?)，D是 擴 散 系 數(shù)，并 且D=γsi-γsw-γiw(取D～ 0.33 J·m-2)。因此，發(fā)生相變的真實孔隙半徑的大小是要略大于ri的(如圖1)，并且有以下關系

圖1 孔喉半徑與溫度壓力間關系。r(T)：考慮液膜厚度的Gibbs-Thomson方程表示的溫度和冰晶侵入孔徑關系；r(P)：壓汞法中Washburn方程表示的汞壓與孔徑關系Fig.1 Relation between the smallest pore access radius and temperature.r(T): relation between temperature and pore radius by Gibbs-Thomson equation considering the premelting film (G-T+e); r(P): relation between the pore radius and Hg pressure in mercury intrusion porosimetry (MIP).

所以，如果已知巖石的孔徑分布，含冰量與溫度間的關系就可以由方程(4)～(6)給出。在實驗方法中 ，孔徑分布數(shù)據(jù)一般是由壓汞法取得，壓汞法的原理就是Washburn方程

其中p是汞壓(MPa)，σHg是汞的表面張力(0.48 mN·m-1)，αHg是汞與固體壁面的接觸角(取cosαHg=0.765)。

結(jié)合方程(4)～(7)，可以得到一定孔徑下相變溫度與汞壓間的關系，

由上式可以看出，含冰量隨溫度的變化可以用壓汞法的數(shù)據(jù)來代替，即通過方程(8)將汞壓轉(zhuǎn)換成溫度。然而，由于液膜的存在，壓汞法求得的孔隙體積分數(shù)(θHg=VHg/V)要比同一孔隙中的含冰量(θi=Vi/V)要大。假設孔隙都是圓柱體孔隙，對應一定孔徑的含冰量與真實孔隙體積之間的關系為:

應該注意的是，以上關于冰含量與溫度間關系的推導是基于在每個孔隙中凍結(jié)和融化現(xiàn)象都獨立存在這樣一個假設。然而，由于結(jié)冰過程是在水中形成冰核，在此過程中產(chǎn)生冰水界面需要消耗表面能。因此，結(jié)冰過程往往是由外部孔隙向內(nèi)逐步推進的，而不是在每個孔道中獨立發(fā)生。換句話說，結(jié)冰過程依賴于孔隙的連接方式，一個孔隙的凍結(jié)依賴于它的入口孔隙半徑大小?？紤]一個入口孔徑較小的大孔，溫度降低時大孔不會先凍結(jié)，而需要等到入口孔喉結(jié)冰后，大孔才會結(jié)冰。相反，由于液膜的存在，融化過程不涉及新表面的生成，因此可以在每個孔隙中獨立發(fā)生。因此，上述模型僅適用于不受孔隙連通性影響的融化過程。

2 儲層的溫度場模擬

在水合物開采過程中，井筒內(nèi)流體溫度相對較高，會對周圍儲層的溫度場產(chǎn)生擾動。針對這一問題，建立如圖2所示平面幾何模型。模型中心的井眼直徑設為200 mm，地層半徑為10 m且為均質(zhì)地層。注意到水合物開采過程中，井筒內(nèi)的流體通過熱對流方式將熱量傳遞給井筒壁面，其后再以熱傳導方式在儲層內(nèi)傳遞，所以外邊界設為自由邊界，內(nèi)邊界設為對流邊界。儲層初始溫度設為-4 °C，在使用隔熱管的情況下，井內(nèi)流體大部分熱量被隔絕，那么設等效的流體溫度為5 °C。

圖2 有限元網(wǎng)格模型Fig.2 Finite element mesh model

利用實驗室得到的一組砂巖的壓汞實驗數(shù)據(jù)求得其孔徑分布曲線(如圖3a所示，這里用巖石A來表示，峰值孔徑約為10 μm)。同時，為了考慮其孔徑分布對傳熱過程的影響，將巖石A的孔徑分布曲線向右平移得到巖石B的孔徑分布圖(峰值孔徑約為50 μm)。再通過方程(8)和(9)得到巖樣A和B孔隙內(nèi)含冰量與溫度間的關系(如圖3b所示)，巖石A和B的含冰量基本都在0 °C附近變化較快，推斷在此區(qū)間內(nèi)的相變反應劇烈，而巖樣B含冰量大于巖石A則是由于巖樣B的孔隙度較大。結(jié)合表1給出的各組分物質(zhì)的物性參數(shù)以及傳熱過程的相關系數(shù)，得到等效熱容值并代入，對地層的傳熱過程進行模擬。

表1 物性參數(shù)和傳熱系數(shù)[21]Table 1 Physical parameters and thermal coefficient in this work[21]

圖3 a.巖樣A和B的孔徑分布曲線；b.孔隙內(nèi)含冰量與溫度間關系Fig.3 a.The pore distribution curve of sample A and B; b.ice content vs temperature in pores

在模擬結(jié)果中，基于模型軸對稱的特點，以軸心為原點，選取任一半徑為X軸。在1小時到30天之間的不同時刻，模擬得到地層溫度場的演化云圖(圖4a)。取其沿X軸方向的地層巖石溫度分布，得到圖4b所示結(jié)果。在不同的t時刻，溫度都呈梯狀分布，沿X軸方向逐漸降低。同時，由于同一熱源，初始相鄰區(qū)域溫差較大，溫度的傳播更快，其溫度梯度更陡；隨著熱量向外的擴散，傳播距離逐漸增大，同時溫差逐漸減小，地層溫度分布梯度逐漸變緩，可以預測一定地層的溫度分布曲線是斜率逐漸變緩的過程。

圖4 a.溫度場演化；b.不同時刻地層溫度分布Fig.4 a.the evolution of temperature field; b.temperature distribution of permafrost at different moments

對于巖樣A和B，取井筒處同一點的溫度變化曲線(如圖5所示)。地層溫度隨時間逐漸升高，當達到零度附近時A，B曲線均出現(xiàn)比較平緩的一段，隨后溫度繼續(xù)上升。B曲線在A曲線上方，表明其溫度上升速率比A快。同時，假設溫度隨時間變化的速率小于0.01以視為溫度曲線的平段區(qū)域(dT/dt＜0.01)，圖中曲線B的平段要比曲線A長。

圖5 巖石A和B的溫度隨時間的變化關系Fig.5 Temperature vs time for sample A and B

3 分析與討論

在模擬結(jié)果中，相變潛熱對地層的溫度變化有極大的影響。在地層的傳熱過程中，孔隙內(nèi)冰不斷發(fā)生相變，平段區(qū)間的出現(xiàn)正是由于其孔隙中相變過程較為明顯，大量的熱量轉(zhuǎn)化為物質(zhì)的潛熱，導致溫度的上升極其緩慢。從控制方程(3)中也可以看出，在其他參數(shù)的一定的條件下，當含冰量隨溫度劇烈變化時，即潛熱項的突然增大，等同于增大了整體的熱容值，系統(tǒng)吸收相同的熱量時其溫度變化遠沒有之前明顯。對于巖樣A和B，在相似的孔徑分布下，巖樣B的孔徑較大，其所對應的相變溫度就較大，即在相同的熱源下，巖樣A比B更早的發(fā)生相變，即B曲線前期的升溫速率要大于A。同時，在2種巖樣的含冰量隨溫度的變化圖中，零度附近B曲線的含冰量隨溫度的變化曲率大于A曲線，即在較小的溫度區(qū)間內(nèi)，巖樣B孔隙中的含冰量變化等效增加的體系熱容較大，所以升溫更慢，溫度曲線的平段更長。

由上可知，曲線的平段在模型中是由含冰量隨溫度的變化決定的，而含冰量隨溫度變化的關系是由壓汞實驗的孔徑分布數(shù)據(jù)獲得的，那么，通過測量材料在低溫下的溫度變化曲線，根據(jù)溫度曲線平段的位置和長度就可以得到物質(zhì)的孔徑分布。值得注意的是雖然壓汞法是實驗室中一種較為常見的測量孔隙分布的方法，它的測量范圍廣，操作簡單，但在理論模型中，利用壓汞法來反求出巖石孔隙中的含冰量與溫度間關系進行了兩次變換，有液膜和溫度壓力間這兩個過程轉(zhuǎn)換誤差的影響。因此，在條件允許的情況下，用更加精確的方法測量巖石孔徑分布或者直接測出巖石孔隙中含冰量隨溫度的變化關系都可以對模型進行優(yōu)化。

4 結(jié)論

針對水合物開發(fā)過程由于相變現(xiàn)象吸收大量的熱，從而會造成儲層溫度場的擾動這一現(xiàn)象，類比冰水相變進行了儲層溫度場的擾動模型研究，并研究了不同孔徑分布對儲層溫度隨時間變化的規(guī)律的影響。主要得到了以下結(jié)論：

(1) 建立了含相變過程的傳熱控制方程。其中，將潛熱和比熱2種能量轉(zhuǎn)化形式合并為等效熱容，從而使控制方程簡化為類似于不含相變過程的經(jīng)典Poisson方程。另外，通過壓汞數(shù)據(jù)可以預測含冰量與溫度的關系。

(2) 當儲層溫度升高時，其地層溫度變化曲線會在0 °C附近出現(xiàn)一段緩慢上升的階段。這主要是因為在這一階段相變劇烈，需要吸收大量熱量，從而使溫度變化并不明顯。

(3) 儲層的溫度演化受巖石本身的孔隙結(jié)構(gòu)影響巨大。當孔隙度較大時，由于骨架的缺失，巖石的整體熱容會相應減少，從而使溫度變化速率升高。同時，孔徑分布也對溫度平段的位置和長度產(chǎn)生。反之，通過測量儲層的溫度變化曲線，也可以反演得到巖石的孔徑分布。