Search Documents

By A. R. Gregory, M. R. J. Wyllie

The literature reveals that scant attention has been paid to the systematic experimental determination of the formation factors of unconsolidated porous media. No experiments appear to have been made on the effect of increasing cementation on the formation factor of an initially unconsolidated porous medium. Measurements have accordingly been made of formation factors as a function of porosity for aggregates of spheres in the porobity range 12-56 per cent and of cubes, cylinders, discs and triangular prisms in the porosity range 30-45 per cent. The results are examined in the light of the theoretical equations of Clerk Maxwell. Fricke and Slawinski. Aggregates of unconsolidated spheres and beach sands have been artificially cemented with silica and the formation factor-porosity relationships determined. A theory is outlined which seeks to explain the results obtained and which postulates that formation factor and porosity data for consolidated sandatones may be used to determine the original porosity of the unconsolidated sands from which the consolidated sandstones were derived. It follows also that for consolidated sandstones F' = C where F is formation factor, fractional porosity and C and k are constants. INTRODUCTION All schemes of electric log interpretation implicitly or explicitly require the use of formation resistivity factors or? as they are now generally called. formation factors. 'The originator both of the concept of formation factor and of the use of formation factor data in electric log interpretation was Archie.' More recently it has been suggested that formation factors are of value in describing fluid flow in porous media. A full discussion of this problem together with adequate reference.: has been given by Wyllie and Spangler.' In spite of the theoretical and practical importance of the subject, there appear to be remarkably few basic formation factor data in the petroleum literature. A wider search reveals that this situation is nut confined to the literature of the to industry. at various times over the last century, attempts have been made to compute the theoretical relationship between the formation factor and porosity of simple systems of unconsolidated particle.. However. the formulae derived are not all identical and the experimental evidence adduced in support of each i; not wholly convincitng. Most of the experimental work has been carried out on system. of unconsolidated spheres, although some attention has been given to systems of spheroids and cylinder. Other simple shapes have received neither theoretical nor experimental investigation. The writers feel that a better understanding of the formation factors of complex systems of practical importance requires a firm foundation of reliable experimental data applicable to simple systems. Accordingly, they have reviewed the relevant literature and have determined experimentally the relationship between the formation factor and porosity of unconsolidated aggregates of particles of regular geometrical shapes. Some experiments have also been carried out to investigate the effect of cementation on the formation factor-porosity relationship. PREVIOUS WORK Although the term formation factor was not used. the relationship between the resistivity of a solution containing dielectric solids and the resistivity of a similar solution not containing solids is an old and important problem. The formation factor of a system of dispersed spheres as a function of porosity was considered by Clerk Maxwell. Maxwell arrived theoretically at the following expression: where F is formation factor and fractional porosity. Lord Rayleigh' also examined the problem of the formation factor of systems of spheres and cylinders and concluded that for spheres Clerk Maxwell'.; expression gave a value of F whirl]. at low porosities particularly, was too small. However,

Jan 1, 1953

By B. K. Larkin, H. R. Bailey

One of the most complicated and potentially one of the most promising secondary recovery methods is that of underground combustion. A number of field tests1,2,3 have been performed recently, apparently indicating that the underground combustion process is technically feasible. Relatively little has been published concerning laboratory or theoretical investigations of the underground combustion process. A few results of laboratory tube runs have been published4,5 and these data are valuable in understanding the process. The tube run data also emphasize the complexity of the thermal recovery process. Under various conditions, problems in heat transfer, multiphase fluid flow and combustion rates arise. For example, at low temperatures the combustion process proceeds so slowly that the front cannot be propagated. It is a problem in heat transfer to predict what injection rate of air is necessary to sustain combustion in a particular reservoir. Certain fluid flow aspects of underground combustion have been studied by Wilson, et a13. In the development of some of the flow equations, assumptions of adiabatic combustion were made. Vogel and Krueger1 have devised an analog model of the heat conduction process for a constant temperature cylindrical source moving radially outward in a homogeneous conducting medium. Jenkins and Ramey7, in a discussion of the Vogel and Krueger paper, present an analytical solution of a heat conduction problem related to problems considered by Vogel and Krueger. In the present paper an idealized model of the heat conduction problem is considered, which is more general than the model of Vogel and Krueger or Jenkins and Ramey. In particular, initial well heating, vertical heat losses and arbitrary frontal velocities are included. Equations for the temperature are obtained in terms of integrals. These integrals have been evaluated numerically for the case of frontal movement corresponding to a constant air injection rate and the case of constant frontal velocity. These results are presented graphically. From these graphs one may estimate conditions under which nearly adiabatic temperatures are obtained. The application of these results requires some knowledge of reservoir thermal properties, combustion rates and residual fuel content. Assuming no vertical losses it is possible to obtain explicit evaluations of these integrals in the linear con- stant velocity case and in the radial case corresponding to a constant air injection rate. PRESENTATION OF RESULTS The geometry of the radial system is illustrated in Fig. 1 where the reservoir is represented as a slab of thickness h. By making an energy balance for an infinitesimal volume, one obtains the following equation. T= x T+F(z.r).........(1) where F (z,r) =HWv1/pc (r-r1), if z h/2 = 0, if z >h/2 . In Eq. 1, V- represents the Laplacian operator and 8 is the Dirac delta function. Eq. 1 is the well-known diff-sion equation with a source term and may be applied to a linear system. For example, if the source is moving in the direction of the x axis, then the argument of the delta function is x - 1,. In deriving Eq. 1 the following approximations were made concerning the actual combustion process: (1) all heat transfer is by conduction in a homogeneous medium, (2) the formations bounding the reservoir have the same thermal properties as the reservoir, (3) the combustion fuel is consumed instantaneously as the front arrives, (4) there is a constant fuel concentration

By B. E. Eakin, A. L. Lee

