A Numerical Mode-Matching Method Based on Multiple-Knot Cubic B-Splines Applied to the Analysis of Well-Logging Tools

Abstract

The numerical mode-matching (NMM) technique is a very efficient approach to analyze well-logging tools used for hydrocarbon exploration. This problem is typically modeled as a cylindrically layered medium (subterranean formation) including both horizontal and vertical stratifications. In this paper, we present a new NMM formulation based on the use of B-splines expansion functions with variable knot multiplicity to represent the dependency of fields on the vertical direction, and cylindrical Bessel functions to represent the fields along the radial (horizontal) direction. We show that this type of expansion allows to better capture the spatial variation of the electromagnetic field near vertical layer interfaces. Illustrative examples elucidate the advantage of knot multiplicity in the modeling of earth formations with anisotropic and lossy layers. A case study on compensated dual resistivity (CDR) well-logging tool sensors is also presented based on the proposed method.

Index Terms
Numerical mode-matching; cubic B-splines; well-logging tools; compensated dual resistivity

I. Introduction

LOGGING-WHILE-DRILLING (LWD) tools are electromagnetic sensors routinely used to estimate constitutive parameters (permittivity and conductivity) of subterranean formations [¹]. This information is useful is assessing the potential for hydrocarbon bearing zones.

Continual advances in computing power have enabled the use of brute-force numerical techniques such as finite-differences (FD) [²-⁵] and finite volumes (FV) [⁷] to model LWD tools operating in complex geophysical environments. However, it is still impractical to employ brute-force methods in inversion efforts where a large number of forwards simulation is required for a parameter space sweep. On the other hand, methods that combine both numerical and analytical steps [⁸-¹³], such as the numerical mode-matching (NMM) [¹⁴-¹⁷], are able to model complex geophysical formations while requiring only a modest amount of computational resources compared to brute-force methods.

In this paper, we present a new NMM formulation based on the use of cubic B-splines with variable knot multiplicity to expand the fields along the vertical direction. The knot multiplicity of these basis functions allow a better representation of the field behavior across interfaces while using a coarser computational grid, leading to a lesser demand on computational resources.

The remaining of this work is organized as follows. Section II summarizes the proposed formulation of NMM using cubic B- splines. Section III shows numerical results to validate the technique, drawing attention to some advantages of the proposed modal expansion based on multiple knots. Section IV presents a case study on compensated dual resistivity (CDR) well-logging tools. Finally, Section V provides some concluding remarks.

II. Numerical Mode-Matching Formulation

The geometry of this problem can be adequately described in terms of a cylindrically layered medium composed of both radial (horizontal) and vertical stratifications [¹⁷]. In the NMM formulation presented in [¹⁵], the electromagnetic field is expanded in terms of basis functions that depend on the vertical coordinate z, while reflection and transmission matrices (obtained by enforcing the appropriate boundary conditions across interfaces) are used to represent the dependency along the radial coordinate ρ in terms of Bessel and Hankel functions.

A. Radial NMM

Assuming a time-harmonic dependence for the fields and sources in the form exp(-iωt), we obtain the two vector wave equations for the electric and magnetic fields, respectively, as [⁶]

(1) $${\mu }_h\nabla \times \left({\mu }^{-1}\nabla \times \overrightarrow{E}\right)-{\omega }^2{\mu }_h\overline{\overline{\epsilon }}\overrightarrow{E}=i\omega {\mu }_h\overrightarrow{J},\text{and}$$

(2) $${\epsilon }_h\nabla \times \left({\overline{\overline{\epsilon }}}^{-1}\nabla \times \overrightarrow{H}\right)-{\omega }^2{\epsilon }_h\overline{\overline{\mu }}\overrightarrow{H}={\epsilon }_h\nabla \times {\overline{\overline{\epsilon }}}^{-1}\overrightarrow{J},$$

