*Article* **Pore Fractal Characteristics of Lacustrine Shale of Upper Cretaceous Nenjiang Formation from the Songliao Basin, NE China**

**Tian Dong 1, Lei Kang 2,\*, Yifan Zhang <sup>1</sup> and Yuan Gao 1,3**


**Abstract:** The evaluation of pore characteristics in lacustrine shales is of significance for shale oil and gas exploration in terrestrial sedimentary basins. However, due to the complex depositional environments, there have been few reports on the quantitative evaluation of lacustrine pores. In this study, we carry out experiments, e.g., scanning electron microscopy (SEM), nitrogen adsorption, X-ray diffraction (XRD), and rock pyrolysis, on sixteen lacustrine shale core samples of the Late Cretaceous Nenjiang Formation (K2n1+2) in the Songliao Basin. At the same time, we bring in the fractal dimensions to quantitatively characterize the pores of shale samples in K2n1+2. The kerogen type of shales in the Nenjiang Formation is mainly type I. The shales of K2n1+2 develop a large number of mineral intergranular pores and a small number of organic pores and mineral intragranular pores. The primary minerals of the samples are clay minerals and quartz. The total organic carbon (TOC) has a correlation with the pore physical parameters, whereas the content of minerals does not correlate significantly with pore development. The average pore fractal dimensions obtained from the analysis of nitrogen adsorption data are D1 = 2.5308 and D2 = 2.6557. There is a significant negative correlation between the fractal dimensions (D1 and D2) and the TOC content: the lower the fractal dimensions, the higher the TOC, and the more the free oil. In low maturity shales with higher TOC content, due to hydrocarbon generation, larger pores with a regular shape (small fractal parameters) can be produced, which can store more free oil. The shale reservoir evaluation of K2n1+2 in the Songliao Basin should be focused on shales with higher TOC and lower pore fractal parameters. Our work provides a reference for the shale oil evaluation of K2n1+2 in the Songliao Basin and complements studies on lacustrine pore characteristics.

**Keywords:** shale oil; pore structure; fractal dimension; lacustrine; Songliao Basin

## **1. Introduction**

With the increasing demand for oil and gas resources, shale oil, as a kind of unconventional resource, has become a research focus [1–5]. Different from conventional oil and gas reservoirs, shale reservoirs tend to have a more complex pore structure, diverse pore morphology, and stronger heterogeneity [6–9]. The study of shale pore characteristics and heterogeneity is of great significance for oil and gas storage, migration, and exploration [1–5].

Previous studies on pore characteristics have been focused on marine shales, such as the Barnet Shale in the Fort Worth Basin, the Eagle Ford Shale in the Western Gulf Basin, and the Bakken Shale in the Williston Basin in the United States [10–13]. China has also made breakthroughs in the exploitation of shale oil and gas resources, especially with the recent development in the marine shale Ordovician Longmaxi Formation and Cambrian Qiongzhusi Formation in the Sichuan Basin [14–21]. Compared to marine shales,

**Citation:** Dong, T.; Kang, L.; Zhang, Y.; Gao, Y. Pore Fractal Characteristics of Lacustrine Shale of Upper Cretaceous Nenjiang Formation from the Songliao Basin, NE China. *Appl. Sci.* **2023**, *13*, 4295. https://doi.org/ 10.3390/app13074295

Academic Editor: Eduardo Ferreira da Silva

Received: 22 February 2023 Revised: 23 March 2023 Accepted: 27 March 2023 Published: 28 March 2023

**Copyright:** © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

the pore systems of lacustrine shales are less known due to more complex depositional environments and lithofacies in lacustrine shales, which may cause different pore structure characteristics and different influences on shale oil/gas enrichment and production. Therefore, the quantitative evaluation of lacustrine pores is significant but has been seldom reported [18–22].

Shale contains a large number of nano-micron pores and has a very strong heterogeneity. The traditional technical methods for pore structure studies have limitations and cannot effectively describe the pore structure and surface morphology. In this regard, the characterization of shale pore structure requires a multi-technology, multi-scale, and multi-faceted method [2,3,8,23–32]. Mandelbrot proposed the fractal geometry theory to characterize the special structure of many objects with complex shapes in nature; the geometry of shale pores is a typical example of such shapes [33–35]. Pfeifer et al. were the first to demonstrate the fractal characteristics of reservoir pores by molecular adsorption [36]. Pape et al. used SEM to observe rock samples and found that pores in the range of 0.2 to 5.0 μm in sandstones, shales, and carbonates have good fractal properties [35]. Katz and Thompson also supported the theory by confirming the fractal characteristics of several sandstones [37]. With the research on the pore structures of shale, fractal theory has become an effective method to characterize the complexity of pore structures in shale [38–40].

