全文快速搜索:   高级搜索

  中国石油大学学报(自然科学版)  2018, Vol. 42 Issue (4): 41-49  DOI:10.3969/j.issn.1673-5005.2018.04.005
0

引用本文 [复制中英文]

袁超, 李潮流, 周灿灿, 等. 基于空间响应分布函数的水平井补偿密度测井快速正演模拟[J]. 中国石油大学学报(自然科学版), 2018, 42(4): 41-49. DOI: 10.3969/j.issn.1673-5005.2018.04.005.
[复制中文]
YUAN Chao, LI Chaoliu, ZHOU Cancan, et al. Fast forward simulation of compensated density logging in horizontal wells based on spatial response distribution function[J]. Journal of China University of Petroleum (Edition of Natural Science), 2018, 42(4): 41-49. DOI: 10.3969/j.issn.1673-5005.2018.04.005.
[复制英文]

基金项目

中国石油天然气集团公司科学研究与技术开发项目(2016A-3601)

作者简介

袁超(1987-), 男, 博士, 研究方向为放射性测井方法、测井资料处理与综合应用。E-mail:yuanc@petrochina.com.cn

文章历史

收稿日期:2017-10-01
基于空间响应分布函数的水平井补偿密度测井快速正演模拟
袁超1 , 李潮流1 , 周灿灿1 , 葛云龙2     
1. 中国石油勘探开发研究院, 北京 100083;
2. 中国石油大学(北京)地球物理与信息工程学院, 北京 102249
摘要: 水平井中补偿密度测井的正演模拟主要采用蒙特卡罗方法,计算时间较长,难以满足生产实时评价的需求。提出种补偿密度测井快速正演模拟的新方法,引入空间响应分布函数表征伽马场分布及地层不同体积元对探测器计数贡献。采用蒙特卡罗程序建立快速正演模拟所需数据库,利用空间响应分布函数计算不同地层区域对探测器计数贡献,实现补偿密度测井的快速正演模拟。数值模拟结果表明,新方法的计算精度与蒙特卡罗方法致,但计算速度提高十万倍以上,能够满足水平井井眼轨迹与地层关系确定中快速正演的需求。
关键词: 补偿密度测井    井眼轨迹与地层关系    快速正演模拟    空间响应分布函数    蒙特卡罗方法    
Fast forward simulation of compensated density logging in horizontal wells based on spatial response distribution function
YUAN Chao1 , LI Chaoliu1 , ZHOU Cancan1 , GE Yunlong2     
1. Research Institute of Petroleum Exploration and Development, PetroChina, Beijing 100083, China;
2. College of Geophysics and Information Ergineering in China University of Petroleum(Beijing), Beijing 102249, China
Abstract: Monte Carlo method is commonly used in forward simulations of compensated density logging in horizontal wells, but is computationally time-consuming and cannot meet the needs for real-time production evaluation. In this paper, a fast forward simulation method of compensated density logging is developed, and a spatial response distribution function is introduced to describe the contribution of gamma ray field distribution and different formation volume elements to the detector counts. Monte Carlo simulation is then employed to establish the database of fast forward simulation in advance. The spatial response distribution function is used to calculate the contribution of different formation regions to the detector counts, and the fast forward simulation of density logging is realized. The numerical simulation results show that the new method can provide results with the same accuracy as regular Monte Carlo simulation method, but the computation speed is increased by over two orders of magnitude, which can meet the needs for determining the relation between well trajectory and drilled formation in the horizontal well.
Keywords: compensated density logging    relation between well trajectory and formation    fast forward simulation    spatial response distribution function    Monte Carlo method    

