Structural Optimization and Stability Analysis for Supercavitating Projectiles
https://doi.org/10.1007/s11804023003528

Abstract
Stability is the key issue for kineticenergy supercavitating projectiles. Our previous work established a six degrees of freedom (DOF) dynamic model for supercavitating projectiles. However, the projectile's structure did not meet our current design specifications (its sailing distance could reach 100 m at an initial speed of 500 m/s). The emphasis of this study lies in optimizing the projectile's configuration. Therefore, a program was developed to optimize the projectile's structure to achieve an optimal design or the largest sailing distance. The program is a working optimal method based on the genetic algorithm (GA). Additionally, the convergence standard and population producing strategy were improved, which greatly elevated the calculation speed and precision. To meet design specifications, the improved GA was combined with the 6DOF model, which establishes a dynamic optimization problem. The new projectile's structure was obtained by solving this problem. Then, the new structures' dynamic features were compared with the ideals proposed in this paper. The criterion of stability, which is called weakened selfstability, was redefined based on the results. The weakened selfstability is the optimal stability for an actual kinetic projectile motion, and it is instructive for the design of supercavitating projectiles in the future.Article Highlights● Derived a universal parameterized model for supercavitating projectiles.● Established the structural global optimization problem of supercavitating projectile.● Revealed the self stabilization mechanism of supercavitating projectiles. 
1 Introduction
For a long time, many scholars and engineers have carried out extensive research on drag reduction and cavities. Underwater vehicles are developing in the direction of long sailing distances, high speed, and intelligence (Ashley, 2001; Savchenko, 2001). To consider a vehicle's index of long sailing distance and high speed, researchers usually optimize the design to meet the requirements. The supercavitating technology represents a successful result of such optimization. Supercavitating vehicles use a special hydrodynamic layout, which makes full use of supercavity features. The drag is greatly reduced as if flying in the air. The kinetic energy supercavitating projectile is efficient in drag reduce. It has a high or ultrahigh speed and can achieve precise and fast damage for targets in water.
The former Soviet Union used cavitation to increase the torpedo velocity to nearly 100 m/s (Zhang et al., 2014). This supercavity drag reduction technology is a milestone for underwater vehicles. However, the singlephase flow becomes a multiphase flow, which has many disturbance factors and is difficult to analyze and calculate in theory (Wei et al., 2013). When a vehicle moves with supercavity, most of the forces that act on it will disappear. The force of the vehicle becomes highly nonlinear, which often leads to deflection and loss of control. These difficulties are an unprecedented challenge for mathematical modeling and stability research of vehicles (Ruzzene and Soranna, 2004; Zhang et al., 2010; Wei et al., 2012; He et al., 2013; Mirzaei et al., 2015).
The stability analysis method mainly uses mathematical modeling or experiments. The experiments of kinetic energy supercavitating vehicles usually require substantial human and financial resources, which is inconvenient.
Theoretical modeling is highly imperative to determine a vehicle's motion features. The model can reveal the mechanism of a vehicle's motion in theory. Rand et al. (1997) proposed the simplistic supercavitating vehicle model, which mainly focuses on the frequency of tail slaps. The body coordinate origin is located at the cavitator, and it supposes that the vehicle rotates around the cavitator on the plane. The results indicate that the frequency of tail slaps will decrease as speed decreases. In addition, the frequency of tail slaps highly depends on the initial condition. Richard Rand's model is great in simplification, and it is deemed instructive for the later supercavitating vehicle research. He et al. (2013) combined Richard Rand's theory and finite element method (FEM) and established a procedure for supercavitating vehicles. The results revealed not only the dynamic features but also the loads of tail slaps. The findings also highly accorded with those of Richard Rand. Zhao et al. (2018, 2019) proposed the hypothesis of mirror tail slap, which is a tail slap model based on Richard Rand's work. The tail slaps were simulated in the cavity cross section at the vehicle's tail position. The results also showed the trajectory of this vehicle. This previous work is preliminary research for stability.
The above works were based on Richard Rand's theory. However, Richard Rand's model is based on the single degree of freedom (DOF) vibration theory. A limitation of this model is that it cannot reveal the mechanism of stability in depth.
Therefore, Kulkarni and Pratap (2000) proposed a 3DOF model in the longitudinal plane as early as 2000. In this model, the motions of tail slaps and those in supercavity are considered separately. In addition, the model is highly nonlinear, which was verified by the computational fluid dynamics (CFD) method. Two kinds of vehicles are simulated by this model. The results indicated that although tail slapping occurred, the trajectory of the vehicle remained linear. This work also proposed an equation about the tail slap angular velocity. Their proposed model is the foundation of many works on stability. To gain a qualitative understanding of stability, Savchenko (2001) proposed four criteria for judging the stability of supercavitating vehicles. These criteria covered the velocity from 50 to 1 000 m/s, which was based on numerous experimental results, and included doublecavity, sliding, tail slapping, and aerodynamichydrodynamic stabilities. Savchenko's results were very institutive for future works. Sail S. Kulkarni and Savchenko provided pioneering works for the study of stability.
The 3DOF dynamic model of Sail S. Kulkarni is considerably better than Richard Rand's work. However, this 3DOF model excludes the influence of attack angle. Zhou et al. (2016) added the influence of attack angle and many uncertainty influences during motion. Moreover, the uncertainty of trajectory was analyzed based on these works. The results were verified by experiments in the literature. This previous work was considered instructive for engineering. Nguyen Thai et al. (2018) established an accurate 3DOF model and evaluated it by experiment. This model is based on exterior ballistics, and it is more scientific compared with previously proposed models. The error was nearly 1.1%. Their work verified for the first time a mathematical model using experiments. Wang et al. (2020) proposed a 6DOF model for kinetic supercavitating vehicles. This model was verified by three kinds of vehicles. In addition, the dynamic features were explored and showed a special dynamic behavior in space. Finally, a new stability evaluation method was proposed, and it was meaningful for theoretical supercavitating research.
Multiple analysis methods have been proposed for years. As a result, the 6DOF dynamics model is more useful than others in exploring the features of kinetics and dynamics and can be used to analyze the mechanism of instability in theory. For kinetic energy supercavitating vehicles, the stability is influenced by the vehicle structure, hydrodynamic environment, and initial conditions. A fluctuating initial condition in a small sailing distance has a minimal impact on stability unless other unpredictable disturbances occur when the vehicle is launched. These factors directly determine whether the target can be destroyed. The vehicle's structure is the most important factor affecting its performance. On the other hand, stability is a key index to measuring the hydrodynamic layout of supercavitating vehicles and is the basis of the following studies. However, quantitative analysis of vehicle stability in the past was difficult. The stability of supercavitating vehicles by experiments is very inconvenient and costs a huge amount of manpower, material, and financial resources.
Most scholars believe that large supercavitating vehicles (LSVs) and small supercavitating vehicles (SSVs or projectiles) have the same dynamic model or behavior. However, the dynamic behavior and model were dramatically different from recent experiments on water. These test data were from authors' previous works and relative literature. From the tail slap angular velocity view, the angular velocity of SSVs is ten times more than that of LSVs. The difference was considered in depth. The established dynamic model also differed. For instance, the supercavity calculation strategy must consider the physical dimension. The memory effect and cavity drift features are considered imperative for LSVs but not for SSVs. These differences will make the dynamic models diverse, which leads to different optimization strategies. The early optimization problem mainly focuses on cavitators because they are important for supercavity configuration. Scholars have conducted many works on the cavitator area, and the research methods used varied. Choi et al. (2005) used the potential theory to optimize the cavitator to minimize the drag, and the cavity shape was calculated. The author's work is based on a gradient optimization algorithm, which can simultaneously optimize the cavitator and cavity. Jiang et al. (2009) converted the freeboundary value problem of supercavities into an equivalent shape optimal problem. The objection function was defined by the integral square error of pressure difference. Then, the optimization of the cavitator and cavity boundary calculation were merged into a multiobjective optimization problem. However, the author designed the cavitator without considering the vehicle's structure layout. Ahn et al. (2010) optimized an LSV with a control system that included configuration, mass distribution, and manipulator. A limitation of this work was the vehicle velocity (170–200 m/s). Shafaghat et al. (2011) defined a multiobjective optimization problem to optimize the cavitator shape. Based on the flow characteristic and cavitator model, the design parameters and constraint conditions were obtained. Then, Nondominated Sorting Genetic AlgorithmsⅡ was employed to optimize this issue. The results were compared with those of the classical optimal method, which showed that the coneshaped cavitator was optimal in the current condition. Mirzaei and Taghvaei (2019) optimized the configuration of a ventilated supercavitating vehicle. The layout was compared between a natural supercavitating vehicle and its ventilated counterpart. The results showed that the ventilated case had a smaller cavity radius than the natural one. The supercavitating vehicle design is aimed at the cavitator structure. The method is mainly based on potential theory and CFD. The dynamic method is rarely used. In addition, some researchers combined the optimal algorithm with FEM, which involves complicated computations; thus, the method is inconvenient to use in engineering.
In our previous work (Wang et al., 2020), we put forward the quantitative criterion for the stability of supercavitating projectiles, and it can be used to predict a projectile's motion by the numerical method. Thus, the hydrodynamic layout of this projectile can be changed quickly, and the resources can be saved. In this work, we mainly focused on solving the optimization problem of highspeed projectiles. A complete program about the projectile's optimal condition, which is also based on genetic algorithm (GA), was explored. The program has a good generality and can meet the accuracy of engineering. An optimal hydrodynamic layout was obtained by the 6DOF model based on the improved GA. For guiding the engineering design, the dynamic features of this optimized projectile were analyzed by the 6DOF model. The results showed that this projectile has a selfstability tendency, which is a special motion. The selfstability tendency is closest to the ideal motion and the optimal design results. This structure is instructive for the design of kinetic energy supercavitating projectile.
2 Mathematical model of supercavitating projectiles
A complete mathematical model was proposed in Wang et al. (2020). Therefore, all the reference coordinate system definitions, coordinate transformation, parameters, and formulation derivation can be referred to in Wang et al. (2020). This work is a continuation of previous work. The supercavitating projectile's research needs not only experiments but also a mathematical model. The quantitative description of the projectile's velocity and attitude is important to study the various phenomena during motion, and the model can reveal the dynamic features in nature. In addition, a reasonable model is required for further research about the projectile's motion. To simplify this model, we analyzed the body coordinate system (Ox_{1}y_{1}z_{1}). The motion equations are listed as follows:
a. Dynamic equations:
$$ \begin{aligned} & m\left(\frac{\mathrm{d} V_{x 1}}{\mathrm{d} t}+\omega_{y 1} V_{z 1}\omega_{z 1} V_{y 1}\right)= \\ & \quad\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) \cos \alpha \cos \beta \\ & \quad\frac{1}{4} C_{x f} \rho V_{x 1}^2 S_w\left[\tanh \left(k_n h_{\mathrm{ef}}\right)+1\right]m g \sin \psi \end{aligned} $$ (1) $$ \begin{aligned} & m\left(\frac{\mathrm{d} V_{y 1}}{\mathrm{d} t}+\omega_{z 1} V_{x 1}\omega_{x 1} V_{z 1}\right) \\ & =\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) \sin \alpha \cos \beta+ \end{aligned}\\ \tanh \left(k_n \psi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) . \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\left[\tanh \left(k_n h_{\mathrm{ef}}\right)+1\right]\right.} \end{array}\right] \\ \sin \eta \cos \theta+ \\ \tanh \left(k_n \varphi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) \cdot \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\left[\tanh \left(k_n h_{\mathrm{ef}}\right)+1\right]\right.} \end{array}\right] \\ \cos \eta \sin \theta+m g \cos \psi \cos \theta $$ (2) $$ \begin{aligned} & m\left(\frac{\mathrm{d} V_{z 1}}{\mathrm{d} t}+\omega_{x 1} V_{y 1}\omega_{y 1} V_{x 1}\right) \\ & =\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) \sin \beta \end{aligned} \\ \tanh \left(k_n \psi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) \cdot \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\left[\tanh \left(k_n h_{\mathrm{ef}}\right)+1\right]\right.} \end{array}\right] \\ \sin \eta \sin \theta \\ +\tanh \left(k_n \varphi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) . \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\left[\tanh \left(k_n h_{\mathrm{ef}}\right)+1\right]\right.} \end{array}\right] \\ \cos \eta \cos \theta+m g \cos \psi \sin \theta $$ (3) $$ \begin{gathered} I_{x 1} \frac{\mathrm{d} \omega_{x 1}}{\mathrm{d} t}+\omega_{y 1} \omega_{z 1}\left(I_{z 1}I_{y 1}\right)= \\ R \tanh \left[k_r \mu_c \ln \left[\exp \left(\frac{\omega_{x 1}}{\mu_c}\right)+\exp \left(\frac{\omega_{x 1}}{\mu_c}\right)\right]\right] \cdot \\ {\left[\frac{1}{2}\left[\tanh \left(k_n h_{\mathrm{ef}}\right)+1\right] K_2+K_1\right]} \end{gathered} $$ (4) $$ K_1=\frac{1}{2} C_{x 0} \rho \omega_{x 1}^2 R_n^4 \pi $$ (4a) $$ K_2=\frac{1}{2} C_{x \mathrm{fr}} \rho \omega_{x 1}^2 R^2 S_w $$ (4b) $$ \begin{aligned} & I_{y 1} \frac{\mathrm{d} \omega_{y 1}}{\mathrm{d} t}+\omega_{x 1} \omega_{z 1}\left(I_{x 1}I_{z 1}\right) \\ & \quad=\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) x_{\mathrm{cm}} \sin \beta+ \\ & \quad\left[\tan h\left(k_n \psi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) \cdot \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\right]\left[\tan h\left(k_n h_{\mathrm{ef}}\right)+1\right]} \end{array}\right]\right. \\ & \sin \eta \sin \theta+ \\ & \tan h\left(k_n \varphi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) \cdot \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\right]\left[\tan h\left(k_n h_{\mathrm{ef}}\right)+1\right]} \end{array}\right] \\ & \cos \eta \cos \theta]\left(Lx_{\mathrm{cm}}\right) \end{aligned} $$ (5) $$ \begin{aligned} & I_{z 1} \frac{\mathrm{d} \omega_{z 1}}{\mathrm{d} t}+\omega_{x 1} \omega_{y 1}\left(I_{y 1}I_{x 1}\right) \\ & =\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) x_{\mathrm{cm}} \sin \alpha \cos \beta \\ & {\left[\tan h\left(k_n \psi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) . \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\right]\left[\tan h\left(k_n h_{\mathrm{ef}}\right)+1\right]} \end{array}\right]\right.} \\ & \end{aligned} \\ \begin{aligned} & \sin \eta \cos \theta+ \\ & \tan h\left(k_n \varphi\right)\left[\begin{array}{l} \frac{1}{2} \pi R^2(t) \rho V^2 \sin \chi \cos \chi\left(\frac{R+h_{\mathrm{ef}}}{R+2 h_{\mathrm{ef}}}\right) . \\ {\left[1\left(\frac{R(t)R}{h_{\mathrm{ef}}+R(t)R}\right)^2\right]\left[\tan h\left(k_n h_{\mathrm{ef}}\right)+1\right]} \end{array}\right] \\ & \cos \eta \sin \theta]\left(Lx_{\mathrm{cm}}\right) \end{aligned} $$ (6) b. Kinetic equations:
$$ \omega_{x 1}=\frac{\mathrm{d} \theta}{\mathrm{d} t}+\frac{\mathrm{d} \varphi}{\mathrm{d} t} \sin \psi $$ (7) $$ \omega_{y 1}=\frac{\mathrm{d} \varphi}{\mathrm{d} t} \cos \psi \cos \theta+\frac{\mathrm{d} \psi}{\mathrm{d} t} \sin \theta $$ (8) $$ \omega_{z 1}=\frac{\mathrm{d} \psi}{\mathrm{d} t} \cos \theta\frac{\mathrm{d} \varphi}{\mathrm{d} t} \cos \psi \sin \theta $$ (9) $$ \begin{aligned} \frac{\mathrm{d} x}{\mathrm{d} t}= & V_{x 1} \cos \psi \cos \varphi+ \\ & V_{y 1}(\sin \psi \cos \varphi \cos \theta+\sin \varphi \sin \theta)+ \\ & V_{z 1}(\sin \psi \cos \varphi \sin \theta+\sin \varphi \cos \theta) \end{aligned} $$ (10) $$ \frac{\mathrm{d} y}{\mathrm{d} t}=V_{x 1} \sin \psi+V_{y 1} \cos \psi \cos \thetaV_{z 1} \cos \psi \sin \theta $$ (11) $$ \begin{aligned} & \frac{\mathrm{d} z}{\mathrm{d} t}=V_{x 1} \cos \psi \sin \varphi+ \\ & \quad V_{y 1}(\sin \psi \sin \varphi \cos \theta+\cos \varphi \sin \theta) \\ & \quad+V_{z 1}(\sin \psi \sin \varphi \sin \theta+\cos \varphi \cos \theta) \end{aligned} $$ (12) c. Geometric equations:
$$ \alpha=\arctan \left(\frac{V_{y 1}}{V_{x 1}}\right) $$ (13) $$ \beta=\arcsin \left(\frac{V_{z 1}}{V}\right) $$ (14) d. Supercavity equations:
$$ \left\{ \begin{array}{l} \begin{aligned} & \sigma(t)=\frac{2\left[P_aP_c+\rho g(hy)\right]}{\rho V^2} \\ & R(t)=R_c(t) \sqrt{1\left(1\frac{R_n^2}{R_c^2(t)}\right)\left(1\frac{2 x_1}{L_c(t)}\right)^{2 / \eta_c}} \\ & L_c(t)=2 R_n \frac{\sqrt{C_0 \ln \sigma(t)}}{\sigma(t)} \end{aligned}\\ \begin{aligned} & R_c(t)=R_n \sqrt{\frac{C_0}{k \sigma(t)}} \\ & k=1 /\left[1+2 \frac{\ln (2 / \sqrt{e})}{\ln (5 / \sigma(t))}\right] \\ & L_p=L \cos \chi=L \cos \varphi \cos \psi \end{aligned} \end{array} \right. $$ (15) e. Hydrodynamics coefficient equations:
$$ C_0=C_{x 0}(1+\sigma) $$ (16) f. Angle equations:
$$ \cos \delta=\cos \alpha \cos \beta $$ (17) g. Motion converting equations:
$$ h_{\mathrm{ef}}=R+L \sin \chiR(t) $$ (18) Eqs. (1)–(18) show the scalar equations used to describe the motion of supercavitating projectiles. These equations are strongly nonlinear and ordinary differential equations. In the process of establishing the above equations, given the limitations of the problems studied, the applicable range of velocity was between 100 and 900 m/s. Moreover, considering the low gas density in the supercavity, the aerodynamic force of the projectile in the cavity was ignored. The motion equation was continued by a hyperbolic function for the convenience of analysis.
Wang et al. (2020) validated the 6DOF model using four methods.
In the first verification case, we employed an experiment on kinetic energy supercavitating projectiles. The supercavity model was validated by the test data. The maximum test sailing distance served as a preliminary validation for this model.
In the second verification case, we compared the model results with those of the classical model, which is a single DOF model. The 6DOF supercavitating model is a precise model that considers many factors. However, the classical model only considers the drag of the cavitator. The displacement of the 6DOF model is less than that of the classical model's counterpart in theory. The 6DOF model greatly rectifies the classical model. Their velocity curves have high similarity. Therefore, the 6DOF model correctly predicts the projectile's motion underwater.
In the third and fourth verification cases, we used two groups test results and projectile data recorded in Zhao et al. (2014) and Nguyen Thai et al. (2018). The results of the test and numerical calculation were consistent, which proves the accuracy of the 6DOF model.
3 Stability of supercavitating projectile
3.1 Stability criterion
The stability of the kinetic energy supercavitating projectile is a key factor in evaluating its ballistics capability. Quantitative research is a significant method for the stability of supercavitating projectiles. When h_{ef} is equal to 0, tail slapping occurs exactly. Therefore, the limitation angle of tail slaps χ_{lim} is defined as follows (Wang et al., 2020):
$$ \chi_{\lim }=\arcsin \left[\frac{R(t)R}{L}\right] $$ (19) The limitation angle of tail slaps χ_{lim} reflects the relative position of the projectile and the supercavity at the projectile's tail location. When χ, which is the angle between the supercavitating body axis and the horizontal line that passes through the mass center, is greater than the χ_{lim}, tail slaps will occur, and vice versa. Therefore, depicting χ and χ_{lim} on the same figure is meaningful. The χ−χ_{lim} curve has the following features (Figure 1):
1) The χ−χ_{lim} curve can reflect the stability of the supercavitating projectile. If the projectile has stability with tail slaps, then the χ curve will sustain oscillation and pass through the χ_{lim} curve.
2) The χ_{lim} was used for the boundary. The proportion of χ below or above the χ_{lim} reflects the degree of motion stability.
When the upper part of χ, which is above χ_{lim}, is greater than the counterpart below, serious tail slaps of the projectile occur. This finding indicates that the projectile's body wetting would increase, and the drag is also increase. These factors will enhance the rate of attitude, reduce the displacement, and deflect the trajectory. In other words, the stability is relatively poor. However, when χ is totally above χ_{lim}, the supercavity cannot reduce the drag, and the projectile loses its stability. However, when the χ is under χ_{lim} completely, the projectile moves into the supercavity without making contact with it, which is an ideal state for supercavitating projectile. Thus, the projectile can balance the force and moment simultaneously due to the attack angle change.
3.2 Ideal stability for supercavitating projectile
Savchenko (2001) believes that the motion of a supercavitating projectile in the cavity varies with speed. Three cavitation schemes can be considered (v < 70 m/s): the stable sliding along the inner surface of the cavity (v = 100–200 m/s), mutual impact with the cavity wall (v < 900 m/s), and aerodynamic interaction with the cavitating steam splashing medium (v > 1 000 m/s). Four different stability mechanisms occur in turn. Theoretically, the ideal motion mode of a supercavitating projectile is the complete movement of the projectile into the cavity during motion. The definition of selfstability is introduced here as follows:
Selfstability: In the underwater motion of a supercavitating projectile, only the cavitator is in contact with water, and the rest is completely in the supercavity.
Section 3.1 presents a special motion mode of supercavitating projectile. In this case, the kinetic energy supercavitating projectile has selfstability. However, the ideal motion is hard to achieve due to the unpredictable disturbance in water, which is mainly caused by turbulence flow and many uncertainty factors.
According to the definition of selfstability, a mathematical model can be established for qualitative analysis. If the ideal selfstability exists, the relation is held as follows:
$$ \chi<\chi_{\lim } $$ (20) Substituting the definition of χ and χ_{lim} into Eq. (20), we have the following:
$$ \arccos (\cos \varphi \cos \psi)<\arcsin \left(\frac{R(t)R}{L}\right) $$ (21) Simplifying Eq. (21), we attain the following equation:
$$ R(t)<L \sqrt{1\cos ^2 \varphi \cos ^2 \psi}+R $$ (22) Substituting the supercavity equation into Eq. (22), we attain the following:
$$ \sqrt{R_n^2 \frac{C_{x 0}(1+\sigma(t))}{k \sigma(t)}\left(R_n^2 \frac{C_{x 0}(1+\sigma(t))}{k \sigma(t)}R_n^2\right) \cdot\left(1\frac{2 L \sigma(t) \cos \varphi \cos \psi}{2 R_n \sqrt{C_{x 0}(1+\sigma(t)) \ln \sigma(t)}}\right)^{2 / \eta_c}} \\ <L \sqrt{1\cos ^2 \varphi \cos ^2 \psi}+R $$ (23) Meanwhile, the selfstability equation also holds:
$$ \begin{aligned} & m\left(\frac{\mathrm{d} V_{x 1}}{\mathrm{d} t}+\omega_{y 1} V_{z 1}\omega_{z 1} V_{y 1}\right)= \\ & \quad\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) \cos \alpha \cos \betam g \sin \psi \end{aligned} $$ (24) $$ \begin{aligned} & m\left(\frac{\mathrm{d} V_{y 1}}{\mathrm{d} t}+\omega_{y 1} V_{z 1}\right)= \\ & \quad\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) \sin \alpha \cos \betam g \cos \psi \end{aligned} $$ (25) $$ m\left(\frac{\mathrm{d} V_{z 1}}{\mathrm{d} t}\omega_{y 1} V_{x 1}\right)=\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) \sin \beta $$ (26) $$ I_{y 1} \frac{\mathrm{d} \omega_{y 1}}{\mathrm{d} t}=\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) x_{\mathrm{cm}} \sin \beta $$ (27) $$ I_{z 1} \frac{\mathrm{d} \omega_{z 1}}{\mathrm{d} t}=\frac{1}{2} C_0 A_c \rho V^2 \sin \left(\frac{\pi}{2}\delta\right) x_{\mathrm{cm}} \sin \alpha \cos \beta $$ (28) The angular velocity ω_{x1} = 0 rad/s and roll angle θ = 0 rad are used in Eqs. (24)–(28), which are constrained by the cavitation number σ(t), pitch angle ψ, and yaw angle φ. The relation between the parameters of Eq. (23) is extremely complicated and strictly defined in mathematics. Meanwhile, the numerical results indicated that the solution is always divergent under a small initial disturbance.
This conclusion was also obtained by the following optimization results. Therefore, designing a kinetic energy supercavitating projectile to keep an ideal selfstability when the projectile moves in water is impossible (the velocity is between 300 and 900 m/s). If the projectile achieves selfstability, then the control strategy must be employed in this projectile, which can cause hydrodynamic forces to constantly shift positions to balance the gravity and disturbances of the projectile.
3.3 Optimizing the projectile's structure
For the above features, the χ−χ_{lim} curve is a key basis for judging the stability of a supercavitating projectile. Figure 1 shows the χ−χ_{lim} curve based on the conditions in Table 1. The curve indicates that this projectile can achieve stable motion under this structure (Figure 2) and the initial conditions.
Number Physical quantity Value 1 Moment of inertia I_{x1} (kg·m^{2}) 5.82×10^{4} 2 Moment of inertia I_{y1} and I_{z1} (kg·m^{2}) 1.59×10^{6} 3 Supercavitating projectile's mass (kg) 0.116 4 Distance between the mass center and cavitator (m) 6.025×10^{2} 5 Projectile's maximum diameter (m) 1.256×10^{2} 6 Projectile's length (m) 1.568 1×10^{2} 7 Cavitator diameter (m) 4.5×10^{3} 8 Fluid density (kg/m^{3}) 1 000 9 Fluid temperature (℃) 20 10 Saturated vapor pressure (Pa) 2 338.8 11 Fluid kinematic viscosity (m^{2}/s) 1.006×10^{6} 12 The atmosphere pressure (Pa) 101 325 13 The initial depth of this projectile motion (m) 2 14 Acceleration of gravity (m/s^{2}) 9.8 15 V_{x1} (m/s) 500 16 V_{y1}, V_{z1} (m/s) 0 17 ω_{x1} (rad/s) 0 18 ω_{y1} (rad/s) 1 19 ω_{z1} (rad/s) 1 20 θ, φ, ψ (rad) 0 21 x, y, z (m) 0 However, the proportion of the χ above χ_{lim} gradually increased with the projectile motion. This condition Eq. (20) suggests the increasingly serious strengthening of the tail slaps. Tail slapping is an important factor in velocity attenuation. The calculation results showed that the displacement of this projectile was 25 m. Although this supercavitating projectile was stable, its sailing distance was short and could not meet the design requirements. Therefore, a reasonable structure design is imperative.
4 Optimization of the structure of supercavitating projectiles
The stability of supercavitating projectiles is affected by many factors, which include initial disturbances, motion disturbances, flow velocity, and fluid structure. The initial disturbances, motion disturbances, and flow velocity are external factors, and the fluid structure is an internal factor. From mathematical point of view, the projectile's configuration is the parameter of dynamic system, which determines the evolution direction of the nonlinear dynamic system. Therefore, the projectile's structure is significant for stability. The design objective meets the design requirements, which allows the projectile to achieve exceptional sailing distance and damage effect.
4.1 Simplification of the projectile structure and calculation
The supercavitating projectile's shape is complicated. Mechanical parameters, such as the moment of inertia, mass, and volume, are usually obtained by the software SolidWorks. The complex structure parameters will increase the number of optimal objective functions. The number of objective functions will greatly increase the difficulty of optimization and reduce the solving speed. Meanwhile, when the supercavity totally covers the projectile, the local structure, such as the fin, which is the transition part between the fin and cylinder, has a minimal effect on the projectile's motion. One of the most notable problems is the fins. For kinetic supercavitating projectiles, the tail is not designed to generate aerodynamic force inside the cavity but to gain an amphibious property. The gas inside the supercavity is thin, and the aerodynamic force will not affect the motion of the supercavitating projectile under a nonhypervelocity condition (V < 900 m/s).
Kinetic energy supercavitating projectiles feature small volume and light mass. The action mode of tail slaps is essentially different from that of large supercavitating projectiles.
The cavity diameter of the tail of the kinetic supercavitating projectile is relatively large and approximately 2–3 times larger than the projectile's tail diameter. In this paper, the cavity diameter was around 2.3 times larger and completely wrapped the projectile. When tail slaps occur, its tail is in unilateral wet. Given the high speed of kinetic energy supercavitating projectiles, secondary cavitation usually occurs at the fin. Secondary cavities on different fins affect each other and merge. Such an effect is equivalent to that of treating the tail as a cylinder. In addition, given the small size of the supercavitating projectile, the secondary cavitation effect of the fin has minimal influence on the cavity during tail slaps. Figure 3 shows the test results of the fin (Wang et al., 2018b). The fin height H was 24 mm. Figure 3 also reveals the test results of the fin piercing cavities with diameters of 18, 12, and 6 mm. The influence of fin on the cavity is very significant in the cases of a and b but very weak in the case of c. In this work, the height of the fin H was 3.28 mm, and the influence on the cavity can be ignored. Therefore, the influence of the fin on the tail slap force was ignored.
Most of the designs of kinetic energy supercavitating projectiles lack fins. Given that almost no aerodynamic stability mechanism is available for kinetic energy supercavitating projectiles when the flight speed is less than 1 000 m/s, the tails are nonfunctional underwter. The details can be found in Zhang et al. (2014).
However, large supercavitating projectiles are different. Large supercavitating projectiles have lower speeds than kinetic supercavitating projectiles and control and propulsion systems, and the cavity diameter of the projectile tail is usually 1.5–2 times that of the projectile tail. The deflection angle is generated by the fin's deflection, which can generate hydrodynamic force and control the attitude of the projectile. Therefore, the fins of supercavitating vehicles usually pierce the supercavity. Thus, the influence of the fins on supercavity cannot be ignored.
The largescale and small supercavitating vehicles mainly involves the following three differences:
1) From the launch point of view, the launch of supercavitating vehicles mainly depends on the launcher and launch system. Its diameter is close to that of a missile or torpedo. However, supercavitating projectiles rely on highdensity and highspeed smallcaliber artillery and machine guns to launch.
2) From the perspective of motion mode, supercavitating vehicles are usually equipped with a control system with a navigation speed of 50–200 m/s, a large sailing distance, and a long navigation time. However, the supercavitating projectile has a sailing speed between 200–1 500 m/s; its motion completely depends on its own kinetic energy, and its sailing distance and sailing time are small.
3) From the perspective of the generation mechanism of supercavity, the supercavity generation of supercavitating vehicles mainly occurs through ventilation. By changing the ventilation and thermodynamic state of the gas, the shape of the supercavity is controlled to change the navigation state of the vehicle. However, supercavitating projectiles rely entirely on their cavitators to produce a supercavity. The change in their velocity will modify their supercavity morphology and eventually cause their gradual evolution toward instability.
A simplified structure was used to derive the formulation of mass center, mass, and moment of inertia. The formulation's result (I_{x}_{1}=1.538×10^{−4} kg·m^{2}, I_{z}_{1}=I_{y}_{1}=1.436×10^{−6} kg·m^{2}, m=0.095 7 kg, x_{cm}=5.99×10^{−2} m) is extremely close with that obtained by the software (I_{x}_{1}=5.82×10^{−4} kg·m^{2}, I_{z}_{1}=I_{y}_{1}= 1.59×10^{−6} kg·m^{2}, m=0.116 kg, x_{cm}=6.025×10^{−2} m). The mass center error was less than 1%. However, the resulting I_{i}_{1} (I_{y}_{1}, I_{z}_{1}) was relatively large. This outcome was caused by the simplification of calculations. The moment of inertia was mainly used to calculate attitude angles. However, the attitude angles change was small for the projectile with stable motion. Meanwhile, the order of magnitude of I_{i}_{1} (I_{y}_{1}, I_{z}_{1}) was consistent with those obtained with formulations and software. Thus, the moment of inertia could not make large errors. For I_{x}_{1}, the fin stability method was employed in the amphibious supercavitating projectiles without a rotating stability. Therefore, I_{x}_{1} had a slight influence on this structure optimization problem. In summary, we ignored the local structures to accelerate the calculation and improve the efficiency of design. Figures 4 and 5 show the original and simplified structures, respectively.
The simplified structure has four parameters: the length of cylinder part L_{1}, diameter of cylinder part D, diameter of cavitator d, and length of frustum of cone L_{2}.
For optimal calculation, the formula of the mass center and moment of inertia must be derived based on the simplified model. We set the coordinate origin at the cavitator center. Figure 6 shows the projectile's axis, that is, the z axis. The material of frustum of cone L_{2} was tungsten alloy (density: 18.75 g/cm^{3}), and that of cylinder part L_{1} was aluminum (density: 2.7 g/cm^{3}). According to the coordinate of the cavitator, the mass center position of x and y were 0 and 0, respectively. The value of z can be calculated as follows:
$$ z_0=\frac{\iiint\limits_{v_1+v_2} z \rho(x, y, z) \mathrm{d} v}{\iiint\limits_{v_1+v_2} \rho(x, y, z) \mathrm{d} v} $$ (29) Eq. (29) can be calculated in two steps as follows:
Step 1: Numerator:
$$ \iiint\limits_{v_1+v_2} z \rho(x, y, z) \mathrm{d} v \\ =\iiint\limits_{v_1} z \rho(x, y, z) \mathrm{d} v+\rho_2 \iiint\limits_{v_2} z \rho(x, y, z) \mathrm{d} v \\ =\rho_1 \iiint\limits_{v_1} z \mathrm{d} v+\rho_2 \iiint\limits_{v_2} z \mathrm{d} v \\ \begin{aligned} = & \rho_1 \int_{L_2}^{L_1+L_2} z \mathrm{d} z \int_0^{2 \pi} \mathrm{d} \theta \int_0^{\frac{D}{2}} r \mathrm{d} r+ \\ & \rho_2 \int_0^{L_2} z \mathrm{d} z \int_0^{2 \pi} \mathrm{d} \theta \int_0^{\frac{z}{2 L_2}[Dd]+\frac{d}{2}} r \mathrm{d} r \\ = & \frac{L_1 D^2 \pi \rho_1}{8}\left(L_1+2 L_2\right)+\frac{L_2^2 \pi \rho_2}{48}\left(d^2+2 d D+3 D^2\right) \end{aligned} $$ (30) Step 2: Denominator:
$$ \iiint\limits_{v_1+v_2} \rho(x, y, z) \mathrm{d} v=\iiint\limits_{v_1} \rho(x, y, z) \mathrm{d} v+\iiint\limits_{v_2} \rho(x, y, z) \mathrm{d} v \\ =\rho_1 \iiint\limits_{v_1} \mathrm{d} v+\rho_2 \iiint\limits_{v_2} \mathrm{d} v \\ =\frac{1}{4} \rho_1 \pi D^2 L_1+\frac{1}{12} \rho_2 \pi L_2\left(D^2+d^2+D d\right) $$ (31) $$ z_0=\frac{\frac{L_1 D^2 \pi \rho_1}{8}\left(L_1+2 L_2\right)+\frac{L_2^2 \pi \rho_2}{48}\left(d^2+2 d D+3 D^2\right)}{\frac{1}{4} \rho_1 \pi D^2 L_1+\frac{1}{12} \rho_2 \pi L_2\left(D^2+d^2+D d\right)} $$ (32) When the mass center was obtained by the above method, the coordinate origin moved to the mass center. The coordinate is called the principal axis system. The matrix of the moment of inertia degenerated to a diagonal matrix. According to the symmetry, the formulation I_{y} = I_{x} holds. The matrix of the moment of inertia was different due to the mass center position and needed to be divided into two cases for the analysis:
When z_{0} > L_{2} holds, the mass center is located on the cylinder part, and the moment of inertia I_{x} can be calculated as follows (Figure 7):
$$ \begin{aligned} I_x & =\iiint\limits_{v_1+v_2}\left(y^2+z^2\right) \rho(x, y, z) \mathrm{d} v=\iiint\limits_{v_1}\left(y^2+z^2\right) \rho(x, y, z) \mathrm{d} v \\ & +\iiint\limits_{v_2}\left(y^2+z^2\right) \rho(x, y, z) \mathrm{d} v \\ & =\rho_1 \iiint\limits_{v_1}\left(y^2+z^2\right) \mathrm{d} v+\rho_2 \iiint\limits_{v_2}\left(y^2+z^2\right) \mathrm{d} v=\rho_1 \iiint\limits_{v_1} y^2 \mathrm{d} v \\ & +\rho_1 \iiint\limits_{v_1} z^2 \mathrm{d} v+\rho_2 \iiint\limits_{v_2} y^2 \mathrm{d} v+\rho_2 \iiint\limits_{v_2} z^2 \mathrm{d} v \end{aligned} \\ \begin{aligned} & =\rho_1 \int_{\left(z_0L_2\right)}^{L_1\left(z_0L_2\right)} \mathrm{d} z \int_0^{\frac{D}{2}} r^3 \mathrm{d} r \int_0^{2 \pi} \sin ^2 \theta \mathrm{d} \theta \\ & +\rho_1 \int_{\left(z_0L_2\right)}^{L_1\left(z_0L_2\right)} z^2 \mathrm{d} z \int_0^{\frac{D}{2}} r \mathrm{d} r \int_0^{2 \pi} \mathrm{d} \theta \\ & +\rho_2 \int_{z_0}^{\left(z_0L_2\right)} \mathrm{d} z \int_0^{\frac{1}{L_2}\left[L_2\left[z\left(z_0L_2\right)\right]\right]\left(\frac{D}{2}\frac{d}{2}\right)+\frac{d}{2}} r^3 \mathrm{d} r \int_0^{2 \pi} \sin ^2 \theta \mathrm{d} \theta \\ & +\rho_2 \int_{z_0}^{\left(z_0L_2\right)} z^2 \mathrm{d} z \int_0^{\frac{1}{L_2}\left[L_2\left[z\left(z_0L_2\right)\right]\right]\left(\frac{D}{2}\frac{d}{2}\right)+\frac{d}{2}} r \mathrm{d} r \int_0^{2 \pi} \mathrm{d} \theta \end{aligned} \\ \begin{aligned} & =\frac{D^4}{64} L_1 \pi \rho_1+\frac{\rho_1}{12} D^2 L_1 \pi\left[L_1^2+3 L_1\left(L_2z_0\right)+3\left(L_2z_0\right)^2\right] \\ & +\frac{\pi \rho_2}{320} L_2\left(d^4+d^3 D+d^2 D^2+d D^3+D^4\right) \\ & +\frac{\pi \rho_2}{120} L_2\left[L_2^2\left(d^2+3 d D+6 D^2\right)\right. \\ & \left.5\left(d^2+2 d D+3 D^2\right) L_2 z_0+10 z_0^2\left(d^2+d D+D^2\right)\right] \end{aligned} $$ (33) When z_{0} < L_{2} holds, the mass center is located on the frustum of the cone part, and the moment of inertia I_{x} can be calculated as follows (Figure 8):
$$ \begin{aligned} I_x & =\iiint\limits_{v_1+v_2}\left(y^2+z^2\right) \rho(x, y, z) \mathrm{d} v=\iiint\limits_{v_1}\left(y^2+z^2\right) \rho(x, y, z) \mathrm{d} v \\ & +\iiint\limits_{v_2}\left(y^2+z^2\right) \rho(x, y, z) \mathrm{d} v \\ & =\rho_1 \iiint\limits_{v_1}\left(y^2+z^2\right) \mathrm{d} v+\rho_2 \iiint\limits_{v_2}\left(y^2+z^2\right) \mathrm{d} v \\ & =\rho_1 \iiint\limits_{v_1} y^2 \mathrm{d} v+\rho_1 \iiint\limits_{v_1} z^2 \mathrm{d} v+\rho_2 \iiint\limits_{v_2} y^2 \mathrm{d} v+\rho_2 \iiint\limits_{v_2} z^2 \mathrm{d} v \end{aligned} \\ \begin{aligned} & =\rho_1 \int_{L_2z_0}^{L_2z_0+L_1} \mathrm{d} z \int_0^{2 \pi} \sin ^2 \theta \mathrm{d} \theta \int_0^{\frac{D}{2}} r^3 \mathrm{d} r \\ & +\rho_1 \int_{L_2z_0}^{L_2z_0+L_1} z^2 \mathrm{d} z \int_0^{2 \pi} \mathrm{d} \theta \int_0^{\frac{D}{2}} r \mathrm{d} r \\ & +\int_{z_0}^{L_2z_0} z^2 \mathrm{d} z \int_0^{2 \pi} \mathrm{d} \theta \int_0^{\frac{d}{2}+\frac{(Dd)}{2 L_2}} r \mathrm{d} r \\ & +\int_{z_0}^{L_2z_0} \mathrm{d} z \int_0^{2 \pi} \sin ^2 \theta \mathrm{d} \theta \int_0^{\frac{d}{2}+\frac{(Dd)}{2 L_2}} r^3 \mathrm{d} r=\frac{D^4 \rho_1 L_1 \pi}{64} \end{aligned} \\ \begin{aligned} & +\frac{D^2 \rho_1 L_1 \pi\left[L_1^2+3 L_1\left(L_2z_0\right)+3\left(L_2z_0\right)^2\right]}{12} \\ & +\frac{\pi\left[d\left(L_2z_0\right)+D z_0\right]^2\left[L_2^23 L_2 z_0+3 z_0^2\right]}{12 L_2} \\ & +\frac{\pi\left[d\left(L_2z_0\right)+D z_0\right]^4}{64 L_2^3} \end{aligned} $$ (34) According to the parallel axis theorem, the coordinates of the mass center are the principal axes, and therefore, the moment of inertia I_{z} is constant. The formulation is as follows:
$$ I_z=\iiint\limits_{v_1+v_2}\left(x^2+y^2\right) \rho(x, y, z) \mathrm{d} v=\iiint\limits_{v_1}\left(x^2+y^2\right) \rho(x, y, z) \mathrm{d} v \\ +\iiint\limits_{v_2}\left(x^2+y^2\right) \rho(x, y, z) \mathrm{d} v \\ =\rho_1 \iiint\limits_{v_1}\left(x^2+y^2\right) \mathrm{d} v+\rho_2 \iiint\limits_{v_2}\left(x^2+y^2\right) \mathrm{d} v \\ \begin{aligned} & =\rho_1 \int_{L_2}^{L_1+L_2} \mathrm{d} z \int_0^{2 \pi} \mathrm{d} \theta \int_0^{\frac{D}{2}} r^3 \mathrm{d} r+\rho_2 \int_0^{L_2} \mathrm{d} z \int_0^{2 \pi} \mathrm{d} \theta \int_0^{\frac{z}{2 L_2}[Dd]+\frac{d}{2}} r^3 \mathrm{d} r \\ & =\frac{D^4 L_1 \pi \rho_1}{32}+\frac{\rho_2 \pi L_2\left(d^4+d^3 D+d^2 D^2+d D^3+D^4\right)}{160} \end{aligned} $$ (35) 4.2 Optimization of algorithm design
Holland (1992) researched the GA theory and method. GA is a highly parallel, random, and selfadaptive algorithm developed from the mechanism of natural selection and evolution in biology. In short, it uses the population to search for optimal solutions, and the population represents a group of solutions. By applying a series of genetic operations, such as selection, crossover, and mutation, to the current population, a new generation is generated, and the population gradually evolves to contain an approximate optimal solution. This optimal solution is not accurate in mathematics, but the approximate solution meets the basic design requirements (Vaissier et al., 2019; Chan et al., 2020; CostaCarrapiço et al., 2020; RamosFigueroa et al., 2020; Wang and Sobey, 2020).
Therefore, in this section, the approximate optimal solution was found based on the GA (Vaissier et al., 2019), and the algorithm was improved to optimize the offline projectile. The improved algorithm can be applied well to the optimization of projectile parameters, and it can accelerate the calculation speed.
For the optimization of this design problem, an optimization problem with a clear mathematical equation must be constructed first. Thus, the time series {t_{n}}, which is the time corresponding to the intersection point between χ and χ_{lim}, was introduced. This time series is important for stability analysis. If the time series {t_{n}} is long, the frequency of tail slap is high. A high tail slap frequency is disadvantageous to the projectile's stability and can increase the drag of flight. During the design of a supercavitating projectile, the number of tail slaps should be minimized as much as possible. The maximum of the time series {t_{n}} reflects the stable flight distance of the supercavitating projectile. When this projectile flight is 0.1 s, the maximum {t_{n}} is close to 0.1 s, which indicates that the projectile can fly stably for a long time. This condition can be written as follows:
$$ \max \left\{t_n\right\} \geqslant k_p t_p $$ (36) where k_{p} is the scale factor between 0 and 1 and t_{p} is the time of this projectile's movement.
If the projectile has minimized tail slap number and long sailing distance, its structural design is optimized. Therefore, {t_{n}} can be one of the constraint conditions, and the sailing distance can be the objective function of this optimization issue. The projectile's sailing distance and the number of {t_{n}} restrict each other. The most important concern in the design is the sailing distance. If the design sailing distance can be reached, then the number of {t_{n}} is not important for engineering.
In summary, this optimization problem can be written as follows:
$$ \begin{array}{ll} \max & x \\ \text { s.t. } & \max \left\{t_n\right\} \geqslant k_p t_p \end{array} \\ \begin{aligned} & m\left(\frac{\mathrm{d} V_{x 1}}{\mathrm{d} t}+\omega_{y 1} V_{z 1}\omega_{z 1} V_{y 1}\right)=F_{x 1} \\ & m\left(\frac{\mathrm{d} V_{y 1}}{\mathrm{d} t}+\omega_{z 1} V_{x 1}\omega_{x 1} V_{z 1}\right)=F_{y 1} \\ & m\left(\frac{\mathrm{d} V_{z 1}}{\mathrm{d} t}+\omega_{x 1} V_{y 1}\omega_{y 1} V_{x 1}\right)=F_{z 1} \\ & I_{x 1} \frac{\mathrm{d} \omega_{x 1}}{\mathrm{d} t}+\omega_{y 1} \omega_{z 1}\left(I_{z 1}I_{y 1}\right)=M_{x 1} \\ & I_{y 1} \frac{\mathrm{d} \omega_{y 1}}{\mathrm{d} t}+\omega_{x 1} \omega_{z 1}\left(I_{x 1}I_{z 1}\right)=M_{y 1} \\ & I_{z 1} \frac{\mathrm{d} \omega_{z 1}}{\mathrm{d} t}+\omega_{x 1} \omega_{y 1}\left(I_{y 1}I_{x 1}\right)=M_{z 1} \end{aligned} \\ \begin{aligned} \omega_{x 1} & =\frac{\mathrm{d} \theta}{\mathrm{d} t}+\frac{\mathrm{d} \varphi}{\mathrm{d} t} \sin \psi \\ \omega_{y 1} & =\frac{\mathrm{d} \varphi}{\mathrm{d} t} \cos \psi \cos \theta+\frac{\mathrm{d} \psi}{\mathrm{d} t} \sin \theta \\ \omega_{z 1} & =\frac{\mathrm{d} \psi}{\mathrm{d} t} \cos \theta\frac{\mathrm{d} \varphi}{\mathrm{d} t} \cos \psi \sin \theta \\ \frac{\mathrm{d} x}{\mathrm{d} t} & =V_x \\ \frac{\mathrm{d} y}{\mathrm{d} t} & =V_y \\ \frac{\mathrm{d} z}{\mathrm{d} t} & =V_z \end{aligned} $$ (37) 4.2.1 Note for initialization of the population
Traditionally, the initial individuals of the population, which are single variables, are generated randomly in a certain range. Four variables were considered parameters of the projectile structure in this work. These parameters feature synergies and control the motion of the projectile. The individuals of the population are the structural parameters and have four genes (cylinder length L_{1}, diameter of cylinder part D, diameter of cavitator d, and length of frustum of cone L_{2}). To ensure that good genes can be inherited between individuals during genetic operation, we assessed the stability of all generated individuals. The stability criterion was based on the maximum value of {t_{n}}. Supposing that the projectile's moving time is t_{p}, the maximum of {t_{n}} is {t_{n}}_{max}, and k_{p} is a coefficient between 0 and 1. If k_{p}t_{p} < {t_{n}}_{max} holds, then the projectile's motion is stable, and vice versa. When an individual is unstable, it needs to be regenerated and judged again until the initial population becomes stable and is outputted for use.
4.2.2 Note for fitness of the population
A range is used for the fitness evaluation of the population, and the fitness of the initial population is stored. In the first iteration, fitness can be used directly to avoid repeated calculations. For engineering, all individuals do not need to reach convergence completely; that is, the difference between individuals is less than 1%. Such accuracy is not significant for practice. Thus, after its calculation, the fitness must be evaluated to determine whether it can meet the engineering requirements. The fitness was evaluated by its mean value, and the difference vector between fitness and its mean value was solved. If the difference is acceptable, the final optimal results can be outputted.
4.2.3 Note for the initial condition of 6DOF dynamic equations
The initial value is the key to trajectory evolution. Therefore, the influence of the initial value of the 6DOF dynamic equations needs to be considered in the optimization problem to ensure that the optimized projectile can meet the requirements of trajectory dispersion, stability, and sailing distance under the small disturbance of the initial value. The projectile needs a certain generalization capability toward the initial value. Intermediate ballistics is a subject that studies the disturbance law of muzzle flow field to projectiles and mainly focuses on the disturbance phenomenon in the air (Li et al., 2015). In addition, the bore damage of the barrel has an important impact on its exterior ballistics (Shen et al., 2020a, 2020b). These two methods have become important in the study of the initial random disturbance of modern exterior ballistics.
These methods are also based on emission dynamics and the CFD method, but their calculation is large and complex, which is not conducive to engineering applications. Therefore, Monte Carlo simulation is still the main research method in modern external ballistics. In this work, ψ and φ were between − 0.05 and 0.05 rad. The ω_{y} and ω_{z} were between −1 and 1 rad/s.
The projectile was launched by a smoothbore gun. Thus, ω_{x} = 0 rad/s and θ = 0 rad.
4.3 Basic flow of the algorithm
The traditional GA was improved to accelerate the calculation speed. Figure 9 shows the flow of this algorithm. The main structure of this algorithm is as follows:
5 Dynamic features of the optimized projectile
5.1 Optimized results
For the maximum use of the initial kinetic energy of the projectile, the maximum sailing distance and stability should be achieved under the initial velocity. The maximum length of the projectile was not more than 200 mm, and the minimum length was not less than 100 mm. Given that the long rod projectile is beneficial for damaging the target, the slenderness ratio of the long rod projectile is generally greater than 10. The projectile's body consists of two materials: the frustum of the cone is tungsten alloy (density: 18.75 g/cm^{3}), and the cylinder part is aluminum (density: 2.7 g/cm^{3}). The lengths of the frustum of the cone and the cylinder were between 50 and 100 mm. A 30 mmcaliber gun was used to launch the projectile. Considering the sabot design, the cylinder diameter cannot be adjusted in a large range. Therefore, the cylinder diameter was set in the range of [12.5, 12.56] mm, and that of the cavitator was [1, 12.5] mm. To ensure accuracy, we used a 22bit binary code for coding. The initial population was 200, and 500 generations were set in the evolution. The cross probability was 0.7, the selection probability was 0.5, and the mutation probability was 0.001 for a general GA.
Algorithmfor optimized supercavitating projectiles 1 Step1: Start; 2 Step2: Calculation of the initial population. 3 Sub2Step1: In a certain range, a certain number of parameters are randomly generated, namely, cylinder length L_{1}, diameter of cylinder part D, diameter of cavitator d, and length of frustum of cone L_{2}, and they are assembled as different individuals. 4 Sub2Step2: The stability is calculated for all individuals. If the result is unstable, turn to Sub2Step3; otherwise, turn to Sub2Step4. 5 Sub2Step3: A new individual is regenerated, and the stability is also calculated; if it is unstable, turn to Sub2Step3; otherwise, turn to Sub2Step4. 6 Sub2Step4: Output the individuals and assemble them into a population. 7 Step3: Start the iteration. 8 Step4: Calculate the initial fitness of the population. 9 Sub4Step1: If it is the initial population, turn to Sub4Step2; otherwise, turn to Sub4Step3. 10 Sub4Step2: The range of the initial population is assigned to the fitness, and Sub4Step4 is applied. 11 Sub4Step3: The individual fitness is calculated, and Sub4Step4 is applied. 12 Sub4Step4: The fitness is outputted. 13 Step5: The mean value of fitness is calculated. 14 Step6: The difference between the fitness vector and its mean value is calculated. 15 Step7: If the difference meets the design, turn to Step8; otherwise, turn to Step10. 16 Step8: The maximum of this difference is solved. 17 Step9: The optimal results are outputted. 18 Step10: Selection. 19 Step11: Coding. 20 Step12: Crossover. 21 Step13: Mutating. 22 Step14: Decoding. 23 Step15: The population is reorganized, and Step3 is applied. Table 2 shows the optimized projectile parameters. It meets the convergence requirements when it evolves to the 49^{th} generation. The projectile can have a large displacement. For the original structure, when the projectile moved for 0.1 s, its velocity was attenuated to 200 m/s. It reached the termination of this code because the code set 0.1 s as the motion time. However, the improved structure was different. When the projectile flight lasted for 0.1 s, the rest velocity was 437.84 m/s. The change in velocity was 62.16 m/s after 0.1 s. Currently, the rest velocity is considerable, and it not only has a strong penetration capability but can also increase the sailing distance. Therefore, the optimized projectile has a better water ballistic performance. The water ballistic performance is mainly evaluated by displacement. For a certain projectile, if the velocity decreases slowly, then the number of tail slaps is lessened, and a good tendency of selfstability is observed. Thus, the projectile is well designed. This projectile structure is instructive for the design of underwater projectiles with a large sailing distance. Therefore, the nature of stability, which is meaningful for engineering design, must be researched in depth. Figure 10 shows the optimized supercavitating projectile configuration and the conceptual design of amphibious supercavitating projectile. This optimized projectile features a large cone part, which has been found in Yi et al. (2008a, 2008b, 2008c, 2009a, 2009b). Meanwhile, the accuracy of this optimization procedure has also been verified. We provided the supercavitating projectile with a tail to inspire more scholars to further carry out the overall design research of amphibious supercavitating projectiles. The amphibious supercavitating projectile will become an important development direction of underwater ammunition in the future.
Items L_{2} (mm) L_{1} (mm) d (mm) D (mm) Max sailing distance (m) Origin 52.92 103.89 4.50 12.56 25.00 Optimized 97.91 96.62 1.00 12.50 46.75 5.2 Mechanism of stability
Table 3 lists the calculation parameters. However, the matrix of the moment of inertia, mass center position, mass, and projectile length can be obtained by calculating the structural parameters. Thus, this optimized projectile is different from the original structure.
Physical quantity Value The density of the frustum of the cone (kg/m^{3}) 18.75×10^{3} The density of the cylinder part (kg/m^{3}) 2.7×10^{3} The length of frustum of cone (m) 9.79×10^{−2} The length of the cylinder part (m) 9.66×10^{−2} The projectile's maximum diameter (m) 1.25×10^{−2} The cavitator diameter (m) 1.0×10^{−3} The fluid density (kg/m^{3}) 1 000 The fluid temperature (℃) 20 The saturated vapor pressure (Pa) 2 338.8 The fluid kinematic viscosity (m^{2}/s) 1.006×10^{−6} The atmosphere pressure (Pa) 101 325 The initial depth of this projectile motion (m) 2 Acceleration of gravity (m/s^{2}) 9.8 V_{x}_{1} (m/s) 500 V_{y}_{1}, V_{z}_{1} (m/s) 0 ω_{x}_{1} (rad/s) 0 ω_{y}_{1}, ω_{z}_{1} (rad/s) 1 θ, φ, ψ (rad) 0 x, y, z (m) 0 By the initial conditions in Table 2, the trajectory and attitude of the projectile can be calculated. First, for the kinetic energy supercavitating projectile, the most important thing is its stability. The good performance of the supercavitating projectile is based on a wellstructured design. It is instructive for the projectile's design in future works. Therefore, the first task is to deeply explore the stability mechanism of this projectile. In our previous work, we proposed the χ−χ_{lim} curve to judge the stability. Figure 11 shows the χ−χ_{lim} curve of the optimized projectile. According to the stability criterion, the number of tail slaps was 8, which is beneficial for the projectile's stability. Thus, the projectile moved by 46.75 m. Meanwhile, when the χ−χ_{lim} curve was enlarged partly (Figure 12), the difference Δχ between χ and χ_{lim} curves was Δχ << 4×10^{−}^{4} rad=0.017°. This result indicates that the supercavitating projectile belongs has a bodyfitted cavity, which is close to the effect of film drag reduction theoretically. The findings also showed that when the tail slaps were about to occur, given its good selfstability, the tail slap motion did not proceed completely. Only the tail of the projectile touched the cavity, and the projectile showed the dynamic behavior of reverse rotation immediately. The tail slap action was largely weakened, and thus, the drag of this projectile can be reduced greatly. The results also revealed that the mathematical solution may be inconsistent with the dynamic behavior. Therefore, the stability criterion, which was proposed in our previous work, needs to be weakened. This relation can be written as follows:
$$ \left\{\begin{array}{l} \left.\frac{\mathrm{d} \chi}{\mathrm{d} t}\right_{t=t_i}=0 \\ \left.\frac{\mathrm{d}^2 \chi}{\mathrm{d} t^2}\right_{t=t_i}<0 \end{array} \Rightarrow \chi_{t_i}\chi_{\lim _i} \leqslant \varepsilon=0.01745\;\mathrm{rad}\right. $$ (38) where {t_{i}} is maximum time array of χ, {χ_{ti}} and {χ_{lim}_{ti}} are the function value of χ and χ_{lim} curves at {t_{i}}, respectively. ε is the weakening quantity. When all values of the sequence {χ_{ti} − χ_{limti}} are less than the ε, then the motion is weakened stability. The value of ε was set as 1°, that is, 0.017 45 rad.
The discussion above showed that selfstability only exists in theory. Notably, given the uncertainty of launch conditions and the environment, selfstability is difficult to maintain. The kinetic projectile approaches selfstability in a state of minimum energy consumption under current conditions, which is a universal law of nature, that is, the principle of minimum energy. For the actual needs, the weakened selfstability is defined as follows: the motion closest to selfstability.
In Section 3.2, we proposed the ideal motion equations of this projectile based on the stability criterion. Selfstability is a special motion that results in the complete projectile movement in the supercavity. In this condition, the best drag reduction effect can be achieved. This motion is impossible in theory because the gravity component in the vertical direction is hard to balance. Eqs. (24) – (28) represent an ideal case. However, in case of any disturbances, this ideal motion will be lost and transferred to the tail slaps.
Figure 11 shows a motion with a good selfstability tendency. This motion is an optimal one under the design. The experimental results also indicated that when the projectile's velocity was more than 1 000 m/s, that is, hypervelocity, the projectile approximately moved to the supercavity. However, the selfstability tendency, whose velocity is between 300 – 900 m/s, is different from the motion of hypervelocity. For the hypervelocity motion, the projectile will interact with the steam medium produced by the cavitation. The selfstability tendency is a projectile's inherent property. We have observed this phenomenon in experiments (Figure 13). The tail of the projectile first touched the smooth and stable supercavity, and then, the tail immediately left the wall and rotated in reverse. At this point, an instantaneous disturbance occurred on the supercavity, and a disturbance wave was generated on the supercavity and expanded continuously. The supercavitating projectile constantly moved in this mode within 0.1 s. Thus, the drag of the projectile was mainly the pressure difference drag, and the friction drag can be ignored during the tail slaps.
Figure 13 displays the velocity curve of the projectile. A comparison of the results with the structure discussed in previous work (Wang et al., 2020) revealed the following characteristics:
1) Similarities: The mean value of kinetic energy dissipation in V_{y}_{1} and V_{z}_{1} directions was zero.
2) Differences: The velocity amplitude of the optimized structure in the direction of V_{y}_{1} and V_{z}_{1} was far less than that before optimization. This result suggests that the selfstability tendency of this projectile has less energy dissipation in V_{y}_{1} and V_{z}_{1} directions. The projectile's velocity attenuated slowly per second. Figure 14 also reveals the resultant velocity and velocity components V_{x}_{1} in x_{1}. The two curves are very close, which indicates that the attenuation of the projectile in the axial was consistent with the resultant velocity. The dissipation of kinetic energy was mainly used to overcome the pressure difference drag movement and increase the sailing distance.
Figure 15 shows the angular velocity curve. Compared with previous work, the angular velocity curve changed dramatically because of the projectile's good selfstability. When the projectile moved into the supercavity, the hydrodynamics acting on the cavitator consistently produced a restoring moment. The moment caused the projectile to rotate inversely. For the maintenance of stability, the hydrodynamics changed directions continually. This performance is related to the structural design. Moreover, although the angular velocity derivative was relatively large, the angular velocity was smaller than that of the projectile with poor stability. In addition, the angular velocity fluctuation range was between − 2.5 and 2.5 rad/s and stable during movement, which is also an important reason for projectile stable motion. Affected by the angular velocity, the pitch angle ψ, yaw angle φ, attack angle α, and sideslip angle β showed changes in a small range. The curve of the pitch angle ψ and yaw angle φ indicated that the behavior of this projectile was stable in water. The attack angle α, and sideslip angle β suggest that the body is less likely to overturn (Figures 16 and 17, respectively). These indicators all show the welldesigned results of this projectile.
Figure 18 presents the displacement of the projectile. When the projectile moved by approximately 47 m in the x direction, the horizontal and vertical dispersions were 4.17 and less than 1 cm, respectively. Therefore, the selfstability tendency of projectile ballistic dispersion was also small. This feature is beneficial for the kinetic energy supercavitating projectile launched by a set number design. The projectile can impact the fixed target accurately in a short time for rapid damage in water.
To fully show the advantages of the optimized projectile, we further provided the simulation results of the optimized projectile within 0.5 s here (Figure 19). Figure 19 shows that the optimized projectile had a sailing distance of nearly 200 m in 0.5 s and a residual speed of 291.4 m/s. Therefore, the optimized kinetic supercavitating projectile showed good performance in terms of the exterior and terminal ballistics.
6 Conclusions
For kinetic energy supercavitating projectiles, stability and sailing distance are key factors in evaluating the projectile performance. These indicators are highly dependent on the structure. Therefore, the structure must have a reasonable design. To improve the sailing distance and stability of the supercavitating projectile, we optimized the projectile structure based on GA. The main conclusions are as follows:
1) Complete selfstability is difficult to achieve. Weakening selfstability is the closest motion state to selfstability in the actual process of supercavitating projectiles. Affected by the initial launch conditions, launching the projectile without any disturbance is difficult. To enable the optimized projectile to reach a certain generalization capability for the initial value, we considered the influence of the initial value by Monte Carlo simulation in the optimization process. Therefore, the optimized projectile approached selfstability by minimum energy dissipation, that is, weakened selfstability. The weakened selfstability is the reflection of selfstability in practice.
2) The kinetic supercavitating projectile has the tendency to destabilize. Whether optimized or unoptimized, the energy of the kinetic supercavitating projectile is always gradually consumed and reduced in underwater motion. The reduction of kinetic energy weakens the stability of the projectile and eventually loses stability, which is the inevitable trend for kinetic projectile movement. However, the purpose of optimizing the projectile configuration is to delay the evolution time of the projectile to instability. The instability cannot be avoided unless the control method is employed.
The optimization method and parametric expression of the projectile used in this paper are universal for the kinetic energy supercavitating projectile, which is launched directly underwater. However, with the development of combat integration, amphibious ammunition will become an important development trend in the future. The research of amphibious ammunition needs to comprehensively consider airexterior, crossmedium, and water ballistics. Therefore, the optimization strategy and method for amphibious ammunition have become an important direction of future research work.
Nomenclatures
m Mass of supercavitating projectile (kg) V_{x1} Velocity component on body coordinate system (O_{x1}) (m/s) V_{y1} Velocity component on body coordinate system (O_{y1}) (m/s) V_{z1} Velocity component on body coordinate system (O_{z1}) (m/s) V Resultant velocity of supercavitating projectile (m/s) ω_{x1} Angular velocity component on body coordinate system (O_{x1}) (rad/s) ω_{y1} Angular velocity component on body coordinate system (O_{y1}) (rad/s) ω_{z1} Angular velocity component on body coordinate system (O_{z1}) (rad/s) t Time (s) C_{0} Drag coefficient A_{c} Cavitator area (m^{2}) ρ Water density (kg/m^{3}) δ Angle between the velocity vector and O_{x1} (rad) α Attack angle (rad) α Sideslip angle (rad) C_{xf} Coefficient of friction drag S_{w} Surface area of the projectile wettingwater (m^{2}) k_{n} Smooth function coefficient h_{ef} Critical parameter (m) g Acceleration of gravity (m/s^{2}) θ Roll angle (rad) φ Yaw angle (rad) ψ Pitch angle (rad) R(t) Supercavity function (m) R Radius of projectile (m) η Angle of the tailslap force direction (rad) χ Angle between the body and the line through the mass center (rad) φ_{lim} Angle between the supercavity and body (rad) R_{n} Cavitator radius (m) C_{xfr} Coefficient of rotation friction drag I_{x1} Moment of inertia on O_{x1} (kg·m^{2}) I_{y1} Moment of inertia on O_{y1} (kg·m^{2}) I_{z1} Moment of inertia on O_{z1} (kg·m^{2}) x_{cm} Length of mass center to the cavitator (m) σ(t) Cavitation number C_{x0} Drag coefficient (0.82) x_{1} Projectile's length projection on the horizon (L_{p}) (m) η_{c} Cavity coefficient (0.85) L Projectile length (m) L_{1} Length of the cylinder part (m) D Diameter of the cylinder part (m) d Cavitator diameter (m) L_{2} Length of frustum of cone (m) z_{0} Mass center in calculation coordinate (m) I_{x}, I_{y}, I_{z} Moment of inertia in calculation coordinate (kg·m^{2}) v_{1} Volume of the cylinder part (m^{3}) v_{2} Frustum of cone volume (m^{3}) P_{a} Standard atmosphere pressure (Pa) P_{c} Supercavity pressure (Pa) R_{c}(t) Maximum radius of supercavity (m) L_{c}(t) Supercavity length (m) H The fin's length (m) Competing interest The authors have no competing interests to declare that are relevant to the content of this article. 
Figure 1 χ−χ_{lim} curve (Wang et al., 2020)
Figure 2 Structure of supercavitating projectile (Wang et al., 2018a)
Figure 3 Experiment results (Wang, 2018b)
Figure 4 Original structure (Wang et al., 2019)
Table 1 Simulation parameters of supercavitating projectile (Wang et al., 2020)
Number Physical quantity Value 1 Moment of inertia I_{x1} (kg·m^{2}) 5.82×10^{4} 2 Moment of inertia I_{y1} and I_{z1} (kg·m^{2}) 1.59×10^{6} 3 Supercavitating projectile's mass (kg) 0.116 4 Distance between the mass center and cavitator (m) 6.025×10^{2} 5 Projectile's maximum diameter (m) 1.256×10^{2} 6 Projectile's length (m) 1.568 1×10^{2} 7 Cavitator diameter (m) 4.5×10^{3} 8 Fluid density (kg/m^{3}) 1 000 9 Fluid temperature (℃) 20 10 Saturated vapor pressure (Pa) 2 338.8 11 Fluid kinematic viscosity (m^{2}/s) 1.006×10^{6} 12 The atmosphere pressure (Pa) 101 325 13 The initial depth of this projectile motion (m) 2 14 Acceleration of gravity (m/s^{2}) 9.8 15 V_{x1} (m/s) 500 16 V_{y1}, V_{z1} (m/s) 0 17 ω_{x1} (rad/s) 0 18 ω_{y1} (rad/s) 1 19 ω_{z1} (rad/s) 1 20 θ, φ, ψ (rad) 0 21 x, y, z (m) 0 Algorithmfor optimized supercavitating projectiles 1 Step1: Start; 2 Step2: Calculation of the initial population. 3 Sub2Step1: In a certain range, a certain number of parameters are randomly generated, namely, cylinder length L_{1}, diameter of cylinder part D, diameter of cavitator d, and length of frustum of cone L_{2}, and they are assembled as different individuals. 4 Sub2Step2: The stability is calculated for all individuals. If the result is unstable, turn to Sub2Step3; otherwise, turn to Sub2Step4. 5 Sub2Step3: A new individual is regenerated, and the stability is also calculated; if it is unstable, turn to Sub2Step3; otherwise, turn to Sub2Step4. 6 Sub2Step4: Output the individuals and assemble them into a population. 7 Step3: Start the iteration. 8 Step4: Calculate the initial fitness of the population. 9 Sub4Step1: If it is the initial population, turn to Sub4Step2; otherwise, turn to Sub4Step3. 10 Sub4Step2: The range of the initial population is assigned to the fitness, and Sub4Step4 is applied. 11 Sub4Step3: The individual fitness is calculated, and Sub4Step4 is applied. 12 Sub4Step4: The fitness is outputted. 13 Step5: The mean value of fitness is calculated. 14 Step6: The difference between the fitness vector and its mean value is calculated. 15 Step7: If the difference meets the design, turn to Step8; otherwise, turn to Step10. 16 Step8: The maximum of this difference is solved. 17 Step9: The optimal results are outputted. 18 Step10: Selection. 19 Step11: Coding. 20 Step12: Crossover. 21 Step13: Mutating. 22 Step14: Decoding. 23 Step15: The population is reorganized, and Step3 is applied. Table 2 Comparison results before and after optimization
Items L_{2} (mm) L_{1} (mm) d (mm) D (mm) Max sailing distance (m) Origin 52.92 103.89 4.50 12.56 25.00 Optimized 97.91 96.62 1.00 12.50 46.75 Table 3 Simulation parameters of optimized projectile
Physical quantity Value The density of the frustum of the cone (kg/m^{3}) 18.75×10^{3} The density of the cylinder part (kg/m^{3}) 2.7×10^{3} The length of frustum of cone (m) 9.79×10^{−2} The length of the cylinder part (m) 9.66×10^{−2} The projectile's maximum diameter (m) 1.25×10^{−2} The cavitator diameter (m) 1.0×10^{−3} The fluid density (kg/m^{3}) 1 000 The fluid temperature (℃) 20 The saturated vapor pressure (Pa) 2 338.8 The fluid kinematic viscosity (m^{2}/s) 1.006×10^{−6} The atmosphere pressure (Pa) 101 325 The initial depth of this projectile motion (m) 2 Acceleration of gravity (m/s^{2}) 9.8 V_{x}_{1} (m/s) 500 V_{y}_{1}, V_{z}_{1} (m/s) 0 ω_{x}_{1} (rad/s) 0 ω_{y}_{1}, ω_{z}_{1} (rad/s) 1 θ, φ, ψ (rad) 0 x, y, z (m) 0 m Mass of supercavitating projectile (kg) V_{x1} Velocity component on body coordinate system (O_{x1}) (m/s) V_{y1} Velocity component on body coordinate system (O_{y1}) (m/s) V_{z1} Velocity component on body coordinate system (O_{z1}) (m/s) V Resultant velocity of supercavitating projectile (m/s) ω_{x1} Angular velocity component on body coordinate system (O_{x1}) (rad/s) ω_{y1} Angular velocity component on body coordinate system (O_{y1}) (rad/s) ω_{z1} Angular velocity component on body coordinate system (O_{z1}) (rad/s) t Time (s) C_{0} Drag coefficient A_{c} Cavitator area (m^{2}) ρ Water density (kg/m^{3}) δ Angle between the velocity vector and O_{x1} (rad) α Attack angle (rad) α Sideslip angle (rad) C_{xf} Coefficient of friction drag S_{w} Surface area of the projectile wettingwater (m^{2}) k_{n} Smooth function coefficient h_{ef} Critical parameter (m) g Acceleration of gravity (m/s^{2}) θ Roll angle (rad) φ Yaw angle (rad) ψ Pitch angle (rad) R(t) Supercavity function (m) R Radius of projectile (m) η Angle of the tailslap force direction (rad) χ Angle between the body and the line through the mass center (rad) φ_{lim} Angle between the supercavity and body (rad) R_{n} Cavitator radius (m) C_{xfr} Coefficient of rotation friction drag I_{x1} Moment of inertia on O_{x1} (kg·m^{2}) I_{y1} Moment of inertia on O_{y1} (kg·m^{2}) I_{z1} Moment of inertia on O_{z1} (kg·m^{2}) x_{cm} Length of mass center to the cavitator (m) σ(t) Cavitation number C_{x0} Drag coefficient (0.82) x_{1} Projectile's length projection on the horizon (L_{p}) (m) η_{c} Cavity coefficient (0.85) L Projectile length (m) L_{1} Length of the cylinder part (m) D Diameter of the cylinder part (m) d Cavitator diameter (m) L_{2} Length of frustum of cone (m) z_{0} Mass center in calculation coordinate (m) I_{x}, I_{y}, I_{z} Moment of inertia in calculation coordinate (kg·m^{2}) v_{1} Volume of the cylinder part (m^{3}) v_{2} Frustum of cone volume (m^{3}) P_{a} Standard atmosphere pressure (Pa) P_{c} Supercavity pressure (Pa) R_{c}(t) Maximum radius of supercavity (m) L_{c}(t) Supercavity length (m) H The fin's length (m) 
Ahn SS, Ruzzene M, Scorcelletti F (2010) Configuration optimization of supercavitating underwater vehicles with maneuvering constraints. IEEE Journal of Oceanic Engineering 35(3): 647662. DOI: 10.1109/JOE.2010.2043576 Ashley S (2001) Warp drive underwater. Scientific American, 7279. DOI: 10.1038/scientificamerican050170 Chan KY, Kwong CK, Kremer GE (2020) Predicting customer satisfaction based on online reviews and hybrid ensemble genetic programming algorithms. Engineering Applications of Artificial Intelligence 95: 103902. https://doi.org/10.1016/j.engappai.2020.103902 Choi JH, Penmetsa RC, Grandhi RV (2005) Shape optimization of the cavitator for a supercavitating torpedo. Structural and Multidisciplinary Optimization 29(2): 159167. DOI: 10.1007/s0015800404660 CostaCarrapiço I, Raslan R, Gonzalez JN (2020) A systematic review of genetic algorithmbased multiobjective optimisation for building retrofitting strategies towards energy efficiency. Energy and Buildings 210: 141. https://doi.org/10.1016/j.enbuild.2019.109690 He QK, Wei YJ, Wang C, Zhang JZ (2013) Impact dynamics of supercavitating projectile with fluid/structure interaction. Journal of Harbin Institute of Technology 20(1): 101106. (in Chinese) Holland JH (1992) Adaptation in natural and artificial systems: An introductory analysis with applications to biology, control, and artificial intelligence. The MIT Press, 123140. DOI: 10.1063/pt.5.031409 Jiang J, Liu Z, Liu WH, Zhou WP (2009) Shape optimization of the cavitator for a supercavitating vehicle based on genetic algorithm. International Conference on Information Engineering and Computer Science, Wuhan, 14. DOI: 10.1109/ICIECS.2009.5366738 Kulkarni SS, Pratap R (2000) Studies on the dynamics of a supercavitating projectile. Applied Mathematical Modelling 24(2): 113129. DOI: 10.1016/S0307904X(99)000281 Li HZ, Jiang XH, Wang Y, Guo ZQ (2015) Intermediate ballistics. Beijing Institute of Technology Press, Beijing, 110. (in Chinese) Mirzaei M, Alishahi MM, Eghtesad M (2015) Highspeed underwater projectiles modeling: a new empirical approach. Journal of the Brazilian Society of Mechanical Sciences and Engineering 37(2): 613626. DOI: 10.1007/s4043001401907 Mirzaei M, Taghvaei H (2019) A novel configuration optimization for highspeed ventilated supercavitating vehicles. Ocean Engineering 179(2): 1321. DOI: 10.1016/j.oceaneng.2019.03.013 Nguyen Thai D, Horak V, Nguyen Van D, Dao Van D, Nguyen Van H, Duc Do L (2018) Ballistics of supercavitating projectiles. Advances in Military Technology 13(2): 237248. DOI: 10.3849/aimt.01243 RamosFigueroa O, QuirozCatellanos M, MezuraMontes E, Kharel R (2020) Variation operators for grouping genetic algorithms: A review. Swarm and Evolutionary Computation 60(1): 164. https://doi.org/10.1016/j.swevo.2020.100796 Rand, R, Pratap R, Ramani D, Cipolla J, Kirschner I (1997) Impact dynamics of a supercavitating underwater projectile. Proceedings of ASME Design Engineering Technical Conferences (DETC'97), Sacramento, 1417 Ruzzene M, Soranna F (2004) Impact dynamics of elastic stiffened supercavitating underwater vehicles. Journal of Vibration and Control 10(2): 243267. DOI: 10.1177/107754604773732179 Savchenko YN (2001) SupercavitationProblems and perspectives. 4th International Symposium on Cavitation, Pasadena, 18 Shafaghat R, Hosseinalipour SM, Lashgari I, Vahedgermi A (2011) Shape optimization of axisymmetric cavitators in supercavitating flows, using the NSGA Ⅱ algorithm. Applied Ocean Research 33(3): 193198. DOI: 10.1016/j.apor.2011.03.001 Shen C, Zhou KD, He L, Lu Y, Li JS (2020a) Research on muzzle disturbances of gun barrel with bore damage during firing. Acta Armamentarii 41(3): 442450. (in Chinese) DOI: 10.3969/j.issn.10001093.2020.03.004 Shen C, Zhou KD, He L, Lu Y, Li JS (2020b) The influence of bore damage on bullet's exterior ballistics process. Transaction of Beijing Institute of Technology 40(6): 623631. (in Chinese) Vaissier B, Pernot J, Chougrani L, Veron P (2019) Geneticalgorithm based framework for lattice support structure optimization in additive manufacturing. ComputerAided Design 110: 1123. https://doi.org/10.1016/j.cad.2018.12.007 Wang KJ, Rong G, Mu Q, Yi WJ, Guan J (2018a) Experiment and FEM simulation of supercavitation bullet penetrating protective plate. 1st ICDTThe First International Conference on Defence Technology, Beijing, 498503 Wang KJ, Rong G, Mu Q, Yi WJ, Chen ZH (2019) Double slapping effects on a supercavitation projectile. AIP Advances 9: 18. DOI: 10.1063/1.5053143 Wang KJ, Rong G, Yin HQ, Yi WJ (2020) Dynamic features of kinetic energy supercavitating vehicles. Applied Ocean Research 102: 119. https://doi.org/10.1016/j.apor.2020.102304 Wang R, Dang JJ, Yao Z (2018b) Hydrodynamic characteristic of fins of supercavitation projectile. Journal of Harbin Institute of Technology 50: 137142. (in Chinese) DOI: 10.11918/j.issn.03676234.201711068 Wang ZZ, Sobey A (2020) A comparative review between genetic algorithm use in composite optimisation and the stateoftheart in evolutionary computation. Composite Structures 233: 111739. https://doi.org/10.1016/jxompstruct2019.111739 Wei P, Hou J, Chen TF (2012) Analysis on the stability of supercavitation projectil. Applied Mechanics and Materials 157158: 5861. DOI: 10.4028/www.scientific.net/AMM.157158.58 Wei YJ, He QK, Wang C, Cao W, Zhang JZ (2013) Review of study on the tail slap problems of supercavitating projectile. Ship Science and Technology 35(1): 714. (in Chinese) Yi WJ, Wang ZY, Xiong TH, Zhou WP, Qian JS (2008a) Analysis of supercavity shape for underwater projectile with typical cavitator. Journal of Ballistics 20(2): 103106. (in Chinese) Yi WJ, Xiong TH, Wang ZY, Wang CM (2009a) Experimental researches on drag characteristics of supercavitation bodies at small cavitation number. Journal of Hydrodynamics 24(1): 16. (in Chinese) DOI: 10.16076/j.cnki.cjhd.2009.01.009 Yi WJ, Zhang Mu, Xiong TH, Chen ZH (2009b) Investigation on the velocity measurement of a underwater highspeed vehicle. Proceedings of the 9th National Congress on Hydrodynamics, Chengdu, 353358. (in Chinese) Yi WJ, Xiong TH, Liu YX (2008b) Numerical simulated research on charateristic of supercavity form around underwater highspeed projectile. Ship Science and Technology 30(4): 117121. (in Chinese) Yi WJ, Xiong, TH, Wang ZY (2008c) Structure design of hydrodynamic distribution for high speed bodies with fin. Proceedings of the 21st National Symposium on hydrodynamics, Jinan, 651656. (in Chinese) Zhang YW, Yuan XL, Deng F (2014) Fluid dynamics of supercavitating underwater vehicles. National Defense Industry Press, Beijing, 5760. (in Chinese) Zhang JZ, Zhang JS, Wei YJ (2010) Simulation of tailslap loads of supercavitating projectiles. 3rd International Conference on Information and Computing, Wuxi, 163166. DOI: 10.1109/ICIC.2010.135 Zhao C, Wang C, Weo YJ, Ma QP (2014) Analysis of the effect of mass center position on tailslap of supercavitating projectile. Journal of Beijing University of Aeronautics and Astronautics 40(12): 17541760 Zhao XF, Lyu XJ, Da L (2018) Modeling of the tail slap for a supercavitating projectile. 8th International Conference on Underwater System Technology: Theory and Application, Wuhan, 15. DOI: 10.1109/USYS.2018.8779216 Zhao XF, Lv XJ, Da L (2019) Modeling of the tail slap for an underwater projectile within supercavitation. Mathematical Problems in Engineering 2019(4): 110. DOI: 10.1155/2019/1290157 Zhou L, Li YH, Xue ZP, Han JZ, Zhang N (2016) Rigidflexiblecavity coupling trajectory and uncertainty trajectory analysis of supercavitating projectiles. Ocean Engineering 126(1): 138151. DOI: 10.1016/j.oceaneng.2016.08.032