Atmospheric pressure viscosity values For light hydrocarbons were calculated by a simple mixing rule based on the Sutherland constants for pure components. This method was further extended so that the atmospheric pressure viscosity of a mixture may be calculated based on the molecular weight information only. INTRODUCTION A proven method for estimating the gas-phase viscosity of a mixture is the one by Herning and zipperer.l This mixing rule has the following form: One method for calculating the pure component viscosity i is the Sutherland equation,2 which has the following form: The Sutherland constants B and S are determined from experimental data. A set of Sutherland constants for light hydrocarbons based on the data of several investigators3-5 is presented in Table 1. APPLICATION The Herning-Zipperer equation and the Sutherland equation were used to obtain an equation for calculating the gas-phase viscosity of light hydrocarbons. This equation has the form, Eq. 3 was used successfully in the study of methane-n-butane mixtures.6 Further analysis indicated that the Sutherland constants are a function of molecular weight.7 The gas-phase viscosity may be calculated based on molecular weight information only. The following equation was developed: The atmospheric pressure viscosity of ethane, propane, n-butaneand four methane-propane mixtures calculated on the basis of Eq. 6 are compared with the data of other investigators in Tables 2, 3, 4 and 5, respectively. A plot of atmospheric pressure viscosity vs molecular weight based on Eq. 6 is presented in Fig. 1.

Jan 1, 1965

By D. W. Peaceman

In previous papers by this author on numerical calculations of multidimensional miscible displacement, some simplifying assumptions were made in writing the dispersion term of the differential equation. It was assumed that the flow vectors were essentially parallel to the x-axis. However, in most multidimensional miscible displacements, the flow vectors are not so simply oriented with respect to the coordinate axes. This paper derives the correct dispersion term /or the more general case; gives a difference approximation /or the dispersion tern; and derives the stability criterion for the corresponding explicit difference eguation.

Jan 1, 1967

By R. F. Jones, N. E. Truitt, M. Prats

In most published discussions and theories of in situ combustion, the combustion fronts are assumed to be vertical. However, evidence from field tests leaves no doubt that combustion fronts often advance more rapidly along the top than near the bottom of a formation as a result of difference in density between injected air and formation liquids. The approximation proposed in this paper to determine the movement of the resultant tilted combustion surfaces states that the vertical rate of movement of combustion surfaces is proportional to the horizontal oxygen flux. Where self-ignition is possible, the proposed method demands that a secondary combustion surface exist around production wells which produce some oxygen. These secondary combustion surfaces may be formed long before the primary combustion surface can advance from injection to production wells. Heat liberated near production wells at these secondary combustion surfaces can contribute to an early increase in production rate. Results indicate that significant oil recoveries cannot be obtained from the usual flood patterns (five-spots, seven-spots, etc.) without producing large volumes of unused oxygen. Ideally, to increase oxygen- consumption efficiency, well patterns should allow oil production from a first line of production wells and gas production from more distant lines of producers. However, it may be desirable to produce some gas at all wells to support (and benefit from) active secondary combustion surfaces. Results indicate that the well spacing through which combustion can be advanced is larger than that predicted by other methods. A large number of production wells may still be desirable to take quick advantage of gravity drainage. From a comparison with results at South Belridge field, California, it appears that this method adequately describes oxygen concentration and temperature histories and combustion-front shapes. However, this method does not accurately locate the most advanced point of the combustion surface. There is some field evidence to substantiate the actual presence of secondary combustion surfaces at South Belridge. Use of the proposed method appears warranted at this time when lay-over of the combustion surface can be anticipated. INTRODUCTION The assumption of vertical combustion fronts has been embodied in all previous publicationsl-6 which use the movement of combustion fronts away from injection wells to determine the temperature and fluid distributions in the reservoir. The only paper concerned with a mathematical model of the combustion process in which a nonvertical combustion front is used was written by Gottfried.11 Actually, nonvertical combustion fronts have been observed in most in situ combustion field tests for which adequate data are available.7, 8, 6 In practice, the typical vertical extent of the burned zone decreases with distance from the injection well, and this burned zone is at or near the top of the sand body. In some cases, such as at South Belridge field near Taft, Calif., the combustion surface is almost horizontal over a very large area.8 Thus, for some years an obvious and serious gap has existed between theory (vertical fronts) and practice (tilted fronts). This is indicated in Fig. 1. Tilted combustion fronts such as observed at South Belridge8 sometimes result from the natural tendency of injected gases to rise to the top of an oil sand. An initially present gas cap, even if thin, would make apparent the effect of gravity segregation from the first, particularly since formations considered for in situ combustion field operations are generally selected for their

Jan 1, 1969

By R. E. Cook

This paper presents results of a study to determine to what extent errors in estimated free gas saturation affect the results of static pressure calculations from build-up curves in two-phase systems. Use is made of the method of pressure build-up curve analysis developed by Miller, Dyes and Hutchinson1 , with modifications for two-phase flow advanced by Perrine2. Cornhined fluid compressibility of gas-oil mixtures is included, with special attention given to the discontinuity that occurs in fluid compressibility at the saturation pressure and its effect on static pressure calculations from build-up curves. The sensitivity of static pressure calculations to errors in gas saturation is considered over wide ranges of formation and fluid properties and reservoir pressure conditions. It is concluded from this study that for analyses of two-phase pressure build-up curves (I) the use of under-saturated crude compressibility data, even when producing GOR's are at or near solution values, will result in very significant negative errors in calculated static pressrrre, and (2) except for large values of buildup curve .slope and/or very low crude gravities, rather large errors in calculated gas saturation can he tolerated without .significant effect on calculated static pressures. The conc1usions arc qualified as referring only to the effects of gas saturation and excludc possible errors in the nuumerous other factors that enter into the ancllysis of pressure build-up curves. INTRODUCTION As part of the analysis of basic data for several major Eocene reservoirs in the Bolivar Coastal field (Western Venezuela), it was necessary to analyze a large number of pressure build-up curves for oil permeability and static reservoir pressure. From performance data it became apparent that gravity segregation is very active in these reservoirs. As a result of this performance it was recognized that the free gas saturation within the oil leg of the reservoirs would be small but would be subject to certain errors in calculation due to the inaccuracy of relative permeability data in the low gas saturation region. Since all of the various methods of build-up curve analysis were originally derived for sin- gle-phase, constant compressibility systems, and since the presence of a small but rather indeterminate gas saturation would greatly affect the reservoir fluid compressibility, a study was undertaken to investigate the effect of gas saturation on the results of two-phase pressure build-up analyses. TWO-PHASE BUILD-UP CURVE ANALYSIS The method of analysis used in this investigation is basically that originally developed by Miller, Dyes and Hutchinson1, with modifications for two-phase flow advanced by Perrine2. Miller, et al, showed that the original equation for single-phase permeability away from the well bore, held for two-phase flow except that the calculated permeability is the effective permeability to oil away from the wellbore. This finding was acknowledged by Per-rine, who also presented the equation for effective gas permeability* as Since Eqs. 1 and 2 contain no terms influenced by fluid compressibility, no distinct problem should develop in calculating permeability from two-phase build-up curves. In regard to the calculation of static reservoir pressure from build-up curves in single-phase systems, the modified Miller, et al, method as presented by Perrine is briefly as follows. 1. From build-up curve data, PVT data for the reservoir fluid and the well spacing and sand thickness, calculate the climensionless shut-in time corresponding to the highest value of actual shut-in time (t, hrs) on the straight-line portion of a semi-log build-up curve. Dimenrsionless shut-in time is related to actual shut-in time by the relation, t- 0.000264 kt 2. From a graph of dimension less shut-in time vs dimensionless pressure difference (Fig. 1 1 ) find the value of dimensionless pressure difference corresponding to the value of 7 calculated under (1) above. Dimen-

