# Effects of Pore Structure on Sandstone Mechanical Properties Based on Micro-CT Reconstruction Model

Authors: Rui Song; Lifu Zheng; Yao Wang; Jianjun Liu
Journal: Advances in Civil Engineering
 (2020)

---

## Abstract

As porous, heterogeneous, and anisotropic material, the microscopic structure of the rock has a significant influence on its mechanical properties. Rare studies were devoted to this area using pore scale modeling and simulations. In this paper, different types of sandstones are imaged using micro-CT technology. The rock porosity is obtained by filtering, binarization, and threshold segmentation. The texture coefficient (TC) and the tortuosity of the rock skeleton are calculated by open source program, where the tortuosity of the rock skeleton is firstly used to characterize the microscopic structure of the rock. Combining with the rock mechanics parameters obtained in the laboratory, the simulation of uniaxial compression is performed on the reconstructed pore scale rock finite element mesh model by ANSYS software. Young’s modulus, compressive strength, yield strength, shear modulus, and other related parameters obtained by numerical simulation are adopted to determine the optimal representative volume element (RVE) size. Moreover, the effects of microscopic structure characteristics on the mechanical properties of the rock are studied quantitatively. The results indicate that the averaged von Mises stress distribution, displacement field, and plastic strain field of rocks show anisotropy and heterogeneity. The stress concentration and the X-shaped conjugate plastic shear zone are investigated. The samples of S1∼S4 reach the elastic limit and enters the plastic yield state, when the strain is about 0.5%. And the critical yield strain of samples S5300-1∼S5400-2 is about 1%. Then, the quantitative relationships between porosity, TC, tortuosity of rock skeleton and rock mechanics parameters of digital rock samples are established and analyzed. The tortuosity of the rock skeleton is highly correlated with the mechanical parameters of the rock, i.e., Young’s modulus (R2 = 0.95), compressive strength (R2 = 0.94), yield strength (R2 = 0.92), and shear modulus (R2 = 0.94), which is believed to be more feasible to reveal the impacts of the microstructure of the rock on its mechanical properties.


---

## Body

## 1. Introduction

Rock is a natural porous media with anisotropy, heterogeneity, and discontinuity, which is related to many engineering applications, e.g., resource exploitation, civil construction, and tunnel excavation [1, 2]. The microstructure of the rock has been regarded as one of the important factors for the formation of shear zone, stress concentration, and cracking, which will threaten the safety of geotechnical engineering. Consequently, it is of practical significance to investigate the effects of the rock microstructure on its mechanical properties [3].For the nontransparent rock with disorder and complex pore structure, the nondestructive and high-resolution CT technology has been used to acquire the 3D structure of the rock. The internal microstructure of the rock is quantitatively analyzed by the microscopic CT image [4–6]. In addition, the CT technology was also used to image the evolution of the rock microstructure under different loads. Argandoña et al. [7] studied the heterogeneity of the rock structure and pore distribution using the established 3D rock model by CT image processing and analyzed the evolution of the pore structure in the freezing-thawing cycles. Kawakata et al. [8] used the method to analyze the deformation and fracture modes of granite under uniaxial and triaxial conditions. Ge et al. [9] studied the microfracture process of the rock by using CT scanning technology and analyzed the deformation of micropores and microcracks under different stress conditions. Ueta et al. [10] investigated the shear bands in sandstones combing CT scanning and shear experiments. Bésuelle et al. [11] used the technology to analyze the strain localization characteristics and the evolution process of the clay rock under triaxial compression. Wang et al. [12] investigated the relationship between various damage parameters (e.g., crack spacing and area) and mechanical properties of concrete materials. Tian et al. [13] identified and extracted cracks from the micro-CT images of concrete and acquired the length, width, and area of cracks by using statistical principles to characterize the microcosmic damage and fracture mechanism of concrete. Li et al. [14] also used the technology to study the mechanical properties of sandstone soaked in acid solution and quantitatively analyzed the damage degree of sandstone.Furthermore, many scientific efforts have been made to reveal the relationships between the parameters acquired from the microscale images (e.g., porosity) and mechanical properties of the rock. Fjar et al. [15] used quantitative and qualitative methods to analyze the relationships between the rock microstructure and its mechanical properties. Al-Harthi et al. [16] found the uniaxial compressive strength and Young’s modulus of the Basalt and Tuff gradually decreased with the increase of porosity. Li and Wang [17] studied the effects of porosity on the mechanical properties and energy dissipation of dry and saturated porous sandstones by SHPB experiments. The results show that the strength of dry sandstone decreases and the energy dissipation increases with the increase of porosity. Howarth and Rowlands [18] firstly proposed the concept of texture coefficient (TC) to characterize the complexity of grains and put up with the mathematical expression. The TC of sandstone, igneous, and marble samples were calculated by electron microscope slices (10 mm × 10 mm), and a good correlation between TC and the mechanical properties were obtained. Ersoy and Waller [19] obtained the TC value of the artificial sandstone flakes (30 mm × 50 mm) and conducted the polycrystalline diamond tablet drilling tests. The results indicate the TC can be used to evaluate rock drillability, mechanical properties, and wear-resisting performance. Azzoni et al. [20] obtained the texture coefficients of several different lithologies, such as igneous rock, metamorphic rock, and sedimentary rock (80 mm × 160 mm), and found close connections between TC and compressive properties of granite, sandstone, and marble. Although the correlation between gneiss and rhyolite is not obvious, TC can be used as a tool to describe and classify different lithologies. Alber and Kahraman [21] obtained the compressive strength of the rock sample and the TC of the fault breccia slices (d = 7.6 mm–73 mm). They found the compressive strength decreases with the increase of TC (R2 = 0.83) when 0 < TC < 0.9, and the modulus of elasticity decreases as TC increases when 0.25 < TC < 0.5. Tandon and Gupta [22] studied the properties of different rocks in the Himalayas and calculated the texture coefficients of minerals, such as mica, feldspar, tourmaline, hornblende, pyroxene, and dolomite. They found the compressive strength increases with the increase of TC (R2 = 0.81) when 1.5 < TC < 3.5. Ozturk and Nasuf [23] calculated the TC values of rocks by using SEM images of samples (9 mm × 9 mm) and studied the relationship between TC and the uniaxial compressive strength of rocks.Scholars carried out the simulations on the reconstructed mesh model to study the mechanical properties of the rock by using micro-CT images. Chen et al. [24] reconstructed the three-dimensional rock model by using the two-dimensional rock image and studied the mechanical response and failure process under conventional Brazilian tensile test conditions by using FLAC3D. Zhu et al. [25] obtained the mineral compositions and extracted the pore distributions of rock by CT technology and simulated the cracking of the reconstructed rock model in the fluid-solid coupling process by RFPA. Yang et al. [26] and Zhu et al. [27] also used RFPA to simulate the deformation and fracture process of reconstructed rock models in uniaxial compression and Brazil disk test, respectively, and studied the impact of voids and fractures on the rock failure mode. Ju et al. [28] combined the CT imaging and X-ray diffraction with physical tests to obtain the pore structure, mineral composition, and mechanical properties of sandstone and established three-dimensional finite element mesh models with different porosity by using Mimics software. The pore scale deformation failure mechanism and energy dissipation of rock were investigated under stress wave propagation by numerical simulation. Tian et al. [29] and Wang [30] used the same modeling methods to simulate the uniaxial compressive and tensile properties of concrete and rock, respectively. Song et al. [31, 32] and Wang et al. [33] adopted rock microscopic CT imaging and structured the finite element mesh model reconstruction method to establish the finite element model in microscale for sandstone skeleton and pore, and the evolution of pore structure and permeability of the rock was carried out by numerical simulation. Wang [34] reconstructed rock microstructure models from CT images by using MATLAB software and analyzed the evolution process of rock failure by FLAC3D. Cheng et al. [35] conducted uniaxial compression simulations on the two-dimensional reconstructed rock models and found that the effective elasticity modulus decreases exponentially with the increase of porosity. Liu et al. [36] established three-dimensional finite element mesh models of synthetic sandstone rock samples and conducted the numerical simulation of uniaxial compression. The cohesive force and internal friction angle decrease linearly with the increase of porosity when the porosity is between 8% and 28%.In this paper, microscopic images of different types of sandstone samples are obtained by CT scanning technology. Combining with image processing, the porosity of rock samples is obtained. Open source programs are used to calculate the TC and tortuosity of the rock skeleton. The finite element mesh models with different porosities are reconstructed by model reconstruction algorithm. Based on the mechanical parameters of core scale, ANSYS software is used to carry out the numerical simulation of uniaxial compression. The calculated mechanical parameters (e.g., Young’s modulus, the compressive strength, and the yield strength) and microstructure parameters (e.g., TC and the tortuosity) of the rock are used to determine the optimal RVE. Moreover, the averaged von Mises stress distribution, displacement field, and plastic strain field of the rock are analyzed qualitatively. In addition, the relationship between the porosity, TC and the tortuosity of the rock skeleton, and mechanical parameters are used to quantitatively characterize mechanical properties of the rock.


