Two-dimensional near-free-surface underwater explosion cavitation forecasting method

Abstract

The invention discloses a two-dimensional near-free-surface underwater explosion cavitation forecasting method, and relates to the technical field of explosion.The method comprises the steps that a mixed fluid is regarded as a mixture of multiple fluid media, and a multiphase compressible fluid motion control equation related to state parameters of the mixed fluid before the mixed fluid reaches an equilibrium state and the mass fraction of the fluid media of all phases is constructed; then constructing a thermodynamic state equation of any k-phase fluid medium based on an NASG equation, and constructing a cavitation phase change equation after the mixed fluid completes phase change relaxation to reach an equilibrium state by considering the phase change process of a liquid phase and a steam phase when a flow field changes; and solving an equation set formed by a multiphase compressible fluid motion control equation, a thermodynamic state equation of each phase of fluid medium and a cavitation phase change equation in each cycle step to obtain a numerical forecasting result. According to the method, cavitation is regarded as the continuous accumulation process of the steam phase content in the mixed fluid, the physical essence is better met, and therefore the forecasting result is more accurate.

Description

Technical Field

[0001] This application relates to the field of explosion technology, and in particular to a method for predicting cavitation in two-dimensional underwater explosion near a free surface. Background Technique

[0002] Underwater explosion near the water surface poses a great threat to surface ships. During the underwater explosion process, in addition to transmitting shock wave loads to the surrounding fluid, the movement of the explosion bubble also exerts pulsating loads on the surrounding fluid. This causes the surface ship to not only bear the action of the shock wave load during the underwater explosion near the water surface, but also bear the action of the bubble pulsating load. In addition, when the shock wave propagates to the free surface, transmission and reflection phenomena will occur on both sides of the free surface. Among them, the reflected one is a rarefaction wave, which will cause the pressure in the water to decrease and is prone to cavitation in the water. The cavitation effect caused by the existence of the free surface will also have an important impact on the deformation and damage of the ship structure. Therefore, studying the cavitation near the free surface caused by underwater explosion is of great significance for predicting the deformation and damage of the ship structure affected by underwater explosion near the water surface.

[0003] The cavitation in underwater explosion near the free surface has the following several remarkable characteristics: First, the cavitation area is large, the energy accumulated during the collapse process is high, and the destructive ability generated is large; Second, although the evolution process of cavitation is longer than the action time of the early shock wave, it is much smaller than the explosion bubble pulsation period and still belongs to a typical transient action process. The shock dynamics behavior during the cavitation collapse process cannot be ignored; Finally, the process of cavitation near the free surface involves multi-fluid interaction and its complex coupling process with the structure, which is a typical multi-substance and multi-interface motion phenomenon and is very difficult to study.

[0004] Kedrinski, Kleine, Cui, etc. recorded and observed the cavitation phenomenon in underwater explosion near the water surface, providing important experimental data for the study of underwater explosion cavitation. However, the experimental conditions are demanding, the testing difficulty is large, and the obtained data is very limited, making it difficult to comprehensively display the fine characteristics of the cavitation domain. Therefore, numerical simulation research on the cavitation movement in underwater explosion has certain advantages compared with experimental methods.

[0005] The numerical simulation of underwater explosion cavitation has mainly gone through two important development stages. First, in the 1980s and 1990s of the last century, cavitation models based on the acoustic element hypothesis were proposed, such as Felippa, Sprague, etc. These models were mainly based on the potential flow theory model and used linear acoustic fluid elements to describe the motion of cavitation bubbles. However, this model has simple calculations and is only applicable to the underwater weak shock wave environment. In the past two decades, with the rapid development of compressible fluid dynamics, cavitation motion has received increasing attention. It is considered that in the process of underwater explosion, the liquid and vapor phases of the fluid inside the cavitation region and the surrounding fluid environment should all be regarded as compressible fluids, and the Euler equation or the NS equation is used to describe the basic motion of the fluid. Thus, one-fluid models represented by the Cut-off model, the isentropic model, the Schmidt model, and the Modified Schmidt model have been developed. This is the commonly used model for most cavitation research in the field of underwater explosion at present. However, this model regards cavitation as a suddenly occurring phenomenon. When the flow field pressure is lower than the saturation pressure, cavitation is considered to have occurred, and the tracking of cavitation evolution is realized by constructing different state equations satisfied in the cavitation region. This kind of numerical simulation can only be applied to simple ideal working conditions, its application range is greatly limited, it is far from the actual engineering background, and the accuracy of the numerical simulation of underwater explosion cavitation motion is not high. Summary of the Invention

[0006] In view of the above problems and technical requirements, the present application proposes a two-dimensional near-free-surface underwater explosion cavitation prediction method. The technical solution of the present application is as follows:

