Evaluating the Performance of a Fracturing Treatment Design

Article pubs.acs.org/IECR

Evaluating the Performance of a Fracturing Treatment Design Qiuying Gu and Karlene A. Hoo* Chemical Engineering, Texas Tech University, Lubbock, Texas 79409-3121, United States ABSTRACT: This paper describes the combination of a Perkins−Kern−Nordgren (PKN) type model for hydraulic fracture propagation with a method for calculation of proppant transport and settling in order to simulate the dynamic growth, stress induced closure, and final geometry of vertically oriented two-dimensional fractures. A mathematical model is developed to describe the fracture growth, fluid flow, and proppant movement along with proppant settling and bank formation. A particle tracking method which uses the concept of pseudoparticles to represent the proppant phase is used for the computation of solids distribution within the fracture and the proppant bank growth. A technique for periodically combining the elements of the computational grid allows for reduced simulation time. Using the computational model, contraction of fracture dimensions after the end of pumping can be simulated in order to determine the final shape of the propped fracture. A sensitivity analysis was conducted to study the effects of pumping rate, inlet proppant concentration, and proppant particle size on the final fracture condition. To evaluate the efficacy of different treatment designs, the resulting geometry of the propped fracture dimensions and the achieved conductivities were compared. Based on the simulation results obtained, specific recommendations on how to avoid premature tip screen-out and achieve desired fracture conductivity are presented.

1. INTRODUCTION Hydraulic fracturing is a well stimulation technique for enhancing the extraction of underground resources such as oil and natural gas. This technology has made a tremendous contribution to the oil and gas industry since its initial application in the 1940s.1 Currently, hydraulic fracturing has evolved as a standard practice applied to various types of reservoir formations.1 A typical hydraulic fracturing treatment starts with a step referred to as “perforation”, in which staged explosions along the wellbore are used to create initial paths and fracture orientation. This is followed by the introduction of fluid, called “pad”, at a high pressure inside the wellbore to initiate fractures in the formation at the perforated sites. Next, a mixture consisting of fracturing fluids, additives, and proppant is pumped into the wellbore at a sufficiently high pressure and flow rate to break the formation and extend the fractures. The suspended proppant (usually sand) is transported into the fracture by the fluid as the fracture propagates. The time scale for one treatment usually is on the order of hours depending on the local formation constraints, predefined fracture size, and the amount of proppant used. After pumping is halted and the remaining fluid seeps (in a process called “leak-off”) into the void space of the reservoir, the natural formation forces close the fracture, trapping the proppant between the fracture walls. The trapped proppant provides a high conductivity channel inside the reservoir that allows for the oil and gas to flow through to the wellbore. The natural geological conditions constrain both sensor type and placement. Downhole measurements are limited to the downhole temperature and pressure sensors and microseismic and tilemeter fracturing monitors.2,3 It is difficult to estimate accurately the actual fracture geometry and proppant transport in real time only from this minimal data. Further exacerbating this estimate is that the scale and sophistication of hydraulic fracturing treatments also have increased in recent years.4 This motivates applying computational simulations to predict the fracture extension and proppant transmitted as the fracture progresses and also to update the treatment design in parallel. © 2014 American Chemical Society

The hydraulic fracturing process is a highly complex process to model mathematically as it demands the coupling of several subprocesses for an accurate representation of the system.5−7 Many models in the open literature attempt to characterize these complex processes phenomenologically under a variety of different assumptions. A reasonable model should at least include (1) the fracture propagation and surface deformation, (2) fluid movement inside the fracture, and (3) proppant transport within the fracture. 1.1. Propagation Geometry. In general, models for computing the geometry of hydraulic fractures mainly involve the coupling between the pressure distribution established using fluid mechanics and fracture shape determined through fracture mechanics.8 Those models can be categorized based on their assumptions about geometry and degree of complexity. Perkins and Kern9 published a bi-wing model based upon Sneddon’s plane strain crack solution.10 Later, Nordgren11 added the leakoff model of Carter12 to the Perkins and Kern model, forming the classic Perkins−Kern−Nordgren (PKN) model. This model assumes constant height along the fracture length and an elliptical-shaped cross section. Other assumptions are that (i) the fracture geometry exhibits a bi-wing shape which is symmetric to the wellbore, (ii) the ratio between the length and the height is large, and (iii) fluid flows only in the fracture length direction. The other classical model is based on the initial work of Khristianovich and Zheltov13 and was further developed by Geertsma and de Klerk (KGD).14 The KGD model assumptions include those of the PKN model, with the exception that the fracture height is assumed to be greater than its length. As a result of the different geometry assumptions, PKN-type models are more suitable to represent long fractures with limited height growth while the KGD model better describes shorter fractures. Received: Revised: Accepted: Published: 10491

December 6, 2013 May 15, 2014 May 18, 2014 May 19, 2014 dx.doi.org/10.1021/ie404134n | Ind. Eng. Chem. Res. 2014, 53, 10491−10503

Industrial & Engineering Chemistry Research

Article

closed when all of the fluid leaks into the surrounding formation. At the end of pumping, proppant suspended in the fracking fluid either continues to settle onto the proppant bank or else is trapped above the bank by the walls of the closing fracture. To achieve high fracture conductivity, it is crucial to distribute the proppant in a manner that maintains as much of the fracture height and length as possible. One way to achieve this state is to keep the proppant suspended inside the fracture until the fracture is fully closed. However, if the fracture closure time is too long, it is likely that the proppant will settle to the fracture bottom, thereby limiting the final fracture height. Hence, modeling fracture closure is essential, as this process will determine the final fracture conductivity and production estimation. This work includes a model of the fracture contraction and proppant settling during shut-in. The objective of this work is to develop a mathematical model of the major phenomena that occur during the fracturing process to enable accurate prediction of the final fracture conductivity and to optimize the treatment design. As indicated, the PKN model can capture the major phenomena of the hydraulic fracturing process, with low requirements with regard to computation load. On the basis of this and also the wide use of the PKN model, it is the starting point for the development of the model in this work. To this base model structure appropriate proppant transport and proppant bank formation models are coupled together with the PKN equations to simulate the fracturing process and predict the final fracture performance and conductivity. This paper is organized as follows. Section 2 describes the mathematical development of the fracture model. Section 3 introduces the numerical implementation used to solve the fracture propagation, proppant transport, bank formation, and fracture closure equations. Section 4 provides the results of several case studies and a parametric sensitivity study. Lastly, section 5 provides a summary of the work and proposes future work.

