1 Introduction

Shale is an important material in a range of engineering and energy-related applications including \(\hbox {CO}_2\) storage (Armitage et al. 2011), nuclear waste disposal (Guéry et al. 2010), and shale gas production. \(\hbox {CO}_2\) storage in depleted oil reservoirs requires a re-assessment of the mechanical integrity of the shale seal prior to long term storage since changes in stress state induced by production promote compaction of the overburden (Goulty 2003; Obradors-Prats et al. 2017, 2019). This increases the chance of brittle failure through overconsolidation of the formation rocks (Daigle et al. 2017; Ingram and Urai 1999; Nygård et al. 2004, 2006), or reactivation of existing fractures (Davies et al. 2013). The production of gas and oil from shale requires hydraulic fracturing of the formation, with the elastic properties of the reservoir being important in characterising brittleness (Rybacki et al. 2016) and to the development of a complex fracture network as a function of lithology, diagenesis, and burial history (Rickman et al. 2008). Elastic properties can also be linked to the viscoplastic response thus providing insights to fracture closure (Sone and Zoback 2013a, b; Rybacki et al. 2017). Sharp declines in gas production rates following fracturing are typical (Baihly et al. 2010), so a clearer understanding of the controls on the mechanical properties of shale should allow improvements in fracture network optimisation.

Despite its importance, the mechanical properties of shale are poorly constrained. Technical and economic challenges in recovering core limit traditional geomechanical testing (Kumar et al. 2012). The multiscale nature of shale is often overlooked, yet the nano- and micro-scale properties ultimately control both mechanical behaviour and anisotropy (Bobko and Ulm 2008). Since predictive modelling of pore network evolution and fracture growth requires a multiscale approach, characteristic mechanical properties for processes acting at appropriate length scales are needed to calibrate models. The grain size of component minerals in shale range from \(< 1\) to \(100\,\upmu \hbox {m}\), making direct mechanical characterisation of individual components difficult. However, developments in instrumental techniques now allow us to access the mechanical properties of materials at smaller length scales, allowing us to more readily access mineralogical and microstructural contributions to mechanical behaviour that were previously impossible to isolate.

Low-load nanoindentation can be used to assess mechanical properties at length scales of 500 nm–3000 nm. Previous studies have assessed the mechanical response of the composite clay–organic matrix of shale using statistical deconvolution methods on large datasets (Bobko and Ulm 2008; Ulm and Abousleiman 2006). Because clay and organic material (OM) are often intimately associated in shale, limited measurements of OM by nanoindentation have been generated (Zeszotarski et al. 2004; Zargari et al. 2013; Abedi et al. 2016). In thermally mature shales, organic matter is highly porous (Loucks et al. 2009, 2012; Zargari et al. 2015), especially at the scale of the nanoindenter tip, so that it is difficult to isolate the mechanical response of pure OM, since the elastic response is heavily influenced by porosity (Goodarzi et al. 2016). Coal has been used as a proxy (Borodich et al. 2015; Vranjes et al. 2018), but is subject to the same problems. The elastic properties of clay are also poorly constrained, with full elastic tensors only being available from molecular dynamics studies (Militzer et al. 2011; Monfared and Ulm 2016).

Atomic force microscopy (AFM) can extract mechanical properties at nanoscale resolution owing to the small diameter (10–30 nm) of the probe. Coupling with continuum mechanics models (e.g., Derjaguin et al. 1975; Sneddon 1965) makes it possible to measure the mechanical properties of the sample at nanometre resolution (Burnham and Colton 1989). Initially used in biomechanics (Smolyakov et al. 2016) and polymer chemistry (Tomasetti et al. 1998), nanomechanical characterisation modes have only recently been applied to geomechanical problems. Eliyahu et al. (2015), Emmanuel et al. (2016, 2017), Goodarzi et al. (2016), Goodarzi et al. (2017), Li et al. (2018b), Li et al. (2018a) have attempted to extract Young’s modulus values from different shales looking at the effects of maturity and temperature on the mechanical properties of OM. These studies all employ the PeakForce™quantitative nanomechanical mapping (PF-QNM) technique and suggest that the Young’s modulus of organic matter is in the range 1–32 GPa. Recently, Fender et al. (2020) have provided the first application of an alternative form of AFM, the quantitative imaging (QI-AFM), to geological materials for the study of coal macerals. It is therefore required to assess how this technique performs against more established methods to be able to full interpret its results, and to understand how QI-AFM may be extended to other geological materials.

Various studies of OM in shale cite other studies to validate their results. However, Jakob et al. (2019) have shown that intra-sample geochemical composition can vary significantly at the nanoscale. Such intra-sample variability therefore precludes inter-sample comparison unless the repeatability of AFM measurements can be demonstrated. The results of Jakob et al. (2019) show that mechanical differences between OM in shales that are equally mature from the perspective of RockEval may still yield mechanical differences at the nanoscale that are the result of maturity-related heterogeneity. Therefore, demonstrating the repeatability of AFM results is crucial to interpret their mechanics reliably.

Furthermore, Eliyahu et al. (2015) observe that AFM appears to yield different Young’s modulus values for mineral phases such as calcite and quartz compared to values reported in the literature (Mavko et al. 2009). It is still unclear if this represents calibration drift in the measurement (Young et al. 2011), or if it represents a true mechanical response due morphological factors, or crystal impurities which can affect its elastic properties (Fournier et al. 2011, 2014). This ambiguity needs to be addressed to fully understand the range of Young’s moduli to which AFM can be applied, an in order better understand how we can apply nano- and micro-mechanical data in models.