随着钻井技术的发展, 为了提高石油天然气采收率, 水平井的数量逐年增加[1-2]。井眼轨迹与地层关系是水平井测井解释的关键和基础[3], 依据正演模拟曲线与实测曲线对比确定两者关系时需要根据初始地层模型快速正演测井曲线。如果正演曲线与实测曲线吻合, 则说明设定地层模型符合真实地层情况, 否则需要重新调整地层模型[4]。目前, 主要利用电阻率、自然伽马和密度测井曲线确定井眼轨迹与地层关系, 国内外大量学者开展电阻率和自然伽马测井快速正演模拟的研究[5-9]。国外学者对补偿密度测井快速正演做了相关研究, Booth等[10]研究了正演模拟中的重要性估算; Watson[11]利用蒙特卡罗模拟探测器微分灵敏度函数; Case等[12]研究了放射性测井建模中的灵敏度函数; Liu等[13]基于精细网格研究了探测器响应的重要性; Mendoza等[14]利用通量灵敏度函数和格林函数快速模拟核测井响应。虽然国外学者对补偿密度测井快速正演模拟做了一系列的研究, 但是国内学者主要采用蒙特卡罗研究测井响应特性[15-16], 未见关于补偿密度测井快速正演模拟的公开发表文献。蒙特卡罗模拟计算耗时长, 不满足井眼轨迹与地层关系确定中快速正演模拟的需求[17]。笔者提出一种补偿密度测井快速正演模拟方法, 利用蒙特卡罗程序建立快速正演模拟所需数据库, 计算不同地层区域对长、短源距探测器计数的贡献近似值, 进而得到长、短源距探测器计数, 并将探测器计数转换为地层密度, 以快速获取高精度的补偿密度测井正演模拟结果。

1 补偿密度测井快速正演模拟原理

伽马源发射出的伽马射线进入地层后, 在地层中的输送过程可用玻尔兹曼方程来表示, 伽马射线在地层中输送过程示意图如图 1所示。利用玻尔兹曼方程可以得出探测器位置处伽马射线计数为

$ N\left( {{r_{\rm{D}}}} \right) = \int {{\rm{d}}r} \int {{\rm{d}}E} \int {{\rm{d}}\mathit{\Omega }\psi \left( {{r_{\rm{S}}},r,E,\mathit{\Omega }} \right)R\left( {{r_{\rm{D}}},r,E,\mathit{\Omega }} \right)} . $ (1)
图 1 补偿密度测井中伽马射线输送过程示意图 Fig.1 Schematic diagram of gamma ray transportation in compensated density logging

式中, r为伽马射线位置; E为伽马射线能量; Ω为入射伽马射线方向; ψ为伽马射线角通量; rS为伽马放射源位置; rD为探测器位置; R为探测器响应函数。

从公式(1)可以看出, 补偿密度测井仪器在探测器位置处的伽马射线计数与入射伽马参数(包括入射伽马能量和方向)、伽马射线角通量、发生康普顿散射反应的伽马射线位置以及探测器响应函数(与探测器位置、发生康普顿散射反应的伽马射线位置以及入射伽马射线参数有关)。在补偿密度测井中, 入射伽马能量和方向等参数已知, 因此在探测器位置处伽马射线计数与伽马射线角通量、发生康普顿散射的伽马射线位置以及探测器响应函数有关。因此探测器位置处伽马射线计数受3个因素控制, 即伽马射线在地层中的空间通量分布、地层中不同位置处发生康普顿散射反应的伽马射线对探测器计数的贡献以及探测器响应函数有关[18]

定义伽马放射源发射出的伽马射线对探测器计数贡献为空间响应分布函数FS:

$ {F_{\rm{S}}} = \mathit{\Phi }\eta . $ (2)

式中, FS为空间响应分布函数, 无量纲; Φ为伽马射线的空间通量分布; η为伽马射线与空间不同位置体积元发生康普顿散射后对探测器计数的贡献。

空间响应分布函数FS的物理意义为地层中不同位置体积元对补偿密度测井伽马探测器计数的贡献。利用蒙特卡罗方法可以模拟单一密度地层中的长、短源距探测器处的空间响应分布函数FS以及长、短源距探测器的伽马计数。

对于如图 2所示的两种密度组合地层模型中, 可认为探测器计数是区域Z1Z2对探测器计数贡献的线性叠加。快速正演模拟中, 选取密度接近组合地层密度平均值的地层作为参考地层, 在不同密度地层区域Z1Z2, 对参考地层的空间响应分布函数积分得到每个区域对长、短源距探测器计数的贡献近似值, 将不同地层密度区域的贡献近似值归一化后与相应地层单一密度条件下记录的长、短源距探测器计数相乘并累加, 计算得到长、短源距探测器计数, 即

