# Modeling the Viscoelastic Behavior of Quartz and Clay Minerals in Shale by Nanoindentation Creep Tests

Authors: Jianfeng Wang; Yuke Liu; Chao Yang; Wenmin Jiang; Yun Li; Yongqiang Xiong
Journal: Geofluids
 (2022)

---

## Abstract

The viscoelastic behavior of minerals in shales is important in predicting the macroscale creep behavior of heterogeneous bulk shale. In this study, in situ indentation measurements of two major constitutive minerals (i.e., quartz and clay) in Longmaxi Formation shale from the Sichuan Basin, South China, were conducted using a nanoindentation technique and high-resolution optical microscope. Firstly, quartz and clay minerals were identified under an optical microscope based on their morphology, surface features, reflection characteristics, particle shapes, and indentation responses. Three viscoelastic models (i.e., three-element Voigt, Burger’s, and two-dashpot Kelvin models) were then used to fit the creep data for both minerals. Finally, the effects of peak load on the viscoelastic behavior of quartz and clay minerals were investigated. Our results show that the sizes of the residual imprints on clay minerals were larger than that of quartz for a specific peak load. Moreover, the initial creep rates and depths in clay minerals were higher than those in quartz. However, the creep rates of quartz and clay minerals displayed similar trends, which were independent of peak load. In addition, all three viscoelastic models produced good fits to the experimental data. However, due to the poor fit in the initial holding stage of the three-element Voigt model and instability of the two-dashpot Kelvin model, Burger’s model is best in obtaining the regression parameters. The regression results indicate that the viscoelastic parameters obtained by these models are associated with peak load, and that a relatively small peak load is more reliable for the determination of viscoelastic parameters. Furthermore, the regression values for the viscoelastic parameters of clay minerals were lower than those of quartz and the bulk shale, suggesting the former facilitates the viscoelastic deformation of shale. Our study provides a better understanding of the nanoscale viscoelastic properties of shale, which can be used to predict the time-dependent deformation of shale.


---

## Body

## 1. Introduction

The marine Longmaxi shale in South China is characterized by extremely low porosity and permeability, and has been exploited for shale gas by multistage hydraulic fracturing and horizontal drilling techniques [1]. Typically, the successful application of these techniques and stable shale gas production require the injection of large volumes of fracturing fluids and proppants. Rock creep, particularly for clay-rich shale, provides an intrinsic contribution to the time-dependent embedment of proppant [2]. Thus, the creep behavior of reservoir shales must be known for accurate prediction of shale gas production [3, 4].Conventional uniaxial–triaxial compression tests can be used to investigate the macroscale creep properties of shale under complex geological conditions [5, 6]. However, this method has disadvantages, as it requires a large sample size while being time-consuming and expensive [7–10]. Comparatively, nanoindentation only requires a small shale sample, and has been developed as a reliable technique to characterize the microscale creep behavior of shales [8, 11–17].The nanoindentation technique presses a small-scale load in millinewtons (mN) onto a sample surface using a diamond indenter [18, 19]. The creep properties of the indented materials are studied by the recorded creep depth, which reflects a change in the indentation depth under a constant peak load for a period of time. Nanoindentation has been used to investigate the anisotropic creep characteristics of shale [13], creep deformation of hard and soft mineral phases within shale [4, 8], relationship between creep response and porosity, total organic carbon (TOC) contents, and mineralogy of shale [14], and the effect of organic matter on shale creep rates [15]. However, in a nanoindentation creep test, the indentation load is important, because variable indentation loads affect different volumes and produce different mechanical responses [20]. A small indentation peak load (several mN) is usually used to measure the creep behavior of individual constituents within shales, as the volume affected by the indenter is small [15, 21, 22]. Although this assumption is theoretically plausible, there are no experimental data that show the correlation between creep parameters and indentation load for individual constituents within shales. Previous study suggested higher peak loads can increase the viscosity of the silica gel in the carbonated wollastonite matrix using nanoindentation tests [23]. Therefore, whether the creep behavior can be influnenced by peak load for individual constituents within Longmaxi shale or not, is still unclear.Previous studies have investigated the viscoelastic behavior of shales using the creep depth, which is a function of time and can be regressed to calculate viscoelastic parameters using viscoelastic models [2, 4]. Viscoelastic models (i.e., three-element Voigt and Burger’s models) consist of spring and dashpot elements, and can predict material viscoelastic behavior throughout the deformation process, as well as link the basic theory, experimental data, and field application of this phenomenon [2]. However, only a few experiments have used these models to simulate the viscoelastic behavior of shale at the nanoscale [2, 4], and the viscoelastic behavior of single mineral phases in shale is not fully understood. The macroscopic creep properties of shale depend on its micro-structures. Therefore, it is highly desirable and necessary that an appropriate viscoelastic model should be selected to predict viscoelastic behavior of single mineral phases in Longmaxi shale.The Longmaxi Formation shale is a fine-grained sedimentary rock composed mainly of large amounts of brittle minerals, such as quartz and abundant clay minerals (≥75 vol.%; [24–26]). The high contents mean these minerals have a critical role in creep deformation of the Longmaxi Formation shales. Therefore, quartz and clay minerals were investigated in this study to constrain their viscoelastic behavior on the nanoscale. It is difficult to directly conduct in situ nanoindentation experiments on minerals in complex and heterogeneous shales, because of the low resolution of optical microscopy. Therefore, statistical nanoindentation methods have commonly been used to determine the creep behavior of shale composites and their constituents [4, 8, 15, 21, 22]. This method requires a large amount of indentation data to be acquired, and the relevant parameters for each phase are extracted with a Gaussian distribution model [27–29]. In theory, the Gaussian distribution model can determine the creep behavior of a single mineral phase in a shale. However, it requires clear identification of the phase composition of the material, and care is needed during the indentation not to indent contacts between two or more phase, which will affect the test results. Moreover, it can also increase the uncertainty of the results, because low-volume phases are not considered.In view of the limited efforts to characterize the creep behavior of minerals in shales, it is necessary to conduct in situ nanoindentation testing of minerals. In this study, quartz and clay minerals in a shale were identified under a high-resolution optical microscope based on reflected colors, surface features, and particle shapes, in combination with the indentation responses. The effects of peak load on the creep rate, creep depth, and viscoelastic parameters of quartz and clay minerals were then determined. Finally, the applicability of three viscoelastic models was compared and creep time constants were obtained. The results improved our understanding of the intrinsic creep behavior and mechanisms of major constituents in shale at nanoscale, which further facilitates an accurate prediction of time-dependent behavior of shale reservoir during gas production.


## 2. Materials and Experiments

### 2.1. Materials and Sample Preparation

#### 2.1.1. Shale Material

A shale core sample from a depth of 89 m was collected from the lower Silurian Longmaxi Formation, Sichuan Basin, South China. X-ray diffraction results showed that the sample consists of quartz (34.4%), feldspar (13.6%), carbonates (4.5%), pyrite (1.5%), and clay minerals (45.9%) (i.e., the same sample used by Yang et al. [24]). Quartz and clay minerals account for 80.3% of the total mineral content. The TOC content of the sample is 3.4 wt.%, and the calculated equivalent vitrinite reflectance is 2.9%, which suggests the sample is over-mature.


#### 2.1.2. Sample Preparation

The surface of a shale sample for nanoindentation testing must be completely flat to obtain reliable results [18, 30–32]. In this study, the shale was cut perpendicular to the bedding and then cast into epoxy resin. Various silicon carbide abrasive papers from 50 to 2000 grit were used to polish the sample surface. Subsequently, 6 and 1.0 μm diamond suspension polishing fluids and 50 nm colloidal silica suspension were used for further surface polishing. Finally, an IM4000 argon ion mill operated at 4.0 kV and a 5° angle, which was produced by Hitachi High-Tech Global Corporation, was used for polishing for 2 h to achieve a completely flat surface.


### 2.2. Nanoindentation

#### 2.2.1. Theoretical Background

In this study, the nanoindentation experiments were conducted with an Anton Paar TTX NHT3 nanoindenter, equipped with a diamond Berkovich indenter that has a three-sided pyramidal diamond tip with a radius of curvature of 50-100 nm (Figures 1(a) and 1(b)). The stiffness threshold was 500 N/m and its spring compliance was 0.571 mm/N. The load of the instrument is 0.1 to 500 mN. The resolutions of the load and displacement settings were 20 nN and 0.01 nm, respectively. The indentation system is also equipped with an optical microscope with a maximum magnification of 4000×, which is suitable for in situ nanoindentation tests.Figure 1
(a) TTX-NHT3 indentation instrument, (b) the main components of the instrument, (c) schematic diagram of a nanoindentation and (d) typical load–displacement (P–h) curve [14]. P is the applied load, Pmax is the peak load, S is the contact stiffness, h is the indentation displacement, and hmax is the maximum indentation depth.
(a)(b)(c)(d)Prior to the indentation tests, a standard fused silicon specimen was used in quasi-static nanoindentation experiments to calibrate the tip shape. For all the indentation tests, the indenter approached and retracted from the sample surface at a rate of 2000 nm/min within a distance of 2000 nm of the sample surface. Typically, nanoindentation is used to apply a small-scale load (several to hundreds of mN) to a material surface through a diamond indenter. The load and displacement data were continuously monitored during the indentation process and a load–displacement curve (P–h curve) was obtained. The hardness and Young’s modulus were calculated with the Pharr–Oliver equations according to the P–h curve [18].Figure1(c) is a schematic diagram of the nanoindentation [14]. Figure 1(d) shows a typical P–h curve, which is divided into three stages (i.e., loading, holding, and unloading). By using a continuous mechanical model, the hardness H and reduced Young’s modulus Er of the material can be obtained from the P–h curve:
(1)H=PmaxAc(2)Er=Sπ2βAcwhere S=dP/dh is obtained from the initial slope of the unloading curve, which represents the stiffness of the material, β is a constant that depends on the geometry of the indenter (for a common Berkovich indenter, β =1.034), and Ac is the projected contact area, which is calculated from the maximum indentation depth [18, 33].Furthermore, based on the reduced Young’s modulus and Poisson’s ratio, the Young’s modulus of the material can be defined as:(3)1Er=1−v2E+1−vi2Eiwhere Ei and νi are the Young’s modulus and Poisson’s ratio of the diamond indenter, respectively. For the diamond indenter, Ei =1140 GPa and νi =0.07. E and ν are the Young’s modulus and Poisson’s ratio of the samples prepared for testing, respectively. In this study, the Poisson’s ratios of the quartz and clay minerals are 0.06 and 0.3, respectively [34, 35]. If the uncertainty on Poisson’s ratio is 0.1, then it results in a 5% uncertainty on the measured Young’s modulus [36].


