全空间瞬变电磁全区视电阻率优化二分搜索算法

姜国庆,程久龙,孙晓云,陈 丁,王海彦

(1.江苏省地质调查研究院基础地质研究所,江苏南京 210018;2.中国矿业大学(北京)煤炭资源与安全开采国家重点实验室,北京 100083)

全空间瞬变电磁全区视电阻率优化二分搜索算法

姜国庆1,2,程久龙2,孙晓云2,陈 丁2,王海彦2

(1.江苏省地质调查研究院基础地质研究所,江苏南京 210018;2.中国矿业大学(北京)煤炭资源与安全开采国家重点实验室,北京 100083)

将优化二分搜索算法应用于全区视电阻率计算,以期提高计算精度和运算效率。以均匀全空间回线源中心的瞬变响应为基础,探讨全空间磁感应强度垂向分量Bz及其时间变化率∂Bz/∂T核函数曲线特征,推导晚期视电阻率表达式。给出全空间条件下由∂Bz/∂T数值积分计算Bz的方法,并通过模拟计算对转换精度进行分析。对瞬变延时T一定时Bz和∂Bz/∂T相对于电阻率ρ的单调性进行研究,结果表明Bz为ρ的单值函数,满足对ρ二分搜索的条件,而∂Bz/∂T为ρ的双值函数,其二分搜索的前提是辅助参数Z在变量ρ的取值区间内恒满足Z≤1.225,并在此基础上对二分搜索的初始区间和具体流程进行了优化。模型算例和应用实例表明,该优化算法计算速度快、精度高,并可同时获得由Bz和∂Bz/∂T定义的全区视电阻率ρall和ρ′all。

全空间瞬变电磁法;中心回线装置;全区视电阻率;二分搜索算法;优化

全空间瞬变电磁法属于近距离观测方法,具有体积效应相对较小、方向性强、分辨率高、对低阻体敏感、施工效率高等优势,因此,在矿井水文地质勘查中得到了广泛应用[1-5]。但是相对而言全空间瞬变电磁法理论的研究相对滞后,以中心回线和重叠回线装置为例,目前,这两种常用装置的视电阻率计算还多是采用半空间计算公式或经校正后的晚期经验公式[1,6],已有全空间全区视电阻率的研究[7]也仅对感应电动势定义全区视电阻率进行了重点探讨,但是理论和实践[8-9]都证明磁场与感应电动势定义全区视电阻率在数值和解释深度上都存在差异,并且在反映地下结构和解释深度方面磁场要明显优于感应电动势。

二分搜索算法是求非线性方程根的常用方法,因其算法稳定、易于实现而得到广泛的应用。Peter George[10]在磁测数据反演中,为了求得著名的吉洪诺夫曲线中合适的折衷参数β,分别使用了插值法和二分法,结果,二分搜索算法因不要求曲线的曲率而更为稳健;Leonard[11]在研究时间域电磁法探测未爆炸的爆炸物时,使用了二分法来求解超越方程的根qs,该参数控制着均匀场中导电球体的时间域二次场响应的表达式。张成范等[12]考虑到二分法对初始隔根区间的限制低而使用二分搜索算法来计算均匀半空间表面矩形大定源回线瞬变电磁测深数据的全区视电阻率;李文尧等[13]通过核函数变量Z的二分法计算中心回线瞬变电磁全期视电阻率。陈清礼等[14]、刘俊等[15]采用对电阻率ρ的二分搜索计算感应电动势定义全区视电阻率,但未对函数单调区间及二分法初始区间进行深入探讨。笔者将二分搜索算法应用到全空间瞬变电磁法全区视电阻率计算中,同时计算Bz和∂Bz/∂T定义的全区视电阻率,并在深入探讨核函数单调性的基础上对算法流程进行优化,使其运算速度更快,适应性更强。

1 全空间中心回线装置全区视电阻率

均匀全空间介质中水平放置的圆形发射框中心的瞬变响应(磁感应强度垂向分量Bz及其时间变化率∂Bz/∂T)可以表示为

此时,电阻率可以表示为

式(1),(2)可表示为

