Triaxial Deformation of the Goldwyer Gas Shale at In Situ Stress Conditions—Part I: Anisotropy of Elastic and Mechanical Properties

The evolution of shale’s mechanical properties with confining pressure, temperature, and mineral composition directly influences fracture closure besides the effect of in situ stress variation across lithologies. We are the first to perform experimental study to characterize the mechanical properties of the Goldwyer gas shale formation located in the Canning Basin, Western Australia. We have performed constant strain rate multistage triaxial tests at in situ stress condition (confining pressure ≤ 22 MPa) on 15 samples of the Goldwyer gas shales with variable minerology, organic content, and heterogeneity. Deformation tests were conducted at room temperature and in drained conditions on cylindrical samples cored parallel (horizontal) and perpendicular (vertical) to the bedding plane. Both triaxial compressive strength (σTCS) and static young’s modulus E show a strong sensitivity to confining pressure and mineralogy, while only E shows a directional dependency, i.e., Eh > Ev. The internal friction coefficient µi in a plane parallel to the bedding is 0.72 ± 0.12, while it is only 0.58 ± 0.17 in the orthogonal direction. Both σTCS and E are significantly lower when larger fractions of weak mineral constituents are present (clays or organic matter). We observe that the Young’s modulus of most vertical samples is best approximated by Reuss’s bound, whereas that of horizontal samples is best approximated by Hill’s average of Voigt and Reuss bounds. The most prospective G-III unit of the Goldwyer shale formation (depth > 1510 m) is semi-brittle to brittle, making it suitable for future development. Global empirical correlations of unconventional shale’s elastic (Young’s modulus) and mechanical properties (compressive rock strength, internal friction coefficient) are established for continuous property prediction in the field Mechanically weak phase fraction ClayTocPHI which is consisting of clay, porosity, and total organic content controls mechanical strength properties of highly heterogeneous gas shale reservoir rocks. A consistent distinction between the first loading and unloading/reloading static Young’s modulus is observed irrespective of whether the applied stress is above or below the in situ effective stress. We suggest acquiring both parameters to better constrain static reservoir deformation behaviour during field development. Global empirical correlations of unconventional shale’s elastic (Young’s modulus) and mechanical properties (compressive rock strength, internal friction coefficient) are established for continuous property prediction in the field Mechanically weak phase fraction ClayTocPHI which is consisting of clay, porosity, and total organic content controls mechanical strength properties of highly heterogeneous gas shale reservoir rocks. A consistent distinction between the first loading and unloading/reloading static Young’s modulus is observed irrespective of whether the applied stress is above or below the in situ effective stress. We suggest acquiring both parameters to better constrain static reservoir deformation behaviour during field development.

• Global empirical correlations of unconventional shale's elastic (Young's modulus) and mechanical properties (compressive rock strength, internal friction coefficient) are established for continuous property prediction in the field • Mechanically weak phase fraction ClayTocPHI which is consisting of clay, porosity, and total organic content controls mechanical strength properties of highly heterogeneous gas shale reservoir rocks. • A consistent distinction between the first loading and unloading/reloading static Young's modulus is observed irrespective of whether the applied stress is above or below the in situ effective stress. We suggest acquiring both parameters to better constrain static reservoir deformation behaviour during field development.

