High-speed projectile water entry super-cavitation flow numerical value prediction method

Abstract

A high-speed projectile water entry supercavitation flow numerical value prediction method comprises the following steps: dividing a flow field computational domain by adopting an overlapped grid, realizing accurate capture of a motion boundary of a projectile through dynamic cutting, and executing a grid independence test; based on the computational domain, constructing a three-phase control equation based on a VOF model; on the basis of the established control equation, a flow field is solved through a Standard k-epsilon model and a Schnerr-Sauer cavitation model, and flow field pressure is output; predicting a projectile motion track by using a six-degree-of-freedom dynamic grid; according to the system, the equipment and the medium, numerical prediction of the high-speed projectile in-water super-cavitation flow process is realized based on the method; according to the method, the accuracy of the method is verified through grid independence test and experimental data comparison (the initial speed is 75.4-142.7 m / s), and finally, a Standard k-epsilon model is selected as an optimal turbulence model; the coupling problem of phase change, turbulent flow and cross-medium flow in the high-speed water entry process is solved, and a high-precision numerical prediction tool is provided for projectile design.

Description

Technical Field

[0001] This invention belongs to the interdisciplinary field of computational fluid dynamics and high-speed projectile water entry dynamics, specifically relating to a numerical prediction method for hypercavitation flow of high-speed projectiles entering water. Background Technology

[0002] The process of high-speed projectiles entering water involves three-phase transient flow of air, water, and vapor, cavitation effects, and multi-scale turbulence characteristics. Traditional experimental methods are costly and struggle to capture microscopic flow details. Existing numerical prediction methods, such as [Xiao Haiyan, Luo Song, Zhu Zhu, et al. Study on ballistic characteristics of high-speed projectiles entering water at small angles [J]. Journal of Beijing Institute of Technology (Social Sciences Edition), 2019.], [Zhu Bangbang, Shi Honghui, Hou Jian, et al. Numerical simulation of air-carrying capacity during high-speed projectile water entry [J]. Journal of Zhejiang Sci-Tech University, 2017, 37-38(Natural Science 3): 402-408.], and [Hao Bo, Dai Hao, Lü Chao. Influence of water entry angle on the process of high-speed projectile water entry [J]. Journal of Northeastern University (Natural Science Edition), 2021.], have the following problems: 1) Turbulence models have poor adaptability to low Reynolds number conditions; 2) Insufficient accuracy in coupling cavitation phase change and multiphase flow; 3) Mesh motion algorithms struggle to handle high-speed boundary conditions. To address the aforementioned issues, there is an urgent need to develop high-precision and high-efficiency numerical prediction methods. Summary of the Invention

[0003] To overcome the shortcomings of the existing technology, the present invention aims to provide a numerical prediction method for hypercavitation flow of high-speed projectiles entering water. This method utilizes the governing equations involved in the hypercavitation flow process of high-speed projectiles entering water, including basic governing equations, turbulence models, multiphase flow models, cavitation models, and overlapping grid techniques and dynamic grid models to achieve grid movement. This allows for the acquisition of the complex flow process of high-speed projectiles entering water, characterized by phase change, strong turbulence, and cross-medium characteristics. The method is then validated and compared with corresponding projectile water entry experimental results. Finally, a model integrating the Standard k-ε turbulence model, the VOF multiphase flow model, and the Schnerr-Sauer cavitation model is selected as the best performing model.

[0004] To achieve the above objectives, the technical solution adopted by the present invention is as follows:

[0005] A numerical prediction method for supercavitation flow of high-speed projectiles entering water includes the following steps:

[0006] Step 1: Divide the flow field computation domain using overlapping meshes, capture the projectile's motion boundary through dynamic cutting, and perform mesh independence checks;

[0007] Step 2: Based on the computational domain of Step 1, construct the three-phase governing equations based on the VOF model, including the air phase, water phase, and vapor phase;

[0008] Step 3: Based on the governing equations established in Step 2, solve for the flow field and output flow field pressure using the Standard k-ε model and the Schnerr-Sauer cavitation model;

[0009] Step 4: Apply a six-degree-of-freedom moving mesh, i.e., a six-degree-of-freedom equation, to predict the trajectory of the projectile.

[0010] The specific method for step one is as follows:

[0011] 1.1 Constructing the background mesh

[0012] The main grid system of the computational domain is constructed using a Cartesian coordinate system.

[0013] 1.1.1) Computational domain planning: Extend outward from the central axis of the high-speed projectile's trajectory, with a reserved wake development area at the trailing edge, and customize the extension and reservation amount;

[0014] 1.1.2) Grid topology: The planned computational domain is meshed using a structured grid, with the basic grid unit being a regular hexahedron;

[0015] 1.1.3) Spatial discretization: A progressive refinement strategy is adopted in the axial direction (X-axis), that is, the mesh near the water surface is refined, and the projectile activity region Δx is set by a user, which gradually increases towards the far field; the basic resolution Δy = Δz is set in the lateral direction (Y / Z-axis), and hyperbolic tangent function refinement is implemented in the boundary layer region; the total number of meshes is controlled between 5×10^6 and 8×10^6 elements, and the convergence condition criterion number CFL is kept ≤1.2;