## 2. Mathematical Model

### 2.1. Plasticity Criterion

Considering the influence of the microstructure on rock mechanical strength, the von Mises yield criterion has been widely applied in the microscale deformation-fracture process of the rock, which are adopted for determining the plastic deformation of the rock skeleton [37–39].The criterion states that yielding begins when the elastic energy of distortion reaches a critical value [40]. That is, when the second invariant of the stress deviator (J2′) reaches a critical value, it represents plastic status, which can be expressed as follows:(1)fσij′=J2′=C,where J2′ is the second invariant of the stress tensor and C is a constant related to material property and independent of the stress state, which can be measured experimentally.The second invariant of the stress tensor is as follows:(2)J2′=16σx−σy2+σy−σx2+σz−σx2+6τxy2+τyz2+τzx2=C,where σx, σy, and σz are normal stresses in different directions (Pa) and τxy, τxz, and τyz are shear stresses (Pa).Equation (2) can be defined as effective stress as follows:(3)σ¯=12σ1−σ22+σ2−σ32+σ3−σ12=3C,where σ¯ is effective stress (Pa) and σ1 and σ3 are principal stresses (Pa).Thus, the Mises yield condition can be represented by more convenient effective stress. It can be expressed as equation (4) in the form of the principal stress and simplified as shown in the following equation:(4)J2′=16σ1−σ22+σ2−σ32+σ3−σ12=C,(5)σ¯=3J2′=3C.


### 2.2. Multilinear Kinematic Hardening

The stress-strain curve of sandstone samples is obtained by laboratory experiments, as shown in Figure1. The samples S1-S2 are drilled from the same origin core sample with the size of 5 cm ×Φ2.5 cm, so it is with the sample S3-S4. The subscript of sample S5 means the cubic voxel with different subscript are extracted from the micro-CT image and used for RVE analysis. According to the determination of tangential modulus proposed by Fairhurst and Hudson [41], 50% of the uniaxial compressive strength of the rock is selected as the yield point, and the tangential modulus is defined as the slope of the curve between the yield point and the compressive strength point.Figure 1
Multilinear kinematic hardening.As shown in Figure1, points A1, A2, and A3 are the yield points and B1, B2, and B3 are the compressive limits. It should be stated that the mechanical parameters used in this paper are core scale test results. Considering the effect of the pore structure characteristics on the rock mechanical properties, the results are based on the model values/input values to reduce the scale effect of the input parameters on the results. The parameters used in the numerical simulation are shown in Table 1. These input mechanical parameters are calculated by the stress-strain curve (Figure 1) of the uniaxial compression test at the core scale. As stated by Liu and Tang [42], Young’s modulus is the slope of the elastic deformation stage (OA1, OA2, and OA3) and the compressive strength is the value of points B1, B2, and B3. The volume modulus and shear modulus are obtained according to equations (6) and (7), respectively:(6)G=E21+v,(7)Kv=E31−2v,where E is Young’s modulus (GPa), G is shear modulus (GPa), v is Poisson’s ratio, and KV is volume modulus (GPa).Table 1
Input parameters of numerical simulation.
Model
Young’s modulus (GPa)
Compressive strength (MPa)
Yield strength (MPa)
Poisson’s ratio
Bulk Modulus (GPa)
Shear modulus (GPa)
S5100-1∼S5400-2
27.83
112.11
61.97
0.296
22.737
10.737
S1-S2
13.864
105.94
52.93
0.296
11.327
53.488
S3-S4
7.803
61.4
44.5
0.296
6.375
3.01The four-linear model proposed by Lu [43] is adopted in the study, which describes the bilinear elasticity, linear softening, and residual ideal plasticity model of the rock. As shown in Figure 1, the red line OA1, OA2, and OA3 are linear elastic deformation stages. The rock changes from elastic deformation into plastic deformation after the yield points A1, A2, and A3, which behaves as nonlinear deformation. For the convenience of analysis and application, it is simplified into the linear hardening stage. After points B1, B2, and B3, plastic failure occurs in rocks. B1C1, B2C2, and B3C3 are simplified to the ideal plastic stage. At this stage, the compressive strength of the rock will not change with the variation of strain.


## 2.1. Plasticity Criterion

Considering the influence of the microstructure on rock mechanical strength, the von Mises yield criterion has been widely applied in the microscale deformation-fracture process of the rock, which are adopted for determining the plastic deformation of the rock skeleton [37–39].The criterion states that yielding begins when the elastic energy of distortion reaches a critical value [40]. That is, when the second invariant of the stress deviator (J2′) reaches a critical value, it represents plastic status, which can be expressed as follows:(1)fσij′=J2′=C,where J2′ is the second invariant of the stress tensor and C is a constant related to material property and independent of the stress state, which can be measured experimentally.The second invariant of the stress tensor is as follows:(2)J2′=16σx−σy2+σy−σx2+σz−σx2+6τxy2+τyz2+τzx2=C,where σx, σy, and σz are normal stresses in different directions (Pa) and τxy, τxz, and τyz are shear stresses (Pa).Equation (2) can be defined as effective stress as follows:(3)σ¯=12σ1−σ22+σ2−σ32+σ3−σ12=3C,where σ¯ is effective stress (Pa) and σ1 and σ3 are principal stresses (Pa).Thus, the Mises yield condition can be represented by more convenient effective stress. It can be expressed as equation (4) in the form of the principal stress and simplified as shown in the following equation:(4)J2′=16σ1−σ22+σ2−σ32+σ3−σ12=C,(5)σ¯=3J2′=3C.


## 2.2. Multilinear Kinematic Hardening

