地球物理学进展  2015, Vol. 30 Issue (2): 628-635   PDF    
二维核磁共振变参量迭代快速反演方法
李鹏举, 谷宇峰     
东北石油大学地球科学学, 大庆 163318
摘要:核磁共振二维谱应用在储层流体识别与评价中具有快速直观识别流体类型、储层参数估算精度高的特点, 使之能高效率的完成测井资料解释工作, 然而拥有稳健、可靠及计算快速的反演方法是获取高质量核磁共振二维谱的关键所在.通过对二维核磁共振子空间类型反演方法的大量研究, 在奇异值截断公式上进行了分析和改进, 使得反演方程组的奇异值矩阵能保留下更多的奇异值, 为反演结果的保真性奠定基础, 并对其进行分组处理使反演方程组的条件数降低, 从而提高求解的稳定性;在迭代计算中提出变参量迭代算法, 使求解更合理、快速及高效.在数值模拟中, 变参量迭代快速反演方法能够准确地还原30×30、60×80扩散—弛豫二维构造谱, 计算用时都在1分钟内;在油水试验中, 变参量迭代快速反演方法反演出的扩散—弛豫二维谱能够正确识别所测流体类型, 含油饱和度估算结果的相对误差为0.6%, 绝对误差为0.79%.变参量迭代快速反演方法能够快速有效得处理二维核磁共振数据, 反演出的扩散—弛豫二维谱质量高, 二维谱应用于解释中所得到的结论可靠性高, 表明该方法具有一定应用价值.
关键词二维核磁共振     TSVD法     变参量迭代快速反演方法     扩散-弛豫二维谱     流体识别与评价    
Rapid variable parameter iteration method of two-dimensional NMR
LI Peng-ju, GU Yu-feng     
Earth Science Institute, Northeast Petroleum University, Daqing 163318, China
Abstract: NMR two-dimensional spectrum has characteristics of fast visual identification of fluid types, and reservoir parameters estimated high accurately when it applied on fluid identification and evaluation, thus makes the NMR two-dimensional spectrum could finished the explanation work of logging data very efficiently, while having a robust, reliable and fast inversion method is the key of obtaining high quality NMR two-dimension spectrum. Though the extensive researching on subspace inversion methods used in two-dimensional NMR, the improvements of singular value's truncated formula has been made, which retains more singular values of inversion equation and lay the foundation of the results' fidelity, meanwhile the values will be packet processed in order to reduce the conditional number of the inversion equation, then the stability of the computing can be enhanced; in iterated processing, the algorithm of variable parameter iteration is proposed in order to solve equations more reasonably, more quickly and efficiently. In the numerical simulation, rapid variable parameter iteration method can accurately restore 30×30 and 60×80 diffusion-relaxation constructional spectrum, the computing time is less than 1 minutes; in oil-water experiments, diffusion-dimensional two-dimensional spectrum got from rapid variable parameter iteration method can identify the types of measured fluid correctly, the relative error of computing results of oil saturation is 0.6%, the absolute error is 0.79%. Rapid variable parameter iteration method can able to process the two-dimensional NMR data quickly and efficiently, the quality of inversed diffusion-relaxation spectrum is high, and the spectrum applied on NMR interpretation can get reliable conclusions, which shows that the inversion method has applied value.
Key words: two-dimensional NMR     TSVD method     rapid variable parameter iteration method     diffusion-relaxation two-dimensional spectrum     fluid identification and evaluation    
0 引 言