[0016] 1.1.4) Boundary conditions of the background mesh: its sides, top surface and bottom surface are all pressure outlet boundary conditions.

[0017] 1.2 Generating the foreground mesh:

[0018] 1.2.1) Construct a body-fitting motion mesh around the water-entry projectile body: adopt a semi-cylindrical-conical composite topology structure, and customize the axial length and radial radius of the projectile;

[0019] 1.2.2) The near-wall region of the projectile is meshed according to the principle that the height of the first layer mesh satisfies the dimensionless wall distance y+≈1, generating several layers of prism layer meshes, and the mesh type is structured mesh;

[0020] 1.3 Coupling in the overlapping region, i.e., three-layer element interpolation transition

[0021] 1.3.1) Set an interpolation layer on the side of the background grid:

[0022] Extrapolation layer setup: In the potential overlap area between the background and foreground meshes, including the outer boundary of the projectile's trajectory, establish three overlapping extrapolation layers:

[0023] Specifically, the mesh is generated by extrapolating the boundary layer function of the mesh generation software layer by layer along the normal of the background mesh boundary, establishing three extrapolation layers, with the interlayer expansion ratio controlled at 1.2 to 1.5.

[0024] The ratio of the outer interpolation cell size to the maximum cell size of the foreground mesh is constrained to 1:3 to 1:5 to balance the interpolation accuracy;

[0025] 1.3.2) Three receiver layers are set around the projectile surface on the foreground grid side, using a radial buffer layer topology: the interlayer spacing of the receiver layers increases geometrically; the normal of the foreground grid is orthogonal to the curvature of the projectile surface to adapt to the gradient change of the interface between the high-speed water entry shock wave and the cavitation bubble.

[0026] The specific method for step two is as follows:

[0027] Let α be the volume fraction of the air phase, water phase, and vapor phase. g α l α v The three satisfy the following relationship:

[0028] α g +α l +α v =1 (1-1)

[0029] The transient changes at the gas-liquid interface are captured by reconstructing the VOF interface using a multiphase flow model. Based on mass conservation, the continuity equation for multiphase flow is:

[0030]

[0031] In the formula, i = 1, 2, 3, and so on;

[0032] The momentum conservation equation is:

[0033]

[0034] In the formula, u i ρ is the fluid velocity component. m The density of the mixed medium, μ m μ is the dynamic viscosity of the mixed medium. t ρ is the turbulent viscosity coefficient. m and μ m The expressions are as follows:

[0035] ρ m =α g ρ g +α l ρ l +α v ρ v (1-4)

[0036] μ m =α g μ g +α l μ l +α v μ v (1-5).

[0037] The specific method for step three is as follows:

[0038] 3.1 The RANS method is used to solve the turbulence. The unsteady Navier-Stokes equations are treated as time-averaged, i.e., the instantaneous velocity u i Decomposed into time average and pulsation value u' i , time average and pulsation value u' i The sum is the actual flow:

[0039]

[0040] Substituting equation (1-6) into the momentum conservation equation, we get:

[0041]

[0042] in, This is called the Reynolds stress term, or turbulent stress term. Introducing the Boussinesq eddy viscosity assumption, we obtain the relationship between the Reynolds stress term and the velocity gradient:

[0043]

[0044] In the formula, k is the turbulent kinetic energy. μ t For turbulent viscosity, C μ ε is an empirical constant, and ε is the turbulent dissipation rate;

[0045] A two-equation turbulence model (Standard k-ε model) is adopted, namely, establishing the transport equation of turbulent kinetic energy (k equation) and the transport equation of dissipation rate (ε equation); by solving the transport equation of turbulent kinetic energy (k equation) and the transport equation of dissipation rate (ε equation), the governing equation is closed.

[0046] The transport equation (k-equation) for turbulent kinetic energy is:

[0047]

[0048] The transport equation (ε equation) for dissipation rate is:

[0049]

[0050] In the formula, μ t =c μ ρ m k 2 / ε,c μ G is an empirical coefficient; k G represents the turbulent kinetic energy caused by the time-averaged velocity gradient. b Y is the turbulent kinetic energy caused by buoyancy. M C represents the contribution of fluctuating expansion in turbulent flow to the total dissipation rate; 1ε C 2ε C 3ε σ k σ ε S is an empirical constant; k and S ε For source terms;

[0051] 3.2 Configure the VOF (Variable Interest Flow) model using Ansys Fluent software. The specific settings are as follows:

[0052] 3.2.1) Set up a multiphase flow model (VOF), set the number of phases to 3, define the primary phase as liquid water, and set the material property to ρ = 1000 kg / m³. 3 μ = 0.001003 Pa·s; Secondary Phase 1 is water vapor, with a material property of ρ = 0.5542 kg / m³. 3 μ = 1.34e-5 Pa·s; Phase 2 is air, keep the default settings;

[0053] 3.2.2) Configure the turbulence model:

[0054] The turbulence selection model adopted is Standard k-epsilon(2eqn);

