Reservoir simulation

The Merriam-Webster Dictionary defines simulate as assuming the appearance of without the reality. Simulation of petroleum reservoir performance refers to the construction and operation of a model whose behavior assumes the appearance of actual reservoir behavior. A model itself is either physical (for example, a laboratory sandpack) or mathematical. A mathematical model is a set of equations that, subject to certain assumptions, describes the physical processes active in the reservoir. Although the model itself obviously lacks the reality of the reservoir, the behavior of a valid model simulates—assumes the appearance of—the actual reservoir.

The purpose of simulation is estimation of field performance (e.g., oil recovery) under one or more producing schemes. Whereas the field can be produced only once, at considerable expense, a model can be produced or run many times at low expense over a short period of time. Observation of model results that represent different producing conditions aids selection of an optimal set of producing conditions for the reservoir.

Tools

The tools of reservoir simulation range from the intuition and judgment of the engineer to complex mathematical models requiring use of digital computers. The question is not whether to simulate, but rather which tool or method to use. The Reservoir Simulation chapter in the 1987 edition of the Petroleum Engineering Handbook included a general description of reservoir simulation models, a discussion related to how and why they are used, choice of different types of models for different-reservoir problems, and reliability of simulation results in the face of model assumptions and uncertainty in reservoir-fluid and rock-description parameters. That material is omitted in this page; instead it summarizes current practices and trends related to development and application of reservoir simulation models.

Models have been referred to by type, such as

• Black-oil
• Compositional
• Thermal
• Generalized
• IMPES (implicit, sequential, adaptive implicit)
• Single-porosity or dual-porosity
• and more.

These types provide a confusing basis for discussing models;

• Some refer to the application (e.g., thermal)
• Others to the model formulation (e.g., implicit)
• Yet others to an attribute of the reservoir formation (e.g., dual-porosity)

The historical trend, though irregular, has been and is toward the generalized model (outlined below), which incorporates all the previously mentioned types and more. The generalized model, which represents most models in use and under development today. Current model capabilities, recent developments, and trends will then be discussed in relation to this generalized model.

Comparative solution project

The SPE Comparative Solution Project was a collaborative effort to enable the calibration and comparison of different simulation approaches using a common data set and common scenarios. The 10 SPE Comparative Solution Project problems, SPE1 through SPE10, are used for some examples below. Table 1 gives a brief description of those problems.

• 

  Table 1

Generalized model

Any reservoir simulator consists of n + m equations for each of N active gridblocks comprising the reservoir. These equations represent conservation of mass of each of n components in each gridblock over a timestep Δt from tⁿ to tⁿ⁺¹. The first n (primary) equations simply express conservation of mass for each of n components such as oil, gas, methane, CO₂, and water.

The components listed above are denoted by subscript I = 1,2,…, n. In the thermal case, one of the "components" is energy and its equation expresses conservation of energy. An additional m (secondary or constraint) equations express constraints such as equal fugacities of each component in all phases where it is present, and the volume balance S_w + S_o + S_g + S_solid = 1.0, where S solid represents any immobile phase such as precipitated solid salt or coke.

There must be n + m variables (unknowns) corresponding to these n + m equations. For example, consider the isothermal, three-phase, compositional case with all components present in all three phases. There are m = 2n + 1 constraint equations consisting of the volume balance and the 2n equations expressing equal fugacities of each component in all three phases, for a total of n + m = 3n + 1 equations. There are 3n + 1 unknowns: p, S_w, S_o, S_g, and the 3(n – 1) independent mol fractions x_ij , where i = 1,2,..., n – 1; j = 1,2,3 denotes the three phases oil, gas, and water. For other cases, such as thermal, dual-porosity, and so on, the m constraint equations, the n + m variables, and equal numbers of equations and unknowns can be defined for each gridblock.

