Research on the influence of rock micro-parameter in particle flow code (PFC) on the rock fracture morphology

3.1. Fracture morphology in uniaxial compression test

According to the simulation results of the uniaxial compression test, all of the 1517 series of results are divided into three groups. The schematic figures and simulation results are shown in Figure 3

Figure 3:

3.2. Fracture morphology in Brazilian splitting test

Similar with the uniaxial compression test, the results of the Brazilian splitting test are also shown in three kinds of fracture morphologies, as shown in Figure 4.

Figure 4:

3.3. Fracture morphology in direct shear test

According to the simulation results, the appeared fracture morphologies can be divided into two groups, as show in Figure 5.

Figure 5:

3.4. Constraints of the fracture morphologies

To investigate the influence of the parameters of the numerical model on the fractures morphologies, the simulation results should be classified according to their fracture morphologies. Before it, criteria of different morphologies are proposed according to the main fracture distribution.

To quantitatively identify the different fracture morphology, the distribution area of the 80% simulated fractures is used as the quantitatively criterion. In detail, as Figure 6 shows, each fracture corresponds to a break of connect between a pair of particles. Accordingly, in PFC 2D, the number of fractures can be counted automatically.

Figure 6:

Further, a proportion k is defined as the threshold of the proposed criteria. If the proportion of fracture located in specific area in the total fractures number are higher than k, the model will be judged as the corresponding classification. The criteria area corresponding to each fracture morphology classification is shown as Figure 7.

Figure 7:

In Figure 7, w and h represent the width and height of the rock numerical model used in the uniaxial compression test, and d represent the diameter of the model used in the Brazilian splitting test, and l represent the length of the model used in the direct shear test. After simulation of the testing, if more than k of the total fractures is located in a criteria area, the model is judged as the corresponding fracture morphology. Take the localized fracture morphology as an example, assuming more than k fractures are located in the criteria area as Figure 7 (a) shows, also the four corners of the numerical model, the model is judged to generate localized fractures morphology after the uniaxial compression test. Specially, due to the large-scale fracture morphologies generated after the uniaxial compression test, the Brazilian splitting test, and the direct shear test are complex, and the generated fractures might be theoretically located in anywhere inside the numerical model, theoretically, there are no criteria area of the three kinds of morphologies. As alternatives, if there is more than (1-k) of the total generated fractures exceeding the criteria area shown in Figure 7, the model will be judged as the large-scale fracture morphology. For example, if less than k of the total fractures of a model generated by the uniaxial compression test are located in the area in both Figure 7(a) and (b), it is judged as the large-scale fracture morphology in uniaxial compressive test. Similarly, if less than k of the fractures generated in the Brazilian splitting test are located in Figure (c) and (d), it is also be judged as the large-scale fracture morphology. In addition, there are three kind of fracture morphologies of each rock numerical model established with the same micro-parameter combination, which is generated by the three tests, respectively, and the three morphologies are always independent.

4.1. Distribution rules of the micro-parameters

Using the micro-parameter combinations listed in Table 1, totally 1517 rock numerical samples are established, and the corresponding tests are simulated. Take one of the numerical models as an example, with the micro-parameters are listed in Table 2, the simulated results are shown as Figure 8.

Figure 8:

To determine the fracture morphology of each numerical model clearly, the criteria k should be selected as a definite value. In this research, it is selected as 0.8. As Figure 8(a) shows, in the uniaxial compression test, there are totally 689 features generated. Among them, 632 features are located in the area shown in Figure 8(a), which accounts about 0.92. Therefore, the fracture morphology of the model generated in uniaxial compression test is recorded as the localized one. Similarly, there are totally 245 and 589 fractures generated in the Brazilian splitting test and the direct shear test, and 207 and 494 fractures are located in the area shown in Figure 8 (b) and (c), and accounts 0.85 and 0.84. Accordingly, it is judged as dumbbell and linearly penetrated fracture morphology. By analyzing the number and position of the generated fractures, the fracture morphologies are recorded. The number and the proportion of each morphology is listed in Table 3.

According to the statistic results, the declining morphology and the linear-penetrated morphology is the most frequently appeared morphology in uniaxial compression test, Brazilian splitting test, and the direct shear test, and the proportion of them is about 60.3%, 74.0%, and 84.8%, respectively. It is also consistent with regular test results. On this basis, the distribution rule of the micro-parameters of each morphology is further analyzed.

