Measuring the Effect of Pack Shape on Gravel’s Pore Characteristics and Permeability Using X-ray Diffraction Computed Tomography

Particle shape is one of the critical parameter factors that affect gravel’s pore structure and permeability. However, few studies have considered its effects on engineering applications due to the difficulty of conducting laboratory tests. To overcome these difficulties, new methods of estimating the gravel pack shape that involve manual work and measuring the surface area of particles and pores based on support vector machine segmentation and the reconstruction of X-ray diffraction computed tomography (CT) images were proposed. Under the same conditions, CT tests were carried out on gravel packs and two other regular-shaped particle packs to investigate the influence of particle shape on the fractal dimension of gravel’s pore–particle interface and the specific surface area of the pore network. Additionally, permeability tests were performed to study the effect of particle shape on gravel’s hydraulic conductivity. The results showed that a gravel pack with a larger aspect ratio and a smaller roundness had a larger specific pore network surface area and a more complex pore structure, leading to lower permeability. This kind of gravel had a more significant length, quantity, and tortuosity of the seepage path when seepage occurred in a two-dimensional seepage field simulation. Therefore, we suggest that the filter materials of hydraulic projects should preferably use blasting gravel with a larger aspect ratio and smaller roundness to achieve better anti-seepage properties. In addition, projects can increase pores’ specific surface area using our method as a control factor in filter construction.


Introduction
Permeability is one of the critical parameters of concern in hydraulic engineering, as seepage is the leading cause of earth rock dam failure, accounting for approximately 25% of such failures worldwide [1]. Gravel is used as the filter material for earth rock-fill dams to ensure seepage stability. According to the specifications for earth rock dams in the United States and China, the particle size distribution and dry density (or porosity) are the main construction control factors in filter design. Due to the influence of field test conditions, the design and construction of water conservancy projects require a formula to accurately predict the permeability of gravel based on the particle size distribution, dry density, and other parameters that can be obtained conveniently from the project site. Particle shape is an essential factor that affects the pore structure and permeability of gravel [2]. The hydraulic conductivity of natural gravel is more than five times that of blasted gravel under the same conditions due to their different particle shapes [3]. However, there are few studies on the influence of blasting gravel shape on its permeability in the context of engineering applications. pack usually includes thousands of particles, and it is challenging to complete th urement of the shape of each particle because of the enormous workload involved fore, it is necessary to estimate the shape of the gravel pack through sampling m ment. The rationality of the sampling method design and whether there is a corr between the shape of gravel bags estimated by sampling and the shape-related pro need to be further discussed In light of the above facts, this work aimed to study the influence of particl on the pore structure and permeability of gravel by combining laboratory tests an lations to provide suggestions for the design and construction of filter materials i conservancy projects and the development of permeability-predicting formulas. fore, X-ray CT scanning and constant head tests were conducted on gravel and tw regular-shaped particles. Meanwhile, a new method was proposed for measuring t cific surface area of the gravel pore network based on the segmentation and reconst of CT images using SVM to ensure that the particle accumulation was unchanged tionally, the relationship between the shape factor and pore characteristics of actua in terms of the fractal dimension, specific surface area, and permeability is discuss results of this study provide suggestions for the design and construction of the ear dams' filter and may help to increase the application of the predicting formula re the shape factor and surface area of gravel in engineering.

Materials
The experimental materials were blasted gravels from the Shuibuya Dam China, which blasts limestone fragments. They were sieved and washed to remove rities and dust before being studied. Glass balls and plastic octahedrons were sele comparison to reduce the workload when artificially measuring the physical dime Since there is no seepage deformation involved in low-head permeability tests, the ent materials did not affect the test results. The sample particle sizes ranged from to 20 mm to diminish the size effect on the test results. The materials are shown in 1.  Figure 2 depicts the experiment flow of this study. First, we manually samp gravel and measured the shape. The proctor compaction test examined the ma maximum and minimum dry density. Then, we designed the particle size distribut porosity of laboratory test samples according to the particle size distribution of th pled gravel and the dry density range. Table 1 shows nine test samples, and the sha the only variable factor. Test samples were prepared as described in Table 1 an used to fill a permeameter in three layers of the same thickness for homogenizati  Figure 2 depicts the experiment flow of this study. First, we manually sampled the gravel and measured the shape. The proctor compaction test examined the materials' maximum and minimum dry density. Then, we designed the particle size distribution and porosity of laboratory test samples according to the particle size distribution of the sampled gravel and the dry density range. Table 1 shows nine test samples, and the shape was the only variable factor. Test samples were prepared as described in Table 1 and were used to fill a permeameter in three layers of the same thickness for homogenization. The permeameter was a plastic cylinder with an inner diameter of 10 cm. The lower water inlet chamber was fitted with a connecting pipe for the pressure measuring tube, and the top was fitted with an overflow pipe. The filling height of the sample was 7 cm, which was more than three times the maximum particle size, to eliminate the effect of size on the results. Next, the CT scanning test was conducted on the samples to obtain their pore structures, and the permeability test was carried out on the samples to study their hydraulic conductivities. We segmented the CT scanning test results to extract the samples' pores for fractal analysis and three-dimensional model reconstruction. Finally, the permeability test was simulated two-dimensionally based on the CT images. permeameter was a plastic cylinder with an inner diameter of 10 cm. The lower water inlet chamber was fitted with a connecting pipe for the pressure measuring tube, and the top was fitted with an overflow pipe. The filling height of the sample was 7 cm, which was more than three times the maximum particle size, to eliminate the effect of size on the results. Next, the CT scanning test was conducted on the samples to obtain their pore structures, and the permeability test was carried out on the samples to study their hydraulic conductivities. We segmented the CT scanning test results to extract the samples' pores for fractal analysis and three-dimensional model reconstruction. Finally, the permeability test was simulated two-dimensionally based on the CT images.  We defined the aspect ratio and roundness as parameters to describe the particle shape. The aspect ratio (α) is the length-to-width ratio, and the roundness (S) is defined as the ratio of the circumference of the equivalent area of the projected area to the projected contour's actual circumference, which is as follows:

=
(1)  We defined the aspect ratio and roundness as parameters to describe the particle shape. The aspect ratio (α) is the length-to-width ratio, and the roundness (S) is defined as the ratio of the circumference of the equivalent area of the projected area to the projected contour's actual circumference, which is as follows: where I is the maximum distance between the projected outline points of a particle and L is the short axis of the equal-area ellipse when the long axis of the particle is I. where A is the projection area of a particle and P is the actual perimeter of the particle's projected contour.
The gravels were sieved into three particle size groups: 2-5 mm, 5-10 mm, and 10-20 mm. A total of 1000 particles in each particle size group were randomly sieved as experimental samples for the CT and permeability tests. Then, 100 out of the 1000 particles in each particle size group were randomly selected and measured. The sampling number was 10% of the total gravel, which conformed to the statistical sampling principle. Therefore, the mean value of the shape parameter of the 100 sample particles was used as the shape parameter of the same particle group.
We designed a fixture to hold and rotate an individual particle. The measurement process was as follows: First, the fixture was fixed onto a gravel particle's longest axis and we assumed that the particle was set at an angle of 0 • . We then rotated the gravel to 60 • and 120 • . We used a laser scanner to capture its outer contour at each of three angles and imported these contours to ImageJ. The geometric dimensions of the outer contours required for Equations (1) and (2) were measured using ImageJ and calculated to obtain an individual particle's aspect ratio and roundness at each of the three angles. The average values of the particle's aspect ratio and roundness at the three angles were taken as α and S and had three-dimensional physical significance.
The design process of the test samples' particle size distributions was as follows: The particle size distribution of 100 particles selected from the 2-5 mm, 5-10 mm, and 10-20 mm size groups was obtained by counting their profiles geometrically and were labeled D I , D II , and D III . Then, D I , D II , and D III were multiplied by the percentage in Table 2 and combined to obtain three new particle size distributions of 2-20 mm, which were labeled D1, D2, and D3, and were used as the test samples' particle size distributions.