#### 2.2.2. Experimental Protocols

A constant loading strain rate test was used to investigate the viscoelastic behavior of quartz and clay minerals, which means that a loading force was applied to the sample surface and then held for a period of time before being unloaded. To understand the effects of the indentation load, different values of peak load (5, 10, 30, and 50 mN) were used at a constant loading strain rate of 0.1 s–1, and then held for a period of 120 s. Finally, the indentation test was unloaded over 10 s. During the holding time stage, the applied load remained unchanged, while the displacement increased due to creep. In general, load–displacement curves with pop-in behavior, due to yielding or dislocation of the material [37], and/or data with a relatively large thermal drift must be discarded. To ensure experimental repeatability and acquire statistically significant results, at least six tests were performed for a specific peak load. Three representative tests with relatively stable data or good load–displacement curves were selected for parameter fitting.


## 2.1. Materials and Sample Preparation

### 2.1.1. Shale Material

A shale core sample from a depth of 89 m was collected from the lower Silurian Longmaxi Formation, Sichuan Basin, South China. X-ray diffraction results showed that the sample consists of quartz (34.4%), feldspar (13.6%), carbonates (4.5%), pyrite (1.5%), and clay minerals (45.9%) (i.e., the same sample used by Yang et al. [24]). Quartz and clay minerals account for 80.3% of the total mineral content. The TOC content of the sample is 3.4 wt.%, and the calculated equivalent vitrinite reflectance is 2.9%, which suggests the sample is over-mature.


### 2.1.2. Sample Preparation

The surface of a shale sample for nanoindentation testing must be completely flat to obtain reliable results [18, 30–32]. In this study, the shale was cut perpendicular to the bedding and then cast into epoxy resin. Various silicon carbide abrasive papers from 50 to 2000 grit were used to polish the sample surface. Subsequently, 6 and 1.0 μm diamond suspension polishing fluids and 50 nm colloidal silica suspension were used for further surface polishing. Finally, an IM4000 argon ion mill operated at 4.0 kV and a 5° angle, which was produced by Hitachi High-Tech Global Corporation, was used for polishing for 2 h to achieve a completely flat surface.


## 2.1.1. Shale Material

A shale core sample from a depth of 89 m was collected from the lower Silurian Longmaxi Formation, Sichuan Basin, South China. X-ray diffraction results showed that the sample consists of quartz (34.4%), feldspar (13.6%), carbonates (4.5%), pyrite (1.5%), and clay minerals (45.9%) (i.e., the same sample used by Yang et al. [24]). Quartz and clay minerals account for 80.3% of the total mineral content. The TOC content of the sample is 3.4 wt.%, and the calculated equivalent vitrinite reflectance is 2.9%, which suggests the sample is over-mature.


## 2.1.2. Sample Preparation

The surface of a shale sample for nanoindentation testing must be completely flat to obtain reliable results [18, 30–32]. In this study, the shale was cut perpendicular to the bedding and then cast into epoxy resin. Various silicon carbide abrasive papers from 50 to 2000 grit were used to polish the sample surface. Subsequently, 6 and 1.0 μm diamond suspension polishing fluids and 50 nm colloidal silica suspension were used for further surface polishing. Finally, an IM4000 argon ion mill operated at 4.0 kV and a 5° angle, which was produced by Hitachi High-Tech Global Corporation, was used for polishing for 2 h to achieve a completely flat surface.


## 2.2. Nanoindentation

### 2.2.1. Theoretical Background

In this study, the nanoindentation experiments were conducted with an Anton Paar TTX NHT3 nanoindenter, equipped with a diamond Berkovich indenter that has a three-sided pyramidal diamond tip with a radius of curvature of 50-100 nm (Figures 1(a) and 1(b)). The stiffness threshold was 500 N/m and its spring compliance was 0.571 mm/N. The load of the instrument is 0.1 to 500 mN. The resolutions of the load and displacement settings were 20 nN and 0.01 nm, respectively. The indentation system is also equipped with an optical microscope with a maximum magnification of 4000×, which is suitable for in situ nanoindentation tests.Figure 1
(a) TTX-NHT3 indentation instrument, (b) the main components of the instrument, (c) schematic diagram of a nanoindentation and (d) typical load–displacement (P–h) curve [14]. P is the applied load, Pmax is the peak load, S is the contact stiffness, h is the indentation displacement, and hmax is the maximum indentation depth.
(a)(b)(c)(d)Prior to the indentation tests, a standard fused silicon specimen was used in quasi-static nanoindentation experiments to calibrate the tip shape. For all the indentation tests, the indenter approached and retracted from the sample surface at a rate of 2000 nm/min within a distance of 2000 nm of the sample surface. Typically, nanoindentation is used to apply a small-scale load (several to hundreds of mN) to a material surface through a diamond indenter. The load and displacement data were continuously monitored during the indentation process and a load–displacement curve (P–h curve) was obtained. The hardness and Young’s modulus were calculated with the Pharr–Oliver equations according to the P–h curve [18].Figure1(c) is a schematic diagram of the nanoindentation [14]. Figure 1(d) shows a typical P–h curve, which is divided into three stages (i.e., loading, holding, and unloading). By using a continuous mechanical model, the hardness H and reduced Young’s modulus Er of the material can be obtained from the P–h curve:
(1)H=PmaxAc(2)Er=Sπ2βAcwhere S=dP/dh is obtained from the initial slope of the unloading curve, which represents the stiffness of the material, β is a constant that depends on the geometry of the indenter (for a common Berkovich indenter, β =1.034), and Ac is the projected contact area, which is calculated from the maximum indentation depth [18, 33].Furthermore, based on the reduced Young’s modulus and Poisson’s ratio, the Young’s modulus of the material can be defined as:(3)1Er=1−v2E+1−vi2Eiwhere Ei and νi are the Young’s modulus and Poisson’s ratio of the diamond indenter, respectively. For the diamond indenter, Ei =1140 GPa and νi =0.07. E and ν are the Young’s modulus and Poisson’s ratio of the samples prepared for testing, respectively. In this study, the Poisson’s ratios of the quartz and clay minerals are 0.06 and 0.3, respectively [34, 35]. If the uncertainty on Poisson’s ratio is 0.1, then it results in a 5% uncertainty on the measured Young’s modulus [36].


### 2.2.2. Experimental Protocols

A constant loading strain rate test was used to investigate the viscoelastic behavior of quartz and clay minerals, which means that a loading force was applied to the sample surface and then held for a period of time before being unloaded. To understand the effects of the indentation load, different values of peak load (5, 10, 30, and 50 mN) were used at a constant loading strain rate of 0.1 s–1, and then held for a period of 120 s. Finally, the indentation test was unloaded over 10 s. During the holding time stage, the applied load remained unchanged, while the displacement increased due to creep. In general, load–displacement curves with pop-in behavior, due to yielding or dislocation of the material [37], and/or data with a relatively large thermal drift must be discarded. To ensure experimental repeatability and acquire statistically significant results, at least six tests were performed for a specific peak load. Three representative tests with relatively stable data or good load–displacement curves were selected for parameter fitting.


## 2.2.1. Theoretical Background

In this study, the nanoindentation experiments were conducted with an Anton Paar TTX NHT3 nanoindenter, equipped with a diamond Berkovich indenter that has a three-sided pyramidal diamond tip with a radius of curvature of 50-100 nm (Figures 1(a) and 1(b)). The stiffness threshold was 500 N/m and its spring compliance was 0.571 mm/N. The load of the instrument is 0.1 to 500 mN. The resolutions of the load and displacement settings were 20 nN and 0.01 nm, respectively. The indentation system is also equipped with an optical microscope with a maximum magnification of 4000×, which is suitable for in situ nanoindentation tests.Figure 1
(a) TTX-NHT3 indentation instrument, (b) the main components of the instrument, (c) schematic diagram of a nanoindentation and (d) typical load–displacement (P–h) curve [14]. P is the applied load, Pmax is the peak load, S is the contact stiffness, h is the indentation displacement, and hmax is the maximum indentation depth.
(a)(b)(c)(d)Prior to the indentation tests, a standard fused silicon specimen was used in quasi-static nanoindentation experiments to calibrate the tip shape. For all the indentation tests, the indenter approached and retracted from the sample surface at a rate of 2000 nm/min within a distance of 2000 nm of the sample surface. Typically, nanoindentation is used to apply a small-scale load (several to hundreds of mN) to a material surface through a diamond indenter. The load and displacement data were continuously monitored during the indentation process and a load–displacement curve (P–h curve) was obtained. The hardness and Young’s modulus were calculated with the Pharr–Oliver equations according to the P–h curve [18].Figure1(c) is a schematic diagram of the nanoindentation [14]. Figure 1(d) shows a typical P–h curve, which is divided into three stages (i.e., loading, holding, and unloading). By using a continuous mechanical model, the hardness H and reduced Young’s modulus Er of the material can be obtained from the P–h curve:
(1)H=PmaxAc(2)Er=Sπ2βAcwhere S=dP/dh is obtained from the initial slope of the unloading curve, which represents the stiffness of the material, β is a constant that depends on the geometry of the indenter (for a common Berkovich indenter, β =1.034), and Ac is the projected contact area, which is calculated from the maximum indentation depth [18, 33].Furthermore, based on the reduced Young’s modulus and Poisson’s ratio, the Young’s modulus of the material can be defined as:(3)1Er=1−v2E+1−vi2Eiwhere Ei and νi are the Young’s modulus and Poisson’s ratio of the diamond indenter, respectively. For the diamond indenter, Ei =1140 GPa and νi =0.07. E and ν are the Young’s modulus and Poisson’s ratio of the samples prepared for testing, respectively. In this study, the Poisson’s ratios of the quartz and clay minerals are 0.06 and 0.3, respectively [34, 35]. If the uncertainty on Poisson’s ratio is 0.1, then it results in a 5% uncertainty on the measured Young’s modulus [36].