Figure 9:

In Figure 9, the longitudinal coordinates represent the micro-parameters used by the numerical model, the horizontal coordinates represent the corresponding fracture morphologies. The color of the cubes represents the proportion of the samples number with corresponding ranges micro-parameters to the number of samples with the same morphology. In other words, the sum of the value corresponding to cubes in one column equals 1. As shown in Figure 9, it is found that the distribution of most fracture morphologies in the micro-parameter ranges are different. Take the tensile strength as an example (Figure 9 (a)), the darkest red cube shows that the tensile strength between particles of the most of localized fracture morphologies models are [10, 20] MPa, and the corresponding proportion is higher than 70%. Meanwhile, there is no localized fracture morphology generated with tensile strength between particles higher than 30 MPa. Similarly, the large-scale fracture morphology generated in the Brazilian splitting test is mainly occurred with low tensile strength between particles, mainly in the range of [10, 30] MPa. Compared with above two kinds of morphology, the distribution of the other morphologies in different tensile strength is relatively uniform, such as the linearly penetrated morphology generated in the Brazilian splitting test, as the corresponding frequency is close. It indicates that, in terms of the linearly morphology generated in the Brazilian splitting test, tensile strength between particles is not the primary influencing factor. According to the statistic results shown in Figure 9, the mainly influencing factors on each fracture morphology are listed in Table 4.

4.2. Evaluating method of fracture morphology based on K-Nearest Neighbors

Based on the above analysis, there are significant impact of the micro-parameters of the numerical model on the fracture morphology. Accordingly, the fracture morphology of numerical models can be evaluated by their parameters before conducting simulations. During establishing numerical model, the micro-parameters always have multiple solutions, and the different micro-parameters always lead to different simulation results, such as the fracture morphology and failure loading. Therefore, evaluation of the fracture morphology is helpful for improving the accuracy of the numerical model.

Due to the large number of numerical samples and high sample density, a simple classification method should be selected to save calculation time consumption. Therefore, K-Nearest Neighbors is used to evaluate the fracture morphology with their micro-parameters as input. In detail, the seven micro-parameters are nondimensionalized to the range from 0 to 1 by 0–1 normalization (Wang et al. 2023; Meng et al. 2025). On this basis, Euclidean distance is used to evaluate the similarity of two samples, as shown in Eq.1.

• i and j - the order number of two samples,

• D_ij - the distance between i^th and j^th sample,

• k – the order number

• n - the total number of the model micro-parameter.

In this paper, n is selected as 7. In addition, y_ik and y_jk is the k^th micro-parameters of the i^th and j^th samples. The more the D_ij, the larger distance between i^th and j^th sample, the more difference between the two samples. According to the Eq.1, the distance between a new unknown sample and any known sample can be calculated. To evaluate the fracture morphology of a new unknown numerical model, the Euclidean distances between it and all the 1517 known samples can be calculated, and the most frequent fracture morphology appeared in the nearest k samples are taken as the evaluated results of the unknown sample. Especially, if there are multiple kind of fracture morphology with equal occurrence times, the selected sample number k will be increased by 1. In this research, the number k is set as 10.

According to Eq.1, the nearest 10 samples can be distinguished. However, the simulated dataset applied on K-Nearest Neighbors is imbalanced, which may have negative influence on the evaluating results. Therefore, the number of nearest samples is modified by Eq.2 to balance the differences between the number of samples with different fracture morphology.

• $n_i^{\prime }$ and n_i - the number of nearest samples with the i^th fracture morphology,

• w_i - the weight of samples with the i^th fracture morphology,

• N and N_i – the total number of the dataset and sample with the i^th fracture morphology,

• C - the number of the fracture morphologies.

Especially, Eq.2 is independently calculated for each kind of rock test. Take the uniaxial compression test as an example, the N and the C is set as 1517 and 3, and the N_i of Localized, Declining, and Large-scale is set as 163, 914, and 440, respectively. Assuming that there is 2 Localized, 6 Declining, and 2 Large-scale samples is the nearest 10 samples, the weighted sample number $n_i^{\prime }$ of the three kind of sample is about 6.2, 3.3, and 2.3. Therefore, the fracture morphology of the assuming sample is evaluated as Localized instead of Declining. The weighted way can effectively avoid the shortcoming of samples being incorrectly evaluated as the majority fracture morphology due to dataset imbalance.

