Quantitative two/three-dimensional spatial characterization and fluid transport prediction of macro/micropores in Gaomiaozi bentonite

Jingfeng Liu,Shiji M,Hongyng Ni,Hi Pu,Xiozho Li,Shojie Chen

a State Key Laboratory for Geomechanics and Deep Underground Engineering,and School of Mechanics and Civil Engineering,China University of Mining and Technology,Xuzhou,221116,China

b CAEA Innovation Center on Geological Disposal of High Level Radioactive Waste,Beijing Research Institute of Uranium Geology,Beijing,100029,China

c Key Laboratory of Mining Disaster Prevention and Control,Shandong University of Science and Technology,Qingdao,266590,China

Keywords:Gaomiaozi(GMZ)bentonite Pore structure Permeability prediction Focused ion beam-scanning electron microscopy(FIB-SEM)

A B S T R A C T The sealing performance of a bentonite barrier is highly dependent on its seepage characteristics,which are directly related to the characteristics of its pore structure.Based on scanning electron microscopy(SEM)and focused ion beam-SEM(FIB-SEM),the pore structure of bentonite was characterized at different scales.First,a reasonable gray threshold was determined through back analysis,and the image was binarized based on the threshold.In addition,binary images were used to analyze bentonite’s pore structure(porosity and pore size distribution).Furthermore,the effects of different algorithms on the pore structure characterization were evaluated.Then,permeability calculations were performed based on the previous pore structure characteristics and a modified permeability prediction model.For permeability prediction based on the three-dimensional model,the effect of pore tortuosity was also considered.Finally,the accuracy of numerical calculations was verified by conducting macroscopic gas and alcohol permeability experiments.This approach provides a better understanding of the microscale mechanism of gas transport in bentonite and the importance of pore structures at different scales in determining its seepage characteristics.

1.Introduction

How to safely store high-level radioactive waste(HLRW)is a global problem.At present,most countries adopt deep geological disposal schemes.The sealing performance of engineering barriers in the geological disposal of HLRW is critical to its safety assessment.In China,granite geological formations are the most important candidates for the storage of HLRW(Zhao et al.,2013;Chen et al.,2014,2017).Bentonite is usually placed around a disposal tank for buffering,backfilling,and sealing(Chen et al.,2014;Yoon et al.,2021),which is traditionally compacted into blocks in situ and placed between the host rock and waste tank(Fig.1).The compacted bentonite is inevitably exposed to water and gas transport during long-term storage.The generation and migration of gas will affect the saturation and sealing characteristics of the bentonite barrier.Therefore,fluid transport(including water and gas)is a crucial parameter for evaluating the sealing efficiency of bentonite buffer.However,fluid transport in rock and soil materials is directly related to pore structure,such as porosity,pore size distribution(PSD),and connectivity of the pore networks.

Currently,several methods are available to analyze the pore information of geomaterials:(i)nitrogen adsorption(NA),pressure-controlled porosimetry(PCP)method,rate-controlled porosimetry(RCP),and nuclear magnetic resonance(NMR);(ii)scanning electron microscopy(SEM)and transmission electron microscopy(TEM);and(iii)X-ray computer tomography(CT)and focused ion beam(FIB)-SEM.Each of these technical methods has its advantages.For example,the NA method can characterize the pore shape by the shape of the hysteresis loop of the isothermal adsorption curve.However,NA can only identify a limited range of PSD(less than 0.1 μm)(Zhang et al.,2017).Compared with the NA method,the PCP method can overcome this shortcoming and characterize larger pores(Yao et al.,2009).NMR technology is characterized by its fast and nondestructive testing to compensate for these deficiencies.In addition,the T2 spectrum information of the sample can be used to invert the pore structure characteristics.

Fig.1.Preliminary concept map of China’s HLRW repository(modified from Chen et al.,2014).

