基于核磁共振T2譜的頁巖巖心孔隙分布量化表征方法

靳軍， 劉偉洲， 王子強， 郭慧英， 李瓊， 李震， 白寧晨， 劉金玉*

(1.中國石油新疆油田分公司實驗檢測研究院， 克拉瑪依 834000； 2.新疆頁巖油勘探開發重點實驗室， 克拉瑪依 834000； 3.中國石油大學(華東)理學院， 青島 266580)

核磁共振(nuclear magnetic resonance，NMR)技術作為一種無損檢測方法已被廣泛應用于表征儲層巖石的孔隙分布[1-2]，其孔隙分布表征原理是基于巖心中油水流體的核磁共振信號，而流體飽和頁巖中除油水流體外，還含有同樣能產生核磁共振信號的干酪根和黏土礦物結構水(類固體)，即頁巖所測的核磁共振信號是巖心中所有發生核磁共振的物質信號的疊加，因此直接利用流體飽和頁巖的核磁共振T2(T2為橫向弛豫時間)譜(即NMR-OS方法)來表征頁巖孔隙分布存在誤差[3- 4]，為了減少誤差就需要將類固體的干擾信號去除，目前主要通過實驗方案改進和T2譜擬合兩種方式去干擾。Li等[3]通過改進實驗方案去除干擾信號，其在巖心洗油烘干前后進行兩次核磁共振測量，利用兩條曲線相減去除干擾，再采用NMR-O法分析頁巖孔隙的分布，這種方法在一定程度上能夠提高孔隙分布測量的精度，但需要洗油、烘干等一系列煩瑣的處理過程，又由于洗油、烘干不能完全去除頁巖納米孔隙中的流體，且會破壞原有的孔隙結構和分布，導致孔隙分布表征仍然不夠精確。在T2譜擬合研究方面，研究者將每一種物質的核磁共振信號用一個峰模型來表示，基于該假設，可將實驗所測流體飽和頁巖的重疊譜分解為多個峰模型。文獻[5-7]中已經建立了T2譜的幾個峰模型，這些模型選用對稱分布的正態函數，對稱峰模型適用于描述純水等單組分流體的T2分布，然而在頁巖中含有多種組分，由于不同組分之間的化學或物理相互作用，會產生不對稱峰，再使用對稱峰模型會產生較大誤差，另外正態分布還存在峰值衰減慢，導致基線遠離峰值[8]，因此選用對稱峰模型不適用于頁巖的T2譜分析。

鑒于以上分析，從峰模型改進入手，結合求導分析和最小二乘法，設計一種新的頁巖孔隙分布分析方法，實現從流體飽和頁巖的核磁共振重疊譜中分離出類固體信號，進而做到快速、準確地表征頁巖孔隙分布，并對所提出的方法進行驗證和開展實例應用。

1 核磁共振-反褶積方法

由核磁共振基本理論可知，原始T2譜A(T2)可視為n個獨立分量的T2分布的線性疊加[9]。假設第k個獨立分量的T2分布用解析峰函數fk(T2)表示，則原始T2譜A(T2)可表示為

(1)

在確定峰函數fk(T2)后，可通過式(1)反演計算得到原始T2譜A(T2)，這個過程稱為“反褶積”。通過反褶積分析，核磁共振T2譜對頁巖中不同相的表征可以轉化為對這些峰函數參數的確定。核磁共振-反褶積(NMR-deconvolution,ND)方法的分析流程如圖1所示。

圖1 核磁共振-反褶積方法的分析流程

1.1 峰模型選擇

為了避免對稱峰模型的缺點，提出了適用于頁巖T2譜的Beta峰模型f(x,f0,a0,μ,s1,s2,w)[10]，即

f(x,f0,a0,μ,w,s1,s2)=

(2)

式(2)中：f0為基線(在T2譜中通常f0=0);a0為峰高；μ為峰值中心；s1和s2為形狀因素，s1和s2均需大于1；w為峰寬；x1和x2分別為左、右端點，其表達式分別為

