Scientific engineering as the basis of modeling processes in field development

M. M. Khasanov; A. N. Sitnikov; A. A. Pustovskikh; A. P. Roshchektayev; N. S. Ismagilov; G. V. Paderin; E. V. Shel

conference proceedings Doi: https://doi.Org/10.18599/grs.2018.3.142-148

Scientific engineering as the basis of modeling processes

in field development

M.M. Khasanov, A.N. Sitnikov, A.A. Pustovskikh*, A.P. Roshchektayev, N.S. Ismagilov,

G.V. Paderin, E.V. Shel

Gazprom Neft Science and Technology Center, Saint-Petersburg, Russian Federation

Abstract. Three characteristic examples of the use of scientific engineering approaches for managing the technological processes of reservoir modeling at different hierarchical levels are presented in the article.

The first example demonstrates the application of the spectral approach for modeling geophysical fields -the available log data is decomposed from the spectrum of Legendre polynomials, after which a stochastic field of the expansion coefficients is constructed. The results obtained by this method of realizing geophysical fields correspond to real data in a wider area of modeling than in classical methods. The speed of building models also increases due to the convenience of parallelization.

The second example demonstrates the use of the source method to optimize the transfer of wells into injection. A flow rate of each well is found by simulating the wells in the development system as linear sources or sinks and recording the resulting system of equations for the flows at each time. According to the optimum of discounted extraction, there is an economically efficient time for well development.

The third example demonstrates the application of the theory of dimensions to the problem of hydraulic fracturing modeling to determine the significance of certain parameters for the design of fracturing. By measuring the dependence of the fracture length on the injection volume, we obtain an empirical formula for the fracture length from its parameters, which determines the level of their significance.

Keywords: geological modeling; geostatistics; geomechanics; hydraulic fracturing; hydrodynamic modeling; source method

Recommended citation: Khasanov M.M., Sitnikov A.N., Pustovskikh A.A., Roshchektayev A.P., Ismagilov N.S., Paderin G.V., Shel E.V. (2018). Scientific engineering as the basis of modeling processes in field development. Georesursy = Georesources, 20(3), Part 1, pp. 142-148. DOI: https://doi.org/10.18599/ grs.2018.3.142-148

The tasks of monitoring and controlling technological processes in field development often lead to the need to simulate geophysical fields, processes associated with filtration in the reservoir and movement through pipes of multiphase mixtures, geomechanical tasks during drilling and during hydraulic fracturing, and many others. Traditionally, the description of oil and gas production processes is carried out on the basis of differential equations of motion of liquids and gases in porous media and pipes. However, this approach does not allow one to describe many of the essential properties of the formation. Like any large system, oil and gas objects require the use of a whole hierarchy of models - from differential to integral, from deterministic to adaptive, capable of describing not only different levels of organization of systems, but also the interaction between these levels. The solution of all these tasks

*Corresponding author: Alexey A. Pustovskikh E-mail: [email protected]

in hierarchical modeling of the processes of the oil and gas industry is engaged in scientific engineering (SciencEngineering) - the field of science, located at the junction of system engineering, basic sciences and cybernetics (Fig. 1).

In this paper, we will consider examples of the approaches of scientific engineering to solving problems that have either not been solved previously in the oil

Fig. 1. Scientific Engineering (SE) - at the junction of areas of various hikes

and gas industry or have been solved in other fields of science for other types of problems. In other words, it will be shown how examples of solving problems of fundamental physics, chemistry, and other natural sciences can be applied to solving problems of oil and gas production.

Geological modeling: application of the spectral approach for modeling geophysical fields

The development of a spectral approach to the problems of geological modeling was first presented in (Baykov et al., 2010), in which the authors proposed a method based on decomposing the logging into coefficients in an orthonormal basis in the space of functions integrable with the L2 square. The method was further developed in works (Baykov et al., 2012; Khasanov et al., 2015), in which a more detailed study of the theoretical foundations of the spectral method was carried out, as well as some results of calculations based on the spectral method of geological modeling were presented. It is also worth mentioning that the spectral analysis of logging has been used for lithofacies analysis (Khasanov et al., 2014).