The PKN and KGD models and variations of these models have been applied widely in treatment designs.15−17 Both the PKN and KGD models assume that fracture propagation proceeds in one direction with uniform height. This assumption is appropriate for the case where fracture growth occurs in a single rock formation layer that has bounding layers above and below. In the absence of this constraint and to make the model applicable to fractures spanning multiple layers within a reservoir, a more complicated model, called the “pseudothree-dimensional (3D) model”,1 was developed. The formulation is similar to the two-dimensional (2D) models but without the assumption of constant height. In addition to lateral propagation, fractures can grow vertically through rock layers containing different stress properties. The use of the pseudo-3D model is used when well log data cannot provide accurate predictions of where or whether the fracture height will be contained. There are certain types of formation conditions in which pseudo-3D models are not suitable because the fractures created may have arbitrary shape and orientation. An alternative to the pseudo-3D model is the planar 3D model,18−20 which assumes a planar fracture shape inside a linear elastic formation. The planar 3D model has been used more widely as computational resources have increased. Recently, some full (the geometry is not preset) 3D hydraulic fracturing models21,22 have been proposed. Due to the excessive computational burden, this model type has seen limited use.6 1.2. Proppant Transport Model. Often, the main performance goal for fracturing treatment design is to optimize the final fracture conductivity.23 To keep the fracture walls open and provide a conductive channel for oil and gas to flow requires the proppant to be distributed properly inside the fracture. It follows that, to predict and optimize proppant placement, proppant transport is a necessary part of the hydraulic fracture model. The proppant undergoes both horizontal convective transport and vertical transport (due to gravity effects). A model describing proppant transport is complex. Indeed the model must describe concentration changes due to fluid leak-off, changes in proppant settling velocity due to concentration dependence, and the phenomenon of tip-bridging, which is the condition when too high of a proppant concentration at the fracture tip halts fracture growth. In this work, a two-dimensional mathematical model is developed that describes proppant transport within a fracture. 1.3. Proppant Bank Formation. As proppant settles from the suspension to the fracture bottom, the process of the proppant bed formation begins with a simultaneous increase in the bank height. Proppant settling is thought to reduce the uniformity of proppant distribution within the fracture with potentially negative effects on the final conductivity. Experimental studies on proppant settling and bank formation have been conducted, allowing for the creation of empirical correlations.24,25 The proppant bank formation at the fracture bottom may physically exhibit complicated behavior, as discussed in the experimental results of Schols and Visser’s work,25 where three consecutive phases were observed. Early work by Kern and Perkins26 showed that the settled bank also may wash out when the fluid velocity increases due to the flow restriction imposed by the bank formation. In this work, a simple model that tracks bank formation is developed; however, the wash-out effect is not included. 1.4. Fracture Closure. Once pumping of the fracturing fluid is halted, the fracture will begin to close (in a process called “shut-in”) due to the stresses within the formation. The closure time is a function of the fluid leak-off rate. The created fracture will be fully

2. MODEL FORMULATION The model assumptions are the following: (i) The vertical fracture is confined within a single horizontal rock layer propagating along a straight line (as in the PKN model). (ii) The layers above and below have sufficiently large stresses to prevent the fracture from growing external to the rock layer, hence the vertical fracture height has a constant value. (iii) The fracture is surrounded by an isotropic homogeneous elastic material. (v) The fracture length is much larger than the fracture width. (iv) The fluid pressure is constant in all vertical cross sections. The model development includes the following phenomena: 1. fracture propagation 2. fracturing fluid flow within the fracture and leak-off into the formation 3. vertical and horizontal transport of proppant within the fracture 4. proppant distribution within the fracture 5. proppant bank shape change with time and fracture geometry The fundamental physical principles applied are continuity of mass and momentum, as well as linear elastic fracture mechanics. 2.1. Governing Equations. 2.1.1. Fluid Momentum. Based on the assumptions of fracture length being much larger than width and negligible fluid inertia, the fluid flow within the fracture 10492

dx.doi.org/10.1021/ie404134n | Ind. Eng. Chem. Res. 2014, 53, 10491−10503


Article

2.1.4. Leak-Off Model. Fluid leak-off during fracture propagation can be described by12

can be approximated by the lubrication equation for a Newtonian fluid between parallel plates:27,28 12μqx dP =− dx hf W 3

U= (1)

2Phf (1 − ν 2) E

(2)

where ν is the Poisson ratio of the formation, and E is the Young’s modulus of the formation. 2.1.3. Local Fluid Mass Balance. Nordgren’s work11 contained a simple mass balance for the fluid inside the fracture, taking account of fracture volume changes and fluid leak-off into the formation, which was based on the assumption of a pure fluid without any settling or bank formation. In the current work, both proppant transport and bank formation due to proppant settling are considered. However, the simple form of Nordgren’s mass continuity equation can be retained, as follows: H

∂q ∂W + x + hf U = 0 ∂x ∂t

t − τ (x )

(4)

where Cl is the overall fluid leak-off coefficient, t is the elapsed time since fracturing was initiated, and τ represents the initial time at which the current location inside the fracture became open. 2.1.5. Proppant Transport and Bank Formation. The proppant bank is assumed to increase within the constant height vertical fracture due to proppant settling. The main assumptions are the following: 1. Both the fluid and proppant are incompressible. 2. The only interactions between the proppant and the flowing fluid is gravity-driven settling. Other processes such as selfdiffusion or particle clustering are neglected. 3. All proppant particles are assumed to be uniform in size. 4. The proppant forms a uniformly distributed suspension inside each finite volume element. 5. Both the proppant settling velocity and concentration are assumed to vary along the lateral and vertical directions, but are uniform across the fracture width. In the current model, it is assumed that a slurry consisting of the fracturing fluids and proppant is injected at a constant rate across the open fracture height at the wellbore by way of equally distributed perforation holes. Possible slurry losses (called “spurt losses”) and other minor effects that accompany injection of the slurry through the perforation holes are not included in the current model development. The suspended proppant is assumed to have a settling velocity relative to the carrier fluid. The proppant bank begins to form and grow once the proppant reaches the fracture bottom. Following initial bank formation, any other proppant that settles as a part of the formed bank contributes to the bank growth. In the vicinity of the fracture tip there is inflow of the proppant but no outflow. As a consequence, the concentration at the fracture tip can be greater than anywhere else. This phenomenon is referred to as “tip screen-out” wherein the proppant concentration becomes so large that the proppant forms “bridges” across the fracture tip preventing further propagation. To model this phenomenon, a maximum proppant concentration value, Cmax, is defined according to the packing density of the proppant particles. When the proppant concentration reaches Cmax at the tip, screen-out is determined to have occurred and propagation of the fracture length ceases. The settling of suspended proppant can be described by

where P is the difference between the fluid pressure inside the fracture and the minimum in situ formation stress, μ is the fluid viscosity, qx is the local flow rate in the x-direction, hf is the local height of the free area fluid flow (difference between fracture height and bank height), and W is the fracture width. Eq 1 implies that a gradient in the fluid pressure develops along the propagation direction as a function of the flow rate and the resistance due to the narrowness of the fracture width. 2.1.2. Pressure−Width Relationship. It has been shown that some types of rock exhibit elastic behavior when subjected to certain ranges of stresses such as those encountered during hydraulic fracturing.29 In these instances, the deformation of a given formation can be described by elasticity theory. Using the assumption of the PKN model development that one dimension is much larger than the other, the pressure drop inside the fracture varies only along the lateral direction; thus plane-strain deformations can be assumed. Large pressure drops that accompany large fluid injection rates will result in wider fractures. Based on the plane-strain solution with the assumption of a constant vertical fracture, Sneddon and Elloit have provided a useful relationship between the maximum fracture width and the net pressure:10 W=