(3)

(4)

圖2給出了一個對數尺度上的Beta峰模型。當s1=s2時，Beta峰為對稱峰；當s1≠s2時，Beta峰不再對稱。因此，Beta峰模型可以同時擬合對稱峰和非對稱峰，且峰的形狀范圍很廣，具有很高的靈活性。另外，該模型是一個分段函數，它由一個Beta分布和一個常數函數y=0組成，這使得峰的前緣或后緣的衰減速率更快，可有效避免基線遠離峰情況的發生。

圖2 對數尺度上的Beta峰模型

1.2 峰檢測

峰檢測的目標是從頁巖T2譜中識別出潛在的峰，可通過導數分析對T2譜進行峰檢測。對T2譜進行平滑六階求導后，直接將局部極值識別為峰。與其他峰檢測方法[10-14]相比，該方法不需要關于圖譜的先驗知識，進而實現簡單、快速的峰檢測。一個典型的純水峰及其平滑六階導數如圖3(a)所示，單峰被成功地識別為峰值范圍內平滑六階導數的局部最小值。由兩個單峰嚴重重疊構成的T2譜及其平滑六階導數如圖3(b)所示，在峰值范圍內六階導數中出現的兩個局部極小值，表明存在兩個峰值。對于更多數量峰重疊的T2譜，可以采用同樣的方法進行峰識別。適當的平滑是成功實現峰檢測的先決條件，選擇Savitsky-Golay平滑方法，該方法通過多項式擬合指定窗口內的數據來進行平滑，它可以通過改變窗口大小實現不同的平滑效果，本研究中窗口大小統一設置為48。

圖3 T2譜峰檢測的導數分析方法

1.3 峰參數確定

通過峰檢測確定峰的個數后，下一步確定每個峰的參數。在原始T2譜峰值范圍內的平滑六階導數中，由于峰被識別為局部極小值，因此峰值中心μ可以估計為

(5)

約束條件為

(6)

a0的估計可表示為

(7)

約束條件為

(8)

平滑六階導數局部極大值一般對應峰值的結束[圖3(a)]，所以峰寬w可以估計為

(9)

由于單個峰的寬度不應超過全譜wspectrum的寬度，因此w的約束條件為

0≤w≤wspectrum

(10)

形狀因子s1、s2在擬合過程中可靈活調節，但要求值應大于1，作為初始估計，可設置為4。

在給定峰參數的初始估計和約束條件后，最后采用最小二乘法確定峰參數的具體值。

(11)

式(11)中：Y為原始T2譜譜與擬合T2譜誤差的絕對值；A為原始T2譜；fk為單峰函數，其中不受約束的峰值可以變為負振幅，或者減小到零寬度。

2 核磁共振-反褶積方法驗證

為了驗證提出的核磁共振-反褶積方法，從準噶爾盆地吉木薩爾凹陷二疊系蘆草溝組選取2塊頁巖巖心進行洗油烘干(EX)和飽和油(OS)，并測量洗油烘干和飽和油后頁巖的核磁共振T2譜，核磁共振實驗采用拉莫爾LMR-1600型低場核磁共振分析儀，儀器的中心頻率為(7±0.1) MHz，樣品腔的直徑為30 mm，樣品測量高度為20 mm。為了更好驗證該方法，又從Li等[3]的研究成果中選取了4塊頁巖巖樣的核磁共振T2譜，6個頁巖樣品的核磁共振T2譜如圖4所示。

圖4 頁巖樣品的T2譜

在樣品2中，EX的T2譜呈單峰分布，弛豫時間小于1 ms。在其他樣品中，EX的T2譜大于1 ms時出現了后續峰，其中一些峰相互疊加。其中弛豫時間小于1 ms的信號主要由類固體(干酪根和結構水)產生，大于1 ms的峰為洗油烘干后殘余流體的信號。由于類固體的弛豫行為受到固體或半固體分子內偶極耦合的嚴重影響，可以用Bloembergen-Purcell-Pound理論描述，即

