Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture Permeability Variation: An Analytical Approach

Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture... Shale gas has recently gained significant attention as one of the most important unconventional gas resources. Shales are fine-grained rocks formed from the compaction of silt- and clay-sized particles and are characterised by their fissured texture and very low permeability. Gas exists in an adsorbed state on the surface of the organic content of the rock and is freely available within the primary and secondary porosity. Geomechanical studies have indicated that, depending on the clay content of the rock, shales can exhibit a brittle failure mechanism. Brittle failure leads to the reduced strength of the plastic zone around a wellbore, which can potentially result in wellbore instability problems. Desorption of gas during production can cause shrinkage of the organic content of the rock. This becomes more important when considering the use of shales for CO sequestration purposes, where C O adsorption-induced swelling can play an important role. These phenomena lead 2 2 to changes in the stress state within the rock mass, which then influence the permeability of the reservoir. Thus, rigorous simulation of material failure within coupled hydro-mechanical analyses is needed to achieve a more systematic and accurate representation of the wellbore. Despite numerous modelling efforts related to permeability, an adequate representation of the geomechanical behaviour of shale and its impact on permeability and gas production has not been achieved. In order to achieve this aim, novel coupled poro-elastoplastic analytical solutions are developed in this paper which take into account the sorption-induced swelling and the brittle failure mechanism. These models employ linear elasticity and a Mohr–Cou- lomb failure criterion in a plane-strain condition with boundary conditions corresponding to both open-hole and cased-hole completions. The post-failure brittle behaviour of the rock is defined using residual strength parameters and a non-associated flow rule. Swelling and shrinkage are considered to be elastic and are defined using a Langmuir-like curve, which is directly related to the reservoir pressure. The models are used to evaluate the stress distribution and the induced change in perme- ability within a reservoir. Results show that development of a plastic zone near the wellbore can significantly impact fracture permeability and gas production. The capabilities and limitations of the models are discussed and potential future develop- ments related to modelling of permeability in brittle shales under elastoplastic deformations are identified. Keywords Coupled geomechanics · Fractured reservoir · Shale gas · Permeability modelling · Poro-elasto-plasticity · Analytical solution List of symbols  In situ horizontal stress under initial reser- voir conditions Applied Mechanics Δ Change in mean horizontal stress from the ,  ,  Stress components in cylindrical coordinates rr  zz initial in situ conditions Mean horizontal stress (i.e. average of ,  ,  Strain components in cylindrical coordinates rr  zz ,  ) rr  e,m The elastic part of mechanically induced strain e,s The sorption-induced volumetric strain * Mohsen S. Masoudian (elastic only) [email protected] The plastic part of strain e,s 1  The in situ swelling strain under initial res- Nottingham Centre for Geomechanics, University ,0 ervoir conditions of Nottingham, University Park, Nottingham NG7 2RD, UK Vol.:(0123456789) 1 3 1566 M. S. Masoudian et al. Biot’s coefficient ZDCS Acronym for zero displacement condition E Elastic modulus on wellbore boundary and constant stress Poisson’s ratio condition on outer boundary c,  Cohesion and friction angle, respectively ZDZD Acronym for zero displacement condition on c ,  Post-failure cohesion and friction angle, both boundaries r r respectively f (r) A special function defined from distribu- , Y Strength parameters in Mohr–Coulomb tions of pore pressure and sorption-induced failure criterion defined using cohesion and strain friction angle F (r) Integral of pore pressure function over the , Y Post-failure strength parameters in Mohr– domain of the reservoir r r Coulomb failure criterion defined using F (r) Integral of sorption-induced strain function post-failure cohesion and friction angle over the domain of the reservoir Dilation angle F(r) A special function defined from F (r) and A function of dilation angle F (r) u Radial displacement at the elastic–plastic ep interface Geometry 1 Introduction r, , z Directions in cylindrical coordinates (the origin is the centre of the wellbore) Unconventional gas reservoirs are increasingly being con- R Radius of the disc-shaped reservoir sidered to provide a relatively clean energy source to meet R Radius of the wellbore burgeoning global demands. Unconventional gas resources, R Radius of plastic zone developed around the ep such as shale, represent around 40% of the remaining wellbore resource of technically recoverable natural gas; hence, they play an important role in the future global energy market Flow and Sorption Parameters P Pore pressure (fracture porosity) (McGlade et al. 2013). In addition, unconventional gas res- ervoirs are being considered as a host rock for CO geo- p Uniform pore pressure in the reservoir under 0 2 initial condition sequestration, which may also provide an enhanced gas recovery (EGR) technique. Shales are fine-grained rocks p Wellbore pressure (constant during gas production) formed from the compaction of silt- and clay-sized particles and are characterised by their fissured texture and very low Q Coefficient of pore pressure reduction due −18 −23 2 to gas production with the steady-state loga-permeability (10 –10  m ) (Alexander et al. 2011). Gas shales can be considered as a dual porosity rock consisting rithmic pressure distribution Q Dimensionless production rate of fractures and matrix blocks (Fig. 1), with a porosity of up to 15% (Wang and Reed 2009). Fracture sets can be subverti- k Permeability of the reservoir k Permeability of the reservoir under initial cal or parallel to bedding, may be filled with mineral gauges, have openings ranging in scale from micrometres to centi- conditions c Fracture compressibility metres, and spacing ranging in scale from a few centimetres to several metres (Gale et al. 2014). The matrix usually con- V Volume of adsorbed gas V , b Maximum monolayer sorption capacity and tains an organic content of 2–10%, while those with higher organic contents are usually considered too immature for Langmuir isotherm constant, respectively Maximum swelling strain production development (Alexander et al. 2011). The matrix contains porosity in the scale of µm to nm in size within Solution‑Dependent Definitions different organic (fibrous) and nonorganic (clastic) struc- A, B, K Constants defined in stress and displacement tures (Mehmani et al. 2013). Natural gas can be free within solutions both the fracture (primary porosity) and clastic structure of C , C , C Integration constants for stress solutions 1 2 3 the matrix (secondary porosity), or in an adsorbed state on CSCS Acronym for constant stress condition on the surface of the porous structure of organic content. The both boundaries adsorption of gas in organic geomaterials can be described CSZD Acronym for constant stress condition on using an isotherm, among which Langmuir’s is most widely wellbore boundary and zero displacement used (Gensterblum et al. 2015). For example, the gas content condition on outer boundary of five reservoirs in the USA, with in situ pressures between 2 and 28 MPa (600–8500 m deep), ranged from 1 to 10 m / 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1567 Fig. 1 a Fractures in shale core samples (Gale et al. 2014), b a subvertical calcite-cemented fracture (Soeder 1988), c simplified dual porosity concept for shale ton, of which 20–85% was adsorbed in the organic content of gas shales to changes in the stress level (e.g. Amann- (Curtis 2002; Hill and Nelson 2000). Hildenbrand et al. 2012; Chen et al. 2015a; Cho et al. 2013; Due to the multi-scale and multiple physical phenomena Li et al. 2016; Spencer 1989; Zhou et al. 2016). However, associated with shale gas production, the mechanism of gas these studies have only considered the elastic behaviour flow is a combination of viscous (Darcy and non-Darcy) of the rock. Plastic deformations within the reservoir rock flow (Huang et al. 2016), and gas-slippage (Klinkenberg) can have significant implications for production, injectiv - effect (Mehmani et al. 2013), as well as continuum (Fickian), ity, and stability of the wellbore. Studies have demonstrated free-molecule (Knudsen), and surface (Yuan et al. 2014) dif- the significance of both elastic and plastic deformations on fusions. On the other hand, the gas sorption (adsorption and/ wellbore producibility and/or injectivity as well as perme- or desorption) also plays an important role in the gas trans- ability prediction in other types of reservoirs (e.g. Cui et al. port mechanism in shale reservoirs, given that the amount of 2007; Han and Dusseault 2003; Masoudian et al. 2016b). adsorbed gas in the organic content is comparable with the In addition, the stability of wellbores is an important issue free gas. Despite the large number of studies on primary and in any oil and gas production project. Therefore, improved CO -EGR, many aspects of shale gas development are still models are needed to properly estimate the deformation and not well understood, mainly due to the complex interactions stress distributions around the wellbore. Experimental stud- of fluid flow with different thermal, chemical, and mechani- ies have also revealed the brittle mechanism of failure in cal phenomena. For instance, withdrawal of water and gas shale (e.g. Amann et al. 2011; Hull et al. 2015; Sone and from a subsurface formation results in significant reservoir Zoback 2013), which implies the need for improved models pressure alteration, which in turn leads to a redistribution of that consider this failure mechanism. The effect of gas and/ stress and strain within the shale formation. Gas and water or water sorption-induced swelling should also be taken into sorption-induced swelling/shrinkage (Lyu et al. 2015; Yuan account when developing permeability models. et al. 2014), thermal phenomena (e.g. during thermal stimu- The most common geomechanical assumption when pre- lation (Carroll et al. 2011)), and chemical reactions of water dicting permeability in gas reservoirs is of one-dimensional or CO (in case of CO -EGR) with shale (Carroll et al. 2011) elastic deformation (vertical deformation). However, this 2 2 can additionally impact the stress–strain relations. Redis- assumption cannot provide an accurate estimation, especially tribution of stress and strain within any reservoir can sub- near the wellbore where plastic deformations may occur. sequently lead to other changes, such as reservoir perme- On the other hand, using fully numerical hydro-mechanical ability alterations, subsidence/uplift of the ground surface, models leads to substantially more expensive computations. fault-reactivation mechanisms, and crack development in the This is especially important in shale gas reservoirs where caprock (Masoudian et al. 2016a). complex multiple porosity reservoir models are utilised. Prediction of permeability has received considerable This study aims to extend previous work (Masoudian and attention due to its significant effect on gas production from Hashemi 2016) by further developing brittle elastoplas- reservoirs. Numerous studies have related the permeability tic models for reservoirs coupled with sorption-induced 1 3 1568 M. S. Masoudian et al. swelling/shrinkage effects. In this paper, models are devel- approach is adopted, where the rock is considered as a con- oped with boundary conditions that simulate different well tinuum at the field scale, while the fractured nature of the completions and geological settings (e.g. open-hole and rock is captured by the effect of stress on fracture porosity. cased-hole completions). The models are then used to esti- Figure 2a illustrates the geometry of the model and the multi- mate the stress distribution within the reservoir, with which continuum approach adopted in this study. Figure 2b depicts the changes in fracture permeability within both elastic and the stress–strain relations within the brittle elastoplastic frame- plastic zones are predicted. This is an important aspect of the work used in this paper. Note that the mechanical behaviour work that has been largely neglected in other studies, mainly of brittle rock is idealised with a sudden drop of strength from due to lack of computational resources for field scale fully peak to residual value. Similar to perfectly plastic behaviour, coupled hydro-elastoplastic simulations. This work aims to there is no further change in stress beyond the yield point as highlight the role of plastic deformations in the prediction plastic strains continue to accumulate. The governing equa- of permeability and its impact on gas production. tions employed in this paper are explained below. Starting with the classical approach in plasticity, strains can e p be decomposed into elastic (  ) and plastic (  ) components, 2 Theoretical Framework and Governing where the elastic strain can be further divided into elastic e,m Equations mechanical strain, (  ), and elastic swelling (or shrinkage) e,s strain, (  ). Note that adsorption induces swelling while des- The reservoir is assumed to be a disc-shaped homogenous iso- orption leads to shrinkage of the matrix blocks. Assuming tropic continuum. A vertical wellbore is assumed to be drilled axial symmetry of the problem and considering a cylindrical at the centre of the disc, through which natural gas can be pro- system of coordinates (r ,  , z for radial, tangential and vertical duced. It should be noted that these simplifying assumptions directions, respectively), the strain compatibility equation is are employed in order to develop a fully analytical framework rr in this paper. Clearly, these assumptions cannot fully represent = (1) r r the heterogeneous nature of shale; however, they do not pre- vent evaluation of the effect of elastoplastic deformations on where r is the radial distance from the centre of a disc- stress and permeability distributions within reservoirs, which shaped reservoir, and subscripts rr and  represent the is the focus of this paper. Here, a simplified multi-continuum radial and tangential coordinates, respectively. In the light of Fig. 2 Schematics of a the geometry of the model, the multi-continuum approach, and the swelling/shrinkage effect on matrix size, and b the stress and strain relationships in the adopted elastoplastic framework 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1569 paper, a simplified pressure solution is sufficient. To arrive Biot’s definition of effective stress (  =  − P ), the equa- tion of equilibrium can be written as at a simplified pressure solution, the reservoir is initially assumed to be under a uniform pore pressure regime. Gas � � � rr  rr production from the wellbore with a fixed pressure ( p ) (2) w = − r r r leads to the development of a non-uniform pore pressure distribution, with pressure remaining equal to the initial where P is pore pressure,  is Biot’s coefficient, and  is pressure at the outer boundary only. Thus, a steady-state effective stress. Assuming that swelling strain is analogous solution of radial flow can be found in the following loga- to thermal strain, the elastic stress–strain constitutive rela- rithmic form (Cui et al. 2007) tions for an isotropic rock can be written as e,s P = p + Q ln (8) =  −   +  − (1 − 2)P + rr rr  zz E 3 e,s where p and R are the initial reservoir pressure and radius =  −   +  − (1 − 2)P + (3) 0 0 zz rr E 3 of the disc-shaped reservoir, respectively. The term Q is e,s defined as =  −   +  − 1 − 2 P + ( ) zz zz rr E 3 p − p w 0 Q = (9) ln R ∕R w 0 e,s where E ,  and  are the elastic modulus, Poisson’s ratio, where R and p are the radius and pressure of the wellbore, and volumetric swelling strain, respectively. The linear w w respectively. Mohr-Coulomb yield criterion can be written as The volume of adsorbed gas, V , can be estimated using � � =  + Y (4) rr Langmuir’s isotherm as a function of pressure ( P) where  and  are the tangential and radial effective rr V b P stresses, respectively. The parameters  and Y are defined V = (10) 1 + b P using cohesion, c , and friction angle,  , of the rock as where V is the maximum sorption capacity and b is the 1 + sin Langmuir isotherm constant. Initially, the matrix of the res- 1 − sin (5) ervoir rock contains a specified amount of adsorbed gas at 2c cos Y = the initial pressure ( p ) which changes as the reservoir pres- 1 − sin sure decreases due to production. As a result, the volumetric e,s e,s swelling strain changes from the initial state of  to  and ,0 Note that by substituting the residual cohesion, c , and is considered to have a linear relationship with adsorbed gas friction angle,  , into these equations, the residual param- volume, which gives eters  and Y can be defined for failed brittle rock. Consid- r r ering a non-associated flow rule, the plastic components of  b P e,s L (11) radial and tangential strains are related as (Park and Kim 1 + b P 2006) where  is the maximum swelling strain that can only take p L +  = 0 rr  place when the gas within the matrix is equal to the sorption e,s (6) u u e e capacity. Note that  is used to account for the in situ equi- +  =  +  ,0 rr r r librium state of the reservoir and can be estimated directly using the equation above, replacing the pressure with p . where u is the radial displacement, and  is a function of There are a number of fracture porosity/permeability dilation angle,  , as models developed for reservoirs where adsorption-induced 1 + sin  swelling or shrinkage is associated with gas injection or (7) production (e.g. Gensterblum et al. 2015; Pan and Connell 1 − sin 2012). However, these models adopt two main approaches There are a large number of studies that have investigated for porosity/permeability modelling: dependency on stress or modelling approaches for flow of gas (methane or carbon strain. Both of these approaches are closely related and may dioxide) in shale reservoirs, where the complex multi-porous lead to very similar results under certain conditions (Palmer nature of shale is carefully considered and different mass 2009), as they both consider a cubic law to link permeabil- transport mechanisms (i.e. viscous and diffusive) are taken ity to fracture porosity (consistent with a Kozeny–Carman into account. However, for the comparative analyses in this type relationship). In this paper, results are obtained with 1 3 1570 M. S. Masoudian et al. the most widely used approach where the permeability is where C and C are the integration constants that can be 1 2 related to the change in horizontal effective stress using an determined from the boundary conditions, and exponential equation as (Shi and Durucan 2005) (1 − 2) E e,s f (r) = P + (1 − ) 3(1 − ) (1 − 2) = exp − 3c Δ E (12) k h F(r) = F (r) + F (r) (1 − ) 3(1 − ) 2 2 r r where k is permeability, the subscript 0 denotes an initial F r = P d =− Q + P ( ) 4 2 � 0 value, c is fracture compressibility, and Δ is the change 2 R 1 + b p 1 + b P 0 L 0 L e,s s in horizontal effective stress, which can be defined as the F (r) =    d =  − exp − 2 Ei 2 × v L 2 b Q b Q b Q L L L mean of the change in the radial and tangential components (14) of effective stress. As illustrated in Fig.  2a, any element of where Ei(x) is the exponential integral function. rock (i.e. at a specific value of r ) can be conceptualised as The stress and strain relations within the plastic zone a hybrid of fractures and matrix blocks and any change in near the wellbore are governed by a Mohr–Coulomb failure stress state (due to pore pressure reduction and swelling/ envelope and a non-associated flow rule. Note that the brittle shrinkage effects) results in a change in the aperture size failure mechanism is implemented using residual strength of the fractures within that element. Equation 12 therefore parameters (i.e. cohesion and friction angle); perfectly plastic provides an efficient method for relating the response of frac- behaviour could be modelled by using peak strength values. tures and their permeability at fracture scale to the changes The general solution in the plastic zone can then be written as in effective stress at field scale. The use of mean horizontal effective stress implies a matchstick matrix model, which � 𝛾 −1 𝜎 = C r − K rr r<R is consistent with the plane-strain assumption. Closed-form ep � 𝛾 −1 solutions of permeability and porosity can then be obtained 𝜎  = Y + 𝛾 C r − K r r 3 𝜃𝜃 r<R ep by substituting the stress and strain equations for different 𝛽 +𝛾 r 𝛽 +𝛾 Eu R − r (1 + 𝜈 ) ep ep −𝛽 𝛽 loading and boundary conditions into the above equations. u = r R − C A r<R ep 3 ep E 1 + 𝜈 𝛽 + 𝛾 It should be noted, however, that experimental results have 𝛽 +1 𝛽 +1 shown that the relationship between post-failure permeabil- R − r ep 1 + 𝛽 − BY − AK + E F R − F (r) r 𝜀 ,𝛽 ep 𝜀 ,𝛽 ity and the stress–strain state can be much more complex 𝛽 + 1 3 than what these models offer (e.g. Carey et al. 2015). How - (15) ever, due to the lack of dedicated post-failure models, this where widely used simplified approach was adopted in this paper. Limitations of the model are discussed later in the paper. A = 1 −  −  +   −  − r r r B =  −  − Y − Q 3 Analytical Solutions for Elastic–Brittle– K = − 1 Plastic Rock (16) 1+ 1+ 1 + b p r 0 L 0 F (r) =  − exp −(1 + ) Following the approach employed by Masoudian and 1 +  b Q b Q L L Hashemi (2016), analytical solutions for stress and strain 1 + b P around a wellbore can be found by assuming two distinct × Ei (1 + ) b Q concentric zones; a plastic zone near the wellbore and an where Y and  are obtained with residual values of cohe- elastic zone towards the outer boundary. The derivation and r r sion and friction angle in Eq.  5, and u is the displace- development of these solutions are partially presented in ep ment at the elastic–plastic interface, which can be found “Appendix 1”. The interface of elastic and plastic zones is at considering continuity of displacement at the interface, i.e. a radial distance of R from the wellbore centre. The general ep by replacing the value of R in the elastic solution of dis- solution within the elastic zone ( r > R ) can be written as ep ep placement (Eq. 13) once the value R is found. To obtain ep F(r) − C R , the integration constants C , C , and C are found by 𝜎 = C + ep 1 2 3 rr 1 r>R ep substituting different sets of boundary conditions into the F(r) − C general solution, as presented in Table 1. Each set of bound- 𝜎 = C − + f (r) 𝜃𝜃 1 r>R ep ary conditions describes a particular geological and well (1 + 𝜈 ) completion condition. The constant stress at the wellbore u = (1 − 2𝜈 )r C − 𝜎 + 𝛼 p r>R 1 0 0 ep represents an uncased hole; hence, the total stress applied F(r) 2 e,s to the inner boundary is equal to the bottomhole pressure (13) + − + r𝜀 v,0 r r 3 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1571 ( p ); i.e. zero effective stress on wellbore wall. On the other outer boundary. As such, case CSCS represents constant hand, zero displacement on the wellbore wall represents an stresses (CS) on both the wellbore wall and outer boundary, ideally cased hole that does not allow the wellbore wall to case ZDCS represents zero displacement (ZD) of the well- converge. For the outer boundary, the constant stress rep- bore wall and constant stress on the outer boundary, case resents an ideal scenario where the reservoir is completely CSZD represents constant stress on the wellbore wall and separated from the adjacent formation. Zero displacement at zero displacement of the outer boundary, and case ZDZD the outer boundary represents a scenario where the reservoir represents zero displacement on both the wellbore wall and is large enough to prevent deformation of the rock mass at outer boundary. the far-field. For example, if the reservoir was cut by low- The interface between elastic and plastic zones is found angle faults at the outer boundary, the reservoir may be considering the continuity of radial stress at r = R and a ep best simulated by applying a constant stress, whereas a zero trial-and-error scheme. To achieve this, an equality is con- horizontal displacement boundary would be more suitable structed using the solution of radial stress within both elastic if the reservoir was cut by steep faults and was in contact and plastic zones, and then its root for R can be found. If this ep with mechanically strong geological strata. However, due equality does not have a root in the R , R domain, it means w 0 to the geological complexity of reservoirs, actual boundary that the rock mass does not include a plastic zone and hence a condition will fall somewhere between these two condi- fully elastic solution is used by replacing R by R in Eq. 13 ep w tions; hence, different boundary conditions should be tested. with the corresponding elastic integration constants presented The boundary condition acronyms consist of four letters, in “Appendix 2”. These formulations and the introduced meth- with the first pair describing the condition of the wellbore odology were implemented in a MATLAB code and used to wall, and the second pair describing the condition of the obtain the results discussed in the following section. Table 1 Integration constants for different sets of boundary conditions Boundary conditions (acronym) Relationships 2 2 R Y−f R −(−1)P R +(+1) F R −F R +R Constant stress at R ( ) ( ) ( ) ( ) (17) ep ep ep ep 0 0 0 C = 1 2 2 R (+1)−R (−1) Constant stress at R 0 ep (CSCS) C = C −  R + F R 2 1 0 0 C = 3  −1 2 2 R Y−f R −(−1)P R +(+1) F R −F R +R Zero displacement at R ( ) ( ) ( ) ( ) (18) ep ep ep ep 0 0 0 C = 1 2 R (+1)−R (−1) Constant stress at R ep 0 0 (ZDCS) C = C −  R + F R 2 1 0 0 +1 +1 Eu R −R +  1+ ep ep w C =   × R − BY − AK + E F R − F R ep 3 + r , ep , w r + r 1+ +1 3 A R −r ep Constant stress at R (19) w 1 C = × R Y − f R − ( − 1)P R 1 2 2 ep ep ep Zero displacement at R R (1−)−R (1+)(1−2) 0 ep 0 (CSZD) E e,s +( + 1) F R − F R + R  − (1 − 2)  − p ep 0 0 0 0 v,0 E e,s C = F R − R (1 − 2) C −  + p + 2 0 1 0 0 0 v,0 C = −1 3 r Zero displacement at R (20) w 1 C = × R Y − f R − ( − 1)P R 1 2 ep ep 2 ep Zero displacement at R R (1−)−R (1+)(1−2) 0 ep (ZDZD) E e,s +( + 1) F R − F R + R  − (1 − 2)  − P ep 0 0 0 v,0 2 e,s C = F R − R (1 − 2) C −  + p + 2 0 1 0 0 0 v,0 +1 +1 Eu R −R +  1+ ep ep w C = × R − BY − AK + E F R − F R ep 3 + r , ep , w + 1+ +1 3 A R −r r ep 1 3 1572 M. S. Masoudian et al. ZDZD scenarios (representing cased-hole completions), the 4 Results and Discussion radial stress significantly increases near the wellbore (mov - ing towards the wellbore from right to left in Fig. 3), since In order to examine the stress change and permeability evo- the rock is restricted on the wellbore wall. On the other hand, lution around a gas production well, a set of values were the tangential stress in these cases monotonically decreases chosen for the model input parameters, as listed in Table 2. as you approach the wellbore. This is because the support These values mainly represent gas shale reservoirs (suitable provided by the casing does not allow radial deformations for hydraulic fracturing due to their brittleness) and were and the release of radial stress within the rock mass. As a chosen to be within the reported ranges from the correspond- result, the change in radial and tangential effective stresses ing references. Studies have shown that hydro-mechanical (which is due to the interaction between pore pressure in the properties of near-wellbore shales are influenced by the pres- fractures and shrinkage of the matrix of the rock) is not large ence of natural and hydraulically generated fractures (Liang enough to satisfy the yield criterion (due to large minimum et al. 2014; Song et al. 2014), and the exposure of rock to principal stress levels) and the rock remains elastic, i.e. there drilling mud (Ma and Chen 2015). However, these effects is no plastic zone around the wellbore. It should also be are not considered in this paper. It should also be noted that noted that the case with constant stress at the outer bound- the sorption and swelling parameters correspond to those of ary (ZDCS) exhibits higher stress levels compared to the organic rich shales. case with a fixed outer boundary (ZDZD). This is because With production of gas from the reservoir, the pore pres- in ZDCS, the applied stress at the far-field compresses the sure decreases following Eq. 8, and the uniform initial pore rock mass towards the wellbore and results in a larger stress pressure is replaced by a radial distribution with maximum applied on the casing of the wellbore. pressure at the outer boundary (initial reservoir pressure,p ) In the CSCS and CSZD cases (representing open-hole and minimum value on the wellbore wall (wellbore pressure, completions), as the wellbore is radially deformed, the radial p ). During production, the gas is desorbed from the matrix stress decreases towards the wellbore pressure of 5 MPa in of the rock which leads to matrix shrinkage. The change the approach to the wellbore. The tangential stress in these in pore pressure combined with sorption-induced shrinkage cases, however, increases to some point and then suddenly leads to a redistribution of stress within the rock mass. The decreases in the approach to the wellbore. This is because distributions of radial and tangential total stresses for all sets the circular opening allows some deformation of the rock of boundary conditions are presented in Fig. 3. Note that the towards the wellbore which is also enhanced by the produc- effective stresses can be evaluated by subtracting pore pres- tion-induced reduction of pore pressure and the desorption- sure from the total stresses (considering  = 1 ); this will be induced shrinkage of the rock matrix. In these cases, the discussed later in the paper. It can be seen that in ZDCS and combination of stress levels satisfies the Mohr–Coulomb Table 2 Input parameters for Parameter Symbol Value Unit References the base simulation case Young’s modulus E 20 GPa Rutqvist et al. (2013) Peak cohesion c 6.0 MPa Ma and Chen, (2015) Peak friction angle  32.23 degrees Ma and Chen (2015) Residual cohesion c 4.39 MPa Ma and Chen (2015) Residual friction angle  25.17 degrees Ma and Chen (2015) Dilation angle  10 degrees Islam and Skalle (2013) Maximum sorption-induced swelling  0.134 % Chen et al. (2015b) Fracture compressibility c 0.02 1/MPa Chen et al. (2015a) Poisson’s ratio  0.22 – Ma and Chen (2015) Langmuir’s isotherm constant b 0.25 1/MPa Chen et al. (2015b) In situ horizontal stress  31.784 MPa Corresponding to 1500 m depth following Rutqvist et al. (2013) Initial reservoir pressure p 24.63 MPa Rutqvist et al. (2013) wellbore pressure p 5 MPa Wellbore radius R 0.216 m Ma and Chen (2015) Biot’s coefficient  1.0 – Chen et al. (2015a) Outer radius R 100.0 m – 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1573 Fig. 3 Profiles of stress within the reservoir: a radial stress, b tangen- Fig. 4 Distribution of mean horizontal stress for all cases: a total tial stress stress, b effective stress failure criterion and a plastic (failure) zone is developed discussed in detail later. The mean horizontal total stress around the wellbore; the maxima in tangential stress pro- decreases as it closely follows the reduction of pore pres- files coincide with the elastic–plastic interface. Note that sure from the far-field towards the wellbore (Fig.  4a). In the sharp drop of the tangential stress profile at the elas - cased-hole analyses (ZDCS and ZDZD), the rock remains tic–plastic interface is due to the brittle failure mechanism fully elastic and therefore the mean horizontal total stress of the rock. In other words, whereas the strength of the rock decreases monotonically from the far-field towards the immediately at the elastic–plastic interface is given by peak wellbore, with the minimum value occurring at the well- properties (peak cohesion and friction angle), the strength bore wall. For CSCS and CSZD, however, the mean total of the rock within the plastic zone is given by the residual stress decreases sharply at the elastic–plastic interface, fol- values of cohesion and friction angle. The radius of the plas- lowing the trends of radial and tangential stresses shown tic zone is larger in the CSCS case (0.324 m) compared to in Fig. 3. The mean horizontal effective stress is shown in that in CSZD (0.304 m) because of lower stress levels in the Fig. 4b, where a maximum can be observed at the interface latter case. of the elastic–plastic boundary in CSCS and CSZD cases. The mean horizontal total and effective stresses, defined It can also be seen that within the elastic zone, the solutions as  =  +  ∕2 and  =  − p , respectively (Cui with similar outer boundary conditions show similar trends m rr  m et  al. 2007), are depicted in Fig.  4. These definitions are (CSCS vs. ZDCS and CSZD vs ZDZD), even though their important in the prediction of permeability, which will be corresponding  and  values are different outside the rr 1 3 1574 M. S. Masoudian et al. open-hole conditions. For example, the radial stress at the elastic zone interface. In other words, for  and  within wall of the cased hole is much larger than that in the open the elastic zone, the effect of wellbore boundary condition hole, while the tangential stress decreases due to the lack diminishes and the outer boundary has the dominant effect of any radial displacement. With an open hole, on the other on the distribution of mean horizontal stress. hand, the reduction in radial stress is compensated by the Since the main purpose of this study is to understand the increase in the tangential stress. Thus, following the discus- effect of plasticity, the profiles of mean horizontal stress for sion of Figs. 4 and 5 shows that the effect of brittle-plastic the two open-hole cases need further discussion. These two behaviour of rock is only significant within the failure zone cases were also modelled using their corresponding fully near the wellbore. elastic solutions, as depicted in Fig. 5a, b. Figure 5a shows The estimated permeability is shown in Fig. 6a, which that the mean horizontal stress estimated using the elastic as expected inversely follows the mean horizontal effective solutions for CSCS and CSZD are almost identical to those stress. For elastoplastic cases, the minima are at the inter- calculated for the cased-hole cases ZDCS and ZDZD (shown face of the elastic and plastic zones. This corresponds to in Fig. 4a), respectively, although the distributions of radial the location where the rock is on the verge of failure and and tangential stresses are very different (Fig.  5c, d). This therefore can sustain the maximum stress levels; hence, the is discussed in more detail in “Appendix 3”, but it is mainly fractures are at their minimum aperture. On the other hand, because the distribution of radial and tangential stresses in in the elastic solutions, the mean effective stress increases elastic models depends on the condition of the wellbore monotonically towards the wellbore and therefore the rock wall; however, they change in a way that their summation near the wellbore is under the maximum compression, (  +  ) remains the same between the cased-hole and rr Fig. 5 Comparison of elastic and plastic solutions for open-hole cases (ZDCS and ZDZD): a mean total stress, b effective mean stress, c (CSCS and CSZD) and the comparison of elastic solutions of open- radial stress, d tangential stress hole cases (CSCS elastic and CSZD elastic) with the cased-hole cases 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1575 Fig. 6 Permeability of the reservoir with open-hole completion: a comparison between elastic and plastic solutions, b comparison between brit- tle–plastic and perfectly plastic solutions, c comparison of the solution with and without swelling effect which leads to the lowest values of permeability. When the and reduction of permeability. On the other hand, the reduc- constant radial stress condition is used for the outer bound- tion of pore pressure leads to desorption of gas from the ary (CSCS), the rock mass is pushed towards the wellbore matrix blocks of the rock. Consequently, the matrix blocks which in turn leads a higher degree of fracture closure and shrink and the width of fracture openings increases. There- lower permeability compared to the CSZD models (zero fore, depending on the poro-mechanical and swelling param- displacement at the outer boundary). The effect of brittle eters of the reservoir, the production of gas can lead to either failure on rock permeability is shown in Fig. 6b. Compared a reduction or an increase in permeability. In the example to perfectly plastic behaviour (with peak strength values), analysed here, the maximum volumetric strain induced by the brittle-plastic behaviour leads to a larger failure zone swelling/shrinkage was chosen to be 0.13% (an order of and lower permeability. This is because, within the plastic magnitude smaller than that of coal). Therefore, the effec- zone, the brittle model sustains lower stress levels than the tive stress change induced by pore pressure reduction had a peak strength perfectly plastic model, which leads to reduced more dominant effect. Figure  6c shows how the desorption- permeability in the brittle model. induced shrinkage can influence the fracture permeability It should be noted that permeability of fractured rock of a reservoir. The cases including swelling effects resulted (and the stress level within the rock mass) depends on the in larger permeability throughout the reservoir due to the opposing effects of two coupled phenomena: pore pressure shrinkage of matrix blocks. Another interesting observation change and sorption-induced swelling/shrinkage. First, the is that the radius of the plastic zone is influenced by the reduction of pressure in fractures leads to an increase in the swelling strain despite the fact that it is considered to be effective stress level, which results in closure of fractures an elastic component. This is because the shrinkage of the 1 3 1576 M. S. Masoudian et al. matrix block results in relaxation of stress; therefore, the towards lower residual values. This is because a stronger compressive stress applied to the rock mass reduces and the failed rock (higher residual strength) can sustain larger failed zone radius reduces. stresses and hence provides better support to the rock mass. Figure 7 depicts the radius of the plastic zone with respect The plastic behaviour of rock has been shown to affect to change in maximum sorption-induced swelling (  ) and the permeability of a reservoir only in the area very near residual strength parameters ( c ,  ). As sorption-induced the wellbore (i.e. within the plastic zone); however, it may r r volumetric change increases from zero to 1%, the matrix have a significant impact on the production flowrate. In order blocks shrink and hence the stress level decreases within to estimate the effect of near-wellbore permeability on pro- the rock mass. As a result, the radius of the plastic zone duction rate, a simplified approach is adopted in which the decreases. No plastic zone developed when using a maxi- equivalent permeability of the whole reservoir is defined by mum swelling strain larger than 1.09% in these examples. considering parallel ring-shaped layers having different per - In other words, the shrinkage of the matrix can help the meability values. In other words, by considering equal radial mechanical integrity of the wellbore. Increasing the values flowrate in each ring-shaped element of the reservoir, the of residual cohesion and friction angle also leads to a smaller equivalent permeability can be found. Hence, a dimension- plastic zone; the effect is proportionally greater at lower val - less production flowrate (somewhat similar to productivity ues of residual strength, i.e. increasing slope of the curves index) can be defined, simply as the ratio of production rate of the reservoir with a spatially variable permeability to that of a reservoir with constant initial permeability, which gives k dr Q = ln R − R (21) D 0 w rk(r) This parameter can be used to illustrate how changes in permeability influence the rate of gas production from a res- ervoir. As shown Fig. 8a, an increase in the maximum swell- ing strain leads to a reduction in the production rate. This is because, although the shrinkage of the shale matrix generally increases the permeability, it restricts the development of the plastic zone (the zone with significantly higher permeability), as demonstrated in Fig. 7, and results in a reduced equivalent permeability of the reservoir. Figure 8b shows how the small plastic zone around the wellbore can significantly affect the production rate. First, note that having a plastic zone devel- oped in elastoplastic solutions with enhanced permeability can lead to higher production rates compared to their elastic counterparts. The production rate also increases with lower residual strength parameters, due to enlargement of the plastic zone around the wellbore, as previously shown in Fig. 7. Such a significant effect of the plastic zone on the equivalent perme- ability of the whole reservoir and its production rate is due to the fact that the wellbore and its surrounding zone is the bot- tleneck of the production. This has also been addressed in the context of formation damage where the adverse effect of drill- ing mud on permeability in the vicinity of the wellbore is stud- ied. Thus, development of a plastic zone with higher perme- ability around the wellbore can have significant implications on production of gas from shale reservoirs. This highlights the fact that elastic solutions are not able to properly predict reservoir permeability and gas production rates. On the other hand, the use of casings for wellbore stability can significantly Fig. 7 Radius of the elastic–plastic interface with respect to change in a maximum sorption-induced swelling, b ratio of residual strength reduce the producibility of wells. Obviously, the mechanical parameters 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1577 Fig. 9 Assessment of the developed methodology for prediction of Fig. 8 Dimensionless production rate with respect to change in: a permeability: a exponential permeability model for different shales maximum sorption-induced swelling, b ratio of residual strength (data from Bustin et al. 2008; Dong et al. 2010; Gutierrez et al. 2000; parameters Pathi 2008), b field data of permeability around a tunnel and model predictions integrity of wellbores is of great importance; however these the permeability for many shales at small scale (core scale), results suggest that the development of novel wellbore casing confirming the validity of the exponential relationship for solutions (e.g. casings that allow a plastic zone to develop but fractured shale reservoirs. Unfortunately, field data present- prevent wellbore caving) could provide significant benefits to ing the effect of elastoplastic deformations on permeability the producibility of a wellbore. of shale reservoirs are not available. However, Souley et al. (2001) reported in situ measurements of permeability around a tunnel excavated in granite at Atomic Energy Canada Lim- 5 Limitations, Verification, and Future ited (AECL) Underground Research Laboratory in Canada Development of the Model to address issues related to the disposal of nuclear wastes. There are obviously some important differences between this As a measure of verification of the model, results are com- case and a shale gas reservoir, but these data are still valu- pared against experimental and field data in Fig.  9. First, the able as they provide insight into how plastic deformations exponential permeability model used in this study (i.e. Eq. 12) can influence permeability around a circular hole. In order is fitted to experimental datasets available in the literature to get a good match to in situ measurements, the proposed (usually from triaxial isotropic compression tests). Figure 9a model was used with no sorption-induced swelling, a fracture −1 shows that the use of Eq. 12 provides a good estimation of compressibility of 0.24 MPa , and a residual strength ratio of 1 3 1578 M. S. Masoudian et al. 0.83. Note that all other parameters were the same as those in well understood and cannot be properly predicted, especially Rutqvist et al. (2009) and Souley et al. (2001), where numeri- for highly anisotropic rocks with sparse organic contents such cal models were employed to replicate these measurements. It as shales (Kumar et al. 2016). Such studies can be even more should be noted that the estimated fracture compressibility in important when considering the CO injection for sequestra- this paper is within the range of values reported by Rutqvist tion and production enhancement purposes. The effect of shear et al. (2009). Although this case is simplified here and many displacement on fracture permeability has been discussed in a processes are neglected, the good agreement between the limited number of studies (Gutierrez et al. 2000; Rutqvist 2015); model and field data provides a verification of the accuracy this is another important topic that requires further investigation. of the proposed model. It can also be seen that considering the brittle behaviour of the rock is imperative to obtain better 6 Conclusions predictions of permeability. Several developments of the model could be considered to In this paper, a series of analytical solutions were presented in address the current limitations. Improvements could be made a coupled poro-elastoplastic framework, taking into account the by developing incremental forms of the solutions such that sorption-induced swelling/shrinkage of the rock matrix and the the stress–strain relations can be implemented within a time- brittle failure mechanisms of the rock mass. These solutions were dependent reservoir modelling framework (e.g. commercial used to model the complex inter-relations between different phe- reservoir simulators). However, care should be taken as the nomena during gas production from a gas shale reservoir. The model presented in this paper ignores the fact that severe results showed that the rock immediately adjacent to the well- softening behaviour may lead to localised failure, wherein bore may undergo plastic deformations, and depending on the the assumption of axisymmetric deformation becomes residual strength of the rock, this plastic zone may be enlarged, invalid (Lubliner 1990). Studies on the localisation effect while the increasing sorption-induced matrix shrinkage leads to on ground reaction curves for circular tunnels (e.g. Alonso smaller plastic zones. The results also showed that the perme- et al. 2003; Varas et al. 2005) found that results of radially ability of shale within the plastic zone can increase significantly symmetric solutions were valid from a practical engineering compared to its initial value, and while the developed plastic zone perspective. Thus, similar studies should be undertaken to is generally very small, its implications on producibility of the evaluate the extent to which the solutions developed in this wellbore can be significant. The advantages and limitations of paper are able to represent the actual changes in permeability the developed model were critically discussed. The model was of shale reservoirs. In addition, unconventional reservoirs verified using data relating to the permeability of a damaged zone such as gas shales are usually stimulated with hydraulic frac- around a tunnel. In the light of this discussion, and considering turing, which can change the mechanical behaviour of the relevant findings from the literature, a number of directions for rock in the vicinity of the wellbore (neglected in this paper), future research were proposed to address the lack of understand- which could be considered for future model developments. ing on the effect of elastoplastic deformations on permeability of Studies have proposed various models for permeability using shales. In addition, suggestions were made for further improve- different forms of the exponential relationship, but these have ment in the models by implementing their incremental forms in mainly concerned with the effect of isotropic effective stress on an advanced reservoir model and by considering more complex the permeability of shale and a correction of fracture compress- hydro-mechanical boundary conditions. These developments ibility (e.g. variation of c with effective stress). The perme - could provide more accurate predictions within the computa- ability model used in this paper assumes that the total sorption tionally efficient analytical framework presented in the paper. strain of the unconfined rock matrix contributes to permeability variation since the sorption strain incorporated in this model was Open Access This article is distributed under the terms of the Creative based on measurements made on an unconfined matrix. Zang Commons Attribution 4.0 International License (http ://crea tive comm ons.org/licenses /b y/4.0/), which permits unrestricted use, distribution, et al. (2015) introduced the concepts of ‘internal swelling’ and and reproduction in any medium, provided you give appropriate credit ‘internal swelling ratio’ for coal samples to evaluate the effect of to the original author(s) and the source, provide a link to the Creative sorption-induced swelling on permeability under confined con- Commons license, and indicate if changes were made. ditions. Such an approach could benefit future developments of permeability models. In addition, recent studies showed that the Appendix 1: Derivation of Stress behaviour of shales can be much more complex, especially when and Displacement Solutions considering shale reservoirs for C O -EGR. Due to the signifi- cance of micro-cracks in shale reservoirs, these findings imply The derivation of the solution of the elastic zone is presented an immediate need for better understanding of the mechanisms in the literature (e.g. Masoudian et al. 2016b), and only the of compression-induced fracture initiation and growth in shales. solution of the plastic zone is presented here. Substituting the Although studies related to hydraulic fracturing cover this topic post-failure yield criterion (an equation in the same form as to some extent, the relevant micro-scale phenomena are not 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1579 Eq. 4 with residual cohesion and friction angle replacing the radial stress will be postulated, which means that both elas- peak values) and the steady-state pressure distribution (Eq. 8) tic and plastic solutions must be equal at r = R . In other ep into the equation of equilibrium (Eq. 2) gives words, the elastic solution needs to satisfy the peak yield criterion (Eq. 4), because the rock at r = R is at its peak � � ep Y Q rr rr r strength. Therefore, substituting Eq. 13 into Eq. 4 gives a (22) + 1 −  = + r r r r solution for stresses r = R ep This is a first-order non-homogenous linear differential = 2C +  −  P rr 1 r=R r=R equation, which can be solved for  , using the integrating ep ep 1 +  1 − rr factor method. The value of  can be found by substituting e,s (27) +   − Y the solution into the post-failure yield criterion. The general r=R 3(1 − ) ep solution for  can be written as rr The effective radial stress at the elastic–plastic interface can Y C r 3 also be found from Eq. 23 and, by including the pore pressure �  −1 = − + = C r − K (23) rr 1− 1 −  1 −  r effect, the total radial stress can be written as r r where C is the integration constant that can be found using r  −1 3 r = − + C R +  P (28) rr 3 ep r=R r=R ep ep the boundary condition applied to the wellbore wall. For 1 −  1 − r r open-hole completions, the bottomhole pressure is equal to The total stress at r = R must be unique (i.e. radial stress ep the total stress applied to the boundary and hence the effec- continuity) and therefore an equality can be written as tive stress is zero. This gives 1  E e,s 2C +  −  P +  − Y 1 r=R C = ep v 3 r=R (24) 1 +  1 −  3(1 − ) −1 ep r  −1 = − + C R +  P 3 ep r=R ep When the wellbore is cased, the displacement at the well- (29) 1 −  1 − r r bore wall is assumed to be zero and thus the solution of Ideally, this equality is solved analytically to give a closed- displacement is needed. The first-order differential equation form expression of R , but this can be difficult, if not impossi- ep of displacement in Eq. 6 can be solved to give ble. Thus, a trial-and-error scheme can be used to find the root. The root of R is then admissible only if it is within the model ep ep ep 1 u = u −  g()d (25) ep r r Note that u is the displacement at the elastic–plastic ep interface and can be found assuming continuity of displace- ment; its value is equal to that found by the elastic solution at r = R . Following the approach suggested by Park and Kim ep (2006) and employed by Masoudian and Hashemi (2016), g(r) can be found, in the light of the first relationship in Eq. 6, as a function of elastic parts of radial and tangential strain by replacing the solution of stress in the plastic zone into elastic stress–strain relationships. Combining these equations and solving the integral term gives r + Eu R − r (1 + ) ep ep u = r R − C A ep 3 E 1 +   + +1 +1 R − r ep 1 + − BY − AK + E F R − F (r) r , ep , + 1 3 (26) Therefore, for a cased-hole completion, C can be found by solving this equation for u = 0 . It is evident from these equations that the solution needs the elastic–plastic interface ( R ) to be known. In order to find R , the continuity of Fig. 10 Verification of CSCS case for: a stress solution, b displace- ep ep ment solution 1 3 1580 M. S. Masoudian et al. Appendix 3: Comparative Analysis of  domain R , R ; i.e. if no root is found or it is smaller than the w 0 for Different Cases wellbore radius or larger than the reservoir radius, the reservoir remains elastic and no failure zone is developed around the It is useful to note how different boundary conditions on the wellbore. It should also be noted that the analytical solutions of wellbore affect the value of  within the elastic zone. As stress and displacement were compared to numerical solutions shown in the paper, CSCS and CSZD cases result in very simi- to verify their performance with perfectly plastic behaviour. lar profiles for  in the elastic zone as ZDCS and ZDZD, Such a verification for the CSCS case is depicted in Fig.  10. respectively. To better explain and justify this, a comparison of all the terms involved in  can be carried out.  is defined Appendix 2: Integration Constants for Elastic m m as the mean of radial and tangential stress as follows: Cases +  f (r) rr See Table 3.  = = C + (34) m 1 2 2 Table 3 Integration constants Boundary conditions (acronym) Relationships for different sets of boundary 2 2 conditions in elastic cases R −p R +F R −F R Constant stress at R ( ) ( ) (30) 0 0 w w w 0 C = 1 2 R −R Constant stress at R 0 w (CSCS) C = C −  R + F R 2 1 0 0 2 2 e,s p R −ER  ∕3−F R +F R ( ) ( ) Zero displacement at R 0 0 w (31) 0 w v,0 C =  − p + 1 0 0 2 (1−2)R +R Constant stress at R w (ZDCS) C = C −  R + F R 2 1 0 0 2 2 e,s 2 p R −ER  ∕3+F R −F R +(1−2)R  −p ( ) ( ) ( ) Constant stress at R w w 0 w 0 0 (32) 0 v,0 0 C = 1 2 (1−2)R +R Zero displacement at R w (CSZD) C = C − p R + F R 2 1 w w e,s F R −F R E Zero displacement at R ( ) ( ) (33) 0 w v,0 C =  − p + − 1 0 0 2 2 (1−2)(R −R ) 1−2 Zero displacement at R 0 w (ZDZD) E e,s C = F R − R (1 − 2) C −  + p + 2 0 1 0 0 0 v,0 Table 4 Integration constants in solution of mean horizontal stress within the elastic zone Boundary conditions (acronym) Relationships F R −F R ( ep) ( 0) CSCS and ZDCS (35) Y−f R −(−1)P R +(+1) + ( ep) ( ep) 0 F(R ) C = ≅  − 1 0 2 (+1)−(−1) R F(R )−F(R ) w 0 CSCS elastic (36) −p + 0 w R F R ( ) 0 0 C = ≅  − 1 0 1− R e,s F(R )−F(R ) 0 w ZDCS elastic (37) p −E ∕3− 0 2 v,0 R F R ( ) 0 0 C =  − p + ≅  − 1 0 0 0 2 (1−2)+1 R CSZD and ZDZD      (38) C =  Y − f R − ( − 1)P R + ( + 1) 1 ep ep (1−)−(1+)(1−2) F R −F R F R ( ep) ( 0) E e,s 1 ( ) E e,s +  − (1 − 2)  − p ≅  − p + − 2 0 0 0 0 2 R 3 v,0 1−2 R 3 v,0 0 0 F(R )−F(R ) e,s 0 w CSZD elastic   (39) p −E ∕3+ +(1−2)( −p ) w 0 0 v,0 2 R F R ( ) 1 0 E e,s C = ≅  − p + − 1 0 0 2 v,0 (1−2)+ 1−2 R 3 e,s F(R )−F(R ) F(R ) ZDZD elastic 0 w v,0 1 0 E e,s (40) C =  − p + − ≅  − p + − 1 0 0 2 0 0 2 v,0 (1−2)R (1−) 1−2 1−2 R 3 0 0 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1581 Dong J-J, Hsu J-Y, Wu W-J, Shimamoto T, Hung J-H, Yeh E-C, Wu Since f (r) is the same for all boundary conditions, the Y-H, Sone H (2010) Stress-dependence of the permeability and differences must relate to C . On the other hand, the plastic porosity of sandstone and shale from TCDP Hole-A. Int J Rock zone is always located very close to the wellbore. Thus, both Mech Min Sci 47(7):1141–1157. http s://doi.org/10.1016 /j.ijrm R and R can be considered very small compared to R . ms.2010 .06.019 w ep 0 2 2 2 2 −6 Gale JFW, Laubach SE, Olson JE, Eichhubl P, Fall A (2014) Natural Consequently, a ratio given by  = R ∕R ≈ R ∕R ≈ 10 w 0 ep 0 fractures in shale: a review and new observations. AAPG Bull can be implemented in all equations for C and neglected. 98(11):2165–2216. http s://doi.org/10.1306 /0812 1413 151 Gensterblum Y, Ghanizadeh A, Cuss RJ, Amann-Hildenbrand A, Furthermore, F R ≪ F R and can also be neglected. ep 0 Krooss BM, Clarkson CR, Harrington JF, Zoback MD (2015) Gas Considering these aspects for all four boundary condition transport and storage capacity in shale gas reservoirs—a review. cases results in the integration constants presented in Part A: transport processes. J Unconv Oil Gas Resour 12:87–122. Table 4. http s://doi.org/10.1016 /j.juog r.2015 .08.001 Gutierrez M, Øino LE, Nygård R (2000) Stress-dependent permeability of a de-mineralised fracture in shale. Mar Petrol Geol 17(8):895– 907. http s://doi.org/10.1016 /S026 4-8172 (00)0002 7-1 References Han G, Dusseault MB (2003) Description of fluid flow around a well- bore with stress-dependent porosity and permeability. J Petrol Sci Hull KL, Abousleiman, YN, Han Y, Al-Muntasheri GA, Hosemann P, Eng 40(1–2):1–16. http s://doi.org/10.1016 /S092 0-4105 (03)0004 Parker SS, Howard CB (2015) New insights on the mechanical 7-0 characterization of Kerogen-rich shale, KRS. In: Abu Dhabi inter- Hill DG, Nelson CR (2000) Gas productive fractured shales: an over- national petroleum exhibition and conference. Society of Petro- view and update. Gas TIPS 6(2):4–13 leum Engineers, Abu Dhabi, UAE, p 25. http://dx.doi.or g/10.2118 Huang X, Bandilla KW, Celia MA (2016) Multi-physics pore-network /1776 28-MS modeling of two-phase shale matrix flows. Transp Porous Media Alexander T, Baihly J, Boyer C, Clark B, Waters G, Jochen V, Le 111(1):123–141. http s://doi.org/10.1007 /s112 42-015-0584 -8 Calvez J, Lewis R, Miller CK, Thaeler J, Toelle BE (2011) Shale Islam MA, Skalle P (2013) Experimentally evaluating shale dilation gas revolution. Oilfield Rev 23(3):40–57 behavior. In: Chatellier J-Y, Jarvie DM (eds) Critical assessment Alonso E, Alejano LR, Varas F, Fdez-Manin G, Carranza-Torres of shale resource plays. American Association of Petroleum Geol- C (2003) Ground response curves for rock masses exhibiting ogists, Tulsa, p 205 strain-softening behaviour. Int J Numer Anal Meth Geomech Kumar H, Elsworth D, Mathews JP, Marone C (2016) Permeabil- 27(13):1153–1185. http s://doi.org/10.1002 /nag.315 ity evolution in sorbing media: analogies between organic-rich Amann F, Button EA, Evans KF, Gischig VS, Blümel M (2011) Experi- shale and coal. Geofluids 16(1):43–55. http s://doi.org/10.1111 mental study of the brittle behavior of clay shale in rapid uncon- /gfl.1213 5 fined compression. Rock Mech Rock Eng 44(4):415–430. https:// Li M, Yin G, Xu J, Cao J, Song Z (2016) Permeability evolution of doi.org/10.1007 /s006 03-011-0156 -3 shale under anisotropic true triaxial stress conditions. Int J Coal Amann-Hildenbrand A, Ghanizadeh A, Krooss BM (2012) Trans- Geol 165:142–148. http s://doi.org/10.1016 /j.coal .2016 .08.017 port properties of unconventional gas systems. Mar Petrol Geol Liang C, Chen M, Jin Y, Lu Y (2014) Wellbore stability model for 31(1):90–99. http s://doi.org/10.1016 /j.marp etge o.2011 .11.009 shale gas reservoir considering the coupling of multi-weakness Bustin RM, Bustin AMM, Cui A, Ross D, Pathi VM (2008) Impact of planes and porous flow. J Nat Gas Sci Eng 21:364–378. https:// shale properties on pore structure and storage characteristics. Soc doi.org/10.1016 /j.jngs e.2014 .08.025 Petrol Eng. http s://doi.org/10.2118 /1198 92-MS Lubliner J (1990) Plasticity theory. Macmillan, New York, p 50. Carey JW, Lei Z, Rougier E, Mori H, Viswanathan H (2015) Fract-per- ISBN 0486462900 meability behavior of shale. J Unconv Oil Gas Resour 11:27–43. Lyu Q, Ranjith PG, Long X, Kang Y, Huang M (2015) A review http s://doi.org/10.1016 /j.juog r.2015 .04.003 of shale swelling by water adsorption. J Nat Gas Sci Eng Part Carroll SA, McNab WW, Torres SC (2011) Experimental study of 27(3):1421–1431. http s://doi.org/10.1016 /j.jngs e.2015 .10.004 cement—sandstone/shale—brine—CO interactions. Geochem Ma T, Chen P (2015) A wellbore stability analysis model with chem- Trans 12(1):9. http s://doi.org/10.1186 /1467 -4866 -12-9 ical-mechanical coupling for shale gas reservoirs. J Nat Gas Chen D, Pan Z, Ye Z (2015a) Dependence of gas shale fracture perme- Sci Eng 26:72–98. https://doi.or g/10.1016/j.jngs e.2015 .05.028 ability on ee ff ctive stress and reservoir pressure: model match and Masoudian MS, Hashemi MA (2016) Analytical solution of a cir- insights. Fuel 139:383–392. http s://doi.org/10.1016 /j.fuel .2014 . cular opening in an axisymmetric elastic-brittle-plastic swell- 09.018 ing rock. J Nat Gas Sci Eng Part A 35:483–496. http s://doi. Chen T, Feng X-T, Pan Z (2015b) Experimental study of swelling of org/10.1016 /j.jngs e.2016 .08.076 organic rich shale in methane. Int J Coal Geol 150–151:64–73. Masoudian MS, Airey DW, El-Zein A (2016a) The role of coal seam http s://doi.org/10.1016 /j.coal .2015 .08.001 properties on coupled processes during C O sequestration: a Cho Y, Ozkan E, Apaydin OG (2013) Pressure-dependent natural- parametric study. Greenh Gases Sci Technol (in Press). http :// fracture permeability in shale and its effect on shale-gas well dx.doi.org/10.1002 /ghg.1575 production. SPE Reserv Eval Eng 16(2):216–228. h t t p s : / / d o i . Masoudian MS, El-Zein A, Airey DW (2016b) Modelling stress and org/10.2118 /1598 01-PA strain in coal seams during CO injection incorporating the Cui X, Bustin RM, Chikatamarla L (2007) Adsorption-induced coal rock–fluid interactions. Comput Geotech 76:51–60. https://doi. swelling and stress: Implications for methane production and acid org/10.1016 /j.comp geo.2016 .02.010 gas sequestration into coal seams. J Geophys Res Solid Earth McGlade C, Speirs J, Sorrell S (2013) Unconventional gas—a review 112:B10. http s://doi.org/10.1029 /2004 jb00 3482 of regional and global resource estimates. Energy 55:571–584. Curtis JB (2002) Fractured shale-gas systems. AAPG Bull http s://doi.org/10.1016 /j.ener gy.2013 .01.048 86(11):1921–1938. http s://doi.or g/10.1306 /61EE DDBE -173E Mehmani A, Prodanović M, Javadpour F (2013) Multiscale, mul- -11D7 -8645 0001 02C1 865D tiphysics network modeling of shale matrix gas flows. Transp 1 3 1582 M. S. Masoudian et al. Porous Media 99(2):377–390. http s://doi.or g/10.1007 /s112 Song F, Warpinski N, Toksöz M (2014) Full-waveform based micro- 42-013-0191 -5 seismic source mechanism studies in the Barnett Shale: link- Palmer I (2009) Permeability changes in coal: analytical modeling. ing microseismicity to reservoir geomechanics. Geophysics Int J Coal Geol 77(1–2):119–126. https://doi.or g/10.1016/j.coal 79(2):KS13–KS30. http s://doi.org/10.1190 /geo2 013-0094 .1 .2008 .09.006 Souley M, Homand F, Pepa S, Hoxha D (2001) Damage-induced Pan Z, Connell LD (2012) Modelling permeability for coal reservoirs: permeability changes in granite: a case example at the URL in a review of analytical models and testing data. Int J Coal Geol Canada. Int J Rock Mech Min Sci 38(2):297–310. http s://doi. 92:1–44. http s://doi.org/10.1016 /j.coal .2011 .12.009org/10.1016 /S136 5-1609 (01)0000 2-8 Park K-H, Kim Y-J (2006) Analytical solution for a circular open- Spencer CW (1989) Review of characteristics of low-permeability gas ing in an elastic–brittle–plastic rock. Int J Rock Mech Min Sci reservoirs in Western United States. AAPG Bull 73(5):613–629 43(4):616–622. http s://doi.org/10.1016 /j.ijrm ms.2005 .11.004 Varas F, Alonso E, Alejano LR, Manín GF (2005) Study of bifurca- Pathi VSM (2008) Factors affecting the permeability of gas shales. PhD tion in the problem of unloading a circular excavation in a strain- Thesis, The University of British Columbia, Vancouver, Canada softening material. Tunn Undergr Space Technol 20(4):311–322. Rutqvist J (2015) Fractured rock stress-permeability relationships from http s://doi.org/10.1016 /j.tust .2004 .12.003 in situ data and effects of temperature and chemical-mechanical Wang FP, Reed RM (2009) Pore networks and fluid flow in gas shales. couplings. Geofluids 15(1–2):48–66. http s://doi.or g/10.1111 / In: SPE annual technical conference and exhibition. Society of gfl.1208 9 Petroleum Engineers, New Orleans, p 3. http://dx.doi.or g/10.2118 Rutqvist J, Börgesson L, Chijimatsu M, Hernelind J, Jing L, Kobayashi /1242 53-MS A, Nguyen S (2009) Modeling of damage, permeability changes Yuan W, Pan Z, Li X, Yang Y, Zhao C, Connell LD, Li S, He J (2014) and pressure responses during excavation of the TSX tunnel in Experimental study and modelling of methane adsorption and dif- granitic rock at URL, Canada. Environ Geol 57(6):1263–1274. fusion in shale. Fuel Part A 117:509–519. https://doi.or g/10.1016 http s://doi.org/10.1007 /s002 54-008-1515 -6/j.fuel .2013 .09.046 Rutqvist J, Rinaldi AP, Cappa F, Moridis GJ (2013) Modeling of fault Zang J, Wang K, Zhao Y (2015) Evaluation of gas sorption-induced reactivation and induced seismicity during hydraulic fracturing internal swelling in coal. Fuel 143(Supplement C):165–172. http of shale-gas reservoirs. J Petrol Sci Eng 107:31–44. http s://doi.s://doi.org/10.1016 /j.fuel .2014 .11.007 org/10.1016 /j.petr ol.2013 .04.023 Zhou T, Zhang S, Feng Y, Shuai Y, Zou Y, Li N (2016) Experimental Shi JQ, Durucan S (2005) A model for changes in coalbed permeability study of permeability characteristics for the cemented natural frac- during primary and enhanced methane recovery. SPE Reserv Eval tures of the shale gas formation. J Nat Gas Sci Eng 29:345–354. Eng. http s://doi.org/10.2118 /8723 0-PAhttp s://doi.org/10.1016 /j.jngs e.2016 .01.005 Soeder DJ (1988) Porosity and permeability of Eastern Devonian gas shale. SPE Form Eval 3(1):116–124. https://doi.or g/10.2118/1521 3-PA Sone H, Zoback MD (2013) Mechanical properties of shale-gas reser- voir rocks—part 2: ductile creep, brittle strength, and their relation to the elastic modulus. Geophysics 78(5):D393–D402. https ://doi. org/10.1190 /geo2 013-0051 .1 1 3 http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Rock Mechanics and Rock Engineering Springer Journals

Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture Permeability Variation: An Analytical Approach

Download PDF
18 pages

Loading next page...
 
/lp/springer-journals/elastic-brittle-plastic-behaviour-of-shale-reservoirs-and-its-E3MDTmroVX

References (54)

Publisher
Springer Journals
Copyright
Copyright © 2018 by The Author(s)
Subject
Earth Sciences; Geophysics/Geodesy; Civil Engineering
ISSN
0723-2632
eISSN
1434-453X
DOI
10.1007/s00603-017-1392-y
Publisher site
See Article on Publisher Site

Abstract

Shale gas has recently gained significant attention as one of the most important unconventional gas resources. Shales are fine-grained rocks formed from the compaction of silt- and clay-sized particles and are characterised by their fissured texture and very low permeability. Gas exists in an adsorbed state on the surface of the organic content of the rock and is freely available within the primary and secondary porosity. Geomechanical studies have indicated that, depending on the clay content of the rock, shales can exhibit a brittle failure mechanism. Brittle failure leads to the reduced strength of the plastic zone around a wellbore, which can potentially result in wellbore instability problems. Desorption of gas during production can cause shrinkage of the organic content of the rock. This becomes more important when considering the use of shales for CO sequestration purposes, where C O adsorption-induced swelling can play an important role. These phenomena lead 2 2 to changes in the stress state within the rock mass, which then influence the permeability of the reservoir. Thus, rigorous simulation of material failure within coupled hydro-mechanical analyses is needed to achieve a more systematic and accurate representation of the wellbore. Despite numerous modelling efforts related to permeability, an adequate representation of the geomechanical behaviour of shale and its impact on permeability and gas production has not been achieved. In order to achieve this aim, novel coupled poro-elastoplastic analytical solutions are developed in this paper which take into account the sorption-induced swelling and the brittle failure mechanism. These models employ linear elasticity and a Mohr–Cou- lomb failure criterion in a plane-strain condition with boundary conditions corresponding to both open-hole and cased-hole completions. The post-failure brittle behaviour of the rock is defined using residual strength parameters and a non-associated flow rule. Swelling and shrinkage are considered to be elastic and are defined using a Langmuir-like curve, which is directly related to the reservoir pressure. The models are used to evaluate the stress distribution and the induced change in perme- ability within a reservoir. Results show that development of a plastic zone near the wellbore can significantly impact fracture permeability and gas production. The capabilities and limitations of the models are discussed and potential future develop- ments related to modelling of permeability in brittle shales under elastoplastic deformations are identified. Keywords Coupled geomechanics · Fractured reservoir · Shale gas · Permeability modelling · Poro-elasto-plasticity · Analytical solution List of symbols  In situ horizontal stress under initial reser- voir conditions Applied Mechanics Δ Change in mean horizontal stress from the ,  ,  Stress components in cylindrical coordinates rr  zz initial in situ conditions Mean horizontal stress (i.e. average of ,  ,  Strain components in cylindrical coordinates rr  zz ,  ) rr  e,m The elastic part of mechanically induced strain e,s The sorption-induced volumetric strain * Mohsen S. Masoudian (elastic only) [email protected] The plastic part of strain e,s 1  The in situ swelling strain under initial res- Nottingham Centre for Geomechanics, University ,0 ervoir conditions of Nottingham, University Park, Nottingham NG7 2RD, UK Vol.:(0123456789) 1 3 1566 M. S. Masoudian et al. Biot’s coefficient ZDCS Acronym for zero displacement condition E Elastic modulus on wellbore boundary and constant stress Poisson’s ratio condition on outer boundary c,  Cohesion and friction angle, respectively ZDZD Acronym for zero displacement condition on c ,  Post-failure cohesion and friction angle, both boundaries r r respectively f (r) A special function defined from distribu- , Y Strength parameters in Mohr–Coulomb tions of pore pressure and sorption-induced failure criterion defined using cohesion and strain friction angle F (r) Integral of pore pressure function over the , Y Post-failure strength parameters in Mohr– domain of the reservoir r r Coulomb failure criterion defined using F (r) Integral of sorption-induced strain function post-failure cohesion and friction angle over the domain of the reservoir Dilation angle F(r) A special function defined from F (r) and A function of dilation angle F (r) u Radial displacement at the elastic–plastic ep interface Geometry 1 Introduction r, , z Directions in cylindrical coordinates (the origin is the centre of the wellbore) Unconventional gas reservoirs are increasingly being con- R Radius of the disc-shaped reservoir sidered to provide a relatively clean energy source to meet R Radius of the wellbore burgeoning global demands. Unconventional gas resources, R Radius of plastic zone developed around the ep such as shale, represent around 40% of the remaining wellbore resource of technically recoverable natural gas; hence, they play an important role in the future global energy market Flow and Sorption Parameters P Pore pressure (fracture porosity) (McGlade et al. 2013). In addition, unconventional gas res- ervoirs are being considered as a host rock for CO geo- p Uniform pore pressure in the reservoir under 0 2 initial condition sequestration, which may also provide an enhanced gas recovery (EGR) technique. Shales are fine-grained rocks p Wellbore pressure (constant during gas production) formed from the compaction of silt- and clay-sized particles and are characterised by their fissured texture and very low Q Coefficient of pore pressure reduction due −18 −23 2 to gas production with the steady-state loga-permeability (10 –10  m ) (Alexander et al. 2011). Gas shales can be considered as a dual porosity rock consisting rithmic pressure distribution Q Dimensionless production rate of fractures and matrix blocks (Fig. 1), with a porosity of up to 15% (Wang and Reed 2009). Fracture sets can be subverti- k Permeability of the reservoir k Permeability of the reservoir under initial cal or parallel to bedding, may be filled with mineral gauges, have openings ranging in scale from micrometres to centi- conditions c Fracture compressibility metres, and spacing ranging in scale from a few centimetres to several metres (Gale et al. 2014). The matrix usually con- V Volume of adsorbed gas V , b Maximum monolayer sorption capacity and tains an organic content of 2–10%, while those with higher organic contents are usually considered too immature for Langmuir isotherm constant, respectively Maximum swelling strain production development (Alexander et al. 2011). The matrix contains porosity in the scale of µm to nm in size within Solution‑Dependent Definitions different organic (fibrous) and nonorganic (clastic) struc- A, B, K Constants defined in stress and displacement tures (Mehmani et al. 2013). Natural gas can be free within solutions both the fracture (primary porosity) and clastic structure of C , C , C Integration constants for stress solutions 1 2 3 the matrix (secondary porosity), or in an adsorbed state on CSCS Acronym for constant stress condition on the surface of the porous structure of organic content. The both boundaries adsorption of gas in organic geomaterials can be described CSZD Acronym for constant stress condition on using an isotherm, among which Langmuir’s is most widely wellbore boundary and zero displacement used (Gensterblum et al. 2015). For example, the gas content condition on outer boundary of five reservoirs in the USA, with in situ pressures between 2 and 28 MPa (600–8500 m deep), ranged from 1 to 10 m / 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1567 Fig. 1 a Fractures in shale core samples (Gale et al. 2014), b a subvertical calcite-cemented fracture (Soeder 1988), c simplified dual porosity concept for shale ton, of which 20–85% was adsorbed in the organic content of gas shales to changes in the stress level (e.g. Amann- (Curtis 2002; Hill and Nelson 2000). Hildenbrand et al. 2012; Chen et al. 2015a; Cho et al. 2013; Due to the multi-scale and multiple physical phenomena Li et al. 2016; Spencer 1989; Zhou et al. 2016). However, associated with shale gas production, the mechanism of gas these studies have only considered the elastic behaviour flow is a combination of viscous (Darcy and non-Darcy) of the rock. Plastic deformations within the reservoir rock flow (Huang et al. 2016), and gas-slippage (Klinkenberg) can have significant implications for production, injectiv - effect (Mehmani et al. 2013), as well as continuum (Fickian), ity, and stability of the wellbore. Studies have demonstrated free-molecule (Knudsen), and surface (Yuan et al. 2014) dif- the significance of both elastic and plastic deformations on fusions. On the other hand, the gas sorption (adsorption and/ wellbore producibility and/or injectivity as well as perme- or desorption) also plays an important role in the gas trans- ability prediction in other types of reservoirs (e.g. Cui et al. port mechanism in shale reservoirs, given that the amount of 2007; Han and Dusseault 2003; Masoudian et al. 2016b). adsorbed gas in the organic content is comparable with the In addition, the stability of wellbores is an important issue free gas. Despite the large number of studies on primary and in any oil and gas production project. Therefore, improved CO -EGR, many aspects of shale gas development are still models are needed to properly estimate the deformation and not well understood, mainly due to the complex interactions stress distributions around the wellbore. Experimental stud- of fluid flow with different thermal, chemical, and mechani- ies have also revealed the brittle mechanism of failure in cal phenomena. For instance, withdrawal of water and gas shale (e.g. Amann et al. 2011; Hull et al. 2015; Sone and from a subsurface formation results in significant reservoir Zoback 2013), which implies the need for improved models pressure alteration, which in turn leads to a redistribution of that consider this failure mechanism. The effect of gas and/ stress and strain within the shale formation. Gas and water or water sorption-induced swelling should also be taken into sorption-induced swelling/shrinkage (Lyu et al. 2015; Yuan account when developing permeability models. et al. 2014), thermal phenomena (e.g. during thermal stimu- The most common geomechanical assumption when pre- lation (Carroll et al. 2011)), and chemical reactions of water dicting permeability in gas reservoirs is of one-dimensional or CO (in case of CO -EGR) with shale (Carroll et al. 2011) elastic deformation (vertical deformation). However, this 2 2 can additionally impact the stress–strain relations. Redis- assumption cannot provide an accurate estimation, especially tribution of stress and strain within any reservoir can sub- near the wellbore where plastic deformations may occur. sequently lead to other changes, such as reservoir perme- On the other hand, using fully numerical hydro-mechanical ability alterations, subsidence/uplift of the ground surface, models leads to substantially more expensive computations. fault-reactivation mechanisms, and crack development in the This is especially important in shale gas reservoirs where caprock (Masoudian et al. 2016a). complex multiple porosity reservoir models are utilised. Prediction of permeability has received considerable This study aims to extend previous work (Masoudian and attention due to its significant effect on gas production from Hashemi 2016) by further developing brittle elastoplas- reservoirs. Numerous studies have related the permeability tic models for reservoirs coupled with sorption-induced 1 3 1568 M. S. Masoudian et al. swelling/shrinkage effects. In this paper, models are devel- approach is adopted, where the rock is considered as a con- oped with boundary conditions that simulate different well tinuum at the field scale, while the fractured nature of the completions and geological settings (e.g. open-hole and rock is captured by the effect of stress on fracture porosity. cased-hole completions). The models are then used to esti- Figure 2a illustrates the geometry of the model and the multi- mate the stress distribution within the reservoir, with which continuum approach adopted in this study. Figure 2b depicts the changes in fracture permeability within both elastic and the stress–strain relations within the brittle elastoplastic frame- plastic zones are predicted. This is an important aspect of the work used in this paper. Note that the mechanical behaviour work that has been largely neglected in other studies, mainly of brittle rock is idealised with a sudden drop of strength from due to lack of computational resources for field scale fully peak to residual value. Similar to perfectly plastic behaviour, coupled hydro-elastoplastic simulations. This work aims to there is no further change in stress beyond the yield point as highlight the role of plastic deformations in the prediction plastic strains continue to accumulate. The governing equa- of permeability and its impact on gas production. tions employed in this paper are explained below. Starting with the classical approach in plasticity, strains can e p be decomposed into elastic (  ) and plastic (  ) components, 2 Theoretical Framework and Governing where the elastic strain can be further divided into elastic e,m Equations mechanical strain, (  ), and elastic swelling (or shrinkage) e,s strain, (  ). Note that adsorption induces swelling while des- The reservoir is assumed to be a disc-shaped homogenous iso- orption leads to shrinkage of the matrix blocks. Assuming tropic continuum. A vertical wellbore is assumed to be drilled axial symmetry of the problem and considering a cylindrical at the centre of the disc, through which natural gas can be pro- system of coordinates (r ,  , z for radial, tangential and vertical duced. It should be noted that these simplifying assumptions directions, respectively), the strain compatibility equation is are employed in order to develop a fully analytical framework rr in this paper. Clearly, these assumptions cannot fully represent = (1) r r the heterogeneous nature of shale; however, they do not pre- vent evaluation of the effect of elastoplastic deformations on where r is the radial distance from the centre of a disc- stress and permeability distributions within reservoirs, which shaped reservoir, and subscripts rr and  represent the is the focus of this paper. Here, a simplified multi-continuum radial and tangential coordinates, respectively. In the light of Fig. 2 Schematics of a the geometry of the model, the multi-continuum approach, and the swelling/shrinkage effect on matrix size, and b the stress and strain relationships in the adopted elastoplastic framework 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1569 paper, a simplified pressure solution is sufficient. To arrive Biot’s definition of effective stress (  =  − P ), the equa- tion of equilibrium can be written as at a simplified pressure solution, the reservoir is initially assumed to be under a uniform pore pressure regime. Gas � � � rr  rr production from the wellbore with a fixed pressure ( p ) (2) w = − r r r leads to the development of a non-uniform pore pressure distribution, with pressure remaining equal to the initial where P is pore pressure,  is Biot’s coefficient, and  is pressure at the outer boundary only. Thus, a steady-state effective stress. Assuming that swelling strain is analogous solution of radial flow can be found in the following loga- to thermal strain, the elastic stress–strain constitutive rela- rithmic form (Cui et al. 2007) tions for an isotropic rock can be written as e,s P = p + Q ln (8) =  −   +  − (1 − 2)P + rr rr  zz E 3 e,s where p and R are the initial reservoir pressure and radius =  −   +  − (1 − 2)P + (3) 0 0 zz rr E 3 of the disc-shaped reservoir, respectively. The term Q is e,s defined as =  −   +  − 1 − 2 P + ( ) zz zz rr E 3 p − p w 0 Q = (9) ln R ∕R w 0 e,s where E ,  and  are the elastic modulus, Poisson’s ratio, where R and p are the radius and pressure of the wellbore, and volumetric swelling strain, respectively. The linear w w respectively. Mohr-Coulomb yield criterion can be written as The volume of adsorbed gas, V , can be estimated using � � =  + Y (4) rr Langmuir’s isotherm as a function of pressure ( P) where  and  are the tangential and radial effective rr V b P stresses, respectively. The parameters  and Y are defined V = (10) 1 + b P using cohesion, c , and friction angle,  , of the rock as where V is the maximum sorption capacity and b is the 1 + sin Langmuir isotherm constant. Initially, the matrix of the res- 1 − sin (5) ervoir rock contains a specified amount of adsorbed gas at 2c cos Y = the initial pressure ( p ) which changes as the reservoir pres- 1 − sin sure decreases due to production. As a result, the volumetric e,s e,s swelling strain changes from the initial state of  to  and ,0 Note that by substituting the residual cohesion, c , and is considered to have a linear relationship with adsorbed gas friction angle,  , into these equations, the residual param- volume, which gives eters  and Y can be defined for failed brittle rock. Consid- r r ering a non-associated flow rule, the plastic components of  b P e,s L (11) radial and tangential strains are related as (Park and Kim 1 + b P 2006) where  is the maximum swelling strain that can only take p L +  = 0 rr  place when the gas within the matrix is equal to the sorption e,s (6) u u e e capacity. Note that  is used to account for the in situ equi- +  =  +  ,0 rr r r librium state of the reservoir and can be estimated directly using the equation above, replacing the pressure with p . where u is the radial displacement, and  is a function of There are a number of fracture porosity/permeability dilation angle,  , as models developed for reservoirs where adsorption-induced 1 + sin  swelling or shrinkage is associated with gas injection or (7) production (e.g. Gensterblum et al. 2015; Pan and Connell 1 − sin 2012). However, these models adopt two main approaches There are a large number of studies that have investigated for porosity/permeability modelling: dependency on stress or modelling approaches for flow of gas (methane or carbon strain. Both of these approaches are closely related and may dioxide) in shale reservoirs, where the complex multi-porous lead to very similar results under certain conditions (Palmer nature of shale is carefully considered and different mass 2009), as they both consider a cubic law to link permeabil- transport mechanisms (i.e. viscous and diffusive) are taken ity to fracture porosity (consistent with a Kozeny–Carman into account. However, for the comparative analyses in this type relationship). In this paper, results are obtained with 1 3 1570 M. S. Masoudian et al. the most widely used approach where the permeability is where C and C are the integration constants that can be 1 2 related to the change in horizontal effective stress using an determined from the boundary conditions, and exponential equation as (Shi and Durucan 2005) (1 − 2) E e,s f (r) = P + (1 − ) 3(1 − ) (1 − 2) = exp − 3c Δ E (12) k h F(r) = F (r) + F (r) (1 − ) 3(1 − ) 2 2 r r where k is permeability, the subscript 0 denotes an initial F r = P d =− Q + P ( ) 4 2 � 0 value, c is fracture compressibility, and Δ is the change 2 R 1 + b p 1 + b P 0 L 0 L e,s s in horizontal effective stress, which can be defined as the F (r) =    d =  − exp − 2 Ei 2 × v L 2 b Q b Q b Q L L L mean of the change in the radial and tangential components (14) of effective stress. As illustrated in Fig.  2a, any element of where Ei(x) is the exponential integral function. rock (i.e. at a specific value of r ) can be conceptualised as The stress and strain relations within the plastic zone a hybrid of fractures and matrix blocks and any change in near the wellbore are governed by a Mohr–Coulomb failure stress state (due to pore pressure reduction and swelling/ envelope and a non-associated flow rule. Note that the brittle shrinkage effects) results in a change in the aperture size failure mechanism is implemented using residual strength of the fractures within that element. Equation 12 therefore parameters (i.e. cohesion and friction angle); perfectly plastic provides an efficient method for relating the response of frac- behaviour could be modelled by using peak strength values. tures and their permeability at fracture scale to the changes The general solution in the plastic zone can then be written as in effective stress at field scale. The use of mean horizontal effective stress implies a matchstick matrix model, which � 𝛾 −1 𝜎 = C r − K rr r<R is consistent with the plane-strain assumption. Closed-form ep � 𝛾 −1 solutions of permeability and porosity can then be obtained 𝜎  = Y + 𝛾 C r − K r r 3 𝜃𝜃 r<R ep by substituting the stress and strain equations for different 𝛽 +𝛾 r 𝛽 +𝛾 Eu R − r (1 + 𝜈 ) ep ep −𝛽 𝛽 loading and boundary conditions into the above equations. u = r R − C A r<R ep 3 ep E 1 + 𝜈 𝛽 + 𝛾 It should be noted, however, that experimental results have 𝛽 +1 𝛽 +1 shown that the relationship between post-failure permeabil- R − r ep 1 + 𝛽 − BY − AK + E F R − F (r) r 𝜀 ,𝛽 ep 𝜀 ,𝛽 ity and the stress–strain state can be much more complex 𝛽 + 1 3 than what these models offer (e.g. Carey et al. 2015). How - (15) ever, due to the lack of dedicated post-failure models, this where widely used simplified approach was adopted in this paper. Limitations of the model are discussed later in the paper. A = 1 −  −  +   −  − r r r B =  −  − Y − Q 3 Analytical Solutions for Elastic–Brittle– K = − 1 Plastic Rock (16) 1+ 1+ 1 + b p r 0 L 0 F (r) =  − exp −(1 + ) Following the approach employed by Masoudian and 1 +  b Q b Q L L Hashemi (2016), analytical solutions for stress and strain 1 + b P around a wellbore can be found by assuming two distinct × Ei (1 + ) b Q concentric zones; a plastic zone near the wellbore and an where Y and  are obtained with residual values of cohe- elastic zone towards the outer boundary. The derivation and r r sion and friction angle in Eq.  5, and u is the displace- development of these solutions are partially presented in ep ment at the elastic–plastic interface, which can be found “Appendix 1”. The interface of elastic and plastic zones is at considering continuity of displacement at the interface, i.e. a radial distance of R from the wellbore centre. The general ep by replacing the value of R in the elastic solution of dis- solution within the elastic zone ( r > R ) can be written as ep ep placement (Eq. 13) once the value R is found. To obtain ep F(r) − C R , the integration constants C , C , and C are found by 𝜎 = C + ep 1 2 3 rr 1 r>R ep substituting different sets of boundary conditions into the F(r) − C general solution, as presented in Table 1. Each set of bound- 𝜎 = C − + f (r) 𝜃𝜃 1 r>R ep ary conditions describes a particular geological and well (1 + 𝜈 ) completion condition. The constant stress at the wellbore u = (1 − 2𝜈 )r C − 𝜎 + 𝛼 p r>R 1 0 0 ep represents an uncased hole; hence, the total stress applied F(r) 2 e,s to the inner boundary is equal to the bottomhole pressure (13) + − + r𝜀 v,0 r r 3 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1571 ( p ); i.e. zero effective stress on wellbore wall. On the other outer boundary. As such, case CSCS represents constant hand, zero displacement on the wellbore wall represents an stresses (CS) on both the wellbore wall and outer boundary, ideally cased hole that does not allow the wellbore wall to case ZDCS represents zero displacement (ZD) of the well- converge. For the outer boundary, the constant stress rep- bore wall and constant stress on the outer boundary, case resents an ideal scenario where the reservoir is completely CSZD represents constant stress on the wellbore wall and separated from the adjacent formation. Zero displacement at zero displacement of the outer boundary, and case ZDZD the outer boundary represents a scenario where the reservoir represents zero displacement on both the wellbore wall and is large enough to prevent deformation of the rock mass at outer boundary. the far-field. For example, if the reservoir was cut by low- The interface between elastic and plastic zones is found angle faults at the outer boundary, the reservoir may be considering the continuity of radial stress at r = R and a ep best simulated by applying a constant stress, whereas a zero trial-and-error scheme. To achieve this, an equality is con- horizontal displacement boundary would be more suitable structed using the solution of radial stress within both elastic if the reservoir was cut by steep faults and was in contact and plastic zones, and then its root for R can be found. If this ep with mechanically strong geological strata. However, due equality does not have a root in the R , R domain, it means w 0 to the geological complexity of reservoirs, actual boundary that the rock mass does not include a plastic zone and hence a condition will fall somewhere between these two condi- fully elastic solution is used by replacing R by R in Eq. 13 ep w tions; hence, different boundary conditions should be tested. with the corresponding elastic integration constants presented The boundary condition acronyms consist of four letters, in “Appendix 2”. These formulations and the introduced meth- with the first pair describing the condition of the wellbore odology were implemented in a MATLAB code and used to wall, and the second pair describing the condition of the obtain the results discussed in the following section. Table 1 Integration constants for different sets of boundary conditions Boundary conditions (acronym) Relationships 2 2 R Y−f R −(−1)P R +(+1) F R −F R +R Constant stress at R ( ) ( ) ( ) ( ) (17) ep ep ep ep 0 0 0 C = 1 2 2 R (+1)−R (−1) Constant stress at R 0 ep (CSCS) C = C −  R + F R 2 1 0 0 C = 3  −1 2 2 R Y−f R −(−1)P R +(+1) F R −F R +R Zero displacement at R ( ) ( ) ( ) ( ) (18) ep ep ep ep 0 0 0 C = 1 2 R (+1)−R (−1) Constant stress at R ep 0 0 (ZDCS) C = C −  R + F R 2 1 0 0 +1 +1 Eu R −R +  1+ ep ep w C =   × R − BY − AK + E F R − F R ep 3 + r , ep , w r + r 1+ +1 3 A R −r ep Constant stress at R (19) w 1 C = × R Y − f R − ( − 1)P R 1 2 2 ep ep ep Zero displacement at R R (1−)−R (1+)(1−2) 0 ep 0 (CSZD) E e,s +( + 1) F R − F R + R  − (1 − 2)  − p ep 0 0 0 0 v,0 E e,s C = F R − R (1 − 2) C −  + p + 2 0 1 0 0 0 v,0 C = −1 3 r Zero displacement at R (20) w 1 C = × R Y − f R − ( − 1)P R 1 2 ep ep 2 ep Zero displacement at R R (1−)−R (1+)(1−2) 0 ep (ZDZD) E e,s +( + 1) F R − F R + R  − (1 − 2)  − P ep 0 0 0 v,0 2 e,s C = F R − R (1 − 2) C −  + p + 2 0 1 0 0 0 v,0 +1 +1 Eu R −R +  1+ ep ep w C = × R − BY − AK + E F R − F R ep 3 + r , ep , w + 1+ +1 3 A R −r r ep 1 3 1572 M. S. Masoudian et al. ZDZD scenarios (representing cased-hole completions), the 4 Results and Discussion radial stress significantly increases near the wellbore (mov - ing towards the wellbore from right to left in Fig. 3), since In order to examine the stress change and permeability evo- the rock is restricted on the wellbore wall. On the other hand, lution around a gas production well, a set of values were the tangential stress in these cases monotonically decreases chosen for the model input parameters, as listed in Table 2. as you approach the wellbore. This is because the support These values mainly represent gas shale reservoirs (suitable provided by the casing does not allow radial deformations for hydraulic fracturing due to their brittleness) and were and the release of radial stress within the rock mass. As a chosen to be within the reported ranges from the correspond- result, the change in radial and tangential effective stresses ing references. Studies have shown that hydro-mechanical (which is due to the interaction between pore pressure in the properties of near-wellbore shales are influenced by the pres- fractures and shrinkage of the matrix of the rock) is not large ence of natural and hydraulically generated fractures (Liang enough to satisfy the yield criterion (due to large minimum et al. 2014; Song et al. 2014), and the exposure of rock to principal stress levels) and the rock remains elastic, i.e. there drilling mud (Ma and Chen 2015). However, these effects is no plastic zone around the wellbore. It should also be are not considered in this paper. It should also be noted that noted that the case with constant stress at the outer bound- the sorption and swelling parameters correspond to those of ary (ZDCS) exhibits higher stress levels compared to the organic rich shales. case with a fixed outer boundary (ZDZD). This is because With production of gas from the reservoir, the pore pres- in ZDCS, the applied stress at the far-field compresses the sure decreases following Eq. 8, and the uniform initial pore rock mass towards the wellbore and results in a larger stress pressure is replaced by a radial distribution with maximum applied on the casing of the wellbore. pressure at the outer boundary (initial reservoir pressure,p ) In the CSCS and CSZD cases (representing open-hole and minimum value on the wellbore wall (wellbore pressure, completions), as the wellbore is radially deformed, the radial p ). During production, the gas is desorbed from the matrix stress decreases towards the wellbore pressure of 5 MPa in of the rock which leads to matrix shrinkage. The change the approach to the wellbore. The tangential stress in these in pore pressure combined with sorption-induced shrinkage cases, however, increases to some point and then suddenly leads to a redistribution of stress within the rock mass. The decreases in the approach to the wellbore. This is because distributions of radial and tangential total stresses for all sets the circular opening allows some deformation of the rock of boundary conditions are presented in Fig. 3. Note that the towards the wellbore which is also enhanced by the produc- effective stresses can be evaluated by subtracting pore pres- tion-induced reduction of pore pressure and the desorption- sure from the total stresses (considering  = 1 ); this will be induced shrinkage of the rock matrix. In these cases, the discussed later in the paper. It can be seen that in ZDCS and combination of stress levels satisfies the Mohr–Coulomb Table 2 Input parameters for Parameter Symbol Value Unit References the base simulation case Young’s modulus E 20 GPa Rutqvist et al. (2013) Peak cohesion c 6.0 MPa Ma and Chen, (2015) Peak friction angle  32.23 degrees Ma and Chen (2015) Residual cohesion c 4.39 MPa Ma and Chen (2015) Residual friction angle  25.17 degrees Ma and Chen (2015) Dilation angle  10 degrees Islam and Skalle (2013) Maximum sorption-induced swelling  0.134 % Chen et al. (2015b) Fracture compressibility c 0.02 1/MPa Chen et al. (2015a) Poisson’s ratio  0.22 – Ma and Chen (2015) Langmuir’s isotherm constant b 0.25 1/MPa Chen et al. (2015b) In situ horizontal stress  31.784 MPa Corresponding to 1500 m depth following Rutqvist et al. (2013) Initial reservoir pressure p 24.63 MPa Rutqvist et al. (2013) wellbore pressure p 5 MPa Wellbore radius R 0.216 m Ma and Chen (2015) Biot’s coefficient  1.0 – Chen et al. (2015a) Outer radius R 100.0 m – 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1573 Fig. 3 Profiles of stress within the reservoir: a radial stress, b tangen- Fig. 4 Distribution of mean horizontal stress for all cases: a total tial stress stress, b effective stress failure criterion and a plastic (failure) zone is developed discussed in detail later. The mean horizontal total stress around the wellbore; the maxima in tangential stress pro- decreases as it closely follows the reduction of pore pres- files coincide with the elastic–plastic interface. Note that sure from the far-field towards the wellbore (Fig.  4a). In the sharp drop of the tangential stress profile at the elas - cased-hole analyses (ZDCS and ZDZD), the rock remains tic–plastic interface is due to the brittle failure mechanism fully elastic and therefore the mean horizontal total stress of the rock. In other words, whereas the strength of the rock decreases monotonically from the far-field towards the immediately at the elastic–plastic interface is given by peak wellbore, with the minimum value occurring at the well- properties (peak cohesion and friction angle), the strength bore wall. For CSCS and CSZD, however, the mean total of the rock within the plastic zone is given by the residual stress decreases sharply at the elastic–plastic interface, fol- values of cohesion and friction angle. The radius of the plas- lowing the trends of radial and tangential stresses shown tic zone is larger in the CSCS case (0.324 m) compared to in Fig. 3. The mean horizontal effective stress is shown in that in CSZD (0.304 m) because of lower stress levels in the Fig. 4b, where a maximum can be observed at the interface latter case. of the elastic–plastic boundary in CSCS and CSZD cases. The mean horizontal total and effective stresses, defined It can also be seen that within the elastic zone, the solutions as  =  +  ∕2 and  =  − p , respectively (Cui with similar outer boundary conditions show similar trends m rr  m et  al. 2007), are depicted in Fig.  4. These definitions are (CSCS vs. ZDCS and CSZD vs ZDZD), even though their important in the prediction of permeability, which will be corresponding  and  values are different outside the rr 1 3 1574 M. S. Masoudian et al. open-hole conditions. For example, the radial stress at the elastic zone interface. In other words, for  and  within wall of the cased hole is much larger than that in the open the elastic zone, the effect of wellbore boundary condition hole, while the tangential stress decreases due to the lack diminishes and the outer boundary has the dominant effect of any radial displacement. With an open hole, on the other on the distribution of mean horizontal stress. hand, the reduction in radial stress is compensated by the Since the main purpose of this study is to understand the increase in the tangential stress. Thus, following the discus- effect of plasticity, the profiles of mean horizontal stress for sion of Figs. 4 and 5 shows that the effect of brittle-plastic the two open-hole cases need further discussion. These two behaviour of rock is only significant within the failure zone cases were also modelled using their corresponding fully near the wellbore. elastic solutions, as depicted in Fig. 5a, b. Figure 5a shows The estimated permeability is shown in Fig. 6a, which that the mean horizontal stress estimated using the elastic as expected inversely follows the mean horizontal effective solutions for CSCS and CSZD are almost identical to those stress. For elastoplastic cases, the minima are at the inter- calculated for the cased-hole cases ZDCS and ZDZD (shown face of the elastic and plastic zones. This corresponds to in Fig. 4a), respectively, although the distributions of radial the location where the rock is on the verge of failure and and tangential stresses are very different (Fig.  5c, d). This therefore can sustain the maximum stress levels; hence, the is discussed in more detail in “Appendix 3”, but it is mainly fractures are at their minimum aperture. On the other hand, because the distribution of radial and tangential stresses in in the elastic solutions, the mean effective stress increases elastic models depends on the condition of the wellbore monotonically towards the wellbore and therefore the rock wall; however, they change in a way that their summation near the wellbore is under the maximum compression, (  +  ) remains the same between the cased-hole and rr Fig. 5 Comparison of elastic and plastic solutions for open-hole cases (ZDCS and ZDZD): a mean total stress, b effective mean stress, c (CSCS and CSZD) and the comparison of elastic solutions of open- radial stress, d tangential stress hole cases (CSCS elastic and CSZD elastic) with the cased-hole cases 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1575 Fig. 6 Permeability of the reservoir with open-hole completion: a comparison between elastic and plastic solutions, b comparison between brit- tle–plastic and perfectly plastic solutions, c comparison of the solution with and without swelling effect which leads to the lowest values of permeability. When the and reduction of permeability. On the other hand, the reduc- constant radial stress condition is used for the outer bound- tion of pore pressure leads to desorption of gas from the ary (CSCS), the rock mass is pushed towards the wellbore matrix blocks of the rock. Consequently, the matrix blocks which in turn leads a higher degree of fracture closure and shrink and the width of fracture openings increases. There- lower permeability compared to the CSZD models (zero fore, depending on the poro-mechanical and swelling param- displacement at the outer boundary). The effect of brittle eters of the reservoir, the production of gas can lead to either failure on rock permeability is shown in Fig. 6b. Compared a reduction or an increase in permeability. In the example to perfectly plastic behaviour (with peak strength values), analysed here, the maximum volumetric strain induced by the brittle-plastic behaviour leads to a larger failure zone swelling/shrinkage was chosen to be 0.13% (an order of and lower permeability. This is because, within the plastic magnitude smaller than that of coal). Therefore, the effec- zone, the brittle model sustains lower stress levels than the tive stress change induced by pore pressure reduction had a peak strength perfectly plastic model, which leads to reduced more dominant effect. Figure  6c shows how the desorption- permeability in the brittle model. induced shrinkage can influence the fracture permeability It should be noted that permeability of fractured rock of a reservoir. The cases including swelling effects resulted (and the stress level within the rock mass) depends on the in larger permeability throughout the reservoir due to the opposing effects of two coupled phenomena: pore pressure shrinkage of matrix blocks. Another interesting observation change and sorption-induced swelling/shrinkage. First, the is that the radius of the plastic zone is influenced by the reduction of pressure in fractures leads to an increase in the swelling strain despite the fact that it is considered to be effective stress level, which results in closure of fractures an elastic component. This is because the shrinkage of the 1 3 1576 M. S. Masoudian et al. matrix block results in relaxation of stress; therefore, the towards lower residual values. This is because a stronger compressive stress applied to the rock mass reduces and the failed rock (higher residual strength) can sustain larger failed zone radius reduces. stresses and hence provides better support to the rock mass. Figure 7 depicts the radius of the plastic zone with respect The plastic behaviour of rock has been shown to affect to change in maximum sorption-induced swelling (  ) and the permeability of a reservoir only in the area very near residual strength parameters ( c ,  ). As sorption-induced the wellbore (i.e. within the plastic zone); however, it may r r volumetric change increases from zero to 1%, the matrix have a significant impact on the production flowrate. In order blocks shrink and hence the stress level decreases within to estimate the effect of near-wellbore permeability on pro- the rock mass. As a result, the radius of the plastic zone duction rate, a simplified approach is adopted in which the decreases. No plastic zone developed when using a maxi- equivalent permeability of the whole reservoir is defined by mum swelling strain larger than 1.09% in these examples. considering parallel ring-shaped layers having different per - In other words, the shrinkage of the matrix can help the meability values. In other words, by considering equal radial mechanical integrity of the wellbore. Increasing the values flowrate in each ring-shaped element of the reservoir, the of residual cohesion and friction angle also leads to a smaller equivalent permeability can be found. Hence, a dimension- plastic zone; the effect is proportionally greater at lower val - less production flowrate (somewhat similar to productivity ues of residual strength, i.e. increasing slope of the curves index) can be defined, simply as the ratio of production rate of the reservoir with a spatially variable permeability to that of a reservoir with constant initial permeability, which gives k dr Q = ln R − R (21) D 0 w rk(r) This parameter can be used to illustrate how changes in permeability influence the rate of gas production from a res- ervoir. As shown Fig. 8a, an increase in the maximum swell- ing strain leads to a reduction in the production rate. This is because, although the shrinkage of the shale matrix generally increases the permeability, it restricts the development of the plastic zone (the zone with significantly higher permeability), as demonstrated in Fig. 7, and results in a reduced equivalent permeability of the reservoir. Figure 8b shows how the small plastic zone around the wellbore can significantly affect the production rate. First, note that having a plastic zone devel- oped in elastoplastic solutions with enhanced permeability can lead to higher production rates compared to their elastic counterparts. The production rate also increases with lower residual strength parameters, due to enlargement of the plastic zone around the wellbore, as previously shown in Fig. 7. Such a significant effect of the plastic zone on the equivalent perme- ability of the whole reservoir and its production rate is due to the fact that the wellbore and its surrounding zone is the bot- tleneck of the production. This has also been addressed in the context of formation damage where the adverse effect of drill- ing mud on permeability in the vicinity of the wellbore is stud- ied. Thus, development of a plastic zone with higher perme- ability around the wellbore can have significant implications on production of gas from shale reservoirs. This highlights the fact that elastic solutions are not able to properly predict reservoir permeability and gas production rates. On the other hand, the use of casings for wellbore stability can significantly Fig. 7 Radius of the elastic–plastic interface with respect to change in a maximum sorption-induced swelling, b ratio of residual strength reduce the producibility of wells. Obviously, the mechanical parameters 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1577 Fig. 9 Assessment of the developed methodology for prediction of Fig. 8 Dimensionless production rate with respect to change in: a permeability: a exponential permeability model for different shales maximum sorption-induced swelling, b ratio of residual strength (data from Bustin et al. 2008; Dong et al. 2010; Gutierrez et al. 2000; parameters Pathi 2008), b field data of permeability around a tunnel and model predictions integrity of wellbores is of great importance; however these the permeability for many shales at small scale (core scale), results suggest that the development of novel wellbore casing confirming the validity of the exponential relationship for solutions (e.g. casings that allow a plastic zone to develop but fractured shale reservoirs. Unfortunately, field data present- prevent wellbore caving) could provide significant benefits to ing the effect of elastoplastic deformations on permeability the producibility of a wellbore. of shale reservoirs are not available. However, Souley et al. (2001) reported in situ measurements of permeability around a tunnel excavated in granite at Atomic Energy Canada Lim- 5 Limitations, Verification, and Future ited (AECL) Underground Research Laboratory in Canada Development of the Model to address issues related to the disposal of nuclear wastes. There are obviously some important differences between this As a measure of verification of the model, results are com- case and a shale gas reservoir, but these data are still valu- pared against experimental and field data in Fig.  9. First, the able as they provide insight into how plastic deformations exponential permeability model used in this study (i.e. Eq. 12) can influence permeability around a circular hole. In order is fitted to experimental datasets available in the literature to get a good match to in situ measurements, the proposed (usually from triaxial isotropic compression tests). Figure 9a model was used with no sorption-induced swelling, a fracture −1 shows that the use of Eq. 12 provides a good estimation of compressibility of 0.24 MPa , and a residual strength ratio of 1 3 1578 M. S. Masoudian et al. 0.83. Note that all other parameters were the same as those in well understood and cannot be properly predicted, especially Rutqvist et al. (2009) and Souley et al. (2001), where numeri- for highly anisotropic rocks with sparse organic contents such cal models were employed to replicate these measurements. It as shales (Kumar et al. 2016). Such studies can be even more should be noted that the estimated fracture compressibility in important when considering the CO injection for sequestra- this paper is within the range of values reported by Rutqvist tion and production enhancement purposes. The effect of shear et al. (2009). Although this case is simplified here and many displacement on fracture permeability has been discussed in a processes are neglected, the good agreement between the limited number of studies (Gutierrez et al. 2000; Rutqvist 2015); model and field data provides a verification of the accuracy this is another important topic that requires further investigation. of the proposed model. It can also be seen that considering the brittle behaviour of the rock is imperative to obtain better 6 Conclusions predictions of permeability. Several developments of the model could be considered to In this paper, a series of analytical solutions were presented in address the current limitations. Improvements could be made a coupled poro-elastoplastic framework, taking into account the by developing incremental forms of the solutions such that sorption-induced swelling/shrinkage of the rock matrix and the the stress–strain relations can be implemented within a time- brittle failure mechanisms of the rock mass. These solutions were dependent reservoir modelling framework (e.g. commercial used to model the complex inter-relations between different phe- reservoir simulators). However, care should be taken as the nomena during gas production from a gas shale reservoir. The model presented in this paper ignores the fact that severe results showed that the rock immediately adjacent to the well- softening behaviour may lead to localised failure, wherein bore may undergo plastic deformations, and depending on the the assumption of axisymmetric deformation becomes residual strength of the rock, this plastic zone may be enlarged, invalid (Lubliner 1990). Studies on the localisation effect while the increasing sorption-induced matrix shrinkage leads to on ground reaction curves for circular tunnels (e.g. Alonso smaller plastic zones. The results also showed that the perme- et al. 2003; Varas et al. 2005) found that results of radially ability of shale within the plastic zone can increase significantly symmetric solutions were valid from a practical engineering compared to its initial value, and while the developed plastic zone perspective. Thus, similar studies should be undertaken to is generally very small, its implications on producibility of the evaluate the extent to which the solutions developed in this wellbore can be significant. The advantages and limitations of paper are able to represent the actual changes in permeability the developed model were critically discussed. The model was of shale reservoirs. In addition, unconventional reservoirs verified using data relating to the permeability of a damaged zone such as gas shales are usually stimulated with hydraulic frac- around a tunnel. In the light of this discussion, and considering turing, which can change the mechanical behaviour of the relevant findings from the literature, a number of directions for rock in the vicinity of the wellbore (neglected in this paper), future research were proposed to address the lack of understand- which could be considered for future model developments. ing on the effect of elastoplastic deformations on permeability of Studies have proposed various models for permeability using shales. In addition, suggestions were made for further improve- different forms of the exponential relationship, but these have ment in the models by implementing their incremental forms in mainly concerned with the effect of isotropic effective stress on an advanced reservoir model and by considering more complex the permeability of shale and a correction of fracture compress- hydro-mechanical boundary conditions. These developments ibility (e.g. variation of c with effective stress). The perme - could provide more accurate predictions within the computa- ability model used in this paper assumes that the total sorption tionally efficient analytical framework presented in the paper. strain of the unconfined rock matrix contributes to permeability variation since the sorption strain incorporated in this model was Open Access This article is distributed under the terms of the Creative based on measurements made on an unconfined matrix. Zang Commons Attribution 4.0 International License (http ://crea tive comm ons.org/licenses /b y/4.0/), which permits unrestricted use, distribution, et al. (2015) introduced the concepts of ‘internal swelling’ and and reproduction in any medium, provided you give appropriate credit ‘internal swelling ratio’ for coal samples to evaluate the effect of to the original author(s) and the source, provide a link to the Creative sorption-induced swelling on permeability under confined con- Commons license, and indicate if changes were made. ditions. Such an approach could benefit future developments of permeability models. In addition, recent studies showed that the Appendix 1: Derivation of Stress behaviour of shales can be much more complex, especially when and Displacement Solutions considering shale reservoirs for C O -EGR. Due to the signifi- cance of micro-cracks in shale reservoirs, these findings imply The derivation of the solution of the elastic zone is presented an immediate need for better understanding of the mechanisms in the literature (e.g. Masoudian et al. 2016b), and only the of compression-induced fracture initiation and growth in shales. solution of the plastic zone is presented here. Substituting the Although studies related to hydraulic fracturing cover this topic post-failure yield criterion (an equation in the same form as to some extent, the relevant micro-scale phenomena are not 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1579 Eq. 4 with residual cohesion and friction angle replacing the radial stress will be postulated, which means that both elas- peak values) and the steady-state pressure distribution (Eq. 8) tic and plastic solutions must be equal at r = R . In other ep into the equation of equilibrium (Eq. 2) gives words, the elastic solution needs to satisfy the peak yield criterion (Eq. 4), because the rock at r = R is at its peak � � ep Y Q rr rr r strength. Therefore, substituting Eq. 13 into Eq. 4 gives a (22) + 1 −  = + r r r r solution for stresses r = R ep This is a first-order non-homogenous linear differential = 2C +  −  P rr 1 r=R r=R equation, which can be solved for  , using the integrating ep ep 1 +  1 − rr factor method. The value of  can be found by substituting e,s (27) +   − Y the solution into the post-failure yield criterion. The general r=R 3(1 − ) ep solution for  can be written as rr The effective radial stress at the elastic–plastic interface can Y C r 3 also be found from Eq. 23 and, by including the pore pressure �  −1 = − + = C r − K (23) rr 1− 1 −  1 −  r effect, the total radial stress can be written as r r where C is the integration constant that can be found using r  −1 3 r = − + C R +  P (28) rr 3 ep r=R r=R ep ep the boundary condition applied to the wellbore wall. For 1 −  1 − r r open-hole completions, the bottomhole pressure is equal to The total stress at r = R must be unique (i.e. radial stress ep the total stress applied to the boundary and hence the effec- continuity) and therefore an equality can be written as tive stress is zero. This gives 1  E e,s 2C +  −  P +  − Y 1 r=R C = ep v 3 r=R (24) 1 +  1 −  3(1 − ) −1 ep r  −1 = − + C R +  P 3 ep r=R ep When the wellbore is cased, the displacement at the well- (29) 1 −  1 − r r bore wall is assumed to be zero and thus the solution of Ideally, this equality is solved analytically to give a closed- displacement is needed. The first-order differential equation form expression of R , but this can be difficult, if not impossi- ep of displacement in Eq. 6 can be solved to give ble. Thus, a trial-and-error scheme can be used to find the root. The root of R is then admissible only if it is within the model ep ep ep 1 u = u −  g()d (25) ep r r Note that u is the displacement at the elastic–plastic ep interface and can be found assuming continuity of displace- ment; its value is equal to that found by the elastic solution at r = R . Following the approach suggested by Park and Kim ep (2006) and employed by Masoudian and Hashemi (2016), g(r) can be found, in the light of the first relationship in Eq. 6, as a function of elastic parts of radial and tangential strain by replacing the solution of stress in the plastic zone into elastic stress–strain relationships. Combining these equations and solving the integral term gives r + Eu R − r (1 + ) ep ep u = r R − C A ep 3 E 1 +   + +1 +1 R − r ep 1 + − BY − AK + E F R − F (r) r , ep , + 1 3 (26) Therefore, for a cased-hole completion, C can be found by solving this equation for u = 0 . It is evident from these equations that the solution needs the elastic–plastic interface ( R ) to be known. In order to find R , the continuity of Fig. 10 Verification of CSCS case for: a stress solution, b displace- ep ep ment solution 1 3 1580 M. S. Masoudian et al. Appendix 3: Comparative Analysis of  domain R , R ; i.e. if no root is found or it is smaller than the w 0 for Different Cases wellbore radius or larger than the reservoir radius, the reservoir remains elastic and no failure zone is developed around the It is useful to note how different boundary conditions on the wellbore. It should also be noted that the analytical solutions of wellbore affect the value of  within the elastic zone. As stress and displacement were compared to numerical solutions shown in the paper, CSCS and CSZD cases result in very simi- to verify their performance with perfectly plastic behaviour. lar profiles for  in the elastic zone as ZDCS and ZDZD, Such a verification for the CSCS case is depicted in Fig.  10. respectively. To better explain and justify this, a comparison of all the terms involved in  can be carried out.  is defined Appendix 2: Integration Constants for Elastic m m as the mean of radial and tangential stress as follows: Cases +  f (r) rr See Table 3.  = = C + (34) m 1 2 2 Table 3 Integration constants Boundary conditions (acronym) Relationships for different sets of boundary 2 2 conditions in elastic cases R −p R +F R −F R Constant stress at R ( ) ( ) (30) 0 0 w w w 0 C = 1 2 R −R Constant stress at R 0 w (CSCS) C = C −  R + F R 2 1 0 0 2 2 e,s p R −ER  ∕3−F R +F R ( ) ( ) Zero displacement at R 0 0 w (31) 0 w v,0 C =  − p + 1 0 0 2 (1−2)R +R Constant stress at R w (ZDCS) C = C −  R + F R 2 1 0 0 2 2 e,s 2 p R −ER  ∕3+F R −F R +(1−2)R  −p ( ) ( ) ( ) Constant stress at R w w 0 w 0 0 (32) 0 v,0 0 C = 1 2 (1−2)R +R Zero displacement at R w (CSZD) C = C − p R + F R 2 1 w w e,s F R −F R E Zero displacement at R ( ) ( ) (33) 0 w v,0 C =  − p + − 1 0 0 2 2 (1−2)(R −R ) 1−2 Zero displacement at R 0 w (ZDZD) E e,s C = F R − R (1 − 2) C −  + p + 2 0 1 0 0 0 v,0 Table 4 Integration constants in solution of mean horizontal stress within the elastic zone Boundary conditions (acronym) Relationships F R −F R ( ep) ( 0) CSCS and ZDCS (35) Y−f R −(−1)P R +(+1) + ( ep) ( ep) 0 F(R ) C = ≅  − 1 0 2 (+1)−(−1) R F(R )−F(R ) w 0 CSCS elastic (36) −p + 0 w R F R ( ) 0 0 C = ≅  − 1 0 1− R e,s F(R )−F(R ) 0 w ZDCS elastic (37) p −E ∕3− 0 2 v,0 R F R ( ) 0 0 C =  − p + ≅  − 1 0 0 0 2 (1−2)+1 R CSZD and ZDZD      (38) C =  Y − f R − ( − 1)P R + ( + 1) 1 ep ep (1−)−(1+)(1−2) F R −F R F R ( ep) ( 0) E e,s 1 ( ) E e,s +  − (1 − 2)  − p ≅  − p + − 2 0 0 0 0 2 R 3 v,0 1−2 R 3 v,0 0 0 F(R )−F(R ) e,s 0 w CSZD elastic   (39) p −E ∕3+ +(1−2)( −p ) w 0 0 v,0 2 R F R ( ) 1 0 E e,s C = ≅  − p + − 1 0 0 2 v,0 (1−2)+ 1−2 R 3 e,s F(R )−F(R ) F(R ) ZDZD elastic 0 w v,0 1 0 E e,s (40) C =  − p + − ≅  − p + − 1 0 0 2 0 0 2 v,0 (1−2)R (1−) 1−2 1−2 R 3 0 0 1 3 Elastic–Brittle–Plastic Behaviour of Shale Reservoirs and Its Implications on Fracture… 1581 Dong J-J, Hsu J-Y, Wu W-J, Shimamoto T, Hung J-H, Yeh E-C, Wu Since f (r) is the same for all boundary conditions, the Y-H, Sone H (2010) Stress-dependence of the permeability and differences must relate to C . On the other hand, the plastic porosity of sandstone and shale from TCDP Hole-A. Int J Rock zone is always located very close to the wellbore. Thus, both Mech Min Sci 47(7):1141–1157. http s://doi.org/10.1016 /j.ijrm R and R can be considered very small compared to R . ms.2010 .06.019 w ep 0 2 2 2 2 −6 Gale JFW, Laubach SE, Olson JE, Eichhubl P, Fall A (2014) Natural Consequently, a ratio given by  = R ∕R ≈ R ∕R ≈ 10 w 0 ep 0 fractures in shale: a review and new observations. AAPG Bull can be implemented in all equations for C and neglected. 98(11):2165–2216. http s://doi.org/10.1306 /0812 1413 151 Gensterblum Y, Ghanizadeh A, Cuss RJ, Amann-Hildenbrand A, Furthermore, F R ≪ F R and can also be neglected. ep 0 Krooss BM, Clarkson CR, Harrington JF, Zoback MD (2015) Gas Considering these aspects for all four boundary condition transport and storage capacity in shale gas reservoirs—a review. cases results in the integration constants presented in Part A: transport processes. J Unconv Oil Gas Resour 12:87–122. Table 4. http s://doi.org/10.1016 /j.juog r.2015 .08.001 Gutierrez M, Øino LE, Nygård R (2000) Stress-dependent permeability of a de-mineralised fracture in shale. Mar Petrol Geol 17(8):895– 907. http s://doi.org/10.1016 /S026 4-8172 (00)0002 7-1 References Han G, Dusseault MB (2003) Description of fluid flow around a well- bore with stress-dependent porosity and permeability. J Petrol Sci Hull KL, Abousleiman, YN, Han Y, Al-Muntasheri GA, Hosemann P, Eng 40(1–2):1–16. http s://doi.org/10.1016 /S092 0-4105 (03)0004 Parker SS, Howard CB (2015) New insights on the mechanical 7-0 characterization of Kerogen-rich shale, KRS. In: Abu Dhabi inter- Hill DG, Nelson CR (2000) Gas productive fractured shales: an over- national petroleum exhibition and conference. Society of Petro- view and update. Gas TIPS 6(2):4–13 leum Engineers, Abu Dhabi, UAE, p 25. http://dx.doi.or g/10.2118 Huang X, Bandilla KW, Celia MA (2016) Multi-physics pore-network /1776 28-MS modeling of two-phase shale matrix flows. Transp Porous Media Alexander T, Baihly J, Boyer C, Clark B, Waters G, Jochen V, Le 111(1):123–141. http s://doi.org/10.1007 /s112 42-015-0584 -8 Calvez J, Lewis R, Miller CK, Thaeler J, Toelle BE (2011) Shale Islam MA, Skalle P (2013) Experimentally evaluating shale dilation gas revolution. Oilfield Rev 23(3):40–57 behavior. In: Chatellier J-Y, Jarvie DM (eds) Critical assessment Alonso E, Alejano LR, Varas F, Fdez-Manin G, Carranza-Torres of shale resource plays. American Association of Petroleum Geol- C (2003) Ground response curves for rock masses exhibiting ogists, Tulsa, p 205 strain-softening behaviour. Int J Numer Anal Meth Geomech Kumar H, Elsworth D, Mathews JP, Marone C (2016) Permeabil- 27(13):1153–1185. http s://doi.org/10.1002 /nag.315 ity evolution in sorbing media: analogies between organic-rich Amann F, Button EA, Evans KF, Gischig VS, Blümel M (2011) Experi- shale and coal. Geofluids 16(1):43–55. http s://doi.org/10.1111 mental study of the brittle behavior of clay shale in rapid uncon- /gfl.1213 5 fined compression. Rock Mech Rock Eng 44(4):415–430. https:// Li M, Yin G, Xu J, Cao J, Song Z (2016) Permeability evolution of doi.org/10.1007 /s006 03-011-0156 -3 shale under anisotropic true triaxial stress conditions. Int J Coal Amann-Hildenbrand A, Ghanizadeh A, Krooss BM (2012) Trans- Geol 165:142–148. http s://doi.org/10.1016 /j.coal .2016 .08.017 port properties of unconventional gas systems. Mar Petrol Geol Liang C, Chen M, Jin Y, Lu Y (2014) Wellbore stability model for 31(1):90–99. http s://doi.org/10.1016 /j.marp etge o.2011 .11.009 shale gas reservoir considering the coupling of multi-weakness Bustin RM, Bustin AMM, Cui A, Ross D, Pathi VM (2008) Impact of planes and porous flow. J Nat Gas Sci Eng 21:364–378. https:// shale properties on pore structure and storage characteristics. Soc doi.org/10.1016 /j.jngs e.2014 .08.025 Petrol Eng. http s://doi.org/10.2118 /1198 92-MS Lubliner J (1990) Plasticity theory. Macmillan, New York, p 50. Carey JW, Lei Z, Rougier E, Mori H, Viswanathan H (2015) Fract-per- ISBN 0486462900 meability behavior of shale. J Unconv Oil Gas Resour 11:27–43. Lyu Q, Ranjith PG, Long X, Kang Y, Huang M (2015) A review http s://doi.org/10.1016 /j.juog r.2015 .04.003 of shale swelling by water adsorption. J Nat Gas Sci Eng Part Carroll SA, McNab WW, Torres SC (2011) Experimental study of 27(3):1421–1431. http s://doi.org/10.1016 /j.jngs e.2015 .10.004 cement—sandstone/shale—brine—CO interactions. Geochem Ma T, Chen P (2015) A wellbore stability analysis model with chem- Trans 12(1):9. http s://doi.org/10.1186 /1467 -4866 -12-9 ical-mechanical coupling for shale gas reservoirs. J Nat Gas Chen D, Pan Z, Ye Z (2015a) Dependence of gas shale fracture perme- Sci Eng 26:72–98. https://doi.or g/10.1016/j.jngs e.2015 .05.028 ability on ee ff ctive stress and reservoir pressure: model match and Masoudian MS, Hashemi MA (2016) Analytical solution of a cir- insights. Fuel 139:383–392. http s://doi.org/10.1016 /j.fuel .2014 . cular opening in an axisymmetric elastic-brittle-plastic swell- 09.018 ing rock. J Nat Gas Sci Eng Part A 35:483–496. http s://doi. Chen T, Feng X-T, Pan Z (2015b) Experimental study of swelling of org/10.1016 /j.jngs e.2016 .08.076 organic rich shale in methane. Int J Coal Geol 150–151:64–73. Masoudian MS, Airey DW, El-Zein A (2016a) The role of coal seam http s://doi.org/10.1016 /j.coal .2015 .08.001 properties on coupled processes during C O sequestration: a Cho Y, Ozkan E, Apaydin OG (2013) Pressure-dependent natural- parametric study. Greenh Gases Sci Technol (in Press). http :// fracture permeability in shale and its effect on shale-gas well dx.doi.org/10.1002 /ghg.1575 production. SPE Reserv Eval Eng 16(2):216–228. h t t p s : / / d o i . Masoudian MS, El-Zein A, Airey DW (2016b) Modelling stress and org/10.2118 /1598 01-PA strain in coal seams during CO injection incorporating the Cui X, Bustin RM, Chikatamarla L (2007) Adsorption-induced coal rock–fluid interactions. Comput Geotech 76:51–60. https://doi. swelling and stress: Implications for methane production and acid org/10.1016 /j.comp geo.2016 .02.010 gas sequestration into coal seams. J Geophys Res Solid Earth McGlade C, Speirs J, Sorrell S (2013) Unconventional gas—a review 112:B10. http s://doi.org/10.1029 /2004 jb00 3482 of regional and global resource estimates. Energy 55:571–584. Curtis JB (2002) Fractured shale-gas systems. AAPG Bull http s://doi.org/10.1016 /j.ener gy.2013 .01.048 86(11):1921–1938. http s://doi.or g/10.1306 /61EE DDBE -173E Mehmani A, Prodanović M, Javadpour F (2013) Multiscale, mul- -11D7 -8645 0001 02C1 865D tiphysics network modeling of shale matrix gas flows. Transp 1 3 1582 M. S. Masoudian et al. Porous Media 99(2):377–390. http s://doi.or g/10.1007 /s112 Song F, Warpinski N, Toksöz M (2014) Full-waveform based micro- 42-013-0191 -5 seismic source mechanism studies in the Barnett Shale: link- Palmer I (2009) Permeability changes in coal: analytical modeling. ing microseismicity to reservoir geomechanics. Geophysics Int J Coal Geol 77(1–2):119–126. https://doi.or g/10.1016/j.coal 79(2):KS13–KS30. http s://doi.org/10.1190 /geo2 013-0094 .1 .2008 .09.006 Souley M, Homand F, Pepa S, Hoxha D (2001) Damage-induced Pan Z, Connell LD (2012) Modelling permeability for coal reservoirs: permeability changes in granite: a case example at the URL in a review of analytical models and testing data. Int J Coal Geol Canada. Int J Rock Mech Min Sci 38(2):297–310. http s://doi. 92:1–44. http s://doi.org/10.1016 /j.coal .2011 .12.009org/10.1016 /S136 5-1609 (01)0000 2-8 Park K-H, Kim Y-J (2006) Analytical solution for a circular open- Spencer CW (1989) Review of characteristics of low-permeability gas ing in an elastic–brittle–plastic rock. Int J Rock Mech Min Sci reservoirs in Western United States. AAPG Bull 73(5):613–629 43(4):616–622. http s://doi.org/10.1016 /j.ijrm ms.2005 .11.004 Varas F, Alonso E, Alejano LR, Manín GF (2005) Study of bifurca- Pathi VSM (2008) Factors affecting the permeability of gas shales. PhD tion in the problem of unloading a circular excavation in a strain- Thesis, The University of British Columbia, Vancouver, Canada softening material. Tunn Undergr Space Technol 20(4):311–322. Rutqvist J (2015) Fractured rock stress-permeability relationships from http s://doi.org/10.1016 /j.tust .2004 .12.003 in situ data and effects of temperature and chemical-mechanical Wang FP, Reed RM (2009) Pore networks and fluid flow in gas shales. couplings. Geofluids 15(1–2):48–66. http s://doi.or g/10.1111 / In: SPE annual technical conference and exhibition. Society of gfl.1208 9 Petroleum Engineers, New Orleans, p 3. http://dx.doi.or g/10.2118 Rutqvist J, Börgesson L, Chijimatsu M, Hernelind J, Jing L, Kobayashi /1242 53-MS A, Nguyen S (2009) Modeling of damage, permeability changes Yuan W, Pan Z, Li X, Yang Y, Zhao C, Connell LD, Li S, He J (2014) and pressure responses during excavation of the TSX tunnel in Experimental study and modelling of methane adsorption and dif- granitic rock at URL, Canada. Environ Geol 57(6):1263–1274. fusion in shale. Fuel Part A 117:509–519. https://doi.or g/10.1016 http s://doi.org/10.1007 /s002 54-008-1515 -6/j.fuel .2013 .09.046 Rutqvist J, Rinaldi AP, Cappa F, Moridis GJ (2013) Modeling of fault Zang J, Wang K, Zhao Y (2015) Evaluation of gas sorption-induced reactivation and induced seismicity during hydraulic fracturing internal swelling in coal. Fuel 143(Supplement C):165–172. http of shale-gas reservoirs. J Petrol Sci Eng 107:31–44. http s://doi.s://doi.org/10.1016 /j.fuel .2014 .11.007 org/10.1016 /j.petr ol.2013 .04.023 Zhou T, Zhang S, Feng Y, Shuai Y, Zou Y, Li N (2016) Experimental Shi JQ, Durucan S (2005) A model for changes in coalbed permeability study of permeability characteristics for the cemented natural frac- during primary and enhanced methane recovery. SPE Reserv Eval tures of the shale gas formation. J Nat Gas Sci Eng 29:345–354. Eng. http s://doi.org/10.2118 /8723 0-PAhttp s://doi.org/10.1016 /j.jngs e.2016 .01.005 Soeder DJ (1988) Porosity and permeability of Eastern Devonian gas shale. SPE Form Eval 3(1):116–124. https://doi.or g/10.2118/1521 3-PA Sone H, Zoback MD (2013) Mechanical properties of shale-gas reser- voir rocks—part 2: ductile creep, brittle strength, and their relation to the elastic modulus. Geophysics 78(5):D393–D402. https ://doi. org/10.1190 /geo2 013-0051 .1 1 3

Journal

Rock Mechanics and Rock Engineering – Springer Journals

Published: Jan 19, 2018

There are no references for this article.