其中,Y(Z)和Y′(Z)分别为全空间中心回线装置Bz和∂Bz/∂T的核函数。图1为根据式(6),(7)所绘制的核函数曲线,从图1可看出,Bz的核函数Y(Z)为单值函数,其相应函数曲线为单调上升曲线,渐进值Y(Z→∞)≈1。由图1可知当实测磁场值Bz小于一次场值(μI/2a)时,由Bz总能求得单值的全区视电阻率值。∂Bz/∂T的核函数Y′(Z)为双值函数,其函数曲线存在一个转折点Z=1.225,该转折点处Y′(Z)取得极大值Y′max=0.41。由图1可知当实测归一化Y′(Z)＞0.41时,全区视电阻率ρ无定义;当Y′(Z)= 0.41时,全区视电阻率ρ存在惟一解;Y′(Z)＜0.41时,全区视电阻率ρ存在分别对应早期和晚期部分的两个解。

晚期条件下Z→0,通过对误差函数和e-Z2进行泰勒级数展开,并省略Z的高次项可以得到

由式(3),(4),(8)可得到由Bz计算的晚期视电阻率表达式

图1 全空间中心回线装置核函数曲线

Fig.1 Kernel function curves of full space central loop unit

由式(3),(5),(9)可得到由∂Bz/∂T计算的晚期视电阻率表达式

2 磁感应强度Bz的计算

瞬变电磁系统一般以固定的采样频率进行信号采集(如TerraTEM系统f=500 kHz),记录数据为时窗范围内采样值的算数平均值。假设某瞬变电磁系统具有N个时窗,对应N+1个时窗边界,定义每个时窗中心所对应延迟时间分别为T1,T2,…,TN,时窗边界所对应延迟时间分别为TB1,TB2,…,TBN+1。根据电磁学原理Bz和∂Bz/∂T的转换关系可以表示为

其中,Bz(TBN+1)为TBN+1时刻的磁场值,磁场计算的关键就是根据观测∂Bz/∂T变化趋势准确估计TBN+1时刻的磁场值。由式(5),(9)可以推导近区时均匀全空间条件下∂Bz/∂T具有与均匀半空间相同的衰减规律,即:∂Bz/∂T∝T-5/2。通过计算可知全空间瞬变电磁法采用对数域线性回归方法[8]计算Bz(TBN+1)的精度较高,这里结合具体算例对转换精度进行分析。

假设均匀全空间电阻率为100 Ω·m,发射线圈为边长2 m的方框(等效半径a=2/π),发射电流为1 A,采样频率f=500 kHz,记录时窗为Terra TEM系统的中间时间序列(intermediate time series)的前100个窗口(Windows),舍弃延迟时间小于10 μs的时窗。图2中∂Bz/∂T曲线为式(5)理论值按f采样所得,采用对数域线性回归方法计算Bz(TBN+1)值为3.431 95×10-18T,由式(4)计算理论值为3.438 94× 10-18T,二者相对误差仅为0.203%。由图2中Bz积分转换值与理论值对比可以看出,二者基本一致。同时,由相对误差曲线可知所有时窗积分转换的相对误差均在0.6%以下,并且表现为两端高、中间低的特点,分析晚期相对误差高值是由Bz(TBN+1)的估算误差引起,早期相对误差高值是由时窗窄、采样相对不足引起。

图2 磁场Bz计算值与理论值对比曲线Fig.2 Contrastive curves of Bzcalculatedvalue and theoretical value

3 函数单调性探讨

由式(1),(2)可知磁感应强度垂向分量Bz(T, ρ)及其对时间变化率∂Bz(T,ρ)/∂T均为电阻率ρ的非线性复杂函数,要通过对ρ的二分搜索来计算全区视电阻率就必须对两函数相对于变量ρ的单调性进行探讨。

3.1 Bz(T,ρ)的单调性探讨

