鮑祥生 侯 斌

(1.廣東石油化工學院 廣東茂名 525000；2.油氣藏地質及開發(fā)工程國家重點實驗室(成都理工大學)四川成都 610059；3.江蘇油田 江蘇揚州 225009)

天然氣水合物是由水和天然氣在低溫高壓條件下形成的一種籠形化合物，是一種類冰狀的物質[1],該物質沒有固定的化學分子式，是一種非化學計量的混合物[2]，主要分布在陸上的永久凍土帶和海域的海底沉積物中[3]。水合物燃燒的產物是二氧化碳和水，因此屬于一種高效潔凈的非常規(guī)能源，是未來化石能源的重要接替能源[4]。水合物大部分分布于海域，全世界范圍內海域水合物分布占到了98%，研究表明水合物主要分布在受到溫壓條件控制的穩(wěn)定帶內[5-11]。但穩(wěn)定帶只是決定水合物存在的空間范圍，并不對水合物富集程度起決定性影響。水合物在穩(wěn)定帶某個地區(qū)富集還受到氣源條件、氣體運移條件、儲集條件等因素制約。

目前預測水合物富集有許多地震方法[1,5-6,8-10,12-13]，如BSR法、振幅空白帶法、屬性預測法、波阻抗預測法、水合物飽和度預測法等。依據(jù)量化特征可將這些方法分為3類：①定性預測方法，如BSR法、振幅空白帶法，該類方法具有使用簡單直觀、易于理解等優(yōu)點，但定量化程度不夠，不利于判斷不同位置富集程度差異；②半定量預測方法，如屬性預測法，該類方法具有一定的量化特征，部分方法已能初步反映出富集程度差異性，但富集程度差異也是相對的；③定量預測方法，如波阻抗預測法、水合物飽和度預測法，該類方法能夠絕對反映水合物的富集程度差異。在這3類方法中，定量預測方法是實際應用中預測水合物富集最受歡迎的一類方法，其中水合物飽和度預測法與水合物資源量計算最相關[14-15]，所以該方法對一個地區(qū)的水合物商業(yè)化進程具有重要意義，目前針對該方法的研究已成為水合物定量預測的熱點課題。由文獻[14-19]可知，Timur 時間平均方程法(簡稱Timur法)、Wood 修正方程法(簡稱Wood法)、Gassmann法是目前被用來預測水合物飽和度最經典的3種方法，但有各自的適用條件[20]。

水合物在海域沉積物中以3種主要賦存方式[19-21]存在，即懸浮模式、顆粒接觸模式及膠結模式。依據(jù)胡高偉 等[22]研究成果可知，我國神狐海域富集水合物[23]是以懸浮模式賦存的。疊后地震反演可直接反演出波阻抗，不同飽和度水合物會導致波阻抗發(fā)生變化，盡管目前Wood法被認為可以較好地預測懸浮模式的水合物飽和度[20]，但并沒有給出不同飽和度水合物與波阻抗之間的關系，因此很難直接將其應用到實際的地震資料中來進行水合物飽和度預測。針對這一問題，本文在對Wood法進行推導分析簡化的基礎上，形成了一種用于懸浮模式水合物飽和度預測的新方法，即Wood阻抗法。

1 基于Wood阻抗法的水合物飽和度計算原理

1.1 Wood法

求取水合物飽和度的Wood法由2個方程構成[17,19]，即

(1)

ρb=(1-Sh)φρw+φShρh+(1-φ)ρm

(2)

式(1)(2)中:Vb、Vpw、Vph、Vpm分別表示含水合物沉積物、海水、純水合物、基質的縱波速度；φ表示孔隙度；Sh表示水合物在孔隙空間中所占的比例;ρb、ρw、ρh、ρm分別表示含水合物沉積物、海水、純水合物、基質的密度。

1.2 基質密度和縱波速度計算

一般來說,構成海底沉積物的基質由多種物質組成，其密度表示為

(3)