With the advancement of micro and meso imaging technologies,SEM and TEM can directly observe the pore structure characteristics.Moreover,the resolution of these two technologies is higher,up to 3 nm for SEM and 0.2 nm for TEM(Tang et al.,2008;Song et al.,2015,2020).Also,the pore structure can be characterized by SEM and TEM in three dimensions,with significant limitations.In recent years,many studies have adopted micro-CT(μCT)and FIBSEM to characterize the pore structure of geological materials,e.g.pore shape,PSD,porosity,connectivity,tortuosity,and fluid transport properties(Tomioka et al.,2010;Liu et al.,2017;Saxena et al.,2017;Tang et al.,2019;Karimpouli et al.,2020).However,μCT has a low three-dimensional(3D)spatial resolution,such as the Carl Zeiss Xradia 510 Versa model with a resolution of 55-0.5 μm.For clay-based pore networks,it is challenging to observe connected pore networks at this scale(Robinet et al.,2012;Song et al.,2015).FIB-SEM overcomes the limitations of CT’s low resolution and the deficiencies of SEM and TEM,which only provide two-dimensional(2D)images.The FIB can continuously cut the sample in thez-direction at intervals up to the nanometer scale.

For the microstructure of bentonite,some researchers observed the initial state of the bentonite powder/pellet mixture by μCT and SEM.The observation of bentonite pellets by μCT can be an excellent supplement to mercury intrusion porosimetry(MIP)and SEM observations.For pores with sizes of less than 5.5 nm,a NA method is needed to characterize the pore structure in the future(Guerra et al.,2017).Once in contact with water,bentonite will absorb water and swell,and the process of bentonite swelling and pore filling can be observed by μCT(Guerra et al.,2020).Several researchers analyzed the interface between steel liners and FEBEX bentonite using SEM images after a long-term(18-year)in situ test at the Grimsel Test Site,which gave essential microscale data on the effect of metal corrosion on bentonite(Hadi et al.,2019).However,the quantitative characterization of bentonite pore structure based on digital images still needs further research.The pore structure of bentonite has a direct effect on its seepage properties(sealing ability),and the permeability of bentonite may vary by up to 6 orders of magnitude(from 10-15m2to 10-21m2)depending on the water content and stress conditions(Liu et al.,2020a).What is the microscopic mechanism behind these permeability changes?Can the macroscopic permeability be inverted through a micro digital core test?Is it possible to establish a bridge between the pore structure and the permeability?This research was conducted to solve the problems mentioned above.

In this study,we focus on the influence of different image magnifications on pore structure characterization and permeability prediction and the effect of different algorithms on the characterization of PSD.Meanwhile,“penetration contribution rate”is introduced and applied to the permeability prediction model combined with data about the PSD.Finally,the results of 2D/3D digital model characterization and sample-scale experimental tests are compared and analyzed to evaluate the compatibility of micro digital and macro experimental tests.

2.Sample information and image processing methods

2.1.Sample information

Gaomiaozi(GMZ)bentonite,obtained from Xinghe County,Inner Mongolia,China,was used in this study.An appropriate amount of bentonite powder was placed in a mold and compacted into a cylindrical sample.Table 1 illustrates the mineral composition of GMZ bentonite powder.The sample’s initial water content,dry density,and porosity are 10.28%,1.7 g/cm3,and 30.42%,respectively.The samples were dried in an oven at 105°C before microscopic observation.One sample was put in the triaxial cell to perform a gas permeability test with the confining pressure of 1 MPa.The measured permeability was 1×10-15m2.Meanwhile,another sample was processed into small pieces with a side length of less than 10 mm following SEM(FEI Quanta TM 250)and FIBSEM(FEI Helios Nanolab 600i)test requirements.More information about the materials can be found in the study reported by Liu et al.(2020a,b).

Table 1The mineral composition of GMZ bentonite powder(Liu et al.,2018).

2.2.Image processing

2.2.1.Image binarization