根据Bz的核函数曲线特征,Y(Z)为Z的单调增函数。针对特定测点的瞬变电磁观测过程,发射线圈半径a和发射电流I一定时,由式(4)可知Bz亦为Z的单调增函数。同时,对于某一瞬变延迟时间T,Z与ρ成反比,综合以上可以判断磁感应强度Bz为电阻率ρ的单调减函数,满足对ρ的二分搜索来计算全区视电阻率的条件。如图3(a)为磁场Bz随电阻率ρ的变化曲线,由图中可以看出对应任意瞬变延迟时间T,Bz均为ρ的单调减函数。

3.2 ∂Bz(T,ρ)/∂T的单调性探讨

根据∂Bz/∂T的核函数曲线特征,Y′(Z)为Z的双值函数,因此,针对特定测点观测的某一瞬变延迟时间T(即参数a,I和T均为常数时),∂Bz/∂T与Y′(Z)相对变量Z具有相同的函数单调性,即均为Z的双值函数,分区间单调。同时,由Z与ρ的反比关系,可以判定在上述条件下∂Bz/∂T亦为电阻率ρ的双值函数。为了满足二分搜索算法对函数单调性的要求,下面对函数∂Bz/∂T的单调区间进行分析。

图3 磁场Bz和磁场对时间导数∂Bz/∂T随电阻率ρ的变化曲线Fig.3 Curves of Bzand∂Bz/∂T changing with resistivity ρ

由式(2),(5)可知∂Bz/∂T是电阻率ρ的非线性复杂函数,通过前面对核函数Y′(Z)的探讨,可以从理论上对函数单调特性的实质进行分析。在参数a, I和T一定的条件下,∂Bz/∂T与Y′(Z)相对变量Z具有相同的函数单调性,因此,二者曲线具有相同的转折点,即Z=1.225,以转折点为界,Z＜1.225对应瞬变延时的晚期部分,随着Z值的减小∂Bz/∂T单调下降;Z＞1.225对应瞬变延时早期部分,随着Z值的增大∂Bz/∂T单调下降。由此,可以得到∂Bz/∂T相对于电阻率ρ的函数单调性的实质是:在特定装置(a为常数)、选定瞬变延时(T为常数)条件下,Z在变量ρ的取值区间内是否恒满足Z≤1.225或者Z≥1.225的条件,在实际勘探中笔者更关注瞬变延时晚期部分,即Z≤1.225的情况。

如果假设方形发射线框边长L=2 m,全空间视电阻率ρ取值区间定为[10-2,104]Ω·m,由Z= aμ/4ρT 可以计算要使得Z在ρ的取值区间内恒满足Z≤1.225的条件,就必须使瞬变延时T≥26.66 μs。物理模拟实验[16]和现场实测结果均表明多匝重叠小回线装置的关断时间一般为n×102μs,且关断时间终点附近观测数据会出现明显震荡,真正有效观测数据的起始瞬变延时T远大于上述临界值,因此,二分搜索算法亦适用于∂Bz/∂T定义全区视电阻率的计算。图3(b)为∂Bz/∂T随电阻率ρ的变化曲线,可以看出,当瞬变延迟时间T=0.001和0.01 ms时,∂Bz/∂T不满足区间内单调的条件,存在极大值点;当瞬变延迟时间T≥0.1 ms时,∂Bz/∂T在ρ的取值区间内单调下降;瞬变延迟时间T= 0.026 66 ms为临界状态。

4 优化二分搜索算法的实现

二分搜索算法为数值分析中二分法的一种应用,它的主要思想就是“取中”和“比较”。为了实现对电阻率ρ的二分搜索,首先必须确定其取值区间[ρa, ρb],文献[14-15]采用固定区间[0.1,10 000],由文献[17]可知常见岩、矿石的电阻率值变化区间为[10-n,10n],跨度达几个数量级,因此,采用固定区间一方面可能导致视电阻率超出区间范围使计算无解,另一方面会增加二分搜索运算量,因此,有必要对二分区间进行优化。