Because the m constraint equations for a block involve unknowns only in the given block, they can be used to eliminate the m secondary variables from the block’s n primary or conservation equations. Thus, in each block, only n primary equations in n unknowns need be considered in discussions of model formulation and the linear solver. The n unknowns are denoted by P_i1, P_i2 ,…, P_in, where P_in is chosen as pressure p_i with no loss of generality. These primary variables may be chosen as any n independent variables from the many available variables: phase and overall mol fractions, mol numbers, saturations, p, and so on. Different authors choose different variables. Any sensible choice of variables and ordering of the primary equations gives for each gridblock a set of n equations in n unknowns which is susceptible to normal Gaussian elimination without pivoting. The (Newton-Raphson) convergence rate for the model’s timestep calculation is independent of the variable choice; the model speed (CPU time) is essentially independent of variable choice.

The Ith primary or conservation equation for block i is

....................(1)

where:

• M_iI is mass of component I in gridblock i
• q_ijI is the interblock flow rate of component I from neighbor block j to block i
• q_iI is a well term

With transposition, this equation is represented by f_iI = 0, the Ith equation of gridblock i. All n equations f_iI = 0 for the block can be expressed as the vector equation F_i = 0 where f_iI is the Ith element of the vector F_i. Finally, the vector equation

....................(2)

represents the entire model, where the ith element of the vector F is F_i. F is a function of the N vector unknowns P_i, where the Ith scalar element of P_i is P_iI. Application of the Newton-Raphson method gives

....................(3)

where δP is P^{l +1}–P^l and the N × N matrix A represents the Jacobian ∂F/∂P. The element A_ij of A is itself an n × n matrix ∂F_i/∂P_j with scalar elements a_rs = ∂f_ir /∂P_js, r and s each = 1,2,..., n. Eq. 3 is solved by the model’s linear solver. The matrix A is very sparse because A_ij is 0 unless block j is a neighbor of block i.

The calculations for a timestep consist of a number of Newton (nonlinear or outer) iterations terminated by satisfaction of specified convergence criteria. Each Newton iteration requires:

• Linearization of the constraint equations and conservation Eq. 1.
• Linear algebra to generate the A matrix coefficients.
• Iterative solution of Eq. 3 (inner or linear iterations).
• Use of the new iterate P^l+1 to obtain from Eq. 1 the moles of each component in the gridblock.
• A flash to give phase compositions, densities, and saturations which allow generation of the A matrix coefficients for the next Newton iteration.

Model formulations

A major portion of the model’s total CPU time is often spent in the linear solver solution of Eq. 3. This CPU time in turn reflects the many multiply operations required. The model formulation has a large effect on the nature and expense of those multiplies.

Implicit vs. explicit

The interblock flow term in Eq. 1,

....................(4)

uses phase mobilities, densities, and mol fractions evaluated at the upstream blocks. A gridblock is implicit in, say, the variable S_g if the new time level value S_gⁿ⁺¹ is used to evaluate interblock flow terms dependent upon it. The block is explicit in S_g if the old time level value S_gⁿ is used.

The implicit formulation

The implicit formulation expresses interblock flow terms using implicit (new time level) values of all variables in all gridblocks. As a consequence, all nonzero A_ij elements of the A matrix of Eq. 3 are full n × n matrices. The resulting multiplies in the linear solver are then either matrix-matrix or matrix-vector multiplies, requiring work (number of scalar multiplies) of order n³ or n², respectively.

The IMPES formulation

Early paper presented the basis of the IMPES (implicit pressure, explicit saturations) formulation for the black-oil case: take all variables in the interblock flow terms explicit, except for pressure, and eliminate all nonpressure variables from the linearized expressions for M_iIⁿ⁺¹ in Eq. 1. The obvious extension to any type model with any number of components was presented later, and numerous IMPES-type compositional models have been published.

The model Eq. 3 can be written as:

....................(5)