## 2.2.2. Experimental Protocols

A constant loading strain rate test was used to investigate the viscoelastic behavior of quartz and clay minerals, which means that a loading force was applied to the sample surface and then held for a period of time before being unloaded. To understand the effects of the indentation load, different values of peak load (5, 10, 30, and 50 mN) were used at a constant loading strain rate of 0.1 s–1, and then held for a period of 120 s. Finally, the indentation test was unloaded over 10 s. During the holding time stage, the applied load remained unchanged, while the displacement increased due to creep. In general, load–displacement curves with pop-in behavior, due to yielding or dislocation of the material [37], and/or data with a relatively large thermal drift must be discarded. To ensure experimental repeatability and acquire statistically significant results, at least six tests were performed for a specific peak load. Three representative tests with relatively stable data or good load–displacement curves were selected for parameter fitting.


## 3. Results and Discussion

### 3.1. Identification of Quartz and Clay Minerals under the Optical Microscope

To conduct the in situ nanoindentation experiments on the minerals in the shale, quartz and clay minerals must firstly be identified under an optical microscope. These minerals can be clearly distinguished based on their different morphological and optical characteristics, combined with their reflection features and particle shapes.In the Longmaxi Formation shale, quartz grains are translucent with small impurities, and have dark gray colors under reflected light (Figures2(a) and 2(b)). They have irregular shapes with straight edges. The grain size is mostly between 15–40 μm. The quartz grains subjected to mechanical polishing have a smoother and cleaner surface than those subjected to ion polishing, which is the most distinct morphological feature (Figure 2(a)–2(c)). Thus, quartz grains subjected to mechanical polishing were chosen for the in situ nanoindentation experiments. The quartz can be divided into authigenic and detrital quartz [26]. However, the authigenic quartz (mostly <2 μm in size) has a similar color as the other dark gray minerals [26], which makes it difficult to identify.Figure 2
Photomicrographs of quartz and clay minerals. The surface morphology of quartz and clay minerals is shown after (a–b) mechanical polishing and (c–f) ion polishing. Quartz is dark gray, translucent, contains small impurities, and has straight and sharp edges without a regular shape. The clay minerals are opaque and black after mechanical polishing, and located in pit areas as lumps or isolated grains along with organic matter. Magnifications: (a) and (c) 20× objective lens; (b) and (d) 50× objective lens; (e–f) 100× objective lens (the eyepiece magnification was 40× for all images).After mechanical polishing, the clay minerals in the shale are opaque and black under the optical microscope. The clay minerals infill pitted areas as lumps or dispersed grains (Figures2(a) and 2(b)), and may have no obvious crystal form after ion polishing (Figure 2(e)). This texture is formed by the dense stacking of flaky clay particles with sizes of tens of nanometers [38, 39]. In general, the clay minerals are continuously distributed and form the mineral framework of the shale. The surfaces of the clay minerals are not flat and cannot be focused on after mechanical polishing, whereas the surface height between quartz and clay minerals is nearly the same after ion polishing (Figure 2(a)–2(d)), which was reported in our previous study [24]. Therefore, unlike quartz, the nanoindentation tests on the clay minerals were carried out on ion-polished shale slices. It should be noted that organic matter is interspersed within the clay grains (Figure 2(d)–2(f)) and has an off-white color, which was also observed in our previous study [24]. Clay minerals mixed with organic matter occur as aggregates in the shale matrix [12] (Figure 2(c)–2(f)). Given that the clay grains are very small, it is not possible to determine the independent creep response of individual clay particles. Thus, the organic–clay matrix was analyzed as a whole. Similarly, clay minerals in the Longmaxi Formation shale also have different origins, including terrigenous clastic, authigenic, and secondary types [40], which are difficult to distinguish under an optical microscope.Although our method is unable to identify the origins of the quartz and clay minerals, it is suitable for the in situ nanoindentation tests. In addition to quartz and clay minerals, the Longmaxi Formation shale also contains other minerals, such as pyrite, carbonates, and feldspar. These minerals are low in content, and have very small grain sizes, especially for pyrite (Figure2(e)–2(f)), thus were excluded from discussion in this study.


### 3.2. Indentation Responses of Quartz and Clay Minerals

In theory, a larger mineral crystal is more suitable for the nanoindentation tests, because it is not readily affected by adjacent phases or the substrate. Due to the small size of the minerals in the shale, representative quartz grains larger than 20μm and an area of clay minerals as large as possible were chosen for the in situ nanoindentation experiments. Figure 3 shows optical microscopy images of the quartz and clay minerals before and after indentation at peak loads of 5–50 mN. The residual imprints indicate that the larger peak loads produce a greater indentation imprint on both the quartz and clay minerals, and the indentations are unclear at a peak load of 5 mN. The reason for this is that a larger peak load causes a larger contact area on the mineral surface. In addition, the sizes of the residual imprints on the clay minerals is larger than those on quartz under the same peak load conditions, which is consistent with the higher strength and stiffness of quartz.Figure 3
Photomicrographs of a quartz grain (a) before and (b) after nanoindentation tests, and clay minerals (c) before and (d–f) after nanoindentation tests.The corresponding load–displacement curves for the quartz and clay minerals are shown in Figure4, which exhibit similar patterns for different peak loads (5, 10, 30, and 50 mN). The unloading slope of quartz was low and the irrecoverable deformation was small, whereas the unloading slope of the clay minerals was high and showed significant irrecoverable deformation. This indicates that the indentation process of quartz and clay minerals was dominated by elastic and plastic work, respectively [24]. The indentations in quartz had a lower maximum depth, while the indentations in the clay minerals had larger maximum depths at the same peak load. The creep displacement during holding period was larger for the clay minerals than quartz (Figure 4). Based on the Oliver–Pharr method, the average Young’s modulus and hardness of quartz were calculated to be 84.08–93.15 and 11.40–12.37 GPa under different peak loads, respectively (Table 1). In contrast, those of the clay minerals are much lower and 25.07–33.97 and 0.67–0.92 GPa, respectively (Table 1). The calculated mechanical parameters are slightly lower than those obtained in previous studies [24, 41]. A possible reason for this is that mechanical parameters are commonly measured using a few seconds of holding time, but a holding time of 120 s was used in this study. A longer holding time could lead to a decrease in the determined mechanical parameters, as verified in previous studies [2, 4]. The Young’s modulus of the clay minerals was more variable than that of quartz, because clay minerals are strongly heterogeneous, anisotropic [42–44] and commonly mixed with variable amounts of organic matter in shale [12, 24]. In general, the quartz and clay minerals in the shale can also be identified by their load–displacement curves and measured mechanical parameters.Figure 4
Typical load–displacement curves for quartz and clay minerals in shale at various peak loads.Table 1
Indentation data related to quartz and minerals in shale.
MineralsPmax (mN)h0 (nm)hmax (nm)Creep depth (nm)Hardness (GPa)Young’s modulus (GPa)Ave.Std.Ave.Std.Ave.Std.Ave.Std.Ave.Std.Clay minerals5555.361.2610.556.155.15.60.670.1225.644.4110714.827.1789.130.974.311.50.820.0526.847.79301192.787.21299.082.1106.311.60.920.1725.072.14501580.271.11720.662.9140.429.40.800.0633.970.62Quartz5156.83.6170.54.813.71.311.720.3893.154.0510240.51.2255.51.815.11.111.400.0686.081.6430428.11.8455.92.927.82.011.930.4789.481.8450574.63.4604.95.930.32.712.370.2584.081.42Pmax: peak load; h0: the indentation depth at the beginning of hold time; hmax: maximum indentation depth.


### 3.3. Creep Rate Change and Depth

The loading–unloading curves (Figure4) can be used to examine the creep of the quartz and clay minerals, but cannot directly provide any quantitative information. To quantitatively assess the creep effect, the creep rates and depths were calculated from the indentations performed with 120 s holding times at various peak loads. Creep rates were calculated as:
(4)h′t=ht+Δt‐ht‐Δt2Δtwhere h’(t) is the creep rate at time t, Δt is the time step, h(t + Δt) is the indentation depth at time t + Δt, and h(t–Δt) is the indentation depth at time t–Δt. A time step of 120 s was used to calculate the creep rate in this study.Figure5 shows representative creep depth and creep rate curves for the quartz and clay minerals during a holding time of 120 s under different peak loads. The creep depth of the clay minerals is larger than that of quartz for a specific peak load. Moreover, a larger peak load resulted in a higher initial creep rate for both these constituents in the shale, while the initial creep rate of the clay minerals was also found to be higher than that of quartz for a specific peak load (e.g., 90 and 30 nm/s for clay minerals and quartz at 30 mN, respectively). However, the trends of the creep rate were similar, and the creep rate decreased rapidly in the first ~20 s, and then it started to maintain a relatively steady value, which was independent of peak load.Figure 5
Examples of creep rate and depth for quartz and clay minerals at different peak loads.The creep of materials can be divided into three stages (i.e., primary, secondary, and tertiary), based on the creep rates. The primary creep stage has a high creep rate that decreases over time, while the secondary creep stage has a relatively steady creep rate. In the tertiary creep stage, the creep deformation accelerates to failure. A holding time of 120 s may not be long enough to allow the creep to enter the tertiary creep stage, and so only two creep stages (i.e., primary creep [first ~20 s] and secondary creep [the following ~100 s]) were observed (Figure5).The variations in creep depths with peak loads for the quartz and clay minerals are shown in Figure6. The creep depths of the clay minerals are obviously higher than those of quartz for each specific peak load, which is consistent with previous data obtained using the deconvolution method [4, 14]. In addition, both the depth–load curves show an increase of creep depth with peak load (Figure 6), from 55.1–140.4 nm for clay minerals and 13.7–30.3 nm for quartz at 5–50 mN (Table 1). The larger and faster growth of creep depth in clay minerals demonstrates that they are more sensitive to peak loads. The larger error bars for the creep depth of the clay minerals at a specific peak load are also consistent with their strongly heterogeneous nature.Figure 6
Creep depths of quartz and clay minerals at various peak loads. Error bars are one standard deviation.