[0007] A two-dimensional near-free-surface underwater explosion cavitation prediction method, the two-dimensional near-free-surface underwater explosion cavitation prediction method includes:

[0008] Construct a computational domain for the two-dimensional near-free-surface underwater explosion scenario and perform discrete grid division on the computational domain. The computational domain includes an air domain above the free surface, a water domain below the free surface, and an explosion bubble domain at the near-free-surface position in the water domain;

[0009] Construct the control equations of multi-phase compressible fluid motion for the state parameters of the mixed fluid in the computational domain before reaching the equilibrium state and the mass fraction Y k of any k-th phase fluid medium in the mixed fluid. The control equations of multi-phase compressible fluid motion include the mass conservation equation, the momentum conservation equation, the energy conservation equation, and the phase change equation. k = 1 represents the liquid phase, and k = 2 represents the vapor phase;

[0010] Based on the NASG equation, construct the thermodynamic state equation of any k-th phase fluid medium. The state parameters of any k-th phase fluid medium include pressure p k 、temperature T k 、Gibbs free energy g kThe sound speed c k and is related to the specific volume υ k and the internal energy per unit mass e k by coupling;

[0011] According to the characteristic that g1 = g2 when the mixed fluid is isothermal and isobaric and reaches the equilibrium state, a cavitation phase change equation after the mixed fluid completes the phase change relaxation and reaches the equilibrium state is constructed;

[0012] The numerical prediction results are obtained by solving the system of equations composed of the motion control equation of the multiphase compressible fluid, the thermodynamic state equation of each phase fluid medium, and the cavitation phase change equation at each cycle step.

[0013] A further technical solution thereof is that the cavitation phase change equation after the mixed fluid completes the phase change relaxation and reaches the equilibrium state is:

[0014]

[0015] where N is the number of components of the fluid medium contained in the mixed fluid, any fluid medium represented by k≥3 is a gas-phase medium and does not participate in the phase change, e is the internal energy per unit mass of the mixed fluid before reaching the equilibrium state, υ is the specific volume of the mixed fluid before reaching the equilibrium state, p * is the pressure of the mixed fluid after reaching the equilibrium state, T * is the temperature of the mixed fluid after reaching the equilibrium state, W k is the molar mass of the k-th phase fluid medium, Y k is the mass fraction of the k-th phase fluid medium in the mixed fluid before reaching the equilibrium state, is the mass fraction of the k-th phase fluid medium in the mixed fluid after reaching the equilibrium state, υ k (T * , p * ) is the specific volume of the k-th phase fluid medium in the mixed fluid after reaching the equilibrium state and is related to p * and T * related, e k (T * , p * ) is the internal energy per unit mass of the k-th phase fluid medium in the mixed fluid after reaching the equilibrium state and is related to p * and T * related, υ k (T * , p * ) and e k (T * , p * )'s relationship with p * and T * is determined according to the thermodynamic state equation of the k-th phase fluid medium, e k (T * , p *) is the internal energy per unit mass of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state and is related to p * and T * where p partial represents the partial pressure of the vapor phase after the mixed fluid reaches the equilibrium state, and p sat (T * ) represents the saturation pressure corresponding to T * .

[0016] A further technical solution thereof is that the thermodynamic state equation of any k-th phase fluid medium includes:

[0017]

[0018] wherein, the pressure p k (υ k , e k ) is related to υ k and e k , the temperature T k (p k , υ k ) is related to p k and υ k , the Gibbs free energy g k (p k , T k ) is related to p k and T k , the speed of sound c k (p k , υ k ) is related to p k and υ k ; γ k , q k , b k , q′ k , C v,k are all state equation parameters obtained by pre-fitting.

[0019] A further technical solution thereof is that the state equation parameters of the liquid phase represented by k = 1 are γ1 = 1.19, q1 = -1177788 J / kg, b1 = 6.61×10 -4 m 3 / kg, q1′ = 0 J / kg / K, C v,1 = 3610 J / kg / K, W1 = 18 g / mol;

[0020] The state equation parameters of the vapor phase represented by k = 2 are γ2 = 1.47, q2 = 2077616 J / kg, b2 = 0 m 3 / kg, q2′ = 2077616 J / kg / K, C v,2 = 955 J / kg / K, W2 = 18 g / mol;

[0021] The state equation parameters of air represented by k = 3 are γ3 = 1.4, q3 = 0 J / kg, b3 = 0 m 3 / kg, q3′ = 0 J / kg / K, C v,3 = 719 J / kg / K, W3 = 29 g / mol;

[0022] The state equation parameters of the explosive gas represented by k = 4 are γ4 = 1.8, q4 = 0 J / kg, b4 = 0 m 3 / kg, q4′ = 0 J / kg / K, C v,4 = 695 J / kg / K, W4 = 290 g / mol.

