Coal has been playing a critical role as the main energy source for many years in the world (Xie et al. 2019; Nandi et al. 2022). Methane gas, stored in the pore-fracture system of coal seams, is an associated product in the process of coal formation (Liang et al. 2021; Liu et al. 2021a, b, c, d). While posing risks to safe mining in underground coal mines, coalbed methane (CBM) is also a high-quality clean energy with wide occurrence, which can be used as fuel for power generation, industry, and residential applications (Liu et al. 2021a, b, c, d; Wang et al. 2022; Wang et al. 2021a, b, c, d). Therefore, CBM development will not only realize the acquisition of low-carbon energy but also promote the control of greenhouse gas ~ methane while reducing gas disaster risk in coal mining (Liu et al. 2017; Li et al. 2018a; Wang et al. 2021d; Qin et al. 2022). Figure 1 illustrates the annual CBM production of the United States, China, and Australia.

Annual CBM Production in the United States, China and Australia (Data from Ma et al. 2022)

The fractures and pores of coal serve as the place for methane storage and the channel for methane gas transport and production (Zhao et al. 2019; Wang et al. 2021b; Li et al. 2023; Xu et al. 2022). The pore-fracture characteristics, including the surface, volume, shape, size, connectivity and tortuosity of systems, will significantly influence the fluid deliverability of coal reservoirs, and are therefore of fundamental importance for the methane sorption, diffusion, seepage and production performance (Akhondzadeh et al. 2020; Li et al. 2021a, b, c; Liu et al. 2021a, b, c, d; Wang et al. 2019a, b; Yu et al. 2020). Thus, effectively identifying and characterizing the pore-fracture structure in coal has always been a prevailing topic in the theoretical study of CBM development (Li et al. 2021a; Ni et al. 2022). At present, the existing methods and analyses for characterizing the microscopic pore-fracture structure of coal are grouped into two main categories: fluid intrusion method and photoelectric radiation technology (Li et al. 2018a, b; Nie et al. 2015; Sun et al. 2020; Wang et al. 2018).

Fluid intrusion method The fluid intrusion method involves injecting a non-wetting fluid (such as mercury) or adsorbent gas (such as CO₂ and N₂) under pressure into the sample, recording the injection volume and pressure data of the fluid or gas, and applying various theoretical methods to calculate the structural characteristic parameters of pores, such as pore size distribution (PSD), volume, and specific surface area (Ross and Marc Bustin 2009; Okolo et al. 2015; Yu et al. 2017; Davudov et al. 2020). The most commonly used fluid intrusion test methods mainly include the low-temperature nitrogen adsorption method (Qi et al. 2017; Wang et al. 2019b), low-temperature carbon dioxide adsorption method (Nie et al. 2020), and mercury intrusion porosimetry (MIP) (Guo et al. 2014; Zhou et al. 2017; Zhang et al. 2021). Although the mercury intrusion method provides a simple and fast operation process for characterizing pore structure, due to the elasticity and compressibility of coal, during the high-pressure mercury injection process, the coal structure will be damaged to a certain extent, which will bring negative impact on the test results (Li et al. 2021a, b, c, 2020a). Regarding the low-temperature N₂/CO₂ adsorption method, while the experimental test ensures the integrity of the sample structure, the experimental environment temperature is required to be very low, and the test cannot be applied to aqueous samples (Xin et al. 2020). Additionally, due to the small molecular dynamics diameters of N₂ (0.364 nm) and CO₂ (0.33 nm), such methods can only identify micropores and mesopores with pore size less than 100 nm (Li et al. 2015a, b).