By R. W. Nelson

The theoretical basis is presented for energy dissipation methods of measuring permeability in saturated heterogeneous media. Analysis starts with the equation for single-phase flow in heterogeneous porous media which is a first-order partial differential equation in the unknown permeability. This equation is reduced to a system of characteristic equations useful in deriving a differential expression for the permeability as a function of the known potential distribution. For steady flow systems the expression is integrable, giving the permeability integral for direct calculation of the distribution along successive streamlines. The boundary condition in permeability is presented and includes a discussion of the special requirement to assure uniqueness. Theoretical results and necessary computations are incorporated into two general computer programs capable of determining the permeability distribution in two- and three-dimensional steady flow systems. Results are then presented for an energy dissipation analysis using the methods and soltware to determine a field permeability thickness distribution. INTRODUCTION Theoretical work leading to the description of fluid flow in heterogeneous porous media has developed slowly and somewhat sporadically, although the areas of engineering applications predominantly involve flow in nonhomogeneous media. The lack of accurate and economical methods for measuring field permeability distributions has been a major limitation to realistic description of natural flow systems. The economical measurement of permeability using energy dissipation methods may overcome this limitation. Among other things, any permeability measurement method involves consideration of energy dissipation; as used here, "energy dissipation method" denotes the following general approach. The areal energy distribution is measured on an existing flow system in porous media. Then through appropriate analysis of the dissipation of energy in the flow system, the permeability distribution is determined throughout the region using the minimum boundary conditions required to assure a mathematically unique result. This paper presents the theoretical basis for such a method, as well as the techniques and computer software needed to carry out the analysis. The method and techniques are then used in a preliminary field application. EQUATIONS DESCRIBING FLOW IN HETEROGENEOUS MEDLA Heterogeneous as used here implies differing in kind, having unlike quantities or possessing different characteristics. Accordingly, a heterogeneous porous medium has unlike quantities or differing characteristics at different locations. More precisely, if a porous medium is heterogeneous with respect to some property, then that property is functionally dependent on the spatial location. (Media properties are considered here as macroscopic.) Specifically considering permeability k,** then k = G(x, y, z) ,............(1) or the permeability is a scalar function G of the location (x, y, z). Spatial variation of several other properties, of which porosity, dispersivity and medium compressibility are but a few, is implied by the broad title of heterogeneous porous media. The term will be used hereafter to imply only heterogeneity with respect to permeability. The porous medium is considered isotropic as implied by G being only a scalar function. Anisotropy

Jan 1, 1969

By C. K. Eilerts

Finite difference equations were programmed and used to integrate the second-order, second-degree, partial differential equation with variable coefficients that represents the transient linear flow of gas-condensate fluids. Effect was given to the change with pressure of the compressibility factor, the viscosity, and the effective permeability and to change of the absolute permeability with distance. Integrations used as illustrations include recovery of fluid from a reservoir at a constant production rate followed by recovery at a declining rate calculated to maintain a constant pressure at the producing boundary. The time required to attain such a limiting pressure and the fraction of the reserve recovered in that time vary markedly with properties of the fluid represented by the coefficients. Fluid also is returned to the reservoir at a constant rate, until initial formation pressure is attained at the input boundary, and then at a calculated rate that will maintain but not exceed the limiting pressure. The computing programs were used to calculate the results that would be obtained in a series of back-pressure tests made at selected intervals of reservoir depletion. If effect is given to the variations in properties of the fluid that actually occur, then a series of back- pressure curves — one for each stage of reserve depletion — is required to indicate open-flow capacity and related flow characteristics dependably. The isochronal performnce method for determining flow characteristics of a well was simulated by computation. INTRODUCTION The back-pressure test procedurel is based on a derivation of the equation for steady-state radial flow of a gas, the properties of which are of necessity assumed to remain unchanged in applying the test results. The properties of most natural gases being recovered from reservoirs change as the reserve is depleted and pressures decline, and the results of an early back-pressure test may not be dependable for estimating the future delivery capacity of a well. The compressibility factor of a fluid under an initial pressure of 10,000 psia can change 45 per cent and the viscosity can change 70 per cent during the productive life of the reservoir. There are indications that the effective permeability to the flowing fluid can become 50 per cent of the original absolute permeability before enough liquid collects in the structure about a well as pressure declines to permit flow of liquid into the well. The advent of programmed electronic computing made practicable the integration of nonlinear, second-order, partial differentia1 equations2.3 pertaining to flow in reservoirs. Aronofsky and Porter4 represented the compressibility factor and the viscosity by a linear relationship, and integrated the equation for radial flow of gas for pressures up to 1,200 psi. Bruce, Peaceman, Rachford and Rice5 integrated the partial differential equations for both linear and radial unsteady-state flow of ideal gas in porous media. Their published results were a substantial guide in this study of integration of the partial differential equation of linear flow with coefficients of the equation variable. The computing program was developed to treat effective permeability as being both distance-dependent and pressure-dependent. In this study all examples of effective permeability are assumptions designed primarily to aid in developing programs for giving effect to this and other variable coefficients. The accumulation of data for expressing the pressure dependency of the effective permeability is the objective of a concurrent investigation. FUNDAMENTAL DIFFERENTIAL EQUATIONS DIFFERENTIAL EQUATION OF FLOW Linear transient fluid flow is represented by the partial differential equation,

Jan 1, 1965

By E. F. Sumner, C. K. Eilerts