Spectral method of geological modeling - a method of modeling three-dimensional cubes of geophysical properties based on well data. The mathematical model adopted in the spectral method of geological modeling is represented by the simulated region D c R3 by the stochastic field G(x, y, h) defined on this region and defined on a certain probability space. For such a model, well data for vertical wells with a fixed lateral coordinate (x*, y*) is a random process, parameterized by a variable characterizing the depth G(x*,y*, h) =f(h). If the simulation involves N wells with coordinates (x., y.) at which the simulated property is given by functions f(h), then these functions can be considered as the known values of a certain realization of the stochastic field G(x, y, h).

The general principle of the spectral method consists in the sequential implementation of several steps: expansion ofthe functions f(h) according to some orthonormal basis, modeling of expansion coefficients in the interwell space, recovery of the simulated stochastic field using these coefficients at each point of the region D.

It is known that in the Hilbert space of square-integrable functions (L2) ortho-normal bases exist that allow one to decompose any function defined in this space into a series in terms of the expansion coefficients determined uniquely. The functions defining the simulated property as functions with finite energy belong to this space. Let us choose in the L2 space a basis of the Legendre polynomials P. (h) which we use to decompose the functions f(h):

fi(!0 = 'Zjcjpi(hl (1)

where the coefficients cj are defined as a scalar product in the space L2. The set of expansion coefficients cji for each level of decomposition j is defined at points (x,y), which are well coordinates, and are known values of some realization of the stochastic field c(x,y .).

To simulate the coefficients cj for all (x, y), the spectral modeling method is used (Prigarin, Mikhailov, 2005), based on the theorem on the integral representation of random processes and fields.

The construction of the conditioned implementation of the simulated c (x , y) field with known c. values is generally accepted for two-step simulation methods -the difference between kriging based on known data and data obtained as a result is constructed (Dyubryul, 2002).

As a result of the implementation of the above steps, we obtain a set of realizations of the stochastic fields of the expansion coefficients c (x ,y.), which are caused by the downhole values of cji. To restore the implementation ofthe entire field G(x,y, h) at each point ofthe simulated region D, summation ofthe basis functions is carried out according to the simulated coefficients:

G(x,y, h) = I; Cj(x,y)pi(h). (2)

By virtue of construction, the field implementation is caused by well data and reproducing the statistical characteristics of the source data.

The spectral modeling method has significant advantages over classical modeling methods, such as, for example, sequential Gaussian modeling. In particular, this method is devoid of the main limitation of classical methods - the hypothesis of stationarity of the modeled property. It is known that most geological processes occurring in nature are non-stationary. A significant advantage of the method is a nonparametric periodogram statistical analysis, instead of the classically used parametric variogram analysis, which does not allow a full assessment of the whole range of variability of the modeled property. The spectral method, by virtue of the algorithm of its construction, is fundamentally parallelizable, which makes it easy to scale calculations. In addition, the method implements building a model without reference to a grid (so-called grid-free simulation), which allows modeling on grids of any configuration and complexity, grinding part of the grid and refining the model in this area and, finally, combining generating implementations of the simulated property in various areas of the model. The performance gain of spectral modeling at a typical workstation of a geologist-modeler is 2 to 4 times as compared with classical methods, depending on the size of the generated model.

The application of the spectral method in practice demonstrated its ability to reproduce well the geophysical fields in the interwell space and in undrilled areas, especially when compared with traditional modeling methods (Fig. 2) (Khasanov et al., 2014; Prigarin, Mikhailov, 2005).

SCIENTIFIC AND TECHNICAL JOURNAL

Fig. 2. Comparison of real and synthetic gamma logging curves. DGK, red curve - real logging, SP_DGK, blue curve - synthetic logging modeled by the spectral method, SGS DGK, black curve - synthetic logging modeled by the SGS method.

Hydrodynamic modeling: the method of sources in the optimization of technological calculations

There are three main groups of methods for solving the differential equation of piezoconductivity (and the heat and diffusion equations that coincide with it) for filtering a fluid in a porous medium: analytical, numerical, and mathematical modeling methods. The analytical methods include: the classical direct integration method, the operational methods (the Laplace integral transform method) and the source method.

The classical method of direct integration consists in finding a set of partial solutions that satisfy the differential equation of piezoconductivity, and then their imposition (superposition), to find the function of interest. In this case, each of the particular solutions is sought, as a rule, in the form of a product of two functions, one of which depends on time, and the other on coordinates. In the underground hydrodynamics, the classical method has found very limited application due to purely mathematical difficulties arising from its use even for relatively simple technological schemes.

Operational methods (in particular, the Laplace integral transformation method) are widely used in petroleum engineering to solve the problem of filtering liquid and gas in a reservoir. In this case, it is not the function itself (the original) that is studied, but its modification (image), obtained by multiplying the original by the exponential function, and the image is integrated within certain limits. After solving the problem in the images, finding the original, i.e. the function describing the pressure field is produced by inverse transformation.

The third method, namely the source/drain method, is a flexible mathematical tool, convenient in engineering applications. It is relatively simple to use

the source method to write an integral that satisfies the piezoconductivity differential equation and boundary conditions; the next task is to calculate its value. In the classical and operational methods, the main task lies in finding the integral itself, which is much more difficult.

The method of sources can be used to solve two-and three-dimensional problems of unsteady filtering both to wells with simple geometry, and complex, with fracture and horizontal wells. In our terminology, a source is a point, line, surface, or volume from which fluid is drawn from a formation (or injected into a formation).

In the present paper, we consider an example of using the source/drain method to solve the problem of choosing the optimal time for injection wells (Sitnikov et al., 2015). The development time is the period of operation of the well in the production mode before putting it into injection.

The pressure distribution from a vertical well, for which the fluid dynamics is given, is written in the form of the Duhamel integral:

(3)

where is the coefficient of piezoconductivity,

P0 is the initial reservoir pressure, ^ is the oil viscosity, K is the permeability, h is the formation thickness, C is the overall compressibility, m is the porosity, r is the distance from the source, t is the calculation time, t is an integration variable.

The problem of determining the flow rate of a well fluid, acting at a constant bottomhole pressure, is reduced to finding the function from the integral equation (3).

For a piecewise constant flow rate of a fluid, expression (3) transforms into equation (4):

where qn is a piecewise constant flow rate at each step, Pcn is the bottomhole pressure, rc is the well radius, n is

the number of time steps, Ei(-x) = f — dx is the integral

function, At is the time step.

Note that, in deriving equation 4, a linear source solution was used, which, in turn, was obtained under the assumption that the wellbore radius is small. In the overwhelming majority of practical situations, this approximation is justified.

Using the principle of superposition (which is ensured by the linearity of the piezoconductivity equation), we can calculate the dynamics of unsteady flow rates for a system of N wells:

j = 1..N.

(5)

For fractures of finite conductivity, the system of equations (5) must be supplemented by a system of equations for calculating bottomhole pressures at sources simulating fractures.

The task of optimizing the time of the development of injection wells in production is to determine the duration of production of fluid from injection wells, at which the cumulative discounted production from the development element will be maximum:

maxfeJ^TV)}

