范宜仁, 李炜,3, 李虎, 王磊, 李栗

(1.中国石油大学(华东)地球科学与技术学院, 山东 青岛 266580; 2.中国石油大学CNPC测井重点实验室, 山东 青岛 266580; 3.中海油田服务股份有限公司物探事业部, 天津 300451)

0 引 言

随钻电磁波测井具有实时性、探测范围大等特征,被广泛应用于随钻地质导向和水平井储层评价等方向[1]。Tsili Wang等使用时域有限差分(FDTD)模拟了倾斜各向同性地层及仪器偏心条件下的测井响应[2],Hwa Ok Lee等使用柱坐标系下的FDTD和PML吸收边界模拟了各向异性地层的随钻电磁响应[3]。董建党采用FDTD模拟了三维各向异性地层条件下随钻电磁波测井响应[4],其地层模型较为简单。范宜仁等研究了仪器偏心条件下三维倾斜地层的测井响应[5]。使用FDTD方法进行传统电测井数值模拟时,网格剖分往往采用阶梯近似方法,具有易于实现、计算速度快等特点,但在处理三维复杂地层边界、仪器偏心、线圈、井眼等情况时,为保证计算精度,需不断加密网格,计算速度与存储量急剧增加。

柴玫、耿轩等研究了二维时域有限差分亚网格技术,并将其应用于带涂层的局部机翼、有刻槽的金属长方柱等模拟[6-7],大大节约了计算内存和时间。张玉、曾碧能等采用共形网格技术或改进的共形网格技术对PBG仿真或涂层三角面元的散射等进行了模拟,极大地提高了模拟精度[8-9]。但中国国内对于将亚网格和共形网格技术应用于随钻电磁波测井响应的数值模拟研究很少。本文把时域有限差分亚网格和共形网格技术应用于复杂地层随钻电磁波测井响应的数值模拟,并与单纯使用时域有限差分算法的模拟结果进行对比,提高了其计算速度和模拟精度。

1 三维圆柱坐标系下FDTD差分格式

模拟的电性各向同性地层为有耗介质,其本构方程为

D=ε·E,J=σ·E

(1)

地层磁性为

B=μ·H,Jm=σm·H

(2)

此时,在无源区Maxwell方程微分形式为

(3)

(4)

式中,ε为地层介电常数;μ为地层磁导率;σ为地层电导率;σm为地层导磁率。

将圆柱坐标系Yee元胞中的E、H的某一分量设为f(r,φ,z,t),表示为

f(r,φ,z,t)=f(iΔr,jΔφ,kΔz,nΔt)=fn(i,j,k)

(5)

式中,Δr、Δφ、Δz分别为空间步长;Δt为时间步长。采用二阶中心差分近似进行离散

(6)

按照Yee元胞可以得出FDTD迭代格式(以Er分量为例)

(7)

2 亚网格技术

FDTD作为模拟随钻电磁波传播的方法其模拟精度主要取决于剖分网格的大小。网格尺寸δ的取值范围必须满足色散条件,即δ≤λ/10,其中λ为FDTD模拟区域内介质尺寸最小的波长。δ的大小还需考虑目标区域的某些细微结构,比如随钻电磁波测井仪器的线圈以及井眼,这时需要对这些细微区域采用较小的网格,否则会影响模拟结果的准确性。对于这种情况,可以在线圈、井眼等细微结构处采用亚网格技术,其他FDTD模拟区域采用粗网格技术。这样,可以实现在保证计算精度的前提下大量节省计算机的运行内存和时间。

2.1 亚网格算法

以轴对称圆柱坐标系中r—φ平面为例,研究三维轴对称坐标系时域有限差分(FDTD)亚网格技术。采用粗网格和细网格同时计算的方式,将粗网格和细网格在时间步上交替计算。粗细网格分界面在切向电场面上,如图1(a)、(b)所示,粗网格和细网格的尺寸为奇数比,这里设为3∶1,则粗细网格比nf=3。因为粗网格和细网格的比例为奇数,所以对于此r—φ平面布局,粗网格的磁场和细网格的磁场在时间和空间上均重合,不需要采用插值。