A partial differential equation was developed and programmed to compute the cumulation and flow of the liquid phase of a gas-condensate fluid. It was found that such cumulation near the wellbore may attain saturations permitting liquid flow within periods of the order of 1 day. The radius at which this saturation obtains advances toward the external reservoir boundary with continued recovery, although it may move only 2 to 5 percent of that distance during the productive life of the well. This relatively high saturation can cause a serious dimunition in delivery capacity. The pressure dilferential required for flow of gas through a formation with such saturation may be as much as three times greater than the corresponding differential in a formation where saturation is less than 2 percent. Because of the radial geometry, total reservoir pressure differential required for flow with a zone of relatively high liquid saturation about the well was indicated to be two to three times greater than without the saturation. INTRODUCTION A computing program developed by the authors 1 to integrate the partial differential equation for transient radial flow of gas-condensate fluids required a prevalued function, W(Po) = kD(p0)/M(p0)z0(Po).......(1) related to cumulation of condensate in pores of the reservoir as pressure declines. Methods were available for relating mobility to saturation, saturation to liquid-vapor volume ratio, liquid-vapor volume ratio to pressure and, hence, mobility to pressure to obtain W(pD). However, certain assumptions had to be made in advance of program execution concerning movement of the liquid phase and change in radial distribution of liquid saturation as recovery progressed. Because reservoir fluid mobility and related pressure gradients are markedly dependent on liquid saturation, it was desirable to eliminate the necessity for these assumptions and in this respect improve the reliability of a result obtained by integration. Experimental research was begun to determine the relationship of variables that control the reservoir flow of low liquid-vapor ratio fluids. With the porous medium partly saturated, permeability to a two-phase flowing fluid could not be measured because the viscosity was unknown. Viscosity could not be measured without knowing permeability of the liquid-containing medium; but mobility, the ratio of permeability to viscosity, was measurable. This property, and liquid-vapor volume ratio of the flowing phases, have values determined by nature of the fluid and porous medium and, in addition, were found to be pressure, velocity and saturation dependent. A program developed by Kniazeff and Naville 2 for computing multicom-ponent, two-phase flow used polynomials to represent saturation dependence of permeabilities for gas and for liquid phases, and an effluent viscosity constant under reservoir conditions. The authors developed a program based on a partial differential equation of similar mass balance origin but designed to represent distribution of an unlimited number of components between gas and liquid phases utilizing the experimental information available as to variable dependence. Attention was given to liquid-vapor volume ratio of the flowing phases, a property which is essential as a liquid saturation boundary condition. A program for computing implicitly saturation increase at uniform intervals along the drainage radius was combined with an existing program for calculating pressures, and the two were used alternately in an iterative solution for both saturations and pressures. The function W(pD) needed for calculating pressures was replaced by the ratio of mobility to compressibility factor. BASIC EQUATIONS The partial differential equation for calculating

By L. F. Elkins

In 1946 Littlefield, Gray and Godbold published a thorough geologic description of the West Edmond Hunton Lime reservoir, located in Central Oklahoma, and discussion of its early performance.l They stressed implications of the nature of the reservoir rock on performance to be expected with production both by natural means and by pressure maintenance with gas injection. The essence of their analysis was that presence of an extensive interconnected system of fractures and solution channels constituting about 10 percent of the total reservoir void space in otherwise very low permeability reservoir rock would result in severe channeling of naturally encroaching water or injected gas through the fractures with little or no benefit to ultimate recovery of oil. Considerations were in progress to unitize the reservoir for pressure maintenance by gas injection. Some other engineers involved in the project had considerably different opinions at the time regarding continuity of fractures and regarding effects of such fractures on performance. The field was unitized Oct. 1, 1947. Subsequently, pilot gas injection tests were conducted in four areas of the field; large scale waterflooding was conducted in two parts of the field in the Bois d'Arc or Upper Hunton; and pilot waterflooding was conducted in the Chimney Hill or Lower Hunton. Full scale water-flooding is now in progress in the Chimney Hill section where it is oil productive. The primary purpose of this paper is to divine some of the secrets of this reservoir through analysis of all aspects of its performance, with the hope that it will be of help to engineers and operators in developing and operating other low permeability fractured reservoirs. In part, this considers variance of actual performance from some widely used conventional reservoir engineering assumptions. For brevity, other aspects of unit operation and its benefits to the owners are not discussed. Oil and Gas in Place For unit participation, net pay in the Bois #Arc defined by the gross SP anomaly on electric logs totalled 2,100,000 acre-ft in the 30,000-acre reservoir. Development of West Edmond predated modem logging techniques that now would be used to determine "net" or "productive" pay. In Fig. 1 this log pick is compared with detailed core analysis and core descripfion by Littlefield et a1.l for one well identified as A on the index map in Fig. 2. Not all of the gross Bois #Arc section had oil staining, but all cores with measureable intergranular permeability were stained. (The committee selection of the Bois d'Arc SP interval also includes the Frisco section as identified from cores by Littlefield.) The Bois d'Arc is oil productive throughout the entire 600-ft monoclinal accumulation from the eastern upstructure erosional pinchout to the west-side oil-water contact. (See Littlefield et al.l and McGee and Jenkins for geologic detail.)

Jan 1, 1970

By H. S. Price, D. A. T. Donohue

The system of equations describing displacement of a hydrocarbon liquid by a hydrocarbon vapor in a porous medium where mass transfer takes place between the phases is solved numerically for a variety of gas injection processes. Even though the method of solution is quite general, only systems with three hydrocarbon components are considered. Computer simulations of displacement processes wherein mass transfer between phases is both considered and neglected are compared, and it is shown that neglecting mass transfer can give pessimistic displacement efficiencies. INTRODUCTION The role of the gas displacement process in the recovery of petroleum has been subjected to a series of detailed analyses; as a result, a number of predictive models have been published in the literature. However, because of major simplifying assumptions, most of these models do not completely represent the physical system. As a result, the effect of making the simplifying assumptions is unknown. Therefore, a complete representation of this process — one without major simplifying assumptions —should lead to a fuller understanding of the process, and perhaps to methods of improving it. The general method of developing a model for two-phase fluid flow in a porous medium is to solve simultaneously the continuity equation, the energy equations and the equation-of-state for each phase under the prescribed initial and boundary conditions. For an isothermal system, the energy equations reduce to the momentum equation. Darcy's law. However, since natural gas is the vapor state of the reservoir liquid, interphase mass transfer may take place with concomitant changes in both the intensive and extensive thermodynamic properties of each phase. It is this phenomenon that has often been omitted in previous mathematical models. An additional relation, then, which accounts for mass transfer between the phases, must be included with the other equations to specify a complete model. Completely formulating the equations to be solved is not a difficult task1 but obtaining their solution has been intractable up to now. Availability of large-memory, high-speed digital computers now makes an attack on this formidable problem possible. This paper presents a preliminary study of the problem. Since this investigation is intended to be exploratory, it is restricted to the linear, horizontal, isothermal, two-phase viscous flow of oil and gas in an oil reservoir. In the early development of predictive models of this process, the reservoir system was considered as a unit and various forms of the material balance equation were proposed.2-4 Pressure and saturation gradients were then added in the Buckley-Leverett model.5 The Buckley-Leverett formulation considered the fluids to be incompressible; thus, the mathematical model reduces to a steady-state system. In the 1950's, studies incorporating numerical techniques were being published.6-9 These mathematical models differed in the efficiency of finite difference techniques, the inclusion or exclusion of capillarity or the number of space dimensions considered. To solve these nonlinear, partial differential equations, each phase was considered to be homogeneous with time; therefore, mass transfer between phases was neglected. The effect of mass transfer on the gas displacement process was first reported by Attra. 10 He simulated the one-dimension flow system by a series of cells in each of which the fluids were equilibrated during a time step. In addition, the pressure throughout the system during each time step was predetermined and constant-phase velocities were calculated according to the Buckley-Leverett incompressible fluid flow model. Welge et al.11 developed a model for the displacement of oil by an enriched gas where composition is considered to be a dependent variable. Fluid

