In this paper we develop a new mathematical formulation for chemical compositional reservoir simulation, and provide a comparison of its results on alkaline-surfactant-polymer ooding with those of UTCHEM simulator. Our research has found that the existing chemical compositional model estimates the adsorption effect on the transport of a component reasonably well but it does not satisfy the principle of mass conservation. Since the total mass conservation equation follows from summing the species-conservation equations over all components, the obtained equation violates the principle of total mass conservation as well. With these partial differential equations as governing equations, several simulators have been developed. In this work, we propose an approach to model the change in pore volume due to adsorption that satis es the mass conservation law, and allows applying a sequential solution approach.
Chemical ooding is one of the most promising and broadly applied enhanced oil recovery (EOR) processes. Chemical ooding can be further subdivided into alkaline ooding, surfactant ooding, polymer ooding, and alkaline-surfactant-polymer (ASP) ooding. Alkali reduces adsorption of the surfactant on the rock surfaces and reacts with acids in the oil to create natural surfactant. Surfactants are chemicals that used to reduce the interfacial tension between the involved uids, increasing oil mobility. ASP ooding is a form of chemical EOR method that can allow operators to extend reservoir pool life and extract incremental reserves currently inaccessible by conventional methods. While ASP ooding has a high efficiency it is technical, costly, and risky. Model studies can assist in this evaluation.
Most multiphase compositional models reported in the literature [1], [2], [3], [4] and [5] are limited in their applicability in one way or another (single species, equilibrium mass transfer, and lack of miscibility modeling etc). The mathematical formulation ⋆ Corresponding author Copyright ⃝c by the paper's authors. Copying permitted for private and academic purposes. In: A. Kononov et al. (eds.): DOOR 2016, Vladivostok, Russia, published at http://ceur-ws.org developed in this work is extended from the UTCHEM model formulation for use in chemical ooding studies that does not have these common limitations.
Sequential schemes are very suitable for chemical compositional ow problems which include a large number of chemical components. Only the implicit pressure and explicit composition (IMPEC) formulation was used for chemical compositional reservoir simulation so far, but there is no obvious reason why the sequential formulation can't be used as well. Because of explicit solution of compositions, the size of time steps is limited to stabilize the general procedure.
Chen et al. [6] presented a numerical approach that solves both pressure and compositions implicitly. Though the approach was claimed to be sequential and extended from the IMPEC approach used in UTCHEM model [7], the mathematical formulations for the governing equations did not undergo any change in their model.
The basic equations used in UTCHEM model that describe multiphase, multicomponent ow in permeable media are the species-conservation, pressure (an overall mass-continuity), and energy conservation equations. Accumulation terms in the species-conservation equations used in UTCHEM model account for the reduction in pore volume caused by adsorption.
During the process of this research, it was revealed that this commonly used approach estimates the adsorption effect on the transport of a component reasonably well but it does not satisfy the species-conservation equation. Since the total mass conservation equation follows from summing the species-conservation equations over all components, the obtained equation violates the principle of total mass conservation as well. In recent years with use of these governing equations several simulators were developed for simulation of the chemical ooding processes [6], [7], [8] and [9].
In this work we introduce a new approach to model the reduction in pore volume due to adsorption that satis es the conservation equations. In certain situations, such as signi cant change in the effective pore size due to adsorption, these enhancements are essential to properly model the physical phenomena occurring in petroleum reservoirs. In addition, this new approach for modeling the adsorption effect on the transport of a component makes it possible to develop a new mathematical formulation for the sequential chemical compositional reservoir simulation. 2
Consider a bulk volume Vb at some point within a porous medium domain. Let us assume that this representative elementary volume (REV) is made up of np + 1 phases (np uid phases and a solid phase consisting of rock grains or soil) with nc chemical species. Conceivably, at least, each species can exist in any phase and can transfer between phases via evaporation, condensation, dissolution, adsorption and so forth.
In our model formalism, each pair (i; ), with i chosen from the species indices and chosen from the phases, is a constituent. Each constituent (i; ) has its own intrinsic mass density i , measured as mass of i per unit volume of , and its own average velocity !ui . Each phase has its own volume fraction ϕ . The volume fraction of phase , ϕ , is the volume of phase divided by the bulk volume Vb.
If the index nc represents the species and the index np + 1 represents the phases
making up the solid phase, then in terms of the above de ned mechanical variables the
mass balance for each constituent (i; ) is
From left to right in Eq. (
The source term Ri accounts for the rate of mass generation (Ri > 0) and consumption (Ri < 0) of component i in phase , either through chemical or biological reactions. There is no general function for Ri . An example of a rst-order reaction rate for radioactive decay or biodegradation is
Ri = kiϕ !i ; {i = 1; :::; nc; = 1; :::; np + 1;} where ki is the decay constant or reaction rate coefficient in units of inverse time.
The second source term rmi expresses the rate of mass transfer of component i from or into the phase owing to vaporization or condensation. Adsorption is described through isotherms.
The last source term in Eq. (
The phase saturation S is de ned as the fraction of the pore volume occupied by
phase , that is, volume of phase V divided by the pore volume Vp. The saturation
of uid phase can also be de ned as S = ϕ =ϕ. For uid phases such as liquids and
vapors, ϕ = ϕS , = 1; :::; np, where ϕS also called the uid content. For the solid
(s) phase ϕs = 1 ϕ, which is the grain volume divided by the bulk volume Vb. We can
rewrite equation (
(
!
+ u′ i ;
fi = 1; :::; nc;
= 1; :::; np:g
Since phase velocities are typically more accessible to measurement than species
velocities, it is convenient to rewrite the constituent mass balance equation (
So far, the mathematical formulation of the mass conservation equations developed above is essentially the same as the standard formulation described in [7]; where it differs is in the treatment of average velocity in the governing equations. Here we start to deviate from the standard formulation.
The uxes of component i in phase with respect to volume-averaged velocity !JDi = ϕS Ki ∇( !i ) and mass-averaged velocity !JDi = ϕS Ki ∇(!i ) owing to hydrodynamic dispersion alone were presented in [10]. The ux with respect to bulk volume-averaged velocity is proposed in this work: ! J Di = (Kxx)i = (Kxy)i = where the subscript l refers to the spatial coordinate in the direction parallel, or longitudinal, to bulk ow, and t is any direction perpendicular, or transverse, to l. Di is the effective binary diffusion coefficient of component i in phase [12], l and t are the longitudinal and transverse dispersivities, and is the permeable medium tortuosity.
A general set of partial differential equations (
Ki
∇(ϕS
!i ) =
kiϕS
!i + rmi + qi ;
fi = 1; :::; nc;
= 1; :::; np; g
(
The term rmi is difficult to calculate without detailed analysis of the transport
occurring within the phases. One typically simpli es the equations by using overall
compositional balance equations. Overall compositional balance equations can be obtained
by summing Eqs. (
In equation (
= 1; :::; np} (1 ϕ) s!is = ϕ ic^i;
fi = 1; :::; nc; g ∇(ϕS i )] = ki i + (1 ϕ) is fi = 1; :::; nc; g =1 ∇ np ∑[Ki =1 np ∑(ϕS =1 [ np ∑ ϕS =1 nc ( np ∑ ∑ S ci + c^i
) i=1 =1
= 1; nc np ∑ ∑ S ci ̸= 1; i=1 =1 unlike existing model. To account for the reduction in pore volume caused by adsorption, the coefficient (1 ∑in=cv1 c^i) is introduced into the overall compositional balance ( np ∑ S ci + c^i )]
+ ∇ ∇ ∑np [Ki ∇(ϕ iS ci )] =
np ϕ i ∑(S ci !u )
=1
kiϕ i
where i is the component mass density in units of mass of i per unit volume of i, ci
is the component concentration in units of volume of i in phase per unit volume of ,
and c^i is the adsorbed component concentration, measured as volume of i in phase
per unit pore volume. The linear, Freundlich and Langmuir adsorption isotherm models
are applied to calculate the adsorbed concentrations c^i. In our de nition
we get
but
(
In the present work we introduce a new approach to model the reduction in pore volume due to adsorption that satis es the continuity equation. Let us denote the modi ed volume fraction of phase due to adsorption by ϕ^ . Porosity ϕ^ is de ned as the fraction of the bulk permeable medium that is pore space remaining after adsorption. This porosity is related to the original porosity ϕ as follows: The saturation of uid phase is de ned as S procedure as carried out above, we obtain = ϕ^ =ϕ^. Now using the same derivation ϕ^ i where k is the permeability tensor, kr is the relative permeability of uid phase , is the dynamic viscosity of uid phase , p is the pressure in uid phase , is the speci c weight for uid phase , and z represents depth.
Variation of pore volume with pore pressure p can be taken into account by the pressure dependence of porosity. The porosity depends on pressure due to rock compressibility, which is often assumed to be constant and can be de ned as ϕ = ϕR[1 + cf (p1 ps)]; where ϕR is the porosity at a speci c pressure ps, p1 is the water phase pressure, and cf is the pore compressibility at ps.
A slightly compressible uid has a small but constant compressibility. For a slightly compressible uid, the component density i can be written as: i = iR[1 + ci0(p1 pR)]; fi = 1; :::; nc; g (17) (18) (19) (20) where iR is the density of component i at the standard pressure pR, a constant value, ci0 is the compressibility of component i.
Now since reference density iR is constant for each component we can divide through both sides of Eq. (18) by iR. In terms of the dimensionless density i = i= iR Eq. (18) can be written as: ϕ^ i ∇ ( np ∑ S ci + c^i =1 np ∑[
Ki =1 )]
np + ∇ ϕ^ i ∑(S ci !u )
=1 ] ∇(ϕ^ iS ci ) = {
np nc } ϕ^ ∑(S !u ∑ ici ) = ∑nc Qi ; =1 i=1 (net dispersive ux in a phase is zero), and according to the total reaction de nition np c~i = ∑ S ci + c^i; =1 ncv ∑ c^i i=1 ) nc ∑ ci0c~i
] i=1 @ ( ncv )
∑ c^i : @t i=1
} pR)] ∑nc (ci0c~i) ; i=1 The total compressibility, ct, is ct = ( 1 ncv ){ ∑ c^i i=1
Ft(c~i; c^i) = (p1 and and by de nition We de ne the overall concentration c~i as
The pressure equation is developed by substituting Darcy's law (Eq. (19)) for the phase ux term of Eq. (23), using the de nition of capillary pressure pc 1 = p p1; = 2; :::; np. The pressure equation in terms of the reference phase (phase 1) pressure is ∇ (
k rT c∇p1) = ∇ ∇
np (k ∑
=1 r c = r r c ∇z
) nc ∑ ici ; i=1
np (k ∑ =2 r c∇pc 1
) ϕFt(c~i; c^i) + ∑nc Qi
; ∇(ϕ^ iS ci ) = K~i ∑np ∇(ϕ^ iS ci ); =1 where !u~ i and K~i can be de ned as some averages. Since differentiation and summation are interchangeable operations in this system, the sum of the gradients can be calculated as the gradient of the sum ~ Ki ∑np ∇(ϕ^ iS ci ) =1 = K~i ∇ (
np ϕ^ i ∑ S ci )
; =1 and
Equation (22) can be written using Eqs. (33), (34), and (35) as below: =1 ( np ∑ S ci + c^i )] + ∇ (
np ϕ^ i !u~ i ∑ S ci
) =1 ∇ [
K~i ∇(ϕ^ i ∑np S ci )] =
This new mathematical formulation of species conservation equations makes it possible to apply a sequential solution approach to solve these equations implicitly for the total concentration ∑np=1 S ci of each component. A ash calculation is then performed to obtain the phase saturations and the concentrations of components in each phase. Numerical values of !u~ i and K~i can be most simply calculated as the weighted averages ! u~ i = np ∑S ci !u =1 np ∑S ci =1 ;
K~i = np ∑Ki =1 ∇(ϕ^ i ∇(ϕ^ iS ci ) ∑np S ci ) =1 ; ASP ooding is the most promising EOR solution for one of the greatest challenges facing the oil industry worldwide: after conventional water ooding the residual oil (drops trapped by capillary forces) in reservoirs around the world is likely to be around 70% of the original oil in place. The mathematical formulation is evaluated in the modeling of a eld scale ASP EOR process.
As illustrated in Fig. 1, the ASP ooding pilot has 4 injection wells and 9 production wells in an inverted ve-spot well pattern. The ASP process was conducted in a 4-slug sequence: pre- ush polymer ood, alkaline/surfactant slug, alkaline/surfactant/polymer slug, and a polymer drive. Total simulation time is 551 days. Reservoir properties include heterogeneous permeability and initial water saturation elds. The reservoir is at a depth of 4150 ft., has an average initial pressure of 1770 psi, and the porosity is assumed to be constant throughout the reservoir and equal to 0.3. Grid dimensions are 19 19 3. The OOIP is 395,427 bbls, the crude oil viscosity is 40 cp, the initial brine salinity is 0.0583 meq/ml and the initial brine divalent cation concentration is 0.0025 meq/ml.
We use S3GRAF software, developed and licensed by Sciencesoft Ltd., for postprocessing the output data.
Three owing phases and eleven components are considered in the numerical simulations. The phases are water, oil and microemulsion, while the components are water, oil, surfactant, polymer, chloride anions, divalent cations (Ca++, Mg++), carbonate, sodium, hydrogen ion, and oil acid. The ASP interactions are modeled using the reactions: in situ generated surfactant, precipitation and dissolution of minerals, cation exchange with clay and micelle, and chemical adsorption. Note the detailed chemical reaction modeling, and the heterogeneous and multiphase petroleum reservoir under consideration.
A comparison with UTCHEM has also been performed. The matches between old and new formulations' numerical results for the matched variables are shown in Figs. 25 for the injected pore volume in the range 0 - 1.0 PV. Comparative studies show that the results obtained from IMPEC implementation of the newly proposed formulation are in a good agreement with that of UTCHEM simulator. In the scope of this research work, through its application to the above-mentioned numerical experiment and comparisons with UTCHEM model results, the newly developed formulation has proven to be reliable, practical, and accurate. The mathematical model and numerical simulation developed in this work can also be used to study the transport of contaminants and remediation of contaminated aquifers surfactants. 4
In the scope of this research work, a new mathematical model formulation for multicomponent, multiphase ow in porous media has been developed. During the process of this research, it was revealed that commonly used approach estimates the adsorption effect on the transport of a component reasonably well but it does not satisfy the mass conservation or continuity equation. In the present work we introduce a new approach to model the reduction in pore volume due to adsorption that satis es the continuity equation. The mathematical formulation developed in the scope of this work is extended from the UTCHEM model formulation for use in chemical ooding studies. A comparison with UTCHEM has also been performed. Comparative studies show that the results obtained from IMPEC implementation of the newly proposed formulation are in a good agreement with that of UTCHEM simulator. The implementation of a sequential solution approach for chemical compositional reservoir simulation based on the formulation described in this paper is scheduled for the future.
Acknowledgments. This paper was supported by the Ministry of Science and Education of the Republic of Kazakhstan under grants No. 1735/GF4 and No. 0128/GF4.