Research on Aerodynamic Characteristics of Trans-Media Vehicles Entering and Exiting the Water in Still Water and Wave Environments

Abstract

The problem of aircraft entering and exiting water is a complex, nonlinear, strongly disturbed, and multi-coupled multiphase flow problem, which involves the precise capture of the air/water interface and the multi-coupling interaction between aircraft, water, and air. Moreover, due to the large difference in medium properties during the crossing, the load on the body will suddenly change. In this paper, the VOF (volume of fluid) algorithm is used to capture the liquid surface at the air/water interface, and since body movement is involved in this process, the overset grid technology is used to avoid the traditional dynamic grid deformation problem. In the process of this numerical simulation prediction, the effects of different water-entry angles and different water-entry heights on the body load and attitude of the trans-medium aircraft, as well as the cavitation evolution law of the body water entry are analyzed. On this basis, to simulate the authenticity and complexity of the water-entry environment, numerical wave-making technology was introduced to analyze the water-entry load, posture, and cavitation evolution law of the body under different wave environments. The numerical parameters under the condition of wave and no wave are compared, and the difference in water-entry performance under the condition of wave and no wave is analyzed.

1. Introduction

With the rapid development of small aircraft, ships, underwater submersibles, and other intelligent mobile platforms in search and rescue, detection, communications, combat reconnaissance, combat, and other types of military activities [1], civilian fields are widely used. However, the three-dimensional process of modern warfare is accelerating, and the environment in which combat equipment is located is gradually changing from the original single type to the multi-phase type [2]. To better adapt to the operational needs of modern warfare, the study of trans-media vehicles has become crucial.

The dynamic model of the trans-medium aircraft [3] is the basis for analyzing its motion characteristics and studying flight stability control. As the core key link of the dynamic model, the numerical prediction of the outlet/inlet section plays a vital role in the accurate response of the dynamic model. The water-entry process is a complex fluid phenomenon with strong nonlinear and unsteady characteristics [4]. From a mechanical point of view, the water-entry process of an object passing through the water surface at a relative speed mainly includes three stages [5]: the slamming when it contacts the free surface, the open cavity with large deformation of the free surface, and the underwater movement after the open cavity closes [6]. The first slamming stage is mainly a matter of impact dynamics, while the last two water-entry stages involve fluid–structure interaction of hydrodynamics and structural dynamics coupling [7]. In the initial stage of water entry, due to the large difference in the characteristics of the liquid and gas media, the slamming process will instantly generate a large load on the airframe, causing the part of the trans-medium aircraft in contact with the free surface to suffer a short-term, high-intensity pressure shock wave [8]. This propagates inside the body, causing the body to produce plastic deformation and have an oscillation effect. After entering the water, due to the large deformation of the free surface and the fluid–solid coupling between the hydrodynamic force and the body [9], cavitation evolution will occur on the surface of key components such as the body rotor, and the unsteady evolution and collapse of the cavitation will produce corresponding damage on the surface of the structure. Pressure pulsations cause significant dynamic responses in the body structure.

A crucial step in the numerical prediction of water entry is the simulation of large deformation of the free interface and the boundary conditions of the dynamic mesh [10]. The simulation of large deformation of free interface mainly includes two types: one is the interface tracking technology, and the other is the interface capture technology [11]. Typical free interface simulation methods include the level-set and VOF methods. The level-set method is a grid-based interface tracking method proposed by Osher and Sethian [12]. This method expresses the free interface as a continuous function, which is convenient for directly using mathematical tools The geometric characteristic parameters and surface tension of the interface are solved, and this method can also be extended to higher-order spaces [13]. The applicable scenarios of the level-set method are mainly distributed in the simulation of the movement of water droplets, the fusion of water droplets on the air-water free interface, and the evolution of bubbles in water [14,15,16,17]. Since this method does not consider the function of the continuity equation, it does not have the property of volume conservation [18]. Therefore, it often faces the problem of mass non-conservation when the numerical simulation is used to track free interfaces [19]. To solve the problem of the non-conservation of mass, Hirt, and Nichols [20] proposed the VOF method. The VOF (volume of fluid) method is an interface tracking method based on the Euler grid. In this method, the mutually incompatible fluid components share a set of momentum equations, and the phase interface tracking in the calculation domain is realized by introducing the variable of phase volume fraction ζ [21,22]. ζ represents the volume of one phase as a percentage of the mesh volume in which it is located. The VOF method can construct the interface by calculating the phase fraction of each grid unit in the entire computational domain [23]. It has very good and accurate numerical prediction results for dam breaks, surf waves, and Rayleigh Taylor instability. Ashgriz and Poo [24] proposed the FLAIR (flux line-segment model for advection and interface reconstruction) technology based on the traditional VOF method. The basic principle is to construct a straight-line segment with an inclination angle to approximately replace the physical free interface. Then calculate the fluid volume flow through the mesh approximation surface per unit of time and use it as a correction value of the numerical flux as a function of fluid volume. The PLIC (piecewise linear interface calculation) reconstruction technology proposed by Youngs [25] is also based on the traditional VOF method. The principle of the PLIC reconstruction technology is similar to that of the FLAIR reconstruction technology. This method is also widely used in liquid surface capture [26,27]. Sussman and Puckett [28] combined the superiority of level-set to calculate the curvature and the conservation of VOF to propose the CLSVOF (coupled level-set volume of fluid) method, which is widely used in numerical simulations of other gas–liquid two-phase interfaces, such as oceans, water conservancy, hydropower, liquid jets, etc. [29,30,31,32].

2. Airframe Design Scheme for Trans-Media Aircraft

Because of its ability to span multiple media, the trans-medium aircraft must adapt to the differences in physical properties of multiple media at the same time. Such prerequisites bring inevitable contradictions to its airframe structure design, but in general, the design concept of the cross-medium aircraft is to balance the difficulties that need to be overcome when the multi-medium crosses under the premise of satisfying the stable operation of the single-medium function, such as (1). However, in the underwater navigation stage, due to the nonlinear disturbance of water flow, to obtain sufficient steering stability, the aircraft should ensure that the aircraft body has considerable weight and structural strength. (2). In terms of the shape of the fuselage, based on the similarity principle of Reynolds number, the speed of the same fuselage shape in water and air media is quite different, which causes instability in the medium transition process of the trans-medium aircraft, which leads to the failure of crossing. (3). In terms of wing layout, when the aircraft is flying in the air, the pressure difference between the upper and lower surfaces of the wing generates lift to overcome gravity. When the flight speed is constant, the lift is related to the wing area. Influenced by airframe drag, the wing surface has a relatively short span and narrow chord. (4). In terms of pressure and tightness, underwater navigation has high requirements for pressure and tightness, while in air flight, especially in low-altitude flight, sealing and pressure performance are hardly considered. Even the pressure-resistant hulls used in high-altitude flights are no match for cabins used in underwater voyages.

Through the analysis of the above difficult problems, this paper proposes a water–air dual-power ducted structure scheme based on the multi-rotor aircraft architecture based on the basic principles of the trans-medium aircraft, as shown in Figure 1. By taking into account the design of the structure and power system, the detailed design of the airframe of the trans-medium aircraft is carried out.

As shown in the figure above, the overall structure of the trans-medium aircraft is based on the “X” layout of the four-rotor aircraft. The air flight and the vertical medium transition between water and air are realized through the air rotors. It has strong operational stability. The structural design of the cross-section duct realizes the reduction in water-entry slamming during the medium transition process and the duct-collection effect on the power of the underwater propeller. The underwater propellers are arranged in pairs in the center of the duct body in a “cross” shape to provide power for underwater multi-directional movement.