By R. H. Barham, G. W. Nabor

Linear aquifers, either limited or essentially infinite, may be encountered in reservoir engineering practice. In areas where faulting fixes reservoir boundaries, the fault block reservoir may have an aquifer of limited extent whose geometry is best approximated as linear. An infinite linear aquifer can occur as a regional feature whenever water movement through the aquifer member is constrained to one direction. Such constraints can arise from major faults, facies changes or pinchout of the member. Miller* pointed out that linear aquifers have received only meager attention in the past. He analyzed the performance of finite and infinite aquifers, developed working equations and curves, and presented examples. While Miller's curves may be used fairly easily, a separate one is required for each size of aquifer. In this paper, Miller's equations have been used as a starting point. By modifying them slightly, they can be reduced to a form which yields a single working curve, applicable to any size of aquifer. Thus. interpolation between curves is eliminated and accuracy is improved. Miller's results for finite aquifers covered only the boundary condition of no flow across the outer aquifer boundary. This paper also includes the case of constant pressure at the outer aquifer boundary.

Jan 1, 1965

By H. Kazemi

A pressure fall-off test on the injection well of a forward combustion project may permit us to calculate the distance to the burning front. In the mathematical description and analysis of test data, we have considered the effect of temperature on the thermodynamic properties of reservoir fluids. The temperature distribution predicted by Chu' for a combustion wave has been used in analyzing a field test. Numerical solutions were obtained by simultaneous solution of the flow equation, the continuity quation and the equation of state. Although this analysis conccntrated on isotropic systems, the results can be extended to systems having directional permeability. INTRODUCTION Conventional pressure fall-off tests have been used extensively to determine reservoir characteristics in gas- or water-injection projects.- such tests consist of injecting fluid into the formation at a constant rate until steady state is nearly achieved, then discontinuing injection and measuring the bottom-hole pressure as a function of time. A plot of pressure vs the logarithm of time (measured from the instant of shut-in) yields a smooth curve characterized by three distinct regions (Fig. 1): Region I reflects the permeability alterations in the vicinity of the wellbore; Region II reflects the reservoir characteristics further away from the wellbore; and Region 111, the characteristics in the displaced fluid bank. The curves in Regions II and III are straight lines! the slopes of which are proportional to the reciprocal of the flow capacity of the formation behind and ahead of the interface as shown by Eq. 1, respectively:

Jan 1, 1967

By W. N. Hiatt

The mathematical theory of the flow of two-phase, incompressible fluid through porous media is clarified, and the development of a general fluid displacement equation for two-phase, incompressible vertical flow through porous media is outlined. The relationship between the three components of the equation — the "fluid drive effect", the "gravity effect", and the "capillary effect'' — is discussed, as well as implications of the analysis in connection with the study of certain types of petroleum reservoirs. INTRODUCTION The purpose of this paper is to clarify and summarize the mathematical basis of the theory of two-phase incompressible fluid flow through porous media, particularly as it pertains to vertical flow, and to present a qualitative* exposition of the controlling factors involved in circumstances where this type of flow may exist. Although the generally accepted theory of heterogeneous fluid flow through porous media was developed a number of years ago by Muskat,1 and Buckley and Leverett2, 3 using the basic premise of Darcy, (and since then expounded upon by many others) there is a need for an over-all summary and clarification of this theory, especially regarding certain types of applications in connection with so-called "gravity-drainage-type reservoirs". The source and derivation is outlined for what is referred to here as the "General Fluid Displacement Equation for Incompressible Vertical Two-Phase Flow through Porous Media''. This equation would apply under the condition where the pressure drop along the path of flow is negligible in comparison with the pressure level on the system, or where compressibility and solubility effects are negligible. These conditions exist ordinarily in gravity-drainage-type reservoirs. Specifically, it will be shown that this general equation may be developed directly from the basic equations, which, according to present theory, govern the general case of viscous heterogeneous flow through porous media. This general equation has appeared in the literature previously, but in most cases its development from basic theory is not indicated clearly, or the equation is presented in terms of groups of constants and variables involving reduced or dimensionless quantities. This presentation is expedient in cases of experimental or computer applications; but it may make it difficult for the practicing petroleum engineer to identify all the factors involved or to understand the relationships between them. Therefore, in the following development, only ordinary mathematical parlance will be used, and no terms involving groups of other equations or reduced or dimensionless quantities will be used. DERIVATION The following equations for two-phase oil and gas flow result from analogy with Muskat's general equations,** if solubility and compressibility effects are neglected, and if all volumes are considered in terms of reservoir conditions, and the water phase, if present, is considered to be immobile (or effectively a part of the porous medium.) In this case, the sum of the oil and gas saturations may be taken as unity, and the saturation of each phase is in terms of the total hydrocarbon saturation. Note that in Eqs. 1 and 2, and in the derivation following, z is taken as positive in the vertically downward direction,

Jan 1, 1969

By G. T. Davis, C. C. Mattax, M. O. Denekas