If all variables but pressure are explicit in the interblock flow terms, then all entries but those in the last column of the n × n A_ij (j ≠ i) matrix are zero (recall, the n th variable in each gridblock, P_in, is pressure p_i). This allows elimination of all nonpressure variables and reduction of the vector Eq. 5 to the scalar equation in pressure only:

....................(6)

or

....................(7)

where A is now a scalar N × N matrix and the P and F vectors have N scalar elements p_i and f_i , respectively. The multiplications required in solution of the IMPES pressure Eq. 7 are scalar multiplications, requiring a small fraction of the work of the matrix-matrix and matrix-vector multiplications of the implicit formulation. Thus, the model CPU time per gridblock per Newton iteration for moderate or large n is much less for the IMPES formulation than for the implicit formulation.

The sequential formulation

The stability of the IMPES formulation for the two-phase water/oil case was improved by following the IMPES pressure equation solution with solution of a water saturation equation using implicit saturations (mobilities). This concept was extended to the three-phase case and called the sequential formulation. For each Newton iteration, this method requires solution of the IMPES pressure Eq. 7, followed by solution for two saturations from a similar equation where the A_ij elements of A are 2 × 2 matrices.

A sequential compositional model was described and mentioned the desirability of a sequential implicit treatment of mol fractions in addition to saturations.

The adaptive implicit formulation

The Adaptive Implicit Method (AIM) uses different levels of implicitness in different blocks. In each gridblock, each of the n variables may be chosen explicit or implicit, independent of the choices in other gridblocks. The choices may change from one timestep to the next. This results in the same equation AδP = – F^l as the Implicit formulation except that the elements A_ij of the A matrix are rectangular matrices of variable size. The numbers of rows and columns in A_ij equal the numbers of implicit variables in blocks i and j, respectively; all A_ii are square matrices. The CPU expense per Newton iteration of an AIM model lies between those of IMPES and Implicit models, tending toward the former as more blocks are taken implicit in pressure only.

Choice of formulation

For a given problem, the previous four formulations generally give widely different CPU times. Generalizations regarding the best formulation have many exceptions. Arguably, the trend is or should be toward sole use of the AIM formulation. Current simulation studies use all of these formulations. The Implicit formulation is generally faster than IMPES for single-well coning studies, and for thermal and naturally fractured reservoir problems. For other problems, IMPES is generally faster than Implicit for moderate or large n (say, n > 4). Most participants used IMPES for SPE Comparative Solution Project problems:

• SPE1
• SPE3
• SPE5
• SPE10

All participants used the Implicit formulation for:

• SPE2
• SPE4
• SPE6
• SPE9

No participants in SPE1 through SPE10 used a Sequential model, and, with few exceptions, none used AIM.

A frequently stated generalization is that numerical dispersion error is significantly larger for Implicit than for IMPES formulations. Truncation error analysis shows this error to be proportional to Δx + uΔt for Implicit and Δx – uΔt for IMPES. Real problem nonlinearities and heterogeneity render the analysis approximate and the generalization of limited merit. For example, Figs. 1 through 3 show virtually identical Implicit and IMPES results for the black-oil 9,000-block SPE9 and 2,000-block gas/oil SPE10 problems. For SPE9 (SPE10), the average timestep was 67 (9.7) times larger for Implicit than for IMPES. The percentage of total CPU time spent in the linear solver for IMPES (Implicit) was 23.7 (57.3) for SPE9 and 35.4 (73.4) for SPE10.

• 

  Fig. 1 – IMPES vs. implicit oil rate and cumulative oil production, SPE9

• 

  Fig. 2 – IMPES vs. implicit GOR and water cut, SPE9.

• 

  Fig. 3 – IMPES vs. implicit oil rate, cumulative oil production, and GOR, SPE10 Model 1.

Implementation is an important factor in the relative efficiencies of different formulations. For a given problem, different models using the same formulation can give widely different CPU times. For example, the IMPES CPU times reported by different participants in the six-component compositional SPE5 problem differed by a factor of over 50.

Advances in model formulations

