Abstract
Hypersonic vehicle has gained increasing attention due to its high cruise speed and long voyage. In this paper, an enhanced Sequential Approximate Optimization method is proposed for aerodynamic optimization of a hypersonic vehicle. In this method, a constrained design of experiment method is adopted to generate the initial sampling set with adequate number of feasible points. A field metamodel is proposed to surrogate the aerodynamic parameters distributions in the flow field obtained by the aerodynamic simulation. A hybrid metamodel combing radial basis functions and polynomial chaos expansion is applied to construct the field metamodel with good approximate performance. A robust mesh morphing scheme based on radial basis functions is developed to generate high-quality meshes for the sequential optimization scheme. The hypersonic vehicle aerodynamic optimization problem is performed using the proposed optimization framework and satisfactory results are obtained with limited computational budgets. Results show that the proposed field metamodel-enhanced Sequential Approximate Optimization method possesses powerful optimization performance and promising prospects in the field of hypersonic vehicle optimization design.
1. Introduction
Over recent decades, hypersonic vehicles have generally gained increasing attention because of potentially significant increases in range capabilities provided by their high lift-to-drag ratio and superior break-defense capabilities provided by their high cruise speeds [1]. The overall performance of hypersonic vehicles, especially the aerodynamic properties, largely depends on their geometric configuration.
Various design methods have been studied to obtain satisfactory hypersonic configurations. In the early stage, just like other aircraft design problems, cut-and-try methods were widely studied in the design problems of hypersonic vehicles. In 1965, Küchemann [2] proposed the aerodynamic design problem of hypersonic vehicle firstly. A conical-derived design-based approach is introduced by Rasmussen [3] in 1990. A few years later, in 1999, Sobieczky et al. [4] developed a characteristic inverse-based method using osculating cones. The traditional cut-and-try approaches are based on the experience of the designers, which, however, are not always available for newly developed and unconventional aircrafts such as hypersonic vehicles. Moreover, the cut-and-try methods are only capable of identifying some feasible candidates rather than locating the optimal solution.
For these circumstances, recent efforts in the field of hypersonic vehicles focus on methods that are aimed at getting the optimum. Many researchers have chosen to combine evolutionary or gradient-based optimization algorithms with semiempirical estimation models to optimize hypersonic vehicles. Kinney utilized the Newtonian Tangent Cone and Tangent Wedge method and a conjugate gradient approach to search for an optimal hypersonic design [5]. Theisinger and Braun employed Newtonian flow theory to rapidly assess the hypersonic aerodynamic performance of candidate aeroshells and utilized the refined Nondominated Sorting Genetic Algorithm (NSGA-II [6]) to conduct a multiobjective optimization [7]. Su et al. utilized local-inclination methods for aerodynamic estimation analysis and SNOPT [8] as an optimization algorithm solver [9]. Zhang et al. further extended Su et al.’s research by combining local-inclination methods with a panel method [10]. In Shen et al.’s documented research, a panel method was applied to numerically investigate the aerodynamic characteristics and combined with the Multi-Island Genetic Algorithm (MIGA) in the aerodynamic optimization of EXPERT configuration [11].
In modern aerodynamic design optimization applications, the high-fidelity computational fluid dynamic (CFD) simulation models have been widely used to improve the confidence and accuracy of analysis [12]. The modelling quality and solution accuracy are improved via the high-fidelity CFD models, which usually require meshes composed of a large number of cells and numerous convergence iterations [13]. Thus, the conventional numerical optimization techniques such as evolutionary algorithms, which normally require thousands of model evaluations to yield an optimized design [14], are practically unacceptable in high-fidelity CFD simulation model-driven aerodynamic optimization problems, especially for hypersonic vehicle optimization design problems, since that the CFD simulations of high-Mach number are particularly time-consuming [15, 16].
To address the challenge of computation-intensive optimization problems with high-fidelity simulation models, metamodels, also known as surrogate methods, have gain intensive attention in the past decade [17–21]. Metamodel-based optimization uses computationally cheap surrogate models as substitutes to the high-fidelity simulation models to reduce the number of evaluations of computationally more expensive true models, thereby achieving the computation-intensive optimization problems with high-fidelity simulation models at much reduced computational cost [22, 23]. Based on the metamodels, the Sequential Approximate Optimization (SAO) [24, 25] method is proposed and has been regarded as one of the most practical methods for engineering optimization problems and successfully used in various high-fidelity CFD simulation model-driven optimization problems [26–28]. A common SAO includes three stages: design of experiment (DoE), metamodeling, and infilling. DoE is applied to generate a limited number of space-filling sampling points, which are simulated using the high-fidelity simulation model. The metamodels are constructed using existing sampling points by surrogating the relationships between the design variables and the simulations outputs, using specific metamodels, such as the radial basis functions (RBF) [29], the Kriging method [30], and the Response Surface Method (RSM) [31]. The metamodels are then gradually updated through a certain adaptive infilling mechanism during the optimization process to improve the efficiency and convergence performance.
The same as evolutionary algorithms, SAO can be used in any formulated optimization problem for the sake of its high generality and flexibility, where the metamodels are often treated as external modules and the simulation models are often treated as “black-box” functions. The metamodels are only coupled with the input and the output of the simulation, which is the most superficial information of the black-box computation-intensive simulation model [32, 33]. This crude metamodel fails to capture the essence of the true model and cannot make full use of the simulation results sufficiently, which limits the further improvement of optimization efficiency. Apparently, the effectiveness and efficiency of SAO can be further improved in the sense that the coupling between the metamodel and the high-fidelity model can be further enhanced.
In this present paper, a field metamodel-enhanced SAO method is proposed for a hypersonic vehicle aerodynamic optimization design problem. This paper is organized as follows: in Section 2, a conceptual lifting-body hypersonic vehicle is developed and the corresponding aerodynamic optimization problem is formulated. Section 3 elaborates on the key techniques and the overall framework of the proposed optimization method. In Section 4, the optimization design of the hypersonic vehicle is performed using the proposed optimization approach and the results are analyzed and discussed in detail. Finally, a brief conclusion of this research study is summarized in Section 5.
2. Aerodynamic Optimization Design Problem Formulation
The parameterization of the hypersonic vehicle considered in this paper is based on the shape of the after-body section generated using the B-spline method [34]. The preliminary configuration of the baseline is shown in Figure 1(a), and the parameterization of the after-body section is demonstrated in Figure 1(b). The geometric configuration of the considered hypersonic vehicle is constructed by an upper part and a lower part, which are both parameterized using three parameters: the height of the symmetrical plane , the curvature radius at the symmetry axis , and the slope at the connection of the curve and the fillet . The fillet radius of the nose and edge is fixed during the optimization design process in consideration of practical feasibility of structural strength and thermal protection [35].
(a)
(b)
The nonlinear aerodynamic optimization design problem can be formulated as follows: where is the design variable, whose variable type and range are determined by the parameterization method; is the objective function of the aerodynamic optimization problem, e.g., the drag coefficient, heat flux, lift coefficient, and lift-to-drag ratio; and refers to the computationally cheap constraints, which can be obtained without time-consuming simulation calculations, such as the airfoil thickness and volume of vehicle, whereas refers to the computationally expensive constraints, which must be obtained through running simulations, such as heat flux, stagnation temperature, and pitching moment.
The objective of the hypersonic vehicle optimization design problem considered in this paper is to maximize the lift-to-drag ratio in light of the long-range flight capability of hypersonic vehicles. An aerodynamic thermal problem is another concern for hypersonic vehicle performance. Stagnation temperature and heat flux are the two most commonly used reference quantities in the study of aerodynamic thermal problems. In the optimization problem considered in this paper, the stagnation temperature is set as a constraint. A reasonable volumetric efficiency should be maintained for the payload, which should be considered a constraint. Therefore, the hypersonic vehicle optimization design problem can be formulated as follows: where and are the lift and drag coefficients, respectively. Thus, is the lift-to-drag ratio. and are the stagnation temperature and volume, respectively. The subscript refers to the baseline.
3. Field Metamodel-Enhanced Sequential Approximate Optimization
In this section, a field metamodel-enhanced SAO method is proposed and the optimization framework is elaborated in detail. The general framework of the field metamodel-enhanced SAO method is shown in Figure 2.
3.1. Constrained OLHD
SAO applies DoE to get a limited number of space-filling sampling points that are performed using the high-fidelity simulation model to obtain the corresponding outputs. The sampling points with corresponding outputs are then employed to initialize the surrogate models. For optimization problems with constraints, the basic idea of DoE is to sample uniformly in the feasible area of a design space. Optimal Latin Hypercube design (OLHD) is widely applied to generate the sampling set in design space without constraints. However, constraints commonly exist in most real-world design problems like aerodynamic optimization design, which divide the initial -dimensional hypercube design space into a feasible area and an infeasible area. Inevitably, infeasible points would exist when constraints are imposed on the design space. If the infeasible points are regarded as the feasible points equally, lots of unnecessary and invalid computational cost would be wasted, while simply removing the infeasible points will result in insufficient and ununiform sampling points. An effective approach to overcome this shortcoming is to generate enough points directly in the feasible region.
In this paper, a modified OLHD for constrained regions proposed in Wu et al.’s research [36] is adopted in the aerodynamic optimization design of the hypersonic vehicle. The adopted method is based on the common OLHD, in which the sampling points are generated by minimizing the -criterion [37]. The modified OLHD for constrained regions distinguishes feasible points and infeasible points by weighting to change the uniformity index of all sampling points to that of feasible points. Moreover, to generate an adequate number of feasible points, the number of feasible points is added to the optimality criterion of modified OLHD. To further ensure that a specified number of feasible points can be generated, a sequential adjustment for the size of the sampling set is developed.
It is worth to mention that the modified OLHD for constrained regions works only for the computationally cheap constraints, i.e., the volume constraint in the optimization problem of this paper. Although the computationally expensive constraints that need to be obtained through high-fidelity CFD simulations are not in consideration, the method still improves the quality of sampling points by eliminating the infeasible points. In this manner, a large number of infeasible points that do not meet the volume constraint can be eliminated in advance, and the corresponding invalid high-precision simulation calculations can be avoided.
3.2. Field Metamodel
Instead of only focusing on the metamodels of objectives and constraints in general SAO, the field metamodel-enhanced SAO focuses on the metamodel of the distribution of flow field parameters. The flow field distributions are obtained in the aerodynamic analysis stage using the high-fidelity CFD simulation, and then the objective and constraints values can be calculated using the flow field distributions in a postprocessing stage. Apparently, the flow field distribution contains much more detailed information of the high-fidelity CFD simulation models compared with only the objectives and the constraints.
A CFD simulation process incorporates a mesh-generating stage, numerical-solving stage, and postprocessing stage. Generally speaking, a numerical solving stage and postprocessing stage are collectively called the aerodynamic analysis process. The postprocessing stage can be separated from the aerodynamic analysis process analytically; thus, the complexity of the numerical analysis that needs to be surrogated by the metamodel can be reduced. Since the coverage of the metamodel on the black-box simulation model is reduced and more CFD simulation results can be stored and utilized in the approximate process, the prediction accuracy of the metamodel would be increased and the optimization convergence would be sped up. Evidently, the optimization efficiency of the enhanced SAO can be increased through using the field metamodel.
The surface mesh of the hypersonic vehicle is displayed in Figure 3, and local surface mesh nodes are enlarged to better display. The numerical stage obtains the flow field distribution including wall pressure , shear stress , and temperature on the mesh nodes, which are utilized to describe the flow field as shown in Figure 3. In the hypersonic vehicle optimization problem considered in this paper, the flow field is described by the wall pressure , shear stress , and temperature on the surface mesh nodes and can be represented by the field vector as follows: where , , and are the wall pressure, the shear stress, and the temperature of the th representation point, respectively, and is the vector of design variables that control the curvature of the upper and lower surfaces of the hypersonic vehicle.
3.3. Hybrid Metamodel of RBF and PCE
RBF was originally proposed by Hardy [38] to fit irregular topographic contours of geographical data. Since proposed, the RBF metamodel has been utilized for various regression purposes, but it has also been proved inefficient for fitting low-dimensional linear models. Polynomial chaos expansion (PCE) proposed by Wiener provides a framework to approximate the solution of a stochastic system by projecting it onto a basis of polynomials of random inputs [39]. PCE is not ideal for the prediction of the nonlinear problems with high degree of nonlinearity, and it is prone to overfitting when polynomial orders are higher. To give more accurate predictions for the complex hypersonic aerodynamic model, a hybrid metamodel combining the RBF basis and PCE basis proposed by Wu et al. [40] is adopted in this paper. The hybrid metamodel first utilizes PCE for global fitting and then uses RBF for local interpolation. The mathematical formulation is showed as follows: where the first item is used to capture the local details by RBF and the second item provides a low-order global prediction by PCE. It is obvious that Equation (4) has unknown coefficients and .
The undetermined coefficients and can be solved using interpolation conditions and orthogonal conditions.
3.4. RBF Based Mesh Morphing Scheme
In the iterative process of an optimization design, the numerical simulation for a new configuration requires a new mesh in every design iteration. It is practically infeasible to regenerate new meshes automatically during each iteration for the reason that the mesh generating process usually needs to be artificially adjusted to ensure adequate mesh quality for the numerical simulation. Alternatively, the mesh morphing technique can morph the existing mesh to match the new geometry without any human intervention. There are many mature mesh morphing deformation techniques that can ensure reasonable mesh continuity and sensitivity, which can meet the requirements of simulation solution. RBF was firstly applied in the mesh morphing field by Boer et al. in 2007 [41] and has received wide attention and application since then. When the RBF-based mesh morphing techniques are applied in 3D configurations, a considerable number of control points bring expensive computational costs and high memory footprint. In this section, the RBF-based mesh morphing approach is modified by reducing the number of RBF control points using a greedy algorithm.
The mesh of the baseline on the symmetric plane is shown in Figure 4, and the local boundary layer meshes are enlarged. The whole flow field domain of the aerodynamic simulation is defined as the deformable domain . There are two types of mesh points in : the volume mesh points and the surface mesh points . The numbers of the volume and the surface mesh points are defined as and , respectively. Because only a small size of surface points are chosen to be control points, the set of RBF control points .
The RBF interpolating function can be described as where the RBF interpolated vector refers to the displacement for any mesh point . is the radial basis function with respect to the design variables , which is defined as follows: in which and refers to the support radius. The support radius is usually set to be much larger than the largest mesh point displacement in order to ensure the mesh continuity and sensitivity.
A relatively small size of control points is screened from all the surface mesh points via a greedy algorithm, which selected the mesh points with the big geometric displacements. The number of the control points is not fixed and will increase until the Root Mean Square Error (RMSE) between deformed mesh and geometry decreases to a certain criterion. Each control point corresponds to a radial basis function, whose undetermined coefficients are solved by interpolating the control points. The interpolating schemes can be formulated considering all the control points as follows: where refers to the given displacements of the control points, and
The displacements to be solved, including all volume mesh points and the remaining surface mesh points except the control points, can be expressed as follows: in which
and represent and . is the basis function between surface points and RBF control points , and is the basis function between volume points and RBF control points .
4. Results and Discussions
The aerodynamic optimization design problem of the hypersonic vehicle is carried out using the field metamodel-enhanced SAO introduced in this paper. The CFD simulation is performed using the ICEM 16.0 for initial meshing and the Spalart-Allmaras (S-A) turbulence model included in FLUENT 16.0 for numerical solving. The simulation conditions are Mach number 6, a turbulent viscosity ratio of 10, and an angle of attack of 5.71°. The pressure-far-field boundary condition is adopted using the atmospheric parameters of 20 km altitude, where the static pressure is 5529.31 Pa and the background temperature is 216.65 K. The boundary condition of the vehicle surface is modelled as the nonslip and adiabatic wall. Since the vehicle is symmetric about plane because of 0° side slip angle, half of the flow field is modelled to reduce the computational cost. It is worth mentioning that, since there is no project background in this paper, the lifting body optimization design problems, including the lifting body shape, the objective function, and the constraints, are all built by the authors to verify the field metamodel and its enhanced SAO optimization framework proposed in this paper. The setting of working condition and objective function of the optimization case refers to the similar lifting body cases in other papers, and the values of the constraints are determined by the results of unconstrained DoE.
4.1. DoE Results
DoE without considering any constraint is first run using ordinary OLHD, and 20 sampling points evenly distributed in the whole design space are obtained. All 20 sampling points are performed using high-fidelity CFD simulation, and the results are displayed in Figure 5 as circles, whose coordinates are volume, temperature, and lift-to-drag ratio. The volume and the temperature constraints are presented as two translucent faces to distinguish feasible and infeasible points. The volume constraint is set to not less than 0.01, so that half of the ordinary OLHD sampling points are on the left side of and cannot meet the volume constraint according to the results of ordinary OLHD. As for stagnation temperature constraint, the average stagnation temperature of all 20 sampling points, 1825.91 K, is set as an upper bound. As shown in Figure 5, nine sampling points on the right of K do not meet the stagnation temperature constraint, while the remaining 11 sampling points meet the stagnation temperature constraint. As a result, only 8 of the 20 sampling points generated by ordinary OLHD satisfy both volume and stagnation temperature constraints, which makes them feasible sampling points and marked as green circles. The configuration showed in Figure 6 is a typical configuration with a high lift-to-drag ratio of 2.9516 but not satisfying both volume and stagnation temperature constraints with 0.008189 volume and 2040.45 K stagnation temperature. The after-body section of the baseline is also displayed in the black line for comparison, and one can observe that the area of the typical configuration’s after-body section is much smaller than that of the baseline. Aside from the small after-body section, the concave upper surface is also one of the reasons for violating the volume constraint.
(a)
(b)
As a comparison, the sampling points generated by the constrained OLHD are also displayed in Figure 5 as diamonds. As demonstrated in Section 3.1, the constrained OLHD locates the sampling points through minimizing -criterion and maximizing the number of feasible points at the same time. To further ensure that 20 feasible points can be generated, sequential adjustments for size of the sampling set are performed, which is shown in Figure 7. In the first run of the constrained OLHD, only 12 of all 20 sampling points are feasible considering volume constraints, which is still more than 10 volume-feasible points generated by ordinary OLHD, because the number of feasible sampling points is added to the sampling criteria of the constrained OLHD. After 10 iterations of adjustment, 20 feasible sampling points are obtained among all 30 sampling points. In this case, an adequate number of feasible sampling points concerning volume constraint are generated and the computational costs of all infeasible points that violate the volume constraint are saved.
4.2. Optimization Iterative Process
With the 20 feasible sampling points generated using the constrained OLHD, the optimization design of the hypersonic vehicle is carried out considering the lift-to-drag ratio as the objective and volume and stagnation temperature as constraints. Based on the metamodel constructed using the sampling points of DoE, the loop of enhanced SAO iterates to locate the optimum. 30 independent runs are performed using the field metamodel-enhanced SAO method to eliminate randomness. The iterative process of the lift-to-drag ratio is presented in Figure 8. The lift-to-drag ratio of the first candidate in the optimization process is 2.724, which is close to the optimal value of DoE with a lift-to-drag ratio of 2.751. The optimization converged after 30 iterations, and the whole process of the modified SAO costs 50 times of high-fidelity CFD simulation calls. The located optimum has a lift-to-drag ratio of 2.874, which is improved by 27.17% compared with the baseline of 2.260. In the process of iterative optimization, only 2 candidates (indicated by hollow diamonds) violate the stagnation temperature constraint, and no candidate violates the volume constraint, which is attributed to preprocessing of the computationally cheap constraint and the reasonable DoE sampling point set obtained by the constrained OLHD. The 6 design variables of the obtained optimum are listed in Table 1.
4.3. Results Analysis
The geometry shape of the optimum is displayed in Figure 9, and the after-body section of the baseline is also showed in black line for comparison. It can be observed that there are great differences between the optimized geometry and the baseline geometry, and the same situation also occurs between the optimized and the typical configurations, whose after-body section is showed in the red line. The upper half of the baseline is much larger than the lower half, while the optimum is completely opposite. And conversely, the height of lower half of the typical configuration is much bigger than its height of upper half . The short upper part height results in the volume constraint being violated, and the tall lower part results in the stagnation temperature constraint not being satisfied. The heights of the upper and lower parts of the baseline and the optimum are completely different, which is the main reason for their differences in aerodynamic performance.
(a)
(b)
The wall pressure coefficient distributions of the baseline and the optimum are exhibited below in Figures 10(a) and 10(b) for comparison. The baseline and the optimum are presented from different angles to better show their respective surfaces with more complex pressure coefficient distributions. As can be seen from Figure 10, the pressure coefficient distributed on the upper surface of the optimum is far less than that of the baseline, while the distribution on the lower surface of the optimum is a bit higher than that of the baseline. Since less pressure acts on the upper surface of the optimum and more pressure acts on its lower surface, the lift force of the optimum is much higher compared with the baseline.
(a)
(b)
The wall shear stress distributions of the baseline and the optimum are exhibited in Figure 11. It can be seen that the head has very obvious stress concentration in both the baseline and the optimum, which brings extremely high stagnation point temperature and further poses great challenges to the head material. The wall shear stress distribution on the upper surface of the optimum is relatively smaller than that of the baseline, while the wall shear stress distribution on the lower surface shows the opposite performance, which can also be explained by the pressure coefficient distribution displayed in Figure 10. The higher lift of the optimum results from the higher pressure distribution on the lower surface, which also results in higher shear stress distribution.
(a)
(b)
The wall temperature distributions of the baseline and the optimum are displayed in Figure 12. During the flight, the stagnation point will appear in the aircraft head area, which can be seen in the temperature distribution of both the optimum and the baseline. The stagnation temperatures of both the optimum and the baseline are both higher than 1700 K but lower than 1825.91 K; thus, both meet the stagnation temperature constraint. One can observe that the temperature distribution of the optimum is less uniform than that of the baseline, and the wall temperature on the after-body section of the optimum is higher than that of the baseline.
(a)
(b)
The Mach number distribution on the symmetric plane of the optimum and the baseline is exhibited in Figure 13. The gradient in the region near the upper surface has little difference between the baseline and the optimum, but it still can be seen from the local enlarged view of the vehicle tail that the upper surface of the basic configuration compresses the nearby flow field more than the optimum. It can also be seen from the figure that the flow field area affected by the lower surface of the optimum is larger than that of the basic shape. The differences between the Mach number distributions of the baseline and the optimum correspond well to the differences of their pressure coefficient distributions.
(a)
(b)
5. Conclusions
In this paper, a field metamodel-enhanced sequence approximation optimization method is proposed and an optimization framework for hypersonic vehicles is built. A specific hypersonic lift body is constructed and parameterized to demonstrate and verified the enhanced SAO method and optimization framework proposed in this paper.
The corresponding optimization problem is formulated firstly in Section 2. The lift-to-drag ratio is chosen as the objection function, and volume and stagnation temperature are set as constraints. An OLHD method for the constrained design space is adopted, which is available for computationally cheap constraints, such as volume constraint in the optimization problem in this paper. With the constrained OLHD, the sampling point that does not meet the volume constraint is discarded before its high-fidelity CFD simulation is performed, so that the computational costs of those points are saved and the overall optimization efficiency is improved. The field metamodel is proposed to surrogate the flow field distribution instead of only objectives and constraints of the optimization problem so that the complexity of the surrogated part of CFD simulation is decreased. Consequently, the difficulty of metamodeling is decreased, which means that a more accurate metamodel can be built with the same number of sample points and thus the efficiency is improved. A hybrid metamodel combing RBF and PCE is applied to construct the field metamodel, which can provide good approximate performance for different types of engineering models when the method proposed in this paper is applied to other optimization problems. An RBF-based mesh morphing scheme is utilized to morph the existing mesh to match the new geometries during the optimization process.
In Section 4, the optimization results of the hypersonic lift body are shown in detail and discussed in depth. The process data of the experimental design is displayed, which can help to better understand the implementation process of the constrained OLHD. The values of volume and stagnation temperature constraints are set according to the results of the constrained OLHD. Enhanced by the field metamodel, the performance of SAO is greatly improved and the optimum is located with very limited CFD evaluations. The optimization results are carefully analyzed, and the obtained optimum is compared with the baseline. The lift-to-drag ratio of the optimum is improved by 27.17% compared with the baseline. The geometries and the flow field distributions of both the optimum and the baseline are exhibited and discussed, whose wall pressure coefficient distribution, wall shear stress distribution, wall temperature distribution on the surface, and the Mach number distribution on the symmetric plane are all displayed in detail. The differences between the flow field distributions of the baseline and the optimum are analyzed in depth.
The optimization design framework of the field metamodel-enhanced SAO introduced in this paper proves to be powerful and flexible and can be modified to the optimization design of other hypersonic vehicles for drag reduction and thermal protection, which can be helpful to researchers in this field.
Data Availability
The simulation and optimization data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that there is no conflict of interest regarding the publication of this paper.
Acknowledgments
This research was sponsored by the National Natural Science Foundation of China (Grant No. 52005502), the Science and Technology Innovation Program of Hunan Province (2020RC2035), and the Research Project of National University of Defense Technology (Project No. ZK19-11).