Photoelectric radiation technology Photoelectric radiation technology mainly includes photoelectric imaging technology and radiation detection technology. The photoelectric imaging techniques mainly include scanning electron microscopy (SEM) (Mastalerz et al. 2012; Tiwari et al. 2013), field emission scanning electron microscopy (FESEM) (Intrusion et al. 2014; Pan et al. 2016), and transmission electron microscopy (TEM) (Pan et al. 2015; Li et al. 2020a). The radiation detection techniques mainly include micro-CT scanning technology (μCT) (Li et al. 2020a, b, c; Wang and Liu 2021), nuclear magnetic resonance (NMR) (Bai et al. 2020; Zhao et al. 2021), focused ion beam scanning electron microscopy (FIB-SEM) (Fang et al. 2019), and small angle X-ray scattering (SAXS) (Luo et al. 2016; Zhao et al. 2020). Among them, the photoelectric imaging technology, such as SEM and FESEM, can provide two-dimensional image information of the surface topography of the structure, which is commonly applied to qualitatively characterize pore structure (Cardott and Curtis 2018; Roslin et al. 2019). Meanwhile, low-field NMR and μCT techniques in radiation detection technology have been proven to be effective for accurately measuring pores of various scales due to their rapidity and non-destructiveness. These techniques are widely used in the microstructural probing for rock and coal (Al-Yaseri et al. 2015; Liu et al. 2021a, b, c, d; Ni et al. 2017). Among them, the CT image-based technology can realize the visualization of the two-dimensional topography and three-dimensional spatial of the porous structure. It also provides quantitative structural information, including pore surface area, volume, size distribution, and connectivity, etc., which is increasingly used to probe the multi-scale structure of coal (Zhou et al. 2018; Wang et al. 2020b, 2021a).

In summary, in the existing studies, CT scanning technology is mostly used for the characterization of coal pores/fractures, and there are fewer studies on the detection of heterogeneous structures in coal, qualitative analysis of pore structure visualization, and quantitative assessment of gas flow capacity in pore space. In this study, two anthracite coals from the Qinshui Basin were probed and characterized in terms of pore structure using MIP, μ-CT, and image processing techniques. Three-dimensional pore structure was reconstructed from the sequence images, and the 3D spatial distribution, connectivity, and tortuosity of the pores were quantitatively analyzed. Subsequently, the permeability and its anisotropic characteristics were estimated with the model of pore network extracted from the 3D structure. The results of this work have an important impetus for the visualization and quantitative analysis of the three-dimensional complex structure of coal.

2.1 Coal sample preparation

Two samples of anthracite coal were collected respectively from the southern and northeastern parts of Qinshui Basin, Shanxi Province, China, and the basic characteristics of these two coal samples are listed in Table 1. The large blocks of coal were extracted from the active mining face, and the coals were transported to the laboratory in sealed packages for sample preparation and laboratory testing. The large blocks of coal were crushed into particles of 5–8 mm in size and cored into cylindrical coals of 25 mm in diameter and 50 mm in height, respectively, for microstructural characterization. The cylindrical cores were drilled perpendicular to the bedding surface of the coal seam. In this work, X-ray CT scanning and MIP were used to qualitatively and quantitatively study the pore and fracture characteristics of coal. The particle coals were used for mercury intrusion testing, and the cylindrical coals were used for the X-ray CT scanning. The experimental procedure of coal sample processing and characterization is shown in Fig. 2.

Coal sample processing and pore characterization testing

2.2 Mercury intrusion tests

Mercury intrusion, a method for measuring the pore structure characteristics of various porous materials, has been widely used in the characterization of the coal pore structure (Zhou et al. 2017). Mercury intrusion method probes the pores by applying controlled mercury pressure to a sample immersed in mercury. External mercury pressure is required to penetrate into pores in the porous media due to the non-wetting characteristics of mercury (Hao et al. 2013; Okolo et al. 2015). The required pressure is inversely proportional to the pore size, and the pore volume and the pore size can be evaluated with the mercury intrusion pressure data using the Washburn equation (Eq. (1)).

where: R indicates the pore radius (m); p indicates the intrusion pressure of mercury (Pa); σ is the mercury surface tension (N/m); θ indicates the contact angle (°) between the mercury and the pore surface.

