Research Article  Open Access
Deflection Mechanism and Safety Analysis of Coal Mine Shaft in Deep Soil Strata
Abstract
The main shaft and auxiliary shaft in the Guotun Coal Mine underwent large deflections, with deflection values of 359 mm and 322 mm, respectively. These two deflections represent the first occurrence of such large vertical shaft deviations in the soil strata in China. The deflection problem has seriously affected the hoisting safety and lining safety and has become a serious impediment to the sustainable production of mines. Therefore, the deflection mechanism must be determined. For this purpose, based on mining subsidence theory, the spatial probability integral method and a more accurate time function were used to establish a model, called 3D dynamic prediction model, for predicting the shaft movement. The formulas for calculating the lining stress caused by coal mining were based on established models. With measured shaft deflection data, the prediction parameters for deep soil strata were calculated on the basis of an inversion analysis. A comparative analysis of measured and calculated deflection values revealed that the reason for shaft deflection in Guotun Coal Mine is the insufficient size of the protection coal pillar (PCP); namely, the design parameters of the PCP in current codes are not applicable to the deep soil strata. As a result, under the asymmetric mining conditions, mining causes the shaft to deflect without damage and under the symmetric mining conditions, mining causes the lining to fracture. The results have an extremely important significance for the prevention and control of shaft deflection, for the rational design of PCP, and for the sustainability of mine production.
1. Introduction
Since 2002, China has built 71 shafts with soil thickness over 400 m, and 86% of these shafts have a soil thickness greater than half of their depths (Figure 1). The shaft is regarded as the throat of a mine, and its safety concerns the survival and sustainability of the entire mine. Since underground mining activities may have a certain effect on the shaft, the protection coal pillar (PCP) must be reserved. Unexpectedly, the shafts located in Guotun Coal Mine have experienced serious deflection disasters (Table 1 [1]). These deflections have caused the cage girders to deform and fall off. Considering that shaft deflection first occurred in China and is rare in the history of mining worldwide, the deflection mechanism must be clarified promptly. As is known to all, the shaft PCPs in European coal mines [2], as well as in China [3, 4], were mainly designed by considering the boundary angle or the movement angle. However, the angle values are similar under different mining conditions (Table 2). Combined with the actual situation of Guotun Coal Mine, the authors’ preliminary analysis indicated that the main reason for shaft deflection might be the insufficient size of PCP; that is, coal mining caused different displacements of the shaft at different depths. Based on the above analysis, mining subsidence theory was used to clarify the shaft deflection mechanism in deep soil strata. To study the problem of surface movement caused by mining, Litwiniszyn [5] proposed the stochastic medium theory. Later, Liu and Liao [6] further developed this theory into a probability integral method to predict mining subsidence. By combining the probability integral method and differential interferometric syntheticaperture radar (SAR) technique, Fan et al. [7] established a model for extracting the large deformation mining subsidence; Diao et al. [8] proposed a new monitoring method that could retrieve highly accurate 3D displacement of mining subsidence; and Yang et al. [9] obtained an approach for costeffective and accurate prediction of 3D mininginduced displacement under different extraction conditions. Using the probability integral method, Li et al. [10] calculated and simulated the influence of pillar mining subsidence on shaft safety. The key to the probability integral method is the solution of prediction parameters. For this purpose, Guo et al. [11] obtained more precise results by applying the theory of artificial neural network to this method, Zhang et al. [12] provided a highly effective method based on least squares support vector machine theory, and Han et al. [13] calculated the parameters that are suitable for deep soil strata by a back analysis. Regarding research on measuring and monitoring of subsidence, Zheng et al. [14] carried out monitoring and analysis of mining 3D deformation by multiplatform SAR images with the probability integral method. Considering that the above methods can predict only the final deformation, Knothe [15] built a time function model that can predict the dynamic surface subsidence based on the Mitscherlich growth law, which has widely been used and further improved [16–23]. Additionally, time functions were established by Liu et al. [24] based on the Harris model and by Xu and Li [25] based on the logistic growth model. However, due to many parameters, their models are relatively complex to solve.


