© M.F. El-Amin et al., published by IFP Energies nouvelles, 2020

Nomenclature

b_d : Factor of slippage

k_0,d: Intrinsic permeability

: Porosity of the matrix blocks

: Porosity of the fractures

: Mass density of the gas in the matrix blocks

: Mass density of the gas in the fractures

: Mass density of the core

: Gas molar weight

: Langmuir pressure

: Gas pressure in matrix blocks

: Gas pressure in fractures

: Critical pressure

: Average fractures pressure around the well

: Langmuir volume

: Mole volume under standard conditions

: General constant of gases

: Temperature

: Critical temperature

: Gas compressibility factor

: Well-radius

: Drainage-radius

: Set of normal fractures

: Mean stress in matrix blocks

: Mean stress in fractures

: Effective stress in matrix blocks

: Effective stress in fractures

: Biot’s effective parameter

: Fracture spacing in

: Fracture spacing in

: Fracture spacing in

Superscripts and subscripts

: Fractures

: Matrix

: or

1 Introduction

Finite Element Methods (FEMs) are effective numerical techniques for solving the complex engineering problems. FEMs appeared in the middle of the last century and were widely used in solid mechanical applications [1]. The FEM has been extensively investigated and developed in the last five decades [2], and now is used for highly nonlinear problems [3]. However, the standard FEM faced certain instabilities for solving elliptic problem that describe the flow in porous media [4, 5]. On the other hand, the Mixed Finite Element Methods (MFEMs) have succeeded in eliminating such instabilities [6, 7] as it may be extended to higher-order approximations as well as it is a locally conservative method [8]. Different physical variables can be treated differently in the MFEM from other variables. For example, the enhanced linear Brezzi-Douglas-Marini (BDM1) space is used for speed approximation, which requires more precision, while the piece-wise constant space is used for pressure approximation. The models of shale-gas transport were mainly developed by adapting the traditional models of flow in fractured porous media. Warren and Root [9] have developed an idealized geometric model to investigate the behavioral flow in fractured porous media. The dual porosity model including stress of fracture was used to investigate the flow of gas in the shale [10, 11], while, the dual-continuum model was used to describe gas flow Kerogen shale [12, 13]. Several versions of the dual-continuum model have been developed to include, for instance, adsorption, heat, deformation and Knudsen diffusion [14–17]. Several authors presented research work in shale gas reservoirs with geomechanical effects, such as Wang et al. [18], Lin et al. [19], and Yang et al. [20]. Wang et al. [18] developed a general model including embedded discrete fractures, multiple interacting continua and geomechanics in shale gas reservoirs with multi-scale fractures. In [19], authors provide insights on fracture propagation using reservoir flow simulation integrated with geomechanics in unconventional tight reservoirs. Effects of multicomponent adsorption and geomechanics are investigated in a shale condensate reservoir [20]. Girault et al. [21] have introduced a priori error estimates for a discretized poro-elastic-elastic system, in which the flow pressure equation is discretized by either a continuous Galerkin scheme or a mixed scheme, while, elastic displacement equations are discretized by a continuous Galerkin scheme. El-Amin et al. [22] have used the MFEM with stability analysis to simulate the problem of natural gas transport in a low-permeability reservoir without considering fractures. The authors [23] extended their work to cover fractured porous media and the rock stress-sensitivity with considered stability analysis of the MFEM. In this work, we present a theoretical basis with proofs of the stability analysis of the MFEM (in Ref. [23]) including the necessary lemmas and theorem.

2 Modeling and formulation

In this section, the mathematical model of the problem under consideration is developed. The Dual Porosity Dual Permeability (DPDP) model is employed to describe the gas transport in fractured porous media which consists of matrix blocks and fractures. The absorbed gas and the free gas coexist in the matrix blocks, however, the free gas exists in fractures only. The mass of free gas accumulation is , and the adsorbed gas accumulation (Langmuir isotherm model) is [24, 25],(1)

The mass density of the gas can be written as,(2)

The gas compressibility factor, which is given by Peng-Robinson equation of state [26],(3)(4)(5)

In low-permeability formations, Klinkenberg effect, which is described by the apparent permeability, takes place and written as,(6)

The DPDP model of gas transport in fracture shale strata is represented as,(7)and(8)where(9)and(10)such that and is the matrix and fractures porosities, respectively,(11)

Moreover, based on the following definitions, and are treated as constants [17],(12)and(13)

