Rapid modeling of space structure of rock mass and the method of equivalent assignment

Accurately grasping the spatial distribution characteristics of fractures is an important basis for computation and analysis of rock engineering. On the basis of field investigation, probability statistics, and numerical simulation, the fast modeling method of fracture with special probability and statistical distribution is proposed. The server connection between COMSOL and MATLAB is established. In the COMSOL with MATLAB programming environment, the fast modeling program of fissure space structure is developed. Aiming at the difficulty of mesh generation due to high fracture aspect ratio, programming language was used to interpolate the material parameters of fissure space structure and construct the equivalent assignment method of fissure space structure. Based on the field fracture investigation in an underground mine, probabilistic data analysis of fracture data is carried out, and the fast modeling and equivalent assignment method is applied to realize the efficient modeling and assignment for the field investigation of fissure spatial structure.


| INTRODUCTION
Because of the geological and engineering transformation, the underground rock mass contains a large number of cracks. The existence of cracks not only affects the mechanical properties of rock mass, but also affects the permeability and thermodynamics of the rock mass. 1 Damage, development, and penetration of fractures are one of the main reasons for deformation and failure of rock mass. 2 They are also important reasons leading to complex characteristics such as heterogeneity and anisotropy. It is very important for the calculation and analysis of rock mass engineering to accurately grasp the distribution characteristics of the spatial structure of the cracks inside the rock mass.
Scholars have done a lot of research on rock fracture modeling. Priest and Hudson 3 used distribution of trace lengths produced by the intersection of planar discontinuities with a planar rock face to determine the distribution of trace lengths, the distribution of semi-trace lengths, and the distribution of censored semi-trace lengths intersected by a randomly located scanline. Kulatilake  Markov chain geostatistical approach to synthesize the discrete permeability structure of moderately fractured rock. A numerical model for rock is proposed in which the rock is represented by a dense packing of non-uniform-sized circular or spherical particles that are bonded together at their contact points. Cundall and Strack 6 introduced measures that were used in a future constitutive model. The main measures, a partitioned stress tensor, and a constraint ratio, which is related to the stability of the assembly, are introduced and illustrated for some numerical experiments. Dershowitz et al 7 illustrates how fracture data obtained in boreholes can be evaluated to determine which type of spatial model best describes the fracture sets present. The mathematical formulation for each type of analysis is presented. Snow 8 found that if the lognormal shape prevails also in undisturbed rocks, then the mean and standard deviation computed from water tests may also disclose fracture sizes at any specified depth in a rock body. Turanboy 9 proposed new geometrical classifications of two discontinuities as a construction method according to the spatial orientations of the discontinuities and their locations relative to each other. Wang et al 10 put forward the idea of "feature point," and the fracture space points are classified to form a description and structural feature analysis method suitable for a variety of pore and fracture spaces. All these studies have made positive progress. However, the difficulty of fracture modeling is that the number of cracks is huge, and the distribution is random. As a structure with high aspect ratio, cracks often cause difficulties in mesh generation. 11,12 Therefore, according to the fracture modeling and meshing problems, a fast modeling and equivalent assignment method, a rock fracture spatial structure, was proposed, and the practical application shows that the rock mass spatial structure of fast modeling and equivalent assignment method has obvious advantages than conventional methods in both mesh efficiency and mesh quality, so as to provide a new and effective method for numerical analysis and engineering calculation of fractured rock mass.

| Modeling ideas
How to generate the fracture network which is closer to the actual situation is still a difficult problem. 13 Generally, the fracture occurrence in the same area will concentrate on several main directions and form several fracture groups. The formation of fissures is mainly influenced by geological history and geostress. For the characterization of concrete crack, the feasible approach is shown in Figure 1.
Although the large-scale fracture data analysis is derived from small-scale fracture data analysis, it cannot reflect the real situation of the internal fissure in the rock mass, but under the current research level and condition, 14 it is a widely accepted and feasible method.

| Sample numbers
The probability and statistical analysis of fissure data requires that the number of fractured samples should not be too small. According to the theory of probability and statistics, the minimum sample size is mainly affected by error and confidence level.
where n min is minimum number of samples, Z a/2 is statistics under confidence level(1-a/2), and m is the ratio of admissible error to standard deviation.

| Probability distribution function
The range of data samples of fracture data is divided into k interval, and the number of samples in each interval is n i (i = 1,2,…,k), and we calculate the sample frequency f i in each interval, that is: When the number of samples is enough, f i can approximate the probability of the random variable in the i interval.
Taking the segmented sample interval as the abscissa axis and the longitudinal coordinate as the sample frequency, we draw the frequency histogram of the fractured data samples, and then, we can analyze the probability function distribution of the fracture parameters.