The stress-strain curve of sandstone samples is obtained by laboratory experiments, as shown in Figure1. The samples S1-S2 are drilled from the same origin core sample with the size of 5 cm ×Φ2.5 cm, so it is with the sample S3-S4. The subscript of sample S5 means the cubic voxel with different subscript are extracted from the micro-CT image and used for RVE analysis. According to the determination of tangential modulus proposed by Fairhurst and Hudson [41], 50% of the uniaxial compressive strength of the rock is selected as the yield point, and the tangential modulus is defined as the slope of the curve between the yield point and the compressive strength point.Figure 1
Multilinear kinematic hardening.As shown in Figure1, points A1, A2, and A3 are the yield points and B1, B2, and B3 are the compressive limits. It should be stated that the mechanical parameters used in this paper are core scale test results. Considering the effect of the pore structure characteristics on the rock mechanical properties, the results are based on the model values/input values to reduce the scale effect of the input parameters on the results. The parameters used in the numerical simulation are shown in Table 1. These input mechanical parameters are calculated by the stress-strain curve (Figure 1) of the uniaxial compression test at the core scale. As stated by Liu and Tang [42], Young’s modulus is the slope of the elastic deformation stage (OA1, OA2, and OA3) and the compressive strength is the value of points B1, B2, and B3. The volume modulus and shear modulus are obtained according to equations (6) and (7), respectively:(6)G=E21+v,(7)Kv=E31−2v,where E is Young’s modulus (GPa), G is shear modulus (GPa), v is Poisson’s ratio, and KV is volume modulus (GPa).Table 1
Input parameters of numerical simulation.
Model
Young’s modulus (GPa)
Compressive strength (MPa)
Yield strength (MPa)
Poisson’s ratio
Bulk Modulus (GPa)
Shear modulus (GPa)
S5100-1∼S5400-2
27.83
112.11
61.97
0.296
22.737
10.737
S1-S2
13.864
105.94
52.93
0.296
11.327
53.488
S3-S4
7.803
61.4
44.5
0.296
6.375
3.01The four-linear model proposed by Lu [43] is adopted in the study, which describes the bilinear elasticity, linear softening, and residual ideal plasticity model of the rock. As shown in Figure 1, the red line OA1, OA2, and OA3 are linear elastic deformation stages. The rock changes from elastic deformation into plastic deformation after the yield points A1, A2, and A3, which behaves as nonlinear deformation. For the convenience of analysis and application, it is simplified into the linear hardening stage. After points B1, B2, and B3, plastic failure occurs in rocks. B1C1, B2C2, and B3C3 are simplified to the ideal plastic stage. At this stage, the compressive strength of the rock will not change with the variation of strain.


## 3. Digital Core RVE Analysis

### 3.1. Digital Core Model

In this paper, micro-CT images are obtained by micro-CT scanning. After filtering, threshold segmentation, and converting into digital images, the finite element mesh model of the rock pore scale skeleton is reconstructed. 16 models with different porosity are obtained. S1 is taken as an example, as shown in Figure2.Figure 2
Rock skeleton finite element mesh model.


### 3.2. Characterized Parameters of Rock Microstructure

The porosity, TC, and tortuosity of the rock skeleton are used as the characterized parameters of the rock microstructure in this paper.


#### 3.2.1. Texture Coefficient

Texture coefficient is a dimensionless parameter, which is used to quantitatively describe the shape, orientation, interlocking degree, and relative properties of grain and matrix [18]. The texture factor takes into account the contour (or grain shape) of the grain, the degree of orientation, and interlocking and the area, length, width, etc. of the grain by image analysis on the thin section of the rock, which can be expressed as follows:(8)TC=AWN0N0+N1×1FF0+N1N0+N1×AR1×AF1,where TC is the texture coefficient; AW is grain area weighting; N0 and N1 represent the number of grain with aspect ratio less than 2 and greater than 2, respectively; FF0 is the arithmetic mean of shape factor of N0 grains, and the formula for each grain shape factor is 4πS/L2, S is the area of grains and L is the perimeter of the grains; AR1 is the arithmetic mean of aspect ratio of N1 grains; AF1 is angle factor orientation of N1 grains, its calculation formula is(9)AF1=15∑i=19XiNN−1/2,where N is the number of elongated grains, Xi is the number of angular differences in each class, and i is the weighting factor and class number.The pores are regarded as grains with the mechanical strength of 0. The images ofX, Y, and Z directions are generated by the VRML of AVIZO software. The pore-related parameters of each section are calculated using the ImageJ software [22], such as area, perimeter, length, width, and alignment direction, as shown in Figure 3. The relevant sizes of the two-dimensional section are shown in Table 2. Moreover, the average TC of all slices of each model are calculated, as shown in Table 3. The TC values of the same sample vary with the different directions, indicating the anisotropy of this microstructure parameter.Figure 3
Schematic diagram of solution of pore TC.Table 2
Size of all slices.
Model
Model size (mm)
Model
Model size (mm)
S1
1.3 × 1.3
S5200-2
0.42 × 0.42
S2
1 × 1
S5200-1
0.42 × 0.42
S3
0.75 × 0.75
S5250-1
0.53 × 0.53
S4
1 × 1
S5250-2
0.53 × 0.53
S5100-1
0.21 × 0.21
S5300-1
0.63 × 0.63
S5100-2
0.21 × 0.21
S5300-2
0.63 × 0.63
S5150-1
0.32 × 0.32
S5400-1
0.84 × 0.84
S5150-2
0.32 × 0.32
S5400-2
0.84 × 0.84Table 3
Average TC values of the model.
Model
Direction
X
Y
Z
S1
0.124
0.132
0.197
S2
0.133
0.163
0.164
S3
0.48
0.46
0.599
S4
0.482
0.47
0.396
S5100-1
0.151
0.191
0.187
S5100-2
0.146
0.153
0.142
S5150-1
0.069
0.06
0.161
S5150-2
0.058
0.045
0.176
S5200-1
0.132
0.095
0.122
S5200-2
0.127
0.123
0.163
S5250-1
0.098
0.089
0.102
S5250-2
0.103
0.105
0.119
S5300-1
0.093
0.086
0.100
S5300-2
0.094
0.091
0.109
S5400-1
0.097
0.083
0.095
S5400-2
0.089
0.089
0.105


#### 3.2.2. Tortuosity of Rock Skeleton

The tortuosity is defined as the ratio of the actual length of the percolation channel to the apparent length (macroscopic distance) through the percolation medium, as shown in Figure4. The tortuosity firstly proposed by Epstein [44] was used to describe the seepage channel. It can be described as(10)τ=LtL0,where Lt is the length of the curved line (mm) and L0 is the length of the straight line of the medium (mm).Figure 4
Schematic diagram of fluid passing through a porous medium.However, the curved lineLt is difficult to calculate directly. Thus, the TauFactor codes programmed by Cooper et al. [45] using MATLAB software are adopted to calculate the tortuosity factor. The codes simulate the reduction in diffusive transport caused by convolution in the geometry of heterogeneous media using the voxels directly as mesh elements, which can be expressed as follows:(11)Deff=Dετ,where ε is the volume fraction of the conductive phase; D is the intrinsic diffusivity of the conductive phase; Deff is the effective diffusivity through a porous volume; and τ is the tortuosity factor. The codes are based on the theory of pixel-incorporated array convolution and Dirichlet boundary of digital rock images.The tortuosity is usually used to describe the bending ratio of a fluid through a rock pore channel. Here, the tortuosity of the rock skeleton is adopted and calculated to quantitatively describe the microstructure of the rock matrix. Based on the acquired microscopic image data, the tortuosity in different directions is obtained, then the skeleton tortuosity of 16 kinds of microscopic images is calculated, as shown in Table4. The tortuosity values of the same sample vary from the different directions, indicating the anisotropy of microstructure parameter of the rock matrix.Table 4
Model’s tortuosity parameters.
Model
Direction
X
Y
Z
S1
1.22
1.2
1.28
S2
1.09
1.53
1.59
S3
2.2
2.36
2.9
S4
2.02
2.22
1.83
S5100-1
1.14
1.13
1.15
S5100-2
1.13
1.4
1.14
S5150-1
1.1
1.12
1.16
S5150-2
1.15
1.12
1.15
S5200-1
1.13
1.13
1.15
S5200-2
1.12
1.15
1.16
S5250-1
1.12
1.12
1.14
S5250-2
1.14
1.14
1.15
S5300-1
1.12
1.12
1.14
S5300-2
1.14
1.14
1.15
S5400-1
1.12
1.12
1.14
S5400-2
1.14
1.14
1.15


### 3.3. RVE Analysis of the Digital Rock Models