已知Bz为ρ的单值函数,由Bz求全区视电阻率的区间最小值ρa可以取为0;∂Bz/∂T为ρ的双值函数,ρa的选取必须满足Z≤1.225的单调条件,因此,对于瞬变延时T,由Zmax=1.225可以计算其区间最小值ρZ=1.225。通过晚期核函数YL(Z)和Y′L(Z)的推导过程发现,对于相同的Z值晚期核函数YL(Z)和Y′L(Z)均大于精确核函数Y(Z)和Y′(Z),从而由式(10),(11)所计算的晚期视电阻率值(ρL,ρ′L)总大于精确计算的全区视电阻率值(ρ,ρ′),并且随着延迟时间T的增大,晚期值逐渐接近于全区值。据此,考虑采用晚期视电阻率(ρL,ρ′L)作为区间上限,从而极大地缩小搜索区间,提高运算速度。笔者将全空间Bz和∂Bz/∂T二分搜索的取值区间分别优化为[0, ρL]和[ρZ=1.225,ρ′L]。以∂Bz/∂T计算全区视电阻率为例,经优化后的计算步骤如下:

(1)对T(i)时刻观测值∂Bz(i)/∂T进行关断时间校正。并由式(5),(7)计算Z=1.225时磁感应强度变化率的理论极大值∂Bzmax(i)/∂T。

(2)比较∂Bz(i)/∂T与∂Bzmax(i)/∂T的相对大小。如果∂Bz(i)/∂T=∂Bzmax(i)/∂T,则全区视电阻率ρ′all(i)=ρ′Z=1.225;如果∂Bz(i)/∂T＞∂Bzmax(i)/∂T,则ρ′all(i)无解;如果∂Bz(i)/∂T＜∂Bzmax(i)/∂T,则进行下一步。

(3)由式(3)计算ρ′Z=1.225,由式(11)计算晚期视电阻率ρ′L,从而得到二分搜索初始取值区间ρa= ρ′Z=1.225,ρb=ρ′L。

(4)计算区间[ρa,ρb]的中点ρm=(ρa+ρb)/2,并将ρm带入式(5),(7)计算中值电阻率所对应的理论值∂Bzm(i)/∂T。

(5)计算观测值∂Bz(i)/∂T与理论值∂Bzm(i)/∂T之间的相对误差ε,如果相对误差ε小于给定误差限ε0,则ρ′all(i)=ρm,计算结束。否则,进行下一步。

(6)如果∂Bz(i)/∂T＞∂Bzm(i)/∂T则ρb=ρm,反之,ρa=ρm。计算搜索区间长度δ=ρb-ρa,如果δ小于给定最小区间长度δ0,则ρ′all(i)=ρm,计算结束。否则,返回步骤(4)继续执行。

5 模型算例与应用实例

5.1 模型算例

均匀全空间条件下磁感应强度垂向分量Bz及其时间变化率∂Bz/∂T具有解析解,并且求取的全区视电阻率具有标准值,因此最能检验计算结果的正确性和精度。这里选用均匀全空间模型对算法的精度及其运算效率进行分析。

假设均匀全空间电阻率为ρ=1 Ω·m,瞬变延迟时间T在对数域近似等间隔分布,共100个时窗。由式(10),(11)计算晚期视电阻率ρL和ρ′L;由二分搜索算法计算全区视电阻率ρall和ρ′all,误差限ε0取为10-10。图4为计算结果,由图4可以看出,二分搜索算法计算的全区视电阻率ρall和ρ′all与模型电阻率ρ吻合较好,不同时窗计算全区视电阻率的相对误差呈锯齿状波动,各时窗相对误差最大值εmax仅为6.7× 10-8%。而晚期近似公式计算的ρL和ρ′L均大于全空间电阻率ρ,并且随着T的增加ρL和ρ′L逐渐接近于全空间电阻率ρ,计算值与模型电阻率相对误差在对数域呈近直线下降,最大相对误差εmax为1%。同时,由于晚期视电阻率值与模型电阻率差异不大,采用晚期视电阻率ρL和ρ′L作为二分搜索区间上限可以有效减少搜索次数,该模型计算结果表明,采用[0,ρL]和[ρZ=1.225,ρ′L]作为搜索区间与采用[0.01,100 000]相比分别可以减少二分搜索次数达32.1%和31.0%。

