Flow towards a guarded sampling probe: Modelling of a 2D flow cell

Journal of Petroleum Science and Engineering 66 (2009) 133–142

Contents lists available at ScienceDirect

Journal of Petroleum Science and Engineering j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / p e t r o l

Flow towards a guarded sampling probe: Modelling of a 2D flow cell J.D. Sherwood Schlumberger Cambridge Research, High Cross, Madingley Road, Cambridge CB3 0EL, UK

a r t i c l e

i n f o

Article history: Received 11 January 2007 Accepted 22 February 2009 Keywords: fluid sampling Darcy flow

a b s t r a c t Analytic and numerical methods are used to predict fluid flow in a thin rectangular sheet of porous material when fluid is extracted through a port at the midpoint of one of the edges of the sheet. The aim is to model (in 2 dimensions) flow towards a sampling probe used to extract reservoir fluids from the rock surrounding a newly-drilled well. Such samples are usually contaminated with drilling fluid filtrate that saturates the rock adjacent to the wellbore, but a guarded probe can eliminate this contamination. Complex variable techniques are used to determine the flow of incompressible fluids when the interface between the two fluids is sharp, and filtrate and original pore fluid have identical viscosities. Results are compared against the predictions of a numerical reservoir simulator, and agreement is reasonable once numerical dispersion has been accounted for. Numerical results are also obtained for the case when the viscosities of the two fluids differ. The results give confidence in the use of 3-dimensional numerical computations to predict flow towards a guarded probe in a real wellbore, and thereby to optimize the design of such probes. © 2009 Elsevier B.V. All rights reserved.