(12)

式(12)中：C為常數；ω為拉莫爾頻率；τ為相關時間，由類固體本身的性質決定。

對于特定環境中的類固體，其τ值一般集中于某一個定值附近[15]。結合樣品2和樣品3洗油烘干后的T2譜可知，類固體的T2分布呈現單峰分布。

與EX的核磁共振信號相比，OS的核磁共振信號強度明顯增加。基于上面分析，可以認為OS的T2譜中大于1 ms的峰信號全部由油貢獻，應計入孔隙分布；小于1 ms的峰信號則由類固體和油共同貢獻。在小于1 ms的T2譜峰中，弛豫時間最短的第一個峰(T2=0.22 ms)一般是類固體信號，不能用于評價孔隙分布；而其他峰信號可能來自于吸附油的貢獻，可以計入孔隙分布表征。因此，想要通過核磁共振獲得精確的孔隙分布信息，需要從OS頁巖T2譜中剝離出類固體信號。

2.1 類固體信號分離

采用ND方法從油飽和頁巖樣品的T2譜中分離出類固體產生的信號。如圖5所示，所有樣品T2譜平滑六階導數的最左側均存在一個局部最小值，表明此處存在一個峰，如前文所述，此峰(T2=0.22 ms)為類固體信號，使用式(2)描述的單個Beta峰來表示該類固體的信號峰。T2譜中其他峰信號全部來自于油的貢獻，用識別出的其他峰的線性疊加來表示。利用式(5)～式(10)給出相應的峰參估計和約束，然后通過式(11)進行曲線擬合，最終確定各峰函數的參數。6個頁巖樣品的T2譜擬合結果如圖5和表1所示，可以看出，擬合T2譜與原始T2譜非常吻合(R2均大于0.999 9，擬合標準差小于0.1)。經ND方法剝離類固體信號后，0.01～1 ms弛豫時間內的核磁共振信號強度明顯減小，而1 ms后的信號強度基本不變，構成了純孔隙流體T2譜(圖4)。

圖5 核磁共振-反褶積方法分離頁巖樣品T2譜

表1 頁巖樣品T2譜曲線擬合結果

2.2 孔隙分布表征

在流體飽和頁巖樣品中，T2值和孔隙半徑r之間的換算關系可表示為[3,16]

r=FρT2

(13)

式(13)中：F為形狀因子，對于球形孔取值3，對于圓柱形孔取值2；ρ為表面弛豫系數。

孔隙體積dV可通過計算孔隙中流體的體積得到，其計算公式為

(14)

式(14)中：dAtotal為給定弛豫速率1/T2時的信號幅值；μfluid為流體密度；AI為振幅指數，其定義為

(15)

式(15)中：Atotal為總振幅，m為流體質量。

與dV相比，dV/dlgD可以更好地揭示不同孔徑范圍的孔隙對總孔隙體積的貢獻，其中D為孔隙等效直徑[17]。經過換算，純孔隙流體T2分布可以轉換為圖6所示的孔隙分布。作為對比，由LTNA、NMR-OS和NMR-O方法得到的孔隙分布也繪制在圖6中。其中LTNA方法主要用于精確表征納米級孔隙(<100 nm)。對于更大的孔隙，氮氣吸附所占的孔隙體積比例較小，使用LTNA方法無法得到準確的孔隙分布[18]。

由圖6可知，對于孔徑小于10 nm的孔隙，由ND、LTNA和NMR-O方法測得的孔隙分布具有比較好的一致性，而NMR-OS方法得到的孔隙分布在小孔徑范圍內出現了一個峰，其dV/dlgD值明顯大于LTNA和NMR-O方法的測試結果，這正是由于NMR-OS方法沒有去除類固體信號導致的。對于孔徑大于10 nm的孔隙，由ND和NMR-OS方法得到的孔隙分布幾乎相同。由于受到殘余流體的影響，NMR-O方法將一部分流體信號從OS譜中減去，往往導致孔隙分布偏小。