The mercury porosimeter AutoPore IV 9620 produced by Micromeritics (as shown in Fig. 2) was used in this study, which can measure the pores ranging in size from 3 nm to 450 μm. Two coal samples were tested in this study: Sample 1 from the Northeast Qinshui coal and Sample 2 from the south Qinshui coal. Before the mercury intrusion tests, the coal particles should be dried in a vacuum environment to remove the moisture in the coal. During the tests, the coal particles were initially placed in the glass penetrometer and evacuated, after which the mercury was intruded into the penetrometer and filled the sample. The mercury intrusion began at 0 Pa, and the mercury began to intrude into the sample due to the existence of the external pressure. When the external pressure was large enough, all the pores were invaded, and no further mercury could be intruded. Mercury retracts spontaneously from the pores after the pressure is relieved, while, the hysteresis related to pore structure characteristics may appear during mercury extraction. The data of pressure and volume of intruded mercury were recorded throughout the tests. The pore structure parameters of the tested coal can be calculated due to these data.

2.3 CT image-based fracture characterization

Micro computed tomography (micro-CT), one of the non-destructive and non-invasive imaging techniques, is able to clearly detect the both cross-sectional and three-dimensional internal structure of the object without causing any damage (Ni et al. 2017; Wang et al. 2020a). In this work, a micro-CT system modeled as the mirco-CT XRM-520 was selected, featuring a maximum X-ray source power of 10 W, a maximum X-ray voltage of 140 kV, and a resolution ranging from 0.7 to 50 μm, as shown in Fig. 2. The imaging system of the equipment is mainly composed of four parts: an X-ray source, an precision sample stage, an X-ray detector, and an X-ray to image conversion system (Fig. 2). X-rays are emitted from the X-ray tube while the sample is rotated on the sample table. Through the photoelectric effect, the X-rays are attenuated as they pass through the sample. The attenuated X-ray beam is collected by a flat panel detector, and its cumulative attenuation intensity is measured. After receiving the signal input from the flat panel detector, the signal conversion system generates slice sequence grayscale images with an imaging reconstruction algorithm, and finally the 3D spatial structure of the sample will be reconstructed. Details of the X-ray CT scan resolution for each sample are listed in Table 2.

3.1 Mercury intrusion test results

3.1.1 Mercury intrusion data for coal

The internal pores in coal include effective pores and isolated pores. Effective pores refer to the pores that are interconnected, including both open pores and semi-closed pores as two basic types (Chen et al. 2018; Ren et al. 2021). Effective pores form the migration channels for CBM in coal seams and are of great significance to CBM well production (Li et al. 2015a, b). Because of the poor connectivity of the pore system, part of the mercury is trapped in the pore structure, resulting in the hysteresis loop of mercury removal in mercury intrusion experiment in coal. The shape and connectivity of the coal pores can be initially evaluated based on mercury intrusion and withdrawal curves, as well as the mercury removal rate. The mercury intrusion and withdrawal curves of the tested samples are shown in Fig. 3.

Mercury intrusion and withdrawal curves of coal

It shows in Fig. 3 that the maximum mercury injection amount of coal sample 2 is significantly higher than that of coal sample 1, which indicates that coal sample 2 has larger porosity. The mercury injection curve of coal sample 1 presents an inverse S shape, and the main section of mercury withdrawal curve decreases linearly. From the release of mercury pressure, the difference between mercury injection and withdrawal corresponding to the same pressure point gradually increases, forming an obvious hysteresis loop and a low mercury withdrawal ratio. This indicates that the pore system of coal sample 1 has a certain degree of connectivity, but its internal connectivity of pore system is weak, which is not conducive to the migration of CBM. For coal sample 2, most sections of the mercury injection and withdrawal curves vary linearly, and the volume difference between mercury injection and withdrawal corresponding to the same pressure point is small, showing a small hysteresis loop. This shows that the interconnectivity between pores in coal sample 2 is better, which can provide more flow channels for CBM transport in coal.

3.1.2 Characteristics of pore size development

With Hodot's classification, pore size in coal can be classified into five types: submicropore < 5 nm, 5 nm < micropore < 10 nm, 10 nm < transition pore < 100 nm, 100 nm < mesopore < 1000 nm, and macropore > 1000 nm (Ye et al. 2019). Pore size distribution (PSD) can reflect the probability distribution of different pore sizes in coal, further characterizing the development of pore sizes. The PSD of the tested coal samples and the cumulative amount of mercury injected into the pores at different scales were obtained from the mercury intrusion data and Washburn's equation, as shown in Fig. 4.