where the (anisotropic) constitutive properties of the medium are represented by the permittivity and permeability tensors, $\overline{\overline{\epsilon }}$ and $\overline{\overline{\mu }}\left[12,15\right]$. The impressed electrical current $\overrightarrow{J}$ is produce by horizontal coil antennas wrapped around a metallic mandrel with radius ρ_T. In addition, we assume the antennal to be axially located at z_T. Using a cylindrical coordinate system (ρ,ϕ, z) aligned with the tool axis, the fields given in (1) and (2) can be represented by means of the Bessel transforms [¹⁸] as

(3) $$F\left(\rho ,\varphi ,z\right)=\frac{1}{\sqrt{2\pi }}{\displaystyle \sum _{v=-\infty }^{\infty }{e}^{iv\varphi }{f}_v\left(\rho ,z\right),}$$

(4) $$f_v\left(\rho ,z\right)=\frac{1}{\sqrt{2\pi }}{\displaystyle \underset{0}{\overset{\infty }{\int }}{\tilde{f}}_v\left(k_{\rho },z\right)J_v\left(k_{\rho }\rho \right)k_{\rho }d{k}_{\rho },}$$

where J_v(k_ρρ) is a Bessel function of order v, and ${\tilde{f}}_v\left(k_{\rho },z\right)$ is the Hankel transformation of f_v(ρ,z). Using (1), (3) and (4), we can rewrite (1) in terms of the z component of the electric field for a source-free region as

(5) $$\left(-\frac{1}{{\epsilon }_v}{k}_{\rho }^2+\frac{\partial }{\partial z}\frac{1}{{\epsilon }_v}\frac{\partial }{\partial z}+{\omega }^2{\mu }_h\right){\tilde{d}}_{zv}\left(k_{\rho },z\right)=0$$

where ${\tilde{d}}_z={\epsilon }_v{\tilde{e}}_z$ is the electric flux density, which is expanded in a set of B-splines functions S_n(z):

(6) $${\tilde{d}}_{zv}\left(k_{\rho },z\right)\to {\tilde{d}}_z\left(z\right)={\displaystyle \sum _{n=1}^{\infty }{a}_n{S}_n(z).}$$

The substitution of (6) in (5) results in

(7) $${\displaystyle \sum _{n=1}^{\infty }{a}_n\left[\frac{\partial }{\partial z}\frac{1}{{\epsilon }_v}\frac{\partial }{\partial z}+w^2{\mu }_h\right]S_n\left(z\right)=k_{\rho }^2{\displaystyle \sum _{n=1}^{\infty }{a}_n\frac{1}{{\epsilon }_v}{S}_n\left(z\right),}}$$

where in practice, the series are truncated using a sufficient number of terms for a given grid resolution. The resulting problem can be arranged in a matrix form as

(8) $${\overline{\overline{L}}}_{\epsilon }\cdot {\overline{a}}_{\epsilon n}=k_{\rho }^2{\overline{\overline{p}}}_{\epsilon }\cdot {\overline{a}}_{\epsilon n}.$$

Equation (8)) represents a generalized eigenvalue problem where $k_{\rho }^2$ are the eigenvalues and $\overline{a}$ are the eigenvectors. Note that that matrices ${\overline{\overline{L}}}_{\epsilon }$ and ${\overline{\overline{p}}}_{\epsilon }$ are symmetric in reciprocal media (i.e., with symmetric $\overline{\overline{\epsilon }}$ and $\overline{\overline{\mu }}$).

In regions with sources, the z-component of the fields obey the following equation:

(9) $$\begin{array}{l}\left(\frac{1}{\rho }\frac{\partial }{\partial \rho }\rho \frac{\partial }{\partial \rho }-\frac{v^2}{{\rho }^2}+{\epsilon }_v\frac{\partial }{\partial z}{\epsilon }_v^{-1}\frac{\partial }{\partial z}+{\omega }^2{\mu }_v{\epsilon }_v\right)d_z\left(\rho ,z\right)=\\ i\omega {\mu }_h{\epsilon }_v{j}_z\left(\rho ,z\right),\end{array}$$

where the electric flux density is given by the expansion:

(10) $$d_z\left(\rho ,z\right)={\displaystyle \sum _{q=1}^N{b}_q\left(\rho \right){\psi }_q\left(z\right),}$$