图4 均匀全空间模型计算视电阻率及相对误差曲线Fig.4 Apparent resistivity and relative error curves of homogeneous full space model

5.2 应用实例

为查明某矿工作面截水巷道底板岩层富水性及含水构造分布情况,开展了全空间瞬变电磁法探测研究,测线长度360 m,点距10 m,共37个测点。图5为0号测点视电阻率曲线,从图中可以看出全区视电阻率ρall和ρ′all的曲线类型一致,但是电阻率值存在差异,具体表现为在相对早期ρ′all明显大于ρall,而在晚期阶段二者重合或较为接近。同时可以看出ρ′all曲线变化较为剧烈,而ρall曲线变化则较为平缓。

理论上讲,测点附近全空间介质的电性结构是确定的,无论是磁感应强度垂向分量Bz还是其时间变化率∂Bz/∂T所定义的全区视电阻率应该是一致的,但是实测结果表明全区视电阻率ρall和ρ′all存在差异,特别是早期阶段,差异较为明显。分析可知差异的产生主要有两方面的原因:一是磁感应强度垂向分量Bz及其时间变化率∂Bz/∂T在地下不同介质界面处的表现不同。ρ′all对地下界面的反映比较敏感,但容易引起假薄层,ρall对地下界面反映则相对平缓。二是磁感应强度垂向分量Bz是∂Bz/∂T经过数值积分得到,存在转换误差,因此,由Bz和∂Bz/∂T计算的全区视电阻率ρall和ρ′all也存在差异[18-20]。

图5 全空间瞬变电磁法0号测点视电阻率曲线Fig.5 Apparent resistivity curves of survey station number 0

完成该测线37个测点ρall和ρ′all计算总耗时不到2 min,可见该算法运算速度较快。从图6可以看出ρall和ρ′all的断面特征基本一致,断面上以150,220, 330 m桩号为中心分别存在①,②和③3个低阻异常区,其中以①号异常区的范围和低阻特征最明显,推断该异常对应范围岩层富水性相对较强或存在含水构造;②和③低阻异常区的范围较小、视电阻率略低,推断②和③异常对应范围岩层为弱富水,对工作面回采不构成威胁。同时,可以看出ρall断面上低阻异常区的顶界埋深均在20 m左右,而ρ′all断面上埋深约为30 m,说明ρall的解释深度要相对小于ρ′all。为了验证物探成果的准确性,矿方在巷道145 m处布置了一个底板探测孔,该钻孔在23 m深度时开始出水,但水量不大,钻孔最大出水量为3 m3/h,物探低阻异常区与底板岩层富水区位置吻合较好,并且ρall解释深度与实际情况更接近。

图6 全空间瞬变电磁法视电阻率断面Fig.6 Apparent resistivity sections of full space TEM

6 结 论

(1)对于任意瞬变延时T,Bz均为ρ的单调减函数,均满足二分搜索算法的条件,而∂Bz/∂T为电阻率ρ的双值函数,二分搜索算法成立的条件为Z在变量ρ的取值区间内恒满足Z≤1.225;

(2)将Bz和∂Bz/∂T二分搜索的初始区间分别优化为[0,ρL]和[ρZ=1.225,ρ′L],既可以解决岩、矿石电阻率值区间跨度大的问题,又可以减少运算量。

(3)全空间条件下磁感应强度Bz可以通过∂Bz/∂T的数值积分求得,文中采用的对数域线性回归方法计算精度相对较高。由Bz和∂Bz/∂T计算的全区视电阻率ρall和ρ′all的曲线类型一致,但在电阻率值和解释深度存在差异。实际探测中,同时计算ρall和ρ′all,分别反演获得拟断面图,并结合实际地质资料对比分析,可以获得更好的探测结果。

[1] 于景邨.矿井瞬变电磁法勘探[M].徐州:中国矿业大学出版社, 2007.

[2] 郭 纯,刘白宙,白登海.地下全空间瞬变电磁技术在煤矿巷道掘进头的连续跟踪超前探测[J].地震地质,2006,28(3):456-462.