The IMPES formulation was improved by concepts of relaxed volume, better choice of variables, 13 and "adaptive" flash calculations.

Relaxed volume

The relaxed volume concept relates to the timestep calculation Steps (d) and (e) given previously.

• Step (d) gives the mass of each component in the gridblock, , which in turn gives overall composition {z_I}^l+1.
• The Step (e) flash then gives phase amounts and densities which in turn give new iterate S_w, S_o, and S_g values.

These saturations do not sum to 1.0 because of the nonlinear nature of the conservation Eq. 1. If the saturations are altered (e.g., divide each by their sum) to exactly satisfy the volume balance Σ_J S_J = 1, then an incremental (timestep) mass-balance error occurs. If the saturations are not altered, then mass is conserved but there is a volume-balance error Σ_J S_J – 1. The authors chose to preserve mass and carry forward the volume balance error from iterate to iterate and step to step. The volume balance going into iterate . This in effect conserves both mass and volume because there is no permanent or accumulating volume error—only that of the given timestep. Equally important, there is no need to iterate out the volume error to a "tight" tolerance, and Newton iterations and model CPU are reduced. In contrast, the previous or historical IMPES procedure reset saturations to preserve volume and iterated out the mass-balance error. Because the latter error was not carried forward, more Newton iteration (and CPU time) was required to keep the permanent, accumulating mass balance error tolerably low. This use of relaxed volume with carryover also reduces Newton iterations and CPU time in the Implicit formulation.

This discussion implies some fundamental advantage of preserving mass and iterating out volume error as opposed to preserving volume and iterating out mass error. In the writer’s opinion, that is not true provided the error is carried forward in both cases. The Newton iteration requirement and CPU time should be similar if "equivalent" mass and volume error tolerances are used as convergence criteria.

Variable choice

The linear algebra required to reduce the gridblock’s n conservation equations to the IMPES pressure equation is influenced by the choice of variables. The influence is absent for black oil, moderate for "moderate" n and up to a factor of three for large n (say, > 15). The choices of p and mol fractions {z_I} or mol numbers are better than the choice of p, saturations, and phase mol fractions for large n. The effect of this variable choice on total CPU time is often small because the affected work is often a small part of total CPU time. This IMPES reduction is absent in the Implicit formulation and the last of the above variable choices is arguably preferable.

Adaptive flash calculations

The work of equations of state (EOS) flash calculations, including the generation of fugacities and their derivatives, can significantly affect model efficiency when the linear solver does not dominate total CPU time. There may be little need to perform (most of) that work in a gridblock when p and composition are changing slowly. Use of internal, intelligent criteria dictating when that work is needed can significantly reduce the total-run flash calculation CPU time. This is similar in principle to the AIM selection of explicit variables for gridblocks which are quiescent in respect to throughput ratio.

Stable timestep and switching criteria

This topic relates to the observation that lower run turnaround time can increase benefits from a reservoir study allotted a budgeted time period. As a corollary, time spent in repeated runs fighting model instabilities or time-stepping is counterproductive. While many factors affect this run time, it always equals the product (CPU time/step) × (number of timesteps). The first factor is "large" and the second "small" for the Implicit formulation, and conversely for the IMPES formulation. IMPES is a conditionally stable formulation requiring that Δt < Δt* to prevent oscillations and error growth, where Δt* is maximum stable timestep. The conditional stability stems from the explicit treatment of nonpressure variables in the interblock flow terms. Mathematicians performed stability analyses for constant-coefficient difference equations bearing some resemblance to IMPES. Authors in our industry extended and applied their results to derive expressions for Δt*, in particular,

....................(8)

for the black-oil 3D case of gas/oil flow. This shows that stable step Δt* is dependent upon the following, which of course vary with time and from one gridblock to another:

• Flow rates
• Phase mobility
• Capillary pressure derivatives

Thus, at a given timestep, there are block-dependent stable step values Δt*_i , where