Based on the Knothe time function, a more accurate time function was proposed. Combined with the probability integral method, a 3D dynamic prediction model was established to reveal the shaft deflection mechanism and to analyse the lining safety in deep soil strata. By comparing and analysing the theoretical and measured horizontal displacement values, the feasibility of the model was proved and the shaft deflection mechanism was clarified. The results in this paper are of great significance to the sustainability of mine production.
2. Models and Methods
2.1. Guotun Coal Mine
Guotun Coal Mine is located in Juye Coalfield in Shandong Province, with an annual production capacity of 2.4 million tons. All the shafts are located at the industrial square. The square is arranged in the coalfree zone. As of July 2018, the entire first mining area and part of the fourth area were already mined (Figure 2). The start times of the first and fourth mining areas were January 2010 and October 2016, respectively. Based on the comprehensive analysis of the distribution, sizes, and mining times of the areas, it was concluded that the shafts in Guotun Coal Mine were mainly affected by the mining activities of the first mining area. The first mining area is approximately 6300 m long in the north and south directions and 3700 m wide in the east and west directions, with an average mining depth of 770 m and average coal thickness of 3.97 m. The coal seam is nearly horizontal, with a general dip angle of 5∼8°. According to the criterion [3], the square PCPs were designed by the vertical section method, with the soil and rock movement angles of 45° and 70°, respectively. The distance from the square boundary to the pillar boundary is approximately 718 m. The main shaft is 140 m from the square boundary in the west.
2.2. Time Function Solution
2.2.1. Function Form
Assuming that the subsidence velocity at a certain moment is directly proportional to the difference between the dynamic subsidence value at this moment and the maximum subsidence value, Knothe [15] obtained a differential function as follows:where is the instantaneous subsidence value of a point on the surface at time t, is the maximum subsidence value of the point, and c is a time coefficient related to lithology.
With the initial condition , equation (1) can be derived as
The Knothe time function can be expressed as . With the characteristics of the velocity and acceleration of the surface subsidence, the subsidence process should be a roughly Sshaped curve with time. Then, the time function can be improved by hypothesis towhere λ is an analogous coefficient to c.
With proper values of c and λ, Figure 3 shows that equation (3) is more accurate than the Knothe function in predicting subsidence.
2.2.2. Parameters Solution
The 1308 working face was first mined in July 2011, and the partial observation points arranged along the strike are shown in Figure 4. The first survey was on July 29, 2011, and the second was on August 15, 2011: the maximum surface subsidence for two times was 6 mm and 340 mm (i.e., point Z_{16}), respectively. To facilitate easy calculation, assume λ = 1. When taking the first survey time as the initial time, based on the measured subsidence values at points Z_{1} and Z_{2}, the fitting curves and the corresponding parameter values are shown in Figure 5. The results showed that the values of c were feasible and it was reasonable to predict the subsidence of a point with λ = 1. If it is assumed that λ is a fixed value for other working faces, then c should be calculated according to the actual mining conditions.
Suppose H is the mining depth, β is the main influence angle, r_{0} is the main influence radius at the surface, s_{0} is the displacement of the inflection point, is the average mining rate, and l_{c} is the goaf critical size, then r_{0} = H/tanβ, and l_{c} = 2r_{0} + 2s_{0} (Figure 6). When the maximum subsidence is equal to , the opening has reached its critical size [26], and the mining time . By combining equations (2) and (3), c can be derived as
2.3. Model Solution for Stratal Movement
With coincident horizontal projection, a surface coordinate system and a mining coordinate system were established separately for each working face. When the parameters of the ith working face are denoted by subscript i, the coordinates of mining unit B and point A are (ξ_{i}, η_{i}, H_{i}) and (x_{i}, y_{i}, z_{i}), respectively (Figure 7). Point A′ and point B′ are the horizontal projections of point A and unit B at the surface. When the strike length of the ith working face is l_{i} and the tendency length is s_{i}, based on the probability integral method and equation (3), the subsidence value of point A caused by mining the ith working face at a certain moment is given aswhere t_{i} is the time of the ith working face from the initial mining to a certain time, is the maximum subsidence value, and r_{zi} is the main influence radius at depth z, which can be written aswhere m_{i} and α_{i} are the thickness and dip angle of the ith working face, respectively, q_{i} is the coefficient of stratal subsidence, and a is a constant.
According to the superposition principle, the subsidence of point A caused by mining n working faces at a certain moment can be given as
Equation (7) can be integrated aswhere erf is the probability integral function, which can be calculated by
As shown in Figure 7, φ_{i} is the rotation angle from the x_{i}axis to the y_{i}axis in the positive direction and then to the specified direction at point A′. Assume χ_{A} is the inclination of point A in the φ_{i} direction, and then
According to reference [27], the horizontal displacement u_{A} can be expressed aswhere b is the horizontal movement factor.
Equation (11), a 3D dynamic prediction function, describes the horizontal movement of the strata caused by mining multiface. The detailed expression can be seen from equation (12), and the prediction parameters can be inversed by measured data. For those indeterminate parameters, the empirical values can be referred to [3].where
2.4. Models Solution for Shaft Movement
Suppose the shaft depth is H_{0}, the initial coordinates of the point at the position of shaft axis are (x_{0i}, y_{0i}, z_{0i}) in the ith surface coordinate system. Substituting the coordinates into equations (8) and (12), the vertical and horizontal displacements can be calculated. The maximum value of z_{0i} should be smaller than the values of H_{0} and H_{i}. Considering that the shaft has a large slenderness ratio and a small lateral bending resistance, the displacement of the shaft is approximately equal to that of the soil at the same position. Taking the centre of the shaft head as the origin, then z = z_{0i} must be true for any point on the shaft axis, and the coordinate system was established to describe the shaft horizontal bending (Figure 8). Using and u(z) to represent the horizontal and vertical displacements of those points on the shaft axis at different depths, equations (8) and (12) can be rewritten as
Under the mining of multiple faces, equations (14) and (15) are the 3D dynamic prediction models for the shaft movement in the vertical and horizontal directions.
2.5. Analysis of Shaft Lining Stress
When the coal was asymmetrically mined relative to the shaft, the strata produced vertical and horizontal displacements due to mining effects, wherein the vertical displacement produced downward additional stress on shaft lining, while the horizontal movement produced bending stress. Assuming that the lining is a homogeneous and continuous elastic material, based on the principle of superposition, the vertical and horizontal lining stresses can be analysed separately. Since the most likely occurrence position of lining rupture was at the inner edge [28], the vertical stress at this position was first analysed. Obviously, the vertical lining stress σ_{z} includes the selfweight stress σ_{0z}, vertical additional stress σ_{c}(z, t), and bending stress σ_{f}(z, t). By assuming that the compressive stress is positive, σ_{z} can be expressed as
2.5.1. SelfWeight Stress
The selfweight stress σ_{0z} can be determined as follows:where γ_{c} is the unit weight of the lining.
2.5.2. Vertical Additional Stress
First, the vertical compressive strain ε_{c}(z, t) of the lining can be calculated by equation (14), which can be derived as
Then, the vertical additional stress σ_{c}(z, t) can be obtained with equation (18) aswhere E_{c} is the elastic modulus of the lining and “−” is due to the positive compressive stress.
2.5.3. Bending Stress
When the shaft can be regarded as an elastic beam, equation (15) can be expressed as the flexural formula of the beam. Based on the theory of the elastic beam, the bending moment of the lining can be obtained aswhere I_{z} is the moment of inertia and d and are the net diameter of the shaft and the outer diameter of the outer lining, respectively.
Therefore, the bending stress σ_{f}(z, t) at the inner edge of the lining can be obtained based on equations (20) and (21) as follows:
Compared to asymmetric mining, the strata produced double vertical displacement to the shaft under the symmetrical mining. The vertical lining stress only included the selfweight stress σ_{0z} and the vertical additional stress , here . Then, can be expressed as
2.6. Analysis of Lining Safety
According to the theory of a thick cylindrical lining, the circumferential stress σ_{θ} and radial stress σ_{r} at the inner edge of the lining can be calculated aswhere is the hydrostatic pressure and D_{n} is the outer diameter of the inner shaft lining.
According to Code [29], when equation (25) is satisfied, the inner edge of the lining is in an unsafe state under the influence of coal mining:where f_{c} is the uniaxial compressive design strength of concrete, is the compressive design strength of reinforcement, η is the improvement coefficient of the concrete strength under multiaxial stress, η is related to σ_{θ} and σ_{r}, and in general, η = 1.2. Finally, ρ_{min} is the minimum steel content.
3. Analysis of Deflection Mechanism and Lining Safety
3.1. Deflection Mechanism Analysis
3.1.1. Prediction Parameters Solution
The prediction parameters in Guotun Coal Mine were obtained as follows: s_{0i} = 0 [3] and q = 1.0 [30], and other parameters were taken based on an inversion analysis. The horizontal displacements of the main shaft were expressed by u_{m}. An inversion analysis was carried out by the measured values of u_{m} in June 2015 (curve I in Figure 9), and the inversion results were basically consistent with the measured values and were suitable for equation (15). Using the inversion parameters, the values of u_{m} in 2017 were calculated as seen from curve II in Figure 9. Compared with the measured values, their average difference was only 2.1 mm and smaller than the measurement error (8 mm). With high precision of anastomosis, the horizontal shaft displacement at different depths could be inverted by the dynamic prediction model. According to the calculation results, the prediction parameters obtained by inversion were a = 0.7, b = 0.66, and tan β = 0.50. The calculation parameters of the working faces and shaft coordinates are shown in Table 3.