In this study,the magnifications are 1000×(image resolution:0.29 μm),2000×(image resolution:0.14 μm),and 5000×(image resolution:0.05 μm),as shown in Fig.2.Segmentation(binarization)of the images is necessary to quantify the pore structure.Binarization is mostly composed of two stages.First,the SEM image is converted into a gray image,and the result is analyzed to obtain a grayscale value matrix.Then,an appropriate threshold was selected based on the obtained gray matrix.The gray threshold(T)can be defined as follows:

wheref(x,y)is the gray value of a pixel point(x,y)in the image,andTis the gray threshold.The matrix value greater than the threshold is set to 1,representing the white part of the image,which is the soil skeleton.The matrix value smaller than this threshold is set to 0,indicating a black image part,i.e.the void portion.

2.2.2.Gray threshold determination

For the gray thresholdT,at the same resolution,when the threshold is lower,the pores are not completely extracted,and the calculated porosity is lower(and vice versa).There are 16 algorithms for determining the threshold(Song et al.,2015),such as Huang algorithm based on fuzzy set theory(Huang and Wang,1995),Otsu(2007)’s threshold clustering algorithm,MaxEntropy algorithm based on the maximum entropy of gray distribution histogram(Kapur et al.,1985),Yen algorithm(Yen et al.,1995)based on the difference between the threshold and original images,and the number of bits required to represent the threshold image.More information can be obtained from the studies of Song et al.(2015,2020)and Liu et al.(2020c).

Fig.2.SEM image of GMZ bentonite with magnifications of 1000×,2000×,and 5000×.

Fig.3.Schematic diagram of DPSD and CPSD algorithms(modified from Münch and Holzer,2008).

Fig.4.Capillary bundle model of porous medium and a modified permeability prediction model.

2.3.Quantitative analysis of pore structure

2.3.1.Surface porosity

Surface porosity φ is defined as the ratio of pore area in an image divided by the total area of the image.Because the SEM image is 2D,the calculated porosity is actually“surface porosity”,and it can be expressed as follows:

whereSφandSare the pore area and the total area of an image,respectively.

2.3.2.PSD

PSD is an essential parameter for the quantitative characterization of pore structure.In this study,the discrete PSD(DPSD)and continuous PSD(CPSD)algorithms were used to analyze the distribution of pores in the image(Münch and Holzer,2008).

(1)DPSD

The principle of DPSD algorithm is to convert an irregular pore structure extracted from the SEM image into a regular circular pore structure with the same area size.As shown in Fig.3,an irregular pore corresponds to a circle of radiusrD.Thus,the PSD of the sample is equivalent toncircles with different radii.

(2)CPSD

Unlike DPSD,CPSD treats the entire pore structure as a continuum.The core idea is to fill the whole pore structure with circles of different sizes,count the circles with different radii,and convert them into the corresponding PSD.As shown in Fig.3,an irregular pore structure can be filled withncircles with different radii,where the maximum radius isbeyond which the circle will exceed the pore boundary.The remaining pore portions are sequentially filled with circles with radii of

2.4.Permeability prediction model

2.4.1.A modified permeability model based on digital image

The Hagen-Poiseuille(H-P)and Darcy’s laws were combined to establish a predictive model for the permeability of porous media(Koz?owski et al.,2011).The H-P equation describes the relationship between fluid flow rate,pressure difference,cylindrical tube size(radius and length),and hydrodynamic viscosity when a fluid is transported through a long cylindrical tube with a constant crosssectional area.It can be expressed as follows:

whereQiis the fluid flow rate,Riis the radius of tubei,ΔPis the pressure loss,μ is the viscosity coefficient of fluid,andLiis the length of tubei.

In Section 2.3.2,the pore network is equivalent to a circle of unequal radius using CPSD and DPSD methods.Therefore,the total flow of fluid flowing through the sample can be equal to the sum of the flows through each tiny pore(Fig.4):

whereQis the total flow in the capillary,andmis the number of capillaries.At the same time,the flow of fluid also follows Darcy’s law,and it can be expressed as follows:

wherekis the permeability,His the sample length in the direction of pressure drop,andAis the cross-sectional area of the sample.

Combining Eqs.(3)-(5),we can obtain

With further simplification,we have

where τi=Li/His the tortuosity.For 2D SEM images,it was assumed that all the capillaries in the sample are straight,i.e.τi=1.Therefore,Eq.(7)can be simplified as follows:

In Eq.(8),we consider that the contribution of different pore sizes to the permeability is the same.The degree of difficulty for fluid migration in various sizes of pores is different.Therefore,we introduce a parameter of permeability contribution rate∈i，as shown in Fig.4:

whereSiis the area of pores with the radiusRi.Therefore,Eq.(8)can be rewritten as below:

2.4.2.Kozeny-Carman and Katz-Thompson equations

In addition to the permeability calculation model based on the H-P and Darcy’s laws,the permeability(kKC)calculation model based on the Kozeny-Carman(K-C)equation(Torquato,1991;Zheng and Li,2015)is expressed as follows:

wheredis the pore diameter,and τ is the tortuosity.

Katz and Thompson(1986)reported that fluid transport in porous media has a close relationship with pore throat,porosity,and tortuosity.The Katz-Thompson(K-T)equation is similar to the K-C equation;however,their basic principles are different.The-T equation describes the fluid transport in porous media with a wide pore throat distribution based on the percolation theory.

The seepage characteristics are controlled by a group of pores larger than the breakthrough pore diameter,corresponding to the peak value of the derivative of the mercury injection curve.The value ofdin the K-C equation mainly depends on the contact degree of particles and pore volume.Recently,the value ofdhas been primarily obtained from the PSD curve:

3.Experimental results and analysis

3.1.Quantitative characterization of pore structure(2D)

3.1.1.Determination of a suitable algorithm and calculation of porosity

The gray threshold of an image is calculated using different segmentation algorithms.Moreover,the porosity was calculated with the threshold.The calculated porosity ranges from 25% to 93.25%(Fig.5).The calculation results of different algorithms varied significantly.According to the MIP and NMR results,the sample porosity ranged from 30% to 35%.Through porosity inversion and visual observation,the calculated porosity obtained by the Yen algorithm was closer to the MIP and NMR results.Therefore,the Yen algorithm was adopted to perform image binarization,as shown in Fig.6.In addition to the threshold algorithm,the image’s brightness and contrast will also impact the gray distribution curve,but this effect is limited relative to the threshold algorithm.

The porosity of different magnifications is different even with the same algorithm(Fig.6).For the three magnifications of 1000×,2000×,and 5000×,the porosity calculated using the Yen algorithm is 25.09%,26.97%,and 25.92%,respectively,with a mean square error of 0.94,indicating that the dispersion is small.Besides,the original and binary images show that the pores are relatively developed,mainly micron-sized pores,and there are a few large pores around the matrix.

Fig.5.Porosity of SEM images at different magnifications obtained by different algorithms.

3.1.2.PSD

According to the binary map shown in Section 3.1.1,the pores in the SEM image were extracted.Further,DPSD and CPSD were used to calculate PSD,as shown in Fig.7a and b.The difference in the morphology of DPSD and CPSD algorithms is slightly significant.The PSD obtained by the DPSD algorithm ranges from 0.01 μm to 10 μm,and most pores are concentrated in 0.01-0.1 μm.The CPSD algorithm achieved a PSD within 0.05-2.4 μm,and most of the pore size is concentrated in 0.14-0.4 μm.According to our previous study,the PCP results show that most pore sizes are concentrated in 0.1-0.4 μm,closer to the CPSD results.The NMR test results are relatively close to the DPSD results.Also,the observed minimum pore is related to magnification.The smallest pores extracted from SEM images at three magnifications(1000×,2000×,and 500×)are 0.1 μm,0.05 μm,and 0.01 μm,and those obtained by the CPSD method are 0.29 μm,0.14 μm and 0.05 μm,respectively.According to Liu et al.(2020a),the MIP(PCP)method can identify nanoscale pores(4 nm),and some pores cannot be observed by the SEM technique.The PSD obtained based on the digital image mainly depends on the image’s resolution and related algorithms.