核磁共振二维谱在储层流体识别与评价上具有快速识别流体性质,储层参数估算精度高的优点(Sun et al., 2003),而获得高质量的核磁二维谱的关键之处是拥有稳健、可靠及计算快速的反演方法(肖立志,2010早先的二维核磁反演应用技术主要以BRD和TSVD两大方法为主,但这些方法在处理数据时存在很多问题,导致反演出的结果时常不可靠(周小龙等,2013).

BRD方法一般是先将求解方程Tikhonov正则化(Kazufumi and Bangti, 2011),再加入罚函数方程的解进行限制以致最终解在约束条件范围内(Borgia et al., 1988).由于正则化因子能够控制罚函数作用的强度(Butler et al., 1981Venkataramanan et al., 2002),因此解的质量高低取决于所选取的正则化因子值是否合适.正则化因子一般用L曲线等方法来确定(Honerkamp and Weese, 1990).如果处理的数据信噪比很低,这时L曲线等方法确定出的正则化因子常常是错误的,导致BRD方法在反演时不收敛或得出没有意义的结果.然而二维核磁共振采集数据的信噪比本就很低,这让BRD方法的适用性受到限制(Borgin et al., 1981;Honerkamp and Weese, 1990Venkataramanan et al., 2002宋公仆等,2013周小龙等,2013).

TSVD方法是通过对方程组系数矩阵进行奇异值分解后对奇异值采取截断处理,一般做法是舍去小的奇异值来降低方程组条件数以保证求解具有稳定性,之后用迭代手段计算出合理解(顾兆斌和刘卫,2007谢然红和肖立志,2009谭茂金和邹友龙,2012周小龙等, 2013ab).然而奇异值截断后所保留的数量非常少,这让采集数据中的很多微小的有效信息都被舍弃掉(顾兆斌和刘卫,2007谢然红和肖立志,2009),导致反演出的结果精度不够,二维谱严重失真(顾兆斌和刘卫,2007谢然红和肖立志,2009谭茂金和邹友龙,2012),使得TSVD方法变得不够可靠.

因此,针对 BRD和TSVD 两个方法中的不足,必须找到一种稳健、可靠的反演方法来处理二维核磁共振数据,为二维核磁共振解释工作奠定坚实基础(肖立志,2010;周小龙等,2013).

1 二维核磁共振扩散-弛豫模型的建立

二维核磁共振一般根据测量采用的脉冲序列的作用规律来建立其物理数学模型(Freedman et al., 2002Song et al., 2002顾兆斌和刘卫,2007谢然红和肖立志,2009Chouzenoux et al., 2010胡法龙等,2012谭茂金和邹友龙,2012周小龙等, 2013ab).当采用扩散编辑(DE)脉冲序列时(Freedman et al., 2002顾兆斌和刘卫,2007谢然红和肖立志,2009胡法龙等,2012谭茂金和邹友龙,2012),其模型为

式中:exp(-t/T2)为弛豫作用核函数(Venkataramanan et al., 2002李鹏举,2010),t是脉冲序列中的回波峰与所施加的90°脉冲峰之间的作用时长,ms;T2为横向弛豫时间常数,ms;P函数为扩散作用核函数,D为扩散系数,m2/s;te,l为长回波间隔,ms;te,s为短回波间隔,ms,f(D,T2)为扩散-弛豫二维谱,S(te,l,te,s)为回波串幅度函数.

P函数的具体形式取决于te,l与te,s之间的大小关系(Freedman et al., 2002Song et al.,2002Chouzenoux et al., 2010),当te,l大于te,s时,为

当te,l等于te,s时,为

式(2)与式(3)中的ad与as为权重系数,与所施加的射频脉冲有关(Freedman et al., 2002);γ为旋磁比,rad·T-1·s-1;g为磁场梯度,T/m.

将式(1)展成离散形式(Freedman et al., 2002胡法龙等,2012谭茂金和邹友龙,2012),为

式中:bk(i)为脉冲序列第k个回波串的第i个回波幅度;k为回波串序数;i为回波序数;p为扩散系数布点个数;m为横向弛豫时间布点个数;TElk为第k个长回波间隔,ms;TEs为短回波间隔,ms;T2j为设定的第j个横向弛豫时间常数;Dl为设定的第l个扩散系数; f lj为扩散-弛豫二维谱.

为简化表达,将式(4)变为

式中: B ik为回波串幅值向量,K ik,lj为弛豫核函数与扩散作用核函数的直积矩阵,ε ik是噪声项,加入此项是考虑了回波串生成的实际情况.

将式(4)进行简化的目的是为了将离散的扩散-弛豫模型形式转化成类似于一维核磁弛豫模型的形式,有便于研究其反演方法(Honerkamp and Weese, 1990Song et al., 2002Venkataramanan et al., 2002谢然红和肖立志,2009Chouzenoux et al., 2010胡法龙等,2012谭茂金和邹友龙,2012周小龙等, 2013ab).

2 变参量迭代快速反演方法研究

2.1 传统TSVD反演方法的缺陷

子空间类型反演方法在二维核磁共振反演方法领域中占有重要的地位,而TSVD方法在子空间类型方法中应用最广.TSVD方法能够很好的解决一维核磁反演问题(Venkataramanan et al., 2002江玉龙,2009陈华等,2009李鹏举,2010李辉等,2013),但在处理二维核磁共振数据时,方法固有的缺点就暴漏了出来(Song et al., 2002周小龙等,2013).

一方面,奇异值矩阵经过截断后所保留的奇异值数量很少(Honerkamp and Weese, 1990Venkataramanan et al., 2002周小龙等,2013),例如,对于一个20×20的奇异值矩阵,所含的奇异值为400个,但经过TSVD方法截断后奇异值个数最多不超多20个.利用这样的奇异值矩阵反演得出的二维谱将会严重失真,原因在于很多有用的信息细节被舍弃掉,不能在图谱中展示出来,从而导致流体评价的结果不可靠.

另一方面,在迭代环节中,TSVD方法的迭代速度与解的个数有关,解的个数越多,迭代速度就越慢(周小龙等, 2013ab).例如,反演一个20×20的二维谱,需求解400个未知数,利用TSVD方法反演时其迭代次数要达到上万次或几十万次,甚至有时不收敛,说明了TSVD方法处理二维核磁共振反演问题的效率极低.

2.2 奇异值截断、分组处理算法

奇异值截断的做法是先对式(5)中的系数矩阵进行奇异值分解,即 K=UWV T,然后用截断公式对 W 做处理(Venkataramanan et al., 2002陈华等,2009李辉等,2013).TSVD方法中的截断公式基本有两种(江玉龙,2008;陈华等,2009李辉等,2013),一种是 SNR1<1/ωi,(i=1,…k),另一种是将上述公式中的SNR1变为(a·SNR+b)/ω1,其中a,b是经验系数.满足截断公式的奇异值由0值替代,但由于奇异值的衰减速度特别快(胡法龙等,2012),通常在第20个左右的值就小于1了,而SNR1或者(a·SNR+b)/ω1一般小于1,因此通过上述的截断公式处理后所留下的奇异值将会非常少.

针对TSVD方法中截断公式的不足,提出一种包含新的截断标准,并对保留的奇异值进行分组处理的算法.

首先,考虑到对奇异值截断的目的就是为了降低方程组的条件数,所以在截断标准上以方程组的条件数作为依据,即先设定一个条件数,然后利用最大的奇异值与后续的奇异值逐个做比直到比值刚小于设定的条件数为止,则这两个值及其之间的所有奇异值就确定为保留的奇异值.当然,为了保留较多的奇异值,条件数要设定的大些,一般在15000到20000之间.

之后,对保留下的奇异值进行分组处理.其目的是因为上述设定的条件数过大,使得方程组高度病态,但如果将保留下的较小的奇异值变大,并保持大的奇异值不变,则条件数就会降低,从而就能求出稳定解.基于这种设想,提出a·ωi+c/ωib处理公式,其中a,b,c为变参量.当处理较小的奇异值时,公式中的幂函数项将会得出较大的值,这样较小的奇异值就达到变大的目的.而当处理较大的奇异值时,幂函数项将会得到很小的值,这样处理后得到的结果与原值相比就不会有太大的变化.处理后,最大的奇异值基本不变,而较小的奇异值却变大了,这就达到了降低条件数的目的.由于较大和较小奇异值要做不同效果的处理,因此所用的变参量是不同的,这也是为什么要进行分组处理.一般是将大于1的奇异值作为大奇异值组,剩下的奇异值作为小奇异值组.有时可根据情况需求调动分组的界限值.

截断、分组处理算法虽然保证了在降低方程组条件数下保留了更多的奇异值,但为了使处理后的值与原值不会有太大变化,公式中的变参量要合理确定.一般,对于大奇异值组,a为1,b和c为0,以保证较大的奇异值只有微小的变化,但有时根据需要,b和c可以都设为1或者稍大的值;对于小奇异值组,a设定在0到10之间,而b为1,c设定在1到2之间,以保证较小的奇异值有一定幅度的增长,不会变化太大.

W 经奇异值截断,分组算法处理后,再结合式(5)得出的求解的公式为

式中,r为保留的奇异值的个数,h为大奇异值组中奇异值的个数,a1,a2,b1,b2,c1,c2为变参量.

2.3 变参量迭代算法

TSVD方法中的迭代计算步骤如下(Venkataramanan et al., 2002;江玉龙,2008;陈华等,2009李鹏举,2010Lin et al., 2011李辉等,2013周小龙等,2013ab):先设定一个初始向量f0,代入到式(5)中,此时不考虑噪声项,得到 B0,将B与B 0做差得到Δ B 0,再代入到式(5)中得到Δ f0,令f 0Δ f0相加得到f1,如果此时f1中的全部分量为非负数,则停止迭代,f1即为所求的解;否则,将f1中的负数分量变为0,再作为下一次迭代的初矢量f 0,如此循环迭过程直到求出符合约束条件的解为止.

TSVD的迭代算法有两个问题(Honerkamp and Weese, 1990冯立新,2011周小龙等, 2013ab),一是差量Δ f 0,它在迭代中是一个逐渐变小的值,而变化程度与解的个数有关,解的个数越多则Δ f0在每次迭代中的变化程度就越小,导致求出的f1以极慢的速度趋向真值边界,这让方法常常出现迭代次数巨多或者不收敛状况.二是每次迭代完求出的解f1,它受到非负条件的限制,因此将其负数分量归零处理,但关键之处在于f1的非负分量,这些值在迭代中具有一种累加的效应.对于f1>1的部分,每次迭代完后其值都会有小幅度的增长,而01<1的部分其值的增长幅度较大,但因为Δ f0的变化程度太小,致使f1中的非负分量的累加速度很慢,f 1的变化就不大,所以迭代计算效率低.

针对上述TSVD迭代中的两个问题,提出一种变参量迭代算法.一是对于Δ f 0,提出公式Δ f i_new=aΔ·Δ f i,aΔ为放缩系数,i为迭代次数.这样做的目的是将每次得到的Δ f0进行扩大化,让求出的f1尽快逼近真值边界.但当f 1超出真值边界时,此时则需要对Δ f0进行缩小,使得f1回退到真值边界.二是对于f1,提出公式f =a· f i_new+c/ f bi_new,其中a,b,c为变参量,i为迭代次数,f i_newΔ f i_new与每次迭代的初矢量相加后的结果.公式本身主要针对 f1中的非负分量.当代入f 1>1的值时,令a为1,则公式中的前部分可保持原值不变,而后部分会得到比原值小得多的值.这样通过这小的值不仅能够加快原值的累加速度,还符合了此部分值累加慢的特点.同样,代入0< f 1<1的值,令a为1,则公式中的前部分可保持原值不变,而后部分则会得到比原值大得多的值,这样通过这大的值也达到了即加快原值的累加速度,又满足此部分值累加快的功效.当然对于负数部分,此公式也是适用的,只需将所有变参量的值取0就可达到归零处理.

通过对Δ f0和f1做上述变化后,其迭代速度将会明显提升,但迭代结果的精度也要同时得到保障,因此还需对两者做进一步处理.在对f1做变参量迭代算法处理时,其实是对它的f1>1、01<1、f1<0三个部分而做的,每个部分所用的变参量值不同,所以f 1要进行分组处理.对于Δ f 0而言,分组处理也同样要做,而且要划分地更细,即对Δ f 0>1、0<Δ f 0<1、-1<Δ f 0<0、Δ f 0<-1四部分做处理.这样做的原因有两点,一是对于负数部分,它不受条件约束,因此不能像 f 1做归零处理,否则每次迭代都会缺少一部分有效信息,致使最终解错误.二是四个部分有各自变化的特点,所用的变参量是不同的.Δ f 0>1的部分代表变化大的差量,它的aΔ值一般为1,如果太大会导致解 f 1很快超出真值边界;0<Δ f 0<1代表变化小的差量,它的aΔ值一般根据解 f 1离真值边界的远近程度来确定,如果离边界远,则aΔ值大些,离得近则aΔ值小些;-1<Δ f 0<0和Δ f 0<-1分别代表变化小的负差量和变化大的负差量,其作用是将解中的部分值变小,这部分值一般属于无效值,在最初的几次迭代时很大,为了尽快消除此部分值对求解的影响,就需要对-1<Δ f 0<0和Δ f 0<-1进行扩大,一般作用于Δ f 0<-1的aΔ值在1到100之间,aΔ越大,对此部分值的消减作用越明显.当这部分值衰减到接近0值边界时,这时主要通过-1<Δ f 0<0部分的值来调节其衰减速度,确定出合适的aΔ值就可让此部分值有效的到达0值边界.

根据上述的奇异值截断、分组处理算法以及变参量迭代算法,整个变参量迭代快速反演方法的计算流程如图 1所示.

图 1 变参量迭代快速反演方法计算流程图 Fig. 1 Processed flow chart of rapid variable parameter iteration method

3 变参量迭代快速反演方法数值模拟验证及分析

进行数值模拟的目的是检验反演方法在理论上是否具有可行性,并在反演精度和计算速度两方面上对方法的反演能力进行考察.

3.1 数值模拟正演

构造30×30和60×80两个扩散-弛豫二维谱(Flaum et al., 2004),如图 2所示,图谱横坐标设为扩散系数,布点范围为1×10-10 m2/s~1×10-2 m2/s,以对数形式分布;纵坐 标设为横向弛豫时间常数,布点范围为 1×10-5 s~1×101 s,以对数形式分布.

图 2 扩散-弛豫二维构造谱
(a)30×30扩散-弛豫二维构造谱;(b)60×80扩散-弛豫二维构造谱.
Fig. 2 Diffusion-relaxation two-dimensional constructional spectrum

采用DE脉冲序列来形成模拟回波串,其中长回波间隔设置16组,分别为0.2 ms、0.5 ms、1 ms、1.5 ms、2 ms、2.5 ms、3 ms、3.5 ms、4 ms、4.5 ms、5 ms、5.5 ms、6 ms、6.5 ms、7 ms、7.5 ms,短回波间隔为0.4 ms,磁旋比γ=2.68×108 rad·T-1·s-1,磁场梯度g=0.1 T/m,式(2)或式(3)中的ad和as分别设为0.5和0.3,由30×30二维谱模拟得到的回波串每组采集1024个回波,60×80二维谱模拟得到的回波串每组采集512个回波.将所有设置参数代入式(4)中计算出回波串的理论值,再加入信噪比为30的高斯白噪声,由此得到模拟回波串,如图 3所示.

图 3 扩散-弛豫二维构造谱的模拟回波串 Fig. 3 Simulated echo strings produced from diffusion-relaxation two-dimensional constructional spectrum
3.2 数值模拟反演

在反演30×30扩散-弛豫二维谱时,为突出变参量快速迭代反演方法的反演能力,引入TSVD方法作对比,并对前面提到该方法的两种截断情况都进行反演以增强对比说服力.为了区分清楚,对截断公式为SNR1<1/ωi的TSVD方法命名为第一种TSVD方法,截断公式为(SNR+b)/ω1的TSVD方法为第二种TSVD方法.

第一种TSVD方法经过截断后奇异值保留9个,在迭代中迭代次数设为3000,计算耗时共用877 s,反演得到的扩散-弛豫二维谱如图 4a所示.图中只显示出一个巨大的峰,与原谱的峰相比无论是位置还是形状都相差甚远,而且峰的幅值只还原到100多,表明第一种TSVD方法经过3000次迭代后计算出的结果没有任何实际意义,方法不适用于二维核磁反演.

图 4 扩散-弛豫二维反演谱图组
(a)第一种TSVD方法反演的30×30扩散-弛豫二维谱;(b)第二种TSVD方法反演的30×30扩散-弛豫二维谱; (c)变参量迭代快速反演方法反演的30×30扩散-弛豫二维谱;(d)变参量迭代快速反演方法反演的60×80扩散-弛豫二维谱.
Fig. 4 Figure set of diffusion-relaxation two-dimensional inversed spectrum

在第二种TSVD方法的截断公式中,a1,b为-10,截断后保留17个奇异值,在迭代中迭代次数设为3000,计算耗时共用875 s,反演得出的扩散-弛豫二维谱如图 4b所示.图中中心位置显示出一个大的峰,与构造谱的峰位置基本一致,可峰的形状没有还原出来,且峰的幅度还原到接近500,表明第二种TSVD方法计算出的结果比第一种TSVD方法要准确些,但仍不可靠.

在变参量迭代快速反演方法中,设定截断的条件数为17000,截断后奇异值保留56个,分组界线值设为1,大奇异值组份中有27个奇异值,系数a,b,c均为1,小奇异值组份中有29个奇异值,系数a为15,b为2,c为1.在迭代中,迭代次数设为200次,Δf0f1的变参量值如表 1所示,计算耗时共用33 s,反演得出的扩散-弛豫二维谱如图 4c所示.从图中可以看出,反演的峰位置正确,峰的形状还原程度高,而且幅值达到600多,非常接近原谱幅值,表明变参量迭代快速反演方法算出的结果精度高.另外,反演耗时在一分钟之内,表明变参量迭代快速反演方法计算效率高.

表 1 变参量迭代快速反演方法在反演30×30 二维谱时Δ f0与f 1的变参量值 Table 1 Variable parameters of Δ f 0 and f 1 when using rapid variable parameter iteration method to inverse 30×30 two-dimensional spectrum.

通过三种方法反演出的图谱对比可知,无论在计算结果的精度上还是计算速度上,变参量迭代快速反演方法都远高于第一种和第二种TSVD方法.

在60×80扩散-弛豫二维谱反演中,只采用变参量迭代快速反演方法.设定截断条件数为20000,截断后奇异值保留53个,分组界线值设为1,大奇异值组份有32个奇异值,系数a,b,c均为1,小奇异值组份有21个奇异值,系数a为13、b为2、c为1.在迭代中,迭代次数设为200次,Δ f0和f 1的变参量值如表 2所示,计算耗时共用51 s,反演出的扩散-弛豫二维谱如图 4d所示.从图中可看出反演的峰位置正确,中间部分的峰还原度高,但右下方的峰还原度较差,表明变参量迭代快速反演方法算出的结果具有一定精度,但准确程度没有反演30×30二维谱的结果高.另外,计算耗时约为1分钟,反演速度快,但慢于反演30×30二维谱的速度.

表 2 变参量迭代快速反演方法在反演60×80 二维谱时Δ f0与f 1的变参量值 Table 2 Variable parameters of Δ f 0 and f 1 when using rapid variable parameter iteration method to inverse 60×80 two-dimensional spectrum

通过30×30和60×80扩散-弛豫二维谱的反演情况来看,变参量迭代快速反演方法在计算结果精度和速度上会 受到谱维数(求解的个数)的影响,反演图谱的维数越高,反演结果精度越差,计算速度越慢,但整体的反演效率仍然很高.

4 变参量迭代快速反演方法油水试验验证及分析

4.1 试验设计

本次试验采用PE-02型核磁共振仪,白油作为模拟油,含有5%NaCl的水作为模拟地层水(Sune and Dunn, 2002).将油和水按31比例配比,其混合溶液作为测量对象,算得实际含油饱和度为76.4%.

在仪器测量设置上,脉冲序列选DE模式,长回波间隔设为16组,分别为0.2 ms、1 ms、2 ms、3 ms、4 ms、5 ms、6 ms、7 ms、8 ms、9 ms、10 ms、11 ms、12 ms、13 ms、14 ms、15 ms,短回波间隔设为0.4 ms,磁场梯度g=0.069063 T/m,扩散函数中的ad和as分别为0.54338和0.38215,每组回波串采集回波1024个;在反演二维谱布置中,横坐标设为扩散系数,布点范围从1×10-12 m2/s~1×10-6 m2/s,布点个数100个,按对数均匀布点,纵坐标设为横向弛豫时间,布点范围从1×10-4s~1×101/s,布点个数50个,按对数均匀布点.

4.2 反演结果验证及分析

试验测得的回波串如图 5a所示,由变参量迭代法确定出所需的变参量值后进行反演,其中条件数设为18000,分组界限值设为1,迭代次数设为200次.计算耗时共用58 s,反演出的扩散-弛豫二维谱如图 5b所示.根据水的扩散系数为2.5×10-9 m2/s,在图中可以确定出中上部的峰为水峰,而另一个峰经确定为油峰,这表明所测流体含有油水两相,与实际情况相符,证明了此二维谱能正确识别出所测流体性质.

图 5 二维核磁共振油水试验测量回波串及二维反演谱
(a)二维核磁共振油水试验测量回波串;(b)油水试验扩散-弛豫二维反演谱.
Fig. 5 Echo strings measured from two-dimensional NMR oil-water experiment and the inversed two-dimensional spectrum

在计算流体饱和度时,要对油和水再各测定一次,目的是要根据反演出的二维谱算出单位体积油和水的信号贡献值.此试验中油相单位体积反演信号贡献值经计算为215.58,水相为218.05.图 5b中油峰信号总和为5379.912,算得油的体积值为24.96,水峰信号总和为1737.28,算得水的体积值为7.97,由此得知测量的混合溶液的含油饱和度为75.8%,其相对误差为0.6%,绝对误差为0.79%,误差值都非常小,证明了此二维谱在估算所测流体饱和度时结果精度高,具有可靠性.

通过二维核磁共振油水试验可知,变参量迭代快速反演方法在一分钟之内处理二维核磁数据所得到的扩散-弛豫二维谱能正确识别所测流体类型,估算地流体饱和度值精度高,表明该方法能够快速、有效的处理二维核磁数据,具有一定应用价值.

5 结 论

5.1     奇异值截断、分组处理算法能够在保留很多奇异值的情况下通过分组处理来降低反演方程组的条件数,为反演计算的稳定性提供了保障.

5.2      变参量迭代算法通过合理调节Δ f0和f 1的大小使迭代计算在短时间内就能得出可靠的结果,达到了快速反演的目的.

5.3      变参量迭代快速反演方法能够快速、有效地处理二维核磁数据,反演出的二维谱能准确识别所测流体类型,估算地流体饱和度结果精度高,表明方法具有一定应用价值,对二维核磁共振快速反演新方法的研究也具有指导意义.

致 谢 感谢审稿专家的宝贵修改意见和编辑部老师的帮助.
参考文献
[1] Borgia G C, Brown R J S, Fantazzini P. 1998. Uniform-penalty inversion of multiexponential decay data[J]. J. Magn. Reson., 132(1): 65-77.
[2] Butler J P, Reeds J A, Dawson S V. 1981. Estimating solutions of first kind integral equations with nonnegative constraints and optimal smoothing[J]. SIAM J. Numer. Anal., 18(3): 381-397.
[3] Chouzenoux E, Moussaoui S, Idier J, et al. 2010. Optimization of a maximum entropy criterion for 2D nuclear magnetic resonance reconstruction[C]. 2010 IEEE International Conference on Acoustics Speech and Signal Processing (ICASSP). Dallas, TX: IEEE, 4154-4157.
[4] Chen H, Pan K J, Tan Y J. 2009. A new method for multi-exponential Inversion of NMR relaxation signal[J]. Well Logging Technology (in Chinese), 33(1): 37-40.
[5] Feng L X. 2011. The computational methods and applications of inverse problems (in Chinese)[M]. Harbin: Harbin Industrial University Press.
[6] Flaum M, Chen J, Hirasaki G. 2004. NMR diffusion editing for D-T2 maps: Application to recognition of wettability change[C].//SPWLA 45th Annual Logging Symposium. Noordwijk, Netherlands: Society of Petrophysicists and Well-Log Analysts.
[7] Freedman R, Heaton N, Flaum M, et al. 2002. Wettability, saturation, and viscosity using the magnetic resonance fluid characterization method and new diffusion-editing pulse sequences[C].//SPE Annual Technical Conference and Exhibition. San Antonio, Texas: Society of Petroleum Engineers.
[8] Gu Z B, Liu W. 2007. The inversion of Two-dimensional spectrum of NMR map[J]. Chinese Journal of Magnetic Resonance (in Chinese), 24(3): 311-319.
[9] Honerkamp J, Weese J. 1990. Tikhonovs regularization method for ill-posed problems. A comparison of different methods for the determination of the regularization parameter[J]. Continuum Mech. Thermodyn., 2(1): 17-30.
[10] Hu F L, Zhou C C, Li C L, et al. 2012. Fluid identification method based on 2D diffusion-relaxation nuclear magnetic resonance (NMR)[J]. Petroleum Exploration and Development (in Chinese), 39(5): 552-558.
[11] Jiang Y L. 2009. Study on data interpretation and experiments of MRIL-Prime logging (in Chinese)[D]. Changchun: Jilin University.
[12] Kazufumi I, Bangti J. 2011. A new approach to nonlinear constrained Tikhonov regularization[J]. Inverse. Probl., 27(10): 105005-105027.
[13] Li H, Xiao Z X, Liu X J. 2013. Study on T2 spectrum inversion of nuclear magnetic resonance based on singular value decomposition[J]. Journal of East China Institute of Technology(Natural Science Edition) (in Chinese), 36(1): 65-68.
[14] Li P J. 2010. Study on T2 spectrum inversion and fluid identification and evaluation of NMR (in Chinese)[D]. Daqing: Northeast Petroleum University.
[15] Lin F, Wang Z W, Li J Y, et al. 2011. Study on algorithms of low SNR inversion of T2 spectrum in NMR[J]. Appl. Geophys., 8(3): 233-238.
[16] Song G P, Wu L, Cheng J J. 2013. Study of sequence and inversion algorithm in 2D NMR well logging[J]. Foreign Electronic Measurement Technology (in Chinese), 32(9): 38-41.
[17] Song Y Q, Venkataramanan L, Hürlimann M D, et al. 2002. T2-T1 correlation spectra obtained using a fast two-dimensional Laplace inversion[J]. Journal of Magnetic Resonance, 154(2): 261-268.
[18] Sun B, Dunn K J. 2002. Core analysis with two dimensional NMR[C]. 2002 International Symposium of the Society of Core Analysts, 22-25.
[19] Sun B Q, Dunn K J, Lopes J, et al. 2003. T1MAS 2D NMR technique provides valuable information on produced oil[C]. SPWLA 44th Annual Logging Symposium. Galveston, Texas: Society of Petrophysicists and Well-Log Analysts .
[20] Tan M J, Zou Y L. 2012. A hybrid inversion method of (T2, D) 2D NMR logging and observation parameters effects[J]. Chinese Journal of Geophysics (in Chinese), 55(2): 683-692, doi: 10.3969/j.issn.0001-5733.2012.02.032.
[21] Venkataramanan L, Song Y Q, Hurliman M D. 2002a. NMR measurements and methods of analyzing NMR data: US, 6462542 B1[P]. 10, 08.
[22] Venkataramanan L, Song Y Q, Hurlimann M D. 2002b. Solving fredholm integrals of the first kind with tensor product structure in 2 and 2. 5 dimensions[J]. IEEE Transactions on Signal Processing, 50(5): 1017-1026.
[23] Xiao L Z. 2007. Some important issues for NMR logging applications in China[J]. Well Logging Technology (in Chinese), 31(5): 401-407.
[24] Xie R H, Xiao L Z. 2009. The (T2, D) NMR logging method for fluids characterization[J]. Chinese Journal of Geophysics (in Chinese), 52(9): 2410-2418.
[25] Zhou X L, Nie S D, Wang Y J, et al. 2013a. A review on the inversion methods in 2D NMR[J]. Chinese Journal of Magnetic Resonance (in Chinese), 30(2): 293-305.
[26] Zhou X L, Nie S D, Wang Y J, et al. 2013b. An iterative truncated singular value decomposition (TSVD) -based inversion methods for 2D NMR[J]. Chinese Journal of Magnetic Resonance (in Chinese), 30(4): 541-551.
[27] 陈华, 潘克家, 谭永基. 2009. 核磁共振弛豫信号多指数反演新方法[J]. 测井技术, 33(1): 37-40.
[28] 冯立新. 2011. 反问题的计算方法及应用[M]. 哈尔滨: 哈尔滨工业大学出版社.
[29] 顾兆斌, 刘卫. 2007. 核磁共振二维谱反演[J]. 波谱学杂志, 24(3): 311-319.
[30] 胡法龙, 周灿灿, 李潮流,等. 2012. 基于弛豫-扩散的二维核磁共振流体识别方法[J]. 石油勘探与开发, 39(5): 552-558.
[31] 江玉龙. 2009. P型核磁测井解释处理方法研究[博士论文]. 长春: 吉林大学.
[32] 李辉, 肖忠祥, 刘晓娟. 2013. 基于奇异值分解的核磁共振T2谱反演研究[J]. 华东理工大学学报(自然科学版), 36(1): 65-68.
[33] 李鹏举. 2010. 核磁共振T2谱反演及流体识别评价方法研究[博士论文]. 大庆: 东北石油大学.
[34] 宋公仆, 吴磊, 程晶晶. 2013. 二维核磁共振测井时序与反演算法实验研究[J]. 国外电子测量技术, 32(9): 38-41.
[36] 谭茂金, 邹友龙. 2012. (T2, D)二维核磁共振测井混合反演方法与参数影响分析[J]. 地球物理学报, 55(2): 683-692, doi: 10.3969/j.issn.0001-5733.2012.02.032.
[37] 肖立志. 2007. 我国核磁共振测井应用中的若干重要问题[J]. 测井技术, 31(5): 401-407.
[38] 谢然红, 肖立志. 2009. (T2, D)二维核磁共振测井识别储层流体的方法[J]. 地球物理学报, 52(9): 2410-2418.
[39] 周小龙, 聂生东, 王远军,等. 2013a. 核磁共振二维谱反演技术综述[J]. 波谱学杂志, 30(2): 293-305.
[40] 周小龙, 聂生东, 王远军,等. 2013b. 基于迭代TSVD的NMR二维谱反演算法[J]. 波谱学杂志, 30(4): 541-551.