$ {N_{\rm{S}}} \cong \frac{{\int_{{Z_1}} {{F_{{\rm{SRS}}}}} }}{{\int_{{Z_1}} {{F_{{\rm{SRS}}}}} + \int_{{Z_2}} {{F_{{\rm{SRS}}}}} }}{N_{\rho {\rm{S1}}}} + \frac{{\int_{{Z_2}} {{F_{{\rm{SRS}}}}} }}{{\int_{{Z_1}} {{F_{{\rm{SRS}}}}} + \int_{{Z_2}} {{F_{{\rm{SRS}}}}} }}{N_{\rho {\rm{S2}}}}, $ (3)
$ {N_{\rm{L}}} \cong \frac{{\int_{{Z_1}} {{F_{{\rm{SRL}}}}} }}{{\int_{{Z_1}} {{F_{{\rm{SRL}}}}} + \int_{{Z_2}} {{F_{{\rm{SRL}}}}} }}{N_{\rho {\rm{L1}}}} + \frac{{\int_{{Z_2}} {{F_{{\rm{SRL}}}}} }}{{\int_{{Z_1}} {{F_{{\rm{SRL}}}}} + \int_{{Z_2}} {{F_{{\rm{SRL}}}}} }}{N_{\rho {\rm{L2}}}}. $ (4)
图 2 两种密度组合地层模型横截面示意图 Fig.2 Cross section of formation model with two density zones

式中, NSNL分别为快速正演模拟的短源距和长源距探测器计数; FSRSFSRL分别为短源距和长源距探测器位置处利用蒙特卡罗方法模拟计算的参考地层的空间响应分布函数; NρS1NρS2分别为区域Z1Z2的单一密度地层中短源距探测器计数; NρL1NρL2分别为区域Z1Z2的单一密度地层中长源距探测器计数。

公式(3)中$\frac{{\int {_{{Z_1}}{F_{{\rm{SRS}}}}} }}{{\int {_{{Z_1}}{F_{{\rm{SRS}}}}} + \int {_{{Z_2}}{F_{{\rm{SRS}}}}} }}$$\frac{{\int {_{{Z_2}}{F_{{\rm{SRS}}}}} }}{{\int {_{{Z_1}}{F_{{\rm{SRS}}}}} + \int {_{{Z_2}}{F_{{\rm{SRS}}}}} }}$分别为区域Z1Z2对短源距探测器计数的贡献; 公式(4)中$\frac{{\int {_{{Z_1}}{F_{{\rm{SRL}}}}} }}{{\int {_{{Z_1}}{F_{{\rm{SRL}}}}} + \int {_{{Z_2}}{F_{{\rm{SRL}}}}} }}$$\frac{{\int {_{{Z_2}}{F_{{\rm{SRL}}}}} }}{{\int {_{{Z_1}}{F_{{\rm{SRL}}}}} + \int {_{{Z_2}}{F_{{\rm{SRL}}}}} }}$分别为区域Z1Z2对长源距探测器计数的贡献。这类贡献数值也就是区域Z1Z2对探测器计数贡献的线性叠加系数, 类似于电法测井中的几何因子, 令

$ {G_{{\rm{rS}}}} = \frac{{\int_{{Z_1}} {{F_{{\rm{SRS}}}}} }}{{\int_{{Z_1}} {{F_{{\rm{SRS}}}}} + \int_{{Z_2}} {{F_{{\rm{SRS}}}}} }}, $ (5)
$ {G_{{\rm{rL}}}} = \frac{{\int_{{Z_1}} {{F_{{\rm{SRL}}}}} }}{{\int_{{Z_1}} {{F_{{\rm{SRL}}}}} + \int_{{Z_2}} {{F_{{\rm{SRL}}}}} }}. $ (6)

则公式(3)和(4)变为

$ {N_{\rm{S}}} \cong {G_{{\rm{rS}}}}{N_{\rho {\rm{S1}}}} + \left( {1 - {G_{{\rm{rS}}}}} \right){N_{\rho {\rm{S2}}}}, $ (7)
$ {N_{\rm{L}}} \cong {G_{{\rm{rL}}}}{N_{\rho {\rm{L1}}}} + \left( {1 - {G_{{\rm{rL}}}}} \right){N_{\rho {\rm{L2}}}}. $ (8)

模拟计算得到探测器计数后, 通过短源距和长源距探测器计数与地层密度转换关系将探测器计数转换为地层密度。

2 快速正演模拟所需数据库建立