Proctor Compacion Test
The proctor compaction test examines gravel's maximum and minimum dry densities. The equipment for this experiment included a compaction cylinder (whose volume was 2103.0 cm 3 ), a hammer (whose mass was 2.5 kg), and a guide cylinder. The falling height of the hammer was 457 mm. The maximum porosity of gravel is the porosity when it is in natural accumulation; the minimum porosity of gravel is the porosity when it is in the densest accumulation. The porosity can be calculated from the dry density as follows [45]: where ϕ is the porosity of the gravel, ρ d is the dry density of the gravel, and ρ G is the density of the gravel.

Computed Tomography Scanning Test
The CT scanning test was performed on particle packs contained in the permeameter, which was a high-spatial-resolution Siemens 40 CT machine from Changjiang River Scientific Research Institute, to observe the internal pore structure. The image reconstruction matrix was 512 × 512, and the minimum spatial resolution was 0.29 mm.

Constant-Head Permeability Test
The constant-head permeability test was performed on particle packs after scanning. The test used boiled purified water to eliminate the influence of bubbles. The sample packs were soaked for more than 2 h to saturate entirely before conducting the permeability test. The test flow direction was from the bottom to the top, and the overflow nozzle controlled the downstream waterhead ( Figure 3). The waterhead of the water tank was gradually raised to provide a 0.05 hydraulic gradient with a waiting period of 10-20 min for each step. The upstream and downstream waterhead values were recorded, and the flow of the overflow pipe was measured over time with a measuring cylinder. The following hydraulic gradient of the permeability test was conducted as long as the flow measured the same values three consecutive times. In each test, the pack had no seepage deformation, and the sedimentation value measured at the top of the sample was maintained at 0.
The CT scanning test was performed on particle packs contained in the permeameter, which was a high-spatial-resolution Siemens 40 CT machine from Changjiang River Scientific Research Institute, to observe the internal pore structure. The image reconstruction matrix was 512 × 512, and the minimum spatial resolution was 0.29 mm.

Constant-Head Permeability Test
The constant-head permeability test was performed on particle packs after scanning. The test used boiled purified water to eliminate the influence of bubbles. The sample packs were soaked for more than 2 h to saturate entirely before conducting the permeability test. The test flow direction was from the bottom to the top, and the overflow nozzle controlled the downstream waterhead ( Figure 3). The waterhead of the water tank was gradually raised to provide a 0.05 hydraulic gradient with a waiting period of 10-20 min for each step. The upstream and downstream waterhead values were recorded, and the flow of the overflow pipe was measured over time with a measuring cylinder. The following hydraulic gradient of the permeability test was conducted as long as the flow measured the same values three consecutive times. In each test, the pack had no seepage deformation, and the sedimentation value measured at the top of the sample was maintained at 0.

The New Method of Surface Area Measurement Based on CT Image Segmentation Using SVM
Before analyzing the pore characteristic of particle packs in CT images, the particles, pores, and containers in the image were segmented separately using SVM. Classifying pixels is the essence of image segmentation. SVM is used for pattern classification and nonlinear regression in multilayer perceptron and principal radial function networks by building a classification hyperplane as the decision surface, which maximizes the isolation edge of positive and negative examples [46]. SVM is the approximate realization of structural risk minimization based on statistics. The learning machine's generalization error rate on test data is bounded by the sum of the training error rate and an item depending on the Vapnik-Chervonenkis dimension [47]. The SVM's value for the former item is zero in a separable mode, and the second item is minimized. Therefore, SVM can provide better generalization performance when it comes to pattern classification, which is unique. The discriminant equation of the SVM model is [47]: where X denotes an eigenvector of an arbitrary instance input and xi is a concrete feature in an eigenvector in which = ( 1 , 2 , … , ). The model is trained with all positive instances of labels for which Y = −1 to pursue the appropriate values for the parameters ω