The Songliao Basin in northeast China is a large Cretaceous terrestrial basin and a petroleum-bearing sedimentary basin [41–45]. There are two sets of shales as source rocks and shale oil reservoirs: the Qingshankou Formation (K2q) and the Nenjiang Formation (K2n) [42,44,46,47]. Multiple methods have been used to characterize the shale pore structure in K2q and K2n of the Songliao Basin, including the application of fractal dimensions. Zhang et al. analyzed the geochemical, pore structural, and fractal characteristics of the K2n samples from an outcrop section located at the junction zone of the Central Depression and Southeastern Uplift and showed that the fractal dimension of the K2n shale is positively correlated with the organic matter, clay content, total pore volume, and specific surface area [47]. Liao et al. measured the content and proportion of retained hydrocarbons at different stages of maturity for outcrop samples from K2n and found that the K2n shale in Songliao Basin has the strongest retention capacity but the weakest expulsion capacity compared to shales in other basins [48]. He et al. showed variation in the permeability of the low-maturity shale of K2n with increasing temperatures [49]. Similar to the previous study, Zhang et al. analyzed the shale pore structures of samples from multiple boreholes and inferred a positive correlation between the fractal dimension and the organic matter content [50]. On the contrary, Wang et al. analyzed the pore structure parameters and the fractal dimension of the K2q shale in the Songliao Basin and found a U-shaped relationship between the fractal dimensions and the TOC content [51]. The different findings on the relationship between the pore structure and shale composition may be due to different studied intervals and discontinuous or weathered samples from multiple cores and outcrop sections.

The SK-1 scientific drilling project, under the International Continental Scientific Drilling Program framework, has recovered ~165 m continuous cores of lacustrine shales from the Nenjiang Formation in Songliao Basin. In this study, we collected 16 samples of shale in the first and second members of the Nenjiang Formation (K2n1+2) of the SK-1(S) well and conducted scanning electron microscopy, nitrogen adsorption, X-ray diffraction (XRD), and rock pyrolysis analysis. The primary objectives of this study were: (1) to study the pore morphology and pore structure of lacustrine shales, (2) to analyze the relationship between the mineral content and pore space, and (3) to introduce fractal dimensions to characterize the pore structure and analyze its relationship with TOC. Our research complements the vertical variation in the pore characteristics of continuous shale samples and provides a reference for shale oil evaluation of the Nengjiang Formation in the Songliao Basin.

## **2. Geological Setting**

The Songliao Basin, which covers roughly 26 × 104 km2, is a large terrestrial basin in Northeastern China [42–44] (Figure 1). Three tectonic stages are recognized in the Songliao Basin: the syn-rift stage (Late Jurassic), the post-rift stage (Early to Middle Cretaceous), and the structural inversion stage (Late Cretaceous) [44,45]. After the structural inversion stage, according to the topographical characteristics, the basin was divided into six structural units: the northern dipping area, the northeast uplift area, the southeast uplift area, the southwest uplift area, the western slope area, and the central depression area [46,47]. K2n was deposited in the post-rift stage, with a development of gray-black lacustrine mudstone, oil shale, marl, and tuff [45,46]. At the beginning of the K2n deposition, the deposition rate accelerated, and the lake expanded rapidly to cover the whole basin. In the later stage, thick gray and black mudstones developed in the deep lacustrine environment until the basin was gradually uplifted, resulting in a shallow lake and the deposition of lighter-color shales [45,46].

**Figure 1.** (**left**): (**a**) Location of the Songliao Basin in China; (**b**) the main structural divisions of Songliao Basin and drilling sites of SK-1(S) and SK-1(N). (**right**): Lithology and sampling depth of K2n1+2 in SK-1(S).

The SK-1 scientific drilling site is located in the central depression area (Figure 1) and composed of two boreholes—the south borehole (SK1(S)) and the north borehole (SK1(N)) [52,53]. The coring ratio of SK-1 is as high as 96.46% [46]. The two cores are correlated by a basin-wide oil shale layer of the Nenjiang Formation (Fm) [46]. The coring layer of the south borehole of SK-1 is from the top of the third member of the Quantou Formation to the bottom of the second member of K2n. This paper focuses on the first and second members of the Nenjiang Formation in the south borehole of SK-1. The oil shale at the bottom of the second member of K2n represents a lake-level peak during the evolution of the Songliao Basin, with a sedimentary range far beyond the present basin boundary. The lithology of K2n1+2 is mainly black and gray shale and mudstone (Figure 1) [41–46].

## **3. Samples and Methods**

Sixteen lacustrine shale samples were obtained from K2n1+2 with a depth range of 972.26–1128.17 m in the south borehole of SK-1 (Figure 1), and scanning electron micro-scopy, nitrogen adsorption, X-ray diffraction (XRD), and rock pyrolysis analysis were conducted.

Field emission scanning electron microscope (FE-SEM) was used to observe the 2D structure of the pores. Before the samples were observed under the microscope, argon ion polishing technology (Ar-BIB) was used to etch the surface of the samples (Leica RES102 Ion Milling System of Leica company, Wetzlar, Germany). The cross-sections of the samples were bombarded by three independent argon ion beams under a voltage of 7.0 kV and current of 2.6 mA. The polishing time was 5–7 h for each sample. The polished surface was crescent-shaped with an area of about 2 mm2. For SEM analysis (Apreo FE-SEM of FEI company, Hillsboro, USA, 5.0 kV, 0.40 nA), the samples were treated by carbon plating to make the surface conductive. Then the samples were placed in the electron microscope instrument, the distance between the SEM lens and the top of the sample table was adjusted to about 4 mm, and the backscattering mode (BSE) was used to observe the rock surface. Using a high-resolution electron microscope, the distribution of the organic matter, pore size, shape, contact relationship between the pores and organic matter was clearly observed. At the same time, the plane distribution of inorganic minerals, such as pyrite, and their contact relationship with the organic matter and pores were also observed [54–56] (Figure 2).

