張家成 張 宮* 黃若坤 覃瑩瑤

1(長江大學(xué)油氣資源與勘探技術(shù)教育部重點實驗室 湖北 武漢 430100)2(塔里木油田分公司勘探開發(fā)研究院 新疆 庫爾勒 841000)

0 引 言

核磁共振(NMR)測井通過觀測儲層中含氫核流體的NMR信號來獲得用于評價儲層參數(shù)的信息[1-2]。自20世紀(jì)50年代被引入到石油測井領(lǐng)域以來，核磁共振技術(shù)在飽和度、滲透率計算、流體識別、孔隙度、孔隙結(jié)構(gòu)評價等方面均得到了成功應(yīng)用[3-5]。2002年Sun等[6]和Hurlimann等[7]在石油測井領(lǐng)域引入了核磁共振波譜學(xué)中的二維核磁共振概念，進(jìn)一步拓展了核磁共振測井的使用范圍，提高了儲層測井解釋評價的準(zhǔn)確率[8]。核磁反演方法也從對單個回波串?dāng)?shù)據(jù)進(jìn)行反演發(fā)展到同時對多個回波串?dāng)?shù)據(jù)進(jìn)行反演[9]，但大量的數(shù)據(jù)加重了反演過程中的運(yùn)算負(fù)擔(dān)。一直以來，國內(nèi)外學(xué)者為滿足反演精度和反演速度上的要求，對二維核磁反演方法進(jìn)行了不斷的探索。顧兆斌等[10]利用改進(jìn)的奇異值分解法對二維核磁共振進(jìn)行反演；謝然紅等[9]針對二維的多回波串?dāng)?shù)據(jù)處理問題提出了多回波串聯(lián)合反演方法；李鵬舉等[11]針對BRD和TSVD兩種反演方法的不足提出了變參量迭代的快速核磁反演方法。在二維核磁反演方法研究過程中，繁瑣的數(shù)據(jù)構(gòu)造和結(jié)果評測過程增加了研究的時間成本。

本文從精度和效率出發(fā)，研發(fā)了一套二維核磁共振正反演模擬軟件，實現(xiàn)了核磁共振二維弛豫及反演的全過程，降低了二維核磁反演方法的研究難度。

1 原 理

1.1 核磁共振弛豫測量

核磁共振測井通過觀測流體中氫原子核的NMR信號，處理得到地層孔隙中流體的含量及其性質(zhì)。巖石物理信息，如滲透率、孔徑分布、孔隙度和束縛水，也可以通過核磁共振弛豫測量得到。核磁共振弛豫過程有兩種，其中縱向弛豫計算如下：

(1)

式中：Pi是第i個弛豫信號占總的縱向弛豫信號的百分比；t表示時間；T1表示縱向弛豫時間，反映了自旋系統(tǒng)從外界吸收或者釋放能量的效率，主要通過反轉(zhuǎn)恢復(fù)脈沖序列或飽和恢復(fù)法脈沖序列進(jìn)行測量[13]。

橫向弛豫計算如下：

(2)

式中：M(t)表示核磁共振測量得到的橫向弛豫信號；Pi是第i個弛豫信號占總的橫向弛豫信號的百分比；t表示時間；T2表示橫向弛豫時間，可以通過CPMG(Carr-Purcell-Meiboom-Gill)脈沖序列進(jìn)行測量[13]。在一系列CPMG序列之間，等待一段時間確保所有的自旋完全恢復(fù)，再準(zhǔn)備下一個激勵，這段時間稱為 “等待時間”。實際測量中，不同地層除了橫向弛豫時間T2不同外，縱向弛豫時間T1、擴(kuò)散系數(shù)D、內(nèi)部梯度G等均存在差異，會直接影響回波信號的衰減，用一組回波串進(jìn)行反演往往不能充分地挖掘流體之間的差異。二維核磁共振采用多組回波串組合測量，通過聯(lián)合反演，得到二維分布譜，使不同流體在處理成果中有更好的區(qū)分度。

1.2 二維核磁弛豫構(gòu)造

構(gòu)造二維弛豫譜是采用對數(shù)坐標(biāo)下的二維高斯分布公式進(jìn)行計算[14],具體實現(xiàn)如下：

(3)

式中：j=1,2,…,m，m是構(gòu)造譜的布點數(shù)目；Aj是第j組分的幅度；F是孔隙度刻度因子；T2j是第j組分的T2弛豫時間；T2mid是T2弛豫時間峰值的位置；Dj是第j組分的擴(kuò)散系數(shù)；Dmid是擴(kuò)散系數(shù)峰值的位置；σ1、σ2可以控制分布的展布寬度；ρ是T2j和Dj的相關(guān)系數(shù)。