(6)

Discounted production from a development element is understood as total discounted production from one injection well and the number of production wells corresponding to the selected development system.

The proposed approach allows one to determine the optimal time for testing injection wells with given reservoir properties of the reservoir, parameters of the development system and process parameters. This method formed the basis of the calculation module, which allows to calculate the dynamics of production of a production and injection well, as well as the dependence of accumulated discounted production from a development unit on the time of injection wells.

The reliability of the calculations by the proposed method was checked by a comparative analysis of the dynamics of fluid production with the results of calculations in a commercial software product and analytical dependencies obtained in the article (Khasanov et al., 2013). Testing was carried out on the basis of data on the structure of one of the fields of PJSC Gazpromneft. The obtained dependence of the dimensionless accumulated discounted production on the time of the injection well testing is shown in Fig. 3. It can be seen that the largest discounted production can be obtained subject to the use of injection wells in production for about nine months. An earlier or late transfer of wells would be less cost effective. This result was not only confirmed by a series of calculations on a full-scale hydrodynamic model, but also agrees with the characteristic time of well transfer, which was obtained experimentally during the operation ofthe field (Sitnikov et al., 2015).

The considered example shows that analytical and, as it is often said now, semi-analytical methods have not lost their relevance in solving engineering problems, even despite the rapid growth in the productivity of computing equipment and, therefore, numerical methods.