Guo Chun,Liu Baizhou,Bai Denghai.Predication of water disasters ahead of tunneling in coal mine using continuous detection by UTEM [J].Seismology and Geology,2006,28(3):456-462.

[3] 程久龙,邱 浩,叶云涛,等.矿井瞬变电磁法波场变换与数据处理方法研究[J].煤炭学报,2013,38(9):1646-1650.

Cheng Jiulong,Qiu Hao,Ye Yuntao.Research on wave-field transformation and data processing of the mine transient electromagenetic method[J].Journal of China Coal Society,2013,38(9):1646-1650.

[4] 于景邨,刘志新,刘树才,等.深部采场突水构造矿井瞬变电磁法探查理论及应用[J].煤炭学报,2007,32(8):818-821.

Yu Jingcun,Liu Zhixin,Liu Shucai,et al.Theoretical analysis of mine transient electromagnetic method and its application in detec-ting water burst structures in deep coal stope[J].Journal of China Coal Society,2007,32(8):818-821.

[5] 姜国庆.矿井瞬变电磁法数据处理系统研究[D].青岛:山东科技大学,2011.

[6] 于景邨,刘树才,王扬州.巷道内金属体瞬变电磁响应特征及处理技术[J].煤炭学报,2008,33(12):1403-1407.

Yu Jingcun,Liu Shucai,Wang Yangzhou.Respounse characteristic of transient electromagnetic to metallic facilities in coal mines and the disposal technology[J].Journal of China Coal Society,2008,33 (12):1403-1407.

[7] 杨海燕,邓居智,张 华,等.矿井瞬变电磁法全空间视电阻率解释方法研究[J].地球物理学报,2010,53(3):651-656.

Yang Haiyan,Deng Juzhi,Zhang Hua,et al.Research on full-space apparent resistivity interpretation technique in mine transient electromagnetic method[J].Chinese Journal of Geophysics,2010,53(3): 651-656.

[8] 许建荣,李爱勇,杨 生.TEM中心回线法瞬变磁场求取和全区视电阻率计算[J].地质与勘探,2008,44(6):62-68.

Xu Jianrong,Li Aiyong,Yang Sheng.Calculation of transient magnetic field and all time apparent resistivity based on central TEM loops method[J].Geology and Prospecting,2008,44(6):62-68.

[9] 白登海,Maxwell A Meju,卢 健,等.时间域瞬变电磁法中心方式全程视电阻率的数值计算[J].地球物理学报,2003,46(5): 697-704.

Bai Denghai,Maxwell A Meju,Lu Jian,et al.Numerical calculation of all-time apparent resistivity for the central loop transient electromagnetic method[J].Chinese Journal of Geophysics,2003,46(5): 697-704.

[10] Peter G L.Forward modeling and inversion of geophysical magnetic data[D].Columbia:The University of British Columbia,2003: 121-126.

[11] Leonard R P.Detecting unexploded ordnance with time domain electromagnetic induction[D].Vancouver:Simon Fraser University,1996:15.

[12] 张成范,翁爱华,孙世栋,等.计算矩形大定源回线瞬变电磁测深全区视电阻率[J].吉林大学学报(地球科学版),2009,39 (4):755-758.

Zhang Chengfan,Weng Aihua,Sun Shidong,et al.Computation of whole-time apparent resistivity of large rectangular loop[J].Journal of Jilin University(Earth Science Edition),2009,39(4):755-758.

[13] 李文尧,晏冲为.中心回线瞬变电磁法全期视电阻率的二分法求解[J].昆明理工大学学报(自然科学版),2013,38(2):26-33.

Li Wenyao,Yan Chongwei.Bisection algorithm of central loop time domain transient electromagnetism all time apparent resistivity[J].Journal of Kunming University of Science and Technology(Natural Science Edition),2013,38(2):26-33.

[14] 陈清礼,严良俊,付志红,等.长偏移距瞬变电磁法全区视电阻率的二分搜索数值算法[J].石油地球物理勘探,2009,44(6): 779-782.

Chen Qingli,Yan Liangjun,Fu Zhihong,et al.Searching algorithm for full time apparent in long-offset transient electromagnetic sounding[J].Oil Geophysical Prospecting,2009,44(6):779-782.