為了更加真實地模擬出實際二維譜，二維構(gòu)造譜必須能夠反映多種不同性質(zhì)流體的疊加效果，所以式(3)需要進(jìn)一步改進(jìn)為：

(4)

式中：j=1,2,…,m，m是構(gòu)造譜的布點數(shù)目；設(shè)具有4類不同性質(zhì)的流體，F(xiàn)g是第g種流體的孔隙度分量；HIg是第g種流體的含氫指數(shù)；T2g,mid是第g種流體的T2弛豫中心值；Dg,mid是第g種流體的擴(kuò)散系數(shù)峰值；σ1g、σ2g分別是第g種流體T2峰值的中心展布寬度和擴(kuò)散系數(shù)峰值的展布寬度；ρ是T2j和Dj的相關(guān)系數(shù)。

1.3 常見的二維核磁共振脈沖序列

針對二維核磁共振測井?dāng)?shù)據(jù)采集的特點，國內(nèi)外學(xué)者均提出了不同的CMPG脈沖序列。Hurlimann等[7]提出了擴(kuò)散編輯脈沖序列，該序列中時序采集分為兩個窗口，第一個窗口內(nèi)僅采集兩個回波，通過調(diào)整這兩個回波的回波間隔來采集孔隙中流體的擴(kuò)散弛豫信息，而第二個窗口則與改良式CPMG脈沖序列一樣，都是采用最小的回波間隔來采集孔隙流體的橫向弛豫信息，如圖1所示。

圖1 擴(kuò)散編輯T2-D測量序列

在回波采集過程中，回波的幅度衰減計算如下：

(5)

式中：bik表示回波間隔為TElk時第一個窗口內(nèi)的第i個回波的幅度；γ表示氫核的旋磁比；G表示磁場梯度；TE表示第二個窗口內(nèi)的回波間隔；當(dāng)擴(kuò)散系數(shù)為Dp，橫向弛豫時間為T2j時的孔隙度分量用f(Dp,T2j)表示。該脈沖序列通過調(diào)整第一個窗口內(nèi)的回波間隔，分別采集多道CPMG回波信號，從而實現(xiàn)T2-D模式下的二維數(shù)據(jù)采集，利用式(5)將采集獲得的回波數(shù)據(jù)進(jìn)行聯(lián)合反演，便可得到地層孔隙中流體的T2-D分布。

謝然紅等[14]則利用回波間隔不同的多條CPMG回波串來完成T2-D二維數(shù)據(jù)的采集，相對于前一種方法更加簡便，利用現(xiàn)有的一維核磁共振測井儀器便可以完成該類數(shù)據(jù)的采集：

(6)

式中：bik是回波間隔為TEk時第i個回波的幅度[9]。

除此之外，還有通過調(diào)整第一個窗口內(nèi)的延時來采集多條CPMG回波信號的T1編輯序列[4]以及測量多組不同等待時間的多TW序列等其他常見的脈沖序列[15]，這里不再贅述。

1.4 二維譜反演方法

目前常見的二維譜反演方法主要有奇異值分解法(SVD)、非負(fù)約束正則化法(BRD)和截斷奇異值分解法(TSVD)等。

SVD法是通過反演測量誤差的核函數(shù)矩陣中的最小奇異值，以實現(xiàn)正則化，求取合理解。顧兆斌等[10]針對SVD法容易導(dǎo)致譜圖不連續(xù)的問題，提出一種適用于對信噪比較高的數(shù)據(jù)進(jìn)行反演的奇異值分解法。譚茂金等[15]提出了同時考慮TSVD的分辨率和精度的混合反演算法。

BRD法是先將求解方程Tikhonov正則化，再加入罰函數(shù)方程的解進(jìn)行限制，使最終的解在約束條件范圍內(nèi)。但是BRD法反演結(jié)果的優(yōu)劣取決于選用的正則化因子值是否合適[11]。一般利用L曲線等方法來確定正則化因子，但這些方法在處理低信噪比的數(shù)據(jù)時常常出現(xiàn)結(jié)果沒有意義或在反演時不收斂的問題。所以國內(nèi)外學(xué)者提出了很多解決辦法，如：文獻(xiàn)[16]利用L曲線的斜率來選擇正則化參數(shù)；文獻(xiàn)[17]提出一種相對誤差較小的基于SVD和BRD的正則化反演算法等。

2 軟件設(shè)計與實現(xiàn)

本文模擬軟件基于.NET4. 6框架，使用C#語言在Visual Studio 2017環(huán)境下進(jìn)行開發(fā)。