The Knudsen diffusion coefficient is defined as,(14)Here, is the transfer parameter which is connecting of the matrix-fracture domains, and defined as,(15)

One may provide the definition of the source of the production well by,(16)

Given the constant, , the drainage-radius is represented as,(17)

The coefficient of crossflow between the matrix and fracture domains is defined as [9],(18)such that is given by,(19)

2.1 Geomechanical effects

Invoking the rock stress-sensitivity [27], the porosity can be expressed as a function of the mean effective-stress,(20)(21)From (20) and (21), we can obtain,(22)

As the porosity increases the difference becomes positive and vice versa. Also, the coefficient becomes small based on the change in the porosity, and it becomes positive if porosity increases and negative if porosity decreases.

2.2 Initial and boundary conditions

The initial conditions can be represented by,(23)

The pressure on the production well (boundary condition) is given as,(24)

The no-flow boundary conditions on matrix are given as,(25)while the no-flow boundary conditions for fracture domain are,(26)

3 MFEM spaces

The MFEM contains two spaces for a scalar variable and it’s flux. The MFEM was developed to approximate both variables simultaneously and to give a higher order approximation for the flux. A compatibility condition must hold for the two spaces to insure stability, consistency and convergence of the mixed method and by considering more constraints to the numerical discretization.

Now, let us define the inner product in as,(27)and the inner product on as,(28)Given,(29)is the largest Hilbert space, such that and are well defined, therefore, . It is needed that and . The Hilbert space with norm given by,(30)(31)

The seeking solution is,(32)such that and are smooth. For the approximate numerical solution, the two spaces become,(33)and,(34)

Thus, the normal component is continuous across the inter-element boundaries.

Moreover, it is useful to present the elements that are designed to approximate [7], which satisfies,(35)and,(36)where is discontinuous piecewise constant and is piecewise linear.

4 Mixed finite element approximation

Let and are, respectively, the matrix and fracture in a polygonal/polyhedral Lipschitz domain (on which we define the standard space such that ), with the boundaries, and . One may write the above dual-porosity model in the following general form,(37)(38)(39)(40)where(41)and(42)

The functions and are moved to the left hand side to avoid discontinuity when we integrate and by parts. Selecting any and , the mixed finite element weak formulation can be written in the following form,(43)(44)(45)(46)

Now, let the approximating subspace duality and be the -th order () Raviart–Thomas space (RT_r) on the partition . The mixed finite element formulations are stated as below: find and such that,(47)(48)(49)(50)for any and .

In order to obtain an explicit formulation for the flux, we employ a quadrature rule along with the MFE method to [8]. The total time interval is divided into time steps with length . The superscript denotes for the current time step, while denotes for the previous one. We use a backward Euler semi-implicit discretization for the time derivative terms. The following scheme has been developed,(51)(52)(53)(54)Given and , the numerical procedure to calculate pressure and velocity is presented here,

1. Update the thermodynamical variables explicitly.

2. Solve the equations (51) and (52) to obtain and .

3. Solve the equations (53) and (54) to obtain and .

5 Stability analysis

In this section, we carry out the stability analysis of the proposed MFE method, which ensures that the discrete solutions are bounded in a physically reasonable range. The key issue encountered in the stability analysis is the nonlinearity of the matrix and fracture pressures. In order to resolve this issue, we need to define some auxiliary functions and analyze their boundedness. Now let us define the following functions,(55)(56)

Therefore, one may(57)(58)(59)and(60)

is a special function called the exponent integral function. As stated above in the model formulation that the coefficient is positive in the case of increasing porosity and negative in the case of decreasing porosity.

Lemma 5.1

For the case of increasing matrix porosity, , and sufficiently large , there exists a positive constant,(61)such that,(62)

Proof. Note that the value of the exponent integral function is negative and (small for big values of ), therefore the quantity, has always a positive value which has lower and upper bounds, namely,(63)where . On the other hand, in shale reservoir, it is well known that the initial pressure has the maximum value, i.e., , while, the pressure of the well has the minimum value, i.e., . Therefore, we have,(64)(65)and,(66)where . Therefore, for the case of increasing matrix porosity, , and sufficiently large , and holding (63)–(66), the following coefficient is positive, i.e.,(67)(68)Therefore,(69)and this completes Lemma 5.1 proof.

Lemma 5.2

For the case of decreasing matrix porosity, , assuming,(70)there exists a positive constant , such that,(71)

