Abstract

Cavability assessment is an important subject during the feasibility stages before determining whether to use block caving mining. This paper provides a fuzzy comprehensive assessment (FCA) approach based on the cavability assessment approaches and its influencing factors, which are all fuzzy. This approach combines the cavability influencing factors with engineering empirical approaches by fuzzy mathematics, which improves the applicability of the cavability assessment results. This approach is applied to assess the cavability via cores in the Luoboling copper molybdenum mine. The spatial distribution of the rock mass cavability at different depths of the borehole is obtained. The cavability ranks of various rocks are determined in different locations. These assessment results can provide a basis for demonstrating the feasibility of block caving mining in the Luoboling copper molybdenum mine. The study can also provide a basis for the design of mining engineering.

1. Introduction

Block caving mining refers to all mining operations in which the ore body caves naturally after undercutting and the caved material is recovered through drawpoints [1]. As a special mechanism in mechanics and technology, an important subject is to assess the rock mass cavability during the feasibility stages before determining whether to use this mining method [2]. The rock mass cavability has a major influence on the mining block heights, production rates, undercutting orientations, undercutting areas, draw controls, mining rates, preconditioning engineering, and so on and is an important guarantee for a mine to achieve the expected economic benefits. Cavability assessment is conducted by classifying the rock mass cavability and determining the ranks according to a given geological environment. And the mine determines whether to use block caving mining at the current industrial level. Cavability assessment is also a multi-index and nonlinear complex system engineering of the rock mass.

In the rock mass cavability field, most approaches of cavability assessment are based on numerical modelling [3, 4], mathematical models [2, 57], and geomechanical classifications [1, 816]. In addition, geomechanical classifications have been widely used in cavability assessment [13, 17]. Geomechanical classifications seem to be in tune with cavability assessment in describing the same problems of a rock mass [17]. Because geomechanical classifications use engineering empirical assessments of the rock mass strength in relation to the existing stresses and measures of the rock structure, the classifications are fuzzy, such as rock quality designation (RQD) [12], rock mass rating (RMR) [8, 9, 16], mining rock mass rating (MRMR) [1, 11], rock mass quality Q-classification (Q) [10, 14], and rock mass basic quality (BQ) [15]. The selection of the influencing factors and the determination of the rock mass ratings are both fuzzy because different approaches use different influencing factors as indicators. In addition, the influencing factors of cavability are interrelated with each other and present great complexities, which lead to different cavability assessment results. Therefore, the cavability assessment and its influencing factors are both fuzzy, and fuzzy mathematics can accurately describe and address the fuzzy phenomena. Therefore, we will adopt fuzzy mathematics to assess the cavability of the rock mass.

Fuzzy mathematics has been applied to predict petrophysical rock parameters [18] and mechanical rock parameters [1922] and to analyse various properties [2, 2327] and phenomena [28, 29] of the rock mass. In a rock mass cavability study, Rafiee et al. [2] designed a fuzzy expert semiquantitative coding methodology to assess the cavability of the rock mass, and Rafiee et al. [7] applied the fuzzy rock engineering systems method to account for the intricate interactions that exist among parameters in real projects. Shaoyong et al. [6] combined fuzzy mathematics and the matter element analysis method and established a model of cavability of the rock mass in terms of complex fuzzy matter element analysis. Although many researchers have applied fuzzy mathematics to assess the cavability of a rock mass, the researchers conducted cavability modelling based on influencing factors. The influencing factors of cavability are not very clear in the current understanding. Under these circumstances, engineering empirical approaches are still important.

In this paper, we combine influencing factors with engineering empirical approaches by fuzzy mathematics and carry out a fuzzy comprehensive assessment (FCA). This approach improves the applicability of the assessment results of cavability. The approach is applied to assess the cavability of cores in a mine, and we obtain the spatial distribution of the rock mass cavability at different depths of the borehole. The cavability ranks of the various rocks in the hanging wall, ore body, and rocks in the ore body and footwall are determined. The assessment results provide a reference and basis to decide whether to adopt the block caving mining method and determine the mining engineering design.

2. Fuzzy Assessment (FA) of Rock Mass Cavability Based on Influencing Factors

In the rock mass cavability field, the approaches of cavability assessment were dependent on influencing factors. It was necessary to take certain influencing factors into consideration in the FA of rock mass cavability.

2.1. Determination of the Influencing Factors and Assessment Ranks

To date, studies on the influencing factors of cavability have been presented in the literature [117]. These influencing factors can be summarized as rock strength, discontinuity properties, water, and in situ stress. In the assessment of rock mass cavability, the appropriate selection of influencing factors was critical to the reliability of the assessment results. When relatively few influencing factors were selected, these factors could not fully reflect the rock mass cavability and even lead to incorrect results of the cavability assessment. When too many influencing factors were selected and these factors connected to each other, this situation might exaggerate the influence of a certain factor on the rock mass cavability and lead to incorrect results. These incorrect results were due to the influencing factors being interrelated with each other and subjectivity in determining the factors. Therefore, we analysed the relationship and difference among the influencing factors of cavability and determined the influencing factors based on the present studies.

In the approaches of cavability assessment, the indices of the rock strength were the uniaxial compressive strength (UCS) or point load strength index (Is(50)). Determining the Is(50) was a more practical, time-saving, and economical method compared to determining the UCS [30]. A large number of studies [3134] have shown that the Is(50) has a good correlation with the UCS. Therefore, the Is(50) was chosen to represent the rock strength in this paper. It was clear that the cavability of a rock mass decreases when the rock strength increases [7]. To quantify the Is(50), the Is(50) of intact rock, based on the RMR classification, was subdivided into five ranks (the results of the ranks are listed in Table 1). The discontinuity properties were some of the most important influencing factors on the cavability of the rock mass. The most important of these properties used to describe the discontinuities were the RQD, joint spacing, intactness index of the rock mass, volumetric joint count of the rock mass, joint roughness (), joint aperture (), and joint filling (). It was clear that the RQD had a good correlation with the joint spacing, intactness index of the rock mass, and volumetric joint count of the rock mass. At the same time, the RQD was most commonly used in cavability assessments. Therefore, these factors were chosen to represent the discontinuity properties, including the RQD, , , and . These factors were also subdivided into five ranks based on the RMR classification (the results of the ranks are listed in Table 1). Water was usually described qualitatively. The water condition () in this study, as in the RMR classification, was divided into five ranks (the results of the ranks are listed in Table 1). Lastly, the ratio of the UCS value to the in situ stress value was chosen to represent the in situ stress [7]. This ratio could be translated into the ratio of the Is(50) value to the in situ stress value according to the relationship between the Is(50) and the UCS (UCS=22.8 Is(50) [15]). The ratio of Is(50) value to in situ stress value is abbreviated to Iss. The results of the Iss ranks are listed in Table 1. The influencing factor set of cavability was determined as U=, u2, u3, u4, u5, u6, =, RQD, , , , , .

The cavability assessment ranks were generally divided into five ranks, including extremely difficult caving I, difficult caving II, fair caving III, easy caving IV, and extremely easy caving V. For the convenience of establishing the membership function during a follow-up operation, the ranks of cavability assessment were quantified (including quantitative value and quantitative range ). The quantitative results of the ranks are listed in Table 1. The rank set of cavability assessment was determined as V, v2, v3, v4, =, II, III, IV, .

2.2. Fuzzy Assessment Matrix and Its Membership Function

It was important for the fuzzy assessment matrix to determine the membership of each influencing factor. These factors can be divided into qualitative indices and quantitative indices according to Table 1. The memberships of the qualitative indices could be determined by counting the assessment frequencies from several surveyors. The memberships of the quantitative indices could be determined by membership functions.

Before establishing the membership functions of the quantitative indices, it was necessary to convert the measured value into a value within the quantitative range for each quantitative index. This method was convenient for establishing the membership function and assessing the rock mass cavability. The converted map function is established as follows:where is the measured value of the quantitative index. In (1), the cavability of the rock mass decreases with increasing measured value, and i is 1, 2, or 7. In (2), the cavability of the rock mass increases with increasing measured value, and i is 4. and are the upper range value and the lower range value of the classification range based on measured value , respectively. When tends to infinity at the boundary, the value should be limited according to the measured value and empirical value. and are the upper range value and the lower range value of the quantitative range , respectively, and the values correspond to . The converted map function uses linear transformation and does not change the assessment results.

After establishing the converted map function, it was necessary to establish the membership functions of the quantitative indices. In fuzzy set theory, the membership function of an index might contain some uncertainty, so the membership is expressed as a degree of belonging to a set [29]. Different people might establish different membership functions for the same fuzzy set because of the limitations of human understanding. However, Yonghua et al. [35] proved that different membership functions had equivalent characteristics in rock mass engineering. Triangular and trapezoidal shapes were the most common types of membership functions in rock or rock mass engineering [18, 22, 23, 26, 28, 29].

Therefore, we adopted the inference method of fuzzy set characteristics to establish the membership function and combined the triangular and trapezoidal shapes of the membership function. The inference method entailed looking for the special elements in the quantitative range, such as the value of membership being equal to 0, 0.5, or 1. The total value of membership was 1 for each influencing factor in the five ranks. In rock or rock mass engineering [6, 18, 22, 23, 26, 28, 29], the membership function usually adopted an intermediate type for each quantitative index. That is, the value of membership was 0.5 at the endpoint for each quantitative range. The value of membership was 1 in the middle range for each quantitative range. In the middle range for each quantitative range, the value of membership was 0 for the neighbourhood range. Eventually, the membership function Aji=Aj(f(ui)) is established based on the above principles as follows:where δ is the neighbourhood value centred on the midpoint of each quantitative range. The default value of δ is 0.05 in this paper.

According to (1) to (7), the memberships of the quantitative indices can be calculated. Combined with the memberships of the qualitative indices, the fuzzy assessment matrix R can be established as follows:where Aji is the value of membership and Aji means that influencing factor ui has a membership of rank vj in the cavability assessment.

2.3. Determination of the Fuzzy Relative Weight Based on the Analytic Hierarchy Process

Because the extent of influence was different for each influencing factor in the cavability assessment, it was necessary to determine the fuzzy relative weight of the factor. Among the approaches for determining the weight, the analytic hierarchy process has been widely applied due to its simplicity, scalability, and pairwise and easy comparison of variables by assigning weights [36]. The analytic hierarchy process has been extensively used in complex decision making with fuzzy mathematics [3639]. The analytic hierarchy process refers to a multicriteria decision-making approach in which factors are arranged in a hierarchic structure [40]. The analytic hierarchy process can combine qualitative analysis with quantitative analysis in the process of cavability assessment. The determination process is as follows.

The first step was the structuring of the rock mass cavability as a hierarchy. The judgement factors that contribute to the cavability were determined. That is, the judgement factors were the influencing factor set U=, u2, u3, u4, u5, u6, =, RQD, , , , , .

The second step was the elicitation of the pairwise comparison judgements and establishing a judgement matrix. The elements were arranged into a matrix and judgements were elicited from the people who had difficulties about the relative importance of the elements with respect to the rock mass cavability. The scale to use in making the judgements was 1~9 and the reciprocal [40]. In evaluating the judgement factors relative to the rock mass cavability, the evaluation was conducted according to the present research findings on the relative importance of the judgement factors with respect to the rock mass cavability. This paper mainly combines the RMR with the Q-classification method, and the judgement matrix P can be established as follows:

The third step was to calculate the order of the relative importance. When calculating the maximum eigenvalue of P with =7.51, the eigenvector is as follows:

The eigenvector is the order of the relative importance. The fuzzy relative weight coefficient of each index was obtained by normalizing the eigenvector . The weight vectors are =[, , , , , , ]=[0.11, 0.34, 0.10, 0.15, 0.17, 0.09, 0.04].

The fourth step was the consistency check. The weight coefficient of each factor was obtained. It was necessary to check whether the distribution of the weight coefficients was reasonable. The formula of the consistency index (CI) is as follows:where is the number of judgement factors and =7.

The CI is compared with the average random consistency index RI (=7, RI7=1.35). The consistency ratio, CR, can be obtained as follows:

According to (12), the judgement matrix meets the consistency check. That is, the distribution of the fuzzy relative weight coefficients is reasonable.

2.4. Fuzzy Mapping and Fuzzy Assessment Based on the Influencing Factors

After determining the weight vector and the fuzzy assessment matrix R, fuzzy subset B can be obtained by fuzzy linear variation. The fuzzy subset B is as follows:where “” is a synthetic operator. The weighted average model is adopted because all kinds of influencing factors affect the rock mass cavability.

The FA method was based on the influencing factors involved in calculating the quantitative FA value. The calculation method considered that membership of quantitative value was a weight coefficient, and the weighted average value of each was taken as a quantitative value of the FA. The quantitative FA value was compared with the quantitative range , and the rank of cavability was obtained. The FA value is calculated as follows:

3. Fuzzy Comprehensive Assessment (FCA) of the Rock Mass Cavability

At present, there are many assessment approaches for the rock mass cavability. However, these approaches were proposed based on certain geological conditions or on a given engineering background. For example, RMR was based on experience gained in numerous visits to construction sites abroad and in South Africa [8, 9, 16], MRMR was combined RMR with mining engineering [1, 11], and Q was originally developed to assist in the empirical design of tunnel and cavern reinforcement and support [10, 14]. There might be certain limitations or inadaptability in applying these approaches directly. If two or more approaches were applied at the same time, the results were often different. However, these approaches were based on a large number of engineering practices and engineering experiences. The approaches had a strong reference value in specific practical projects. Therefore, it was necessary to synthesize these assessment approaches according to the specific mining geology. Furthermore, more objective and reasonable results of the cavability were obtained. The results provided a strong reference and basis to decide whether to adopt the block caving mining and determine the mining engineering design. Fuzzy mathematics provided the method for synthesizing these different assessment approaches. The method was fuzzy mapping and fuzzy comprehensive assessment, which was based on fuzzy comprehensive assessment matrix and fuzzy relative weight.

3.1. Selection of the Assessment Approaches and Assessment Ranks

Among the assessment approaches, the most widely used approaches include rock quality designation RQD, rock mass rating RMR, mining rock mass rating MRMR, and rock mass quality Q-classification Q. The MRMR was proposed for mining but was most affected by engineering experience in the assessment process. Because the mine was in the feasibility stage, it was not put into production and no rock mass excavation engineering occurred. There was hardly any engineering experience that could be referenced, and it was impossible to revise the parameters in the MRMR. Therefore, the MRMR was not selected as an index for the FCA approach. At the same time, the FA of the rock mass cavability based on the influencing factors was introduced. The final selection of the assessment approach set was determined as Z=, z2, z3, z4, =, RMR, Q, BQ, . The rock cavability classifications based on the assessment approaches are listed in Table 2.

The ranks of cavability assessment were also divided into five ranks, including extremely difficult caving I, difficult caving II, fair caving III, easy caving IV, and extremely easy caving V. The rank set of the cavability assessment was determined as V, v2, v3, v4, =, II, III, IV, . The ranks of cavability are listed in Table 2.

3.2. Fuzzy Comprehensive Assessment Matrix and Its Membership Function

As shown in Table 2, the indices (assessment approaches) are quantitative indices in the FCA. The memberships of the quantitative indices could be determined by membership functions. Before establishing the membership function, it was necessary to convert the measured value into a value within the quantitative range for each quantitative index. This method was also convenient for establishing the membership function and assessing the rock mass cavability. The converted map function f(zi) is established as follows:where zi is the calculated value of the quantitative index and i is 1, 2, 3, 4, 5. Here, and are the upper range value and the lower range value of the classification range based on calculated value zi, respectively. When tends to infinity at the boundary, the value should be limited according to the calculated value and empirical value. and are the upper range value and the lower range value of the quantitative range , respectively, and correspond to zi.

After establishing the converted map function, it was necessary to establish the membership functions of the quantitative indices. The inference method was also adopted to establish the membership function. The established method was similar to the FA method. The membership function Zji=Zj(f(zi)) is established as follows:where δ is the neighbourhood value centred on the midpoint of each quantitative range. The default value of δ is 0.05 in this paper.

According to (15) to (20), the memberships of the quantitative indices can be calculated. The fuzzy comprehensive assessment matrix can be established:where Zji is the value of membership and Zji means that assessment approach zi has a membership of rank vj of the cavability assessment.

3.3. Determination of the Fuzzy Relative Weight Based on the Analytic Hierarchy Process

Because the extent of the influence was different for each assessment approach in the cavability assessment, it was necessary to determine the fuzzy relative weight. The analytic hierarchy process was also adopted. The determination process was as follows.

The first step was the structuring of the rock mass cavability as a hierarchy. The judgement factors that contribute to the cavability were determined. That is, the assessment approach set constituted the judgement factors, Z=, z2, z3, z4, =, RMR, Q, BQ, .

The second step was the elicitation of pairwise comparison judgements and establishing a judgement matrix. The established method was similar to the FA method. In making the judgements of the assessment approaches relative to the rock mass cavability, the process was done according to the applicable conditions, engineering backgrounds, and application statuses of these assessment approaches. This paper mainly relies on the specific mining geology of the Luoboling copper-molybdenum mine. The judgement matrix can be established as follows:

The third step was to calculate the order of the relative importance. When calculating the maximum eigenvalue of with =5.003, the eigenvector is as follows:

The eigenvector was the order of the relative importance. The fuzzy relative weight coefficient of each assessment approach was obtained by normalizing the eigenvector . The weight vector is =[0.05, 0.26, 0.17, 0.13, 0.39].

The fourth step was the consistency check. The weight coefficient of each factor was obtained. It was necessary to check whether the distribution of the weight coefficients was reasonable. The formula of the CI is as follows:where is the number of judgement factors and =5.

Compare the CI with the average RI (=5, RI5=1.12). The consistency ratio, , can be obtained as follows:

According to (25), the judgement matrix meets the consistency check. That is, the distribution of the fuzzy relative weight coefficients is reasonable.

3.4. Fuzzy Mapping and Fuzzy Comprehensive Assessment

After determining the weight vector and fuzzy assessment matrix , the fuzzy subset can be obtained by fuzzy linear variation. The fuzzy subset is as follows:where “” is a synthetic operator. The weighted average model is adopted because all kinds of assessment approaches affect the rock mass cavability.

The FCA method was used to calculate the quantitative FCA value. The calculation method considered that membership of quantitative value was a weight coefficient, and the weighted average value of each was taken as a quantitative value of the FCA. The quantitative value of the FCA was compared with the quantitative range , and the rank of cavability was obtained. The FCA value is calculated as follows:

4. Practical Application in the Luoboling Copper-Molybdenum Mine

The Luoboling copper-molybdenum mine belongs to a porphyry deposit. The characteristics of the ore body include deep burial, large distribution area, large thickness, large dip change, low grade, large reserves, and complex shape. According to the characteristics of the ore body, block caving mining was determined during the feasibility stages. Therefore, it was a crucial step to assess the rock mass cavability and obtain the spatial distribution maps of the cavability. Cavability assessment was beneficial for the engineering layout, the stope structure parameter selection, and determining whether to use block caving mining.

However, there was no mining excavation engineering during the feasibility stages for the Luoboling copper-molybdenum mine, and the characteristic parameters of the rock mass could not be obtained. This situation was also the same for other mines during the feasibility stages. However, a total of 176 boreholes were completed in the exploration stage of the Luoboling copper-molybdenum mine. The total footage was 130761.23 m, and the controlled area was 6.77 km2. A large number of cores were retained. Therefore, this paper determined the rock mass cavability through cores.

First, 25 boreholes were determined from the 176 boreholes according to the spatial position relationship between the borehole and the ore body (as shown in Figure 1). Other considerations included the spacing of borehole, the volume of work, and the shape of ore body. After that, each borehole core was divided into several groups in the vertical direction according to lithology and RQD value. If the lithology was consistent and the RQD value was close within a vertical distance, the cores were divided into a group. The rock mass cavability was assessed according to the groups.

4.1. Fuzzy Assessment (FA) of the Rock Mass Cavability

According to the groups, the measured values of the seven influencing factors were determined. The methods of measurement included mechanics experiments, field surveys, and measurements, which referred to the hydrogeology and engineering geology of the boreholes. The FA approach was adopted to assess the rock mass cavability. The result of the assessment is shown in Figure 2. As can be seen from Figure 2, the general trend of cavability is from easy caving IV to difficult caving II in the vertical direction of the boreholes.

4.2. Fuzzy Comprehensive Assessment (FCA) of the Rock Mass Cavability

The rock mass cavability was assessed by the approaches of RQD, RMR, Q, and BQ. The results of the assessment are shown in Figures 36.

It can be seen from Figures 26 that the trend of cavability ranged from extremely easy caving V to extremely difficult caving I in the vertical direction of the boreholes. However, the results of the assessment were quite different for the different approaches. Therefore, the FCA approach was also necessary for the cavability assessment.

The comprehensive approach of the FCA was adopted to assess the rock mass cavability. The result of the assessment is shown in Figure 7. As is seen from Figure 7, the general trend of cavability is from easy caving IV to difficult caving II in the vertical direction of the boreholes. In the boreholes as a whole, the upper part belongs to easy caving IV, and the lower part belongs to fair caving III and difficult caving II.

4.3. The Assessment Results of the Rock Mass Cavability

We obtained the spatial distribution of the rock mass cavability at different depths of the borehole. The advantage of the assessment cavability for cores was that we combined the cavability with the geological information from the boreholes (such as the lithology and location). We could count the lengths of the cores that had the same ranks of cavability and same locations. In this way, it was beneficial for analysing the rock mass cavability in different locations. We obtained the cavability ranks of the various rocks in the hanging wall, ore body, and rocks in the ore body and footwall (as shown in Figure 8), which were based on the spatial distribution of the FCA value in the vertical direction of the boreholes (Figure 7).

It can be seen from (a) of Figure 8 that the rock mass cavability of the rock in the hanging wall is mainly easy caving IV, a small portion of the rock is fair caving III, and a very small amount of the rock is difficult caving II or extremely easy caving V. It can be seen from (b) of Figure 8 that the rock mass cavability of the ore in the ore body is difficult caving II, fair caving III, and easy caving IV, and a very small amount of the rock is extremely easy caving V; fair caving III and easy caving IV account for approximately 70%. It can be seen from (c) of Figure 8 that the rock mass cavability of the rock in the ore body is difficult caving II, fair caving III, and easy caving IV; a very small amount of the rock is extremely easy caving V; and fair caving III and easy caving IV account for approximately 65%. It can be seen from (d) of Figure 8 that the rock mass cavability of the rock in the footwall is mainly fair caving III, a small amount of the rock is difficult caving II, and a small amount of the rock is easy caving IV. On the whole, the mine has hardly any rock mass of extremely difficult caving I and extremely easy caving V. And the rock mass cavability is mainly fair caving III and easy caving IV. The borehole cores can verify these conclusions during the feasibility stages. As shown in Figure 9, some typical cores of mine are presented. It can be seen from (a) of Figure 9 that a part of the rock is fragmented at the top of the borehole. As the depth increases the cores are relatively intact (as shown in (b) and (c) of Figure 9), but the cores have many discontinuities and help the rock mass naturally cave. Therefore, the assessment results of mine are correct. These assessment results can provide a basis for demonstrating the feasibility of block caving mining in the Luoboling copper-molybdenum mine. The study can also provide a basis for designing the mining engineering next.

5. Conclusions

In this study, a fuzzy comprehensive assessment (FCA) approach was provided (as shown in Figure 10) that was based on cavability assessment and its influencing factors, which were fuzzy. For this purpose, we determined the influencing factors and assessment approaches of cavability, established the converted map functions and membership functions, and adopted the analytic hierarchy process to determine the fuzzy relative weights. This approach combined the cavability influencing factors with engineering empirical approaches by fuzzy mathematics. The method improved the applicability of the assessment results of cavability. Finally, the FCA approach was applied to assess the cavability of cores in the Luoboling copper–molybdenum mine, which, according to the mine stage, has no mining excavation engineering at present. The spatial distribution of the rock mass cavability at different depths of the borehole was obtained. We combined cavability with the geological information from the boreholes and determined the cavability ranks of the various rocks in the hanging wall, ore body, and rocks in the ore body and footwall. The assessment results provided a reference and basis to decide whether to adopt the block caving mining method and to design the mining engineering.

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Authors’ Contributions

Rongxing He and Huan Liu contributed to the formulation of the overarching research goals and aims and conducted the FCA; Rongxing He and Fengyu Ren determined the influencing factors and assessment approaches of cavability; Fengyu Ren determined the measured boreholes in the Luoboling copper-molybdenum mine; Huan Liu carried out the FCA and wrote the paper; and Huan Liu, Guanghui Li, and Jing Zhang completed the mechanics experiments, field surveys, and measurements and referred to the hydrogeology and engineering geology of the boreholes in the Luoboling copper,molybdenum mine.

Acknowledgments

The study is jointly supported by grants from the National Key Research and Development Program of China (Grant no. 2016YFC0801604) and the Key Program of the National Natural Science Foundation of China (Grant no. 51534003). The authors are grateful for the support.