2C l

(3)

where H is the predefined fracture height, and U is the fluid leakoff per unit height accounting for both walls of the fracture. Since both proppant and fluid are assumed to be incompressible, eq 3 is in the form of a volume balance. This formulation of the continuity equation was created by considering a control volume that has infinitesimal extent in the x-direction while covering the entire fracture height, and does not differentiate between fluid and proppant whether suspended or settled. Therefore, the first term describing the control volume change due to fracture expansion or contraction is independent of the bank height. The bank height affects the area open to flow or leak-off, so the terms including flow qx and height dependent leak-off, hfU, will include the effect of restricted flow area due to bank formation. Equation 3 requires two boundary conditions and an initial condition. At the wellbore, the flow rate is specified by qx(x,t) = Q0(t), where Q0 is the rate of slurry injection. At the fracture tip, x = L(t), the commonly chosen boundary condition11,30 is that the net pressure is zero, which leads to the condition W(L(t),t) = 0, via eq 2. Initially, the fracture is closed, providing W(x, 0) = 0.

CVsW d(H − hf )W = dt 1−ϕ

(5)

where Vs is the proppant settling velocity, hf is the height of the slurry, the difference H − hf is the bank height, and ϕ is the proppant bank porosity. 2.1.6. Proppant Settling Velocity. The proppant entering at the wellbore is either convectively moving at the carrier fluid velocity or settling due to the effects of gravity. Assuming the proppant particles are spherical in shape, the drag force that results from a viscous fluid onto a single particle falling in an infinite laminar flow can be estimated by Stokes’ law.31 Then, the settling velocity is given by

Vs = 10493

(ρsd − ρf )gd 2 18μ

(6)



Article

where ρsd is the proppant particle density, ρf is the pure fluid density, g is the acceleration due to gravity, d is the proppant particle diameter, and μ is the fracture fluid viscosity. However, eq 6 is only suitable for the free motion of a single particle. In the hydraulic fracturing process, proppant particles are suspended in a slurry. Their settling velocity has to be modified to account for the slurry rheological properties and their changes, such as that due to the effect of slurry viscosity on the proppant settling velocity. Conversely, variations in proppant concentration that result from proppant settling also will affect the fluid viscosity. Previous studies have been concerned with modeling the particle and particle fluid interactions to adjust the slurry viscosity.32−34 Alternatively, the observations can be lumped into empirical formulas, which is a common modeling practice:32 −α ⎛ C ⎞ μa = μ0 ⎜1 − ⎟ Cmax ⎠ ⎝

In addition to the modified mass balance eq 13, other equations to be solved include eqs 5 and 8. The latter requires the concentration C(x,y,t). Local concentrations can be found by considering the proppant as a finite number of representative particles and determining where these particles travel with time. Each particle is assumed to be moving with the local fluid velocity in the direction of fracture propagation and takes on the local settling settling velocity Vs (see eq 8). 2.2. Fracture Closure. In this work, the closure (called “shutin”) of the vertical fracture is described by the same governing equations but the slurry pumping rate is set to zero. The following assumptions are made: 1. The fracture stops propagating after pumping stops. 2. The formed proppant bank is assumed to be tightly packed and will not deform during the fracture closure process. 3. If the local concentration reaches the packed proppant concentration because of the contraction of the local fracture dimensions due to leak-off, the fracture is considered to be closed (and propped) in that area. 4. If the local fracture width becomes smaller than the proppant particle diameter before the propping condition is satisfied, then the fracture is considered to be fully closed in this area with no propping achieved. 5. The fracturing fluid does not experience reverse flow toward the wellbore during shut-in. With these assumptions, the final bank geometry after the fracture is fully closed can be determined.

(7)

where μ0 is the pure fluid viscosity, μa is the modified fluid viscosity corrected for proppant concentration effect, C is the proppant concentration, Cmax is the maximum packed proppant concentration, α is the viscosity-slurry exponent in the range 1.2−1.8, depending on proppant and fluid type. The corrected settling velocity8 can be described by Vs =

2 (1 − C)2 3gd ρsd − ρf 101.82C 108 72μa

(8)

2.1.7. Coupling the Governing Equations. From eq 2, the gradient of fluid pressure along the propagation direction can be expressed as ⎛ 1 dW dP E W dh f ⎞ ⎜ ⎟ = − 2 2 dx 2(1 − ν ) ⎝ hf dx h f dx ⎠

3. NUMERICAL IMPLEMENTATION The mathematical model introduced above is simulated to investigate the behavior of fracture propagation, proppant transport, and bank formation when subjected to realistic values of volumetric pumping rates and fracturing fluid and proppant properties. Fluid density, proppant particle strength, and proppant density are held constant. The governing partial differential equations are complex and highly nonlinear. To solve them accurately requires a suitable numerical algorithm to update the following seven variables at each time step: qx(x,t), hf(x,t), W(x,t), U(x,t), C(x,y,t), Vs(x,y,t), and L(t). Cartesian coordinates are used in all simulations. 3.1. Meshing Strategy. A two-dimensional (2D) grid system is generated to represent the created fracture geometry. Under the assumption of constant vertical height, the y-direction mesh is evenly distributed and fixed spatially and temporally. Two special features of the fracturing process need to be noted when designing a numerical solution procedure. First, in the xdirection (direction of fracture propagation), as the fracture tip advances, the number of equations required to solve for all of the variables at the growing number of cells will increase, making the whole simulation excessively burdensome to compute. Second, several important quantities including pressure, width, and flow rate vary slowly along the length of the fracture, except near the tip area, where relatively large changes occur. In previous hydraulic fracturing simulation models, two commonly implemented numerical methods for dealing with the growing geometry are (i) having a moving coordinate frame that changes position with time and (ii) periodic remeshing of the fixed coordinate system.6 Detournay30 introduced a method of the first type, in which a fixed number of elements are scaled by the fracture length, thus limiting the number of elements while still retaining a high element density near the fracture tip. On the other hand, obtaining a high number of elements near the tip within a fixed coordinate system causes the number of

(9)

2

The term W/hf dhf/dx can be ignored under the assumption that fracture width W is much smaller than the slurry height hf. Substitution of eq 9 into eq 1 provides an expression that relates the flow rate, qx, and fracture width, W: qx = −

hf W 3 ∂W 12μMc ∂x

(10)

where Mc the fracture compliance which is given by Mc =

2hf (1 − ν 2) E

(11)

Hence, the gradient of the flow rate, ∂qx/∂x, can be expressed as ∂qx ∂x

=−

hf ∂ 2W 4 48μMc ∂x 2

(12)

and eq 3 as H