Proof. It is clear that, , so, we always have, . Holding the assumption (70), one can easily prove this lemma in a similar way to Lemma 5.1.

Lemma 5.3

For the case of increasing matrix porosity, , and for sufficiently large , there exist two positive constants,(72)and(73)such that,(74)

Proof. Holding (63)–(66) with for sufficiently large , and integrate (58) over for the case of increasing porosity, , it is easy to prove (74), which completes the lemma proof.

Lemma 5.4

For the case of decreasing matrix porosity, , with assuming that,(75)and(76)one can prove,(77)

Proof. It is clear that, , so,we always have, . Holding the assumptions (75) and (76), therefore we find, , and integrate (58) over for the case of decreasing porosity, , which proves (74), and this completes the proof.

Lemma 5.5

For both cases of increasing or decreasing fracture porosity, namely, or , and sufficiently large , there exists a positive constant,(78)such that,(79)

Proof. Again, we have,(80)where . Therefore, for the case of increasing fracture porosity, , and sufficiently large , Therefore,(81)

As stated above, for the case of the decreasing porosity, , the coefficient remains positive, then the above inequality holds. This completes the lemma proof.

Lemma 5.6

For both cases of increasing or decreasing fracture porosity, namely, or , and sufficiently large , there exists a positive constant,(82)such that,(83)

Proof. Similar to the previous proofs.

Lemma 5.7

Assuming that,(84)then, is bounded.

Proof. Since is bounded by and , i.e., , then, , such that is positive number. Holding the conditions (84), one may find that,(85)

Theorem 8

Holding the results in the above lemmas ( 5.1 – 5.6 ), namely, (86) and (87) with assuming that, (88) and (89)

Moreover, suppose that . Then,(90)

Proof. Let in (47) and in (48), and sum the two equations, we get,(91)

It follows from the definitions of and that,(92)

Integrating (92) from 0 to () yields,(93)

It is similar to obtain,(94)

There exists a positive constant such that,(95)

Thus, it is obtained from, (93)–(95) that,(96)Finally, (90) is obtained by 1Gronwall’s lemma.

Remark 5.1. In Lemma 5.6, if we do not use the coefficient, , and just use the coefficient, , then the coefficient of the term, in will be moved from into . Thus, the definitions of both and will be slightly changed without loss of the stability concept.

6 Conclusion

In this work, stability analysis of the MFE solution of the gas transport in a low-permeability fractured reservoir with stress-sensitivity effect has been considered. The dual-porosity dual permeability model with the slippage effect and the apparent permeability have been used to describe the flow in a low-permeability fractured porous media. Stability analysis of the MFEM is presented theoretically and numerically. We proved a seven lemmas and a theorem on the stability of the MFEM. Stability conditions are stated and estimated. Variation in the porosity is correlated to the stress-sensitivity effect and depends on the values of the corresponding physical parameters. For example, the coefficient is relatively small depending on the variation of the porosity. It has a positive value in the case of porosity increasing while it has a negative value in the case of porosity decreasing. Lemmas 5.1 and 5.3 discuss the positive case, while Lemmas 5.2 and 5.4 discuss the positive case. Lemmas 5.5 and 5.6 discuss the two cases.

References