[0023] Its further technical solution is to calculate T * The corresponding saturation pressure p sat (T * ) The method includes:

[0024] Taking T * As the saturation temperature T sat Substitute it into the following relational formula to calculate the corresponding saturation pressure p sat As T * The corresponding saturation pressure p sat (T * ):

[0025]

[0026] Among them, A, B, C, and D are calculated according to the state equation parameters of the liquid phase and the vapor phase.

[0027] Its further technical solution is,[[]] C p,1 = 4285 J / kg / K, C p,2 = 1401 J / kg / K.

[0028] Its further technical solution is that the numerical prediction results obtained in each cycle step include:

[0029] Combining the thermodynamic state equations of the fluid media in each phase of the mixed fluid to solve the motion control equations of the multiphase compressible fluid to obtain the specific volume υ, the internal energy per unit mass e, the pressure p, the temperature T, and the mass fraction Y of the fluid media in each phase before the mixed fluid reaches the equilibrium state kAnd the state parameters of each phase of the fluid medium before the mixed fluid reaches the equilibrium state are substituted into the cavitation phase change equation and solved by an iterative method to obtain the pressure p after the mixed fluid reaches the equilibrium state * , temperature T * and the mass fractions of each phase of the fluid medium

[0030] A further technical solution thereof is that the numerical prediction results obtained by solving in each cycle step further include:

[0031] According to Calculate the volume fraction z of the k-th phase fluid medium in the mixed fluid k , when the volume fraction z2 of the vapor phase in the mixed fluid reaches a predetermined threshold, the corresponding grid cell is classified as a cavitation grid, and a cavitation domain composed of all cavitation grids in each cycle step is obtained; wherein, ρ is the density of the mixed fluid before reaching the equilibrium state, is the density of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state and

[0032] A further technical solution thereof is that solving the multi-phase compressible fluid motion control equation in combination with the thermodynamic state equations of each phase of the fluid medium in the mixed fluid includes:

[0033] Use the thermodynamic state equation of the k-th phase fluid medium to solve the sound speed c of the k-th phase fluid medium k (p k , υ k ), and use the second-order MUSCL-Hancock method and the HLLC approximate Riemann solver to solve the multi-phase compressible fluid motion control equation in combination with the obtained sound speed.

[0034] A further technical solution thereof is that the multi-phase compressible fluid motion control equation is:

[0035]

[0036] Among them, ρ is the density of the mixed fluid before the equilibrium state and υ = 1 / ρ, u is the velocity vector of the mixed fluid, E0 is the total energy per unit mass of the mixed fluid and E0 = e + 0.5u×u, I is the unit matrix, is the gradient operator, and t represents time.

[0037] The beneficial technical effects of this application are:

[0038] The present application discloses a method for predicting cavitation in two-dimensional underwater explosions near a free surface. In this prediction method, the mixed fluid is regarded as a mixture of multiple fluid media. When the flow field changes, the liquid phase and the vapor phase will have convective motion and a phase change process accompanied by mass and heat exchange, and finally reach equilibrium. Cavitation is regarded as a process of continuous accumulation of the vapor phase content in the mixed fluid, so as to realize the numerical prediction of cavitation in the two-dimensional underwater explosion scenario near the free surface. The cavitation evolution process described by the numerical prediction model constructed by this method is more in line with the physical essence. Therefore, the prediction result of cavitation in the two-dimensional underwater explosion near the free surface is more accurate, and the prediction efficiency is also higher, which can be well applied to engineering applications. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] Figure 1 is a schematic flow chart of the method for predicting cavitation in two-dimensional underwater explosions near a free surface according to an embodiment of the present application.

[0040] Figure 2 is a schematic diagram of the computational domain in an example.

[0041] Figure 3 is an evolution diagram of the prediction result of the density of the mixed fluid obtained by using the prediction method of the present application in an example.

[0042] Figure 4 is Figure 3 an evolution diagram of the prediction result of the pressure of the mixed fluid obtained by using the prediction method of the present application in an example.

[0043] Figure 5 is Figure 3 an evolution diagram of the prediction result of the volume fraction of the vapor phase in the mixed fluid obtained by using the prediction method of the present application in an example.

[0044] Figure 6 is Figure 3 an evolution diagram of the prediction result of the volume fraction of the cavitation domain obtained by using the prediction method of the present application in an example. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0045] The following further describes the specific embodiments of the present application with reference to the drawings.

[0046] The present application discloses a method for predicting cavitation in two-dimensional underwater explosions near a free surface. Please refer to Figure 1 the flow chart shown. This method for predicting cavitation in two-dimensional underwater explosions near a free surface includes:

[0047] Step S1, constructing a computational domain for the two-dimensional underwater explosion scenario near the free surface and performing discrete grid division on the computational domain.

[0048] In the two-dimensional underwater explosion scenario near the free surface, the constructed computational domain is as shown in Figure 2As shown in the figure, the computational domain includes an air domain 210 above the free surface, a water domain 220 below the free surface, and an explosive bubble domain 230 at a position near the free surface in the water domain. The position near the free surface here means that the explosive bubble domain 230 is located below the free surface and the distance from the free surface does not exceed a distance threshold.

[0049] The fluid in the computational domain is a mixed fluid, which is obtained by mixing N different fluid media. In the two-dimensional underwater explosion scenario near the free surface, the mixed fluid includes four fluid media. In this application, for any k-th phase fluid medium in the mixed fluid, k = 1 represents the liquid phase, k = 2 represents the vapor phase, k = 3 represents air (non-condensable), k = 4 represents the explosive gas. Among them, the liquid phase and the vapor phase will undergo phase changes, while air and the explosive gas do not undergo phase changes. There are also other fluid media with smaller volume fractions in the mixed fluid that do not participate in phase changes and are negligible within the error range.

[0050] Under the initial conditions: (1) The mixed fluid in the air domain 210 is composed of air and the vapor phase, without the liquid phase and the explosive gas. The volume fraction of air in the mixed fluid in the air domain 210 is much larger than that of the vapor phase. That is, the mixed fluid in the air domain contains mostly air and a small amount of vapor phase medium, and the mass fractions of the liquid phase and the explosive gas are both 0. The mixed fluid in the air domain 210 is in a superheated state. (2) The mixed fluid in the water domain 220 is composed of the liquid phase, air, and the vapor phase, without the explosive gas. The volume fraction of the liquid phase in the mixed fluid in the water domain 220 is much larger than that of air and the vapor phase. That is, the mixed fluid in the water domain contains mostly liquid phase medium and a small amount of air and vapor phase. The mixed fluid in the water domain 220 is in a saturated state. (3) The mixed fluid in the explosive bubble domain 230 is composed of the explosive gas and the vapor phase, without the liquid phase and air. The volume fraction of the explosive gas in the mixed fluid in the explosive bubble domain 230 is much larger than that of the vapor phase. That is, the mixed fluid in the explosive bubble domain 230 contains mostly the explosive gas and a small amount of vapor phase. The mixed fluid in the explosive bubble domain 230 is in a superheated state. When giving the initial conditions, it is necessary to give the mass fraction / volume fraction of the fluid media contained in the mixed fluid in the air domain 210, the water domain 220, and the explosive bubble domain 230 respectively, as well as the initial state parameters of the mixed fluid.

[0051] Step S2: Construct the multi-phase compressible fluid motion control equation for the mixed fluid in the computational domain. This multi-phase compressible fluid motion control equation characterizes the relationship between the state parameters of the mixed fluid before reaching the equilibrium state and the mass fraction Y k of any k-th phase fluid medium in the mixed fluid.

[0052] The state parameters of the mixed fluid before reaching the equilibrium state include the density ρ, velocity vector u, pressure p, and total energy per unit mass E0 of the mixed fluid before reaching the equilibrium state. The established governing equations of the multiphase compressible fluid motion include the mass conservation equation, momentum conservation equation, energy conservation equation, and phase change equation, which can be expressed as:

[0053]

[0054] Among them, the total energy per unit mass E0 = e + 0.5u×u, where e is the internal energy per unit mass of the mixed fluid, I is the identity matrix, is the gradient operator, and t represents time.

[0055] Step S3: Based on the NASG equation, establish the thermodynamic state equation of any k-th phase fluid medium.

[0056] The NASG equation comprehensively considers the thermal vibration effects of molecules and atoms in the fluid medium, as well as the repulsive and attractive effects between molecules. The state parameters of any k-th phase fluid medium include pressure p k , temperature T k , Gibbs free energy g k , and sound speed c k . Therefore, the thermodynamic state equation can be expressed as:

[0057]

[0058] It can be seen from the above equation (2) that the pressure p k (υ k , e k ) of the k-th phase fluid medium is related to the specific volume υ k and the internal energy per unit mass e k . The temperature T k (p k , υ k ) of the k-th phase fluid medium is related to the pressure p k and the specific volume υ k . The Gibbs free energy g k (p k , T k ) of the k-th phase fluid medium is related to the pressure p k and the temperature T k . The sound speed c k (p k , υ k ) of the k-th phase fluid medium is related to the pressure p k and the specific volume υ k . That is, the state parameters of the k-th phase fluid medium are coupled and related to the specific volume υ k and the internal energy per unit mass e k of the k-th phase fluid medium. The specific volume υk = 1 / ρ k , where ρ k is the density of the k-th phase fluid medium.