利用蒙特卡罗方法在单一密度地层中建立补偿密度测井快速正演模拟所需要的数据库, 包括单一密度地层与伽马射线发生康普顿散射反应对探测器计数贡献的空间响应分布函数、探测器计数以及探测器计数与地层密度的转换关系。

2.1 蒙特卡罗计算模型

利用蒙特卡罗方法[19-20]建立仪器-井眼-地层的三维柱状模型, 计算模型设定参数如下:

(1) 地层模型高度和半径分别为60和35 cm, 对模拟地层进行网格化以提高计算精度, 每个网格尺寸为0.5 cm×0.5 cm。

(2) 井眼半径为10 cm, 井眼内充满淡水。

(3) 放射源采用Cs137, 发射伽马射线能量为662 keV。

(4) 短源距探测器和长源距探测器的长度分别为5和10 cm, 距离放射源的距离分别为18和42 cm。

(5) 放射源和短源距探测器之间放置钨屏蔽体, 屏蔽体长度为3 cm。

2.2 空间响应分布函数

利用如图 3所示的蒙特卡罗计算模型, 改变地层孔隙度为0、5%、10%、15%、20%、25%、30%、35%和40%, 设定不同的地层密度为2.65、2.567 5、2.485、2.402 5、2.32、2.237 5、2.155、2.072 5和1.99 g/cm3, 模拟伽马射线的空间分布以及伽马射线与空间不同位置体积元发生康普顿散射反应对伽马探测器计数的贡献, 两者乘积为空间响应分布函数。对于未模拟地层密度的空间响应分布函数需要利用插值的方法获取。地层密度设定为2.485 g/cm3时, 短源距和长源距探测器的空间响应分布函数如图 4所示。

图 3 蒙特卡罗计算模型横截面示意图 Fig.3 Cross section of Monte Carlo simulation model
图 4 空间响应分布函数 Fig.4 Spatial response distribution function

图 4可以看出, 伽马探测器计数的贡献主要来自于放射源以及探测器附近地层, 其他地层区域对伽马探测器计数影响较小。另外, 当计算模型发生变化时, 空间响应分布函数也会发生相应变化, 因此需要根据测井仪器类型以及不同井眼和地层条件, 建立补偿密度测井快速正演模拟所需要的空间响应分布函数的数据库。

2.3 地层密度刻度关系

在模拟空间响应分布函数的同时记录短源距和长源距探测器计数, 建立不同地层密度条件下短源距和长源距探测器计数的数据库, 并建立短源距和长源距探测器计数与地层密度的转换关系数据库。由于实际仪器测量结果为补偿密度, 尤其是当测量结果受围岩影响时, 长、短源距密度与实际仪器响应有较大误差, 通过脊肋图的方法获得补偿密度, 通过模拟地层模型与实际仪器响应对比分析, 获得地层真实密度值。在设定的计算模型条件下, 短源距和长源距探测器计数与地层密度的转换关系如图 5所示。

图 5 地层密度刻度关系 Fig.5 Formation density calibration relationship

图 5所示的地层密度转换关系的脊肋图, 在设定的计算模型条件下, 脊线和肋线方程分别为

$ \ln {N_{\rm{L}}} = 3.2615\ln {N_{\rm{S}}} - 2.453, $ (9)
$ \ln {N_{\rm{L}}} = 0.9662\ln {N_{\rm{S}}} + b. $ (10)

把长、短源距探测器计数NLNS代入公式(9)肋线方程, 得到截距b, 与脊线方程(10)联立, 得到交点的长、短源距探测器计数NLCNSC, 密度校正量为

$ \Delta \rho = k\left( {\ln {N_{{\rm{LC}}}} - \ln {N_{\rm{L}}}} \right). $

式中, k为长源距密度斜率。则补偿密度为

$ \rho = {\rho_{\rm{L}}} + \Delta \rho . $
3 补偿密度测井快速正演模拟及误差分析

利用建立的数据库, 开展补偿密度测井快速正演模拟, 并对比快速正演模拟结果与蒙特卡罗模拟结果, 分析快速正演模拟结果的误差。

3.1 直井模型

在如图 2所示的地层模型中, 区域Z1Z2的地层密度分别为2.567和2.155 g/cm3, 井轴与Z轴的夹角(相对地层倾角)为0°, 测井仪器不断上提模拟在直井中的测井过程。选取区域Z2为参考地层, 根据在密度为2.155 g/cm3的单一密度地层中模拟的空间响应分布函数, 在每个深度点近似计算区域Z1Z2对探测器计数的贡献:

$ {C_{{\rm{S1}}}} = \int_{{Z_1}} {{F_{{\rm{SRS}}}}} , $ (11)
$ {C_{{\rm{S2}}}} = \int_{{Z_2}} {{F_{{\rm{SRS}}}}} , $ (12)
$ {C_{{\rm{L1}}}} = \int_{{Z_1}} {{F_{{\rm{SRL}}}}} , $ (13)
$ {C_{{\rm{L2}}}} = \int_{{Z_2}} {{F_{{\rm{SRL}}}}} . $ (14)

式中, CS1CS2分别为地层区域Z1Z2对短源距探测器计数贡献近似值; CL1CL2分别为地层区域Z1Z2对长源距探测器计数贡献近似值; FSRSFSRL分别为数据库中密度为2.155 g/cm3的参考地层与伽马射线发生康普顿散射反应对短源距和长源距探测器计数贡献的空间响应分布函数。

根据上述计算的地层区域Z1Z2对短源距和长源距探测器计数贡献的近似值, 以及数据库中在地层区域Z1Z2的单一密度地层中记录的探测器计数, 可得到在该地层模型中短源距和长源距探测器计数:

$ {N_{\rm{S}}} = \frac{{{C_{{\rm{S1}}}}}}{{{C_{{\rm{S1}}}} + {C_{{\rm{S2}}}}}}{N_{{\rm{S1}}}} + \frac{{{C_{{\rm{S2}}}}}}{{{C_{{\rm{S1}}}} + {C_{{\rm{S2}}}}}}{N_{{\rm{S2}}}}, $ (15)
$ {N_{\rm{L}}} = \frac{{{C_{{\rm{L1}}}}}}{{{C_{{\rm{L1}}}} + {C_{{\rm{L2}}}}}}{N_{{\rm{L1}}}} + \frac{{{C_{{\rm{L2}}}}}}{{{C_{{\rm{L1}}}} + {C_{{\rm{L2}}}}}}{N_{{\rm{L2}}}}. $ (16)

式中, NSNL分别为快速正演模拟计算的短源距和长源距探测器计数; NS1NS2分别为数据库中区域Z1Z2的单一密度地层中记录的短源距探测器计数; NL1NL2分别为数据库中在区域Z1Z2的单一密度地层中记录的长源距探测器计数。

在相同地层模型条件下, 选取10个位置点, 利用蒙特卡罗方法模拟短源距和长源距探测器计数。根据如图 5建立的短源距和长源距探测器计数与地层密度的转换关系, 将快速正演模拟和蒙特卡罗模拟的探测器计数转换为地层密度, 结果如图 6所示。

图 6 直井中快速正演模拟与蒙特卡罗模拟结果对比 Fig.6 Comparison between fast forward simulation and Monte Carlo simulation results in vertical wells
3.2 斜井模型

在如图 2所示的地层模型中, 区域Z1Z2的地层密度分别为2.567和2.155 g/cm3, 井轴与Z轴的夹角(相对地层倾角)为70°, 测井仪器不断上提模拟在斜井中的测井过程。利用与在直井中快速正演模拟的方法, 模拟计算每个深度点的短源距和长源距探测器计数。在相同地层模型条件下, 选取10个位置点, 利用蒙特卡罗方法模拟短源距和长源距探测器计数。根据如图 5建立的短源距和长源距探测器计数与地层密度的转换关系, 将快速正演模拟和蒙特卡罗模拟的探测器计数转换为地层密度, 结果如图 7所示。

图 7 斜井中快速正演模拟与蒙特卡罗模拟结果对比 Fig.7 Comparison between fast forward simulation and Monte Carlo simulation results in deviated wells

密度成像测井是在测井过程中记录不同方位的地层密度值, 实现井周地层密度的方位测量, 并将地层密度值按方位更加直观的显示[21]。本文中还实现了密度成像测井的快速正演模拟, 如图 8(a)所示, 倾斜地层界面与井轴的交点为b, 但在不同测量方位上地层界面与井轴的交点不同, 利用柱坐标下的地层界面与井筒方程可确定不同测量方位下地层界面与井轴的交点。按如图 8(b)所示的记录扇区示意图, 在每个采样点模拟计算16个扇区方位的地层密度测井值, 得到的密度成像如图 9所示。

