NOMENCLATURE
- c, γ, n
prescribed parameter
- ERBF
extended radial basis function
- LHS
Latin hypercube sampling
- np
the number of sampling points
- RBF
radial basis function
- RSM
response surface methodology
- wi
weight factor
- x
the coordinates of the sampling point
- xi
the coordinates of the sampling point xi
- xj
point x along the j − th dimension
- xij
point xj along the j − th dimension
- σi,αLij,αRij,βLij
unknown factor
1.0 INTRODUCTION
Researchers make heavy use of computer simulation codes to replace expensive physical experiments and improve the quality and performance of engineered products and devices in many science and engineering fields. Such simulation activities are collectively referred to as computational science or engineering. Unfortunately, while allowing scientists more flexibility to study the phenomena under controlled conditions, computer simulations require a substantial investment on computation time. One simulation may take many minutes, hours, days or even weeks, quickly rendering parameter studies impractical(Reference Forrester, Sóbester and Keane1-Reference Braun, Kleditzsch and Scharler3). To meet the challenge of increasing model complexity, the process of building approximate models, or surrogate model, has gained wide acceptance from the design community(Reference Braun and Riedel4-Reference Queipo, Haftka and Shyy6). Such surrogate models have been successfully applied to optimisation problems(Reference Rosenow, Lindner and Fricke7,Reference Conway8) . Among the available surrogate model techniques, the Response Surface Methodology (RSM) approach was introduced by Box and Wilson in 1951(Reference Box and Wilson9). This approach was used to optimise fire performance of ultra-low density fiberboards polynomial by Wu(Reference Wu, Huang and Wang10) in 2017. Kriging approach was named by French mathematician Matheron in 1963, after the South African mining engineer Krige, as it is still known in spatial statistics today. The Kriging approach is used to optimise the structure of water axial piston pump and cavitation of plunger cavity in 2016(Reference Sun, Xiao and Xu11). Radial basis function (RBF) methods are effective multi-dimensional approximation approaches. The performances of RBF methods are independent of the dimensionality to an extent(Reference Akhtar and Shoemaker12). An approach of extended radial basis function (ERBF) is proposed by Mullur and Messac(Reference Mullur and Messac13) in 2005. The ERBF approaches are more flexible compare with RSM, Kriging and RBF. Since the ERBF approach is used to build the surrogate model in this research. However, these surrogate model approaches are all have no adaptive capability to the sampling points chosen. In the field of adaptive sampling techniques, Chen(Reference Chen, Qiu and Gao14) used a local adaptive sampling to enhance the efficiency of constructing Kriging models for reliability-based design optimisation problems in 2014. But this approach is not applicable if sampling points are close to the limit state boundaries. An adaptive importance sampling techniques based on stochastic Newton recursions was used to accurately predict the power penalty induced(Reference Remondo, Srinivasan and Nicola15). Li, Wuand Chuang(Reference Li, Wu and Chuang16) apply Stein's Unbiased Risk Estimator (SURE) to adaptive sampling and reconstruction to reduce noise in Monte Carlo rendering. Aerodynamic shape of high-speed train nose is optimised using the adaptive surrogate model(Reference Vytla, Huang and Penmetsa17). An approach for constructing adaptive surrogate models with application in production optimisation problem is proposed(Reference Golzari, Sefat and Jamshidi18). However, these approaches are not applicable, considering that increases new sampling points in the blank region after initial training points are selected using Latin hypercube sampling (LHS). A stochastic process model is proposed by Jones, Schonlau and Welch(Reference Jones, Schonlau and Welch19) based on response surface methodology in 1998. This stochastic process model is widely used in global optimisation. This approach can be used to construct an efficient global optimisation algorithm with a credible stopping rule. However, the aim of this research is to find sample points with the lowest computational cost under specified accuracy. An adaptive sampling approach that new sampling points are placed in the blank area and all the sampling points are uniformly distributed in the design region is used. The Euclidean distances between new sampling points and selected sampling points are evaluated to locate where the new points should be using Multi-Island GA. This adaptive sampling approach combines with ERBF surrogate model is used to optimise aerodynamic and stealthy performance for a cruise missile head shape.
The remainder of the paper is organised as follows. In Section 2, the introduction to methods of optimisation and ERBF is provided. The adaptive surrogate model is built using a novel adaptive sampling approach and ERBF. In Section 3, a cruise missile head shape is optimised using the adaptive surrogate model, and some results are described. In Section 4, our contributions are briefly summarised.
2.0 RESEARCH METHOD
The process of aerodynamic and stealthy optimisation for cruise missile head is divided into two steps in this paper. The first step is to build an adaptive surrogate model. The second step is to combine the adaptive surrogate model with Multi-Island genetic algorithm to determine the optimum shape of the cruise missile head. The first step will be described in this section.
2.1 Optimisation methods
There are several meaningful optimisation methods in the optimised research field but one of the most used algorithms in engineering is the Genetic Algorithm (GA). Genetic algorithms are classical stochastic optimisation algorithms inspired by evolutionary analogy. Because of their robustness and ease of application, genetic algorithms are used for machine learning, automatic control, and so on. Instead of the traditional genetic algorithm, Multi-Island GA is employed for optimisation. In Multi-Island GA, the population is divided into several sub-populations staying on isolated “islands”, whereas traditional genetic algorithm operations are performed on each sub-population separately. A certain number of individuals between the islands migrate after a certain number of generations. Thus, Multi-Island GA can prevent the problem of “premature” by maintaining the diversity of the population(Reference Zhang, Xu and Gao20). In addition, the calculation speed of Multi-Island GA can be greater than that of traditional genetic algorithms.
2.2 Extended radial basis function models
The computer simulations to study and analyse designs are usually very expensive. In order to achieve the result, many computational resources and lots of time are needed. Since the problem of simulation cost becomes more severe. ERBF surrogate model-based design optimisation helps in reducing the number of real computer simulations necessary to solve this problem.
ERBF approaches are the extension of RBF. RBF is expressed according to the Euclidean distance (r = ‖x − xi‖) of a generic point x from a given point data xi which can be mathematically defined as
where c is a prescribed parameter. The radial basis function is a linear combination equation, as described by
where σi is the unknown factor to be solved and np represents the number of sampling points. The preceding equation is expressed in matrix form as follows:
where Aik, σ, F are written as
The interpolation result of RBF for the generic point x is expressed as
The typical RBF approaches provide only an interpolative solving method to the surrogate model problem, but they do not provide patterns for the designer to deliver desirable performance for the meta-models. Mullur and Messac proposed a surrogate model of extended radial basis function(Reference Mullur and Messac13). They defined a coordinate vector as ξij = xj − xij, which is the coordinate of any point x in the design domain relative to the sampling point xi along the j − th dimension, and defined non-radial basis function as
with the functions ϕL, ϕR and ϕβ as described in Table 1.
where γ and n are prescribed parameters, extended radial basis function method is a surrogate model approach that combines radial with non-radial basis functions. It can be expressed as
and also
define:
where the k − th row of B is expressed as
where BLk is expressed as
define
The vector $\overline{\bm \alpha } $ is defined by solving the Equation (20). The Euclidean distances and non-radial basis values that sampling point x relative to all data points (x1, x2, . . ., xnp) into
The interpolation result of ERBF for the generic point x can be expressed as
The ERBF approaches provide the designer with significant flexibility and freedom compared with conventional RBFs in the surrogate model construction process, and their research results show that the ERBF approaches are more accurate in some research field(Reference Mullur and Messac13).
To measure the surrogate model accuracy of the results, the standard error measure is used: Normalised root-mean-squared error (NRMSE). The error measure is defined as follows:
The descriptive error measure represents the error level of the adaptive surrogate model. A tiny value of NRMSE indicates a good fit, whereas a high value of NRMSE indicates a poor fit. This standard error measure as convergence criteria to judge whether the iteration is converged or not.
2.3 An adaptive surrogate model for multi-objective optimisation problem
The adaptive surrogate model is built on an adaptive sampling approach and an extended radial basis function (ERBF) in our recent study(Reference Guo, Ang and Cai21). We get the same surrogate model accuracy with minimal time and resources using this approach. The adaptive sampling is an approach in which new sampling points are placed in the blank area and all the sampling points are uniformly distributed in the design region using Multi-Island GA. The flow chart of the adaptive surrogate model approach for the aerodynamic and stealthy multidisciplinary optimisation is shown in Fig. 1. The original sampling points for the ERBF surrogate model building is picked using the LHS technique. These original sampling points are used for aerodynamic and stealthy numerical evaluation, and the results of numerical evaluation are used to fit the surrogate model of ERBF. Standard error measure as convergence criteria to judge whether the iteration is converged or not. A new sampling point is searched using the Multi-Island GA technique. This new sampling point is added to the numerical evaluation if the convergence condition is not satisfied. The Euclidean distances between the new sampling point and the selected sampling points are evaluated during the process of searching new sampling point. The new sampling point is added to the numerical evaluation if the iteration is convergent in the process of searching new sampling point. Otherwise, update the Multi-Island GA population. Black box of surrogate model is fitted if the convergence condition is satisfied. The black box of the surrogate model is used to find the response in the multi-objective optimisation. For this multi-objective optimisation problem, the linear weighted sum method is employed as follows:
Go to finished if the iteration is converged; otherwise, update the Multi-Island GA population.
The adaptive sampling approach is illustrated by a sampling case of two design variables. The circular points are the initial sampling points, and the square points from 1 to 11 are ready to add to the numeral calculations in Fig. 2. The square points are selected by the adaptive sampling approach. Points 1, 2, . . ., 11 are selected one by one using the Multi-Island GA. Point 1 is selected first, because the minimum Euclidean distance of Point 1 to all initial sampling points is maximum compared to Points 2, 3,. . ., 11. New sampling points are placed in the blank area and all the sampling points are uniformly distributed in the design region using the adaptive sampling approach. So the most appropriate sample points are chosen to improve the precision of the surrogate model.
3.0 AERODYNAMIC AND STEALTHY OPTIMISATION FOR A CRUISE MISSILE HEAD USING THE ADAPTIVE SURROGATE MODEL
In the aerodynamic and stealthy design process for cruise missile head, the aerodynamic drag and the RCS value should be reduced. This is a multidisciplinary optimisation problem. In this section, the aerodynamic characteristics and the radar target characteristics of the cruise missile head are considered, and aerodynamic and stealthy optimisation for the cruise missile head is studied using the adaptive surrogate model and Multi-Island GA algorithm.
3.1 CAD model of cruise missile head and numeral calculations
The shape of the cruise missile head is shown in Fig. 3. The shape is controlled by the control point. In order to ensure the control point is not too low or high, the maximum value of z is equal to 270 mm. The minimum value of z is equal to 0 mm. For the purpose of the cruise missile head size is not too long or too short. The maximum of x is equal to 1500 mm, and the minimum of x is equal to 500 mm. The shape of the cruise missile head ensures that the curve and the surface are smooth.
The velocity of the cruise missile is equal to 0.65 M. The flying height is equal to 6000 m. In order to ensure that the subsonic CFD numerical calculation is correct, a complete aircraft fuselage is used as shown in Fig. 4. However, we only monitor the drag of the cruise missile head in the CFD numerical calculation.
According to the windward area of the cruise missile head and the air conditions, the value of the feature size is equal to 0.773 m, and the value of the Reynolds number is equal to 6.4×106. In the grid generation process, the thickness of the first layer of the boundary layer will have an impact on the drag. The value of y+is equal to 30. According to the boundary-layer empirical formula, the value of first layer thickness is equal to 0.0001 m. The grid of the cruise missile head is shown in Fig. 5. The number of grid cells is equal to 3 m. Velocity-inlet and pressure-outlet as the boundary conditions of the CFD calculation.
Because of the different power and wavelength of radar, the detection distance is different. The radar detection distance is up to 4000 km or more. In order to ensure that the cruise missile is not detected by the ground radar at cruise, the RCS of the missile head in the range of 0° to 5° should be small. So the pitch angle ranges from 0° to 5° in the RCS numerical calculation as shown in Fig. 6. A radar-wave frequency of 9 GHz as the input wave in the RCS numerical calculation. The calculation fields in the plane-wave incident direction only be considered. Horizontal polarization as the polarization direction. The RCS average value, which is calculated once every 0.1 degrees in the range from 0° to 5°.
3.2 Problem set-up
The Pareto front can be found using the Pareto ranking method in this research. However, the main issues we are concerned with are the drag coefficient and the RCS value, with the absolute value of the RCS value divided by the drag coefficient. Concerning the drag coefficient and the RCS value, the linear-weighted sum method is employed as follows:
where F(X) is the global objective function and f 1(X) and f 2(X) are the drag coefficient and the RCS value, respectively. The w 1 and w2 are the weight factors. In order to more clearly distinguish between aerodynamic performance and stealth performance, two kinds of weight factors are only considered in this paper. One is that w 1 = 1 and w 2 = 0. Another one is that w 1 = 0 and w 2 = 1. The two kinds of weight factors consider aerodynamic performance and stealth performance, respectively. Concerning the absolute value of the RCS value divided by the drag coefficient, there is a global objective function as follows:
The f 2(X) is a negative number. The smaller the value of f 2(X), the stronger the stealthy capability. The f 1(X) is a positive number. The smaller the value of f 1(X), the better the aerodynamic performance. Minimising the F(X) for Equation 26, and maximising the F(X) for Equation 27 is the direction in this optimisation.
The geometry constraints as follows:
3.3 Results
The initial sampling points required to build the ERBF surrogate model are selected using the LHS approach. The whole number of 18 training points is picked to build the original surrogate model as shown in Fig. 7. The number of original training points may be fewer or more. The 18 training points selected here are artificially defined as examples. The adaptive algorithm is able to check the accuracy of the surrogate model by comparing the value of NRMSE and the request for additional sampling points whenever necessary. At 1.5% confidence bounds, the adaptive algorithm requested a supernumerary 13 sampling points to improve the accuracy of the ERBF surrogate model. The confidence bound of 1.5% is satisfied using a total of 31 sampling points.
The selection of the weights for the multi-objective functions will have a greater effect on the final optimum drag and RCS, since the weights will decide the locations of the optimum design variables(Reference Peng, Wu and Yi22). Aerodynamic and stealthy optimisations are studied using different weights. The optimisation history of the multi-objective function with the first weights (w 1 = 1, w 2 = 0) for the first objective function (Equation 26) is illustrated in Fig. 8. The aerodynamic optimisation history shows the decrease of F(X). The number of islands is equal to 10 in the Multi-Islands GA algorithm, since there are ten locally optimal solutions. But the globally optimal solution has only one, the value of F(X) equal to 0.037. The values of x and z are equal to 1337.8 mm and 30.9 mm, respectively. The location of the control point is marked using a red circle in Fig. 8.
The optimisation history of the multi-objective function with the second weights (w 1 = 0, w 2 = 1) for the first objective function (Equation 26) is illustrated in Fig. 9. The stealthy optimisation history shows the decrease of F(X). The number of islands is equal to 10 in the Multi-Islands GA algorithm, since there are ten locally optimal solutions. But the globally optimal solution has only one, the value of F(X) equal to 35.8 dBm2. The values of x and z are equal to 640.2 mm and 0.6 mm, respectively. The location of the control point is marked using a red circle in Fig. 9.
The optimisation history of the multi-objective function for the second objective function (Equation 27) is illustrated in Fig. 10. The aerodynamic and stealthy optimisation history shows the increase of F(X). The number of islands is equal to 10 in the Multi-Islands GA algorithm, since there are ten locally optimal solutions. But the globally optimal solution has only one, the value of F(X) equal to 880.6. The values of x and z are equal to 1339.3 mm and 9.6 mm, respectively. The location of control point is marked using a red circle in Fig. 10.
Stealthy and aerodynamic performance are researched between the original model and the optimised model in this section. In the original model, the values of x and z are equal to 800 mm and 200 mm, respectively. The value of w 1 equal to 1 and the value of w 2 equal to 0 in Equation 26 means that aerodynamic performance is considered only in the optimisation process. The cruise-missile-head pressure distributions of the original model and the optimised model are shown in Fig. 11. The left one is the pressure distribution of the optimised model (x = 1337.8 mm, z = 30.9 mm), and the right one is the pressure distribution of the original model (x = 800 mm, z = 200 mm). The high-pressure area of the original model is significantly larger than that of the optimised model. The drag coefficient of the optimised model is equal to 0.036 and the drag coefficient of the original model is equal to 0.064.
The cruise-missile-head streamlines of the original model and the optimised model are shown in Fig. 12. The left one is the streamline of the optimised model (x = 1337.8 mm, z = 30.9 mm), and the right one is the streamline of the original model (x = 800 mm, z = 200 mm). The airflow at the head of the cruise missile has been separated in the original model. The airflow at the head of the cruise missile is transferred smoothly in the optimised model.
The value of w 1 equal to 0 and the value of w 2 equal to 1 in Equation 26 means that stealth performance is considered only in the optimisation process. The RCS values of the original model and the optimised model are shown in Fig. 13. The blue line is the RCS value of the optimised model (x = 640.2 mm, z = 0.6 mm) and the black line is the RCS value of the original model (x = 800 mm, z = 200 mm). The blue line is under the black line completely. The RCS average value of the optimised model is equal to 35.8 dBm2 and the RCS average value of the original model is equal to 31.7 dBm2.
The global objective function (Equation 27) means that stealthy and aerodynamic performance are all considered in the optimisation process. The pressure distribution three-dimensional views of the original model and the optimised model are shown in Fig. 14. The left side is the pressure distribution three-dimensional view of the original model (x = 800 mm, z = 200 mm). The right side is the pressure distribution three-dimensional view of the optimised model (x = 1339.3 mm, z = 9.6 mm). The top view of the pressure distribution indicates that the airflow is separated in the original model and the airflow is smooth in the optimised model. The front view of the pressure distribution indicates that the high-pressure area of the original model is significantly larger than that of the optimised model. The drag coefficient of the optimised model is equal to 0.038.
The RCS values of the original model and the optimised model are shown in Fig. 15. The red line is the RCS value of the optimised model (x = 1339.3 mm, z = 9.6 mm) and the black line is the RCS value of the original model (x = 800 mm, z = 200 mm). The RCS average value of the optimised model is equal to −32.3 dBm2.
The above results show that the choice of the weights for the multi-objective functions has a major effect on the final optimum results. The shape of optimised cruise missile head is shown in Fig. 16. When the optimised shape of the first weight factor (Equation 26, w 1 = 1, w 1 = 0) and the optimised shape of global objective function (Equation 27) are compared, we find that the two optimised shapes are similar. The value of x is equal to 1339.3 mm and the value of z is equal to 9.6 mm for the global objective function (Equation 27), and the value of x is equal to 1337.8 mm and the value of z is equal to 30.9 mm for the global objective function (Equation 26, w 1 = 1, w 1 = 0). But optimised shape of second weight factor (Equation 26, w 1 = 0, w 1 = 1) is different from the above two optimised shapes. The value of x is equal to 640.2 mm and the value of z is equal to 0.6 mm for the global objective function (Equation 26, w 1 = 0, w 1 = 1).
4.0 CONCLUSIONS
An adaptive sampling approach in which new sampling points are placed in the blank area and all the sampling points are uniformly distributed in the design region is used. Multi-objective optimisation of the cruise missile head shape by considering two objectives: drag coefficient and RCS value.
Minimising the drag coefficient and RCS value is performed using Multi-Island GA algorithm. The executable of optimising using the adaptive surrogate model in combination with Multi-Island GA optimisation approach is demonstrated successfully.
It is observed that a total of 18 LHS sampling points are picked to calculate the original ERBF surrogate model. The optimisation presented in this paper required 13 sampling points for 1.5% confidence bound.
The effect of the selection of different weights on the optimum cruise missile head shape is researched. Considering drag coefficient only, the drag coefficient is optimised to 0.036. Considering RCS only, the value of RCS is optimised to −35.8 dBm2. Considering drag coefficient and RCS all, the ratio of drag coefficient to the absolute value of RCS is optimised to 850, and the drag coefficient is optimised to 0.038, the value of RCS is optimised to −32.3 dBm2.
ACKNOWLEDGEMENTS
Project supported by the Natural Science Foundation of Jiangsu Province, China (Grant No. BK20160817) and the Fundametnal Research Funds for the Central Universities (Grant No. 30915118807, 30917011302).