Mechanical data from AFM are normally validated by comparison to nanoindentation; for example, Young et al. (2011) show good agreement between PF-QNM and nanoindentation results for stiff polymers up to \(E\approx 3\,\mathrm{GPa}\). However, for organic-rich shales, this method of validation is problematic, since isolating intimately mixed organic matter and clay is difficult. Furthermore, both OM and clay minerals undergo a complex series of chemical reactions during burial, resulting in compositional changes (e.g., Thyberg and Jahren 2011; Jakob et al. 2019) making it difficult to create samples in an equivalent diagenetic state to the shale.

This study aims to address issues surrounding repeatability of AFM measurements combining the use of two independent AFM approaches on a single sample. This comparison also serves to assess the performance of QI-AFM relative to PF-QNM as this will show the range of geological materials this method can be applied to. For inorganic, non-clay phases, AFM will be compared to low-load nanoindentation to address the issue of the low Young’s moduli for calcite and quartz reported by Eliyahu et al. (2015).                                            

The organic-rich shales used in this study are highly calcareous. Existing micro-mechanical data for well-characterised shales focus on the siliciclastic end-member. Despite this, the Niobrara, Eagle Ford, Green River, Tarfaya, and parts of the Horn River shales are all examples of notable organic-rich, calcareous shales that contribute significantly to hydrocarbon production. The results presented in this study will provide a robust micro-mechanical dataset that is more appropriate for modelling the macro-scale elasticity of calcareous shales.

2 Sample Characterisation and Preparation

2.1 Sample Characterisation

This study focuses on two shale samples from the dry gas window (Mullen 2010). Sample mineralogy was characterised by X-ray powder diffraction (Table 1). Both samples are carbonate–clay mudstones with clay present as predominantly mixed layer illite–smectite. In addition to bulk rock mineralogy, local compositional mapping was carried out using energy-dispersive X-ray spectroscopy (EDX), with simultaneous acquisition of backscattered electron images (BSE) (Figs. 1, 2). Samples were also analysed for total organic carbon (TOC).

Table 1 Sample mineralogy and total organic carbon (vol %)
Full size table
Fig. 1
figure 1

BSE image of CC8 with associated EDX maps used for mineralogical interpretation. The red box outlines the area targeted for AFM mechanical mapping. OM organic matter, Cal calcite, Clay clay mineral, Pyr pyrite

Full size image
Fig. 2
figure 2

BSE image of TS10 with associated EDX maps used for mineralogical interpretation. The red box outlines the area targeted for AFM mechanical mapping. OM organic matter, Cal calcite, Qtz quartz, Fsp feldspar

Full size image

2.2 Sample Preparation

AFM requires a flat sample surface. To achieve this, samples were first prepared as thick sections cut perpendicular to the bedding planes. Pucks of 2.5 mm diameter were then cut from the thick sections and further polished by broad ion beam milling to remove plucking artefacts (Amirmajdi et al. 2009; Loucks et al. 2009) which contribute to significant short wavelength changes in topography. This procedure results in topographic variations on the order of 100 nm and low surface angles (Fig. 3a). Changes in surface angle are noted where there is a significant contrast in stiffness between adjacent phases (compare Fig. 2 with Fig. 3a), since softer phases are less resistant to milling; however, over 90% of the surface is inclined at less than \(10^{\circ }\) (Fig. 3b), sufficient for AFM analysis.

Fig. 3
figure 3

a Topography map of TS10. b Probability distribution of surface angle for region in (a). The modal angle is \(\approx \,4^{\circ }\) with 90% of the surface with inclination of less than \(10^{\circ }\)

Full size image

3 Analytical Methods

3.1 Atomic Force Microscopy

Two different types of nanomechanical characterisation AFM modes were employed in this study. First, PeakForce quantitative nanomechanical mapping (PF-QNM) was deployed using a Bruker multimode AFM equipped with a Nanoscope V controller. In addition, and for comparison, we use the quantitative imaging (QI) mode developed by JPK, and referred to here as QI-AFM. This uses the nanowizard™III platform. In force spectroscopy mode, a piezoelectric scanner moves the sample/probe vertically. As the tip makes contact with the sample surface, the deflection in the cantilever is measured by the position of the laser on a quadruple photo-diode detector (Fig. 4a). As it approaches, electrostatic forces pull the tip to the surface (adhesive regime). The point at which adhesion is at its maximum value is defined as the tip-sample contact point (\(z_\mathrm{c}\) in Fig. 4b). Beyond this, the tip enters the repulsive contact regime where force increases until the user-defined setpoint (defined in calibration, see Sect. 3.1.2) is reached. This point corresponds to the position of maximum z-piezo extension, after which the tip is retracted. The portion of the curve between \(z_c\) and the setpoint is used for fitting of the contact model from which mechanical properties are calculated (Sect. 3.1.2).

Fig. 4
figure 4

a Schematic representation of AFM operation. b Force–distance curves during tip approach–withdraw cycle. The adhesive component is exaggerated for illustrative purposes, especially on the approach curve

Full size image

The operation of the AFM outlined above is common to PF-QNM and QI-AFM, with the main difference being the motion z-piezo as a function of time. PF-QNM uses sinusoidal motion in the z direction, whereas QI-AFM uses a sawtooth motion. The sinusoidal motion of the z-piezo in PF-QNM results in velocity and acceleration profiles that are also smooth and continuous. The contact model employed by the two systems also differs (see Sect. 3.1.2 for more detail). Another important difference is the oscillation speed of the probe, which varies between 0.25 and 2 kHz for PF-QNM and between 0.5 and 100 Hz for QI-AFM.