hf ∂W ∂ 2W 4 − + hf U = 0 48μMc ∂x 2 ∂t

(13)

Under the assumptions of constant viscosity and total fracture height, eq 13 now contains only two unknowns, the slurry flow height (hf) and fracture width (W). The boundary condition at the tip remains W(L(t),t) = 0, while the condition qx(x,t) = Q0(t) becomes ⎡ ∂W 4 ⎤ 48μMc Q0 =− ⎢ ⎥ hf ⎣ ∂x ⎦x = 0 10494



Article

indices are renumbered with the combined cell marked as cell 1. Also illustrated is combining larger size cells at propagation step 145 since the number of cells at this larger size has reached Ne. As the treatment progresses, a reasonable number of computational cells are maintained even as the number of iterations increases. Without cell combination, the number of equations to solve would increase by one each iteration, leading to the total number of equations solved to be approximately N2/2, where N is the total number of time steps. With cell combination, a maximum of M, M ≪ N, computational cells can be maintained for the length of the simulation, leading to a maximum of M × N times the equations were solved, so that a linear growth rate exists between the number of solutions achieved and the total number of iterations, instead of a quadratic growth rate. For example, in this work, less than 100 cells are required through each of the 10 000 iterations. During remeshing, the fracture and fluid properties are adapted to the new mesh. The updated width is taken to be the average of the cells before they were combined, in order to maintain mass conservation. Flow rates into and out of the new larger cell are taken from the entrance flow rate of the first cell and exit flow rate of the last cell that were combined. Because particle tracking is used to determine the proppant location, recalculation of the concentration within the new mesh is straightforward. 3.2. Numerical Solution Procedure. The steps of the numerical algorithm are shown in Figure 3. k−1 1. At time step tk, the values of the variables {qk−1 xi ,hfi ,k−1 k−1 k−1 k−1 } are known, where subscript i = 1, ..., Wk−1 i ,Ui ,Ci,j ,Vsi,j ,L

computation elements to become too large. Therefore, periodic remeshing is useful,6 to reduce the total number of elements while concentrating elements near the fracture tip, which is the method applied in this work. This is illustrated in Figure 1, where

Figure 1. Schematic of the computational mesh.

a finer mesh is applied near the tip in the x-direction. Nx is the number of cells in the x-direction, which varies as the fracture propagates, and Ny is the number of cells in the y-direction, which remains constant throughout the simulation. At each time step, the fracture is assumed to extend a predefined length, ΔL, with the time step size Δt iteratively solved for. As the treatment time increases and the number of cells increases, in order to mitigate the computational burden, a simple remeshing strategy is proposed, in which only a fixed number of cells, Ne, of size ΔL are retained. When the number of cells exceeds the threshold Ne, a fixed number of the cells, Nc, of size ΔL are combined to form a larger cell of size NcΔL. As the treatment progresses, if the number of cells at this larger size NcΔL exceeds Ne, they also are combined into even larger size cells; this process repeats for the length of the fracture propagation. The process of element combination begins with cells furthest from the wellbore in order to preserve smaller cells close to the fracture tip where larger variations in the solution exist. By performing cell remeshing, the number of cells required to solve all of the equations can be managed within a relatively small range, which helps to avoid numerical problems. The process of cell remeshing is illustrated in Figure 2, where Ne and Nc are 25 and 5, respectively. On the 25th iteration, the number of cells of the same length reaches the threshold value Ne; thus Nc of these cells, 1−5, located on the wellbore (left-hand) side are combined into a larger cell and the cell

Figure 3. Numerical solution procedure for fracture propagation.

Figure 2. Illustration showing combination of computational elements. 10495



Article

Nx is the index of the cell in the x-direction and subscript j = 1, ..., Ny is the index of the cell in the y-direction. Indices such as i + 1/2 are introduced in the x-direction to indicate the flow rate qx and velocity Vx at the exit of cell i. 2. The coupled mass balance eq 13 is discretized with an implicit Euler difference in time and upwind differentiation scheme in space: h fki − 1 ⎛ ∂ 2W 4 ⎞k ΔW ik k−1 k − H ⎟ + h fi Ui = 0 ⎜ Δt 48μMc ⎝ ∂x 2 ⎠i

fracture is updated at each time step. A rectangular cross section is assumed during particle tracking. With the available updated fracture mesh data and the updated lateral and vertical positions of the pseudoparticles, the index of the containing cell can be located for each pseudoparticle. By dividing the sum of the mass of all the pseudoparticles within a cell by the cell volume, the proppant concentration in each cell is updated. 6. Using the bank height at time step k − 1 and the particle settling velocities, the number and locations of particles that settle at time step k can be determined. Thus, the proppant bank size at each node i in the x-direction is updated. 7. Using the updated proppant concentration data, the proppant particle settling velocity is updated based on eq 8. 8. Lastly, the combining of cells is done to limit the number of x-direction elements under a certain limit to maintain numerical stability. With the exception of the tip area where a finer mesh is necessary for the duration of the treatment process, whenever the number of cells of one size reaches a predefined value, these cells are combined into a cell of larger size and the associated variables are updated for this new cell. The numerical method presented is capable of solving the model that describes the contraction of the fracture geometry and the settling of the suspended proppant as the fracture width and length recede simultaneously. This results in determination of the final proppant bank shape and an estimate of the fracture conductivity. 3.2.1. Fracture Closure. The procedure to solve the model that describes fracture closure (shut-in) is shown in Figure 4.

(14)

where ΔWki ≡ Wki − Wk−1 is the increase in the width at time step i k at cell i, and the derivative (∂2W4/∂x2)ki is computed using a spatial central difference in space: ⎡ 4 k 4 k ⎛ ∂ 2W 4 ⎞k 1 ⎢ (W )i + 1 − (W )i = ⎟ ⎜ Δxi ⎢⎣ 1 (Δxi + 1 + Δxi) ⎝ ∂x 2 ⎠i 2 −

⎤ (W 4)ik − (W 4)ik− 1 ⎥ 1 (Δxi + Δxi − 1) ⎥⎦ 2

(15)

Equation 15 is nonlinear with respect to Wki . As suggested k4 previously,11 since Wk−1 is known, W can be approximated from i i a Taylor series expansion: (W ik)4 ≈ (W ik − 1)4 + 4(W ik − 1)3 ΔW ik + higher order terms (16)

Substitution of this equation into eq 14 yields a linear tridiagonal system in ΔWki ⎛ k−1 ⎧ k−1 3 1 ⎪ 4(W i ) ⎜ H + h fi ⎨1 ⎪ ⎜ Δt 48μMc Δxi ⎩ (Δxi + Δxi + 1) ⎝ 2 ⎞ ⎫ ⎪ ⎬⎟ΔW ik + 1 ⎪⎟ ( ) Δ + Δ x x 1 i − i ⎭⎠ 2 4(W ik − 1)3

