Numerical Simulation of Fracture Permeability Change in Production of Pressure-sensitive Reservoirs with In-situ Stress Field

In pressure sensitive reservoirs, interaction effects among the porous media flow field, the fracture field and the stress field can cause some specific flow characteristics entirely different from those in conventional reservoirs. Dynamic fracture behavior is one of them, which generates a change in the value of fracture aperture and even a variation in the anisotropy of permeability. In this paper, we focus on the dynamic behavior of fractures and some affecting factors, including driving pressure and in-situ stress. Numerical discrete fracture network (DFN) models are built and solved by the finite element method to investigate what the range-ability the fracture presents and what impact these affecting factors have. In these mathematical models, both dynamic fractures and the fluid-solid coupling are taken into account, and a stress-strain model, a flow field model and a fluid-solid coupling model are included. Based on the models, the variation of fracture aperture in pressure sensitive reservoirs is studied and the results show that a different direction and connectivity of fractures lead fracture dilation to varying degrees as pressure changes so that the idea of anisotropic fracture porosity is proposed for reservoir scale simulation. The study also indicates that the drop of formation pressure determines the conductivity of fractures and anisotropy of permeability but just has a slight impact on the direction of principal permeability. Finally, the study shows the interaction of the in-situ stress pressure and the fracture field.


INTRODUCTION
Low permeable oil reservoirs represent more than 50% of newly discovered oil reservoirs around the globe, and a considerable number of them contain fractures.In this type of reservoirs, matrix blocks are tight and poorly permeable so they can only function as reservoir space.The fractures can serve as channels or baffles for fluids to flow, which is crucial for efficient and rational development of these reservoirs.
This type of reservoirs always presents two features: The first one is anisotropic permeability, especially common in fractured reservoirs.These reservoirs have been found in abundance, including the Spraberry Frend oil reservoir in Texas [1] and the Lisburne reservoir in Alaska [2], and the Jingan oil reservoirs [3] in Ordos Basin.In particular, China is a good representative in possessing anisotropic oil reservoirs.The main difficulty in such reservoirs is gathering insitu fracture attributes from limited information and assessing what the permeability is in these reservoirs [4].The study in this area has been performed by many researchers in recent years.Guo et al. [5] proposed a method for estimating directional fracture permeability from pressure data recorded during single well tests.Khamis et al. [6] and Chen et al. [7] both used interference and pulse well testing to analyze anisotropic permeability in full tensor.Liu et al. [8] proposed a comprehensive method to measure the anisotropic permeability of a core sample from formation.Graziano et al. [9] used a combination of image logs, WFTs and Mini DSTs for integrated formation evaluation in an anisotropic reservoir located in Indonesia.Kadet et al. [10] proposed a laboratory method for anisotropic reservoir properties determination.
The other feature is a fluid-solid coupling effect in pressure sensitive formation.This feature brings more challenges to reservoir engineers, and some researchers have also contributed in this area.Wang et al. [11] using an experimental method found that reservoir pore pressure fluctuates to some extent during a CO 2 flood, causing a change in effective confining pressure.The result is rock deformation and a reduction in permeability (particularly, with a reduction in fracture permeability), causing increased flow resistance in the fracture space.Dunayevsky et al. [12] used a numerical method to investigate the effect of water flooding to a stress field.Other researchers showed that fracture aperture and fracture permeability are not only related to flow pressure [13,14] but also are influenced by the in-situ stress field [15][16][17], even by temperature in non-isothermal reservoirs [18].
From the above analysis, a low permeable fractured reservoir is a complex system affected by flow in porous media, fracture shape and distribution, in situ stress field, and even a production approach, such as injection pressure and thermal production.To solve this multi-physics field coupled problem, it is necessary to include as many factors as possible.Only in this way can research be consistent essentially with the reservoir prototype in nature.Even faced with lots of hardship, researchers have made some progress.Min et al. [19] used a numerical method to study the stress-dependent permeability of fractured rock in consideration of stress fields and fracture networks.Wong et al. [20,21] used a simulator to investigate the behavior of dynamic permeability in reservoirs.
In this paper, we focus on understanding the impact of far-field stress and injection pressure on the dynamic anisotropic permeability behavior.In order to investigate this problem, we build numerical DFN (discrete fracture network) models by using the finite element method to solve multiphysics coupled problems and simulate the behavior of dynamic anisotropic permeability in anisotropic fractured reservoirs.A tensor theory is applied to analyze the results, with some discussions followed.

