- Email: [email protected]
3D simulation of transient electromagnetic field for geosteering horizontal wells
Russian Geology and Geophysics 52 (2011) 725–729 www.elsevier.com/locate/rgg
3D simulation of transient electromagnetic field for geosteering horizontal wells E.V. Onegova a,*, M.I. Epov a,b a
b
Baker Hughes Russian Science Center, ul. Kutateladze 4a, Novosibirsk, 630090, Russia A.A. Trofimuk Institute of Petroleum Geology and Geophysics, Siberian Branch of the Russian Academy of Sciences, pr. Akademika Koptyuga 3, Novosibirsk, 630090, Russia Received 13 September 2010; accepted 7 December 2010
Abstract The article discusses numerical simulation of the transient electromagnetic field. The field source is an induction coil. We consider the situation when a logging tool is in a horizontal well in a medium with horizontal and vertical boundaries. The specific features of this problem are the metallic mandrel of the tool, 3D geometry, and distant boundaries. The method of separate computation of normal and anomalous fields is proposed. The finite element method is used for spatial approximation of the field, and the implicit finite-difference scheme is used for approximation of the field in time. Correctness and advantages of the method are shown. Some numerical results are demonstrated. The method proposed can be used when designing tools for geosteering. © 2011, V.S. Sobolev IGM, Siberian Branch of the RAS. Published by Elsevier B.V. All rights reserved. Keywords: geosteering; transient electromagnetic method; finite element method
Introduction As more deviated and horizontal wells are drilled today, geosteering has become a relevant issue. Geosteering is used for correcting the well trajectory based on the measurements performed while drilling, the goal being to achieve the optimal trajectory and thus, to maximize well productivity. When drilling deviated and horizontal wells, it is essential to determine the distances to the top and base of the reservoir, and to its inner boundaries (water-oil and gas-oil contact, clay layers), as early as possible. In this case, the drillstring can be redirected in time to avoid opening the water-bearing bed or the gas cap. The deeper the logging while drilling for the data one may obtain, the more efficient are the solutions which could be available. It is well known that in the transient electromagnetic method (TEM), sensitivity to remote areas of the formation tends to enlarge with the increase of the signal registration time, while closely-located objects seem to lose their influence. The main limitation for using TEM while drilling is the highly-conductive drillstring, as its signal can exceed that of the medium by several orders of magnitude.
First publications on inductive logging using TEM appeared in the 1970s (Kaufman and Sokolov, 1972; Plyusnin and Vilge, 1969). These papers considered the signal from a vertical magnetic dipole in the homogeneous medium (Plyusnin and Vilge, 1969), on the well axis and in beds of finite thickness (Kaufman and Sokolov, 1972). In (Anderson and Chew, 1989), numerical modeling of a well tool with an insulating mandrel in an axially symmetrical medium was described. A number of articles focused on TEM-based well flaw detection (Epov et al., 2002; Potapov and Kneller, 2000; Sidorov, 1996). In this case, the medium can be adequately described by a 2D model, calculations here must take into account both electric conductivity and magnetic permeability of the casing. The paper offers numerical schemes of modeling the situation typical of geosteering. The finite element method is used for modeling the transient electric field in the medium penetrated by the horizontal well, with horizontal layers and vertical boundary, taking into account electromagnetic properties of the drillstring (Fig. 1).
* Corresponding author. E-mail address: [email protected] (E.V. Onegova) 1068-7971/$ - see front matter D 201 1, V . S. S o bolev IGM, Siberian Branch of the RAS. Published by Elsevier B.V. All rights reserved. doi:10.1016/j.rgg.2011.06+.005
726
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
that the medium corresponding to the normal field can contain any axisymmetrical objects, including parts of the tool. The normal field E0 can be described by the following equation: ∂E0 1 ∂J =− , rot rot E0 + σ0 µ ∂t ∂t
Fig. 1. Model of the medium. 1, reservoir cap; 2, drillstring with transmitter and receiver coils; 3, well; 4, reservoir; 5, underlying medium; 6, vertical boundary.
Mathematical model Without taking into account displacement currents, the transient electric field can be described by the following equation: 1 ∂E ∂J =− , rot rot E + σ µ ∂t ∂t
(1)
where E is the vector of the electric field, µ is magnetic permeability, σ is electric conductivity, t is time, and J is density of the extraneous current. Let us set the z axis of the Cartesian coordinates in the line of well (Fig. 1). Because a coil coaxial with the pipe is used as the source, density of the extraneous current J will have only one nonzero component, Jϕ. Changes in density of the extraneous current with time are described by the Heaviside step function. In accordance to the Faraday law of induction, after the current in the source is turned off, eddy currents in the medium tend to be located next to the source. Then, in the process of evolution, the currents begin to diffuse into the medium. Before the currents reach horizontal boundaries, the electric field remains axisymmetrical just as in the coil, and has only one nonzero component, Eϕ. Let us use this fact and, for early times, solve a 2D problem in cylindrical (r, z) coordinates, rather than 3D. The presence of the metal pipe next to the source, as well as distant boundaries whose influence becomes noticeable only at later times, make the problem quite complex and difficult to compute. Because the signal from the currents in the pipe dominates in the total signal, let us divide the total field into the normal field due to the pipe in the axisymmetrical domain E0 = (0, Eϕ (r, z), 0)T and the anomalous field due to the rest of the medium E+ (Onegova, 2010; Tabarovsky and Sokolov, 1982). These fields will be computed separately. In this case, fine grids will be needed only in the 2D problem, and thus, computation time will be greatly reduced. It should be noted
(2)
where σ0 is electric conductivity of the axisymmetrical medium. Therefore, σ0 differs from σ only in those parts of the computational domain that are not axisymmetrical. For the model shown in Fig. 1, σ0 differs from σ in the upper and lower horizontal layers and in the layer ahead the well, because electric conductivity σ0 coincides with that of the middle layer. In this case, the normal field is the field of the pipe in the well in the homogeneous medium. Subtracting Eq. 2 from Eq. 1, we arrive at the equation for the anomalous field E+: ∂E+ ∂E0 1 = (σ0 − σ) rot rot E+ + σ . ∂t ∂t µ
(3)
Thus, the computation process can be reduced to the following: 1) From zero to a certain value of t0, we will be solving the axisymmetrical problem. 2) Starting with the moment t0, we will be each time layer be solving the axisymmetrical problem and the 3D problem. In both problems, we have homogeneous initial conditions: = 0, E0 t = 0
(4)
= 0. E+ t = t0
(5)
Approximate homogeneous Dirichlet boundary conditions are used here: = 0, E0 × n ∂Ω0
(6)
= 0. E+ × n ∂Ω
(7)
where ∂Ω0, ∂Ω are the boundaries of the axisymmetrical and 3D computation problems, n denotes the external normal to the corresponding boundary.
Method of solution To solve initial boundary problems (2), (4), and (6) and (3), (5), and (7) numerically, the finite elements method is used here. The variational formulation of the 2D problem has the following form: 1
∫ µ rot E0 ⋅ rot ψ dΩ + ∫ σ0
Ω0
Ω0
∂E0 ∂t
⋅ ψ dΩ =
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
−∫
Ω0
727
∂J ⋅ ψ dΩ ∀ ψ ∈ H0 (rot, Ω0), ∂t
where H0 (rot, Ω0) = υ ∈ (L2(Ω0))3 rot υ ∈ (L2(Ω0))3, υ × n|∂Ω = 0
0
is the space of test functions, ψ = (0, ψϕ, 0)T. The variational equation for the 3D problem has the following form: 1
∫ µ rot E+ ⋅ rot ψ dΩ + ∫ σ
Ω
Ω
(σ0 − σ) ∫
Ω
∂E0 ∂t
∂E+ ∂t
⋅ ψ dΩ =
Fig. 2. The electromotive force in the homogeneous medium. 1, with pipe; 2, no pipe.
⋅ ψ dΩ ∀ ψ ∈ H0 (rot, Ω),
where ψ = (ψx, ψy, ψz)T. To approximate the time derivative of the solution, the completely implicit three-layer scheme is used (Marchuk 1980) ∂E ≈ γ0 Ej − γ1 Ej−1 + γ2 Ej−2, ∂t t = t j
where γ0 =
∆01 + ∆02
, γ1 =
∆02
, γ12 =
∆01
, ∆12 = ∆01 ∆02 ∆12 ∆01 ∆12 ∆02 tj−1 − tj−2, ∆02 = tj − tj−2, ∆01 = tj − tj−1, ti is the time layer, and Ei is the solution at the ith layer. For the finite-element discretization, rectangles with piecewise-bilinear basis functions are used for the 2D problem, and hexahedrons with edge basis functions of the first order, for the 3D problem (Nedelec, 1980).
Testing The results of 3D modeling using the bedded medium were compared to known solutions (Tabarovsky and Sokolov, 1982). The magnetic dipole was used as the field source. In the 3D formulation, a coil with the radius of 0.01 m was used as the field source. Relative discrepancy between solutions did not exceed 1% in the time range between 10–7 and 10–2 s. Besides, the results of solving the 2D problem were compared with computations obtained using FEMAX (Bespalov, 2002) and COMSOL (Pryor, 2009). The transient process in the homogeneous medium with a metal pipe of a finite length was modeled. A coil coaxial with the pipe was used as the field source. Pipes with diverse electric conductivity and magnetic permeability (magnetic steel, nonmagnetic steel, and copper) and diverse containing media were considered. Average relative discrepancy between solutions obtained using three different methods was 2–3% in the time range between 10–7 and 10–2 s. To understand the degree of the impact the pipe might have upon the signal and show the possibilities of the computational
scheme developed, let us consider some computations. Figure 2 shows the electromotive force induced in the receiver coil as a function of time. Two results are presented: with pipe and without pipe. Formation resistivity was 100 Ohm⋅m; steel pipe resistivity was 7.14 × 10–7 Ohm⋅m; its magnetic permeability was 100µ0 (µ0 = 4π × 10−7H/m); the outer and inner radii of the pipe were 0.04 and 0.07 m, respectively; coil radius was 0.085 m; the distance between the transmitter coil and receiver coil was 5 m; the electric current was 5 A; the number of loops in the transmitter coil was 100; the number of loops in the receiver coil was 10. In the graph, the pipe affects the transient process greatly, both quantitatively and qualitatively. In the range between 5 × 10−4 and 5 × 10−2 s, the signal from the pipe decreases quite slowly, and one can even observe some increase at the time of 2 × 10−2 s. After the time of 0.1 s, the signal from the pipe decreases greatly, and after the time of 2 s, it attenuates as t−5 / 2. Physically, this means there are no more currents in the pipe, and the other currents are far from the initial source that the situation becomes equivalent to that of a magnetic dipole in the homogeneous medium. On the other hand, the signal level at these times is much lower than the measured values. Generally, it should be noted that the computational scheme can model the signal in the dynamic range around 200 dB.
Numerical experiments Let us consider a geoelectrical model typical of the west Siberian oil and gas province. The horizontal interval of the well is in the oil-saturated reservoir (ρ1 = 15 Ohm⋅m, Fig. 3). The reservoir has a clay cap (ρ2 = 4 Ohm⋅m). Under the reservoir, there is a water-saturated interval (ρ3 = 8 Ohm⋅m). In front of the tool, orthogonally with respect to the well, there is a clay layer (ρ4 = 2.0–3.5 Ohm⋅m). The drilling mud with a high clay content has electric resistivity of 2 Ohm⋅m. Resistivity of the pipe made of nonmagnetic steel is 7.14 × 10−7 Ohm⋅m. The outer and inner radii of the pipe were 0.04 and 0.07 m, respectively; coil radius was 0.085 m for
728
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
Fig. 3. Problem domain and model parameters. Fig. 4. Cross section via the xy plane with the hexahedron grid.
both transmitter coil (T) and receiver coil (R); the well radius was 0.108 m. The electric current in the transmitter coil was 5 A; the number of loops in the transmitter coil was 100; the number of loops in the receiver coil was 10. The object of measurements was the electromotive force induced in the receiver coil. Figure 4 shows a grid of hexahedrons for the computational domain. Inside the well, the grid is radial in the xy plane, and in the neighborhood of the horizontal boundaries, it becomes rectangular. Because the computational domain is symmetrical with respect to the axis y = 0, the problem was solved in a half of the domain only (y ≥ 0), and the homogeneous Dirichlet boundary conditions were applied at the line of symmetry. As it was mentioned, the total field was divided into the normal and anomalous fields, and these components were computed separately. The medium for the normal component consisted of the pipe and the well penetrating the homogeneous formation (ρ1 = 15 Ohm⋅m). The normal field was computed with the relative error of 0.1%. This estimate was made on the basis of comparing two solutions obtained at two nested grids. The relative error of computing the anomalous field, on the average, was about 2%. It should be noted that the anomalous component does not exceed 10% of the total field. Thus, the error of computing the total field does not exceed
Fig. 5. Anomalous component of emf. 1, top layer only (ρ4 = ρ1, ρ3 = ρ1, h1 = 4 m); 2, bottom layer only (ρ4 = ρ1, ρ2 = ρ1, h2 = 9 m); 3, horizontal layers only (ρ4 = ρ1, h1 = 4 m, h2 = 9 m); 4, with all boundaries (ρ4 = 3.5 Ohm⋅m, ∆ = 11 m); 5, with all boundaries (ρ4 = 2.0 Ohm⋅m, ∆ = 11 m)
0.3% when the grid of 606 thousand finite elements is used. Without dividing the total field into the normal field and the anomalous field, this error reached 2% even though the grid comprised 1035 thousand finite elements. Time dependence of the modulus of the anomalous component of the electromotive force for different resistivities of layers ρ2, ρ3, and ρ4 are shown in Fig. 5. The time interval shown is the one where the curves are especially dissimilar. We consider anomalous fields, because the relative difference between the total fields in different models does not exceed 9%. The magenta curve corresponds to the model with the top layer only (ρ4 = ρ1, ρ3 = ρ1); green, with the bottom layer only; black, with horizontal layers (no vertical interlayering); red and blue are the curves corresponding to the models with all boundaries with ρ4 = 3.5 Ohm⋅m and 2.0 Ohm⋅m, respectively. The top horizontal layer gives a greater contribution to the signal than the bottom one, as it is closer to the source and has higher conductivity. All transient curves differ with respect to the moment of time when the sign change takes place. Figure 6 shows the ratio between the modules of the anomalous component of the electromotive force in the model with and without the vertical boundary, in the presence of the horizontal boundaries. Results are given for different distances
Fig. 6. Anomalous component of emf relative to the signal in model with no vertical boundary for different distances to the vertical boundary, ρ4 = 2.0 Ohm⋅m.
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
729
Conclusions A computational scheme for simulation of the transient electromagnetic field in the 3D medium with highly conductive objects present has been developed. Adequacy of the scheme developed and advantages of separate computation of normal and anomalous fields have been shown. Examples of computations have been presented. Sensitivity of the signal measured to the horizontal boundaries and the boundary ahead of the tool has been analyzed. The scheme developed may be used when designing tools for geosteering.
Fig. 7. Anomalous component of emf relative to the signal in model with no horizontal boundary for different distances to the horizontal boundaries, ρ4 = 2.0 Ohm⋅m.
∆. Two time intervals can be set apart, based on the curve behavior: less and greater than 10–5 s. For the first interval, the relative value of the electromotive force reaches its maximum, which corresponds to the sign switch of its anomalous component. For the second interval, the electromotive force in the model with a vertical boundary differs from the corresponding model with no vertical boundary by the factor of 1.5 at ∆ = 11 m, by the factor of 1.2 at ∆ = 15 m, and by the factor of less than 1.1 at ∆ = 20 m. That is, sensitivity to the boundary is observed when it is located 11–15 m away from the source. The vertical boundary at 20 m does not seem to have any effect on the signal. Thus, there are limitations on determining the distance to a remote boundary ahead of the tool. Let us consider sensitivity of the signal to the horizontal boundaries. Figure 7 shows the ratio between the modules of the anomalous component of the electromotive force in the model with and without the horizontal boundaries, in the presence of the vertical boundary. As one can see, the red and black curves, corresponding to a fixed distance to the top boundary h1 and different distances to the bottom boundary h2, agree quite well. This means that sensitivity to the lower boundary is weak. At the same time, the tool retains sensitivity to the top boundary located at h1 = 20 m.
References Anderson, B., Chew, W.C., 1989. Transient response of some borehole mandrel tools. Geophys. 54 (2), 216–224. Bespalov, A., 2002. FEMAX—Software for simulation of magnetic induction tools in vertical wells. SEG Expanded Abstracts 21, 708–711. Epov, M.I., Morozova G.M., Antonov, E.Yu., 2002. Electromagnetic Method of Finding Faults in Casing of Oil and Gas Wellbores (Theoretical Basis and Method) [in Russian]. NITs OIGGM, SB RAS Publishing House, Novosibirsk. Kaufman, A.A., Sokolov, V.P., 1972. Theory of Inductive Logging Using Transient Electromagnetic Method [in Russian]. Nauka, Siberian Branch, Novosibirsk. Marchuk, G.I., 1980. Methods of Computational Mathematics [in Russian]. Nauka, Moscow. Nedelec, J.C., 1980. Mixed finite elements in R3. Numer. Math. 35, 315–341. Onegova, E.V., 2010. Effect of multicoil electromagnetic tool eccentricity on measured signals. Russian Geology and Geophysics (Geologiya i Geofizika) 51 (4), 423–427 (540–545). Plyusnin, M.I., Vilge, B.I., 1969. Foundations of inductive logging using transient electromagnetic method. Izvestiya Vysshikh Uchebnykh Zavedenii. Geologiya i Razvedka, No. 5, 158–165. Potapov, A.P., Kneller, L.E., 2000. Mathematical modeling and interpretation of borehole pulsed electromagnetic wall thickness measurements. Geofizika, No. 5, 27–30. Pryor, R., 2009. Multiphysics Modeling Using COMSOL: A First Principles Approach. Jones & Bartlett Learning, Sudbury. Sidorov, V.A., 1996. Borehole fault analyzers and wall thickness measurement tools for logging multistring boreholes. NTV Karotazhnik, No. 24, pp. 83–94. Tabarovsky, L.A., Sokolov, V.P., 1982. Code for computing transient field of dipole sources in horizontally layered medium (ALEKS), in: Elektromagnetic Methods in Geophysics [in Russian]. IGiG SO AN SSSR, Novosibirsk, pp. 57–77.
Editorial responsibility: A.D. Duchkov
3D simulation of transient electromagnetic field for geosteering horizontal wells E.V. Onegova a,*, M.I. Epov a,b a
b
Baker Hughes Russian Science Center, ul. Kutateladze 4a, Novosibirsk, 630090, Russia A.A. Trofimuk Institute of Petroleum Geology and Geophysics, Siberian Branch of the Russian Academy of Sciences, pr. Akademika Koptyuga 3, Novosibirsk, 630090, Russia Received 13 September 2010; accepted 7 December 2010
Abstract The article discusses numerical simulation of the transient electromagnetic field. The field source is an induction coil. We consider the situation when a logging tool is in a horizontal well in a medium with horizontal and vertical boundaries. The specific features of this problem are the metallic mandrel of the tool, 3D geometry, and distant boundaries. The method of separate computation of normal and anomalous fields is proposed. The finite element method is used for spatial approximation of the field, and the implicit finite-difference scheme is used for approximation of the field in time. Correctness and advantages of the method are shown. Some numerical results are demonstrated. The method proposed can be used when designing tools for geosteering. © 2011, V.S. Sobolev IGM, Siberian Branch of the RAS. Published by Elsevier B.V. All rights reserved. Keywords: geosteering; transient electromagnetic method; finite element method
Introduction As more deviated and horizontal wells are drilled today, geosteering has become a relevant issue. Geosteering is used for correcting the well trajectory based on the measurements performed while drilling, the goal being to achieve the optimal trajectory and thus, to maximize well productivity. When drilling deviated and horizontal wells, it is essential to determine the distances to the top and base of the reservoir, and to its inner boundaries (water-oil and gas-oil contact, clay layers), as early as possible. In this case, the drillstring can be redirected in time to avoid opening the water-bearing bed or the gas cap. The deeper the logging while drilling for the data one may obtain, the more efficient are the solutions which could be available. It is well known that in the transient electromagnetic method (TEM), sensitivity to remote areas of the formation tends to enlarge with the increase of the signal registration time, while closely-located objects seem to lose their influence. The main limitation for using TEM while drilling is the highly-conductive drillstring, as its signal can exceed that of the medium by several orders of magnitude.
First publications on inductive logging using TEM appeared in the 1970s (Kaufman and Sokolov, 1972; Plyusnin and Vilge, 1969). These papers considered the signal from a vertical magnetic dipole in the homogeneous medium (Plyusnin and Vilge, 1969), on the well axis and in beds of finite thickness (Kaufman and Sokolov, 1972). In (Anderson and Chew, 1989), numerical modeling of a well tool with an insulating mandrel in an axially symmetrical medium was described. A number of articles focused on TEM-based well flaw detection (Epov et al., 2002; Potapov and Kneller, 2000; Sidorov, 1996). In this case, the medium can be adequately described by a 2D model, calculations here must take into account both electric conductivity and magnetic permeability of the casing. The paper offers numerical schemes of modeling the situation typical of geosteering. The finite element method is used for modeling the transient electric field in the medium penetrated by the horizontal well, with horizontal layers and vertical boundary, taking into account electromagnetic properties of the drillstring (Fig. 1).
* Corresponding author. E-mail address: [email protected] (E.V. Onegova) 1068-7971/$ - see front matter D 201 1, V . S. S o bolev IGM, Siberian Branch of the RAS. Published by Elsevier B.V. All rights reserved. doi:10.1016/j.rgg.2011.06+.005
726
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
that the medium corresponding to the normal field can contain any axisymmetrical objects, including parts of the tool. The normal field E0 can be described by the following equation: ∂E0 1 ∂J =− , rot rot E0 + σ0 µ ∂t ∂t
Fig. 1. Model of the medium. 1, reservoir cap; 2, drillstring with transmitter and receiver coils; 3, well; 4, reservoir; 5, underlying medium; 6, vertical boundary.
Mathematical model Without taking into account displacement currents, the transient electric field can be described by the following equation: 1 ∂E ∂J =− , rot rot E + σ µ ∂t ∂t
(1)
where E is the vector of the electric field, µ is magnetic permeability, σ is electric conductivity, t is time, and J is density of the extraneous current. Let us set the z axis of the Cartesian coordinates in the line of well (Fig. 1). Because a coil coaxial with the pipe is used as the source, density of the extraneous current J will have only one nonzero component, Jϕ. Changes in density of the extraneous current with time are described by the Heaviside step function. In accordance to the Faraday law of induction, after the current in the source is turned off, eddy currents in the medium tend to be located next to the source. Then, in the process of evolution, the currents begin to diffuse into the medium. Before the currents reach horizontal boundaries, the electric field remains axisymmetrical just as in the coil, and has only one nonzero component, Eϕ. Let us use this fact and, for early times, solve a 2D problem in cylindrical (r, z) coordinates, rather than 3D. The presence of the metal pipe next to the source, as well as distant boundaries whose influence becomes noticeable only at later times, make the problem quite complex and difficult to compute. Because the signal from the currents in the pipe dominates in the total signal, let us divide the total field into the normal field due to the pipe in the axisymmetrical domain E0 = (0, Eϕ (r, z), 0)T and the anomalous field due to the rest of the medium E+ (Onegova, 2010; Tabarovsky and Sokolov, 1982). These fields will be computed separately. In this case, fine grids will be needed only in the 2D problem, and thus, computation time will be greatly reduced. It should be noted
(2)
where σ0 is electric conductivity of the axisymmetrical medium. Therefore, σ0 differs from σ only in those parts of the computational domain that are not axisymmetrical. For the model shown in Fig. 1, σ0 differs from σ in the upper and lower horizontal layers and in the layer ahead the well, because electric conductivity σ0 coincides with that of the middle layer. In this case, the normal field is the field of the pipe in the well in the homogeneous medium. Subtracting Eq. 2 from Eq. 1, we arrive at the equation for the anomalous field E+: ∂E+ ∂E0 1 = (σ0 − σ) rot rot E+ + σ . ∂t ∂t µ
(3)
Thus, the computation process can be reduced to the following: 1) From zero to a certain value of t0, we will be solving the axisymmetrical problem. 2) Starting with the moment t0, we will be each time layer be solving the axisymmetrical problem and the 3D problem. In both problems, we have homogeneous initial conditions: = 0, E0 t = 0
(4)
= 0. E+ t = t0
(5)
Approximate homogeneous Dirichlet boundary conditions are used here: = 0, E0 × n ∂Ω0
(6)
= 0. E+ × n ∂Ω
(7)
where ∂Ω0, ∂Ω are the boundaries of the axisymmetrical and 3D computation problems, n denotes the external normal to the corresponding boundary.
Method of solution To solve initial boundary problems (2), (4), and (6) and (3), (5), and (7) numerically, the finite elements method is used here. The variational formulation of the 2D problem has the following form: 1
∫ µ rot E0 ⋅ rot ψ dΩ + ∫ σ0
Ω0
Ω0
∂E0 ∂t
⋅ ψ dΩ =
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
−∫
Ω0
727
∂J ⋅ ψ dΩ ∀ ψ ∈ H0 (rot, Ω0), ∂t
where H0 (rot, Ω0) = υ ∈ (L2(Ω0))3 rot υ ∈ (L2(Ω0))3, υ × n|∂Ω = 0
0
is the space of test functions, ψ = (0, ψϕ, 0)T. The variational equation for the 3D problem has the following form: 1
∫ µ rot E+ ⋅ rot ψ dΩ + ∫ σ
Ω
Ω
(σ0 − σ) ∫
Ω
∂E0 ∂t
∂E+ ∂t
⋅ ψ dΩ =
Fig. 2. The electromotive force in the homogeneous medium. 1, with pipe; 2, no pipe.
⋅ ψ dΩ ∀ ψ ∈ H0 (rot, Ω),
where ψ = (ψx, ψy, ψz)T. To approximate the time derivative of the solution, the completely implicit three-layer scheme is used (Marchuk 1980) ∂E ≈ γ0 Ej − γ1 Ej−1 + γ2 Ej−2, ∂t t = t j
where γ0 =
∆01 + ∆02
, γ1 =
∆02
, γ12 =
∆01
, ∆12 = ∆01 ∆02 ∆12 ∆01 ∆12 ∆02 tj−1 − tj−2, ∆02 = tj − tj−2, ∆01 = tj − tj−1, ti is the time layer, and Ei is the solution at the ith layer. For the finite-element discretization, rectangles with piecewise-bilinear basis functions are used for the 2D problem, and hexahedrons with edge basis functions of the first order, for the 3D problem (Nedelec, 1980).
Testing The results of 3D modeling using the bedded medium were compared to known solutions (Tabarovsky and Sokolov, 1982). The magnetic dipole was used as the field source. In the 3D formulation, a coil with the radius of 0.01 m was used as the field source. Relative discrepancy between solutions did not exceed 1% in the time range between 10–7 and 10–2 s. Besides, the results of solving the 2D problem were compared with computations obtained using FEMAX (Bespalov, 2002) and COMSOL (Pryor, 2009). The transient process in the homogeneous medium with a metal pipe of a finite length was modeled. A coil coaxial with the pipe was used as the field source. Pipes with diverse electric conductivity and magnetic permeability (magnetic steel, nonmagnetic steel, and copper) and diverse containing media were considered. Average relative discrepancy between solutions obtained using three different methods was 2–3% in the time range between 10–7 and 10–2 s. To understand the degree of the impact the pipe might have upon the signal and show the possibilities of the computational
scheme developed, let us consider some computations. Figure 2 shows the electromotive force induced in the receiver coil as a function of time. Two results are presented: with pipe and without pipe. Formation resistivity was 100 Ohm⋅m; steel pipe resistivity was 7.14 × 10–7 Ohm⋅m; its magnetic permeability was 100µ0 (µ0 = 4π × 10−7H/m); the outer and inner radii of the pipe were 0.04 and 0.07 m, respectively; coil radius was 0.085 m; the distance between the transmitter coil and receiver coil was 5 m; the electric current was 5 A; the number of loops in the transmitter coil was 100; the number of loops in the receiver coil was 10. In the graph, the pipe affects the transient process greatly, both quantitatively and qualitatively. In the range between 5 × 10−4 and 5 × 10−2 s, the signal from the pipe decreases quite slowly, and one can even observe some increase at the time of 2 × 10−2 s. After the time of 0.1 s, the signal from the pipe decreases greatly, and after the time of 2 s, it attenuates as t−5 / 2. Physically, this means there are no more currents in the pipe, and the other currents are far from the initial source that the situation becomes equivalent to that of a magnetic dipole in the homogeneous medium. On the other hand, the signal level at these times is much lower than the measured values. Generally, it should be noted that the computational scheme can model the signal in the dynamic range around 200 dB.
Numerical experiments Let us consider a geoelectrical model typical of the west Siberian oil and gas province. The horizontal interval of the well is in the oil-saturated reservoir (ρ1 = 15 Ohm⋅m, Fig. 3). The reservoir has a clay cap (ρ2 = 4 Ohm⋅m). Under the reservoir, there is a water-saturated interval (ρ3 = 8 Ohm⋅m). In front of the tool, orthogonally with respect to the well, there is a clay layer (ρ4 = 2.0–3.5 Ohm⋅m). The drilling mud with a high clay content has electric resistivity of 2 Ohm⋅m. Resistivity of the pipe made of nonmagnetic steel is 7.14 × 10−7 Ohm⋅m. The outer and inner radii of the pipe were 0.04 and 0.07 m, respectively; coil radius was 0.085 m for
728
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
Fig. 3. Problem domain and model parameters. Fig. 4. Cross section via the xy plane with the hexahedron grid.
both transmitter coil (T) and receiver coil (R); the well radius was 0.108 m. The electric current in the transmitter coil was 5 A; the number of loops in the transmitter coil was 100; the number of loops in the receiver coil was 10. The object of measurements was the electromotive force induced in the receiver coil. Figure 4 shows a grid of hexahedrons for the computational domain. Inside the well, the grid is radial in the xy plane, and in the neighborhood of the horizontal boundaries, it becomes rectangular. Because the computational domain is symmetrical with respect to the axis y = 0, the problem was solved in a half of the domain only (y ≥ 0), and the homogeneous Dirichlet boundary conditions were applied at the line of symmetry. As it was mentioned, the total field was divided into the normal and anomalous fields, and these components were computed separately. The medium for the normal component consisted of the pipe and the well penetrating the homogeneous formation (ρ1 = 15 Ohm⋅m). The normal field was computed with the relative error of 0.1%. This estimate was made on the basis of comparing two solutions obtained at two nested grids. The relative error of computing the anomalous field, on the average, was about 2%. It should be noted that the anomalous component does not exceed 10% of the total field. Thus, the error of computing the total field does not exceed
Fig. 5. Anomalous component of emf. 1, top layer only (ρ4 = ρ1, ρ3 = ρ1, h1 = 4 m); 2, bottom layer only (ρ4 = ρ1, ρ2 = ρ1, h2 = 9 m); 3, horizontal layers only (ρ4 = ρ1, h1 = 4 m, h2 = 9 m); 4, with all boundaries (ρ4 = 3.5 Ohm⋅m, ∆ = 11 m); 5, with all boundaries (ρ4 = 2.0 Ohm⋅m, ∆ = 11 m)
0.3% when the grid of 606 thousand finite elements is used. Without dividing the total field into the normal field and the anomalous field, this error reached 2% even though the grid comprised 1035 thousand finite elements. Time dependence of the modulus of the anomalous component of the electromotive force for different resistivities of layers ρ2, ρ3, and ρ4 are shown in Fig. 5. The time interval shown is the one where the curves are especially dissimilar. We consider anomalous fields, because the relative difference between the total fields in different models does not exceed 9%. The magenta curve corresponds to the model with the top layer only (ρ4 = ρ1, ρ3 = ρ1); green, with the bottom layer only; black, with horizontal layers (no vertical interlayering); red and blue are the curves corresponding to the models with all boundaries with ρ4 = 3.5 Ohm⋅m and 2.0 Ohm⋅m, respectively. The top horizontal layer gives a greater contribution to the signal than the bottom one, as it is closer to the source and has higher conductivity. All transient curves differ with respect to the moment of time when the sign change takes place. Figure 6 shows the ratio between the modules of the anomalous component of the electromotive force in the model with and without the vertical boundary, in the presence of the horizontal boundaries. Results are given for different distances
Fig. 6. Anomalous component of emf relative to the signal in model with no vertical boundary for different distances to the vertical boundary, ρ4 = 2.0 Ohm⋅m.
E.V. Onegova and M.I. Epov / Russian Geology and Geophysics 52 (2011) 725–729
729
Conclusions A computational scheme for simulation of the transient electromagnetic field in the 3D medium with highly conductive objects present has been developed. Adequacy of the scheme developed and advantages of separate computation of normal and anomalous fields have been shown. Examples of computations have been presented. Sensitivity of the signal measured to the horizontal boundaries and the boundary ahead of the tool has been analyzed. The scheme developed may be used when designing tools for geosteering.
Fig. 7. Anomalous component of emf relative to the signal in model with no horizontal boundary for different distances to the horizontal boundaries, ρ4 = 2.0 Ohm⋅m.
∆. Two time intervals can be set apart, based on the curve behavior: less and greater than 10–5 s. For the first interval, the relative value of the electromotive force reaches its maximum, which corresponds to the sign switch of its anomalous component. For the second interval, the electromotive force in the model with a vertical boundary differs from the corresponding model with no vertical boundary by the factor of 1.5 at ∆ = 11 m, by the factor of 1.2 at ∆ = 15 m, and by the factor of less than 1.1 at ∆ = 20 m. That is, sensitivity to the boundary is observed when it is located 11–15 m away from the source. The vertical boundary at 20 m does not seem to have any effect on the signal. Thus, there are limitations on determining the distance to a remote boundary ahead of the tool. Let us consider sensitivity of the signal to the horizontal boundaries. Figure 7 shows the ratio between the modules of the anomalous component of the electromotive force in the model with and without the horizontal boundaries, in the presence of the vertical boundary. As one can see, the red and black curves, corresponding to a fixed distance to the top boundary h1 and different distances to the bottom boundary h2, agree quite well. This means that sensitivity to the lower boundary is weak. At the same time, the tool retains sensitivity to the top boundary located at h1 = 20 m.
References Anderson, B., Chew, W.C., 1989. Transient response of some borehole mandrel tools. Geophys. 54 (2), 216–224. Bespalov, A., 2002. FEMAX—Software for simulation of magnetic induction tools in vertical wells. SEG Expanded Abstracts 21, 708–711. Epov, M.I., Morozova G.M., Antonov, E.Yu., 2002. Electromagnetic Method of Finding Faults in Casing of Oil and Gas Wellbores (Theoretical Basis and Method) [in Russian]. NITs OIGGM, SB RAS Publishing House, Novosibirsk. Kaufman, A.A., Sokolov, V.P., 1972. Theory of Inductive Logging Using Transient Electromagnetic Method [in Russian]. Nauka, Siberian Branch, Novosibirsk. Marchuk, G.I., 1980. Methods of Computational Mathematics [in Russian]. Nauka, Moscow. Nedelec, J.C., 1980. Mixed finite elements in R3. Numer. Math. 35, 315–341. Onegova, E.V., 2010. Effect of multicoil electromagnetic tool eccentricity on measured signals. Russian Geology and Geophysics (Geologiya i Geofizika) 51 (4), 423–427 (540–545). Plyusnin, M.I., Vilge, B.I., 1969. Foundations of inductive logging using transient electromagnetic method. Izvestiya Vysshikh Uchebnykh Zavedenii. Geologiya i Razvedka, No. 5, 158–165. Potapov, A.P., Kneller, L.E., 2000. Mathematical modeling and interpretation of borehole pulsed electromagnetic wall thickness measurements. Geofizika, No. 5, 27–30. Pryor, R., 2009. Multiphysics Modeling Using COMSOL: A First Principles Approach. Jones & Bartlett Learning, Sudbury. Sidorov, V.A., 1996. Borehole fault analyzers and wall thickness measurement tools for logging multistring boreholes. NTV Karotazhnik, No. 24, pp. 83–94. Tabarovsky, L.A., Sokolov, V.P., 1982. Code for computing transient field of dipole sources in horizontally layered medium (ALEKS), in: Elektromagnetic Methods in Geophysics [in Russian]. IGiG SO AN SSSR, Novosibirsk, pp. 57–77.
Editorial responsibility: A.D. Duchkov