Take a numerical model with randomly generated micro-parameters as an example, its fracture morphology under uniaxial compression test, the Brazilian splitting test, and the direct shear test is unknown, and its micro-parameters are shown as the 1^st row in Table 5. After normalization, the distances between the example model and the 1517 models are calculated, and partial the micro-parameters of the 10 closest models are also listed in Table 5, and their simulated fracture morphologies are listed in Table 6.

According to the K-nearest neighbors, the fracture morphology of the example model in the uniaxial compression test is evaluated as declining fracture morphology, since there are 9 declining fracture morphology appeared in the 10 closest model. Similarly, the fracture morphology in the Brazilian splitting test and the direct shear test is judged as the Dumbbell and the Linearly penetration one.

To verify the evaluating method, a numerical rock model is established according to the example micro-parameter combination, and the rock tests are conducted, the image of the fracture distribution is shown in Figure 10.

Figure 10:

In the simulation process of the uniaxial compression test, totally 840 fractures are generated, and 728 of them are located in the area shown in Figure 10 (a). Therefore, the fracture morphology is judged as declining one. In the Brazilian splitting test and direct shear test, totally 489 and 702 fractures are generated, and 441 and 552 fractures are located in the area shown in Figure 10 (b) and (c). In terms of the Brazilian splitting test results, the occupation of the generated fractures located in Figure 10 (b) is 0.90, which is higher than the criteria, and the fracture morphology is recorded as dumbbell, which is consistent with the K-nearest neighbors. However, in the direct shear test, the occupation of fractures in Figure (c) is only 0.79, so the fracture morphology is regarded as the large-scale one, which is different with the K-nearest neighbors method. Although the evaluating results differ with the simulation results, the occupation is extremely close to the criteria, which proves that the morphology of the example model has similarity of the linearly penetrated fracture morphology. The effectiveness of the evaluating method should be further verified by more numerical models with different micro-parameter combinations.

5.1. Verification results of the evaluation method

To verify the classification and evaluating method of the fracture morphology of the rock numerical model. Totally 40 combinations of the micro-parameter are uniformly and randomly generated in the range listed in Table 1, partial of which is shown in Table 7. In evaluating process, the k number of K-nearest neighbors is selected as 10, and the partial evaluating results of the 40 numerical model are listed in Table 8.

40 numerical models are established by the combinations of the micro-parameters listed in Table 7, and the simulation of the test is conducted. Among them, partial simulation results are shown in Figure 11–13.

Figure 11:

Figure 12:

Figure 13:

As Figure 11(a) shows, during the uniaxial compression test by the 4^th model, there are totally 1573 fractures generated, and 1290 fractures are located in the orange area. Therefore, the fracture morphology generated in uniaxial compression test of the model is recorded as declining morphology. Similarly, the fracture morphologies of all models are recorded and shown in Figure 11–13. These simulated morphologies are compared with the evaluated results, as Figure 14 shows.

Figure 14:

In Figure 14 (a)–(c), the horizontal and vertical axis represent the evaluated and simulated fracture morphology, and the color of the cubes represent the frequency of each kind of sample. In these figures, the cube located on the diagonal line from the top left to the bottom right represents the frequency of the sample with corrected evaluated results. There are totally 33, 34, and 34 corrected samples appeared in the three kinds of rock test, and the corresponding accurate rate reaches 82.5%, 85%, and 85%. Further, the precision and recall of the totally 8 kinds of fracture morphologies are listed in Table 9.

As Table 9 shows, the evaluated precision and recall is relatively acceptable. Except the Large-scale morphology in Brazilian splitting test and direct shear test, which might be affected by the small sample size, nearly all the precision and recall are higher than 60%. In terms of the main morphology with large sample size, such as the Declining morphology in uniaxial compression test, the Dumbbell morphology in Brazilian splitting test, and the Linearly morphology in direct shear test, the precision and recall are higher than 80%, and the corresponding F1-score is 0.89, 0.91, and 0.91. The average F1-score reaches 0.73. It proves that the evaluation results balance precision and recall, and the evaluation performance of the method is relatively good.