2.1. Layout and Principle of Airframe Propulsion System

This design combines the advantages of the layout of the four-rotor unmanned aerial vehicle. The air rotor power propeller adopts the “X” layout scheme, and the underwater propeller power adopts the “cross” layout scheme. The propeller types and their layout are shown in Figure 2, the red triangle represents the forward direction of the aircraft, and the front and back blades of the rotor are distinguished by different colors.

Through the description of the propeller types and their arrangement above, a brief introduction will be given to the working principles of the different missions of the trans-medium aircraft described in this paper, as shown in Figure 3. The four aerodynamic propellers are mainly responsible for the air flight mission of the aircraft. The basic propulsion principle is consistent with the drive mode of the quadrotor aircraft, and the purpose of flying in six degrees of freedom in space is achieved by changing the speed input of the four motors. Arrows of different sizes in the figure represent the differences in motor speed resulting in differences in the lift, thereby achieving the purpose of changing the aircraft’s motion attitude.

2.2. Internal Structure Design of Trans-Media Aircraft

2.2.1. Conduit Shape Design

For the layout design of the underwater propeller, the presence of the culvert body can effectively increase the efficiency of underwater propulsion. The duct body of the water–air dual-power ducted system mainly plays the role of protecting the propeller (underwater turbulence impact), supporting the body structure, and storing the integrated electronic control system inside the aircraft. Due to the existence of the duct body, the impact noise of the propeller tip is reduced, and the induced resistance is small so that the operating efficiency of the propeller is improved and greater thrust will be generated under the same power consumption condition, which greatly improves the aerodynamic performance of the ducted system. In the basic shape of the duct body design, a part of the cylindrical duct body section resembles a curved airfoil. The performance enhancement of ducted propellers has three geometrical features: a properly designed air inlet accelerates the flow of the propeller rotating medium, the gap between the propeller blades and the duct wall (tip clearance) improves the blade tip flow, and the diffuser partially suppresses the shrinkage of the flowing medium [33].

As shown in Figure 4, the duct is generally divided into three parts, the inlet section, the middle section, and the expanded outlet section. Its main shape and structure parameters are the inlet curvature R, the inner diameter D_c, the expanded outlet diameter D_e, the duct expansion angle Φ, and the culvert road height H. The above geometric variables can be dimensionless through the inner diameter or propeller diameter.

By placing the propeller inside the duct, the trans-medium aircraft has unparalleled safety for the operator and the surrounding environment, and through its VOTL capability, it can take off and land and operate in a narrow environment. With the support of the designed duct, the duct wall can effectively convert the propeller slipstream into thrust, thereby generating additional thrust.

2.2.2. Frame Structure Design

The frame structure is the bearing platform of the aircraft. According to the core, internal components of the trans-media aircraft and the overall layout structure under different media conditions, the air four-rotor support frame, and the underwater four-propeller thrust support frame are respectively designed. The air-supported frame only uses a wheelbase of 500 mm and is directly connected to the duct body in a rigid body. The frame layout adopts an “X” shape, as shown in Figure 5a, and the material is mainly carbon fiber, which greatly reduces the weight of the body and enhances the water entry of the aircraft, its anti-slamming ability. The main supporting component of the underwater support frame is the propeller drive motor. To realize the internal layout of the culvert body of the propeller, a cross-shaped frame body is used, and its material is also carbon fiber. The propeller support base and the cross-shaped frame are fixed to realize the culvert body. The layout of the four propellers in the “ten” shape inside the channel body is shown in Figure 5b,c. The propeller support base and the cross-shaped frame are fixed by pipe clamps. The design form of the overall aircraft frame is shown in Figure 5d. The frame bears the load of the underwater thrust device and is driven by independent power systems to complete single-medium and multi-medium crossing and multi-modal movements.

3. Based on OpenFOAM Numerical Calculation Theory

3.1. Numerical Water/Gas Two-Phase Flow Mathematical Model

For the numerical simulation of water–air two-phase flow, its solution models are mainly divided into two categories: high-phase fraction models and interface capture models [34]. In this study, the trans-medium vehicle is used to cross the water–air two-phase flow medium. Such a fluid–fluid problem with a clear boundary can be well-applied to the interface capture model. The two key issues of the interface capture idea are: (1) the interface should be as thin as possible; (2) the algorithm should be as stable as possible when the interface is unstable [35]. The interFOAM solver based on the VOF method is used to solve the interface of water–gas two-phase flow in OpenFOAM [36]. The interFOAM solver is based on incompressible, adiabatic, insoluble, non-mixing (mixing refers to the distribution of one fluid in another fluid, such as the movement of bubbles in liquid) two-phase interface capture to solve. The volume fraction method (VOF) is used to capture the free interface of two-phase flow, while the finite volume method [37] (FVM) is used for spatial discretization.

The complex flow conditions in nature are controlled by three basic physical principles, namely the law of conservation of mass, Newton’s second law, and the law of conservation of energy. These three basic physical principles correspond to three governing equations, that is, the governing equations of fluid mechanics (continuity equation, momentum equation, and energy equation), and these three equations are the mathematical description of the corresponding physical principles. According to the above three basic physical principles, the governing equation of the incompressible Newtonian fluid is as follows:

(1)$\left\{\begin{array}{c}\nabla ·\mathit{U}=0\\ \frac{\partial \mathit{U}}{\partial t}+\nabla ·(\mathit{U}\mathit{U})=-\nabla \frac{p}{\rho }+\nabla ·(\mathit{\upsilon }\nabla \mathit{U})\end{array}\right.$

where ρ is the density, U is the velocity of the control body, ρ is the pressure, and υ is the viscosity.

Through the above basic incompressible N-S equation, the momentum equation of the insoluble multiphase system can be obtained under the condition of gravity and source term:

(2)$\frac{\partial \rho \mathit{U}}{\partial t}+\nabla ·(\rho \mathit{U}\mathit{U})-\nabla ·\mathit{\tau }=-\nabla p+\rho g+\mathit{F}$

where F is the surface tension. The continuity equation can be expressed as:

(3)$\frac{\partial \rho }{\partial t}+\nabla ·(\rho \mathit{U})=0$

If two kinds of fluids are considered to be incompressible, that is to say, the fluid element of a certain fluid motion is considered, and its material derivative is 0, that is:

(4)$\frac{D\rho }{Dt}=\frac{\partial \rho }{\partial t}+\mathit{U}·\nabla p=0$

The density in Equation (4) can be solved for the pure convection equation when the velocity U is known. In the volume of fluid (VOF) method, $\zeta$ is defined to represent the phase fraction of the fluid. Consider a gas–liquid two-phase system of a certain grid unit, if the grid unit is filled with fluid, then $\zeta$ = 1; if the grid unit is filled with gas, then $\zeta$ = 0. If the value of $\zeta$ is between 0 and 1, there is gas–liquid mixing in this grid cell.

As shown in Figure 6, if the phase fraction of the entire grid is calculated to construct the interface, it will be found that the volume fraction is a step function in space and is discontinuous in space, so the reconstructed phase interface is discontinuous. The interface between two adjacent grids is discontinuous, and the physical quantity is also discontinuous when passing through the interface. This phenomenon is called parasitic current. The main work of the VOF method at present is to alleviate the parasitic flow phenomenon caused by the numerical method.

Next, the surface tension equation in the momentum equation of the insoluble multiphase system is deduced. The most important feature of the surface tension is that there is a sharp pressure drop Δp at the interface. Take the fluid microelement of the two-phase interface, and let p₁ be the downward pressure exerted by the air above the element; p₂ is the upward pressure exerted by the liquid below the surface microelement, so Δp = p₁ + p₂ can be obtained, and the surface tension of the curved surface microelement is defined as F. As shown, Figure 7 is a schematic diagram of the force on the surface micro-element at the two-phase interface.

In the figure, $s_1$ and $s_2$ are the lengths of the sides of the surface micro-element, so the area of the surface micro-element is $d{s}_{1}d{s}_2$. $d\alpha$ and $d\beta$ are the curvature angles of different curves on the surface micro-element, and $R_1$ and $R_2$ are the radii of curvature of different curves on the surface micro-element, respectively. The magnitude of the combined force on the surface micro-element is:

(5)${\mathit{F}}_p=(p_1-p_2)d{s}_{1}d{s}_2$

The definition of surface tension shows that it acts on the four edges of the surface micro-element and that the combined force of tension is vertical and balanced with ${\mathit{F}}_p$ above. Defining the coefficient of surface tension as $\sigma$, the magnitude of the tension on the $s_1$-edge is $\sigma d{s}_1$, and the magnitude of its component in the vertical direction is $\sigma d{s}_1\sin (d\beta /2)$. In the case of small $\beta$, $\sigma d{s}_1\sin (d\beta /2)$ ≈ $\sigma d{s}_1(d\beta /2)$. Similarly, the tension acting on the $s_2$-edge is $\sigma d{s}_2(d\alpha /2)$. Thus, the surface tension acting on the surface micro-element is:

(6)${\mathit{F}}_{\sigma }=\sigma (d{s}_{1}d\beta +d{s}_2\alpha )=\sigma (\frac{d\beta }{d{s}_1}+\frac{d\alpha }{d{s}_2})d{s}_{1}d{s}_2$

Now denote by $\mathit{n}$ the unit normal to the surface differential element, which points from a higher phase fraction to a lower phase fraction, modulo 1. Then we have:

(7)$\frac{\partial {\mathit{n}}_x}{\partial x}=\frac{d\alpha }{dx}=\frac{1}{R_1}$

(8)$\frac{\partial {\mathit{n}}_y}{\partial y}=\frac{d\beta }{dy}=\frac{1}{R_2}$

(9)$\frac{\partial {\mathit{n}}_z}{\partial z}=0$

So the surface tension expression can be translated as:

(10)${\mathit{F}}_{\sigma }=\sigma (\frac{1}{R_1}+\frac{1}{R_2})d{s}_{1}d{s}_2$

In addition because:

(11)$\nabla ·\mathit{n}=\frac{\partial \mathit{n}}{\partial x}+\frac{\partial \mathit{n}}{\partial y}+\frac{\partial \mathit{n}}{\partial z}=\frac{1}{R_1}+\frac{1}{R_2}$

Therefore:

(12)${\mathit{F}}_{\sigma }=\sigma \nabla ·\mathit{n}d{s}_{1}d{s}_2$

With the surface micro-element in force equilibrium, then we have:

(13)$\Delta p=p_1-p_2=\sigma (\frac{1}{R_1}+\frac{1}{R_2})=\sigma \nabla ·\mathit{n}$

Next, the surface tension is modeled according to the continuum surface force model to obtain:

(14)$\nabla p=\sigma \kappa \nabla \zeta$

where K represents the curvature at the interface.

Next, to solve the pressure term and body force term in the momentum conservation equation, define p_rgh in OpenFOAM as follows:

(15)$p_{rgh}=p-\rho \mathit{g}·\mathit{h}$

where the p_rgh is the reference pressure, which is obtained by splitting the pressure $p$. It has very good stability during the solution calculation; therefore, in general form, the reference pressure is used as the variable for the iterative calculation.

Since the two-phase flow problem is involved, the density in the basic governing equation is related to the physical properties of water and gas, so for the unification of the equations, ρ is defined as follows:

(16)$\rho =\zeta {\rho }_1+(1-\zeta ){\rho }_2$

where ${\rho }_1$ and ${\rho }_2$ are different phase densities, respectively.

Therefore, by substituting the previously derived expressions for each variable into the incompressible N-S Equation (1), the fluid micro-element control equation under the VOF method is obtained as follows:

(17)$\frac{\partial \zeta }{\partial t}+\mathit{U}·\nabla \zeta =0$

(18)$\frac{\partial \rho \mathit{U}}{\partial t}+\nabla ·(\rho \mathit{U}\mathit{U})-\nabla ·(\mathit{\upsilon }\nabla \mathit{U})-\nabla \mathit{U}·\nabla \mathit{\upsilon }=-\nabla p_{rgh}-g·h\nabla p+\sigma \kappa \nabla \zeta$

(19)$\nabla ·\mathit{U}=0$

It can be seen from the above governing equations that there is no essential difference between the VOF equation and the single-phase flow equation, both of which contain a momentum equation and a continuity equation. The continuity equations are exactly the same. The momentum equation, the VOF model is supplemented with a space-varying viscosity and an additional surface tension term. In addition, a phase transport equation is included. To sharpen the free interface, OpenFOAM uses the method proposed by Waller [38] to add artificial convection terms to the phase equation to squeeze the phase fraction near the phase interface to counterbalance the solution. The equation causes the ambiguity of the phase boundary caused by numerical dissipation to ensure the clarity of the interface. According to the idea of adding an artificial convection term, the VOF model can be expressed as:

(20)$\frac{\partial \zeta }{\partial t}+\nabla ·(\zeta \mathit{U})+\nabla ·(\zeta (1-\zeta ){\mathit{U}}_c)=\zeta \nabla ·\mathit{U}$

The third item in Equation (20) is an artificially added compressible item. According to the mathematical properties of this item, when $\zeta$ = 0 or 1, the calculated value is 0; that is, its value is 0 in pure phase (non-interface). A non-zero value exists when $\zeta$ is between 0 and 1. U_c is the compression velocity, the direction is the same as the phase gradient, so:

(21)${\mathit{U}}_c=c\left|\mathit{U}\right|\frac{\nabla \zeta }{\left|\nabla \zeta \right|}$

where c represents a controllable compression factor. When c = 0, there is no compression effect. The larger c is, the faster and more obvious the compression effect will be. The final phase equation embedded in the OpenFOAM solver is as follows:

(22)$\frac{\partial \zeta }{\partial t}+\nabla ·(\zeta \mathit{U})+\nabla ·(\zeta (1-\zeta )c\left|\mathit{U}\right|\frac{\nabla \zeta }{\left|\nabla \zeta \right|})=\zeta \nabla ·\mathit{U}$

3.2. Trans-Medium Aircraft Multi-Medium Spanning Moving Mesh Technology

It is very convenient to use dynamic overset grids for complex structures with relative motion. The relative motion between subdomains does not require grid deformation, let alone regenerate the grid, and only needs to define the law of grid motion under the condition that the subdomains overlap each other.

The advantage of the overset grid method is that it allows complex integrated or separated geometries to be processed using multiple overlapping meshes. Each component grid can handle a part of the domain and can be moved freely or as specified. Before solving the flow-governing equations, all component meshes need to be assembled to determine the information connectivity between different domains. This process includes the following key steps: (1) cutting holes, (2) overlapping boundary search, and (3) contribution unit search. The specific schematic diagram is shown in Figure 8.

3.3. Numerical Wave Making

The application scenario of modern numerical wave-making technology is the simulation of buildings in the marine environment. In this paper, to verify the ability of the medium-crossing medium of the trans-medium aircraft in the wave environment and compare the load impact on the body with the data under the condition of still water, the aerodynamic performance of the trans-media aircraft entering and exiting water under different external environments is obtained.

The numerical wave-making technology in this paper is still the numerical solution of the two-phase flow, and the solver used is the same as that of the aircraft entering the water under still water conditions, which is overInterDyFoam, which uses finite volume discretization and the volume fraction of fluid (VOF) method to solve the 3D Reynolds-averaged Navier–Stokes (RANS) equations for two incompressible phases and handle dynamic meshes. The calculation domain of numerical wave-making is shown in Figure 9.

4. Numerical Prediction of Aircraft Entering and Exiting the Water

4.1. Mesh Independence Test

To verify the grid independence, the trans-medium aircraft is selected as the research object, and different levels of grids are calculated under the condition of the air rotor of the aircraft at 2000 r/min. The air rotor model is APC-1047SF.

Figure 10a–c are the rotor thrust curves, the tip speed curve when R/r = 1.1, and the tip speed curve when R/r = 1.2 under different grid levels, respectively. R is the air rotor radius and r is the sampling radius.

From the above schematic diagram of the curve change, it can be seen that better computational accuracy can also be obtained when the grid is controlled at a medium level. Therefore, in this study, the computational grid level for the numerical prediction of aircraft entry and exit is controlled at a medium level, which can effectively save computational power and obtain more accurate results with limited computational resources.

To verify the convergence of the solution iteration process, the residual curves of the turbulence parameters involved in the RANS modeling method are extracted by sampling in this paper, and the curves are shown in Figure 11.

By observing the above curves, it can be found that each parameter gradually converges to a constant value as the solution time increases. By such a trend of functional variation, it can be judged that the numerical simulation converges during the computational solution of the grid-independent numerical experimental validation. The simulated debugging parameters are shown in Table 1, Table 2 and Table 3.

4.2. Mesh Generation Verification Analysis

Numerical prediction of trans-medium aircraft entering and exiting water involves body movement, so there is a dynamic grid problem. During the solution process, the movement and deformation of the grid will have a huge impact on the solution accuracy and solution convergence, so the overlapping grid method will be used in the pre-processing grid generation process for the numerical simulation of aircraft entering and exiting the water. The overlapping grid method needs to establish two sets of grids, in which the calculation domain of the external flow field of the trans-media aircraft is used as the background grid, and the aircraft itself is used as a moving overlapping grid. The role of moving in the grid, the numerical communication between the two sets of grids adopts the Lagrange interpolation shape function method because this process does not involve the deformation of the moving grid, thus ensuring the stability and convergence of the solution process. However, the numerical communication in the solution process leads to a longer calculation time for this kind of mesh method than that of traditional deformable dynamic mesh. Figure 12 is a schematic diagram of the grid generation of the aircraft entering and exiting the water.

4.3. Research on the Water-Entry and -Exit Characteristics of the Aircraft in the Static Water State

4.3.1. Water-Entry Characteristics at Different Angles

In this paper, the initial conditions for the water-entry characteristics of the trans-media aircraft at different angles are set as follows: the water-entry height is 0.4 m, the initial velocity of water entry is 0, the free fall hits the water surface, and the body enters the water body under the action of gravity. To verify the water-entry performance of the aircraft at different angles, this paper sets its water-entry angles to 0°, 10°, 20°, 30°, and 40°, respectively. The angle of entry into the water is defined as the angle between the plane of the air rotor blades and the water surface. The vehicle is in the process of entering the water and is to turn off all power devices at a certain angle and height of free fall into the water.

To reveal the influence of the static water flow disturbance on the speed and impact load of the aircraft during the water-entry process, this paper takes the aircraft at a height of 0.4 m and an inclination angle of 40° as an example and obtains the phase change of the water flow field and the motion of the aircraft through numerical simulation calculations. The streamlined diagram is shown in Figure 13.

It can also be seen from the above water flow phase change diagram and body motion streamline diagram that at the moment of water entry and slamming, the high-speed zone is near the water/air interface, and as the water-entry process proceeds, the water flow buffers the body speed to gradually decrease; the kinetic energy of the body movement is transformed into the kinetic energy of water flow disturbance cavitation sputtering.

Next, to analyze the vortex structure of the flow field during the water-entry process of the trans-medium aircraft, it is necessary to be able to identify and detect the vortex in the flow field through an effective visualization method. In this paper, the Q criterion [39] method is used for the above analysis.

The theory of the Q criterion comes from the velocity gradient tensor. The velocity gradient tensor expression is as follows:

(23)$\frac{\partial u_i}{\partial x_j}=0.5[\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}]+0.5[\frac{\partial u_i}{\partial x_j}-\frac{\partial u_j}{\partial x_i}]$