The RVE is the smallest subvolume of the heterogeneous material to be used to determine the corresponding effective properties for the macroscopic model. The RVE should be large enough to contain sufficient information about the microstructure in order to be representative; however, it should be much smaller than the macroscopic body [46]. RVE is considered the basic concept for estimating average or overall properties of heterogeneous solids. It can not only represent the effective response of the microstructure but also minimize the computational effort. As shown in Figure 5(a), a cubic of 100–400 pixels is extracted from the center of the original three-dimensional cylindrical rock sample image, which is gradually expanded to the surrounding area. The microscale deformation analysis of the digital core requires an accurate physical model, which can represent the microstructure of the core and meet the needs of computer simulation.Table 5
Detailed parameters of RVEs of digital cores.
Model
Mesh model volume (m3)
Model mesh number
Porosity (%)
Side length (mm)
REV edge length (pixel)
S5300-1
2.29 × 10−10
1408315
8.36
0.63
300
S5300-2
2.27 × 10−10
1744435
9.07
0.63
S5400-1
5.44 × 10−10
3572342
8.30
0.84
400
S5400-2
5.39 × 10−10
3989611
8.39
0.84
S1
1.81 × 10−9
2787997
17.6
1.3
300
S2
7.97 × 10−10
797051
20.3
1
S3
2.6 × 10−10
2072375
38.6
0.75
400
S4
5.96 × 10−10
596233
40.4
1Figure 5
RVE analysis of mechanical properties and microstructure: (a) the model selection for RVE analysis and (b)∼(h) the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, ratio of shear modulus, ratio of Poisson’s ratio, tortuosity, and TC against the model size for the RVE analysis.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)The displacement loading is used to simulate the pore scale uniaxial compression in theX, Y, and Z directions. The fixed constraint is set at the bottom of the model, and the top is the displacement loading condition, which only allows the core to be displaced in the vertical direction and the horizontal direction is unconstrained. The mechanical parameters (e.g., Young’s modulus, Poisson’s ratio, compressive strength, yield strength, and shear modulus) are obtained by using ANSYS software. Thus, the RVE analysis of the rock skeleton on the microstructure parameters (e.g., tortuosity and TC) and mechanical parameters against the model size are conducted; the RVE analysis results of sample S5 are shown in Figure 5.In order to reduce the size effect of the model, when the mechanical parameters are used for RVE analysis, the ratio of the output parameters to the input parameters is taken. As shown in Figure5, all the parameters fluctuate greatly when the length of the model is between 100 and 250 pixels and approach steady when the length is beyond 300 pixels. Thus, the RVE of sample S5 is 300. Through the same methods, the optimal RVE of different rock samples is determined. Based on this, the finite element mesh model of each rock sample is established. The RVE parameters of the models S1, S2, S3, S4, and S5 are shown in Table 5.


## 3.1. Digital Core Model

In this paper, micro-CT images are obtained by micro-CT scanning. After filtering, threshold segmentation, and converting into digital images, the finite element mesh model of the rock pore scale skeleton is reconstructed. 16 models with different porosity are obtained. S1 is taken as an example, as shown in Figure2.Figure 2
Rock skeleton finite element mesh model.


## 3.2. Characterized Parameters of Rock Microstructure

The porosity, TC, and tortuosity of the rock skeleton are used as the characterized parameters of the rock microstructure in this paper.


### 3.2.1. Texture Coefficient

Texture coefficient is a dimensionless parameter, which is used to quantitatively describe the shape, orientation, interlocking degree, and relative properties of grain and matrix [18]. The texture factor takes into account the contour (or grain shape) of the grain, the degree of orientation, and interlocking and the area, length, width, etc. of the grain by image analysis on the thin section of the rock, which can be expressed as follows:(8)TC=AWN0N0+N1×1FF0+N1N0+N1×AR1×AF1,where TC is the texture coefficient; AW is grain area weighting; N0 and N1 represent the number of grain with aspect ratio less than 2 and greater than 2, respectively; FF0 is the arithmetic mean of shape factor of N0 grains, and the formula for each grain shape factor is 4πS/L2, S is the area of grains and L is the perimeter of the grains; AR1 is the arithmetic mean of aspect ratio of N1 grains; AF1 is angle factor orientation of N1 grains, its calculation formula is(9)AF1=15∑i=19XiNN−1/2,where N is the number of elongated grains, Xi is the number of angular differences in each class, and i is the weighting factor and class number.The pores are regarded as grains with the mechanical strength of 0. The images ofX, Y, and Z directions are generated by the VRML of AVIZO software. The pore-related parameters of each section are calculated using the ImageJ software [22], such as area, perimeter, length, width, and alignment direction, as shown in Figure 3. The relevant sizes of the two-dimensional section are shown in Table 2. Moreover, the average TC of all slices of each model are calculated, as shown in Table 3. The TC values of the same sample vary with the different directions, indicating the anisotropy of this microstructure parameter.Figure 3
Schematic diagram of solution of pore TC.Table 2
Size of all slices.
Model
Model size (mm)
Model
Model size (mm)
S1
1.3 × 1.3
S5200-2
0.42 × 0.42
S2
1 × 1
S5200-1
0.42 × 0.42
S3
0.75 × 0.75
S5250-1
0.53 × 0.53
S4
1 × 1
S5250-2
0.53 × 0.53
S5100-1
0.21 × 0.21
S5300-1
0.63 × 0.63
S5100-2
0.21 × 0.21
S5300-2
0.63 × 0.63
S5150-1
0.32 × 0.32
S5400-1
0.84 × 0.84
S5150-2
0.32 × 0.32
S5400-2
0.84 × 0.84Table 3
Average TC values of the model.
Model
Direction
X
Y
Z
S1
0.124
0.132
0.197
S2
0.133
0.163
0.164
S3
0.48
0.46
0.599
S4
0.482
0.47
0.396
S5100-1
0.151
0.191
0.187
S5100-2
0.146
0.153
0.142
S5150-1
0.069
0.06
0.161
S5150-2
0.058
0.045
0.176
S5200-1
0.132
0.095
0.122
S5200-2
0.127
0.123
0.163
S5250-1
0.098
0.089
0.102
S5250-2
0.103
0.105
0.119
S5300-1
0.093
0.086
0.100
S5300-2
0.094
0.091
0.109
S5400-1
0.097
0.083
0.095
S5400-2
0.089
0.089
0.105


### 3.2.2. Tortuosity of Rock Skeleton

The tortuosity is defined as the ratio of the actual length of the percolation channel to the apparent length (macroscopic distance) through the percolation medium, as shown in Figure4. The tortuosity firstly proposed by Epstein [44] was used to describe the seepage channel. It can be described as(10)τ=LtL0,where Lt is the length of the curved line (mm) and L0 is the length of the straight line of the medium (mm).Figure 4
Schematic diagram of fluid passing through a porous medium.However, the curved lineLt is difficult to calculate directly. Thus, the TauFactor codes programmed by Cooper et al. [45] using MATLAB software are adopted to calculate the tortuosity factor. The codes simulate the reduction in diffusive transport caused by convolution in the geometry of heterogeneous media using the voxels directly as mesh elements, which can be expressed as follows:(11)Deff=Dετ,where ε is the volume fraction of the conductive phase; D is the intrinsic diffusivity of the conductive phase; Deff is the effective diffusivity through a porous volume; and τ is the tortuosity factor. The codes are based on the theory of pixel-incorporated array convolution and Dirichlet boundary of digital rock images.The tortuosity is usually used to describe the bending ratio of a fluid through a rock pore channel. Here, the tortuosity of the rock skeleton is adopted and calculated to quantitatively describe the microstructure of the rock matrix. Based on the acquired microscopic image data, the tortuosity in different directions is obtained, then the skeleton tortuosity of 16 kinds of microscopic images is calculated, as shown in Table4. The tortuosity values of the same sample vary from the different directions, indicating the anisotropy of microstructure parameter of the rock matrix.Table 4
Model’s tortuosity parameters.
Model
Direction
X
Y
Z
S1
1.22
1.2
1.28
S2
1.09
1.53
1.59
S3
2.2
2.36
2.9
S4
2.02
2.22
1.83
S5100-1
1.14
1.13
1.15
S5100-2
1.13
1.4
1.14
S5150-1
1.1
1.12
1.16
S5150-2
1.15
1.12
1.15
S5200-1
1.13
1.13
1.15
S5200-2
1.12
1.15
1.16
S5250-1
1.12
1.12
1.14
S5250-2
1.14
1.14
1.15
S5300-1
1.12
1.12
1.14
S5300-2
1.14
1.14
1.15
S5400-1
1.12
1.12
1.14
S5400-2
1.14
1.14
1.15