⎛ h k−1 ⎞ 4(W ik−−11)3 ⎟ 1 f k − ⎜⎜ i 1 ⎟ΔW i − 1 48 M x μ Δ ( ) x x Δ + Δ c i 1 i i − ⎝ ⎠ 2 ⎛ h k−1 ⎞ 4(W ik+−11)3 ⎟ 1 f − ⎜⎜ i ΔW ik+ 1 1 ⎟ 48 M x μ Δ c i 2 (Δxi + Δxi + 1) ⎠ ⎝ =−

⎛ k−1 4 k−1 4 1 ⎜ (W i ) − (W i + 1 ) 48μMc Δxi ⎜⎝ 1 (Δxi + Δxi + 1) 2

+

⎞ (W ik−−11)4 − (W ik − 1)4 ⎟ k−1 k 1 ⎟ + h fi Ui ( x x ) Δ + Δ i i−1 ⎠ 2

h fki − 1

Figure 4. Numerical solution procedure for fracture closure.

(17)

that can be solved at each cell i (i = 1, ..., Nkx). This updates the width at every cell at the current time step. 3. The time interval Δtk is determined iteratively by repeating the procedure until the mass balance error approaches zero. 4. The flow rate q and fluid velocity Vx are updated using the value of Wki . 5. Calculation of the local proppant concentration at cells in both the x- and y-directions proceeds. For computational efficiency, at each time step, it is assumed that the proppant particles entering the well can be grouped and thus represented by pseudoparticles. The pseudoparticle position inside the

During fracture closure, no mesh reconfiguration occurs. Once fracture closure is initiated, at each predefined time step, the solution procedure is as follows: 1. Update the average fracture width and flow rate at each mesh position along the fracture length using eq 13. 2. Based on the local fracture width and proppant concentration, it is determined if the fracture walls are closed in each cell at the current time step. 3. The position of the suspended proppant particles are tracked and the local proppant concentrations are updated along with the height of the settled proppant bank. 10496



Article

4. The fracture length is reduced if all of the cells at any x-position become closed. 3.2.2. Fracture Conductivity Calculation. During the hydraulic fracturing treatment, by placing layers of proppant inside the created fracture, the conductive area for oil and gas to flow from the reservoir to the wellbore is greatly enhanced. One of the most important criteria in evaluating productivity improvement is the measure of fracture conductivity, Cf = KfW̅ , where Kf is the proppant bank permeability, and W̅ is the average fracture width. This conductivity, along with the rock permeability (Kr) and bank length (Lf), is used to calculate the dimensionless fracture conductivity:

CfD =

WK ̅ f Lf K r

Table 1. Model Parameters and Treatment Schedule Used in the Case Study30,35 definition

symbol

proppant diam acceleration due to gravity fluid-loss coeff max proppant concn Young’s modulus of formation proppant permeability rock permeability vertical fracture height proppant particle density pure fluid density pure fracturing fluid viscosity Poisson ratio of formation proppant bank porosity proppant concn vol (kg/m3) (m3)

(18)

In this work, fracture conductivity is calculated as follows for the proposed fracturing model. First, as illustrated in Figure 5, at

0 (pad) 240 480 720

50 20 25 10

d g Cl Cmax E Kf Kr H ρsd ρf μ ν ϕ proppant diam (μm) 0 450 450 800

value 1000 μm 9.8 m/s2 6.3 × 10−5 m/s1/2 0.64 104 MPa 500 darcy 0.01 mD 10 m 2648 kg/m3 1000 kg/m3 0.56 Pa·s 0.2 0.12 pumping rate (m3/s) 0.04 0.04 0.04 0.06

Figure 5. Fracture cross section at one position along the fracture length.

each grid point i along the fracture length, the fracture width is defined as an averaged value along the fracture height N

Wi̅ =

∑ j =y 1 ΔhjWj N

∑ j =y 1 Δhj

(19)

or mean this quantity is averaged along the bank length to give the final fracture bank width N

W̅ =

∑i = 1 Wi̅ Δxi N

∑i = 1 Δxi

Figure 6. Fracture width at the wellbore. (20)

In section 4, as in industry practice, fracture conductivity serves as an important criterion to evaluate the fracture treatment design.

propagates and then subsequently decreases as the fracture closes. The progression in time of the fracture width along the length of the fracture is shown in Figure 7. The changes in the fracture width in the tip region support the use of a finer mesh in this region. Figure 8 shows the fracture propagation and proppant distribution at the end of the pumping schedule (just prior to shut-in). It is seen that the proppant bank height decreases from the wellbore to the fracture tip, and the amount of suspended proppant near the wellbore is also much greater than the amount suspended near the fracture tip. Although not shown, proppants with a larger diameter are suspended at a higher concentration near the wellbore at the end of pumping. In contrast, proppants with a smaller diameter have traveled further from the wellbore before shut-in and have settled in the region around the tip. This type of pumping schedule can prevent early tip screen-out while maximizing conductivity in the vicinity of the wellbore. For the treatment listed in Table 1 and the assumptions made, approximately 40% of the total injected proppant has settled by

4. RESULTS AND DISCUSSION OF RESULTS 4.1. Case Study. The model’s parameters and the pumping schedule are listed in Table 1. The fracturing fluid is scheduled to be pumped at different rates along with varying proppant types and concentrations. Following the pumping of pure fluid, proppants of two different diameters are pumped sequentially. First, the smaller diameter proppant is pumped at relatively low concentration and pumping rate, and next the larger diameter proppant is pumped at a higher flow rate and concentration. The pumping duration is fixed at 34 min followed by a shut-in duration of 23 min. The fracture height is assumed to be constant and is set at a value of 10 m. The fracture propagation is assumed to take place within a homogeneous reservoir environment. Figure 6 shows how the fracture width changes at the wellbore as a function of time. It is observed that the fracture width near the wellbore at first increases monotonically as the fracture 10497



Article

Figure 7. Evolution of the fracture width along the fracture length at several snapshots in time.

Figure 9. Suspended proppant concentration at the end of the pumping schedule.

Figure 10. Final proppant bank surface geometry after the fracture is closed. Figure 8. Suspended proppant distribution and proppant bank height (black line) at the end of the pumping schedule.

the time pumping is stopped (34 min). The furthest distance traveled within the fracture by the proppant is about 240 m, 86% of the fracture length of 280 m; the settled proppant covers 5% of the created fracture area. Figure 9 shows the distribution of the suspended proppant, which exhibits higher concentrations in the vicinity near the wellbore and above the proppant bank versus near the tip and away from the wellbore. After the designed amount of proppant and fluid are pumped according to the treatment schedule in Table 1, fluid injection is stopped. Leak-off and proppant bank formation continue as close-in begins. When the fracture is fully closed, the suspended proppant either becomes a part of the proppant bank or is trapped by the fracture walls. The proppant bank geometry after closure is shown in Figure 10. The former is a surface plot of the proppant bank after completion of fluid leak-off, at which time all of the suspended proppant has either settled or become trapped between the fracture walls. Figure 11 compares the deposited bank height along the fracture length before and after shut-in. It is observed that the proppant bank appears to have a triangular-like cross section. Analysis reveals that