**Figure 2.** Pore types of shales by SEM in K2n1+2, Songliao Basin. (**a**) Mudstone sample at 1082.4 m; (**b**) shale sample at 978.2 m; (**c**) mudstone sample at 1052.1 m; (**d**) mudstone sample at 1127.8 m.

Low-pressure N2 adsorption is an effective method to analyze the nanopore of shale [55,57–62]. The samples were ground to a 40–80 mesh powder (~178–425 μm particle size) and dried at 110 ◦C for 8 h to remove any adsorbed moisture and volatiles. An N2 adsorption experiment was carried out on the ultra-high performance automatic gas adsorption system ASAP 2020 produced by the Micromeritics Company (Georgia State, USA). The range of the pore size distribution was 1.7 to 300 nm. The equivalent surface areas were determined with the Brunauer–Emmett–Teller (BET) method. The pore volumes were calculated using the density functional theory (DFT) model. The average pore diameters were calculated using the adsorption branch of the Barrette–Joyner–Halenda (BJH) model [54,57–66].

The whole rock mineral composition of shale was analyzed by an X-ray diffractometer produced by Rigaku Motor Co., Ltd., Tokyo, Japan. The test conditions were carried out with reference to industry standard SY/T 5163-2010. The whole rock mineral analysis was carried out on samples of 300 mesh powder. With reference to the K value of the international standard sample, the mass percentage of each mineral was calculated using MDI JADE software (6.5, Central South University, Hunan, China).

The analysis of the total organic carbon content (TOC) was carried out in accordance with the national standard GB/T 19145-2003. The Leco CS230 carbon and sulfur instrument was used to analyze the TOC. In the first step, a sample powder of ~10 mg was weighed by electronic balance and put into the porous porcelain crucible (which was heated in a muffle furnace to 1000 ◦C for 2 h). In the second step, a sufficient amount of 12.5% HCl was added to the electric heating plate at 60 ◦C for 2 h until the reaction was complete. In the third step, the crucible was put into the filter container and washed with distilled water every half an hour to one hour for three days. In the fourth step, the crucible was dried in a furnace at 60 ◦C and cooled, and the TOC was determined.

A Rock-Eval 6 pyrolyzer, manufactured by Vinci Technologies, was used for pyrolysis. The samples were crushed into 0.07–0.15 mm, and then weighed and loaded onto sample crucibles, which were put into a pyrolysis oven. The pyrolyzed samples were put into an oxidation oven. Determination parameters: free hydrocarbon (S1) was analyzed at 300 ◦C for 3 min to detect the free hydrocarbon content in the sample (mg/g); pyrolysis hydrocarbon (S2) was detected at a temperature rate of 25 ◦C/min between 300 and 650 ◦C to detect the pyrolysis hydrocarbon content in the sample (mg/g). Tmax is the pyrolysis temperature corresponding to the highest point of the S2 peak (◦C). S3 is the amount of CO2 produced by the pyrolysis of organic matter (mg CO2/g Rock), which is generally used to calculate the oxygen index. The pyrolyzable carbon (PC; from the pyrolysis stage), residual carbon (RC; from the oxidation stage), total organic (TOC; sum of RC and PC), hydrogen index [HI; (S2/TOC) × 100], oxygen index [OI; (S3/TOC) × 100], oil saturation index [OSI; (S1/TOC) × 100] were determined and calculated. The analysis and calculation procedures were carried out in accordance with the national standard GB/T 18602-2012.

Correlations among the variables were evaluated by the R2, the F-test, the T-test, and the *p*-value. In general, the degree of the variables' statistical correlation (R2) was judged by the magnitude of the correlation coefficient (1.0–0.8: very strong correlation, 0.8–0.6: strong correlation, 0.6–0.4: moderate correlation, 0.4–0.2: weak correlation, 0.2–0.0: very weak or no correlation); the F-test was used for the overall regression equation significance test; the T-test was used to determine the significance of the coefficient of variation; and the *p*-value was used as a measure of the T-test and F-test [67–70]. In simple linear regression, the square of the T-statistic of the T-test for the explanatory variables is equal to the F-test for the overall significance of the equation of the F-statistic, and the two tests are equivalent [63–68]. We measured the correlation and significance of the regression equation using the T-statistic, *p*-value, and R2: T-absolute value > 1.65, *p* < 0.10 regression coefficient is significant at the 10% level; T-absolute value > 1.96, *p* < 0.05 regression coefficient is significant at the 5% level; T-absolute value > 2.58, *p* < 0.01 regression coefficient is significant at the 1% level [67–70].

#### **4. Results**

*4.1. 2D Pore Structure*

4.1.1. SEM Analysis

Under the scanning electron microscope, the types of pores observed in the SK1-S-2, SK1-S-4, SK1-S-5, SK1-S-6, SK1-S-9, and SK1-S-11 samples were intergranular pores and intragranular pores, but no organic pores were observed. The pores were mainly intergranular pores, including mineral intergranular pores, organic-mineral intergranular pores, and clay mineral intergranular pores. Most of the mineral intergranular pores and organic–mineral intergranular pores were narrow or acicular, and a few were round or sub-circular. The pore diameters observed under the microscope were about 0.1–1 μm. The intergranular pores of the clay minerals were mostly in the shape of long strips, and the pore size was about 0.1 μm (Figure 2a).

