Generalized law of multiphase filtration and new effects of surface phenomena at two-phase flows in a porous medium

The generalized law of multiphase filtration, which explicitly takes into account the viscous interaction between fluids, capillary and gravitational effects, is analyzed. A special case is considered – a two-phase filtration flow as the basis of the technology of oil reservoirs flooding. A method for determining new material functions – «coupling» phase permeabilities, and their quantitative contribution to the field development indicators is given. A new effect of the dependence of capillary forces on the direction of flow in anisotropic media has been revealed.


Introduction
Designing the development of hydrocarbon fields relies largely on the theory of multiphase filtration.The latter includes the basic conservation laws and the Darcy Law fairness postulate for the w i α filtration rate of each phase as an equation of motion (Basniyev et al., 1993): where α = 1.2 (in the two-phase option); k ij is the absolute permeability tensor; p α -pressure in phases; latin indices are summed up.
Relation (1) contain two fundamental hypotheses: for a given rock and a pair of filtered fluids, capillary pressure p c and relative phase permeability (RPD) k α are universal, equilibrium, stationary, scalar functions of local saturation s of one of the phases.
Studies in recent decades show that in many cases these hypotheses are not fulfilled.It is often necessary to take into account the effects of nonequilibrium and nonstationarity in the measurement of material functions.In anisotropic media, such as real reservoirs, the RPD and capillary pressure depend not only on saturation, but also on the direction of flow, and have a tensor nature.Consideration of the wettability of the rock, which varies in the process of oil and gas production, is of great importance.The index of wettability is the integral of capillary pressure.At the same time, the derivative of capillary pressure enters the Rappoport-Lis equations, which only makes it possible to control the change in the saturation gradient, for which these equations are asymptotically transformed into the Buckley-Leverett model at too large values.
The motion description of multiphase multicomponent media, accompanied by interconnected physicochemical, various relaxation and other processes, reaches sufficient completeness and logical harmony when using the methods of non-equilibrium thermodynamics.Its main provisions on which the continual modeling of specific systems is based are associated with the names of I. Prigogin, S. de Grot, P. Glensdorf, P. Mazur, N. Dyarmati, R. Defei and other researchers.
The formalism of linear nonequilibrium thermodynamics was successfully applied to the study of complex filtration flows (Nikolaevsky, 1984;Nigmatulin, 1987, etc.).In the listed and other well-known works, when describing multiphase media, they proceed from the balance equations for each phase (or components).At the same time, it is necessary to set heuristically the laws of interfacial and intercomponent interactions.
The specificity of the approach presented in this paper is that, starting from the conservation laws for the entire continuum as a whole and using the methods of nonequilibrium thermodynamics, we can derive the equations of phase motion and obtain the laws of interfacial interaction.
For generalization and development of the theory of multiphase filtration, some results are considered of applying the methods of nonequilibrium thermodynamics to processes in heterogeneous environments.

the generalized equation of motion of multiphase mixtures in a porous medium
In the works (Maksimov, 1994;Kolesnichenko, Maksimov, 1999) using the thermodynamic methods of non-equilibrium processes, a rather general theory of a multi-phase multicomponent chemically active continuum was developed, containing a number of thermodynamic, physicochemical, thermal and other effects.The author of this paper is limited to the analysis of the motion equations for multiphase flows of immiscible liquids, which generalize the Darcy law.For the case of isothermal processes in an isotropic reservoir, these equations have the form: where the index ф corresponds to the formation skeleton, p αβ is the multiphase RPD coefficients, w β and ρ α are the filtration rate and the true density of the α-phase, respectively.Equations ( 2) are used to determine the phase filtration rates and express their dependence on pressure gradients in all phases and gravitational forces in the phase.The last term on the right-hand side of (2) is the Rakhmatulin force, which is associated with the force effect of the system on the selected phase due to the pressure mismatch in the individual phases (capillary effects, strength effects, phase inertia in their small-scale movement).In relation (2), thermal diffusion relations and the contribution of thermophoretic forces, which are a second-order effect, are neglected.
An alternative to equations (2) form of representing the equations of motion in the form of generalized Stefan-Maxwell relations for a heterogeneous medium is shown in (Kolesnichenko, Maksimov, 1999), which includes a symmetric matrix of binary resistance coefficients R αβ .
It is shown (Maksimov et al., 1994) that these coefficients are related to the multiphase RPD coefficients k αβ by the system of algebraic equations: (3) In the particular case of a two-phase flow describing the process of oil displacement by water, the problem of determining the material functions in equations ( 2) and (3) is simplified.Matrices of resistance coefficients R αβ and relative phase permeabilities k αβ are given by three saturation functions.
In this case, the relation (3) takes the form (α = 3 corresponds to the skeleton of the rock): where R αβ -is the binary friction coefficients determined in standard experiments.
The representation of the motion equations in the form permitted with respect to the acting forces has a double advantage.
First, the determination of binary coefficients of resistance R αβ using standard experiments on the joint motion of two phases in a porous medium makes it possible in principle to calculate multiphase RPD k αβ by solving the system of equations ( 3).It should be noted that the methods of direct experimental determination of relative phase permeabilities are absent even for threephase systems.
Secondly, the equations of filtration in this form, being substituted into the equations of continuity, lead to equations that are not resolved with respect to the highest derivative.This creates certain difficulties in solving problems in numerical implementation.