## 3.2.1. Texture Coefficient

Texture coefficient is a dimensionless parameter, which is used to quantitatively describe the shape, orientation, interlocking degree, and relative properties of grain and matrix [18]. The texture factor takes into account the contour (or grain shape) of the grain, the degree of orientation, and interlocking and the area, length, width, etc. of the grain by image analysis on the thin section of the rock, which can be expressed as follows:(8)TC=AWN0N0+N1×1FF0+N1N0+N1×AR1×AF1,where TC is the texture coefficient; AW is grain area weighting; N0 and N1 represent the number of grain with aspect ratio less than 2 and greater than 2, respectively; FF0 is the arithmetic mean of shape factor of N0 grains, and the formula for each grain shape factor is 4πS/L2, S is the area of grains and L is the perimeter of the grains; AR1 is the arithmetic mean of aspect ratio of N1 grains; AF1 is angle factor orientation of N1 grains, its calculation formula is(9)AF1=15∑i=19XiNN−1/2,where N is the number of elongated grains, Xi is the number of angular differences in each class, and i is the weighting factor and class number.The pores are regarded as grains with the mechanical strength of 0. The images ofX, Y, and Z directions are generated by the VRML of AVIZO software. The pore-related parameters of each section are calculated using the ImageJ software [22], such as area, perimeter, length, width, and alignment direction, as shown in Figure 3. The relevant sizes of the two-dimensional section are shown in Table 2. Moreover, the average TC of all slices of each model are calculated, as shown in Table 3. The TC values of the same sample vary with the different directions, indicating the anisotropy of this microstructure parameter.Figure 3
Schematic diagram of solution of pore TC.Table 2
Size of all slices.
Model
Model size (mm)
Model
Model size (mm)
S1
1.3 × 1.3
S5200-2
0.42 × 0.42
S2
1 × 1
S5200-1
0.42 × 0.42
S3
0.75 × 0.75
S5250-1
0.53 × 0.53
S4
1 × 1
S5250-2
0.53 × 0.53
S5100-1
0.21 × 0.21
S5300-1
0.63 × 0.63
S5100-2
0.21 × 0.21
S5300-2
0.63 × 0.63
S5150-1
0.32 × 0.32
S5400-1
0.84 × 0.84
S5150-2
0.32 × 0.32
S5400-2
0.84 × 0.84Table 3
Average TC values of the model.
Model
Direction
X
Y
Z
S1
0.124
0.132
0.197
S2
0.133
0.163
0.164
S3
0.48
0.46
0.599
S4
0.482
0.47
0.396
S5100-1
0.151
0.191
0.187
S5100-2
0.146
0.153
0.142
S5150-1
0.069
0.06
0.161
S5150-2
0.058
0.045
0.176
S5200-1
0.132
0.095
0.122
S5200-2
0.127
0.123
0.163
S5250-1
0.098
0.089
0.102
S5250-2
0.103
0.105
0.119
S5300-1
0.093
0.086
0.100
S5300-2
0.094
0.091
0.109
S5400-1
0.097
0.083
0.095
S5400-2
0.089
0.089
0.105


## 3.2.2. Tortuosity of Rock Skeleton

The tortuosity is defined as the ratio of the actual length of the percolation channel to the apparent length (macroscopic distance) through the percolation medium, as shown in Figure4. The tortuosity firstly proposed by Epstein [44] was used to describe the seepage channel. It can be described as(10)τ=LtL0,where Lt is the length of the curved line (mm) and L0 is the length of the straight line of the medium (mm).Figure 4
Schematic diagram of fluid passing through a porous medium.However, the curved lineLt is difficult to calculate directly. Thus, the TauFactor codes programmed by Cooper et al. [45] using MATLAB software are adopted to calculate the tortuosity factor. The codes simulate the reduction in diffusive transport caused by convolution in the geometry of heterogeneous media using the voxels directly as mesh elements, which can be expressed as follows:(11)Deff=Dετ,where ε is the volume fraction of the conductive phase; D is the intrinsic diffusivity of the conductive phase; Deff is the effective diffusivity through a porous volume; and τ is the tortuosity factor. The codes are based on the theory of pixel-incorporated array convolution and Dirichlet boundary of digital rock images.The tortuosity is usually used to describe the bending ratio of a fluid through a rock pore channel. Here, the tortuosity of the rock skeleton is adopted and calculated to quantitatively describe the microstructure of the rock matrix. Based on the acquired microscopic image data, the tortuosity in different directions is obtained, then the skeleton tortuosity of 16 kinds of microscopic images is calculated, as shown in Table4. The tortuosity values of the same sample vary from the different directions, indicating the anisotropy of microstructure parameter of the rock matrix.Table 4
Model’s tortuosity parameters.
Model
Direction
X
Y
Z
S1
1.22
1.2
1.28
S2
1.09
1.53
1.59
S3
2.2
2.36
2.9
S4
2.02
2.22
1.83
S5100-1
1.14
1.13
1.15
S5100-2
1.13
1.4
1.14
S5150-1
1.1
1.12
1.16
S5150-2
1.15
1.12
1.15
S5200-1
1.13
1.13
1.15
S5200-2
1.12
1.15
1.16
S5250-1
1.12
1.12
1.14
S5250-2
1.14
1.14
1.15
S5300-1
1.12
1.12
1.14
S5300-2
1.14
1.14
1.15
S5400-1
1.12
1.12
1.14
S5400-2
1.14
1.14
1.15


## 3.3. RVE Analysis of the Digital Rock Models

The RVE is the smallest subvolume of the heterogeneous material to be used to determine the corresponding effective properties for the macroscopic model. The RVE should be large enough to contain sufficient information about the microstructure in order to be representative; however, it should be much smaller than the macroscopic body [46]. RVE is considered the basic concept for estimating average or overall properties of heterogeneous solids. It can not only represent the effective response of the microstructure but also minimize the computational effort. As shown in Figure 5(a), a cubic of 100–400 pixels is extracted from the center of the original three-dimensional cylindrical rock sample image, which is gradually expanded to the surrounding area. The microscale deformation analysis of the digital core requires an accurate physical model, which can represent the microstructure of the core and meet the needs of computer simulation.Table 5
Detailed parameters of RVEs of digital cores.
Model
Mesh model volume (m3)
Model mesh number
Porosity (%)
Side length (mm)
REV edge length (pixel)
S5300-1
2.29 × 10−10
1408315
8.36
0.63
300
S5300-2
2.27 × 10−10
1744435
9.07
0.63
S5400-1
5.44 × 10−10
3572342
8.30
0.84
400
S5400-2
5.39 × 10−10
3989611
8.39
0.84
S1
1.81 × 10−9
2787997
17.6
1.3
300
S2
7.97 × 10−10
797051
20.3
1
S3
2.6 × 10−10
2072375
38.6
0.75
400
S4
5.96 × 10−10
596233
40.4
1Figure 5
RVE analysis of mechanical properties and microstructure: (a) the model selection for RVE analysis and (b)∼(h) the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, ratio of shear modulus, ratio of Poisson’s ratio, tortuosity, and TC against the model size for the RVE analysis.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)The displacement loading is used to simulate the pore scale uniaxial compression in theX, Y, and Z directions. The fixed constraint is set at the bottom of the model, and the top is the displacement loading condition, which only allows the core to be displaced in the vertical direction and the horizontal direction is unconstrained. The mechanical parameters (e.g., Young’s modulus, Poisson’s ratio, compressive strength, yield strength, and shear modulus) are obtained by using ANSYS software. Thus, the RVE analysis of the rock skeleton on the microstructure parameters (e.g., tortuosity and TC) and mechanical parameters against the model size are conducted; the RVE analysis results of sample S5 are shown in Figure 5.In order to reduce the size effect of the model, when the mechanical parameters are used for RVE analysis, the ratio of the output parameters to the input parameters is taken. As shown in Figure5, all the parameters fluctuate greatly when the length of the model is between 100 and 250 pixels and approach steady when the length is beyond 300 pixels. Thus, the RVE of sample S5 is 300. Through the same methods, the optimal RVE of different rock samples is determined. Based on this, the finite element mesh model of each rock sample is established. The RVE parameters of the models S1, S2, S3, S4, and S5 are shown in Table 5.