PSD and cumulative amount of mercury injected in the tested coal

The PSD curve of coal sample 1 shows multiple peaks, indicating that the micropore, transition pores, mesopores, and macropores in the sample are relatively developed. Moreover, the cumulative mercury injection amount of pores of different scales has a small difference, showing that the pore system in Sample 1 is well-developed across multiple scales. The PSD curve of coal sample 2 presents a double-peak feature, and the peaks appear between 6 nm–50 nm and 20,000 nm–120,000 nm, indicating that the micropores, transition pores, and macropores are relatively developed. The amount of mercury injected at the macropore scale is about 0.1006 mL/g, accounting for 80.34% of the total pore volume, indicating that the pores in Sample 2 are mainly macropores. In addition, the mercury injection data of the tested coal samples show that the calculated porosities of coal sample 1 and coal sample 2 are 8.3637% and 14.4349%, respectively.

3.2 CT-based characterization

The original sequence images of coal samples can be obtained based on the CT scanning. Subsequently, the 3D digital core of coal samples can be reconstructed by Avizo software, and the fracture structure parameters of coal samples can be quantitatively extracted. This allows for the establishment of a pore network model of coal samples, which can be used to analyze and evaluate the connectivity and permeability of the internal pore structure system of coal samples.

3.2.1 Original CT image

Each sample generates a series of 1800 × 1800 pixels 2D grayscale images, as shown in Fig. 5. The white areas indicate high-density minerals, the black indicates pores, and the gray indicates the coal matrix. We scanned the sample at two resolutions: 1 μm and 12 μm, and with the scanning resolution increased the microscopic pore structure became clearer in the original CT image. It can be seen from the image at a resolution of 1 μm that the material in the view field is mainly the coal matrix, with pores densely distributed within it. The extensional morphology of coal fractures is displayed in the images with 12 μm resolution. There are few white areas in the image of Sample 1, indicating that the coal cores 1-1 and 1-2 contain fewer minerals. The image of Sample 2-1 (12 μm) contains many concentrated white areas, indicating that coal core 2-1 contains a significant amount of high-density minerals. In the image of sample 2-2 (12 μm), there are many widely distributed white areas, suggesting that coal core 2-2 contains a large number of minerals mixed with the coal matrix and distributed throughout Sample 2-2.

Images of structural composition of coal based on CT scanning

3.2.2 Three-dimensional digital core reconstruction

The series of grayscale images obtained from CT scanning are noise-containing. Perform noise reduction processing on grayscale images and determine the spatial distribution of different-phase substances, such as pores, coal matrix, and minerals using the grayscale value information contained in the images. The above process involves the image filtering, denoising, and image segmentation. Based on the results of image segmentation, a reconstructed 3D digital core can be obtained. The detailed reconstruction process of the 3D digital core is shown in Fig. 6.

Reconstruction process of 3D digital coal core

With image segmentation, the detailed structural information of the coal sample can be extracted from the original image. It also allows for distinguish the pore space, coal matrix, and mineral distribution from other areas. In this work, the structural information of pores is extracted from the CT grayscale images with the help of a marker-based watershed segmentation algorithm. It is able to build up a 3D digital core of the coal sample by stacking multiple slices after segmentation. The real 3D pore structure of the coal sample is shown in Fig. 7.

3D coal pore structure

3.2.3 Quantitative analysis of coal structure

In this work, the reconstructed digital coal core can be divided into pore and solid structure components, which can be expressed by the phase function in the following way:

where: \(\overrightarrow{r}\) means the vector of a certain length and direction starting from the origin and pointing towards the pore or solid structure. Porosity is one of the most basic but important structural parameters of porous media, which represents the ratio of the volume of the pore structure to the total volume of the coal sample:

Similarly, the porosity of a 3D digital core can be calculated with the following expressions:

