MODELING THE PROCESS OF OIL DISPLACEMENT BY A HEAT CARRIER CONSIDERING THE CAPILLARY EFFECT

The development of oil fields is accompanied by a series of factors that affect the efficiency of extraction. Complexity of the oil production process necessitates the development of new and the modification of existing methods that forecast the impact of different procedures at their intensification. This applies to both the geological characteristics of an oil-bearing bed and its technical parameters [1]. In particular, flooded layers in which the temperature regime is disrupted demonstrate the so-called hydrodynamic and capillary effects. It is known [2] that a decrease in the temperature of oil-bearing beds increases the viscosity of a fluid. This, in turn, leads to the saturation of collectors, which cause additional resistance to fluid flow in a porous medium. Over the past decades, extensive experience has been gained in using such intensive systems to develop deposits as rock hydraulic fracturing (RHF) [3‒6]. The result is that the cracks formed expand the area of influence at operating wells, creating a relationship with zones of high permeability. In the case of heavy oil, for example, bitumen, effect of the hydraulic fracturing technique is strengthened by feeding different working reagents to an oil-bearing bed under pressure: hot water, steam, a mixture of different surface-active substances. The combination of the described factors that influence an oil-bearing bed makes it possible to reduce the areas of an immovable fluid, the so-called stagnation zones, and increase output from operating wells. It is relevant, while difficult, to account for each of them in order to build a mathematical model that is adapted to modern problematic oil fields. As well as, consequently, an appropriate mathematical apparatus to solve boundary value problems with the subsequent optimization of their parameters.


Introduction
The development of oil fields is accompanied by a series of factors that affect the efficiency of extraction. Complexity of the oil production process necessitates the development of new and the modification of existing methods that forecast the impact of different procedures at their intensification. This applies to both the geological characteristics of an oil-bearing bed and its technical parameters [1]. In particular, flooded layers in which the temperature regime is disrupted demonstrate the so-called hydrodynamic and capillary effects. It is known [2] that a decrease in the temperature of oil-bearing beds increases the viscosity of a fluid. This, in turn, leads to the saturation of collectors, which cause additional resistance to fluid flow in a porous medium. Over the past decades, extensive experience has been gained in using such intensive systems to develop deposits as rock hydraulic fracturing (RHF) [3][4][5][6]. The result is that the cracks formed expand the area of influence at operating wells, creating a relationship with zones of high permeability. In the case of heavy oil, for example, bitumen, effect of the hydraulic fracturing technique is strengthened by feeding different working reagents to an oil-bearing bed under pressure: hot water, steam, a mixture of different surface-active substances.
The combination of the described factors that influence an oil-bearing bed makes it possible to reduce the areas of an immovable fluid, the so-called stagnation zones, and increase output from operating wells. It is relevant, while difficult, to account for each of them in order to build a mathematical model that is adapted to modern problematic oil fields. As well as, consequently, an appropriate mathematical apparatus to solve boundary value problems with the subsequent optimization of their parameters.

Literature review and problem statement
Paper [7] gives an analytical solution to the boundary value problem in which cracks are represented in the form of a section of zero thickness and finite conductivity. However, such a model does not reflect the actual filtration properties of the displacement process. One of the first studies whose author investigated the process of filtering through an elliptical crack is reported in [8]. Such a model takes into consideration only the perfect process of hydraulic fracturing. It is known that crack formation is affected by a non-uniform environment around the well, as well as the existing cracks from hydraulic fracturing. Specifically, work [9] proposed a study into the mutual influence on the displacement process exerted by the arrangement of multiple cracks from RHF at a single operating well. In this case, the authors calculated filtration flow rate at oil-extracting well, but the issue that remains to be solved is finding a saturation field that would predict the flooding rate of operating wells. That would also make it possible to define patterns in the operation of any deposit under the condition of the projected location of wells and RHF cracks in them.