圖6 不同方法測定的頁巖樣品的孔隙分布

為了進一步測試ND方法的性能，在T2<1 ms的范圍內比較了4種方法得到的孔隙體積，結果如表2所示。可以看到，ND方法得出的孔隙體積相對誤差絕對值(<13%)明顯小于NMR-OS方法的相對誤差絕對值(高達134.76%)。這是由于傳統的NMR-OS方法將油和類固體信號都計入孔隙流體信號，因此這種方法有著巨大誤差。另外NMR-O方法的精度受洗油烘干的效果影響大，這一點在頁巖樣品3和樣品4的結果上體現比較明顯。其T2譜中連續多峰的存在[圖4(c)和圖4(d)]表明小孔隙中的流體沒有完全洗出，而殘余流體信號被算作類固體信號，因此經轉換后得到的孔隙體積偏小(相對誤差分別為-13.04%和-13.57%)。為此，需要大量的萃取劑和較長的干燥時間盡量完全去除頁巖中的流體，但這又會破壞孔隙的原始結構和分布。相比之下，ND方法不需要煩瑣耗時的洗油和干燥處理，在保證相同精度的前提下，具有簡單、快速、無損等優點。

表2 不同方法得到的T2<1 ms的范圍內孔隙體積及相對誤差

3 實例應用

從準噶爾盆地吉木薩爾凹陷二疊系蘆草溝組中的云屑砂巖、砂屑云巖、巖屑砂巖三類巖性中各選取1塊新鮮巖心，不做洗油烘干處理，直接進行核磁共振實驗，測量其流體飽和T2譜。利用已驗證過的ND方法分離出流體信號，再結合T2值和孔隙半徑r之間的換算關系，對3塊新鮮巖心的孔隙分布進行量化表征，3塊巖心的流體飽和T2譜、分離出的流體信號及半徑分布如圖7所示。

由圖7可知，巖屑砂巖新鮮巖心的T2譜呈現為三峰，表明巖心包含三類孔隙，包含納米級孔隙、微納孔隙和微米級孔隙，孔隙半徑分布范圍廣，最高峰T2弛豫時間為11.8 ms，對應孔隙半徑為201.17 nm，最大孔隙半徑為4 262.23 nm，孔隙半徑相對較大，表明其孔隙度和滲透率較高，物性較好。砂屑云巖新鮮巖心的T2譜基本呈現為單峰，表明巖心主要包含一類孔隙，其均質性好，但其孔隙半徑分布范圍窄，不包含微米級孔，納米級孔隙占比高，最高峰T2弛豫時間為4.0 ms，對應孔隙半徑為38.04 nm，最大孔隙半徑僅為610.63 nm，其孔隙半徑相對較小，孔隙度和滲透率較低，物性差。云屑砂巖新鮮巖心的T2譜呈現為雙峰，表明巖心包含兩類孔隙，最高峰T2弛豫時間為9.4 ms，對應孔隙半徑為142.54 nm，最大孔隙半徑為1 853.46 nm，其孔物性介于砂屑云巖、巖屑砂巖之間。

圖7 新鮮巖心流體飽和T2譜、分離出的流體信號及半徑分布綜合圖

4 結論

(1)ND方法無需煩瑣耗時的洗油和干燥處理過程，可從飽和流體頁巖樣品T2譜中分離出類固體信號，進而實現簡單、快速、準確的分析頁巖孔隙分布。

(2)對于孔徑小于10 nm的孔隙，由ND、LTNA和NMR-O方法測得的孔隙半徑的起始點和孔隙分布的變化趨勢具有可比性；對于孔徑大于10 nm的孔隙，由ND方法和NMR-OS方法得到的孔隙分布幾乎完全相同，而由NMR-O方法測得的孔隙分布往往偏小。

(3)當T2<1 ms時，ND方法測得的孔隙分布相對誤差絕對值(<15%)明顯小于通過NMR-OS方法得到的孔隙分布結果(>135%)。