3.1.2. Prediction Effect Verification
To further verify the feasibility of the prediction parameters, the auxiliary shaft was taken as a case. The horizontal displacements of the auxiliary shaft were expressed by u_{a}. The theoretical values and the measured values in 2015 (curve I) and 2017 (curve II) are shown in Figure 10. The errors in two times were only approximately 7% and 2% of the measured values at the shaft head, respectively, which indicated that the calculation results were still highly reliable considering the measurement error and the difference between the main shaft and auxiliary shaft.
3.1.3. Analysis of Results
The prediction model can accurately predict the shaft deflection in deep soil strata and calculate the horizontal displacement magnitudes at different depths under conditions of multiface mining. The calculated results agreed with the measured values. Based on the prediction parameters of inversion, it could be calculated that when H > 700 m, r_{0} = H/tan β > 1400 m. Therefore, the reason for the shaft deflection in Guotun Coal Mine was the insufficient size of PCP.
3.2. Lining Safety Analysis
Taking the main shaft as an example, the vertical theoretical displacement was calculated using the 3D dynamic prediction model (Figure 11). The calculation results showed that the maximum vertical displacement at two times was 37.3% and 35.4% of the maximum horizontal displacement, indicating that the mining effect was greater on the horizontal displacement of the shaft than on the vertical displacement in deep soil strata.
The azimuth angles of the main shaft at different depths in 2015 were calculated by the measured and theoretical displacement values, respectively. The azimuth angles did not change with depth in the soil section, and the two calculation results were basically consistent (Figure 12). When using the 3D dynamic prediction model to calculate the force of the shaft in the soil section, it was reasonable to consider the main shaft as a plane bending deformation.
To comprehensively analyse the influence of mining on lining safety, the vertical stress was analysed under symmetrical and asymmetric mining conditions in Guotun Coal Mine. According to the actual parameters of the lining structure of the main shaft in the soil section (Table 4), the vertical stress in 2015 and 2017 was calculated. The main shaft lining adopted an HRB335 steel bar, for which . When ρ_{min} = 0.2%, the stress values at different times are shown in Figure 13.

