Abstract
Diffuser’s aerodynamic performance is crucial for the centrifugal compressors, while at present the universal principles for the optimization design of the vaned diffusers are still not available. In this paper, three vaned diffusers with different inlet Mach numbers were numerically studied in order to explore new design principles for the centrifugal compressor vaned diffusers. It proved that there are practical and effective design principles for the vaned diffuser optimizations, the performance of the vaned diffuser can be improved by carefully control of two aerodynamic parameter distributions: Tangential velocity (Vt) and Meridional velocity (Vm). The vaned diffusers with subsonic, transonic and supersonic inlet conditions were optimized with the new design principles, and the peak efficiencies were increased by 4.23%, 2.15% and 2.59%, respectively. The stage pressure ratios were increased by 3.36%, 1.39% and 6.49%, respectively, and their surge margins were also improved substantially. Finally, since the Vt and Vm could affect each other during the optimization process, an interactive optimization design procedure was also presented in this paper in order to accelerate the optimization process.
1. Introduction
Centrifugal compressors have the advantages of compact structure, high pressure ratio, good reliability, and wide operating range. They are widely used in small and medium-sized aircraft engines and gas turbines. Better aerodynamic performances and wider operating range are desired for the centrifugal compressors. Researchers have been conducting related studies all the time. HA Klassen [1] experimentally studied a pressure ratio of 6: 1 centrifugal compressor by testing a swept-back impeller. C Osbosne [2] performed a mechanical design and aerodynamic analysis of an 8: 1 centrifugal compressor. H Krain [3] carried out the numerical and experimental researches on the transonic flow fields of the centrifugal compressors, identified areas of large flow loss in centrifugal compressors, and improved compressor performance by modifying the impeller thickness and blade angle. M Voges [4] used PIV to conduct experimental measurement and analysis of the transonic flow field. With the deepening understands of the internal flow field of the compressor, the design method of the impeller blade has been developed from two dimensions to three dimensions. The centrifugal compressor impeller efficiency can be above 90%. While, the stage performance of the centrifugal compressor has not been greatly improved due to the unsatisfying performance of the diffusers [5].
The diffuser could greatly affect the stage performance of the centrifugal compressor. The impeller outlet airflow velocity is very high and unsteady, the height of the centrifugal impeller outlet is very small, and the distance between the diffuser inlet and the impeller outlet is very small, there is strong unsteady interaction between the centrifugal impeller and the diffuser, the working conditions for the diffuser are harsh. All these factors will greatly reduce the efficiency and stability of the centrifugal compressor. The design of compact and efficient diffuser is very challenging, and it is also a major technical obstacle for the high-pressure-ratio centrifugal compressors [6, 7]. At present, there are many researches on the design of the diffusers, such as low solidity vanes, half guide vanes, pipe diffuser and leading edge slotted design, these could improve the performance of the diffusers [8–12].
Blade optimization design is time-consuming work, and it is even more difficult without a practical and universal design principle. Nowadays, a few articles proposed diffuser geometry design methods but few articles proposed universal design principles. To investigate the sensitivity of how compressor performances change as the vaned diffuser geometry is varied, Hong-Won Kim [13] studied and compared three kinds of vaned diffusers. It found that the NACA 65 diffuser geometry satisfied a wider operating range and had a better performance. For diffuser channel design, Reneau [14] proposed that excessive divergence angle would cause a rapid change of flow velocity, and there would be flow separation in the diffuser channel, resulting in failure of the expansion capacity. In order to make the design method more general, Zangeneh M [15] proposed an inverse design method. Although it is difficult to give universal design principles, it is feasible to design the diffuser by using the blade load distribution. Wang [16] and Engeda [17] both use intelligent methods to optimize three-dimensional blades, while it is somewhat difficult for the current industrial production and application.
For the conventional diffuser optimization design, direct geometry modifications are the common choice, and the geometric parameters are usually modified without directions. Therefore, it is necessary to explore new design principles for higher performance and more convenience. In this paper, by controlling the aerodynamic parameter distributions of Vt and Vm along the flow passage, the universal design principles of the vaned diffusers were investigated and verified for three centrifugal compressors, and significant improvements of the stage performance have been achieved.
2. Case Descriptions
In this paper, for exploring the universal design principles of the vaned diffusers with different Mach numbers, 3 diffusers of subsonic, transonic and supersonic inlet conditions were investigated and optimized in the centrifugal compressors, respectively.
For diffuser of transonic inlet condition, NASA DDA 404-III high-pressure ratio centrifugal compressors [18] were studied, which are mainly used for advanced gas turbine power generation devices. The flow path of the compressor and its 3D model are illustrated in Figure 1. When a vaneless diffuser was equipped in the compressor, the pressure and temperature were measured at station A, which is 1.18 times of the impeller tip radius. When a wedge diffuser was used, outlet flow parameters were measured at station B. The impeller and diffuser geometries were provided by McKain and Holbrook [19], as shown in Table 1.
For diffuser of supersonic inlet condition, a pressure ratio of 8: 1 centrifugal compressor from our laboratory was used, which consists of a centrifugal impeller, a radial diffuser, and an axial diffuser. Table 2 shows the design parameters of this compressor.
For diffuser of subsonic inlet condition, a pressure ratio of 2.2: 1 centrifugal compressor from our laboratory was used, including a centrifugal impeller and a radial diffuser. This centrifugal compressor is used in a micro gas turbine, Table 3 shows the design parameters of the compressor.
3. Numerical Methods and Verification
3.1. Numerical Methods
CFD software NUMECA was used for the calculation. Spatial discretization is achieved with the cell-centered Jameson finite volume explicit scheme and temporal discretization is performed via a fourth-order Runge-Kutta scheme. Compared with the SST and κ-ε models, One-equation Spalart-Allamaras (SA) has better simulation of stage performances and is selected as the turbulence model. A multi-grid scheme is employed to reduce the computation time.
The topology of the impeller and diffuser is generated by AutoGrid, the grids were refined near the blade surfaces and the endwalls. To ensure boundary layer resolution, the y + is set around 2 at the wall surfaces. The compressor stage is of about 1.5× 106 grids. The computational grid of NASA DDA 404-III high-pressure ratio centrifugal compressor is shown in Figure 2.
A real gas model is used. The total pressure, the total temperature and the airflow angle (axial intake) are specified at the domain inlet. When the parameter change is less than 0.1% and the difference between the inlet and outlet mass flow is less than 0.1%, the calculation is considered to be convergent. Non-slip and adiabatic wall boundary conditions are used. The mixing plane is applied between the rotary and stationary domains.
3.2. Grid Independence Study
Grid independence was checked for all the simulations. As for the NASA DDA 404-III high-pressure ratio centrifugal compressor, grid number varies from 1.4 to 2.5 million. The calculation results at the design mass flow are shown in Table 4. Figure 3 shows the trend of calculation results for all grids. When the grid number surpasses 2.0 million, CFD simulations have consistent results.
3.3. CFD Result Verifications
Comparisons for the CFD prediction and experimental results of NASA DDA 404-III baseline compressor performance (80%, 90%, and 100% of design speed) are shown in Figure 4. Both the predicted choke mass flow rate and the total pressure ratio are just slightly higher than the experimental results, it can be concluded that the calculations can be used to predict the performance of the compressors.
For the pressure ratio 8 : 1 centrifugal compressor, Figure 5 shows the performance comparisons of the compressor at 100%, 95%, and 90% rotating speeds. It can be found that stage performances of the numerical calculation is close to the experimental results, and the biggest mass flow rate difference is about 2.7%.
For the pressure ratio of 2.2 : 1 centrifugal compressor, we have also verified. For the 100% rotating speed, compared to the experiment, the results of the numerical calculation are as follows: the total pressure ratio is 2.17% higher, the peak efficiency is 0.65% higher, and the choke mass flow rate is 2.20% higher.
To sum up the above results, it is noticed that the satisfactory agreement is obtained at low speeds. And the predicted choke flow rate and the pressure ratio are both slightly higher than that of the experiment results at design speed. The discrepancies not only occur in the current research, but also were reported by former researchers [20, 21]. And the numerical performance line distributions are well-matched with the experiment results. Therefore, we can conclude that the numerical simulations are good enough for performance prediction of the three centrifugal compressors.
4. Exploration of the Optimization Design Principles
For the vaned diffusers, the geometry variations are the common method to perform optimization design, while which geometry parameter to be optimized and what is the best value are difficult to tell, so far there are no universal principles for the optimization design of the vaned diffusers.
From the prospect of aerodynamic view, the universal optimization design principles for the vaned diffusers were studied in this paper. We determine the most important aerodynamic parameters of the vaned diffusers to be optimized, and try to find out their best distribution laws.
In this section, we first explore the optimal Vt, Vm distributions of transonic compressor diffuser by modifying the diffuser geometries, then the preliminary optimization design principles were concluded. Finally, we applied the design principles to the vaned diffuses of the transonic, supersonic and subsonic compressors, and the design principles were further investigated and verified for different conditions.
4.1. Main Aerodynamic Parameters for Optimization Design
For the performance evaluation of the diffuser, the static pressure recovery coefficient Cp and the total pressure loss coefficient Cpt are defined as [22]:
These two coefficients are defined uniformly as follow:
The subscripts 3 and 4 are the cross-sections of the diffuser inlet and outlet, respectively. For a diffuser, the higher the static pressure recovery coefficient is, and the lower the total pressure loss coefficient is, the δ is higher, and it means that the vaned diffuser could convert the kinetic energy of the impeller out flow into the pressure energy more efficiently.
As shown in Figure 6, radial velocity Vr, tangential velocity Vt and Meridional velocity Vm are defined as follows:
For radial diffusers, the velocity along the z-direction Vz is very small, so Vm is approximately equal to Vr.
For the stage performance evaluation, the surge margin is defined as:
The subscripts s and d represent the surge point and design point parameters, respectively.
4.2. Exploration of the Design Principles and Procedure
4.2.1. Geometry Parameterization
Blade parameterization was first performed for the optimization process. Bezier curves are used to control the blade profile. Given n + 1 = 6 control points Pi, forming n = 5th order Bezier curve polygon [23]. For given coordinates and order of the control points, the Bezier curve can be quickly fitted to the blade profile, and then the blade geometry can be controlled by adjusting the control points. A specific geometry sample is shown in Figure 7.
4.2.2. Studies for Vaned Diffuser with the Transonic Inlet Condition
The flow control equations of the vaneless diffuser were established with radial and tangential velocity [24]:
Where Cf is the friction coefficient, b is the blade height, and α is the airflow angle. It can be seen from formula (10) that the mass flow rate is only related to Vm and ρ when the blade height was unchanged; from formulas (8) and (9), the internal velocity of the diffuser could be expressed by two components, Vt and Vm.
Generally, for a vaneless diffuser, Vt decreases slowly with a typical spiral flow pattern. The long streamlines with high friction losses resulting in lower performance. One way to improve the performance is by a larger reduction of the Vt and the flow path by means of vanes turning the flow more radially. And Vm is defined as a constant value except for the variations of density or channel width [25]. While for a vaned diffuser, Vt decreases drastically and more kinetic energy is converted into pressure energy. Vm reflects the velocity of through flow and is related to the flow-through area, the Vm should not change drastically to avoid unnecessary flow losses. As for the radial diffuser, with the growth of radius and gas density, the Vm should decrease slowly to avoid larger flow separation [26] and flow losses.
In summary, the velocities of Vt and Vm are two important aerodynamic parameters, Vt is the largest velocity of all the velocity components, and represents the capability of pressure recovery; Vm is related to the flow-through area, and represents the capability of through flow. These two velocities are extremely important for the diffuser performance.
With the same meridional channel contour, the NASA DDA 404-III centrifugal compressor is equipped with a vaned diffuser and a vaneless diffuser, respectively. The velocity distributions along the meridional channel from the experimental datum are shown in Figure 8. For the vaneless diffuser, the Vt is reduced by about 1/6 from inlet to outlet, while for the vaned diffuser, Vt is reduced by 4/5, so the Vt decreases drastically in the vaned diffuser. Vm of the vaneless diffuser decreases gradually, but for this vaned diffuser, it first accelerated and then decelerated.
So the tangential velocity Vt and the meridional velocity Vm are also the critical indicators for the performance of the vaned diffuser. If their best distributions can be determined, then the universal design principles for the vaned diffusers can be found.
Geometry variations of the Vaned diffuser were performed in order to explore the best distributions of Vt and Vm for the transonic compressor diffuser. As shown in Figure 9, vaned diffusers with different blade angle and blade thickness distributions were studied. Case1-5 just change the blade angle distributions, and Case6-Case10 just change the blade thickness distributions.
For Case 1-5 (the thickness distribution unchanged) and Case 6-10 (the blade angle distribution unchanged), CFD calculations are performed with the same boundary conditions which taken from the peak efficiency point of the baseline case. The performance parameters for different cases are shown in Tables 5 and 6. Note that both Case 2 and Case 8 are the same baseline case (represented by the red lines in Figure 9).
As shown in Table 5, Case1, Case2 and Case3 are better in terms of δ and relative mass flow. As shown in Table 6, Case7, Case8 and Case10 are better in terms of δ and relative mass flow.
Figure 10 shows the Vt distributions along the meridional channel for different diffusers. Based on the diffuser performance with Tables 5 and 6, it shows that when Vt decreases faster at the beginning then slower along the meridional channel, the diffuser has better performance. But Vt should not drop too fast at the leading edge to avoid boundary layer separation. Therefore, better diffuser performance can be achieved by careful control of Vt distribution.
As for the Vm distribution, it should decrease gradually. There is a sudden increase in Vm at the leading edge of the baseline diffuser in Figure 8, it is undesirable. Case 11 was studied, which keeps the blade angle consistent with the baseline diffuser, but slightly change the thickness near the leading edge, as shown in Figure 11.
Table 7 shows the performances of Case 11 and the baseline diffuser. It is clear that Case 11 has better performance, and the relative mass flow is also improved. Figure 12 shows the distributions of Vt and Vm (Case11 and Baseline) along the meridional channel. Slight thickness variation near the diffuser leading edge will greatly change the distributions of Vt and Vm, and affect the performances of the diffuser.
Studies indicated that when Vt decreases fast at the leading edge and then decreases slowly in the diffusion channel, the kinetic energy will convert into pressure energy more efficiently; For the Vm, a steady decrease throughout the flow passage could avoid the flow losses caused by the sudden change in the meridional velocity.
Both the blade angle and thickness distributions determine the diffuser geometry, which could change the flow direction greatly, then affect the Vt distributions. Flow-through area will affect the Vm greatly, and the flow-through section of vaned diffuser is consist of the pressure surface, the suction surface and the hub/shroud endwalls, any change of these endwalls could affect the Vm distributions greatly. Generally, for the vaned diffusers, the blade angle distributions could greatly affect the Vt distributions, and the blade thickness distributions could change the diffuser profile and flow-through area, and affect the distributions of Vt and Vm at the same time, and we can also use the hub/shroud endwalls to control the Vm distributions.
Based on the above studies, the design principles for the vaned diffuser can be summarized as follows: 1) The Vt and Vm distributions of the vaned diffusers should be optimized for better diffuser performance; 2) Geometry modifications should be performed according to the desired Vt and Vm distributions; 3) The Vm should decrease gradually through the whole flow passage; 4) The Vt should decrease smoothly, and should decrease fast at the leading edge and then slowly along the flow passage for the transonic compressor diffuser. Note that the laws of Vt and Vm distributions may be different for vanned diffusers with different inlet Mach numbers.
Using the above design principles, the vaned diffuser of the transonic centrifugal compressor was optimized again, the velocity distributions and stage performance are shown in Figures 13 and 14. After optimizations of the diffuser, the compressor peak efficiency was increased by 2.15%, the pressure ratio was increased by 1.39%, and the surge margin was increased by 38.45%.
4.3. Optimizations for the Diffusers of the Subsonic and Supersonic Compressors
After the study of the diffuser for the transonic inlet condition, we further explore the distribution laws of Vt and Vm for the diffusers of supersonic and subsonic inlet conditions.
For the diffuser of supersonic inlet condition, optimizations were performed for the compressor from Table 2, and the result of a wedge diffuser (baseline diffuser is revised from this wedge diffuser) is also provided for comparisons, the profiles of the wedge diffuser and optimized case are illustrated in Figure 15. The velocity distributions and stage performances are shown in Figures 16 and 17. Compared with the result of the wedge diffuser, the peak efficiency of the optimized compressor was increased by 2.59%, the pressure ratio was increased by 6.49%, and the surge margin was increased by 15.46%. From Figures 16 and 17, it can be seen that a proper decrease rate in Vt is beneficial to the performance improvement, but we found that excessive Vt decrease at the inlet of the diffuser resulting in reduced performance.
For the diffuser of subsonic inlet condition, optimizations were performed for the compressor from Table 3, and the velocity distributions and stage performances are shown in Figures 18 and 19. After optimizing the velocity distributions of the diffuser, compressor peak efficiency was increased by 4.23%, the pressure ratio was increased by 3.36%, and the surge margin was increased by 28.94%.
For all cases, the compressor efficiencies and pressure ratios were improved as a result of the better performances of the diffusers. Specially, the optimal Vm distributions certainly have enhanced the flow-through capabilities of the diffusers and the surge margins as well.
4.4. Discussion of the Distribution Laws for Different Compressor Conditions
Three optimization results for the diffusers of different compressor conditions have shown that the design principles by controlling the distributions of Vt and Vm along the flow passage are feasible and effective.
Generally, Vm should decrease gradually and slowly for all inlet conditions, but the specific Vt distribution laws for diffusers of different inlet conditions are different. A general Vt distribution map based on a large number of optimizations is provided in Figure 20, both the Vt and the meridional channel are normalized.
There are 3 dash lines in Figure 20, the upper and lower dash lines represent the “subsonic boundary” and “transonic & supersonic boundary”, which are the extreme Vt distributions. And the middle dash line separate the map into a “lower section” and “upper section”, the positions and shapes of all the three dash lines are the results of approximation. Three solid lines with square, circle and triangle symbol are the actual diffuser Vt distributions from the above optimizations for the subsonic, transonic and supersonic diffusers, respectively.
Apparently, the Vt should decrease smoothly for all conditions. For the transonic and supersonic inlet conditions, their Vt distributions are similar, Vt decreases fast at the inlet then decreases slowly in the flow passage, their Vt distributions are of curved lines in the lower section. During the optimization process, we found that if the Vt decreased too fast at the inlet will result in flow separation. For the subsonic inlet conditions, the Vt distributions are flatter curved lines in the upper section. Note that the Mach number of the subsonic case in this paper is pretty low, and its Vt distribution is close to a linear distribution. Though there aren’t enough cases to fill the upper section of this map, it’s reasonable to speculate that with the growth of the Mach number, the Vt distribution curve should be steeper and more curved in this map for the subsonic compressor conditions.
The specific law of the Vt distribution for a given diffuser could be found in the related sections of Figure 20 according to the inlet conditions of the vaned diffusers, and the very best Vt distribution could be explored in the adjacent areas of this map with further optimizations.
4.5. The Optimization Design Procedure
Since the Vt and Vm will affect each other during the optimization process, an optimization design procedure is presented in Figure 21 to accelerate this process.
The main steps of the procedure are as follows: (1)Parameterization: Baseline diffuser is parameterized for the optimization process. Note that R4/R3, R3/R2 and the inlet/outlet blade angles should not change too much for the given diffuser working condition(2)Simulation: Vt, Vm, m and Ma are provided by the numerical calculations(3)Optimization: Adjust the diffuser geometric parameters to meet the desired Vt and Vm distributions along the meridional channel with all the limitations in check. Iterations are performed to meet the desired distribution laws of Vt and Vm
Note that the Vt is very sensitive to the blade angle. Blade thickness could affect the Vt and Vm at the same time. Furthermore, hub/shroud endwall geometry could also be changed to obtain the desired Vm distributions. Usually, the desired velocity distributions and satisfying diffuser performance could be achieved after couple iterations.
5. Summary
In this section, the design principles of vaned diffuser are preliminarily determined for the transonic compressor condition, then this design principles were applied for the transonic, supersonic and subsonic compressor conditions to explore the specific distribution laws of Vt and Vm. The universal design principles of vaned diffuser are investigated and verified for all the typical compressor inlet conditions. A general Vt distribution map is provided and an optimization design procedure is also presented.
6. Analysis and Discussion
After optimizations of the diffusers, the overall performances of the compressors are improved substantially, and the details of the diffuser flow field have changed significantly. Figure 15 shows the Mach number contours of the diffusers at 50% span for all cases. With the decreases of the Vt and Vm at the inlet, a lower inlet Mach number was found for all cases, which could greatly improve the overall performance of the diffusers. We found that the overall gradients of the Mach number was lowered for all cases after optimizations, and more even distribution of the Mach number along the flow passage were observed at the same time.
The better Mach number contours certainly improved the entropy contours as shown in Figure 22. Both the highest entropy and the area of high entropy were reduced substantially, which indicate that less flow losses and better overall performances of the optimized vaned diffusers. Specially, we found great entropy improvements near the pressure side surface of the diffusers for all cases, which could be the results of lower Mach numbers and less flow separations.
As shown in Figure 23, since the Mach number decreases faster at the inlet, and more kinetic energy was converted into pressure energy, higher pressure gradients at the inlet are observed for all cases.
The total pressure recovery coefficients were improved substantially for all cases as shown in Figure 24. The total pressure recovery coefficients have not been increased much at the inlet because the fast velocity decreases, but have been increased a lot through the whole passage for better distributions of Vt and Vm.
As shown in Figure 25, the overall distributions of Vm were improved greatly for all cases. The highest Vm and the Vm gradients are lower, and more even Vm distributions were observed, especially for the subsonic case. Better Vm distributions certainly could improve the efficiencies and flow-through capabilities of the diffusers, and we do observed great improvements of the surge margins for all cases.
After optimizations with the new design principles of the vaned diffusers, not only the overall performances of the diffuser were improved greatly, but also the details of the flow fields, the distributions of Mach number, entropy, static pressure and Vm are better, and the flow separations are restricted as well.
7. Conclusions
In this paper, numerical simulations and optimizations were performed on the vaned diffusers with different inlet Mach numbers. The new vaned diffuser design principles were investigated and verified, and the main conclusions are as follows: (a)Vaned diffuser design principles: (1)The performance of the vaned diffuser can be improved by carefully control of two aerodynamic parameter distributions: Vt and Vm, which represent the capability of pressure recovery and through flow, respectively(2)The Vm should decrease gradually and slowly along the flow passage for all inlet conditions to avoid the flow losses caused by the sudden change in the meridional velocity, and improve the flow-through capability and surge margin(3)The diffuser Vt distribution should vary from the inlet conditions in order to ensure kinetic energy converted into pressure energy more efficiently. Generally the Vt should decrease smoothly, and the Vt distribution is of a tilted concave curve, which means that the Vt decreases faster at the inlet then slower along the flow passage. A general Vt distribution map is provided, but the specific Vt distribution for a given vaned diffuser should be further optimized accordingly(b)Diffuser geometry modifications: (4)The geometry of the vaned diffuser should be adjusted to guarantee the desired Vt and Vm distributions(5)Blade angle could change Vt directly, blade thickness could change blade profile and flow-through area, then affects the Vt and Vm at the same time, and hub/shroud endwall geometry could be changed also to achieve the desired Vm distributions(c)Optimization results for vaned diffusers of different conditions: (6)For the subsonic, transonic and the supersonic vaned diffusers, the stage peak efficiencies were increased by 4.23%, 2.15% and 2.59%, respectively, the stage pressure ratios increased by 3.36%, 1.39% and 6.49%, respectively, and the surge margins were improved by 28.94%, 38.45% and 15.46%, respectively(d)Optimized design procedure: (7)Vt and Vm will affect each other during the optimization process, an optimization design procedure is presented to accelerate the process in this paper
Nomenclature
| Cp: | Static pressure recovery coefficient |
| Cp,t: | Total pressure loss coefficient |
| CFD: | Computational fluid dynamics |
| Ma: | Mach number |
| Ma0: | Mach number at diffuser inlet |
| p: | Static pressure |
| pt: | Total pressure |
| PS: | Pressure surface |
| SS: | Suction surface |
| R: | Radius |
| V: | Velocity |
| Vt: | Tangential velocity |
| Vm: | Meridional velocity |
| Vr: | Radial velocity |
| η: | Total-to-total isentropic efficiency |
| π: | Total pressure ratio |
| m: | Mass flow rate |
| δ: | Diffuser performance parameter |
| x, y: | x and y components of rectangular coordinate system |
Data Availability
All data generated or analyzed during this study are included in this article.
Conflicts of Interest
No conflict of interest exists in the submission of this manuscript, and manuscript is approved by all authors for publication. I would like to declare on behalf of my co-authors that the work described was original research that has not been published previously, and not under consideration for publication elsewhere, in whole or in part.
Acknowledgments
The authors would like to acknowledge the financial support of the National Major Science and Technology Project of China (No. 2017-II-0002-0014) and the National Natural Science Foundation of China (Project No. 51706223) for this project.