Abstract
The model pump is a high-speed, high-power pump designed to achieve rapid mine flooding rescue. This study conducted experiments to investigate pump performance curves, including head, efficiency, and power for the following six different rotation speeds: 3000, 3600, 4200, 4800, 5400, and 6000 rpm. Then, the numerical simulation method based on computational fluid dynamics commercial code Ansys was used to present the internal flow of the pump for the six different rotation speeds through steady and unsteady analyses. Results show that the numerical results agree well with experimental data. The designs of outlet and inlet angles of the impeller match each other well at high rotation speeds. The pressure pulsation coefficient Cp in the impeller and the diffuser channel remain constant at the same monitor point under different rotation speed conditions. The varying trend of the pressure-augmented coefficient ΔP indicates that, with the increase in rotation speed, the effect on pressure rise induced by the back part of the impeller channel is more evident than that by the front part. Also, the main frequency components of ΔP are concentrated on the region with low frequency. Moreover, the rotation speed has no significant effect on ΔP in the diffuser region. This study provides effective guidance and valuable reference for the design of high-speed, high-power pumps.
1. Introduction
Today, high-power, high-head centrifugal pumps worldwide work under the rotation speed of 3000 rpm. This low speed limits the use of these pumps in various fields. Especially in mine flooding accidents, the speed of rescue pumps has become the key technology to achieving rapid rescue to save fortune and lives [1–3]. With the development of frequency conversion technology, the high-speed feature is becoming one of the developing trends of high-power centrifugal pumps. This study presents a new rescue pump that fills the gap between rescue pumps and the field of high-speed, high-power centrifugal pumps.
Owing to the complex conditions of mine flooding, such as changes in flooding water level and variations in pipe diameter, the flexible operation of rescue pumps requires running at different flows and speeds to achieve multioperation requirements. The pump speed has a considerable influence on pump performance, which is a key index of pump standard. In recent years, many researchers have been inspired to conduct studies on the effect of pump speed on the characteristics of pumps. Draghici et al. [4] presented experimental investigations performed for a centrifugal pump which operated at two speeds. Georgescu et al. [5] studied the characteristic pump curves with two speeds based on EPANET. Rakibuzzaman et al. [6] researched multistage centrifugal pump performance characteristics based on inverter-controlled variable speed. Sha et al. [7] obtained the internal flow of axial flow pump at variable speeds from numerical simulations and experiments.
In addition, the application of computational fluid dynamics (CFD) on pump design has been proven to be a very useful tool in the flow analysis inside pumps, which can be used in numerical simulations to obtain flow field performance inside the pump. With the advantages of CFD in terms of time and cost, an increasing number of researchers have carried out CFD for the design and analysis of pumps [8–15]. For example, in the field of pressure pulsation, Zhang et al. [16] investigated unsteady flow in a centrifugal pump with a special slope volute under various conditions through simulation and experimental methods. Spence et al. [17, 18] suggested that monitoring pulsations at the top center of a centrifugal pump could be used as an indication of internal flow and found that the time variation analysis of pressure in a complete double-entry, double-volute centrifugal pump could help increase component life and reduce vibration. Luo et al. [19] studied the relationship among pressure pulsations, operation mode, and water head for a low head bidirectional tidal bulb turbine.
However, most of the preceding studies were conducted on low- or normal-speed pumps, whereas studies on the effect of rotation speed on high-speed and high-power pumps remain scarce. Consequently, this study investigates a high-speed and high-power mine rescue pump with six different frequencies ranging from 3000 rpm to 6000 rpm using experimental and numerical simulation methods. Pressure distribution and pressure pulsation in the pump are also analyzed. The results of this study may provide guidance for the development of a series of high-speed and high-power pumps.
2. Hydraulic Design
The GPQ 200-300 is one type in the high-speed rescue pump series presented in this study. The main design parameters are as follows: flow = 200 m3/h, head = 300 m, specific speed = 71.6, and speed = 6000 rpm. As main flow passage components, the impeller and diffuser have significant effects on pump performance. The parameters of the impeller and diffuser vanes are shown in Table 1, and the models are shown in Figure 1.
(a)
(b)
(c)
(d)
3. Experiment System
Hydraulic testing of the unit was executed at Hefei Hengda Group, China (Figure 2). Testing was carried out with an output range of 3000–6000 rpm from the frequency converter. The flow rate was altered by the control valve. The output pressure was obtained by the pressure collection at the pump outlet through the WT200 intelligent pressure sensor with 0.1% precision. The flow rate was measured at the outlet pipe by means of a LWGY-200A turbine flow meter with a measurement error of 0.5%. Experimental measurements were performed for the following six speeds: 3000, 3600, 4200, 4800, 5400, and 6000 rpm with 50, 60, 70, 80, 90, and 100 Hz, respectively. The flow rates ranged from 50% to 120% of the nominal discharge based on pump affinity law corresponding to respective frequencies.
(a) Design diagram
(b) Test site
Figure 3 shows that the trends of performance curves follow the same change law. With the increase in output frequency of the frequency converter, the crew power and head at monitor flow points demonstrate a regular increase. Simultaneously, the crew efficiency curves indicate regular deviations. Regarding head curves, the head declines with the growth of flow before exceeding the design flow point. However, flow growth demonstrates an evident decrease in the head after the design point due to the increasing mismatch in characteristics of flow passage components under heavy flow conditions. Furthermore, a small frequency leads to an evident decrease when the flow rate exceeds the design point, which may be attributed to the following reason: the flow components are designed for the working condition under 6000 rpm, which may not perfectly match the flow characteristic with the reduction of rotation speed. In particular, impact loss and flow blockage may be large under big flow points with low rotation speed, thereby leading to the decrease in head curves. For the crew efficiency curves, the best crew efficiency points (BEP) all appear near the design point with six different frequencies. However, the BEP values with sundry rotating velocities show an irregular characteristic. The BEP values are minimum for 6000 rpm and 4800 rpm and maximum for 3000 rpm and 5400 rpm. Considering that the efficiency in Figure 3 is crew efficiency, which includes motor and frequency converter efficiency, this phenomenon may result in motor and frequency converter characteristics with different frequencies.
(a) Head curves
(b) Efficiency and power curves
4. Numerical Modeling
Figure 4 shows the flowchart of numerical simulation in this paper. Especially in the validation of the simulation results, the head, the efficiency, and the power obtained from the experiment and the simulation results, based on different turbulence models, are used for comparison to get the best models. In addition, the efficiency and power obtained from the experiment are crew efficiency () and crew power (). However, the corresponding parameters obtained from the numerical simulation are hydraulic efficiency () and hydraulic power (). The relationship between and could be seen in (1), and the relationship between Pc and could be seen in (2). Besides, in the following, all the efficiency and power mentioned, which appear without special note, are hydraulic efficiency () and hydraulic power ().
where is the conversion factor, is the frequency converter efficiency, is the motor efficiency, is the mechanical efficiency, is the volumetric efficiency, and is the hydraulic efficiency. Furthermore, in the paper, the head and the efficiency are defined as follows:
where is the outlet pressure, is the inlet pressure, is density, is gravity acceleration, is flow rate, and is head.
4.1. Computational Domain
The inlet, impeller, diffuser vane, gaps of front and back covers, and outlet form the computation domain, which benefits from a high-quality fully structured grid. In addition, considering flow stability and simulation accuracy, the inlet model is extended to five times the diameter of the inlet pipe. Advantages of the structured grid include small truncation error and improved convergence to simulate flow in the rotation machinery. In addition, the 3D model and the structured grid are obtained by commercial codes Creo 2.0 and Icem 16.0, respectively (Figures 5 and 6).
4.2. Mesh Independence
A sensitivity study considering five levels of mesh refinement is conducted to assess the influence of mesh density on the simulation results. Table 2 shows five groups of the entire computation domain with different mesh numbers, and Figure 7 presents the head and efficiency computation for different mesh densities. Simulation accuracy and calculation time are also assessed to assure grid-independent solutions. Although the total number exceeds 1769 thousands, the simulation results appear to remain constant. Then, the third method in Table 2 is selected for the final simulation. The resulting computational domain has values below 50, thereby implying adequate coverage of the critical regions.
4.3. Boundary Condition
The simulation in this study is carried out based on the commercial code Ansys CFX 16.0. For steady simulation, the pressure is set to inlet condition with a value of 1 atm, and the mass flow is set to outlet condition with a value below the current flow. The impeller domain is also set to rotation with a value corresponding to the current power frequency, and the other components are set to stationary. In addition, the frozen rotor is set to the interference between the rotating and stationary domains. Referring to the machining precision standards in China, wall roughness is set to 25 µm to obtain exact results [20]. For transient numerical simulation, the steady numerical result obtained by the steady simulation is used for initial conditions. For time duration, time step is set to a series of 3° changes in the impeller circumference corresponding to different rotation speeds. Total time is set to 10 full impeller revolutions, which means 1200 times of . Last, the maximum residual of the convergence criterion is set to 1 × 10−5 in accordance with the high-resolution scheme for time discretization.
4.4. Turbulence Model
Considering the reliability and feasibility of dynamic calculations, the most widely used two-equation models in turbomachine numerical simulations include the k-ω model, k-ε model, and SST k-ω model. Simulation results based on the three different turbulence models are selected for comparison with experimental data with 6000 rpm to determine the most suitable turbulence model for the numerical calculation.
However, hydraulic efficiency and hydraulic power acquired from the simulation differ from crew efficiency and crew power mentioned in the experimental data. Crew efficiency could be transformed into hydraulic efficiency from (1). Meanwhile, crew power could be obtained from hydraulic power from (2). Therefore, efficiency and power indicator are applied to validate simulation accuracy with the head indicator.
Figure 8 shows the comparative results of the numerical simulation and the experimental data. Despite a few limitations at certain points, the performance curves of the simulation present a tendency consistent with that of experimental data. These limitations may be attributed to the following reasons. (1) The flow near the wall may be affected by the rough surface, and the resulting machining precision may fail to achieve the standard. (2) An ideal condition in the numerical simulation is when the flow in the inlet is considered uniform, which may differ from the actual flow condition. Hence, this finding indicates that the SST k-ω model has been successfully applied to predict the performance of the tested pump, thereby obtaining results closest to the experimental data among the three turbulence models. Furthermore, considering wall function and turbulent vortex, the SST k-ω model is suitable for flows, such as complex shear flows and high-speed rotation flows near the wall [21, 22]. Thus, adopting the SST k-ω model with standard wall function in this study is appropriate.
(a) Head curves
(b) Efficiency curves
(c) Power curves
5. Result Analysis
All numerical simulation results are acquired from the nominal flow with corresponding rotation speeds based on the affinity law to analyze and compare the performance and flow field of the pump under different rotation speeds.
5.1. Steady Analysis
Figure 9 is the velocity streamline distribution on the middle section plan under the design flow corresponding to the different rotation speeds. Based on Figure 9, the streamline is well distributed and no evident flow separation exists in the impeller channel for all six cases. In addition, the maximum velocity vector appears in the regions with interference of the impeller and diffuser and the location near the diffuser inlet. Furthermore, the range of regions expands with the increase in rotation speed. This observation demonstrates that the outlet and inlet angle designs of the impeller match each other well for high rotation speed, thereby reducing the impact loss and back flow near these regions. The streamline in the diffuser channel has asymmetric characteristics due to the effect of dynamic coupling calculation between the impeller and diffuser. In addition, although a small region of vortex and back flow appear in the diffuser channel, this region does not have obvious influence on the overall diffuser performance.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
Figure 10 shows the pressure distribution on the blade-to-blade surface of the impeller with six different rotation speeds. This figure also shows that the maximum pressure gradually increases with the increase in rotation speed, which could be expressed for the reason that the head increases with the rotation speed. In terms of minimum pressure, all six rotation speeds follow the rule that high speed leads to small pressure, which may be attributed to the low-pressure area for the high rotation speed. Especially in terms of the pressure counter, the minimum pressure appears in the back of the blade inlet. When the liquid begins to flow into the rotation domain, the flow begins to rotate with a high speed, which leads to the instant increase in kinetic energy corresponding to the instantaneous pressure reduction. In addition, in the back of impeller, flow separation as a response to the instant flow condition change of liquid near the inlet region is prone to appear, thus leading to the minimum pressure. The discrepancy of two adjacent rotation speeds with maximal pressure magnitude mainly becomes large when the rotation speed increases. On the contrary, the minimal pressure disparity is small. The impeller is designed for the 6000 rpm work condition. Therefore, proximity to the 6000 rpm leads to improved impeller characteristics.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
5.2. Unsteady Analysis
Owing to transient initial characteristics, the transient data used to monitor the flow in this study are all from the 6th to 10th period at the design flow () with different rotation speeds. During calculations, pressure data are saved at recording points in the middle streamline in the impeller and diffuser flow passages. As shown in Figure 11, these points are evenly distributed according to the same radial spacing. Among these recording points, yl-1, yl-2, yl-3, yl-4, yl-5, and yl-6 represent monitors in the impeller. Similarly, dy-1, dy-2, dy-3, dy-4, dy-5, and dy-6 represent the monitors in the diffuser.
5.2.1. Pressure Coefficient
Pressure pulsation is caused by comprehensive effects, such as rotor-stator interaction, turbulence, and recirculation. The pressure pulsation energy is evaluated alongside the effects of amplitude errors to intuitively and clearly compare pressure pulsation. Subsequently, a pressure pulsation coefficient is defined to (5) as follows:where represents static pressure at the monitor point, represents mean pressure, is density, and is velocity at the impeller outlet.
Figure 12 shows time histories of the pressure pulsation coefficient in the impeller with different rotation speeds. From the comparison of Cp at different monitor points in the impeller under the same rotation speed conditions, the mean magnitude of Cp is found to increase with the increase in radial location. This phenomenon could be attributed to the relationship between pulsation strength and tangential speed. Generally, with the increase in tangential speed due to radial growth, liquid in the impeller flow channel acquires more kinetic energy near the outlet than near the inlet. Furthermore, the peak magnitude of Cp appearing in the impeller channel allows for the possible speculation that impact loss manifests in the region near the impeller inlet and evolves along the entire impeller passage. By contrast, when rotation speed varies, identifying the difference between peak value and distribution trend of Cp at a certain monitor point in the impeller is difficult. This difficulty may be related to the flow characteristics in the impeller channel.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
Figure 13 provides time histories of pressure pulsation coefficient in the diffuser with different rotation speeds. Evidently, the peak magnitude of Cp only appears at the monitor point dy-1. When the liquid enters the diffuser channel, the flow has significant instability due to the rotor-stator interaction. At the dy-1 point, the peak value simultaneously appears corresponding to the impeller monitor points. Therefore, the flow condition at the inlet region of the diffuser is evidently affected by the outlet liquid status of the impeller tip. Subsequently, the relationship between the peak magnitude of Cp and location of monitor points under the same rotation speed condition could be described as follows: with the increase in radial distance from the axis, the peak and mean magnitude of CP decrease step by step. When liquid passes the diffuser passage, the flow status becomes stable and smooth due to the diffuser diffusion effect. No significant difference exists in the peak magnitude of Cp at the same monitor point from 3000 rpm to 6000 rpm. The variation of Cp indicates the possibility of a similar flow structure in the diffuser.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
5.2.2. Pressure-Augment Coefficient
A dimensionless parameter and pressure-augmented coefficient is introduced to clarify the effect of each impeller part on flow characteristics during the entire work process, as shown in (7):
where and present pressure at the monitor points. For example, ΔPyl-23 indicates the difference between Pyl-3 and Pyl-2 divided by Pyl-2.
Time histories of pressure-augmented coefficient in the impeller with different rotation speeds are presented in Figure 14. The pressure of monitor Pyl-1 is not discussed in this section due to the complex flow results in the transition process from stator to rotor. Regarding the mean value of ΔP, ΔPyl-45 and ΔPyl-56 present similar variation range with different rotation speeds because the liquid experiences similar pressure increase in the back region of the impeller channel, which does not depend on the rotation speed. Furthermore, their mean magnitude is smaller than that of ΔPyl-23 and ΔPyl-34. Therefore, the effect on the pressure increase induced by the back part of the impeller channel is more evident than that by the front part. Also, with the increase in radial distance, the liquid acquires large kinetic energy due to the large tangential speed. Consequently, this result leads to a smaller pressure increase gradient compared with that in the region near the inlet.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
Notably, when the rotation speed increases, the trends of ΔPyl-23 and ΔPyl-34 significantly change. First, the amplitudes and averages of ΔPyl-23 and ΔPyl-34 evidently increase. This finding can be attributed to the pressure increase in the front region of the impeller channel, which is related to the initial liquid velocity speed as the working fluid enters the impeller region. In addition, the magnitude increase in ΔPyl-23 is larger than that of ΔPyl-34 with the increase in rotation speed. This finding may be attributed to the relationship between the location in the impeller front part and the effect on the result. Proximity to the impeller inlet possibly affects the rotor-stator interaction near the impeller inlet.
Finally, when the rotation speed exceeds 3600 rpm, the mean value of ΔPyl-23 becomes larger than that of ΔPyl-34. Therefore, for ΔP in the impeller channel, different regions of the impeller front part play various roles in increasing pressure with different rotation speeds. Moreover, the contribution to the pressure rise near the impeller region of the inlet is evident with the increase in rotation speed.
Figure 15 shows the frequency spectrum of pressure-augmented coefficient in impeller to comprehensively understand the role of each impeller part in the progress of pressure-augmented pulsation under different rotation speed conditions. For each rotating speed, the axis frequency is nearly equal to the output frequency from the frequency converter. Hence, corresponds to 50, 60, 70, 80, 90, and 100 Hz.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
The Fast Fourier transform (FFT) method is adopted to transform the time domain signals into frequency domain signals to analyze the pressure spectra. Figure 15 evidently shows that the main frequency components of in the low-frequency area are concentrated in the low-frequency range (, corresponds to the diffuser passing frequency). Furthermore, in terms of largest frequency, the value of the remaining regions is , except for ΔPyl-23. This phenomenon is due to the appearance of back flow and impact loss in the gap of impeller outlet and diffuser inlet. The larger amplitude of ΔPyl-56 at than that of ΔPyl-45 could be attributed to the region near the impeller outlet, which is more likely affected by the back flow appearance at the gap. However, the amplitudes at of ΔPyl-34 and ΔPyl-23 are larger than that of ΔPyl-56. Furthermore, the main frequency of ΔPyl-23 is . These results could be attributed to the pressure pulsation in the region between monitor points yl-2 and yl-3, which are more likely affected by the axis rotation. The large amplitude at of ΔPyl-23 and ΔPyl-34 could be clarified by the mutual gain effect between the vortex flow emerging in the middle area of the impeller passage and the inlet turbulent flow due to the rotating axis.
Figure 16 presents time histories of the pressure-augmented coefficient in diffuser with different rotation speeds to investigate the change of ΔP in each region of the diffuser. ΔPdy12 has similar apex to ΔPdy-23 among the five regions for all rotation speeds. Simultaneously, the maximal magnitude of ΔPdy-34 is larger than that of ΔPdy-45, and ΔPdy-56 has a minimal value. Therefore, pressure augmentation for diffuser channel flow mainly occurs in the front region near the diffuser inlet. This finding is also in accordance with the rapid decrease in velocity due to the impact loss at the diffuser inlet. When the liquid passes through the front diffuser region, the flow becomes smooth and the velocity loss is reduced. Consequently, the pressure increase is not evident in the back region of diffuser. Furthermore, the trends of each ΔP are nearly uniform under different rotation speed conditions. The effect on the flow status due to the diffuser channel structure has no obvious relationship with rotation speed.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
Figure 17 is the frequency spectrum of pressure-augmented coefficient in diffuser with different rotation speeds. After FFT, the main frequency for each rotation speed is focused on the low-frequency region. In this region, the peak value appears at 3fn location, which means that the influence of rotating axis is gradually highlighted as the fluid enters the vane. Furthermore, among the five regions, ΔPdy-23 has the largest peak amplitude due to the vortex and separated flow produced in this region from the weakened bondage capability of the fluid on the diffuser blade surface. In high-frequency regions, the order of each blade passing frequency is dominant. Therefore, this phenomenon can be attributed to the blade wake when the impeller blade sweeps across the stationary diffuser blade.
(a) 3000 rpm
(b) 3600 rpm
(c) 4200 rpm
(d) 4800 rpm
(e) 5400 rpm
(f) 6000 rpm
6. Conclusion
The influence of rotation speed on performance and flow characteristics of high-speed, high-power rescue pumps is studied in this paper via experiments and CFD. The agreement between the experimental data and numerical simulation results provides the validation of the CFD method. The results show that the design and structure of the main flow component (impeller and diffuser) are suitable for six different rotation speeds varying from 3000 rpm to 6000 rpm.
In terms of the steady flow field, the velocity field and pressure distribution at the design flow with different rotation speeds show that the outlet and inlet angle designs of impeller match each other well for the high rotation speed by reducing the impact loss and back flow occurrence near these regions.
Furthermore, several interesting features of pressure pulsation inside the impeller and diffuser flow channel are captured well based on the simulations. Under a certain rotation speed, the magnitude of Cp presents an increasing trend from the inlet to the outlet of impeller. On the contrary, the values present a decreasing trend from the diffuser inlet to the outlet. However, the rotation speed has no significant impact on the distribution of Cp for the same monitor points.
In addition, a pressure-augmented coefficient ΔP is introduced to compare the pressure increase ratio at each region of the impeller and diffuser with different rotation speeds. With the increase in rotation speed, the variation of ΔP at the same impeller region indicates that the effect on the pressure increase applied by the back part of the impeller channel is larger than that of the front part, and the main components of ΔP are concentrated on the ranges with low frequency . However, ΔP in the diffuser channel has no evident relation with rotation speed.
The present work not only provides significant guidance in the design of high-speed, high-powered rescue pumps but also contributes to an enhanced understanding of flow characteristics in the new types of pumps, such as velocity field, pressure field, and pressure pulsation due to rotation speeds.
Appendix
| : | Specific speed: |
| : | Throat height of diffuser (mm) |
| : | Outlet height of diffuser (mm) |
| : | Conversion factor |
| : | Outlet width of impeller (mm) |
| : | Pressure pulsation coefficient |
| : | Hub diameter (mm) |
| : | Inlet diameter of impeller (mm) |
| : | Outlet diameter of impeller (mm) |
| : | Base diameter of diffuser (mm) |
| : | Outer diameter of diffuser (mm) |
| : | Outlet diameter of diffuser (mm) |
| : | Axis frequency (Hz) |
| : | Gravity acceleration (m/s2) |
| : | Head (m) |
| : | Design head (m) |
| : | Rotational speed (rpm) |
| : | Pressure (Pa) |
| : | Inlet pressure (Pa) |
| : | Outlet pressure (Pa) |
| : | Hydraulic power (kW) |
| : | Crew power (kW) |
| : | Average pressure fluctuation (Pa) |
| : | Pressure-augmented coefficient |
| : | Flow rate (m3/s) |
| : | Design flow rate (m3/s) |
| : | Velocity of impeller outlet (m/s) |
| : | Blade number of impeller |
| : | Blade number of diffuser |
| : | Inlet angle of diffuser (°) |
| : | Outlet blade angle of impeller (°) |
| : | Wrap angle (°) |
| : | Hydraulic efficiency (%) |
| : | Crew efficiency (%) |
| : | Frequency converter efficiency (%) |
| : | Mechanical efficiency (%) |
| : | Motor efficiency (%) |
| : | Volume efficiency (%) |
| : | Density (kg/m3). |
Data Availability
The authors declare that materials described in the article, including all relevant raw data, will be freely available to any scientist wishing to use them for noncommercial purposes, without breaching participant confidentiality.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The authors disclose receipt of the following financial support for the research, authorship, and/or publication of this article: this work was partly supported by the grants of National Natural Science Foundation of China (no. 11602097, no. 51509111, and 2013BAK06B02), Top-notch Academic Programs Project of Jiangsu Higher Education Institutions (PPZY2015A086), Industry-University-Research Cooperation Project of Jiangsu Province (BY2016023-04), Jiangsu Higher Education Excellent Scientific and Technological Innovation Team Project (New Material Forming Processing Technology), and the Priority Academic Program Development (PAPD) of Jiangsu Higher Education Institutions.