式(3)中:ρm為基質的密度；fi為第i種礦物組分在基質中的體積百分比；ρi為組成基質中第i種物質的密度;n表示組成基質所包含物質的數(shù)量。

基質的縱波速度可表示為[24]

(4)

式(4)中：Km、Gm分別為基質的體積模量和剪切模量。

基質的體積模量和剪切模量可依據(jù) Hill平均公式來求取[25],即

(5)

(6)

式(5)中：Ki為第i種礦物組分的體積模量；Km為基質的體積模量；Gi為第i種礦物組分的剪切模量；Gm為基質的剪切模量。

1.3 純水合物及海水的縱橫波速度計算

依據(jù)文獻[24]中有關介質彈性參數(shù)之間的關系，對于純水合物的縱橫波速度可直接依據(jù)彈性參數(shù)并利用式(7)、(8)來進行計算。

(7)

(8)

式(7)、(8)中：E為楊氏模量；ρ為物質的密度；σ為泊松比；Vph為縱波速度；Vs為橫波速度。

對于楊氏模量和泊松比，可利用物質的體積模量和剪切模量來進行求取,即

(9)

(10)

式(9)、(10)中：K、G分別為物質的體積模量和剪切模量。

對于海水的縱波速度，一般取1 480 m/s;由于海水為液體，橫波在其中無法傳播，所以海水中的橫波速度為0。

1.4 Wood阻抗法推導

將式(1)兩邊同時乘以1/ρb，則有

(11)

令含水合物沉積物的縱波阻抗Zb=ρbVb,海水的縱波阻抗Zw=ρwVpw,純水合物的縱波阻抗Zh=ρhVph；基質的縱波阻抗Zm=ρmVpm,則式(11)可表示為

(12)

由文獻[11]、[26]可知,海水的密度一般可取1.03 g/cm3,純水合物的密度一般可取0.9 g/cm3。由式(2)可知

ρb/ρw=(1-Sh)φ+φShρh/ρw+(1-φ)ρm/ρw

≈(1-Sh)φ+0.87φSh+0.97(1-φ)ρm

=0.97(1-φ)ρm+φ-0.13φSh

(13)

ρb/ρh=(1-Sh)φρw/ρh+φSh+(1-φ)ρm/ρh

≈1.15(1-Sh)φ+φSh+1.11(1-φ)ρm

=1.11(1-φ)ρm+1.15φ+0.15φSh

(14)

ρb/ρm=(1-Sh)φρw/ρm+φShρh/ρm+(1-φ)

≈(1-φ)+1.03φ/ρm-0.13φSh/ρm

(15)

令Cbw=ρb/ρw,Cbh=ρb/ρh,Cbm=ρb/ρm,則有Cbh>Cbw>1>Cbm。由于ρb一般大于1.5 g/cm3，而基質密度比ρb大，但最大一般不超過3 g/cm3，所以最小的系數(shù)Cbm都大于0.5;而φSh一般在0.1左右，式(13)～(15)中φSh所在項往往在0.01左右，因此相對于大于0.5的數(shù)值而言，該數(shù)值可忽略不計。鑒于φSh所在項可以忽略，則Cbw、Cbh、Cbm可近似表示為

Cbw≈0.97(1-φ)ρm+φ

(16)

Cbh≈1.11(1-φ)ρm+1.15φ

(17)

Cbm≈(1-φ)+1.03φ/ρm

(18)

由公式(16)～(18)可知，Cbw、Cbh、Cbm可被認為是與孔隙度和基質密度有關的系數(shù)。

將系數(shù)Cbw、Cbh、Cbm代入式(12)，可得

(19)

式(19)即為本文提出的Wood阻抗法求取水合物飽和度公式。

由此可見，在已知一個地區(qū)含水合物沉積物的縱波阻抗、海水的縱波阻抗、純水合物的縱波阻抗、基質的縱波阻抗、孔隙度以及基質密度的情況下，則可求出水合物飽和度。其中，含水合物沉積物的縱波阻抗可以通過縱波阻抗反演求得，海水的縱波阻抗、純水合物的縱波阻抗則依據(jù)彈性參數(shù)測量可以獲得，基質的縱波阻抗、孔隙度以及基質密度一般依據(jù)一個地區(qū)鉆井資料可以獲得。