Loucks et al.(2012)defined the PSD criteria,ranging from 1 nm to 1 μm for nanopores,1-62.51 μm for micropores,and 62.51 μm to 4 mm for mesopores.According to the standard,the DPSD results show that the proportion of nanopores is 32.63%(1000×),29.23%(2000×),and 50.03%(5000×),and the remaining pores are micron-sized.The results at magnifications of 1000×and 2000×have a minor difference but are quite different from those at magnification of 5000×because the pores within the range of 0.01-0.05 μm cannot be observed at magnifications of 1000×and 2000×,and the proportion of pores in this range is about 16%.The CPSD results show that the proportion of nanopores is 89.24%(1000×),92.87%(2000×),and 100%(5000×),quite different from the DPSD results.The MIP tests show that the nanopores proportion is about 72%,which is closer to the CPSD results.

Fig.6.Binary images of different magnifications based on the Yen algorithm.

Fig.7.PSD distribution based on(a)DPSD and(b)CPSD algorithms.

3.2.Permeability prediction(2D)

After quantitatively characterizing the pore characteristics of the image,a permeability prediction model(i.e.Eq.(8))was established according to the H-P and Darcy’s laws.To improve the credibility,five different areas(upper left,upper right,middle,lower left,and lower right)were selected for processing and calculation for each SEM image(Fig.8).Then,each area(defined as the frame size)was expanded by 50×50,100×100,…,800×800,850×850 pixels until the image boundary was reached.Similarly,the permeability calculation was also based on DPSD and CPSD algorithms.

3.2.1.Permeability prediction based on the DPSD algorithm

The results of the permeability calculation for SEM images at three different magnifications are shown in Fig.9a-c.It can be observed that the results are very discrete,and the curves fluctuate considerably for smaller boxes.When the box length is 50-100 pixels,the difference in permeability results in different regions is up to three orders of magnitude.Besides,a large gray value is observed in local areas of the image(e.g.upper right corner of the image at a magnification of 1000×).Thus,when the image is binarized,more pores are extracted,resulting in a higher permeability in this region.The difference between the images with different sizes can be attributed to the heterogeneity of the pore distribution of material.For bentonite,the pore size varies broadly.Therefore,the calculated permeability is very large when the box covers only the pore accumulation area.

When the box size increases gradually,the calculation results of permeability in different regions tend to converge.Overall,when the box size exceeds 850 pixels(1000×),850 pixels(2000×),and 800 pixels(5000×),the calculated permeability is closer to the full-scale image calculations,namely,3.56×10-13m2(1000×),1.93×10-13m2(2000×),and 1.26×10-13m2(5000×).Therefore,the representative elementary areas(REA)calculated according to the box overlay method are 0.2095 mm2,0.2095 mm2,and 0.1469 mm2,respectively,while the image size is 0.2613 mm2.

Fig.8.Five different locations for permeability calculation.

Fig.9.Permeability prediction results based on the DPSD algorithm at magnifications of(a)1000×,(b)2000×,and(c)5000×.

Fig.10.Permeability prediction results based on the CPSD algorithm at magnifications of(a)1000×,(b)2000×,and(c)5000×.

3.2.2.Permeability prediction based on the CPSD algorithm

As shown in Fig.10,the permeability calculated by the CPSD algorithm has a similar trend to that of the DPSD algorithm.Similar to the trend of DPSD results,when the box size is small,the permeability results fluctuate significantly.Moreover,the higher the resolution,the greater the dispersion of results.Therefore,when the cover size of the box is small,the heterogeneity of the image has a significant effect on the calculation result.