图1 柱坐标系下的亚网格示意图

设粗网格的时间步长为Δtf,尺寸大小为Δr和Δφ,细网格的时间步Δtf=Δtc/nf。不论是粗网格还是亚网格,在网格内部,传统FDTD的计算都适用。因此,亚网格算法的重点是粗细网格边界上E和H的计算,采用空间上与时间上的内插插值法进行计算。对于粗网格直接采取FDTD计算,局部细网格采取蛙跳方式进行计算。

(8)

式中

(9)

(10)

式中,E0、E1、E2由上步得到。粗细网格边界处其他与粗网格不重合点的细网格电场值亦可通过式(10)得到。为了提高计算时电场的稳定性,对于粗细网格边界内一层的电场进行修正,采用加权平均的方式。针对图1所示,采用公式

(11)

(12)

图2 亚网格效果图

测井仪器、井眼、发射线圈、接收线圈等处需要精细刻画,以保证电磁波传播模拟更加精准,所以需要在这些地方采用细网格,在地层等其他地方采用粗网格,其中粗网格与细网格尺寸比例为3∶1。具体横截面和径向的网格剖分见图2(a)、(b)。图2(a)中,在r—φ平面上,对线圈、井眼处采用了细网格,地层其他地方采用了粗网格;图2(b)中,在r—z平面上,对发射线圈和接收线圈处,采用了细网格,地层其他地方采用了粗网格。

2.2 计算速度/精度对比与实例

通过在均匀地层条件下模拟单发双收随钻电磁波仪器的测井响应,对比分析了传统FDTD算法和亚网格算法的计算时间和精度。取发射频率为2 MHz,线圈距为[24 30] in*非法定计量单位,1 ft=12 in=0.304 8 m,下同,模拟地层电阻率分别取1、10 Ω·m,2种算法的对比结果见表1、表2。

表1 R=1 Ω·m模拟结果精度与时间

表2 R=10 Ω·m模拟结果精度与时间

由传统FDTD和亚网格技术的计算时间和模拟精度对比可以发现,亚网格技术计算时间仅为传统FDTD的19%,亚网格技术模拟相对误差均控制在0.8%之内。图3(a)和(b)为传统FDTD和亚网格技术计算的发射频率2 MHz、线圈距[24 30] in条件下的相位差、幅度比刻度图版,可以看出两者吻合很好,验证了亚网格技术的正确性。

图3 利用亚网格技术计算刻度图版与FDTD对比

3 共形网格技术

图4 倾斜地层切割Yee元胞示意图

图5 仪器偏心横截面示意图

使用FDTD模拟三维复杂地层随钻电磁波测井响应时,对于复杂地层比如倾斜地层(见图4)或仪器偏心(见图5)等会产生阶梯误差。剖分的网格尺寸越小,模拟过程中产生的阶梯误差越小。但是减少网格尺寸势必造成网格总数的增加和时间步长的减少,导致计算机计算时间和存储空间的急剧增加。对于复杂地层,这种减少网格尺寸的做法并不可取。这一难题可以用共形网格技术予以解决。

3.1 共形网格技术原理

图6 r—z平面共形网格中地层参数及其等效

设地层1的电磁参数为ε1、σ1、μ1、σm1;地层2的电磁参数为ε2、σ2、μ2、σm2。电场分量位于元胞棱边的中点,根据Maxwell方程积分形式,电磁参数等效值可以由加权平均得到,因此,图6中A、C处4个电场分量中的等效ε和等效σ为

(13)

(14)

同理,F处的磁场分量可由不同地层所占的元胞面积加权平均得到,F处磁场分量中的等效μ和等效σm为

(15)