2 Wood阻抗法實驗誤差分析

為驗證Wood阻抗法的可靠性，這里以一組數(shù)據(jù)為例來進行誤差分析。在應用方法前先做一些基本的數(shù)據(jù)假設：

1)假設海域含水合物沉積物的基質由粉砂和黏土2種物質組成，它們在基質中的比例分別是75%、25%；

2)假設懸浮模式水合物飽和度為30%；

3)假設水合物沉積物的孔隙度為40%。

表1列出了一些研究所涉及物質的彈性參數(shù)[11,26]。由表1和式(3)可算出基質的密度約為2.63 g/cm3；由表1和式(5)、(6)可算出基質的體積模量和剪切模量分別約為33.94、19.32 GPa；由基質的密度、體積模量和剪切模量，利用式(4)可算出基質的縱波速度約為4762.34 m/s；由表1和式(9)、(10)可算出水合物的楊氏模量和泊松比分別約為6.3 GPa、0.31；由表1和式(7)以及計算的楊氏模量和泊松比可算出水合物的縱波速度為3 126.94 m/s；由式(16)～(18)可算出Cbw、Cbh、Cbm分別約為1.93、2.21、0.76；由式(2)可算出含水合物沉積物的密度約為1.97 g/cm3；由式(1)可算出含水合物沉積物的縱波速度約為1 855.96 m/s；結合已算出的數(shù)據(jù)可知含水合物沉積物的縱波阻抗Zb、流體的縱波阻抗Zw、純水合物的縱波阻抗Zp、基質的縱波阻抗Zm約為3 651.61、1 524.4、2 814.25、12 536.88(m·s-1)(g·cm-3)。由式(19)可算出預測的水合物飽和度是28.5%，與實際已知值30%非常接近，說明本文提出的Wood阻抗法具有一定的適用性。

表1 物質的彈性參數(shù)Table 1 Elastic parameters of matter

3 應用效果分析

3.1 基礎數(shù)據(jù)準備

以神狐海域數(shù)據(jù)為參照建立二維模型。陳芳 等[27]指出,神狐海域水合物樣品附近沉積物的組成按照粒度大小分為砂、粉砂、黏土等3種;沉積物中最主要的顆粒類型是粉砂，其次是黏土，再次是砂，其中粉砂的平均含量大致占75%，黏土的平均含量大致占20%，砂的平均含量大致占5%。本次應用中，參照上述沉積物粒度組成分析成果，設置用于構成沉積物基質的各物質比例是：粉砂75%、黏土20%、砂5%。

圖1、2分別是神狐海域鉆獲水合物的S 井的縱波速度和密度曲線。從圖1可以看出，該井縱波速度曲線可大致分為3段:第1段為低速層，深度介于100～198 m，縱波速度一般低于1 900 m/s；第2段為高速層，深度介于 198～218 m，縱波速度明顯增加，最高達到 2 300 m/s，已被確認為含水合物層[11]；第3段為中低速度層，深度介于198～245 m，縱波速度相對含水合物層出現(xiàn)了明顯的減小，大致在1900 m/s左右。從圖2可以看出，該井密度曲線也可大致分為3段：第1段為中低密度段，深度介于100～200 m，密度主要在1.98 g/cm3左右；第2段為低密度段，深度介于200～210 m，密度主要在1.90 g/cm3左右；第3段為高密度段，深度介于210～245 m，密度主要在2.03 g/cm3左右。從圖1、2的對比可以看出，水合物高速層并不完全表現(xiàn)為低密度，210 m之上的水合物層密度主要表現(xiàn)為低值，但大于210 m時主要表現(xiàn)為高值。陳芳 等[27]對S井開展了沉積物組成成分分析，指出水合物層下段高密度處的砂巖及粗砂巖含量明顯高于上部低密度處的砂巖及粗砂巖含量，水合物層下段高密度處的黏土含量要低于上部低密度處的黏土含量，因此推測水合物層表現(xiàn)不同的密度差異性可能主要是巖性差異造成的。盡管水合物層不同位置的密度存在一定的差異性，但一般來說，密度與縱波速度相比對縱波阻抗的貢獻要小得多，因此時常可在研究時忽略密度變化對阻抗的影響。

