本发明涉及油田开发技术领域，特别是涉及到一种水平井复杂环境随钻电磁波电阻率测井快速仿真方法。
背景技术：
随着水平井和大斜度井钻井技术广泛应用，电缆测井应用受到一定制约，测井数据采集广泛采用随钻测量方式。随钻电阻率由于钻开时间短，测量点受泥浆侵入影响较小，且能够提供不同探测深度的测井曲线。随钻电阻率数据能反映地层真实的参数，因此随钻测井仪器有着广阔的应用前景。随着复杂油气田勘探开发的不断深入，大斜度井、水平井等复杂工艺井的广泛应用，随钻测井技术研究与随钻测井仪器研发得到足够的重视，并得到了快速的发展。
随钻电磁波电阻率类仪器电阻率测量与电缆侧向类仪器有很大不同，侧向类仪器利用直流电模式进行测量，通过一定的仪器常数进行刻度，即可把测量采集的信号转换成地层电阻率信号；随钻电磁波电阻率仪器是利用距发射源不同距离的两个接收线圈获得感应电磁波的相位移和幅度衰减变化来描述地层电阻率。相位移和幅度衰减与地层电阻率是一种非线性关系，没有固定的刻度系数，通常利用相位移-电阻率转换关系将相位移转换成相位移电阻率，利用幅度衰减-电阻率转换关系将幅度衰减转换成幅度衰减电阻率。不同仪器、不同源距、不同工作频率所对应的电阻率转换关系不同。
传统的随钻电磁波测井仪器发射线圈与接收线圈共轴，测量得到的地层信号为地层信息的平均值；近年来三大石油测井服务公司相继推出了具有方位探测能力的随钻方位电磁波测井仪器。2005年斯伦贝谢公司推出了PeriScope方位电阻率测量仪，2006年贝克休斯公司推出了随钻方位电磁波测井仪APR，2007年哈里伯顿公司推出了方位深电阻率测量仪ADR。方位电磁波测量仪器均采用轴向，倾斜或横向线圈混合，能够更好地提供有关地层方位的信息，指示地层的各向异性，识别地层边界。
目前，随钻电磁波电阻率仪器正演仿真主要有有限元方法、模式匹配法、伪解析解方法、积分方程方法等，以上的所有方法在二维模型中应用较为成熟，而对于复杂的三维地层模型，伪解析解、积分方程法、模式匹配法等已经不再适用，而有限元方法虽然可以模拟复杂三维结构，但是面临计算速度慢、计算要求计算机内存空间大等问题，严重制约了三维数值仿真的发展和应用。
为此，我们从随钻电磁波电阻率测量信号源头出发，将不存在解析解的复杂的三维地层模型按矢量叠加原理进行拆分，获得纵向成层、径向成层两套简单结构模型，纵向成层、径向成层地层模型电磁波数值仿真存在解析解，利用解析解分别计算，并将纵向成层模型仿真计算结果叠加到径向成层模型中，最终对径向成层模型数值仿真耦合了三维地层模型全参数，对最终径向成层介质仿真结果进行环境校正处理，得到最终的数值仿真结果。本方法解决了传统随钻测井有限元三维仿真计算速度慢、计算量大从而导致在工程应用中受限，同时也解决了三维复杂地层环境中随钻电磁波仿真不存在解析解，导致解析解无法应用在三维空间中随钻电磁仿真的问题。本技术将传统有限元仿真一个计算单元耗时大于30分钟，提高到仿真一个计算单元耗时小于10秒，大大提高了计算速度，解决了三维仿真工程应用技术瓶颈问题。
技术实现要素：
本发明的目的是提供一种水平井复杂环境随钻电磁波电阻率测井快速仿真方法。
本发明的目的可通过如下技术措施来实现：
步骤1，依照随钻电磁波电阻率仪器工作原理，设定工作频率和收发线圈源距，将发射线圈等价为点源，根据单位磁偶极子发射源随时间变化关系，建立单位磁偶极子电磁场分量关系式；
步骤2，将水平井地层模型等价拆分包含井斜、层界面、原状地层水平电阻率、各向异性系数、地层介电常数、边界距参数的纵向成层介质模型和包含仪器结构、井眼直径、泥浆电阻率、井眼偏心距、偏心方位、侵入深度、侵入带电阻率、原状地层电阻率的径向成层介质模型；
步骤3，将步骤2纵向成层介质模型利用解析方法进行包括数值模拟、信号补偿、电阻率转换的操作，获得仪器纵向成层介质中的测井响应；
步骤4，将步骤3中的测井响应，替代步骤2中径向成层介质原状地层电阻率，调用径向成层解析解数值仿真方法，计算测井响应，对计算得到的信号进行包括对称补偿、电阻率工程转换的操作获得电阻率曲线；
步骤5，对计算后得到的电阻率曲线进行包括井眼环境校正、偏心校正的环境因素校正操作，最后获得输出电阻率仿真结果。
在步骤1中包含：
步骤11：单位磁偶极子源随时间的变化关系为eiωt，其中ω为角频率，得到方向单位磁偶极子在均匀各向异性介质中产生的Hertz势；
步骤12：将所述Hertz势经过转化处理后表示为Sommerfeld积分形式；
步骤13：由电磁场与Hertz势之间的关系式可以得到沿三个方向单位磁偶极子产生的电场和磁场各分量的解析式，并将解析式表示为Sommerfeld积分形式，其中电场和磁场z分量的Sommerfeld积分形式分别表示成波模积分的形式；
步骤14：根据Maxwell方程组，电场和磁场波模的切向分量与纵向分量之间的关系以及所述电场和磁场z分量的Sommerfeld积分形式，分别计算得到电场和磁场x、y分量的Sommerfeld积分表达式。
在步骤2中包含：
步骤21：按照三维空间地层电性分布和结构参数，充分考虑层、侵入、倾斜、井眼变化，建立三维原始模型；
步骤22：按照原始模型的层数、层界面位置、层内水平电阻率、垂直电阻率、地层倾角、地层方位角、井斜角、采样深度、采样间隔建立新纵向成层介质模型，纵向成层介质模型不考虑模型参数径向上变化；
步骤23：按照原始模型在径向上变化，逐采样点建立径向成层介质模型，模型参数与原始模型在对应采样点处井径、泥浆电阻率、侵入深度、侵入带电阻率、仪器偏心距一致，原状地层电阻率不予赋值。
步骤3包括：
步骤31：纵向共有n+1层介质中，各层层界面位置、层厚和地层参数不尽相同，各层编号为k＝0,1,…,n，磁偶极子源在第j层，方向单位磁偶极子产生TE波，方向单位磁偶极子既产生TM波又产生TE波，方向单位磁偶极子也是既产生TM波又产生TE波。
步骤32：根据所述电场和磁场z分量的Sommerfeld积分形式以及层界面处电场和磁场的切向量连续，得到z分量的电场和磁场强度表达式。
步骤33：根据所述z分量的电场和磁场强度表达式以及Maxwell方程组电场和磁场波模的切向分量与纵向分量之间的关系，得到方向单位磁偶极子在x，y分量上的电场和磁场强度表达式。
步骤34：对推导得到的电场和磁场强度表达式进行积分计算，即可得到各个方向磁偶极子x，y，z分量的电场和磁场值。
步骤35：根据接收线圈处的电场强度或磁场强度计算得到接收线圈处的感应电动势，分别计算出两个接收线圈的相位和幅度。
步骤36：用两个接收线圈得到的相位做差，幅度做比值，得到相位移和幅度衰减。
步骤4包括：
步骤41：径向成层介质模型采用逐源距和逐个频率仿真，利用步骤3中纵向成层介质正演结果对应源距和频率下的电阻率值，赋值径向成层模型原状地层电阻率，将径向成层介质模型补充完整；
步骤42：根据环形电流源随时间的变化关系为exp(iωt)，其中ω为角频率，有ω＝2πf，f为交变电流频率；并假设在地层直角坐标系中，水平面为xy面，源点的位置坐标为rt＝(xt,yt,zt)，场点的位置坐标为r＝(x,y,z)，则方向单位磁偶极子在均匀各向异性介质中产生的Hertz势可表示为：
式中，μb为均匀介质磁导率，σhb为均匀各向异性介质的水平复电导率；
步骤43：建立电磁场与Hertz势之间的关系式：
其中E为电场强度，Π为Hertz势，i为虚部单位，为均匀各向异性介质的电导率张量；
步骤44：根据Maxwell方程组，建立电场和磁场波模的切向分量与纵向分量之间方程组：
式中表示磁场波模切向分量，表示磁场波模纵向分量，表示电场切向分量，表示电场纵向分量，μ为磁导率，εh为介电常数，λ为积分变量，为各向异性系数，σvb为均匀各向异性介质的垂向复电导率；
步骤45：求解上式得到线圈处的电场强度或磁场强度计算得到接收线圈处的感应电动势，进而计算得到两个接收线圈之间的相位差和幅度比值。
步骤46：将相位差值和幅度比值通过电阻率转换链表转换成电阻率值。
本发明有益效果：随钻电磁波电阻率仪器结构复杂，三维正演仿真难度大，在水平井复杂地层环境中，现有随钻电磁波电阻率正演仿真要么不适用，要么计算速度慢，制约数值仿真在仪器开发领域的应用效率，同时更严重制约正演仿真在实时地质导向、井场资料处理等领域的应用。本发明可有效提高正演仿真计算效率，将三维地质模型按照矢量合成方式拆分，通过纵向-径向数值仿真算法联合，解决了三维复杂地层环境中随钻电磁波仿真不存在解析解，导致解析法无法应用在三维空间中随钻电磁仿真的问题，且开发出新的水平井复杂环境下三维快速正演仿真算法，将传统有限元仿真一个计算单元耗时大于30分钟，提高到仿真一个计算单元耗时小于10秒，大大提高了计算效率，解决传统三维数值仿真在井场资料处理应用中效率低、速度慢的技术瓶颈问题。
附图说明
图1为本发明的一种水平井复杂环境随钻电磁波电阻率测井快速仿真方法具体实施例的流程图；
图2为纵向成层介质地层模型；
图3为本发明的幅度比-电阻率转换链表；
图4为本发明的相位差-电阻率转换链表；
图5为本发明的复杂模型正演结果示例。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
本发明的水平井复杂环境随钻电磁波电阻率测井快速仿真方法，包括：
步骤1，依照随钻电磁波电阻率仪器工作原理，设定工作频率和收发线圈源距，将发射线圈等价为点源，根据单位磁偶极子发射源随时间变化关系，建立单位磁偶极子电磁场分量关系式，包含几个步骤：
1)设单位磁偶极子源随时间的变化关系为exp(iωt)，其中ω为角频率，并假设在地层直角坐标系中(水平面为xy面)源点的位置坐标为rt＝(xt,yt,zt)，场点的位置坐标为r＝(x,y,z)，则方向单位磁偶极子在均匀各向异性介质中产生的Hertz势可表示为：
式中，μb为均匀介质磁导率，σhb为均匀各向异性介质的水平复电导率。
2)经过转化处理，(11)式可以表示为如下Sommerfeld积分形式：
式中，Jv为v阶Bessel函数，λ为积分变量。
方向单位磁偶极子在均匀各向异性介质中产生的Hertz势可以表示为：
式中，为各向异性系数，i为复数单位，σhb、σvb分别为均匀各向异性介质的水平、垂向复电导率。经推导(13)、(14)式可以分别表示为如下Sommerfeld积分形式：
式中，
3)由电磁场与Hertz势之间的关系式
其中为均匀各向异性介质的电导率张量。将(2)、(5)、(6)式带入(7)式，可以得到沿三个方向单位磁偶极子产生的电场和磁场各分量的解析式及Sommerfeld积分形式，其中电场和磁场z分量的Sommerfeld积分形式分别表示为：
上述各分量均被表示成了波模积分的形式。例如
式中为某一λ对应的波模。
4)根据Maxwell方程组，电场和磁场波模的切向分量与纵向分量之间的关系可表示为：
步骤2，如图2所示，将水平井地层模型等价拆分包含井斜、层界面、原状地层水平电阻率、各向异性系数、地层介电常数、边界距参数的纵向成层介质模型和包含仪器结构、井眼直径、泥浆电阻率、井眼偏心距、偏心方位、侵入深度、侵入带电阻率、原状地层电阻率的径向成层介质模型，包含几个步骤：
1)按照三维空间地层电性分布和结构参数，充分考虑层、侵入、倾斜、井眼等变化，建立三维初始模型；
2)按照原始三维模型的层数、层界面位置、层内水平电阻率、垂直电阻率、地层倾角、地层方位角、井斜角、采样深度、采样间隔建立新纵向成层介质模型，纵向成层介质模型不考虑模型参数径向上变化；
3)按照原始模型在径向上变化，逐采样点建立径向成层介质模型，模型参数与原始模型在对应采样点处井径、泥浆电阻率、侵入深度、侵入带电阻率、仪器偏心距一致，原状地层电阻率未赋值。
步骤3，将步骤2纵向成层介质模型利用解析解方法进行数值模拟、信号补偿、电阻率转换等操作，获得仪器纵向成层介质中的测井响应，包括几个步骤：
1)设纵向成层各向异性介质共有n+1层，各层编号为k＝0,1,…,n，磁偶极子源在第j层，各层层界面位置、层厚和地层参数如图2所示。方向单位磁偶极只产生TE波，方向单位磁偶极子既产生TM波又产生TE波，方向单位磁偶极子也是既产生TM波又产生TE波。
2)根据式(18a)-(18f)可以得到z分量的电场和磁场强度表达式为：
上式中l＝0,1,…,n，μl表示第l层磁导率，σhl表示第l层水平电导率，σvl表示第l层垂直电导率，当l＝j时δlj＝1，否则δlj＝0；Al(λ)，Bl(λ)，Cl(λ)，Dl(λ)，El(λ)，Fl(λ)为待定系数，由层界面处，电场和磁场的切向量连续得到。
3)将式(20a)-(23d)中对应的变量代入式(19a)-(19d)中，即可以计算得到方向单位磁偶极子在x，y分量上的电场和磁场强度。
4)对推导得到的电场和磁场强度表达式进行积分计算，即可得到各个方向磁偶极子x，y，z分量的电场和磁场值。
5)然后根据接收线圈处的电场强度或磁场强度计算得到接收线圈处的感应电动势，利用距发射线圈不同距离接受线圈感应电动势，可以计算出两接收线圈相位(式23)和幅度(式24)：
式中PS表示相位，AT表示幅度，βi为方位角，为方位角为βi视接收线圈电压。
6)用两个接收线圈得到的相位做差，两个接收线圈得到的幅度做比值，即可得到相位移和幅度衰减。
步骤4，将步骤3中的测井响应，替代步骤2中径向成层介质原状地层电阻率，调用径向成层解析解数值仿真方法，计算测井响应，对计算得到的信号进行对称补偿、电阻率工程转换等获得电阻率曲线，包括以下步骤：
1)径向成层介质模型采用逐源距和频率仿真，利用步骤3中纵向成层介质正演结果对应源距和频率下的电阻率值，赋值径向成层模型原状地层电阻率，将径向成层介质模型补充完整；
2)根据环形电流源随时间的变化关系为exp(iωt)，其中ω为角频率，有ω＝2πf，f为交变电流频率；并假设在地层直角坐标系中，水平面为xy面，源点的位置坐标为rt＝(xt,yt,zt)，场点的位置坐标为r＝(x,y,z)，则方向单位磁偶极子在均匀各向异性介质中产生的Hertz势可表示为：
式中，μb为均匀介质磁导率，σhb为均匀各向异性介质的水平复电导率；
3)建立电磁场与Hertz势之间的关系式：
其中E为电场强度，Π为Hertz势，i为虚部单位，为均匀各向异性介质的电导率张量；
4)根据Maxwell方程组，建立电场和磁场波模的切向分量与纵向分量之间方程组：
式中表示磁场波模切向分量，表示磁场波模纵向分量，表示电场切向分量，表示电场纵向分量，μ为磁导率，εh为介电常数，λ为积分变量，为各向异性系数，σvb为均匀各向异性介质的垂向复电导率；
5)求解上式线圈处的电场强度或磁场强度，计算得到接收线圈处的感应电动势，进而计算得到两个接收线圈之间的相位差和幅度比值。
6)将相位差值和幅度比值通过电阻率转换链表如图3和图4所示，转换成电阻率值。
步骤5，对计算后得到的电阻率曲线进行井眼环境校正、偏心校正等环境因素校正，最后获得输出电阻率仿真结果，如图5所示，图中第一道为电阻率剖面道，第二道为侵入深度反演结果道，第三道为高频相位差电阻率道，第四道为高频幅度比电阻率道。
本发明可以解决现有数值仿真方法在复杂地层环境下随钻电磁波电阻率仪器数值仿真计算速度慢的问题，大大提高计算效率，节约计算时间。
以上，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本方法领域的方法人员在本发明揭露的方法范围内，根据本发明的方法方案及其发明构思加以等同替换或改变，都应涵盖在本发明的保护范围之内。