1. Zienkiewicz O.C., Taylor R.L. (2000) The finite element method (5th edn.), vol. 1 – The basis, Butterworth-Heinemann, Oxford.
2. Bathe K.J. (1996) Finite element procedures, Prentice Hall, Englewood Cliffs, NJ.
3. Zienkiewicz O.C., Taylor R.L. (2000) The finite element method (5th edn), vol. 3 – fluid dynamics, Butterworth-Heinemann, Oxford.
4. Brezzi F., Fortin V. (1991) Mixed and hybrid finite element methods, Springer-Verlag, New York.
5. Chen Z. (2010) Finite element methods and their applications, Ch. 3. Springer-Verlag Berlin and Heidelberg GmbH & Co. KG.
6. Brezzi F., Douglas J., Marini L.D. (1985) Two families of mixed finite elements for second order elliptic problems, Numer. Math. 47, 217–235.
7. Raviart P.A., Thomas J.M. (1977) A mixed finite element method for 2nd order elliptic problems. Mathematical aspects of finite element methods, in: Proc. Conf., Consiglio Naz. delle Ricerche (C.N.R.), Rome, 1975, Lect. Notes Math., Vol. 606, Springer, Berlin, pp. 292–315.
8. Arbogast T., Wheeler M.F., Yotov I. (1997) Mixed finite elements for elliptic problems with tensor coefficients as cell-centered finite differences, SIAM J. Num. Anal. 34, 2, 828–852.
9. Warren J.E., Root P.J. (1963) The behavior of naturally fractured reservoirs, Soc. Petrol. Eng. J. 3, 245–255.
10. Ozkan E., Raghavan R., Apaydin O. (2010) Modeling of fluid transfer from shale matrix to fracture network, in: Annual Technical Conference and Exhibition, Lima, Peru, 1–3 December. SPE-134830-MS.
11. Ertekin T., King G.R., Schwerer F.C. (1986) Dynamic gas slippage: A unique dual-mechanism approach to the flow of gas in tight formations, SPE Form. Evalu. 1,1, 43–52.
12. Bustin A., Bustin R., Cui X. (2008) Importance of fabric on the production of gas shales, in: Unconventional Reservoirs Conference in Colorado, USA, 10–12 February. SPE-114167-MS.
13. Moridis G., Blasingame T., Freeman C. (2010) Analysis of mechanisms of flow in fractured tight-gas and shale-gas reservoirs, in: Latin American and Caribbean Petroleum Engineering Conference, 1–3 December. SPE-139250-MS.
14. Shu W.Y., Fakcharoenphol P. (2011) A unified mathematical model for unconventional reservoir simulation, in: EUROPEC/EAGE Annual Conference and Exhibition in Vienna, Austria, 23–26 May. SPE-142884-MS.
15. Kazemi H. (1969) Pressure transient analysis of naturally fractured reservoirs with uniform fracture distribution, Soc. Pet. Eng. 9, 4, 451–462.
16. Guo C, Wei M, Chen H, He X, Bai B (2014) Improved numerical simulation for shale gas reservoirs, in: Offshore Technology Conference-Asia, 25–28 March 2014, Kuala Lumpur, Malaysia.
17. Javadpour F. (2009) Nanopores and apparent permeability of gas flow in mudrocks (shales and siltstone), J. Can. Pet. Technol. 48, 8, 16–21.
18. Wang K., Liu H., Luo J., Wu K., Chen Z. (2017) A comprehensive model coupling embedded discrete fractures, multiple interacting continua and geomechanics in shale gas reservoirs with multi-scale fractures, Energy Fuel 31, 7758–7776. doi: 10.1021/acs.energyfuels.7b00394.
19. Lin M., Chen S., Mbia S., Chen Z. (2018) Application of reservoir flow simulation integrated with geomechanics in unconventional tight play, Rock Mech. Rock Eng. 51, 315–328. doi: 10.1007/s00603-017-1304-1.
20. Yang S., Wu K., Xu J., Chen Z. (2018) Roles of multicomponent adsorption and geomechanics in the development of an Eagle Ford shale condensate reservoir, Fuel 242, 710–718. doi: 10.1016/j.fuel.2019.01.016.
21. Girault V., Wheeler M.F., Almani T., Dana S. (2019) A priori error estimates for a discretized poro-elastic–elastic system solved by a fixed-stress algorithm, Oil Gas Sci. Technol. - Rev. IFP Energies nouvelles 74, 24. doi: 10.2516/ogst/2018071.
22. El-Amin M.F., Kou J., Sun S. (2018) Mixed finite element simulation with stability analysis for gas transport in low-permeability reservoirs, Energies 111, 208–226.
23. El-Amin M.F., Kou J., Sun S. (2018) Numerical modeling and simulation of shale-gas transport with geomechanical effect, Trans. Porous Med. 126, 779–806. doi: 10.1007/s11242-018-1206-z.
24. Cui X., Bustin A.M.M., Bustin R.M. (2009) Measurements of gas permeability and diffusivity of tight reservoir rocks: Different approaches and their applications, Geofluids 9, 3, 208–223.
25. Civan F., Rai C.S., Sondergeld C.H. (2011) Shale-gas permeability and diffusivity inferred by improved formulation of relevant retention and transport mechanisms, Trans. Porous Med. 86, 3, 925–944.
26. Firoozabadi A. (2015) Thermodynamics and applications in hydrocarbon reservoirs and production. McGraw-Hill Education – Europe, United States.
27. Terzaghi K. (1936) The shearing resistance of saturated soils and the angle between the planes of shear, in: Proceedings of International Conference on Soil Mechanics and Foundation Engineering, Harvard University Press, Cambridge, MA., Vol. 1, pp. 54–56.