• 1 < i < N
• the IMPES stable step is Min(i) Δt*_i

An IMPES model using this internally determined stable step will run stably but may suffer from the weakest-link principle. As an extreme example, consider a 500,000-gridblock problem where, over a 100-day period, the Δt*_i value is:

• 0.01 day for one block
• > 30 days for all other blocks

The IMPES model will require 10,000 timesteps over the 100-day period.

In the AIM formulation, the stable step Δt*_i depends upon the number and identities of variables chosen explicit in block i; theoretically, Δt*_i = ∞ if all block i variables are chosen implicit. In the previous example, all nonpressure variables could be chosen implicit in the block where Δt*_i = 0.01 and explicit in all other blocks. The AIM model would then require CPU time/step essentially no greater than the IMPES model but would require only three timesteps for the 100-day period.

Numerous papers address the problem of determining expressions for the Δt*_i for use internally as switching criteria to select block variables as explicit or implicit in the AIM model. The stability analyses involved are complex and may be impractically complex when allowing the implicit vs. explicit variable choice to include all permutations (in number and identity) of the n variables. The most reliable and efficient AIM models in the future will stem from continuing research leading to the following: (a) Δt*_i estimates which are "accurate," and (b) implicit vs. explicit variable choices, block by block, which are near-optimal and minimize total CPU time, (CPU time/step) × (number of steps).

The linear solver

Preconditioned Orthomin is the most widely used method for iterative solution of Eqs. 3 or7. Nested Factorization (NF) and incomplete LU factorization [ILU(n)] are the two most widely used preconditioners. (See Reservoir simulation linear equation solver for more information.)

• The term "LU factorization" refers to the factoring of the matrix A into the product of a lower triangular matrix L and an upper triangular matrix U. That is an expensive operation but is straightforward, involving only Gaussian elimination.
• The term "ILU(n)" denotes incomplete LU factorization, where only limited fill-in is allowed and n is the "order of fill."

NF performs exceptionally well when transmissibilities associated with a particular direction (in a structured grid) dominate those in other directions uniformly throughout the grid. In general, ILU(n) or red-black ILU(n) [RBILU(n)] is less sensitive than NF to ordering of the blocks and spatial variation of the direction of dominant transmissibilities. In addition, RBILU(n) or ILU(n) have the parameter n (order of allowed infill) which can be increased as needed to solve problems of any difficulty.

A literature search and discussions with numerous developers and users have failed to establish consensus on whether NF or ILU preconditioning is better. Some are strong advocates of one method and others are just as adamantly supportive of the other. But many find, like this writer, that the better method is problem-dependent and it is difficult to find a reliable a priori indicator for making an up-front choice. In the writer’s experience, (a) when NF works well, it is faster than ILU methods, (b) RBILU(0) with no residual constraint is frequently the best of the ILU variants and a good default choice, and (c) in some cases, global residual constraint with the ILU or RBILU method is beneficial.

Cartesian grids and reservoir definition

For many years, simulation used orthogonal Cartesian grids. In the past 15 years, numerous papers have described local grid refinement and various non-Cartesian grids. These papers show that non-Cartesian grids can reduce grid-orientation effects and provide definition and accuracy near wells, faults, highly heterogeneous areas, and so on more efficiently than Cartesian grids.

The premise that Cartesian grids cannot provide required accuracy efficiently in these respects has come to be accepted as a fact. In addition, advances in geophysics have led to geostatistical description of permeability and porosity on a fine scale once unimaginable. Increasingly, our papers include examples using thousands of gridblocks for two- or few-well "patterns," in part to reflect these geostatistical descriptions. The purpose of this section is to show, using a few examples, that Cartesian grids can provide adequate accuracy and reservoir and near-well definition efficiently in some cases, even without local grid refinement. No generalizations from the examples used are intended. For the most part, the examples are taken from the literature.