Generalized two-phase filtration model
For an isothermal two-phase flow, an alternative (4) approach is possible, based on a series of well-known experiments (Maximov, 1994).
In this case, the generalized law (2) of multiphase filtration is represented as: Here Rakhmatullin's force is neglected, given the capillary effects of conventional dependence of capillary pressure in the form of Leverett (Basniev et al., 1993).
As can be seen from system (5), the coefficients λ αβ have the meaning of generalized mobility of the phases.In this case, the off-diagonal coefficients of the matrix λ αβ (α ≠ β) are due to the viscous interaction between fluids and capillary forces, and the diagonal coefficients represent the contribution of both phases to the total flow, as if each of them moved independently in the pore space, modified the presence of another phase.
The work (Maksimov, 2015) shows cases in which the classical and generalized two-phase filtration formalisms are equivalent to each other.
Assuming the incompressibility of both phases and the reservoir skeleton, for the one-dimensional displacement of oil by water in an isotropic inclined formation after standard transformations from system (5), taking into account the continuity equation, a determining equation is derived for saturation s of the wetting phase, similar in structure to the Rappoport-Lis equation of the classical model (Basniev et al, 1993): Here, dimensionless variables τ = wt/mL, ξ = x/L are introduced (w is the total phase filtration rate, L is the characteristic length of the reservoir).
For the generalized model: For the Rappoport-Lis equation: In equations ( 7)-( 8), the symbol «ʹ» means differentiation according saturation (for example, etc.); λ α and f(s) are the phase mobility and the Buckley-Leverett function in the classical theory (Basniev et al., 1993).
Comparison of the generalized and classical models leads to the conclusion that the corresponding equations are identical for two independent processes: for countercurrent capillary-gravity impregnation (w = 0) and unidirectional displacement without taking capillary forces into account.This leads to the following relationships: In formulas (9), the symbol «ʹ» corresponds to the impregnation process.From relations (9), the expressions for the coefficients of the matrix λ αβ are found through the standard mobilities λ α ʹ and λ α for the two indicated displacement processes in the form: To determine the standard mobilities λ α ʹ and λ α , two independent experiments are sufficient: unidirectional displacement and counter-current impregnation.
There are a large number of laboratory experiments to determine the RPD for impregnation and drainage processes in the same type of samples.Let us consider one such experiment (Jerauld, Salter, 1990).
As a result, the following dependencies for RPD were obtained.
These results and the dependences of matrix of generalized motions, defined similarly to the previous one, are presented in Fig. 1, 2.
In the works (Maksimov, 2015;2016), a quantitative assessment of the contribution of the surface interaction of fluids to displacement indicators is given.In the calculations, we used the numerical values of parameters and coefficients of the generalized motion matrix corresponding to the experiment (Jerauld, Salter, 1990) (formulas ( 11)-( 12)).
In solving the boundary value problem for the Rappoport-Lis equation, the RPD and capillary pressure functions were used, which correspond to the impregnation cycle (an increase in the saturation of the wetting phase).In this case, the capillary According to calculations for the classical model, watering begins earlier (τ = 241) than for the generalized model (τ = 342).After water breakthrough, the pattern of oil displacement remains different.But asymptotically, with large times and volumes of water injection, the amount of oil displaced, calculated on both models, will be almost the same (Maksimov, 2016).
Thus, a comparison of the generalized filtering law, which takes into account the «cross» phase permeability, with the classical theory showed a different picture of changes in the saturation profiles.The type of displacement in the generalized model is slower; the start time of flooding occurs later, and the displacement rate is lower.P k = -0.276s+ 0.068/s + 0.051.Figures 3 and 4 show the water saturation profiles for different points of time for the classical and generalized models, as well as the dynamics of the displacement coefficient η (Fig. 5).The latter is defined as ratio of the volume of oil displaced by water to the total volume of pores occupied by oil before the start of displacement.After the water breakthrough through the production gallery, the graphs show the dimensionless water breakthrough time τ and the anhydrous oil recovery coefficient η.
The graphs show that in the case of a generalized model, the process of displacement is slower, and the moment of breaking through the water front occurs later.This is due to the consideration of the surface inter-viscous interaction between the phases in the generalized model.At the interface between the media, the more viscous fluid «slows down» the less viscous and, in accordance with this, creates additional resistance in the promotion of the latter (compared to the classical model).The vertical dashed lines correspond to the start of watering for each of the models.