Solid Modeling
Here we have a 2D plane model.Displacement of the solid matrix is designated by the components u x , u y .The stress component in each direction satisfies According to constitutive equation, the normal strain is In the above model, the shear modulus is With no respect to body force, the constitutive equation set can be written as Navier-Stokes equation in porous media is = ( l + s (1 ))g (7) According to the Biot's theory [22], then we have a 2D stress field model from the above equations

Fluid Flow Modeling
In matrix, the fluid equation from conservation of mass is In fractures, it is The velocity from Darcy's law in matrix is In fractures, it is Furthermore, the storage model is Thus the pressure formulation in matrix is In fractures, it is

Fluid-solid Coupling Model
From the fluid pressure, the stress-strain equation is According to the Biot's theory [22], the fluid pore pressure is: In a matrix-fracture system, fractures and matrix have different flow characteristics.The storage coefficient in matrix is But in fractures, according to the data of experiments in [23], we know that a nonlinear relationship exists between permeability and the fluid pressure.Using the concept of a DFN model, fluid is full of the fracture space so the storage coefficient in fractures is The aperture of fractures is [24,25] The effective stress is Finally, the fracture flow is

Initial and Boundary Conditions
The solution of the governing equations for the above model requires a description of the initial state of the reservoir and appropriate far-field boundary conditions.The model can be depicted as in Fig. (1).

Procedure of Computing
The above model of full partial differential equations describing the reservoir behavior is solved using the finite element method, and the calculation process is described as in Fig. (2).( 1) The flow equation, including the matrix flow and fracture flow, is firstly solved to find the velocity and pressure fields; (2) then the stress-strain equation is solved to obtain the variation of stresses field; (3) new deformation of rock is next updated through solving the fluid-solid coupled model; (4) reservoir properties, including facture aperture and average permeability, can be calculated; (5) the flow model is solved to get the new pressure distribution according to new fracture permeability; (6) finally, the whole process of permeability variation can be realized.The parameters used for solution are listed in Table 1, and all values are fixed unless otherwise stated.

Variation of Fracture Aperture and Porosity
The geometry of the simulated fractures is built up in a model as in Fig. (3).In this model, four numbered fractures form a radicalized pattern with one cross-point and each of them is placed every 45 degrees.The rock with a surrounding rock pressure is drained from the left side to the right side so the variation of fracture permeability can be detected through numerical modeling.The model is calculated five times, and the calculation time is 360 days till the state stage reached.The inlet pressure is set as 12, 14, 16, 18, and 20 MPa, respectively.The outlet pressure is 10Mpa and the in-situ stress is 10MPa; other parameters are given in Table 1.After calculation, the pressure distribution is shown in Fig. (4), and the fracture aperture with different directions is displayed in Fig. (5).
From the calculations, the aperture of every fracture increases with the pressure rising on the same tendency, and the more the pressure rises, the more the fractures open, which can be seen from Fig. (5).However, the change of average aperture of each fracture is to an apparently varying     From the analysis of Fracture 1 and Fracture 2, this phenomenon is not hard to understand.The first part has one side (inlet of Fracture 3) to link with the matrix so the pressure consumes less in this fracture but more in matrix.The second part has one side to link with the connected fracture network, and the other side link with outlet so the pressure has no matrix to consume and the pressure gradient is naturally larger.Fracture 4 parallels to inlet and outlet and no pressure gradient exits so the aperture along this fracture increases on the same level as the pressure rises.
Through the above analysis, we know that not just pressure can influence the fracture aperture in pressure sensitive formation; the distribution of a fracture plays a significant role in this area as well.Depending on the fractures in pressure sensitive formation with a different direction or conductivity, the performance of aperture change can be wholly different.Because of the aperture linear correlation with fracture porosity [26], the direction and conductivity of a fracture is different so the change of fracture porosity is no longer in the same way.The variation of porosity of an anisotropic fracture can be investigated from a DFN model.Yet, due to the lack of precise fracture data and huge amount of computation, the DFN model cannot support a wide range of applications for large-scale reservoir simulation; therefore, a fracture network model from seismic interpretation and geo-statistical analysis still needs to be transformed to a dual-pore or dual-perm model for simulation.Some researchers used Eq.23 to establish the relationship between fracture permeability and porosity [27] k f = C f μ f (23) If we want to use the above model to simulate an anisotropic permeability reservoir and take this relationship between the permeability and porosity into consideration, then Eq. 23 can be written in a nine-point scheme of a 2D grid with three sets of fractures (as Fig. 6 shows) as follows: But if the formation is pressure-sensitive, Eq. 24 with a bulk fracture volume is no longer suitable.Because the previous studies have already illustrated that a different location of a fracture leads to different performance in an anisotropic reservoir.The technically available way is using a stressporosity coupling method, rather than a stress-aperture coupling model, to simulate a pressure-sensitive fractured reservoir.However, using the bulk fracture porosity for simulation cannot characterize this kind of performance.Therefore, it is necessary to introduce a new concept, anisotropic fracture porosity, for simulation of this type of reservoirs.The permeability in an anisotropic reservoir is in tensor form, while the fracture porosity in the same reservoir should be in tensor form as well.This means that every fracture in a different direction in Fig. (6) has its own porosity so the equation can be written as In recent years, porous media with anisotropic porosity have already become a hot area of research [28][29][30].As for the fractured porous media, the fracture porosity plays a more important role in seepage flow compared with matrix porosity, and even sometimes the anisotropy of the matrix porosity can be ignored in highly fractured media, e.g., in a dual-porosity and single-permeability model.Therefore, in a process of simulating a pressure sensitive fractured reservoir, using an anisotropic fracture porosity model would obtain more reliable results.