S . Y a r o s h c h a k
Great interest has been paid to the process of a twophase filtration, for example, in work [10]. However, accounting for overflows between fluids has remained an unresolved issue until now. Studying the interaction among several fluids and extracted oil raises a question regarding the impact of capillary forces on the bottom-hole zone of wells. In particular, paper [11] attempted to form a basis for the further study into the flow of a fluid in the process of thermal displacement. In addition, article [12] has in more detail investigated the three-phase non-isothermal oil filtering for the case of vapor-gravity drainage. As regards the functions of capillary pressure p ow and p go , which depend on the saturation and temperature of respective phases (steam, condensed water, and oil), they are considered to be known. Authors of work [13] take into consideration the pore sizes of the environment where filtering is performed. It was found that the adverse ratio of hydrodynamic to capillary forces leads to forming stagnant zones not covered by flooding. However, the unresolved issue is the distribution of an oil-bearing bed saturation and a velocity field. It is useful to consider a change in the viscosity of oil when increasing the share of water at different stages of oil fields development, as explored in paper [14]. Thus, there is a need to solve the wider classes of problems in order to account for a capillary effect at non-isothermal multi-phase filtering in the inhomogeneous oil-bearing beds. In this case, it is necessary to calculate a hydrodynamic grid in order to find the fields of velocities, distribution of temperature and saturation at any node of the grid. As well as the coordinates for critical points and values for a filtration flow rate in an element of symmetry, for example a five-point system of flooding.

The aim and objectives of the study
The aim of this study is to model mathematically the processes of fluids displacement in oil-bearing beds considering the capillary effect and influence of cracks from hydraulic fracturing by using numerical methods of quasiconformal mappings.
To accomplish the aim, the following tasks have been set: -to adapt a mathematical model of the process of oil displacement from porous inhomogeneous sedimentary rocks for the case of non-isothermal multiphase filtering considering a capillary effect; -to advance a procedure and to build the appropriate numerical algorithm for solving respective boundary value problems; -to perform numeric calculations and analyze the results obtained.

Modeling the impact of a capillary effect on the displacement process in oil-bearing beds
Consider the process of non-isothermal displacement of oil by a heat-carrier (such as water) in a horizontal bed, induced by the difference in pressure at injection and operating wells. Fig. 1 shows the element of a bed whose outer contour L is impenetrable or is a feed path. This bed contains * 1, n α = forcing (•) and * 1, n b = operational (○) wells whose contours are, respectively, L α . Pumping and operating wells along with an external contour define the area z G of fluids filtering.
An oil-bearing bed is considered to be non-uniform, which in turn makes it possible to take into consideration the existence of cracks from rock hydraulic fracturing (RHF). In addition, it is taken into consideration that the dynamic viscosity of phases changes with a change in temperature. In the process of displacement, the fluid movement is slow and proceeds without phase transitions. The functions of relative phase permeabilities and capillary pressure are known and represent unambiguous saturation functions. Taking into consideration the capillary forces, law of filtration [15] for two phases is recorded in the form of continuity equation and motion equation where p l is the pressure in phases, µ l (T) is the dynamic viscosity, T is temperature, l υ is the filtration rate of the l-th phase, s l is the saturation of bed with the l-th phase; σ, k=k(x, y) are the coefficients of porosity and absolute permeability of soil; are the relative phase permeability (s=s w ). Indexes "w" and "o" denote the magnitudes that characterize water and oil. The difference of pressure in phases is considered to be equal to the capillary pressure: θ s χ is the coefficient of surface tension, θ is the boundary wetting angle, J(s) is the dimensionless Leverett function. Hence, considering the total velocity w o υ = υ + υ of the filtration flow, we have: where φ=φ(x ,y ,t) is the filtering rate quasi-potential, To determine pressure p o via φ, perform the following transformations: For the case of accounting for the cracks from hydraulic fracturing, similar to [16,17], we shall require that the condition for the continuity of flow and pressure should be satisfied at their border. The respective coefficients of absolute permeability of soil are represented in the form where D α is some area of the bed, which corresponds to a crack with index α. Initial distributions of water saturation in a bed and its value for pumping wells is denoted, respec- To describe the process of redistribution of heat between the phases and the skeleton, we accept a single-temperature model, according to which there is an instantaneous transfer of heat from a fluid to the skeleton and in the opposite direction. Thus, to calculate a thermal field, we use the following equation [2]: where