Introduction
Significant progress has been made over the past decade to extract gas from ultra-low permeable unconventional shale reservoirs (Guo et al. 2015;Jarvie et al. 2007;Rezaee 2015;Sone 2012), a cleaner and readily available energy resource to accelerate the transition towards a lower carbon economy. Rapid progress in technology, notably hydraulic fracturing, and horizontal drilling, made shale gas economically producible (Herrmann et al. 2018;Rybacki et al. 2017;Singh et al. 2019;Sone and Zoback 2013a;Yang et al. 2015;Yang and Zoback 2014). The success of North America's massive shale gas revolution allows other countries like China, Australia, Argentina, or Poland to seriously consider extraction of natural gas from the subsurface where economic and recoverable reserves are identified (EIA 2013). The petrophysical and mechanical properties of these gas shale formations are, therefore, in high demand to assess the potential producing zones. However, the generalization of North America's organic-rich shales to other shale formations worldwide is not straightforward due to significant differences in minerology, geomechanical properties, fracability, and the presence of 100-1000 s of parasequence units within the reservoir. Over the past decade, Sone (2012), Sone and Zoback (2013a), Villamor Lora et al. (2016), Delle , Rybacki et al. (2015), Herrmann et al. (2018), Yang et al. (2015) and Guo et al. (2015) conducted very detailed deformation experiments on low permeable unconventional shale samples from different reservoirs worldwide. These experiments were complemented by a study of the micro-structure, orientation of bedding plane, mineral composition, organic matter, depositional history, in-situ stress state, brittleness, and creep behaviour. Although hydraulic fracturing (HF) reopens and/or creates fractures at various scales (Norris et al. 2016), and proppant is placed to keep them open, the rapid decline in production observed in such fractured reservoirs is mostly attributed to progressive (time-dependent) fracture closure (Al-Rbeawi 2018; Wang 2016). Fracture closure is controlled mainly by in-situ conditions such as confining pressure, temperature, stress orientation Islam and Skalle 2013;Masri et al. 2014;Niandou et al. 1997;Rybacki et al. 2015Rybacki et al. , 2016Sone and Zoback 2013a, b;Villamor Lora et al. 2016), compositional stress layering (Ma and Zoback 2017;Mandal et al. 2021;Singh et al. 2019;Xu et al. 2019;Yang et al. 2015), time-dependent deformation (Herrmann et al. 2020;Rassouli and Zoback 2018;Rybacki et al. 2017;Sone and Zoback 2014;Xu et al. 2019), and the mechanical/petrophysical properties of the shale, e.g., porosity, minerology, brittleness, friction (Cerasi et al. 2017;Kohli and Zoback 2013;Rybacki et al. 2015Rybacki et al. , 2016Sone and Zoback 2013a, b). Therefore, the detailed mechanical characterization of gas shales is essential, not only for the selection of favourable intervals for hydraulic fracturing stimulation, but also to better predict and mitigate post-stimulation fracture closure.
Here we report the mechanical properties (compressive strength, Young's modulus) of dry samples of the Goldwyer shale formation, onshore Canning Basin in Western Australia, with varying mineralogy from two distinct stratigraphic units (G-I, and G-III). The reported data are derived from triaxial deformation experiments conducted at constant strain rate, at room temperature, in drained conditions, and at different confining pressures up to 22 MPa, including the insitu mean stress value. We also establish empirical relations between mechanical properties (compressive strength, internal friction coefficient, peak axial strain at failure) with rock composition and static Young's modulus, which presents the advantage of being directly applicable by field engineers to build continuous depth profiles of rock properties from wireline logs. Further, we quantitatively relate the anisotropy of the static Young's modulus (horizontal-to-vertical ratio E h /E v ) to the fraction of the mechanically weaker phase, i.e., the combination of clay minerals, organic matter, and porosity. In addition, effective medium theory is used to estimate the theoretical bounds of the elastic modulus from their weak phase fraction (Voight iso-stress and Reuss iso-strain). The deformation response of the Goldwyer shale (brittle or ductile in terms of brittleness index) is investigated and correlated with axial strain at failure. Moreover, these geomechanical data support the estimation of the in-situ horizontal principal stress and building continuous mechanical strength profile and, therefore, assist in the determination of the most likely faulting regime, the failure envelope, identification of best operational zone at depth and the design and optimization of hydraulic fracture operations in the field.
Note that time-dependent deformation under constant stress (creep), frictional failure characteristics and viscoelastic stress relaxation modelling are necessary to predict the long-term geomechanical implications of fracture closure and lithological variation of the least principal stress magnitude. To limit the scope of this manuscript, these aspects will be reported in a subsequent paper.

Materials and Methods
The rock samples characterised and tested in the laboratory originate from the middle Ordovician Goldwyer formation (~ 400 m thick within the interval 1188-1592 m), crossed by the vertical exploration well Theia-1. This well is located approximately 155 km southwest of the Broome platform. The vertical trajectory of the well cuts across the sub-horizontal to nearly horizontal bedding of the Goldwyer formation. The investigated cores (63 mm in diameter) cover two stratigraphic units G-I (sample Th1), and G-III (samples Th2 to Th10). Among them black shale samples are mature for hydrocarbon generation. A preliminary naked-eye visualization of the cores suggests that the G-III unit is dominated by black mudstones with the presence of thin calcareous laminations, concretions and thin veins filled with a white material. In contrast, the samples from the G-I units are a mixture of clay-and carbonate-rich material with interbedded mudstone (Fig. 1). Earlier work on the Goldwyer formation indicated that this shale is strongly laminated, and the presence of fractures partially filled with density-contrasting minerals was also noted , which is confirmed by our observations.
Powdered sample groups were subjected to X-ray diffraction (XRD) and geochemical analysis (Rock-Eval pyrolysis) to determine mineralogy and total organic carbon (TOC in wt%). Porosity is evaluated from the dry bulk density and average grain density (Helium Pycnometer on crushed samples oven-dried at 95 °C for 24 h). The remainder of this section will detail the sample selection strategy, the minerology, and the micro-structure of the selected samples.

Sample Selection Strategy
The well was drilled in 2015 and the recovered cores were stored unpreserved (room-dry) in the core library at the West Australian Department of Mines, Industry, Regulation and Safety (DMIRS 2019) since. The recovered cores from the Theia-1 well within the Goldwyer formation interval were cut in 1 m-long sections. Samples were exposed to room humidity conditions at the core library (~ 50%), until they were shipped to us for triaxial testing. Therefore, we tested all the samples in "as-received" conditions. It is expected that the cores were in equilibrium with the atmospheric conditions within the core library so that free water must have evaporated from the shale, whereas clay-bound water should not have. The water saturation of the shale as estimated form drying it at 105 °C is about 30-40%. The selection of samples for this geomechanical study was conducted at the WA DMIRS Core Library, where the cores were exposed for inspection. This involved (i) a physical inspection of the core trays to assess the mechanical integrity of the various core sections available, which is necessary to allow sampling for further characterisation and testing in the laboratory; and (ii) a detailed visual inspection of the various lithologies present and their heterogeneity, i.e., colour variations, sedimentary features like calcite veins, limestone nodules, laminations, fractures, etc.
Thereafter, ten target depths were selected to extract both horizontal and vertical plugs (a total of 15 samples, see Table 1) to allow for an analysis of (i) the static elastic properties accounting for anisotropy, (ii) the stress orientation and magnitude, and (iii) the brittleness index (Mandal et al. 2020c). Two criteria were used for the selection: (i) Heterogeneity: from simple to complex, i.e., from homogeneous carbonate-rich shale to clayey and heterogeneous shale (ii) Organic richness: from light grey to black organicrich mudstone.
Core sections approximately 120 mm in length were cut off the main cores at the ten target depths. These sections were then imaged with a medical X-ray CT to generate 3D density maps (voxel size = 100 μm). These maps were used to guide where to extract the samples to be further characterised and tested in the laboratory, i.e., identify laminations, nodules, and other heterogeneities, while avoiding incipient/ open fractures and increasing the success rate of the sampling program. In fact, the same criteria were used for selecting the location of these samples as for selecting the parent core sections (see above). As an example, Fig. 2 shows two mutually orthogonal slice images of a vertical core section extracted from its 3D density map. One can see a vertical calcite vein in the XY view, and few horizontal contrasting composition layers in the YZ view.
A total of 15 cylindrical samples with a length-to-diameter ratio of about 2 were extracted and their end-faces Fig. 1 Core sections from the Theia-1 well selected for this laboratory geomechanical study. Sample depths are indicated in the Gamma ray (GR) log with red triangles along the Goldwyer formation. One core section is selected from the G-I unit, and the other nine from the G-III unit. The representative set of samples selected covers the range of lithologies visually observed on the cores, i.e., qualitative mineralogy and organic content (based on colour), sedimentary features and heterogeneities such as laminations, veins, nodules, fractures, and bedding orientation (vertical and horizontal plugs). (Colour figure online)  Table 1).

Mineralogy, Density, Porosity, and Organic Matter
The X-Ray Diffraction (XRD) method is used to determine the mineral composition in weight percentage (wt%) for representative powdered samples from the ten target depths. The dry bulk density of the samples ranges between 2.44 and 2.67 g/cm 3 , and the total porosity φ ranges between 4 and 13%. Following Jarvie's approach, thermal maturity is calculated and expressed in terms of vitrinite reflectance from T max (Jarvie et al. 2007), and ranges between 0.22 and 0.96% VR 0 , with an average of 0.71%. The Total Organic Content (TOC) exhibits a broad distribution from 0.1 to 8 vol%. We have reported here three constituent groups with contrasting mechanical properties: (i) QFP refers to the strong phases comprising Quartz, Feldspar, and Pyrite; (ii) Cb refers to the intermediately strong phases comprising Calcite and Dolomite; and (iii) ClayTocPHI refers to the weak phases comprising all clay minerals, organic matter, and porosity. Note that the Fig. 2 Two mutually orthogonal slice images of a vertical core section extracted from the 3D density map acquired with a Medical X-ray CT (voxel size = 100 microns): a a horizontal slice image along the XY plane; and b a vertical slice image along the YZ plane calculated total clay phase includes mica since clay minerals are generally formed from Mica (Churchman 1980). Among the investigated samples, Th2, Th4, and Th8 are considered predominantly calcareous, whereas the other samples are clay-dominated (see Fig. 3 or Table 2). Petrographic and depositional analysis reveal the presence of five distinct types of lithofacies, namely thinly laminated, black, heterolithic, calcareous, and banded concretionary mudstone. This analysis shows the degree of heterogeneity of the Goldwyer formation and hints at the challenges associated with the characterisation of the formation's mechanical behaviour. Cross-plots of the volume fractions of TOC, clay, and porosity are presented in Fig. 4. We observe a weak correlation between them. The qualitative increase in porosity with either clay or TOC in the Goldwyer formation is consistent with the presence of organic, intraparticle, and interparticle pores (Iqbal et al. 2021). It has also been reported for other unconventional shales that pores tend to reside in solid organic components and/or in between clay platelets (Curtis et al. 2010). However, the relatively weak 1-to-1 correlations between TOC, clay, and porosity also suggest that these three variables are reasonably independent. Note that the newly introduced volumetric fraction of the weak phase ClayTocPHI is defined as the sum of the volume fractions of clay, TOC, and porosity.