The porosities of the six reconstructed 3D digital cores were calculated with Eq. (4), and the porosities were obtained as shown in Table 3.

From the total pore structure, the connected pore structure can be extracted and analyzed based on the Axis Connectivity function in Avizo software, as shown in Fig. 8.

3D connected pore structure

From the comparison of the connected pore structure of samples, it is not difficult to find that the connected pore volume of Sample 1 is relatively smaller than that of Sample 2. Connectivity is commonly used to quantitatively describe the degree of connectivity of the pore structure within a porous medium and the expression is as follows:

where: ϕ′ represents the porosity of the connected pore structure, \(\phi\) represents the total porosity of the 3D pore structure.

The data of total porosity, connected porosity, and connectivity of the six 3D pore structures were calculated and shown in Fig. 9. It is clear that the connected porosity of coal sample 1 is significantly smaller than total porosity, while the connected porosity of coal sample 2 is only slightly smaller than the total porosity. Therefore, it can be inferred that the pores in coal sample 1 are mostly isolated pores, while there are a large number of connected pores in coal sample 2, which is in good agreement with the results of the mercury intrusion tests.

Porosity and connectivity of 3D pore structures

It can also be seen from Fig. 9 that the connectivity is influenced by the type of coal samples, and the connectivity of coal sample 1 and coal sample 2 varies within a certain range of intervals, respectively. Furthermore, the connectivity of coal sample 1 is lower while the connectivity of coal sample 2 is higher. Besides, the connectivity of the same kind of samples also changes to some extent with sampling. Compared with the variation of pore connectivity induced by different types of coal samples, the variation induced by different sampling methods shows a less significant degree.

Parameters such as total porosity and connectivity of the 3D pore structure reveal the basic pore structure information in coal. In order to further analyze the complex pore structure within the coal, the surface porosity is introduced here. When the surface porosity is 0 or approximately 0, it indicates that the pore structure within this surface is either not connected or poorly connected along this direction. Based on the descriptive statistics method, the surface porosity of the 3D connected pore structure was calculated and analyzed. The results can be seen in Fig. 10.

Statistical plots of surface porosity in x (red), y (green) and z (blue) directions for the 3D pore structures

As shown in Fig. 10, the surface porosity of coal sample 2 was relatively larger than that of coal sample 1, which presented a good agreement with the calculated porosity and pore connectivity. For the digital cores built with 1 μm resolution images (S1-1-1 μm and S2-1-1 μm), the difference of the surface porosity distribution among the three directions of the 3D pore structures is minimal, whereas the surface porosity distribution in the 3D direction differs significantly for the 12 μm digital cores. It indicates that the anisotropic characteristics of the pore structure at different scales will be different, and the large-scale complex pore structures tend to have more significant heterogeneity.

3.2.4 Fractal analysis of connected structures

Hausdorff proposed the theory of fractal dimension in 1919, which can be applied to characterize the complexity and irregularity of porous media (Mart́nez-López et al. 2001). Because of the simplicity of the box counting algorithm and its applicability to fractal digital images, the box counting method has been applied to estimate the fractal dimension of the connected pore structure. The basic principle of the box counting method is to divide a 3D space of size M × N × P into cubic boxes of size D × D × D. The number of boxes N that span the edges of the geometric structure is counted, which means the number of boxes whose values inside the boxes contain both 0 and 1. For geometric structures with fractal dimension, the relationship between N and D can be expressed as:

Therefore, the fractal dimension F can be calculated as follows:

With a program in Python, the relationship between logN(D) and log(D) are calculated and plotted, as shown in Fig. 11. The absolute value of the slope is the fractal dimension of the connected pore structure.

Fractal dimension of 3D connected pore structure

As shown in Fig. 11, the fractal dimensions of the connected pore structures of the six digital cores are 2.0971, 2.1244, 2.0558, 2.3952, 2.3914 and 2.4940, respectively. Overall, the fractal dimension of coal sample 1 is lower than that of coal sample 2. It indicates that the pore structure of coal sample 2 is more complex and the pore surfaces are rougher, which can be verified from Fig. 8. The pore structure of coal sample 1 is relatively less developed, and there are a large number of areas in the structure without pore distribution. The pore structure inside coal sample 2 performs well-developed, distributed throughout the porous media structure space, and the connection of pores is also more complicated.

