CN111400968A - Method for constructing compressible two-phase flow interface condition under curve coordinate system - Google Patents

Abstract

The invention relates to a construction method of a compressible two-phase flow interface condition under a curve coordinate system, which is used for popularizing an interface flow numerical simulation technology under a traditional Cartesian coordinate system and is suitable for flow field simulation of different fluid interface evolutions under a general curve coordinate system. The method comprises the following steps: converting a level set equation and a reinitialization equation representing interface position information into a general curve coordinate system for solving through coordinate transformation; the virtual fluid method reflecting the two-phase flow effect is popularized to a general curve coordinate system for interface condition definition. The advantages are that: the improved level set method under the curve coordinate system is suitable for representing the dynamic position information of the interface under any structural grid; the improved virtual fluid method under the curve coordinate system can be applied to conventional CFD software, so that the CFD software has the capacity of simulating two-phase flow; the interface condition expression form designed by adopting a coordinate transformation method is suitable for solving the problem of the compressible two-phase flow interface based on the curve grid under the complex shape.

Description

A method for constructing interface conditions for compressible two-phase flow in curvilinear coordinate system

technical field

The invention belongs to the technical field of computational fluid mechanics, and in particular relates to a method for constructing interface conditions of a compressible two-phase flow in a curvilinear coordinate system.

Background technique

The dynamic characteristics of compressible two-phase flow interfacial motion are of great significance to the study of high-speed flow fields containing droplets, bubbles, and elastic-plastic walls. The role of high-confidence numerical simulation is increasing day by day, and it can completely replace some expensive, dangerous, and even difficult to implement experiments, greatly reducing the development cost and shortening the development cycle, and is playing an increasingly important role.

Since the physical parameters of different fluids on both sides of the interface are very different, if the numerical format for simulating single-medium fluids is directly applied to the simulation of multi-medium fluids, numerical instability may occur. Especially when highly nonlinear waves such as shock waves interact with the interface, non-physical oscillations are likely to occur near the interface, which even makes the calculation difficult. The virtual fluid method treats the interface as a special kind of boundary. Based on the solution of the two-phase flow Riemann problem, the boundary conditions are constructed in the virtual flow field near the interface, which can avoid numerical instability near the interface and is easy to generalize to multi-dimensional.

However, the above techniques are mainly applied to Cartesian grids and cannot be used in general curvilinear coordinate systems at present.

SUMMARY OF THE INVENTION

The present invention aims to propose a method for constructing compressible two-phase flow interface conditions in a curvilinear coordinate system. Based on coordinate transformation, the level set method and virtual fluid method in the traditional Cartesian coordinate system are extended to make them suitable for general curves. The flow field simulation of the evolution of different fluid interfaces in the coordinate system has important application value.

The technical solution of the present invention is:

A method for constructing compressible two-phase flow interface conditions in a curvilinear coordinate system, comprising the following steps:

1) According to the initial value of the signed distance function φ of each physical grid node, the virtual fluid area outside each fluid interface is respectively defined; the media properties of two adjacent fluids are different and the two adjacent fluids are incompatible with each other; the medium The different properties are as follows: at least one of the velocity, pressure or density of two adjacent fluids is different; the fluids are all compressible fluids; the virtual fluid region contains 3-5 physical mesh nodes of another fluid near the interface ; The physical grid is obtained by using grid generation software to mesh the computational flow field area;

2) utilize the level set equation under the general curvilinear coordinate system, update the initial value of the described symbolic distance function φ of step 1), obtain the updated symbolic distance function φ;

3) According to the updated signed distance function φ described in step 2), solve the reinitialization equation of the signed distance function under the general curvilinear coordinate system, until the stable solution of the reinitialized equation is obtained, and the stable solution is used as the revised signed distance function φ. ;

4) According to the modified signed distance function φ described in step 3), update the fluid interface position information;

5) According to the virtual fluid area outside the fluid interface in step 1), for each virtual fluid node P in the virtual fluid area, find the real fluid node P' of the real fluid area on the other side of the interface, so that the interface of the virtual fluid node P is The angle θ formed by the unit normal vector n and the interface unit normal vector n' of the real fluid node P' is the smallest; the interface normal direction of each node is determined according to the initial value of the signed distance function φ in step 1); wherein, θ=arccos(n·n′);

6) Projecting the velocity field of the virtual fluid node P and its corresponding real fluid node P' in step 5) to the interface normal, respectively, to obtain the interface normal velocity u _n of the virtual fluid node P and the real fluid node P';