[0055] The near-wall treatment model uses StandardWall Functions;

[0056] 3.2.3) Implanting the Schnerr-Sauer cavitation model:

[0057] The mass transfer mechanism is set as a gas-liquid-vapor multiphase continuous flow, and the output is the flow field pressure.

[0058] The specific method for step four is as follows:

[0059] 4.1) Based on the multiphase flow model VOF from step 3.2, define the rigid body motion properties:

[0060] 4.1.1) Characteristic length: Input projectile diameter D (unit: m);

[0061] 4.1.2) Center of Moment: Input the coordinates of the projectile's center of mass (Xc, Yc, Zc);

[0062] 4.1.3) Reference density: set as water density ρ = 1000 kg / m³ 3 ;

[0063] 4.2) Activate the six-degree-of-freedom solver; set the rigid body mass properties:

[0064] 4.2.1) Input the mass m (kg) of the projectile;

[0065] 4.2.2) Input the projectile inertia tensor Ixx / Iyy / Izz (kg·m 2 );

[0066] 4.2.3) Set the gravitational acceleration g = 9.81 m / s² 2 (Negative Z-axis direction);

[0067] 4.3) Define the projectile motion region: Select the projectile surface boundary, and select "Rigid Body" as the type; associate 6DOF parameters, and select the compiled 6DOF_Control macro; activate the translational and rotational degrees of freedom to predict and output the projectile motion trajectory.

[0068] The prediction system based on the above numerical prediction method for supercavitation flow of high-speed projectiles entering water includes:

[0069] Preprocessing module: Used in step one, it establishes a parametric geometric model of the projectile head configuration using unstructured mesh topology technology; implements an overlapping mesh partitioning strategy based on dynamic boundary conditions to achieve independent discretization of the moving body and the flow field; and ensures the analytical accuracy of near-wall turbulence by satisfying the principle that the first layer height of the near-wall mesh meets y+≈1.

[0070] Solver module: Used in steps two and three, integrating the Standard k-ε turbulence model to solve the Reynolds time-averaged NS equations, configuring enhanced wall functions to handle high-shear flows; using the VOF multiphase flow model to track the evolution of the gas-liquid interface; embedding the Schnerr-Sauer cavitation mass transport model to capture cavitation phenomena during projectile entry into water in real time;

[0071] Post-processing module: Used in step four in Ansys Fluent to set rigid body motion properties, activate the six-degree-of-freedom solver, define the projectile motion region, and predict the projectile trajectory.

[0072] A device for predicting the numerical values ​​of supercavitation flow when a high-speed projectile enters water, comprising:

[0073] Memory, used to store computer programs;

[0074] A processor is used to execute the computer program to implement the numerical prediction method for the supercavitation flow process of a high-speed projectile entering water, as described in steps 1 to 4.

[0075] A computer-readable storage medium storing a computer program that, when executed by a processor, enables numerical prediction of the hypercavitation flow process of a high-speed projectile entering water based on the method described in steps 1 to 4.

[0076] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0077] 1. Significantly improved multiphysics coupling accuracy

[0078] By constructing a homogeneous equilibrium multiphase flow model including air (αg), water (αl), and vapor (αv), and combining it with VOF interface tracking technology, the dynamic evolution of the gas-liquid-vapor three-phase interface during projectile entry into water was accurately captured. Figure 4 The prediction error of cavitation morphology is reduced by 35% compared to the traditional Mixture model. Using the Schnerr-Sauer cavitation model to describe the phase transition process, and incorporating a temperature-corrected saturated vapor pressure formula, the prediction error of the peak pressure at the instant of cavitation closure is <8%. Figure 5 This model improves accuracy by 20% compared to the traditional Singhal model.

[0079] 2. Robustness optimization of turbulence model

[0080] By comparing four turbulence models—Standard k-ε, Realizable k-ε, and SST k-ω—the Standard k-ε model was selected as the core turbulence model. Its velocity decay curve in the velocity range of 75.4–142.7 m / s showed the best agreement with experimental data, especially improving the prediction stability by 40% for low Reynolds number conditions (late-stage water entry velocity <30 m / s). The use of an eddy viscosity constraint function suppressed the overestimation of turbulent viscosity, reducing the calculation error of the projectile surface shear stress from 12% to less than 5%.

[0081] 3. Breakthrough in dynamic grid computing efficiency

[0082] Based on overlapping mesh technology Figure 1-3This method decouples the local motion domain (foreground mesh) of the projectile from the global flow field domain (background mesh), reducing mesh reconstruction time. For the same example, the computation time is reduced by 40% compared to traditional moving mesh methods. By constraining the 6DOF moving mesh to 3 degrees of freedom (translational x / y, rotation θz) using a 1 / 2 symmetric model, memory usage is reduced by 50% while maintaining projectile dynamics accuracy.

[0083] 4. Enhanced engineering applicability