Fig. 3. The dependence of the dimensionless discounted production Q/max (QdcJ on the transfer time of the well to the injection T (Sitnikov et al., 2015)

Geomechanics: comparison of parameters of hydraulic fracturing in dimensionless variables of design and hydrodynamic studies

Due to the frequent use of hydraulic fracturing (GF) technology in oil and gas fields, there is a large amount of statistical information on conducted operations. Based on the results of processing this information, it is possible to conclude about the effectiveness of the hydrulic fracturing carried out, which will make it possible to make adjustments to the design and develop further recommendations. However, this task is complicated by the fact that, when conducting hydrodynamic studies in wells with hydraulic fracturing, it is found that the values of fracture parameters, in particular, half the length, differ significantly from those planned by design.

Below are the results of a study of the possible causes of this discrepancy using dimensionless variables, the introduction of which allows for the analysis of information on the performed fracturing operations.

From the graph shown in Fig. 4, it can be seen that in the general case there is no obvious relationship between the fracture length and the injection volume. The reasons that led to this are fairly obvious, and the fact is that, in addition to the injection volume of the cross-linked

%

* *

Fig. 4. Dependence of the fracture length in design and well testing on the volume of injected fluid

SCIENTIFIC AND TECHNICAL JOURNAL

gel (and the associated mass of the proppant), the formation parameters such as the Young's modulus for the section, the thickness of the layers compressing stresses perpendicular to the fracture, fracture toughness coefficients, as well as the technological parameters of the design of the hydraulic fracturing - fluid rheology, injection rate and proppant concentration.

To analyze the obtained data, it is required to find a sample of wells that would be completely identical in all parameters of the fracture except for one, after which this procedure should be repeated for all parameters of the fracturing to determine the degree of influence of each of them. Having determined these degrees of influence, it would be possible to establish the most significant parameters for the fracture geometry. This would suggest the parameters, the error in which could lead to a significant systematic discrepancy between the length of the fracture, planned for the design, and the length of the fracture, obtained by hydrodynamic studies of the well. In practice, however, it is not possible to select such samples of fractures that would not differ in all parameters except one.

All cases are selected so that the injected fluids of the hydraulic fracturing are identical, including the concentration of the polymer. Then we can talk about the similar rheology of fluid fracturing for these cases. Thus, at least for this parameter, the data of the fracturing operation are identical.

For the method proposed in this article, this is a prerequisite.

The design fluid injection rates also vary slightly in design, so that this factor could not have a decisive impact on the large scatter and formation of a "cloud" in Fig. 4.

Thus, mostly various geomechanical parameters remain that can vary along the reservoir quite significantly (especially the thickness of the reservoir). A successful statistical analysis requires a method that reduces the dimensionality of the problem according to these parameters, and leads to fractures in reservoirs with different geomechanical properties to the same "denominator". Such a method is the introduction of dimensionless parameters of the problem of the development of a fracturing.

The fracture propagation according to the Planar3D model in the case of the injection of non-Newtonian fluid is described by three laws (Khasanov et al., 2017):

- Hooke's law;

- The law of viscous friction;

- The law of conservation of mass.

We will conduct the de-dimensioning of the equations in the same way as described in (José I. Adachi et al., 2010) for the case of non-proppant fractures; we will also assume that there are no leaks in the reservoir, Young's modulus is homogeneous, and lithology is three-layer

and symmetric. As a result, we obtain the following list of dimensionless parameters:

k'E'2n+V

Y =

V = V

fj 371^0-271+2 <2n+2 kE Qn

H

3n+3A_2n+3

L = L

C =

,2n+l kE Qn ^371+1^271+2

3n 1 1 Qtf 2 +2A(j"+2

Ï ;TT~

<77 ,71 —I—

k2E Q2r2

K

Hü—

~ -\j H La

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

(7)

(8)

(9)

(10)

(11)

where E = —is the flat strain modulus, E is Young's

modulus, v is Poisson's ratio, n is the fluid behavior, k' is the flow density coefficient, Q is the fluid flow, H is the formation thickness, Ac is the stress contrast, L - half-length of a fracture, V - volume of injected fluid, C-leakage coefficient according to Carter, K - coefficient of fracture resistance.

Comparing the obtained dimensionless variables with those given in the article (José I. Adachi et al., 2010) for the Pseudo3D model, one can see that they are identical up to constants. The only important difference is the g parameter. It is clear from the work that if we divide the reservoir thickness by the length scale factor from this work, then the dimensionless parameter g is obtained, up to a constant. Its physical meaning is the ratio of the thickness of the reservoir to the fracture length attained.

As a result, the non-dimensioning decreases the dimension of the problem by 5, since the Young modulus, the Poisson's ratio, stress contrast, reservoir thickness and fluid viscosity are replaced by one dimensionless parameter g (7). Obtaining dimensionless parameters makes it possible to carry out an analysis of the work done by the hydraulic fracturing, namely, to identify some regularities when comparing the half-lengths obtained in design and well testing in comparison with the total injection volume. Thus, the necessary data for analysis in addition to the values of half-lengths of fractures are: reservoir thickness, Young's modulus and Poisson's ratio (included in the flat strain module), stress contrast in the reservoir, rheology of the injected fluid, fluid flow and total injection volume of the formation fluid. Strength properties of rocks in this analysis are not taken into account. Figures 5-7 show the dependence of the dimensionless length on the dimensionless volume for three types of injected fluid. From the graphs it can be seen that the predicted dependence of the dimensionless parameters is preserved, only the degree depending on the rheology of the injected fluid differs. In this case, the degree dependence is preserved both for the data

Fig. 5. Dependence of the dimensionless fracture length on the dimensionless volume: liquid 1, m/p A; a = 0.69

Fig. 6. Dependence of the dimensionless fracture length on the dimensionless volume: liquid 2, m/p A; a = 0.63

Fig. 7. Dependence of the dimensionless fracture length on the dimensionless volume: liquid 3, m/p A; a = 0.66

obtained by design and hydrodynamic studies.

The above results prove that the dependence of the dimensionless fracture length on the dimensionless volume is a degree dependence:

L = A-Va (12)

Then substituting expressions (7) and (8) into (12),

we can obtain a formula reflecting the dependence of the dimensional length of a fracture on the injection volume for an arbitrary degree a:

fj3n+l-(3n+3}a^a2n+2-(2n+3)a L = A k,i-aQn(1_a)E,2n+l-(2n+2)a V

(13)

From the graphs (Fig. 5-7), it can be concluded that the degree of dependence of the dimensionless fracture length on the injection volume a is in the range of 0.6-0.7. By taking the degree of dependence equal to 0.6, and the indicator of the behavior of a fluid n equal to 0.5, you can get an empirical formula for calculating the half-length of a fracture:

L = A

3

Act S

12 1 HSk'5QSE'S

3

TF5

(14)

This formula allows you to evaluate the effect of error (geomechanical parameters) or changes in input parameters (process parameters) on the change in the half-length of the fracture (Table 1).

Parameter ^parameter Geomechanical parameters Process parameters

H E' Ka k' Q

+10% - 2% - 2% 6% - 4% - 2%

+20% - 4% - 4% 12% - 8% - 4%

+30% - 6% - 6% 18% - 12% - 6%

Table 1. The effect of changing parameters on the change of fracture length

In this paper we analyzed the problem of the divergence of half-lengths of fractures in design and in hydrodynamic studies of the well. For the analysis we applied the method of dimensionless variables, developed on the basis of parametrization of the fundamental equations of hydraulic fracturing. The dependence of the half-length of the fracture on the injection volume, geomechanical parameters and rheology of the fluid is analyzed. This method allowed us to reduce the dimension of the problem and obtain a fairly universal empirical degree dependence of the dimensionless length on the dimensionless volume, which in its dimensional form gives a simple empirical formula for estimating the fracture length. It is concluded that the geomechanical parameters have a weak effect on the fracture length. The validity of the equations used in modeling hydraulic fracturing has been confirmed.

References

Baykov V. A., Bakirov N.K., Yakovlev A.A. (2010). New approaches to the theory of geostatistical modeling. Vestnik UGATU, 37(2), pp. 209-215. (In Russ.)

Baykov V.A., Bakirov N.K., Yakovlev A.A. (2012). Mathematical Geology. V. 1. Introduction in geostatistics. Moscow-Izhevsk: IKI, Biblioteka neftyanogo inzhiniringa, 228 p. (In Russ.)

Dyubryul O. (2002). Use of geostatistics to include data in the geological

SCIENTIFIC AND TECHNICAL JOURNAL

model. EAGE: SEG. 295 p. (In Russ.)

José I. Adachi, Emmanuel Detournay, Anthony P. Peirce. (2010). Analysis of the classical pseudo-3D model for hydraulic fracture with equilibrium height growth across stress barriers. International Journal of Rock Mechanics and Mining Sciences, 47(4), pp. 625-639.

Khasanov M.M., Belozerov B.V., Bochkov A.S., Fuks O.M., Tengelidi D.I. (2015). Automated lithologic facies analysis based on the spectral theory. Neftyanoe khozyaistvo = Oil Industry, 12, pp. 48-51. (In Russ.)

Khasanov M.M., Belozerov B.V., Bochkov A.S., Ushmaev O.S., Fuks O.M. (2014). Application of the spectral theory to the analysis and modelling of the rock properties of the reservoir. Neftyanoe khozyaistvo = Oil Industry, 12, pp. 60-64. (In Russ.)

Khasanov M.M., Melchaeva O.Yu., Sitnikov A.N., Roshchektaev A.P. (2013). Dynamics of hydraulically fractured wells production for economically optimal development systems. Neftyanoe khozyaistvo = Oil Industry, 12, pp. 36-39. (In Russ.)

Khasanov M.M., Paderin G.V., Shel E.V., Yakovlev A.A., Pustovskikh A.A. (2017). Approaches to modeling hydraulic fracturing and their development. Neftyanoe khozyaistvo = Oil Industry, 12, pp. 37-41. (In Russ.)

Prigarin S.M., Mikhailov G.A. (2005). Metods of numerical modelling of random processes and fields. Novosibirsk: IVMIMG SO RAN, 259 p. (In Russ.)

Sitnikov A.N., Pustovskikh A.A., Roshchektaev A.P., Andzhukaev T.V. (2015). A method to determine optimal switching time to injection mode for field development system. Neftyanoe khozyaistvo = Oil Industry, 3, pp. 84-87. (In Russ.)

About the Authors

M.M. Khasanov - DSc (Engineering), Director General

Gazprom Neft Science and Technology Center