趙振偉，尚新春

(1.北京科技大學(xué)土木與環(huán)境工程學(xué)院，北京 100083；2.北京科技大學(xué)應(yīng)用科學(xué)學(xué)院，北京 100083)

天然氣水合物是由天然氣和水在一定溫度和壓力條件下形成的一種冰狀晶體，其廣泛存在于永久凍土帶和海底沉積物孔隙當中。天然氣水合物的含碳量大約是地球上其他所有常規(guī)能源含碳量的兩倍，被認為是未來一種理想能源[1-2]。天然氣水合的分解，主要受壓力和溫度控制，研究水合物開采過程中溫度場和壓力場的分布，能夠在實際生產(chǎn)中確定水合物的分解范圍，預(yù)測天然氣的產(chǎn)量。通過對天然氣水合物分解面位置的分析，可以找出影響分解面移動速度的因素，為提高天然氣的產(chǎn)出效率提供可行的方法。對天然氣水合物開采過程中溫度場和壓力場的分析，是近20年多來的研究熱點。

作者最近研究了水合物開采過程中的壓力場和溫度場，在求解壓力場時只考慮了水的滲流，于是壓力場控制方程為線性方程[9]。本文在文獻[1]的基礎(chǔ)上，進一步考慮了水合物相變潛熱的影響，在分解面上加入能量守恒方程，對水合物開采過程中溫壓場的軸對稱模型進行了改進，并通過波爾茲曼變換方法，得到溫度場壓力場耦合模型的半解析-半數(shù)值解。在確定分解面位置時，將分解面上的質(zhì)量守恒方程轉(zhuǎn)化為一常微分方程來求解。與文獻[2]不同的是，本文考慮了氣體的流動，壓力場為非線性方程，并在分解面上考慮了氣體密度隨壓力和溫度的變化。

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

考慮水合物區(qū)一豎直井，通過降壓法開采水合物。假設(shè)壓力場和溫度場以井筒垂直軸線為中心成軸對稱分布且與深度無關(guān)。引入柱坐標系(r,θ,z)。本文將水合物的開采過程看作是一個移動邊界問題，在水合物開采過程中，礦區(qū)被分為兩個區(qū)域：分解區(qū)r0≤r≤R和未分解區(qū)R≤r≤+∞，r=R(t)表示移動邊界的位置。記Pn(r,t)為壓力場和Tn(r,t)為溫度場。本文各個力學(xué)量的下標n=1、2，分別代表分解區(qū)和未分解區(qū)。

Goodrz Ahmadi(2007)給出了如下壓力場的非線性控制方程[3]：

(1)

(2)

因沉積物的導(dǎo)熱率較低，這里忽略了熱傳導(dǎo)效應(yīng)，僅考慮對流效應(yīng)，其溫度場控制方程為[3]：

(3)

式中，常數(shù)δ和η分別為孔隙中氣體的節(jié)流系數(shù)和絕熱系數(shù)；cv和c分別為氣體和沉積物的熱容。

假設(shè)在井口以及無窮遠處壓力和溫度均為常數(shù)，則邊界條件如下：

(4)

(5)

并且初始條件如下：

(6)

(7)

在分解面上水合物因分解而產(chǎn)生氣體，其質(zhì)量守恒方程為[5]：

(8)

其中氣體體密度滿足如下關(guān)系[8]：

(9)

T0=273.15K；P0=1.01×105Pa；ρ0是氣體在T0、P0條件下的密度；z=0.88，是氣體的壓縮系數(shù)；TD、PD是水合物相變時的溫度和壓力。

假設(shè)氣體流速滿足滲流的達西定律：

(10)

將式(9)和式(10)帶入質(zhì)量守恒方程可得到：

(11)

其次，在分解面上，由于水合物分解會吸收很大部分熱量，所以在分解面上的能量平衡方程時，應(yīng)考慮水合物相變潛熱的影響。而沉積物的導(dǎo)熱率相對較低，在此忽略熱傳導(dǎo)的影響，可以推導(dǎo)出考慮水合物相變影響的的能量守恒方程如下：

(12)

式中：L是水合物的相變潛熱；ρh是水合物的密度。

將式(9)和式(10)代入式(12)可得：

(13)

最后，在分解面上，水合物的溫度壓力滿足水合物的相平衡條件[2]：

(14)

式中：a=0.0342K-1；b=0.0005K-2；c=6.4804。

在分解面上需滿足壓力和溫度的連續(xù)性條件：

(15)

(16)

至此，確定水合物分解過程的壓力場Pn(r,t)和溫度場Tn(r,t)(n=1,2)在數(shù)學(xué)上被歸結(jié)為：在分解區(qū)r0≤r≤R(t)和未分解區(qū)R(t)≤r≤+∞上, 求解偏微分方程組(2)和(3)的解, 使其滿足初-邊值條件(4)～(7), 同時滿足移動邊界R=Rc(t)上的質(zhì)量和能量守恒方程(11)和(13) 以及連續(xù)性條件(15)～(16)。

2 問題的求解

引入波爾茲曼變換[4]

(17)

將式(17)代入方程(2)可得到壓力場的常微分方程表達形式：

(18)

(19)

利用式(17)，邊界條件(4)、初始條件(6)和界面連續(xù)性條件(15)可轉(zhuǎn)換為：

(20)

(21)

(22)

對方程(18)、(19)進行兩次積分，再利用邊界條件(20)～(22)，可得到壓力場的解析表達式：

(23)

(24)

同樣，將式(17)以及壓力場解析表達式(23)、(24)代入方程(3)，可得溫度場的常微分方程表達形式：

(25)