[0084] Propose a grid independence test standard ( Figure 4-5 When the grid number reaches 478,800 (grid 2), the errors in cavitation morphology and axial force are both <5%, providing a basis for grid density selection for similar problems. A modular numerical prediction system was developed to support full automation of the preprocessing (geometry / mesh), solving (model parameterization), and post-processing (cone plots / curves).

[0085] 5. Interdisciplinary technological integration and innovation

[0086] For the first time, the shear stress transport mechanism of the SST k-ω model is combined with the overlapping mesh interpolation algorithm to solve the vortex shedding distortion problem caused by high-speed moving boundaries, with a vortex shedding frequency prediction error of <5%.

[0087] Real-time coupling of 6DOF motion and flow field solutions is achieved through compiled UDFs, with a time step that can be stabilized to 1e-6s, which is 30% faster than conventional explicit coupling methods.

[0088] In summary, this invention, through multi-level coupled modeling of multiphase flow, turbulence, cavitation, and dynamic mesh, has the advantage of overcoming the challenges of phase change, cross-medium, and strong transient characteristics in high-speed water ingress processes. Attached Figure Description

[0089] Figure 1 This is a schematic diagram of the computational domain of the present invention; wherein, Figure 1 (a) is a schematic diagram of the foreground mesh geometry. Figure 1 (b) is a schematic diagram of the background mesh geometry. Both are divided using structured meshes. For the foreground mesh, the projectile surface is set as a wall boundary condition, and the outer surface mesh is set as an overset boundary condition. For the background mesh, since the corresponding experiment involves launching the projectile into a water tank, the top surface is a pressure outlet boundary condition, and the sides and bottom surfaces are both wall boundary conditions. The cross-section of the semi-cylinder of the foreground mesh and the plane at the overlapping position of the background mesh are both symmetric boundary conditions.

[0090] Figure 2 This refers to the foreground mesh and boundary conditions of this invention.

[0091] Figure 3 This is the background mesh and boundary conditions of the present invention.

[0092] Figure 4 This invention provides a comparison of cavitation cloud images during water entry processes with different grid densities; among which, Figure 4 (a) is the cavitation cloud diagram when the water entry time is 0.5 ms. Figure 4 (b) is the cavitation cloud diagram when the water entry time is 1.5 ms. Figure 4 (c) is the cavitation cloud diagram when the water entry time is 2.5 ms. Figure 4 (d) is the cavitation cloud diagram when the water entry time is 3.5 ms. Figure 4 (e) is the cavitation cloud diagram when the water entry time is 4.5 ms.

[0093] Figure 5 This invention compares the axial pressure during water entry processes with different mesh densities. Figure 6 This is a comparison of experimental and numerical simulation results of projectile velocity decay curves under different initial velocity conditions. Figure 7 This is a comparison of the projectile velocity decay curves under different turbulence model examples with an initial velocity of 75.4 m / s. Detailed Implementation

[0094] The specific embodiments of the present invention will be described in detail below with reference to the accompanying drawings. It should be understood that the specific embodiments described herein are for illustration and explanation only and are not intended to limit the scope of the present disclosure.

[0095] A numerical prediction method for supercavitation flow of high-speed projectiles entering water.

[0096] A numerical prediction method for supercavitation flow of high-speed projectiles entering water includes the following steps:

[0097] Step 1: Divide the flow field computation domain using overlapping meshes, achieve accurate capture of the projectile's motion boundary through dynamic cutting, and perform mesh independence checks;

[0098] 1.1 Constructing the background mesh

[0099] The main grid system of the computational domain is constructed using a Cartesian coordinate system.

[0100] 1.1.1) Computational domain planning: Extend outward from the central axis of the high-speed projectile's trajectory, with a reserved wake development area at the trailing edge, and customize the extension and reservation amount;

[0101] 1.1.2) Grid topology: The planned computational domain is meshed using a structured grid, with the basic grid unit being a regular hexahedron;

[0102] 1.1.3) Spatial discretization: A progressive refinement strategy is adopted in the axial direction (X-axis), that is, the mesh near the water surface is refined, and the projectile activity region Δx is set by a user, which gradually increases towards the far field; the basic resolution Δy = Δz is set in the lateral direction (Y / Z-axis), and hyperbolic tangent function refinement is implemented in the boundary layer region; the total number of meshes is controlled between 5×10^6 and 8×10^6 elements, and the convergence condition criterion number CFL is kept ≤1.2;

[0103] 1.1.4) Boundary conditions of the background mesh: its sides, top surface and bottom surface are all pressure outlet boundary conditions.

[0104] 1.2 Generating the foreground mesh:

[0105] 1.2.1) Construct a body-fitting motion mesh around the water-entry projectile body: adopt a semi-cylindrical-conical composite topology structure, and customize the axial length and radial radius of the projectile;

[0106] 1.2.2) The near-wall region of the projectile is meshed according to the principle that the height of the first layer mesh satisfies the dimensionless wall distance y+≈1, generating several layers of prism layer meshes, and the mesh type is structured mesh;

[0107] 1.3 Coupling in the overlapping region, i.e., three-layer element interpolation transition

[0108] 1.3.1) Set an interpolation layer on the side of the background grid:

[0109] Extrapolation layer setup: In the potential overlap area between the background and foreground meshes, including the outer boundary of the projectile's trajectory, establish three overlapping extrapolation layers:

[0110] Specifically, the mesh is generated by extrapolating the boundary layer function of the mesh generation software layer by layer along the normal of the background mesh boundary, establishing three extrapolation layers, with the interlayer expansion ratio controlled at 1.2 to 1.5.

[0111] The ratio of the outer interpolation cell size to the maximum cell size of the foreground mesh is constrained to 1:3 to 1:5 to balance the interpolation accuracy;

[0112] 1.3.2) Three receiver layers are set around the projectile surface on the foreground grid side, using a radial buffer layer topology: the interlayer spacing of the receiver layers increases geometrically; the normal of the foreground grid is orthogonal to the curvature of the projectile surface to adapt to the gradient change of the interface between the high-speed water entry shock wave and the cavitation bubble.

[0113] Step 2: Based on the computational domain of Step 1, construct the three-phase governing equations based on the VOF model, including the air phase, water phase, and vapor phase;

[0114] Let α be the volume fraction of the air phase, water phase, and vapor phase. g α l α v The three satisfy the following relationship:

[0115] α g +α l +α v =1 (1-1)

[0116] The transient changes at the gas-liquid interface are captured by reconstructing the VOF interface using a multiphase flow model. Based on mass conservation, the continuity equation for multiphase flow is:

[0117]

[0118] In the formula, i = 1, 2, 3, and so on;

[0119] The momentum conservation equation is:

[0120]

[0121] In the formula, u i ρ is the fluid velocity component. m The density of the mixed medium, μ m μ is the dynamic viscosity of the mixed medium. t ρ is the turbulent viscosity coefficient. m and μ m The expressions are as follows:

[0122] ρ m =α g ρ g +α l ρ l +α v ρ v (1-4)

[0123] μ m =α g μ g +α l μ l +α v μ v (1-5).

[0124] Step 3: Based on the governing equations established in Step 2, solve for the flow field and output flow field pressure using the Standard k-ε model and the Schnerr-Sauer cavitation model;

[0125] 3.1 The RANS method is used to solve the turbulence. The unsteady Navier-Stokes equations are treated as time-averaged, i.e., the instantaneous velocity u i Decomposed into time average and pulsation value u' i , time average and pulsation value u' i The sum is the actual flow:

[0126]

[0127] Substituting equation (1-6) into the momentum conservation equation, we get:

[0128]

[0129] in, This is called the Reynolds stress term, or turbulent stress term. Introducing the Boussinesq eddy viscosity assumption, we obtain the relationship between the Reynolds stress term and the velocity gradient:

[0130]

[0131] In the formula, k is the turbulent kinetic energy. μ t For turbulent viscosity, C μ ε is an empirical constant, and ε is the turbulent dissipation rate;

[0132] A two-equation turbulence model (Standard k-ε model) is adopted, namely, establishing the transport equation of turbulent kinetic energy (k equation) and the transport equation of dissipation rate (ε equation); by solving the transport equation of turbulent kinetic energy (k equation) and the transport equation of dissipation rate (ε equation), the governing equation is closed.

[0133] The transport equation (k-equation) for turbulent kinetic energy is:

[0134]

[0135] The transport equation (ε equation) for dissipation rate is:

[0136]

[0137] In the formula, μ t =c μ ρ m k 2 / ε,c μ G is an empirical coefficient; k G represents the turbulent kinetic energy caused by the time-averaged velocity gradient. b Y is the turbulent kinetic energy caused by buoyancy. M C represents the contribution of fluctuating expansion in turbulent flow to the total dissipation rate; 1ε C 2ε C 3ε σ k σ εS is an empirical constant; k and S ε For source terms;

[0138] 3.2 Configure the VOF (Variable Interest Flow) model using Ansys Fluent software. The specific settings are as follows:

[0139] 3.2.1) Set up a multiphase flow model (VOF), set the number of phases to 3, define the primary phase as liquid water, and set the material property to ρ = 1000 kg / m³. 3 μ = 0.001003 Pa·s; Secondary Phase 1 is water vapor, with a material property of ρ = 0.5542 kg / m³. 3 μ = 1.34e-5 Pa·s; Phase 2 is air, keep the default settings;

[0140] 3.2.2) Configure the turbulence model:

[0141] The turbulence selection model adopted is Standard k-epsilon(2eqn);

[0142] The near-wall treatment model uses StandardWall Functions;

[0143] 3.2.3) Implanting the Schnerr-Sauer cavitation model:

[0144] The mass transfer mechanism is set as a gas-liquid-vapor multiphase continuous flow, and the output is the flow field pressure. Step 4: Apply a six-degree-of-freedom moving mesh, i.e., six-degree-of-freedom equations, to predict the projectile trajectory:

[0145] 4.1) Based on the multiphase flow model VOF from step 3.2, define the rigid body motion properties:

[0146] 4.1.1) Characteristic length: Input projectile diameter D (unit: m);

[0147] 4.1.2) Center of Moment: Input the coordinates of the projectile's center of mass (Xc, Yc, Zc);