7) According to the interface normal velocity u _n of the virtual fluid node P and the real fluid node P' described in step 6) and the pressure and density corresponding to the virtual fluid node P and the real fluid node P', construct and solve the two equations in the interface normal direction. Phase flow Riemann problem, obtain the interface pressure p _I , the interface normal velocity u _n,I and the density on both sides of the interface ρ _IL ,ρ _IR ;

8) According to the interface pressure p _I , the interface normal velocity u _n,I and the density ρ _IL , ρ _IR on both sides of the interface described in step 7), define the interface condition of the virtual fluid node P, that is, the virtual fluid state;

9) Repeat steps 5) to 8) until the interface conditions of all virtual fluid nodes P of each fluid virtual fluid region are obtained.

10) According to step 9), the interface condition is the interface condition at the current calculation time, and steps 1) to 9) are repeated for the construction of the interface condition at the next calculation time.

Step 2) the level set equation under the general curvilinear coordinate system, specifically:

φ _t +Uφ _ξ +Vφ _η +Wφ _ζ =0;

U=ξ _x u+ξ _y v+ξ _z w;

V=η _x u+η _y v+η _z w;

W= _ζxu + _ζyv + _ζzw ;

Among them, φ _t represents the partial derivative of the signed distance function φ with respect to time t, t represents time, U, V, W are the fluid reversal velocities of each node; u, v, w are the fluid velocities of each node; (x , y, z) represents the physical grid coordinates of the flow field, and the calculation flow field area can be meshed by grid generation software to obtain the x, y, z coordinate values of each node; (ξ, η, ζ) Represents the computational grid coordinates obtained _by transforming the physical grid _coordinates from the _{general curvilinear coordinate system} to the _{Cartesian coordinate} system; _z is the metric coefficient of the calculation grid, which represents the partial derivatives of the calculation grid coordinates ξ, η, ζ to the physical grid coordinates x, y, z, usually calculated by the central difference; φ _ξ , φ _η and φ _ζ represent the signed distance The partial derivatives of the function φ with respect to the calculated grid coordinates ξ, η, ζ can be differentially discretized according to the selected format accuracy.

Step 3) The reinitialization equation of the signed distance function under the general curvilinear coordinate system, specifically:

a=ξ _x η _x +ξ _y η _y +ξ _z η _z , b=ξ _x ζ _x +ξ _y ζ _y +ξ _z ζ _z , c=η _x ζ _x +η _y ζ _y +η _z ζ _z ;

Among them, τ represents the pseudo-time, ε is the side length of the physical grid where each node is located; φ ₀ represents the value of the signed distance function before reinitialization.

Step 5) The unit normal vector n of the nodes near the interface under the general curvilinear coordinate system, specifically:

Step 6) The interface normal velocity u _n of the nodes near the interface in the general curvilinear coordinate system, specifically:

Step 8) The method for defining the interface conditions of the virtual fluid node P based on the obtained interface pressure p _I , the interface normal velocity u _n,I and the density on both sides of the interface ρ _IL , ρ _IR , specifically:

If the signed distance function φ corresponding to the virtual fluid node P is less than 0, then

If the signed distance function φ corresponding to the virtual fluid node P is greater than 0, then

where the superscript * denotes virtual fluid. ρ ^* represents the density of the virtual fluid, p ^* represents the pressure of the virtual fluid; u ^* , v ^* and w ^* represent the virtual fluid velocity; L represents the node where the fluid signed distance function φ is less than 0; R represents the fluid signed distance function φ greater than 0 node; let H=L or R, u _H , v _H , w _H and U _H , V _H , W _H be the fluid velocity and fluid reversal velocity of the real fluid node P', respectively.

Compared with the prior art of the present invention, its advantages and innovations are mainly reflected in the following aspects:

1) The present invention extends the level set equation and its re-initialization equation under the traditional Cartesian coordinate system to the general curvilinear coordinate system, and is suitable for characterizing the interface dynamic position information under any structural grid in principle;

2) The present invention extends the process of constructing interface conditions by the virtual fluid method under the Cartesian coordinate system to the general curvilinear coordinate system, which can be applied to conventional CFD software in principle, so that it has the ability to simulate two-phase flow;

3) The present invention adopts the coordinate transformation method to design the expression form of the interface condition, which is suitable for solving the interface problem of compressible two-phase flow based on the curve grid under the complex shape.

Description of drawings

Figure 1 is a schematic diagram of the coordinate transformation from physical space to computational space for the interface problem of two-phase flow field;

Figure 2 shows the construction process of the interface conditions of the two-phase flow field in the general curvilinear coordinate system;

Fig. 3 is the distribution of the virtual fluid area and the corresponding relationship of the nodes under the general curvilinear coordinate system;

Figure 4 shows the mesh division of the simulated underwater explosion problem in the general curvilinear coordinate system, and the comparison of the numerical schlieren at 1.0ms time with the results in the Cartesian coordinate system;

Figure 5 shows the comparison of the density distribution on the x=0 axis at the time of 0.2ms and 1.2ms to simulate the underwater explosion problem in the general curvilinear coordinate system and the results of the self-developed program and the literature in the Cartesian coordinate system.

Detailed ways

Research problems such as gas-liquid mixing in engine combustion chambers and underwater explosions all involve interface evolution processes such as the deformation and collapse of bubbles/droplets, and the interaction between shock waves and interfaces. Under normal circumstances, the calculation area or shape is complex, and the mesh is usually divided according to the body. The mesh is generally irregular and not strictly orthogonal. It is necessary to establish the interface condition model in the general curvilinear coordinate system, and use the interface condition model to calculate the two-phase flow field, so as to obtain the simulation of the interface evolution and instability process.

The present invention uses the coordinate transformation method to transform the state of matter and related equations near the interface of the two-phase flow field from the physical space (x, y, z) to the computational space (ξ, η, ζ). The two-dimensional schematic diagram of the coordinate transformation is shown in Figure 1 shown. In the general curvilinear coordinate system, the construction process of the interface conditions of the two-phase flow field is shown in Figure 2, which mainly includes the acquisition and update of the interface position information, and the definition of the interface conditions. Include the following steps:

1) According to the initial value of the signed distance function φ of each physical grid node, the virtual fluid area outside each fluid interface is respectively defined; the media properties of two adjacent fluids are different and the two adjacent fluids are incompatible with each other; the medium The different properties are as follows: at least one of the velocity, pressure or density of two adjacent fluids is different; the fluids are all compressible fluids; the virtual fluid region contains 3-5 physical mesh nodes of another fluid near the interface ; The physical grid is obtained by using grid generation software to mesh the computational flow field area;

2) utilize the level set equation under the general curvilinear coordinate system, update the initial value of the described symbolic distance function φ of step 1), obtain the updated symbolic distance function φ;

3) According to the updated signed distance function φ described in step 2), solve the reinitialization equation of the signed distance function under the general curvilinear coordinate system, until the stable solution of the reinitialized equation is obtained, and the stable solution is used as the revised signed distance function φ. ;

4) According to the modified signed distance function φ described in step 3), update the fluid interface position information;

5) According to the virtual fluid area outside the fluid interface in step 1), for each virtual fluid node P in the virtual fluid area, find the real fluid node P' of the real fluid area on the other side of the interface, so that the interface of the virtual fluid node P is The angle θ formed by the unit normal vector n and the interface unit normal vector n' of the real fluid node P' is the smallest; the interface normal direction of each node is determined according to the initial value of the signed distance function φ in step 1); wherein, θ=arccos(n·n′);

6) Projecting the velocity field of the virtual fluid node P and its corresponding real fluid node P' described in step 5) to the interface normal direction, respectively, to obtain the interface normal velocity u _n of the virtual fluid node P and the real fluid node P';

7) According to the interface normal velocity u _n of the virtual fluid node P and the real fluid node P' described in step 6) and the pressure and density corresponding to the virtual fluid node P and the real fluid node P', construct and solve the two equations in the interface normal direction. Phase flow Riemann problem, obtain the interface pressure p _I , the interface normal velocity u _n,I and the density on both sides of the interface ρ _IL ,ρ _IR ;

8) According to the interface pressure p _I , the interface normal velocity u _n,I and the densities ρ _IL , ρ _IR on both sides of the interface described in step 7), define the interface condition of the virtual fluid node P, that is, the virtual fluid state;

9) Repeat steps 5) to 8) until the interface conditions of all virtual fluid nodes P in each fluid virtual fluid region are obtained;

10) According to step 9), the interface condition is the interface condition at the current calculation time, and steps 1) to 9) are repeated for the construction of the interface condition at the next calculation time.

The specific implementation process is as follows:

(I) Acquisition and update of interface location information