(26)

利用式(17)，邊界條件(5)初始條件(7)及界面連續(xù)性條件(16)可轉(zhuǎn)換為：

(27)

(28)

(29)

結(jié)合邊界條件(27)～(29)，求解常微分方程(25)、(26)，可得到溫度場的解析表達式：

(30)

(31)

上述壓力場和溫度場的表達式中，PD、TD及與分解面位置有關(guān)的uc為不確定量，可通過界面質(zhì)量守恒方程(11)和能量守恒方程(13)及水合物的相平衡方程(14)來確定。在確定上述三個不定量時，首先將方程(11)和(13)相減可得到：

(32)

將方程(14)代入到式(32)，可得關(guān)于TD的非線性方程，利用數(shù)值方法求解出TD，再由方程(14)確定出PD。將壓力場Pn的表達式(23)、(24)代入能量守恒方程(13)，得到如下關(guān)于Rc(t)常微分方程的初值問題：

(33)

Rc(t0)=r0

通過數(shù)值方法求解上述常微分方程的初值問題(33)，即可確定t時刻的分解面位置Rc(t)。

3 數(shù)值算例

下面給出溫度場和壓力場分布以及分解面位置的數(shù)值計算結(jié)果，計算參數(shù)如下：

c=2000 J/kg,cv=3000 J/kg,r0=0.1 m,

α=0.08,β=0.1,L=43500 J/kg,T0=273.15K,ε=0.129,

P0=1.01×105Pa,ρ0=0.706 kg/m3,Ρh=910 kg/m3,δ=8.0×10-7K/Pa,η=3.2×10-6K/Pa,K1=5.2md,K2=0.4md,Φ= 0.2。

圖1給出了井筒壓力為4MPa、分解時間分別為30d、60d和90d時，壓力和溫度沿水平徑向的分布。由圖1中可以看出，由于水合物開采時井筒壓力急劇降低，井筒處的壓力變化梯度較大；在未分解區(qū)，分解面附近的壓力最小，向無窮遠處逐漸增大到地層的初始壓力，并且壓力梯度越來越小。在同一水平位置，壓力隨著時間逐漸降低，但是對時間變化率卻在減小。由圖2可以看出，由于這里考慮了水合物相變潛熱的影響，水合物分解會吸收部分熱量，溫度的最低值在分解面附近，向兩邊逐漸過渡到地層的初始溫度，分解面附近溫度梯度較大。

圖1 不同時間壓力分布

圖2 不同時間溫度分布

圖3給出了不同的井筒壓力條件下，水合物的分解面位置距井筒中心的水平距離與時間的關(guān)系。從圖中可以看出，在井筒壓力一定的條件下，分解距離隨時間逐步增大，但對時間的變化率卻在不斷減小。在開始時刻，分解前緣移動速度較快，隨著分解時間的進行，分解前緣的移動速度會逐漸變慢。通過對不同井筒壓力下分解距離的比較可以看出，井筒壓力越低，分解速度越快。當時間為90d時,井筒壓力從5MPa降低到3MPa，分解距離從7m增加到24m。由此可見，水合物的分解速度對井筒壓力比較敏感，所以在開采過程中，應(yīng)盡量降低井筒壓力以加快水合物的分解。

圖3 分解面距井筒中心的水平距離與時間的關(guān)系

通過壓力場的解析表達式，單位高度上，天然氣水合物開采時的產(chǎn)氣速率，可通過井筒處氣體的流速確定：

(34)

通過氣體產(chǎn)生速率的表達式(34)，可以得到降壓法開采水合物時產(chǎn)氣速率與時間的關(guān)系。圖4給出了井筒壓力為3MPa、4MPa、5MPa時，氣體產(chǎn)生速率隨時間的變化曲線。從圖4中可以看出，產(chǎn)氣速率隨時間的增加在減小，開始時減小幅度較大，但是經(jīng)過一段時間，產(chǎn)氣速率逐步趨于一穩(wěn)定值。隨著井口壓力的降低，氣體產(chǎn)氣速率會明顯提高。

圖4 產(chǎn)氣速率與時間的關(guān)系

4 結(jié) 論

(1)通過水合物分解過程中的壓力場和溫度場分布情況看，在井口處壓力梯度最大，壓力值最?。浑S著徑向坐標的增大，壓力逐漸變大并趨于地層的初始壓力；受水合物分解的影響，分解面附近溫度較低，溫度梯度較大。

(2)水合物開始分解時，天然氣的產(chǎn)氣速率隨時間在逐步遞減，當經(jīng)過一段時間后，產(chǎn)氣率趨于一穩(wěn)定值。

(3)天然氣的產(chǎn)氣速率受井口壓力的影響明顯，通過降低井口的壓力，可以明顯提高天然氣的產(chǎn)出速率。在天然氣水合物的實際開采中，可以通過降低井口壓力來提高生產(chǎn)效率。

[1] Ahmadi G, Ji C, Smith DH. Production of natural gas from methane hydrate by a constant down-hole pressure well[J]. Energy Conversion and Management，2007，48(7)： 2053-2068.

[2] Zhenwei Zhao，Xinchun Shang. Analysis for temperature and pressure fields in process of hydrate dissociation by depressurization[J]. International Journal for Numerical and Analytical Methods in Geomechanics，2010(3).

[3] Ji C, Ahmadi G, Smith DH. Natural gas production from hydrate decomposition by depressurization [J]. Chemical Engineering Science，2001，56(20)：5801-5814.

[4] M.Janz. Moisture diffusivities evaluated at high moisture levels from a series of water absorption tests [J]. Materials and Structures，2007(35)：141-148.