[0148] 4.1.3) Reference density: set as water density ρ = 1000 kg / m³ 3 ;

[0149] 4.2) Activate the six-degree-of-freedom solver; set the rigid body mass properties:

[0150] 4.2.1) Input the mass m (kg) of the projectile;

[0151] 4.2.2) Input the projectile inertia tensor Ixx / Iyy / Izz (kg·m 2 );

[0152] 4.2.3) Set the gravitational acceleration g = 9.81 m / s² 2 (Negative Z-axis direction);

[0153] 4.3) Define the projectile motion region:

[0154] After selecting the projectile surface boundary, choose "RigidBody" as the type; associate the 6DOF parameters and select the compiled 6DOF_Control macro; activate the translational and rotational degrees of freedom to predict the projectile's trajectory.

[0155] Simulation experiment: Taking the initial velocity of 95.4 m / s as an example:

[0156] Introduction to Geometric Models and Meshes

[0157] like Figure 1 As shown, the geometry is a 1 / 2 symmetric model. The projectile is a solid cylinder (1 / 2 is a semi-cylinder). The projectile diameter D = 12.65 mm, length L = 25.4 mm, mass m = 23.5 g, and density ρ = 7850 kg / m³. 3 The foreground grid region is a semi-cylinder with a diameter of 50 mm and a height of 60 mm; the background grid region is a cuboid area with a length of 360 mm, a width of 180 mm, and a height of 650 mm. The initial velocities of the projectiles entering the water are 75.4 m / s, 95.4 m / s, 118.8 m / s, and 142.7 m / s, respectively, and the water surface height is set to 600 mm.

[0158] Figure 2 and 3 These are schematic diagrams of the foreground and background meshes of an overlapping mesh, along with their corresponding boundary conditions. As shown, both are meshed using structured meshing. For the foreground mesh, the projectile surface is set as a wall boundary condition, and the outer surface mesh is set as an overset boundary condition. For the background mesh, since the experiment involves launching the projectile into a water tank, the top surface is a pressure outlet boundary condition, while the sides and bottom are wall boundary conditions. The cross-section of the semi-cylindrical foreground mesh and the plane where it overlaps with the background mesh are both symmetrical boundary conditions. Based on overset meshing technology... Figure 1-3 This method decouples the local motion domain (foreground mesh) of the projectile from the global flow field domain (background mesh), reducing mesh reconstruction time. For the same example, the computation time is reduced by 40% compared to traditional moving mesh methods. By constraining the 6DOF moving mesh to 3 degrees of freedom (translational x / y, rotation θz) using a 1 / 2 symmetric model, memory usage is reduced by 50% while maintaining projectile dynamics accuracy.

[0159] Mesh independence test

[0160] To eliminate the influence of grid density on the numerical prediction results, this paper refines the grid for the projectile surface, trajectory, and water surface regions, dividing the data into three sets of grids with grid numbers of 236,900 (grid 1), 478,800 (grid 2), and 752,300 (grid 3), respectively. With an initial water entry velocity of 95.4 m / s, by comparing the cavitation contour maps of the three regions at corresponding underwater moments and the magnitude of the axial pressure on the projectile surface, a suitable grid is ultimately selected to obtain a grid-independent solution.

[0161] like Figure 4 As shown, observing the cavitation morphology at five moments from Δt = 0.5 ms to 4.5 ms after the projectile begins its motion, it can be seen that from the initial expulsion stage to the cavitation convergence stage and then to the cavitation closure stage, the difference between grid 2 and grid 3 is small, while grid 1 differs significantly from the former two. Figure 5 As shown, by comparing the magnitude of the axial force on the projectile throughout the entire water entry process, it can be seen that grid 1 differs significantly from grids 2 and 3 in the water entry stage and the final cavitation collapse stage. Therefore, grid 2 (478,800 grids) was ultimately chosen for subsequent calculations.

[0162] Since the high-speed projectile's entry into water is a physical process with high turbulence intensity, the choice of turbulence model will significantly affect the calculation results when using the RANS method for numerical prediction. Furthermore, considering that the Reynolds number decreases after the projectile velocity decays to a certain level, the chosen turbulence model should ideally account for this flow state. In this section, using experimentally set initial projectile velocities of 75.4 m / s, 118.8 m / s, and 142.7 m / s as references, four turbulence models—Standard k-ε, Realizable k-ε, SST k-ω, and SST k-ω with low Reynolds number correction—are combined and compared to determine the final complete numerical model.

[0163] First, the velocity decay curves of the experiment were compared with those of the Standard k-ε model and the Realizable k-ε model under three initial water entry velocities. Figure 6 As shown, the Standard k-ε model agrees well with experimental results under different initial velocities; however, the Realizable k-ε model deviates from experimental results from the moment of entry into the water at an initial velocity of 75.4 m / s. This indicates that for high Reynolds number models, the focus of comparison should be on lower velocity conditions. In other words, the turbulence model that performs best under low-velocity conditions should be selected. Figure 7As shown, by comparing the four turbulence models under an initial velocity of 75.4 m / s, it can be seen that the Standard k-ε model has the smallest discrepancy with the experimental results. In summary, the Standard k-ε turbulence model is selected for this numerical model. The velocity decay curve in the velocity range of 75.4–142.7 m / s shows the best agreement with the experimental data, especially improving the prediction stability by 40% for low Reynolds number conditions (late-stage water entry velocity <30 m / s). The use of an eddy viscosity constraint function suppresses the overestimation of turbulent viscosity, reducing the calculation error of the projectile surface shear stress from 12% to less than 5%.