| Hypothesis test
Due to the unknown probability distribution type of fracture parameters, it is necessary to use nonparametric test to test probability density function of fracture parameters when analyzing probability distribution. In this paper, the Pearson test method is adopted, and the principle is as follows: Assuming that the fracture parameter sample (total number: n) obey a specific probability and statistical distribution, the number of theoretical samples is obtained by using the probability density function, and the sample data are combined with the actual sample data of the interval.
For a given level of significance α, making: The criteria for Pearson test are as follows 15 : The value of the critical value 2 is less than the level of the Pearson 2 (through the look-up table), and we accept the probability density function of the assumed distribution fitting.

| Fast modeling method of fracture structure
According to the fracture characterization method proposed by ISRM 16 and actual demand, the fracture modeling parameters include crack attitude (orientation, orientation, angle), crack density (per unit length or volume within the fracture number), and length and width (fissure surface vertical distance and other parameters). The fracture sample data are obtained by field investigation, and the statistical analysis of probability distribution is carried out. The shape of the crack can be ellipsoid, cylindrical, and polygon, and it can also be a combination of three shapes.
The number of cracks in the natural rock mass is huge, it is not realistic to model the fracture individually, and the fracture may obey a certain probability distribution. 17,18 The traditional method is that the main faults and joints can only be modeled by manual assignment. For a large number of cracks, it is impossible to complete the modeling. Therefore, it is necessary to use the programming method to realize fast modeling of fissure obeying specific probability distribution. The numerical modeling platform used in this paper is the finite element software COMSOL Multiphysics, which provides abundant geometric types and has powerful programming functions. Therefore, by establishing the server connection between COMSOL and MATLAB ( Figure 2), with MATLAB's powerful programming function, we can realize the programming of fracture geometry structure modeling.
In the programming environment shown in Figures 2 and 3, a fast modeling of fractured space structure that obeys a specific probability distribution is realized through a custom programming language. Here, to illustrate the method of custom programming, we assumes that the fracture parameters obey uniform distribution of rock mass, the unit size is 1 m × 1 m × 1 m, the fissure number is 30, the trace length In the programming environment of COMSOL connecting to MATLAB, we run the crack modeling code written above and introduce the crack space structure model into the COMSOL software through the server connection between COMSOL and MATLAB. Through the above custom programming method, fast modeling of fracture can be realized according to the specific probability and statistical distribution. As shown in Figure 3, the ellipsoidal, cylindrical, and polygonal fissure spatial structure model of the ellipsoidal, cylindrical, and polygonal is set in turn by the custom programming method.

ASSIGNMENT METHOD
The length of fracture trace direction is usually larger than that of the other two directions (the general opening is 0.1~2 mm magnitude, and the trace length is m magnitude), and the crack size is relatively small relative to the engineering rock mass. If we regard the fracture as another kind of different materials to rock matrix, the mesh generation will be extremely difficult because of the high aspect ratio of fracture, which leads to difficulty that dissecting is huge, the computing demand is very high, and it cannot even be dissecting. 19,20 Therefore, a method of fracture equivalent assignment is put forward. The core idea is to substitute the equivalent interpolation of fractured space location to represent the spatial structure model of cracks, so as to reduce the number and difficulty of mesh generation.
The COMSOL software default method is used to mesh the fissure model established in the 2.3 section, and the results are shown in Table 1. From Table 1, we can see that the smallest element size of crack is 0.005 m, while the minimum unit size of rock mass is 0.18 m, and the size difference is more than 30 times. So the number of generated cells is about 3 000 000. In spite of such small size division, the average element mass of the mesh is only 0.51 (the optimal Interpolation equivalent assignment value is 1), because the gap between the fracture and the rock mass is too large. The number of networks will directly affect the efficiency of numerical computation. The algorithm used by COMSOL software is the finite element method. [21][22][23] Similarly, most algorithm-based finite element method is discrete (the computational domain is discretized into a finite element calculation, and the number of meshes directly affects the computational complexity of the solution).
Here, there are only three independent field variables for the analysis of the simple linear elastic mechanical problem in rock mechanics. Therefore, to carry out the analysis on fracture mesh model of linear elastic mechanics, the calculation of degree of freedom will reach about 12 000 000. Figure 4 is the correspondence between the degree of freedom in COMSOL and the calculated memory. It is not difficult to see that when the degree of freedom is 4 million, the memory required by the computer reaches 16GB, and when the degree of freedom is more than 10 million, the computer needs at least 48GB of memory. The analysis is based on linear elastic problems in rock mechanics for relatively simple assumptions, and when the analysis involves rock nonlinear, anisotropic, multifield coupling, and other complex problems, the degree of freedom will be far more than this number and the computing configuration required is very high.
Therefore, based on the feasibility and cost of the calculation, an equivalent assignment method for the fracture space structure is proposed, which is shown in Figure 5.
In the proposed equivalent assignment method, the derivation of the parameter value of the fissure material and the interpolation function need to be realized by programming. Similarly, in the programming language environment of the COMSOL connection to MATLAB, the derivation of the parameter value of the fissure material using the MATLAB M language is as follows: The implementation statements of the interpolation function are as follows: According to the COMSOL software default method and the proposed equivalent method, the fracture mesh is divided   and assigned. The comparison results are shown in Table 2. It is not difficult to see that the equivalent method can achieve the result of material parameter assignment that is consistent with the default method. Compared with the about 3 000 000 mesh numbers formed by the default partition, the number of mesh dissected by the equivalent method is less than 100 thousand, which greatly reduces the number of mesh generation. At the same time, the average quality of the grid cell divided by the default method is 0.51, while the average mass of the grid cell with the equivalent method is 0.77, which is also significantly higher than the default method. Therefore, the fracture equivalent assignment method is obviously superior to the conventional method, whether it is the efficiency of mesh generation and the quality of mesh generation.