Under the scanning electron microscope, the pore types of samples SK1-S-1, SK1-S-8 and SK1-S-12 were observed as mineral intergranular pores, organic–mineral intergranular pores, mineral intragranular pores, and organic pores. Mineral intergranular pores were mainly observed in SK1-S-12 with a few organic pores. The pores of SK1-S-1 and SK1-S-8 were mostly circular or subcircular, and their sizes were between 0.1 and 0.5 μm (Figure 2b).

The pore types observed in SK1-S-3, SK1-S-13, and SK1-S-7 were mainly intergranular pores and mineral intragranular pores, but no organic pores were observed. Intergranular pores of pyrites have developed, which were in the form of microspheres with a pore size of less than 0.1 μm (Figure 2c).

Intragranular dissolution pores are those formed by dissolution inside of the particles of feldspar, carbonate, and other soluble minerals. They are formed by port-like dissolution of the particles around the pores. Intergranular dissolution pores are larger than intragrain dissolution pores. The intragranular pores observed in the SK1-S-16 sample were developed in the calcite (Figure 2d).

In summary, the shale of K2n1+2 in SK-1(S) has mainly developed intergranular pores, including mineral intergranular pores, organic–mineral intergranular pores, and clay mineral intergranular pores. A small number of mineral intragranular pores have developed, including calcite intragranular dissolution pores. The distribution of the organic pores was uneven, and the number of organic pores was small, which shows that K2n1+2 in Songliao Basin has oil generation potential.

#### 4.1.2. PCAS Analysis

The pore (particle) and fracture image identification and analysis system (PCAS) is professional software for the quantitative analysis of the pore system and fracture system. The software can automatically obtain geometric and statistical parameters and calculate the depth, area, and fractal dimensions of individual pores. The pores are identified by importing the SEM pictures of the core samples into the software and adjusting the grayscale threshold of the pictures. In the SEM pictures, the areas with pores have a darker color compared to other areas, and the adjustment process can make the darker parts of the pictures appear darker and the brighter parts appear brighter, thus identifying the physical parameters of the pores [63].

In this study, the probability entropy and shape factor obtained were mainly calculated by PCAS software [63]. The probability entropy is mainly used to describe the directionality of the 2D structure of the pores, which is calculated as follows (Equation (1)):

$$H = -\sum\_{i=1}^{n} P\_i \log\_n P\_i \tag{1}$$

where *H* is the probability entropy; *Pi* denotes the percentage of pores in a specific range, for example, between 0◦ and 10◦ when *i* = 1, in the 2D direction ranging from 0◦ to 180◦, and divided into eighteen parts, i.e., *n* = 18, the probability entropy has values between 0 (all pores are arranged in the same orientation) and 1 (pores have random orientation).

The shape factor (*ff*) is a descriptor of the characteristic shape and reflects the roundness and roughness of the hole edges. It has a maximum value of 1.0 (for a circle) or 0.785 (for a square). Its calculation formula is as follows (Equation (2)):

$$ff = 4\pi \text{S/C}^2\tag{2}$$

where *S* is the area of the pore, and *C* is the perimeter of the pore.

The scanning electron microscope images of the shale samples from K2n1+2 were imported into PCAS software, and the organic and inorganic pores could be identified by adjusting the gray level of the images. Taking SEM images of sample SK1-S-1 as an example, the gray level of the SEM images was adjusted to identify the organic pores. Each organic pore was colored with different colors, and the other inorganic pores were eliminated. Then some parameters of each organic pore with different colors were obtained by further calculation in the software (Figure 3) [63].

**Figure 3.** Original (**a**) and colorized (**b**) images of organic pores in sample SK1-S-1 using PCAS software processing.

Using PCAS software, the scanning electron microscope images were divided into organic and inorganic pores, the probability entropy and shape factors of different types of pores were analyzed separately, and the distribution of the pore size were analyzed at the same time (Table 1). By processing 31 scanned images of K2n1+2 with PCAS software, and rejecting the unqualified images and obtaining each parameter of the pores, the following results were obtained by comprehensive analysis: (i) The shape factor of the organic pores was mainly concentrated between 0.3 and 0.7. The average shape factor was 0.56. The average probability entropy of the organic pores was ~0.6 (Table 2). (ii) The shape factor of the inorganic pores was mainly concentrated between 0.5 and 0.6. The average shape factor was 0.54. The average probability entropy of the inorganic pores was ~0.9 (Table 2). (iii) As a whole, the number of inorganic pores was much larger than the number of organic pores.

The use of backscattered electron imaging is beneficial to observe the shape, size, and distribution pattern of nanoscale pores on the smooth and flat surface of the sample after argon ion polishing [10,15–18]. SEM analysis and PCAS analysis showed the 2D pore structure characteristics. As image-processing tests on the number of pores have some limitations, these were specific analyses of typical micro-regions, aiming to complement the information on pore characteristics. To identify the pore types and quantitatively characterize the pore structures, nitrogen adsorption and desorption experiments were applied and are described in the following section.