With increasing box size,the predicted permeability results tend to stabilize.The REA values calculated according to the box overlay method are 0.1044 mm2(1000×),0.1044 mm2(2000×),and 0.1469 mm2(5000×),respectively.The permeabilities of the full-size binary image at three different magnifications are 1.12×10-14m2(1000×),4.92×10-15m2(2000×),and 1.29×10-15m2(5000×),respectively.

3.2.3.Comparison of the permeability of images at different magnifications

As shown in Fig.11,the results of full-scale image calculations at different magnifications are still different.For the DPSD algorithm,the permeability calculation results at the three magnifications are 6.64×10-12m2(1000×),4.25×10-13m2(2000×),and 2.26×10-13m2(5000×).The difference in permeability calculations for different magnifications can be nearly ten times.For the results based on the CPSD algorithm,the calculation results at the three magnifications are 3.63×10-15m2(1000×),9.96×10-16m2(2000×),and 1.48×10-16m2(5000×),respectively.The calculated values at magnifications of 1000×and 5000×differ byabout 29.38(DPSD)and 24.53(CPSD)times,respectively.

The difference in the calculation results of permeability at different magnifications is mainly due to the heterogeneity of the pore distribution of bentonite.Different magnifications cover other areas of the sample,resulting in different calculation results for permeability.According to Fig.6a and b,at the magnifications of 1000×,2000×,and 5000×,the maximum pores observed are 10.17 μm,6.83 μm and 3.72 μm,while those obtained by the CPSD method are 2.32 μm,1.89 μm and 0.98 μm,respectively.According to Eq.(8),the permeability prediction results mainly depend on the large pore size and distribution frequency.

3.2.4.Comparison of permeability between DPSD and CPSD algorithms and its experimental verification

Fig.11.Comparison of the permeability of images at different magnifications obtained by(a)DPSD and(b)CPSD algorithms.

As shown in Fig.12,the permeability calculated by the DPSD method is substantially 2-3 orders of magnitude larger than that calculated by CPSD because the maximum aperture predicted by the CPSD algorithm is much lower than that of the DPSD method.Moreover,the number of large pores extracted by the DPSD algorithm is also relatively large.

Fig.12.Comparison of calculation(Cal.)results based on DPSD and CPSD methods(with different magnifications)with experimental(exp)results(gas permeability).

We selected a sample for the gas permeability test to verify the accuracy of the predicted result,and the calculated permeability is 10-15m2.The results of the CPSD algorithm are close to the observed values,especially at higher magnifications.For the image with a magnification of 2000×,the prediction results and laboratory test values are 9.96×10-16m2and 1×10-15m2,respectively,which are very close.The image resolution at the magnification of 1000×is relatively low,resulting in parts of the pores that cannot be observed.The resolution of the image at the magnification of 5000×is higher,and the observation area is too small,which causes the calculation result to be significantly affected by the selected region.

Moreover,the permeability calculated by the CPSD algorithm is slightly lower than that obtained in the laboratory.In Eq.(8),it is assumed that all pores are involved in fluid migration.However,only the interconnected pores participate in fluid migration.The sophisticated pore network is simplified into many separate connected tubes in the H-P and Darcy’s law prediction models.The interconnection between pore networks is intricate.The connectivity and tortuosity of the pore network not only exist in the 2D plane but also in the 3D space.Because the SEM image is 2D,the connectivity and tortuosity of the 3D pore network cannot be considered.Therefore,the predicted value is slightly oversized.However,compared with the conventional H-P law,the improved permeability prediction model can evaluate the contribution of pores of different sizes to the seepage characteristics,and the prediction results are relatively accurate.

Fig.13.(a)3D reconstruction model of GMZ bentonite sample;and(b)Comparison of calculation results of different algorithms.

Fig.14.3D binary images of pore structure:(a)Matrix(gray)+pore(blue+red);and(b)Connected pore(blue)+isolated pore(red).

Fig.15.(a)Connected and(b)isolated pores.

3.3.Quantitative characterization of pore structure(3D)

3.3.1.3D reconstruction and porosity calculation