Figure 11. Comparison between the bank height along the fracture length at the end of pumping (dash line) and after fracture closure (solid line).

the bank height increased ∼2 m near the wellbore and ∼1 m elsewhere, with the additional bank volume contributed by the 10498



Article

length) is detected and only ∼30% of the designed proppant treatment is actually injected. When the proppant concentration increases to C2, similar results are found. Screen-out also occurs at pumping rates Q3 and Q4 with only ∼70% of the designed proppant treatment injected. The suspected reason is that treatments with higher pumping rates allow the proppant to travel at a higher velocity leading to a high proppant concentration at the tip. In addition, for a fixed amount of treatment fluid and proppant, pumping proppant at a lower concentration requires a longer pumping time and larger slurry volume, reducing the available pad volume and increasing the likelihood of tip screen-out. It is observed that, at higher concentrations, C3 and C4, tip screen-out is absent. Therefore, it appears that, to avoid tip screen-out at a given pumping rate, increasing the pad volume is more effective than decreasing the slurry proppant concentration. At higher concentrations, C3 and C4, as pumping rates increase the average fracture width during fracture propagation increases while the final proppant bank’s average width and the amount of settled proppant at the end of pumping decrease. The final fracture conductivity shows a decreasing trend with increasing pumping rate. From this result, it appears that proppant that settled prior to fracture closure contributes more to the final fracture conductivity when compared to suspended proppant that settled later or is trapped during shut-in. Ultimately, the final proppant bank appears to be more uniformly distributed at higher pumping rates. 4.2.2. Proppant Concentration. For flow rate Q1, the results show that, as the proppant concentration increases, the fraction of proppant settled before closure decreases; thus the final fracture conductivity increases. The proppant bank formed at higher concentrations exhibits a more triangular-like bank shape. This shape is more obvious at higher pumping rates (C3Q3, C4Q3 and C3Q4, C4Q4). Also at these higher pumping rates, there exists the presence of a tiny bank tail near the tip that is disconnected from the main bank, which does not contribute much to the final conductivity. Comparing all cases, the most uniform proppant placement is obtained at C3Q3. On the basis of the results, it appears that higher inlet proppant concentrations give higher fracture conductivities and more uniform-like proppant bank distributions. In contrast, low proppant concentrations result in low conductivities and nonuniform fractures after closure. 4.2.3. Proppant Size. The proppant plays a critical role in preventing complete fracture closure, and its final distribution directly affects the flow of oil and gas to the wellbore. Several physical properties, such as proppant strength, grain size, sphericity, quality, and proppant density affect the proppant conductivity. To illustrate the effect of proppant size on its transport, bank formation, and fracture conductivity, simulations at three different proppant diameters with the other properties held constant are performed. The proppant particle diameters chosen are 1700, 800, and 450 μm; these are typical sizes used in practice.35 The pumping rate is set at Q2 = 0.04 m3/s and proppant concentration is set at C3 = 720 kg/m3 for each experiment. Figure 14 shows the proppant distributions at the end of pumping and closure. It is observed that larger diameter proppant tends to settle near the wellbore. The bank height for 1700 μm diameter proppant reaches 6 m compared to the bank heights of approximately 1 m or less achieved with the other two proppant diameters. The final bank formed by the smaller diameter proppant is located closer to the fracture tip, while the

suspended proppant after pumping is stopped. The fracture width shrinks while the local bank height increases as the fracture closes. At the tip cells where the smallest fracture width maintains, during fracture closure, a height of ∼2 m was formed by the suspended proppant. Any proppant suspended at the fracture closure starting time also contributes to the fracture treatment success, and thus closure is an important phenomenon to model. For the selected treatment schedule, it has been shown that proppant bank formation on the bottom of a vertical fracture is affected by the pumping schedule and the injected proppant properties. The computational model predicts that using a smaller diameter particle at the start of proppant pumping is effective in preventing screen-out (as is industry practice). Larger proppant can then be added toward the end of treatment. The analysis also reveals that to achieve a desired proppant placement these factors should be considered when designing the treatment schedule. 4.2. Sensitivity Analysis. As demonstrated above, some of the design operating conditions in the fracture pumping schedule have a strong impact on the success of the fracture treatment. In this section, an analysis of the sensitivity of the fracture model to several input parameters is presented with the objective of evaluating their effects on the simulation model’s predictions of fracture propagation and proppant bank formation. Three main treatment parameters are studied: proppant injection concentration, pumping rate, and proppant particle size (diameter). Values of the proppant concentration and fluid pumping rate are listed in Table 2. The other input parameters remain unchanged from the values in Table 1. Table 2. Proppant Concentration and Pumping Rate Values proppant concn (kg/m3) pumping rate (m3/s)

C1 = 240

C2 = 480

C3 = 720

C4 = 960

Q1 = 0.02

Q2 = 0.04

Q3 = 0.08

Q4 = 0.10