### 3.4. Application of Viscoelastic Models for Quartz and Clay Minerals

The viscoelastic properties of rocks can be expressed by the creep function. Thus, numerous creep element models have been developed to represent the viscoelastic behavior of rocks. In these creep models, different viscoelastic properties are accounted by connecting spring and dashpot elements in parallel and/or in series. A spring is considered to be an elastic element (i.e., obeys Hooke’s law), while the dashpot is a perfectly viscous element (i.e., obeys Newton’s law). The strain of a spring element is completely recovered after the release of a load, while the strain of a dashpot element is permanent [2, 23]. The traditional linear viscoelastic models, which consist of spring and dashpot elements, such as the Maxwell, Kelvin, and Burger models, have been widely used to represent the creep behavior of materials. Therefore, it is necessary to select the appropriate viscoelastic model to describe the actual material properties.In this study, three viscoelastic models (i.e., the three-element Voigt, Burger’s, and two-dashpot Kelvin models) were chosen to fit the creep curves during the holding time, which can recover the creep time constants.


#### 3.4.1. Three-Element Voigt Model

The three-element Voigt model is also known as the standard solid model, and consists of a Kelvin Voigt element and another spring in series (Figure7(a)). Its constitutive equation is:
(5)εt=σ01E1+1E21‐exp‐E2t/η1where ε(t) is the time-dependent strain, σ0 is a constant stress, E1 and E2 are the moduli that represent the elastic behavior, and η1 is a creep time constant that represents the viscosity of the material. Previous studies have combined the three-element Voigt model with the instrumented indentation technique using an observational approach [2, 23, 45, 46], and the indentation depth of this model can be represented by:
(6)h2t=π2P0cotα1E1+1E21‐exp‐E2t/η1where P0 is the peak load, h(t) is the indentation depth, and α is the equivalent cone semi-angle (70.3° for a Berkovich indenter). By regressing the depth–time data during the holding time, we can calculate the parameters E1, E2, and η1. Equation Equation (6) was simplified to group the linear terms into the minimal number of fitting parameters and to only reflect the constant load part of the nanoindentation test [2, 23, 45, 46], which leads to:
(7)h2=B0+B11‐exp‐t/τwhere B0 is the elastic deformation after the loading stage, B1 is the primary creep coefficient, and τ is the creep time constant.Figure 7
Viscoelastic models constructed from spring and dashpot elements to account for the viscoelastic behavior of quartz and clay minerals: (a) three-element Voigt model, (b) Burger’s model, and (c) two-dashpot Kelvin model.
(a)(b)(c)Figure8 shows the variations in the regression parameters of the three-element Voigt model with peak load. This model provided a good match for the creep behavior of both the quartz and clay minerals, with R2>0.940 (Table 2). The instantaneous elastic modulus E1 for quartz decreases with increasing peak load, suggesting quartz behaved relatively less elastically at the beginning of the holding time, whereas E2 and the viscosity η1 show the opposite trend. In contrast, the regression parameters appeared to be less affected by the variations in peak load for clay minerals and exhibit a slight increasing trend with increasing peak load, similar to the silica gel phase of carbonated wollastonite [23]. In general, the R2 of all the viscoelastic parameters (E1, E2, and η1) for quartz are larger than those of the clay minerals. E1, E2, and η1 values of quartz are also much larger than those of the clay minerals.Figure 8
Linear regressions for employing the three-element Voigt model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, and (c) η1. Error bars are one standard deviation.
(a)(b)(c)Table 2
Fitting parameters obtained by different viscoelastic models of typical creep curves at 30 mN.
ModelThree element VoigtBurgersTwo dashpot kelvinParameterE1E2η1R2E1E2η1η2R2E1E2E3η1η2R2(GPa)(GPa)(GPa·s)(GPa)(GPa)(GPa·s)(GPa·s)(GPa)(GPa)(GPa)(GPa·s)(GPa·s)Quartz88402988120.9899718871272077970.9961924322936241060870.998Clay1113250350.94212162251129480.9841611391132988530.995


#### 3.4.2. Four-Element Burger’s Model

Burger’s model is a widely used viscoelastic model for a variety of materials, including bone [47], carbonated cement matrixes [23], cement pastes [48], and shale [2, 4]. It consists of four linear elements (Figure 7(b)) and its constitutive equation is:
(8)εt=σ01E1+1E21‐exp‐E2t/η2+1η1twhere E1 and E2 are the moduli, η1 is the long-term creep viscosity, and η2 is the creep time constant. Similar to the three-element Voigt model, the indentation depth of this model can be represented by:
(9)h2t=π2P0cotα1E1+1E21‐exp‐E2t/η2+tη1.Equation (9) was simplified to group the linear terms into the minimal number of parameters, obtaining:
(10)h2=B0+B11‐exp‐t/τ+B2t.Equation (10) was used to obtain the regression parameters of Burger’s model from the depth–time data obtained by nanoindentation. The Burger’s model fitted the depth–time data better than the three-element Voigt model. The average R2 values of the quartz and clay minerals for the Burger’s models are 0.996 and 0.984, respectively (Table 2). Figure 9 shows the variations in the regression parameters obtained from Burger’s models with peak loads. The elastic moduli (E1 and E2) decrease with increasing peak loads for quartz, although the decrease of E2 is more pronounced. The data exhibit a relatively good linear correlation (R2>0.730), which suggests the quartz became relatively less elastic with increasing peak load. The long-term creep viscosity η1 increased with increasing peak load for quartz (Figure 9(d)). Similar to the three-element Voigt model, all the regression parameters exhibit a positive relationship with peak load for the clay minerals, E2 and η2 have relatively good linear correlation (R2>0.820). In general, all the viscoelastic parameters for quartz are much larger than those of the clay minerals.Figure 9
Linear regressions for employing Burger’s model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, (c) η1, and (d) η2. Error bars are one standard deviation.
(a)(b)(c)(d)


#### 3.4.3. Two-Dashpot Kelvin Model

The two-dashpot Kelvin model consists of two viscous and three elastic elements, and it has five constants (Figure7(c)). Its constitutive equation is:
(11)εt=σ01E3+1E11‐exp‐E1t/η1+1E21‐exp‐E2t/η2where E1, E2, and E3 are the moduli, and η1 and η2 are the creep time constants. Like the first two models, the indentation depth of this model can be represented by:
(12)h2t=π2P0cotα1E3+1E11‐exp‐E1t/η1+1E21‐exp‐E2t/η2.Similar to the first two models, Eq. (12) can be simplified as:
(13)h2=B0+B11‐exp‐t/τ1+B21‐exp‐t/τ2.Equation (13) was fitted to the depth–time data to determine the regression parameters. Compared with the first two models, the two-dashpot Kelvin model has the best fit to the experimental data. The average R2 values for the quartz and clay minerals are 0.998 and 0.995, respectively (Table 2). Figure 10 shows the variations of the regression parameters from this model with the peak loads. Similar to Burger’s model, E1 and E3 decrease with increasing peak load and exhibit good linear relationships for quartz (R2>0.770), whereas E2 has the opposite trend. Like the first two models, the regression parameters except for η2 have a positive relationship with peak load for the clay minerals. Of them, E1 and η1 have the best linear correlations (R2>0.920). The regression parameters for quartz are also larger than those of the clay minerals. The fitting parameters are variably affected by the peak load.Figure 10
Linear regressions for employing the two-dashpot Kelvin model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, (c) E3, (d) η1, and (e) η2. Error bars are one standard deviation.
(a)(b)(c)(d)(e)


#### 3.4.4. Comparison of the Viscoelastic Models