Further,we attempted to characterize the 3D pore structure of bentonite using FIB-SEM,as shown in Fig.13a.First of all,the influence of several algorithms on the gray threshold is compared.The porosity calculated with the Shanbhag algorithm(29.73%)was close to the test results with the MIP and NMR(25%-35%),while those obtained by the Yen algorithm were better for the SEM image,as shown in Fig.13b.Therefore,the Shanbhag algorithm was used for FIB-SEM image binarization.The porosity(29.73%)obtained based on the Shanbhag algorithm was calculated using ImageJ.According to this algorithm,each slice has a gray threshold.However,in the 3D reconstruction,a unified threshold was used,and it was impossible to apply a threshold to each image individually for binarization.

Fig.14a shows a 3D binary image of the pore structure containing both matrix and pores of GMZ bentonite based on the Shanbhag algorithm.Extraction of the pore network from the matrix allowed for further evaluation of its percolation characteristics(Fig.14b).It was found that the pores of bentonite are relatively developed,verifying that the value of permeability measured in the laboratory test is relatively high(with an order of magnitude of 10-15m2).The samples used for the FIB-SEM and gas permeation tests are all in dry state.

The seepage characteristics mainly originate from the connected pores,while the isolated pores do not contribute to the seepage characteristics.There are typically 6-,18-,and 26-neighborhood connectivities regarding the pore connectivity.In this study,26-neighborhood connectivities were used since 6-and 18-neighborhood connectivities resulted in smaller interconnected pores that differed significantly from the test results.Fig.15a and b shows the distribution of interconnected and isolated pores.The pore connectivity of bentonite is relatively good.The connectivity(Lc)was 78%,and the effective porosity(φe-3D)was 0.21.These two indicators can be calculated using the following formulas:

wherenc-pis the pixel number of the connected pore,andntis the total pixel number of pore and matrix.

wherent-pis the pixel number of the pore.

3.3.2.PSD

The image obtained by the FIB-SEM scan has an actual size of 20.2541 nm per pixel in thex-andy-directions(pxandpy),and the interval between every two layers is 20 nm in thez-direction(pz).Because the electron gun is not perpendicular to the imaging plane(at an angle of 52°),the actual size represented by each pixel in they-direction should be divided by sin52°to reflect the original size.Because the image resolution obtained by the FIB-SEM scan is inconsistent in thex-,y-,andz-directions,it is inconvenient for subsequent image size calculation.The voxel points were transformed in a volume-equivalent manner to facilitate data processing.The converted resolution(

p)is the cube root of the original voxel volume:

The calculated PSD is shown in Fig.16.It can be found from the figure that the densest region of PSD is located near the image resolution(21.84 nm).This means that image resolution has a significant effect on PSD characterization.The proportion of pores of around 21.84 nm is about 42%,while the proportion of pores larger than 400 nm(0.4 μm)is small.The predicted pore radii for the SEM image with the magnification of 5000×(continuous algorithm)are 0.02 μm(with a proporation of 42%)and 0.4 μm(3%),and the results based on FIB-SEM are 0.05 μm(54%)and 0.4 μm(0.2%).Again,it is shown that nanopores occupy a significant part of GMZ bentonite.

3.3.3.Tortuosity

Path tortuosity(τ0-n)is an essential indicator for evaluating the connectivity of pore networks and a key parameter for predicting permeability.The path curvature is defined as the ratio of path length to the straight line distance between the two ends(H0-n)along thez-direction(L0-n)(Keller et al.,2011)(Fig.17):

In our model,the path length along thez-direction(L0-n)is

The calculated tortuosity is 2.794,indicating that the pore network is quite tortuous.In the mathematical sense,path tortuosity is a geometrical parameter.This value(2.794)is similar to that used to predict permeability in the preceding sections.

Fig.16.PSD of GMZ bentonite based on FIB-SEM(continuous algorithm).

Fig.17.Schematic diagram for calculation of tortuosity.

3.4.Permeability prediction(3D)