(I-1) According to the value of the signed distance function φ, define the virtual fluid area outside the interface for each fluid. A node with a negative sign (φ<0) represents a fluid on the left side of the interface, a node with a positive sign (φ>0) represents another fluid on the right side of the interface, and a node with a zero sign (φ=0) The location represents the interface. The virtual fluid region outside the interface generally spans at least two grids of the interface, depending on the numerical format used.

(I-2) Solve the level set equation in the general curvilinear coordinate system:

φ _t +Uφ _ξ +Vφ _η +Wφ _ζ =0

U=ξ _x u+ξ _y v+ξ _z w;

V=η _x u+η _y v+η _z w;

W= _ζxu + _ζyv + _ζzw ;

和

然后取：Among them, φ _t represents the partial derivative of the signed distance function φ with respect to time t, t represents time, U, V, W are the fluid reversal velocities of each node; u, v, w are the fluid velocities of each node; (x , y, z) represents the physical grid coordinates of the flow field, and the calculation flow field area can be meshed with grid generation software to obtain the x, y, z coordinate values of each node; (ξ, η, ζ) Represents the computational grid coordinates obtained _by converting the physical grid _coordinates from the _{general curvilinear coordinate system} to the _{Cartesian coordinate} system; _z is the metric coefficient of the calculation grid, which represents the partial derivative of the calculation grid coordinates ξ, η, ζ to the physical grid coordinates x, y, z, usually calculated by the central difference; φ _ξ , φ _η and φ _ζ represent the signed distance The partial derivatives of the function φ with respect to the calculated grid coordinates ξ, η, ζ can be differentially discretized according to the selected format accuracy. For example, using the WENO format to construct the discrete derivative values on the left and right sides

and

Then take:

Do similar processing for φ _η and φ _ζ .

(I-3) Reinitialize the signed distance function in the general curvilinear coordinate system, that is, solve:

a=ξ _x η _x +ξ _y η _y +ξ _z η _z , b=ξ _x ζ _x +ξ _y ζ _y +ξ _z ζ _z , c=η _x ζ _x +η _y ζ _y +η _z ζ _z ;

和

然后取Among them, τ represents pseudo-time, ε is a small quantity, which can correspond to the side length of the physical grid where each node is located; φ ₀ represents the value of the signed distance function before re-initialization; φ _ξ , φ _η and φ _ζ represent the signed distance The partial derivatives of the function φ with respect to the calculated grid coordinates ξ, η, ζ can be differentially discretized according to the selected format accuracy. For example, using the WENO format to construct the discrete derivative values on the left and right sides

and

then take

对φ_η和φ_ζ做类似的处理。in,

Do similar processing for φ _η and φ _ζ .

(I-4) Obtain the interface position information at the new time, and the position of φ=0 represents the interface position at the new time.

(II) Definition of Interface Conditions

For each virtual fluid node P in the virtual fluid region near the interface:

(II-1) Find the real fluid node P' in the real fluid region on the other side of the interface, so that the angle θ formed by the interface unit normal vector n of the virtual fluid node P and the interface unit normal vector n' of the real fluid node P' is the smallest , see Figure 3. Among them, θ=arccos(n·n′), the interface unit normal vector is expressed as

(II-2) Project the velocity fields of the virtual fluid node P and the real fluid node P' to the interface normal respectively, and establish the interface normal velocities u _n of the virtual fluid node P and the real fluid node P',

Among them, the subscript n represents the interface normal.

(II-3) Based on the obtained interface normal velocity u _n of the virtual fluid node P and the real fluid node P' and other fluid states, construct and solve the two-phase flow Riemann problem in the interface normal direction

分别表示界面法向上的守恒变量和通量。ρ是流体密度；p是流体压力；u_n是流体法向速度；

是流体法向上的总能量，其中e为流体单位质量的内能。U_n,L和U_n,R是由位于界面位置n_I处的物质界面分开的两种流体常值状态，下标L和R表示界面两侧位于虚拟流体节点P和真实流体节点P’的流体。L表示流场界面左侧(φ<0的节点)，R表示流场界面右侧(φ>0的节点)。Obtain the interface pressure p _I , the interface normal velocity u _n,I and the density ρ _IL ,ρ _IR on both sides of the interface. Here, U _n =[ρ,ρu _n ,E _n ] ^T and

are the conserved variables and fluxes in the interface normal, respectively. ρ is the fluid density; p is the fluid pressure; u _n is the fluid normal velocity;