The regression curves at different peak loads are similar. Typical regressions of depth–time data for a holding time of 120 s at a peak load of 30 mN are shown in Figure11. The R2 values are more than 0.940 for all three models for the quartz and clay minerals. The three-element Voigt model has the worst fit in the initial holding stage, particularly for clay minerals (Figure 11). This phenomenon can also be observed on the fitting to the creep data of bulk shale [2], the silica gel phase of carbonated wollastonite [23] and bone [46] in nanoindentation tests. The possible reason is that in the initial holding stage, there are much elastic energy near the indenter that cannot be released in time, which leads to the large creep displacement and obvious creep deformation. The two-dashpot Kelvin model has the best fit of the three models. Compared with the other two models, it consists of five fitting parameters, and has the highest calculation cost for the regression of the nanoindentation data. However, some of the depth–time curves could not be fitted or the obtained regression parameters were unreasonable for this model, similar to the result of extended viscoelastic model with two Kelvin models in series that was used to fit quartz nanoindentation creep data [49]. Consequently, it appears that this model is unstable. Burger’s model provides the second-best fitting results and reasonable regression parameters. This suggests that this model is not only suitable for a complex material like shale [2, 4], but also quartz and clay minerals in shale.Figure 11
Regression curves obtained from different viscoelastic models at 30 mN for (a) quartz and (b) clay minerals.
(a)(b)Regardless of the suitability of these models, we can obtain some specific rules from Burger’s model of the viscoelasticity:(1)
The instantaneous elastic modulus (E1), the modulus E2 that represents the long-term reversible creep, and creep time constant (η2) decreased with increasing peak loads for quartz, whereas they exhibited an opposite trend for clay minerals(2)
The long-term creep viscosity (η1) increased with increasing peak loads for both the quartz and clay minerals(3)
Quartz and clay minerals in shale have different viscoelastic properties. The viscoelastic parameters of quartz are higher than those of clay minerals, and independent of peak loadThe viscoelastic parameters (E1, E2, and η2) decrease with increasing peak loads for quartz. This is similar to the mechanical parameters (i.e., Young’s modulus) that exhibit a slight decrease at higher indentation peak loads for quartz [24, 50]. This downward trend is similar to the results obtained from a hard film on a soft substrate in the thin film-substrate systems (Figure 12) [50–52]. In general, no more than 10% of hmax of the film thickness is used to avoid the effect of the underlying substrate on the determination of thin film properties. However, this rule is not necessarily correct for measuring the elastic modulus of certain materials or composites [50, 52, 53]. Moreover, for the highly heterogeneous shale, creep deformation is not just from the deformation of the mineral directly in contact with the indenter, but also from the time-dependent strains that develop in the area around the indenter [54]. As such, even though the indent acts on a single phase, with increasing peak load, the elastic zone expands from a harder phase (quartz) into a softer phase (clay minerals) (Figure 12), the overall mechanical response will become soft, leading to a decrease in these viscoelastic parameters.Figure 12
Schematic diagram of the nanoindentation performed on quartz and clay minerals in shale [50].In contrast, these viscoelastic parameters increased with increasing peak loads for clay minerals, which is similar to the increase in Young’s modulus of clay matrix with increasing depth or peak load [55]. One possible reason may be due to grain compaction and rearrangement [56, 57] and grain motion by frictional sliding [4] of clay minerals under larger peak loads, which enhances the creep deformation. Another possible explanation is that the indent locating on clay minerals in shale is similar to that of a soft film on a hard substrate in the thin film-substrate systems, which is opposite from that of quartz. When the elastic zone expands from clay minerals into quartz or other hard phase (Figure 12), the overall mechanical response will be stiff, which may lead to an increase of these viscoelastic parameters. The η1 value increased with increasing peak loads for both the quartz and clay minerals, although the increase for quartz is more pronounced, suggesting the long-term creep behavior of them behaved relatively more viscously during the holding time. The reason for this is likely to be the different crystal structures between quartz and clay minerals.The viscoelastic parameters (E1 and E2) of the clay minerals are less than those of quartz. They are calculated by nonlinear curve fitting, and the initial values in Eqs (7), (10), and (13) can be estimated from the initial (h0) and final depth (hmax) [46]. For example, B0 should be equal to the indentation depth (h0) at the beginning of the holding time, and B1 should be equal to the creep depth (Δh). Because of the compact structure and high stiffness of quartz, h0 and Δh of quartz are both much smaller than for clay minerals under the same peak load (Table 1), which will result in the larger E1 and E2 values for quartz. The creep time constant τ (E2/η2) is also much larger for the clay minerals than quartz under the same conditions (i.e., τ is around 6 s for clay minerals and 0.5 s for quartz at 10 mN from Burger’s model), which demonstrates that the clay minerals are more viscous [4]. However, the time constants (η1 and η2) are not easily estimated, and the relatively small time constants for the clay minerals were possibly controlled by the creep mechanisms, which may be due to contact interactions between single plate-like clay grains [58]. Compared to a previous study [2], the viscoelastic parameters of shale are between those of clay minerals and quartz (i.e., E1 values of shale, quartz, and clay minerals are ~25, ~100, and~10 GPa), indicating that clay minerals probably facilitate viscoelastic deformation, while quartz may inhibit viscoelastic deformation.


## 3.1. Identification of Quartz and Clay Minerals under the Optical Microscope

To conduct the in situ nanoindentation experiments on the minerals in the shale, quartz and clay minerals must firstly be identified under an optical microscope. These minerals can be clearly distinguished based on their different morphological and optical characteristics, combined with their reflection features and particle shapes.In the Longmaxi Formation shale, quartz grains are translucent with small impurities, and have dark gray colors under reflected light (Figures2(a) and 2(b)). They have irregular shapes with straight edges. The grain size is mostly between 15–40 μm. The quartz grains subjected to mechanical polishing have a smoother and cleaner surface than those subjected to ion polishing, which is the most distinct morphological feature (Figure 2(a)–2(c)). Thus, quartz grains subjected to mechanical polishing were chosen for the in situ nanoindentation experiments. The quartz can be divided into authigenic and detrital quartz [26]. However, the authigenic quartz (mostly <2 μm in size) has a similar color as the other dark gray minerals [26], which makes it difficult to identify.Figure 2
Photomicrographs of quartz and clay minerals. The surface morphology of quartz and clay minerals is shown after (a–b) mechanical polishing and (c–f) ion polishing. Quartz is dark gray, translucent, contains small impurities, and has straight and sharp edges without a regular shape. The clay minerals are opaque and black after mechanical polishing, and located in pit areas as lumps or isolated grains along with organic matter. Magnifications: (a) and (c) 20× objective lens; (b) and (d) 50× objective lens; (e–f) 100× objective lens (the eyepiece magnification was 40× for all images).After mechanical polishing, the clay minerals in the shale are opaque and black under the optical microscope. The clay minerals infill pitted areas as lumps or dispersed grains (Figures2(a) and 2(b)), and may have no obvious crystal form after ion polishing (Figure 2(e)). This texture is formed by the dense stacking of flaky clay particles with sizes of tens of nanometers [38, 39]. In general, the clay minerals are continuously distributed and form the mineral framework of the shale. The surfaces of the clay minerals are not flat and cannot be focused on after mechanical polishing, whereas the surface height between quartz and clay minerals is nearly the same after ion polishing (Figure 2(a)–2(d)), which was reported in our previous study [24]. Therefore, unlike quartz, the nanoindentation tests on the clay minerals were carried out on ion-polished shale slices. It should be noted that organic matter is interspersed within the clay grains (Figure 2(d)–2(f)) and has an off-white color, which was also observed in our previous study [24]. Clay minerals mixed with organic matter occur as aggregates in the shale matrix [12] (Figure 2(c)–2(f)). Given that the clay grains are very small, it is not possible to determine the independent creep response of individual clay particles. Thus, the organic–clay matrix was analyzed as a whole. Similarly, clay minerals in the Longmaxi Formation shale also have different origins, including terrigenous clastic, authigenic, and secondary types [40], which are difficult to distinguish under an optical microscope.Although our method is unable to identify the origins of the quartz and clay minerals, it is suitable for the in situ nanoindentation tests. In addition to quartz and clay minerals, the Longmaxi Formation shale also contains other minerals, such as pyrite, carbonates, and feldspar. These minerals are low in content, and have very small grain sizes, especially for pyrite (Figure2(e)–2(f)), thus were excluded from discussion in this study.


## 3.2. Indentation Responses of Quartz and Clay Minerals

In theory, a larger mineral crystal is more suitable for the nanoindentation tests, because it is not readily affected by adjacent phases or the substrate. Due to the small size of the minerals in the shale, representative quartz grains larger than 20μm and an area of clay minerals as large as possible were chosen for the in situ nanoindentation experiments. Figure 3 shows optical microscopy images of the quartz and clay minerals before and after indentation at peak loads of 5–50 mN. The residual imprints indicate that the larger peak loads produce a greater indentation imprint on both the quartz and clay minerals, and the indentations are unclear at a peak load of 5 mN. The reason for this is that a larger peak load causes a larger contact area on the mineral surface. In addition, the sizes of the residual imprints on the clay minerals is larger than those on quartz under the same peak load conditions, which is consistent with the higher strength and stiffness of quartz.Figure 3
Photomicrographs of a quartz grain (a) before and (b) after nanoindentation tests, and clay minerals (c) before and (d–f) after nanoindentation tests.The corresponding load–displacement curves for the quartz and clay minerals are shown in Figure4, which exhibit similar patterns for different peak loads (5, 10, 30, and 50 mN). The unloading slope of quartz was low and the irrecoverable deformation was small, whereas the unloading slope of the clay minerals was high and showed significant irrecoverable deformation. This indicates that the indentation process of quartz and clay minerals was dominated by elastic and plastic work, respectively [24]. The indentations in quartz had a lower maximum depth, while the indentations in the clay minerals had larger maximum depths at the same peak load. The creep displacement during holding period was larger for the clay minerals than quartz (Figure 4). Based on the Oliver–Pharr method, the average Young’s modulus and hardness of quartz were calculated to be 84.08–93.15 and 11.40–12.37 GPa under different peak loads, respectively (Table 1). In contrast, those of the clay minerals are much lower and 25.07–33.97 and 0.67–0.92 GPa, respectively (Table 1). The calculated mechanical parameters are slightly lower than those obtained in previous studies [24, 41]. A possible reason for this is that mechanical parameters are commonly measured using a few seconds of holding time, but a holding time of 120 s was used in this study. A longer holding time could lead to a decrease in the determined mechanical parameters, as verified in previous studies [2, 4]. The Young’s modulus of the clay minerals was more variable than that of quartz, because clay minerals are strongly heterogeneous, anisotropic [42–44] and commonly mixed with variable amounts of organic matter in shale [12, 24]. In general, the quartz and clay minerals in the shale can also be identified by their load–displacement curves and measured mechanical parameters.Figure 4
Typical load–displacement curves for quartz and clay minerals in shale at various peak loads.Table 1
Indentation data related to quartz and minerals in shale.
MineralsPmax (mN)h0 (nm)hmax (nm)Creep depth (nm)Hardness (GPa)Young’s modulus (GPa)Ave.Std.Ave.Std.Ave.Std.Ave.Std.Ave.Std.Clay minerals5555.361.2610.556.155.15.60.670.1225.644.4110714.827.1789.130.974.311.50.820.0526.847.79301192.787.21299.082.1106.311.60.920.1725.072.14501580.271.11720.662.9140.429.40.800.0633.970.62Quartz5156.83.6170.54.813.71.311.720.3893.154.0510240.51.2255.51.815.11.111.400.0686.081.6430428.11.8455.92.927.82.011.930.4789.481.8450574.63.4604.95.930.32.712.370.2584.081.42Pmax: peak load; h0: the indentation depth at the beginning of hold time; hmax: maximum indentation depth.


## 3.3. Creep Rate Change and Depth