**Table 1.** Parameters of each organic pore obtained from sample SK1-S-1 after PCAS software processing. Length: longest diameter; width: shortest diameter; average pore size: average diameter (longest+ shortest)/2.

**Table 2.** Average pore parameters of organic pores and inorganic pores obtained from K2n1+2 by PCAS software.


## *4.2. N2 Adsorption and Desorption Isotherms*

According to IUPAC classification, there are four types of N2 adsorption isotherms [61,62]. In this study, there were obvious hysteresis loops in the samples, and all of the hysteresis loops were closed. The existence of a hysteresis loop shows that evaporation from the pores is obviously different from condensation in the pores, and that capillary condensation occurs in mesoporous pores [71]. Useful information about the pore structure can be obtained from the shape of the hysteresis loop. The hysteresis loops of 16 samples were analyzed and it was found that there was only one type of hysteresis loop: type H3 (Figure 4). H3 hysteresis occurs with the accumulation of plate-like particles that produce slit-shaped pores, which corresponds to the slit-shaped and acicular intergranular pores of clay minerals observed under a scanning electron microscope (SEM) (Figure 2).

**Figure 4.** Nitrogen adsorption–desorption isotherms of K2n1+2 samples.

The sample crush size may have an influence on the behavior of low-pressure gas adsorption and desorption. Through experimental analysis, Hazra et al. showed that the macro-pore volume in N2 adsorption and desorption increased if the particle-crush size was too small [71]. In this study, the crushed size of the samples was in the range of 178–425 μm, which was chosen to take into account the shale matrix and the particle size distribution within the mineral. This size range falls into the medium particle size range in the study of Hazra et al. (2018), which showed no significant influence on the N2 adsorption and desorption results [71]. At the same time, the surface chemistry of the sample also affects the results of the pore characterization by nitrogen adsorption [72]. The presence of a catalyst on the sample surface attracts the adsorption of nitrogen molecules, resulting in a decrease in the ability of gas molecules to enter the pores and a shift in the measured pore size distribution curve toward smaller pore sizes [72]. The adsorption behavior of nitrogen molecules is also affected by the hydrophilicity or hydrophobicity of the sample surface. The hydrophilic surface attracts nitrogen molecules, resulting in a high measured specific surface area, while the hydrophobic surface repels nitrogen molecules, resulting in a low measured specific surface area [72]. This effect was minimized in the pre-treatments of the experiments to accurately characterize the pore structure.

Some useful information about the surface area, pore volume, and pore size distribution (PSD) can be provided by N2 adsorption experiments. The equivalent surface area was calculated by the multi-point BET equation (equation describing the theory of adsorption of multilayers based on the Langmuir equation, proposed by Brunauer, Emmet and Teller), and the BET repeat calculation was corrected using "BET Surface Identification" (BETSI) software [23]. The range of the BET surface area of the samples was 9.69–65.91 m2/g. The pore volume and average pore size were estimated by the Kelvin equation and the Barret–Joyner–Halenda (BJH) model. The pore volume range was 0.0444–0.1187 cm3/g (average 0.0768 cm3/g), and the average pore diameter was 8.22–15.40 nm (average 10.10 nm), which is slightly larger than that of K2q in the Songliao Basin [51] (Table 3).


**Table 3.** Pore parameters of K2n1+2 by low-pressure N2 adsorption.

## *4.3. Fractal Dimensions*

At present, several methods to calculate the fractal dimension based on the gas adsorption method have been proposed, including the BET model method, the PSD model method, and the Frenkel–Halsey–Hill (FHH) model method [64,65], but the FHH model qA the most widely used and the most effective method [39,65–67]. In this study, the fractal dimensions of the 16 shale samples of K2n1+2 in the Songliao Basin were calculated by the FHH model. The calculation formula (Equation (3)) is as follows:

$$\ln\left(\frac{v\_n}{v\_{n0}}\right) = A\left[\ln\left(\ln\left(\frac{P\_0}{P}\right)\right)\right] + \beta\tag{3}$$

where *vn* is the volume of adsorbed gas at pressure equilibrium *P*; *vn*<sup>0</sup> is the volume of the saturated adsorption of the monolayer; *A* is a power-law exponent related to the fractal dimension; *β* is the intercept; *P*<sup>0</sup> is the saturated vapor pressure of the gas. The value *A* can be calculated by plotting the *ln vn vn*<sup>0</sup> and ln ln *<sup>P</sup>*<sup>0</sup> *P* of the gas adsorption isotherm data, and the slope of the straight line is equal to *A*. The method for calculating the fractal dimension *Dn* (Equation (4)) is as follows:

$$D\_n = A + \mathfrak{Z} \tag{4}$$

Figure 5 shows the FHH diagram of the sample SK1-S-2. At the relative pressures *P <sup>P</sup>*<sup>0</sup> of 0–0.5 (R2 = 0.9999) and 0.5–1 (R<sup>2</sup> = 0.9971), there are two different linear segments because the gas adsorption mechanism was different at different relative pressures (Appendix A) [63,72–77]. The characteristics of the pores were also different. Considering the difference in the gas adsorption behavior between *<sup>P</sup> <sup>P</sup>*<sup>0</sup> < 0.5 and *<sup>P</sup> <sup>P</sup>*<sup>0</sup> > 0.5, the fractal dimensions of the samples were calculated in these two regions, which are defined as D1 and D2. D1 is monolayer–multilayer adsorption controlled by van der Waals forces, and D2 represents the adsorption behavior controlled by capillary condensation at a higher *P <sup>P</sup>*<sup>0</sup> [73–77].