Sixteen simulation trials are conducted with each possible combination of proppant concentration and pumping rate. As for the notation used: “C1Q1” represents that proppant was pumped at concentration C1 and rate Q1. For each experiment, the pumped pad volume, slurry volume, and pumping time varied with changes in the pumping rate. Higher pumping rates require less pumping time, while greater proppant concentrations correspond to less slurry pumping time and larger pad volumes. The fracture experiments consisted of pumping 140 m3 of fluid when slurry concentration C1 = 240 kg/m3 is used and 100 m3 of fluid at other concentrations. The total amount of proppant is fixed at 10 m3 for each experiment. The predicted fracture shape and bank profile are shown in Figures 12 and 13, respectively. Table 3 compares the results of the different experiments. 4.2.1. Pumping Rate. To evaluate the effect of the pumping rate on treatment performance, consider the lowest pumping rate Q1 as the base case. Starting with C1Q1, around 90% of the injected proppant settles to form the proppant bank before fracture shut-in. This occurs since, at a low pumping rate, the proppant travels slowly into the fracture and has more time to settle. At an increased pumping rate, Q2, the amount of proppant settled at the end of pumping decreases to 87% and the formed proppant bank has a more uniform distribution and increased length. At even higher pumping rates, the proppant travels further along the fracture length. In both the cases of C1Q3 and C1Q4, the condition of screen-out (proppant bank completely blocking the fluid flow at some locations along the fracture 10499



Article

Figure 12. Proppant bank formation as a function of proppant concentration and pumping rate.

larger diameter proppant formed a bank closer to the wellbore. Note that the proppant bank porosity (and therefore fracture conductivity) decreases with proppant diameter. This is because the smaller diameter proppant packing is more dense, leaving less void space to conduct fluid flow. However, the larger diameter proppant may be at risk of being crushed during fracture closure.35 There also exists a trade-off between the larger porosity (and area available to flow) achieved with proppants of larger diameters and a proppants of smaller diameters that can be transported further into the fracture. The sequential pumping of different proppant diameters is the usual means to try to take advantage of these different characteristics. Altogether, the results presented in this section show the importance of simulating the treatment design in order to reduce the probability of screen-out and to obtain a reasonable prediction of fracture performance. 4.3. Potential Applications in Model-Based Control. The foregoing results demonstrate that significant variation in fracture geometry and final conductivity can be obtained by manipulating the pumping rate and fluid makeup (including proppant size and concentration). The simulation model developed in this work could potentially be used in design of a fracture treatment schedule, provided that the physical

Table 3. Calculated Fracture Conductivity, Average Fracture Width, and Settled Proppant When Pumping Stops for the 16 Experiments conductivity

C1

Q1 16.44 Q2 15.38 Q3 8.00 Q4 7.65 av fracture width (mm)

C2

C3

21.62 15.67 7.89 7.76 C1

C2

C4

22.74 15.70 11.96 8.30 C3

12.21 12.15 5.12 4.53

28.8 18.98 14.72 12.73 C4

Q1 13.72 Q2 14.41 Q3 5.42 Q4 4.76 settled proppant at end of pumping (%)

12.10 11.40 9.52 6.30

12.00 11.15 8.35 8.64

C1

C2

C3

C4

Q1 Q2 Q3 Q4

92.59 87.14 18.10 13.00

87.48 74.03 27.76 19.02

83.11 56.19 25.22 18.97

78.24 44.81 20.24 14.68

parameters of the formation are obtained such as through a mini-frac test. This is similar to how existing fracture simulators 10500



Article

Figure 13. Proppant bank height as a function of proppant concentration and pumping rate.

are used in practice.1 However, it is desired to develop techniques for the online adjustment of the fracture treatment schedule in real time based on the available sensing information, in order to correct for uncertainties and deficiencies in the simulator model or formation parameters that cause the fracture propagation to deviate from the predicted trajectory. Therefore, current work is focused on the application of existing process control techniques to the hydraulic fracturing system. This goal is complicated by the presence of distributed and nonlinear dynamics as well as a lack of available sensing information. In fact, downhole data of the hydraulic fracturing process are limited to pressure, temperature, and microseismic measurements. Downhole pressure can be used to estimate the fracture width near the wellbore, such as via eq 2, while the microseismic readings can be used to estimate changes in the fracture length and height over time. This leaves many important system states, such as local fracture width, proppant concentration, proppant bank shape, and so on, that need to be estimated for control purposes. The simulation model developed

in this work can be applied to the problems of state estimation and forecasting that is required within a model-based control framework. Solutions for the feedback control problem of many types of partial differential equation systems have been formulated, including for systems with nonlinearities and time-dependent spatial domains.36 However, the current simulation model falls outside of the class of problems having an exact solution for feedback control. A more appropriate approach may be the use of empirical eigenfunctions to decompose the model into simpler dynamics which capture the most important system behavior such as work by Christofides and co-workers36,37 and by Hoo and co-workers.38−40 Using reduced order models, the implementation of control techniques originating in the process industries can be explored, such as DMC (dynamic matrix control) and other model-based control techniques. The fracture simulator and associated reduced order models can also provide the state estimates of the variables to be controlled. With particle filtering or ensemble Kalman filtering, estimation of model parameters, 10501



Article

Figure 14. Effect on proppant bank geometry of proppant particle size.

proppant concentration results in higher conductivities and for a given inlet proppant concentration the fracture conductivity decreases as pumping rate increases. The final proppant bank shows more uniform distribution at higher pumping rates and higher proppant inlet concentrations. Sequential pumping of different sizes of proppant will most likely result in a more uniform bank distribution. These findings about how the parameters affect the placement of proppant are valuable in the design of a successful fracture treatment. A discussion of future research directions in real-time model-based control of fracture propagation and proppant placement has also been included.

such as the leak-off coefficient, can be achieved. In this framework, the possibility of adaptive model-based control of fracture propagation and proppant placement could be realized, in order to incorporate field measurements in real time for the optimization of predicted final fracture performance.

5. SUMMARY This work described the development of a hydraulic fracturing computational model that allows for evaluating the performance of a designed fracture treatment. The model includes the following features: 1. Determination of the fracture geometry shape change with timethe fracture propagation component. 2. Dynamic and spatial tracking of the proppant location inside the fracture and computation of the distribution of the injected proppant during both fracture propagation and closurethe proppant transport component. 3. Determination of the proppant bank formation and shape during shut-inthe proppant bank formation component. 4. Characterization of fracture shut-in phenomena including fluid leak-off and the final propped fracture shapethe fracture closure component. 5. Analysis of fracture conductivity as a function of the designed fracture treatment. A sensitivity analysis was carried out with respect to three main design parameters: proppant concentration, pumping rate, and proppant size (diameter). The results of the fracturing model revealed that the proppant bank formation was sensitive to these parameters, displaying different size, shape, and final conductivity based on the combination of proppant concentration, pumping rate, and proppant size chosen. Guidelines for successful fracture treatment designs to achieve the desired performance objective based on this parameter sensitivity analysis were provided. For example, it appears that, for a given pumping rate, a higher inlet

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS Q.G. acknowledges the financial support of the TTU Process Control & Optimization Consortium.

10502

NOMENCLATURE d = proppant diameter, μm g = acceleration due to gravity, m/s2 hf(x,t) = local slurry flow height, m qx(x,t) = fracturing fluid flow rate, m3/s x = lateral fracture propagating direction, m t = lapse time since start of pumping, s C(x,y,t) = proppant concentration (v/v) C1 = fluid-loss coefficient, m/s1/2 Cmax = maximum proppant concentration CfD = fracture conductivity E = Young’s modulus of formation dx.doi.org/10.1021/ie404134n | Ind. Eng. Chem. Res. 2014, 53, 10491−10503


Article

(19) Naceur, K.; Thiercelin, M.; Touboul, E. Simulation of fluid flow in hydraulic fracturing: implications for 3D propagation. SPE Prod. Eng. 1990, 5, 133−141. (20) Vandamme, L.; Curran, J. A three-dimensional hydraulic fracturing simulator. Int. J. Numer. Methods Eng. 1989, 28, 909−927. (21) Carter, B.; Desroches, J.; Ingraffea, A.; Wawrzynek, P. Simulating fully 3D hydraulic fracturing. In Modeling in Geomechanics; Zaman, M., Gioda, G., Booker, J., Eds.; Wiley: Hoboken, NJ, 2000; pp 525−557. (22) Hsu, Y.; Dang, X.; Chilton, W.; Chang, P.; Stelin, I.; Gusain, D.; Northington, N.; De Pater, H. New physics-based 3D hydraulic fracture model. SPE Hydraulic Fracturing Technology Conference; Society of Petroleum Engineers: Richardson, TX, 2012. (23) Davies, D.; Kuiper, T. Fracture conductivity in hydraulic fracture stimulation. J. Pet. Technol. 1988, 40, 550−552. (24) Babcock, R. E.; Prokop, C.; Kehle, R. Distribution of propping agent in vertical fractures. Drilling and Production Practice; American Petroleum Institute: Washington, DC, 1967. (25) Schols, R.; Visser, W. Proppant bank buildup in a vertical fracture without fluid loss. SPE European Spring Meeting; Society of Petroleum Engineers: Richardson, TX, 1974. (26) Kern, L.; Perkins, T.; Wyant, R. The mechanics of sand movement in fracturing. J. Pet. Technol. 1959, 11, 55−57. (27) Pinkus, O.; Sternlicht, B. Theory of Hydrodynamic Lubrication; McGraw-Hill: New York, 1961. (28) Economides, M. J.; Nolte, K. G.; Ahmed, U. Reservoir Stimulation; Wiley: Chichester, U.K., 2000; Chapter 6. (29) Wuerker, R. G. Annotated tables of strength and elastic properties of rocks. Petroleum Transactions Report No. 6; AIME: Englewood, CO, 1969; pp 23−45. (30) Detournay, E.; Cheng, A.-D.; McLennan, J. A poroelastic PKN hydraulic fracture model based on an explicit moving mesh algorithm. J. Energy Resour. Technol. 1990, 112, 224−230. (31) Batchelor, G. K. An Introduction to Fluid Dynamics; Cambridge University Press: Cambridge, U.K., 2000. (32) Barree, R.; Conway, M. Experimental and numerical modeling of convective proppant transport. J. Pet. Technol. 1995, 47, 216−222. (33) Roodhart, L. Proppant settling in non-Newtonian fracturing fluids. SPE/DOE Low Permeability Gas Reservoirs Symposium; Society of Petroleum Engineers: Richardson, TX, 1985. (34) Riddle, M. J.; Narvaez, C.; Bird, R. B. Interactions between two spheres falling along their line of centers in a viscoelastic fluid. J. NonNewtonian Fluid Mech. 1977, 2, 23−35. (35) Gidley, J. L.; Holditch, S. A.; et al. Recent Advances in Hydraulic Fracturing; Society of Petroleum Engineers: Richardson, TX, 1989. (36) Christofides, P. Nonlinear and Robust Control of PDE Systems: Methods and Applications to Transport-Reaction Processes; Springer: Berlin, 2001. (37) Armaou, A.; Christofides, P. D. Finite-dimensional control of nonlinear parabolic PDE systems with time-dependent spatial domains using empirical eigenfunctions. Appl. Math. Comput. Sci. 2001, 11, 287− 318. (38) Hoo, K. A.; Zheng, D. Low-order control-relevant models for a class of distributed parameter systems. Chem. Eng. Sci. 2001, 56, 6683− 6710. (39) Zheng, D.; Hoo, K. System identification and model-based control for distributed parameter systems. Comput. Chem. Eng. 2004, 28, 1361−1375. (40) Zheng, D.; Hoo, K.; Piovoso, M. Low-order model identification of distributed parameter systems by a combination of singular value decomposition and the Karhunen-Loève expansion. Ind. Eng. Chem. Res. 2002, 41, 1545−1556.

Kf = proppant permeability Kr = rock permeability L(t) = fracture length (one wing), m H = vertical fracture height, m P(x,t) = net pressure, MPa U(x,t) = fluid leak-off velocity, m/s Vs(x,y,t) = proppant settling velocity, m/s W(x,t) = fracture width, m ρsd = proppant particle density, kg/m3 ρf = fluid density, kg/m3 μ = fluid viscosity, Pa·s μa(x,y,t) = fluid viscosity accounted for concentration effects, Pa·s ν = Poisson ratio of the formation τ(x) = fracture opening time at position x, s ϕ = proppant bank porosity

■

REFERENCES

(1) Modern Fracturing: Enhancing Natural Gas Production; Economides, M. J., Martin, T., Eds.; Energy Tribune Publishing: Houston, TX, 2007. (2) Warpinski, N.; Griffin, L.; Davis, E.; Grant, T. Improving hydraulic fracture diagnostics by joint inversion of downhole microseismic and tiltmeter data. SPE Annual Technical Conference and Exhibition, San Antonio, TX; Society of Petroleum Engineers: Richardson, TX, 2006. (3) Warpinski, N. Microseismic monitoring: Inside and out. J. Pet. Technol. 2009, 61, 80−85. (4) Differentiating between hydraulic fracturing applications. Presented at the GSA Annual Meeting, Charlotte, NC, 2012. (5) Veatch, R. Overview of Current Hydraulic Fracturing Design and Treatment TechnologyPart 1. J. Pet. Technol. 1983, 35, 677−687. (6) Adachi, J.; Siebrits, E.; Peirce, A.; Desroches, J. Computer simulation of hydraulic fractures. Int. J. Rock Mech. Min. Sci. 2007, 44, 739−757. (7) Mendelsohn, D. A review of hydraulic fracture modeling. I: General concepts, 2D models, motivation for 3D modeling. J. Energy Resour. Technol. 1984, 106, 369−376. (8) Daneshy, A. Numerical solution of sand transport in hydraulic fracturing. J. Pet. Technol. 1978, 30, 132−140. (9) Perkins, T.; Kern, L. Widths of hydraulic fractures. J. Pet. Technol. 1961, 13, 937−949. (10) Sneddon, I.; Elliot, H. The opening of a Griffith crack under internal pressure. Q. Appl. Math. 1946, 4, 262−267. (11) Nordgren, R. Propagation of a vertical hydraulic fracture. Soc. Pet. Eng. J. 1972, 12, 306−314. (12) Howard, G. C.; Fast, C. R. Drilling and Production Practice; American Petroleum Institute: Washington, DC, 1957; Appendix I: Optimum Fluid Characteristics for Fracture Extension, by R. Carter, pp 261−269. (13) Khristianovich, S.; Zheltov, Y. P. Formation of vertical fractures by means of highly viscous liquid. Proc. 4th World Pet. Congr. 1955, 2, 579−586. (14) Geertsma, J.; De Klerk, F. A rapid method of predicting width and extent of hydraulically induced fractures. J. Pet. Technol. 1969, 21, 1571− 1581. (15) Daneshy, A. On the design of vertical hydraulic fractures. J. Pet. Technol. 1973, 25, 83−97. (16) Spence, D.; Sharp, P. Self-similar solutions for elastohydrodynamic cavity flow. Proc. R. Soc. London, A: Math. Phys. Sci. 1985, 400, 289−313. (17) Howard, G. C.; Fast, C. R. Hydraulic Fracturing; Society of Petroleum Engineers: Richardson, TX, 1970; p 210. (18) Advani, S.; Lee, T.; Lee, J. Three-dimensional modeling of hydraulic fractures in layered media. I, Finite element formulations. J. Energy Resour. Technol. 1990, 112, 1−9. 10503


Evaluating the Performance of a Fracturing Treatment Design

Recommend Documents