共形网格中处的FDTD格式与传统FDTD相同,只需将等效电磁参数替代普通地层电磁参数,代入FDTD差分格式即可。共形网格技术的这种简单近似,应用到介质突变分界面等情况,可以很好地减少阶梯近似产生的误差,提高模拟的精度,且不会增加计算时间。图7(a)、(b)为3层45 °倾斜地层共形网格和阶梯近似对比图,可以看出倾斜地层的层界面处的电阻率为上下地层电阻率的加权平均。

图7 45 °的3层倾斜地层电阻率

3.2 计算速度/精度对比与实例

模拟单发双收随钻电磁波仪器偏心时的测井响应,对比分析了传统FDTD阶梯近似和共形网格算法的计算时间和精度。取发射频率为2 MHz,线圈距为[24 30] in,模拟地层电阻率为0.1 Ω·m,并与阶梯近似模拟结果、Hue一维半数值解进行对比(见表3,以一维半数值解为参照)。图8为共形网格、阶梯近似、一维半数值解3种方法模拟的仪器偏心条件下随钻电磁波测井仪器响应特征对比。

图8 利用亚网格技术计算刻度图版与FDTD对比

由表3和图8可知,仪器偏心条件下随钻电磁波测井响应数值模拟,共形网格技术在保证计算时间的前提下,减少了相对误差,精确度较高。

表3 仪器偏心条件下模拟结果精度与时间

4 结 论

(1) 研究了三维柱坐标系下亚网格技术与共形网格技术,将其应用于复杂地层条件随钻电磁波测井响应的数值模拟。亚网格技术解决了对井眼、线圈等精细结构网格剖分过小、导致CPU计算时间过多的难题;精细结构处用细网格剖分,其他结构采用粗网格,其计算时间远远小于传统FDTD,且相对误差控制在0.8%以内。

(2) 共形网格技术解决了复杂地层边界或仪器偏心等情况下采用阶梯近似产生较大误差的难题,在保证计算时间的条件下大大减少了阶梯误差,提高了精度。

(3) 这2种技术可用于模拟三维复杂地层随钻电磁波测井响应,为大斜度井/水平井评价和随钻快速反演提供理论指导和技术支持。

参考文献:

[1] 杨震. 非均匀复杂地层随钻电磁波测井响应研究 [D]. 青岛: 中国石油大学(华东), 2009.

[2] Wang T, Signorelli J. Finite-difference Modeling of Electromagnetic Tool Response for Logging While Drilling

[J]. Geophysics, 2004, 69(1): 152-160.

[3] Hue Y K, Teixeiral F L, San L E, et al. Modeling of EM Logging Tools in Arbitrary 3-D Borehole Geometries Using PML FDTD [J]. IEEE Geosci Remote Sens Lett, 2005, 2(1): 78-81.

[4] 董建党. 各向异性地层电磁测井响应的时域有限差分方法研究 [D]. 青岛: 中国石油大学(华东), 2008.

[5] 范宜仁, 胡云云, 李虎, 等. 随钻电磁波测井仪器偏心条件下响应模拟与分析 [J]. 中国石油大学学报: 自然科学版, 2014, 38(2): 59-66.

[6] 柴玫, 阎玉波, 葛德彪. 一种可穿越介质的二维时域有限差分亚网格技术 [J]. 西安电子科技大学学报: 自然科学版, 2000, 27(5): 581-584.

[7] 耿轩, 徐金平. R-FDTD与亚网格相结合技术及其应用 [J]. 东南大学学报: 自然科学版, 2004, 34(2): 161-165.

[8] 张玉, 宋健, 梁昌洪. 并行共形FDTD算法及其在PBG结构仿真中的应用 [J]. 电子学报, 2003, 12(A): 2142-2144.

[9] 曾碧能. 基于FDTD算法的非均匀网格及共形网格技术的实现 [D]. 成都: 电子科技大学, 2007.