3.2.5 Pore tortuosity

The geometric tortuosity of the pore structure can measure the degree to which the flow path of the fluid in the porous media deviates from the straight line, which is defined as the ratio of the effective length L_h of the flow path in the macroscopic flow direction to the linear distance L of the macroscopic flow direction (Jarrar et al. 2021), expressed as:

The transport flow path in porous media is formed by the pore network, so the shortest pore channel L_s is usually used as the effective flow path L_h to calculate the geometric tortuosity, and the Eq. (8) can be rearranged as:

The pore center method was applied to calculate the geometric tortuosity of the 3D connected pore structures in three different directions of x, y and z. The basic principle of the pore center method is that the geometric tortuosity can be estimated by calculating the average change in the position of the pore centroid among adjacent 2D slices of the 3D microscopic pore structure, as shown in Fig. 12.

Schematic diagram of the geometrical tortuosity of the centroid method (Blue indicates the pore, black indicates the solid)

The centroid of each slice is C(i), and in three-dimensional space, the centroid coordinate is (x_i, y_i, z_i). The distance li between adjacent slice i and slice i+1 can be calculated by the Euclidean distance formula (Faizah et al.  2020), and thus the geometric tortuosity is calculated as:

where: i represents the ith slice; n represents the total number of slices. The calculated geometric tortuosity of the pore structure of different cores in the x, y and z directions is shown in Table 4. The results show that the pore tortuosity in the three directions of each 3D pore structure is significantly different. Comprehensively, the tortuosity of the pore structure of coal sample 1 is higher than that of coal sample 2, which will cause the permeability of coal sample 1 to be lower than that of coal sample 2. Furthermore, digital cores S2-1-12 μm and S2-2-12 μm have a relatively high tortuosity in the z direction. According to Fig. 10, the surface porosity distribution interval in the z direction of the two digital cores is wide, and the minimal surface porosity is close to zero, which will lead to very narrow pore channels in a certain layer, forming a 'funnel'-shaped pore structure. Consequently, this increases the effective transport path length and results in greater pore tortuosity.

4.1 Pore network model generation

The pore network model (PNM), a simplified model, can clearly and intuitively reflect the complex and irregular pore network structure inside the porous media. PNM is not only an important approach to investigate the fluid migration characteristics in microscopic pores and fractures, but  it also establishes a bridge connecting microstructural characteristics and macroscopic seepage behaviors (Lv et al. 2020; Lanetc et al. 2022). PNM is a major tool to investigate the pore structure characteristics inside porous media. The study divides the pore space of complex porous media into units, with narrow pore spaces defined as pore throats and connections among pore throats defined as equivalent pores. The pores and throats are given corresponding characteristic parameters, thus forming a topology based on the original pore structure (TIMUR A 1969; Zheng et al. 2018). In this model, a ball represents the equivalent pore, with its radius corresponding to the equivalent radius of the pore; a circular tube represents the pore throat, with its radius and length corresponding to the width and length of the pore channel, as shown in Fig. 13.

Schematic of pore network model extraction

As shown in Fig. 8, the 3D pore structures of the two coal samples selected in this work were established for each of them. Their equivalent pore network models could be extracted from the 3D pore structure by the central axis method, as shown in Fig. 14. It can be seen from Fig. 14 that compared with coal sample 1, the number of equivalent pores and throats of coal sample 2 are obviously more. This quantitative difference shows that the pore space structure inside coal sample 2 is more developed and the connectivity of the pore structure of coal sample 2 is higher.

Equivalent pore network model for the 3D connected pore structure

4.2 PNM-based quantitative analysis of pore structure