**Figure 5.** Fractal calculation results (intercept: *A* of Equation (4)) from FHH model (Equation (3)) of N2 adsorption at low temperature for sample SK1-S-12 from K2n1+2 of the Songliao Basin, NE China.

The general fractal dimension was between 2 and 3, which was determined by the geometric irregularity and roughness of the surface. The larger value of D represents a more complex and irregular surface (with a value close to 2, the surface is smoother and flatter; with a value close to 3, the surface is rougher and the heterogeneity is stronger) [35–37]. The fractal dimension D1 indicates the irregularity of the sample surface (the larger the fractal dimension D1, the more irregular and rougher surface) [65]. The fractal dimension D2 indicates the irregularity of the pore structure (the larger the fractal dimension D2, the higher the heterogeneity of the pore distribution, the smaller the pore size, the higher the liquid/gas surface tension, and the lower the gas adsorption capacity) [65]. Fractal dimensions D1 and D2 represent the pore surface fractal dimensions and pore structure fractal dimensions, respectively (Figure 5, Appendix A).

The range of the fractal dimension D1 was 2.4406 to 2.6174, and the average value was 2.5308. The fractal dimension D2 ranged from 2.5036 to 2.7197, and the average value was 2.6556. The fractal dimension D values were all higher than 2.5 (Table 4). These results indicate that the shale pores have good fractal properties in both their internal structure and surface. The D2 values were larger than the D1 values of the K2n1+2 shale samples, demonstrating that the pore internal structure has stronger complexity and heterogeneity than the pore surface, while the smaller pores have stronger complexity and heterogeneity

than the larger pores. There was a positive correlation (R<sup>2</sup> = 0.7372, T = 6.56, *p* = 0.00001) between D1 and D2 of K2n1+2 in the Songliao Basin, indicating that these two fractal dimensions can be used to characterize the pore structure and pore surface of the shale [51] (Figure 6).


**Table 4.** Fractal dimensions calculated by the FHH model.

**Figure 6.** Relationship between fractal dimensions D1 and D2 of K2n1+2 shale.

## *4.4. Mineralogy*

All samples contained clay, quartz, feldspar, calcite, pyrite, and other minerals. Clay minerals were the most abundant, with a content between 40.2 wt.% and 56.1 wt.% and an average of 48.92 wt.%. The second most abundant mineral was quartz, which accounted for 17 wt.% to 35.6 wt.% of the samples (average 27.79 wt.%). The content of carbonate (calcite and dolomite) was between 3.4 and 24.3 wt.%, and the content of feldspar was between 6.5 wt.% and 15.1 wt.%. The content of pyrite was low, not more than 5.5 wt.% (Table 5).


**Table 5.** Mineral composition of K2n1+2 by XRD (unit: wt.%).

## *4.5. Organic Geochemistry*

The hydrocarbon generation potential of shale can be evaluated according to the TOC content and hydrogen enrichment. Although the crushed particle size has an effect on the rock analysis results, the rock pyrolysis and the total organic carbon content (TOC) were analyzed more than twice in parallel strictly according to the national standards GB/T 18602-2012 and GB/T 19145-2003 [78] so as to guarantee analytical precision. Flame ionization detector (FID) saturation during pyrolysis is a common source of error found in the analysis of organic-rich rocks, especially those that generate large quantities of hydrocarbons from kerogen during pyrolysis (large S2 yields). FID saturation usually occurs when analyzing type II kerogen with thermal evolution in the immature-to-medium oil window range or with high TOC samples. Considering the TOC content of the samples and the type of kerogen, the experimental results would not have FID saturation [79–82].

The S2 value of the K2n1+2 shale in the Songliao Basin varied greatly, with an average of 24.43 mg HC/g rock, ranging from 0.39 to 75.41 mg HC/g rock. The content of TOC ranged from 0.3 to 8.74 wt.%, with an average of 3.21 wt.% (Table 4). The S2-TOC diagram shows that there is type I kerogen in the shale of K2n1+2 [45] (Figure 7). S2 and the TOC had a good linear correlation, and the correlation coefficient was R2 = 0.9781 (Figure 7). The slope (S2/TOC) indicates that the trend-based HI (hydrogen index) had a value of about 980 mg HC/g TOC (Slope × 100 = HI = 980 mg HC/g TOC), and the sample was of type I kerogen oil-prone.

The value of free hydrocarbon (S1) was between 0.04 and 3.38 mg HC/g rock, with an average of 0.79 mg HC/g rock. S1 had a positive correlation (R<sup>2</sup> = 0.8973, T > 2.58, *p* < 0.01) with the TOC (Figure 7). The range of HI was 131 to 975 mg HC/g TOC, and the average was 600.8 mg HC/g TOC. The value of the oil saturation index (OSI) ranged from 7.87 to 43.9 mg HC/g TOC, and the average was 18.62 mgHC/g TOC. The Tmax value represents the temperature at which the S2 peak reached its maximum, ranging from 434 ◦C to 443 ◦C, with an average of 439 ◦C (Table 6).