5.2. Comparison between multiple algorithm

On the basis of the proposed constraints of the fracture morphology, a fracture morphology evaluation method is developed by KNN, and the evaluation results are proved to be accurate and acceptable. In order to further verify the reliability of KNN on the dataset, the other classifiers, including support vector machine (SVM) and artificial neural networks (ANN) are tried to solve the same dataset. The two classifiers are achieved in Weka 3.6.1. Toally 1517 samples introduced in Section 2 is used as the training set, and the 40 simulation samples introduced in Section 5.1 is used as the testing set. By error- and trail, a series of hyper-parameters of SVM and ANN is selected, and the testing results is shown as Figure 15, and the corresponding precision and recall are listed in Table 10.

Figure 15:

As Figure 15 shows, the evaluated results by SVM and ANN is also acceptable, as the darkest color appeared on the diagonal line from the top left to the bottom right, which means the correct evaluated results are given by SVM and ANN for most of samples. In terms of the three kinds of tests, there are 32, 33, and 36 correct results are given by SVM, and 32, 34, and 35 correct results are given by ANN. On the whole, the evaluation accuracies of SVM and ANN are closed to the one of KNN.

As Table 10 shows, the average F1-scores of SVM and ANN are 0.71 and 0.70, respectively, which is also closed to the level of KNN. However, there is a significant imbalance of the precision and recall of SVM and ANN. There are 6 and 5 precision or recall values lower than 0.60 appear in the results of SVM and ANN, which is only 3 values of KNN, while the lowest precision or recall of SVM and ANN is 40.0% and 33.0%, which is 50.0% of KNN. According to Figure 15, the evaluated results by SVM and ANN tend to cluster towards the majority of samples, which is also the reason of the imbalance of SVM and ANN results. Therefore, although the accuracy of SVM and ANN keeps similar level of KNN, considering the imbalance of the evaluated results, KNN is the most suitable algorithm to evaluate the generated fracture morphology.

The reason for the results is speculated to be that the SVM and ANN are suitable for solving the strongly nonlinear and complex problem by small size of samples. In this research, the generated fracture morphology is closely related to the micro-parameters, and the generated fracture morphologies of the numerical models constructed by similar micro-parameter combinations are hardly to exhibit significant differences. In addition, there is only two or three groups are set for each kind of rock test in this research, which also results in the nonlinearity and complexity of this problem is slightly weaker than conventional classification problem in the field of engineering geology. In other words, simple and direct classification method, such as KNN, is more suitable for evaluating generated fracture morphology than methods with high complexity.

5.3. Discussion

This section mainly discusses the potential and limitation of this research.

1. Potential

This research is mainly aiming at the fracture morphologies generated in the numerical simulation process of the different rock test, by using numerical models with different micro-parameter combinations. In geotechnical engineering, the macro-parameters, such as the uniaxial compressive strength and tensile strength are used to describe the rock properties. However, in numerical simulation, the rock materials are defined by micro-parameters, part of which are listed in Table 1, instead of the macro-parameters. Therefore, during building numerical models, an accurate relationship between micro-parameters and macro-parameters should be established. In terms of a numerical model with known macro-parameters, the selection of micro-parameters is always not unique, which means that there are usually more than one micro-parameter combinations meet the requirement of the model strength, and they should be screened. Especially, the fracture morphology generated in the rock test is an available filter to suppress the non-uniqueness of the micro-parameter combinations.

Take a field rock sample came from a bridge project in Southwest China as an example, its uniaxial compressive strength tested by the uniaxial compression test is 84.0 MPa, and its broken photo is shown as Figure 16. According to its compressive strength, 10 micro-parameter combinations with the uniaxial compressive strength error lower than 5%, also the uniaxial compressive strength ranges from 79.8MPa to 88.2MPa are randomly selected among the 1517 samples. The errors between the actual strength and the simulated strength of each selected sample are lower than 5%. In terms of the compressive strength, all of the combinations listed in Table 11 meet the requirement. In fact, part of the simulated fracture morphologies are not matched the actual test results.

Figure 16:

Table 11:

The micro-parameters of the selected samples

	σmi [MPa]	cmi [MPa]	λmi	$E_{pmi}^{*}[GPa]$	ϕmi [°]	μmi	pmi	Morphology
1	29.04	75.33	0.97	18.82	51.82	0.44	0.13	Declining
2	61.17	30.16	0.87	77.71	52.14	0.80	0.13	Large-scale
3	63.29	33.45	0.78	60.39	41.95	0.23	0.12	Declining
4	78.89	35.12	0.92	70.39	58.92	0.31	0.12	Declining
5	65.17	28.22	0.92	85.61	48.28	0.29	0.12	Large-scale
6	58.47	67.00	0.53	24.55	50.73	0.40	0.10	Declining
7	52.55	29.80	0.68	26.13	41.89	0.67	0.11	Large-scale
8	26.12	33.79	0.97	73.89	51.69	0.65	0.11	Declining
9	72.77	36.09	0.90	34.68	59.95	0.66	0.10	Declining
10	45.59	29.63	0.94	35.86	54.30	0.78	0.13	Large-scale

According to the proposed classification of the morphologies, the morphology of the actual sample is declining one, and only 1#, 3#, 4#, 6#, 8#, and 9# meet the requirement. The rest of the listed combinations cannot accurately reflect the rock behavior in the test, and they should be excluded from the selection for establishing the numerical model. The fracture morphologies generated in the Brazilian splitting test and the direct shear test can also be used as the basis to reduce the selection area of the micro-parameters to suppress their non-uniqueness.

2. Limitation

This paper investigated the effect of the micro-parameters of the numerical model on the fracture morphologies generated in several common rock test, and method for evaluating the fracture morphology is proposed. However, there is still several limitations of the method.

1. This research is conducted on the basis of simulation results of the totally 1517 numerical models with different micro-parameters. However, the ranges of the micro-parameters are limited by the limited computing resources. In fact, there are significant differences of the micro-parameters of different rock material. In terms of the rock with the micro-parameters outside the ranges shown in Table 1, the generated fracture morphologies may not be accurately evaluated by the proposed method. More simulations with extended micro-parameters ranges should be conducted in the futural studies.

2. In this research, the fractures generated by rock test need to be classified, and it need to determine the spatial distribution area of each kind of morphology and a criteria coefficient k. In this paper, the area is given as Figure 7, and the coefficient is selected as 0.8. However, these selections are subjective and empirical and can be used as a reference. During the actual simulation process, the area and coefficient should be further adjusted.

3. According to the previous research and the simulation process, there are multiple factors beyond micro-parameter affecting the generated fracture morphology. In fact, the stoppage time node of the simulation has significant influence of the fracture morphology. Take one of the 1517 models with the uniaxial compressive strength of 84MPa as an example, the generated fractures in different time is shown as Figure 17.

Figure 17:

As Figure 17 shows, the fracture number increases with the increasing of the iterations. In the 1^st and 2^nd stage corresponding to the Figure 17 (b) and (c), the generated fractures are not penetrated the rock sample, and there are no fragments peeled off the model. With the conducting of iteration, in the 3^rd stage corresponding to the Figure 17(d), numbers of fractures are generated in the lower left and lower right corner, and some particles are peeled off the model. Correspondingly, there is a significant fluctuation appeared on the load curve near the 3^rd stage, shown in Figure 17(a). Also as shown in Figure 17 (d), if the simulation is ended near this stage, the morphologies of this sample might be judged as the localized one. In fact, the load is continuously increasing, and numbers of fractures appeared in the center area of the model in the 4^th stage, as shown in Figure 17 (e). Accordingly, the morphology of this model in the uniaxial compression test are recorded as the declining one.

The example indicates that with the different ending criterion, the recorded fracture morphologies might be different. Especially, a universal and automatic ending criterion is used in this research, because there are numbers of models are used in simulation. In this research, the axial load is used as the ending criterion, i.e. when the current load decreased to be lower than the maximum of the historical load, the simulation should be ended. This criterion achieved acceptable results in the used 1517 models. However, the criterion is still subjective and empirical, and may not be applied to every rock numerical model. It also results in the fracture morphology results of some original models, which are not matched with the actual results. By conducting more rock test simulation to proposed a more reasonable ending criterion, is one of significant works in the futural research.