3.1.1 AFM Probe Selection

Successful measurement of the elastic properties of a material by AFM requires an appropriate tip/cantilever assembly. The cantilever needs an adequate spring constant to prevent excessive deformation during contact. The tip must be sufficiently stiff to avoid significant deformation itself. Pittenger et al. (2014, 2011) provide the working ranges of commercially available cantilevers and calibration standards. Two tips are suitable (Fig. 5a) for organic matter (OM), but only the DNISP probe is capable of measuring stiffer inorganic phases based on results from existing studies (Eliyahu et al. 2015; Emmanuel et al. 2016, 2017; Goodarzi et al. 2016; Li et al. 2018a). The DNISP probe has a diamond corner-cube tip on a stainless steel cantilever (Fig. 5b), whereas the TAP525 made from silicon nitride. DNISP probes are superior in that they theoretically allow all material components in shale to be characterised without changing tip, but are significantly more expensive.

Fig. 5
figure 5

a The operational limits of AFM probes, modified from Pittenger et al. (2011). In addition to the working ranges of the tips, commercially available calibration standards are also shown. Expected ranges of minerals are shown for comparison. O.M. organic material, C.M. clay minerals. DNISP probe is a diamond tip on stainless steel cantilever; TAP525 is silicon nitride. b SEM secondary electron image of a DNISP cantilever after scanning. The presence of material removed from the sample ate tip apex should be noted

Full size image

3.1.2 Contact Model

In line with previous studies that use PF-QNM (Eliyahu et al. 2015; Emmanuel et al. 2016; Goodarzi et al. 2016; Emmanuel et al. 2017; Li et al. 2018b), we make use of the Dejaguin–Muller–Toporov (DMT) contact model (Derjaguin et al. 1975) as the basis for our calculation of Young’s modulus. This model is an extension of the classical Hertz–Sneddon contact model for the indentation of an elastic half space (Sneddon 1965). The Hertz–Sneddon model is derived for a spherical indenter and is used in this case to approximate the tip of the AFM probes. In this scheme, the so-called reduced Young’s modulus is given by:

$$\begin{aligned} F_{\mathrm{hertz}} = \frac{4E^*\sqrt{R(d-d_0)^3}}{3}, \end{aligned}$$
(1)

where F is the force acting on the sample/tip, \(E^*\) is the reduced Young’s modulus, and \((d-d_0)\) is the tip-sample separation distance. Reduced Young’s modulus can be converted to absolute Young’s modulus using the following relationship:

$$\begin{aligned} E^* = \left[ \frac{1-\nu _t^2}{E_t} + \frac{1-\nu _s^2}{E_s}\right] ^{-1}, \end{aligned}$$
(2)

where \(\nu _t\) and \(\nu _s\) are the Poisson’s ratio of the tip and sample, respectively, and \(E_t\) is the Young’s modulus of the tip. For a DNISP probe, \(E_t\) is on the order of 1.1TPa, so \((1-\nu _t^2)/E_t\rightarrow 0\). Throughout the rest of this paper, we only consider the reduced modulus value, but drop the use of reduced for brevity. The DMT model additionally accounts for adhesion between the tip and the sample via an additional linear term:

$$\begin{aligned} F_{\mathrm{int}} = F_{\mathrm{hertz}} + F_{\mathrm{adh}}, \end{aligned}$$
(3)

where \(F_{\mathrm{adh}}\) is the adhesive component of the total interaction force \(F_{\mathrm{int}}\). PF-QNM fits Eq. 3 to the retract curve (Fig 4b), (Bruker 2011a, b). QI-AFM fits Eq. 1 to the approach curve (JPK 2017), since the adhesive component is negligible compared to the retract curve, and also negligible relative to the peak interaction force/setpoint.

3.1.3 Instrument Calibration

The deflection sensitivity \(S [\mathrm {nm/V}]\) of the cantilever, and its spring constant \(k [{\mathrm {N m}}^{-1}]\) need to be determined. Cantilever deflection, D, is initially in volts, and is converted to length through the deflection sensitivity constant. The spring constant allows applied force to be calculated by Hooke’s law.

To determine S, the tip is ramped against a sapphire standard of comparable stiffness to the tip to force negligible deformation of the sample. Representative deflection sensitivities found during calibration for both DNISP and TAP525 tips using QI-AFM and PF-QNM are given in Table 2. The QI system is less sensitive to cantilever deflection than the PF-QNM system. This results from the different photo-diode detectors are used in each system. For quartzofeldspathic and carbonate phases, the expected elastic moduli are high and the indentation of the sample surface is expected to be small. For a fixed instrumental precision, high sensitivity may lead to lower measured sample indentation and overestimation of the elastic modulus. In the case of the DNISP probes, k is determined by the manufacturer. For the TAP525, the spring constant was determined by the thermal tune method (Bruker 2011a).