Fig. (6).
A schematic representation of fracture in formation of a nine-point scheme grid.

Effect of Production Process
During practical production of a really low permeable reservoir with fractures, some variation of fracture permeability can be detected through a well logging analysis, e.g., in a geothermal field in Dixie Valley, Nevada, US [31].It is evident that a production process has an effect on formation productivity so the effect of the development process needs to be discussed.Another model is used here as in Fig. (7).
The model in Fig. (7) is different from the early model: Fracture 4 is deleted so the permeability can be more anisotropic.Here we use different inlet and outlet pressures to simulate the two different production methods in a low permeable reservoir: inlet with 20MPa and outlet with 10MPa (20-10MPa drain) to simulate leading water injection; inlet with 10MPa and outlet with 0MPa (10-0MPa drain) for delayed water injection.The main distinction of these two methods is the formation pressure difference, which can be realized by adjusting the inlet pressure and outlet pressure of the model.All other parameters are constant, and the calculation results are shown in Fig. (8 and 9).
From the calculation results, we can clearly see that as the field pressure drops, the aperture of fractures decreases simultaneously.In addition, the two obvious distinctions can be detected: one difference between 10-0MPa drain and 20-10MPa drain is that all fractures have a linear pressure variation along the fracture, and the pressure gradient in the fractures is nearly the same as that in matrix, as Fig. (8) shows; the other difference is that the aperture of all parts of all fractures decreases to a nearly same value in the 10-0MPa drain process, as Fig. (9) shows.The reason is that as the pressure drops fractures closes more with less aperture and lose their conductivity, and permeability of the closed fractures approaches the permeability of matrix so the pressure gradient in both fractures and matrix are nearly the same.
It is also at the point where the anisotropy of permeability changes as the fracture closes and permeability declines.Here we use the method of fracture permeability in tensor form by Cartesian coordinates [32] to calculate the permeability, which is By ignorance of matrix permeability, Eq. 26 can be written in short form as follows: Then the principal permeability value and the principal permeability direction can be calculated, and the calculation results are given in Fig. (10).Fig. (10) shows that the principal permeability propagates some distance in the 20-10MPa drain process and the first principal propagates more compared with other slightly changed ones so the permeability of formation becomes more anisotropic, and the principal permeability directions rotates 8.4 degrees in the clockwise direction.Moreover, the first principal permeability is nearly twice as the secondary principal permeability in the 20-10MPa drain process, and, as a result, the formation presents a strong anisotropy.However, as the field pressure drops, both values of the first and second principal permeability decline to a nearly same value so the anisotropic permeability becomes approximately   isotropic permeability, and the principal permeability direction rotates back by 5.1 degrees.Therefore, the leading water injection raises the field pressure and opens the fracture more, causing more anisotropic permeability and changing the permeability direction.While the delayed water injection reduces the field pressure and closes the fracture, it results in more isotropic permeability and changes more permeability direction.From the above analysis, the leading water injection seems more suitable for fracture permeability since it is capable to increase the conductivity of fractures, but the oversized fractures may lead to earlier water break though.Hence it is a case by case issue to determine a reasonable field pressure.