[0059] In the above formula (2), for the k-th phase fluid medium, γ k , q k , b k , q′ k , C v,k are all state equation parameters obtained by pre-fitting. These state equation parameters can be obtained by fitting data such as the saturation pressure related to temperature, the specific volume of each phase, and the enthalpy of each phase in the fluid thermodynamics experiment. The state equation parameters corresponding to various fluid media are shown in Table 1 below:

[0060] Table 1 State Equation Parameters of Each Phase Fluid Medium

[0061]

[0062] Phase transformation occurs between the liquid phase and the vapor phase. If there are both the liquid phase and the vapor phase in the mixed fluid after the phase transformation relaxation, then the saturation temperature T sat and the saturation pressure p sat also need to satisfy the vapor-liquid two-phase saturation curve relationship, and the vapor-liquid two-phase saturation curve relationship can also be constructed using the state equation parameters corresponding to various fluid media, expressed as:

[0063]

[0064] Among them, A, B, C, and D are calculated based on the state equation parameters of the liquid phase and the vapor phase, and specifically:

[0065]

[0066] Among them, C p,k represents the specific heat capacity at constant pressure of the k-th phase fluid medium. The corresponding values for the liquid phase and the vapor phase can refer to the state equation parameters shown in the above table.

[0067] Step S4: According to the characteristic that g1 = g2 when the mixed fluid reaches equilibrium at constant temperature and constant pressure, construct the cavitation phase transformation equation after the mixed fluid completes phase transformation relaxation and reaches equilibrium.

[0068] For the case where mass and energy conversion occur between the liquid phase and its corresponding vapor phase medium in the fluid, it is considered that the Gibbs free energy between the liquid phase and its vapor phase is equal after reaching equilibrium. At the same time, combined with the assumption of the multi-phase compressible fluid motion control equation, it can be considered that the mixed fluid satisfies the following after reaching equilibrium:

[0069]

[0070] In Equation (4), g1 = g2 means that the Gibbs free energies of the liquid phase and the vapor phase of the mixed fluid are equal after reaching the state, which means that the mixed fluid satisfies the isothermal and isobaric conditions after reaching the equilibrium state, is the mass fraction of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state, is the specific volume of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state, is the internal energy per unit mass of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state.

[0071] Since phase transitions occur between the liquid phase and the vapor phase, the mass fractions of the liquid phase and the vapor phase will change before and after reaching the equilibrium state. However, neither air nor explosive gas undergoes phase transitions, so their mass fractions remain unchanged, that is, for k ≥ 3, there is In addition, according to the thermodynamic state equation (2) of any k-th phase fluid medium, it can be known that the specific volume υ k of the k-th phase fluid medium is related to the pressure p k and the temperature T k of the k-th phase fluid medium. The internal energy per unit mass e k of the k-th phase fluid medium is also related to the pressure p k and the temperature T k of the k-th phase fluid medium. Therefore, after determining the relationship between υ k and p k and T k , and determining the relationship between e k and p k and T k , combined with the isothermal and isobaric assumptions, when using p * to represent the pressure of the mixed fluid after reaching the equilibrium state and using T * to represent the temperature of the mixed fluid after reaching the equilibrium state, the specific volume of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state can be expressed as υ k (T * , p * ), and the internal energy per unit mass of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state can be expressed as e k (T * , p * ). Thus, the cavitation phase transition equation after the mixed fluid completes the phase transition relaxation and reaches the equilibrium state can be obtained as:

[0072]

[0073] In Equation (5), W kis the molar mass of the k-th phase fluid medium, and the molar masses of the phase fluid media can be referred to as shown in Table 1 above. p partial represents the partial pressure of the vapor phase after the mixed fluid reaches the equilibrium state, and it can be known from the formula that it is proportional to the molar mass fraction of the vapor phase. p sat (T * ) represents T * corresponding saturation pressure. The saturation pressure p sat (T * ) is obtained by substituting T * as the saturation temperature T sat into the vapor-liquid two-phase saturation curve relationship (3) to determine the saturation pressure p sat .

[0074] Step S5, solve the system of equations composed of the multiphase compressible fluid motion control equation, the thermodynamic state equations of each phase fluid medium, and the cavitation phase change equation in each cycle step to obtain the numerical prediction result.

[0075] The multiphase compressible fluid motion control equation (1), the thermodynamic state equations of each phase fluid medium (2), and the cavitation phase change equation (5) together constitute the system of equations of the multiphase compressible fluid system considering the phase change effect, and can be solved distributively using the relaxation method, including:

[0076] First, combine the thermodynamic state equations of each phase fluid medium in the mixed fluid (2) to solve the multiphase compressible fluid motion control equation (1) to obtain the specific volume υ, internal energy per unit mass e, pressure p, temperature T, and the mass fraction Y of each phase fluid medium before the mixed fluid reaches the equilibrium state. k . Specifically, use the second-order MUSCL-Hancock method and the HLLC approximate Riemann solver to solve the multiphase compressible fluid motion control equation (1). The HLLC approximate Riemann solver needs to use the mixed sound speed to calculate the flux. Then, this step first uses the thermodynamic state equation of each phase fluid medium (2) to solve the sound speed c of the k-th phase fluid medium k (p k , υ k ), and then calculates the mixed sound speed using the Allaire calculation formula and applies it to the HLLC approximate Riemann solver. The specific solution method and principle can refer to the prior art and will not be elaborated here.

[0077] After the above solution, the state parameters of the mixed fluid before reaching the equilibrium state and the mass fraction Y of each phase fluid medium are obtained kAnd the state parameters of each phase of fluid medium before the mixed fluid reaches the equilibrium state. Further, the specific volume υ = 1 / ρ of the mixed fluid before reaching the equilibrium state can be calculated using the density ρ of the mixed fluid before reaching the equilibrium state. The specific internal energy e of the mixed fluid before reaching the equilibrium state can be calculated using the total energy per unit mass E0 and the velocity vector u of the mixed fluid before reaching the equilibrium state. Additionally, by combining the isothermal and isobaric assumptions of the mixed fluid, the pressure p and temperature T of the mixed fluid before reaching the equilibrium state can be obtained.

[0078] Although the above solution results satisfy the isothermal and isobaric conditions of each phase, the liquid phase and the vapor phase are in a sub-equilibrium state and mass transfer will occur. They will reach the equilibrium state only after the phase change relaxation of the liquid phase and the vapor phase is completed. And the conversion process is based on criterion (4), and actually numerical discrete solution is carried out by solving (5). Therefore, substituting the above solution results into the cavitation phase change equation (5) and using the iterative method to solve, the pressure p of the mixed fluid after reaching the equilibrium state can be obtained. * 、temperature T * and the mass fraction of each phase of fluid medium The iterative method used is, for example, the Newton-Raphson iterative method.

[0079] Furthermore, according to the volume fraction z of the k-th phase fluid medium in the mixed fluid can be calculated. k When the volume fraction z2 of the vapor phase in the mixed fluid reaches a predetermined threshold, the corresponding grid cell can be classified as a cavitation grid. Summarizing all the cavitation grids in the current cycle step can form the cavitation domain of the current cycle step, thereby visually showing the evolution process of the cavitation domain. Among them, ρ is the density of the mixed fluid before reaching the equilibrium state, is the density of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state and

[0080] In an example, take the computational domain as Figure 2 shown. The explosion bubble domain 230 is generated by the explosion of 50 g of equivalent TNT charge at a depth of 0.5 m underwater. The computational domain is discretely divided into 1000*1600 grid cells. The initial conditions set for the computational domain are:

[0081] In the mixed fluid in the air domain 210 above the free surface, the mass fraction Y2 of the vapor phase = 10 -9 、the volume fraction is z2 = 1.56×10 -9 , and in the mixed fluid in the air domain 210, except for the vapor phase, it is all air, and the mass fractions of the liquid phase and the explosive gas are both 0, that is, Y1 = Y4 = 0. The density ρ of the mixed fluid in the air domain 210 is initialized to be 1.18 kg / m 3 、the pressure p = 10 5Pa, temperature T = 295K, and the mixed fluid in the air domain 210 is in a superheated state.

[0082] In the mixed fluid in the water domain 220 below the free surface, the mass fraction Y2 of the vapor phase is 10 -8 , and the volume fraction is z2 = 1.39×10 -5 . In the mixed fluid in the water domain 220, the mass fraction Y3 of air is 4.73×10 -7 . In the mixed fluid in the water domain 220, except for the vapor phase and air, the rest is the liquid phase, and the mass fraction of the explosive gas is 0, that is, Y4 = 0. The density ρ of the mixed fluid in the water domain 220 is initialized to be 1054 kg / m 3 , pressure p = 10 5 Pa, temperature T = 295K, and the mixed fluid in the water domain 220 is in a saturated state.

[0083] In the mixed fluid in the explosion bubble domain 230, the volume fraction of the vapor phase is z2 = 8.07×10 -7 . In the mixed fluid in the explosion bubble domain 230, except for the vapor phase, the rest is the explosive gas, and the mass fractions of the liquid phase and air are both 0. The density ρ of the mixed fluid in the explosion bubble domain 230 is initialized to be 1606 kg / m 3 , pressure p = 10 9 Pa, temperature T = 1120K, and the mixed fluid in the explosion bubble domain 230 is in a superheated state.