图 8 斜井地层界面与井筒关系及不同记录扇区示意图 Fig.8 Relationship between formation interface and wellbore in deviated wells, and schematic diagram of different recording sectors
图 9 密度成像测井模拟结果 Fig.9 Results of density imaging logging simulation

图 9为测井仪器从密度为2.567 g/cm3的地层上提到密度为2.155 g/cm3的地层中模拟测量的方位密度成像图, 模型的地层界面相对倾角α为70°。利用模拟的方位密度成像图计算得出的地层界面相对倾角为68.5°, 与地层真实相对倾角的相对误差为2.1%。

3.3 模拟结果误差分析

蒙特卡罗方法是放射性测井中常用的数值模拟方法, 可以精确设定计算模型, 得到高精度的数值模拟结果, 所以将蒙特卡罗方法模拟结果作为理论精确值。在直井和斜井条件下的快速正演模拟和蒙特卡罗模拟结果的交会图如图 10所示, 图中蓝色线为对角线, 计算结果越靠近蓝色线表明精度越高。

图 10 快速正演模拟与蒙特卡罗模拟结果交会图 Fig.10 Cross plot of fast forward simulation and Monte Carlo simulation results

为验证快速正演模拟结果的精确性, 利用如下公式计算快速正演模拟与蒙特卡罗模拟结果的绝对误差平均值EAA、相对误差平均值EAR及皮尔逊相关系数Pr:

$ {E_{{\rm{AA}}}} = \frac{{\sum\limits_{i = 1}^{i = n} {\left| {{S_{{\rm{F}}i}} - {S_{{\rm{M}}i}}} \right|} }}{n}, $ (17)
$ {E_{{\rm{AR}}}} = \frac{{\sum\limits_{i = 1}^{i = n} {\left( {\left| {{S_{{\rm{F}}i}} - {S_{{\rm{M}}i}}} \right|/{S_{{\rm{M}}i}}} \right)} }}{n}, $ (18)
$ {P_{\rm{r}}} = \frac{{n\sum\limits_{i = 1}^n {{S_{{\rm{F}}i}}{S_{{\rm{M}}i}}} - \sum\limits_{i = 1}^n {{S_{{\rm{F}}i}}\sum\limits_{i = 1}^n {{S_{{\rm{M}}i}}} } }}{{\sqrt {n\sum\limits_{i = 1}^n {S_{{\rm{F}}i}^2} - {{\left( {\sum\limits_{i = 1}^n {{S_{{\rm{F}}i}}} } \right)}^2}} \sqrt {n\sum\limits_{i = 1}^n {S_{{\rm{M}}i}^2} - {{\left( {\sum\limits_{i = 1}^n {{S_{{\rm{M}}i}}} } \right)}^2}} }}. $ (19)

式中, n为采样数据点的个数; SFi为在第i采样点的快速正演模拟结果; SMi为在第i采样点的蒙特卡罗方法模拟的结果。

快速正演模拟和蒙特卡罗模拟结果的绝对误差平均值、相对误差平均值及相关系数计算值以及计算速度如表 1所示。

表 1 快速正演模拟和蒙特卡罗模拟结果误差分析和计算速度 Table 1 Error analysis and calculation speed of fast simulation and Monte Carlo simulation results

根据图 10所示的交会图以及表 1所示误差分析和计算速度结果, 可以看出快速正演模拟与蒙特卡罗模拟结果基本一致, 快速正演模拟可提供与蒙特卡罗模拟精度一致的补偿密度测井正演模拟结果。但是, 本文快速正演模拟方法的计算速度非常快, 利用快速正演模拟方法计算本文中直井或斜井实例的密度测井曲线可在2 s内完成, 而利用蒙特卡罗方法则需要大约60 h, 计算速度提高了十万倍以上。本文方法可在极短时间内得到与蒙特卡罗方法精度基本一致的补偿密度测井正演模拟结果, 可以满足水平井井眼轨迹与地层关系确定补偿密度测井实时正演模拟计算的需求。

4 快速正演模拟应用实例