The loading–unloading curves (Figure4) can be used to examine the creep of the quartz and clay minerals, but cannot directly provide any quantitative information. To quantitatively assess the creep effect, the creep rates and depths were calculated from the indentations performed with 120 s holding times at various peak loads. Creep rates were calculated as:
(4)h′t=ht+Δt‐ht‐Δt2Δtwhere h’(t) is the creep rate at time t, Δt is the time step, h(t + Δt) is the indentation depth at time t + Δt, and h(t–Δt) is the indentation depth at time t–Δt. A time step of 120 s was used to calculate the creep rate in this study.Figure5 shows representative creep depth and creep rate curves for the quartz and clay minerals during a holding time of 120 s under different peak loads. The creep depth of the clay minerals is larger than that of quartz for a specific peak load. Moreover, a larger peak load resulted in a higher initial creep rate for both these constituents in the shale, while the initial creep rate of the clay minerals was also found to be higher than that of quartz for a specific peak load (e.g., 90 and 30 nm/s for clay minerals and quartz at 30 mN, respectively). However, the trends of the creep rate were similar, and the creep rate decreased rapidly in the first ~20 s, and then it started to maintain a relatively steady value, which was independent of peak load.Figure 5
Examples of creep rate and depth for quartz and clay minerals at different peak loads.The creep of materials can be divided into three stages (i.e., primary, secondary, and tertiary), based on the creep rates. The primary creep stage has a high creep rate that decreases over time, while the secondary creep stage has a relatively steady creep rate. In the tertiary creep stage, the creep deformation accelerates to failure. A holding time of 120 s may not be long enough to allow the creep to enter the tertiary creep stage, and so only two creep stages (i.e., primary creep [first ~20 s] and secondary creep [the following ~100 s]) were observed (Figure5).The variations in creep depths with peak loads for the quartz and clay minerals are shown in Figure6. The creep depths of the clay minerals are obviously higher than those of quartz for each specific peak load, which is consistent with previous data obtained using the deconvolution method [4, 14]. In addition, both the depth–load curves show an increase of creep depth with peak load (Figure 6), from 55.1–140.4 nm for clay minerals and 13.7–30.3 nm for quartz at 5–50 mN (Table 1). The larger and faster growth of creep depth in clay minerals demonstrates that they are more sensitive to peak loads. The larger error bars for the creep depth of the clay minerals at a specific peak load are also consistent with their strongly heterogeneous nature.Figure 6
Creep depths of quartz and clay minerals at various peak loads. Error bars are one standard deviation.


## 3.4. Application of Viscoelastic Models for Quartz and Clay Minerals

The viscoelastic properties of rocks can be expressed by the creep function. Thus, numerous creep element models have been developed to represent the viscoelastic behavior of rocks. In these creep models, different viscoelastic properties are accounted by connecting spring and dashpot elements in parallel and/or in series. A spring is considered to be an elastic element (i.e., obeys Hooke’s law), while the dashpot is a perfectly viscous element (i.e., obeys Newton’s law). The strain of a spring element is completely recovered after the release of a load, while the strain of a dashpot element is permanent [2, 23]. The traditional linear viscoelastic models, which consist of spring and dashpot elements, such as the Maxwell, Kelvin, and Burger models, have been widely used to represent the creep behavior of materials. Therefore, it is necessary to select the appropriate viscoelastic model to describe the actual material properties.In this study, three viscoelastic models (i.e., the three-element Voigt, Burger’s, and two-dashpot Kelvin models) were chosen to fit the creep curves during the holding time, which can recover the creep time constants.


### 3.4.1. Three-Element Voigt Model

The three-element Voigt model is also known as the standard solid model, and consists of a Kelvin Voigt element and another spring in series (Figure7(a)). Its constitutive equation is:
(5)εt=σ01E1+1E21‐exp‐E2t/η1where ε(t) is the time-dependent strain, σ0 is a constant stress, E1 and E2 are the moduli that represent the elastic behavior, and η1 is a creep time constant that represents the viscosity of the material. Previous studies have combined the three-element Voigt model with the instrumented indentation technique using an observational approach [2, 23, 45, 46], and the indentation depth of this model can be represented by:
(6)h2t=π2P0cotα1E1+1E21‐exp‐E2t/η1where P0 is the peak load, h(t) is the indentation depth, and α is the equivalent cone semi-angle (70.3° for a Berkovich indenter). By regressing the depth–time data during the holding time, we can calculate the parameters E1, E2, and η1. Equation Equation (6) was simplified to group the linear terms into the minimal number of fitting parameters and to only reflect the constant load part of the nanoindentation test [2, 23, 45, 46], which leads to:
(7)h2=B0+B11‐exp‐t/τwhere B0 is the elastic deformation after the loading stage, B1 is the primary creep coefficient, and τ is the creep time constant.Figure 7
Viscoelastic models constructed from spring and dashpot elements to account for the viscoelastic behavior of quartz and clay minerals: (a) three-element Voigt model, (b) Burger’s model, and (c) two-dashpot Kelvin model.
(a)(b)(c)Figure8 shows the variations in the regression parameters of the three-element Voigt model with peak load. This model provided a good match for the creep behavior of both the quartz and clay minerals, with R2>0.940 (Table 2). The instantaneous elastic modulus E1 for quartz decreases with increasing peak load, suggesting quartz behaved relatively less elastically at the beginning of the holding time, whereas E2 and the viscosity η1 show the opposite trend. In contrast, the regression parameters appeared to be less affected by the variations in peak load for clay minerals and exhibit a slight increasing trend with increasing peak load, similar to the silica gel phase of carbonated wollastonite [23]. In general, the R2 of all the viscoelastic parameters (E1, E2, and η1) for quartz are larger than those of the clay minerals. E1, E2, and η1 values of quartz are also much larger than those of the clay minerals.Figure 8
Linear regressions for employing the three-element Voigt model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, and (c) η1. Error bars are one standard deviation.
(a)(b)(c)Table 2
Fitting parameters obtained by different viscoelastic models of typical creep curves at 30 mN.
ModelThree element VoigtBurgersTwo dashpot kelvinParameterE1E2η1R2E1E2η1η2R2E1E2E3η1η2R2(GPa)(GPa)(GPa·s)(GPa)(GPa)(GPa·s)(GPa·s)(GPa)(GPa)(GPa)(GPa·s)(GPa·s)Quartz88402988120.9899718871272077970.9961924322936241060870.998Clay1113250350.94212162251129480.9841611391132988530.995


### 3.4.2. Four-Element Burger’s Model

Burger’s model is a widely used viscoelastic model for a variety of materials, including bone [47], carbonated cement matrixes [23], cement pastes [48], and shale [2, 4]. It consists of four linear elements (Figure 7(b)) and its constitutive equation is:
(8)εt=σ01E1+1E21‐exp‐E2t/η2+1η1twhere E1 and E2 are the moduli, η1 is the long-term creep viscosity, and η2 is the creep time constant. Similar to the three-element Voigt model, the indentation depth of this model can be represented by:
(9)h2t=π2P0cotα1E1+1E21‐exp‐E2t/η2+tη1.Equation (9) was simplified to group the linear terms into the minimal number of parameters, obtaining:
(10)h2=B0+B11‐exp‐t/τ+B2t.Equation (10) was used to obtain the regression parameters of Burger’s model from the depth–time data obtained by nanoindentation. The Burger’s model fitted the depth–time data better than the three-element Voigt model. The average R2 values of the quartz and clay minerals for the Burger’s models are 0.996 and 0.984, respectively (Table 2). Figure 9 shows the variations in the regression parameters obtained from Burger’s models with peak loads. The elastic moduli (E1 and E2) decrease with increasing peak loads for quartz, although the decrease of E2 is more pronounced. The data exhibit a relatively good linear correlation (R2>0.730), which suggests the quartz became relatively less elastic with increasing peak load. The long-term creep viscosity η1 increased with increasing peak load for quartz (Figure 9(d)). Similar to the three-element Voigt model, all the regression parameters exhibit a positive relationship with peak load for the clay minerals, E2 and η2 have relatively good linear correlation (R2>0.820). In general, all the viscoelastic parameters for quartz are much larger than those of the clay minerals.Figure 9
Linear regressions for employing Burger’s model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, (c) η1, and (d) η2. Error bars are one standard deviation.
(a)(b)(c)(d)


### 3.4.3. Two-Dashpot Kelvin Model

The two-dashpot Kelvin model consists of two viscous and three elastic elements, and it has five constants (Figure7(c)). Its constitutive equation is:
(11)εt=σ01E3+1E11‐exp‐E1t/η1+1E21‐exp‐E2t/η2where E1, E2, and E3 are the moduli, and η1 and η2 are the creep time constants. Like the first two models, the indentation depth of this model can be represented by:
(12)h2t=π2P0cotα1E3+1E11‐exp‐E1t/η1+1E21‐exp‐E2t/η2.Similar to the first two models, Eq. (12) can be simplified as:
(13)h2=B0+B11‐exp‐t/τ1+B21‐exp‐t/τ2.Equation (13) was fitted to the depth–time data to determine the regression parameters. Compared with the first two models, the two-dashpot Kelvin model has the best fit to the experimental data. The average R2 values for the quartz and clay minerals are 0.998 and 0.995, respectively (Table 2). Figure 10 shows the variations of the regression parameters from this model with the peak loads. Similar to Burger’s model, E1 and E3 decrease with increasing peak load and exhibit good linear relationships for quartz (R2>0.770), whereas E2 has the opposite trend. Like the first two models, the regression parameters except for η2 have a positive relationship with peak load for the clay minerals. Of them, E1 and η1 have the best linear correlations (R2>0.920). The regression parameters for quartz are also larger than those of the clay minerals. The fitting parameters are variably affected by the peak load.Figure 10
Linear regressions for employing the two-dashpot Kelvin model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, (c) E3, (d) η1, and (e) η2. Error bars are one standard deviation.
(a)(b)(c)(d)(e)


### 3.4.4. Comparison of the Viscoelastic Models