[0084] The boundary conditions of the computational domain are set as follows: The upper, lower, and right boundaries of the air domain 210 are set to non-reflective conditions, and the CFL is set to 0.5.

[0085] Based on the above initial conditions and boundary conditions, the evolution process of the density ρ of the mixed fluid at multiple typical moments obtained by numerical prediction using the two-dimensional near-free-surface underwater explosion cavitation prediction method of this application is as Figure 3 shown, the evolution process of the pressure p of the mixed fluid is as Figure 4 shown, the evolution process of the volume fraction z2 of the vapor phase of the mixed fluid is as Figure 5 shown, and the evolution process of the cavitation domain (the red marked area in the figure) is as Figure 6 shown. From Figures 3 to 6It can be seen that the initial shock wave caused by the explosion caused cavitation underwater after being reflected by the free surface. From the relatively high-pressure flow field near the explosion bubble and the wave propagation speed, the rarefaction wave reflected from the free surface has propagated to the upper wall of the explosion bubble and reflected the compression wave. At the same time, from the evolution cloud diagram of the volume fraction z2 of the vapor phase, local cavitation also occurred near the free surface above the water surface. This is because the initial shock wave transmits the compression wave into the air after encountering the free surface, causing the liquid water and steam near the water surface to convect to the air directly above the water surface under high pressure and high speed.

[0086] The above is only a preferred embodiment of the present application, and the present application is not limited to the above embodiments. It is understood that other improvements and changes directly derived or associated by those skilled in the art without departing from the spirit and concept of the present application should be considered to be included in the protection scope of the present application.

Claims

1. A two-dimensional near-free surface underwater explosion cavitation prediction method, characterized in that, The two-dimensional near-free surface underwater explosion cavitation prediction method includes: Constructing a computational domain for the two-dimensional near-free surface underwater explosion scenario and performing discrete grid division on the computational domain. The computational domain includes an air domain above the free surface, a water domain below the free surface, and an explosion bubble domain at the near-free surface position in the water domain; Construct the state parameters of the mixed fluid in the computational domain before reaching the equilibrium state and the mass fraction Y of any k-th phase fluid medium in the mixed fluid k of the multiphase compressible fluid motion control equations. The multiphase compressible fluid motion control equations include the mass conservation equation, the momentum conservation equation, the energy conservation equation, and the phase change equation. k = 1 represents the liquid phase, and k = 2 represents the vapor phase; Based on the NASG equation, a thermodynamic state equation for any k-th phase fluid medium is constructed. The state parameters of any k-th phase fluid medium include pressure p k , temperature T k , Gibbs free energy g k and sound speed c k and are coupled with specific volume υ k and internal energy e per unit mass k ; According to the characteristic that the mixed fluid has equal temperature and pressure and reaches an equilibrium state with g1 = g2, constructing a cavitation phase change equation after the mixed fluid completes phase change relaxation and reaches an equilibrium state; Solving a system of equations composed of the multi-phase compressible fluid motion control equation, the thermodynamic state equation of each phase fluid medium, and the cavitation phase change equation in each cycle step to obtain a numerical prediction result.

2. The two-dimensional near-free surface underwater explosion cavitation prediction method according to claim 1, characterized in that, The constructed cavitation phase change equation after the mixed fluid completes phase change relaxation and reaches an equilibrium state is: Where N is the number of components of the fluid medium contained in the mixed fluid, any fluid medium represented by k≥3 is a gas-phase medium and does not participate in the phase change, e is the internal energy per unit mass of the mixed fluid before reaching the equilibrium state, υ is the specific volume of the mixed fluid before reaching the equilibrium state, p * is the pressure of the mixed fluid after reaching the equilibrium state, T * is the temperature of the mixed fluid after reaching the equilibrium state, W k is the molar mass of the k-th phase fluid medium, Y k is the mass fraction of the k-th phase fluid medium in the mixed fluid before reaching the equilibrium state, is the mass fraction of the k-th phase fluid medium in the mixed fluid after reaching the equilibrium state, υ k (T * , p * ) is the specific volume of the k-th phase fluid medium in the mixed fluid after reaching the equilibrium state and is related to p * and T * , e k (T * , p * ) is the internal energy per unit mass of the k-th phase fluid medium in the mixed fluid after reaching the equilibrium state and is related to p * and T * , υ k (T * , p * ) and e k (T * , p * )'s relationship with p * and T * is determined according to the thermodynamic state equation of the k-th phase fluid medium, e k (T * , p * ) is the internal energy per unit mass of the k-th phase fluid medium in the mixed fluid after reaching the equilibrium state and is related to p * and T * , p partial represents the partial pressure of the vapor phase of the mixed fluid after reaching the equilibrium state, p sat (T * ) represents the saturation state pressure corresponding to T * .