SPE7 is an x-direction horizontal well problem with a 9 × 9 × 6 Cartesian grid representing a 2,700 × 2,700 × 160-ft reservoir section. The specified block Δy values decrease from 620 to 60 ft at the well, presumably to increase near-well definition and accuracy of results. The Δx are uniformly 300 ft. Fig. 4 compares Case 1A results for the SPE7 grid with results using uniform areal spacing Δx = Δy = 300 ft. The near-well y-direction refinement of the specified grid has no effect and is not necessary in this problem.

• 

  Fig. 4 – Effect of near-well grid refinement, SPE7 Case 1A.

SPE8 is a gas/oil problem with one gas injector and two producers on the corners of a 5,000 × 5,000 × 325-ft square reservoir. A 10 × 10 × 4 Cartesian grid with uniform Δx = Δy = 500 ft is specified. Five participants compared their results for that grid with results from their (areally) locally refined or unstructured grids. They showed good agreement for grids having approximately four times fewer blocks than the 10 × 10 × 4 grid. Fig. 5 shows equally good agreement for a 5 × 5 × 4 (Δx = Δy = 1,000 ft) Cartesian grid with no local refinement.

• 

  Fig. 5 – Effects of coarser grids on producing GOR.

SPE10 (Model 1) is a 2D cross-sectional gas/oil problem with a geostatistical permeability distribution given on a 100 × 1 × 20 Cartesian fine grid. Coarse-grid submittals included results using upscaling and local grid refinement. A homogeneous 5 × 1 × 5 Cartesian grid with no alteration of relative permeability matched the 100 × 1 × 20 results nearly exactly.

SPE10 (Model 2) is a 3D water/oil problem with a 1.122 million-cell geostatistical grid. Some coarse grid submittals included sophisticated upscaling and gridding techniques with no pseudoization of relative permeability and grids from 4,810 to 70,224 blocks. Others used simple flow-based upscaling to 75- to 2,000-block Cartesian grids with moderate k_r changes. In general, the latter submittals showed the best agreement with the fine-grid solution.

Numerous papers show that non-Cartesian grids can significantly reduce the grid-orientation effects of Cartesian grids. However, most of the examples used to study those effects are highly adverse mobility ratio displacements in homogeneous, horizontal reservoirs. In reservoirs with more normal fluid mobilities, areal fluid movement is more strongly affected by heterogeneity and/or gravity forces associated with reservoir structure (variable dip), and grid-orientation effects tend toward a second-order effect. As an example, the SPE10 (Model 2) water/oil problem reservoir is highly heterogeneous. Fig. 6 compares five-point and nine-point field results for an upscaled 28 × 55 × 85 Cartesian grid. The close agreement indicates an absence of grid-orientation effects even though the unfavorable oil/water viscosity ratio is 10 and there is no dip.

• 

  Fig. 6 – 9-point vs. 5-point, SPE10 Model 2.

Example 1

Table 1 gives data for Example 1, a ¼ five-spot, vertical-well problem. Fig. 7 shows two block-centered grids (a) and (b) used for this type of problem. The four-fold smaller well blocks of grid (b) provide finer well definition and presumably increase the accuracy of results. Fig. 8 shows the identical results for 10 × 10 grid (a) and 11 × 11 grid (b). Results are nearly identical for the 5 × 5 grid (a) and 6 × 6 grid (b), and Fig 9 shows insignificant difference between 3 × 3 grid (a) and 4 × 4 grid (b) results. The grid (b) doubles the grid (a) IMPES run CPU time but contributes no greater accuracy. Well-index effects are not important here. When they are, a single one-layer single-phase run can be made to determine the index correction factor for grid (a) wells located in the corners of their gridblocks.

• 

  Table 1 - Example 1 data.

• 

  Fig. 7 – Block-centered 3 x 3 and 4 x 4 Cartesian grids used in Example 1

• 

  Fig. 8 – Effect of (a) 10 x 10 vs. (b) 11 x 11 grid on Example 1 results.

