CN114740532B - 一种水泥声阻抗、泥浆声阻抗的计算方法 - Google Patents

一种水泥声阻抗、泥浆声阻抗的计算方法 Download PDF

Info

Publication number
CN114740532B
CN114740532B CN202210272669.0A CN202210272669A CN114740532B CN 114740532 B CN114740532 B CN 114740532B CN 202210272669 A CN202210272669 A CN 202210272669A CN 114740532 B CN114740532 B CN 114740532B
Authority
CN
China
Prior art keywords
acoustic impedance
waveform
slurry
cement
time
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.)
Active
Application number
CN202210272669.0A
Other languages
English (en)
Other versions
CN114740532A (zh
Inventor
王�华
杨耿骁
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN202210272669.0A priority Critical patent/CN114740532B/zh
Publication of CN114740532A publication Critical patent/CN114740532A/zh
Application granted granted Critical
Publication of CN114740532B publication Critical patent/CN114740532B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/40—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging
    • G01V1/44—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging using generators and receivers in the same well
    • G01V1/48—Processing data
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Abstract

本发明公开了一种水泥声阻抗、泥浆声阻抗的计算方法，方法通过预先建立好的正演模型研究脉冲回波的特征参数与水泥声阻抗、泥浆声阻抗的变化规律，得到迭代方案。通过实测信号归一化幅度包络的对数在共振波区间内的一次拟合直线的斜率和截距，信号的频谱的最大分量对应的频率，套管厚度。将得到的信号的频谱的最大分量对应的频率，套管厚度以及水泥声阻抗、泥浆声阻抗初值带入正演模型中计算模拟波形并计算模拟波形的斜率和截距。先比较模拟波形与实际波形的截距，迭代得到泥浆声阻抗，在将得到的泥浆声阻抗带入迭代计算水泥声阻抗。本发明的方案能够在泥浆声阻抗未知的情况下计算水泥声阻抗和泥浆声阻抗。

Description

一种水泥声阻抗、泥浆声阻抗的计算方法
技术领域
本发明涉及固井质量评价领域，尤其涉及一种水泥声阻抗、泥浆声阻抗计算方法。
背景技术
在地下油气资源的开采过程中需要先进行钻井工作，之后需要在钻好的裸眼井中下放套管，并在套管和井壁之间填充水泥以固定套管和封隔油、气、水层。随着时间的推移，填充的水泥可能或出现裂缝从而导致套管受到液体腐蚀以及不同储层之间的窜通。不管是哪一种情况都会给油气资源的生产带来巨大的麻烦，严重的会导致重大的生产事故造成严重的人员伤亡、巨大的经济损失以及难以挽回的环境破坏。因此固井质量评价是油气资源勘探开发过程中极为重要的一个环节。
目前固井质量评价一般是将超声波测井仪器居中下放到套管中，使超声波测井仪器中搭载的可360度旋转的超声换能器发射出的超声波垂直入射到套管中并激发套管共振，由接收器接收回波信号。脉冲回波信号反映了套管和地层之间环空区域材料的声阻抗特性，一般气体的声阻抗比较小接近于0MRay液体声阻抗在0-3MRay之间固体声阻抗在2MRay以上，通过声阻抗可以判断环空区域材料的种类并进而判断固井质量。利用专门的声阻抗反演算法可以从脉冲回波波形中计算得到环空区域的声阻抗。
声阻抗反演算法随着测井仪器的更新有过多次发展更新，目前常用的水泥声阻抗反演算法存在如下的缺点：需要提前已知泥浆声阻抗的值，然而提前测量套管井中的泥浆声阻抗是不太现实的，因此给出的解决方案是假设套管井中的泥浆声阻抗为一固定不变的值。但是实际情况下泥浆会因为重力作用不断在井底沉积，导致泥浆声阻抗随着深度的增加而增大。这样的假设导致反演得到的水泥声阻抗值必然存在误差。
发明内容
本发明的目的在于克服当前常用声阻抗反演算法的不足，提供一种同时反演水泥声阻抗、泥浆声阻抗的计算方法，为固井质量评价提供有利的保障。
本发明的目的是通过以下技术方案来实现的：
一种水泥声阻抗、泥浆声阻抗计算方法，方法包括如下步骤：
Step1：设置水泥声阻抗、泥浆声阻抗、套管厚度的值，并利用已有的脉冲回波正演模型计算模拟波形；
Step2：计算模块1，计算原始波形最大值对应的时间以及共振波的索引值；
Step3：计算模块2，计算斜率k和截距b；
Step4：根据步骤一和步骤三得到特征参数k、b与水泥声阻抗、泥浆声阻抗的数据表，根据数据表分析待反演参数(泥浆声阻抗、水泥声阻抗)与特征参数之间的变化规律，并根据变化规律得到待反演参数的迭代方案；
Step5：计算模块3，计算波形中声源主频f0、套管厚度Th、波形最大值对应的时间、共振波的索引范围；
Step6：计算模块2；
Step7：泥浆声阻抗迭代计算模块，利用步骤四得到的迭代方案得到实际的泥浆声阻抗；
Step8：水泥声阻抗计算模块，利用步骤四得到的迭代方案以及步骤七得到的泥浆声阻抗得到实际的水泥声阻抗。
前述的水泥声阻抗、泥浆声阻抗反演方法，其特征在于，在Step1中，利用已有的正演模型计算模拟波形过程包括：模拟波形的频域计算的原理如下：
wave(ω)＝R(ω)F(ω)
其中wave(ω)为模拟波形的频域表现形式，R(ω)为模型在频域上的系统函数，F(ω)为模拟声源，选择雷克子波作为正演模型的声源，声源的频域表达式如下：
Figure BDA0003554338060000021
其中f是频率，范围为0～8MHz，f0是声源主频，主频设置为250kHz，i 是虚数单位，ts＝1.5/f0。其中正演模型的频域系统函数为：
Figure BDA0003554338060000022
其中Zm为泥浆声阻抗、Zc为水泥声阻抗、Zc为套管声阻抗、ω为角频率且ω＝ 2πf，Th为套管厚度，i为虚数单位，水泥声阻抗的设置的变化范围为0-10MRay，步长为0.5MRay，泥浆声阻抗的设置的变化范围为1-3MRay，步长为0.1MRay，对模拟波形的频域表现形式做傅里叶逆变换可得时域模拟波形：
wave(t)＝F-1(wave(ω))
其中wave(t)为模拟波形的时域表现形式，F-1为傅里叶逆变换。
具体的，在Step2中，计算模块1的计算过程包括：
找到原始波形最大的三个极大值点以及对应的时间，对这三个极大值点做二次拟合，并确定二次拟合曲线的最大值点对应的时间，该时间即是波形最大值点对应的时间time_start，计算原始波形的傅里叶变换得到其频谱，确定频谱最大分量对应的频率，并得到该频率的倒数为周期period，根据经验设置共振波存在的时间区间为[time_start+3*period，time_start+8*period]，根据时间区间得到共振波的索引index_resonance。
具体的，在Step3中，计算模块2计算过程包括：
对原始波形做归一化，将原始波形的幅度与原始波形最大绝对值作比值得到归一化波形，对归一化波形计算希尔伯特变化得到归一化波形包络，对归一化波形包络取自然对数，并在共振波区间内对包络做一次线性拟合，拟合直线的斜率为k、拟合直线在t＝time_start处的取值定义为截距b，斜率和截距即为波形的特征参数。
具体的，在Step4中，步骤四中根据步骤一和步骤三得到特征参数k、b与水泥声阻抗、泥浆声阻抗的数据表，根据数据表分析待反演参数(泥浆声阻抗、水泥声阻抗)与特征参数之间的变化规律，并根据变化规律得到待反演参数的迭代方案过程包括：
可视化水泥声阻抗与斜率Zc-k、水泥声阻抗与截距Zc-b、泥浆声阻抗与斜率Zm-k、泥浆声阻抗与截距Zm-b，得到如下规律：①k随水泥声阻抗的增大而减小，并且呈现很强的线性关系；②b几乎不随水泥声阻抗的变化而变化；③k随泥浆声阻抗的增大而减小，并且呈现很强的线性关系；④b随着泥浆声阻抗的增大而增大。据此得到的迭代方案如下：设置预期的水泥声阻抗与泥浆声阻抗变化范围，利用二分法迭代，水泥声阻抗不影响截距，因此先给出大致的水泥声阻抗初值并通过迭代得到泥浆声阻抗值，得到泥浆声阻抗值之后利用已得到的水泥声阻抗和二分法迭代得到水泥声阻抗。
具体的，在Step5中，计算模块3的计算过程包括：
找到原始波形最大的三个极大值点以及对应的时间，对这三个极大值点做二次拟合，并确定二次拟合曲线的最大值点对应的时间，该时间即是波形最大值点对应的时间time_start，对原始波形计算傅里叶变换得到其频谱，确定频谱最大分量对映的频率即为声源主频f0，计算声源主频的倒数得到周期period，设置时间窗[time_start-1.5*period，time_start+1.5*period]用来提取原始波形的反射波部分，并认为反射波部分为原始声源，计算原始声源的频谱，将原始声源频谱与原始波形频谱做差并求差值的最大值对应的频谱为共振频率 fresonance由公式：
Th＝Vs/(2*fresonance)
其中Th为计算得到的套管厚度，Vs为声波在套管中的速度，据经验设置共振波存在的时间区间为[time_start+3*period，time_start+8*period]，根据时间区间得到共振波的索引index_resonance。
具体的，在Step6中，计算模块2计算过程包括：
对原始波形做归一化，将原始波形的幅度与原始波形最大绝对值作比值得到归一化波形，对归一化波形计算希尔伯特变化得到归一化波形包络，对归一化波形包络取自然对数，并在共振波区间内对包络做一次线性拟合，拟合直线的斜率为k、拟合直线在t＝time_start处的取值定义为截距b，斜率和截距即为波形的特征参数
具体的，在Step7中，泥浆声阻抗迭代模块的计算过程包括：
预先设置泥浆声阻抗值的迭代区间1-3MRay，取一个水泥声阻抗的初值以及在迭代区间内的泥浆声阻抗初值，以步骤五计算得到的声源主频作为模拟声源的主频，步骤五计算得到的套管厚度作为正演模型的套管厚度，其中正演模型中声波在套管中的速度以及套管声阻抗为已知固定值且Vs＝5860m/s、Zs＝ 46MRay，通过这些参数计算模拟波形，根据精度要求预先设置泥浆声阻抗迭代次数为5(具体迭代次数这里不做限制)，计算模拟波形的最大值对应的时间、共振波对应的索引，计算模拟波形的斜率以及截距，比较步骤六得到的截距与模拟波形的截距，并根据二分法迭代泥浆声阻抗，将迭代得到的泥浆声阻抗再次用于正演模拟波形直到迭代次数到达预设的迭代次数，取最后一次迭代得到的泥浆声阻抗作为实际泥浆声阻抗。
具体的，在Step8中，水泥声阻抗迭代模块计算过程包括：
预先设置水泥声阻抗的迭代区间0-10MRay，在区间内取一水泥声阻抗初值，泥浆声阻抗的值由步骤七确定，步骤五计算得到的声源主频作为模拟声源的主频，步骤五计算得到的套管厚度作为正演模型的套管厚度，其中正演模型中声波在套管中的速度以及套管声阻抗为已知固定值，通过这些参数计算模拟波形，根据精度要求提前设置好水泥声阻抗迭代次数6(具体迭代次数这里不做限制)，计算模拟波形的最大值对应的时间、共振波对应的索引，计算模拟波形的斜率以及截距，比较步骤六得到的斜率与模拟波形的斜率，并根据二分法迭代水泥声阻抗，将迭代得到的水泥声阻抗再次用于正演模拟波形直到迭代次数到达预设的迭代次数，取最后一次迭代得到的水泥声阻抗作为实际水泥声阻抗。
本发明的有益效果：
本发明通过匹配模拟波形与实际波形在时域上变换后共振波区间的一次拟合直线的斜率与截距同时反演水泥声阻抗与泥浆声阻抗，为固井质量评价提供数据支撑。
附图说明
图1是本发明提供的方法流程图；
图2是是正演模型示意图；
图3是计算模块1的流程图
图4是计算模块2的流程图；
图5是计算模块3的流程图；
图6是水泥声阻抗计算模块流程图；
图7是泥浆声阻抗计算模块流程图；
图8是Zm-k、Zm-b、Zc-k、Zc-b的结果图；
图9是模拟数据的计算结果示意图；
图10是实测数据的计算结果示意图。
具体实施方式
为了对本发明的技术特征、目的和有益效果有更加清楚的理解，现对本发明的技术方案精选以下详细说明。显然，所描述的实施案例是本发明一部分实施例，而不是全部实施例，不能理解为对本发明可实施范围的限定。基于本发明的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的其他所有实施例，都属于本发明的保护范围。
实施例一：
本实施例中，如图1所示，一种水泥声阻抗、泥浆声阻抗计算方法包括以下步骤：
Step1：设置水泥声阻抗、泥浆声阻抗、套管厚度值，并利用已有的脉冲回波正演模型计算模拟波形。本实施例中，如图2所示，正演模型包含频域系统函数R(ω)、频域正演波形wave(t)、时域正演波形wave(t)和声源模型F(ω)，其中频域系统函数R(ω)分析处理泥浆声阻抗(Zm)、水泥声阻抗(Zc)、套管声阻抗(Zs)、套管厚度(Th)和套管中的声速(Vs)。声源模型F(ω)分析处理声源主频、声源频谱范围和声源幅度。
模拟波形的频域计算的原理如下：
wave(ω)＝R(ω)F(ω)
其中wave(ω)为模拟波形的频域表现形式，R(ω)为模型在频域上的系统函数，F(ω)为模拟声源，选择雷克子波作为正演模型的声源，声源的频域表达式如下：
Figure BDA0003554338060000061
/>
其中f是频率，范围为0～8MHz，f0是声源主频，主频设置为250kHz，i 是虚数单位，ts＝1.5/f0。其中正演模型的频域系统函数为：
Figure BDA0003554338060000062
其中Zm为泥浆声阻抗、Zc为水泥声阻抗、Zs为套管声阻抗、ω为角频率且ω＝2πf，Th为套管厚度，i为虚数单位，水泥声阻抗的设置的变化范围为 0-10MRay，步长为0.5MRay，泥浆声阻抗的设置的变化范围为1-3MRay，步长为0.1MRay，对模拟波形的频域表现形式做傅里叶逆变换可得时域模拟波形：
wave(t)＝F-1(wave(ω))
其中wave(t)为模拟波形的时域表现形式，F-1为傅里叶逆变换。
Step2：计算模块1
如图3所示，找到原始波形最大的三个极大值点以及对应的时间，对这三个极大值点做二次拟合，并确定二次拟合曲线的最大值点对应的时间，该时间即是波形最大值点对应的时间time_start，计算原始波形的傅里叶变换得到其频谱，确定频谱最大分量对应的频率，并得到该频率的倒数为周期period，根据经验设置共振波存在的时间区间为[time_start+3*period，time_start+8* period]，根据时间区间得到共振波的索引index_resonance。
Step3：计算模块2
如图4所示，对原始波形做归一化，将原始波形的幅度与原始波形最大绝对值作比值得到归一化波形，对归一化波形计算希尔伯特变化得到归一化波形包络，对归一化波形包络取自然对数，并在共振波区间内对包络做一次线性拟合，拟合直线的斜率为k、拟合直线在t＝time_start处的取值定义为截距b，斜率和截距即为波形的特征参数；
Step4：根据步骤一和步骤三得到特征参数k、b与水泥声阻抗、泥浆声阻抗的数据表，根据数据表分析待反演参数(泥浆声阻抗、水泥声阻抗)与特征参数之间的变化规律，并根据变化规律得到待反演参数的迭代方案。
如图8所示，可视化水泥声阻抗与斜率Zc-k、水泥声阻抗与截距Zc-b、泥浆声阻抗与斜率Zm-k、泥浆声阻抗与截距Zm-b，得到如下规律：①k随水泥声阻抗的增大而减小，并且呈现很强的线性关系；②b几乎不随水泥声阻抗的变化而变化；③k随泥浆声阻抗的增大而减小，并且呈现很强的线性关系；④b 随着泥浆声阻抗的增大而增大。如图6所示，据此得到的迭代方案如下：设置预期的水泥声阻抗与泥浆声阻抗变化范围，利用二分法迭代，水泥声阻抗不影响截距，因此先给出大致的水泥声阻抗初值并通过迭代得到泥浆声阻抗值，得到泥浆声阻抗值之后利用已得到的水泥声阻抗和二分法迭代得到水泥声阻抗。
Step5：计算模块3
如图5所示，找到原始波形最大的三个极大值点以及对应的时间，对这三个极大值点做二次拟合，并确定二次拟合曲线的最大值点对应的时间，该时间即是波形最大值点对应的时间time_start，对原始波形计算傅里叶变换得到其频谱，确定频谱最大分量对映的频率即为声源主频f0，计算声源主频的倒数得到周期period，设置时间窗[time_start-1.5*period，time_start+1.5*period]用来提取原始波形的反射波部分，并认为反射波部分为原始声源，计算原始声源的频谱，将原始声源频谱与原始波形频谱做差并求差值的最大值对应的频谱为共振频率fresonance由公式：
Th＝Vs/(2*fresonance)
其中Th为计算得到的套管厚度，Vs为声波在套管中的速度，据经验设置共振波存在的时间区间为[time_start+3*period，time_start+8*period]，根据时间区间得到共振波的索引index_resonance；
Step6：计算模块2
对原始波形做归一化，将原始波形的幅度与原始波形最大绝对值作比值得到归一化波形，对归一化波形计算希尔伯特变化得到归一化波形包络，对归一化波形包络取自然对数，并在共振波区间内对包络做一次线性拟合，拟合直线的斜率为k、拟合直线在t＝time_start处的取值定义为截距b，斜率和截距即为波形的特征参数；
Step7：泥浆声阻抗迭代模块
如图7所示，预先设置泥浆声阻抗值的迭代区间1-3MRay，取一个水泥声阻抗的初值以及在迭代区间内的泥浆声阻抗初值，以步骤五计算得到的声源主频作为模拟声源的主频，步骤五计算得到的套管厚度作为正演模型的套管厚度，其中正演模型中声波在套管中的速度以及套管声阻抗为已知固定值且Vs＝ 5860m/s、Zs＝46Mray，通过这些参数计算模拟波形，根据精度要求预先设置泥浆声阻抗迭代次数为5(具体迭代次数这里不做限制)，计算模拟波形的最大值对应的时间、共振波对应的索引，计算模拟波形的斜率以及截距，比较步骤六得到的截距与模拟波形的截距，并根据二分法迭代泥浆声阻抗，将迭代得到的泥浆声阻抗再次用于正演模拟波形直到迭代次数到达预设的迭代次数，取最后一次迭代得到的泥浆声阻抗作为实际泥浆声阻抗。
Step8：水泥声阻抗迭代模块
预先设置水泥声阻抗的迭代区间0-10MRay，在区间内取一水泥声阻抗初值，泥浆声阻抗的值由步骤七确定，步骤五计算得到的声源主频作为模拟声源的主频，步骤五计算得到的套管厚度作为正演模型的套管厚度，其中正演模型中声波在套管中的速度以及套管声阻抗为已知固定值，通过这些参数计算模拟波形，根据精度要求提前设置好水泥声阻抗迭代次数6(具体迭代次数这里不做限制)，计算模拟波形的最大值对应的时间、共振波对应的索引，计算模拟波形的斜率以及截距，比较步骤六得到的斜率与模拟波形的斜率，并根据二分法迭代水泥声阻抗，将迭代得到的水泥声阻抗再次用于正演模拟波形直到迭代次数到达预设的迭代次数，取最后一次迭代得到的水泥声阻抗作为实际水泥声阻抗。
本实施例中，模块1、模块2以及模块3之间并不存在直接的关系，只是为了方便的描述该方法而集成的模块。
为了验证本实施例提供的水泥声阻抗、泥浆声阻抗计算方法具有良好的水泥声阻抗、泥浆声阻抗反演效果，我们分别对模拟脉冲回波数据与实测脉冲回波数据进行试算，分别为实施例二和实施例三。
实施例二：
本实施例为模拟脉冲回波数据实施案例，具体试算流程如下：
(1)设置水泥声阻抗、泥浆声阻抗、套管厚度值，并利用已有的脉冲回波正演模型计算模拟波形；
(2)计算模块1，计算原始波形最大值对应的时间以及共振波的索引值；
(3)计算模块2，计算斜率k和截距b；
(4)根据步骤一和步骤三得到特征参数k、b与水泥声阻抗、泥浆声阻抗的数据表，根据数据表分析待反演参数(泥浆声阻抗、水泥声阻抗)与特征参数之间的变化规律，并根据变化规律得到待反演参数的迭代方案；
(5)计算模块3，计算波形中声源主频f0、套管厚度Th、波形最大值对应的时间、共振波的索引范围；
(6)计算模块2；
(7)泥浆声阻抗迭代计算模块，利用步骤四得到的迭代方案得到实际的泥浆声阻抗；
(8)水泥声阻抗计算模块，利用步骤四得到的迭代方案以及步骤七得到的泥浆声阻抗反演得到实际的水泥声阻抗。
图9为模拟脉冲回波数据处理结果图，从图中可以看到计算得到的水泥声阻抗以及泥浆声阻抗结果与理论值能够很好的吻合，验证了本算法的有效性。
实施例三：
本实施例为实测脉冲回波数据实施案例，具体计算流程如下：
(1)设置水泥声阻抗、泥浆声阻抗、套管厚度值，并利用已有的脉冲回波正演模型计算模拟波形；
(2)计算模块1，计算原始波形最大值对应的时间以及共振波的索引值；
(3)计算模块2，计算斜率k和截距b；
(4)根据步骤一和步骤三得到特征参数k、b与水泥声阻抗、泥浆声阻抗的数据表，根据数据表分析待反演参数(泥浆声阻抗、水泥声阻抗)与特征参数之间的变化规律，并根据变化规律得到待反演参数的迭代方案；
(5)计算模块3，计算波形中声源主频f0、套管厚度Th、波形最大值对应的时间、共振波的索引范围；
(6)计算模块2；
(7)泥浆声阻抗迭代计算模块，利用步骤四得到的迭代方案得到实际的泥浆声阻抗；
(8)水泥声阻抗计算模块，利用步骤四得到的迭代方案以及步骤七得到的泥浆声阻抗反演得到实际的水泥声阻抗，该结果如图10所示。
由此可见，本实施例提供的水泥声阻抗、泥浆声阻抗方法，为固井质量评价提供有力的保障。
以上显示和描述了本发明的基本原理、主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护的范围由所附的权利要求书及其等效物界定。

Claims (1)

1.一种水泥声阻抗、泥浆声阻抗的计算方法，其特征在于，包括：
步骤S1：设置水泥声阻抗、泥浆声阻抗和套管厚度值，并利用已有的脉冲回波正演模型计算模拟波形，包括：模拟波形的频域计算的原理如下：
wave(ω)＝R(ω)F(ω)
其中wave(ω)为模拟波形的频域表现形式，R(ω)为模型在频域上的系统函数，F(ω)为模拟声源，选择雷克子波作为正演模型的声源，声源的频域表达式如下：
Figure FDA0004186751940000011
其中f是频率，范围为0～8MHz，f0是声源主频，主频设置为250kHz，i是虚数单位，ts＝1.5/f0；其中正演模型的频域系统函数为：
Figure FDA0004186751940000012
其中Zm为泥浆声阻抗、Zc为水泥声阻抗、Zs为套管声阻抗、ω为角频率且ω＝2πf，Th为套管厚度，i为虚数单位，水泥声阻抗的设置的变化范围为0-10MRay，步长为0.5MRay，泥浆声阻抗的设置的变化范围为1-3MRay，步长为0.1MRay，对模拟波形的频域表现形式做傅里叶逆变换可得时域模拟波形：
wave(t)＝F-1(wave(ω))
其中wave(t)为模拟波形的时域表现形式，F-1为傅里叶逆变换；
步骤S2：计算模块1，计算原始波形最大值对应的时间以及共振波的索引值计算过程包括：找到原始波形最大的三个极大值点以及对应的时间，对这三个极大值点做二次拟合，并确定二次拟合曲线的最大值点对应的时间，该时间即是波形最大值点对应的时间time_start，计算原始波形的傅里叶变换得到其频谱，确定频谱最大分量对应的频率，并得到该频率的倒数为周期period，根据经验设置共振波存在的时间区间为[time_start+3*period，time_start+8*period]，根据时间区间得到共振波的索引index_resonance；
步骤S3：计算模块2，计算斜率k和截距b，计算过程包括：对原始波形做归一化，将原始波形的幅度与原始波形最大绝对值作比值得到归一化波形，对归一化波形计算希尔伯特变化得到归一化波形包络，对归一化波形包络取自然对数，并在共振波区间内对包络做一次线性拟合，拟合直线的斜率为k、拟合直线在t＝time_start处的取值定义为截距b，斜率k和截距b即为波形的特征参数；
步骤S4：根据步骤S1和步骤S3得到特征参数k、b与水泥声阻抗、泥浆声阻抗的数据表，根据数据表分析待反演参数与特征参数之间的变化规律，并根据变化规律得到待反演参数的迭代方案，包括：可视化水泥声阻抗与斜率Zc-k、水泥声阻抗与截距Zc-b、泥浆声阻抗与斜率Zm-k、泥浆声阻抗与截距Zm-b，得到如下规律：①k随水泥声阻抗的增大而减小，并且呈现很强的线性关系；②b几乎不随水泥声阻抗的变化而变化；③k随泥浆声阻抗的增大而减小，并且呈现很强的线性关系；④b随着泥浆声阻抗的增大而增大；据此得到的迭代方案如下：设置预期的水泥声阻抗与泥浆声阻抗变化范围，利用二分法迭代，水泥声阻抗不影响截距，先给出水泥声阻抗初值并通过迭代得到泥浆声阻抗值，得到泥浆声阻抗值之后利用已得到的水泥声阻抗和二分法迭代得到水泥声阻抗；
步骤S5：计算模块3，计算波形中声源主频f0、套管厚度Th、波形最大值对应的时间、共振波的索引范围，计算过程包括：找到原始波形最大的三个极大值点以及对应的时间，对这三个极大值点做二次拟合，并确定二次拟合曲线的最大值点对应的时间，该时间即是波形最大值点对应的时间time_start，对原始波形计算傅里叶变换得到其频谱，确定频谱最大分量对映的频率即为声源主频f0，计算声源主频的倒数得到周期period，设置时间窗[time_start-1.5*period，time_start+1.5*period]用来提取原始波形的反射波部分，并认为反射波部分为原始声源，计算原始声源的频谱，将原始声源频谱与原始波形频谱做差并求差值的最大值对应的频谱为共振频率fresonance由公式：
Th＝Vs/(2*fresonance)
其中Th为计算得到的套管厚度，Vs为声波在套管中的速度，据经验设置共振波存在的时间区间为[time_start+3*period，time_start+8*period]，根据时间区间得到共振波的索引index_resonance；
步骤S6：计算模块2，包括：对原始波形做归一化，将原始波形的幅度与原始波形最大绝对值作比值得到归一化波形，对归一化波形计算希尔伯特变化得到归一化波形包络，对归一化波形包络取自然对数，并在共振波区间内对包络做一次线性拟合，拟合直线的斜率为k、拟合直线在t＝time_start处的取值定义为截距b，斜率和截距即为波形的特征参数；
步骤S7：泥浆声阻抗迭代计算模块，利用步骤S4得到的迭代方案得到实际的泥浆声阻抗，计算过程包括：预先设置泥浆声阻抗值的迭代区间1-3MRay，取一个水泥声阻抗的初值以及在迭代区间内的泥浆声阻抗初值，以步骤五计算得到的声源主频作为模拟声源的主频，步骤五计算得到的套管厚度作为正演模型的套管厚度，其中正演模型中声波在套管中的速度以及套管声阻抗为已知固定值且Vs＝5860m/s、Zs＝46MRay，通过这些参数计算模拟波形，根据精度要求预先设置泥浆声阻抗迭代次数为5，计算模拟波形的最大值对应的时间、共振波对应的索引，计算模拟波形的斜率以及截距，比较步骤六得到的截距与模拟波形的截距，并根据二分法迭代泥浆声阻抗，将迭代得到的泥浆声阻抗再次用于正演模拟波形直到迭代次数到达预设的迭代次数，取最后一次迭代得到的泥浆声阻抗作为实际泥浆声阻抗；
步骤S8：水泥声阻抗计算模块，利用步骤S4得到的迭代方案以及步骤S7得到的泥浆声阻抗反演得到实际的水泥声阻抗，计算过程包括：预先设置水泥声阻抗的迭代区间0-10MRay，在区间内取一水泥声阻抗初值，泥浆声阻抗的值由步骤七确定，步骤五计算得到的声源主频作为模拟声源的主频，步骤五计算得到的套管厚度作为正演模型的套管厚度，其中正演模型中声波在套管中的速度以及套管声阻抗为已知固定值，通过这些参数计算模拟波形，根据精度要求提前设置好水泥声阻抗迭代次数6，计算模拟波形的最大值对应的时间、共振波对应的索引，计算模拟波形的斜率以及截距，比较步骤六得到的斜率与模拟波形的斜率，并根据二分法迭代水泥声阻抗，将迭代得到的水泥声阻抗再次用于正演模拟波形直到迭代次数到达预设的迭代次数，取最后一次迭代得到的水泥声阻抗作为实际水泥声阻抗。
CN202210272669.0A 2022-03-18 2022-03-18 一种水泥声阻抗、泥浆声阻抗的计算方法 Active CN114740532B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210272669.0A CN114740532B (zh) 2022-03-18 2022-03-18 一种水泥声阻抗、泥浆声阻抗的计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210272669.0A CN114740532B (zh) 2022-03-18 2022-03-18 一种水泥声阻抗、泥浆声阻抗的计算方法

Publications (2)

Publication Number Publication Date
CN114740532A CN114740532A (zh) 2022-07-12
CN114740532B true CN114740532B (zh) 2023-05-30

Family

ID=82277225

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210272669.0A Active CN114740532B (zh) 2022-03-18 2022-03-18 一种水泥声阻抗、泥浆声阻抗的计算方法

Country Status (1)

Country Link
CN (1) CN114740532B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116660997B (zh) * 2023-08-02 2023-09-29 中海油田服务股份有限公司 套管内外介质声阻抗反演方法、装置及电子设备

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP0314573A2 (en) * 1987-10-30 1989-05-03 Schlumberger Limited Well logging apparatus and method
EP1635193A1 (en) * 2004-09-13 2006-03-15 Services Petroliers Schlumberger Determining impedance of material behind a casing in a borehole
EP2803816A1 (en) * 2013-05-16 2014-11-19 Services Pétroliers Schlumberger Systems and Methods for Cement Evaluation
CN104747178A (zh) * 2015-02-10 2015-07-01 中国海洋石油总公司 一种获取声阻抗的方法和装置
CN106383173A (zh) * 2016-09-13 2017-02-08 中国海洋石油总公司 一种水泥声阻抗计算方法和装置
CN110159253A (zh) * 2019-06-20 2019-08-23 太平洋远景石油技术(北京)有限公司 一种超声波成像测井方法
WO2019200242A1 (en) * 2018-04-12 2019-10-17 Schlumberger Technology Corporation Evaluating casing cement using automated detection of clinging compression wave (p) arrivals
CN110806601A (zh) * 2019-11-08 2020-02-18 中国海洋石油集团有限公司 套内外物质阻抗计算、关系确定、固井评价方法和设备

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20140052376A1 (en) * 2012-08-15 2014-02-20 Pingjun Guo Method for Cement Evaluation with Acoustic and Nuclear Density Logs
US11397081B2 (en) * 2019-07-02 2022-07-26 Baker Hughes Oilfield Operations Llc Method and apparatus for determining a tubular thickness using a pulse echo waveform signal

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP0314573A2 (en) * 1987-10-30 1989-05-03 Schlumberger Limited Well logging apparatus and method
EP1635193A1 (en) * 2004-09-13 2006-03-15 Services Petroliers Schlumberger Determining impedance of material behind a casing in a borehole
EP2803816A1 (en) * 2013-05-16 2014-11-19 Services Pétroliers Schlumberger Systems and Methods for Cement Evaluation
CN104747178A (zh) * 2015-02-10 2015-07-01 中国海洋石油总公司 一种获取声阻抗的方法和装置
CN106383173A (zh) * 2016-09-13 2017-02-08 中国海洋石油总公司 一种水泥声阻抗计算方法和装置
WO2019200242A1 (en) * 2018-04-12 2019-10-17 Schlumberger Technology Corporation Evaluating casing cement using automated detection of clinging compression wave (p) arrivals
CN110159253A (zh) * 2019-06-20 2019-08-23 太平洋远景石油技术(北京)有限公司 一种超声波成像测井方法
CN110806601A (zh) * 2019-11-08 2020-02-18 中国海洋石油集团有限公司 套内外物质阻抗计算、关系确定、固井评价方法和设备

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
套管中模式波的响应特征;陈雪莲;唐晓明;刘临政;;应用声学(01);第2-8页 *
水泥混凝土电化学进展──交流阻抗谱理论;许仲梓;硅酸盐学报(02);第173-180页 *

Also Published As

Publication number Publication date
CN114740532A (zh) 2022-07-12

Similar Documents

Publication Publication Date Title
US10345465B2 (en) Resonance-based inversion of acoustic impedance of annulus behind casing
US9829597B2 (en) Model based inversion of acoustic impedance of annulus behind casing
US10012749B2 (en) Fast model based inversion of acoustic impedance of annulus behind casing
CA2490936C (en) Method and apparatus for determining earth formation shear-wave transverse isotropy from borehole stoneley-wave measurements
CN110348135B (zh) 一种随钻声波测井评价地层渗透率的方法
CN110318740B (zh) 一种随钻声波测井评价地层各向异性的方法
GB2421079A (en) Acoustic waveform processing using dispersion curves based on a borehole-formation model having a set of borehole-formation parameters
MX2014007352A (es) Uso de armonicas de orden superior para determinar la impedancia acustica de material por detras de la tuberia de revestimiento de los agujeros.
CN106595538B (zh) 一种检测套管外水泥与套管间裂隙宽度的系统及方法
CN111980676B (zh) 一种阵列声波测井评价固井质量的方法及处理装置
CN101285381B (zh) 一种泄漏模式波反演软地层横波速度的方法
CN114740532B (zh) 一种水泥声阻抗、泥浆声阻抗的计算方法
CN115992691A (zh) 基于超声兰姆波的固井质量检测方法及装置
CN115324564B (zh) 固井质量检测方法、装置、计算设备及存储介质
CN114233276B (zh) 基于套管井响应的阵列声波测井固井质量评价解释图版
CN106383173B (zh) 一种水泥声阻抗计算方法和装置
CN114017014B (zh) 声波测井波动声学方法
CN110967742A (zh) 一种孔隙度反演方法及系统
CN114139584A (zh) 一种随钻声波测井套管波衰减反演方法和装置
CN111335888A (zh) 一种确定地质层中井筒性质的方法
CN113356834B (zh) 井眼声波时差的校正方法、校正装置和电子设备
CN116241239B (zh) 基于远近单极的固井评价的方法、装置、设备及存储介质
US11662496B2 (en) Through tubing acoustic measurements
CN112228046B (zh) 一种随钻超声井径数据校正方法
Yang et al. A parameterized-model-based radial profiling for formation shear slowness in cased boreholes

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