Micro-structure
The mechanical properties of shales are largely controlled by their microstructure and mineralogy Dewhurst et al. 2015;Josh et al. 2012;Rybacki et al. 2015;Sarout and Guéguen 2008a, b). Due to their fine-grained and heterogeneous nature, gas shales require a combination of visualization techniques in order to characterise their fabric and texture at various scales, up to the Fig. 3 Ternary diagram representing the mineralogy at the ten selected depths in the Theia-1 well (numbered as Th1 to Th10). Each data point is colour-coded according to its TOC. Calcareous (Th2, 4) and clay-dominated shale (Th6, 9, 10) samples are highlighted with a light blue and a red circle, respectively  Organic-rich laminated mudstone with abundant laminar clay lay-ers and organic matter; c Argillaceous mudstone with abundant calcitic matrix and organic-rich detrital clay. Qtz, Py, and Cal stand for quartz, pyrite, and calcite, respectively (Colour figure online). Modified from unpublished report (DMIRS 2019) illite and mica-dominated detrital clays (< 5 µm), calcitic matrix (CM), thin lens to angular shaped organic matter, and elongated, lenticular to rounded clasts of quartz, feldspar, and calcite (10-100 µm). Figure 5b and c displays clay-and organic-rich, and carbonate-dominated samples, respectively. The textural fabric development associated with the preferred orientation of clay platelets, elongated organic matter, and detrital fossils depicting the bedding plane are clearly visible in Fig. 5b. The mineralogy of the samples combined with thin section imagery suggests that fabric development becomes pronounced when the clay fraction reaches approximately 40 vol%. With decreasing clay fraction, the clastic fraction becomes the main loadbearing structure, i.e., calcitic matrix [see Fig. 5a, and (Kohli and Zoback 2013;Sone and Zoback 2013a;Villamor Lora et al. 2016)]. The distribution of organic matter appears discontinuous, which is consistent with previous reports [organic-rich Bakken shale in Vernik and Nur (1992)]. The pore space remains essentially unresolved with the available imaging techniques, i.e., nanometer scale in the Goldwyer formation according to Yuan et al. 2018.

Triaxial Deformation Experiments
Because the available cores were not preserved since recovery at the well site (room-dry), and because of the prohibitive cost and duration of geomechanical tests with shales under controlled pore pressure conditions, the selected samples were tested dry. Considering that this is a gas shale, it is not likely to be fully water-saturated in situ, in contrast with sealing shales/caprocks, which are often fully water-saturated and require a full saturation of the samples when testing in the laboratory. In this context, the pore pressure during triaxial deformation is expected to remain negligible (compressible air), and the total stress applied in the laboratory corresponds to the effective stress estimated in the field. The samples are subjected to the in situ stress prior to loading them with a deviatoric stress, with the aim to minimise the impact of desiccation damage by closing potential micro-cracks. The triaxiality of the stress field was reproduced in the laboratory within the limitations associated with a conventional axi-symmetric triaxial stress vessel, i.e., (i) the two radial stresses must be equal, i.e., S hmin = S Hmax = S h ; and (ii) the vertical stress must be greater than or equal to the horizontal stress S v ≥ S h .

Fig. 6
Representative backscatter SEM images of the Goldwyer formation in the investigated depth interval (1188.5-1593 m, approx. 400 m). a Argillaceous mudstone. The red filled black arrows point to detrital clay flakes, while the yellow arrows highlight micropores/ nanopores. Distribution of detrital clay is highlighted by red dashed contour. Some of the micropores and nanopores are filled with organic matter. b Organic-rich laminated mudstone. Detrital clay and organic matter are outlined with a white dashed line. Clay is the controlling grain framework. Qtz, Fes, Mi, OM, Py, Cal, Dol stand for quartz, feldspar, mica, organic matter, pyrite, calcite, and dolomite, respectively. Modified from unpublished report (DMIRS 2019) Using a combination of well logs (including cross-dipole sonic) and laboratory triaxial deformation data, while accounting for shale anisotropy, Mandal et al. (2020c) estimated the three principal stress magnitudes (vertical stress S v , minimum horizontal stress S hmin , and maximum horizontal stress S Hmax ) and pore pressure (P p ) along the Theia-1 well. This analysis suggests a hybrid stress state, i.e., normal-to-strike-slip faulting in the deeper G-III unit, whereas in the shallower G-I unit a dominant strike-slip faulting stress regime is estimated. Predicted stress anisotropy has been attenuated due to the observed overpressure in these intervals. In other words, the effective stress field (external stress minus pore pressure) is closer to isotropy than the "external" stress field alone (S v , S hmin , S Hmax ).
Therefore, for triaxial testing in the laboratory, the selected dry shale samples (large vertical or smaller horizontal samples) are subjected to a confining pressure commensurate with the mean effective stress usual at the sample's recovery depth point (p c = (S v + S hmin + S Hmax )/3 − P p ). In fact, due to the scarcity of shale samples and their variability, each sample is subjected to a so-called multistage triaxial test MST, (Kim and Ko 1979;Youn and Tonon 2010), i.e., triaxial loading at constant confining pressure, repeated at multiple confining pressures, including the estimated in situ mean effective stress. In practice, five confining pressure stages were defined with respect to the estimated mean effective stress p c at depth: 0.25 × p c , 0.5 × p c , 0.75 × p c , 1 × p c , and 1.25 × p c (Fig. 7). The MST allows for the determination of multiple mechanical properties for each sample separately, i.e., Young's modulus and Poisson's ratio at each confining pressure (or depth), and Mohr-Coulomb failure envelope spanning the interval of confining pressures tested.
In practice, CSIRO's autonomous triaxial cell (ATC) with axial and radial ultrasonic P-wave transducers was used for the tests on the larger vertical plugs (Sarout et al. 2014), whereas for the smaller horizontal plugs only axial ultrasonic P-wave transducers were available (Fig. 7). The ATC consists of a high-stiffness load cell, a pressure vessel and independent stepping motor pumps for independently controlling confining pressure, pore pressure and axial load up to 70 MPa, 70 MPa and 400 MPa (for 38 mm diameter specimen), respectively. The sample assembly includes the following: i. The cylindrical sample is set between a top and a bottom steel platens; ii. A Viton sleeve isolates the sample from the confining oil, with provision for multiple radial P-wave transducers to be attached to the large vertical samples for acquiring multi-directional ultrasonic wave velocity data; iii. The two steel platens are fitted with ultrasonic P and S-wave transducers (with provision for pore fluid injection and pore pressure control, not used here); iv. Two linear differential displacement transducers (LVDTs) attached to the platens at diametrically opposite positions to measure the average axial displacement; v. A C-shaped cantilever radial strain transducer (for large vertical samples), or two mutually orthogonal radial strain transducers (for small horizontal samples); vi. An axial load cell located underneath the bottom platen.
The stability of the target confining and pore pressures can be maintained within ± 5 kPa over a duration of few months. Temperature of the test is maintained at 22 ± 0.5 °C.
Previous researchers focusing on triaxial testing of shales used various axial strain rates, ranging between 5 × 10 -6 and 10 -4 s −1 . No significant impact of the strain rate on the mechanical behaviour of the shale was observed in this range (Herrmann et al. 2018;Rybacki et al. 2017;Sarout et al. 2017). Here, we used a rate of 5 × 10 -7 s −1 during the first four stages of each MST, and a slightly faster rate of 1 × 10 -5 s −1 during the final stage (up to failure). A faster rate was selected for the last stage to promote micro-seismic activity prior and during failure. Experimental procedures are summarised in Fig. 7b and Table 3.
Note that during the first four stages of each MST, the axial loading is effectively stopped before the actual failure of the rock once the axial stress-strain curve departs from linearity (slope change by 50%). Only at the end of the fifth and last stage of the MST, the sample is loaded until actual mechanical failure. Consequently, the first four Mohr circles are expected to be smaller than if a separate (but identical) sample was brought to failure at each of these stages. The fifth Mohr circle is more representative of the actual failure of the sample, although by the time the fifth loading stage is carried out, the sample might have accumulated some mechanical damage during the first four stages. The Mohr-Coulomb failure parameters, cohesion, and friction angle are estimated from the first four stages of the MST only to (i) avoid the influence of the different strain rate used during the fifth and final stage of the MST (actual rock failure), and (ii) account for the difference between the first four stages (proxy to mechanical failure), and the fifth stage (actual failure). This approach leads to an overestimation of the internal friction coefficient (by 0.02-0.03), and an underestimation of the cohesion (by less than ~ 1 MPa), leading to an error in the determination of the strength of the rock of about 5 to 10% (Dewhurst et al. 2011;Dewhurst et al. 2015). On the other hand, accounting for measurement errors (error propagation method), the uncertainty associated with the determination of the mechanical strength, the strain at failure and Young's modulus is estimated to be approximately ± 2%, ± 3%, ± 2%, respectively. The maximum strength error combining both uncertainty sources vary between 2 and 12%. The (tangent) Young's modulus E and Poisson's ratio ν are estimated from the axial stress-strain and axial-radial strain curves in the range 40-60% of the differential peak stress value.