[0164] The preferred embodiments of the present invention have been described in detail above with reference to the accompanying drawings. However, the present invention is not limited to the specific details of the above embodiments. Within the scope of the technical concept of the present invention, various simple modifications can be made to the technical solution of the present invention, and these simple modifications all fall within the protection scope of the present invention.

Claims

1. A numerical prediction method for supercavitation flow of high-speed projectiles entering water, specifically including the following steps: Step 1: Divide the flow field computation domain using overlapping meshes, capture the projectile's motion boundary through dynamic cutting, and perform mesh independence checks; Step 2: Based on the computational domain of Step 1, construct the three-phase governing equations based on the VOF model, including the air phase, water phase, and vapor phase; Step 3: Based on the governing equations established in Step 2, solve for the flow field and output flow field pressure using the Standard k-ε model and the Schnerr-Sauer cavitation model; Step 4: Apply a six-degree-of-freedom moving mesh, i.e., a six-degree-of-freedom equation, to predict the trajectory of the projectile.

2. The numerical prediction method for supercavitation flow of high-speed projectiles entering water according to claim 1, characterized in that, The specific method for step one is as follows: 1.1 Constructing the background mesh The main grid system of the computational domain is constructed using a Cartesian coordinate system. 1.1.1) Computational domain planning: Extend outward from the central axis of the high-speed projectile's trajectory, with a reserved wake development area at the trailing edge, and customize the extension and reservation amount; 1.1.2) Grid topology: The planned computational domain is meshed using a structured grid, with the basic grid unit being a regular hexahedron; 1.1.3) Spatial discretization: A progressive refinement strategy is adopted in the axial direction (X-axis), that is, the mesh near the water surface is refined, and the projectile activity region Δx is set by a user, which gradually increases towards the far field; the basic resolution Δy = Δz is set in the lateral direction (Y / Z-axis), and hyperbolic tangent function refinement is implemented in the boundary layer region; the total number of meshes is controlled between 5×10^6 and 8×10^6 elements, and the convergence condition criterion number CFL is kept ≤1.2; 1.1.4) Boundary conditions of the background mesh: its sides, top surface and bottom surface are all pressure outlet boundary conditions. 1.2 Generating the foreground mesh: 1.2.1) Construct a body-fitting motion mesh around the water-entry projectile body: adopt a semi-cylindrical-conical composite topology structure, and customize the axial length and radial radius of the projectile; 1.2.2) The near-wall region of the projectile is meshed according to the principle that the height of the first layer mesh satisfies the dimensionless wall distance y+≈1, generating several layers of prism layer meshes, and the mesh type is structured mesh; 1.3 Coupling in the overlapping region, i.e., three-layer element interpolation transition 1.3.1) Set an interpolation layer on the side of the background grid: Extrapolation layer setup: In the potential overlap area between the background and foreground meshes, including the outer boundary of the projectile's trajectory, establish three overlapping extrapolation layers: Specifically, the mesh is generated by extrapolating the boundary layer function of the mesh generation software layer by layer along the normal of the background mesh boundary, establishing three extrapolation layers, with the interlayer expansion ratio controlled at 1.2 to 1.

5. The ratio of the outer interpolation cell size to the maximum cell size of the foreground mesh is constrained to 1:3 to 1:5 to balance the interpolation accuracy; 1.3.2) Three receiver layers are set around the projectile surface on the foreground grid side, using a radial buffer layer topology: the interlayer spacing of the receiver layers increases geometrically; the normal of the foreground grid is orthogonal to the curvature of the projectile surface to adapt to the gradient change of the interface between the high-speed water entry shock wave and the cavitation bubble.

3. The numerical prediction method for supercavitation flow of high-speed projectiles entering water according to claim 1, characterized in that, The specific method for step two is as follows: Let α be the volume fraction of the air phase, water phase, and vapor phase. g α l α v The three satisfy the following relationship: α g +α l +α v ＝1 (1-1) The transient changes at the gas-liquid interface are captured by reconstructing the VOF interface using a multiphase flow model. Based on mass conservation, the continuity equation for multiphase flow is: In the formula, i = 1, 2, 3, and so on; The momentum conservation equation is: In the formula, u i ρ is the fluid velocity component. m The density of the mixed medium, μ m μ is the dynamic viscosity of the mixed medium. t ρ is the turbulent viscosity coefficient. m and μ m The expressions are as follows: r m =a g r g +a l r l +a v r v (1-4) m m =a g m g +a l m l +a v m v (1-5).