[15] 刘 俊,杨海燕,林 云.地下瞬变电磁法全区视电阻率二分搜索算法[J].工程地球物理学报,2012,9(5):526-530.

Liu Jun,Yang Haiyan,Lin Yun.Binary search algorithm of all-time apparent resistivity in underground transient electromagnetic method[J].Chinese Journal of Engineering Geophysics,2012,9(5): 526-530.

[16] 姜志海,岳建华,刘树才.多匝重叠小回线装置的矿井瞬变电磁观测系统[J].煤炭学报,2007,32(11):1152-1156.

Jiang Zhihai,Yue Jianhua,Liu Shucai.Mine transient electromagnetic observation system of small multi-turn coincident configuration [J].Journal of China Coal Society,2007,32(11):1152-1156.

[17] 李金铭.地电场与电法勘探[M].北京:地质出版社,2005.

[18] Fullagar P K.Generation of conductivity-depth pseudo-section from coincident and in-loop TEM data[J].Exploration Geophysics, 1995,26:43-55.

[19] Raiche A P.The use and misuse of apparent resistivity in electromagnetic methods[J].Geophysics,1984,15:37-41.

[20] 牛之琏.时间域电磁法原理[M].长沙:中南大学出版社,2007.

Optimized binary search algorithm of full space transient electromagnetic method all-time apparent resistivity

JIANG Guo-qing1,2,CHENG Jiu-long2,SUN Xiao-yun2,CHEN Ding2,WANG Hai-yan2

(1.Basic Geological Research Institute,Geological Survey of Jiangsu Province,Nanjing 210018,China;2.State Key Laboratory of Coal Resources and Safe Mining,China University of Mining and Technology(Beijing),Beijing 100083,China)

For the purpose of higher precision and efficiency,the optimized binary search algorithm was used for the calculation of all-time apparent resistivity.Based on the transient response of loop center in homogeneous full space, the kernel function curve feature of full space was discussed and the expression of late period apparent resistivity was derived.Method for calculating Bzfrom∂Bz/∂T by numerical integration in full space was given,and the transition precision was analyzed through simulation.The function monotonicity of Bzand∂Bz/∂T according to resistivity ρ at certain transient delay T was studied,and the results indicate that Bzis a single value function of ρ and it constantly satisfies the condition of binary search for ρ,but∂Bz/∂T is a double value function of ρ and the precondition for binary search is instrumental variable Z constantly meeting the formula Z≤1.225 during the value interval of variable ρ,on the basis of systematic study both the initial value interval and detailed flow chart of binary search algorithm were optimized.Model calculation and applied case indicate that the optimized binary search algorithm proposed in this paper has the characteristics of fast calculation speed,high calculation precision and can obtain all-time apparent resistivity ρalland ρ′allat the same time respectively from Bzand∂Bz/∂T,it has a better practicability and better detection effect in fullspace transient electromagnetic method.

full space transient electromagnetic method;central loop unit;all-time apparent resistivity;binary search algorithm;optimized

P631.3

A

0253-9993(2014)12-2482-07

2014-04-11 责任编辑:韩晋平

国家自然科学基金资助项目(51034003,51174210);高等学校博士学科点专项科研基金资助项目(20120023110014)

姜国庆(1986—),男,山东济宁人,工程师,硕士研究生。E-mail:jiangguoqing220@126.com。通讯作者:程久龙(1965—),男,安徽安庆人,教授,博士生导师。E-mail:jlcheng@126.com

姜国庆,程久龙,孙晓云,等.全空间瞬变电磁全区视电阻率优化二分搜索算法[J].煤炭学报,2014,39(12):2482-2488.

10.13225/j.cnki.jccs.2014.0477

Jiang Guoqing,Cheng Jiulong,Sun Xiaoyun,et al.Optimized binary search algorithm of full space transient electromagnetic method all-time apparent resistivity[J].Journal of China Coal Society,2014,39(12):2482-2488.doi:10.13225/j.cnki.jccs.2014.0477