Results
Multistage triaxial deformation experiments were conducted at room temperature on 15 vertical and horizontal shale samples from the Goldwyer formation. The first four stages were loaded at a constant strain rates of έ = 5 × 10 −7 s −1 , over a range of confining pressures comprised between 0.25 × p c and 1 × p c , where p c is the estimated in situ mean stress. A

Mechanical and Elastic Properties
For each MST experiment, the post-failure residual stress evolution (fifth stage) was recorded until the axial strain reached ~ 4%. The stress-strain curves recorded for the vertical and horizontal samples during this final loading-to-failure stage are displayed in Fig. 8. The triaxial compressive strength σ TCS , the static Young's modulus E, and the axial strain at failure ε max from the deformation experiment at the last loading-to-failure stage (at 1.25 × p c ) are presented in Table 4 for the Goldwyer shale formation. In a similar way the static Young's modulus and Poisson's ratio have been calculated during loading and unloading/reloading cycles at each stage of the MST. The recorded stress-strain curves at the failure stage typically show relatively small amounts of inelastic deformation until failure, which occurs at stress values beyond 100 MPa, and reflect the brittleness of this shale formation (see Fig. 8). Samples deformed parallel to the bedding generally reached failure at a slightly higher axial stress compared to the vertical samples by 5-10 MPa, except few samples which showed no difference. Static Young's modulus of the horizontal shale plug is 1.5-to 4-times larger than that of the vertical plugs (Table 4). The axial strain value at failure ranges between 0.5 and 1%, with larger values displayed by the vertical samples (Fig. 8c).
All tested samples exhibit a semi-brittle to brittle response. For brittle deformation, the rock experiences small amounts of inelastic strain prior to an abrupt failure, while semi-brittle rocks exhibit a more non-linear hardening prior to the peak stress, followed by a progressive strain weakening toward a stable residual stress. The organic-rich and clay-dominated samples with higher porosity (Th3, Th6, Th9, Th10) showed a semi-brittle failure, while the carbonate-dominated samples (Th1, Th2) displayed a more brittle failure. Overall, samples with a high clay content and TOC are mechanically weaker (σ TCS ~ 100-130 MPa), whereas calcareous samples are stronger (σ TCS ~ 165-205 MPa) (see Table 4 and Fig. 8). In addition, the volumetric strain is larger for the organic-rich vertical samples (Fig. 9a) than for their horizontal counterparts (Fig. 9b). Relatively small amounts of volumetric strain are accumulated by calcareous samples prior to failure.

Rock Strength Evolution with Confining Pressure (Depth)
To understand the role of confining pressure (depth) on the mechanical properties of the Goldwyer shale, multi-stage triaxial deformation experiments are performed at a constant strain rate έ = 5 × 10 −7 s −1 for the first four stages, covering a confining pressure range of 4-17 MPa. The stress-strain curves reported in Fig. 10 show that at low confining pressure all samples exhibit a non-linear axial strain response at the early stages of the loading; a more linear response is recorded in this early stage for higher confining pressures. A higher peak stress value is generally observed at higher confining pressures (see Fig. 10a, b). Confining pressure also controls the deformation regime of the organic-rich samples: transition from a brittle failure at low confining pressure to a mixed brittle/ductile yield at higher confining pressure (see Figs. 8 and 10). However, the axial strain achieved at peak stress remains similar for all confining pressure stages after the first stage (at 4 MPa). Figure 10c and d shows that with elevated confining pressure, the static Young's modulus exhibits an incremental increase by 5-30% while random dependency observed for static Poisson's ratio measurement (Fig. 10e, f)  Note that a clear shear failure plane is observed for most of the samples after multi-stage triaxial testing. A more detailed analysis of the failure modes, shear plane angles and associated structural changes is out of the scope of the present article and will be reported in a subsequent publication.

