CN107561112A - 一种获取岩石渗透率剖面的核磁共振方法 - Google Patents
一种获取岩石渗透率剖面的核磁共振方法 Download PDFInfo
- Publication number
- CN107561112A CN107561112A CN201710869691.2A CN201710869691A CN107561112A CN 107561112 A CN107561112 A CN 107561112A CN 201710869691 A CN201710869691 A CN 201710869691A CN 107561112 A CN107561112 A CN 107561112A
- Authority
- CN
- China
- Prior art keywords
- mrow
- msubsup
- passages
- msup
- rock
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Landscapes
- Magnetic Resonance Imaging Apparatus (AREA)
Abstract
本发明涉及一种获取岩石渗透率剖面的核磁共振方法。该方法依托核磁共振加权成像方法获取岩石的原位孔隙度及孔隙尺寸信息。通过相邻层位处孔隙尺寸分布的偶合关系定义局部连通因子,进一步考虑了相邻位置的孔隙连通性,不仅提高了渗透率的评价效果,还能提供岩心在某一轴向方向上的连续渗透率剖面,得到的渗透信息更加详细。在石油地质勘探及提高采收率方面有重要的应用空间。
Description
技术领域
本发明涉及地球物理领域,具体为一种获取岩石渗透率剖面的核磁共振方法。
技术背景
岩石渗透率是一个重要的岩石物理参数,用于评价岩石内饱含流体在孔隙中的流动难易程度,因此通常可用于油藏储层建模,进而优化实际油气开发方案。就本质而言,岩石渗透率是一个内部饱含流体在孔隙结构中的流动函数,可等效考虑为有效孔隙流通管道的横截面积。目前在油田现场开发或实验室中,有两种主要的方法对岩石渗透率进行评价。第一种为直接流动测量,测量过程可以在岩石某横截面上通气或水获得该横截面上的有效渗透率。这种测量模式依托的为Darcy定律,该方法的缺点是需要依托复杂的测量设备获取该参数;另外一种为间接测量方式,这些测量方式通过一些物理测量手段如声波、电阻率和核磁共振等来获取岩石的孔隙度及孔隙尺寸,随后利用这些参数对岩石渗透率进行评价。
作为一种已被广泛应用在岩石物理和油气勘探过程中的方法,核磁共振成像和核磁共振弛豫能够提供岩石的宏观流体含量和微观孔隙结构等信息。常见的技术手段如弛豫时间T1/T2测量能够给出岩石的孔隙尺寸和孔隙结构等信息。通过大量的实际研究工作发现,核磁共振弛豫测量结果可以通过一些经验公式和孔隙网络模型与渗透率参数相关联。第一种评价渗透率的核磁共振模型公式为:
其中,κ为岩石渗透率;c为评价经验因子,与岩石岩性相关;φ为岩石孔隙度;T1为纵向弛豫时间。另外一种评价渗透率的核磁共振模型公式为:
其中,κ为岩石渗透率;c为评价经验因子,与岩石岩性相关;φ为岩石孔隙度;FFI为岩石孔隙内束缚水含量,BVI为岩石孔隙内自由流体含量。因此为了得到以上两种流体的相对含量,需要对岩石进行离心实验,分别进行离心前后的核磁共振测量,得到离心前后的孔隙含水率,以此来评价自身渗透率。
在岩石测量实际情况中,尤其是针对非均质性较强的岩石,单一渗透率值已经无法准确评价岩石内部流体的流动特性。举例说明如下:假定岩石左侧T1弛豫时间为0.5s,孔隙度为0.15,这样的区域代表了岩石中孔隙相对较大的溶洞;岩石右侧T1弛豫时间为0.05s,孔隙度为0.45,这样的区域代表了岩石中孔隙空间较小的微管道。如果按照第一种核磁共振模型公式对该岩石渗透率进行评价的话,这两个区间的渗透率几乎相同,流体可较易流动。但是对于整体岩石的真实情况来说,这两个区域由于孔隙尺寸的差异性,内部流体很难流动,整体的渗透率相对较低。因此需要引入新方法对岩石渗透率进行准确评价。
发明内容
针对现有技术中的上述缺陷,本发明提出一种获取岩石渗透率剖面的核磁共振方法,具体为。
步骤1、在TRS通道上向被测样品的含氢质子自旋系统施加90°脉冲将所述宏观磁化强度矢量M0扳转至横向平面;
步骤2、等待极短时间T后,在TRS通道上向被测样品自旋系统施加180°脉冲,重聚散相后的横向平面磁化矢量;
步骤3、再次等待极短时间T后,在TRS通道上向被测样品自旋系统施加90°脉冲,将重聚后的横向平面磁化矢量扳转90°至纵向轴向上(与静磁场方向一致);
步骤4、在GRD通道上向被测样品施加一个恒定幅度的梯度,该梯度用于消除目前自旋系统中仍然残留在横向平面上的磁化矢量;
步骤5、接下来在TRS通道上继续施加一小角度α射频脉冲,之后在GRD通道上施加相位编码梯度脉冲,该脉冲宽度为δ,高度从-gmax到+gmax变化m步。
步骤6、在GRD通道上的梯度施加完毕后,在TRS通道上重复施加180°射频脉冲,一定时间之后在ACQ通道采集到一个自旋回波信号。实际中,只记录该回波信号的峰值即可;
步骤7、在ACQ通道上的信号采集完毕后,在TRS通道上再次施加一个180°射频脉冲,将残留在纵向方向上的磁化矢量进行翻转。之后,在GRD通道上再次施加一个恒定幅度的梯度,用于消除目前自旋系统中仍然残留在横向平面上的磁化矢量;
步骤8、从TRS通道上的第一个小角度α射频脉冲开始至GRD通道上的最后一个恒定幅度的梯度为止,整个时序持续时间为Δ。采集系统持续循环这一部分时序N次,在ACQ通道中一次将采集得到N个自旋回波峰值。改变GRD梯度中的相位编码梯度幅度值m次,在ACQ通道中最终将采集得到m*N个自旋回波峰值。对采集到的所述回波串信号进行核磁共振数据处理,可得到我们所需的快速核磁共振加权成像结果。
本发明首先给出一种具体的快速核磁共振加权成像方法,该方法与常规的核磁共振弛豫时间测量相比,融合了核磁共振成像技术和弛豫技术,因此能够给出更为准确的岩心测量结果。随后,在渗透率计算过程中,根据相邻层位岩石的孔隙尺寸分布定义内部连通因子,将该连通因子带入到最终渗透率剖面的计算过程中,通过与直接测量方法比对发现,本发明阐述的方法优化并提高了对岩石渗透率的评判效果,对油气行业中的开发和多次采收方案有积极贡献。
附图说明
图1为本发明实施例一所提供的岩石核磁共振加权成像脉冲序列图;
图2为本发明实施例二所提供的岩石核磁共振加权成像脉冲序列图;
图3为本发明实施例所提供的岩石核磁共振加权成像图谱;
图4为本发明实施例所提供的局部连通因子定义方法详述;
图5为本发明实施例所提供的岩石渗透率剖面对比图。
其中,TRS代表核磁共振系统的脉冲发射通道,GRD代表梯度脉冲发射通道,ACQ代表核磁共振系统的信号接受通道。
具体实施例
结合说明书附图说明本发明的具体实施方式。在此需要说明的是,对于这些实施例方式的说明用于帮助理解本发明,但并不构成对本发明的限定。
首先对必要的核磁共振的基本原理概念和理论进行介绍。
静磁场B0。静磁场由磁体提供,决定核磁共振信号的信噪比。被测样品置于静磁场中,自旋系统内发生能级分裂,沿着静磁场方向会产生一个宏观磁化矢量M0。M0由静磁场强度B0,温度等参数决定。磁体材料通常有永磁体和超导体两种。永磁体基本用于低场核磁共振测量;超导体通常用于医学成像和实验室内高场仪器化学谱分析中,需要使用液氦和液氮保持磁体温度恒定。
射频磁场B1与脉冲。射频脉冲为电磁信号,通常由线圈产生。射频脉冲产生的磁场为射频磁场。射频磁场的方向与静磁场方向垂直,实现对在静磁场中形成的磁化矢量的扳转操作,扳转角度为:θ=γB1tp。其中,γ为质子的磁旋比,B1为射频磁场强度,tp为射频脉冲的持续时间。因此可通过控制射频脉冲的幅值或持续时间达到改变扳转角的目的。核磁共振脉冲序列由不同数量和频率属性的射频脉冲按照设定时序组成。通过调整脉冲间时间间隔,脉冲角度及脉冲的频率选择性,实现对自旋系统的弛豫、扩散等测量。
磁场梯度与成像。脉冲磁场梯度由梯度线圈产生,通常在施加过程中考虑脉冲梯度线圈与射频线圈的涡流效应,注意屏蔽效果。通过空间磁场强度与梯度值之间的关系,可对被测样品进行相应的空间相位、频率和选层编码,实现不同维度上的空间成像。对于空间位置的某一方向,以z为例,通过在该方向上施加一幅度为g的梯度脉冲后,不同空间位置的质子Larmor频率为:
ω(z)=γB0+γgz
因此,采集得到的回波信号与空间所求的成像质子密度之间的相互关系为:
M(k)=∫ε(z)ei2πkzdz
ε(z)=∫M(k)e-i2πkzdk
其中,k为定义的波函数,与梯度脉冲的参数有关。当系统采用频率编码模式进行成像实验时,k=γgmaxδ/2π;当采用相位编码时,k=γgmaxδ/mπ。由上述公式可知,M(k)与ε(z)为Fourier变换对。因此对采集得到的回波信号进行Fourier变换即可得到成像结果。通过配合使用频率梯度编码和相位梯度编码模式和方法,即可得到高维度下的核磁共振成像结果.
自旋回波。自旋回波为核磁共振测量最常见的一种信号。首先对被测样品施加90°脉冲,将磁化矢量M0扳转至垂直于静磁场方向的横向平面上。由于分子的扩散及静磁场的空间非均匀性等原因,磁化矢量M0发生散相。这一段时间如果打开信号采集通道对信号进行采集,得到自由衰减信号。经历一定时间τ之后,施加180°脉冲。散相后的磁化矢量会在同等时间τ之后实现重聚,形成一个回波信号。改回波信号称之为自旋回波信号。自旋回波在核磁共振应用中主要有以下三个方面:(1)通过施加一连串的180°射频脉冲,反复形成自旋回波,记录回波串信号,此脉冲序列为CPMG脉冲序列。该信号对于研究孔隙介质的横向弛豫特性极为重要,在一定条件下可以得到孔隙尺寸相关信息;(2)在梯度磁场下通过改变梯度幅值或梯度持续时间,记录自旋回波幅值的变化,可以得到流体分子的自扩散系数;(3)通过施加成对的频率编码或相位编码梯度,解析被测样品的空间自旋密度信息,实现核磁共振成像。
弛豫。自旋系统从共振状态恢复至热平衡状态的过程。该过程在不同方向上由纵向弛豫时间T1或横向弛豫时间T2表征。T1又称为自旋-晶格弛豫时间,反映自旋系统与外部环境的能量交换。T2又称为自旋-自旋弛豫时间,反映自旋系统内部能量损耗。自旋系统弛豫过程可由Bloch方程进行描述。纵向弛豫时间T1可采用饱和恢复脉冲序列进行测量。通过改变两个脉冲之间的时间间隔TW,记录信号幅值,反映纵向磁化矢量在不同编辑时间下的演化过程:
以上方法所需时间较长。每一步TW下都需要质子自旋系统等待较长时间并达到热平衡状态下才能进行下一步实验,因此采集过程极慢。
本申请的T1成像技术基于一种快速T1测量方法,采用的是小角度α射频脉冲串采集得到最终测量结果。通过施加这一串含有N个小角度α射频脉冲的脉冲串,其被测样品本身自旋系统磁化矢量的分量可表述为:
采用相关适应的射频脉冲相位循环并对采集信号幅度,可得到每一个小角度α射频脉冲下的磁化矢量为:
其中,N为小角度射频脉冲的个数,Δ为相邻两个小角度脉冲之间的时间间隔。相对于常规的T1测量方法,由于对纵向磁化矢量只进行一次操作,因此上述方法能够在较短时间内完成T1测量。
弛豫成像。在实际测量应用中发现,仅仅得到被测样品的质子密度信息,即成像,对于分析样品的微观信息是远远不够的。因此,如果能够将弛豫信息的获得与成像技术结合可完美的实现对被测样品宏观和微观的观测。
不同生物组织或样本间的纵向弛豫时间T1差异相对明显,因此通常被选为加权信息与成像技术结合。如果采用常规的T1测量方法与成像结合,实际测量时间相较较长,不利于快速动态观测样品本身信息。因此在本发明中将快速纵向弛豫时间T1测量方法与相关成像技术融合,在优化调整相关参数的步骤下,给出快速纵向弛豫时间T1成像技术的可行方案。
根据测量被测样品纵向弛豫时间T1和空间梯度编码的相互关系,设计如图1和图2所示的两种快速核磁共振一维T1成像脉冲序列,通过采集信号可得以下响应公式:
M(kz,NΔ)=∫∫F(z,T1)·K1·K2dzdT1
其中两个核函数K1,K2的具体形式为:
K1=exp(i2πkzz)
对于图1中所示的相位编码模式,kz=γgmaxδ/mπ。当采用相位编码时,通过改变梯度幅度步数m,T1编辑中的小角度脉冲个数N获取二维数据。采用后续数据反演步骤对所得数据进行处理,得到被测样品的核磁共振加权成像结果F(z,T1)。
实施例一
关于该发明中涉及的新型核磁共振加权成像技术,如图1所示,其特征在于,所述方法包括:
步骤1、在TRS通道上向被测样品的含氢质子自旋系统施加90°射频脉冲将所述宏观磁化强度矢量M0扳转至横向平面;
步骤2、等待极短时间T后,在TRS通道上向所述自旋系统施加180°射频脉冲,重聚散相后的横向平面磁化强度矢量;所述极短时间T通常为几十微秒。
步骤3、再次等待极短时间T后,在TRS通道上向所述自旋系统施加90°射频脉冲,将重聚后的横向平面磁化强度矢量扳转90°至纵向轴向上(与静磁场方向一致);
步骤4、在GRD通道上向被测样品施加一个恒定幅度的梯度脉冲,所述梯度脉冲用于消除所述自旋系统中仍然残留在横向平面上的磁化强度矢量;
步骤5、接下来在TRS通道上继续施加一小角度α射频脉冲,之后在GRD通道上施加相位编码梯度脉冲,该脉冲宽度为δ,高度从-gmax到+gmax变化m步;所述小角度α通常在15°以下。
步骤6、在GRD通道上的梯度施加完毕后,在TRS通道上重复施加180°射频脉冲,一定时间之后在ACQ通道采集到一个自旋回波信号。实际中,只记录该回波信号的峰值即可;
步骤7、在ACQ通道上的信号采集完毕后,在TRS通道上再次施加一个180°射频脉冲,将残留在纵向方向上的磁化强度矢量进行翻转。之后,在GRD通道上再次施加一个恒定幅度的梯度,用于消除所述自旋系统中仍然残留在横向平面上的磁化矢量;
步骤8、从TRS通道上的第一个小角度α射频脉冲开始至GRD通道上的最后一个恒定幅度的梯度为止,整个时序持续时间为Δ。采集系统持续循环所述时序N次,在ACQ通道中一次将采集得到N个自旋回波峰值。改变GRD梯度中的相位编码梯度幅度值m次,在ACQ通道中最终将采集得到m*N个自旋回波峰值。对采集到的所述回波串信号进行核磁共振数据处理,可得到我们所需的快速核磁共振加权成像结果。
关于该发明中涉及的岩心局部连通因子定义方法,如图4所示,其特征在于,可通过相邻层位的T1分布函数的耦合程度进行定义。该定义基于一个合理假设,即如果相邻层位的T1分布函数耦合程度较高,则该局部区域的孔隙连通性相对较好,岩石渗透率较高;反之亦然。因此该局部连通因子的定义方法如下:
其中,Pi,i+1为第i层岩石与第i+1层岩石之间的连通因子;代表第i层岩石中的第j个T1分布对应的幅度值;a b c d分别代表相邻层位T1分布函数的与逻辑运算和并逻辑运算的区间值。由上述连通因子的定义可知,该因子在相邻层位的T1分布完全相同的情况下为1,完全不同的情况下为0。
本发明结合图3、4、5对给出的快速核磁共振加权成像的数据处理及后续渗透率剖面评价步骤进行详细阐述:
步骤1:首先对采集得到的M(k,NΔ)数据进行Fourier变换,在成像维度上对数据进行解编。Fourier变换为线性变换,为非病态问题,因此此处不再赘述。
步骤2:对得到的解编数据进行一维Inverse Laplace反演,得到最终不同空间位置处的T1分布。由于Inverse Laplace反演为病态问题,因此本处对Inverse Laplace变换进行简述。此处将引入正则化项对此数据矩阵进行反演。为了得到稳定准确的解F,通常采用Tikhonov正则化方法,引入平滑项来求解该问题:
其中,s是正则化因子,与采集数据的信噪比相关,||·||项代表矩阵的Frobenius范数。引入的正则化项决定求解结果的稳定性与准确性。正则化因子选取过大,尽管求解得到的分布越稳定,但是解的准确性越差,即所谓的过平滑;正则化因子选取过小,解的求取越准确,但是解的稳定性降低,出现的伪信号越多,即欠平滑。因此,综合考虑解的真实性和解的稳定性,使用合理的正则化因子,是该方法的重点。
步骤3:通过非负约束步骤可以得到特定正则化因子s下的非负约束解F(z,T1)。得到的结果如图3所示。图3的二维图谱中,纵轴代表岩石位置,横轴代表T1分布数值。顶部投影代表的是岩石整体的T1分布,通过标准样品刻度,变为区间孔隙度分布曲线。投影图中虚线代表独立的T1测量实验结果。与此实验结果比对发现二者吻合情况极好。右侧投影为含流体剖面曲线,代表岩石不同位置处的含流体孔隙度。竖条虚线为独立测量的孔隙度值。
步骤4:由于加权成像F(z,T1)函数给出了不同层位岩石的局部孔隙度和孔隙尺寸分布信息,因此可以依据公式计算求出每个层位处的局部渗透率,进而得到岩石的渗透率剖面。其中为第i层位处T1分布的对数平均值,计算方法为:
步骤5:上述步骤直接采用不同层位岩石的局部孔隙度和孔隙尺寸分布信息对岩石的局部渗透率进行评价。在进一步考虑本发明提出的局部连通因子Pi,i+1后,该评价公式可校正为:
通过比对由步骤4和步骤5得到的岩石渗透率剖面,并与实际渗透率直接测量方法的结果相比对,验证本发明方法在准确获取岩石渗透率的效果。对比结果如图5所示。该图中,曲线a代表考虑连通因子计算后得到的岩石渗透率剖面,竖虚线c为曲线a的几何平均值;曲线b代表未考虑连通因子计算后得到的岩石渗透率剖面,竖虚线d为曲线b的几何平均值;竖虚线e代表直接测量流体驱动得到的岩石渗透率值。通过三者比对可以发现,考虑了连通因子能够更准确的获取岩石的渗透特性。
实施例二
在实施例一的基础上,数据采集的方法可以替换为另外一种实现形式,其他步骤不做变化,其脉冲序列图如图2所示,所述方法包括:
步骤1、在TRS通道上向被测样品的含氢质子自旋系统施加90°射频脉冲将所述宏观磁化强度矢量M0扳转至横向平面;
步骤2、等待极短时间T后,在TRS通道上向被测样品自旋系统施加180°射频脉冲,重聚散相后的横向平面磁化矢量;所述极短时间T通常为几十微秒;
步骤3、再次等待极短时间T后,在TRS通道上向被测样品自旋系统施加90°脉冲,将重聚后的横向平面磁化矢量扳转90°至纵向轴向上(与静磁场方向一致);
步骤4、在GRD通道上向被测样品施加一个恒定幅度的梯度,该梯度用于消除目前自旋系统中仍然残留在横向平面上的磁化矢量;
步骤5、接下来在TRS通道上继续施加一小角度α射频脉冲,之后在GRD通道上施加频率编码梯度脉冲,该脉冲宽度为δ,高度固定为gmax;所述小角度α通常在15°以下。
步骤6、在GRD通道上的梯度施加完毕后,在TRS通道上重复施加180°射频脉冲,一定时间之后在ACQ通道采集到一个完整离散的自旋回波信号。在开始采集此自旋回波信号时,在GRD通道上施加频率解码梯度脉冲。该梯度脉冲与之前的频率编码梯度脉冲幅度相同,持续时间成2倍关系;
步骤7、在ACQ通道上的信号采集完毕后,在TRS通道上再次施加一个180°射频脉冲,将残留在纵向方向上的磁化矢量进行翻转。之后,在GRD通道上再次施加一个恒定幅度的梯度,用于消除目前自旋系统中仍然残留在横向平面上的磁化矢量;
步骤8、从TRS通道上的第一个小角度α射频脉冲开始至GRD通道上的最后一个恒定幅度的梯度为止,整个时序持续时间为Δ。持续循环所述时序N次,在ACQ通道中依次将采集得到N个自旋回波。在ACQ通道中最终将采集得到N个自旋回波。对采集到的所述回波串信号进行核磁共振数据处理,可得到我们所需的快速核磁共振一维T1成像结果。
实施例三
在实施例一到二的基础上,关于该发明中涉及的岩心局部连通因子定义方法,如图4所示,其特征在于,可通过相邻层位的T1分布函数的耦合程度进行定义。该定义基于一个合理假设,即如果相邻层位的T1分布函数耦合程度较高,则该局部区域的孔隙连通性相对较好,岩石渗透率较高;反之亦然。因此该局部连通因子的定义方法如下:
其中,pi,i+1为第i层岩石与第i+1层岩石之间的连通因子;代表第i层岩石中的第j个T1分布对应的幅度值;abcd分别代表相邻层位T1分布函数的与逻辑运算和并逻辑运算的区间值。由上述连通因子的定义可知,该因子在相邻层位的T1分布完全相同的情况下为1,完全不同的情况下为0。
如上所述,可较好地实现本发明。在不脱离本发明的原理和精神的情况下对这些实施例进行变化、修改、替换、整合和变型仍落入本发明的保护范围内。
Claims (6)
1.一种获取岩石渗透率剖面的核磁共振方法,其特征在于,所述方法包括,数据采集、数据处理和数据解释。
2.根据权利要求1所述的一种获取岩石渗透率剖面的核磁共振方法,其特征在于,所述数据采集具体包括:
步骤1.1、在TRS通道上向被测样品的含氢质子自旋系统施加90°射频脉冲将宏观磁化强度矢量M0扳转至横向平面;
步骤1.2、等待极短时间τ后,在TRS通道上向所述自旋系统施加180°射频脉冲,重聚散相后的横向平面磁化强度矢量;
步骤1.3、再次等待极短时间τ后,在TRS通道上向所述自旋系统施加90°射频脉冲,将重聚后的横向平面磁化强度矢量扳转90°至纵向轴向上;
步骤1.4、在GRD通道上向被测样品施加一个恒定幅度的梯度脉冲,所述梯度脉冲用于消除所述自旋系统中仍然残留在横向平面上的磁化强度矢量;
步骤1.5、在TRS通道上继续施加一小角度α射频脉冲,之后在GRD通道上施加相位编码梯度脉冲;其中,所述相位编码梯度脉冲宽度为δ,幅度从-gmax到+gmax变化m步;
步骤1.6、在TRS通道上重复施加180°射频脉冲,一定时间之后在ACQ通道采集到一个自旋回波峰值;
步骤1.7、在TRS通道上再次施加一个180°射频脉冲,将残留在纵向方向上的磁化强度矢量进行翻转;之后,在GRD通道上再次施加一个恒定幅度的梯度脉冲,用于消除所述自旋系统中仍然残留在横向平面上的磁化强度矢量;
步骤1.8、从所述步骤1.5中在TRS通道上施加一小角度α射频脉冲开始至GRD通道上施加的最后一个恒定幅度的梯度脉冲为止,整个时序持续时间为Δ;持续循环所述时序N次,在ACQ通道中依次将采集得到N个自旋回波峰值;同时,改变GRD通道上施加的相位编码梯度脉冲的幅度值m次,在ACQ通道中最终将采集得到m*N个自旋回波峰值,得到采集信号的响应矩阵为M(k,NΔ),k为波函数。
3.根据权利要求2所述的一种获取岩石渗透率剖面的核磁共振方法,其特征在于,所述数据处理步骤具体包括:
步骤2.1、对M(k,NΔ)数据进行Fourier变换,在成像维度上对数据进行解编;
步骤2.2、对得到的解编数据进行一维Inverse Laplace反演,得到最终不同空间位置处的T1分布;
步骤2.3、通过非负约束步骤可以得到特定正则化因子s下的非负约束解F(z,T1),F(z,T1)即为加权成像函数。
4.根据权利要求3所述的一种获取岩石渗透率剖面的核磁共振方法,其特征在于,所述数据解释的步骤具体包括:
所述加权成像函数F(z,T1)给出了不同层位岩石的局部孔隙度和孔隙尺寸分布信息,通过评价公式可计算求出每个层位处的局部渗透率κ,进而得到岩石的渗透率剖面;其中,所述评价公式为:
<mrow>
<msubsup>
<mi>&kappa;</mi>
<mrow>
<mi>l</mi>
<mi>o</mi>
<mi>c</mi>
<mi>a</mi>
<mi>l</mi>
</mrow>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mo>=</mo>
<mi>c</mi>
<mo>&CenterDot;</mo>
<msup>
<msub>
<mi>&phi;</mi>
<mi>i</mi>
</msub>
<mn>2</mn>
</msup>
<mo>&CenterDot;</mo>
<msup>
<msub>
<mi>&phi;</mi>
<mrow>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
<mn>2</mn>
</msup>
<mo>&CenterDot;</mo>
<msubsup>
<mi>T</mi>
<mrow>
<mn>1</mn>
<mi>l</mi>
<mi>g</mi>
<mi>m</mi>
</mrow>
<mi>i</mi>
</msubsup>
<mo>&CenterDot;</mo>
<msubsup>
<mi>T</mi>
<mrow>
<mn>1</mn>
<mi>lg</mi>
<mi>m</mi>
</mrow>
<mrow>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mo>&CenterDot;</mo>
<msup>
<mi>P</mi>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msup>
</mrow>
其中,为第i层岩石与第i+1层岩石之间的局部渗透率;c为评价经验因子;φi为第i层岩石的岩石孔隙度;Pi,i+1为第i层岩石与第i+1层岩石之间的连通因子;为第i层位处T1分布的对数平均值,两者的计算公式分别为:
<mrow>
<msup>
<mi>P</mi>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msup>
<mo>=</mo>
<mfrac>
<mrow>
<msubsup>
<mi>&Sigma;</mi>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mi>a</mi>
</mrow>
<mi>b</mi>
</msubsup>
<mfrac>
<msup>
<mrow>
<mo>&lsqb;</mo>
<mi>M</mi>
<mi>i</mi>
<mi>n</mi>
<mrow>
<mo>(</mo>
<msubsup>
<mi>A</mi>
<mi>j</mi>
<mi>i</mi>
</msubsup>
<mo>,</mo>
<msubsup>
<mi>A</mi>
<mi>j</mi>
<mrow>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
<mn>2</mn>
</msup>
<mrow>
<msubsup>
<mi>A</mi>
<mi>j</mi>
<mi>i</mi>
</msubsup>
<mo>&CenterDot;</mo>
<msubsup>
<mi>A</mi>
<mi>j</mi>
<mrow>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msubsup>
</mrow>
</mfrac>
</mrow>
<mrow>
<mi>L</mi>
<mi>e</mi>
<mi>w</mi>
<mi>t</mi>
<mi>h</mi>
<mo>&lsqb;</mo>
<mi>U</mi>
<mi>n</mi>
<mi>i</mi>
<mi>o</mi>
<mi>n</mi>
<mrow>
<mo>(</mo>
<msup>
<mi>A</mi>
<mi>i</mi>
</msup>
<mo>,</mo>
<msup>
<mi>A</mi>
<mrow>
<mi>i</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msup>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
</mfrac>
</mrow>
<mrow>
<msubsup>
<mi>T</mi>
<mrow>
<mn>1</mn>
<mi>lg</mi>
<mi>m</mi>
</mrow>
<mi>i</mi>
</msubsup>
<mo>=</mo>
<msup>
<mn>10</mn>
<mrow>
<mo>&lsqb;</mo>
<mo>-</mo>
<mfrac>
<mrow>
<msubsup>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>n</mi>
</msubsup>
<msub>
<mi>A</mi>
<mi>j</mi>
</msub>
<mo>&CenterDot;</mo>
<mi>log</mi>
<msub>
<mrow>
<mo>(</mo>
<msub>
<mi>T</mi>
<mn>1</mn>
</msub>
<mo>)</mo>
</mrow>
<mi>j</mi>
</msub>
</mrow>
<mrow>
<msubsup>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>N</mi>
</msubsup>
<msub>
<mi>A</mi>
<mi>j</mi>
</msub>
</mrow>
</mfrac>
<mo>&rsqb;</mo>
</mrow>
</msup>
</mrow>
其中,代表第i层岩石中的第j个T1分布对应的幅度值;a、b分别代表相邻层位T1分布函数的与逻辑运算和并逻辑运算的区间值;n代表层位T1分布函数的布点个数。
5.根据权利要求2所述的一种获取岩石渗透率剖面的核磁共振方法,其特征在于,所述步骤1.3中的所述纵向轴向与静磁场方向一致。
6.根据权利要求3所述的一种获取岩石渗透率剖面的核磁共振方法,其特征在于,所述步骤2.2中的反演采用Tikhonov正则化方法,正则化项为:
<mrow>
<mi>arg</mi>
<munder>
<mi>min</mi>
<mrow>
<mi>F</mi>
<mo>&GreaterEqual;</mo>
<mn>0</mn>
</mrow>
</munder>
<mo>|</mo>
<mo>|</mo>
<mi>M</mi>
<mo>-</mo>
<mi>K</mi>
<mi>F</mi>
<mo>|</mo>
<msup>
<mo>|</mo>
<mn>2</mn>
</msup>
<mo>+</mo>
<mi>s</mi>
<mo>|</mo>
<mo>|</mo>
<mi>F</mi>
<mo>|</mo>
<msup>
<mo>|</mo>
<mn>2</mn>
</msup>
</mrow>
其中,s是正则化因子,与采集数据的信噪比相关,||·||项代表矩阵的Frobenius范数,M为采集信号的响应矩阵,F为加权成像函数,K为数据处理中定义的核函数。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710869691.2A CN107561112B (zh) | 2017-09-23 | 2017-09-23 | 一种获取岩石渗透率剖面的核磁共振方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710869691.2A CN107561112B (zh) | 2017-09-23 | 2017-09-23 | 一种获取岩石渗透率剖面的核磁共振方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN107561112A true CN107561112A (zh) | 2018-01-09 |
CN107561112B CN107561112B (zh) | 2019-04-30 |
Family
ID=60982685
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201710869691.2A Active CN107561112B (zh) | 2017-09-23 | 2017-09-23 | 一种获取岩石渗透率剖面的核磁共振方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN107561112B (zh) |
Cited By (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108982567A (zh) * | 2018-06-04 | 2018-12-11 | 中国科学院电工研究所无锡分所 | 一种核磁共振岩心自旋回波成像方法 |
CN109375283A (zh) * | 2018-11-07 | 2019-02-22 | 中国石油大学(华东) | 一种砂岩储层3d渗透率演化史的分析方法 |
CN109613041A (zh) * | 2018-12-11 | 2019-04-12 | 西安理工大学 | 一种基于核磁共振的冻土渗透系数测定系统及方法 |
CN109932301A (zh) * | 2019-04-10 | 2019-06-25 | 西南石油大学 | 一种计算致密储层自发渗吸两相流体相对渗透率的方法 |
CN110702581A (zh) * | 2019-10-23 | 2020-01-17 | 山东省科学院海洋仪器仪表研究所 | 一种强非均质多孔介质的多尺度渗透率计算方法 |
CN112129801A (zh) * | 2020-09-14 | 2020-12-25 | 北京大学 | 一种多孔介质的选层t2弛豫谱快速测试方法 |
CN112414917A (zh) * | 2020-11-03 | 2021-02-26 | 西安石油大学 | 一种页岩油储层有机孔隙和无机孔隙的划分与表征方法 |
CN112461726A (zh) * | 2020-10-06 | 2021-03-09 | 长江大学 | 一种致密气藏砂岩岩石渗吸水驱前缘速度评价方法和系统 |
CN112834542A (zh) * | 2020-02-27 | 2021-05-25 | 苏州纽迈分析仪器股份有限公司 | 一种同时测量岩心分层含水率和孔径分布的方法 |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5387865A (en) * | 1991-09-20 | 1995-02-07 | Exxon Research And Engineering Company | Permeability determination from NMR relaxation measurements for fluids in porous media |
EP1655617A2 (en) * | 2004-10-29 | 2006-05-10 | University of Brunswick | Methods and apparatus for measuring capillary pressure in a sample |
CA2634104A1 (en) * | 2007-06-08 | 2008-12-08 | University Of New Brunswick | Methods suitable for measuring capillary pressure and relative permeability curves of porous rocks |
CN103353462A (zh) * | 2013-06-17 | 2013-10-16 | 中国石油大学(华东) | 一种基于核磁共振成像的岩石非均质性定量评价方法 |
CN103954639A (zh) * | 2014-04-09 | 2014-07-30 | 上海大学 | 一种探测凝胶在微孔道中分布的方法 |
CN104237284A (zh) * | 2014-10-17 | 2014-12-24 | 西安石油大学 | 脆硬性泥页岩微裂缝损伤的核磁共振检测方法 |
CN105241913A (zh) * | 2015-10-10 | 2016-01-13 | 西安石油大学 | 岩石微裂缝损伤变量的核磁共振定量分析方法 |
-
2017
- 2017-09-23 CN CN201710869691.2A patent/CN107561112B/zh active Active
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5387865A (en) * | 1991-09-20 | 1995-02-07 | Exxon Research And Engineering Company | Permeability determination from NMR relaxation measurements for fluids in porous media |
EP1655617A2 (en) * | 2004-10-29 | 2006-05-10 | University of Brunswick | Methods and apparatus for measuring capillary pressure in a sample |
CA2634104A1 (en) * | 2007-06-08 | 2008-12-08 | University Of New Brunswick | Methods suitable for measuring capillary pressure and relative permeability curves of porous rocks |
CN103353462A (zh) * | 2013-06-17 | 2013-10-16 | 中国石油大学(华东) | 一种基于核磁共振成像的岩石非均质性定量评价方法 |
CN103954639A (zh) * | 2014-04-09 | 2014-07-30 | 上海大学 | 一种探测凝胶在微孔道中分布的方法 |
CN104237284A (zh) * | 2014-10-17 | 2014-12-24 | 西安石油大学 | 脆硬性泥页岩微裂缝损伤的核磁共振检测方法 |
CN105241913A (zh) * | 2015-10-10 | 2016-01-13 | 西安石油大学 | 岩石微裂缝损伤变量的核磁共振定量分析方法 |
Non-Patent Citations (1)
Title |
---|
ELIZABETH 等 (夏宏泉 等 译): "一种根据核磁共振成像检验岩石物理数据的模拟方法", 《国外测井技术》 * |
Cited By (13)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108982567A (zh) * | 2018-06-04 | 2018-12-11 | 中国科学院电工研究所无锡分所 | 一种核磁共振岩心自旋回波成像方法 |
CN109375283A (zh) * | 2018-11-07 | 2019-02-22 | 中国石油大学(华东) | 一种砂岩储层3d渗透率演化史的分析方法 |
CN109613041B (zh) * | 2018-12-11 | 2022-04-22 | 西安理工大学 | 一种基于核磁共振的冻土渗透系数测定系统及方法 |
CN109613041A (zh) * | 2018-12-11 | 2019-04-12 | 西安理工大学 | 一种基于核磁共振的冻土渗透系数测定系统及方法 |
CN109932301A (zh) * | 2019-04-10 | 2019-06-25 | 西南石油大学 | 一种计算致密储层自发渗吸两相流体相对渗透率的方法 |
CN110702581A (zh) * | 2019-10-23 | 2020-01-17 | 山东省科学院海洋仪器仪表研究所 | 一种强非均质多孔介质的多尺度渗透率计算方法 |
CN112834542A (zh) * | 2020-02-27 | 2021-05-25 | 苏州纽迈分析仪器股份有限公司 | 一种同时测量岩心分层含水率和孔径分布的方法 |
CN112834542B (zh) * | 2020-02-27 | 2024-03-22 | 苏州纽迈分析仪器股份有限公司 | 一种同时测量岩心分层含水率和孔径分布的方法 |
CN112129801A (zh) * | 2020-09-14 | 2020-12-25 | 北京大学 | 一种多孔介质的选层t2弛豫谱快速测试方法 |
CN112461726A (zh) * | 2020-10-06 | 2021-03-09 | 长江大学 | 一种致密气藏砂岩岩石渗吸水驱前缘速度评价方法和系统 |
CN112461726B (zh) * | 2020-10-06 | 2023-06-30 | 长江大学 | 一种致密气藏砂岩岩石渗吸水驱前缘速度评价方法和系统 |
CN112414917A (zh) * | 2020-11-03 | 2021-02-26 | 西安石油大学 | 一种页岩油储层有机孔隙和无机孔隙的划分与表征方法 |
CN112414917B (zh) * | 2020-11-03 | 2023-09-01 | 西安石油大学 | 一种页岩油储层有机孔隙和无机孔隙的划分与表征方法 |
Also Published As
Publication number | Publication date |
---|---|
CN107561112B (zh) | 2019-04-30 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN107561112B (zh) | 一种获取岩石渗透率剖面的核磁共振方法 | |
CN107102020B (zh) | 多维核磁共振测量方法 | |
CN108956417B (zh) | 分析岩石孔隙无效吸水量的同位素核磁方法 | |
Davies et al. | Pore‐size distributions from nuclear magnetic resonance spin‐lattice relaxation measurements of fluid‐saturated porous solids. I. Theory and simulation | |
Guo et al. | Advances in low-field nuclear magnetic resonance (NMR) technologies applied for characterization of pore space inside rocks: a critical review | |
CN107748126B (zh) | 一种获取岩石孔隙尺寸和孔隙表面弛豫率的核磁共振方法 | |
US9551769B2 (en) | Magnetic resonance examination of porous samples | |
Song | Magnetic resonance of porous media (MRPM): A perspective | |
US6937014B2 (en) | Method for obtaining multi-dimensional proton density distributions from a system of nuclear spins | |
CN104215652B (zh) | 确定油气饱和度的方法和装置 | |
Sun et al. | Developing a new NMR-based permeability model for fractured carbonate gas reservoirs | |
US9541513B2 (en) | Method for nuclear magnetic resonance diffusion measurements | |
WO2019232299A1 (en) | Method for assessment of pore-throat size distribution and permeability in porous media | |
US8532929B2 (en) | Method and apparatus to incorporate internal gradient and restricted diffusion in NMR inversion | |
Ramskill et al. | In situ chemically-selective monitoring of multiphase displacement processes in a carbonate rock using 3D magnetic resonance imaging | |
CN106199474A (zh) | 一种低场核磁共振二维谱反演算法 | |
Zou et al. | A novel method for NMR data compression | |
Yan et al. | Direct measurement of pore size and surface relaxivity with magnetic resonance at variable temperature | |
CN107728088B (zh) | 一种快速核磁共振t1成像方法 | |
Gu et al. | Two-step inversion method for NMR relaxometry data using norm smoothing and artificial fish swarm algorithm | |
CN107727678B (zh) | 一种核磁共振弛豫高低本征模态耦合方法 | |
CN109507221B (zh) | 一种多维核磁共振分子扩散耦合成像方法 | |
CN108120944A (zh) | 一种加权迭代的低场核磁共振t2谱反演算法 | |
CN107991710A (zh) | 一种储层孔径分布获取方法及装置 | |
Zhao et al. | Characterization of single-phase flow through carbonate rocks: quantitative comparison of NMR flow propagator measurements with a realistic pore network model |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |