# A Novel Method to Evaluate the Effect of Slope Blasting under Impact Loading

Authors: Wenhai Zhou; Zhenbao Li; Xin Wu; Shujiang Wang; Yaru Lv; Xiaoming Lou; Chaofei Du; Rui Liang
Journal: Shock and Vibration
 (2020)

---

## Abstract

In this paper, the quantitative evaluation of the explosion effect based on the fuzzy comprehensive evaluation method is proposed to describe the qualitative evaluation results. The selected state characteristic parameters are expressed by two kinds of membership functions, fuzzy normal and triangular distribution membership functions, and preliminary evaluation results are obtained. The validity index of the maximum membership principle is used to assess the accuracy of the evaluation results of two algorithms, and a relevant approaching degree is chosen to optimize the results. The entire evaluation process selects eleven indicators to form an evaluation set, including the boulder yield, root rate, flying distance of flyrock, explosive consumption, postcracking distance, detonator unit consumption, vibration velocity, loose coefficient, cast distance, throw rate, and blasted volume per meter of hole. Part of the indicator parameters are derived from field test monitoring, and another part of the indicator parameters are derived from numerical simulation. The simulation process uses the user-defined material interface function provided by LS-DYNA. And the numerical model of slope blasting is established by embedding the evolution relationship of tensile and compressive damage into the elastoplastic constitutive material. The evaluation method proposed in this paper is used to evaluate the postexplosion effect of Zijin Mountain gold-copper mine slope cast blasting. The results demonstrate that the fuzzy normal distribution membership function can correlate the state characteristic information and evaluation index effectively, and the working condition after explosion can be reflected accurately. Additionally, the influencing factors can be ranked by the importance degrees according to the calculated value of the evaluation index.


---

## Body

## 1. Introduction

Blasting technology is widely used in infrastructure construction, such as mine engineering, tunnel and underground engineering, municipal engineering, civil air defense and national defense engineering, water conservancy and hydropower engineering, and nuclear power construction engineering. In the last decade, most studies focus on optimization of blasting process parameters. Kuzu et al. [1, 2] highlighted that reasonable blasting design parameters and blasting effect can greatly improve production efficiency and achieve significant economic and social benefits. Wang et al. [3] identified the main parameters affecting the vibration spectrum based on a large number of explosion signals from a series of events to detect differences in continuous parameters. They also discussed how the parameters affect the occurrence and evolution of vibration signals. In recent researches, Xu et al. [4] investigated the dynamic response of deep rock masses during excavation of a stable blast. A new safety standard for blasting vibration of deep rock masses was proposed through field testing and numerical simulation and the failure mechanism of deep rock masses. Yan et al. [5] aimed to suggest a revised criterion to evaluate human response to the vibration from blasting excavation of a large-scale rock slope. Chen et al. [6] evaluated the stability of the slope by comprehensively analyzing the distribution of failure factors and vibration velocity characteristics. The parameters were analyzed by controlling the blasting load to optimize the blasting design. Albert et al. [7] measured results of seismic signals generated by airborne near-ground explosions of various types of explosives on the ground which showed that there are two different types of ground vibrations arriving, the earlier arrival time is mainly the underground path, and the explosion sound is close to the source. Through research and analysis, An et al. [8] concluded that there are more than 20 types of indices that affect rock-soil blasting and mainly focused on rock properties, explosive characteristics, and blasting parameters. However, because of designers’ different understanding of the blasting principle and design method, different design results are often obtained, and the blasting effect is clearly different. The increasing maturity of computer technology and the improvement of the evaluation algorithm make it possible to solve these problems. By establishing a mathematical model, not only can the required blasting parameters be obtained but these parameters can also be optimized and the postblasting effect can be evaluated [9, 10]. From the comparison of the development and application of BP network and fuzzy mathematics in risk prediction and assessment, Jong and Lee [11] earlier applied neural network (NN) models to training and testing models on the optimal NN model using collected data based on a series of observations and numerical experiments. The results showed that the neural network model can predict the powder coefficient. Xue [12] established a neural network model in which 14 geological parameters were selected as the network input, and the unit consumption of highway tunnel blasting explosives was calculated and the blasting vibration velocity was predicted. Ding et al. [13] developed digital three-dimensional (3D) processing software for cast blasting, which can be used to design the blasting process and prediction analysis of the postblasting effect. Dadashzadeh et al. [14] proposed an expert system to design presplit blasting in which parameters are calculated based on a three-layer BP (backpropagation) neural network. Additionally, Monjezi et al. [15], Li et al. [16], and Lei et al. [17] used gray correlation and a fuzzy mathematical algorithm to evaluate the blasting effect of step casting.The postblasting effect of slope blasting is jointly determined by a number of evaluation indicators. The dimensions of each index are different, and the importance of the evaluation results is different. Therefore, the evaluation of postblasting effects is a vague problem. Fuzzy comprehensive evaluation is a method based on fuzzy mathematics, which measures factors with unclear boundaries and difficult to quantify. It can evaluate the membership status of the evaluated object from multiple factors. This method can reflect the importance of different evaluation indexes to the evaluation results through the setting of weights and can solve these types of problems similar to the postexplosion effect evaluation.Based on the fuzzy comprehensive evaluation method for the step-cast blasting effect, in this paper, we mainly discuss the difference between the blasting effect evaluation results when membership functions for the factor evaluation set are selected for Gaussian and triangle distributions, and then relevant nearness analysis is chosen to optimize the results. We propose that when nonlinear factors, such as the impact of the blasting effect, are selected to measure the evaluation results of the target layer, the related membership function should be adopted to optimize the evaluation result so that the fuzzy probability obtained for the step-cast blasting effect can reflect engineering practice more accurately.


## 2. Fuzzy Comprehensive Evaluation of Cast Blasting

### 2.1. Fundamentals

(1)
Establishing a set of factors for object of evaluationU: for the subject of evaluation, the evaluation indicators are divided into multiple subsets U = {U1, U2, U3, … , Um}, where Ui is the factor that influences the object of evaluation.
(2)
Establishing a set of evaluationsF: the evaluation set is a collection of the results for the object of evaluation, and F = {f1, f2, … , fn} has n evaluations of the state of each factor.
(3)
Constructing an evaluation matrix: each factor ofUi (i = 1, 2, … , m) is evaluated to obtain the membership degree rij of Ui to fj (j = 1, 2, … , n), and then the single factor evaluation set ri = (ri1, ri2, … , rin) of Ui is obtained. The evaluation matrix that consists of m evaluated factors can be represented as follows:
(1)
R
=
r
i
j
m
×
n
=
r
11
r
12
⋯
r
1
n
r
21
r
22
⋯
r
2
n
⋮
⋮
⋮
⋮
r
m
1
r
m
2
⋯
r
m
n
.
(4)
Comprehensive evaluation: the weights of each indicator are combined with the fuzzy relation matrix to determine the comprehensive evaluation setD = [d1, d2, … , dn], and the final evaluation result is determined by comparing d1, d2, … , dn:
(2)
D
=
W
⋅
R
=
ω
1
ω
2
⋯
ω
m
r
11
r
12
⋯
r
1
n
r
21
r
22
⋯
r
2
n
⋮
⋮
⋮
⋮
r
m
1
r
m
2
⋯
r
m
n
=
d
1
d
2
⋯
d
n
,where W, ωi, and “·” represent the comprehensive weight vector and the fuzzy operator, respectively. Different operator definitions have different fuzzy comprehensive evaluation models. The main factor determinant model is M∧,∨, which corresponds to dj=Vi=1mωi∧rij, and the weighted average model is M⋅,∨, which corresponds to dj=min1,∑i=1mωirij.Typically, the main factor deterministic definition of the Zadeh operator is defined by the membership function, which is obtained by taking the larger, smaller, and remainder points in sequence. Therefore, in the specific process, the result of the operation often deviates from the nature of the problem because too much information is missing. It is better to use triangular modulus operators to reflect the nature of “And logic” and “Or logic.” The influence of the main information, in addition to the impact of nonprimary information on the evaluation results, can be highlighted. Triangular modular operators include the drastic operator, Lukasiewicz operator, Einstein operator, Hamacher operator, Dubois−Prade operator, Yager operator, and Dombi operator, which are commonly used [18, 19]. The Dombi operator is used in the main factor determinant M∧,∨ to define the fuzzy operation, and the specific expressions of “And logic ∧” and “Or logic ∨” are represented as follows:(3)Tdopωi,rij=0,ωi=0orrij=0,11+1/ωi−1ξ+1/rij−1ξ1/ξ,others,(4)Sdosωi,rij=0,ωi=rij=0,1,rij=0orωi=0,1+11/ωi−1−ξ+1/rij−1−ξ−1/ξ,others,where ξ is a parameter, which requires ξ>0, and the value of ξ=1.


### 2.2. Establishing the Membership Degree Function

Typically, the form used to define the membership function can be divided into two types: quantitative and qualitative [20, 21]. According to research and membership function principles, the Gaussian distribution quantitative method and triangle quantitative method have been used to study the effect of cast blasting separately. According to the nature of the selected factors, the factor set can be divided into two categories: smaller is better and larger is better.A preliminary comparison of the Gaussian distribution and triangular function forms has demonstrated that when determining the interval of the parameter value, the Gaussian distribution function can better circumvent information with a low membership degree when the membership degree is small, and the Gaussian distribution function can better collect information with a high membership degree when the membership degree is large, as shown in Figure1.Figure 1
Comparison of membership in normal and triangular distributions.


#### 2.2.1. Membership Degree Function of the Fuzzy Gaussian Distribution

The evaluation units use a seminormal distribution near the extreme value when the Gaussian distribution function is selected for each parameter membership degree. The function can be expressed as follows:(5)γu=e−u−μ2/2σ2,where γu is the value of the membership degree of the evaluated element u, μ is the expected value, and σ represents the width of the membership function.Three comment sets are selected to evaluate the elements; the corresponding membership distribution functions areγnoru, γattu, and γabnu; and the three interval distribution expectation values μ are in the ascending order: μmin, μ0, and μmax, where μmin and μmax are the minimum and maximum values of the evaluated elements and μ0 is the given expectation of the attention state (common). As is known from the characteristics of the Gaussian function, the area of 99.73% under the function curve is within the range of 3 standard deviations (3σ) around the expected value μ. Therefore, 6σ is chosen as the domain of the function in this paper, which is commonly used. In light of the randomness, fuzziness, and interval features of values of the geotechnical parameter, Zhang et al. [22] constructed the triangle and normal distribution membership functions of geotechnical parameters to evaluate the slope stability, respectively. It is concluded that the evaluation result is the best when the range of geotechnical parameters is ±2.5σ and ±3σ, respectively. Juang et al. [23] and Su et al. [24] analyzed the reliability of slope stability when the geotechnical parameters are 1, 2, and 3 times the standard deviation. The results showed that the reliability is the best when 3 times the standard deviation is taken as the limit. The formulas are given as follows, where equations (6)–(8) are of the smaller-better type of the Gaussian function expression and equations (9)–(11) are of the larger-better type of the Gaussian function expression:(6)γnoru=1,u≤μmin,e−u−μmin2/2σnor2,μmin<u≤μmax,σnor=μ0−μmin3,(7)γattu=e−u−μ02/2σatt12,μmin≤u≤μ0,σatt1=μ0−μmin3,e−u−μ02/2σatt22,μ0<u≤μmax,σatt2=μmax−μ03,(8)γabnu=e−umax−μ2/2σabn2,μmin≤u<μmax,1,μ≥μmax,σabn=μmax−μ03,(9)γnoru=0,u≤μmin,e−u−μmin2/2σnor2,μmin<u≤μmax,σnor=μ0−μmin3,(10)γattu=1−e−u−μ02/2σ2att1,μmin≤u≤μ0,σatt1=μ0−μmin3,1−e−u−μ02/2σ2att2,μ0<u≤μmax,σatt2=μmax−μ03,(11)γabnu=1−e−umax−μ2/2σabn2,μmin≤u<μmax,1,μ≥μmax,σabn=μmax−μ03.According to equations (6)–(11), the calculation model of each parameter for the three categories of membership degrees is determined, as shown in Figure 2.Figure 2
Model of membership in the fuzzy normal distribution.