As shown in Fig. 14, the established pore network model with balls and circular tubes exhibits a large amount of quantitative data about the equivalent pore structure, and the results obtained by counting and calculating the pore parameters (pore radius), throat parameters (throat length, throat radius), and coordination numbers in the PNM model are shown in Fig. 15 and listed in Table 5. It can be seen that the number of equivalent pores and throats in the pore network of Sample 2 is significantly larger than that of Sample 1. For example, the number of pores for PNMs S1-1-1 μm, S1-1-12 μm, and S1-2-2 μm are 634, 918, and 970, respectively. The number of pores for PNMs S2-1-1 μm, S2-1-12 μm, and S2-2-2 μm are 1595, 1591, and 2875, respectively. The equivalent pore size distribution of the two samples is similar. The pore size of the 1 μm PNMs is mainly between 2–7.5 μm, and most of the pore size of the 12 μm PNMs is between 20–70 μm.

Statistical distribution of PNM structure

In addition, the throat in the pore network model acts as an important mass transfer channel between equivalent pores, and the number of throats connected with the individual pores is defined as the coordination number (An et al. 2016). The larger the coordination number is, the more throats are connected to the pores, and the more transport channels it has with other pores, corresponding to the higher pore connectivity. The coordination number in PNM of Sample 1 model is 2.33, 2.82, and 2.63, while the number of Sample 2 is 2.69, 3.58, and 4.74. The above data indicate that the pore system of coal sample 2 has a better pore connectivity compared with coal sample 1. This finding is in high agreement with the analytical results of porosity and connectivity of the extracted 3D pore structure, as shown in Fig. 9.

4.3 PNM-based permeability estimation

4.3.1 Permeability estimation based on the PNM simulation

According to the PNM, it is enabled to analyze the transport properties inside the pore structure, such as permeability. Based on Poiseuille's law (Li and Feng 2021), it can be calculated from the pore throat length l_ij, the throat radius r_ij, and the dynamic viscosity μ to obtain the conduction coefficient, which is expressed as:

Due to the existence of pressure gradient in the pore network, it is known from the law of laminar flow that the fluid will flow driven by the pressure gradient, and the flow velocity q_ij is linearly related to the pressure gradient (P_i–P_j), the expression of which is as follows:

For each node in the pore network system, the continuity equation can be established based on the law of conservation of mass:

where: the summation symbol represents all pores ~ j which are connected to pore ~ i. Through applying Eq. (12) to all nodes of the whole pore network, the matrix expression could be given as:

where: G represents the matrix of conduction coefficients with a size of N_p × N_p; N_p represents the number of pore nodes. P is the pressure vector of each pore node; S represents the source vector constrained by the boundary conditions of inlet and outlet pressure. G_ij, an element of matrix ~ G, could be expressed as:

By adding up the flow contribution q_ij through each pore throat, the total flow flux (Q) can be obtained.

After that, it is able to calculate the gas permeability of the pore network model in each of the three directions based on the relationship between the permeability and the flow flux and pressure gradients (as shown in Eq. (16)), and the values are listed in Table 6. The parameters used in the calculation are set as follows: methane is selected as the fluid, the dynamic viscosity coefficient is set to 1.08 × 10⁻⁵ Pa∙s, and the inlet and outlet pressures are selected as 2 bar and 1 bar, respectively.

From the data in Table 6, it is clear that the absolute permeability of PNM at 1 μm resolution is significantly smaller than that at 12 μm, which is due to that the 12 μm digital core contains much larger pores and larger throat radius, which can be found in Fig. 15. The larger pores contribute more significantly to permeability. In addition, Sample 2 presents a higher permeability compared to Sample 1, attributed to its more developed pore structure. And its fractal dimension is significantly larger than that of Sample 1 (shown in Fig. 11), providing more flow channels in the coal. Moreover, the overall pore tortuosity of Sample 2 is smaller (shown in Table 4), leading to Less flow resistance in the pore system and higher gas deliverability. Therefore, it can be seen that the permeability characteristics obtained by PNM simulations are in good agreement with the pore structure characteristics.

4.3.2 Linking the gas deliverability to the pore structure in coal