**Figure 7.** Relationship between TOC content and hydrocarbon generation potential of K2n1+2 shale. (**a**) S2 versus TOC plot showing the quantity and quality of the organic matter, kerogen type, trend HI (HIL), and dead carbon content. (**b**) Relationship between TOC content and S1.


**Table 6.** Rock-Eval data of K2n1+2 shale, Songliao Basin.

## **5. Discussion**

As both the generation and storage of shale oil and gas occur in shales, understanding the pore characteristics, e.g., the distribution of pore spaces, the fractal dimensions of the pores, etc., could assist in the study of shale reservoirs. Below, we discuss the relationships among the mineral composition (quartz, clay minerals, etc.), organic geochemical parameters (TOC content, free oil content S1, etc.), and shale pore parameters (e.g., BET surface area, average pore diameter, and pore volume), analyze the trends of pore evolution and the shale oil exploitation potentials of the K2n1+2 shales in the Songliao Basin, and make comparisons with similar pore shales.

## *5.1. Relationships between Fractal Dimensions and Pore Structure Parameters*

The average shape factor was 0.56 for organic pores and 0.54 for inorganic pores (Table 1). These indicate that both the organic pores and inorganic pores in the shale of K2n have similar near-square pore shapes. The average probability entropy was about 0.6 for the organic pores and about 0.9 for the inorganic pores (Table 2). These indicate that more than half of the pores in the shale of K2n have same directional property. This is consistent with the fracture pores shown by the hysteresis loop of the nitrogen adsorption isotherm.

The fractal dimension D1 ranged from 2.4406 to 2.6174, with a mean value of 2.5308. The fractal dimension D2 ranged from 2.5036 to 2.7197, with a mean value of 2.6556. The fractal dimensions indicate that the shale pores have good fractal properties in both internal structure and surface, the internal structure of the pores has stronger complexity and nonhomogeneity than the surface of the pores, and the smaller pores have stronger complexity and non-homogeneity than the larger ones. The BET specific surface area of the samples ranged from 9.69 to 65.91 m2/g, the pore volume ranged from 0.0444 to 0.1187 cm3/g (average 0.0768 cm3/g), and the average pore size was 8.22–15.40 nm (average 10.10 nm). The average pore diameter of the shale of K2n was negatively correlated with the fractal dimensions D1 (R2 = 0.5948, T > 2.58, *p* < 0.01) and D2 (R2 = 0.9716, T > 2.58, *p* < 0.01), indicating that the fractal dimension increases with the decrease in the average pore diameter (Figure 8). This trend is consistent with the shale data of the Qingshankou Formation recorded by Wang et al. and Cao et al. [51,78]. The lacustrine shale samples with higher D2 values had a more complex pore structure, which may be due to the fact that the shale samples with a smaller average pore size also contained more micropores.

**Figure 8.** Relationship between fractal dimensions and total pore structure parameters of K2n1+2 and K2q shale. The K2q compared to [78].

The fractal dimensions D1 and D2 increased with the increase in the BET surface area. The correlation coefficient between D2 and the BET surface area (R2 = 0.8002, T > 2.58, *p* < 0.01) was slightly higher than that between D1 and the BET surface area (R<sup>2</sup> = 0.7757, T > 2.58, *p* < 0.01, Figure 8). This correlation indicates that the irregularity of the shale pore structure and the roughness of the surface make the pore surface area increase, and the complex pore structure provides a greater contribution to the pore surface area than the pore roughness itself. The fractal dimensions (D1 and D2) also had a good positive

correlation (R2 > 0.5, T > 2.58, *p* < 0.01, R2 D2 > R2 D1) with the total pore volume of the lacustrine shale samples, indicating that the irregularity of the pores further increases the total pore volume (Figure 8). This trend is consistent with the shale data of the Qingshankou Formation recorded by Cao et al. (Figure 8) [78].

## *5.2. Relationships among the Sedimentary Environment, Fractal Dimensions, and Mineral Compositions*

The minerals comprising the K2n<sup>1</sup> shale in the Songliao Basin are clay minerals (50.2%), quartz (35.6%), calcite (2.3%), and dolomite (2.4%); the K2n<sup>2</sup> shale consists of clay minerals (48.7%), quartz (21%), calcite (8%), and dolomite (7.3%) (Tables 1 and 5). Compared to the K2n2 shale, the K2n1 shale consists of more quartz, less calcite, and comparable clay content. During the evolution of the Songliao Basin, K2n1 followed the trend of rapid lake development by the end of the Yaojia Formation, and the semi-deep lacustrine sedimentation covered the central depression of the basin [51]. Strong input from terrestrial sources may have caused the higher quartz content but less authigenic calcite. The basin area was further expanded during the K2n2 deposition period, and the whole basin was in a deep lacustrine environment with rich organic matter [51], consistent with a decrease in the quartz content and increases in the calcite and dolomite contents.