• 

  Fig. 9 – Effect of (a) 3 x 3 vs. (b) 4 x 4 grid on Example 1 results.

Fig. 10 shows a small effect of grid refinement on Example 1 results for grids from 20 × 20 to 3 × 3. The results indicate little need to enhance near-well definition by unstructured grids or by grid refinement (global or local) for grids finer than 3 × 3 for this problem.

• 

  Fig. 10 – Effects of Cartesian grid coarsening, Example 17.1 results.

Example 2

Flexible non-Cartesian grids are shown to significantly reduce the required number of gridblocks. An example provided was water/oil coning in a horizontal well in a 600 × 300 × 230 m homogeneous reservoir. Results were: (a) a 25,823-block 31 × 17 × 49 Cartesian grid was required to obtain a converged solution, and (b) a 3D 2,066-block control-volume finite-element (CVFE) unstructured grid gave correct results. Table 2 gives data for Example 2, a similar problem. Fig 11 compares Example 2 results for 60 × 31 × 48 and 10 × 7 × 9 Cartesian grids. The 630-block coarse Cartesian grid results here agree as well with the Cartesian 60 × 31 × 48 fine-grid results as the reported 2,066-block CVFE results agree with the 31 × 17 × 49 Cartesian fine-grid results.

• 

  Table 2 - Example 2 data.

• 

  Fig. 11 – Fine vs. coarse Cartesian grids, Example 2

Non-Cartesian grids are argued to define irregular reservoir boundaries more efficiently than Cartesian grids. This is not necessarily true. For over 30 years, many models have used active-block coding. While the Cartesian grid extends past boundaries to numerous inactive blocks, those inactive blocks are dropped by the model and require no computer storage or CPU. These numerous inactive blocks pose a problem only for models, if any, that do not use active-block coding.

Nomenclature

fg′	=	∂fg / ∂Sg
qIji	=	interblock flow rate of component I from block j to block i, volume/time
M	=	mass, mobility, molecular weight, or preconditioner matrix depending on use
N	=	number of active gridblocks
Pcgo′	=	dPcgo / dSg
Vp	=	pore volume, bbl
Δp	=	pressure difference, psi
Δt	=	timestep
ΔZ	=	Zj – Zi
γ	=	phase density expressed as psi/ft or bar/m
ρ	=	density, mass/volume or molar density, mols/volume unless noted otherwise
λ	=	phase mobility, kr / μ
ψ	=	λoλg / (λo + λg)
∂	=	normalized difference parameter

Subscripts

I	=	component number, index counter, or initial condition
g	=	gas phase
j	=	index counter
o	=	oil phase
p	=	phase designation
x	=	cell designation
x,y,z	=	Cartesian x,y,z directions

Noteworthy papers in OnePetro

Use this section to list papers in OnePetro that a reader who wants to learn more should definitely read

Online multimedia

Jafarpour, Behnam, Ryosuke Okuno, and Rami M. Younis. 2013. Advances On Sparse Representation of Reservoir Properties, Fluid Characterization Methods for Steam-Solvent Co-Injection, and Full-Resolution Simulations of Multi-Scale Processes. https://webevents.spe.org/products/advances-on-sparse-representation-of-reservoir-properties-fluid-characterization-methods-for-steam-solvent-co-injection-and-full-resolution-simulations-of-multi-scale-processes

Lake, Larry W. 2013. Modeling and Simplicity: Occam’s Razor in the 21st Century. https://webevents.spe.org/products/modeling-and-simplicity-occams-razor-in-the-21st-century

External links

Use this section to provide links to relevant material on websites other than PetroWiki and OnePetro

See also

Reservoir simulation linear equation solver

Gridding in reservoir simulation

Upscaling of grid properties in reservoir simulation

Phase behavior in reservoir simulation

Geomechanics in reservoir simulation

Streamline simulation

Reservoir simulation applications

PEH:Reservoir_Simulation