## 4. Pore Scale Simulation on Uniaxial Compression

The optimal REV model is determined in Section3.3. The uniaxial compression numerical simulation in this section is conducted on the RVE model. The averaged von Mises stress distribution, total displacement field, plastic strain field, and stress-strain curve are obtained. The following analysis is investigated by taking model S2 as an example.Figure6 shows averaged von Mises stress distribution of model S2 at the termination of displacement loading. The red and blue zones indicate the maximum and minimum stresses, respectively. It can be seen that stress concentration occurs on each loading surface, and the shape of the rock skeleton is extremely irregular in the stress concentration region. In the case of applying the same displacement, the stress of the rock sample in the X direction is relatively uniform, and the stress distribution in the Y direction is mainly concentrated in the right half of the rock sample. The compressive performance in the X direction is significantly larger than that of the Y direction. The existence of pores leads to nonuniform distribution of stress, which is shown in detail in Figure 7.Figure 6
Averaged von Mises stress distribution of model S2 (unit: Pa) represents the loading applied (a) inX direction and (b) in Y direction.
(a)
(b)Figure 7
Cross-sections averaged von Mises stress distribution of model S2 (unit: Pa): (a) spatial location of slices at (b)y = 1 mm, (c) y = 0.8 mm, (d) y = 0.6 mm, (e) y = 0.4 mm, and (f) y = 0.2 mm.
(a)
(b)
(c)
(d)
(e)
(f)Figure8 shows the distribution of the rock total displacement. The deformation of the model right part is serious, as shown in Figure 8(a), and it indicates that the compressive capacity of this part is relatively weak, and Figure 8(b) is the opposite. During the deformation process of the rock, the displacement field exhibits the translational characteristics, and arch-shaped deformation regions are formed between the deformation zones. The slice displacement details of different spatial positions are shown in Figure 9.Figure 8
Total displacement distribution of model S2 (unit: m) represents the loading applied (a) inX direction and (b) in Y direction.
(a)
(b)Figure 9
Cross-sections total displacement distribution of model S2 (unit: m): (a) spatial location of slices at (b)y = 1 mm, (c) y = 0.8 mm, (d) y = 0.6 mm, (e) y = 0.4 mm, and (f) y = 0.2 mm.
(a)
(b)
(c)
(d)
(e)
(f)Figure10 is the plastic strain distribution of model S2, which directly reflects the deformation of the rock after yield. When the bearing capacity reaches the limit, plastic failure will occur in the rock. It can be seen from Figure 11 that at the end of displacement loading, plastic deformation passes through the pores, forming various shape shear zones. An obvious X-shaped plastic conjugate shear zone (red circle) is formed, as shown in Figure 10(a). Combining with the averaged von Mises stress distribution, the plastic failure priorly occurs where the stress is relatively high.Figure 10
Plastic stain distribution of model S2 represents the loading applied (a) inX direction and (b) in Y direction.
(a)
(b)Figure 11
Cross-sections plastic strain distribution of model S2: (a) spatial location of slices at (b)y = 1 mm, (c) y = 0.8 mm, (d) y = 0.6 mm, (e) y = 0.4 mm, and (f) y = 0.2 mm.
(a)
(b)
(c)
(d)
(e)
(f)The stress-strain curves of different models are shown in Figure12. In the process of numerical simulation, models S1-S2, S3-S4, and S5300-1∼S5400-2 used three groups of input parameters. It can be seen that the stress-strain curve of the models varies with the loading directions, indicating the anisotropy and heterogeneity of the rock. The anisotropy of the models S1 to S4 is significantly higher than model S5300-1∼S5400-2. In model S1∼S4, when the strain is about 0.5%, the rock reaches the elastic limit and enters the plastic yield state. The yield strength of model S1 and S2 are about 40 MPa.Figure 12
Stress-strain curves of digital rock samples. (a)∼(h) are the stress-strain curves of model S1, S2, S3, S4, S5300-1, S5300-2, S5400-1, and S5400-2. (a) φ = 17.6%. (b) φ = 20.3%. (c) φ = 38.6%. (d) φ = 40.4%. (e) φ = 8.36%. (f) φ = 9.07%. (g) φ = 8.3%. (h) φ = 8.39%.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)In model S5300-1∼S5400-2, there is little difference in porosity, the compressive strength of rock is about 100 MPa, the yield strength is about 60 MPa, and mechanical parameters in each direction are not much different, indicating that the anisotropy of the rock is not obvious. When the strain is 1%, the rock is destroyed and the compressive limit is reached. When the model reaches the compressive limit, the strain increases continuously, and the stress is almost constant. The reason is the plastic skeleton of the rock is assumed to be perfectly plastic when reaching the compressive strength.


## 5. Analysis of Related Parameters of Digital Core

### 5.1. Correlation Analysis between Porosity and Mechanical Parameters

Porosity refers to the ratio of pore volume to the entire volume of the rock sample, which is an important parameter used to characterize the rock microstructure. The mechanical parameters of different rock samples are obtained by numerical simulation and the relationship between porosity and mechanical parameters, as shown in Figure13.Figure 13
Correlation analysis between porosity and rock mechanics parameters. (a)∼(d) are the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, and ratio of shear modulus against the model porosity for correlation analysis.
(a)
(b)
(c)
(d)As shown in Figure13(a), the ratio of Young’s modulus of digital rock samples decreases linearly with the increase of porosity (R2 = 0.87). Moreover, it is found that the effective Young’s modulus of the rock has a good exponential relationship with porosity (R2 = 0.975) by Cheng et al. [35]. In addition, Pabst et al. [47] also found Young’s modulus has a good correlation with porosity from the perspective of the exponential model (R2 = 0.99). The fitting coefficient of the simulation is lower than that of the literature. The reason is that tight sandstone and loose sandstone from different sites are adopted in this study, while the rock samples are drilled from the same sites in the literature. In general, compressive strength decreases with the increase of the porosity. The ratio of compressive strength is negatively correlated with porosity (Figure 13(b) R2 = 0.87), the result is consistent with Zhou and Xiao [48]. The ratio of yield strength is shown in Figure 13(c) and its corresponding coefficient R2 is 0.90. The ratio of shear modulus is deeply related to porosity, and the fitting curve is shown in Figure 13(d) (R2 = 0.84). With the increase of porosity, all mechanical parameters of the rock decrease linearly. The main reason is that the cementation becomes worse and the bearing capacity becomes weaker with the rising of porosity. However, for the same model with the same porosity, there are great differences in Young’s modulus, compressive strength, yield strength, and shear modulus in different directions. It indicates that the porosity is possible to reproduce the anisotropy and heterogeneity of the rock microstructure.


### 5.2. Correlation Analysis between TC and Mechanical Parameters