Among them, the symmetric part is denoted as S, which is usually called the strain rate tensor; the antisymmetric part is denoted as Ω, which is usually called the rotation rate or vorticity tensor.

(24)$S=0.5[\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}]$

(25)$\mathrm{\Omega }=0.5[\frac{\partial u_i}{\partial x_j}-\frac{\partial u_j}{\partial x_i}]$

The Q value is defined as the second invariant of the velocity gradient tensor, namely:

(26)$Q=0.5[|\left|\mathrm{\Omega }\right|{|}_{\mathit{F}}^2-\left|\right|S|{|}_{\mathit{F}}^2]$

It can be seen from the above formula that the positive value of Q indicates the region where the vorticity is dominant in the flow field, and the negative value indicates the region where the strain rate or viscous stress is dominant. As shown in Figure 14, when the aircraft enters the water at 40°, the cloud diagram of the vorticity changes the Q criterion when entering the water.

The flow mechanism of this changing cloud image reveals that the vortex at the boundary and front end of the aircraft body is mainly affected by shearing, and the continuous development of the cavity leads to the complexity and fragmentation of the positive and negative Q value alternating area of the cavity development. This indicates that the formation of the cavitation region is mainly affected by the joint action of rotation and shear.

4.3.2. Water-Entry Characteristics at Different Heights

This paper explores the aerodynamic characteristics of the vertical entry of the aircraft into the water at different altitudes. The initial conditions for the aircraft are set as follows: the aircraft maintains a vertical attitude and enters the water, and the initial heights of entering the water are set to 0.2 m, 0.3 m, and 0.5 m, respectively. Through the data sampling and analysis of the numerical results of the aircraft entering the water at different heights, the body load impact curves and speed curves suffered by the airframe in the process of entering the water at different heights are obtained, as shown in Figure 15 and Figure 16.

It can be seen from the curve of load impact change that, when entering the water at different heights, there is a second extreme value in the peak value of the load impact on the body, indicating that in the process of entering the water, there is a periodical change in the load impact of the body in a part of a very short period. This process will produce alternating loads and stress mutations on the body structure, thus reflecting that the structural stability of the trans-medium aircraft needs to be extremely guaranteed. According to the different water-entry heights, the peak value of the load impact is also different. It can be seen from the curve that as the initial height increases, the water-entry slamming is more severe and the load peak value is also higher. It can be seen from the graph of body speed change that because of the difference in the initial height, the speed at which the aircraft reaches the water–air interface is also different. The impact effect is also more obvious, so the rate of change of its velocity is also greater, and the momentum impact is also stronger.

The Figure 17 shows the cavitation phase transition diagram of the aircraft vertically entering the water at different heights.

Through the above cavitation phase transition diagram, it can be found that when different initial heights enter the water, the degree of fragmentation of the water flow cavitation is different at the moment of water collision. The higher the height, the greater the speed of water entry, and the stronger the degree of cavitation fragmentation. In addition, the bigger it is, the greater the kinetic energy of the body is converted into the kinetic energy of water flow. As mentioned above, the streamlined distribution map of the flow field body is obtained through numerical calculation, as shown in Figure 18.

Through the above-streamlined distribution diagram, it can be found that the velocity field of the overall calculation domain presents a dynamic unsteady change with the evolution of the cavitation phase transition. As the water-entry speed increases, the high-speed area at the tail also gradually moves up. This high-speed area will gradually expand upward and finally form a certain far-field airflow.