3. The two-dimensional near-free surface underwater explosion cavitation prediction method according to claim 2, wherein The thermodynamic state equation of any k-th phase fluid medium includes: where the pressure p k (υ k ,e k ) is related to υ k and e k , the temperature T k (p k ,υ k ) is related to p k and υ k , the Gibbs free energy g k (p k ,T k ) is related to p k and T k , the speed of sound c k (p k ,υ k ) is related to p k and υ k ; γ k 、q k 、b k 、 ′ q k , C v,k are all state equation parameters obtained by pre-fitting.

4. The two-dimensional near-free surface underwater explosion cavitation prediction method according to claim 3, wherein For the liquid phase with k = 1, the state equation parameters are γ1 = 1.19, q1 = -1177788 J / kg, b1 = 6.61×10 -4 m 3 / kg, q1′ = 0 J / kg / K, C v,1 = 3610 J / kg / K, W1 = 18 g / mol; The state equation parameters of the vapor phase represented by k = 2 are γ2 = 1.47, q2 = 2077616 J / kg, b2 = 0 m 3 / kg, q2′ = 2077616 J / kg / K, C v,2 = 955 J / kg / K, W2 = 18 g / mol; The state equation parameters of air represented by k = 3 are γ3 = 1.4, q3 = 0 J / kg, b3 = 0 m 3 / kg, q3′ = 0 J / kg / K, C v,3 = 719 J / kg / K, W3 = 29 g / mol; The state equation parameters of the explosive gas represented by k = 4 are γ4 = 1.8, q4 = 0 J / kg, and b4 = 0 m 3 / kg, q4′ = 0 J / kg / K, C v,4 = 695 J / kg / K, W4 = 290 g / mol.

5. The two-dimensional near-free-surface underwater explosion cavitation prediction method according to claim 3, wherein, Calculate T * The corresponding saturation pressure p sat (T * ) The method includes: Take T * as the saturation temperature T sat and substitute it into the following relational expression to calculate the corresponding saturation pressure p sat Take p * corresponding to T sat (T * ): where A, B, C, and D are calculated based on the state equation parameters of the liquid phase and the vapor phase.

6. The two-dimensional near-free surface underwater explosion cavitation prediction method according to claim 5, wherein C p,1 = 4285 J / kg / K, C p,2 = 1401 J / kg / K.

7. The two-dimensional near-free surface underwater explosion cavitation prediction method according to claim 2, wherein The numerical prediction results obtained by solving in each cycle step include: Solving the motion control equations of multiphase compressible fluids by combining with the thermodynamic state equations of each phase fluid medium in the mixed fluid to obtain the specific volume υ, internal energy per unit mass e, pressure p, temperature T, and mass fraction Y of each phase fluid medium before the mixed fluid reaches the equilibrium state k And the state parameters of each phase fluid medium before the mixed fluid reaches the equilibrium state are substituted into the cavitation phase change equation and solved by an iterative method to obtain the pressure p after the mixed fluid reaches the equilibrium state * , temperature T * And the mass fraction of each phase fluid medium 8. The two-dimensional near-free surface underwater explosion cavitation prediction method according to claim 7, characterized in that The numerical prediction results obtained by solving in each cycle step further include: According to calculate the volume fraction z of the k-th phase fluid medium in the mixed fluid k , when the volume fraction z2 of the vapor phase in the mixed fluid reaches a predetermined threshold, assign the corresponding grid cell to the cavitation grid, and obtain the cavitation domain composed of all cavitation grids in each cycle step; where ρ is the density of the mixed fluid before reaching the equilibrium state is the density of the k-th phase fluid medium after the mixed fluid reaches the equilibrium state and 9. The two-dimensional near-free surface underwater explosion cavitation prediction method according to claim 7, wherein Solving the multi-phase compressible fluid motion control equation in combination with the thermodynamic state equation of each phase fluid medium in the mixed fluid includes: Solving for the speed of sound \(c\) of the \(k\)-th phase fluid medium using the thermodynamic state equation of the \(k\)-th phase fluid medium k (p k , υ k ), and solving the governing equations of multiphase compressible fluid motion using the second-order MUSCL-Hancock method and the HLLC approximate Riemann solver in combination with the obtained speed of sound.

10. The two-dimensional near-free-surface underwater explosion cavitation prediction method according to claim 7, wherein, The multi-phase compressible fluid motion control equation is: where ρ is the density of the mixed fluid before the equilibrium state and υ = 1 / ρ, u is the velocity vector of the mixed fluid, E0 is the total energy per unit mass of the mixed fluid and E0 = e + 0.5u×u, I is the identity matrix, is the gradient operator, and t represents time.

Images