Figure 13 shows that, under the conditions of asymmetric mining, the vertical stress values are less than the yield strength values, and the main shaft remains safe. Shaft deflection mainly affects hoisting safety and has a smaller effect on lining safety. However, if the coal is mined symmetrically in Guotun Coal Mine, the mining will cause the lining to fracture. In conclusion, the problems are mainly caused by unreasonable design size of PCP in deep soil strata.
4. Conclusions
Combining the probability integral method with a new time function, a 3D dynamic prediction model was established to describe the shaft movement. Based on the shaft deflection in Guotun Coal Mine, the prediction parameters suitable for deep soil strata were obtained. By comparing and analysing the theoretical and measured values, the following conclusions were obtained:(1)The shaft deflection mechanism in deep soil strata was as follows: under the conditions of asymmetric mining, the mining caused the soil around the shaft to produce different displacement magnitudes at different depths, resulting in the vertical and horizontal displacements of the shaft within the range of mining influence. According to the current regulations, the shafts in Guotun Coal Mine are not located in the scope of mining influence, while they have experienced severe deflections. These deflections indicate that the stipulated design parameters for the shaft PCPs in current codes are not applicable to the deep soil strata.(2)The new time function was more accurate in predicting subsidence than the Knothe function. With λ = 1, the parameters of the prediction model were obtained by inverting the measured data. The prediction model can accurately predict the shaft movement, which provides an important reference for preventing shaft deflection and designing the shaft PCPs in deep soil strata.(3)Considering the current situation in Guotun Coal Mine, the grouting technology should be adopted to correct the deflection. For similar newly built shafts, the symmetrical or filling mining is recommended according to the actual situation, and lining structures should adapt themselves to a certain bending deformation and be installed with contractible devices.
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The authors would like to thank Guotun Coal Mine, Shandong University of Science and Technology, and China Coal Technology & Engineering Group Nanjing Design and Research Institute for their valuable data and technical support. This research was supported by the National Key Research and Development Program of China (Grant no. 2016YFC0600904).
References
 X. L. Cheng, S. G. Liu, and H. X. Zhang, “Technical summary of monitoring project of main and auxiliary shaft, derrick, head sheave and and winch in guotun coal mine,” Tech. Rep., Shandong Science and Technology University Technology Development Company, Qingdao, China, 2018, Research Report. View at: Google Scholar
 I. H. Kratzsch, “Mining subsidence engineering,” Environmental Geology and Water Sciences, vol. 8, no. 3, pp. 133–136, 1986. View at: Publisher Site  Google Scholar
 China Coal Industry Publishing House, State Bureau of Coal Industry, Criterion for Coal Pillars and Coal Mining in Buildings, Water Bodies, Railways and Main Shafts, China Coal Industry Publishing House, Beijing, China, 2017.
 R. L. Zhang, G. W. He, and Z. Li, Mining Engineering Design Manual (Part I), China Coal Industry Publishing House, Beijing, China, 2003.
 T. Litwiniszyn, “Application of the equation of stochastic processes to mechanics of loose bodies,” Archiwum Mechaniki Stosowanej, vol. 8, no. 4, pp. 393–411, 1956. View at: Google Scholar
 B. C. Liu and G. H. Liao, The Law of Surface Movement in Coal Mine, China Architecture & Building Press, Beijing, China, 1965.
 H.d. Fan, W. Gu, Y. Qin, J.q. Xue, and B.q. Chen, “A model for extracting large deformation mining subsidence using DInSAR technique and probability integral method,” Transactions of Nonferrous Metals Society of China, vol. 24, no. 4, pp. 1242–1247, 2014. View at: Publisher Site  Google Scholar
 X. Diao, K. Wu, D. Hu, L. Li, and D. Zhou, “Combining differential SAR interferometry and the probability integral method for threedimensional deformation monitoring of mining areas,” International Journal of Remote Sensing, vol. 37, no. 21, pp. 5196–5212, 2016. View at: Publisher Site  Google Scholar
 Z. F. Yang, Z. W. Li, J. J. Zhu et al., “An extension of the InSARbased probability integral method and its application for predicting 3D mininginduced displacements under different extraction conditions,” IEEE Transactions on Geoscience and Remote Sensing, vol. 55, no. 7, pp. 3835–3845, 2017. View at: Publisher Site  Google Scholar
 P. Li, Z. Tan, and L. Yan, “A shaft pillar mining subsidence calculation using both probability integral method and numerical simulation,” Computer Modeling in Engineering & Sciences, vol. 117, no. 2, pp. 231–250, 2018. View at: Publisher Site  Google Scholar
 W. B. Guo, K. Z. Deng, and Y. F. Zou, “Artificial neural network model for predicting parameters of probabilityintegral method,” Journal of China University of Mining and Technology, vol. 33, no. 3, pp. 322–326, 2004. View at: Google Scholar
 H. Zhang, Y.j. Wang, and Y.f. Li, “SVM model for estimating the parameters of the probabilityintegral method of predicting mining subsidence,” Mining Science and Technology (China), vol. 19, no. 3, pp. 385–388, 2009. View at: Publisher Site  Google Scholar
 J. Han, J. Zou, C. Hu, and W. Yang, “Study on size design of shaft protection rock/coal pillars in thick soil and thin rock strata,” Energies, vol. 12, no. 13, p. 2553, 2019. View at: Publisher Site  Google Scholar
 M. Zheng, K. Deng, H. Fan, and J. Huang, “Monitoring and analysis of mining 3D deformation by multiplatform SAR images with the probability integral method,” Frontiers of Earth Science, vol. 13, no. 1, pp. 169–179, 2019. View at: Publisher Site  Google Scholar
 S. Knothe, “Effect of time on formation of basin subsidence,” Archives of Mining and Steel Industry, vol. 1, no. 1, pp. 1–7, 1953. View at: Google Scholar
 Q. F. Hu, X.M. Cui, G. Wang, M.R. Wang, Y.X. Ji, and W. Xue, “Key technology of predicting dynamic surface subsidence based on Knothe time function,” Journal of Software, vol. 6, no. 7, pp. 1273–1280, 2011. View at: Publisher Site  Google Scholar
 P. Polanin, “Application of two parameter groups of the Knothe–Budryk theory in subsidence prediction,” Journal of Sustainable Mining, vol. 14, no. 2, pp. 67–75, 2015. View at: Publisher Site  Google Scholar
 B. Zhang, Y. Li, N. Fantuzzi et al., “Investigation of the flow properties of CBM based on stochastic fracture network modeling,” Materials, vol. 12, no. 15, p. 2387, 2019. View at: Publisher Site  Google Scholar
 X. Zhu, G. Guo, J. Zha, T. Chen, Q. Fang, and X. Yang, “Surface dynamic subsidence prediction model of solid backfill mining,” Environmental Earth Sciences, vol. 75, no. 12, p. 1007, 2016. View at: Publisher Site  Google Scholar
 Z.G. Zhang, S.G. Cao, Y. Li, P. Guo, H. Y. Yang, and T. Yang, “Effect of moisture content on methane adsorption and desorption induced deformation of tectonically deformed coal,” Adsorption Science and Technology, vol. 36, no. 910, pp. 1648–1668, 2018. View at: Publisher Site  Google Scholar
 Z. Q. Chang and J. Z. Wang, “Study on time function of surface subsidence—the improved knothe time function,” Chinese Journal of Rock Mechanics and Engineering, vol. 22, no. 9, pp. 1496–1499, 2003. View at: Google Scholar
 Y. F. Gao, J. Y. Jia, B. Li, and Q. S. Zhang, “The attenuation function of surface subsidence and stability analysis due to mining,” Journal of China Coal Society, vol. 34, no. 7, pp. 892–896, 2009. View at: Google Scholar
 Y. Zhao, S. Cao, Y. Li et al., “The occurrence state of moisture in coal and its influence model on pore seepage,” RSC Advances, vol. 8, no. 10, pp. 5420–5432, 2018. View at: Publisher Site  Google Scholar
 X. Liu, J. Wang, J. Guo, H. Yuan, and P. Li, “Time function of surface subsidence based on Harris model in minedout area,” International Journal of Mining Science and Technology, vol. 23, no. 2, pp. 245–248, 2013. View at: Publisher Site  Google Scholar
 H. Z. Xu and X. H. Li, “Time function of surface subsidence based on Logistic growth model,” Rock and Soil Mechanics, vol. 26, pp. 151–153, 2005. View at: Google Scholar
 Q. Hu, X. Deng, R. Feng, C. Li, X. Wang, and T. Jiang, “Model for calculating the parameter of the Knothe time function based on angle of full subsidence,” International Journal of Rock Mechanics and Mining Sciences, vol. 78, pp. 19–26, 2015. View at: Publisher Site  Google Scholar
 G. Q. He, L. Yang, G. D. Ling, C. F. Jia, and D. Hong, Mining Subsidence Theory, China University of Mining and Technology Press, Xuzhou, China, 1991.
 J. Han, J. Zou, W. Yang, and C. Hu, “Mechanism of fracturing in shaft lining caused by highpressure pore water in stable rock strata,” Mathematical Problems in Engineering, vol. 2019, Article ID 5234642, 9 pages, 2019. View at: Publisher Site  Google Scholar
 China Architecture and Building Press, GB 500102010, Code for Design of Concrete Structures, China Architecture and Building Press, Beijing, China, 2015.
 Shandong Luneng Group Heze Coal Power Development Co., Ltd, “Research report on surface movement and deformation law of 1308 working face in Guotun Coal Mine,” Guotun Coal Mine, Heze, China, 2012. View at: Google Scholar
Copyright
Copyright © 2019 Jihuan Han et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.