4.3.3. Aircraft Water-Exit Characteristics

After the aircraft completes the underwater diving mission, it transitions to the air medium to perform the air flight mission. This process must go through the water-exit action from water to air. The water and air media are quite different, so this process will lead to sudden changes in the body load at the moment of water exit. In addition, the airframe undergoes continuous and stable water resistance load until the load almost disappears, which will also lead to the imbalance of the control of the trans-medium aircraft. Therefore, the numerical prediction of the water-outlet characteristics of the aircraft can provide mechanical sensing analysis for the imbalance adjustment of the controller in the water-outlet stage.

Consistent with the previous analysis method, based on the OpenFOAM open-source numerical platform, the numerical aerodynamic characteristics of the water-exit process are calculated using the overlapping grid method. The initial conditions are set as follows: the initial speed of the trans-media aircraft is 2 m/s, and the vertical water-exit attitude is used to complete the water-exit process. In the process of the aircraft out of the water, the underwater propellers’ upward propulsion close to the water’s surface, when the air rotor is exposed to the water’s surface, starts the rotor, the underwater propeller, and the air rotor dual thrust to complete the role of leaving the water surface. In the simulation process, the underwater propeller and air rotor are divided into rotating areas, and the stationariness and motion of the rotating areas are controlled to achieve the above physical process. Figure 19 shows the water load change curve of the body.

It can be seen from the change curve that as the aircraft gradually moves from the deep-water area to the shallow-water area, the body load will decrease. When approaching the free water surface at a certain initial speed, the impact effect still exists, and the body load will reach the peak value in a short time. When the body jumps out of the water, the load of the body is gradually reduced until the body jumps out of the water completely and the load reaches a smaller value. Due to the inertia effect of the body and the viscosity effect of the water flow during the water exit, the body still entrains the water flow, so the load is still at the same level. The value underwater is almost the same. As it continues to rise, the water flow leaves the fuselage, and the load reaches a minimum value. At this time, the load on the aircraft is equal to that in the air flight stage.

Figure 20 shows the change curve of the attitude angular velocity of the trans-medium aircraft during the water exit.

Through the analysis of the above angular velocity change curve, it can be seen that the pitch attitude of the aircraft remains stable during the process of exiting the water, while the roll attitude and yaw attitude have some changes because of the disturbance of the water flow in this process, and the change range of the yaw attitude is larger than that of the roll. The attitude is a little bigger, but it remains stable after changing to the peak value, while the roll attitude is irregular and nonlinear under the wave disturbance because the water flow disturbance wave is consistent with the roll direction during the water-discharge process. According to the analysis of the simulated attitude change of the aircraft entering and exiting the water, the controller can be adjusted reasonably, effectively, and stably, thus providing a certain simulation practice basis for the control design of the trans-medium aircraft.

The generation of the above mechanical properties is related to the cavitation phase transition of the water flow during the water-exit process. Through the numerical calculation of the convection field, the schematic diagram of the cavitation phase transition of the water flow during the water exit of the aircraft is obtained, as shown in Figure 21 and Figure 22.

Through the above longitudinal and top-down views of the water flow cavitation phase transition diagram of the aircraft’s water-exit process, the analogy with the body load and attitude change results obtained above shows that since the aircraft body is close to the free surface of the water, the influence of the free surface on the water-exit process cannot be ignored. When the aircraft passes through the free surface, the water resistance will suddenly decrease, and it is often asymmetrical concerning the body, which will cause the body to deflect, which is consistent with the attitude change results mentioned above. In addition, the water discharge hits the free liquid surface to cause the liquid level to fluctuate, which further increases the randomness of the water-discharge process. The influence of water flow cavitation on the viscosity of the fuselage during water discharge cannot be ignored, just like the load results of the fuselage mentioned above.

4.4. Research on Water-Entry Characteristics in Wave State

In the real physical environment, the trans-medium aircraft will be impacted by strong nonlinear waves during the process of entering and exiting the water. Therefore, to be close to the real physical scene, the water column impact nonlinear wave is used in the numerical prediction of the aircraft entering and exiting the water in the wave environment. The basic method is as follows: set a water column with a certain height and width in the calculation domain, and these water columns will hit the free water surface under the action of gravity, thereby generating strong nonlinear wave disturbance. The schematic diagram of the shock wave generation is shown in Figure 23.

As shown in the figure above, this paper uses a fixed water-column width of 0.7 m and different water-column heights to simulate the wave intensity variable. The idea of variable parameter analysis in this paper has two aspects. On the one hand, under the condition of fixed wave intensity, the aircraft enters the water at different angles. On the other hand, under the condition that the aircraft enters the water at a fixed angle, the water-entry characteristics of the aircraft are analyzed with different wave intensities.

4.4.1. Fixed Wave Intensity (the Height of the Water Column Is a Fixed Value)

Set the height of the water column to 0.2 m, and the aircraft enters the water at 0°, 10°, 20°, 30°, and 40°, respectively. Through numerical calculation, the schematic diagram of the wave height change curve at the central axis position of the calculation domain under different water-entry angles is shown in Figure 24.

From the changing trend of the above curves, it can be judged that the nonlinearity produced by the water-column impact wave-making method is very strong, which is closer to the real physical scene. It can also be seen from the figure that the wave disturbance on the water surface when entering the water vertically is stronger than that when entering the water obliquely.

By visualizing the calculation results, a schematic diagram of the cavitation evolution of the phase volume fraction of the aircraft entering the water at different angles is obtained, as shown in Figure 25.

It can be seen from the schematic diagram of the cavitation evolution of the volume fraction of the water-entry phase from the above-mentioned angles that when the body enters the water vertically, the degree of wave interference is the largest, and the cavitation evolution of the water flow is the most severe; when the body is inclined into the water, the impact of wave disturbance on the body posture is larger. In the case of a small inclination angle, the body posture can naturally achieve self-warming, but as the inclination angle increases, the body posture loses its self-stabilization condition and rolls significantly during the water-entry process because of the drooping impact of waves. Therefore, when it comes to the controller, external interference should be added to make the aircraft perform feedback adjustments under the drooping impact of the wave to achieve stability of the body attitude.

To quantitatively analyze the attitude change of the aircraft, the angular velocity change curves of the aircraft entering the water at different angles are obtained through numerical simulation, as shown in Figure 26.

It can be seen from the above curves that the yaw attitude of the airframe changes significantly in the state of vertical entry into the water, while the pitch and roll attitudes do not change significantly. In the case of entering the water at a small angle, the roll attitude of the body changes slightly because of the impact of waves, and the change in the pitch attitude is a necessary change because of the influence of the restoring moment of the body. After the water-entry process is completed, the body can still achieve stability. When the body enters the water at a relatively large angle, there is a second peak phenomenon in the pitching angular velocity. The generation of the initial peak value is still affected by the restoration torque. With the completion of the water-entry process, the body cannot achieve self-stabilization conditions because of the excessive water-entry angle under the sagging impact of the waves, resulting in the loss of balance of the body posture, roll posture, and yaw. Attitude will change significantly

Figure 27 is a schematic diagram of the change of the mechanical parameter curve during the body’s entry into the water, showing the changing trend of the body’s load impact, body speed, and body movement displacement when entering the water at different angles.