The H-P equation is the most straightforward and most practical way to calculate permeability.If the H-P equation is used to calculate a 3D image,we first need to calculate the PSD of the 3D image.Then,each of the pores is equivalent to a corresponding circular tube.This means that the 3D pore network is ultimately equal to the 2D pore network.Therefore,we adopt the K-C and KT models to calculate permeability.

The formulas for the calculation of permeability are shown in Eqs.(9)and(10).In the formulas,the porosity(φ)is the connected porosity(φe-3D),not the total porosity.The tortuosity(τ)was calculated in the previous section.The pore diameter(d)refers to the breakthrough diameter,and it is usually assumed as the peak value of the PSD curve.

Fig.18.Comparison of calculated permeability based on different models and experimental values(gas/alcohol permeability).

Fig.18 shows the calculated permeability for the two models,i.e.0.2×10-18m2for the K-T model and 0.14×10-17m2for the K-C model.The results are much smaller than the tested gas permeability(1×10-15m2).Song et al.(2015)also indicated that the results obtained with the K-T model are lower(two orders of magnitude)than the measured gas permeability.Since the pore distribution observed by FIB-SEM ranges from 0.02 to 1.02 μm,the pores outside this range cannot be characterized.The PSD range obtained based on MIP is 0.003-53.667 μm,thus the permeability prediction based on the pore information obtained by FIB-SEM ignores some pores smaller than 0.02 μm and larger than 1.02 μm(Liu et al.,2020a).The permeabilities calculated based on the MIP results in Liu et al.(2020a)and the K-C and K-T models were 1.91×10-16m2and 1.35×10-15m2,respectively.The comparative analysis shows that the calculation result based on MIP is more closer to the measured gas permeability than the calculated result based on FIB-SEM.

Further,we tested GMZ samples dried with ethanol.Previous experiments showed that if water is used as the fluid medium,the water molecules will react with bentonite particles,affecting the pore structure.The ethanol molecules do not react with the bentonite particles;therefore,the test results are more accurate.The test result(ethanol permeability)is 4.6×10-17m2,close to the calculation result of FIB-SEM based on the K-C model(0.14×10-17m2).However,the difference between the gas and alcohol permeabilities reaches two orders of magnitude.Therefore,the influence of the fluid medium on the permeability test results is more significant and is mainly related to the viscosity coefficient of the fluid and the size of the liquid molecules.

4.Conclusions

This study mainly discussed the characterization of the pore structure of bentonite based on digital image technology and the relationship between pore structure and seepage characteristics.In this process,the effects of different algorithms on digital image binarization were compared.The results showed that the threshold and binary image obtained based on the Yen algorithm are more consistent with the actual situation.Furthermore,the effect of different algorithms on PSD calculations was evaluated.Both the DPSD and CPSD algorithms have unique advantages within a specific aperture range,and the results of the CPSD method are relatively close to those of PCP.

In terms of permeability prediction,different from the conventional H-P law,the permeability contribution rate is introduced to make the prediction result more accurate.Similarly,the calculated results of permeability based on the CPSD algorithm are close to the measured values.Besides,an REA exists,beyond which the calculation results tend to stabilize(converge).In addition,too low or too high image resolution significantly influences the accuracy of the permeability prediction.Therefore,it is critical to reasonablely choose the image resolution for subsequent quantitative characterization.

Finally,FIB-SEM was used to characterize the pore structure of bentonite.Unlike in 2D digital images,two parameters,i.e.connectivity and tortuosity,were obtained.Tortuosity is very important for subsequent permeability calculations,and it can better reflect the actual situation of pore network.The calculation results of the tortuosity were compared with those based on the K-C and K-T models.It was indicated that the result based on the K-C model is closer to the alcohol permeability and lower than the gas permeability.Furthermore,this method can also be used to analyze other materials such as claystone and shale.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

The authors are grateful for the support of the National Natural Science Foundation of China(Grant Nos.52174133 and 51809263),and China Atomic Energy Authority.