Discussion
As extensively demonstrated in the literature for various shales, triaxial compressive strength, static Young's modulus, internal friction coefficient, cohesion, and the axial strain at failure are primarily controlled by mineralogy, organic content, porosity, bedding orientation, microstructure, and the test conditions such as confining pressure, temperature and strain rate Herrmann et al. 2018;Rybacki et al. 2015;Sarout and Guéguen 2008a, b;Sone and Zoback 2013a;Villamor Lora et al. 2016). We analyse here the influence of mineralogy, confining pressure, and microstructure on the mechanical and elastic properties of the Goldwyer shale formation. A comparison of the Goldwyer shale to USA's unconventional shales is also carried out (data from vertical samples of the Barnett, Haynesville, Eagle ford and Fort St Fig. 9 Volumetric strain recorded during the (last) failure stage of each multi-stage triaxial test conducted on vertical (a) and horizontal (b) samples. Continuous and dashed lines represent vertical and horizontal plugs, respectively. Clay-rich and organic-rich samples (Th3, 6, 9) accommodate a higher volumetric strain prior to failure compared to calcareous samples (Th1, 2, 4, 5) (Colour figure online) John formations reported by Sone (2012). Although 16 piezoceramic ultrasonic transducers are attached to the vertical core samples to record stress-induced seismicity (acoustic emissions) potentially induced by rock damage and cracking, the recording system did not detect any significant or representative activity in the frequency range of the transducers (0.1-1 MHz), neither during confining pressure nor during deviatoric stress changes.

Mohr-Coulomb Envelope and Mechanical Strength
A linearised Mohr-Coulomb failure envelope can be fitted to each MST test data, involving the Unconfined Compressive Strength UCS, the internal friction coefficient µ i , and the cohesion S 0 of the intact shale material through (Zoback 2010) where 1 = TCS + p c is the total peak axial stress at failure, n is the best-fit slope of the σ 1 versus p c data, σ n is effective normal stress, and τ is shear stress. Note that µ i is related to the best-fit slope n through i = (n−1) 2 √ n (see Table 5). Figure 11a shows the corresponding Mohr circles in the plane shear stress versus effective normal stress, and the derived linear Mohr-Coulomb failure envelope for samples Th6-V and Th6-H. Recall that this failure/yield analysis is based on the first four stages of each MST experiment conducted on a given sample at the same strain rate έ = 5 × 10 −7 s −1 . The fifth and final stage at the highest confining pressure is not included in this analysis because it was conducted at a higher loading rate to promote micro-seismic activity, and also because it was the only stage of each MST experiment at which the sample was brought to actual failure (observed peak stress). Figure 11b shows the failure envelope in the σ 1 versus p c plane for the available set of horizontal and vertical sample pairs. The peak stress of the thinly laminated and heterogeneous shale samples Th6-V and Th6-H exhibits a larger sensitivity to the confining pressure/normal stress than the (1) 1 = UCS + np c or = S 0 + i n , black homogeneous shale samples Th9-V and Th9-H. Moreover, the sensitivity of the total peak axial stress to confining pressure is generally higher for vertical samples (axial stress applied orthogonal to the bedding) than for horizontal samples (axial stress applied along the bedding).
The cohesion of the tested samples ranges between 9 and 15 MPa (see Table 5) and seems uncorrelated with their mineralogy. Figure 12 shows how the internal friction coefficient and unconfined compressive strength of the vertical Goldwyer shale samples correlate with their Young's modulus, bulk density, and the volume fraction of the mechanically weak phase, which includes clay minerals, organic matter, and porosity (ClayTocPHI). Implementing the propagation of errors method, the uncertainty on the determination of the internal friction coefficient from the triaxial test data is estimated to range between 0.02 and 0.15. In this figure USA's gas shale dataset from Sone (2012) are included for comparison purposes.
Both UCS and µ i slightly decrease with the increase in the weak phase fraction ClayTocPHI, although the correlation is relatively low (Fig. 12b, d). For the vertical samples of Goldwyer shale, the UCS ranges between 34 and 64 MPa, while µ i varies between 0.32 and 0.91, with an average of 0.58 ± 0.17. When the volume fraction of clay minerals, porosity and TOC is greater than 40% (Fig. 12b), the average internal friction coefficient reaches a value of 0.44, which is comparable to Europe's Posidonia and Bowland black shales (Herrmann et al. 2018).
An empirical linear relationship is established between µ i and the static vertical Young's modulus separately for the Goldwyer and USA gas shales (see Fig. 12a).  Moreover, a positive correlation is found between µ i and the bulk density when the data for the vertical samples of the Goldwyer and USA shales are combined (see Fig. 12c), as follows: The horizontal samples of the Goldwyer shale typically exhibit a larger internal friction coefficient than the vertical samples by about 0.2. Mandal et al. (2020b) showed how change in frictional strength along with viscoelastic behaviour of gas shale reservoirs control lithology dependent least principal stress magnitude. This suggests that larger amounts of fluid may be required to be injected in a horizontal well than in a vertical well during hydraulic fracturing in a horizontally layered shale formation (2) i = 1.44 b − 3.13, with R 2 = 0.43, for b ≥ 2.4 g∕cm 3 , global gas shales.
assuming that natural fracture networks are not favourability orientated in the specific in situ stress conditions. Figure 13a and b compares the static elastic properties (E and ν) obtained from the vertical and horizontal samples of the Goldwyer shale during the unloading part of the successive loading-unloading stages, to the static elastic properties obtained during the first loading stage of each multi-stage triaxial test. The successive unloading Young's modulus are 10-30% larger (10% for horizontal and 30% for vertical samples) than the first loading modulus, regardless of the mineralogy and petrophysical properties of the sample. In contrast, no clear correlation is observed for Poisson's ratio. This difference in Young's modulus raises the question of which ones are more suitable for analyzing reservoir deformation during production and the related stress perturbations. It is driven by strain hysteresis during loading-unloading cycles and has been extensively reported Fig. 11 a Mohr circles derived from five stages triaxial test data on the organic-rich shale samples (Th6-V and Th6-H). Mohr circles for vertical samples are plotted with dashed blue lines, the light blue dashed circle representing the last loading-tofailure stage. Green plain lines represent horizontal samples, while the red circle represents the last loading-to-failure stage. Only Mohr circles for the first four loading stages of the MST experiment are used to compute the linearized Mohr-Coulomb failure envelope, i.e., plain red line for the horizontal plug Th6-H, and plain black line for the vertical sample Th6-V. The intersection of these linear envelopes with the y-axis corresponds to the cohesion of the intact rock. b Influence of the confining pressure on total peak axial stress for all depth intervals for which two mutually orthogonal samples were available (Th3, Th5, Th6, Th9) (Colour figure online) in the literature for various rocks. Rate-dependent strain hysteresis is often associated with inelastic (irreversible) and/or anelastic (reversible with time delay) deformation mechanisms (Zoback 2010).
The effective mean stress value at which this hysteresis phenomenon becomes noticeable reflects the memory of the rock to the maximum stress it has experienced. In our dataset this distinction between loading and unloading Young's modulus is noticed for all confining pressures tested, above and below the in situ effective mean stress, which suggests that our shale samples may have not reached, during the MST experiment, the highest mean effective stress experienced by the formation in the past.
As shown in Fig. 13c, a strong positive correlation between compressive strength and first loading static Young's modulus is derived for the vertical samples of the Goldwyer formation. This correlation is confirmed when USA gas shale data from vertical samples reported by Sone (2012) are incorporated, i.e., Figure 13d shows a reasonable power-law correlation between compressive strength and bulk density for our Goldwyer shale samples (see Eq. 4). However, this local correlation is no longer valid once USA gas shale data are incorporated. However, USA gas shale data exhibit a similar power-law correlation, but with different coefficients (see also Eq. 4). This is possibly because the bulk density of gas shale is generally severely affected by TOC and porosity (low density constituents). Nevertheless, at and above a bulk density value of about 2.4 g/cm 3 , the compressive strength of the Goldwyer shale can be predicted from wireline density curve, after applying the standard  correction factor to remove the impact of organic matter content and gas effect from bulk density logs acquired at the subsurface gas shale reservoirs (Green and Vernik 2021).

Impact of Mineralogy and Petrophysical Properties
The constituents of the shale are ranked according to their mechanical strength and grouped into (i) a weak phase made of clay minerals, organic matter (TOC) and porosity (ClayTocPHI); (ii) an intermediate phase made of carbonates (Cb); and (iii) a strong phase made of quartz, feldspar, and pyrite (QFP). Figure 14 reports the empirical correlations devised for the laboratory data obtained on the vertical samples of Goldwyer and USA shales, i.e., between the triaxial compressive strength (Fig. 14a and  d), the first loading static Young's modulus ( Fig. 15b and e), the axial strain at failure (Fig. 14c and f)   influence on these mechanical properties, which is consistent with observations reported for European black shales (Herrmann et al. 2018 Rybacki et al. (2015) reported a positive correlation of the compressive strength with both the strong and intermediatestrength mineral phases when the distribution of the QFP and Cb phases were within the narrow range 15 to 20% and 30-40 vol%, respectively. Several factors could explain this apparently contrasting result such as the depth of the various shale formations compared here (confining pressure and Fig. 15 Static elastic anisotropy defined as the ratio of the horizontal to vertical static Young's modulus E h ∕E v (≥ 1) , for the Goldwyer and USA shales. a Correlation with the vertical Young's modulus E v , combining all shale data. b Correlation with the weak phase fraction ClayTocPHI, where the Goldwyer and USA shales are taken separately. For Goldwyer shale, elastic anisotropy empirical relationship is valid only when ClayTocPHI is greater than 40 vol% where anisotropy intercept is close to 1 temperature), the strain rate used to triaxially deform the shale, and the difference in microstructure, porosity, and water content. Studies conducted on several shales such as the Longmaxi shale form China (Jia et al. 2021), USA shales (Sone and Zoback 2013a, b), and European black shales (Herrmann et al. 2018;Rybacki et al. 2015) have shown a general negative correlation between compressive strength and the combined fractions of clay minerals and organic matter. However, above a fraction of 40 vol% of weak phase, compressive strength σ TCS remains virtually constant. A similar transition at around ~ 30-40 vol% of weak phase is also reported by Kohli and Zoback (2013). This transition in compressive strength probably reflects a change in shale's grain packing framework from strong-mineral support (strong phase QFP) to weak-mineral support (ClayTocPHI). This behaviour is consistent with the observed variation in compressive strength σ TCS with porosity (see Tables 1 and  3). Since intermediate-to-strong carbonates are the dominant mineral phase in our Goldwyer shale samples, they could well be acting as the primary load-bearing framework when the fraction of weak ClayTocPHI is lower than ~ 40%.
Young's modulus decreases with the increase in the weak phase fraction ClayTocPHI, following a power-law. Gas shale data reported by Sone (2012) on USA shales suggest a similar power-law trend, regardless of the reservoir conditions under consideration. This further confirms that for a given gas shale formation, and in the absence of additional knowledge/data, the combined fraction of clay minerals, organic matter and porosity is a very effective proxy for the static Young's modulus at depth. In agreement with the present study, Herrmann et al. (2018) and found a similar dependency of the static Young's modulus with the weak mineral phase. The empirical correlation reported in Eq. (3) can be indirectly used to predict the depth profile of compressive strength from the profile of dynamic Young's modulus derived from wireline logs after applying a suitable dynamic-to-static conversion.
Because carbonates predominantly control the loadbearing framework, and therefore the strength and elastic properties of the shale when the fraction of weak phase is lower than 40 vol%, their presence is expected to reduce axial strain for carbonate-dominated mudstones (Fig. 14c). Similarly, the maximum strain at failure ε max exhibits a positive correlation with the fraction of weak phases Clay-TocPHI when the latter become the load-bearing framework (Fig. 14f). For convenience, the empirical correlations reported in Fig. 14 are summarised here (5) TCS = 0.94 × Cb + 113.86, with R 2 = 0.52, for the Goldwyer gas shale, E v = 0.25 × Cb + 7.87, with R 2 = 0.77, for the Goldwyer gas shale, ε max = −0.06 × Cb + 10.97, with R 2 = 0.65, for the Goldwyer gas shale, Table 6 Spearman's rank correlation between different variables used for establishing empirical relationship reported in Eqs. (3), (4), (5), (6), (7), (8) and (15) Strong correlation, moderate correlation, uncorrelated

Input
Output Rank correlation-Spearman's rank correlation Since r-squared value (R 2 ) does not indicate the full robustness of the regression model, it is necessary to draw conclusions about the model by analysing r-squared together with the other variables in a statistical model. To improve reliability of the reported regression based empirical relationships (refer to Eqs. 3-8, and 15), a nonparametric rank correlation analysis is done to quantify association between two variables with a non-Gaussian distribution (Corder and Foreman 2009). Correlation coefficient quantifies the relationship between the observed values of two variables between − 1 (perfectly negative) to 1 (perfectly positive). It can be interpretated as trivial (0), weak (− 0.1 to 0.1), moderate (− 0.3 to 0.3), strong (− 0.5 to 0.5) and perfect (− 1 to 1). Here Spearman's rank correlation is introduced, and the outcome is presented in Table 6. From this rank correlation along with r-squared value, it is obvious that weak mineral phase ClayTocPHI is strongly correlated with static Young's modulus E and triaxial compression strength (σ TCS ) for the combined dataset (Goldwyer + USA shales), while intermediate phase Cb is statistically insignificant when correlating with compressive strength σ TCS . Further, it is noticed that the correlation of internal friction coefficient μ i with rock's bulk density ρ b showed moderate correlation and can be considered when building continuous profile at depth.
The above reported empirical correlations (ε max , E, and σ TCS with weak and intermediate-strong phases) can be used along the reservoir intervals to estimate reservoir compaction under given stress conditions. Note that Herrmann et al. (2018) reported similar correlations between the axial strain value at failure and the weak and intermediate-strong phases at room conditions, but not at elevated confining pressure and temperature conditions. In summary, the weak phase ClayTocPHI, and to a lesser extent the intermediate-strong phase Cb are found to play a key control in the deformation and strength of the studied gas shales.
be made for the Unconfined Compressive Strength (UCS). This minor anisotropy of compressive strength could be due to a preferential fracture propagation along the shale's weak planes/bedding. As often reported in the literature, shales loaded at an angle of ~ 45 deg generally exhibit the lowest compressive strength, regardless of the confining pressure/depth Gholami and Rasouli 2014;Holt et al. 2015;Islam and Skalle 2013).
In contrast, Young's modulus is generally lower orthogonal to the bedding than parallel to it (see Tables 3 and 4, and Fig. 14b, e). The closure of pre-existing microcracks aligned with the bedding (planes of weakness) is the most likely source of higher deformability in the direction orthogonal to this bedding. Similarly, a higher axial strain at failure ε max is generally expected for vertical shale samples (Herrmann et al. 2018;Sarout et al. 2014;Sarout and Guéguen 2008a, b). Figure 15a reports the static elastic anisotropy, defined as the ratio of horizontal to vertical static Young's modulus E h ∕E v (≥ 1) , as a function of the vertical Young's modulus for all studied shales, including USA shales. A strong power-law correlation is devised in this case (Eq. 7), which seems independent of the shale's geographical location, hence suggesting a universal correlation.
The degree of anisotropy is also plotted against the weak phase fraction in Fig. 15b and shows a different linear correlation for the Goldwyer and USA shales (Eq. 8), as expected from previous studies (Sone 2012;Sone and Zoback 2013a). The preferred orientation of platy clay minerals, the intrinsic anisotropy of clay minerals (Sarout et al. 2014;Sarout and Guéguen 2008a, b;Sone and Zoback 2013a), and the fraction of organic matter (Vernik and Liu 1997;Vernik and Nur 1992) have been found to control the overall anisotropy often observed in shales. Qualitatively, the increase in mechanical anisotropy can be explained by the presence of anisotropic clay platelets and organic matter since a general correlation already exists between TOC and the fraction of clay minerals for the Goldwyer formation (see Fig. 4a).

Prediction of Young's Modulus from Effective Medium Theory
It is possible to predict the elastic modulus of rocks using the effective medium concept if quantitative information about the fraction of individual mineral constituents, their elastic modulus and the geometrical arrangement of the constituents are available (Dvorkin et al. 2014;Sarout and Guéguen 2008b). In the absence of quantitative microstructural information, effective medium bounds can yet be devised (Dvorkin et al. 2014;Mavko et al. 2009). To this end, we group the mineral constituents into two endmembers based on their elastic stiffness: stiff minerals (Qtz, Cb, Fsp, Py) and soft phases (clays, TOC, porosity). We can then compute the Voigt, the Reuss, and the Hashin-Shtrikman bounds, as well as the simple arithmetic average of the Voight upper bound and Reuss lower bound known as Voigt-Reuss-Hill average, for the composite made of stiff and soft phases (Herrmann et al. 2018;Mandal et al. 2020c;Sone and Zoback 2013a). The upper and lower bounds of Young's modulus E for the homogeneous composite of these two phases are defined as (Dvorkin et al. 2014) where E is the effective Young's modulus of the composite; E stiff and E soft are the Young's modulus of the stiff and soft phase, respectively; f stiff and f soft are the volume fraction of the stiff and soft phase, respectively. The scaling exponent α takes the value − 1 for the lower bound (Reuss iso-stress model), and the value 1 for the upper bound (Voigt isostrain model). The simplest prediction of the actual Young's modulus is obtained by averaging Voigt and Reuss bounds and is known as Voigt-Reuss-Hill average (V-R-H average).
The equivalent Young's modulus E stiff (or E soft ) of the stiff (or soft) phase is the average of the Young's modulus of each of its constituents as repor ted in Table 7 [from Mavko et al. (2009)], where the weight factors are normalised by the constituent's volume fraction as reported in Table 2. A similar procedure can be applied to the bulk and shear modulus if these moduli are available instead of Young's modulus, i.e., Table 7 Elastic properties of individual mineral constituents (Mavko et al. 2009) and their weighting factor (normalized fraction of stiff and soft phase from Table 2), are used to bound and predict the vertical static Young's modulus from effective medium theory as a function of the weak phase ClayTocPHI (see Fig. 16) The weighting factors are defined from XRD analysis covered in Sect. 3.1. The stiff and soft fractions are the normalized weighted average of each individual mean composition from the combination of QFP + Cb (stiff phase) and ClayTocPHI (weak/soft phase), respectively   (10) can be used to compute the equivalent bulk and shear modulus of the stiff and soft phases, then derive the corresponding static Young's modulus for Voigt and Reuss bounds. Hashin-Shtrikman upper and lower bounds (HS ±) constitute tighter bounds than Voigt-Reuss's (Mavko et al. 2009), and can be used to further constrain the elastic properties of the simplified shale composite. Using the same bulk and shear modulus for the soft and stiff phases (Eq. 10), and HS ± equations reported in Mavko et al. (2009), we also derive the corresponding static Young's modulus for the HS ± bounds and for the effective V-R-H average. Figure 16 shows the evolution of the Voigt, the Reuss, and the Hashin-Shtrikman bounds, as well as the Voigt-Reuss-Hill average, as a function of the weak phase fraction ClayTocPHI. The laboratory-derived Young's modulus for the vertical and horizontal samples of the Goldwyer shale are also reported in this figure. The static Young's modulus of the vertical samples E v is confined near the Hashin-Shtrikman lower bound (HS-), whereas that of the horizontal samples E h lies between the Voigt-Reuss-Hill average and the Hashin-Shtrikman's lower bound (HS-). Note, however, that the simplified effective medium modelling (Herrmann et al. 2018) reported here neglects (i) shale anisotropy, and (ii) the impact of shale microstructure (relative spatial arrangement of the constituents).
In the literature, the static Young's modulus of vertical gas shale samples (bedding-perpendicular) from China, Europe, and USA mostly falls between the HS-and the Reuss bounds, generally closer to the Reuss bound. Overall, as intuitively expected, these effective medium models validate our observation that Young's modulus decreases with the increase in the weak phase fraction. Depending on the information available, these models can be used to constrain the Young's modulus of shale composites given their clay and organic matter content and their porosity. Voigt-Reuss bounds are the widest, Hashin-Shtrikman bounds are narrower and therefore more useful in practice, and the (10) Voigt-Reuss-Hill average can be used for direct first-order prediction.

Brittleness Index
The Brittleness Index (BI) of a rock is an empirical parameter characterising a rock's deformation regime under prevailing stress and temperature conditions. It ranges between 0 for ductile and 1 for brittle deformation. Brittle rocks are prone to fracturing and those fractures are more likely to remain open for a longer duration, whereas ductile rocks tend to deform uniformly, without strain localisation features. Recently, Mandal et al. (2020c) compared the various definitions of the Brittleness Index (BI) available in the literature, and analysed their suitability for optimal zonation of the prospective depths for successful hydraulic fracturing operations. Deformation-based brittleness index (BI E ) is devised as follows (Hucka and Das 1974) i.e., in which ε elastic and ε failure stand for elastic and total strain at failure from the failure stage of the stress-strain data. Mandal et al. (2020a) showed that a strong power-law relationship exists between the static Young's modulus and the deformation-based BI E parameter derived from triaxial stress-strain data for the Goldwyer shale formation, i.e., laboratory-derived static Young's modulus E sta Rybacki et al. (2016) derived a mineralogy-based brittleness index BI min based on data from a suite of shales where C x is the volumetric fraction of phase x defined in vol%; W QFP , W Cb , W ClayToc and W φ are the weighting factors statistically derived by Rybacki et al. (2016) from a global suite of shale data, i.e., W QFP = W ClayToc = W φ = 1 and W Cb = 0.5. Table 4 reports these two brittleness indices for all the tested Goldwyer shale samples. Generally, BI E is larger than BI min since QFP is not the dominant phase in our sample suite. Lower values of BI E are obtained as the weak phase dominates the sample's load-bearing framework. Figure 17 shows the correlation between the axial strain at failure ε max and the brittleness indices BI min (Fig. 17a) or BI E (Fig. 17b). The deformation-based brittleness index BI E provides a more discriminative range of values than the mineralogy-based one, i.e., BI E exhibits a clearer threshold value of 0.6 ± 0.05 above which a predominantly brittle deformation is confirmed by the stress-strain curves reported in Fig. 8a. The brittleness of horizontal samples seems 10-20% higher than that of vertical samples, suggesting that if we perform hydraulic fracturing operation in a horizontal well, a more pervasive and multi-directional fracture network will develop, without unnecessary vertical growth under consideration of several assumptions such as perforating layer with lower S hmin magnitude compared to overlying and underlying layers, more orientated natural fractures along the direction of S Hmax , stress shadow effect, etc. The correlation between ε max and BI E for the vertical samples of Goldwyer shale is (13)

Continuous Strength Profiling
Laboratory rock mechanical and elastic data can be used to derive continuous strength profiles from existing wireline logs. Several empirical equations have been devised for conventional overburden shales in multiple geological basins (Chang et al. 2006;Dewhurst et al. 2015;Horsrud 2001). Significant differences in mineralogical composition, structural complexity, maturity, saturation, rock strength, and frictional behaviour preclude the use of these correlations for the Goldwyer gas shale formation to derive profiles of UCS, internal friction coefficient or elastic anisotropy. Using our data for the Goldwyer formation and the parameter introduced earlier as the weak phase fraction ClayTocPHI, we derive the required correlations. In addition, we include in this analysis USA gas shale data because of their compositional similarity and the relative availability of similar data in the literature (Sone 2012;Sone and Zoback 2013a, b). Global empirical Eqs. (2), (3) and (4) have already been introduced that relate the triaxial compressive strength σ TCS and internal friction coefficient µ i to the vertical Young's modulus and to the bulk density, respectively, with the provided underlying assumptions and limitations. It is also possible to relate the UCS perpendicular to the bedding plane to the porosity or to the first loading static Young's modulus E v (see Fig. 18), i.e., In which ϕ is the porosity defined in %. In practice, knowledge of the first loading static Young's modulus allows for the estimation of the UCS, the triaxial compressive strength, while knowledge of the weak phase fraction indirectly provides bounds on the static Young's modulus. Based on the results presented so far for the Goldwyer shale formation, a suitable prospect interval exists below 1510 m depth. For example, the organic-rich mudstone samples with the highest maturity (Th3, Th6, Th9, Th10) located in the lower G-III interval exhibit the highest brittleness in the horizontal direction (along the bedding). Lower fracture closure rates when performing hydraulic fracturing operations in the horizontal direction are expected in this deeper interval of the Goldwyer-III unit. The samples in this interval also contain a significant amount of organic matter (TOC of 3.5-7 vol%), a reasonable porosity (6-8%) and a higher brittleness index BI E , above the threshold of 0.6.

Conclusions
Based on (i) the analysis and interpretation of multistage triaxial deformation experiments conducted on vertical and horizontal samples of Goldwyer shale formation (stress-strain data), (ii) their mineralogical, micro-structural and petrophysical properties, and (ii) the published data for USA gas shales, the following conclusion can be drawn: i. The mechanically weak phase ClayTocPHI mostly controls the mechanical (compressive strength, internal friction coefficient, unconfined compressive strength) and elastic (Young's modulus) properties of the Goldwyer gas shale. The same is valid for USA gas shales when included in the analysis. Both compressive strength and Young's modulus are lower in magnitude when the mechanically weak phase fraction (clay minerals + TOC + porosity) is larger than approximately 40 vol%. ii. The sensitivity of the compressive strength and the static Young's modulus to confining pressure is significant in the Goldwyer gas shale formation. The sensitivity of the compressive strength is higher for vertical than for horizontal samples. iii. The Goldwyer gas shale samples deformed perpendicular to the bedding exhibit a lower internal friction coefficient (0.58 ± 0.17) than the bedding-parallel samples (0.72 ± 0.12). Henceforth, it is anticipated that a higher amount of injection energy is required to generate hydraulic fractures in a horizontal well. iv. A consistent distinction between the first loading and unloading static Young's modulus (hysteresis) is observed at all stages of each MST, i.e., for all the confining pressures applied. This suggests that our shale samples may have not reached, during the MST experiment, the highest mean effective stress experienced by the formation in the past. v. The presence of a clay fabric correlates with the observed anisotropy of static Young's modulus (E h > E v ), which seems to hold for the combined Goldwyer and USA gas shales data. vi. For our Goldwyer gas shale samples, effective medium theory shows that the vertical static Young's modulus follows Reuss's lower bound, whereas the horizontal samples are best approximated by a simple arithmetic average of Voigt's upper and Reuss' lower bounds (the so-called Voigt-Reuss-Hill average). vii. The laboratory deformation experiments conducted on Goldwyer gas shale samples with varying composition and heterogeneity indicate a semi-brittle to brittle behaviour at their respective in situ stress conditions. Bedding-perpendicular deformation is more brittle compared to the observed semi-brittle behaviour of bedding parallel samples in the most prospective G-III unit (> 1510 m depth).

A Appendix: Weight-to-Volume Conversion of Rock Mineralogy
XRD analysis provides semi-quantitative mineral compositions of rock in weight percentage (wt%). For calculating the amount of reserve exists in a proven hydrocarbon field, volumetric percentage (vol%) is required to make reasonable prediction about field development plan and return of investment. In addition, mechanical behaviour of gas shale is also driven by volumetric composition (Herrmann et al. 2018;Rybacki et al. 2015;Sone and Zoback 2013a). Hence, in this study we reported here standard conversion process implemented from weight to volume through individual component's bulk density as in Table 8. First individual mineral components are converted and then normalized to 100 vol% for columns 4-8 as reported in Table 2. Later, when grouped into three phases (Strong, intermediate, and weak), further renormalization is done as can be seen from column 9-11 in Table 2.