It can be seen from the load curve that the load of the body entering the water vertically is greater than the load impact of the oblique entry into the water, and the load under the small angle of the water entry will gradually tend to balance, while the load of the large angle of the water entry will show a nonlinear change that is due to instability. When the body speed enters the water vertically or at a small angle, after the water-entry process is completed, the speed will increase linearly because of the drooping impact of the wave. The large-angle water-entry velocity curve also has a secondary peak phenomenon, which is still caused by the drooping impact of the wave, which causes the body to become unstable, so that the body moves downward in the form of a rollover, which greatly reduces the resistance on the body; thus, this makes the body obtain greater downward speed. It can be seen from the displacement curve that after the body rolls over, the speed increases, resulting in a greater displacement of the body at the same time.

Figure 28 is a binary diagram of the phase volume fraction. Through the binary diagram, we can observe the generation and closure of cavitation under different angles of water entry.

To visualize the impact of wave drooping on the body, this paper performs post-processing on the velocity flow field in the calculation domain to obtain the flow field velocity streamline diagram in the vertical direction, as shown in Figure 29.

Through the visual streamline diagram processing, under the condition of large-angle water entry, the impact of wave drooping impact on the rollover of the body can be seen. When the angle of the body deviates to one side, the drooping force of the wave will form a flow around this place, so that the body that is already inclined cannot maintain its original state. This is the root cause of the body rollover.

4.4.2. Fixed Water-Entry Angle (Changes Wave Intensity)

Enter the water at a fixed angle of 40° and create waves of different intensities by changing the height of the impacting water column. In this paper, the heights of the impacting water columns are set to be 0.2 m, 0.3 m, and 0.4 m, respectively. Under different wave intensities, the wave height curves at the central axis of the calculation domain are obtained, as shown in Figure 30.

Similar to the previous ones, the generated waves have strong nonlinearity, which is close to the scene of real aircraft entering the water.

Figure 31 is a schematic diagram of the cavitation evolution of the phase volume fraction caused by the body entering the water at an angle of 40° under different wave intensities.

From the figure above, it can be seen that in the calculation time of 0.6 s, the aircraft was unstable when the height of the water column was 0.2 m, but it was not unstable when the height of the water column was 0.3 m and 0.4 m. This is because, with the increase in the height of the impacting water column, the wave strength is too large, so that when the water column initially sinks, the weight of the body cannot overcome the influence of the wave force, so that the body obtains an upward moment, which cancels out the rollover moment. Therefore, the aircraft can still realize the self-stabilizing condition in a short time. With the arrival of the next wave-dropping moment, the body will still roll over and become unstable.

The Figure 32 shows the change curve of the mechanical characteristics of the body, which is similar to the previous one. They are the change curve of the body load, the change curve of the body motion speed, and the change curve of the body motion displacement.

It can also be seen from the above mechanical characteristic curve that when the wave intensity is too large, the body will be temporarily given an upward wave moment to maintain the posture balance. Through the speed curve, it can be seen that under the condition of the highest wave intensity, such an upward wave moment will temporarily change the speed direction of the body and make the body move upward.

To show the evolution law of water flow cavitation more intuitively, similar to the above, a binary map of the phase volume fraction of the body entering the water is extracted, as shown in Figure 33.

The above binary diagram clearly shows the evolution diagram of water-entry cavitation under two different wave intensities. As mentioned above, under the action of upward wave torque, the body gradually achieves self-stabilization in a short period. The expansion and closure of the cavity during the water-entry process have also been displayed.

5. Conclusions

In this paper, the VOF method with overset meshes is used to capture the air–water interface for vehicle entry into the water, and the results of mechanical properties are obtained for different angles and different heights under static water conditions and different angles and different wave intensities under wave conditions. Under still water conditions, due to the great disparity of the water–air medium, a restoring moment opposite to the tilt direction will be generated during the tilt angle entering the water, resulting in obvious changes in the pitching attitude of the aircraft. Under the condition of small wave intensity, after the body enters the water at a large angle, due to the influence of the wave-dropping moment, the attitude of the body will be unstable and rollover will occur. Under the condition of large wave intensity, due to the large drooping moment squeezing the fluid around the body, the body generates an upward moment, which cancels out the rollover moment, so that the body can achieve attitude balance in a short time. Through the analysis of the load curve graph of the body into the water and the horizontal projection view of the phase volume fraction, we know that in the horizontal projection, the area of the water bubble evolution of the vertical water entry is larger, resulting in the vertical resistance to the body also being larger, and the submersion speed also being slower. Therefore, the vertical resistance to water entry can be reduced by tilting the water entry at a certain angle. The above research on the hydraulic mechanism of entry and exit provides a theoretical basis for the design of the medium-crossing controller of trans-medium aircraft.

The different stages through which the trans-media vehicle enters the water are analyzed by the schematic diagram of the evolution of the water bubble for different operating conditions as follows:

• (1). At the beginning of the water entry, the wet water area of the vehicle increases rapidly, and if the speed of entry is relatively large, there will be a relatively stable vacuole separation line; the water is separated from the surface of the body, and the entry vacuole begins to form.

• (2). After the vehicle enters the water for a period of time, the vacuole will be connected to the atmosphere; air continues to fill the rear space of the object entering the water, and the vacuole continues to increase, which is partly air and partly steam. As the vacuole increases, the buoyancy force on the vehicle increases.

• (3). The next step in the development of the incoming vacuole is the closure of the vacuole, where the air on the surface of the water no longer enters the vacuole. When the forces that determine vacuole closure (hydrostatic pressure, dynamic pressure of air flowing in the vacuole, and surface tension) dominate, the vacuole begins to narrow, neck, and finally close. When the bubble closes, the waters moving inward collide together and produce upward and downward water jets.

• (4). After the bubble is closed, as the object continues to move, the bubble, due to the water hostage effect, gradually reduces to completely disappear. After the disappearance of the air bubble, the object in the water begins to enter the full wetting motion; at this time, the object is far from the free surface of the water, and the free surface of the flow of the impact can be disregarded. The flow is unbounded.

Author Contributions

Conceptualization, J.W. and Y.-L.C.; methodology, J.W.; software, D.Z.; validation, D.-P.C.; formal analysis, J.W.; investigation, Z.C.; resources, J.W. and Y.-L.C.; data curation, J.W.; writing—original draft preparation, J.W.; writing—review and editing, J.W.; visualization, J.W.; supervision, Y.-L.C.; project administration, X.-Y.H.; funding acquisition, Y.-B.S. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

Not applicable.

Acknowledgments

This work is supported by the Jilin Province Key Science and Technology R&D project under Grant No. 20210203175SF, the Aeronautical Science Foundation of China under Grant No. 2019ZA0R4001, the National Natural Science Foundation of China under Grant No. 51505174, the Foundation of Education Bureau of Jilin Province under Grant No. JJKH20220988KJ, and the Interdisciplinary Integration, Innovation, and Cultivation Project of Jilin University under Grant No. JLUXKJC2020105.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. Schematic diagram of the overall structure of the water–air dual-power ducted trans-medium aircraft.

Figure 2. Schematic diagram of thruster type and its arrangement.

Figure 3. Schematic diagram of the working principle of the power system of the trans-medium aircraft: (a–d) translation, roll, pitch, yaw (air layout); (e–h) translation, roll, pitch, yaw (underwater layout).

Figure 3. Schematic diagram of the working principle of the power system of the trans-medium aircraft: (a–d) translation, roll, pitch, yaw (air layout); (e–h) translation, roll, pitch, yaw (underwater layout).

Figure 4. Schematic diagram of the cross-sectional structure of the duct body.

Figure 5. Schematic diagram of the frame structure design of the trans-media aircraft: (a) air quadcopter support frame; (b) cross frame structure; (c) underwater propeller support rack; (d) aircraft overall rack layout design.

Figure 6. Spatial distribution of phase fraction and its interface reconstruction by VOF method.

Figure 7. Schematic diagram of the force on the surface microelement at the interface of two phases.