4. The numerical prediction method for supercavitation flow of high-speed projectiles entering water according to claim 1, characterized in that, The specific method for step three is as follows: 3.1 The RANS method is used to solve the turbulence. The unsteady Navier-Stokes equations are treated as time-averaged, i.e., the instantaneous velocity u i Decomposed into time average and pulsation value u' i , time average and pulsation value u' i The sum is the actual flow: Substituting equation (1-6) into the momentum conservation equation, we get: in, This is called the Reynolds stress term, or turbulent stress term. Introducing the Boussinesq eddy viscosity assumption, we obtain the relationship between the Reynolds stress term and the velocity gradient: In the formula, k is the turbulent kinetic energy. μ t For turbulent viscosity, C μ ε is an empirical constant, and ε is the turbulent dissipation rate; A two-equation turbulence model (Standard k-ε model) is adopted, namely, establishing the transport equation of turbulent kinetic energy (k equation) and the transport equation of dissipation rate (ε equation); by solving the transport equation of turbulent kinetic energy (k equation) and the transport equation of dissipation rate (ε equation), the governing equation is closed. The transport equation (k-equation) for turbulent kinetic energy is: The transport equation (ε equation) for dissipation rate is: In the formula, μ t =c μ ρ m k 2 / ε,c μ G is an empirical coefficient; k G represents the turbulent kinetic energy caused by the time-averaged velocity gradient. b Y is the turbulent kinetic energy caused by buoyancy. M C represents the contribution of fluctuating expansion in turbulent flow to the total dissipation rate; 1ε C 2ε C 3ε σ k σ ε S is an empirical constant; k and S ε For source terms; 3.2 Configure the VOF (Variable Interest Flow) model using Ansys Fluent software. The specific settings are as follows: 3.2.1) Set up a multiphase flow model (VOF), set the number of phases to 3, define the primary phase as liquid water, and set the material property to ρ = 1000 kg / m³. 3 μ = 0.001003 Pa·s; Secondary Phase 1 is water vapor, with a material property of ρ = 0.5542 kg / m³. 3 μ = 1.34e-5 Pa·s; Phase 2 is air, keep the default settings; 3.2.2) Configure the turbulence model: The turbulence selection model adopted is Standard k-epsilon(2eqn); The near-wall treatment model uses StandardWall Functions; 3.2.3) Implanting the Schnerr-Sauer cavitation model: The mass transfer mechanism is set as a gas-liquid-vapor multiphase continuous flow, and the output is the flow field pressure.

5. The numerical prediction method for supercavitation flow of high-speed projectiles entering water according to claim 1, characterized in that, The specific method for step four is as follows: 4.1) Based on the multiphase flow model VOF from step 3.2, define the rigid body motion properties: 4.1.1) Characteristic length: Input projectile diameter D (unit: m); 4.1.2) Center of Moment: Input the coordinates of the projectile's center of mass (Xc, Yc, Zc); 4.1.3) Reference density: set as water density ρ = 1000 kg / m³ 3 ; 4.2) Activate the six-degree-of-freedom solver; 4.2.1) Input the mass m (kg) of the projectile; 4.2.2) Input the projectile inertia tensor Ixx / Iyy / Izz (kg·m 2 ); 4.2.3) Set the gravitational acceleration g = 9.81 m / s² 2 (Negative Z-axis direction); 4.3) Define the projectile motion region: Select the projectile surface boundary, and select "Rigid Body" as the type; associate 6DOF parameters, and select the compiled 6DOF_Control macro; activate the translational and rotational degrees of freedom to predict the projectile's trajectory.

6. A prediction system based on the numerical prediction method for supercavitation flow of high-speed projectiles entering water according to any one of claims 1 to 5, characterized in that, include: Preprocessing module: Used in step one to establish a parametric geometric model of the projectile head configuration using unstructured mesh topology technology; An overlapping grid partitioning strategy is implemented based on dynamic boundary conditions to achieve independent discretization of the moving body and the flow field; the analytical accuracy of near-wall turbulence is guaranteed by ensuring that the height of the first layer of the near-wall grid satisfies the principle of y+≈1. Solver module: Used in steps two and three, integrating the Standard k-ε turbulence model to solve the Reynolds time-averaged NS equations, configuring enhanced wall functions to handle high-shear flows; using the VOF multiphase flow model to track the evolution of the gas-liquid interface; embedding the Schnerr-Sauer cavitation mass transport model to capture cavitation phenomena during projectile entry into water in real time; Post-processing module: Used in step four in Ansys Fluent to set rigid body motion properties, activate the six-degree-of-freedom solver, define the projectile motion region, and predict the projectile trajectory.

7. A device for numerical prediction of supercavitation flow in water following high-speed projectile entry, characterized in that, include: Memory, used to store computer programs; A processor is used to execute the computer program to implement the numerical prediction method for the supercavitation flow process of a high-speed projectile entering water as described in claims 1 to 6.

8. A computer-readable storage medium storing a computer program, characterized in that, When the computer program is executed by the processor, it can perform numerical prediction of the supercavitation flow process of high-speed projectiles entering water based on the method described in claims 1 to 6.

Images