The New Method of Surface Area Measurement Based on CT Image Segmentation Using SVM
Before analyzing the pore characteristic of particle packs in CT images, the particles, pores, and containers in the image were segmented separately using SVM. Classifying pixels is the essence of image segmentation. SVM is used for pattern classification and nonlinear regression in multilayer perceptron and principal radial function networks by building a classification hyperplane as the decision surface, which maximizes the isolation edge of positive and negative examples [46]. SVM is the approximate realization of structural risk minimization based on statistics. The learning machine's generalization error rate on test data is bounded by the sum of the training error rate and an item depending on the Vapnik-Chervonenkis dimension [47]. The SVM's value for the former item is zero in a separable mode, and the second item is minimized. Therefore, SVM can provide better generalization performance when it comes to pattern classification, which is unique. The discriminant equation of the SVM model is [47]: where X denotes an eigenvector of an arbitrary instance input and x i is a concrete feature in an eigenvector in which X = (x 1 , x 2 , . . . , x m ). The model is trained with all positive instances of labels for which Y = −1 to pursue the appropriate values for the parameters ω and b. Thus, an unknown X i would be classified as a positive case when ω T X i + b > 0, and vice versa. The C-SVC (a type of SVM solution model) is suitable for binary problem judgment. In the sample, set T is the input [48]: where x i ∈ X = R n , y i ∈ Y = {1, -1} (i = 1, 2, . . . , l), and x i is an eigenvector. After selecting a suitable kernel function K(x,x ) and an appropriate parameter C, we constructed and solved the global optimization: Therefore, the optimal solution α* is α* = (α 1 *, . . . , α l *) T . We selected a positive element of vector α* (0 < α j < C) and calculated the threshold b*: Then, we constructed the decision function f(x): We classified the pixel values of different objects in CT images first. A CT image of S1 is given as an example in Figure 4. The pores between particles in the CT image are black, while the particles and the container are white with inhomogeneous saturation levels. The pixel values of particles and pores were set as training samples in the limited area according to the SVM method (Equations (4)-(9)) using a code we compiled in MATLAB. The 8-bit image data picture is a three-dimensional matrix with 484 rows, 484 columns, and 3 pages. and b. Thus, an unknown Xi would be classified as a positive case when + > 0, and vice versa.
The C-SVC (a type of SVM solution model) is suitable for binary problem judgment. In the sample, set T is the input [48]: where xi∈X = R n , yi∈Y = {1, -1} (i = 1, 2, …, l), and xi is an eigenvector. After selecting a suitable kernel function K(x,x ' ) and an appropriate parameter C, we constructed and solved the global optimization: Therefore, the optimal solution α * is α * = (α1 * ,…,αl * ) T . We selected a positive element of vector α * (0 < αj < C) and calculated the threshold b * : Then, we constructed the decision function f(x): We classified the pixel values of different objects in CT images first. A CT image of S1 is given as an example in Figure 4. The pores between particles in the CT image are black, while the particles and the container are white with inhomogeneous saturation levels. The pixel values of particles and pores were set as training samples in the limited area according to the SVM method (Equations (4)-(9)) using a code we compiled in MATLAB. The 8-bit image data picture is a three-dimensional matrix with 484 rows, 484 columns, and 3 pages. The material's interactive medical image control system (MIMICS) was used to reconstruct and measure the 3D model's volume and surface area after segmenting the CT images. The pores and particles were reconstructed separately. The reconstruction The material's interactive medical image control system (MIMICS) was used to reconstruct and measure the 3D model's volume and surface area after segmenting the CT images. The pores and particles were reconstructed separately. The reconstruction method was gray value interpolation [49], which considers the partial volume effect. The advantage of gray value interpolation is that it gives lots of detail and the correct dimensions. The validity of the surface area measurement of the reconstructed pore network model was verified by comparing the porosity of the reconstructed particle pack model with the same test sample using the relative error. The porosity of the reconstructed particle pack model is as follows: where ϕ sim is the porosity of the pore network model, V pore is the volume of the pore network model, and V 0 is the sum of the volume of the pore network model and particle pack model. The relative error of the comparison of the porosity of the reconstructed particle pack model with the same test sample is as follows: where R is the relative error and ϕ Lab is the porosity of the test sample. A smaller R means that the two quantities being compared are closer.

Fractal Analysis Method
The pore-particle interface of porous media has good fractal characteristics, and the fractal dimension can be used to study the pore structure's complexity quantitatively using the box-counting method [50]. Assume that the objects are covered by orthogonal line grids with an increasing lattice constant. The number (N) of those grids (boxes) that contain any part of the structure is calculated by the size of each box (equal to the lattice constant ε) and stored in the data list. The macro increases the box size (ε) in selectable step sizes; for each box size, any boxes that contain at least 1 pixel of the contour line were counted (N). This count depends on the box size ε, box-counting b, length L, and fractal dimension D according to Equation (12) [50]: Thus, for fractal objects, a double-logarithmic plot yields a straight line: where D can be determined as the absolute value of its slope, and the constant c describes the ordinate intercept. Binary image processing is required before calculating the box dimension of the pore-particle interface (D). The binary image is divided by an equivalent grid with a side length ε, where the grid occupied by white pixels is defined as N(ε). Thus, lgN(ε i )/lg(1/ε)→D when ε→0. For a specific decreasing sequence {ε i }, the definition of the fractal dimension approximates its slope. The decreasing sequence {ε i } is usually defined via a dichotomy. Then, the box-counting dimension can be defined as follows: For boxes with different side lengths (ε), the required number of boxes (N(ε)) to cover the pore pixels varies. A gray pixel can be covered only by a box with a specific side length considering the relativity of a binary gray image. For the linear equation, the least squares method was used to fit the data points linearly: where the slope a is dimension D.
We compiled a code in MATLAB to calculate the box dimension of the pore-particle interface in CT images of test samples according to Equations (12)- (15). Each sample had 140 CT images.

Pore Scale Simulation
The water flow in particle packs was simulated in FLUENT. Figure 5 shows a schematic of the computational domain and boundary conditions in the numerical simulations consisting of the CT image. We took four central cross-sectional CT images of each pack for the simulation (Figure 5c). The calculation area was 10 × 7 cm. The CT sections were exported to the DXF format file using MIMICS to build the model and generate regions using Au-toCAD. The pressure inlet was determined using the corresponding value obtained from the permeability test while the outlet was free flow. Both the sidewalls and particles were impermeable.

Pore Scale Simulation
The water flow in particle packs was simulated in FLUENT. Figure 5 shows a schematic of the computational domain and boundary conditions in the numerical simulations consisting of the CT image. We took four central cross-sectional CT images of each pack for the simulation (Figure 5c). The calculation area was 10 × 7 cm. The CT sections were exported to the DXF format file using MIMICS to build the model and generate regions using AutoCAD. The pressure inlet was determined using the corresponding value obtained from the permeability test while the outlet was free flow. Both the sidewalls and particles were impermeable. The control equation of the simulation in the domain area Ω is given below [51]. For velocity v: Ω → R 2 ; and pressure p: Ω → R: where ̇= / , is the gradient operator, ρ is the density of the fluid, μ is the dynamic viscosity of the fluid, g is the body force, n is the unit outward normal on the boundary Γof Ω, ΓD is the part of the boundary that experiences the Dirichlet boundary condition, ̅ is the applied velocity, ΓN is the part of the boundary that experiences the Neumann boundary condition, and ̅ is the applied time. Here, Γ = ΓD ∩ ΓN and ΓD ∩ ΓN = ∅. v0 and p0 are the initial velocities and pressure fields in the fluid in the domain, respectively. The pseudo-stress σ is given as follows: where I is the second-order identity tensor. The code uses an unstructured quadrilateral grid of the finite volume method, which has excellent adaptability. The calculation of the time-dependent physical quantities was transient. The number of grids in the calculation area was more than 50,000; therefore, the pressure field was calculated based on staggered grids using the semi-implicit method for pressure-linked equations (SIMPLE) to solve Equations (16)- (17) in order to reduce the calculation time. The models were saturated, and only the solid and liquid phases were in the calculation area.
The inlet water pressure in the model was consistent with that of the laboratory test. The flow state in the simulation was laminar flow; therefore, the velocity ratio at the outlet of the same model to the hydraulic gradient was invariant under different hydraulic gradients. This was defined as the hydraulic conductivity of the model. The control equation of the simulation in the domain area Ω is given below [51]. For velocity v: Ω → R 2 ; and pressure p: Ω → R: where . v = ∂v/∂t, ∇ is the gradient operator, ρ is the density of the fluid, µ is the dynamic viscosity of the fluid, g is the body force, n is the unit outward normal on the boundary Γ of Ω, Γ D is the part of the boundary that experiences the Dirichlet boundary condition, v is the applied velocity, Γ N is the part of the boundary that experiences the Neumann boundary condition, and t is the applied time. Here, Γ = Γ D ∩ Γ N and Γ D ∩ Γ N = ∅. v 0 and p 0 are the initial velocities and pressure fields in the fluid in the domain, respectively. The pseudo-stress σ is given as follows: where I is the second-order identity tensor. The code uses an unstructured quadrilateral grid of the finite volume method, which has excellent adaptability. The calculation of the time-dependent physical quantities was transient. The number of grids in the calculation area was more than 50,000; therefore, the pressure field was calculated based on staggered grids using the semi-implicit method for pressure-linked equations (SIMPLE) to solve Equations (16)- (17) in order to reduce the calculation time. The models were saturated, and only the solid and liquid phases were in the calculation area.
The inlet water pressure in the model was consistent with that of the laboratory test. The flow state in the simulation was laminar flow; therefore, the velocity ratio at the outlet of the same model to the hydraulic gradient was invariant under different hydraulic gradients. This was defined as the hydraulic conductivity of the model. The validity of the simulation was verified by comparing the hydraulic conductivity of the simulation and the results of the constant-head permeability test. The root-mean-square error (RMSE) [52], Pearson correlation coefficient (PCC) [53], and Nash-Sutcliffe model efficiency coefficient (NSE) [54] were used to evaluate the accuracy of the simulation results: where m is the number of data points; k Lab,i and k Sim,i are the ith experimental and simulated k, respectively; and k Lab and k Sim are the equivalent experimental and simulated mean k, respectively. The RMSE can vary from 0 to +∞. A smaller RMSE indicates a bettersimulated data fit for the experimental data. The PCC varies from −1 to 1, where higher values indicate better data congruence. The NSE varies from −∞ to 1 and can assess hydrological models' predictive power, where a value closer to 1.0 indicates a better match between the observed and modeled values. We used Tecplot to process the calculation results from FLUENT in order to generate flow field diagrams. The velocity field is a vector field that describes a fluid's velocity distribution at several points in space. The velocity field is constant when the fluid flows stably and does not change over time.

Gravel Pack Shape, Particle Size Distribution, and Porosity
Box-plots of the α and S of the 300 sampled particles from the different size groups are shown in Figure 6. The figure indicates that the distribution of the particles' α was positively skewed, while that of S was negatively skewed. The mean value of α and S of 100 samples in each size group was taken as the gravel's shape parameter in this size group, which were recorded as α A and S A , respectively. The D I , D II , and D III particle size distributions of the sampled gravels are shown in Figure 7a-c. The particle size distributions D1, D2, and D3 of the particle packs after the CT scanning and permeability tests, following the method described in Section 2.2.1, are shown in Figure 7d.     Figure 8 depicts the maximum and minimum porosities of the samples with D1, D2, and D3 particle size distributions according to the results of the proctor compaction test. As can be seen from the bar graphs, the porosity range variation was the most significant for gravels and the most minor for glass balls for the same particle size distribution. Therefore, the porosity of packs should be determined using the glass balls' porosity range. The porosities of the samples in D1, D2, and D3 were set to 38.81% (P1), 32.29% (P2), and 31.22% (P3), respectively.
Because the particle size distribution of the test samples was obtained by multiplying the particle content of the sampled particle size distributions of 2-5 mm, 5-10 mm, and 10-20 mm by the percentage in Table 2, the aspect ratio of a particle pack (α P ) was equal to the sum of α A for each particle size distribution multiplied by the percentage in Table 2. Similarly, the roundness of a particle pack (S p ) was equal to the sum of S A for each particle size distribution multiplied by the percentage from Table 2. The particle shape parameters of the glass ball pack and plastic octahedron pack were fixed values that were independent of the size distribution. and D3 particle size distributions according to the results of the proctor compaction test. As can be seen from the bar graphs, the porosity range variation was the most significant for gravels and the most minor for glass balls for the same particle size distribution. Therefore, the porosity of packs should be determined using the glass balls' porosity range. The porosities of the samples in D1, D2, and D3 were set to 38.81% (P1), 32.29% (P2), and 31.22% (P3), respectively. Because the particle size distribution of the test samples was obtained by multiplying the particle content of the sampled particle size distributions of 2-5 mm, 5-10 mm, and 10-20 mm by the percentage in Table 2, the aspect ratio of a particle pack (αP) was equal to the sum of αA for each particle size distribution multiplied by the percentage in Table  2. Similarly, the roundness of a particle pack (Sp) was equal to the sum of SA for each particle size distribution multiplied by the percentage from Table 2. The particle shape parameters of the glass ball pack and plastic octahedron pack were fixed values that were independent of the size distribution.
Information on the test samples is listed in Table 3. Information on the test samples is listed in Table 3. α p is the aspect ratio of the particle pack, S p is the roundness of the particle pack, and ϕ is the porosity of the pack.  Figure 9 shows the segmentation result of a CT image from S1 using SVM. As shown in Figure 7b, the outlines of pores and particles after image segmentation were clear, demonstrating that the image segmentation code was effective. All CT images were processed in batches. Due to the significant number of CT images of all samples, only the segmentation results of one image are shown here. The edges of pores in the other segmented images are also precise.

D3
Glass ball B3 αp is the aspect ratio of the particle pack, Sp is the roundness of the particle pack, and φ is the poros of the pack. Figure 9 shows the segmentation result of a CT image from S1 using SVM. As show in Figure 7b, the outlines of pores and particles after image segmentation were cle demonstrating that the image segmentation code was effective. All CT images were pr cessed in batches. Due to the significant number of CT images of all samples, only t segmentation results of one image are shown here. The edges of pores in the other se mented images are also precise. There is currently no standard for verifying the validity of image segmentation r sults. The outline of the segmented target can only be observed by the human eye to s whether it is clear and complete. The accuracy of calculating the pore-particle interface fractal dimension and the pore surface area depends entirely on the accuracy of the ima segmentation. Therefore, we verified the validity of the SVM by comparing the segme tation results of the same CT image using SVM, the gray-scale morphology method [5 There is currently no standard for verifying the validity of image segmentation results. The outline of the segmented target can only be observed by the human eye to see whether it is clear and complete. The accuracy of calculating the pore-particle interface's fractal dimension and the pore surface area depends entirely on the accuracy of the image segmentation. Therefore, we verified the validity of the SVM by comparing the segmentation results of the same CT image using SVM, the gray-scale morphology method [55], and the histogram method [56]. Gray-scale morphology and histogram segmentation are common segmentation methods. Figure 10 shows the comparison of O1 s CT image segmentation result. Among the three materials, the CT image segmentation of plastic octahedrons was the most difficult because of its non-uniform color in the CT images due to the uneven particle density. The plastic particles had a high density at the edges and a low density in the middle, which could not be avoided when manufacturing due to their pouring process. In addition to the problem of the non-uniform density of plastic particles, the densities of the plastic particles and the resin permeameter were close, which was also an obstacle. In Figure 10, purple represents the particles segmented from the image. The original CT image was used as the background to discriminate the segmentation effect. The gray-scale morphology method could not separate the particles from the permeameter in the CT image of O1, according to Figure 10b. Using the histogram method not only failed to separate the particles from the permeameter but also over-divided the particles in the same image, according to Figure 10c. This showed that the SVM method segmented the CT image of O1 with significantly better results than the other two methods. SVM had advantages when segmenting the same heterogeneous object or different objects with similar densities. gray-scale morphology method could not separate the particles from the permeameter in the CT image of O1, according to Figure 10b. Using the histogram method not only failed to separate the particles from the permeameter but also over-divided the particles in the same image, according to Figure 10c. This showed that the SVM method segmented the CT image of O1 with significantly better results than the other two methods. SVM had advantages when segmenting the same heterogeneous object or different objects with similar densities.

Fractal
Analysis of the Pore-Particle Interface Based on CT Images Figure 11 illustrates the change in the box dimension of the pore-particle interface with different particle pack shapes. According to Figure 11a, the box dimension of the pore-particle interface (D) in each CT image of all the packs ranged from 1.804 to 1.865. All box dimensions of the pore-particle interface were above their topological dimension (the value of which is 1) and less than two dimensions (the value of which is 2), demonstrating that the pore-particle interface of the samples showed prominent fractal characteristics based on the theory of Mandelbrot [30]. The median line of boxes was not located in the middle of the box, indicating that the distribution of D of each sample presented a skewed distribution. The points outside the box in the figure are outliers. The reason was that the CT machine scanned the sample spiral forward for slice scanning, and some particles in the CT image appeared to be suspended without contact. Particles came into contact in three-dimensional space; however, the two-dimensional slice may not cut to the contact point, which is a normal phenomenon. When calculating the average box dimension of each sample's pore-particle interface (D A ), these outliers must be deleted before the calculation. The calculation results are shown in Figure 11b. Unless otherwise specified, later box dimensions mentioned refer to the average box dimensions. Except for the standard deviation of the D A of the octahedron pack for the D2 and D3 particle size distributions being 0.003, the standard deviation of the D A of the other packs was 0.002. As can be seen from the second bar chart, the D A of the ball pack was the smallest, and there was little difference between the D A of the gravel pack and that of the octahedron pack for the same particle size distribution. For the D1 particle size distribution, the D A of the ball pack was 1.810. For the D2 particle size distribution and D3 particle size distribution, the D A values of the glass ball pack were 1.822 and 1.846, respectively. The D A values of the gravel pack and octahedron pack were 0.31% and 0.64% larger than that of the ball pack for the D1 particle size distribution, respectively; the D A values of the gravel pack and octahedron pack were 1.12% and 1.14% larger than that of the ball pack for the D2 particle size distribution, respectively; both the D A values of the gravel pack and the octahedron pack were 0.86% larger than that of the ball pack for the D3 particle size distribution. For the D2 particle size distribution, the content of fine particles (less than 5 mm) in the packs ranked second among the three particle size distributions. However, for the D2 particle size distribution, the difference between the D A of the ball pack and the D A of the other two shaped particle packs was the largest.
sion of each sample's pore-particle interface (DA), these outliers must be deleted before the calculation. The calculation results are shown in Figure 11b. Unless otherwise specified, later box dimensions mentioned refer to the average box dimensions. Except for the standard deviation of the DA of the octahedron pack for the D2 and D3 particle size distributions being 0.003, the standard deviation of the DA of the other packs was 0.002. As can be seen from the second bar chart, the DA of the ball pack was the smallest, and there was little difference between the DA of the gravel pack and that of the octahedron pack for the same particle size distribution. For the D1 particle size distribution, the DA of the ball pack was 1.810. For the D2 particle size distribution and D3 particle size distribution, the DA values of the glass ball pack were 1.822 and 1.846, respectively. The DA values of the gravel pack and octahedron pack were 0.31% and 0.64% larger than that of the ball pack for the D1 particle size distribution, respectively; the DA values of the gravel pack and octahedron pack were 1.12% and 1.14% larger than that of the ball pack for the D2 particle size distribution, respectively; both the DA values of the gravel pack and the octahedron pack were 0.86% larger than that of the ball pack for the D3 particle size distribution. For the D2 particle size distribution, the content of fine particles (less than 5 mm) in the packs ranked second among the three particle size distributions. However, for the D2 particle size distribution, the difference between the DA of the ball pack and the DA of the other two shaped particle packs was the largest.   Figure 12 describes the change in the average box dimension of the pore-particle interface with a change in the shape of the particle pack for three particle size distributions. According to Figure 10a, D A increased with an increase in α p when α p was less than 2.154. D A remained almost stable when α p was larger than 2 and less than 5. As shown in Figure 10b, D A decreased with a decrease in S p when S p was larger than 0.862 and less than 1. D A remained nearly unchanged when S p was larger than 0.72 and less than 0.862. Both figures reflect that the influence of the particle pack shape parameters on D A was the most significant for the D2 particle distribution.  Figure 12 describes the change in the average box dimension of the pore-particle interface with a change in the shape of the particle pack for three particle size distributions. According to Figure 10a, DA increased with an increase in αp when αp was less than 2.154. DA remained almost stable when αp was larger than 2 and less than 5. As shown in Figure  10b, DA decreased with a decrease in Sp when Sp was larger than 0.862 and less than 1. DA remained nearly unchanged when Sp was larger than 0.72 and less than 0.862. Both figures reflect that the influence of the particle pack shape parameters on DA was the most significant for the D2 particle distribution.
(a) (b) Figure 12. The relationship between the particle pack shape and the average box dimension of the pore-particle interfaces. The error lines in the figures were based on the standard deviation. (a) The relationship between the average box dimension of the pore-particle interface and the aspect ratio of the packs. (b) The relationship between the average box dimension of the pore-particle interface and the roundness of the packs. Figure 13 shows the CT image of the middle vertical sections of B1, S1, and T1. In the figure, pink marks the point-to-point contact relationships and blue marks the edgeto-edge contact relationships. The length marked in the figure is not the actual contact length. The format of Figure 13 is TIFF, which is not the standard format of CT images used for display in other papers. In fact, CT images come in a DICOM format, which has a higher resolution and can only be read using special CT image editing software. We annotated and counted the contact relationship of CT images in MIMICS. We defined a Figure 12. The relationship between the particle pack shape and the average box dimension of the pore-particle interfaces. The error lines in the figures were based on the standard deviation. (a) The relationship between the average box dimension of the pore-particle interface and the aspect ratio of the packs. (b) The relationship between the average box dimension of the pore-particle interface and the roundness of the packs. Figure 13 shows the CT image of the middle vertical sections of B1, S1, and T1. In the figure, pink marks the point-to-point contact relationships and blue marks the edge-toedge contact relationships. The length marked in the figure is not the actual contact length. The format of Figure 13 is TIFF, which is not the standard format of CT images used for display in other papers. In fact, CT images come in a DICOM format, which has a higher resolution and can only be read using special CT image editing software. We annotated and counted the contact relationship of CT images in MIMICS. We defined a particle contact length of less than 1 mm as a point-to-point contact relationship and a length greater than 1 mm as an edge-to-edge contact relationship. In Figure 13, B1 has 42 point-to-point contact relationships, S1 has 68 point-to-point contact relationships and 24 edge-to-edge contact relationships, and O1 has 93 point-to-point contact relationships and 54 edgeto-edge contact relationships. CT images can only show the two-dimensional contact of particles. In fact, in three-dimensional space, there were also edge-to-surface contact relationships in the gravel pack and octahedral pack. We calculated the type and quantity of contact relationship of all CT images of nine test samples. There was only point-topoint contact relationship between spheres, which made it the most special. Under the same condition, compared with other shaped particle packs, the number of contact points between spheres was the least. From the perspective of the pore space, the pore structure of spheres was the simplest, and the average box dimension of the pore-particle interface was the smallest. Once the particles deviated from the sphere, the number of contact points soared due to there being more contact relations between particles. This led to the average box dimension of the pore-particle interface and pore structure complexity skyrocketing. When the shape parameters of particle packs reach a critical value, the box dimension of the pore-particle interface reaches a plateau. Unfortunately, we did not find that critical value due to the few kinds of materials used in the laboratory test. We will research the critical value through numerical simulation in future studies.
Li et al. [57] established a new method for evaluating the complexity of digital rock pore structure using the relative value of the box dimension and verified this method's effectiveness by calculating some stones' CT images. They set up a series of rock models with porosities of 0.05 to 0.4 and calculated the three-dimensional box dimension of the corresponding pore-rock interface in the range of 2.2-2.6. Their calculation results showed that when the rock's porosity was 0.05-0.2, the box dimension increased sharply with a rise in porosity. However, when the porosity exceeded 0.2, the box dimension grew slowly. This indicated that there was a critical porosity. When the porosity of particles exceeded the critical porosity, the box dimension increased slightly. The box dimension under the critical porosity can be used as a parameter to judge the complexity of the pore structure. Although our model and method differed from Li's calculation, we obtained the same law regarding the box dimension. As shown in Figure 8, when the three shapes of particle packs with the same particle size distribution had their minimum porosity, the porosity ranking of the samples was as follows: porosity of the sphere pack > porosity of the gravel pack > porosity of the plastic octahedron pack. This indicated that the larger the particle pack's aspect ratio, the smaller the porosity. The relationship between the particle pack's roundness and particles' porosity was the opposite. The average box dimension we calculated was the average of the box dimensions of all CT images of each sample, which is a parameter that was between two-dimensional and three-dimensional. The relationship between the average box dimension of the pore-particle interface and the particle pack shape parameters in Figure 12 showed that there was a critical shape parameter, and the average box dimension corresponding to the critical shape parameter could be used as the basis for judging the complexity of the pore structure and evaluating the particle pack shape of the gravel. Our calculations were based on the CT images of actual particle packs, while Li used numerical models. Our conclusions were consistent because of the relationship between the particle shape parameters and the porosity. Wang et al. [58], Han et al. [59], Xiu et al. [60], and Ari et al. [61] also reported similar conclusions when using simulations.
We illustrated the existence of the critical box dimension from the actual particle pack and verified their conclusions via a laboratory experiment. critical value, the box dimension of the pore-particle interface reaches a plateau. Unfortunately, we did not find that critical value due to the few kinds of materials used in the laboratory test. We will research the critical value through numerical simulation in future studies.
(a) (b) (c) Li et al. [57] established a new method for evaluating the complexity of digital rock pore structure using the relative value of the box dimension and verified this method's effectiveness by calculating some stones' CT images. They set up a series of rock models with porosities of 0.05 to 0.4 and calculated the three-dimensional box dimension of the corresponding pore-rock interface in the range of 2.2-2.6. Their calculation results showed that when the rock's porosity was 0.05-0.2, the box dimension increased sharply with a rise in porosity. However, when the porosity exceeded 0.2, the box dimension grew slowly. This indicated that there was a critical porosity. When the porosity of particles exceeded the critical porosity, the box dimension increased slightly. The box dimension under the critical porosity can be used as a parameter to judge the complexity of the pore structure. Although our model and method differed from Li's calculation, we obtained the same law regarding the box dimension. As shown in Figure 8, when the three shapes of particle packs with the same particle size distribution had their minimum porosity, the porosity ranking of the samples was as follows: porosity of the sphere pack > porosity of

Analysis of the Specific Surface Area of the Reconstructed Pore Network Model
After the image segmentation, the pore images were imported into MIMICS to generate a 3D model (as shown in Figure 14) of the pore network of packs. We used porosity as the determination criterion to verify the accuracy of the pore model according to Equations (10) and (11), as shown in Table 4. In Table 4, all the relative errors between the porosity of the model and the sample were less than 4%, which implied that the relative errors between the model's surface area and the sample's surface area were also less than 4%. This demonstrated that the method used to measure the surface area of porous media in Section 2.3 was appropriate. erate a 3D model (as shown in Figure 14) of the pore network of packs. We used porosity as the determination criterion to verify the accuracy of the pore model according to Equations (10) and (11), as shown in Table 4. In Table 4, all the relative errors between the porosity of the model and the sample were less than 4%, which implied that the relative errors between the model's surface area and the sample's surface area were also less than 4%. This demonstrated that the method used to measure the surface area of porous media in Section 2.3 was appropriate. Figure 14. The pore network model of the packs. All models were cut in the middle to show their internal structure.  V is the volume of the model of the pore network, ϕ sim is the porosity of the model of the pore network, and ϕ Lab is the porosity of the pack. R is the relative error. Figure 15 shows the specific surface area of the pore network (A p ) values of all the packs. As can be seen from this graph, the A p of the ball pack was the smallest for each particle size distribution. For the D1 particle size distribution, the A p of the ball pack was 7.39 cm −1 . For the D2 particle size distribution and D3 particle size distribution, the A p values of the ball packs were 13.37 cm −1 and 16.88 cm −1 , respectively. Moreover, the A p values of the gravel pack and octahedron pack were 23.7% and 49.7% larger than that of the ball pack for the D1 particle size distribution, respectively; the A p values of the gravel pack and octahedron pack were 26.4% and 33.1% larger than that of ball pack for the D2 particle size distribution, respectively; and the A p values of the gravel pack and octahedron pack were 14.0% and 20.1% larger than that of the ball pack for the D3 particle size distribution, respectively. The A p of the octahedron pack was the largest for each particle size distribution. In particular, the difference between the A p of the octahedron pack and that of the ball pack was the largest for the D1 particle size distribution. The A p of the gravel pack was the smallest for the D1 particle size distribution and the largest for the D3 particle size distribution. The A p values of the ball pack and octahedron pack also showed the same trend. dron pack were 14.0% and 20.1% larger than that of the ball pack for the D3 particle size distribution, respectively. The Ap of the octahedron pack was the largest for each particle size distribution. In particular, the difference between the Ap of the octahedron pack and that of the ball pack was the largest for the D1 particle size distribution. The Ap of the gravel pack was the smallest for the D1 particle size distribution and the largest for the D3 particle size distribution. The Ap values of the ball pack and octahedron pack also showed the same trend.  Figure 16 demonstrates that the pore network's specific surface area varied with the particle pack's shape for three particle size distributions. According to Figure 16a, the specific surface area of the pore network increased with an increase in the aspect ratio of the particle pack. Conversely, there was a downward trend in the specific surface area of the pore network with the growth in the roundness of the particle pack based on Figure 16b.  Figure 16 demonstrates that the pore network's specific surface area varied with the particle pack's shape for three particle size distributions. According to Figure 16a, the specific surface area of the pore network increased with an increase in the aspect ratio of the particle pack. Conversely, there was a downward trend in the specific surface area of the pore network with the growth in the roundness of the particle pack based on Figure 16b. The KC formula and other modified forms [6,[13][14][15][16][62][63][64] require the surface area of particles rather than the surface area of the pore network to predict the permeability coefficient. However, the fluid flows along the particle surfaces and the calculation area's side wall. The KC formula and other modified forms ignore the flow along the side wall, resulting in a smaller calculated hydraulic conductivity. Considering the influence of the sample volume on the pore surface area, we suggest using the surface area of the pore network per unit volume, i.e., the specific surface area of the pore network, to predict the permeability of the gravel. Furthermore, those formulas are built through simulation, which is difficult to verify and use directly in laboratory and field tests. The new method we proposed to measure the surface area of the pore network of real gravel was reliable and is expected to increase the practicality of these formulas. Figure 17 illustrates the results of the constant-head permeability tests of the packs. The first line graph shows the relationship between the water velocity (V) and the hydraulic gradient (J) at the test water temperature (which was 14 °C). In the test, the flow velocity The KC formula and other modified forms [6,[13][14][15][16][62][63][64] require the surface area of particles rather than the surface area of the pore network to predict the permeability coefficient. However, the fluid flows along the particle surfaces and the calculation area's side wall. The KC formula and other modified forms ignore the flow along the side wall, resulting in a smaller calculated hydraulic conductivity. Considering the influence of the sample volume on the pore surface area, we suggest using the surface area of the pore network per unit volume, i.e., the specific surface area of the pore network, to predict the permeability of the gravel. Furthermore, those formulas are built through simulation, which is difficult to verify and use directly in laboratory and field tests. The new method we proposed to measure the surface area of the pore network of real gravel was reliable and is expected to increase the practicality of these formulas. Figure 17 illustrates the results of the constant-head permeability tests of the packs. The first line graph shows the relationship between the water velocity (V) and the hydraulic gradient (J) at the test water temperature (which was 14 • C). In the test, the flow velocity of all packs was stable, and the relationship between the flow velocity and the hydraulic gradient was basically linear. The water flow in the particles was laminar, and the hydraulic conductivity (k) could be calculated using the Darcy formula. The particles did not move, and there was no seepage deformation in the test. In Figure 15a, when the hydraulic gradient was less than 0.1, the calculated hydraulic conductivity differed from that calculated for other gradients. When debugging the instrument, we found that this was due to the limited accuracy of the water head measuring instrument when the hydraulic gradient was less than 0.1. Therefore, we eliminated the velocity corresponding to a hydraulic gradient of less than 0.1. The number of velocity groups of the samples was sufficient after applying the criterion that the hydraulic gradient needed to be larger than 0.1. The hydraulic conductivity of the samples is shown in Figure 15b. As can be seen from the bar graph, for the D1 particle size distribution, the hydraulic conductivity of the ball pack was 0.32 times and 2.39 times higher than those of the gravel pack and octahedron pack, respectively. Additionally, for the D2 particle size distribution, the hydraulic conductivity of the ball pack was 1.28 times and 3.12 times higher than those of the gravel pack and octahedron pack, respectively. For the D3 particle size distribution, although the hydraulic conductivity of the ball pack was more than 1.5 times higher than those of the gravel pack and octahedron pack, there was a slight difference between the hydraulic conductivity of the gravel pack and the octahedron pack. Therefore, we believe that when the content of fine particles (particle size less than 5 mm) was in a particular range with the same particle size distribution and porosity, the main factor that affected the hydraulic conductivity of porous media is the particle shape. However, with an increase in the fine particle content, the effect of the particle shape on the permeability of porous media decreased (except for spheres). pack and octahedron pack, there was a slight difference between the hydraulic conductivity of the gravel pack and the octahedron pack. Therefore, we believe that when the content of fine particles (particle size less than 5 mm) was in a particular range with the same particle size distribution and porosity, the main factor that affected the hydraulic conductivity of porous media is the particle shape. However, with an increase in the fine particle content, the effect of the particle shape on the permeability of porous media decreased (except for spheres).  Figure 18 illustrates the relationship between the hydraulic conductivity and the pore characteristics of the packs. According to Figure 18a, there was a downward trend in the hydraulic conductivity with an increase in the aspect ratio of the particle pack for the same particle size distribution. The hydraulic conductivity of the same material pack was the largest for the D1 particle size distribution and the smallest for the D3 particle size distribution, which is consistent with the Hazen formula [4], Kozeny formula [5], Kozeny-Carman formula [6], and Terzaghi formula [7]. For the D1 particle size distribution, the hydraulic conductivity of the packs seemed to decrease uniformly with an increase in the particle pack's aspect ratio. However, for the particle size distribution of D2 and D3, the hydraulic conductivity of the packs decreased gradually when the aspect ratio of the particle pack was larger than 2.1. According to Figure 18b, the relationship between the packs' hydraulic conductivity and roundness was the opposite. Moreover, the hydraulic conduc-  Figure 18 illustrates the relationship between the hydraulic conductivity and the pore characteristics of the packs. According to Figure 18a, there was a downward trend in the hydraulic conductivity with an increase in the aspect ratio of the particle pack for the same particle size distribution. The hydraulic conductivity of the same material pack was the largest for the D1 particle size distribution and the smallest for the D3 particle size distribution, which is consistent with the Hazen formula [4], Kozeny formula [5], Kozeny-Carman formula [6], and Terzaghi formula [7]. For the D1 particle size distribution, the hydraulic conductivity of the packs seemed to decrease uniformly with an increase in the particle pack's aspect ratio. However, for the particle size distribution of D2 and D3, the hydraulic conductivity of the packs decreased gradually when the aspect ratio of the particle pack was larger than 2.1. According to Figure 18b, the relationship between the packs' hydraulic conductivity and roundness was the opposite. Moreover, the hydraulic conductivity of the packs grew significantly when the roundness of the pack was larger than 0.86 for the particle size distribution of D2 and D3. Meanwhile, for the D1 particle size distribution, the hydraulic conductivity of the packs rose uniformly with an increase in the roundness of the pack. Therefore, similar to the conclusion drawn from Figure 18a, when studying the influence of the particle shape parameters on the permeability of the gravel, packs with different fine particle contents differed from one another. Figure 18c depicts the variation in the hydraulic conductivity of the packs with the average box dimension of the pore-particle interface. As can be seen from this line graph, the hydraulic conductivity of the packs decreased with an increase in the average box dimension of the pore-particle interface for the same particle size distribution. For the D2 particle size distribution, the variation in hydraulic conductivity with the average box dimension of the pore-particle interface from 1.822 to 1.842 was less than those of others. We think this was because the spatial distribution of particles was random and uneven. Coarser pore channels were formed in the middle of some particles in the S2 sample, resulting in its hydraulic conductivity deviating from those of the other packs. As shown in Figure 18d, the hydraulic conductivity of the packs for the same particle size distribution decreased with an increase in the specific surface area of the pore network. The three lines in this figure are approximately parallel, which makes them significantly different from the second line in Figure 18c. Therefore, when studying the relationship between pore characteristics and permeability of the porous media, multiple parameters describing pore characteristics should be used to describe its characteristics comprehensively.

Constant-Head Permeability Test Result
Our results are in agreement with those of Zakhari et al. [65], Liu et al. [66], Conzelmann et al. [67], and Li et al. [68]. In addition to the particle pack shape, we also found that the content of fine particles less than 5 mm in the porous media dramatically influenced a shape's sensitivity to permeability. We believe this was because there were different contact relationships between particles with different shapes. Contact relationships can be divided into five types: point-to-point, point-to-face, edge-to-edge, edge-to-face, and face-to-face. The sphere pack only had the point-to-point contact relationship, according to Figure 13. With an increase in the aspect ratio or a decrease in the roundness of the particle packs, the other four contact relationships appeared. Generally, the two relationships of edge-to-edge and edge-to-face appeared more in the frontal body pack and less in naturally formed particle packs, such as the gravel pack. Different contact relationships and the number of each relationship affected the pore structure characteristics of the porous media, resulting in a difference in permeability. Some scholars think that the influence of the shape on the permeability of porous media is negligible [18]. We believe that it depends on the different types and number of contact relationships and the fine particle content. For porous media with less fine particle content, such as the D1 grading in our test, it was evident that the influence of shape on permeability could not be ignored. Our laboratory test results implied that the type and number of contact relations affected the shape sensitivity to permeability. The influence of shape on permeability can be ignored in a particular range, but the shape needs to be considered beyond this range. When the shape factor is considered in the formula for predicting permeability, the conditions for the type and number of contact relations should be given. That is also why some scholars add the shape factor to the formula to predict hydraulic conductivity, but the results are only consistent with their experiments or simulations.
conductivity deviating from those of the other packs. As shown in Figure 18d, the hydraulic conductivity of the packs for the same particle size distribution decreased with an increase in the specific surface area of the pore network. The three lines in this figure are approximately parallel, which makes them significantly different from the second line in Figure 18c. Therefore, when studying the relationship between pore characteristics and permeability of the porous media, multiple parameters describing pore characteristics should be used to describe its characteristics comprehensively. Our results are in agreement with those of Zakhari et al. [65], Liu et al. [66], Conzelmann et al. [67], and Li et al. [68]. In addition to the particle pack shape, we also found that the content of fine particles less than 5 mm in the porous media dramatically influenced a shape's sensitivity to permeability. We believe this was because there were different contact relationships between particles with different shapes. Contact relationships can be divided into five types: point-to-point, point-to-face, edge-to-edge, edge-to-face, and face-to-face. The sphere pack only had the point-to-point contact relationship, according to Figure 13. With an increase in the aspect ratio or a decrease in the roundness of the

Simulation
As was mentioned in Section 3.2, the instrument for measuring the water head was inaccurate when the hydraulic gradient was less than 0.1. Therefore, when using FLUENT to simulate the laboratory test, we only simulated a flow field with a hydraulic gradient ranging from 0.1 to 0.3.
The RMSE, PCC, and NSE values, according to Equations (18)- (20), are listed in Table 5. As can be seen from the table, the maximum RMSE was 0.322 for S1, the minimum PCC was 0.998 for both S1 and B1, and the minimum NSE was 0.839 for S1. The data analysis results showed that the model, grid, and parameter settings used in our numerical simulation were appropriate. Figure 19 shows the flow fields of the models for the D1 particle size distribution with the 0.1 hydraulic gradient. We found that the flow field differences between the models were similar. Due to the limited space in this paper, we only present the model calculation results for the D1 particle size distribution for illustration purposes. The streamline diagrams with a velocity greater than 0.09 cm/s and running from the water inlet to the water outlet with the 0.1 hydraulic gradient are also shown in the same graph. (b) S1 (c) T1 Figure 19. Partial flow field diagrams obtained from the simulation. The left column graphs are the velocity field diagrams with a hydraulic gradient of 0.1, and the right column graphs are the streamline diagrams with a velocity greater than 0.09 cm/s, flowing from the inlet to the outlet with the same hydraulic gradient. Figure 19. Partial flow field diagrams obtained from the simulation. The left column graphs are the velocity field diagrams with a hydraulic gradient of 0.1, and the right column graphs are the streamline diagrams with a velocity greater than 0.09 cm/s, flowing from the inlet to the outlet with the same hydraulic gradient.
The model presented a preferential flow according to the velocity field diagrams in Figure 19. The lighter color bands in the velocity field represent the path of the water particles with a velocity greater than 0.9 cm/s. Under pressure from the bottom to the top, the water did not flow evenly along all waterways but concentrated on a few paths. Macropore flow, bypass flow, and pipe flow existed in the models. These all belonged to typical preferential flow phenomena. As can be seen from the graph, the number and diameter of the preferential flow channels in the sphere pack model were the largest, and those of the octahedron pack model were the smallest. There were more preferential flow channels in the sphere packs from the inlet to the outlet. The preferential flow of the model of the other two shaped packs diverged after encountering the corner, showing a decrease in the flow rate from the macroscopic perspective.
The streamline is a curve that is tangent to the velocity vector at every point in the flow field, while the trace is the curve depicted by the fluid mass as it moves through space. For a constant flow, the streamline and trace coincide. Therefore, as long as the tortuosity of the streamline is calculated, the tortuosity of the trace is obtained. The number of seed points for all simulation results that generated the streamlines was set to 10. The tortuosity (τ) is an important parameter that describes the seepage channels. It is defined as the actual length of the seepage channel and the apparent length through the samples; that is, it is the exact length of the water's motion track in the channel when the permeate water passes through a unit distance of a pack: where L t is the length of the curved line and L o is the length of the line along the direction of the hydraulic gradient. Figure 20 depicts the statistical results of the tortuosity of the traces of all models. According to the box plot, for the same particle size distribution, the τ of the ball pack was the smallest and that of the octahedron pack was the largest. The greater the tortuosity, the longer the actual path of water flow, and the macroscopic results were that the cross-section flow of porous media decreased and the permeability decreased over time. The same material pack had the smallest τ for D1 and the largest τ for D3, which is consistent with the law of hydraulic conductivity. The distribution of τ was not uniform because the two-dimensional simulation of the seepage field had limitations. Since a three-dimensional porous media pore network model is too complex in space to mesh appropriately, we only simulated the laboratory test in two dimensions. Although the τ in Figure 20 does not represent the actual tortuosity of the water path of porous media in the laboratory test, our two-dimensional simulation was based on the CT images of the packs, and the pore structure in the model was the same as the actual pore structure. The characteristic of the flow field is convincing. (c) T1 Figure 19. Partial flow field diagrams obtained from the simulation. The left column graphs are the velocity field diagrams with a hydraulic gradient of 0.1, and the right column graphs are the streamline diagrams with a velocity greater than 0.09 cm/s, flowing from the inlet to the outlet with the same hydraulic gradient.
The streamline is a curve that is tangent to the velocity vector at every point in the flow field, while the trace is the curve depicted by the fluid mass as it moves through space. For a constant flow, the streamline and trace coincide. Therefore, as long as the tortuosity of the streamline is calculated, the tortuosity of the trace is obtained. The number of seed points for all simulation results that generated the streamlines was set to 10. The tortuosity (τ) is an important parameter that describes the seepage channels. It is defined as the actual length of the seepage channel and the apparent length through the samples; that is, it is the exact length of the water's motion track in the channel when the permeate water passes through a unit distance of a pack: where Lt is the length of the curved line and Lo is the length of the line along the direction of the hydraulic gradient. Figure 20 depicts the statistical results of the tortuosity of the traces of all models. According to the box plot, for the same particle size distribution, the τ of the ball pack was the smallest and that of the octahedron pack was the largest. The greater the tortuosity, the longer the actual path of water flow, and the macroscopic results were that the cross-section flow of porous media decreased and the permeability decreased over time. The same material pack had the smallest τ for D1 and the largest τ for D3, which is consistent with the law of hydraulic conductivity. The distribution of τ was not uniform because the two-dimensional simulation of the seepage field had limitations. Since a three-dimensional porous media pore network model is too complex in space to mesh appropriately, we only simulated the laboratory test in two dimensions. Although the τ in Figure 20 does not represent the actual tortuosity of the water path of porous media in the laboratory test, our two-dimensional simulation was based on the CT images of the packs, and the pore structure in the model was the same as the actual pore structure. The characteristic of the flow field is convincing.

Contributions, Applications, and Limitations
In this study, the influence of particle shape on the pore characteristics and permeability of gravel was studied by combining laboratory tests and simulations. The contributions and applications are as follows:

1.
A method for estimating the shape of a gravel pack via manual sampling and measuring the shape factor of gravel with three-dimensional physical meaning was proposed. Previous studies showed that the hydraulic conductivity of the filter material in a natural quarry is higher than that in a blasting quarry under the same conditions [3,22]. However, because the shape of the gravel pack is difficult to quantify through the laboratory test, few studies have considered the particle shape in the design and selection of the filter material. Many dams are built in alpine and canyon areas for higher economic benefits. They use blasted gravel as dam building materials, which is different from the current design specification and does not consider the particle shape. According to our results, for better anti-seepage, the dam building material should preferably be blasted gravel with a large aspect ratio and small roundness. After particle size sieving in the stockyard, gravel of the same particle size group can be sieved again for particle shape using a rectangular or rhombus sieve.

2.
A method for measuring the surface area of gravel and its pore network based on SVM segmentation and the reconstruction of CT images was proposed. This method can self-adjust the parameters through deep learning to measure the surface area of particles with different densities and sizes, which is suitable for engineering applications. In addition, it has few requirements in terms of vessel materials and can be coupled with other hydraulic and mechanical tests. This method increases the practicability of the formula for predicting the hydraulic conductivity of gravel using the specific surface area in engineering applications [11,[14][15][16][62][63][64][69][70][71]. During a dam's construction, a laboratory is set up on the site to test the particle size distribution and dry density of the filled part to control the construction quality. As filling uses rolling technology, the gravel may break during the rolling process, which causes the actual particle size distribution and dry density of the filling material to deviate from the design value. Based on the consensus that material with a more significant specific surface area has low permeability [11,[14][15][16][62][63][64][69][70][71], the specific surface area can be added for construction control in the filter using our measurement method.
A limitation of this study was that the number of research samples was small due to the heavy workload of manually measuring the gravel shape. A new formula for predicting permeability by considering the gravel shape has not been proposed. In addition, due to the difficulty of meshing a three-dimensional reconstructed model using CT images, this study only calculated the seepage field of four CT slices of each sample from the two-dimensionality. We will conduct further research to expand the total number of measurements and overcome the obstacles of three-dimensional model meshing.

Conclusions
This study aimed to investigate the influence of gravel's shape on its pore structure and permeability through CT scanning tests, permeability tests, and simulations to provide suggestions for the design and construction of the earth rock dams' filter and increase the application of predictive formulas related to the shape factor and surface area of gravel in engineering applications. Some valuable conclusions were as follows:

1.
A new method was proposed to estimate the gravel pack shape; this method involved manual sampling and measuring the gravel's aspect ratio and roundness with threedimensional physical significance, which is expected to be popularized for the study of the shape of actual gravel packs and their related hydraulic or mechanical properties. One should pay attention to making the gravel pack's particle size distribution consistent with the sample's particle size distribution and control the particle size to the centimeter level when using this method.

2.
A new method was put forward that uses SVM segmentation and the reconstruction of CT images to measure the surface area of a gravel pack and its pore network. The advantage of the method is that it can be coupled with other hydraulic and mechanical tests and can automatically adjust the parameters according to different testing materials for convenient use in engineering. The specific surface area can be added for construction control of the filter using this method.

3.
A gravel pack with a larger aspect ratio and smaller roundness had a larger box dimension associated with its pore-particle interface and a greater specific surface area of the pore network, which meant it had a more complex pore structure. The content of particles less than 5 mm affected the relationship between the shape factor, pore-particle interface, and specific surface area of the pore network. The influence degree of particle shape was dependent upon the content of fine particles.

4.
A gravel pack with a larger aspect ratio and smaller roundness had a smaller hydraulic conductivity. This was because the CT scanning results showed that the larger the aspect ratio and the smaller the roundness, the more contact points and contact relationships there were between the particles in a gravel pack with a complex pore structure. This would increase the number, length, and tortuosity of the seepage channels when seepage occurs in such gravel packs in a two-dimensional simulation seepage field.

5.
In addition to allowing the particle size distribution and dry density to meet the requirements of the dam design specifications, the filter material should preferentially use blasting gravel with a larger aspect ratio and a smaller roundness for better anti-seepage performance.