Figure 8. Schematic diagram of the assembly principle of overset grid.

Figure 9. Schematic diagram of the numerical wave-making calculation domain.

Figure 10. Schematic diagram of variable curves under different grid levels: (a) rotor thrust curve; (b) the tip speed curve when R/r = 1.1; (c) the tip speed curve when R/r = 1.2.

Figure 11. Residual curve.

Figure 12. Schematic diagram of grid generation for numerical prediction of water inflow and outflow.

Figure 13. Water flow phase change and body motion streamline diagram of the inflow flow field.

Figure 14. Q-value cloud chart when entering the water at 40°.

Figure 14. Q-value cloud chart when entering the water at 40°.

Figure 15. Variation curve of water load impact at different heights.

Figure 16. Velocity curves of the body entering the water at different heights.

Figure 17. Phase transition diagram of water flows cavitation at different heights: (a) 0.2 m high into the water; (b) 0.3 m high into the water; (c) 0.5 m high into the water.

Figure 17. Phase transition diagram of water flows cavitation at different heights: (a) 0.2 m high into the water; (b) 0.3 m high into the water; (c) 0.5 m high into the water.

Figure 18. Velocity flow line diagram of the calculation domain under different altitude conditions.

Figure 19. Body water load change curve.

Figure 20. Angular velocity variation curve of aircraft out of water attitude.

Figure 21. Schematic diagram of cavitation phase transition of water flow in water-outlet process (portrait).

Figure 22. Schematic diagram of cavitation phase transition of water flow in water-outlet process (overlook).

Figure 23. Schematic diagram of shock wave generation.

Figure 24. Schematic diagram of wave height curve changes under different water-entry angles.

Figure 25. Schematic diagram of cavitation evolution of water phase volume fraction at different angles: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 25. Schematic diagram of cavitation evolution of water phase volume fraction at different angles: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 25. Schematic diagram of cavitation evolution of water phase volume fraction at different angles: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 26. Schematic diagram of the angular velocity change curve of entering the water at different angles: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 27. Schematic diagram of the change curve of the mechanical parameters of the body: (a) body load impact; (b) body movement speed; (c) body motion displacement.

Figure 28. Phase volume fraction binary map: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 28. Phase volume fraction binary map: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 29. Schematic diagram of vertical velocity streamlines in the computational domain flow field: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 29. Schematic diagram of vertical velocity streamlines in the computational domain flow field: (a–e) are, respectively, 0°, 10°, 20°, 30°, and 40° into the water.

Figure 30. Variation curves of wave height at the central axis position under different wave intensities.

Figure 31. Schematic diagram of cavitation evolution: (a) 0.3 m; (b) 0.4 m.

Figure 31. Schematic diagram of cavitation evolution: (a) 0.3 m; (b) 0.4 m.

Figure 32. The mechanical characteristic curves of the body: (a) load curve; (b) speed curve; (c) displacement curve.

Figure 32. The mechanical characteristic curves of the body: (a) load curve; (b) speed curve; (c) displacement curve.

Figure 33. Phase volume fraction binary map: (a) 0.3 m; (b) 0.4 m.

Table 1

Mesh convergence analysis.

Mesh Type	Number of Grid Cells	Grid Base Size/m
		Background Area	Overlapping Area
coarse	1,554,953	0.13	0.06
medium	3,663,283	0.10	0.05
fine	4,437,138	0.07	0.04

Table 2

Boundary Conditions.

Turbulence Parameters	Inlet	Outlet	Wall	Blades
k	fixedValue	inletOutlet	kqRWallFunction	kqRWallFunction
nut	calculated	calculated	nutUSpaldingWallFunction	nutUSpaldingWallFunction
omega	fixedValue	inletOutlet	omegaWallFunction	omegaWallFunction
p	fixedFluxPressure	totalPressure	fixedFluxPressure	fixedFluxPressure
U	fixedValue	pressureInletOutletVelocity	fixedValue	movingWallVelocity

Table 3

Simulation of debugging parameters.

Mesh Type	Turbulence Model	Physical Time Step	Maximum y+ Value	Whether to Use Wall Functions	Solving Time
coarse	k-omega two-equation model	$3\times {10}^{-4}$ s	40	YES	0.88 h (64 CPUs)
medium	k-omega two-equation model	$2\times {10}^{-4}$ s	35	YES	2.07 h (64 CPUs)
fine	k-omega two-equation model	$1\times {10}^{-4}$ s	32	YES	2.54 h (64 CPUs)

References