The new effect of capillary forces at the interface in anisotropic media
In the previous sections, the author limited himself to the case of isotropy to reveal the effect of viscous interaction between fluids at the interface.Capillary pressure still remained a scalar function of saturation.Experiments performed on an anisotropic core showed that this is not true.
Initially, a series of complex experiments (Dmitriev et al., 2014) was conducted on selected cylindrical core of cemented layered sandstone 100 mm in height and 100 mm in diameter, which was prepared for research.In accordance with the methodology (Kuznetsov et al., 2010), the type of anisotropy (orthotropic symmetry) and the main directions of the absolute permeability tensor were determined.Next, 4 smaller samples were cut from the initial core, cores with a diameter of 25 mm and a length of 30 mm.Three cores (samples 1-3) were cut along the main directions: samples 2 and 3 -along the X and Y axes in the bedding plane; sample 1 -in the vertical direction Z (Fig. 6).The fourth samplethe control -along the bisector of the angle between the extreme directions in the bedding plane.The last sample was used to test the assumption that the axis of core symmetry coincides with the main direction of the absolute permeability tensor, and that the values obtained as a result of the experiment are indeed components of the tensor.Further, on all oriented samples, the porosity and the main values of absolute permeability were determined for helium filtration: m = 18.6%; k 1 = 668 md, k 2 = 689 md, k 3 = 579 mD.
Using the obtained values of k 2 and k 3 , we can calculate the value of permeability in any direction ( ) , find the theoretical value of k 4 = 638 mD and compare it with the experimental value obtained k 4Э = 644 mD.The difference is less than 2%.
Subsequent studies of oriented core samples on a SkyScan 1172 computer tomograph led to the establishment of the tensor nature of pore distribution density functions over radii and the introduction of the characteristic linear dimensions tensor R ij (effective pore radii).
Further, an additional series of studies was carried out, consisting in determining the residual water saturation for the same oriented samples.To recalculate the laboratory parameters, the Hassler-Bruner technique (Kuznetsov et al., 2010) was used, which allows determining the magnitude of capillary pressure at the outer end of the sample with a high degree of accuracy, estimating its corresponding saturation and building the dependences of capillary pressure on saturation.The research results are shown in Fig. 7 for samples along the X and Y directions in the plane of bedding, and for control sample 4, and in Fig. 8 for sample Z along the vertical direction.The results of a comparison of the theoretical and experimental values of capillary pressure on the control sample suggest the tensor nature of capillary pressure in anisotropic media.
An important consequence of this study is the effect of capillary pressure dependence on the direction of measurement.It follows that capillary pressure is not a universal function of saturation for an anisotropic rock sample, but depends on the direction of flow, on the direction of impact on the reservoir.Confirmation of this fact requires a more in-depth experimental and theoretical study, identifying the physical mechanisms of this effect, the characteristics of the physical and physicochemical interaction of the «skeleton -fluid» system, the nature of the wettability of the formation, and other factors.Repeatability of experiments for various rock samples (terrigenous, carbonate) with a different scale of heterogeneity is required.
Strictly speaking, we should not speak about the tensor nature of capillary pressure, but that the tensor is the surface stress of capillary forces at the interface pijc, whose connection with the tensor p ij c is the inverse tensor of characteristic linear dimensions R ij , can be represented as: ( ) where a c is determined by the interfacial tension coefficient and the wetting angle; n i n j -dyad; J ijkl -rank 4 tensor, symmetric with respect to a pair of indices and their permutation; The tensor r kl is found experimentally (Dmitriev et al., 2014).The structure of the J ijkl tensor in ( 16) requires additional studies for different rock compositions and fluid properties.When processing this experiment, we used the following approximation of this tensor along the main directions of the tensor r kl : where I 1 (r) is the first invariant of the tensor r ij , the parameters а i , ε i are determined experimentally, while а i = J i * (s i * ).The integral of the capillary pressure distribution function is the most informative indicator of wettability, which changes during the filtration of fluids and the motion of the interface.Therefore, the dependence of capillary pressure on the flow direction has a certain physical basis.

Conclusion
U s i n g t h e m e t h o d s o f n o n e q u i l i b r i u m thermodynamics, the motion equations of multiphase filtration (the generalized Darcy law) and the laws of interphase interaction are obtained and investigated.The possibility of direct determination of the coefficients of «multiphase» relative phase permeabilities is shown.For two-phase filtration, the effect of viscous surface interaction between fluids is studied; a method for determining the «cross-phase» permeabilities and their quantitative contribution to the field development indicators is given.
On the basis of modern research methods of anisotropic core material, a new effect of capillary pressure dependence on the measurement direction has been established.This means that capillary pressure, as well as relative phase permeability, is not a universal scalar function of saturation.The physical explanation of this fact may be due to a change in the wettability of the rock during the extraction of oil and gas.It is also known the emergence of «heterogeneous» wettability, when in some areas of the deposit the surface of the rock is hydrophilic, and on others it is hydrophobic.The most informative integral indicator of wettability is the integral of the distribution function of capillary pressure.

Fig. 7 .
Fig. 7. Dependence of capillary pressure on water saturation for X, Y directions and control sample