(11) $${\psi }_q={\displaystyle \sum _{n=1}^N{\left({\overline{a}}_q\right)}_n{S}_n\left(z\right).}$$

Equation (11) can represent the electric flux density or the magnetic flux density. For TE modes, we have b_q(ρ) = b_μq(ρ), which represents the ρ -dependence for each modal field, and similarly ψ_q(z) = ψ_μq(z). The expansion described in (10) and (11) can be used to rearrange (9): by taking the product of both sides of (9) with by ${\psi }_{\mu \rho }\left(z\right){\mu }_v^{-1}$ and then applying an inner product along z, we find

(12) $$\begin{array}{l}\left(\frac{1}{\rho }\frac{\partial }{\partial \rho }\rho \frac{\partial }{\partial \rho }-\frac{v^2}{{\rho }^2}+k_{\rho q}^2\right)b_{\mu q}\left(\rho \right)=\\ -j_{T\mu \upsilon q}\frac{1}{\rho }\frac{\partial }{\partial \rho }\left(\rho \delta \left(\rho -{\rho }_T\right)\right),\end{array}$$

where j_Tμvq represents the spectral contribution of the source for the TE modes. Note that only TE modes are excited by horizontal coil antennas as considered here.

B. Transimpedance Calculation

With the expressions of the fields in the source region, we can now obtain the values of the induced voltages at a given receiver coil due to a current in the transmitter antenna. The transimpedance is defined as

(13) $$Z_{RT}=\frac{V_R}{I_T}=\frac{-1}{I_T}{\displaystyle \int \overrightarrow{E}\cdot d{\overrightarrow{l}}_R,}$$

Where, since the receiver coil antennas are horizontal, we only need to consider the azimuthal electric field component. The latter can be written in the following form [⁶] :

(14) $${\overline{E}}_{R\varphi }\left(\rho ,z\right)=\overline{\overline{S^t}}\left(z\right)\cdot {\overline{\mathrm{\Gamma }}}_{R\varphi }\left(\rho \right),$$

where