1. Ma, Z.; Chen, D.; Li, G.; Jing, X.; Xiao, S. Configuration Design and Trans-Media Control Status of the Hybrid Aerial Underwater Vehicles. Appl. Sci.; 2022; 12, 765. [DOI: https://dx.doi.org/10.3390/app12020765]

2. Zhang, X.; Huang, J.; Huang, Y.; Huang, K.; Yang, L.; Han, Y.; Wang, L.; Liu, H.; Luo, J.; Li, J. Intelligent amphibious ground-aerial vehicles: State of the art technology for future transportation. IEEE Trans. Intell. Veh.; 2022; pp. 1-19. [DOI: https://dx.doi.org/10.1109/TIV.2022.3193418]

3. Bi, Y.; Lu, D.; Zeng, Z.; Lian, L. Dynamics, and control of hybrid aerial underwater vehicle subject to disturbances. Ocean Eng.; 2022; 250, 110933. [DOI: https://dx.doi.org/10.1016/j.oceaneng.2022.110933]

4. Shi, Y.; Pan, G.; Yim, S.C.; Yan, G.; Zhang, D. Numerical investigation of hydroelastic water-entry impact dynamics of AUVs. J. Fluids Struct.; 2019; 91, 102760. [DOI: https://dx.doi.org/10.1016/j.jfluidstructs.2019.102760]

5. Zhou, D.-H.; Shi, H.-H.; Jia, H.-X. Characteristics of the multiphase flow field with super-cavitation induced by successively fired projectiles under-water and cross-medium. J. Mech. Sci. Technol.; 2022; 36, pp. 247-258. [DOI: https://dx.doi.org/10.1007/s12206-021-1223-7]

6. Shi, Y.; Pan, G.; Yan, G.-X.; Yim, S.C.; Jiang, J. Numerical study on the cavity characteristics and impact loads of AUV water entry. Appl. Ocean Res.; 2019; 89, pp. 44-58. [DOI: https://dx.doi.org/10.1016/j.apor.2019.05.012]

7. Huang, L.; Ren, K.; Li, M.; Tuković, Ž.; Cardiff, P.; Thomas, G. Fluid-structure interaction of a large ice sheet in waves. Ocean Eng.; 2019; 182, pp. 102-111. [DOI: https://dx.doi.org/10.1016/j.oceaneng.2019.04.015]

8. Shi, F.; Xin, J.; Jin, Q. A Cartesian grid based multiphase flow model for water impact of an arbitrary complex body. Int. J. Multiph. Flow; 2019; 110, pp. 132-147. [DOI: https://dx.doi.org/10.1016/j.ijmultiphaseflow.2018.09.008]

9. Wu, Q.; Ni, B.; Bai, X.; Cui, B.; Sun, S. Experimental study on large deformation of free surface during water exit of a sphere. Ocean Eng.; 2017; 140, pp. 369-376. [DOI: https://dx.doi.org/10.1016/j.oceaneng.2017.06.009]

10. Peng, B.; Zhou, C. An approach of dynamic mesh adaptation for simulating 3-dimensional unsteady moving-immersed-boundary flows. Int. J. Numer. Methods Fluids; 2018; 87, pp. 180-201. [DOI: https://dx.doi.org/10.1002/fld.4486]

11. Xin, J.; Shi, F.; Fan, S.; Jin, Q.; Chang, X. Parametric studies on the water impact of one and twin free-falling wedges by a Cartesian grid multiphase flow model. Ocean Eng.; 2022; 249, 110854. [DOI: https://dx.doi.org/10.1016/j.oceaneng.2022.110854]

12. Osher, S.; Sethian, J.A. Fronts propagating with curvature-dependent speed: Algorithms based on Hamilton-Jacobi formulations. J. Comput. Phys.; 1988; 79, pp. 12-49. [DOI: https://dx.doi.org/10.1016/0021-9991(88)90002-2]

13. Osher, S.; Fedkiw, R.P. Level set methods: An overview and some recent results. J. Comput. Phys.; 2001; 169, pp. 463-502. [DOI: https://dx.doi.org/10.1006/jcph.2000.6636]

14. Chang, Y.-C.; Hou, T.; Merriman, B.; Osher, S. A level set formulation of Eulerian interface capturing methods for incompressible fluid flows. J. Comput. Phys.; 1996; 124, pp. 449-464. [DOI: https://dx.doi.org/10.1006/jcph.1996.0072]

15. Sussman, M.; Fatemi, E.; Smereka, P.; Osher, S. An improved level set method for incompressible two-phase flows. Comput. Fluids; 1998; 27, pp. 663-680. [DOI: https://dx.doi.org/10.1016/S0045-7930(97)00053-4]

16. Yang, C.; Mao, Z.-S. An improved level set approach to the simulation of drop and bubble motion. Chin. J. Chem. Eng. Engl. Version; 2002; 10, pp. 263-272.

17. Marchandise, E.; Remacle, J.-F. A stabilized finite element method using a discontinuous level set approach for solving two phase incompressible flows. J. Comput. Phys.; 2006; 219, pp. 780-800. [DOI: https://dx.doi.org/10.1016/j.jcp.2006.04.015]

18. Sethian, J.A.; Smereka, P. Level set methods for fluid interfaces. Annu. Rev. Fluid Mech.; 2003; 35, pp. 341-372. [DOI: https://dx.doi.org/10.1146/annurev.fluid.35.101101.161105]

19. Olsson, E.; Kreiss, G. A conservative level set method for two phase flow. J. Comput. Phys.; 2005; 210, pp. 225-246. [DOI: https://dx.doi.org/10.1016/j.jcp.2005.04.007]

20. Hirt, C.W.; Nichols, B.D. Volume of fluid (VOF) method for the dynamics of free boundaries. J. Comput. Phys.; 1981; 39, pp. 201-225. [DOI: https://dx.doi.org/10.1016/0021-9991(81)90145-5]

21. Rider, W.J.; Kothe, D.B. Reconstructing volume tracking. J. Comput. Phys.; 1998; 141, pp. 112-152. [DOI: https://dx.doi.org/10.1006/jcph.1998.5906]

22. Harvie, D.J.; Fletcher, D.F. A new volume of fluid advection algorithm: The stream scheme. J. Comput. Phys.; 2000; 162, pp. 1-32. [DOI: https://dx.doi.org/10.1006/jcph.2000.6510]

23. Cifani, P.; Michalek, W.; Priems, G.; Kuerten, J.G.; van der Geld, C.; Geurts, B.J. A comparison between the surface compression method and an interface reconstruction method for the VOF approach. Comput. Fluids; 2016; 136, pp. 421-435. [DOI: https://dx.doi.org/10.1016/j.compfluid.2016.06.026]

24. Ashgriz, N.; Poo, J. FLAIR: Flux line-segment model for advection and interface reconstruction. J. Comput. Phys.; 1991; 93, pp. 449-468. [DOI: https://dx.doi.org/10.1016/0021-9991(91)90194-P]

25. Youngs, D.L. Time-dependent multi-material flow with large fluid distortion. Numerical Methods in Fluid Dynamics; Academic Press: Cambridge, MA, USA, 1982.

26. Yang, A.-S.; Yang, J.-C.; Hong, M.-C. Droplet ejection study of a Picojet printhead. J. Micromechanics Microengineering; 2005; 16, 180. [DOI: https://dx.doi.org/10.1088/0960-1317/16/1/024]

27. Meier, M.; Yadigaroglu, G.; Smith, B.L. A novel technique for including surface tension in PLIC-VOF methods. Eur. J. Mech. -B/Fluids; 2002; 21, pp. 61-73. [DOI: https://dx.doi.org/10.1016/S0997-7546(01)01161-X]

28. Sussman, M.; Puckett, E.G. A coupled level set and volume-of-fluid method for computing 3D and axisymmetric incompressible two-phase flows. J. Comput. Phys.; 2000; 162, pp. 301-337. [DOI: https://dx.doi.org/10.1006/jcph.2000.6537]

29. Lv, X.; Zou, Q.; Reeve, D. Numerical simulation of overflow at vertical weirs using a hybrid level set/VOF method. Adv. Water Resour.; 2011; 34, pp. 1320-1334. [DOI: https://dx.doi.org/10.1016/j.advwatres.2011.06.009]

30. Wang, Z.; Yang, J.; Koo, B.; Stern, F. A coupled level set and volume-of-fluid method for sharp interface simulation of plunging breaking waves. Int. J. Multiph. Flow; 2009; 35, pp. 227-246. [DOI: https://dx.doi.org/10.1016/j.ijmultiphaseflow.2008.11.004]

31. Liu, D.; Tang, W.; Wang, J.; Xue, H.; Wang, K. Modelling of liquid sloshing using CLSVOF method and very large eddy simulation. Ocean Eng.; 2017; 129, pp. 160-176. [DOI: https://dx.doi.org/10.1016/j.oceaneng.2016.11.027]

32. Sun, D.; Tao, W. A coupled volume-of-fluid and level set (VOSET) method for computing incompressible two-phase flows. Int. J. Heat Mass Transf.; 2010; 53, pp. 645-655. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2009.10.030]

33. Zhao, L.; Shkarayev, S. Characterization of ducted contra-rotating propeller propulsions. Int. J. Micro Air Veh.; 2019; 11, 1756829319837661. [DOI: https://dx.doi.org/10.1177/1756829319837661]

34. Noh, W.F.; Woodward, P. SLIC (simple line interface calculation). Proceedings of the Fifth International Conference on Numerical Methods in Fluid Dynamics, Twente University; Enschede, The Netherlands, 28 June–2 July 1976; pp. 330-340.

35. Lopez, J.; Hernandez, J.; Gomez, P.; Faura, F. An improved PLIC-VOF method for tracking thin fluid structures in incompressible two-phase flows. J. Comput. Phys.; 2005; 208, pp. 51-74. [DOI: https://dx.doi.org/10.1016/j.jcp.2005.01.031]

36. Roenby, J.; Larsen, B.E.; Bredmose, H.; Jasak, H. A new volume-of-fluid method in OpenFOAM. Proceedings of the MARINE VI: Proceedings of the VI International Conference on Computational Methods in Marine Engineering; Nantes, France, 5–17 May 2017; pp. 266-277.

37. Eymard, R.; Gallouët, T.; Herbin, R. Finite volume methods. Handb. Numer. Anal.; 2000; 7, pp. 713-1018.

38. Weller, H. A Code Independent Notation for Finite Volume Algorithms; Technical Report TR/HGW/02 Nabla Ltd.: Port Louis, Mauritius, 2002.

39. Dubief, Y.; Delcayre, F. On coherent-vortex identification in turbulence. J. Turbul.; 2000; 1, N11. [DOI: https://dx.doi.org/10.1088/1468-5248/1/1/011]