2.1 結(jié)構(gòu)設(shè)計

軟件設(shè)計結(jié)構(gòu)如圖2所示，整體分為四大部分：算法、應(yīng)用模塊、數(shù)據(jù)及圖形顯示。

圖2 軟件設(shè)計結(jié)構(gòu)

軟件設(shè)計時采用模塊化編程方式，將所有數(shù)學(xué)算法函數(shù)集中打包放在統(tǒng)一的算法庫中，從而獨立于交互界面，便于算法的持續(xù)改進(jìn)和擴(kuò)展；應(yīng)用模塊包括二維譜構(gòu)造模塊(用于根據(jù)用戶設(shè)置的流體類型及參數(shù)，構(gòu)造出滿足二維高斯分布的構(gòu)造譜)，模擬回波采集模塊(用于根據(jù)用戶設(shè)置的采集模式信息，模擬儀器進(jìn)行回波采集的過程，同時保存采集的回波串?dāng)?shù)據(jù))，二維譜反演模塊(實現(xiàn)了SVD和BRD二維譜反演算法，用于將采集得到的回波串?dāng)?shù)據(jù)反演成二維譜)；數(shù)據(jù)管理模塊和繪圖模塊同樣相對獨立，并設(shè)置有外部平臺接口，用于導(dǎo)入實驗室?guī)r心分析數(shù)據(jù)或?qū)嶋H測井?dāng)?shù)據(jù)。

主界面分為菜單欄、繪圖區(qū)、參數(shù)區(qū)三部分，如圖3所示。菜單欄有文件(NMR File)、參數(shù)設(shè)置(Parameters Set)和關(guān)于(About)。其中：NMR File菜單主要功能包括新建模擬文件、保存模擬文件、打開模擬文件和退出；Parameters Set菜單主要功能包括加載各個應(yīng)用模塊的參數(shù)設(shè)置界面、設(shè)置參數(shù)和導(dǎo)入其他平臺的數(shù)據(jù)；About菜單用于顯示本文軟件的開發(fā)過程及人員信息。

圖3 軟件主界面

2.2 功能實現(xiàn)

(1) 二維譜構(gòu)造模塊。根據(jù)回波串組合序列完成“二維譜構(gòu)造”模塊，該模塊是基于二維高斯分布進(jìn)行設(shè)計的。通過選擇相應(yīng)流體類型、填寫流體T2值及其展布、分量孔隙度(各類流體含量)、含氫指數(shù)和另一維度參數(shù)(縱向弛豫時間T1、擴(kuò)散系數(shù)D、或內(nèi)部梯度G)，以及最后的布點方案：T2及另一維度參數(shù)的布點的起始值、終止值和布點數(shù)目，來進(jìn)行二維譜的構(gòu)造。

對話框中的各個參數(shù)用戶可根據(jù)實際模擬的流體進(jìn)行設(shè)置，如圖4所示，完成后軟件會根據(jù)用戶填寫的結(jié)果自動計算出總孔隙度和核磁孔隙度。參數(shù)設(shè)置分為流體參數(shù)(各類流體的核磁屬性信息及孔隙度分量)和全局參數(shù)(模擬模式構(gòu)造譜分布參數(shù))，所有參數(shù)設(shè)置完成后，可以實時地顯示出二維構(gòu)造譜，如圖5所示。

圖4 二維譜構(gòu)造參數(shù)設(shè)置對話框

圖5 構(gòu)造的T2-T1二維分布譜

(2) 模擬回波串采集模塊。該模塊根據(jù)上一步設(shè)置構(gòu)造的二維譜數(shù)據(jù)，結(jié)合給定的模擬采集指令模擬采集回波串，并存儲回波數(shù)據(jù)。參數(shù)包括：需要采集的回波組數(shù)，各組回波串對應(yīng)的回波間隔(TE)值、等待時間(TW)值、回波采集個數(shù)(NE)、磁場梯度，以及設(shè)置的信噪比(S/N)。

參數(shù)設(shè)置完成后，軟件會自動顯示出模擬采集得到的回波數(shù)據(jù)及對應(yīng)的噪音數(shù)據(jù)，用戶可以根據(jù)圖形實時對參數(shù)進(jìn)行修訂，如圖6所示，模擬多TE二維采集得到的回波串?dāng)?shù)據(jù)如圖7所示。另外，該模塊允許用戶導(dǎo)入實驗室?guī)r心分析得到的二維回波數(shù)據(jù)，以及實際測井得到的二維回波數(shù)據(jù)，并對數(shù)據(jù)進(jìn)行后續(xù)的反演處理和分析。