Through a comparison of the mineral content and the fractal dimensions, it was found that there was no correlation between the quartz and the fractal dimensions. Although the fracturing of quartz could further induce development of the pore space and gain more storage space, quartz may not be the main factor affecting pore development in the samples of this study. There was also no clear correlation between the clay minerals and fractal dimensions, indicating little or no correlation between the clay minerals and fractal dimension D1 (R<sup>2</sup> < 0.2, T > 1.96, *p* > 0.05) and weak correlation between the clay minerals and D2 (R<sup>2</sup> > 0.2, T > 1.96, *p* < 0.05, Figure 9). It is interpreted that the contribution of clay minerals to the pore structure complexity of K2n1+2 shale is limited. This is consistent with the Qingshankou formation results (Figure 9).

**Figure 9.** Relationships between mineral compositions and fractal dimensions of K2n1+2 and K2q shale. The K2q compared to [78].

#### *5.3. Relationships between TOC Content and Fractal Dimensions of Shale Pores*

The fractal dimensions (D1 and D2) show a clear negative correlation with the total organic carbon content, indicating that the fractal dimension decreases significantly with the increase in the TOC content (Figure 10a). This phenomenon is not consistent with the findings of previous studies on the Qingshankou and the Nenjiang Formations [46,50,51,65,78], which are discussed below.

**Figure 10.** (**a**) Relationship between TOC content and fractal dimensions D1 and D2 of K2n1+2 and K2q shale samples. (**b**) Relationship between TOC content and average pore diameter of K2n1+2 and K2q shale samples. The K2q compared to [78].

For the shales of the Qingshankou Formation within a depth range of 500–2300 m and a Tmax range of 427–460 ◦C, Wang et al. found that D1 and D2 decreased rapidly with increasing TOC content when the TOC content was less than 2.5~3.0 wt.% [51]. This is interpreted as the result of hydrocarbon generation and expansion in low-maturity shale to create more uniform mesopores, which further lead to lower D1 and D2 and a higher average pore size [51]. When the TOC value is greater than 2.5 to 3.0%, D1 and D2 increase with the increase in the TOC content [51]. This indicates that the organic pores may collapse due to higher TOC levels, burial, and thermal maturation, leading to more complex and smaller pores.

For the shales of the Nenjiang Formation, previous studies have indicated a positive correlation between the fractal dimensions and the organic matter content [49,51]. However, in the present study, we found a negative correlation between the fractal dimensions and the organic matter content (Figure 10). One possible explanation is that previous studies have selected discontinuous samples from multiple cores and/or outcrop sections. The analysis of our 16 continuously distributed shale samples from the SK-1(S) core indicates that, although the Tmax is low, the limited maturation process may still cause expansion of hydrocarbons to some extent, thus increasing the pore sizes of the K2n1+2 shale and decreasing the fractal dimensions (Figure 10b). The results of the study by Cao et al. support our explanation (Figure 10).

## *5.4. Relationships between Fractal Dimensions and Free Oil Content*

There was an obvious linear negative correlation between S1, OSI, and the fractal dimensions of the shale of K2n1+2. The correlation between S1 and D2 (R2 = 0.9464, T > 2.58, *p* < 0.01) was better than that of S1 and D1 (R<sup>2</sup> = 0.6603, T > 2.58, *p* < 0.01), indicating that D2 can be used to characterize the oil content of shale (Figure 11). This is consistent with the negative correlation shown between the Qingshankou Formation S1 and the fractal dimensions (Figure 11, R<sup>2</sup> D2q > R2 D1q). Shale samples with higher free oil levels always have a lower fractal dimension D2 value, which may be due to the fact that free oil is mainly stored in large pores, and these pores have lower D2 values. Contrary to shale gas adsorption, the molecular diameter of shale oil is much larger than that of shale gas. Therefore, shale oil is intelligently stored in pores with larger diameters, which causes the phenomenon that shale with higher free oil content usually has lower D2 [51]. As a result, combining the influence of the TOC content on the fractal dimensions and the relationship between the free oil and fractal dimensions, we conclude that the best area for exploitation in the low maturity shales of the Nengjiang Formation in the Songliao Basin would be the area with lower pore fractal dimensions and higher TOC content.

**Figure 11.** Relationship between fractal dimensions and oil content of K2n1+2 and K2q shale samples, (**a**) for S1; (**b**) for OSI. The K2q compared to [78].

## **6. Conclusions**

In this study, based on experiments using scanning electron microscopy, XRD, rock pyrolysis, and nitrogen adsorption, the shale pores of K2n1+2 in the Songliao Basin were studied, and the following conclusions are addressed based on the advantages of the fractal dimensions in characterizing pores.


**Author Contributions:** Conceptualization, L.K. and Y.G.; methodology, T.D.; software, L.K.; investigation, T.D. and L.K.; writing—original draft preparation, T.D., L.K. and Y.G.; writing—review and editing, T.D., L.K., Y.G. and Y.Z.; funding acquisition, Y.G. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the China Geological Survey Program (DD20190502).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** All data needed to evaluate the conclusions in the paper are present in the paper. Additional data related to this paper are available from the authors upon request.

**Acknowledgments:** We deeply appreciate the Research Institute of Petroleum Exploration and Development for the analytical equipment and help with the experiments.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Appendix A**

**Figure A1.** *Cont*.

**Figure A1.** Fractal calculation results from FHH model of N<sup>2</sup> adsorption at low temperature for all samples.

## **References**


**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