is the total energy normal to the fluid, where e is the internal energy per unit mass of the fluid. U _n,L and U _n,R are two fluid constant states separated by the material interface located at the interface position n _I , and the subscripts L and R indicate that the two sides of the interface are located at the virtual fluid node P and the real fluid node P' fluid. L represents the left side of the flow field interface (nodes with φ<0), and R represents the right side of the flow field interface (nodes with φ>0).

(II-4) Based on the obtained interface pressure p _I , the interface normal velocity u _n,I and the densities ρ _IL ,ρ _IR on both sides of the interface, define the interface conditions of the virtual fluid node P (that is, the virtual fluid state):

If the signed distance function φ corresponding to the virtual fluid node P is less than 0, that is, the virtual fluid node P is located on the left side of the interface, then:

If the signed distance function φ corresponding to the virtual fluid node P is greater than 0, that is, the virtual fluid node P is located on the right side of the interface, then:

Here, the superscript ^* denotes virtual fluid. ρ ^* represents the density of the virtual fluid, p ^* represents the pressure of the virtual fluid; u ^* , v ^* and w ^* represent the virtual fluid velocity; L represents the node where the fluid sign distance function φ is less than 0, that is, the node is located on the left side of the interface, R represents The node with the fluid sign distance function φ greater than 0, that is, the node is located on the right side of the interface; let H=L or R, u _H , v _H , w _H and U _H , V _H , _WH are the fluid of the real fluid node P' respectively Velocity and Fluid Inversion Velocity.

Figures 4 and 5 show the numerical results of solving the underwater explosion problem in the curvilinear coordinate system by applying the technology of the present invention. In the initial conditions, the explosion gas is p=1GPa, ρ=1250kg/m ³ , the water is p=100kPa, ρ=1000kg/m ³ , and the air is p=100kPa, ρ=1kg/m ³ , all in a static state. The state equation adopts the rigid gas state equation p=(N-1)ρe-NB, for gas N=1.4, B=0, for water N=4.4, B=600MPa. Figure 4 compares the results of grid division, numerical schlieren at 1.0ms and the results of the self-developed program in the Cartesian coordinate system. Literature (Hu et al, J. Comput. Phys. 228(17)(2009) 6572-6589) compares the results in the Cartesian coordinate system. It can be seen that although the grid division and topology in the curvilinear coordinate system are completely different from the Cartesian coordinate system, the calculation results are in good agreement with the calculation results in the Cartesian coordinate system in the self-developed program and the literature.

Aiming at the problem that the prior art cannot be used in the general curvilinear coordinate system, the present invention aims to apply the virtual fluid method to the general curvilinear coordinate system, and construct a method for defining the interface conditions of the compressible two-phase flow in the general curvilinear coordinate system through coordinate transformation, Realize the numerical simulation of compressible two-phase flow under any grid structure. The method can solve practical problems with complex shapes and has important application value.

The above is only the best specific embodiment of the present invention, but the protection scope of the present invention is not limited to this. Substitutions should be covered within the protection scope of the present invention.

The content not described in detail in the present invention is known to those skilled in the art.

Claims (6)

1. A method for constructing a compressible two-phase flow interface condition under a curve coordinate system is characterized by comprising the following steps:

1) respectively defining a virtual fluid area outside each fluid interface according to the initial value of the symbol distance function phi of each physical grid node; the medium properties of two adjacent fluids are different, and the two adjacent fluids are not mutually fused; the media properties are different and specifically: at least one of the velocity, pressure or density of two adjacent fluids is different; the fluids are all compressible fluids; the virtual fluid zone contains 3-5 physical grid nodes of another fluid near the interface; the physical grid is obtained by utilizing grid generation software to perform grid division on a calculation flow field area;

2) updating the initial value of the symbol distance function phi in the step 1) by using a level set equation under a general curve coordinate system to obtain an updated symbol distance function phi;

3) solving a reinitialization equation of the symbol distance function under the general curve coordinate system according to the updated symbol distance function phi in the step 2) until a stable solution of the reinitialization equation is obtained, and taking the stable solution as the corrected symbol distance function phi;

4) updating the position information of the fluid interface according to the corrected symbol distance function phi in the step 3);

5) according to the virtual fluid area outside the fluid interface in the step 1), for each virtual fluid node P in the virtual fluid area, searching a real fluid node P ' of a real fluid area on the other side of the interface, so that an included angle theta formed by an interface unit normal vector n of the virtual fluid node P and an interface unit normal vector n ' of the real fluid node P ' is minimum; the interface normal direction of each node is determined according to the initial value of the symbol distance function phi in the step 1); wherein θ is arccos (n · n');