The permeabilities in the three-dimensional directions show the strong anisotropy, especially for the PNMs at 12 μm resolution, and the permeability difference in different directions can be up to several hundred times, as shown in Fig. 16. The distribution of surface porosity in the 3D directions interprets the anisotropic permeability of the digital structure. Large surface porosity in local structure does not necessarily correspond to high permeability, and the existence of local structure with lower surface porosity greatly reduces the gas transportability of the whole pore structure. Besides, the pore tortuosity in the z-direction of S2-1-12 μm core is 4.4 (Table 4), which is significantly greater than that of other pore structures, indicating that the pore channel route in this direction is very long. This elongation is not conducive to efficient mass transfer, resulting in very low permeability.

Anisotropic permeability and 3D surface porosity

From the calculated permeability values, it can be found that in the x-direction and y-direction the permeability is significantly greater than that in the z-direction. This indicates that the coal seams have greater permeability in the horizontal direction compared to the permeability in the vertical direction. From the data of surface porosity, it can be seen that the minimum surface porosities of pore structures S2-1-12 μm and S2-2-12 μm in the z-direction are nearly zero, indicating that there is almost no pore space in some cross-sections of the pore structures (as shown in Fig. 17). Consequently, gas transport will be hindered by the partial areas with low surface porosity, leading to a relatively low permeability in the z-direction.

Surface porosity of the PNM in the z-direction ~ S2-1-12 μm and S2-2-12 μm

These values show that two anthracite pore structures studied have significant anisotropic characteristics. The pore structure in the horizontal direction is the dominant transportation channel, which can offer a theoretical basis for the design of borehole drilling for CBM extraction and reservoir permeability enhancement.

In this work, MIP and μ-CT image-based techniques were applied to probe the multi-scale pore structure of two anthracite coals from Qinshui Basin. 3D pore structures were reconstructed with processing the sequential CT images, and the structural parameters including pore volume were quantitatively analyzed. Permeability and anisotropy were estimated by extracting PNM from the 3D pore structure. Based on these findings, the following conclusions can be drawn:

1. (1)

   The visualization of the 3D spatial distribution of coal pores, coal matrix, and minerals was realized by reconstructing 3D digital coal cores, enabling quantitative characterization of coal's structural parameters. The results show that the porosity ranges from 7.04% to 8.47% for Sample 1 and from 10.88% to 12.11% for Sample 2. The connectivity ranges from 0.5422 to 0.6852 for Sample 1 and from 0.7948 to 0.9186 for Sample 2.

2. (2)

   The surface porosity of the three-dimensional structures in x, y, and z directions can be used to characterize the anisotropy of pore structure development in different directions. The results show that the surface porosity of the studied samples varies significantly in the z direction, and the minimum surface porosity is almost 0, indicating that there are some poorly developed pores in the z direction, which may hinder gas transport in local regions.

3. (3)

   The results show that the fractal dimensions of the connected pore structures of the six digital cores are 2.0971, 2.1244, 2.0558, 2.3952, 2.3914, and 2.474, respectively, indicating that coal sample 2 has a more complex pore structure and rougher pore surface compared to coal sample 1. The pore tortuosity in all three directions of each 3D pore structure is significantly different. The overall tortuosity of the pore structure of coal sample 1 is higher than that of coal sample 2, which results in lower permeability for coal sample 1 compared to coal sample 2.

4. (4)

   PNM is a simplification and representation of the original pore structure. The coordination number is another quantitative indicator used to characterize connectivity. The coordination numbers of PNMs in Sample 1 are 2.33, 2.82, and 2.63, while those in Sample 2 are 2.69, 3.58, and 4.74. The estimated permeability of coal structures based on PNM shows that permeability is higher in the horizontal direction, indicating that horizontal direction is the dominant flow direction of the coal seam.

The published datasets also indicate that pore size, connectivity, tortuosity, and anisotropic permeability can be fully captured at the micron-scale resolution. In summary, it is concluded that multiscale imaging is powerful and can realize the three-dimensional, visual, and quantitative characterization of complex pore structures. This provides a basis for studying the transport behavior of CBM in real coal pore structures and for the efficient development of CBM.