Interaction Effect with in-situ Stress Field
In-situ stress has a major impact on reservoir development, especially in a pressure sensitive reservoir.If the formation is pressure sensitive, the fracture is stress-dependent as well, and what effect the in-situ stress has is necessary to be investigated.Here the model in Fig. (7) is still used, but the in-situ stress is conducted with 10MPa, 15MPa, and 20MPa, while keeping other parameters constant.The calculation results can be seen in Figs.(11 and 12).
As can be seen from Fig. (11), the more the in-situ stress is, the less deformation the fracture performs in the development process, and the smaller aperture can be detected in all fractures.Different from the above situation, the initial aperture of fractures also declines as the in-situ stress rises.A slight variation in the direction still exists but is hard to be observed as the stress rise; anisotropic permeability remains.It illustrates that the stress has an intensive impact on the value of permeability but little impact on anisotropy and the principal direction of permeability.Actually, the in-situ stress enhances with depth, and it means if the reservoir is deeper, the aperture of fractures becomes smaller and the reservoir should be more isotropic.

Fig. ( 1
Fig. (1).A schematic representation of initial and boundary conditions of dynamic flow in matrix-fracture media.

Fig. ( 5 ).
Fig. (5).Effect of different fracture direction on the extent of fracture aperture vs. pressure difference, (a) fracture 1 (b) fracture 2 (c) fracture 3 and (d) fracture 4.degree.As for Fracture 1, the aperture of the head end (the end of this fracture close to inlet) opens in the maximum level, and the tail end (the end of this fracture close to outlet) presents little difference with the initial condition.The drop of fracture aperture from the head end to the tail end keeps a near-linear relationship.It is because Fracture 1 links inlet with outlet directly so its pressure gradient is maximum among all fractures, which can be seen from Fig.(5).Fracture 2 presents the same linear variation trend, but to a lesser extent.It is because Fracture 2 is not connected to both inlet and outlet.Compared with the matrix, opening fractures have more conductivity so the pressure consumes more in matrix from inlet to the head end and from the tail end to outlet than in the fracture, as Fig.(5) shows.Different from the above two examples, Fracture 3 shows a nonlinear variation trend; it can be divided into two parts, where the first part from the head end to the cross-point has a relatively small pressure gradient and the second part from the crosspoint to the tail end has a relatively large pressure gradient.From the analysis of Fracture 1 and Fracture 2, this phenomenon is not hard to understand.The first part has one side (inlet of Fracture 3) to link with the matrix so the pressure consumes less in this fracture but more in matrix.The second part has one side to link with the connected fracture

Fig. ( 8
Fig. (8).A distribution map of pressure fields in (a) 2D and (b) 2D field elevation when inlet pressure is 10MP and outlet pressure is 0MPa.

Fig. ( 10
Fig. (10).A diagrammatic drawing of change of principle permeability direction in polar coordinates with different drained condition.

Fig. ( 12
Fig. (12).A diagrammatic drawing of change of principle permeability direction in polar coordinates with different in-situ stress field.

Fluid-solid coupling=
is an existing phenomenon in porous media flow, and what influence the stress field has for anisotropic fractures remains as an open research question.Here we simulate two models, one with the fractures as in Fig. (7) and the other with no fracture.The drain condition is 20-10MPa and in-situ stress 10MPa; no other parameters change.The outcomes from the calculations are represented in Fig. (13), from which we can observe that if a fracture is parallel to the drain pressure gradient, no obvious change can be seen in Fig. (13a).If the fracture and drain pressure gradient are unparalleled, the principal stress in the area among the fractures or near the fractures rotates, and this can be seen from Fig. (13b).It indicates a strong effect by distribution of fractures to the stress field.term in fractures (kg m -3 s -1 ) Q m = Mass source term in matrix (kg m -3 s -velocity in fractures (m s 1 ) u xx = Displacement in x direction (m) u yy = Displacement in y direction (m) Porosity of fractures (dimensionless)