The effects of crude oil cornponents on the wellabil-ities of sandstone and limestone were investigated. Fractions containing cornponents differing in molecular weight and molecular structure were obtained from crude oils by distillation, extration and chromatography. Individual fractions were then tested for their effects on rock wettability. Tests indicate that sundstone wetta-bility may he changed by a complex variety of surfactants varying both in molecular structure and molecular weight. Limestone appears to be particularly sensitive to basic, nitrogenous surfactants. INTRODUCTION Investigations in recent years have shown that petroleum reservoir rock wettability can exert a significant influence on the efficiency with which oil can be produced by water flooding. While most reservoirs are presumably water-wet, they niay range in their degree of water-wettability from near-neutral to strongly water-wet.'" Reservoir wettabilities other than strongly water-wet are likely to be induced by adsorption of surface-active components froni the crude oil on the pore walls of reservoir rock.:' Little is known, however, about the nature of the surface-active materials which are likely to be adsorbed by the reservoir rock. Due to the complexity of crude oils. attempts made in the past90 isolate these surface-active components have met with only limited success. It is probable that many different types of surface-active materials arc indigenous to crude oils and that many of these may be adsorbed to varying degrees by reservoir rock. This was cxolored in the studies discussed in this paper. The over-all objective in these studies is to ascertain whether the wettability of a given reservoir can be determined by examining the surfactant content of the reservoir crude. To this end, crude oils were examined to determine the variability of indigeneous surfactants with regard to chemical type and molecular weight. Crude oils were separated by distillation into fractions differing principally in molecular weight, by chroma-tography into fractions containing compounds differing in polarity, and by solvent extraction into nitrogenous and non-nitrogenous fractions. Individual fractions were then tested for their effects on the wettabilities of sandstone or limestone rock samples. EXPERIMENTAL PROCEDURES Fractionation of the Crude Oils Samples of Miocene, Eocene and Jurassic crudes were distilled at temperatures not exceeding 200°C. The final stages of distillation were completed in a molecular still at pressures down to three microns of mercury. Fifteen to 30 fractions were obtained from each crude oil. These cuts were sufficiently broad that separation can be considered to have been effected principally on the basis of the molecular weights of the constituents of the crude oil. A considerable portion (20 to 40 per cent) of the crudes would not distill under these conditions. The residues were recovered and tested with the other fractions. Fractions differing in polarity were separated from a crude of Pennsylvanian age and an extracted sample of Miocene oil by chromatography, using a solid adsorbent. Since surfactants are, for the most part, polar compounds, chromatography should separate many of the surfactants from the crude oil. Such a separation should provide fractions containing compounds differing in molecular structure. Nitrogeneous compounds were extracted from Miocene crude oil with a solution of sulfuric acid in meth-anol. The residual oil was further processed by chonia-tography. Each of the fractions obtained by thesc procedures was dissolved in a non-polar solvent (xylene) and diluted to its original concentration in the crude oil. No attempt was made to maintain an anaerobic atmosphere above the samples while they were being dissolved. These solutions of the fractions were then tested for their effects on the wettability of sandstone and limcstone as discussed in the next section. Measurement of the Effects of Crude Oil Fractions on Rock Wettability No entirely satisfactory method for measuring rock wettability has yet been developed. All methods used are empirical. The imbibition test was used in these studies. This test is based on the tendency of a rock to imbibe the wetting phase spontaneously. For example, if a strongly water-wet rock is first saturated with oil and then placed in water, the water will quickly invade the rock by capillarity and much of the oil will be displaced. If the rock is slightly water-wet, water irnbibition will proceed more slowly and, in many instances, considerably less oil will be displaced. A water-saturated, oil-wet rock will imbibe oil. The initial rate with which water (or oil) imbibition takes place indicates, qualitatively, the degree of water (or oil) wettability of the rock.

By G. E. Perry, J. F. Bruskotter, D. G. Russell, J. H. Goodrich

The depletion performance of gas wells has been investigated by mathematical simulation techniques. The gas well model which was studied consisted of a single well located in the center of a bounded, cylindrical, homogeneous reservoir. Dependency of gas compressibility and viscosity on pressure was considered in studies of well performance at both constant mass flow rate and constant flowing pressure conditions. To carry out the investigation, the nonlinear, second-order, partial differential equation which describes Darcy flow of a nonideal gas through porous media was solved numerically. Some of the previous investigations of gas well performance have been of limited general use, because assumptions were introduced to simplify either the gas properties or the basic differential equation. Other studies have been rigorous in these respects but have presented a very limited set of calculated results. The present study was attempted to present a rigorous theoretical model and sufficient numerical results to permit meaningful conclusions to be drawn. It was found that all terms must be retained in the partial differential equation to make accurate predictions. The neglect of higher-order terms, e.g., terms of the order of the "gradient squared", leads to serious material balance errors at large times and to conservative estimates of gas well performance. The higher the gas flow rate and/or the lower the permeability-thickness product of the formation, the more pronounced are these deviations. For example, in a well draining 640 acres in a 25-md-ft formation (8,120 MMcf gas in place) at a constant rate of 993 Mcf/D, the rigorous solution predicts a bottom-hole pressure decline from 4,000 to 1,000 psia in 8.7 years. If higher-order terms are neglected in the differential equation, this decline in pressure is predicted to occur in 5.3 years. With the results of the numerical solution of the differential equation as a basis, simple, easy-to-use approximations for predicting gas well performance for Darcy flow conditions have been developed. These simple approximations are based on the familiar equations for flow of a single, slightly compressible fluid. The approximate methods possess a high degree of accuracy and enable the prediction of long-term gas well performance to be made quickly and accurately without the use of a digital computer. Both transient and stable flow period approximations were developed. INTRODUCTION In recent years income from the sale of natural gas and associated products has represented an ever-incre as ing fraction of the industry's total revenue from operations. To meet the surge in demand for natural gas, the industry has depended heavily upon established reserves and has actively pursued development of new reserves. The search has progressively led to reservoirs which yesterday were too tight and/or deep to yield the desired return on invested capital. More than ever before, evaluation accuracy is now required to forecast the criteria upon which engineering recommendations and management decisions are based. Considerable effort has been expended by both research and operations personnel on the development and application of methods for analyzing and predicting the performance of gas wells. Fundamentally, the problem is the familiar one of extracting data during the drilling, testing and early production life of a well and applying these data within an accurate simulation model to predict long-term behavior. During the past 30 or more years a voluminous literature dealing specifically with gas field problems has been generated. A recent book' lists a comprehensive bibliography of published material through 1959. Over 1,200 references are cited. Since then 39 additional articles on natural gas technology have been published in Transactions volumes of the Society of Petroleum Engineers of AIME. Most existing theory for predicting gas well performance requires that one or more idealizations (e.g., steady-state flow, ideal gas of constant viscosity, small and constant compressibility and constant-viscosity fluid) be applied. Although existing theory may apply directly or be adapted by various artifices to describe specific gas well and reservoir behavior, no widely applicable method is available, and existing methods appear to be subject to appreciable error unless better limits of applicability are defined. The objectives of this paper are (1) to present numerical solutions to the partial differential equation describing gas flow under conditions of general interest in gas well performance prediction work, (2) to present solutions which possess a high order of accuracy for both transient and stabilized flow periods of a well producing at constant rate or constant pressure, and (3) to develop and present simplified relationships which can be used as high-order approximations to the exact numerical results for fast and accurate predictions of gas well performance at the operating level. Combined, these objectives are designed gen-