C(s)=σρ o c o (1-s)+σρ w c w s+(1-σ)ρ e c e
is the volumetric heat capacity of a porous medium, c o , c w , c e , ρ o , ρ w , ρ e are specific heat capacity and density for oil, water, and a bed's skeleton, respectively. By using the formulae for determining the velocities of movement of oil and water [18,19]: To determine a temperature field, we shall accept that its distribution at initial moment T(x, y, 0)=T°(x ,y), as well as the value for temperature at the region's contours, are assigned.
Algorithm to solve the problem. We shall employ the procedure, proposed in [19], for solving the nonlinear problems on two-phase filtering in oil-bearing beds. The region of the displacement process is multi-connected. In this case, it is necessary to use the method of numerical quasiconformal mapping and basic concepts from a procedure for the stagewise registration of characteristics of the environment and the process. The goal is a transition from a direct problem to inverse, namely from the physical domain to the simply connected region of complex quasi-potential. Thus, we introduce the series of conditional sections of region G z along the flow lines that pass through the critical points.
Upon such a transform of region G z and following the construction of a respective region of the complex quasi-potential, it is necessary to rewrite the problem on non-isothermal multiphase filtering relative to the new estimated region -G ω . By applying transition formulae derived in [21], we obtain: , Equations (7) to (10), depending on the chosen configuration of the filtration region, shall be supplemented with appropriate boundary and initial conditions [21].
The algorithm for solving the corresponding boundary value problem shall be represented in the form of the following sequence of steps: 1) calculate the field of velocity potential based on the current fields of saturation and temperature, considering the geometry of filtration region, initial and boundary conditions; 2) solve the problem on quasiconformal mapping: build a hydrodynamic grid, characteristic lines of the flow separation, find quasi-potential φ, flow rate and other unknown filtering parameters; 3) find a redistribution of saturation (based on the established filtering characteristics, according to (9)); 4) find a redistribution of temperature (based on the established filtering characteristics and the recalculated field of saturation, according to (10)); 5) check the conditions for terminating the algorithm operation, whose violation leads to the repeated clarification of step 1 of this algorithm. One of such conditions for termination could be the condition of exceeding the permissible share of a displacing liquid in the output from an operational well.
The difference analog of the problem and its solving algorithm are constructed similarly to the procedure from [5]. By introducing to the field of the complex quasi-potential of the uniform orthogonal grid with nodes at points (φ i , ψ j ) ( 1, , i n = 1, j m = ), relative to which the approximation of equations (7) to (10) υ is the velocity (found similarly to work [5]). By consequently selecting a time-dependent step, the parameters for splitting the field of complex quasi-potential (position of the grid's nodes (φ i , ψ j )), the initial approximations of coordinates for the boundary and inner nodes at a hydrodynamic grid, we derive values for the quasiconformal invariants. Next, we refine the coordinates for internal nodes at a hydrodynamic grid by solving the respective (8) difference analogs. After this, we adjust the boundary nodes under the conditions for fixing the surrounding boundary and near-boundary conditions, applying the orthogonality conditions, and derive the approximation of flow rate magnitudes (well production rates). The conditions for finalizing the algorithm for constructing a hydrodynamic grid (finding the unknown filtration parameters, in particular a velocity field) at a given iterative stage include: the stabilization of flow rates, the stabilization of boundary nodes, etc. By employing the constructed velocity speed υ i,j , the fields of saturation s i,j and temperature , i j T from the preceding iterative step (taking into consideration the boundary conditions), we find the distribution of saturation , i j s in a bed at a given point in time. According to (11) and using (12), we calculate a temperature field , i j T , and recalculate then the field of velocity and potential.

Numerical calculation of the model problem
We consider a problem on the non-isothermal two-phase filtration in the symmetry element of a five-point flooding system, for the following model's parameters:  Fig. 2 shows the hydrodynamic grid in region G z , where a single operating well is exposed to the action from 4 pumping wells. The field of temperatures when oil is displaced with water for the case of non-isothermal filtering process is shown in Fig. 4.
The respective saturation field considering the impact of capillary effect is shown in Fig. 5.  proposed in paper [6] for the selected element of symmetry under respective conditions.

Discussion of results of studying the processes of displacement in oil-bearing beds in the elements of areal flooding
The mathematical model (1) to (5), constructed in the current study, in a combination with the built algorithm solve the problem on establishing the field of velocity, saturation, and temperature ( Fig. 3-6) under the respective operating conditions at wells. The model and the algorithm imply a possibility to account for the capillary effect, different kinds of inhomogeneous inclusions, such as RHF cracks and overflows between wells. This is the basis for establishing the optimal location of wells and for conducting hydraulic fracturing at them.
When carrying out numerical experiments in the symmetry element of a five-point flooding system it was found that the front of displacement in case of not taking into consideration the capillary effect propagates by 20 % faster than when accounting for it. This, in turn, would lead to an increase in the actual time of water breaking through to the operating well, which would not make it possible to timely perform appropriate waterproofing operations.
The qualitative pattern of filtration flow, the region of low velocities, and the redistribution of flows are demonstrated by the hydrodynamic grid (Fig. 2), which is computed at each point in time. In addition, the grid shows the influence of wells on one another, and on the process of flooding in particular.
Construction of a hydrodynamic grid by using other existing methods, such as finite elements or differences, would require significantly larger computation efforts. In a given case, it is a key element of the built algorithm and is an integral part of the solution to the set problem.
It is also worth noting that the method for solving the problems on multiphase non-isothermal filtering, proposed in the current paper, is to be used for the case when it is possible to introduce a quasi-potential function φ and the corresponding conjugated thereto flow function ψ. Introducing them is not always possible and requires additional research. However, when it becomes possible, this greatly simplifies the calculation of fields of saturation and temperature. In this case, the spatial dimensionality of the respective subproblems reduces; solutions to them are found along the lines of the flow and it is one of the advantages of the method. Such an approach, in contrast to conventional methods, of finite differences and finite elements, does not require a separate recalculation of the field of velocities and employing the interpolating and gradient methods to construct the equipotential lines and the lines of the flow.
It is also worth noting that the devised model does not take into consideration the possibility of inter-phase transition of a fluid during the displacement process. Therefore, the temperature of an environment and a fluid should not exceed its critical value at which inter-phase transitions occur.
In the future, there is an obvious prospect to advance both the mathematical model taking into consideration the inter-phase transitions of a fluid and the procedure for solving respective boundary value problems.

1.
We have improved a mathematical model of oil displacement from porous sedimentary rocks in the presence of cracks from hydraulic fracturing. The model accounts for the reverse effect of process characteristics on the initial characteristics of an environment susceptible to deformation in a horizontal oil-bearing bed, constrained by the lines of the flow and the equipotential lines.
2. We have advanced a procedure for solving the respective boundary value problems provided the capillary effect is considered. Specifically, a formula has been proposed for the introduction of quasi-potential, using which helped reduce the spatial dimensionality of subproblems, to derive the fields of saturation and temperature. Their solution has been obtained along the lines of the flow.
3. We have constructed a numerical algorithm for solving respective boundary value problems. It was established in the process of numerical calculations that the front of displacement in case the capillary effect is not taken into consideration propagates 1.2 times faster than that when it is accounted for. This, in turn, could lead to an increase in the actual time of water breaking through to the operating well, which would not make it possible to timely conduct appropriate waterproofing operations.