The regression curves at different peak loads are similar. Typical regressions of depth–time data for a holding time of 120 s at a peak load of 30 mN are shown in Figure11. The R2 values are more than 0.940 for all three models for the quartz and clay minerals. The three-element Voigt model has the worst fit in the initial holding stage, particularly for clay minerals (Figure 11). This phenomenon can also be observed on the fitting to the creep data of bulk shale [2], the silica gel phase of carbonated wollastonite [23] and bone [46] in nanoindentation tests. The possible reason is that in the initial holding stage, there are much elastic energy near the indenter that cannot be released in time, which leads to the large creep displacement and obvious creep deformation. The two-dashpot Kelvin model has the best fit of the three models. Compared with the other two models, it consists of five fitting parameters, and has the highest calculation cost for the regression of the nanoindentation data. However, some of the depth–time curves could not be fitted or the obtained regression parameters were unreasonable for this model, similar to the result of extended viscoelastic model with two Kelvin models in series that was used to fit quartz nanoindentation creep data [49]. Consequently, it appears that this model is unstable. Burger’s model provides the second-best fitting results and reasonable regression parameters. This suggests that this model is not only suitable for a complex material like shale [2, 4], but also quartz and clay minerals in shale.Figure 11
Regression curves obtained from different viscoelastic models at 30 mN for (a) quartz and (b) clay minerals.
(a)(b)Regardless of the suitability of these models, we can obtain some specific rules from Burger’s model of the viscoelasticity:(1)
The instantaneous elastic modulus (E1), the modulus E2 that represents the long-term reversible creep, and creep time constant (η2) decreased with increasing peak loads for quartz, whereas they exhibited an opposite trend for clay minerals(2)
The long-term creep viscosity (η1) increased with increasing peak loads for both the quartz and clay minerals(3)
Quartz and clay minerals in shale have different viscoelastic properties. The viscoelastic parameters of quartz are higher than those of clay minerals, and independent of peak loadThe viscoelastic parameters (E1, E2, and η2) decrease with increasing peak loads for quartz. This is similar to the mechanical parameters (i.e., Young’s modulus) that exhibit a slight decrease at higher indentation peak loads for quartz [24, 50]. This downward trend is similar to the results obtained from a hard film on a soft substrate in the thin film-substrate systems (Figure 12) [50–52]. In general, no more than 10% of hmax of the film thickness is used to avoid the effect of the underlying substrate on the determination of thin film properties. However, this rule is not necessarily correct for measuring the elastic modulus of certain materials or composites [50, 52, 53]. Moreover, for the highly heterogeneous shale, creep deformation is not just from the deformation of the mineral directly in contact with the indenter, but also from the time-dependent strains that develop in the area around the indenter [54]. As such, even though the indent acts on a single phase, with increasing peak load, the elastic zone expands from a harder phase (quartz) into a softer phase (clay minerals) (Figure 12), the overall mechanical response will become soft, leading to a decrease in these viscoelastic parameters.Figure 12
Schematic diagram of the nanoindentation performed on quartz and clay minerals in shale [50].In contrast, these viscoelastic parameters increased with increasing peak loads for clay minerals, which is similar to the increase in Young’s modulus of clay matrix with increasing depth or peak load [55]. One possible reason may be due to grain compaction and rearrangement [56, 57] and grain motion by frictional sliding [4] of clay minerals under larger peak loads, which enhances the creep deformation. Another possible explanation is that the indent locating on clay minerals in shale is similar to that of a soft film on a hard substrate in the thin film-substrate systems, which is opposite from that of quartz. When the elastic zone expands from clay minerals into quartz or other hard phase (Figure 12), the overall mechanical response will be stiff, which may lead to an increase of these viscoelastic parameters. The η1 value increased with increasing peak loads for both the quartz and clay minerals, although the increase for quartz is more pronounced, suggesting the long-term creep behavior of them behaved relatively more viscously during the holding time. The reason for this is likely to be the different crystal structures between quartz and clay minerals.The viscoelastic parameters (E1 and E2) of the clay minerals are less than those of quartz. They are calculated by nonlinear curve fitting, and the initial values in Eqs (7), (10), and (13) can be estimated from the initial (h0) and final depth (hmax) [46]. For example, B0 should be equal to the indentation depth (h0) at the beginning of the holding time, and B1 should be equal to the creep depth (Δh). Because of the compact structure and high stiffness of quartz, h0 and Δh of quartz are both much smaller than for clay minerals under the same peak load (Table 1), which will result in the larger E1 and E2 values for quartz. The creep time constant τ (E2/η2) is also much larger for the clay minerals than quartz under the same conditions (i.e., τ is around 6 s for clay minerals and 0.5 s for quartz at 10 mN from Burger’s model), which demonstrates that the clay minerals are more viscous [4]. However, the time constants (η1 and η2) are not easily estimated, and the relatively small time constants for the clay minerals were possibly controlled by the creep mechanisms, which may be due to contact interactions between single plate-like clay grains [58]. Compared to a previous study [2], the viscoelastic parameters of shale are between those of clay minerals and quartz (i.e., E1 values of shale, quartz, and clay minerals are ~25, ~100, and~10 GPa), indicating that clay minerals probably facilitate viscoelastic deformation, while quartz may inhibit viscoelastic deformation.


## 3.4.1. Three-Element Voigt Model

The three-element Voigt model is also known as the standard solid model, and consists of a Kelvin Voigt element and another spring in series (Figure7(a)). Its constitutive equation is:
(5)εt=σ01E1+1E21‐exp‐E2t/η1where ε(t) is the time-dependent strain, σ0 is a constant stress, E1 and E2 are the moduli that represent the elastic behavior, and η1 is a creep time constant that represents the viscosity of the material. Previous studies have combined the three-element Voigt model with the instrumented indentation technique using an observational approach [2, 23, 45, 46], and the indentation depth of this model can be represented by:
(6)h2t=π2P0cotα1E1+1E21‐exp‐E2t/η1where P0 is the peak load, h(t) is the indentation depth, and α is the equivalent cone semi-angle (70.3° for a Berkovich indenter). By regressing the depth–time data during the holding time, we can calculate the parameters E1, E2, and η1. Equation Equation (6) was simplified to group the linear terms into the minimal number of fitting parameters and to only reflect the constant load part of the nanoindentation test [2, 23, 45, 46], which leads to:
(7)h2=B0+B11‐exp‐t/τwhere B0 is the elastic deformation after the loading stage, B1 is the primary creep coefficient, and τ is the creep time constant.Figure 7
Viscoelastic models constructed from spring and dashpot elements to account for the viscoelastic behavior of quartz and clay minerals: (a) three-element Voigt model, (b) Burger’s model, and (c) two-dashpot Kelvin model.
(a)(b)(c)Figure8 shows the variations in the regression parameters of the three-element Voigt model with peak load. This model provided a good match for the creep behavior of both the quartz and clay minerals, with R2>0.940 (Table 2). The instantaneous elastic modulus E1 for quartz decreases with increasing peak load, suggesting quartz behaved relatively less elastically at the beginning of the holding time, whereas E2 and the viscosity η1 show the opposite trend. In contrast, the regression parameters appeared to be less affected by the variations in peak load for clay minerals and exhibit a slight increasing trend with increasing peak load, similar to the silica gel phase of carbonated wollastonite [23]. In general, the R2 of all the viscoelastic parameters (E1, E2, and η1) for quartz are larger than those of the clay minerals. E1, E2, and η1 values of quartz are also much larger than those of the clay minerals.Figure 8
Linear regressions for employing the three-element Voigt model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, and (c) η1. Error bars are one standard deviation.
(a)(b)(c)Table 2
Fitting parameters obtained by different viscoelastic models of typical creep curves at 30 mN.
ModelThree element VoigtBurgersTwo dashpot kelvinParameterE1E2η1R2E1E2η1η2R2E1E2E3η1η2R2(GPa)(GPa)(GPa·s)(GPa)(GPa)(GPa·s)(GPa·s)(GPa)(GPa)(GPa)(GPa·s)(GPa·s)Quartz88402988120.9899718871272077970.9961924322936241060870.998Clay1113250350.94212162251129480.9841611391132988530.995


## 3.4.2. Four-Element Burger’s Model

Burger’s model is a widely used viscoelastic model for a variety of materials, including bone [47], carbonated cement matrixes [23], cement pastes [48], and shale [2, 4]. It consists of four linear elements (Figure 7(b)) and its constitutive equation is:
(8)εt=σ01E1+1E21‐exp‐E2t/η2+1η1twhere E1 and E2 are the moduli, η1 is the long-term creep viscosity, and η2 is the creep time constant. Similar to the three-element Voigt model, the indentation depth of this model can be represented by:
(9)h2t=π2P0cotα1E1+1E21‐exp‐E2t/η2+tη1.Equation (9) was simplified to group the linear terms into the minimal number of parameters, obtaining:
(10)h2=B0+B11‐exp‐t/τ+B2t.Equation (10) was used to obtain the regression parameters of Burger’s model from the depth–time data obtained by nanoindentation. The Burger’s model fitted the depth–time data better than the three-element Voigt model. The average R2 values of the quartz and clay minerals for the Burger’s models are 0.996 and 0.984, respectively (Table 2). Figure 9 shows the variations in the regression parameters obtained from Burger’s models with peak loads. The elastic moduli (E1 and E2) decrease with increasing peak loads for quartz, although the decrease of E2 is more pronounced. The data exhibit a relatively good linear correlation (R2>0.730), which suggests the quartz became relatively less elastic with increasing peak load. The long-term creep viscosity η1 increased with increasing peak load for quartz (Figure 9(d)). Similar to the three-element Voigt model, all the regression parameters exhibit a positive relationship with peak load for the clay minerals, E2 and η2 have relatively good linear correlation (R2>0.820). In general, all the viscoelastic parameters for quartz are much larger than those of the clay minerals.Figure 9
Linear regressions for employing Burger’s model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, (c) η1, and (d) η2. Error bars are one standard deviation.
(a)(b)(c)(d)


## 3.4.3. Two-Dashpot Kelvin Model

The two-dashpot Kelvin model consists of two viscous and three elastic elements, and it has five constants (Figure7(c)). Its constitutive equation is:
(11)εt=σ01E3+1E11‐exp‐E1t/η1+1E21‐exp‐E2t/η2where E1, E2, and E3 are the moduli, and η1 and η2 are the creep time constants. Like the first two models, the indentation depth of this model can be represented by:
(12)h2t=π2P0cotα1E3+1E11‐exp‐E1t/η1+1E21‐exp‐E2t/η2.Similar to the first two models, Eq. (12) can be simplified as:
(13)h2=B0+B11‐exp‐t/τ1+B21‐exp‐t/τ2.Equation (13) was fitted to the depth–time data to determine the regression parameters. Compared with the first two models, the two-dashpot Kelvin model has the best fit to the experimental data. The average R2 values for the quartz and clay minerals are 0.998 and 0.995, respectively (Table 2). Figure 10 shows the variations of the regression parameters from this model with the peak loads. Similar to Burger’s model, E1 and E3 decrease with increasing peak load and exhibit good linear relationships for quartz (R2>0.770), whereas E2 has the opposite trend. Like the first two models, the regression parameters except for η2 have a positive relationship with peak load for the clay minerals. Of them, E1 and η1 have the best linear correlations (R2>0.920). The regression parameters for quartz are also larger than those of the clay minerals. The fitting parameters are variably affected by the peak load.Figure 10
Linear regressions for employing the two-dashpot Kelvin model for quartz and clay minerals at various peak loads. Fitting parameters: (a) E1, (b) E2, (c) E3, (d) η1, and (e) η2. Error bars are one standard deviation.
(a)(b)(c)(d)(e)