图 11为在四川页岩气水平井X井中快速正演模拟实例, 第1道为自然伽马曲线, 第2道为深度道, 第3道为相位电阻率曲线, 第4道为密度测井曲线及提取的地层模型, 第5道为实际测量的随钻方位密度成像图像, 第6道为根据密度测井曲线提取的地层模型快速正演模拟计算的方位密度成像图像。该模拟实例中30 m的井段正演模拟计算时间约需要20 s, 快速正演模拟计算的结果与实际测井资料基本一致。

图 11 页岩气X井中快速正演模拟应用 Fig.11 Application of fast forward simulation in shale gas X wells
5 结论

(1) 提出一种补偿密度测井快速正演模拟方法, 引入空间响应分布函数表征伽马场分布及地层不同位置体积元对探测器计数贡献。采用蒙特卡罗方法预先建立快速正演模拟所需数据库, 利用数据库中的空间响应分布函数计算不同地层区域对探测器计数贡献, 并结合数据库中单一密度地层中模拟记录的探测器计数, 得到长、短源距探测器计数, 进而根据数据库中的转换关系将探测器计数转换为地层密度, 可获取高精度的补偿密度测井快速正演模拟结果。

(2) 本文快速正演模拟方法的计算精度与蒙特卡罗方法基本一致, 但计算速度提高十万倍以上, 能够满足水平井井眼轨迹与地层关系确定中快速正演计算的需求。

(3) 补偿密度测井快速正演模拟结果的计算精度和计算速度与实际测井资料基本一致, 验证了新方法的实用性。