1. Introduction Downhole sampling of hydrocarbon reservoir fluids immediately after a well has been drilled has been radically improved by the introduction of guarded probes. Such probes markedly enhance the quality of the fluid sample, reducing the level of contamination by drilling fluid filtrate that has entered the rock pores during drilling (Hrametz et al., 2001; Sherwood et al., 2004; Schlumberger, 2006). The guarded probe is divided into two regions: an inner sampling probe, surrounded by an outer annular guard probe. To obtain a sample of pore fluid, the probe is pushed against the rock surface surrounding the wellbore, and fluid is pumped through the probe from the rock into the sampling device. Fluid arriving at the outer, guard probe is always contaminated by drilling fluid filtrate, and is discarded. Theory (Sherwood, 2005) predicts that the level of contamination in fluid arriving at the central sampling probe should drop to zero in finite time, and this is confirmed by field experience (Weinheber and Vasques, 2006; O'Keefe et al., 2006). Tarvin et al. (2008) discuss experiments that investigate 2dimensional flow towards a guarded downhole sampling probe. The experiments were compared to the predictions of a finite-difference numerical reservoir simulator. Here we show how simple analytic models can be used to understand both the experiments and the numerical predictions. The flow field in the 2-dimensional geometry considered here is much simpler than that in the real wellbore

E-mail address: [email protected]field.slb.com. 0920-4105/$ – see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.petrol.2009.02.002

geometry downhole. This enables us to confirm our understanding of the underlying physics before going on to model the more complicated downhole geometry. In Section 2 we consider a simple model for flow in an unbounded half-space. This is a plane 2-dimensional analogue of the axisymmetric analysis of Sherwood (2005), which was itself an extension of a model for unguarded probes due to Hammond (1991). In Sections 3–5 we show how the important lateral boundaries of the experimental cell may be introduced into the analysis by means of a transformation in the complex plane. Finally in Section 6 we compare analytic results with the predictions of a numerical reservoir simulator. 2. Flow towards a guarded probe in a half-space Fig. 1 shows the plane 2-dimensional geometry we are considering. The region x N 0 is porous rock, and x = 0 represents the (plane) wellbore wall. This is sealed by drilling fluid filtercake, which we regard as impermeable, except over the region |y| b a where the filter cake has been removed by a probe. The rock porosity is ϕ. Initially the rock pores in the region 0 ≤ x ≤ βa are saturated with drilling fluid filtrate of viscosity μf which has displaced original pore fluid, of viscosity μp, which now occupies the region x N βa. We assume that the interface between the two fluids is sharp and remains so when fluid flows within the rock. Thus we are neglecting dispersion, which tends to smear the interface both when filtrate invades the rock and when fluid is pumped towards the probe. We are also neglecting viscous fingering instabilities. The interface will be unstable when filtrate displaces original pore fluid if μf b μp, and will be unstable when fluid is extracted if μf N μp. From now on (with the exception of Section

134

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

In limit V → ∞ we find  1 = 2 2/ F fβa ; πV

ð4Þ

  and F(V) decays more slowly than the equivalent prediction F~V − 2 = 3 for the analogous 3-dimensional problem (Hammond, 1991). The fluid velocity u in the x direction at the probe x = 0 is given by Morse and Feshbach (1953, p. 1196) Q  1 = 2 ; π a2 − y2

uð0; yÞ =

ð5Þ

and the total flow through the probe over the region |y| b y1 is Z y 1 2Q −1 sin ðy1 = aÞ: Q1 ðy1 Þ = 2 uð0; yÞdy = π 0

Fig. 1. Schematic showing the motion of the interface between filtrate and original pore fluid, originally at x = x0 = βa. The fluid particle initially at (x0, y0) reaches the probe at (0, yp).

6.3) we set μf = μp = μ and neglect these effects. We assume that the fluid velocity u = (u, v) obeys Darcy's law k u = − jp μ

ð1Þ

where k is the hydraulic permeability and p the pressure. The fluids are assumed incompressible, so that ∇ · u = 0. Hence p satisfies the Laplace equation 2

j p = 0:

ð2Þ

Fluid is withdrawn by the probe at x = 0, where the region |y| b a is at uniform pressure. The interface between filtrate and original pore fluid is initially at x = βa, shown as the line ABC on Fig. 1. When fluid is withdrawn by the probe at x = 0, the interface moves towards, and eventually reaches, the probe. Fig. 1 shows the moment when fluid originally at B (βa, y0) reaches the probe at E (0, yp) and the interface has moved to ADE and CFG. The interface between filtrate and original pore fluid first reaches the probe on the axis of symmetry y = 0. The points (0, ±yp) at which the interface meets the probe subsequently move away from the centre of the probe, and when yp N αa the central sampling region of the probe, |y| b αa, collects only original pore fluid. At distances r = (x2 + y2)1/2 ≫ a the fluid velocity is approximately radial towards a point sink at the probe P, as shown on Fig. 2. If the total volumetric flow rate through the probe is Q (per unit distance in the z direction perpendicular to the (x, y) plane), with Q b 0 for a sink flow, then at time t the total volume of fluid withdrawn by the probe is Rt V ðt Þ = − 0 Q ðt VÞdt V, and thus fluid initially in the half-space x N 0 within the disc r b rw(t) = (2V / πϕ)1/2 will have been withdrawn. If rw b βa the half-disc r b rw lies entirely within filtrate, so that the fraction of filtrate in the fluid arriving at the probe is F = 1. At later times, the proportion of fluid arriving at the probe that is original pore fluid is 2θ/π, where θ = cos− 1(βa / rw) is shown on Fig. 2, so that F = 1; " #1 = 2; 2 2 2 − 1 β a π/ = sin π 2V

1 2 2 Vb /πβ a 2 2 2

V N 12/πβ a :

ð6Þ

The motion of the interface between filtrate and original pore fluid can be determined by integration of the velocity field around the probe. However, Sherwood (2005) showed that the far-field analysis based upon radial flow gives a good approximation of the fraction of original pore fluid entering the probe for the analogous axisymmetric geometry, and we expect that Eqs. (3) is similarly a good approximation in the plane problem considered here. If we now consider only the central, sampling region of the probe, |y| b αa, the total flow through this region is Q1 ðαaÞ =

2Q −1 sin ðα Þ; π

ð7Þ

of which 2(Q / π)cos− 1(βa /rw) is original pore fluid until the interface between filtrate and original pore fluid meets the probe at y1 =αa, after which the central sampling region collects only original pore fluid. The filtrate fraction of the flow into the central sampling region is therefore 1 2 2 Vb /πβ a 2

F = 1;

= 1−

cos − 1 ðβa = rw Þ ; sin − 1 α

/πβ2 a2 /πβ2 a2  bVb  2 2 1 − α2

ð8aÞ

ð8bÞ

2 2

= 0;

/πβ a  VN  2 1 − α2

ð8cÞ

where rw = (2V /π/)1/2. 3. The fluid velocity field in a bounded geometry The flow geometry used in the experiments of Tarvin et al. (2008) is shown in Fig. 3. The cell has width 2b and hence the total amount of filtrate in the region x b βa is finite. In the absence of boundaries at

ð3aÞ

ð3bÞ Fig. 2. Schematic representing the simple 2-dimensional radial flow analysis.

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

135

Fig. 4. The z plane. The cell is supposed infinitely long, so that the side DEF has gone off to infinity.

Thus Ψ is a stream function, and Φ= −

Fig. 3. Schematic of the plane flow cell. The interface between filtrate and original pore fluid is assumed sharp, and is initially at x = βa = x0.

y ± b, the analysis of Section 2 predicted F ~ V− 1/2 as V → ∞ (Eqs. (4)), RV so that the volume of filtrate collected, 0 F ðV VÞdV V, is unbounded, and cannot correctly represent the finite volume of filtrate in the experimental cell. The walls of the cell therefore play an important role, and the analysis of Section 2 fails when the volume V of fluid withdrawn is such that rw = (2V / πϕ)1/2 is comparable to b. We now use complex variable theory to solve the Laplace Eq. (2) and thereby determine the fluid velocity within the cell. Schwarz– Christoffel transformations of the complex plane might be used to handle the exact geometry of the cell. However, here we follow Dufrêche et al. (2002) and assume that the cell is long compared to its width, i.e. l ≫ b. Fig. 4 shows a schematic of the cell, which occupies the region ABCDEFGHA. The sides CBAHG and DEF are of length 2b. BH represents the (compound) probe of total length 2a through which fluid is withdrawn from the cell, and over which we assume that the pressure p = p0 = 0 is uniform. CD and FG are impermeable walls of length 2l; GH and BC are also impermeable and of length b–a. Fluid enters the cell through the side DEF, which is at uniform pressure p1. We assume that DEF is far from the probe (l ≫ b), and take the limit in which D and F go off to infinity. The fluid velocity u is given by Darcy's law (Eq. (1)) and the pore pressure p satisfies the Laplace Eq. (2). Standard theory (e.g. Morse and Feshbach, 1953) tells us that if G(z) = Φ + iΨ is an analytic function in the complex plane z = x + iy, then AΦ AW = ; Ax Ay

AΦ AW = − : Ay Ax

ð9Þ

dG u − iv = dz

ð10Þ

ð13Þ

is proportional to the pore pressure p. The flow geometry shown in Fig. 4 lies in the z plane, and we transform this to the w plane depicted in Fig. 5 by means of the transformation w = i sinhðπz = 2bÞ:

ð14Þ

The region of interest in the z plane is mapped into the upper half of the w plane. C and G go to w = b1, whereas B and H are mapped to b sinðaπ = 2bÞ = bλ, say. Far from the probe the flow in the w plane looks like that towards a point sink so that p ~ log w at infinity. We now transform to the w1 plane via   2 2 1=2 : ð15Þ w1 = w + w −λ The cut BAH in the w plane, which represents the probe, is mapped onto a circle of radius λ in the w1 plane. The region outside the cut BAH in the w plane is mapped onto the exterior of this circle in the w1 plane. A function which has real part zero (corresponding to p = 0) on the circle |w1| = λ with logarithmic behaviour at infinity is G=

Q lnðw1 = λÞ: π

ð16Þ

The velocity field is u − iv =

dG dG dw1 dw Q coshðπz = 2bÞ = =  1 = 2 : dz dw1 dw dz 2b λ2 + sinh2 ðπz =2bÞ

Q cosð πy = 2bÞ  1 = 2 ; 2b λ2 − sin2 ð πy=2bÞ

AΦ AW = ; Ax Ay

v=

AΦ AW = − Ay Ax

ð11Þ

with Au Av 2 + = j Φ = 0: Ax Ay

ð12Þ

ð18Þ

and we note that in the limit b → ∞ we recover the corresponding result (Eq. (5)) for flow within an unbounded half-space. There is an

where u=

ð17Þ

As z → ∞ we see that the fluid velocity tends to a uniform flow u ~ Q / (2b) where Q b 0 if, as here, the probe acts as a sink. Since the flow cell has width 2b the total volumetric flow rate is Q. Over the probe (x = 0,|y| ≤ a), the Darcy velocity is u=

Hence we may look for a fluid velocity field (u, v) such that

pμ k

Fig. 5. The w plane.

136

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

(integrable) singularity in the velocity at the edges y = ±a of the probe. The cumulative flow out of the probe over the region −y1 b y b y1 is Z Qc ðy1 Þ = 2

y1

u dy =

0

   2Q πy −1 1 sin sin 1 π λ 2b

ð19Þ

and we see that Qc(a) = Q, as required by volume conservation.

4.1. The pressure as x → ∞ The pressure field corresponding to Eq. (16) is, by Eq. (13), μ  ReðGÞ p= − k      zπ   zπ   1 = 2 Qμ = − Re ln i sinh + − sinh2 −λ2 − ln λ : πk 2b 2b ð20Þ As x → ∞ we obtain ð21Þ

which may be thought of as a uniform pressure gradient together with a pressure drop caused by converging flow into the probe. This pressure drop vanishes if a = b. We now investigate the errors created when we use the velocity field (Eq. (17)) which is exact in the limit l → ∞ rather than the velocity field proper to a bounded flow cell with finite length 2l. In the experimental cell the pressure p = p1 is uniform over the side x = 2l = b. For the approximate solution used here, at any given finite value of x there will be a pressure variation as y varies over (− b,b). We examine the variation of pressure along the lines x = 2b and x = b. 4.2. Pressure variation along x = 2b

ð22Þ

which gives p(2b, 0)=−Qμ/k as expected when λ=1, i.e. when a=b and the probe occupies the entire width of the cell. If λ≪sinh π=11.55, then #  " 2 Qμ λ ::: ; lnð2 sinh πÞ − lnλ + + pð2b; 0Þ = − πk 4 sinh2 π

ð23Þ

where ln(2 sinh(π))=3.139. At the cell walls we find pð2b; − bÞ = −

     1 Qμ 2 2 Re ln cosh π + cosh π − λ 2 − ln λ πk ð24Þ

with p(2b, −b)=−Qμ/k in the limit a=b, and #  " Qμ λ2 : : : lnð2 cosh π Þ − ln λ − + pð2b; − bÞ≈ − πk 4 cosh2 π

5. Motion of the filtrate interface We now integrate the velocity field obtained in Section 3 in order to follow the motion of points on the interface between filtrate and original pore fluid, which is situated at x = x0 = βa at time t = 0. In the absence of forces (e.g. capillarity or gravity) other than the applied pressure gradient, the fluid velocity is everywhere proportional to the volumetric flow rate Q, and the fluid displacement depends solely upon the volume of fluid withdrawn, whether or not Q(t) is constant. 5.1. Streamlines close to the centreline, y ≪ a We first consider the streamlines close to the axis y = 0. These are important since the sampling probe is likely to occupy only a small portion |y|≤αa of the total probe, with α ≪ 1. We shall later make the additional assumption that a ≪ x0. The velocity close to the axis of symmetry y = 0 may be obtained by expanding Eq. (17) as   2    Q S 1 − λ πy 3 +O y   4b2 λ2 + S2 3 = 2

 v=

ð26Þ

( " #)  Q C π2 y2 C 2 − λ2 + 2S2 3C 2 S2 u= 1+ − 2   2b λ2 + S2 1 = 2 8b2 λ2 + S 2 λ2 + S 2   4 +O y ; ð27Þ 

On the centreline, at (x, y) = (2b, 0) the pressure is predicted by Eq. (20) to be     Qμ − 1 sinh π pð2b; 0Þ = − sinh πk λ

4.3. Pressure variation along x = b The pressure p(b, 0) is given by replacing sinh π by sinh (π/2) = 2.30 in Eq. (23) whereas p(b, −b) is obtained by replacing cosh π (Eq. (25)) by cosh (π/2) = 2.51, where ln(2 sinh (π/2)) = 1.53 and ln(2 cosh (π/2)) = 1.61. The errors in the pressure across the halfwidth of the cell are of order 6 % and are still sufficiently small to allow us to use Eq. (17) as the velocity field.

4. The pressure field

 h aπi Q μ xπ pf − − ln sin ; πk 2b 2b

therefore closely approximated by Eq. (17), although this is strictly only valid in the limit in which l ≫ b. Note that to evaluate p(2b, b) = p(2b, − b) we must be careful to choose the correct branch of the square root in Eq. (20).

ð25Þ

when λ ≪ cosh π = 11.59 and ln(2 cosh π) = 3.143. Comparing Eqs. (23) and (25) we see that the difference in pressure across the half-width of the cell, p(2b, b) − p(2b, 0), is negligibly small when x = 2b. The velocity field within a square cell, with l = b, is

where C = coshðπx = 2bÞ;

S = sinhðπx = 2bÞ:

ð28Þ

Note that the term O(y2) in the expansion for u is zero when λ = 1, i.e. when a = b and there is uniform flow along streamlines y constant. Streamlines are given (when y ≪ b) by   S 1 − λ2 πy dy dy dt   = = dx dt dx 2bC λ2 + S2

ð29Þ

which integrates to give 2

y = y20

2

2

λ +S 1 + S2

!

2

1 + S0 λ2 + S20

! ð30Þ

where y = y0 when x = x0 and S = S0 = sinh (πx0 / 2b). Thus the streamline that starts at (x0, y0) reaches the probe on the wall x = 0 at 1 + S20 yp = y0 λ 2 λ + S20 ≈ y0 λ

!1 = 2

when

ð31Þ S0 ≫1:

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

If the porosity of the porous medium is ϕ, the velocity of a fluid particle is larger than the Darcy velocity by a factor ϕ− 1. The time t0 taken for a particle initially at (x0, y0) on the interface between filtrate and original pore fluid to reach the probe at x = 0 may be found by integrating Eq. (27) along a streamline Z 0

t0

Q ðt Þ dt = 2b/

Z

0

 1 = 2 λ2 + S 2

dx −

π2 y20 8b2

1 + S20 λ2 + S20

!

In the next section we shall discuss the limit b ≫a, for which Eq. (38) predicts Vα =

 πx20 /  2 1 + α +: : : : 2

ð39Þ

5.2. Large cell: a2 ≪ V ≪ b2

C   # 2 2 3=2 " 2 0 λ + S C − λ2 + 2S2 3C 2 S2 − dx:    2 C3 λ2 + S2 x0 λ2 + S2 x0

137

Z

If the radius rw from which fluid has been withdraw lies in the range a ≪ rw ≪ b, the approximate analysis of Section 2 holds. We see from Eq. (8c) that we expect the sampling probe to be free of filtrate after a time

ð32Þ The left-hand side of Eq. (32) is simply − V0 / (2bϕ), where V0 is the volume of fluid that must be withdrawn for pore fluid initially at (x0, y0) to reach the probe. The first term on the right-hand side of Eq. (32) gives the volume Vf at which the interface first reaches the probe along a streamline y = 0. It may be evaluated exactly, and we find   2  4b λ − 1 ðz − z1 Þð1 − z2 Þ Vf 2b −1 ln ð33Þ sinh ðS0 = λÞ + = ðz − z2 Þð1 − z1 Þ π 2b/ πλ2 ðz1 − z2 Þ where   −1 z = exp 2 sinh ðS0 = λÞ   −2 −4 −2 1=2 + 2 λ −λ z1 = 1 − 2λ   −2 −4 −2 1=2 z2 = 1 − 2λ − 2 λ −λ :

Vα =

/πx20   2 1 − α2

ð40Þ

and this agrees with Eq. (39) to O(α2). Ultimately, however, the influence of the walls at y = ±b will be felt at the probe. In particular, the total amount of filtrate within the cell is finite, as will be the integrated amount of filtrate arriving at the combined guard and sampling probes. This is not so in an unbounded geometry, where the contamination decays as V− 1/2 (Eq. (4)). Thus the finite size of the cell has an important effect on the contamination at long time. However, the stagnation flow in the corners at (0, ±b) prevents the contamination in the guard probe from dropping to zero in finite time. During sampling of pore fluid from the rock surrounding a real wellbore, any layers of impermeable shale can play a role similar to that of the cell boundaries in the laboratory experiment discussed here. 5.3. Numerical integration

The second term on the right-hand side of Eq. (32) is more difficult to integrate. We assume λ ≪ S0 (i.e. a ≪ x0), so that Eq. (32) becomes V0 = 2b/

Z

x0

0

(

S π2 y20 S + C 4b2 C 3

1 + S20 S20

!)

  2 dx + O λ

ð34Þ

where O(λ2) errors come from the integral over the region near x = 0 where S is no longer large compared to λ. Integrating Eq. (34) we find 2 2 2   π yp S0 πV0 2 + +O λ ; = ln C 0 2 2 2 2 4b / 8C0 λ b

ð35Þ

where C0 = cosh(πx0 / 2b), so that C20 − S20 = 1, and yp is given by Eq. (31). However, near the centreline y = 0 it is clearly more accurate to use the exact result (Eq. (33)) for Vf, leading to

V0 = Vf +

πy2p S20 / 2C02 λ2

  2 +O λ :

ð36Þ

We now integrate the fluid velocity (Eq. (17)) numerically in order to determine the displacement of the interface. We scale distances by a and pumped volumes V by Vs = a2ϕ. Fig. 6 shows the position of an interface, initially at x0 = 5a in a flow cell of half-width b = 10a. Only half of the cell is shown. The interface is shown at a series of non-dimensional pumped volumes Vˆ = V / Vs = 0, [1], 9. Fig. 7 shows the pumped volume Vˆa at which a sampling probe occupying the region − αa b y b αa will be free of filtrate. Curve (a) gives the result of full numerical calculations, whereas curve (b) gives the radial flow approximation (Eq. (40)) of Section 5.2 and (c) gives the approximation (Eq. (37)) of Section 5.1. The approximations are only valid for small α. At large V the boundaries play an important role and the radial flow approximation breaks down. The radial flow approximation predicts a velocity that is too low far from the probe, and too high close to the probe. We see that on the centre-line these errors approximately cancel and the time at which original pore fluid first reaches the probe is correctly predicted.

Now suppose the sampling region of the probe occupies the region |y| b αa. In the absence of dispersion, original pore fluid will first reach the probe at its centre. The sampling probe will be completely free of filtrate when the interface intersects the probe at (0, ±αa): Fig. 1 shows a sketch of the geometry of the filtrate interface and probe. We conclude, following Sherwood (2005), that the sampling probe will be free of filtrate when the volume of fluid withdrawn is b2 /πα 2 S20 2C02

ð37Þ

! 4b2 / α 2 S20 ::: : ln C0 + + π 2C02

ð38Þ

Vα = Vf +

≈

Fig. 6. Position of the interface at V̂ = V / Vs = 0,[1],9, with flow from right to left. Interface initially at x0 = 5a, flow cell half-width b = 10a.

138

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

Fig. 7. Non-dimensional pumped volume Vα̂ = Vα / Vs at which a sampling probe of size 2αa becomes filtrate-free. (a) Full numerical computation, (b) radial flow approximation (Eq. (40)), (c) expansion (Eq. (37)) for α ≪ 1, x0 = 5a, b = 10a.

Fig. 8 shows similar results for Vˆ α for the case when the bounding box is moved further away, to b = 100a, while keeping x0 = 5. The radial flow approximation (curve b) is now good at late time, but is poorer for small α. As before, the radial flow approximation assumes a velocity that is too high close to the probe. However, now that the boundary is further away, the radial approximation is better at large distances and under-estimates the far-field velocity less than before. Consequently the first arrival time of original pore fluid is now underestimated. Fig. 9 shows the filtrate fraction F in (a) the combined guard + sampling probes, and in sampling probes with (b) α = 0.2, (c) α = 0.4, (d) α = 0.6. The cell half-width is b = 10a and the interface is initially at x0 = 5a, corresponding to the case considered in Figs. 6 and 7. Note from Eq. (39) that once original pore fluid has arrived at the probe at pumped volume Vf, the additional volume required for the sampling probe to become filtrate-free is small when α ≪ 1 and is approximately α2Vf. 6. Numerical reservoir simulations 6.1. Geometry The complex variable methods discussed above cannot be used to predict 3-dimensional flow around a wellbore. An obvious method to investigate such flow is by means of a numerical reservoir simulator. We shall not present such a 3-dimensional investigation here, but instead we discuss the ability of numerical reservoir simulators to predict the 2-dimensional flow considered in Sections 1–5. The ECLIPSE⁎ simulations were performed by the ECLIPSE FrontSim⁎ streamline simulator (Schlumberger, 2007), using a Cartesian coordinate system, with −25 m b y b 25 m and 0 b x b 25 m. In the simulations the Cartesian grid divided the (x, y) plane into (200 × 400) cells, of which the region (133 × 266) closest to the point of withdrawal had size dx = dy = 0.075 m. Grid sizes (dx, dy) outside this region increased in geometric progression with the ratios dxi + 1/dxi and dyi + 1/dyi constant. The gridblocks in the corners furthest from the probe had dimension dx = dy = 0.49 m. In the z direction, the distance between the plates was represented by a single grid cell of dimension dz = 10 m. The porosity was taken to be uniform, with ϕ = 0.4. Permeabilities were assumed uniform and isotropic for the cases considered here. The filtrate was represented by oil, and initial pore fluid by water. No attempt was made to include connate water. The case of most concern in the field is that in which filtrate from oil base mud contaminates the oil initially within the reservoir, so that filtrate and original pore fluid are miscible. Capillary pressure and residual

⁎

Mark of Schlumberger.

Fig. 8. As for Fig. 7, but with x0 = 5a, b = 100a.

saturations were therefore set to zero, and relative permeabilities were taken to be straight lines between 0 and 1, so that the ECLIPSE simulations effectively treated the fluids as miscible. In all the results presented here, the region x b 9.975 m was initially saturated with filtrate, with x b 9.975 m saturated with original pore fluid. Thus the interface was initially sharp. The effect of gravity has been considered by Sherwood (2005). If the filtrate density differs from that of the original pore fluid by an amount Δρ, buoyancy-induced Darcy velocities of order kgΔρ / μ can be created, where g is the acceleration due to gravity. In a vertical wellbore this can lead to vertical flow of filtrate until an impermable barrier is encountered (Dussan V and Auzerais, 1993). The timescale for drilling is sufficiently long that this can appreciably modify the position of the interface between filtrate and original pore fluid prior to sampling. However, pumping rates during sampling are sufficiently high that velocities due to buoyancy are small compared to those due to pumping in the region of interest around the sampling probe. Gravitational effects were negligible in the experiments of Tarvin et al. (2008): in particular, in a bounded flow cell the velocities due to fluid withdrawal do not decay far from the probe and therefore remain greater than velocities due to buoyancy everywhere within the cell. We therefore take the densities of the two fluids to be that of water, so that gravity plays no role in the simulations: results can therefore be compared directly to results of the complex variable analysis of Sections 3–5 in which gravity was similarly absent. Two production wells were situated at x = 0. One was open over the interval −a1 b y b a1 and represented the sampling port, the other was open over the two intervals a1 + a2 b |y| b a1 + a2 + a3 = a, and represented the (combined) guard ports. At x = 25 m was an injection well, open over the entire interval −25 m b y b 25 m. This allowed connate fluid to enter the simulated flow cell and replace the fluid withdrawn through the sample and guard ports.

Fig. 9. Filtrate fraction F as a function of non-dimensional pumped volume V̂ = V / Vs for the case b = 10a, x0 = 5a, corresponding to Figs. 6 and 7. (a) F in combined guard and sampling probes. Other curves show F in sampling probes of dimension (b) α = 0.2, (c) α = 0.4, (d) α = 0.6.

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

139

Fig. 10. Initial saturation for the ECLIPSE simulations. White represents original pore fluid, black represents filtrate. The probe is at the centre of the lower boundary (see Fig. 11).

Fig. 11. Computed saturations, showing displacement of the interface when the volume withdrawn is V / Vs = 70.1. Equal viscosities μp = μf.

The total probe dimension was taken to be 2a = 3 m, so that b /a =16.66 and β =x0 /a= 6.65. In all cases the gap a2 between the sampling and guard probes was set to zero, so that results are directly comparable to the analysis of Section 5. Three cases were considered, namely a1 = 0.3, 0.6 or 0.9 m, corresponding to α =a1 /a = 0.2, 0.4, 0.6. For the results presented here, the pressure in the injection well was 0.1 bar higher than that in the two producers, and this pressure difference drove the flow. The pressure in the guard ports was set to be identical to that within the sampling port. Computed results for the volume of fluid withdrawn have been scaled by the thickness dz to give the volume V per unit thickness in the z direction. V has then been scaled by Vs = a2ϕ, for direct comparison with the predictions of the complex variable analysis of Sections 3–5. Fig. 10 shows the initial position of the interface between filtrate and original pore fluid, shown in the same orientation as the experiments of Tarvin et al. (2008) i.e. with the x axis vertically upwards and with the probe situated at the midpoint of the lower boundary. Original pore fluid (white) occupies the upper portion of the cell, and filtrate (black) the lower portion. Intermediate saturations are depicted as intermediate shades of grey. The interface is assumed to be sharp initially, as in the experiments of Tarvin et al. (2008). Fig. 11 shows the saturations predicted by ECLIPSE after a volume V / Vs = 70.1 has been withdrawn from the cell, at which point the filtrate contamination F = 0.47. The initially sharp interface has been smeared by the numerical computations.

other, except at early time. The early breakthrough of original pore fluid in the ECLIPSE simulation (e) is caused by numerical diffusion, which can be seen in the saturation plot of Fig. 11. Analytic and numerical predictions for filtrate contamination in the sampling probe are shown in curves (b–d) and (f–h). The change in filtrate contamination in the sampling probe is rapid, and the effects of numerical diffusion are therefore much more apparent in results for the sampling probe than in those for the combined guard + sampling probes. Breakthrough of original pore fluid occurs too early in the ECLIPSE simulations; at later time the contamination predicted by the ECLIPSE simulation is larger than the analytic predictions of Section 5. The effect of physical dispersion is discussed by Sherwood (2005), and we expect numerical dispersion to have a qualitatively similar effect.

6.2. Comparison against complex variable analysis Results for the fraction F of filtrate in the fluid withdrawn are shown in Fig. 12. Curves (a–d) show the analytic predictions of Section 5.3, and curves (e–h) show the predictions of the ECLIPSE FrontSim simulations. Curves (a) and (e) show the fraction of filtrate F in the combined sample + guard probe: these curves agree well with each

Fig. 12. Fraction of filtrate F as a function of the fluid volume V / Vs withdrawn. (a–d) Predictions of analytic model of Section 5, (e–h) Predictions of the ECLIPSE FrontSim reservoir simulator. (a,e) Combined sampling + guard probe; (b,f) sampling probe, α = 0.2; (c,g) sampling probe, α = 0.4; (d,h) sampling probe, α = 0.6. Equal viscosities μp = μf.

140

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

Fig. 13. Fraction of filtrate F as a function of the fluid volume V / Vs withdrawn, as predicted by ECLIPSE. (a,c,e) sampling probe with α = 0.2; (b,d,f) total guard + sampling probes. (a,b) μf = μp, (c,d) μf = 10μp, (e,f) μf = 0.1 μp.

However, even with numerical dispersion, the numerical simulations indicate (in both the 2D simulations presented here and in 3D simulations) that the filtrate contamination in samples collected by a guarded probe will be substantially smaller than that in samples collected by an unguarded probe. The square root singularity in the velocity field at the edge of the probe (Eq. (5)) was poorly captured by the numerical simulations. For the case α = 0.2, the simulations predicted that the flow rate into the sampling probe was a constant fraction 0.196 of the total flow into the combined sampling + guard probes, larger than the fraction 0.13 predicted by Eq. (7). As a result, the analysis for an unbounded 2-dimensional flow predicts that the sampling probe will be filtrate free after a volume V /Vs = 76.5 has been pumped, rather than at V /Vs = 72.4 as predicted by Eq. (8c). This accounts for some of the discrepancy between curves (b) and (f) of Fig. 12 as F decreases to 0. Tarvin et al. (2008) compare the predictions of the complex variable analysis against experimental results (their Figs. 12 and 13). Physical dispersion affects the rate of change of filtrate contamination in the sampling probe in much the same way as numerical diffusion. Dispersion in the experiment was of the same order of magnitude as in the ECLIPSE simulations of Fig. 12. The filtrate contamination collected by the guard probe varied much more slowly than that collected by the sampling probe, and dispersion was therefore of little consequence, as was found in the predictions of the reservoir simulator (curve (e) of Fig. 12). 6.3. The effect of viscosity The complex variable analysis of Sections 3–5 required the original pore fluid and filtrate to have the same viscosity, as in the experiments of Tarvin et al. (2008). Within the ECLIPSE reservoir simulator, on the other hand, it is straightforward to change the viscosities of the two fluids. We present results for the case in which the ratio of sampling

probe area to total probe area was α = 0.2. When viscosities differed, the flow rates changed during the course of the simulation: all results will be presented in terms of the total volume V of fluid withdrawn. We have already seen that numerical dispersion is to be expected in the simulations, leading to a zone in which original pore fluid and filtrate are mixed. The viscosity of the fluid mixture depends upon the properties of the two fluids. Our assumption that relative permeabilities depend linearly on the saturations corresponds to a mixture viscosity that varies linearly with saturation, rather than exponentially, as is more commonly assumed (e.g. Chen and Meiburg, 1998). Viscous fingering is expected when μp b μf. It is known that the effect of sub-grid-scale fingering can be modelled by modified mixture viscosities, as proposed by e.g. Todd and Longstaff (1972) and discussed by Fayers et al. (1992), and this can improve agreement with experiment. However, in the absence of any experiments against which computations of the effect of viscosity ratio could be compared, this modification was not implemented. Fig. 13 shows predictions of filtrate contamination for three values of the viscosity ratio κ = μp / μf. If μp = μf the displacement is marginally stable. Curves (a) and (b) show the filtrate concentration in the sampling probe and in the combined guard + sampling probes, and correspond to curves (b) and (a) of Fig. 12. If μp = 10μf the displacement is stable and the interface remains sharp. The arrival of original pore fluid at the probe is delayed (curves (e) and (f)), since low viscosity filtrate is removed by the probe more easily than high viscosity original pore fluid. This can be seen in Fig. 14, which shows the position of the interface at V / Vs = 137.1, when F = 0.50. A much larger volume of filtrate has been removed than was required to reach a similar value of F when viscosities were equal (Fig. 11). The interface stays sharp, and remains almost plane until it is very close to the probe withdrawing fluid, where velocities and velocity gradients are high. It is easier for low viscosity filtrate to flow from the sides of the cell (far from the probe) than for high viscosity original pore fluid to flow rapidly in the central region close to the probe. Even when original pore fluid has broken through to the sampling probe, filtrate continues to be preferentially removed, and the time taken for filtrate contamination to drop to zero in the sampling probe is long. If μp = 0.1μf the displacement is unstable. This shows up as very early breakthrough of original pore fluid, in both the sampling probe (c) and the combined guard + sampling probes (d). However, the viscous filtrate is removed only slowly, and there is a long tail of contamination, even when V becomes large. Motion of the interface is unstable. This shows up as considerable noise in curves (c) and (d) of Fig. 13, and has led to viscous fingers in Fig. 15, which corresponds to V / Vs = 38.2 and F = 0.13. Fig. 16 shows data from a series of simulations over a wide range of viscosity ratios κ = μp / μf. Curve (a) shows the volume V / Vs at which F = 0.99 in the sampling probe, i.e. the volume at which original pore fluid breaks through to the probe. Curve (b) shows the volume at which F = 0.5, and curve (c) shows the volume at which F = 0.01.

Fig. 14. Computed saturations, showing displacement of the interface when the volume withdrawn is V / Vs = 137.1. Stable interface displacement, μp = 10μf.

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

141

Fig. 15. Computed saturations, showing displacement of the interface when the volume withdrawn is V / Vs = 38.2. Unstable interface displacement, μp = 0.1μf.

There is considerable scatter in curve (c) when κ ≪ 1: results are noisy due to the physical instability of the interface. We see that breakthrough (curve a of Fig. 16) occurs earlier when κ ≪ 1. The effect of viscosity ratio on 2-dimensional viscous fingering towards a point sink was studied in a near-axisymmetric geometry by Ceniceros et al. (1999). When the invading fluid had zero viscosity, the volume at which breakthrough occurred was smaller, by a factor ≈0.7, than when the two fluids had the same viscosity. Chen and Meiburg (1998) performed careful simulations of miscible flow in a 5-spot geometry, and found that breakthrough volumes at highly unstable viscosity ratios could be smaller than those for identical viscosities by a factor approaching 0.4. In Fig. 16 we see that when κ = 0.02 the breakthrough volume is approximately 0.48 that required when μ p = μ f. Ceniceros et al. (1999) also studied the effect of interfacial tension between the two fluids, but since their study concerned flow in a Hele–Shaw cell, the effect of interfacial tension differed from that in porous media (which has in any case been neglected entirely in the miscible studies presented here). Chen and Meiburg (1998) studied the effect of diffusion in detail, using a Peclet number to represent the ratio of convection to diffusion. Sampling pump rates used in the field are sufficiently high that molecular diffusion is likely to be small compared to velocity-dependent dispersion (Marle, 1981). A thorough study of the combined effect of Peclet number and viscosity ratio for flow towards a sampling probe would be a major undertaking. When μp / μf = κ N 1, filtrate is withdrawn in preference to original pore fluid. A crude estimate for unbounded flow might increase filtrate velocities in the y direction (see Fig. 3) by a factor κ− 1 compared to original pore fluid velocities in the x direction, suggesting that the breakthrough volume should increase as κ− 1/2. However, this

Fig. 16. (a) The volume V / Vs at which F = 0.99 in the sampling probe, i.e. the volume at which original pore fluid breaks through to the sampling probe, as a function of viscosity ratio μp / μf, as computed by ECLIPSE. (b) The volume at which F = 0.5; (c) the volume at which F = 0.01.

clearly overestimates the effect of viscosity ratios κ N 1 seen in Fig. 16, and suggests that the presence of the lateral boundaries will be felt when βa / b ≈ κ− 1. For the geometry considered here βa / b ≈ 0.4, and so the effect of the lateral boundaries of the cell is quickly felt when κ N 1, and as a result the breakthrough curve (a) does not increase beyond V / Vs = 89. However, even after breakthrough, filtrate continues to flow towards the probes in preference to original pore fluid, and curves (b) and (c) continue to increase with κ. The initial volume of filtrate in the cell is 2bβaϕ = 222Vs, and the maximum value on curve (c) is at V = 229Vs, indicating that filtrate is entering the sampling probe until almost all filtrate has been removed from the cell. 7. Conclusions The results of numerical reservoir simulations presented here agree reasonably well with the complex variable analysis of Section 5 as long as we take into account the effects of numerical dispersion. Tarvin et al. (2008) show good agreement between reservoir simulation and experiment when the interface is taken to be at the 50% saturation contour, where the saturation of original pore fluid equals that of filtrate. This agreement gives us confidence in using the results of numerical reservoir simulators to predict 3-dimensional fluid flow towards a guarded probe pushed against rock in a cylindrical borehole. Such results can be used to improve the design of probes in order to handle a wide range of depths of filtrate invasion, fluid viscosities and anisotropic rock permeabilities. References Ceniceros, H.D., Hou, T.Y., Si, H., 1999. Numerical study of Hele–Shaw flow with suction. Phys. Fluids 11, 2471–2486. Chen, C.-Y., Meiburg, E., 1998. Miscible porous media displacements in the quarter fivespot configuration. Part 1. The homogeneous case. J. Fluid Mech. 371, 233–268. Dufrêche, J., Prat, M., Schmitz, P., Sherwood, J.D., 2002. On the apparent permeability of a porous layer backed by a perforated plate. Chem. Engng. Sci. 57, 2933–2944. Dussan V, E.B., Auzerais, F.M., 1993. Buoyancy-induced flow in porous media generated near a drilled oil well. Part 1. The accumulation of filtrate at a horizontal impermeable boundary. J. Fluid Mech. 254, 283–311. Fayers, F.J., Blunt, M.J., Christie, M.A., 1992. Comparison of empirical viscous-fingering models and their calibration for heterogeneous problems. SPE Reserv. Eng. 7, 195–203. Hammond, P.S., 1991. One- and two-phase flow during fluid sampling by a wireline tool. Transport in Porous Media 6, 299–330. Hrametz, A.A., Gardner, C.C., Wais, M.C., Proett, M.A., 2001. Focused formation fluid sampling probe. US Patent 6301959. Marle, C.M., 1981. Multiphase Flow in Porous Media. Editions Technip, Paris. Morse, P.M., Feshbach, H., 1953. Methods of Theoretical Physics, vol. II. McGraw Hill, New York. O'Keefe, M., Eriksen, K.O., Williams, S., Stensland, D., Vasques, R., 2006. Focused sampling of reservoir fluids achieves undetectable levels of contamination. Proc. Soc. Petroleum Engineers Asia Pacific Oil & Gas Conference, Adelaide, Australia, 11–13 September 2006. SPE paper 101084. Schlumberger, 2006. Fundamentals of Formation Testing. Schlumberger, Sugar Land TX.

142

J.D. Sherwood / Journal of Petroleum Science and Engineering 66 (2009) 133–142

Schlumberger, 2007. FrontSim User Guide, 2007.1. Schlumberger, Abingdon U.K. Sherwood, J.D., 2005. Optimal probes for withdrawal of uncontaminated fluid samples. Phys. Fluids 17, 083102. Sherwood, J.D., Fitzgerald, J.B., Hill, B.M., 2004. Fluid sampling methods and apparatus for use in boreholes. US Patent 6719049 B2. Tarvin, J.A., Gustavson, G., Balkunas, S., Sherwood, J.D., 2008. Two-dimensional flow towards a guarded downhole sampling probe: an experimental study. J. Pet. Sci. Eng. 61, 75–87.

Todd, M.R., Longstaff, W.J., 1972. The development, testing and application of a numerical simulator for predicting miscible flood performance. J. Pet. Techn., Trans A.I.M.E. 253, 874–882. Weinheber, P., Vasques, R., 2006. New formation tester probe design for low-contamination sampling. Proc. Soc. Petrophysicists and Well Log Analysts 47th Annual Logging Symposium, Veracruz, Mexico, 4–7 June 2006. SPWLA paper 132425Q.