TC is an important parameter to quantitatively describe the internal microstructure of the rock and evaluate the mechanical properties of the rock. As shown in Figure14, with the increase of the TC value, the mechanical parameters of digital rock samples show logarithmic decline. In Figure 14(a), Young’s modulus ratio of rock samples gradually decreases with the increase of TC, and its correlation coefficient R2 is 0.83. Meanwhile, the TC value has a high correlation with Young’s modulus (Figure 14(b) R2 = 0.87), yield strength (Figure 14(c) R2 = 0.89), and shear modulus (Figure 14(d) R2 = 0.85). These fitting results are also consistent with the results of other scholars. Howarth and Rowlands [18] found the TC was highly correlated with the compressive strength (R2 = 0.92). The reported R2 correlations of TC between compressive strength was 0.69 and Young’s modulus was 0.69 by Erosy and Waller [19]. Alber and Kahraman [21] reported that the compressive strength of the prediction model had a very high correlation with TC (R2 = 0.90). The difference between the results of this paper and the literature is that the pores are regarded as grains with a mechanical strength of 0, while the TC of the grains was investigated in the literature. According to the definition of TC, it is determined by the number, shape, and randomness of the orientation of the pores. Přikryl [49] found mechanical properties of rocks which have been mostly affected by variable number, shape, and orientation of grains in some samples. TC increases with the rising of the proportion and aspect ratio of pores [18]. Compared to the porosity, TC varies with directions for the same model, which is able to represent the anisotropy and heterogeneity of the rock microstructure.Figure 14
Correlation analysis between TC and digital core mechanical parameters. (a)∼(d) are the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, and ratio of shear modulus against the model TC for correlation analysis.
(a)
(b)
(c)
(d)


### 5.3. Correlation Analysis between Tortuosity and Mechanical Parameters

The tortuosity has been regarded as an important parameter of the microstructure for porous media, which is related to the porosity and grain radius as well as permeability. Salem and Chilingarian [50] found the variation of tortuosity depends on the internal structure and the composition, such as the orientation of grains, pore-water salinity, degree of compaction, and flow direction. Wang et al. [51] proposed the inverse correlations between tortuosity and permeability. Attia [52] found the tortuosity factor varies with many parameters, e.g., the geometry of porous media, the number of fine grains, and cementation factor. These researchers mainly studied the tortuosity of fluid flowing through pore channels and analyzed the influencing factors of tortuosity. However, the quantitative characterization of rock mechanical properties by tortuosity has not been found in existing studies. In this paper, it is proposed for the first time to analyze the mechanical properties of the rock by using tortuosity of the rock skeleton. The correlations between the tortuosity of the digital rock sample skeleton and the mechanical parameters of the digital cores are shown in Figure 15. It is observed from Figure 15(a) that with the increase of the tortuosity of the rock skeleton, the ratio of Young’s modulus of digital rock samples shows logarithmic decrement (R2 = 0.95). Moreover, the ratio of compressive strength, yield strength, and shear modulus are in good correlation with the tortuosity of the rock skeleton, and the fitting coefficients are 0.94, 0.92, and 0.94, respectively. The corresponding fitting curves are depicted in Figures 15(b)–15(d). It is assumed that the rock is composed of the skeleton and pore. When calculating the tortuosity of the rock skeleton, the skeleton of the rock is taken as the pore channel of fluid. The greater the tortuosity of the rock skeleton is, the longer the distance for the fluid pass through is. It indicates that the more pores there are, the worse the mechanical properties of the rock will be.Figure 15
Correlation analysis between the tortuosity and the mechanical parameters of the digital core. (a)∼(d) are the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, and ratio of shear modulus against the model tortuosity for correlation analysis.
(a)
(b)
(c)
(d)The porosity and TC are well-known parameters describing the features of the rock interior. The tortuosity of the rock skeleton is firstly used to characterize the microstructure of the rock. The correlations between the three parameters and the mechanical properties of rock are presented in Table6. It can be seen from the table that the porosity, TC, and tortuosity all affect the mechanical parameters of the rock, while the correlation between the tortuosity and the rock mechanical parameters is higher than that of the porosity and TC. Thus, the tortuosity of the rock skeleton is believed to be more feasible to reveal the impacts of the microstructure of the rock on its mechanical properties.Table 6
Correlation between porosity, TC, tortuosity and mechanical parameters.
Microstructural parameters
Mechanical parameters
Fitting formula
Correlation coefficientR2
Porosity (φ)
Ratio of Young's modulus (E0/E)
E
0
/
E
=
−
0.02
φ
+
1.00
0.82
Ratio of compressive strength (σC0/σC)
σ
C
0
/
σ
C
=
−
0.02
φ
+
1.00
0.87
Ratio of yield strength (σS0/σS)
σ
S
0
/
σ
S
=
−
0.02
φ
+
0.94
0.90
Ratio of shear modulus (G0/G)
G
0
/
G
=
−
0.02
φ
+
1.04
0.84
Texture coefficient (TC)
Ratio of Young’s modulus (E0/E)
E
0
/
E
=
−
0.31
ln
TC
+
0.13
0.83
Ratio of compressive strength(σC0/σC)
σ
C
0
/
σ
C
=
−
0.33
ln
TC
+
0.07
0.87
Ratio of yield strength (σS0/σS)
σ
S
0
/
σ
S
=
−
0.33
ln
TC
+
0.006
0.89
Ratio of shear modulus (G0/G)
G
0
/
G
=
−
0.32
ln
TC
+
0.14
0.85
Tortuosity (τ)
Ratio of Young’s modulus (E0/E)
E
0
/
E
=
−
0.24
ln
τ
−
0.98
+
0.40
0.95
Ratio of compressive strength(σC0/σC)
σ
C
0
/
σ
C
=
−
0.26
ln
τ
−
0.97
+
0.37
0.94
Ratio of yield strength (σS0/σS)
σ
S
0
/
σ
S
=
−
0.24
ln
τ
−
0.98
+
0.31
0.92
Ratio of shear modulus (G0/G)
G
0
/
G
=
−
0.24
ln
τ
−
0.99
+
0.43
0.94


## 5.1. Correlation Analysis between Porosity and Mechanical Parameters

Porosity refers to the ratio of pore volume to the entire volume of the rock sample, which is an important parameter used to characterize the rock microstructure. The mechanical parameters of different rock samples are obtained by numerical simulation and the relationship between porosity and mechanical parameters, as shown in Figure13.Figure 13
Correlation analysis between porosity and rock mechanics parameters. (a)∼(d) are the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, and ratio of shear modulus against the model porosity for correlation analysis.
(a)
(b)
(c)
(d)As shown in Figure13(a), the ratio of Young’s modulus of digital rock samples decreases linearly with the increase of porosity (R2 = 0.87). Moreover, it is found that the effective Young’s modulus of the rock has a good exponential relationship with porosity (R2 = 0.975) by Cheng et al. [35]. In addition, Pabst et al. [47] also found Young’s modulus has a good correlation with porosity from the perspective of the exponential model (R2 = 0.99). The fitting coefficient of the simulation is lower than that of the literature. The reason is that tight sandstone and loose sandstone from different sites are adopted in this study, while the rock samples are drilled from the same sites in the literature. In general, compressive strength decreases with the increase of the porosity. The ratio of compressive strength is negatively correlated with porosity (Figure 13(b) R2 = 0.87), the result is consistent with Zhou and Xiao [48]. The ratio of yield strength is shown in Figure 13(c) and its corresponding coefficient R2 is 0.90. The ratio of shear modulus is deeply related to porosity, and the fitting curve is shown in Figure 13(d) (R2 = 0.84). With the increase of porosity, all mechanical parameters of the rock decrease linearly. The main reason is that the cementation becomes worse and the bearing capacity becomes weaker with the rising of porosity. However, for the same model with the same porosity, there are great differences in Young’s modulus, compressive strength, yield strength, and shear modulus in different directions. It indicates that the porosity is possible to reproduce the anisotropy and heterogeneity of the rock microstructure.


## 5.2. Correlation Analysis between TC and Mechanical Parameters