(15) $${\overline{\mathrm{\Gamma }}}_{R\varphi }=\{\begin{array}{l}\left[{\overline{\overline{Y}}}_1^{+}\cdot {\overline{\overline{H}}}_1\left(\rho \right)\cdot {\widetilde{\overline{\overline{R}}}}_{1,0}^{+}+{\overline{\overline{Y}}}_1^{-}\cdot {\overline{\overline{J}}}_1\left(\rho \right)\right]\cdot {\overline{d}}_{T1}^{-}\rho \le {\rho }_T\\ \left[{\overline{\overline{Y}}}_1^{+}\cdot {\overline{\overline{H}}}_1\left(\rho \right)+{\overline{\overline{Y}}}_1^{-}\cdot {\overline{\overline{J}}}_1\left(\rho \right)\cdot {\widetilde{\overline{\overline{R}}}}_{1,2}^{-}\right]\cdot {\overline{d}}_{T1}^{+}\rho \ge {\rho }_T\end{array}$$

is a factor depending on the reflection coefficients across layers. See [¹⁵] for more details on the evaluation of (13).

C. Cubic B-Spline Eigenfunctions

In the NMM formulations described in [⁶] and [¹⁹], the electromagnetic fields were approximated by triangular, sinusoidal and quadratic B-splines basis functions. In the present work, we adopt cubic B-splines as base (expansion) functions [²⁰] because of their high-order nature and ability to handle knots of higher multiplicity. The function S_n (z) employed in (6) can be each written as a combination of four polynomial terms:

(16) $$S_n\left(z\right)=Se{g}_1\left(z\right)+Se{g}_2\left(z\right)+Se{g}_3\left(z\right)+Se{g}_4\left(z\right),$$

where

(17) $$Se{g}_1\left(z\right)=\frac{z_{10}{}^2}{z_{30}{z}_{20}}{u}^3,$$

(18) $$\begin{array}{l}Se{g}_2\left(z\right)=\frac{z_{10}{}^2}{z_{30}{z}_{20}}+3\left(\frac{z_{21}{z}_{10}}{z_{30}{z}_{20}}\right)u+3\left(\frac{z_{21}{}^2}{z_{30}{z}_{20}}\right)u^2\\ -\left(\frac{z_{21}{}^2}{z_{30}{z}_{20}}+\frac{z_{21}{}^2}{z_{30}{z}_{31}}+\frac{z_{21}{}^2}{z_{41}{z}_{31}}\right)u^3,\end{array}$$

(19) $$\begin{array}{l}Se{g}_3\left(z\right)=\left(\frac{z_{32}{z}_{20}}{z_{30}{z}_{31}}+\frac{z_{42}{z}_{21}}{z_{41}{z}_{31}}\right)+3\left(\frac{z_{32}{}^2}{z_{30}{z}_{31}}-\frac{z_{21}{z}_{32}}{z_{41}{z}_{31}}\right)u\\ -3\left(\frac{z_{32}{}^2}{z_{30}{z}_{31}}+\frac{z_{32}{}^2}{z_{41}{z}_{31}}\right)u^2+\left(\frac{z_{32}{}^2}{z_{30}{z}_{31}}+\frac{z_{32}{}^2}{z_{41}{z}_{31}}+\frac{z_{32}{}^2}{z_{41}{z}_{42}}\right)u^3,\end{array}$$

(20) $$Se{g}_4\left(z\right)=\frac{z_{43}{}^2}{z_{41}{z}_{42}}{\left(1-u\right)}^3,$$

where

(21) $$\begin{array}{l}{z}_{10}=z_1-z_0\dots z_{40}=z_4-z_0,\\ z_{21}=z_2-z_1\dots z_{41}=z_4-z_1,\\ u=\frac{z-z_{i+1}}{z_{i+2}-z_{i+1}}\end{array}$$

Each B-spline has five knots (nodes). When these knots are equally spaced, we obtain a uniform B-spline, as shown in in Fig. 1.

Fig. 1
The four polynomials that compose a uniform cubic B-spline with five equally space knots.

When two or more knots are collocated in space, we obtain a knot with higher multiplicity. Fig. 2 illustrates the behavior of a cubic B-spline function for different knot multiplicities. The top curve shows a case with single knots only. The mid curve shows a case with one knot with double multiplicity. The bottom curve shows a case with one knot with triple multiplicity. By increasing the multiplicity, one decreases the regularity of the function. At the knot with double multiplicity the cubic B -spline has discontinuous second derivative, whereas at the knot with triple multiplicity the cubic B-spline has discontinuous first derivative. This properties allows for properly capturing the behavior of the electromagnetic fields from the boundary conditions across layer interfaces.

Fig. 2
Cubic B-spline behavior for different knot multiplicities.

III. The Knots-Multiplicity Effects

In order to first examine the effect of the multiplicity of knots of the cubic B-splines functions in a simpler setting, we consider the problem of a parallel plate guide filled by two layers in the z-direction, as shown in Fig. 3. This problem admits an analytical solution as described in [²¹]. The eigenvalues for TM and TE (to z) waves satisfy

(22) $$k_{\rho \epsilon }^2+\frac{{\epsilon }_v}{{\epsilon }_h}{k}_z^2={\epsilon }_v{k}_0^2,\text{and}$$

(23) $$k_{\rho \mu }^2={\epsilon }_h{k}_0^2-k_z^2,$$

respectively. In the examples that follow, we compare the eigenvalues obtained by the NMM using the cubic B-splines described above with the values predicted by (22) and (23).

Fig. 3
Parallel plate waveguide filled with two vertical layers.

A. Isotropic and Lossless Layers

Consider a scenario with d = 1000 m and a = 500 m, as indicated in Fig. 3. The upper layer has ε₁ = ε₀ and μ₁ = μ₀ and the bottom layer has ε₂ = 2.55 ε₀ and μ₂= μ₀. The relative error in the solution of the radial wavenumber k_ρ for the first 9 modes with expansions using multiple knots (MK) or only single knots (SK) are depicted in Fig. 4 and Fig. 5, for the TM and TE fields, respectively. We can observe that lower order modes present a significant reduction by using MK. The differences between the MK and SK results are more pronounced for TM modes. This behavior can be explained by the presence of a double derivative with respect to z in the first term within the brackets in (7). In this case, the MK expansion provides, of course, a more accurate result compared to the SK expansion.

Fig. 4
Relative error on TM eigenvalues, in the case of isotropic and lossless layers.

Fig. 5
Relative error on TE eigenvalues, in the case of isotropic and lossless layers.

B. Isotropic and Lossy Layers

As a second example, consider now a parallel plate waveguide filled with lossy layers. Layer 1 has ε₁ = ε₀ and σ₁ = 1 S/m, and layer 2 has ε₂ = ε₀ and σ₂ = 100 S/m. In addition the vacuum permeability is present everywhere. The distance between the plates is now d = 20 m. Fig. 6 and Fig. 7 present the relative error in the solution of k_ρ for the first 10 modes with MK and SK. We again observe smaller errors by using MK compared to SK, for both TM and TE modes.

Fig. 6
Relative error on TM eigenvalues, in the case of isotropic and lossy layers.

Fig. 7
Relative error on TE eigenvalues, in the case of isotropic and lossy layers.

C. Anisotropic and Lossy Layers

In a third example, consider uniaxially anisotropic and lossy layers. The upper layer has ε_1,v,h = ε₀ and σ_1,ν = 1 S/m and σ_1,h = 100 S/m. The bottom layer has ε_2,v,h = ε₀, σ_2,ν = 100 S/m and σ_{2 h} = 1 S/m. Again, we consider d = 20m. Fig. 8 and Fig. 9 show results for the first ten modes using SK and MK, for both TM and TE modes, respectively. Once more, the advantage of the MK approach is clear compared to SK.

Fig. 8
Relative error on TM eigenvalues, in the case of anisotropic and lossy layers.

Fig. 9
Relative error on TE eigenvalues, in the case of anisotropic and lossy layers.

IV. Case Study: Borehole with Rugosity

Next, we consider a compensated dual-resistivity (CDR) LWD tool traversing a rugose borehole. This problem was considered before in [²², ²³] and is examined here to validate our formulation and illustrate the ability of the proposed algorithm based on MK cubic B-splines to efficiently model this kind of problem. The CDR sensor operates at 2MHz and is composed by two transmitters and two receivers as depicted in Fig. 10. Typically, this sensor yields two apparent resistivities, R_ps and R_da, associated with the phase difference and the amplitude ratio of the voltage inducted at the receiver antennas. The R_ps resistivity provides information of radial regions close to the antennas while the R_da resistivity provides information from regions farther from the antennas.

Fig. 10
A CDR LWD sensor tool.

The transmitters and receivers are modelled as coil antennas with radius of 3.5 in. The conductor mandrel is modeled as a PEC cylinder with radius of 3.25 in. The sensor antennas are placed concentrically around the mandrel, and their axial positions are represented in Fig. 10. The CDR sensor moves along a 2.5-in-rugose borehole wall, as depicted in Fig. 11, inside a 3.25-inborehole with resistivity 0.1 Ω.m, surrounded by a soil formation with resistivity 2 Ω.m. For further details on this geometry, see [²², ²³].

Fig. 11
CDR sensor tool inside a rugose borehole.

To eliminate spurious reflections arising from the truncation of the domain along z, we employed a cylindrical perfectly matched layer with properties similar to that used in [²⁴]. Initially, we use a uniform grid along z with Δz = 0.020 m within the domain of interest and based on SK. The resistivity results due to transmitter 2 excitation are shown in Fig. 12 (note that, in this uncompensated case, transmitter 1 is not active). The abscissa axis represents the midpoint between the two receivers, i.e., (z_r1 + z_r2) / 2. We can clear see that the rugosity in the borehole wall disturbs the response of the sensor resistivity compared to the smooth wall solution. This can lead to misinterpretations when reading the LWD sensor response.

Fig. 12
Amplitude and phase apparent resistivities of a CDR sensor tool due the transmitter 2 alone (uncompensated).

In addition, Fig. 13 shows the compensated CDR response obtained via the proposed NMM with MK cubic B-splines. The compensated response is obtained by using both transmitters individually, and subsequently averaging the phase and amplitude results. The discretization grid along z in this case has Δz=0.03048 m in the washout region, and Δz=0.0217 m out of it. Very good agreement is observed versus the results reported in [²²]. It should be noted that the MK approach requires 55% less discretization points near the vertical interfaces compared to the (SK based) uniform grid.

Fig. 13
Amplitude and phase apparent resistivities of a compensated CDR sensor tool modeled by the NMM employing MK cubic B-splines.

V. Conclusion

We have introduced a new NMM formulation based on the use of cubic B-splines with multiple knots as expansion functions, for the modeling of well-logging tools. Numerical results showed that the use of knots with higher multiplicity allows for a better representation of the fields near interface layers of typical geophysical scenarios. The presented technique is able to accurately model complex media including anisotropic and lossy layers using less grid points than conventional (uniform) cubic B-splines. As such, the proposed technique represents a useful improvement on the modeling of well-logging-tools.

• This work was supported in part by the Brazilian Agencies CAPES and CNPq under doctoral fellowships 2444/11-9 and 150454/2016-3.

References

• ^[1]
  S. Y. Chen, W. C. Chew, V. R. N. Santos, K. Sainath, and F. L. Teixeira, “Electromagnetic subsurface remote sensing,” in: Wiley Encyclopedia of Electrical and Electronics Engineering, John Wiley & Sons, New York, 2016.
• ^[2]
  Y.-K. Hue, F. L. Teixeira, L. San Martin, and M. Bittar, ‘Three-dimensional simulation of eccentric LWD tool response in boreholes through dipping formations,’ IEEE Trans. Geosci. Remote Sens., vol. 43, no. 2, pp. 257-268, 2005.
• ^[3]
  H.O. Lee and F. L. Teixeira, ‘Cylindrical FDTD analysis of LWD tools through anisotropic dipping layered earth media,’ IEEE Trans. Geosci. Remote Sens, vol. 45, no. 2, pp. 383-388, 2007.
• ^[4]
  D. Wu, C. R. Liu, and J. Chen, “An efficient FDTD method for axially symmetric LWD environments,” IEEE Trans. Geosci. Remote Sens., vol. 46, no. 6, pp. 1652-1656, 2008.
• ^[5]
  H. O. Lee, F. L. Teixeira, L. E. San Martin, and M. S. Bittar, ‘Numerical modeling of eccentered LWD borehole sensors in dipping and fully anisotropic Earth formations,’ IEEE Trans. Geosci. Remote Sens., vol. 50, no. 3, pp. 727-735, 2012.
• ^[6]
  M. S. Novo, L. C. Silva, and F. L. Teixeira, ‘Three-dimensional finite-volume analysis of directional resistivity logging sensors,’ IEEE Trans. Geosci. Remote Sens., vol. 48, no. 3, pp. 1151-1158, 2010.
• ^[7]
  M. S. Novo, L. C. da Silva and F. L. Teixeira, “Finite volume modeling of borehole electromagnetic logging in 3-D anisotropic formations using coupled scalar-vector potentials,” IEEE Antennas Wireless Propag. Lett., vol. 6, pp. 549-552, 2007.
• ^[8]
  J. R. Lovell and W. C. Chew, “Response of a point-source in a multicylindrically layered medium,” IEEE Trans. Geosci. Remote Sens., vol. 25, no. 6, pp. 850-858, 1987.
• ^[9]
  J. R. Lovell and W. C. Chew, “Effect of tool eccentricity on some electrical well-logging tools,” IEEE Trans. Geosci. Remote Sens., vol. 28, no. 1, pp. 127-136, 1990.
• ^[10]
  Z.-Q. Zhang and Q.-H. Liu, “Applications of the BCGS-FFT method to 3-D induction well-logging problems,” IEEE Trans. Geosci. Remote Sens., vol. 41, no. 5, pp. 998-1004, 2003.
• ^[11]
  Y.-K. Hue and F. L. Teixeira, ‘Analysis of tilted-coil eccentric borehole antennas in cylindrical multilayered formations for well-logging applications,’ IEEE Trans. Antennas Propagat., vol. 54, no. 4, pp. 1058-1064, 2006.
• ^[12]
  G.-S. Liu, F. L. Teixeira, and G.-J. Zhang, ‘Analysis of directional logging tools in anisotropic and multieccentric cylindrically-layered Earth formations,’ IEEE Trans. Antennas Propagat., vol. 60, no. 1, pp. 318-327, 2012.
• ^[13]
  H. Moon, B. Donderici, and F. L. Teixeira, ‘Stable evaluation of Green's functions in cylindrically stratified regions with uniaxial anisotropic layers,’ J. Comp. Phys., vol. 325, pp. 174-200, 2016.
• ^[14]
  G.-X. Fan, Q. H. Liu, and S. P. Blanchard, “3-D numerical mode-matching (NMM) method for resistivity well-logging tools,” IEEE Trans. Geosci. Remote Sens, vol. 38, no. 10, pp. 1544-1552, 2000.
• ^[15]
  Y.-K. Hue and F. L. Teixeira, ‘Numerical mode-matching method for tilted coil antennas in cylindrically layered anisotropic media with multiple horizontal beds,’ IEEE Trans. Geosci. Remote Sens., vol. 45, no. 8, pp. 2451-2462, 2007.
• ^[16]
  W. C. Chew, Waves and Fields in Inhomogeneous Media. New York, NY, USA: John Wiley & Sons, 1995.
• ^[17]
  J. Li and L. C. Shen, “Vertical eigenstate method for simulation of induction and MWD resistivity sensors,” IEEE Trans. Geosci. Remote Sens., vol. 31, no. 2, pp. 399-406, Mar. 1993.
• ^[18]
  E. Cavanagh and B. Cook, “Numerical evaluation of Hankel transforms via Gaussian-Laguerre polynomial expansions,” IEEE Trans. Acoustics, Speech, Signal Proc., vol. 27, no. 4, pp. 361-366, Aug. 1979.
• ^[19]
  D. M. Pai, “Induction log modeling using vertical eigenstates,” IEEE Trans. Geosci. Remote Sens., vol. 29, no. 2, pp. 209-213, Mar. 1991.
• ^[20]
  H. Yang, W. Yue, Y. He, H. Huang, and H. Xia, “The deduction of coefficient matrix for cubic non-uniform B-Spline curves,” in 1st Int. Workshop Edu. Technol. Comp. Sci., Hubei, China, pp. 607-609, 2009.
• ^[21]
  R. F. Harrington, Time-Harmonic Electromagnetic Fields. New York, NY, USA: McGraw-Hill, 1961.
• ^[22]
  B. Clark., M. G. Luling, J. Jundt, M. O. Ross, and D. L. Best, “A dual depth resistivity measurement for FEWD,” in SPWLA 29th Annual Logging Symp., paper SPWLA-1988-A, San Antonio, TX, USA, 5-8 June, 1988.
• ^[23]
  B. Clark., D. F. Allen, D. L. Best, S. D. Bonner, J. Jundt, M. G. Luling, and M. O. Ross, “Electromagnetic propagation logging while drilling: Theory and experiment,” SPE Format. Eval., vol. 5, no. 3, pp. 263-271, 1990.
• ^[24]
  F. L. Teixeira and W. C. Chew, “Finite-difference computation of transient electromagnetic waves for cylindrical geometries in complex media,” IEEE Trans. Geosci. Remote Sens., vol. 38, no. 4, pp. 1530-1543, 2000.

Publication Dates

• Publication in this collection
  Sept 2017

History

• Received
  17 Jan 2017
• Reviewed
  19 Jan 2017
• Accepted
  07 June 2017