## 3.4.4. Comparison of the Viscoelastic Models

The regression curves at different peak loads are similar. Typical regressions of depth–time data for a holding time of 120 s at a peak load of 30 mN are shown in Figure11. The R2 values are more than 0.940 for all three models for the quartz and clay minerals. The three-element Voigt model has the worst fit in the initial holding stage, particularly for clay minerals (Figure 11). This phenomenon can also be observed on the fitting to the creep data of bulk shale [2], the silica gel phase of carbonated wollastonite [23] and bone [46] in nanoindentation tests. The possible reason is that in the initial holding stage, there are much elastic energy near the indenter that cannot be released in time, which leads to the large creep displacement and obvious creep deformation. The two-dashpot Kelvin model has the best fit of the three models. Compared with the other two models, it consists of five fitting parameters, and has the highest calculation cost for the regression of the nanoindentation data. However, some of the depth–time curves could not be fitted or the obtained regression parameters were unreasonable for this model, similar to the result of extended viscoelastic model with two Kelvin models in series that was used to fit quartz nanoindentation creep data [49]. Consequently, it appears that this model is unstable. Burger’s model provides the second-best fitting results and reasonable regression parameters. This suggests that this model is not only suitable for a complex material like shale [2, 4], but also quartz and clay minerals in shale.Figure 11
Regression curves obtained from different viscoelastic models at 30 mN for (a) quartz and (b) clay minerals.
(a)(b)Regardless of the suitability of these models, we can obtain some specific rules from Burger’s model of the viscoelasticity:(1)
The instantaneous elastic modulus (E1), the modulus E2 that represents the long-term reversible creep, and creep time constant (η2) decreased with increasing peak loads for quartz, whereas they exhibited an opposite trend for clay minerals(2)
The long-term creep viscosity (η1) increased with increasing peak loads for both the quartz and clay minerals(3)
Quartz and clay minerals in shale have different viscoelastic properties. The viscoelastic parameters of quartz are higher than those of clay minerals, and independent of peak loadThe viscoelastic parameters (E1, E2, and η2) decrease with increasing peak loads for quartz. This is similar to the mechanical parameters (i.e., Young’s modulus) that exhibit a slight decrease at higher indentation peak loads for quartz [24, 50]. This downward trend is similar to the results obtained from a hard film on a soft substrate in the thin film-substrate systems (Figure 12) [50–52]. In general, no more than 10% of hmax of the film thickness is used to avoid the effect of the underlying substrate on the determination of thin film properties. However, this rule is not necessarily correct for measuring the elastic modulus of certain materials or composites [50, 52, 53]. Moreover, for the highly heterogeneous shale, creep deformation is not just from the deformation of the mineral directly in contact with the indenter, but also from the time-dependent strains that develop in the area around the indenter [54]. As such, even though the indent acts on a single phase, with increasing peak load, the elastic zone expands from a harder phase (quartz) into a softer phase (clay minerals) (Figure 12), the overall mechanical response will become soft, leading to a decrease in these viscoelastic parameters.Figure 12
Schematic diagram of the nanoindentation performed on quartz and clay minerals in shale [50].In contrast, these viscoelastic parameters increased with increasing peak loads for clay minerals, which is similar to the increase in Young’s modulus of clay matrix with increasing depth or peak load [55]. One possible reason may be due to grain compaction and rearrangement [56, 57] and grain motion by frictional sliding [4] of clay minerals under larger peak loads, which enhances the creep deformation. Another possible explanation is that the indent locating on clay minerals in shale is similar to that of a soft film on a hard substrate in the thin film-substrate systems, which is opposite from that of quartz. When the elastic zone expands from clay minerals into quartz or other hard phase (Figure 12), the overall mechanical response will be stiff, which may lead to an increase of these viscoelastic parameters. The η1 value increased with increasing peak loads for both the quartz and clay minerals, although the increase for quartz is more pronounced, suggesting the long-term creep behavior of them behaved relatively more viscously during the holding time. The reason for this is likely to be the different crystal structures between quartz and clay minerals.The viscoelastic parameters (E1 and E2) of the clay minerals are less than those of quartz. They are calculated by nonlinear curve fitting, and the initial values in Eqs (7), (10), and (13) can be estimated from the initial (h0) and final depth (hmax) [46]. For example, B0 should be equal to the indentation depth (h0) at the beginning of the holding time, and B1 should be equal to the creep depth (Δh). Because of the compact structure and high stiffness of quartz, h0 and Δh of quartz are both much smaller than for clay minerals under the same peak load (Table 1), which will result in the larger E1 and E2 values for quartz. The creep time constant τ (E2/η2) is also much larger for the clay minerals than quartz under the same conditions (i.e., τ is around 6 s for clay minerals and 0.5 s for quartz at 10 mN from Burger’s model), which demonstrates that the clay minerals are more viscous [4]. However, the time constants (η1 and η2) are not easily estimated, and the relatively small time constants for the clay minerals were possibly controlled by the creep mechanisms, which may be due to contact interactions between single plate-like clay grains [58]. Compared to a previous study [2], the viscoelastic parameters of shale are between those of clay minerals and quartz (i.e., E1 values of shale, quartz, and clay minerals are ~25, ~100, and~10 GPa), indicating that clay minerals probably facilitate viscoelastic deformation, while quartz may inhibit viscoelastic deformation.


## 4. Implications and Limitations

Shale is a heterogeneous composite material composed of various minerals. The viscoelastic parameters of bulk shale have been obtained using the three aforementioned viscoelastic models and nanoindentation creep datasets [2], and we have further obtained the viscoelastic parameters of quartz and clay minerals in shale. However, the relationship between the viscoelastic parameters of the mineral components and bulk shale is not fully understood and requires further research. For polymers, the viscoelastic constant measured by nanoindentation is consistent with that measured by macroscopic methods [59]. The ultimate goal of nanoindentation research on the microscale viscoelastic properties of shale was to compare with and predict the macroscopic viscoelastic properties. It is expected that the obtained viscoelastic parameters will not be affected by the peak loads or the loading protocol within a certain range. However, Figure 8–10 show that the viscoelastic parameters were variably influenced by the peak load. In previous studies, different peak loads, such as 4.8 mN [15, 57], 8 mN [21], 12.5 mN [22], and 50 mN [49], were used to characterize the deformation behavior of single mineral phases in shale or sandstone. Consequently, the effect of peak load on creep was not discussed in these studies. Our research indicated that a larger peak load produced a greater indentation imprint or deeper displacement. Therefore, it is likely that the indentation of quartz was affected by the underlying clay minerals or a pressure-bearing framework at a larger peak load. In some cases, clay minerals are not infinite continuous blocks, and appear to be influenced by stiffer mineral phases in shale. The E2 value from the three-element Voigt model and the two-dashpot Kelvin model at 50 mN are much greater than those of the other peak loads (Figure 8(b) and 10(b)). The η2 value from Burger’s model at 50 mN for quartz is approximately equal to that of the clay minerals (Figure 9(d)), which appears to be unreasonable. As such, the relatively small peak load data (at 5 or 10 mN) should be best for the determination of the viscoelastic parameters of the quartz and clay minerals. Furthermore, the relationship between the viscoelastic parameters obtained from nanoindentation and macroscopic methods should be explored. In addition, the influences of fluid and temperature were not considered here. Therefore, the determined viscoelastic parameters of quartz and clay minerals in shale do not fully represent those at the subsurface. It is reported that the friction coefficient between clay minerals and organic matter would increase with increasing temperature owing to the loss of water from the internal structure of the clay minerals, thus altering the whole creep behavior of the shale rock [57]. Moreover, the effect of temperature on micro-properties of rock cannot be ignored [60, 61], which will be more representative of the actual creep behavior in geological situation. Details of in situ shale reservoir conditions (i.e., temperature and fluids) in studying the minerals or bulk shale should be considered in the future.


## 5. Conclusions

To better understand the creep behavior of shale,in situ nanoindentation measurements of quartz and clay minerals in shale were conducted using a nanoindentation technique and high-resolution optical microscope. Our main conclusions are as follows.
(1)
A larger peak load produced greater indentation imprints on both quartz and clay minerals. However, the size of the residual imprints observed on clay minerals was larger than that of quartz for a specific peak load, which is consistent with the higher strength and stiffness of quartz(2)
The initial creep rates and creep depths of the clay minerals were higher than that of quartz. However, the creep rate trends were similar (i.e., the creep rate first fluctuated rapidly in the primary creep stage [first ~20 s], and then it started to reach a relatively steady value in the secondary creep stage [~100 s], which was independent of peak load. A holding time of 120 s was not long enough to identify the tertiary creep stage for the quartz and clay minerals in the Longmaxi shale(3)
The three-element Voigt, Burger’s, and two-dashpot Kelvin models can all produce good fits to the depth–time data for calculating the viscoelastic parameters of quartz and clay minerals. Given the poor fit in the initial holding stage of the three-element Voigt model and instability of the two-dashpot Kelvin model, Burger’s model provided the most reasonable fitting parameters to determine the viscoelastic behavior of quartz and clay minerals(4)
The fitting parameters of these models were affected by peak loads. Given the sizes and shapes of quartz and clay mineral grains, and their spatial arrangement in shale, relatively small peak loads are best for the determination of the viscoelastic parameters. The viscoelastic parameters of clay minerals were lower than those of quartz and the bulk shale, suggesting that clay minerals probably facilitate the viscoelastic deformation of shale



---

*Source: 2860077-2022-01-13.xml*