TC is an important parameter to quantitatively describe the internal microstructure of the rock and evaluate the mechanical properties of the rock. As shown in Figure14, with the increase of the TC value, the mechanical parameters of digital rock samples show logarithmic decline. In Figure 14(a), Young’s modulus ratio of rock samples gradually decreases with the increase of TC, and its correlation coefficient R2 is 0.83. Meanwhile, the TC value has a high correlation with Young’s modulus (Figure 14(b) R2 = 0.87), yield strength (Figure 14(c) R2 = 0.89), and shear modulus (Figure 14(d) R2 = 0.85). These fitting results are also consistent with the results of other scholars. Howarth and Rowlands [18] found the TC was highly correlated with the compressive strength (R2 = 0.92). The reported R2 correlations of TC between compressive strength was 0.69 and Young’s modulus was 0.69 by Erosy and Waller [19]. Alber and Kahraman [21] reported that the compressive strength of the prediction model had a very high correlation with TC (R2 = 0.90). The difference between the results of this paper and the literature is that the pores are regarded as grains with a mechanical strength of 0, while the TC of the grains was investigated in the literature. According to the definition of TC, it is determined by the number, shape, and randomness of the orientation of the pores. Přikryl [49] found mechanical properties of rocks which have been mostly affected by variable number, shape, and orientation of grains in some samples. TC increases with the rising of the proportion and aspect ratio of pores [18]. Compared to the porosity, TC varies with directions for the same model, which is able to represent the anisotropy and heterogeneity of the rock microstructure.Figure 14
Correlation analysis between TC and digital core mechanical parameters. (a)∼(d) are the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, and ratio of shear modulus against the model TC for correlation analysis.
(a)
(b)
(c)
(d)


## 5.3. Correlation Analysis between Tortuosity and Mechanical Parameters

The tortuosity has been regarded as an important parameter of the microstructure for porous media, which is related to the porosity and grain radius as well as permeability. Salem and Chilingarian [50] found the variation of tortuosity depends on the internal structure and the composition, such as the orientation of grains, pore-water salinity, degree of compaction, and flow direction. Wang et al. [51] proposed the inverse correlations between tortuosity and permeability. Attia [52] found the tortuosity factor varies with many parameters, e.g., the geometry of porous media, the number of fine grains, and cementation factor. These researchers mainly studied the tortuosity of fluid flowing through pore channels and analyzed the influencing factors of tortuosity. However, the quantitative characterization of rock mechanical properties by tortuosity has not been found in existing studies. In this paper, it is proposed for the first time to analyze the mechanical properties of the rock by using tortuosity of the rock skeleton. The correlations between the tortuosity of the digital rock sample skeleton and the mechanical parameters of the digital cores are shown in Figure 15. It is observed from Figure 15(a) that with the increase of the tortuosity of the rock skeleton, the ratio of Young’s modulus of digital rock samples shows logarithmic decrement (R2 = 0.95). Moreover, the ratio of compressive strength, yield strength, and shear modulus are in good correlation with the tortuosity of the rock skeleton, and the fitting coefficients are 0.94, 0.92, and 0.94, respectively. The corresponding fitting curves are depicted in Figures 15(b)–15(d). It is assumed that the rock is composed of the skeleton and pore. When calculating the tortuosity of the rock skeleton, the skeleton of the rock is taken as the pore channel of fluid. The greater the tortuosity of the rock skeleton is, the longer the distance for the fluid pass through is. It indicates that the more pores there are, the worse the mechanical properties of the rock will be.Figure 15
Correlation analysis between the tortuosity and the mechanical parameters of the digital core. (a)∼(d) are the ratio of Young’s modulus, ratio of compressive strength, ratio of yield strength, and ratio of shear modulus against the model tortuosity for correlation analysis.
(a)
(b)
(c)
(d)The porosity and TC are well-known parameters describing the features of the rock interior. The tortuosity of the rock skeleton is firstly used to characterize the microstructure of the rock. The correlations between the three parameters and the mechanical properties of rock are presented in Table6. It can be seen from the table that the porosity, TC, and tortuosity all affect the mechanical parameters of the rock, while the correlation between the tortuosity and the rock mechanical parameters is higher than that of the porosity and TC. Thus, the tortuosity of the rock skeleton is believed to be more feasible to reveal the impacts of the microstructure of the rock on its mechanical properties.Table 6
Correlation between porosity, TC, tortuosity and mechanical parameters.
Microstructural parameters
Mechanical parameters
Fitting formula
Correlation coefficientR2
Porosity (φ)
Ratio of Young's modulus (E0/E)
E
0
/
E
=
−
0.02
φ
+
1.00
0.82
Ratio of compressive strength (σC0/σC)
σ
C
0
/
σ
C
=
−
0.02
φ
+
1.00
0.87
Ratio of yield strength (σS0/σS)
σ
S
0
/
σ
S
=
−
0.02
φ
+
0.94
0.90
Ratio of shear modulus (G0/G)
G
0
/
G
=
−
0.02
φ
+
1.04
0.84
Texture coefficient (TC)
Ratio of Young’s modulus (E0/E)
E
0
/
E
=
−
0.31
ln
TC
+
0.13
0.83
Ratio of compressive strength(σC0/σC)
σ
C
0
/
σ
C
=
−
0.33
ln
TC
+
0.07
0.87
Ratio of yield strength (σS0/σS)
σ
S
0
/
σ
S
=
−
0.33
ln
TC
+
0.006
0.89
Ratio of shear modulus (G0/G)
G
0
/
G
=
−
0.32
ln
TC
+
0.14
0.85
Tortuosity (τ)
Ratio of Young’s modulus (E0/E)
E
0
/
E
=
−
0.24
ln
τ
−
0.98
+
0.40
0.95
Ratio of compressive strength(σC0/σC)
σ
C
0
/
σ
C
=
−
0.26
ln
τ
−
0.97
+
0.37
0.94
Ratio of yield strength (σS0/σS)
σ
S
0
/
σ
S
=
−
0.24
ln
τ
−
0.98
+
0.31
0.92
Ratio of shear modulus (G0/G)
G
0
/
G
=
−
0.24
ln
τ
−
0.99
+
0.43
0.94


## 6. Conclusions

Different types of sandstones are imaged using micro-CT technology. The rock porosity, the texture coefficient (TC), and the tortuosity of the rock skeleton are used to characterize the microscopic structure of the rock. Meanwhile, pore scale finite element mesh models of the rock with different pore structures are established. Combining with the rock mechanical parameters obtained in the laboratory, the microscale numerical simulation of uniaxial compression of the rock is carried out by ANSYS software. The following conclusions can be made:(1)
Through the uniaxial compression simulation of 12 digital cores in theX, Y, and Z directions, the optimal REV of different rock samples is determined.(2)
In model S1∼S4, when the strain is about 0.5%, the rock reaches the elastic limit and enters the plastic yield state. In model S5300-1∼S5400-2, there is little difference in porosity, the compressive strength of rock is about 100 MPa, and the yield strength is about 60 MPa. When the strain is about 1%, plastic failure occurs in rocks.(3)
The relationships between the mechanical parameters of rock and the porosity, texture coefficient (TC), and the tortuosity of rock skeleton are fitted quantitatively. With the increase of porosity, mechanical parameters of digital rock samples decrease linearly. Moreover, with the increase of TC, mechanical parameters of the rock are logarithmically decreasing. In addition, the tortuosity of the rock skeleton is also highly correlated with the mechanical parameters of the rock, i.e., Young’s modulus (R2 = 0.95), compressive strength (R2 = 0.94), yield strength (R2 = 0.92), and shear modulus (R2 = 0.94). The tortuosity of the rock skeleton is prior to the other two parameters in reproducing the anisotropy of the rock microstructure with the highest fitted correlation coefficient. Thus, the tortuosity of the rock skeleton, firstly presented in this paper, is believed to be more feasible to reveal the impacts of the microstructure of the rock on its mechanical properties.



---

*Source: 9085045-2020-02-24.xml*