Numerical Simulation of Coning Using Implicit Production Terms

1970 ◽  
Vol 10 (03) ◽  
pp. 257-267 ◽  
Author(s):  
A. Spivak ◽  
K.H. Coats

Abstract This paper describes the use of a multiphase, multidimensional mathematical model to predict two- and three-phase coning behavior. Severe computational instability in the form of saturation oscillations in grid blocks near the wellbore is commonly encountered in the mathematical simulation of coning. This instability is due to the explicit (dated at the beginning of a time step and held constant for that time step) handling of saturation - dependent transmissibilities and production terms in the finite-difference solution of production terms in the finite-difference solution of the flow equations. An analysis of stability with respect to explicit handling of saturation-dependent transmissibilities is presented in this paper. This analysis shows why explicit transmissibilities can result in a severe time-step restriction for coning simulation. The use of implicit production terms in the difference equations to reduce instabilities is discussed and examples are given. These examples show that the implicit handling of production terms alone can result in a fivefold increase and permissible time step for a coning simulation with virtually no increase in computing time per time step. A laboratory water-coning experiment was simulated and excellent agreement was obtained between computed and observed results. A three-phase coning example for a gravity-segregation reservoir is also presented. Introduction Simulation of coning behavior is normally done by numerically solving the flow equations expressed in cylindrical (r, z, theta) coordinates with symmetry in the theta direction. The finite-difference technique of numerical solution of differential equations requires that the portion of the reservoir being simulated be divided into grid blocks as shown in Fig. 1. Since coning is a well phenomenon and not a gross reservoir phenomenon, the grid blocks must necessarily be relatively small in the vicinity of the wellbore because both pressures and saturations vary rapidly in this region. Severe computational instability is commonly encountered in the simulation of coning due to the relatively small grid-block sizes and high flow velocities in the vicinity of the wellbore. During a time step that would be considered normal for most reservoir simulation problems, a block near the wellbore is required to pass a volume of fluid many times its pore volume. SPEJ P. 257

1973 ◽  
Vol 13 (04) ◽  
pp. 200-210 ◽  
Author(s):  
C.A. Chase ◽  
P.M. O'Dell

Abstract Variational principles stated by, Biot have been applied to obtain a two-parameter (approximation for heat losses to cap and base rock from a reservoir undergoing thermal recovery. The approximation predicts heat losses to within a few percent of the predicts heat losses to within a few percent of the exact value when the beat losses result from one-dimensional conduction into cap and base rock in the direction normal to the reservoir boundary surfaces. Conduction in the longitudinal direction is neglected. Therefore, the approximate temperature distribution is valid only when the temperature gradient in this direction is small. But because the Peclet number (ratio of convective to conductive heat transport) is high in most reservoir thermal processes, the horizontal temperature gradient will processes, the horizontal temperature gradient will be small everywhere except in the vicinity of a thermal front, and the approximation will be valid. Comparison with a finite-difference solution in cap and base rock shows that reasonable accuracy is obtained when the Peclet number is 100 or greater. The variation solution has been incorporated into our thermal simulator and yields a considerable sailings in core storage. It is no longer necessary to store grid-block temperatures for cap and base rock nor to solve the finite-difference form of the energy balance in this region. Instead a system of two nonlinear ordinary differential equations must be solved for each grid block at the interface of the reservoir and the cap rock. In addition to savings in core storage, a reduction in computation time is achieved because fewer finite-difference grid blocks are needed. Introduction Heat losses to cap and base rock must be considered in modeling thermal processes in petroleum reservoirs. Since there is no mass petroleum reservoirs. Since there is no mass transport in the cap and base rock, the only mechanism for heat transfer is conduction. One of the most obvious ways of determining heat losses from the reservoir is to solve the energy equation in the cap- and base-rock region by finite differences. To do this, the reservoir finite-difference grid must be extended into the cap- and base-rock region. This can consume a good deal of computer core storage - at a time when all available core storage is needed to adequately model mass and energy transport in the reservoir region. Furthermore, since there is no mass transport in the cap and base rock, one would like to eliminate having to solve the conservation-of-mass equations in this region, but to do so requires a special computer code. Hence, a finite-difference solution can be costly. It does, however, have the advantage of generality in that a minimum of assumptions is involved in formulating the conservation equations. There are ways of calculating heat losses to cap and base rock other than by finite differences. However, for a method to be competitive with the finite-difference method, it must offer some advantage such as accuracy, reduced computer core storage, or lower computation time. One alternative to finite differences is the use of superposition to couple an analytic solution for the cap and base-rock temperature distribution with the finite-difference solution of the reservoir energy balance. But, during the course of the simulation, the superposition principle would necessitate having temperature data for all previous time steps for each grid block adjacent to the cap and base rock. This requires an appreciable amount of computer core storage, perhaps even more than would be required for a complete finite-difference solution. Hence, this method does not seem attractive. The use of variational principles appeared to offer the advantages of both reduced core storage and lower computation time and was therefore considered as a means of treating heat losses to cap and base rock. The advantage of the variational method is that a priori knowledge of the approximate shape of the temperature profile can be used to choose the functional form of the temperature distribution. The chosen functional form will contain several free parameters. SPEJ P. 200


Author(s):  
T. N. Krishnamurti ◽  
H. S. Bedi ◽  
V. M. Hardiker

Atmospheric models generally require the solutions of partial differential equations. In spectral models, the governing partial differential equations reduce to a set of coupled ordinary nonlinear differential equations where the dependent variables contain derivatives with respect to time as well. To march forward in time in numerical weather prediction, one needs to use a time-differencing scheme. Although much sophistication has emerged for the spatial derivatives (i.e., second- and fourth-order differencing), the time derivative has remained constructed mostly around the first- and second-order accurate schemes. Higher-order schemes in time require the specification of more than a single initial state, which has been considered to be rather cumbersome. Therefore, following the current state of the art, we focus on the first- and second-order accurate schemes. However, higher-order schemes, especially for long-term integrations such as climate modeling, deserve examination. We start with the differential equation dF/dt = G, where F = F(t) and G = G(t). If the exact solution of the above equation can be expressed by trigonometric functions, then our problem would be to choose an appropriate time step in order to obtain a solution which behaves properly; that is, it remains bounded with time. This is illustrated in Fig. 3.1. We next show that: (1) if an improper time step is chosen, then the approximate finite difference solution may become unbounded, and (2) if a proper time step is chosen, then the finite difference solution will behave quite similar to the exact solution. The stability or instability of a numerical scheme will be discussed for a single Fourier wave. This would also be valid for a somewhat more general case, since the total solution is a linear combination of sine and cosine functions, which are all bounded. We next define an amplification factor |λ|, the magnitude of which would determine whether a scheme is stable or not.


Sign in / Sign up

Export Citation Format

Share Document