参考文献
[1]
曲占庆, 田雨, 李建雄, 等. 水平井多段分簇压裂裂缝扩展形态数值模拟[J]. 中国石油大学学报(自然科学版), 2017, 41(1): 102-109.
QU Zhanqing, TIAN Yu, LI Jianxiong, et al. Numerical simulation study on fracture extension and morphology of multi-cluster staged fracturing for horizontal wells[J]. Journal of China University of Petroleum (Edition of Natural Science), 2017, 41(1): 102-109.
[2]
何佑伟, 程时清, 胡利民, 等. 多段压裂水平井不均匀产油试井模型[J]. 中国石油大学学报(自然科学版), 2017, 41(4): 116-123.
HE Youwei, CHENG Shiqing, HU Limin, et al. A pressure transient analysis model of multi-fractured horizontal well in consideration of unequal production of each fracture[J]. Journal of China University of Petroleum (Edition of Natural Science), 2017, 41(4): 116-123.
[3]
胡松, 周灿灿, 王昌学, 等. 泥页岩油气藏水平井评价对策与实践[J]. 地球物理学进展, 2013, 28(4): 1877-1885.
HU Song, ZHOU Cancan, WANG Changxue, et al. The strategies and practice of horizontal well evaluation in shale hydrocarbon reservoir[J]. Progress in Geophysics, 2013, 28(4): 1877-1885. DOI:10.6038/pg20130429
[4]
ZHOU Q. Log interpretation in high-deviation wells through user-friendly tool-response processing[C/OL]//SPWLA 49th Annual Logging Symposium, Austin, Texas, USA, 25-28 May, 2008[2017-07-22]. https://www.onepetro.org/conference-paper/SPWLA-2008-AAAA.
[5]
巫振观, 范宜仁, 王磊, 等. 随钻方位电磁波测井仪器偏心响应模拟及分析[J]. 中国石油大学学报(自然科学版), 2017, 41(5): 69-79.
WU Zhenguan, FAN Yiren, WANG Lei, et al. Numerical modeling and analysis of eccentricity effect on borehole response of azimuthal electromagnetic logging while drilling tool[J]. Journal of China University of Petroleum (Edition of Natural Science), 2017, 41(5): 69-79.
[6]
许巍, 柯式镇, 李安宗, 等. 随钻电磁波测井仪器结构影响的三维有限元模拟[J]. 中国石油大学学报(自然科学版), 2016, 40(6): 50-56.
XU Wei, KE Shizhen, LI Anzong, et al. Structural effects analysis of an electromagnetic wave propagation resistivity LWD tool by 3D finite element method[J]. Journal of China University of Petroleum (Edition of Natural Science), 2016, 40(6): 50-56.
[7]
BITTAR M S, KLEIN J D, RANDY B, et al. A new azimuthal deep-reading resistivity tool for geosteering and advanced formation evaluation[J]. SPE Reservoir Evaluation and Engineering, 2009, 12(2): 270-279. DOI:10.2118/109971-PA
[8]
YIN H Z, ZHOU J J, GUO P J, et al. A hybrid solution for fast 3D gamma ray tool modeling in high angle and horizontal wells[C/OL]//SPWLA 49th Annual Logging Symposium, Austin, Texas, USA, 25-28 May, 2008[2017-07-22]. https://www.onepetro.org/conference-paper/SPWLA-2008-YYY.
[9]
邵才瑞, 曹先军, 陈国兴, 等. 随钻伽马测井快速正演算法及地质导向应用[J]. 地球物理学报, 2013, 56(11): 3932-3942.
SHAO Cairui, CAO Xianjun, CHEN Guoxing, et al. A fast forward algorithm for LWD gamma-ray response and its geosteering application[J]. Chinese Journal of Geophysics, 2013, 56(11): 3932-3942. DOI:10.6038/cjg20131135
[10]
BOOTH T E, HENDRICKS J S. Importance estimation in forward Monte Carlo calculation[J]. Nuclear Technology, 1984, 5: 90-100.
[11]
WATSON C C. Monte Carlo computation of differential sensitivity functions[J]. Transactions of American Nuclear Society, 1984, 46: 655-656.
[12]
CASE C, WATSON C C, LAWRENCE R. Sensitivity function technique for modeling nuclear tools: 5334833[P]. 1994-06-14.
[13]
LIU L, GARDNER R P. A geometry independent fine-mesh-based Monte Carlo importance generator[J]. Nuclear Science and Engineering, 1997, 125: 188-195. DOI:10.13182/NSE97-A24265
[14]
MENDOZA A, TORRES-VERDÍN C, PREEG B. Linear iterative refinement for the rapid simulation of borehole nuclear measurements:part Ⅱ-high-angle and horizontal wells[J]. Geophysics, 2010, 75(2): 80-89.
[15]
吴文圣, 黄隆基. 三探测器密度测井的Monte Carlo模拟[J]. 地球物理学报, 2004, 47(1): 164-170.
WU Wensheng, HUANG Longji. Monte Carlo simulation of three-detector density logging[J]. Chinese Journal of Geophysics, 2004, 47(1): 164-170.
[16]
于华伟, 孙建孟, 李晓鹏, 等. 蒙特卡罗模拟研究水平井密度测井中岩屑层的影响[J]. 核技术, 2009, 32(10): 752-755.
YU Huawei, SUN Jianmeng, LI Xiaopeng, et al. Monte Carlo simulation on influence of cutting beds on density log in the horizontal wells[J]. Nuclear Techniques, 2009, 32(10): 752-755.
[17]
袁超, 周灿灿, 张锋, 等. 蒙特卡罗模拟在放射性测井中的应用[J]. 地球物理学进展, 2014, 29(5): 2301-2306.
YUAN Chao, ZHOU Cancan, ZHANG Feng, et al. Application of Monte Carlo simulation in nuclear well logging[J]. Progress in Geophysics, 2014, 29(5): 2301-2306. DOI:10.6038/pg20140546
[18]
葛云龙. 水平井轨迹与地层之间关系的确定以及密度测井校正[D]. 大庆: 东北石油大学, 2017.
GE Yunlong. Determining the relationship between horizontal well trajectory and formation and density log correction[D]. Daqing: Northeast Petroleum University, 2017. https://www.doc88.com/p-9062887891518.html
[19]
BRIESMEISTER J F. MCNP:a general Monte Carlo n-particle transport code[M]. Washington: Los Alamos National Laboratory, 2000.
[20]
张泉滢, 张锋, 王玉伟, 等. 随钻中子伽马密度测井的双源距含氢指数校正方法[J]. 中国石油大学学报(自然科学版), 2017, 41(4): 78-84.
ZHANG Quanying, ZHANG Feng, WANG Yuwei, et al. Dual-spacing hydrogen index correction method for neutron gamma density measuring in LWD[J]. Journal of China University of Petroleum (Edition of Natural Science), 2017, 41(4): 78-84.
[21]
TRUAX J, GUO W J, GOODYEAR G, et al. Response of LWD density logs at a bed boundary[C/OL]//SPWLA 49th Annual Logging Symposium, Austin, Texas, USA, 25-28 May, 2008[2017-07-22]. https://www.onepetro.org/conference-paper/SPWLA-2008-T.