#### 2.2.2. Fuzzy Triangle Membership Degree Function

The specific expressions of the membership functions of the triangle distribution are listed. According to the parameter expectation value and distribution law of each evaluation interval of the fuzzy normal distribution, the formulas for the smaller-better type are equations (12)–(14) and the formulas for the larger-better type are equations (15)–(17):(12)rnorμ=1,μ<μmin,μmax−μμmax−μmin,μmin≤μ≤μmax,0,μ>μmax,(13)rattμ=1,μ<μmin,μ−μminμmax−μmin,μmin<μ<μmax,1−μ−μmaxμmax−μmin,μ>μmax,(14)rabnμ=0,μ<μmin,0,μmin<μ<μmax,μ−μmaxμmax−μmin,μ>μmax,(15)rnorμ=μmin−μμmax−μmin,μ<μmin,0,μmin<μ<μmax,0,x>μmax,(16)rattμ=μμmax−μmin,μ<μmin,μmax−μμmax−μmin,μmin<μ<μmax,0,μ>μmax,(17)rabnμ=0,μ<μmin,μ−μminμmax−μmin,μmin<μ<μmax,1,μ>μmax.According to equations (12)–(17), the calculation model of each parameter for the three categories of membership is determined, as shown in Figure 3.Figure 3
Model of membership in the triangular distribution.The bench cast blasting effect includes many aspects, and it is difficult to take every factor into consideration. Considering the comprehensive evaluation index system of the step-blasting effect in previous research [16, 17, 25], 11 factors are selected as the evaluation indicators in this paper. Then, based on the results of 16 effective blasting operations conducted in [25] and others from 2009 to 2010 in a particular open-pit mine, parameters μi of the evaluation factors of step blasting are determined. The distribution expectations for the three evaluation intervals are shown in Table 1. Among them, the boulder yield, root rate, flying distance of flyrock, explosive consumption, postcracking distance, detonator unit consumption, and vibration velocity are of the smaller-better type. The four factors of the loose coefficient, casting distance, throw rate, and blasted volume per meter of hole are of the larger-better type. According to the expected value in Table 1, the membership calculation model of each parameter can be determined.Table 1
Expected value of the state parameters.
State parameter
Unit
Factor set
Comment set
Better (V1), μmin
Common (V2), α
Poor (V3), μmax
Boulder yield
%
u
1
15
30
60
Root rate
%
u
2
0.1
0.4
0.8
Flying distance of flyrock
m
u
3
50
125
250
Explosive consumption
kg/m3
u
4
0.51
0.74
0.82
Postcracking distance
m
u
5
1.5
1.75
2.5
Detonator unit consumption
One/m3
u
6
0.01
0.04
0.1
Vibration velocity
cm/s
u
7
0.5
1.25
2
Poor (V3)
Common (V2)
Better (V1)
Loose coefficient
%
u
8
1.02
1.15
1.26
Cast distance
m
u
9
197
210
231
Throw rate
%
u
10
28.2
33
38.3
Blasted volume per meter of hole
m3/m
u
11
5
7.5
20


### 2.3. Index Weight Determination

The analytic hierarchy process is used to determine the weight.(1)
Establishing the evaluation matrix: using the expert survey method to analyze the importance of indicators in the state evaluation, the construction judgment matrix is as follows:(18)A=A11A12⋯A1nA21A22⋯A2n⋮⋮⋯⋮Am1Am2⋯Amn,whereAi1i2 is the importance degree of the evaluation index i1 relative to i2, and its value is from 1 to 9.(2)
Calculating the weights of indicators: the weights of the index obtained by the analytic hierarchy process are as follows:(19)ωi1=∏i2=1mAi1i2m∑ii=1m∏i2=1mAi1i2m.(3)
Applying the consistency test to the judgment matrix: the matrix consistency ratioCR is determined as follows:(20)CR=CIRI=λmax−m/m−1RI,whereλmax is the largest eigenvalue of the judgment matrix, λmax=∑i=1mAωi1/nωi1, and RI is the average random consistency indicator. When CR < 0.1, the judgment matrix is considered to have acceptable consistency; otherwise, the judgment matrix needs to be readjusted.


### 2.4. Effectiveness Index of the Maximum Subordination Principle

In the fuzzy comprehensive evaluation, for the comment setB = (b1, b2, … , bn), the relative index of the maximum membership principle of ∑i=1nbi=1 is defined as(21)α=nβ−12γn−1,where β=max1≤i≤nbi and γ=max1≤j≤n,j≠ibj.(1)
Whenα=+∞, the maximum membership principle is fully effective.(2)
When1≤α<+∞, the maximum membership principle is very effective.(3)
When0.5≤α<1, the maximum membership principle is slightly effective.(4)
When0≤α<0.5, the maximum membership principle is the least effective.(5)
Whenα=0, the maximum membership principle has no effect.The validity of theα index enables the explanation of the relative confidence level of the implementation of the maximum membership principle and provides a quantitative description for the rationality of the evaluation results.


### 2.5. Approaching Degree Method

The approaching degree method is divided into symmetric closeness and asymmetric closeness. The former considers the comment level in the same level, whereas the latter considers the comment set in different positions. From the level of the blasting effect, the status of the commentary level is unequal, so in this paper, we choose the asymmetric sticking schedule method.(1)
ClosenessN(Dj, Ej) is the degree of closeness between the judgment object Dj and its ideal subset Ei (0, 0, … , 1, 0, … , 0). Standardize the evaluation object Dj series and rank it as Dj, Dj −1, Dj +1, Dj −2, Dj +2, …. Convert Ei to E1 (E1 = (1, 0, … , 0)) accordingly.(2)
Calculate the closeness as follows:(22)ND,E=1−1nn+1∑k=1nμDPVk−μEPVkk,P=1,2⋯,where μDVk is the membership of the object corresponding to D belonging to Vk and P is a parameter that adjusts the calculation result. P should be determined according to the specific scenario, and it is not easy to take a too large value. Otherwise, when there are more objects to be evaluated, the evaluation result will be detrimental to the grading. In this paper, P=1 (μDVk is the membership of Vk corresponding to D).


## 2.1. Fundamentals

(1)
Establishing a set of factors for object of evaluationU: for the subject of evaluation, the evaluation indicators are divided into multiple subsets U = {U1, U2, U3, … , Um}, where Ui is the factor that influences the object of evaluation.
(2)
Establishing a set of evaluationsF: the evaluation set is a collection of the results for the object of evaluation, and F = {f1, f2, … , fn} has n evaluations of the state of each factor.
(3)
Constructing an evaluation matrix: each factor ofUi (i = 1, 2, … , m) is evaluated to obtain the membership degree rij of Ui to fj (j = 1, 2, … , n), and then the single factor evaluation set ri = (ri1, ri2, … , rin) of Ui is obtained. The evaluation matrix that consists of m evaluated factors can be represented as follows:
(1)
R
=
r
i
j
m
×
n
=
r
11
r
12
⋯
r
1
n
r
21
r
22
⋯
r
2
n
⋮
⋮
⋮
⋮
r
m
1
r
m
2
⋯
r
m
n
.
(4)
Comprehensive evaluation: the weights of each indicator are combined with the fuzzy relation matrix to determine the comprehensive evaluation setD = [d1, d2, … , dn], and the final evaluation result is determined by comparing d1, d2, … , dn:
(2)
D
=
W
⋅
R
=
ω
1
ω
2
⋯
ω
m
r
11
r
12
⋯
r
1
n
r
21
r
22
⋯
r
2
n
⋮
⋮
⋮
⋮
r
m
1
r
m
2
⋯
r
m
n
=
d
1
d
2
⋯
d
n
,where W, ωi, and “·” represent the comprehensive weight vector and the fuzzy operator, respectively. Different operator definitions have different fuzzy comprehensive evaluation models. The main factor determinant model is M∧,∨, which corresponds to dj=Vi=1mωi∧rij, and the weighted average model is M⋅,∨, which corresponds to dj=min1,∑i=1mωirij.Typically, the main factor deterministic definition of the Zadeh operator is defined by the membership function, which is obtained by taking the larger, smaller, and remainder points in sequence. Therefore, in the specific process, the result of the operation often deviates from the nature of the problem because too much information is missing. It is better to use triangular modulus operators to reflect the nature of “And logic” and “Or logic.” The influence of the main information, in addition to the impact of nonprimary information on the evaluation results, can be highlighted. Triangular modular operators include the drastic operator, Lukasiewicz operator, Einstein operator, Hamacher operator, Dubois−Prade operator, Yager operator, and Dombi operator, which are commonly used [18, 19]. The Dombi operator is used in the main factor determinant M∧,∨ to define the fuzzy operation, and the specific expressions of “And logic ∧” and “Or logic ∨” are represented as follows:(3)Tdopωi,rij=0,ωi=0orrij=0,11+1/ωi−1ξ+1/rij−1ξ1/ξ,others,(4)Sdosωi,rij=0,ωi=rij=0,1,rij=0orωi=0,1+11/ωi−1−ξ+1/rij−1−ξ−1/ξ,others,where ξ is a parameter, which requires ξ>0, and the value of ξ=1.


## 2.2. Establishing the Membership Degree Function

Typically, the form used to define the membership function can be divided into two types: quantitative and qualitative [20, 21]. According to research and membership function principles, the Gaussian distribution quantitative method and triangle quantitative method have been used to study the effect of cast blasting separately. According to the nature of the selected factors, the factor set can be divided into two categories: smaller is better and larger is better.A preliminary comparison of the Gaussian distribution and triangular function forms has demonstrated that when determining the interval of the parameter value, the Gaussian distribution function can better circumvent information with a low membership degree when the membership degree is small, and the Gaussian distribution function can better collect information with a high membership degree when the membership degree is large, as shown in Figure1.Figure 1
Comparison of membership in normal and triangular distributions.


### 2.2.1. Membership Degree Function of the Fuzzy Gaussian Distribution

The evaluation units use a seminormal distribution near the extreme value when the Gaussian distribution function is selected for each parameter membership degree. The function can be expressed as follows:(5)γu=e−u−μ2/2σ2,where γu is the value of the membership degree of the evaluated element u, μ is the expected value, and σ represents the width of the membership function.Three comment sets are selected to evaluate the elements; the corresponding membership distribution functions areγnoru, γattu, and γabnu; and the three interval distribution expectation values μ are in the ascending order: μmin, μ0, and μmax, where μmin and μmax are the minimum and maximum values of the evaluated elements and μ0 is the given expectation of the attention state (common). As is known from the characteristics of the Gaussian function, the area of 99.73% under the function curve is within the range of 3 standard deviations (3σ) around the expected value μ. Therefore, 6σ is chosen as the domain of the function in this paper, which is commonly used. In light of the randomness, fuzziness, and interval features of values of the geotechnical parameter, Zhang et al. [22] constructed the triangle and normal distribution membership functions of geotechnical parameters to evaluate the slope stability, respectively. It is concluded that the evaluation result is the best when the range of geotechnical parameters is ±2.5σ and ±3σ, respectively. Juang et al. [23] and Su et al. [24] analyzed the reliability of slope stability when the geotechnical parameters are 1, 2, and 3 times the standard deviation. The results showed that the reliability is the best when 3 times the standard deviation is taken as the limit. The formulas are given as follows, where equations (6)–(8) are of the smaller-better type of the Gaussian function expression and equations (9)–(11) are of the larger-better type of the Gaussian function expression:(6)γnoru=1,u≤μmin,e−u−μmin2/2σnor2,μmin<u≤μmax,σnor=μ0−μmin3,(7)γattu=e−u−μ02/2σatt12,μmin≤u≤μ0,σatt1=μ0−μmin3,e−u−μ02/2σatt22,μ0<u≤μmax,σatt2=μmax−μ03,(8)γabnu=e−umax−μ2/2σabn2,μmin≤u<μmax,1,μ≥μmax,σabn=μmax−μ03,(9)γnoru=0,u≤μmin,e−u−μmin2/2σnor2,μmin<u≤μmax,σnor=μ0−μmin3,(10)γattu=1−e−u−μ02/2σ2att1,μmin≤u≤μ0,σatt1=μ0−μmin3,1−e−u−μ02/2σ2att2,μ0<u≤μmax,σatt2=μmax−μ03,(11)γabnu=1−e−umax−μ2/2σabn2,μmin≤u<μmax,1,μ≥μmax,σabn=μmax−μ03.According to equations (6)–(11), the calculation model of each parameter for the three categories of membership degrees is determined, as shown in Figure 2.Figure 2
Model of membership in the fuzzy normal distribution.


### 2.2.2. Fuzzy Triangle Membership Degree Function

The specific expressions of the membership functions of the triangle distribution are listed. According to the parameter expectation value and distribution law of each evaluation interval of the fuzzy normal distribution, the formulas for the smaller-better type are equations (12)–(14) and the formulas for the larger-better type are equations (15)–(17):(12)rnorμ=1,μ<μmin,μmax−μμmax−μmin,μmin≤μ≤μmax,0,μ>μmax,(13)rattμ=1,μ<μmin,μ−μminμmax−μmin,μmin<μ<μmax,1−μ−μmaxμmax−μmin,μ>μmax,(14)rabnμ=0,μ<μmin,0,μmin<μ<μmax,μ−μmaxμmax−μmin,μ>μmax,(15)rnorμ=μmin−μμmax−μmin,μ<μmin,0,μmin<μ<μmax,0,x>μmax,(16)rattμ=μμmax−μmin,μ<μmin,μmax−μμmax−μmin,μmin<μ<μmax,0,μ>μmax,(17)rabnμ=0,μ<μmin,μ−μminμmax−μmin,μmin<μ<μmax,1,μ>μmax.According to equations (12)–(17), the calculation model of each parameter for the three categories of membership is determined, as shown in Figure 3.Figure 3
Model of membership in the triangular distribution.The bench cast blasting effect includes many aspects, and it is difficult to take every factor into consideration. Considering the comprehensive evaluation index system of the step-blasting effect in previous research [16, 17, 25], 11 factors are selected as the evaluation indicators in this paper. Then, based on the results of 16 effective blasting operations conducted in [25] and others from 2009 to 2010 in a particular open-pit mine, parameters μi of the evaluation factors of step blasting are determined. The distribution expectations for the three evaluation intervals are shown in Table 1. Among them, the boulder yield, root rate, flying distance of flyrock, explosive consumption, postcracking distance, detonator unit consumption, and vibration velocity are of the smaller-better type. The four factors of the loose coefficient, casting distance, throw rate, and blasted volume per meter of hole are of the larger-better type. According to the expected value in Table 1, the membership calculation model of each parameter can be determined.Table 1
Expected value of the state parameters.
State parameter
Unit
Factor set
Comment set
Better (V1), μmin
Common (V2), α
Poor (V3), μmax
Boulder yield
%
u
1
15
30
60
Root rate
%
u
2
0.1
0.4
0.8
Flying distance of flyrock
m
u
3
50
125
250
Explosive consumption
kg/m3
u
4
0.51
0.74
0.82
Postcracking distance
m
u
5
1.5
1.75
2.5
Detonator unit consumption
One/m3
u
6
0.01
0.04
0.1
Vibration velocity
cm/s
u
7
0.5
1.25
2
Poor (V3)
Common (V2)
Better (V1)
Loose coefficient
%
u
8
1.02
1.15
1.26
Cast distance
m
u
9
197
210
231
Throw rate
%
u
10
28.2
33
38.3
Blasted volume per meter of hole
m3/m
u
11
5
7.5
20


## 2.2.1. Membership Degree Function of the Fuzzy Gaussian Distribution

The evaluation units use a seminormal distribution near the extreme value when the Gaussian distribution function is selected for each parameter membership degree. The function can be expressed as follows:(5)γu=e−u−μ2/2σ2,where γu is the value of the membership degree of the evaluated element u, μ is the expected value, and σ represents the width of the membership function.Three comment sets are selected to evaluate the elements; the corresponding membership distribution functions areγnoru, γattu, and γabnu; and the three interval distribution expectation values μ are in the ascending order: μmin, μ0, and μmax, where μmin and μmax are the minimum and maximum values of the evaluated elements and μ0 is the given expectation of the attention state (common). As is known from the characteristics of the Gaussian function, the area of 99.73% under the function curve is within the range of 3 standard deviations (3σ) around the expected value μ. Therefore, 6σ is chosen as the domain of the function in this paper, which is commonly used. In light of the randomness, fuzziness, and interval features of values of the geotechnical parameter, Zhang et al. [22] constructed the triangle and normal distribution membership functions of geotechnical parameters to evaluate the slope stability, respectively. It is concluded that the evaluation result is the best when the range of geotechnical parameters is ±2.5σ and ±3σ, respectively. Juang et al. [23] and Su et al. [24] analyzed the reliability of slope stability when the geotechnical parameters are 1, 2, and 3 times the standard deviation. The results showed that the reliability is the best when 3 times the standard deviation is taken as the limit. The formulas are given as follows, where equations (6)–(8) are of the smaller-better type of the Gaussian function expression and equations (9)–(11) are of the larger-better type of the Gaussian function expression:(6)γnoru=1,u≤μmin,e−u−μmin2/2σnor2,μmin<u≤μmax,σnor=μ0−μmin3,(7)γattu=e−u−μ02/2σatt12,μmin≤u≤μ0,σatt1=μ0−μmin3,e−u−μ02/2σatt22,μ0<u≤μmax,σatt2=μmax−μ03,(8)γabnu=e−umax−μ2/2σabn2,μmin≤u<μmax,1,μ≥μmax,σabn=μmax−μ03,(9)γnoru=0,u≤μmin,e−u−μmin2/2σnor2,μmin<u≤μmax,σnor=μ0−μmin3,(10)γattu=1−e−u−μ02/2σ2att1,μmin≤u≤μ0,σatt1=μ0−μmin3,1−e−u−μ02/2σ2att2,μ0<u≤μmax,σatt2=μmax−μ03,(11)γabnu=1−e−umax−μ2/2σabn2,μmin≤u<μmax,1,μ≥μmax,σabn=μmax−μ03.According to equations (6)–(11), the calculation model of each parameter for the three categories of membership degrees is determined, as shown in Figure 2.Figure 2
Model of membership in the fuzzy normal distribution.


## 2.2.2. Fuzzy Triangle Membership Degree Function

The specific expressions of the membership functions of the triangle distribution are listed. According to the parameter expectation value and distribution law of each evaluation interval of the fuzzy normal distribution, the formulas for the smaller-better type are equations (12)–(14) and the formulas for the larger-better type are equations (15)–(17):(12)rnorμ=1,μ<μmin,μmax−μμmax−μmin,μmin≤μ≤μmax,0,μ>μmax,(13)rattμ=1,μ<μmin,μ−μminμmax−μmin,μmin<μ<μmax,1−μ−μmaxμmax−μmin,μ>μmax,(14)rabnμ=0,μ<μmin,0,μmin<μ<μmax,μ−μmaxμmax−μmin,μ>μmax,(15)rnorμ=μmin−μμmax−μmin,μ<μmin,0,μmin<μ<μmax,0,x>μmax,(16)rattμ=μμmax−μmin,μ<μmin,μmax−μμmax−μmin,μmin<μ<μmax,0,μ>μmax,(17)rabnμ=0,μ<μmin,μ−μminμmax−μmin,μmin<μ<μmax,1,μ>μmax.According to equations (12)–(17), the calculation model of each parameter for the three categories of membership is determined, as shown in Figure 3.Figure 3
Model of membership in the triangular distribution.The bench cast blasting effect includes many aspects, and it is difficult to take every factor into consideration. Considering the comprehensive evaluation index system of the step-blasting effect in previous research [16, 17, 25], 11 factors are selected as the evaluation indicators in this paper. Then, based on the results of 16 effective blasting operations conducted in [25] and others from 2009 to 2010 in a particular open-pit mine, parameters μi of the evaluation factors of step blasting are determined. The distribution expectations for the three evaluation intervals are shown in Table 1. Among them, the boulder yield, root rate, flying distance of flyrock, explosive consumption, postcracking distance, detonator unit consumption, and vibration velocity are of the smaller-better type. The four factors of the loose coefficient, casting distance, throw rate, and blasted volume per meter of hole are of the larger-better type. According to the expected value in Table 1, the membership calculation model of each parameter can be determined.Table 1
Expected value of the state parameters.
State parameter
Unit
Factor set
Comment set
Better (V1), μmin
Common (V2), α
Poor (V3), μmax
Boulder yield
%
u
1
15
30
60
Root rate
%
u
2
0.1
0.4
0.8
Flying distance of flyrock
m
u
3
50
125
250
Explosive consumption
kg/m3
u
4
0.51
0.74
0.82
Postcracking distance
m
u
5
1.5
1.75
2.5
Detonator unit consumption
One/m3
u
6
0.01
0.04
0.1
Vibration velocity
cm/s
u
7
0.5
1.25
2
Poor (V3)
Common (V2)
Better (V1)
Loose coefficient
%
u
8
1.02
1.15
1.26
Cast distance
m
u
9
197
210
231
Throw rate
%
u
10
28.2
33
38.3
Blasted volume per meter of hole
m3/m
u
11
5
7.5
20


## 2.3. Index Weight Determination

The analytic hierarchy process is used to determine the weight.(1)
Establishing the evaluation matrix: using the expert survey method to analyze the importance of indicators in the state evaluation, the construction judgment matrix is as follows:(18)A=A11A12⋯A1nA21A22⋯A2n⋮⋮⋯⋮Am1Am2⋯Amn,whereAi1i2 is the importance degree of the evaluation index i1 relative to i2, and its value is from 1 to 9.(2)
Calculating the weights of indicators: the weights of the index obtained by the analytic hierarchy process are as follows:(19)ωi1=∏i2=1mAi1i2m∑ii=1m∏i2=1mAi1i2m.(3)
Applying the consistency test to the judgment matrix: the matrix consistency ratioCR is determined as follows:(20)CR=CIRI=λmax−m/m−1RI,whereλmax is the largest eigenvalue of the judgment matrix, λmax=∑i=1mAωi1/nωi1, and RI is the average random consistency indicator. When CR < 0.1, the judgment matrix is considered to have acceptable consistency; otherwise, the judgment matrix needs to be readjusted.


## 2.4. Effectiveness Index of the Maximum Subordination Principle

In the fuzzy comprehensive evaluation, for the comment setB = (b1, b2, … , bn), the relative index of the maximum membership principle of ∑i=1nbi=1 is defined as(21)α=nβ−12γn−1,where β=max1≤i≤nbi and γ=max1≤j≤n,j≠ibj.(1)
Whenα=+∞, the maximum membership principle is fully effective.(2)
When1≤α<+∞, the maximum membership principle is very effective.(3)
When0.5≤α<1, the maximum membership principle is slightly effective.(4)
When0≤α<0.5, the maximum membership principle is the least effective.(5)
Whenα=0, the maximum membership principle has no effect.The validity of theα index enables the explanation of the relative confidence level of the implementation of the maximum membership principle and provides a quantitative description for the rationality of the evaluation results.


## 2.5. Approaching Degree Method

The approaching degree method is divided into symmetric closeness and asymmetric closeness. The former considers the comment level in the same level, whereas the latter considers the comment set in different positions. From the level of the blasting effect, the status of the commentary level is unequal, so in this paper, we choose the asymmetric sticking schedule method.(1)
ClosenessN(Dj, Ej) is the degree of closeness between the judgment object Dj and its ideal subset Ei (0, 0, … , 1, 0, … , 0). Standardize the evaluation object Dj series and rank it as Dj, Dj −1, Dj +1, Dj −2, Dj +2, …. Convert Ei to E1 (E1 = (1, 0, … , 0)) accordingly.(2)
Calculate the closeness as follows:(22)ND,E=1−1nn+1∑k=1nμDPVk−μEPVkk,P=1,2⋯,where μDVk is the membership of the object corresponding to D belonging to Vk and P is a parameter that adjusts the calculation result. P should be determined according to the specific scenario, and it is not easy to take a too large value. Otherwise, when there are more objects to be evaluated, the evaluation result will be detrimental to the grading. In this paper, P=1 (μDVk is the membership of Vk corresponding to D).


## 3. Case Analysis

### 3.1. Field Instance

The project of the cast blasting of a particular explosion area in the open pit of the Zijin Mountain gold-copper mine is selected as the evaluation object. The elevation of the slope is 544 m to 784 m, and the relative height difference is 240 m. The height of the explosion step is 12 m, the width of the platform is 12 m, the angle of the slope of the step is 75°, the hole spacing is 6.5 m × 5 m, and the unit consumption of explosives is controlled in the range of 0.61–0.65 kg/m3. A blasting vibrometer is used in this construction, and the vibration test is conducted at a distance of approximately 150 m from the explosion area. The location of the explosion area, the effect after the explosion, and the arrangement of the vibrometer are shown in Figure 4.Figure 4
Explosion site and the effect after blasting. (a) Before the explosion. (b) After the explosion. (c) Mounting vibrometer.
(a)
(b)
(c)The vibration wave velocity propagation law measured at the monitoring point is shown in Figure5, and the vibration-related data measured are shown in Table 2. The three-axis direction is shown in Figure 6.Figure 5
Vibration velocity of the monitoring point by the field test. (a)X-direction vibration of the monitoring point by the field test. (b) Y-direction vibration of the monitoring point by the field test. (c) Z-direction vibration of the monitoring point by the field test.
(a)
(b)
(c)Table 2
Measured vibration-related data.
Channel name
Maximum (cm/s)
Half-wave frequency (Hz)
FFT frequency (Hz)
Range (cm/s)
Sensitivity coefficientV (m/s)
X-direction vibration
1.1150
21.7
16.9
36.193
0.276
Y-direction vibration
1.2861
17.6
20.5
37.133
0.269
Z-direction vibration
0.4754
14.0
20.8
35.714
0.280Figure 6
Illustration of the three-axis direction.Through a field investigation and analysis after blasting, the results demonstrate that the block rate is approximately 42% and the root rate is approximately 0.25%. The blasting rockfill looseness coefficient is approximately 1.2, and the farthest casting distance is up to 215 m. According to Jiang et al. [26], the ground particle vibration velocity is mainly determined by the Y-axis vibration velocity. Then, the vibration velocity of the blasting surface can be conservatively 1.1 cm/s.


### 3.2. Simulation

#### 3.2.1. Setting Up Numerical Model according to Field Size

Based on the stage-cast blasting of an open pit in the Zijin Mountain gold-copper mine, the similarity 3D solid model of single hole detonation is established proportionally. Additionally, the values of the blasting process parameters and material parameters are in accordance with the actual scenario on-site. The overall specifications of the model are determined as follows by referring to the research experience by Ma et al. [27]: The measurements of the platform are a width of 12 m, step slope of 75°, step height of 12 m, and lower step height of 8 m. To save computing time and space, and dynamically display the crack propagation and throwing tendency of the rock mass during blasting, the model thickness is defined as one unit thickness. For the rock material, a typical elastoplastic material is selected as the research object; that is, the MAT-PLASTIC-KINEMATIC model is adopted. The main physical and mechanical parameters of the rock are shown in Table 3, where ρ is the density, E is Young’s modulus, v is Poisson’s ratio, σ0 is the yield stress, G is the shear modulus, K is the bulk modulus, β is the enhancement factor, and KIC is the fracture toughness of the medium.Table 3
Physical and mechanical parameters of rock.
ρ (g/cm3)
E (GPa)
v
σ
0 (MPa)
β
G (GPa)
K(GPa)
K
IC (105N⋅m−3/2)
2.7
18.23
0.23
102
1.0
7.41
11.25
5.32
Note:ρ is the density, E is Young’s modulus, v is Poisson’s ratio, σ0 is the yield stress, G is the shear modulus, K is the bulk modulus, β is the enhancement factor, and KIC is the fracture toughness of the medium.For the explosive, the 2# rock emulsion explosive is selected, its material type is HIGH_EXPLOSIVE_BURN, and the state equation of explosives is selected without considering the JWL equation of the explosive product composition. The values of explosive material parameters are selected according to the simulation experience in [28, 29]. The specific parameters are shown in Table 4, where A, B, R1, R2, and ω are material parameters; E0 is the initial specific energy; p is the pressure; and V0 is the initial specific volume of the explosive before detonation.Table 4
Parameters of the explosive.
ρ (g/cm3)
p (GPa)
A (GPa)
B (GPa)
R
1
R
2
ω
E
0 (GPa)
V
0
1.15
4.63
214
0.182
4.15
0.95
0.15
4.192
1.0
Note:A, B, R1, R2, and ω are material parameters; E0 is the initial specific energy; p is the pressure; and V0 is the initial specific volume of the explosive before detonation.The specific form of the JWL equation of the state when the explosive is exploded [28–30] is(23)p=A1−ωRVe−R1V+B1−ωRVe−R2V+ωEV.


#### 3.2.2. Constitutive Relation of Rock Mass Damage under Tension and Compression

(1) Calculation Formula for Tensile Crack Damage. There are various kinds of cracks in the rock mass itself, also known as original damage, which further expands under the action of tensile stress of blasting load. Budiansky and O’connell [31] calculated the elastic modulus of the damage medium based on the “self-consistent method.” The TCK model proposes that cracks in rock mass will be activated under the tensile stress of explosive load and the cracks will squeeze and collide with each other to form interaction force. The effective bulk modulus of a homogeneous rock mass with original damage can be expressed as(24)KK0=1−169f1vCd=1−D,where K and K0, respectively, represent the bulk modulus of the damage medium and initial bulk modulus; f1v=1−v2/1−2v, in which v is Poisson’s ratio; Cd is the medium crack density; and D is the damage scalar.Xie et al. [32] put forward Poisson’s ratio change formula of the rock mass medium under impact load:(25)v=v01−169Cd,where v0 is the initial Poisson’s ratio of the medium.Grady and Kipp [33] proposed that the number of original damage cracks N in a unit medium under impact tensile stress can be accurately calculated by Weibull function:(26)N=ξp3km,where p is the impact stress and ξ and m represent the distribution parameters of Weibull function.Taylor et al. [34] proposed a formula for the radius of fragments in brittle media such as rock mass under dynamic impact loading:(27)r=1220KICρCε˙v,max2/3,where KIC and ρ, respectively, represent the fracture toughness and density of the medium, C is the longitudinal wave velocity, and ε˙v,max represents the volume strain rate at fracture. The functional relationship between the average radius of fragments r and the crack density Cd can be expressed as(28)Cd=Nϖr3,where ϖ is the proportional coefficient.Substituting formulas (26) and (27) into (28) gives(29)Cd=ξp3kmϖ⋅1820KICρCε˙v,max2=52ξϖp3kmKICρCε˙v,max2.Takingk=ξϖ in equation (29) yields(30)Cd=52kP3KmKICρC2ε˙vmax−2.The formula of tensile crack damage function can be obtained from equation (24) as(31)D=169f1v¯Cd.From the generalized Hook’s law, the stress-strain function relationship is as follows:(32)Δσij=KdδijΔεv+2GdΔeij,where Kd and Gd represent the modulus of the damaged medium and Δεv and Δeij are body strain increments and partial strain increments, Δεv=3εm, in which εm is the average normal stress.(2) Constitutive Relationship of Hardening Effect under Compressive Stress. Zhang et al. [35] pointed out in the report that the elastoplastic constitutive relationship of the rock mass medium obeying isotropic hardening under the action of detonation gas compression (yield equation) can be expressed as(33)Φ=σi2−σy2,where σi=3/2sij. According to the generalized Hooke’s law, the tensor form is used to represent the deviatoric stress Sij as(34)sij=σij−13σkkδij.During the explosion process, the functional relationship between the effective plastic strainεeffP of the rock mass and the yield limit σy of the fracture failure under the dynamic load is(35)εeffP=∫0t32ε˙ijPdt,σy=σ0+βEpεeffP,where Ep is the plastic hardening modulus, usually set as 0.1E, β is the hardening parameter, and ε˙ijp=ε˙ij−ε˙ije. In formula (33), when Φ≤0, the condition is elastic or neutral loading; when Φ>0, the condition is plastic hardening.(3) Pull-Pressure Damage Constitutive Model Embedding Simulation Software. The large-scale finite element analysis software LS-DYNA provides a rich material library for instantaneous power and large deformation damage simulation and also provides users with a convenient model import window. In this study, the calculation model of the cracking damage of the slope rock mass with cracks and the constitutive relation of the hardening effect under compressive stress are coupled; then the simulation software is imported for step-blasting calculation, which can simulate the trajectory of the throwing motion and the fragmentation of the rock mass under the impact of the explosion. The specific process is shown in Figure 7.Figure 7
Flowchart of the damage model.


### 3.3. Simulation Result Analysis

#### 3.3.1. Casting State of Motion

The distribution of the casting state of the rock at different times is shown in Figure8. The model initiates the cast velocity from the foot of the slope at 7 ms. As time elapses, the casting phenomenon evolves from the foot of the slope to the free surface. By monitoring the model casting units, we found that the maximum casting velocity appeared at the middle of the blast-hole pack, with a velocity of 21.6 m/s. Additionally, this value is consistent with the conclusion regarding the maximum casting velocity proposed by Huang et al. [36, 37], which is 18–28 m/s. The average scattering distance is approximately 47 m, the postcracking distance is approximately 1.2 m, and the casting rate is approximately 34%.Figure 8
Variety of damage and deformation of the rock at different times.
(a)
(b)
(c)
(d)
(e)
(f)


#### 3.3.2. Rock Mass Fragmentation Analysis

The bulks are mainly distributed at the foot of the slope, the sides of the contact surface between the explosive and the blockage, and the free surface at the top of the step. The bulks at the bottom of the slope mainly depend on the selection of the chassis resistance line value, and the difference in the resistance line value results in different bulk distributions. As shown for rock No. 5 in Figure9, the diameter is approximately 2.3 m. The bulks appear on both sides of the contact surface of the explosive, and the blockage is mainly caused by the different mechanical properties and physical parameters of the material on both sides of the contact surface. The stress wave is reflected and transmitted at the contact surface, and the reflection stretch acts on the penetrating microcrack, which causes the rock mass to break away from the parent rock under tensile stress to form bulks. As shown in rock block Nos. 1 and 2 and Nos. 3 and 4 in Figure 9, the diameter distribution is in the range of 1.6 to 2.7 m. The bulks appear on the free surface at the top of the step because of the transmission of stress waves at the interface of the explosive and blockage to the top free surface. The transmitted wave also forms a reflection stretching effect, which causes the rock mass in the region to break along the internal microfracture to form bulks. As shown in rock block No. 6 in Figure 9, the diameter is approximately 1.7 m.Figure 9
Distribution of rock blocks.


#### 3.3.3. Analysis of Simulated Vibration Wave Velocity

A particle at the surface location of the simulation model is selected as the vibration speed-monitoring station, the installation position of the vibration-measuring instrument in the monitoring point and the field monitoring test is spatially consistent, and the specific change law is shown in Figure10. The simulated three-axis peak vibration velocity of the monitoring points is Vx = 0.68 cm/s, Vy = 0.74 cm/s, and Vz = 0.59 cm/s, respectively. The triaxial velocity variation law, the peak size, and the time point at which the peak occurs are in good agreement with the data obtained by the vibrometer monitoring, It is reasonable to take the surface vibration velocity Vy = 1.1 cm/s, when the effect evaluation is carried out in the text (the three-axis velocity direction is consistent with the data in Figure 6).Figure 10
Vibration velocity of the monitoring point by numerical simulation. (a)X-direction vibration of the monitoring point by numerical simulation. (b) Y-direction vibration of the monitoring point by numerical simulation. (c) Z-direction vibration of the monitoring point by numerical simulation.
(a)
(b)
(c)


### 3.4. Fuzzy Comprehensive Evaluation

#### 3.4.1. The Basic Steps

According to the actual postblasting scenario and a similar simulation test, the 11 evaluation indicators selected are considered sequentially, and the set of evaluation object factors is established as follows:(36)U=42,0.25,47,0.63,1.2,0.07,1.1,1.2,215,34.1,6.3.Step 1 (determining the membership degree of each index). The factor set U of the evaluation objects is substituted into equations (6)–(11), and the evaluation matrix determined by the Gaussian distribution membership function is obtained after calculation:(37)Rifgaussian=00.4870.1980.32510.294100.0560.7380.9260.94600.32500.35700.3250.0110.1770.2250.1760.64500000.3250.1970000.355T.The factor set U of the evaluation objects is substituted into equations (12)–(17), and the evaluation matrix determined by the triangle membership function is obtained after calculation:(38)Riftriangle=0.40.600.78610.6130.33300.60.7500.5290.5840.0870.21400.3870.6670.3250.40.2500.4710.4160.91300000.32500000T.In the factor set, the smaller-better type of evaluation objects u1, u2, u3, u4, u5, u6, and u7 is substituted into equations (6)–(8) and (12)–(14) for calculation, whereas the larger-better type u8, u9, u10, and u11 is substituted into equations (9)–(11) and (15)–(17) for calculation. The comment set constructed is of the decreasing type in this paper, so the evaluation object is calculated in a smaller-better membership function, and the calculated values are sequentially put into the evaluation matrix Rif according to the order of comments. Conversely, the membership values calculated by substituting the evaluation object into a larger-better membership function need to be sorted back first and then put into the evaluation matrix Rif in turn.Step 2 (determining the weight of each index using the AHP). The judgment matrix is constructed as follows using the expert investigation method according to equation (18):(39)R=1561773463215111532131211314161116211413121515156177234211713121711141413151617121171114131214151334124411211214231343111121316121432121113141335125412311124516523411.The weight vector is obtained after normalizing each index according to equation (19):(40)W=0.2300.0400.0290.1840.0210.0250.0920.0680.0470.1130.151.CR is obtained after calculation according to the formula CR = 0.052 < 0.1, which shows that the judgment matrix is consistent and there is no need to revise it.Step 3 (comprehensive evaluation). The evaluation matrix R obtained after calculation and weight vector W obtained after normalization are substituted into equations (2)–(4) for calculation to obtain the evaluation results, which are presented in Table 5.Table 5
Comparison of fuzzy comprehensive evaluations.
Evaluation object
Membership function distribution
Membership degree
Better (V1), bnor
Common (V2), batt
Poor (V3), rabn
Blasting effect
Gaussian distribution
0.227
0.409
0.364
Triangle distribution
0.493
0.507
0According to the maximum membership degree principle, the evaluation results of this blasting result controlled by the Gaussian distribution membership function and triangle membership function are both good. However, from the numerical distribution of each component of the two evaluation results, it is obvious that there are great differences between the two evaluation results, which show that the maximum membership principle is invalid in this case. This is because, in the application process, the maximum membership degree principle loses the information of nonmajor components in the evaluation results, which leads to the obvious unreasonable evaluation results in Table5.


#### 3.4.2. Analysis of the Effectiveness Index of the Maximum Membership Degree Principle

The results obtained by the two methods of evaluation are substituted into equation (21) for calculation, and two effectiveness indices are obtained as follows: αgauss=0.47 and αtriangle=0.46. The range of the confidence level judged according to the implementation of the principle of maximum membership indicates that the principle of maximum membership is inefficient in this case. Therefore, using the closeness analysis method is considered to process the evaluation results to determine the specific level of the evaluation object. The nearness method is considered to process evaluation results to determine the specific grade of the evaluation objects.


#### 3.4.3. Applying the Nearness Method to Analyze the Evaluation Results

The results of two types of evaluation are input into equation (22) for calculation. Two evaluation result vectors controlled by the Gaussian distribution and triangle membership function are, respectively, denoted as D1 and D2, and the results of the nearness calculation are as follows:(41)ND1,E1=ND11,E1=1−11100.773×1+0.409×2+0.364×3=0.976,ND1,E2=ND12,E1=1−11100.591×1+0.227×2+0.364×3=0.981,ND1,E3=ND13,E1=1−11100.636×1+0.409×2+0.227×3=0.981,ND2,E1=ND21,E1=1−11100.507×1+0.507×2+0×3=0.986,ND2,E2=ND22,E1=1−11100.493×1+0.493×2+0×3=0.987,ND2,E3=ND23,E1=1−11101×1+0.507×2+0.493×3=0.967.After processing the fuzzy evaluation results controlled by two different membership function forms using the nearness degree algorithm, it is found that when the factor set defined by the Gaussian-type membership function is considered the value of the membership degree in the comment set, the degree of the blasting effect evaluation which belongs to the common comment set is consistent with that which belongs to the poor comment set; that is, the final evaluation result of the blasting effect using the Gaussian-type membership function is between the “common” and “poor” grades. Similarly, the degree of the blasting effect evaluation that belongs to “better” is essentially consistent with that which belongs to “common” when using the triangle membership function, which demonstrates that the final evaluation result of the blasting effect using the triangle membership function is between the “better” and “common” grades.Through further observation of the blasting site and postblasting effect, we found that the indices of explosive consumption and detonator consumption in the blasting construction process are slightly beyond the reasonable range; the longest cast distance and blasting output per meter hole are also slightly less than the reasonable values given in [16, 17, 25], and there are obvious toe rocks and bulks; loose rock poured into goaf accounts for approximately 20% of the stripped rock mass, and the accumulation of the big pile is not ideal. Thus, the result of the subjective blasting effect evaluation should be consistent with the final result of the fuzzy comprehensive evaluation under the control of the Gaussian membership function after processing the nearness algorithm. The evaluation results controlled by two different membership functions can all be considered common by the nearness algorithm, which is because the number of comment sets selected is small (three comment sets are selected, which are better, common, and poor, in this paper). If the evaluation set in this paper extends to five or more items, then the results of the fuzzy comprehensive evaluation processed by the nearness algorithm, and the final evaluation results controlled by two different membership functions, will be different.


## 3.1. Field Instance

The project of the cast blasting of a particular explosion area in the open pit of the Zijin Mountain gold-copper mine is selected as the evaluation object. The elevation of the slope is 544 m to 784 m, and the relative height difference is 240 m. The height of the explosion step is 12 m, the width of the platform is 12 m, the angle of the slope of the step is 75°, the hole spacing is 6.5 m × 5 m, and the unit consumption of explosives is controlled in the range of 0.61–0.65 kg/m3. A blasting vibrometer is used in this construction, and the vibration test is conducted at a distance of approximately 150 m from the explosion area. The location of the explosion area, the effect after the explosion, and the arrangement of the vibrometer are shown in Figure 4.Figure 4
Explosion site and the effect after blasting. (a) Before the explosion. (b) After the explosion. (c) Mounting vibrometer.
(a)
(b)
(c)The vibration wave velocity propagation law measured at the monitoring point is shown in Figure5, and the vibration-related data measured are shown in Table 2. The three-axis direction is shown in Figure 6.Figure 5
Vibration velocity of the monitoring point by the field test. (a)X-direction vibration of the monitoring point by the field test. (b) Y-direction vibration of the monitoring point by the field test. (c) Z-direction vibration of the monitoring point by the field test.
(a)
(b)
(c)Table 2
Measured vibration-related data.
Channel name
Maximum (cm/s)
Half-wave frequency (Hz)
FFT frequency (Hz)
Range (cm/s)
Sensitivity coefficientV (m/s)
X-direction vibration
1.1150
21.7
16.9
36.193
0.276
Y-direction vibration
1.2861
17.6
20.5
37.133
0.269
Z-direction vibration
0.4754
14.0
20.8
35.714
0.280Figure 6
Illustration of the three-axis direction.Through a field investigation and analysis after blasting, the results demonstrate that the block rate is approximately 42% and the root rate is approximately 0.25%. The blasting rockfill looseness coefficient is approximately 1.2, and the farthest casting distance is up to 215 m. According to Jiang et al. [26], the ground particle vibration velocity is mainly determined by the Y-axis vibration velocity. Then, the vibration velocity of the blasting surface can be conservatively 1.1 cm/s.


## 3.2. Simulation

### 3.2.1. Setting Up Numerical Model according to Field Size

Based on the stage-cast blasting of an open pit in the Zijin Mountain gold-copper mine, the similarity 3D solid model of single hole detonation is established proportionally. Additionally, the values of the blasting process parameters and material parameters are in accordance with the actual scenario on-site. The overall specifications of the model are determined as follows by referring to the research experience by Ma et al. [27]: The measurements of the platform are a width of 12 m, step slope of 75°, step height of 12 m, and lower step height of 8 m. To save computing time and space, and dynamically display the crack propagation and throwing tendency of the rock mass during blasting, the model thickness is defined as one unit thickness. For the rock material, a typical elastoplastic material is selected as the research object; that is, the MAT-PLASTIC-KINEMATIC model is adopted. The main physical and mechanical parameters of the rock are shown in Table 3, where ρ is the density, E is Young’s modulus, v is Poisson’s ratio, σ0 is the yield stress, G is the shear modulus, K is the bulk modulus, β is the enhancement factor, and KIC is the fracture toughness of the medium.Table 3
Physical and mechanical parameters of rock.
ρ (g/cm3)
E (GPa)
v
σ
0 (MPa)
β
G (GPa)
K(GPa)
K
IC (105N⋅m−3/2)
2.7
18.23
0.23
102
1.0
7.41
11.25
5.32
Note:ρ is the density, E is Young’s modulus, v is Poisson’s ratio, σ0 is the yield stress, G is the shear modulus, K is the bulk modulus, β is the enhancement factor, and KIC is the fracture toughness of the medium.For the explosive, the 2# rock emulsion explosive is selected, its material type is HIGH_EXPLOSIVE_BURN, and the state equation of explosives is selected without considering the JWL equation of the explosive product composition. The values of explosive material parameters are selected according to the simulation experience in [28, 29]. The specific parameters are shown in Table 4, where A, B, R1, R2, and ω are material parameters; E0 is the initial specific energy; p is the pressure; and V0 is the initial specific volume of the explosive before detonation.Table 4
Parameters of the explosive.
ρ (g/cm3)
p (GPa)
A (GPa)
B (GPa)
R
1
R
2
ω
E
0 (GPa)
V
0
1.15
4.63
214
0.182
4.15
0.95
0.15
4.192
1.0
Note:A, B, R1, R2, and ω are material parameters; E0 is the initial specific energy; p is the pressure; and V0 is the initial specific volume of the explosive before detonation.The specific form of the JWL equation of the state when the explosive is exploded [28–30] is(23)p=A1−ωRVe−R1V+B1−ωRVe−R2V+ωEV.


### 3.2.2. Constitutive Relation of Rock Mass Damage under Tension and Compression

(1) Calculation Formula for Tensile Crack Damage. There are various kinds of cracks in the rock mass itself, also known as original damage, which further expands under the action of tensile stress of blasting load. Budiansky and O’connell [31] calculated the elastic modulus of the damage medium based on the “self-consistent method.” The TCK model proposes that cracks in rock mass will be activated under the tensile stress of explosive load and the cracks will squeeze and collide with each other to form interaction force. The effective bulk modulus of a homogeneous rock mass with original damage can be expressed as(24)KK0=1−169f1vCd=1−D,where K and K0, respectively, represent the bulk modulus of the damage medium and initial bulk modulus; f1v=1−v2/1−2v, in which v is Poisson’s ratio; Cd is the medium crack density; and D is the damage scalar.Xie et al. [32] put forward Poisson’s ratio change formula of the rock mass medium under impact load:(25)v=v01−169Cd,where v0 is the initial Poisson’s ratio of the medium.Grady and Kipp [33] proposed that the number of original damage cracks N in a unit medium under impact tensile stress can be accurately calculated by Weibull function:(26)N=ξp3km,where p is the impact stress and ξ and m represent the distribution parameters of Weibull function.Taylor et al. [34] proposed a formula for the radius of fragments in brittle media such as rock mass under dynamic impact loading:(27)r=1220KICρCε˙v,max2/3,where KIC and ρ, respectively, represent the fracture toughness and density of the medium, C is the longitudinal wave velocity, and ε˙v,max represents the volume strain rate at fracture. The functional relationship between the average radius of fragments r and the crack density Cd can be expressed as(28)Cd=Nϖr3,where ϖ is the proportional coefficient.Substituting formulas (26) and (27) into (28) gives(29)Cd=ξp3kmϖ⋅1820KICρCε˙v,max2=52ξϖp3kmKICρCε˙v,max2.Takingk=ξϖ in equation (29) yields(30)Cd=52kP3KmKICρC2ε˙vmax−2.The formula of tensile crack damage function can be obtained from equation (24) as(31)D=169f1v¯Cd.From the generalized Hook’s law, the stress-strain function relationship is as follows:(32)Δσij=KdδijΔεv+2GdΔeij,where Kd and Gd represent the modulus of the damaged medium and Δεv and Δeij are body strain increments and partial strain increments, Δεv=3εm, in which εm is the average normal stress.(2) Constitutive Relationship of Hardening Effect under Compressive Stress. Zhang et al. [35] pointed out in the report that the elastoplastic constitutive relationship of the rock mass medium obeying isotropic hardening under the action of detonation gas compression (yield equation) can be expressed as(33)Φ=σi2−σy2,where σi=3/2sij. According to the generalized Hooke’s law, the tensor form is used to represent the deviatoric stress Sij as(34)sij=σij−13σkkδij.During the explosion process, the functional relationship between the effective plastic strainεeffP of the rock mass and the yield limit σy of the fracture failure under the dynamic load is(35)εeffP=∫0t32ε˙ijPdt,σy=σ0+βEpεeffP,where Ep is the plastic hardening modulus, usually set as 0.1E, β is the hardening parameter, and ε˙ijp=ε˙ij−ε˙ije. In formula (33), when Φ≤0, the condition is elastic or neutral loading; when Φ>0, the condition is plastic hardening.(3) Pull-Pressure Damage Constitutive Model Embedding Simulation Software. The large-scale finite element analysis software LS-DYNA provides a rich material library for instantaneous power and large deformation damage simulation and also provides users with a convenient model import window. In this study, the calculation model of the cracking damage of the slope rock mass with cracks and the constitutive relation of the hardening effect under compressive stress are coupled; then the simulation software is imported for step-blasting calculation, which can simulate the trajectory of the throwing motion and the fragmentation of the rock mass under the impact of the explosion. The specific process is shown in Figure 7.Figure 7
Flowchart of the damage model.


## 3.2.1. Setting Up Numerical Model according to Field Size

Based on the stage-cast blasting of an open pit in the Zijin Mountain gold-copper mine, the similarity 3D solid model of single hole detonation is established proportionally. Additionally, the values of the blasting process parameters and material parameters are in accordance with the actual scenario on-site. The overall specifications of the model are determined as follows by referring to the research experience by Ma et al. [27]: The measurements of the platform are a width of 12 m, step slope of 75°, step height of 12 m, and lower step height of 8 m. To save computing time and space, and dynamically display the crack propagation and throwing tendency of the rock mass during blasting, the model thickness is defined as one unit thickness. For the rock material, a typical elastoplastic material is selected as the research object; that is, the MAT-PLASTIC-KINEMATIC model is adopted. The main physical and mechanical parameters of the rock are shown in Table 3, where ρ is the density, E is Young’s modulus, v is Poisson’s ratio, σ0 is the yield stress, G is the shear modulus, K is the bulk modulus, β is the enhancement factor, and KIC is the fracture toughness of the medium.Table 3
Physical and mechanical parameters of rock.
ρ (g/cm3)
E (GPa)
v
σ
0 (MPa)
β
G (GPa)
K(GPa)
K
IC (105N⋅m−3/2)
2.7
18.23
0.23
102
1.0
7.41
11.25
5.32
Note:ρ is the density, E is Young’s modulus, v is Poisson’s ratio, σ0 is the yield stress, G is the shear modulus, K is the bulk modulus, β is the enhancement factor, and KIC is the fracture toughness of the medium.For the explosive, the 2# rock emulsion explosive is selected, its material type is HIGH_EXPLOSIVE_BURN, and the state equation of explosives is selected without considering the JWL equation of the explosive product composition. The values of explosive material parameters are selected according to the simulation experience in [28, 29]. The specific parameters are shown in Table 4, where A, B, R1, R2, and ω are material parameters; E0 is the initial specific energy; p is the pressure; and V0 is the initial specific volume of the explosive before detonation.Table 4
Parameters of the explosive.
ρ (g/cm3)
p (GPa)
A (GPa)
B (GPa)
R
1
R
2
ω
E
0 (GPa)
V
0
1.15
4.63
214
0.182
4.15
0.95
0.15
4.192
1.0
Note:A, B, R1, R2, and ω are material parameters; E0 is the initial specific energy; p is the pressure; and V0 is the initial specific volume of the explosive before detonation.The specific form of the JWL equation of the state when the explosive is exploded [28–30] is(23)p=A1−ωRVe−R1V+B1−ωRVe−R2V+ωEV.


## 3.2.2. Constitutive Relation of Rock Mass Damage under Tension and Compression

(1) Calculation Formula for Tensile Crack Damage. There are various kinds of cracks in the rock mass itself, also known as original damage, which further expands under the action of tensile stress of blasting load. Budiansky and O’connell [31] calculated the elastic modulus of the damage medium based on the “self-consistent method.” The TCK model proposes that cracks in rock mass will be activated under the tensile stress of explosive load and the cracks will squeeze and collide with each other to form interaction force. The effective bulk modulus of a homogeneous rock mass with original damage can be expressed as(24)KK0=1−169f1vCd=1−D,where K and K0, respectively, represent the bulk modulus of the damage medium and initial bulk modulus; f1v=1−v2/1−2v, in which v is Poisson’s ratio; Cd is the medium crack density; and D is the damage scalar.Xie et al. [32] put forward Poisson’s ratio change formula of the rock mass medium under impact load:(25)v=v01−169Cd,where v0 is the initial Poisson’s ratio of the medium.Grady and Kipp [33] proposed that the number of original damage cracks N in a unit medium under impact tensile stress can be accurately calculated by Weibull function:(26)N=ξp3km,where p is the impact stress and ξ and m represent the distribution parameters of Weibull function.Taylor et al. [34] proposed a formula for the radius of fragments in brittle media such as rock mass under dynamic impact loading:(27)r=1220KICρCε˙v,max2/3,where KIC and ρ, respectively, represent the fracture toughness and density of the medium, C is the longitudinal wave velocity, and ε˙v,max represents the volume strain rate at fracture. The functional relationship between the average radius of fragments r and the crack density Cd can be expressed as(28)Cd=Nϖr3,where ϖ is the proportional coefficient.Substituting formulas (26) and (27) into (28) gives(29)Cd=ξp3kmϖ⋅1820KICρCε˙v,max2=52ξϖp3kmKICρCε˙v,max2.Takingk=ξϖ in equation (29) yields(30)Cd=52kP3KmKICρC2ε˙vmax−2.The formula of tensile crack damage function can be obtained from equation (24) as(31)D=169f1v¯Cd.From the generalized Hook’s law, the stress-strain function relationship is as follows:(32)Δσij=KdδijΔεv+2GdΔeij,where Kd and Gd represent the modulus of the damaged medium and Δεv and Δeij are body strain increments and partial strain increments, Δεv=3εm, in which εm is the average normal stress.(2) Constitutive Relationship of Hardening Effect under Compressive Stress. Zhang et al. [35] pointed out in the report that the elastoplastic constitutive relationship of the rock mass medium obeying isotropic hardening under the action of detonation gas compression (yield equation) can be expressed as(33)Φ=σi2−σy2,where σi=3/2sij. According to the generalized Hooke’s law, the tensor form is used to represent the deviatoric stress Sij as(34)sij=σij−13σkkδij.During the explosion process, the functional relationship between the effective plastic strainεeffP of the rock mass and the yield limit σy of the fracture failure under the dynamic load is(35)εeffP=∫0t32ε˙ijPdt,σy=σ0+βEpεeffP,where Ep is the plastic hardening modulus, usually set as 0.1E, β is the hardening parameter, and ε˙ijp=ε˙ij−ε˙ije. In formula (33), when Φ≤0, the condition is elastic or neutral loading; when Φ>0, the condition is plastic hardening.(3) Pull-Pressure Damage Constitutive Model Embedding Simulation Software. The large-scale finite element analysis software LS-DYNA provides a rich material library for instantaneous power and large deformation damage simulation and also provides users with a convenient model import window. In this study, the calculation model of the cracking damage of the slope rock mass with cracks and the constitutive relation of the hardening effect under compressive stress are coupled; then the simulation software is imported for step-blasting calculation, which can simulate the trajectory of the throwing motion and the fragmentation of the rock mass under the impact of the explosion. The specific process is shown in Figure 7.Figure 7
Flowchart of the damage model.


## 3.3. Simulation Result Analysis

### 3.3.1. Casting State of Motion

The distribution of the casting state of the rock at different times is shown in Figure8. The model initiates the cast velocity from the foot of the slope at 7 ms. As time elapses, the casting phenomenon evolves from the foot of the slope to the free surface. By monitoring the model casting units, we found that the maximum casting velocity appeared at the middle of the blast-hole pack, with a velocity of 21.6 m/s. Additionally, this value is consistent with the conclusion regarding the maximum casting velocity proposed by Huang et al. [36, 37], which is 18–28 m/s. The average scattering distance is approximately 47 m, the postcracking distance is approximately 1.2 m, and the casting rate is approximately 34%.Figure 8
Variety of damage and deformation of the rock at different times.
(a)
(b)
(c)
(d)
(e)
(f)


### 3.3.2. Rock Mass Fragmentation Analysis

The bulks are mainly distributed at the foot of the slope, the sides of the contact surface between the explosive and the blockage, and the free surface at the top of the step. The bulks at the bottom of the slope mainly depend on the selection of the chassis resistance line value, and the difference in the resistance line value results in different bulk distributions. As shown for rock No. 5 in Figure9, the diameter is approximately 2.3 m. The bulks appear on both sides of the contact surface of the explosive, and the blockage is mainly caused by the different mechanical properties and physical parameters of the material on both sides of the contact surface. The stress wave is reflected and transmitted at the contact surface, and the reflection stretch acts on the penetrating microcrack, which causes the rock mass to break away from the parent rock under tensile stress to form bulks. As shown in rock block Nos. 1 and 2 and Nos. 3 and 4 in Figure 9, the diameter distribution is in the range of 1.6 to 2.7 m. The bulks appear on the free surface at the top of the step because of the transmission of stress waves at the interface of the explosive and blockage to the top free surface. The transmitted wave also forms a reflection stretching effect, which causes the rock mass in the region to break along the internal microfracture to form bulks. As shown in rock block No. 6 in Figure 9, the diameter is approximately 1.7 m.Figure 9
Distribution of rock blocks.


### 3.3.3. Analysis of Simulated Vibration Wave Velocity

A particle at the surface location of the simulation model is selected as the vibration speed-monitoring station, the installation position of the vibration-measuring instrument in the monitoring point and the field monitoring test is spatially consistent, and the specific change law is shown in Figure10. The simulated three-axis peak vibration velocity of the monitoring points is Vx = 0.68 cm/s, Vy = 0.74 cm/s, and Vz = 0.59 cm/s, respectively. The triaxial velocity variation law, the peak size, and the time point at which the peak occurs are in good agreement with the data obtained by the vibrometer monitoring, It is reasonable to take the surface vibration velocity Vy = 1.1 cm/s, when the effect evaluation is carried out in the text (the three-axis velocity direction is consistent with the data in Figure 6).Figure 10
Vibration velocity of the monitoring point by numerical simulation. (a)X-direction vibration of the monitoring point by numerical simulation. (b) Y-direction vibration of the monitoring point by numerical simulation. (c) Z-direction vibration of the monitoring point by numerical simulation.
(a)
(b)
(c)


## 3.3.1. Casting State of Motion

The distribution of the casting state of the rock at different times is shown in Figure8. The model initiates the cast velocity from the foot of the slope at 7 ms. As time elapses, the casting phenomenon evolves from the foot of the slope to the free surface. By monitoring the model casting units, we found that the maximum casting velocity appeared at the middle of the blast-hole pack, with a velocity of 21.6 m/s. Additionally, this value is consistent with the conclusion regarding the maximum casting velocity proposed by Huang et al. [36, 37], which is 18–28 m/s. The average scattering distance is approximately 47 m, the postcracking distance is approximately 1.2 m, and the casting rate is approximately 34%.Figure 8
Variety of damage and deformation of the rock at different times.
(a)
(b)
(c)
(d)
(e)
(f)


## 3.3.2. Rock Mass Fragmentation Analysis

The bulks are mainly distributed at the foot of the slope, the sides of the contact surface between the explosive and the blockage, and the free surface at the top of the step. The bulks at the bottom of the slope mainly depend on the selection of the chassis resistance line value, and the difference in the resistance line value results in different bulk distributions. As shown for rock No. 5 in Figure9, the diameter is approximately 2.3 m. The bulks appear on both sides of the contact surface of the explosive, and the blockage is mainly caused by the different mechanical properties and physical parameters of the material on both sides of the contact surface. The stress wave is reflected and transmitted at the contact surface, and the reflection stretch acts on the penetrating microcrack, which causes the rock mass to break away from the parent rock under tensile stress to form bulks. As shown in rock block Nos. 1 and 2 and Nos. 3 and 4 in Figure 9, the diameter distribution is in the range of 1.6 to 2.7 m. The bulks appear on the free surface at the top of the step because of the transmission of stress waves at the interface of the explosive and blockage to the top free surface. The transmitted wave also forms a reflection stretching effect, which causes the rock mass in the region to break along the internal microfracture to form bulks. As shown in rock block No. 6 in Figure 9, the diameter is approximately 1.7 m.Figure 9
Distribution of rock blocks.


## 3.3.3. Analysis of Simulated Vibration Wave Velocity

A particle at the surface location of the simulation model is selected as the vibration speed-monitoring station, the installation position of the vibration-measuring instrument in the monitoring point and the field monitoring test is spatially consistent, and the specific change law is shown in Figure10. The simulated three-axis peak vibration velocity of the monitoring points is Vx = 0.68 cm/s, Vy = 0.74 cm/s, and Vz = 0.59 cm/s, respectively. The triaxial velocity variation law, the peak size, and the time point at which the peak occurs are in good agreement with the data obtained by the vibrometer monitoring, It is reasonable to take the surface vibration velocity Vy = 1.1 cm/s, when the effect evaluation is carried out in the text (the three-axis velocity direction is consistent with the data in Figure 6).Figure 10
Vibration velocity of the monitoring point by numerical simulation. (a)X-direction vibration of the monitoring point by numerical simulation. (b) Y-direction vibration of the monitoring point by numerical simulation. (c) Z-direction vibration of the monitoring point by numerical simulation.
(a)
(b)
(c)


## 3.4. Fuzzy Comprehensive Evaluation

### 3.4.1. The Basic Steps

According to the actual postblasting scenario and a similar simulation test, the 11 evaluation indicators selected are considered sequentially, and the set of evaluation object factors is established as follows:(36)U=42,0.25,47,0.63,1.2,0.07,1.1,1.2,215,34.1,6.3.Step 1 (determining the membership degree of each index). The factor set U of the evaluation objects is substituted into equations (6)–(11), and the evaluation matrix determined by the Gaussian distribution membership function is obtained after calculation:(37)Rifgaussian=00.4870.1980.32510.294100.0560.7380.9260.94600.32500.35700.3250.0110.1770.2250.1760.64500000.3250.1970000.355T.The factor set U of the evaluation objects is substituted into equations (12)–(17), and the evaluation matrix determined by the triangle membership function is obtained after calculation:(38)Riftriangle=0.40.600.78610.6130.33300.60.7500.5290.5840.0870.21400.3870.6670.3250.40.2500.4710.4160.91300000.32500000T.In the factor set, the smaller-better type of evaluation objects u1, u2, u3, u4, u5, u6, and u7 is substituted into equations (6)–(8) and (12)–(14) for calculation, whereas the larger-better type u8, u9, u10, and u11 is substituted into equations (9)–(11) and (15)–(17) for calculation. The comment set constructed is of the decreasing type in this paper, so the evaluation object is calculated in a smaller-better membership function, and the calculated values are sequentially put into the evaluation matrix Rif according to the order of comments. Conversely, the membership values calculated by substituting the evaluation object into a larger-better membership function need to be sorted back first and then put into the evaluation matrix Rif in turn.Step 2 (determining the weight of each index using the AHP). The judgment matrix is constructed as follows using the expert investigation method according to equation (18):(39)R=1561773463215111532131211314161116211413121515156177234211713121711141413151617121171114131214151334124411211214231343111121316121432121113141335125412311124516523411.The weight vector is obtained after normalizing each index according to equation (19):(40)W=0.2300.0400.0290.1840.0210.0250.0920.0680.0470.1130.151.CR is obtained after calculation according to the formula CR = 0.052 < 0.1, which shows that the judgment matrix is consistent and there is no need to revise it.Step 3 (comprehensive evaluation). The evaluation matrix R obtained after calculation and weight vector W obtained after normalization are substituted into equations (2)–(4) for calculation to obtain the evaluation results, which are presented in Table 5.Table 5
Comparison of fuzzy comprehensive evaluations.
Evaluation object
Membership function distribution
Membership degree
Better (V1), bnor
Common (V2), batt
Poor (V3), rabn
Blasting effect
Gaussian distribution
0.227
0.409
0.364
Triangle distribution
0.493
0.507
0According to the maximum membership degree principle, the evaluation results of this blasting result controlled by the Gaussian distribution membership function and triangle membership function are both good. However, from the numerical distribution of each component of the two evaluation results, it is obvious that there are great differences between the two evaluation results, which show that the maximum membership principle is invalid in this case. This is because, in the application process, the maximum membership degree principle loses the information of nonmajor components in the evaluation results, which leads to the obvious unreasonable evaluation results in Table5.


### 3.4.2. Analysis of the Effectiveness Index of the Maximum Membership Degree Principle

The results obtained by the two methods of evaluation are substituted into equation (21) for calculation, and two effectiveness indices are obtained as follows: αgauss=0.47 and αtriangle=0.46. The range of the confidence level judged according to the implementation of the principle of maximum membership indicates that the principle of maximum membership is inefficient in this case. Therefore, using the closeness analysis method is considered to process the evaluation results to determine the specific level of the evaluation object. The nearness method is considered to process evaluation results to determine the specific grade of the evaluation objects.


### 3.4.3. Applying the Nearness Method to Analyze the Evaluation Results

The results of two types of evaluation are input into equation (22) for calculation. Two evaluation result vectors controlled by the Gaussian distribution and triangle membership function are, respectively, denoted as D1 and D2, and the results of the nearness calculation are as follows:(41)ND1,E1=ND11,E1=1−11100.773×1+0.409×2+0.364×3=0.976,ND1,E2=ND12,E1=1−11100.591×1+0.227×2+0.364×3=0.981,ND1,E3=ND13,E1=1−11100.636×1+0.409×2+0.227×3=0.981,ND2,E1=ND21,E1=1−11100.507×1+0.507×2+0×3=0.986,ND2,E2=ND22,E1=1−11100.493×1+0.493×2+0×3=0.987,ND2,E3=ND23,E1=1−11101×1+0.507×2+0.493×3=0.967.After processing the fuzzy evaluation results controlled by two different membership function forms using the nearness degree algorithm, it is found that when the factor set defined by the Gaussian-type membership function is considered the value of the membership degree in the comment set, the degree of the blasting effect evaluation which belongs to the common comment set is consistent with that which belongs to the poor comment set; that is, the final evaluation result of the blasting effect using the Gaussian-type membership function is between the “common” and “poor” grades. Similarly, the degree of the blasting effect evaluation that belongs to “better” is essentially consistent with that which belongs to “common” when using the triangle membership function, which demonstrates that the final evaluation result of the blasting effect using the triangle membership function is between the “better” and “common” grades.Through further observation of the blasting site and postblasting effect, we found that the indices of explosive consumption and detonator consumption in the blasting construction process are slightly beyond the reasonable range; the longest cast distance and blasting output per meter hole are also slightly less than the reasonable values given in [16, 17, 25], and there are obvious toe rocks and bulks; loose rock poured into goaf accounts for approximately 20% of the stripped rock mass, and the accumulation of the big pile is not ideal. Thus, the result of the subjective blasting effect evaluation should be consistent with the final result of the fuzzy comprehensive evaluation under the control of the Gaussian membership function after processing the nearness algorithm. The evaluation results controlled by two different membership functions can all be considered common by the nearness algorithm, which is because the number of comment sets selected is small (three comment sets are selected, which are better, common, and poor, in this paper). If the evaluation set in this paper extends to five or more items, then the results of the fuzzy comprehensive evaluation processed by the nearness algorithm, and the final evaluation results controlled by two different membership functions, will be different.


## 3.4.1. The Basic Steps

According to the actual postblasting scenario and a similar simulation test, the 11 evaluation indicators selected are considered sequentially, and the set of evaluation object factors is established as follows:(36)U=42,0.25,47,0.63,1.2,0.07,1.1,1.2,215,34.1,6.3.Step 1 (determining the membership degree of each index). The factor set U of the evaluation objects is substituted into equations (6)–(11), and the evaluation matrix determined by the Gaussian distribution membership function is obtained after calculation:(37)Rifgaussian=00.4870.1980.32510.294100.0560.7380.9260.94600.32500.35700.3250.0110.1770.2250.1760.64500000.3250.1970000.355T.The factor set U of the evaluation objects is substituted into equations (12)–(17), and the evaluation matrix determined by the triangle membership function is obtained after calculation:(38)Riftriangle=0.40.600.78610.6130.33300.60.7500.5290.5840.0870.21400.3870.6670.3250.40.2500.4710.4160.91300000.32500000T.In the factor set, the smaller-better type of evaluation objects u1, u2, u3, u4, u5, u6, and u7 is substituted into equations (6)–(8) and (12)–(14) for calculation, whereas the larger-better type u8, u9, u10, and u11 is substituted into equations (9)–(11) and (15)–(17) for calculation. The comment set constructed is of the decreasing type in this paper, so the evaluation object is calculated in a smaller-better membership function, and the calculated values are sequentially put into the evaluation matrix Rif according to the order of comments. Conversely, the membership values calculated by substituting the evaluation object into a larger-better membership function need to be sorted back first and then put into the evaluation matrix Rif in turn.Step 2 (determining the weight of each index using the AHP). The judgment matrix is constructed as follows using the expert investigation method according to equation (18):(39)R=1561773463215111532131211314161116211413121515156177234211713121711141413151617121171114131214151334124411211214231343111121316121432121113141335125412311124516523411.The weight vector is obtained after normalizing each index according to equation (19):(40)W=0.2300.0400.0290.1840.0210.0250.0920.0680.0470.1130.151.CR is obtained after calculation according to the formula CR = 0.052 < 0.1, which shows that the judgment matrix is consistent and there is no need to revise it.Step 3 (comprehensive evaluation). The evaluation matrix R obtained after calculation and weight vector W obtained after normalization are substituted into equations (2)–(4) for calculation to obtain the evaluation results, which are presented in Table 5.Table 5
Comparison of fuzzy comprehensive evaluations.
Evaluation object
Membership function distribution
Membership degree
Better (V1), bnor
Common (V2), batt
Poor (V3), rabn
Blasting effect
Gaussian distribution
0.227
0.409
0.364
Triangle distribution
0.493
0.507
0According to the maximum membership degree principle, the evaluation results of this blasting result controlled by the Gaussian distribution membership function and triangle membership function are both good. However, from the numerical distribution of each component of the two evaluation results, it is obvious that there are great differences between the two evaluation results, which show that the maximum membership principle is invalid in this case. This is because, in the application process, the maximum membership degree principle loses the information of nonmajor components in the evaluation results, which leads to the obvious unreasonable evaluation results in Table5.


## 3.4.2. Analysis of the Effectiveness Index of the Maximum Membership Degree Principle

The results obtained by the two methods of evaluation are substituted into equation (21) for calculation, and two effectiveness indices are obtained as follows: αgauss=0.47 and αtriangle=0.46. The range of the confidence level judged according to the implementation of the principle of maximum membership indicates that the principle of maximum membership is inefficient in this case. Therefore, using the closeness analysis method is considered to process the evaluation results to determine the specific level of the evaluation object. The nearness method is considered to process evaluation results to determine the specific grade of the evaluation objects.


## 3.4.3. Applying the Nearness Method to Analyze the Evaluation Results

The results of two types of evaluation are input into equation (22) for calculation. Two evaluation result vectors controlled by the Gaussian distribution and triangle membership function are, respectively, denoted as D1 and D2, and the results of the nearness calculation are as follows:(41)ND1,E1=ND11,E1=1−11100.773×1+0.409×2+0.364×3=0.976,ND1,E2=ND12,E1=1−11100.591×1+0.227×2+0.364×3=0.981,ND1,E3=ND13,E1=1−11100.636×1+0.409×2+0.227×3=0.981,ND2,E1=ND21,E1=1−11100.507×1+0.507×2+0×3=0.986,ND2,E2=ND22,E1=1−11100.493×1+0.493×2+0×3=0.987,ND2,E3=ND23,E1=1−11101×1+0.507×2+0.493×3=0.967.After processing the fuzzy evaluation results controlled by two different membership function forms using the nearness degree algorithm, it is found that when the factor set defined by the Gaussian-type membership function is considered the value of the membership degree in the comment set, the degree of the blasting effect evaluation which belongs to the common comment set is consistent with that which belongs to the poor comment set; that is, the final evaluation result of the blasting effect using the Gaussian-type membership function is between the “common” and “poor” grades. Similarly, the degree of the blasting effect evaluation that belongs to “better” is essentially consistent with that which belongs to “common” when using the triangle membership function, which demonstrates that the final evaluation result of the blasting effect using the triangle membership function is between the “better” and “common” grades.Through further observation of the blasting site and postblasting effect, we found that the indices of explosive consumption and detonator consumption in the blasting construction process are slightly beyond the reasonable range; the longest cast distance and blasting output per meter hole are also slightly less than the reasonable values given in [16, 17, 25], and there are obvious toe rocks and bulks; loose rock poured into goaf accounts for approximately 20% of the stripped rock mass, and the accumulation of the big pile is not ideal. Thus, the result of the subjective blasting effect evaluation should be consistent with the final result of the fuzzy comprehensive evaluation under the control of the Gaussian membership function after processing the nearness algorithm. The evaluation results controlled by two different membership functions can all be considered common by the nearness algorithm, which is because the number of comment sets selected is small (three comment sets are selected, which are better, common, and poor, in this paper). If the evaluation set in this paper extends to five or more items, then the results of the fuzzy comprehensive evaluation processed by the nearness algorithm, and the final evaluation results controlled by two different membership functions, will be different.


## 4. Conclusions

This paper proposes a new evaluation method for the nonlinear relationship between the blasting effect and evaluation factors. This method selects 11 evaluation indicators to evaluate the postexplosion effect. Part of the indicator parameters are derived from field test monitoring, and another part of the indicator parameters are derived from numerical simulation. The main conclusions of this study are as follows:(1)
Two different membership functions are defined for fuzzy evaluation of the blasting effect. Gaussian and triangle membership functions are selected, respectively. The preliminary evaluation results of the two membership functions are common. But the initial results have a large difference in the weight of the comment set, so the initial evaluation results are not credible. It is essential to further optimize the evaluation results by the proximity analysis method. The appraisal result is between common and poor with the Gaussian membership function, while the result is between better and common with the triangle membership function. Through the observation, it is found that the local blasting effect is not ideal, and the evaluation result of the Gaussian membership function is closer to the actual situation. Therefore, the Gaussian function instead of the linear function should be adopted to define the membership function of each index when nonlinear factors such as the blasting effect are selected to measure the evaluation results of the target layer.(2)
Slope blasting is the problem of fracture damage of the rock mass under blasting impact load. The rupture process includes two parts: the rupture effect of tensile cracking and compression hardening. Therefore, the establishment of a reasonable constitutive relationship of tensile and compression damage can accurately simulate the blasting phenomenon of the rock mass. It is assumed that the compression state of the rock mass under the impact load obeys Mises yield failure, the stress-strain relationship satisfies the bilinear elastoplastic constitutive relation, and then the constitutive damage model of hardening effect under compressive stress can be constructed. In contrast, the rock mass satisfies the TCK relationship under the impact load, the material response obeys the elastic Hooke’s law, and the calculation model of the fracture damage can be constructed. By coupling the rock mass fracture damage calculation model and the compressive stress hardening effect constitutive relationship, the constitutive relation of rock mass tensile compression damage can be constructed. And the user-defined material interface function provided by LS-DYNA is used to embed the tensile and compressive damage evolution model into the elastoplastic constitutive material to accurately simulate the blasting effect of the slope.(3)
Based on the actual postexplosion situation and similar simulation tests, the 11 evaluation indicators in this paper are assigned. The strongest indicator of the impact on the blasting effect obtained through fuzzy calculation is boulder yield, and the weakest indicator is postcracking distance. The impact of each indicator varies greatly, and the weight of the strongest indicator is about 10 times that of the weakest. The 11 indicators are arranged in the order of impact from strong to weak as follows: boulder yield (23%), explosive consumption (18.4%), blasted volume per meter of hole (15.7%), throw rate (11.3%), vibration velocity (9.2%), loose coefficient (6.8%), cast distance (4.7%), root rate (4%), flying distance of flyrock (2.9%), detonator unit consumption (2.5%), and postcracking distance (2.1%).



---

*Source: 6029190-2020-04-15.xml*