Jan 1, 1967

By E. L. Anders

The Mile Six pool is located on the La Brea-Parinas Cullcession of International Petroleum Co., Ltd., in northwestern Peru on the west coast of South America. The reservoir pressure in this pool has been maintained within 200 psi of its initial value throughout its history, and gravity drainage has played an important role in the production behavior. It has now produced 95 per cent of its estimated ultimate recovery. It is estimated that this interesting oil pool will ultimately produce 67 per cent of the initial oil in place and that the resulting residual oil saturation may be as low as 19 per cent of the pore volume (29 per cent of the hydrocarbon pore volume). An evaluation of reservoir rock and fluid characteristics and ultimate oil recovery is presented. INTRODUCTION This study of Mile Six pool was made to evaluate its performance according to latest available information. The production performance of this pool has been discussed in various articles in the past. and the reported behavior has been used as an example for application of computation procedures for gravity drainage depletion' and as an illustration of field behavior under gravity drainage or expanding gas cap drive. There have been wide variations in reported values of initial oil in place, reservoir oil volume factor. connate-water sauration, volume of effective sand, and ultimate recovery because of the paucity of reliable basic data. These various factors have been determined as accurately as practicable with the latest available information, and this evaluation is presented herein. The production history of Mile Six is an excellent example of gravity drainage depleion with effective pressure maintenance by gar injection. GENERAL Mile Six pool was discovered by cable-tool drilling in November. 1927. when well 1996 was completed in the Parinas sand. After slow development with cable tools and sporadic production. the pool was opened to continuous pro(iuction in November. 1933. and develpment was completed with rotary rigs. Pressure maintenance was started in December, 1933. by returning gas to upstructure wells. Most of the development was 'completed by 1937, but some additional wells were drilled in the period 1939-1947. and several old wells were deepened. A total of 46 oil and gas wells and 4 dry holes were drilled on approximately 7-acre spacing. Of the producers. 21 are now flowing. 2 are pumping. 44 are gas input wells. 3 are abandoned. I is a gas well shut in. and 15 are shut in because of non-commercial production or high gas-oil ratio. The locations of all wells are shown on the map of Fig. I. Total oil production on Dec. 31. 1952, was 30,867,373 bbl: cumulative gas production was 22,023,777 Mcf; and 26,410,946 Mcf of gas had been returned to the reservoir. These figures do not include oil and gas lost ill a blowout in January. 1940. GEOLOGICAL DESCRIPTION Mile Six pool is located on the northern end of a structural spur projecting from the La Brea-Negritos uplift.' The spur is probably a reflection of a basement structure. It plunges gently to the north, i broken into a complex series of fault blocks. and contain.; the Verdun Alto. Section Sixteen. and Mile Six pools. The Parinas handstone (lower Eocene!. which is the producing formation in Mile Six. occurs at an average depth of 2,200 ft in the pool and dips north and east at from 15' to 20". The pool covers an area of approximately 350 acres. Mile Six is downfaulted about 600 ft from Section Sixteen pool to the soutb. and a major fault forms its western boundarv. The north and east boundaries are formed by the intersection of the sand top with the water-oil contact which occurs at approximately 2.440 ft subsea. An original gas-oil contact probably existed at about 1.875 ft subsea. Fig. 1 presents the latest structural interpretation of the pool. and Fig. 2 is an isopach map showing thickness of the total Parinas formation above the original water-oil contact. The heavy lints of Fig. 1 are contours on the sand top, and the fine lines are contours on the fault planes. This type of straight-line structural map was developed 1)) International's geologists to reflect structural conditions where the bedding planes dip and have no curvature. 'The La Brea-Parinas Concession is highly faulted by normal fault.. The beds are flat wherever exposed. The Parinas formation is approximately 635 ft thick. and it is etimated that 62 percent II the formation is effective sand. The original oil zone was about 565 ft thick. Fig. 211 presents an electric log showing typical Parinas sand development ill Mile Six pool. The Parinas band in Mile Six i. a well-sorted. medium- to-coal-se-grained. cross-bedded sand with minor lenses of shale and small lenses and pockets of pebble conglomerate. The sand grains are subangular to rounded and consist chiefly of quartz with feldspars. biotite hornblende, and augite as accessory minerals. Because of faulting of the Parinas- formation to the east and north of the pool. there is probably l possibility of a significant, natural water drive in Mile six. The faults within the pool. as indicated in Figs. 1 and 2. are of smaller disI,laceInent and seem to act only , partial barriers to fluid movement within the reservoir. RESERVOIR CHARACTERISTICS Core analysis data are available from five wells. The data were obtained from three wells (Nos. 3401. 3586 and 3719) at the time of their completion and from well 1996 when the original liner was sidetracked and the well was deepened ill 1946. Data from well 2779 were obtained in 1943 from old cores taken when the well was deepened in 1934. From these core analyses. the average porosity was estimated to be 22.6 per cent. and the average permeability to dry ail. was estimated to be 780 and Measured productility indices varied from 3.1 to 71.4 B/D per psi differential. Specific productivity indices varied approximately from 0.1 to 0.3 B./D per psi per ft of sand.

Jan 1, 1953

By H. H. Rachford