6) respectively projecting the velocity fields of the virtual fluid node P and the real fluid node P' in the step 5) to the interface normal direction to obtain the interface normal velocity u of the virtual fluid node P and the real fluid node P_n；

7) The interface normal velocities u of the virtual fluid node P and the real fluid node P' according to step 6)_nConstructing and solving a two-phase flow Riemann problem in the interface normal direction according to the pressure and density corresponding to the virtual fluid node P and the real fluid node P', and obtaining the interface pressure P_IInterfacial normal velocity u_n,IAnd density ρ at both sides of the interface_IL,ρ_IR；

8) The interface pressure p according to step 7)_IInterfacial normal velocity u_n,IAnd density ρ at both sides of the interface_IL,ρ_IRDefining the interface condition of the virtual fluid node P;

9) repeating the steps 5) to 8) until the interface conditions of all virtual fluid nodes P of each fluid virtual fluid area are obtained;

10) and according to the step 9), the interface condition is the interface condition of the current computing time, and the steps 1) to 9) are repeated for the construction of the interface condition of the next computing time.

2. The method for constructing the condition of the compressible two-phase flow interface under the curved coordinate system according to claim 1, wherein the level set equation under the general curved coordinate system in the step 2) is specifically as follows:

φ_t+Uφ_ξ+Vφ_η+Wφ_ζ＝0；

U＝ξ_xu+ξ_yv+ξ_zw；

V＝η_xu+η_yv+η_zw；

W＝ζ_xu+ζ_yv+ζ_zw；

wherein phi is_tRepresenting the partial derivative of a symbolic distance function phi to time t, wherein t represents time, and U, V and W are the fluid inversion speed of each node; u, v, w are per nodeFluid velocity (x, y, z) physical grid coordinates of the flow field (ξ, ζ) computational grid coordinates obtained by converting the physical grid coordinates from a generally curvilinear coordinate system to a cartesian coordinate system (ξ)_x,ξ_y,ξ_z,η_x,η_y,η_z,ζ_x,ζ_y,ζ_zIs a metric coefficient of the computational grid representing the partial derivative of the computational grid coordinates ξ, ζ, to the physical grid coordinates x, y, z, #_ξ、φ_ηAnd phi_ζRepresenting the partial derivative of the signed distance function phi with respect to the computed grid coordinates ξ, zeta.

3. The method for constructing the condition of the compressible two-phase flow interface under the curved coordinate system according to claim 2, wherein the step 3) is a reinitialization equation of a symbolic distance function under the general curved coordinate system, and specifically comprises the following steps:

a＝ξ_xη_x+ξ_yη_y+ξ_zη_z，b＝ξ_xζ_x+ξ_yζ_y+ξ_zζ_z，c＝η_xζ_x+η_yζ_y+η_zζ_z；

wherein, tau represents pseudo time and is the side length of the physical grid where each node is located; phi is a₀Representing the value of the symbol distance function before re-initialization.

4. The method for constructing the condition of the compressible two-phase flow interface under the curved coordinate system according to claim 3, wherein the unit normal vector n of the node near the lower interface of the general curved coordinate system in the step 5) is specifically as follows:

5. the method for constructing the condition of the compressible two-phase flow interface under the curved coordinate system according to claim 3, wherein the normal velocity u of the interface of the node near the interface under the general curved coordinate system in the step 6)_nThe method specifically comprises the following steps:

6. the method for constructing the compressible two-phase flow interface condition under the curved coordinate system according to any one of claims 3 to 5, wherein the step 8) is based on the obtained interface pressure p_IInterfacial normal velocity u_n,IAnd density ρ at both sides of the interface_IL,ρ_IRThe method for defining the interface condition of the virtual fluid node P specifically includes:

if the symbol distance function phi corresponding to the virtual fluid node P is less than 0, then:

if the symbol distance function phi corresponding to the virtual fluid node P is greater than 0, then:

where ρ is^*Representing the density, p, of the virtual fluid^*Representing the pressure of the virtual fluid; u. of^*、v^*And w^*Representing virtual fluid velocity, L representing nodes where the fluid-symbol distance function phi is less than 0, R representing nodes where the fluid-symbol distance function phi is greater than 0, H L or R, u_H、v_H、w_HAnd U_H、V_H、W_HThe fluid velocity and the fluid inversion velocity of the true fluid node P' are respectively.

Images