| Field investigation of fissures
In a underground mine, a field survey of joints and cracks is carried out by measuring line method ( Figure 6). The primary problem of the survey is the identification of fractures. From the hydrogeological point of view, it is determined that fractures are permeable cracks in the rock mass, which are able to pass through the flow. In this investigation, joints that are filled with impermeable substances are not identified as fissures. On the other hand, the small fissure (or small fissure with very small gap) in the rock mass, which has only small water storage meaning, is also not the object of measurement. The investigation of fissures in the field should guarantee the number of samples to ensure the representativeness and reliability of the statistical results of the fracture. A total of 350 cracks were measured in the field survey, and the results were discrete. The field survey data are shown in Table 3 and Table 4. According to the dominant direction of fracture can be divided into four groups: 5°, 35°, 303°, and 345°; according to tend to have two groups: 75° with ∠60° and 293° with ∠67°.

| Probability and statistical analysis of fracture data
Take the field survey data of the slit angle as an example, the probability distribution method is used to analyze the probability distribution of the fissure dip angle. First, the frequency histogram of the fracture angle sample is drawn. As shown in Figure 7, the probability distribution function that the dip angle may obey is normal distribution. According to the normal distribution, fitting the inclination data, the average and standard deviation of the normal distribution are 60.92 and 7.18, respectively.
The Pearson test is used to test the probability distribution function of the dip angle, which assumes that the dip angle  obeys the normal distribution, and the results of the test are shown in Table 5. From Table 5, we can see that the Pearson test value of crack inclination obeys normal distribution (test value is 18.02), while the Pearson threshold value of the confidence interval 95% is 21.66, and the test value is less than threshold value. Therefore, it can be determined that the samples of fracture dip data obey normal distribution. Similarly, the frequency histograms of other fracture parameters (orientation, orientation, gap width, etc.) are shown in Figure 8, and the distribution test is done according to the same method. Finally, the distribution data of each fracture parameter sample are shown in Table 6.

| Fracture modeling and equivalent assignment
According to the results of fracture field investigation and probability statistical analysis (2.5.2), the rapid modeling and equivalent assignment of fissure space structure subject to specific probability distribution function are presented. The results are shown in Figures 9 (constructed by MATLAB) and 10 (constructed by COMSOL). The modeling and equivalent assignment method of fractured space structure lay a solid foundation for in-depth consideration of the complex structural characteristics of rock mass and accurate engineering calculation and analysis of rock mass ( Figure 10).

| CONCLUSIONS
Fracture structure is an important factor affecting rock mechanics and permeability. Accurately grasping the spatial distribution characteristics of fractures is an important basis for computation and analysis of rock engineering. Based on the field fracture investigation, a fast modeling and assignment method for fissure space structure is constructed by using probability statistics, programming, and numerical simulation. It is applied to the engineering example, and the conclusions are as follows: 1. On the basis of field investigation, probability and statistics, programming, and numerical simulation, the fast modeling method of fracture with special probability and statistical distribution is proposed. The server connection between COMSOL and MATLAB is established.
In the COMSOL with MATLAB programming environment, the fast modeling program of fissure space structure is written by MATLAB M language.

2.
Aiming at the difficulty of mesh generation due to high fracture aspect ratio (large number of units and high computation cost), we use programming language to interpolate the material parameters of fissure space structure and construct the equivalent assignment method of fissure space structure. The application shows that the equivalent method can achieve the  same assignment effect as the conventional method, and the equivalent method can greatly reduce the number of mesh generation and improve the modeling efficiency compared with the conventional method and significantly improve the quality of the mesh generation unit. 3. Based on the field fracture investigation in an underground mine, probabilistic data analysis of fracture data is carried out, and the fast modeling and equivalent assignment method is applied to realize the efficient modeling and assignment for the field investigation of fissure spatial structure.