圖1 S井縱波速度隨深度的變化Fig.1 Variation of P-wave velocity with depth in Well S

圖2 S井密度隨深度的變化Fig.2 Variation of density with depth in Well S

結合圖1，本文將低速段稱為上覆非水合物層，高速段稱為水合物層，中低速段稱為下伏非水合物層，表2列出了S井不同層段密度和速度的平均取值。

表2 S井不同層段密度和速度的平均值Table 2 Average density and velocity in different sections of Well S

3.2 二維地質模型描述

圖3 是依據(jù)海水彈性參數(shù)和表2數(shù)據(jù)構建的二維地質模型，該模型共包括4層，第1、2、3、4層分別是海水層、上覆非水合物層、水合物層和下伏非水合物層，其縱波速度和密度數(shù)據(jù)見表3。依據(jù)式(1)及表1、2可估算第2層和第4層的非水合物層孔隙度分別約為36%、28%。一般來講，海底沉積物中的孔隙度往往隨深度的增加逐漸減小，模型中水合物層介于上覆非水合物層和下伏非水合物層之間，所以其孔隙度一般應介于這2個非水合物層的孔隙度之間，因此本文取這2個非水合物層的孔隙度平均值32%作為水合物層(即第3層)的孔隙度。

圖3 含水合物層的二維地質模型Fig.3 2D geological model of gas hydrate-bearing sediments

表3 含水合物層的二維地質模型的縱波速度和密度數(shù)據(jù)Table 3 P-wave velocity and density data for 2D gelogical model of gas hydrate-bearing sediments

模擬海上產生零偏移距地震記錄，假設震源子波是雷克子波，地震子波主頻為50 Hz，則圖3 對應的零偏移距記錄如圖4所示。圖5為依據(jù)圖4振幅特征解釋的3個層位，從上到下依次為B1、B2、B3，其中B1代表海底的反射，具有正極性反射特征且振幅最強；B2代表上覆非水合物層與水合物層之間的反射，具有正極性反射特征且振幅相對較強；B3代表水合物與下伏非水合物層的反射，具有負極性反射特征且振幅相對較弱。圖6是反演的縱波阻抗剖面，B2和B3層位之間水合物層具有較高的縱波阻抗，在4 100 (m·s-1)(g·cm-3)左右?；趫D6，再依據(jù)式(19)預測出水合物飽和度(圖7)，預測的水合物飽和度大致介于28%～38%。而Zhang Haiqi等[23]利用氯離子法估算水合物含量富集層段水合物含量主要介于30%～40%，可見本文方法預測結果可信，說明該方法具有一定的適用性。

圖4 二維地質模型對應零偏移距地震記錄Fig.4 Zero-offset seismic records of 2D geological model

圖5 帶有解釋層位的二維地質模型對應零偏移距地震記錄Fig.5 Zero offset seismic records of 2D geological model with interpreted horizons

圖6 反演的縱波阻抗剖面Fig.6 P-wave impedance inversion profile

圖7 反演的水合物飽和度剖面Fig.7 Hydrate saturation inversion profile

4 結論與建議

1)本文提出的Wood阻抗法是一種基于縱波阻抗預測水合物飽和度預測方法，該方法克服了Wood法在實際疊后反演中進行應用難的問題。

2)實驗誤差分析和以神狐海域實際資料應用表明，本文方法具有較好的實用性。

3)孔隙度對于本文方法應用非常重要，因此建議在應用本文方法時要對海底沉積物孔隙度分布規(guī)律進行詳細研究，以有利于取得更好的預測效果。