For PF-QNM, the optimal scanning setpoint (maximum applied force) is determined by scanning a highly ordered pyrolytic graphite (HOPG) standard (\(E^*\approx\) 18 GPa). This is chosen as it lies within the range of clay mineral and OM stiffness (Fig. 5a). Calibration determines the force needed to give a sample deformation of at least \(2\,\mathrm{nm}\) according to the manufacturer’s recommendations Bruker (2011a, (2011b). Young et al. (2011) report that if \(|E_\mathrm{sample} - E_{\mathrm{HOPG}}| \gg 0\), the error in AFM measurement increases. Variation is expected in the \(E^*\) of the standard, and this is accounted for by letting the value of \(E^*\) fall within one standard deviation of the mean measured value (Trtik et al. 2012). The same approach is taken in this study. The results of calibration for both probes are given in Fig. 6. There is good agreement between expected and obtained values. QI-AFM does not require the scanning of a standard, since modulus is not calculated in real time and can only be done by processing the results in the output software. The choice of setpoint is then somewhat arbitrary. For all scans using QI-AFM in this study, a setpoint of 500 nN was chosen.

Fig. 6
figure 6

Histograms derived from data acquired by PF-QNM scanning the HOPG standard using a DNISP probe and b TAP525 probe

Full size image

3.2 Nanoindentation

To access the mechanical response at length scales above \(\sim 10\,\upmu \mathrm{m}\), nanoindentation is employed. A relatively well-established technique in geomechanics compared to AFM (Bobko and Ulm 2008; Liu and Ostadhassan 2017; Liu et al. 2018; Lu et al. 2019; Akono et al. 2019; Abedi et al. 2016), only a brief overview of the technique is provided here. For composite materials, the load applied in indentation testing is important to distinguish between local and bulk mechanical responses (Constantinides et al. 2006). In this study, nanoindentation is used to look at both the local and bulk response of the samples, so two load cases are considered (Fig. 7).

Fig. 7
figure 7

Reduced modulus maps of sample TS10 as measured with a TAP525 cantilever. From the scan conducted with PF-QNM, it can be seen that in inorganic regions, a uniform background reading of \(E^* \sim\) 25 GPa is recorded, with some regions in white where the instrument has not resolved a stiffness value. This suggests that an upper limit on the resolution of the TAP525 is set at 25 GPa

Full size image

The local response was measured by the low-load grid indentation technique (Constantinides et al. 2006). 400 indents were applied to each sample as four randomly positioned \(10\times 10\) measurement grids (100 indents per grid). This was done to capture any contribution of local spatial variation in the sample. Indentation spacing was \(5\,\upmu \mathrm{m}\) within each grid, with an applied load of 5 mN in all cases. Low-load testing was carried out using a Hysitron triboindenter system (Newcastle University, UK). Bulk shale response was characterised by a single grid of \(7\times 7\) indents at 100 \(\upmu\)m separation. as per the methodology of (Goodarzi et al. 2017; Goodarzi 2018). Maximum load was 250 mN. High-load testing was carried out using a NanoTest™ Vantage system.

All nanoindentation curves were analysed using the classical approach of Oliver and Pharr (1992); Pharr et al. (1992), in which the reduced modulus, \(E^*\), is calculated as:

$$\begin{aligned} E^* = \left. \frac{\sqrt{\pi }}{2\sqrt{A_\mathrm{c}}}\frac{\partial P}{\partial h}\right| _{h = h_{\mathrm{max}}}, \end{aligned}$$
(4)

where \(A_\mathrm{c}\) is the contact area between the sample and the indenter, and \(\partial P/\partial h\) is the gradient of the unloading curve at maximum indentation depth \(h_{\mathrm{max}}\). Reduced modulus can be converted to true Young’s modulus by Eq. 2. Only reduced modulus values are considered. Data from low-load tests were subject to statistical deconvolution methods to separate individual chemo-mechanical components. This was carried out using the procedure of Ulm et al. (2007).

3.3 Mechanical Segmentation and Data Analysis

Segmentation of AFM elastic moduli maps into organic and inorganic components was carried out through comparison to backscattered electron scanning electron microscopy (BSE), and EDX images of the same sample regions. The methodology of Camp and Wawak (2013) was applied to BSE images to produce pseudo-greyscale images for thresholding. Four greyscale classes were identified, and comparison with EDX element maps showed that these can be interpreted as organic carbon, carbonate, aluminosilicate \(+\) quartz, and ‘heavy’ phases (pyrite, apatite, etc.). Since thresholding is controlled by pseudo-greyscale number derived from BSE images, separation of aluminosilicate and quartz phases was difficult due to low contrast difference. These phases were manually segmented by comparing AFM maps with SEM and EDX images. AFM data were analysed using NanoScope Wizard software v.1.9 (PF-QNM) and JPKDP software v.6.1 (QI-AFM).

4 Results and Discussion

We now present the results of our comparison of QI-AFM and PF-QNM measurements of elastic responses, along with a comparison between the values of elastic modulus as a function of the tip used to determine them.

4.1 Organic Matter

For sample CC8, the average modal value of \(E^*\) for OM based on all instrument/probe configurations is \(9.8 \pm 1.1\,\hbox { GPa}\). In the case of sample TS10, the average modal value is \(5.2 \pm 1.1\,\hbox { GPa}\). In both cases, the error is the standard deviation of all instrument/probe configurations for a given sample (see Table 4 for the values used in these calculations). For each sample, the resulting organic probability density functions are shown in Fig. 8; these were used to construct confidence intervals for the expected distribution at the 95% confidence interval. For both samples, all instrument-probe configurations produce distributions that fall within the 95% confidence interval, suggesting that a high level of consistency in the mechanical responses of organic material is achieved by both PF-QNM and QI-AFM systems. Additionally, the DNISP tips work as well as the TAP525 tips in the low stiffness range.

Fig. 8
figure 8

Probability density functions (PDF) for the reduced Young’s modulus, \(E^*\), of organic matter in samples CC8 and TS10. The PDFs are determined using the kernel density estimator approach

Full size image

Our measurements are at the lower end of the 5–30 GPa range measured by Emmanuel et al. (2016) in mature kerogen. The different mean moduli between samples probably reflect chemical differences in OM, supported by observations by Jakob et al. (2019) who showed that there is significant chemo-mechanical variation even within a single sample. The effects of OM chemistry on mechanical properties are still an open question (Li et al. 2018b, a; Goodarzi 2018), although some preliminary work has been done by Emmanuel et al. (2016), who suggest that increasing maturity causes stiffening in organic material. Eliyahu et al. (2015) suggested that co-existing of kerogen and bitumen has distinct mechanical properties.

For sample TS10, PF-QNM produces a slightly higher modal reduced modulus for OM when compared to QI-AFM. The higher deflection sensitivity of the cantilevers used on the PF-QNM system may explain this, since this lowers the resolution with which the instrument can detect cantilever bending. This results in underestimation of the true deflection, leading to the measurements being artificially stiff. The opposite is seen for CC8. For sensitivity alone to explain these differences, we should also expect the experiments run with the TAP525 to yield a lower \(E^*\) than the DNISP probes, regardless of instrument, as a TAP525 has lower deflection sensitivity than a DNISP probe for both instruments (Table 2). For TS10, however, TAP525 tips give higher modal values than DNISP probes. These results suggest that overall the differences in sensitivities between probes and instruments do not play a significant role in the measured Young’s modulus distribution, with differences being attributable to experimental error.

Table 2 Cantilever and tip calibrated and nominal values for the types used in this study
Full size table

DNISP probes produce a greater difference in modal \(E^*\) between instruments, 2.2 GPa for CC8 and 1.8 GPa for TS10 compared to differences of 0.9 GPa and 1.4 GPa for CC8 and TS10, respectively using TAP525. This is attributed to the fact that in the case of OM, DNISP probes are working at their lower operational limit. Indeed, Pittenger et al. (2014) put the lower operational limit of DNISP probes at 10 GPa (solid blue line in Fig. 5a), so we are theoretically out of cantilever capability when dealing with OM. The generally comparable results when compared to the TAP525 response for both samples suggest that the values obtained are still valid. Moreover, this suggests that the lower limit of DNISP probes can be taken down to \(\sim \,5 \mathrm{GPa}\).

4.2 Inorganic Phases

In these samples, non-clay inorganic components are generally physically larger, and more volumetrically dominant than both clay minerals and organic material (Table 1). It is thus possible to measure their mechanical response using low-load nanoindentation as well as AFM.

Table 3 Numerical results of deconvolution of nanoindentation results shown in Fig. 9
Full size table

4.2.1 Atomic Force Microscopy and Nanoindentation

Once segmented, AFM data were fitted to a kernel density function to estimate probability density functions. Due to the 25 GPa operational limit of TAP525 tips (see Fig. 7), only DNISP probes are considered for this part of the study. Reduced moduli results from AFM are shown in Fig. 9a, b, with low-load nanoindentation results in Fig. 9c, d.

Fig. 9
figure 9

Probability density functions for the reduced modulus of inorganic phases from a CC8 and b TS10 measured by both AFM systems (upper panel) and histograms of raw low-load nanoindentation results (lower panel) overlain by probability density functions of the reduced modulus of the mechanically active phases based on the results of deconvolution. The numerical results of deconvolution are provided in Table 3. 400 indents were used in the analysis of each sample

Full size image

The fine-scale intergrowth of clay minerals and organic matter means that segmentation of these phases is incomplete, reflected in small secondary peaks in the inorganic probability density functions at \(E^*\approx\) 10 GPa. PF-QNM results for sample CC8 show a peak at \(E^*\approx\) 22 GPa that is not present in sample TS10. This is interpreted as the response of clay minerals, as EDS mapping of CC8 indicates localised packets of clay materials (Figs. 1, 2) which are not present in the region studied within TS10.

A third peak is present in both samples between 40 and 50 GPa. This peak defines the dominant inorganic phase in TS10 and contributes half the total probability density in CC8. We infer that this phase is calcite and note that the modulus values between 40 and 50 GPa are somewhat lower than those previously measured for calcite in similar fine-grained sedimentary rocks. Using PF-QNM and a DNISP probe, Eliyahu et al. (2015) reported Young’s moduli of \(53\pm 6\hbox { GPa}\) for calcite, whilst Emmanuel et al. (2016, (2017) measured values around 55 GPa, again using a PF-QNM setup.

The modulus values corresponding to calcium carbonate are in agreement with other studies of similar rocks. Eliyahu et al. (2015) report calcite to have reduced Young’s modulus of \(53 \pm 6\,\hbox { GPa}\), using PF-QNM and a DNISP probe. Emmanuel et al. (2016, (2017) also note calcite in the range of 55 GPa, again using a PF-QNM setup.

The results of the deconvolution of nanoindentation data are presented in Table 3. For CC8 and TS10, the phase interpreted as calcite has reduced moduli of \(\sim \,40\pm 6.5\,\mathrm{GPa}\) (\(1\sigma\)), very similar to the results from PF-QNM. For TS10 values of 40 GPa are measured by both techniques, with CC8 yielding 50 GPa from AFM and 42 GPa from nanoindentation. The AFM results for CC8 is outside one standard deviation of the nanoindentation result of \(41.9\,\mathrm{GPa}\,\pm 6.3\); however, inspection of the histogram CC8 shows that mean value for calcite should be closer to 45 GPa. This mismatch arises, because the deconvolution algorithm simultaneously solves for material hardness (Table 4).

Table 4 Summary results of reduced Young’s modulus for minerals and organic matter measured by AFM and low-load nanoindentation
Full size table

Moduli of 40–50 GPa are much lower than those reported for calcite crystals obtained by acoustic methods. Chen et al. (2001) report bulk and shear moduli of 76.1 GPa and 32.8 GPa, respectively using Brillouin spectroscopy, corresponding to a Young’s modulus of 86 GPa. Mavko et al. (2009) suggest a value of 60 GPa. The difference in Young’s moduli reflects the morphology of calcite. Experiments used to determine the full set of elastic constants typically use large crystals, whereas calcite in sedimentary rocks typically is of biological origin or is micritic, with the majority of calcite existing as clay like platelets rather than well-defined euhedral crystals. Measured Young’s moduli for micrite areas show areas as low as 1–20 GPa (Zhang et al. 2018) up to 70 GPa (Fournier et al. 2011), with many values in the range 30–50 GPa (Vialle and Lebedev 2015). De Paula et al. (2010) present nanoindentation data for oolitic limestone that has a bimodal distribution with a major peak at 56 GPa.

Bobko and Ulm (2008) carried out an extensive nanoindentation study of clay-rich shales, from which they determined the mechanical properties of ‘solid clay’ by extrapolation to a packing density of unity, i.e., zero porosity. For the bedding perpendicular direction, they report a Young’s modulus of 25 GPa, similar to values (a) obtained through inversion of ultrasonic pulse velocity tests in a micro-mechanics framework (Ortega et al. 2007) and (b) this study \(\sim 22\) GPa by PF-QNM.

Our segmented AFM results for combined plagioclase and quartz give a single peak with a mean \(E^*\) of 55 GPa (Fig. 9, in good agreement with the results of deconvolution from low-load nanoindentation tests. As for calcite, these values are lower than those reported in the literature on larger crystals. Mavko et al. (2009)’s compilation of Young’s moduli of quartz gives a range of 77–96 GPa, and Zhu et al. (2007) carried out nanoindentation quartz crystals in granite yielding \(E^*\) of 101 GPa. The elastic properties of feldspars are strongly dependent on their chemistry. Brown et al. (2016) give a range of Young’s modulus which varies in the range of 87 GPa for albite \(\hbox {An}_0\) to 103 GPa for anorthite \(\hbox {An}_{100}\), with labradorite (An\(_{50{-}70}\)) having a Young’s modulus of 97 GPa. Using nanoindentation, Zhu et al. (2007) measured \(E^* \approx\) 62 GPa for albite in granite.

One explanation for the lower Young’s moduli for quartz and feldspar in our study would be the occurrence of microporosity (Ulm et al. 2007; Goodarzi et al. 2016), which can amount to a few volume percent in feldspars (Worden et al. 1990) as a result of dissolution reactions and also diagenetic transformation to clay minerals (Jin et al. 2011). The small tip size of the DNISP probe used in this study (\(\sim\)30 nm) results in deformations of \(\sim\)2 nm, so that the mechanical properties measured by AFM may be expected to vary depending on the amount of microporosity sampled at each point of contact. On the other hand, the larger tip size of the Berkovich indenter will sample a more representative volume of feldspar microporosity. Since our results from AFM and nanoindentation are similar, we infer that both methods approximate the bulk mechanical response of the mineral grains, and that it is the differing micro-morphology of both feldspars and quartz in fine-grained sedimentary rocks, compared to igneous rocks and more ideal crystals, which is most likely to explain the lower Young’s modulus.

4.2.2 High-Load Nanoindentation

The results of high-load nanoindentation testing are presented in Fig. 10. Figure 10a, c shows the map of reduced modulus over an area \(600\times 600\) microns for CC8 and TS10, respectively, compared to an area of \(95\times 95\) microns for low-load tests, and \(10\times 10\) microns for AFM. The corresponding modulus histograms are shown in Fig. 10b, d.

Fig. 10
figure 10

a Reduced modulus map of sample CC8. b Histogram of raw data, overlain by best-fit normal distribution. c Reduced modulus map of sample TS10. d Histogram and best-fit normal distribution

Full size image

Figure 10b, d, shows a single well-defined modulus peak for both samples, in contrast to the bimodal and trimodal distributions from low-load tests in Fig. 9c, d. Summary statistics for the best-fit Gaussian distributions to the high-load test data are given in Table 5. The high-load bulk responses of both samples are similar, around 42 GPa, and also in agreement with the results of the low-load indentation and AFM tests. Moduli around 42 GPa are consistent with the bulk response being controlled by the major mineral phase present in these samples, calcite (Table 1), with similar amounts of both stiffer (quartz and feldspar) and softer (clay and organic material) in both samples.

Table 5 Mean (\(E^*_{\mathrm{Bulk}}\)) and standard deviation (\(\sigma\)) of bulk Young’s modulus derived from n high-load nanoindentation tests for CC8 and TS10
Full size table

4.3 Implications for Upscalling and Anisotropy

It is useful to consider the general implications of our new results for future micro-mechanical modelling of the elasticity of organic-rich calcareous shales. We do this by computing the Voigt and Reuss bounds for a simplified four-component shale comprising of organic matter, clay, calcite, and quartzofeldspathic phases (Fig. 11). Similar simplified conceptual models are have been used previously, however, our use of a four-component system reflects the availability of a robust micro-mechanical data set. This is less conservative in its predictions than the very general two-component models used by previous authors where shales are divided into generally soft and stiff phases (Mavko et al. 2009; Bayuk et al. 2008; Sone and Zoback 2013a, b; Herrmann et al. 2018). Hill (1963) used an energy balance argument to show that the Voigt and Reuss bounds must be considered as upper and lower bounds, respectively, on the Young’s modulus of a composite material. Therefore, we must observe the high-load nanoindentation composite Young’s modulus falling within these bounds.

The Voigt bound is the arithmetic mean of the individual Young’s moduli of minerals in a shale sample, and the Reuss bound is a harmonic mean based on the same input parameters. They are calculated from information about components in the composite using Eq. 5, where \(f_i\) represents the ith fraction of an N component mixture with a corresponding Young’s modulus \(E_i\). These bounds are calculated as follows:

$$\begin{aligned} \displaystyle {E_{\mathrm{Voigt}} = \sum _{i=1}^N f_iE_i \quad \mathrm{and} \quad E_{\mathrm{Reuss}} = \left[ \sum _{i = 1}^N \frac{f_i}{E_i} \right] ^{-1}}. \end{aligned}$$
(5)

The robustness of the measured data can be demonstrated by computing Voigt and Reuss bounds on the composite material using AFM and low-load nanoindentation data. Results of these calculations are presented in Fig. 11. It can be seen that Voigt and Reuss bounds obtained from mineral moduli using either AFM or low-load nanoindentation provide better constraints on the high-load bulk response. It should be noted that there is a difference between locally measured AFM and low-load data. This is likely due to the fact that volume fractions used for AFM data in Eq. 5 are derived from XRD (Table 1), whereas the bounds computed from low-load indentation data use the fractions extracted from deconvolution. This suggests that it may be better to represent shale in terms of its mechanically active phases rather than its mineralogical ones, allowing consideration of both microstructural (e.g., inter-grain interactions, cementation, etc.) and mechanical variation (e.g., grain orientation and intra-grain chemo-mechanical variation).

In contrast, it can be observed from Fig. 11 that when \(E_{\mathrm{Bulk}}\) is calculated using mineral moduli from the literature, the results for the bounds are higher than the bulk response measured using high-load nanoindentation. Even when the experimental error in the high-load indentation data is taken into account, the bulk responses do not reach the Reuss bound using published moduli for highly crystalline quartz, calcite, and feldspar. The error bars on the composite Young’s modulus in Fig. 11 represent one standard deviation. If we consider the 95% confidence interval on the mean as less conservative measure of error, then it can be seen that previously assumed Young’s modulus values for quartz and calcite do not seem to provide bounds on the composite that adequately reflect the homogenised Young’s modulus.

The use of previously assumed moduli for mineral phases for estimation of the composite elastic responses is still prevalent in the homogenisation literature (Giraud et al. 2007; Ortega et al. 2007; Guéry et al. 2010; Goodarzi et al. 2017; Zhang et al. 2020) despite the potential for significant overestimation of the true value as demonstrated by our results. This emphasises the need to use micro-mechanically appropriate Young’s modulus values for all component phases when estimating the components of the composite stiffness matrix \(C_{ij}\) from mineralogical data, not just the clay matrix as in previous studies (Ulm and Abousleiman 2006; Bobko and Ulm 2008). The assumption that macroscopically derived elastic moduli are appropriate may lead to bounding values on the composite that are energetically inconsistent based on the results of Hill (1963). Indeed, where calibration of a mean-field homogenisation model is carried out with respect to the clay matrix, the interpretation of the obtained constants (e.g., in terms clay mineralogy) may be undermined by the assumed Young’s moduli of the inclusion phases. This is due to the calibrated clay phase which will be acting to correct for other inaccuracies in the model.

The anisotropy of shale is of importance in many field applications relating to wellbore stability and seismic interpretation (Ortega et al. 2007; Sayers et al. 2015). Therefore, it is advantageous to include it in micro-mechanical studies to understand the origins of macro-scale anisotropy. The scaling of the Young’s modulus of shale building blocks over several orders of magnitude conducted in this work allows for advanced micro-mechanical modelling work, where anisotropy is explicitly considered in determining the elastic constants. Figure 11 shows that when XRD mineral fractions are used in conjunction with the values reported in Table 4, the composite Young’s modulus in both the \(x_1\) and \(x_3\) material directions (bedding parallel and normal respectively) falls within the calculated bounds. Provided that a suitable estimated on the transversely anisotropic clay matrix is used, this implies that the values for calcite and quartz-feldspar obtained in this study are valid for use in typical two-level conceptual homogenisation models (Ulm et al. 2005).

Fig. 11
figure 11

Young’s modulus from high-load indentation (\(E_{\mathrm{Bulk}}\)) and the corresponding Voigt, and Reuss bounds on the homogenised Young’s modulus based on the measured component responses using AFM and low-load nanoindentation and literature data. AFM bounds are calculated using the XRD mineral volumes using a simplified mineralogy of Quartz + Feldspar, Calcite, Clay (= Chl + I-S/ML), and TOC. Moduli used in the calculations are listed in Table 4, with the Voigt and Reuss bounds calculated using the mean ± the peak half-width, respectively. Low-load indentation bounds are calculated using the fractions in Table 3. The corresponding reduced moduli are the mean ± one standard deviation for the Voigt and Reuss bounds, respectively. Values based on literature data table values use the fractions from low-load indentation tests as it is felt that allows the bounds based on literature data to reflect component interactions in the shale composite, thus providing better constraints on the value of the bounds. \(E_{\mathrm{Cly+TOC+Pore}} = 25\,\mathrm{GPa}\) (Bobko and Ulm 2008), \(E_{\mathrm{QFP}} = 60\,\mathrm{GPa}\) and \(E_{\mathrm{QFP}} = 83\,\mathrm{GPa}\) (Mavko et al. 2009)

Full size image

The mineral components of shale, with the exception of OM and pyrite are intrinsically anisotropic. Deposition of mineral grains of quartz, calcite, and feldspar on the bedding surface can be considered a random process. Additionally, the low energy of deposition associated with shale can lead to no preferred crystal orientation. This has the effect of averaging out the effects of intrinsic crystal anisotropy as it is manifested in the low-load indention results. Since the grids cover a large region of the sample relative to grain size, this results in the indentations encountering numerous locally varied crystal orientations. This results in the deconvolution procedure producing a singular peak representative of the equivalent isotropic material.

In Fig. 11, if the bounds are calculated using the fractions from deconvolution in Eq. 5, then only the composite Young’s modulus in the \(x_1\) direction falls within them. This being the same direction in which low-load indentation was carried out. Deconvolution fractions are partially a product of inter-particle interactions representing mechanically active phases. Calcite is normally considered an inclusion, and considered mechanically isotropic from a modelling perspective. However, given the volumetric dominance of calcite, it is interesting to note that the bulk Young’s modulus in the \(x_3\) direction does not fall inside the mechanical bounds for either sample when deconvolution fractions are used. This observation poses an interesting question as to whether or not calcite needs to be considered in an anisotropic manner given the comparatively low volume of clay relative to calcite in the studied samples. From a practical perspective, the majority of homogenisation studies continue to use XRD data to inform component volumes. Therefore, there is still practical merit in using the values in Table 4 for isotropic inclusions. Future work should focus on the mechanical character of calcite in calcareous shales. We have already noted the good agreement of the indentation modulus for clay measured by AFM in this study, and that reported by Ulm and Abousleiman (2006); Bobko and Ulm (2008); Ortega (2010). This suggest in part that the elastic constants of Ortega et al. (2007) for a transversely isotropic clay matrix would be suitable starting points for more advanced micro-mechanical modelling of the shales considered in this study in conjunctions with the results in Table 4.

4.4 Measured Modulus and Length Scales

The total area covered by our low-load indentation grids is 9025 \(\upmu \mathrm{m}^2\), compared to \(100\,\upmu \mathrm{m}^2\) for AFM scans. AFM measures variation in Young’s modulus over a single or a few mineral grains, whereas deconvolution of nanoindentation data samples the average stiffness of a phase from a much larger number of individual grains. Despite the difference in the length scale of observation, the spread in the distribution of moduli for both calcite and quartz-feldspar about their respective means is similar for both AFM and low-load indentation (see Fig. 9). The agreement between AFM and low-load nanoindentation is surprisingly good given that the Young’s modulus on inorganic phases is 20–30 GPa larger than the HOPG standard used to calibrate the AFM. In future, the use of AFM to determine mechanical properties of complex natural composite materials may be advantageous since whilst individual indentation measurements can take up to 5 min to obtain, and AFM scan of \(\sim 2.6\times 10^5\) measurements can be acquired in around 30 min.

5 Conclusions

This study has presented the first direct in situ comparison between the mechanical response of rock forming minerals over length scales from \(\sim 10\,\hbox {nm}\) to \(\sim 100\,\upmu m\) by combining AFM, low-load and high-load nanoindentation approaches. We present the first multi-instrument approach to assess the mechanical response of organic matter, and demonstrate a highly consistent response which can be obtained, yielding Young’s modulus in the range \(E^*=\) 5–10 GPa. This implies that AFM is a viable tool to answer more detailed questions about organic matter mechanics with implications for cap rock leakage and petroleum expulsion.

AFM with a DNISP probe and low-load nanoindentation give similar results for the Young’s moduli of clay (\(22\pm 2.0\,\hbox { GPa}\)), calcite (40–50 \(\pm 8.5\,\hbox { GPa}\)), quartz, and feldspar (\(55\pm 10.1\)). For non-clay, inorganic phases, the mechanical responses in terms of Young’s modulus are significantly lower than previous measurements obtained on macro-scale highly crystalline phases, but similar to other in situ studies of individual minerals in sedimentary rocks.

High-load nanoindentation generates a unimodal responses in the range 40—50 GPa for both samples studied here, and are consistent with calcite being the major mechanically active phase, from low-load nanoindentation and AFM results. We have calculated, Voigt, Reuss, and Voigt–Reuss–Hill bounds on the bulk response using low-load nanoindentation, AFM, and literature values for individual mineral phases. Low-load nanoindentation and AFM data provide the best estimates for constraints on the composite response at length scales of \(\sim\) 100–600 \(\upmu \mathrm{m}\) measured from high-load indentation data. Conversely, moduli measured on highly crystalline macro-scale mineral phases do not provide effective bounds on the composite response. This questions the validity of using such macro-scale measurements in homogenisation schemes, and suggests that more emphasis should be placed on the use of nano- and micro-mechanical measurements when using effective medium theories and homogenisation schemes to predict the bulk mechanical response. Homogenising to length scales of \(\sim \,100\,\upmu \mathrm{m}\) may be considered the ‘undamaged’ homogeneous response, thus better prediction of properties at this scale, which may allow for a more realistic inclusion of damage and microcracks into homogenisation schemes and a better prediction of Young’s modulus at triaxial scale.

Since low-load nanoindentation data intrinsically couple mineral mechanical properties with how they interact with their microstructural context, it may be better to conceptualise shale in terms of its mechanically active phases, rather than purely in terms of its mineralogical phases thus allowing us to include mechanical responses arising from microstructural factors into homogenisation schemes without additional mathematical complexity.