Numerical solutions of immiscible flow problems in which dispersive effects of capillarity are dominated by convection require excessively fine grid spacing with attendant high computing costs. The use of coarser spacing reduces cost but often produces oscillation or undue dispersion associated with displacement fronts, A numerical formulation is proposed here which should be applicable to two-dimensional flow problems. it is in part analogous to an approach previously tested for miscible systems. The convective transport is approximated using a change of variables to yield a coordinate system moving approximately with the local characteristic velocity. The capillarity-induced dispersive terms in the differential system describing the process are approximated with respect to a fixed coordinate system by the usual implicit formulation. One-dimensional tests of the procedure yielded results in which the saturation profiles tended smoothly to the zero-capillary pressure solution as the ratio of viscous to capillary forces was successively increased in a sequence of calculations. This contrasted favorably with solutions by other numerical procedures which would require attendant grid refinements to approach the zero capillary pressure results. INTRODUCTION Numerical solution of displacement problems has until recently relied on applying methods developed primarily for transient heat-flow problems. Such problems are classified as parabolic in type, and where the heat transport is purely by diffusion their solutions are characterized by a high degree of smoothness. It is not surprising, therefore, that for approximating these solutions available finite difference methods are quite adequate. In flow problems the transport is partly by diffusion, partly by convection or flow. Although the problem remains of parabolic type because the dispersive effects of capillary forces or diffusion play some role in every displacement, at high flow rates the problem is dominated by convection, and solutions tend toward those of equations of the hyperbolic type. Solutions of hyperbolic problems are characterized by the translation of fronts, or discontinuities, that may progressively increase in sharpness. Numerical methods for treating parabolic problems become less and less satisfactory as displacement rates increase and the role of dispersion due to concentration or capillary pressure gradients becomes small relative to transport due to flow. In computation the difficulty manifests itself as an error associated with the grid size chosen. 1-6 In summary, if the heat-flow type approximations are to include the terms arising due to convection, one of several choices may be made: (1) an upstream (to the direction of flow) approximation for the convection terms may be used; (2) a centered-in-distance (CID) approximation may be used; or (3) a recently developed approximation based on the theory of oscillation matrices may be chosen.6 The last appears to have significant promise for one-dimensional flow problems; its extendibility to two or three dimensions is an open question. In either of the first two approaches, a suitably small ratio of v&/D must be maintained, where v is the velocity, & is the grid spacing and D the effective dispersivity in the direction of flow. In the first choice, the approximation of the convective part is only first-order correct and errors introduced appear as a numerically induced dispersivity of magnitude proportional to v?x. In the CID choice, the approximation can be second-order correct, but the difference formulation fails to satisfy the maximum principle unless a condition on v?x/D is met. Practically, this means that for high flow rates oscillatory solutions may result in the neighborhood of a front unless exceedingly small grid intervals are taken. While the procedure proposed by Stone and Brian4 permits a less severe limitation to be placed on this ratio, ultimately the flow rates increase relative to the dispersivity the oscillation obtains. Further, extensions of their approach to higher dimensional systems may be attended by considerable

Jan 1, 1967

By H. S. Price, J. C. Cavendish, R. S. Varga

A numerical formulation of bigh-order accuracy, based on variational methods, is proposed for the solution of multidimensional diffusion-convection-type equations. Accurate solutions are obtained without the difficulties that standard finite difference approximations present. In addition, tests show that accurate solutions of a one-dimensional problem can be obtained in the neighborhood of a sharp front without the need for a large number of calculations for the entire region of interest. Results using these variational methods are compared with several standard finite difference approximations and with a technique based on the method of characteristics. The variational methods are shown to yield higher accuracies in less computer time. Finally, it is indicated how one can use these attractive features of the variational methods for solving miscible displacement problems in two dimensions. INTRODUCTION The problem of finding suitable numerical approximations for equations describing the transport of heat or mass by diffusion and convection simultaneously has been of interest for some time. Equations of this type, which will be called diffusionconvection equations, arise in describing many diverse physical processes. Of particular interest here is the equation describing the process by which one miscible liquid displaces another liquid in a one-dimensional porous medium. The behavior of such a system is described by the following parabolic partial differential equation:

Jan 1, 1969

By R. W. Snyder, W. H. Guilinger, B. S. Gottfried

Two numerical methods are presented JOT solving the equations Jor one-dimensional, multiphase /low in porous media. The case oF variable physical properties is included in the Formulation, although gravity and capillarity are ignored. Both methods are analyzed mathematically, resulting in upper and lower bounds for the ratio of time step to mesh spacing. The methods are applied to two- and three-phase waterflooding problems in laboratory-size cores, and resulting saturation and pressure distributions and production histories are presented graphically. Results oJ the tuo-phase flow problem are in agreement with the predictions of the Buckley-Leuerett theory. Several three-phase flow problems are presented which consider variations in the water injection rate and changes in the initial oil-and ujater-saturation distributions. The results are different physically from the two- phase case; however, it is shown that the Buckley-Lez~erett theory can accurately predict fluid interlace velocities and displacing-fluid frontal snturations for three-phase flow, providing the correct assumptions are made. The above solutions are used as a basis JOT et~aluuting the numerical methods with respect to machine time requirements and allowable time step for a fixed mesh spacing. INTRODUCTION Considerable progress has been made in recent years in obtaining numerical solutions of the equations for two-phase flow in porous media. Douglas, Blair and Wagner2and McEwenll present different methods for solving the one-dimensional case for incompressible fluids with capillarity (the former using finite differences, the latter with an approach based upon characteristics). Fayers and Sheldon4 and Hovanesian and Fayers8 have extended these studies to include the effects of gravity. West, Garvin and Sheldon, l4 in a pioneer paper, treat linear and radial systems with both capillarity and gravity and they also include the effects of compressibility. Douglas, Peaceman and Rachford3 consider two-dimensional, two-phase, in compressible flow with gravity and capillarity and Blair and Peaceman1 have extended this method to allow for compressible fluids. No one, however, has examined the case of three-phase flow, even for the relatively simple case of one-dimensional flow of incompressible fluids in the absence of gravity and capillarity. In obtaining a numerical technique for simulating forward in situ combustion laboratory experiments, Gottfried5 has developed a method for solving the one-dimensional, compressible flow equations with any number of flowing phases. Gravity and capillarity are not included in the formulation. The method has been used successfully, however, for two- and three-phase problems in a variable-temperature field with sources and sinks. This paper examines the algorithm of Gottfried more critically. Two numerical methods are presented for solving the one-dimensional, multi-phase flow equations with variable physical properties. Both methods are analyzed mathematically, and are used to simulate two- and three-phase waterflooding problems. The numerical solutions are then taken as a basis for comparing the utility of the methods.

Jan 1, 1967