圖6 回波采集參數(shù)設(shè)置對話框

圖7 回波采集結(jié)果

(3) 二維譜反演模塊。該模塊的主要功能是將上一步模擬采集得到的回波串?dāng)?shù)據(jù)進(jìn)行二維反演，得到最終的二維反演譜。需要設(shè)置的參數(shù)包括：T2維度布點參數(shù)(布點的數(shù)量，起始值和終止值)、第二維參數(shù)譜布點的方式(布點的數(shù)量，起始值和終止值)、反演使用的回波串選擇、平滑因子以及反演方法的選擇，參數(shù)設(shè)置對話框如圖8所示。圖9為反演得到的T2-T1二維譜。

圖8 二維譜反演參數(shù)設(shè)置對話框

除了軟件內(nèi)置的二維反演方法(SVD及BRD)外，還可以利用軟件中提供的C#接口(該接口命名為：INMR_InversionMethod)來編寫自己的反演算法，從而研究不同二維譜反演算法的適用性及影響因素。反演完成后，軟件會再次將反演得到的二維譜進(jìn)行正演(擬合)，觀測擬合結(jié)果與采集得到的回波串是否吻合，從而判斷反演效果好壞。二維回波串?dāng)M合效果如圖10所示。

圖10 二維譜回波串?dāng)M合效果

3 可靠性驗證及應(yīng)用實例

3.1 可靠性驗證

軟件編寫完成后，以硫酸銅溶液的核磁物理實驗為例，對正反演模擬軟件的可靠性進(jìn)行驗證，具體實施步驟如下：

第一步配置三種不同濃度的硫酸銅溶液。

第二步分別取三種不同濃度的硫酸銅溶液3 ml放入玻璃質(zhì)的集氣瓶中，然后將三只集氣瓶封口后放入試管內(nèi)進(jìn)行二維核磁共振(T2-T1)測量。實驗結(jié)果顯示：三種濃度的硫酸銅溶液的弛豫時間分別約為2 500、80、2 ms，如圖11所示。

圖11 二維反演譜(物理實驗)

第三步利用本文研發(fā)的正反演模擬軟件采用與物理實驗相同的參數(shù)進(jìn)行二維譜構(gòu)造、回波采集以及二維譜反演，得到構(gòu)造譜和反演譜的數(shù)值模擬結(jié)果分別如圖12、圖13所示。

圖12 二維構(gòu)造譜(軟件模擬)

圖13 二維反演譜(軟件模擬)

將軟件模擬硫酸銅實驗得到的二維構(gòu)造譜和二維反演譜與實際測量得到的二維反演譜進(jìn)行對比。可以看到數(shù)值模擬結(jié)果與物理實驗結(jié)果無論是分布大小還是分布位置都十分接近，驗證了軟件的可靠性。

3.2 應(yīng)用實例

研發(fā)的二維核磁共振正反演模擬軟件可以進(jìn)行多種二維模式的正反演數(shù)值模擬。現(xiàn)以T2-T1二維模式下核磁共振方法在含有束縛水的油層中識別流體的效果為例進(jìn)行模擬，驗證軟件的可靠性。

油水的特征參數(shù)如表1所示，根據(jù)該表參數(shù)構(gòu)造出油水模型的理想分布譜即構(gòu)造譜，如圖14所示，然后設(shè)置不同的等待時間模擬出15組回波串?dāng)?shù)據(jù)：TE均設(shè)為1 ms，TW分別為10、30、50、70、100、300、500、700、1 000、2 000、4 000、8 000、10 000、12 000、15 000 ms。

表1 油水模型參數(shù)表[14]

圖14 二維構(gòu)造譜

同時利用軟件的反演模塊將模擬數(shù)據(jù)進(jìn)行反演處理，結(jié)果如圖15所示。通過對比二維反演譜與構(gòu)造譜，可見油水的分布位置基本重合,表明了軟件對于T2-T1模式下的油水模型的正反演數(shù)值模擬具有較高的可信度。

圖15 二維反演譜

4 結(jié) 語

本文研發(fā)的二維核磁共振測井正反演模擬軟件可以進(jìn)行多種二維模式的正反演數(shù)值模擬，其模擬結(jié)果與實際巖心核磁共振儀器實驗結(jié)果非常接近。該軟件可以同時模擬4種不同的核磁共振弛豫過程，并且能夠?qū)Ω黝惡舜殴舱駥傩赃M(jìn)行單參量因素分析，為二維核磁共振相關(guān)研究工作提供了高效的技術(shù)手段。