Research Article | Open Access
Zhiyuan Zhang, Rui Ma, Chunge Wang, Jun Zhao, "Research on Springback Control in Stretch Bending Based on Iterative Compensation Method", Mathematical Problems in Engineering, vol. 2019, Article ID 2025717, 11 pages, 2019. https://doi.org/10.1155/2019/2025717
Research on Springback Control in Stretch Bending Based on Iterative Compensation Method
An unavoidable problem in the stretch-bending process is springback, which dictates the shape and dimensional accuracy of the product. This problem can be solved by adjusting the geometry of the die or through active process control. This study focuses on the design of the die shape to achieve the target product. Based on the fixed-point iterative method and displacement adjustment (DA) method, this paper proposes an iterative compensation method, which has a higher convergence rate, lower number of iterations, and higher precision compared to the DA method with only one control parameter. In addition, like the DA method, the proposed method does not depend on the material properties or mechanical model, but the difference is that it can quickly and effectively find out the iteration parameter, determine whether the parameter has convergence or not, and has no compensation factor. According to the deviation of iterative parameters between the value after stretch bending and the target value, the iterative compensation method can be used to calculate the compensation magnitude and compensation direction of the iterative parameter. For stretch-bending processes with invariable- and variable-curvature die shapes, the convergence of control parameters is verified mathematically with the convergence theorem of the method, and experiments are conducted to verify the iterative compensation method. The experimental results show that the target products can be obtained with a small number of iterations without knowing the specific material properties.
Stretch-bending processes are widely applied in shipbuilding, automobile, aerospace, and other manufacturing industries because they yield components with high precision and good surface quality. Springback is a key issue in the stretch-bending process, and it dictates the dimensional accuracy of the product [1–3]. In conventional stretch-bending production, several methods have been applied to reduce springback. Most of them focus on adjusting the main process parameters such as axial tension and on optimizing the loading mode [4–7]. These methods are effective and do not require the adjustment of the die shape, but they fail to eliminate springback completely and require time-consuming trial-and-error procedures.
In some methods, instead of reducing springback, the die shape is adjusted to compensate for springback. Consequently, although a large springback would remain, the final product would be close to the requirements because of the modified die shape. These methods are relatively more cost-effective and have the potential to compensate for springback completely, even for complex parts. Lingbeek et al.  adopted the displacement adjustment (DA) method, which is a strictly geometrical method that displaces the die geometry in the direction opposite to the geometrical error. In the DA method, a large amount of data is needed to determine the compensation factor. The method is general in that it is not limited to operations having a particular symmetry, die shape, or magnitude of springback shape change . However, tests of industrial cases show that the effectiveness of the method depends on the material, process, and geometrical parameters [10, 11], and the method involves inefficient trail-and-error procedures. Thus, increasing the convergence rate and consequently reducing the number of iterations are important for the iterative DA method.
To achieve the above objectives, various methods have been proposed based on finite-element simulation. The comprehensive compensation (CC) method  mainly focuses on the large rotation and displacement of advanced high-strength steel and automotive panel stamping. Although this method can have high precision, the compensation direction of its control parameter, which depends on the material properties, should be determined in advance. The alternate hybrid method (HM)  is a combination of the DA and SF methods, whereas the SF method is based on stress state in the forming process, which is multiplied by negative factor to calculate the deformed shape to realize the compensated surface . The accelerated compensation (AC) method  reduces the number of compensation steps by establishing mechanical models. The enhanced displacement adjustment (E-DA) method  applies additional point topology mappings, which establish corresponding interrelations between the discretized point topologies by simulating or testing the forming process many times. The sheet elements compensation (SEC) method  is based on trial and error and a theoretical analysis, because of which the specific material properties should be determined in advance. The smooth displacement adjustment (SDA) method and the surface controlled overbending (SCO) method  both aim at deep drawing and optimize the tool shape through finite-element simulation. Xiong et al.  proposed a springback compensation algorithm for age forming based on a quasi-Newton method, which requires the specific material properties to establish a mechanical model for calculating the iteration factors. Lingbeek et al.  proposed a noniterative variant, also known as the one-step DA method, which depends on the material, process, and geometrical parameters.
Although the above methods can quickly meet the requirements, they need either the specific material properties to establish mechanical models or a large number of finite-element simulations. Ma et al.  established an iterative compensation mechanism based on an iterative solution method, but they only analyzed pure bending under a simple stress state. A clear and stable compensation direction has a great effect on the precision of the results. When the compensation direction coincides with the direction of the desired die, the geometrical error and the precision of results can be improved. On the other hand, a significant deviation from the direction of the desired die will decrease the precision of results and increase the number of iterations. Based on the above analysis, the present study developed an iterative compensation method for stretch bending involving simple geometries. The proposed method has only one control parameter and a clear and stable compensation direction, and it requires only a small number of iterations and experimental data of the control parameter before and after springback.
2. Principle and Method of Iterative Compensation
2.1. Iterative Compensation Principle
Compensating for springback by adjusting the die shape is one of the methods adopted to obtain the required products. Specifically, according to the predicted or measured springback value, the die shape is modified so that the shape of the product after springback approximates the design requirements. Conventionally, this approach is a trial-and-error method that mainly depends on experience and repeated experimentation, which is time-consuming and costly because of its nondirectional and nonquantitative nature. To improve the efficiency of springback compensation, this paper proposes the fixed-point iterative method , which obtains the approximate solution of an equation. Its basic principle is as follows.
Firstly, it is assumed that the product shape after springback is represented by a function with a single control parameter , and the target value of the control parameter is equal to a constant . The equationis transformed into an implicit equation to obtain its root:
To use the fixed-point iterative method, one must examine whether the equation has iterative convergence in the neighborhood of the equation root , which is a necessary condition for this method. The local convergence theorem of the fixed-point iterative method states that if is continuous in some neighborhood of and , the equation has local convergence and can be solved iteratively. The process of determining the is as follows:
When , where is a preset error, it can be considered that . From the above mathematical analysis, the iterative method can make the iterative value gradually converge to the exact solution.
The derivative of the implicit equation is
If , the compensation magnitude becomes excessive, and equation (1) results inincreasing the number of iterations. It can be concluded that if the equation needs to be solved iteratively, is an indispensable condition.
According to the equivalent stress-strain relationship in the elastic-plastic deformation shown in Figure 1 , the springback produced by the larger equivalent strain is distinctly greater than the springback produced by the equivalent strain . Therefore, for the same ordinary metal materials, the larger the deformation, the greater is the springback when maintaining deformation conditions. Let be the value of the control parameter before springback and be the function expressing the control parameter after springback. Then, the springback value can be expressed as . According to the above theoretical basis, both and are monotonically increasing functions; that is, and . Thus, .
Based on the above analysis, an iterative compensation mechanism is established to compensate for the springback in stretch bending.
2.2. Iterative Compensation Method
The iterative compensation method to compensate for the springback in stretch bending is described as follows. For a general stretch-bending process, the relationship between the shape of the profile and the controllable parameter (before and after springback) satisfies , where . Then, when the target value of the controllable parameter after springback is , set . The desired value for the controllable parameter before springback can be obtained according to the flowchart shown in Figure 2.
However, for a general stretch-bending springback problem, is a nonlinear function that is too difficult to solve mathematically. Therefore, the fixed-point method is only used to help in selecting the control variables. In the practical application of the iterative compensation mechanism, the value of is usually measured using experiments.
The iterative compensation method can reduce the error in the iterative parameter and guide the die shape in the correct direction. In addition, each compensation value is defined solely by the difference between the experimental data of the iterative parameter before and after springback.
3. Mechanical Model of Stretch-Bending Process
The object under study is a metal sheet with a rectangular section, and Figure 3 shows the loading method of stretch bending. The tension is applied at the two ends of the sheet, following which the sheet is fitted with the die under the moment . The magnitude of the tension is constant, and its direction is always along the neutral layer of the sheet in stretch bending. The bending radius of the neutral layer of the sheet is denoted by .
3.1. Basic Assumptions
The following basic assumptions are made in the mechanical model. First, the cross section is assumed to be always planar and perpendicular to the neutral layer. Second, the blank is assumed to be in an approximately unidirectional stress-strain state. Third, it is assumed that the neutral stress layer, neutral strain layer, and geometric center layer of the sheet always coincide. The strain is linearly distributed and satisfieswhere is the strain, is the distance from the particle to the neutral layer, and is the bending radius of the neutral layer. Fourth, the material properties are in accordance with the bilinear hardening model:where is the stress, is the strain, is the elastic modulus, is the plastic modulus, and is the initial yield stress.
The cross section of the sheet is rectangular, with a cross-sectional area and moment of inertia of and , respectively, where b and t are the width and thickness of the sheet, respectively.
3.2. Relationship between Curvature of the Neutral Layer before and after Springback
When the sheet is stretched, the stress of the cross section isand the strain of the cross section iswhere is the tensile load and A is the cross-sectional area.
Since only the changes to the curvature with springback need to be analyzed, a part of the specimen is selected as the study subject. The length, width, and thickness directions of the specimen are represented as the x-, y-, and z-axis, respectively, as shown in Figure 4. After the bending load is completed, the stress in the cross section satisfieswhere is the bending radius at the elastic limit.
Because the magnitude of the tension is constant,
The springback equation of stretch bending can be obtained according to the springback equation of small-curvature plane bending  as follows:where and are the radii of the neutral layer before and after springback, respectively.
Substituting and into equation (13) and letting represent the relationship between the curvatures before and after springback, thenwhere and are the curvatures of the neutral layer before and after springback, respectively.
3.3. Iterative Convergence of Curvature Parameter
In order to use the fixed-point iterative method, is taken as the iterative parameter, and its iterative convergence needs to be proved. The derivative of the function defined in equation (14) is
For general materials, it is considered that . In this case, it can be deduced that
According to the convergence theorem of the fixed-point iterative method, equation (16) demonstrates the iterative convergence of the curvature of the neutral layer in the stretch-bending process, but the curvature is not a controllable parameter unless it is invariable everywhere along the neutral layer. Therefore, it can be concluded that, to apply the iterative compensation method, the curvature of the neutral layer can be used as the iterative parameter in the stretch-bending process with an invariable-curvature die shape.
4. Springback of Stretch Bending with Variable-Curvature Die Shape
To apply the iterative compensation method in the stretch-bending process with a variable-curvature die shape, the die shape having the curve with and is considered in this paper.
It is assumed that the curve shape of the inner cross section in the effective deformation zone of the sheet conforms to the equation after stretch bending. Therefore, the radius of curvature of each point after springback satisfies
Since the analytical solution of the equation is difficult to determine, in order to discuss the problem in practical application, the numerical solution is obtained using MATLAB software with and the material properties given. In the process of using MATLAB software, the ode45 function, namely, the fourth-fifth-order Runge–Kutta algorithm, is selected to solve the problem. It provides candidate solutions with the fourth-order method and controls errors with the fifth-order method.
Table 1 lists the material properties of a cold-rolled ST12 steel sheet obtained from a uniaxial tension test. The sample size is considered to be in the table.
Because the sheet should be in the elastoplastic state during the stretch-bending process and the length of sheet is 490 mm, the stress of the cross section is set as , and the value of the abscissa ranges from −200 to 200. For an initial profile with the given as 0.0035, the profile after springback was obtained using equation (18).
As shown in Figure 5, the equation of the fitting curve is , and the coefficient of correlation , which represents the fitting accuracy, is 0.99998. It is considered that the characteristics of the quadratic curve are maintained after springback, and can be described by .
To discuss the local convergence of , various values of are selected, which represent different curvatures, and the corresponding values of are obtained according to the above steps, as listed in Table 2. Subsequently, these data are fitted, and the result is shown in Figure 6. The equation of the fitted curve is , and the coefficient of correlation is 0.9996. It can be concluded that when . Therefore, the in can be considered as the iterative parameter of the iterative compensation method.
Before springback, suppose (, ) is a point in the curve . The curvature at this point is
After springback, the coordinates of point are transformed to (, ), which represent a point in the curve . The curvature at this point is
For the same point, , , , , and . Therefore, it can be deduced that . Hence, it is proved that has iterative convergence in the stretch bending for a die shape defined by a quadratic curve and can be used as the iterative parameter for the iterative compensation method.
In order to discuss whether has iterative convergence when stretch bending is performed with a die shape defined by a cubic curve or higher-degree polynomial curve, the stretch bending for is studied using the same procedure as that applied for a quadratic curve.
Before springback, suppose is a point in the curve with the curvature
After springback, the coordinates of are transformed to , which represent a point in the curve . The curvature at this point is
For point , , , , , and . Therefore, it can be deduced that . Hence, it is proved that has iterative convergence in the stretch bending for a die shape defined by a cubic curve. It can be inferred that, for higher values of in , the above derivation can be repeated to prove the iterative convergence of .
5. Simulation and Discussion
In order to verify the convergence of , the iterative process is simulated for die shapes defined by and . In addition, the iterative process of the DA method is simulated.
5.1. Finite-Element Model
The finite-element model is established using ABAQUS 6.10 software. The mechanical properties and geometric dimensions of the selected sheet are listed in Table 1. The sheet is discretized by 8-node linear hexagonal incompatible mode elements (C3D8I). The die is modeled as a discrete rigid body. The contact between the sheet and punch is set as a purely master-slave and kinematic contact condition, and the frictional coefficient is set as 0.1. In addition, the ABAQUS/Standard solver is employed (Figure 7).
5.2. Result and Discussion
In the process of finite-element simulation, the tension applied at the two ends of the sheet is T = 8000 N, and the compensation accuracy is set at 1% of the control parameters. Under the same conditions and requirements, the compensation of the DA method is also simulated. The compensation factor is chosen as 0.8 and 1.5 due to that the compensation factor of the DA method ranging from 0.7 to 2.5 according the practical experience, and the value of compensation factor is different for each forming process and cannot be predicted effectively .
The simulation results obtained using the proposed iteration method and DA method with the die shape are shown in Figure 8, where is the target value of the control parameter and is the value of the control parameter after springback. It can be seen that, with the proposed iterative method, only one iteration is needed to obtain the target value, and the accuracy is high. With the DA method and h = 0.8, although only one iteration is needed and the accuracy meets the requirement, its accuracy is not as high as that of the proposed iterative compensation method, especially at , which requires two iterations. When the compensation factor h = 1.5, the value of oscillates around zero and the compensation direction keeps changing, which not only increases the number of iterations but also fails in improving the compensation accuracy owing to the inappropriate compensation factor.
Figure 9 shows the results of simulating the iterative compensation process with the die shape . As shown in the figure, with high accuracy, the proposed iterative method needs two iterations at most, and the DA method with h = 0.8 needs the same number of iterations to meet the compensation requirements. Moreover, the DA method requires more iterations to achieve the same accuracy as the iterative compensation method. When the compensation factor h = 1.5 and , the value of oscillates around zero and the compensation direction keeps changing. However, under the condition of h = 1.5 and , the accuracy requirement is satisfied after one iteration, but the accuracy after two iterations has not reached that of the proposed iterative compensation method.
It can be inferred from the simulation that has iterative convergence. Compared with the DA method, the proposed iterative compensation method has a higher convergence speed and higher accuracy.
6. Experiment on Iterative Compensation Method of Stretch Bending
To verify the iterative compensation mechanism in stretch bending, an experimental device was built, which includes a stretch-bending testing machine, a tension-detecting device, and measurement equipment. The stretch-bending testing machine includes a hydraulic power back, an actuator, and an electric control unit. Figure 10 shows the structure of the concrete. The tension-detecting device detects the magnitude of tension during the stretch-bending process. The equipment used to measure the sheet curvature after springback is a 3000iTM series of portable three-coordinate measurement machine (ROMER CimCore Inc., USA) with a measurement accuracy of 0.01 mm.
The stretch-bending process is performed using the following steps. Firstly, the stretch cylinders are rotated such that they are collinear with the jigs, which are fixed on the locators. Next, the stretch cylinders apply a certain tension on the specimen, the ends of which are fixed on the jigs, and the tension is equal to the value required in the experiment process. Then, the bending cylinder moves the bending die forward to finish the bending process. Finally, the loading is stopped, and the springback occurs on the specimen.
Subsequently, a series of data points of the formed specimen contour is obtained using the measurement equipment, and the approximate distance between two adjacent points is 1 mm. These data points are fitted into the required curve form.
The test sheets are a cold-rolled ST12 steel sheet and a cold-rolled H59 brass sheet and their geometry and material properties are listed in Tables 1 and 3, respectively. Since the proposed iterative compensation method does not need the material properties in advance, Table 3 shows only the tensile strength of the H59 sheet, which was provided by the supplier.
In Table 4, Kd is the curvature of the target sheets; T is the value of the tension; K and Kp are the curvature of the sheets before and after springback, respectively; and Knext is the curvature of the next bending process before springback, which is also the curvature of next bending die shape.
Taking the bending curvature as an example, the compensation accuracy is determined to be 0.1%. That is, the error value of curvature is less than . The iterative compensation process is as follows.
First, the ST12 sheet is stretched and bent by the die with a curvature of the target value . The curvature of the sheet is after springback, and the iterative error is calculated to be , which does not meet the accuracy requirements. Therefore, a second compensation step is required, and the curvature of the second compensation step is calculated as .
The second stretch-bending process is performed with a curvature of , and the curvature of the sheet is after springback. The error is , which does not meet the accuracy requirement either. Therefore, a third compensation step is required with a curvature of .
The third stretch-bending process is performed with a curvature of . The curvature is after springback, and the error is , which meets the accuracy requirement. Thus, the iterative compensation ends.
The curvature of the die is determined to be after three iterative compensation steps, and the precision is controlled within 0.1%. According to Table 4, the error decreases rapidly and the iterative parameter approaches the target value rapidly with the increase of compensation steps. Simultaneously, it is verified that the curvature can be used as an iterative parameter to apply the iterative compensation method in stretch bending with an invariable-curvature die shape, and the die shape can be determined using limited iterations to obtain products with the required precision.
Similarly, experiments for stretch bending with a variable-curvature die shape were performed under a stable tension, and the shape of the target products conforms to . The experimental results are presented in Figure 12 and Table 5.
In Table 5, is the target value of the controllable parameter; is the value of the tension; and are the values of the controllable parameter before and after springback, respectively; and is the value of the controllable parameter in the next bending process.
According to Table 5, the target products with the required accuracy are obtained with two or three iterative compensation steps, and the compensation error decreases rapidly with increase of compensation steps. It is suggested that the coefficient can be used as the iterative parameter to apply the iterative compensation method in the stretch-bending process with variable-curvature die shapes of and . In addition, for H59 sheets, the target products are obtained without measuring the material properties or establishing mechanical models.
The experimental results for H59 and ST12 sheets show that the proposed iterative compensation method can effectively compensate for springback in the stretch-bending process with a finite number of compensation steps and does not depend on the material properties or mechanical model. Moreover, it has a stable convergence direction, and the error accuracy can be evaluated. Therefore, the iterative compensation method has high practical value.
In this paper, an iterative compensation method was established to compensate for springback in the stretch-bending process. Based on a theoretical analysis of the springback mechanism, the controllable parameter in the stretch-bending process with invariable- and variable-curvature die shapes was found, and the iterative convergence of the parameter was proved. Physical experiments were designed to verify the iterative compensation method. From the results, the following conclusions are drawn.
The iterative compensation method not only determines the die shape through finite compensation to obtain the target products but also is independent of the material properties and mechanical model. Further, it was theoretically proven that the curvature K in the stretch-bending process with an invariable-curvature die shape has iterative convergence. Therefore, K can be used as the iterative parameter to apply the iterative compensation method. Moreover, through numerical and theoretical analyses, it was proven that the coefficient of the stretch-bending process with the die shape has iterative convergence. Therefore, can be used as the iterative parameter to apply the iterative compensation method.
Iterative compensation experiments on H59 and ST12 sheets were conducted based on curvature control. It was shown that the iterative compensation method has a stable convergence direction, and the compensation error decreases rapidly with the increase of the number of compensation steps. Finally, the error was found to be within acceptable levels.
The experimental 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.
This work was funded and supported by the Natural Science Foundation of Hebei Province in China (Grant no. E2015203244).
- L. Xiaoqiang, C. Zengqiang, W. Junbiao et al., “Progress in the research on profile stretch bending process,” Aeronautical Manufacturing Technology, vol. 4, pp. 80–82, 2004.
- C. F. Jin and J. Wu, “Study on stretch-bending process for aluminium profile,” Mechanical Engineering & Automation, vol. 3, pp. 211-212, 2010.
- X. Yi, L. Luo-xing, L. Guang-yao et al., “Extrusion bending: method, application and recent development,” Journal of Plasticity Engineering, vol. 15, no. 3, pp. 61–67, 2008.
- J. Zhao, R. Zhai, Z. Qian, and R. Ma, “A study on springback of profile plane stretch-bending in the loading method of pretension and moment,” International Journal of Mechanical Sciences, vol. 75, no. 10, pp. 45–54, 2013.
- T. Kuwabara, S. Takahashi, K. Akiyama et al., “2-D springback analysis for stretch-bending processes based on total strain theory,” Transactions of the Society of Automotive Engineers, vol. 5, pp. 504–513, 1995.
- A. El-Domiaty, “Stretch forming of beams of non-uniform section,” Journal of Materials Processing Technology, vol. 22, no. 1, pp. 21–28, 1990.
- T. Liu, Y. Wang, J. Wu et al., “Springback analysis of Z & T-section 2196-T8511 and 2099-T83 Al-Li alloys extrusions in displacement controlled cold stretch bending,” Journal of Materials Processing Technology, vol. 225, pp. 295–309, 2015.
- R. A. Lingbeek, W. Gan, R. H. Wagoner, T. Meinders, and J. Weiher, “Theoretical verification of the displacement adjustment and springforward algorithms for springback compensation,” International Journal of Material Forming, vol. 1, no. 3, pp. 159–168, 2008.
- W. Gan and R. H. Wagoner, “Die design method for sheet springback,” International Journal of Mechanical Sciences, vol. 46, no. 7, pp. 1097–1113, 2004.
- R. A. Lingbeek, J. Huétink, S. Ohnimus et al., “Iterative springback compensation of NUMISHEET benchmark #1,” in Proceedings of the AIP Conference Proceedings, vol. 778, pp. 328–333, Berhampur, India, June 2005.
- R. Lingbeek, T. Meinders, S. Ohnimus et al., “Springback compensation: fundamental topics and practical application,” in Proceedings of the 9th ESAFORM Conference, N. Juster and A. Rosochowski, Eds., pp. 403–406, Glasgow, UK, April 2006.
- X. A. Yang and F. Ruan, “A die design method for springback compensation based on displacement adjustment,” International Journal of Mechanical Sciences, vol. 53, no. 5, pp. 399–406, 2011.
- W. A. Siswanto, A. D. Anggono, B. Omar et al., “An alternate method to springback compensation for sheet metal forming,” The Scientific World Journal, vol. 2014, Article ID 301271, 13 pages, 2014.
- A. P. Karafillis and M. C. Boyce, “Tooling design in sheet metal forming using springback calculations,” International Journal of Mechanical Sciences, vol. 34, no. 2, pp. 113–131, 1992.
- H. S. Cheng, J. Cao, and Z. C. Xia, “An accelerated springback compensation method,” International Journal of Mechanical Sciences, vol. 49, no. 3, pp. 267–279, 2007.
- G. Cafuta, N. Mole, and B. Štok, “An enhanced displacement adjustment method: springback and thinning compensation,” Materials & Design, vol. 40, no. 3, pp. 476–487, 2012.
- Z. Zhang, J. Wu, S. Zhang, M. Wang, R. Guo, and S. Guo, “A new iterative method for springback control based on theory analysis and displacement adjustment,” International Journal of Mechanical Sciences, vol. 105, pp. 330–339, 2016.
- R. Lingbeek, J. Huétink, S. Ohnimus et al., “The development of a finite elements based springback compensation tool for sheet metal products,” Journal of Materials Processing Technology, vol. 169, no. 1, pp. 115–125, 2003.
- W. Xiong, Z. Gan, S. Xiong, and Y. Xia, “Rapid springback compensation for age forming based on quasi Newton method,” Chinese Journal of Mechanical Engineering, vol. 27, no. 3, pp. 551–557, 2014.
- R. Ma, C. Wang, R. Zhai et al., “An iterative compensation algorithm for springback control in plane deformation and its application,” Chinese Journal of Mechanical Engineering, vol. 494, no. 14, pp. 62–69, 2018.
- Q. Li, N. Wang, and D. Yi, Numerical Analysis, Huazhong University of Science & Technology Press, Wuhan, China, 4th edition, 2006.
- J. Zhao, J. Yin, R. Ma, and L. Ma, “Springback equation of small curvature plane bending,” Science China Technological Sciences, vol. 54, no. 9, pp. 2386–2396, 2011.
Copyright © 2019 Zhiyuan Zhang 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.