A curvature-corrected volume function method for simulating propellant injection atomization

Abstract

The application provides a curvature correction volume function method for simulating propellant injection atomization, which comprises the following steps: solving an interface transport equation by using a VOF method, updating a VOF function field; calculating a predicted velocity by using an explicit time advancing method; discretizing a surface tension term based on a virtual fluid method; solving a discrete pressure Poisson equation; correcting a velocity field by using a pressure gradient; and judging whether to end the calculation. The curvature correction volume function method for simulating propellant injection atomization can obtain a HF function increment to correct the VOF function, greatly improve the accuracy and smoothness of local curvature calculation, and thus realize high stability numerical simulation of the propellant injection atomization process.

Description

Technical Field

[0001] This invention relates to the field of computational fluid dynamics, and more specifically, to a method for simulating the curvature-corrected volume function of propellant injection atomization. Background Technology

[0002] The atomization of propellant injection in rocket engines is closely related to combustion stability and propulsion efficiency. Research on propellant injection atomization is time-consuming and costly. With advancements in computer hardware and numerical simulation methods, using gas-liquid two-phase flow numerical models for refined simulation of the injection atomization process has become an effective means of studying related issues. One of the main challenges in refined simulation of gas-liquid two-phase flow is the accurate simulation of surface tension. A core issue in surface tension models is the estimation of interface curvature, and the model's core performance often depends on the accuracy of its local curvature calculations.

[0003] Existing curvature estimation methods mainly include those based on Level Set (LS) functions, Smooth Volume of Fluid (VOF) functions, and Height Function (HF) methods. Among these, the LS method needs to be used in conjunction with the LS interface transport method, while the Smooth VOF and HF methods are generally used in conjunction with the VOF interface transport method. Although the traditional HF curvature calculation method can achieve second-order spatial convergence accuracy, the smoothness of the calculated curvature distribution is severely insufficient, especially for coarse background meshes. Non-smooth curvature easily generates spurious, non-physical capillary waves (surface tension waves). On the one hand, spurious capillary waves at the gas-liquid interface do not truly exist and cannot reflect the actual physical problem, severely interfering with the simulation of the atomization process dominated by surface tension. On the other hand, the propagation speed of spurious capillary waves is very fast, severely limiting the computation time step. Forcing a large time step in the actual simulation of the atomization process easily leads to numerical stability problems. Summary of the Invention

[0004] To address the shortcomings of existing technologies, the present invention aims to provide a curvature-corrected volume function method for simulating propellant injection atomization. This method obtains the incremental correction of the HF function to the VOF function, which significantly improves the accuracy and smoothness of local curvature calculation, thereby achieving highly stable numerical simulation of the propellant injection atomization process.

[0005] To solve the above problems, the technical solution of the present invention is as follows:

[0006] A curvature-corrected volume function method for simulating propellant injection atomization includes the following steps:

[0007] The VOF method is used to solve the interface transport equations and update the VOF function field.

[0008] The estimated velocity is calculated using an explicit time-progression method;

[0009] Discretizing surface tension terms based on virtual fluid method;

[0010] Solve the discrete pressure Poisson equation;

[0011] Correcting the velocity field using pressure gradients;

[0012] Determine whether the calculation has ended.

[0013] Preferably, in the step of using the VOF method to solve the interface transport equation and update the VOF function field, the interface transport equation is: Where Γ(x) is a step function representing the distribution of free surfaces within the grid, and its integral within the grid is the VOF volume fraction C.

[0014] Preferably, in the step of calculating the estimated velocity using an explicit time-progression method, the estimated velocity U... * The calculation formula is: in, and These are the convection term, diffusion term, and external force term, respectively. Where M = ρU is the mass flux of the interface element, F g F represents the acceleration due to gravity. σ This represents the surface tension term.

[0015] Preferably, the step of discretizing the surface tension term based on the virtual fluid method specifically includes: discretizing F based on the surface tension model using the virtual fluid method. σ The surface tension effect is simulated in the gradient operator by constructing pressure discontinuities. in, For local curvature, considering the discontinuous pressure gradient due to surface tension, the following is given:

[0016]

[0017] in,

[0018]

[0019] Virtual point pressure that includes both gradient discontinuities and pressure discontinuities. for:

[0020]

[0021] Preferably, the calculation step of the local curvature specifically includes:

[0022] Step 1: Loop through and check if a cell is a UI element;

[0023] Step 2: Use the HF function method to complete the local curvature prediction;

[0024] Step 3: Calculate the local mean curvature field;

[0025] Step 4: Calculate the HF function correction increment based on VPMCF;

[0026] Step 5: Based on the monotonicity of the VOF function column, superimpose the corresponding HF function increment on the VOF function of the interface unit to complete the correction of the VOF function field.

[0027] Step 6: Determine whether the smoothness requirement has been met.

[0028] Preferably, the step of using the HF function method to complete the local curvature prediction specifically includes:

[0029] According to C i,j and C i±1,j±1 Calculate the interface normal vector n = (n x ,n y );

[0030] Determine the main direction of the interface based on |n|;

[0031] The curvature at (i,j) is calculated based on the obtained HF function, and its discrete form is as follows:

[0032] Preferably, the step of superimposing the corresponding HF function increment onto the VOF function of the interface unit based on the monotonicity of the VOF function series to complete the correction of the VOF function field specifically includes:

[0033] like Then, using the HF function in the corresponding column, search in a monotonically decreasing direction to find the first C. m,j ∈[∈,1.-∈] unit, if Then the new VOF function End; if Then C' m,j =1., let Continue searching for the next C m,j For units ∈[0.,1.], perform the operation of superimposing the increments of the HF function until... ∈ is the cutoff value, and it is generally taken as ∈ = 10. -5～-11 ;

[0034] like Then, using the HF function in the corresponding column, search in a monotonically increasing direction to find the first C. m,j ∈[∈,1.-∈] unit, if Then the new VOF function End; if Then C'm,j =0., let Continue searching for the next C m,j For units ∈[0.,1.], perform the operation of superimposing the increments of the HF function until...

[0035] Preferably, the step of solving the discrete pressure Poisson equation specifically includes: constructing the discrete pressure Poisson equation based on the pressure gradient using the virtual fluid method.

[0036] Preferably, the step of correcting the velocity field using a pressure gradient specifically includes: correcting the velocity field based on a pressure gradient constructed using the virtual fluid method, so that the velocity field U... n+1 Satisfy the continuity condition:

[0037] Compared with existing technologies, this invention addresses the problem of local curvature calculation based on the VOF function. It employs the height function method as a local curvature prediction method, deriving the curvature correction equation based on the principle of volume-preserving mean curvature flow (VPMCF). By solving this equation, the incremental correction of the VOF function using the HF function is obtained, significantly improving the accuracy and smoothness of local curvature calculation. This enables highly stable numerical simulation of the propellant injection atomization process. Furthermore, the method of this invention, based on the principle of volume-preserving mean curvature flow, improves the smoothness of local curvature calculation while ensuring the conservation characteristics of the VOF function field. Attached Figure Description

[0038] Other features, objects, and advantages of the present invention will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings:

[0039] Figure 1 Flowchart of the curvature correction volume function method for simulating propellant injection atomization provided in this embodiment of the invention.

[0040] Figure 2 A flowchart illustrating the method for calculating local smooth curvature based on a volume function, as provided in an embodiment of the present invention.

[0041] Figure 3 This is a schematic diagram illustrating the calculation of local curvature of a two-dimensional curve using the height function method, provided in an embodiment of the present invention.

[0042] Figure 4 HF function increment provided in the embodiments of the present invention Corrected VOF function diagram;

[0043] Figure 5The image shows the result of simulating the atomization of a rectangular fuel injection nozzle jet in a rocket engine, as provided in an embodiment of the present invention.

[0044] Figure 6 The droplet size distribution at different times obtained by simulated jet atomization is shown in the embodiment of the present invention. Detailed Implementation

[0045] The present invention will now be described in detail with reference to specific embodiments. These embodiments will help those skilled in the art to further understand the present invention, but do not limit the invention in any way. It should be noted that those skilled in the art can make several changes and improvements without departing from the concept of the present invention. These all fall within the scope of protection of the present invention.

[0046] This invention can be applied to the numerical simulation of the fuel injection atomization process in rocket engines, simulating typical atomization processes such as fuel liquid film breakup, droplet separation, and coalescence driven by surface tension. It is used to evaluate and optimize fuel injection methods to achieve the optimal propellant-oxidizer mixing ratio and ensure complete fuel combustion. Addressing the shortcomings of existing curvature calculation methods based on the HF function, this invention proposes a smooth curvature calculation method based on Volume Preserving Mean Curvature Flow (VPMCF). VPMCF has important applications in surface evolution research, as the volume of the surface remains constant during the evolution of the volume-preserving mean curvature flow. Therefore, the special properties of VPMCF can be utilized to improve the smoothness of local curvature calculations while ensuring the conservation characteristics of the VOF function field.

[0047] Specifically, such as Figure 1 As shown, the curvature correction volume function method for simulating propellant injection atomization includes the following steps:

[0048] S1: Use the VOF method to solve the interface transport equations and update the VOF function field;

[0049] Specifically, the interface transport equation is:

[0050]

[0051] Where Γ(x) is a step function representing the distribution of free surfaces within the grid, and its integral within the grid is the VOF volume fraction C.

[0052] S2: The estimated velocity is calculated using an explicit time-progression method;

[0053] The estimated speed U * The calculation formula is:

[0054]

[0055] in, and These are the convection term, diffusion term, and external force term, respectively.

[0056]

[0057]

[0058] Where M = ρU is the mass flux of the interface element. g F represents the acceleration due to gravity. σ This represents the surface tension term.

[0059] S3: Discretization of surface tension terms based on virtual fluid method;

[0060] Specifically, the surface tension model based on the virtual fluid method is discretized as F σ The surface tension effect is simulated in the gradient operator by constructing pressure discontinuities. in, For local curvature, considering the discontinuous pressure gradient due to surface tension, the following is given:

[0061]

[0062] in,

[0063]

[0064] Virtual point pressure that includes both gradient discontinuities and pressure discontinuities. for:

[0065]

[0066] like Figure 2 As shown, the local curvature Calculated using steps 1-6:

[0067] Step 1: Loop through and check if a cell is a UI element;

[0068] Specifically, the VOF function value C of cell (i,j) is checked in a loop. i,j Whether it is located in [∈, 1.-∈] (where ∈ is the truncation error, usually chosen as ∈ < 10). -6 If the quantity is small, it is marked as an interface unit, on which the curvature estimation method can be implemented.

[0069] Step 2: Use the HF function method to complete the local curvature prediction;

[0070] Specifically, the HF function method is used to perform local curvature prediction. The specific steps for (i,j) lattice curvature prediction are as follows:

[0071] Step 21: According to C i,j and C i±1,j±1 Calculate the interface normal vector n = (n x ,n y );

[0072] Step 22: Determine the main direction of the interface based on |n|;

[0073] Specifically, such as Figure 3 As shown, if x- is the dominant direction, then template 1 is selected, and the HF function is calculated by integrating in the x- direction. If y- is the dominant direction, then template 2 is selected, and the HF function is calculated by integrating in the y- direction.

[0074]

[0075] Step 23: Calculate the curvature of (i,j) based on the obtained HF function. The discrete form is as follows:

[0076]

[0077] Step 3: Calculate the local mean curvature field;

[0078] Specifically, using the predicted curvature obtained in step 2 The local mean curvature field was calculated. The calculation formula is: Where β is the weighting coefficient and L is the number of non-zero curvatures in the weighted calculation template.

[0079]

[0080] Step 4: Calculate the HF function correction increment based on VPMCF;

[0081] Specifically, the predicted curvature obtained in step 2 is based on the concept of VPMCF. A smoothing correction is performed, and the HF function correction increment is calculated. It's important to note that during curvature correction, local curvature... With VOF function C i,j and C i±1,j±1 As a bidirectional tightly coupled system, not just a unidirectional correction of curvature, the iterative correction formula is designed based on the VPMCF principle as follows:

[0082]

[0083] in, n is the interface normal vector, τ is the virtual time step, and σ is the surface tension coefficient. For the current local curvature estimation, Let Γ be the local average curvature, and Γ be a coefficient. Assuming y is the principal normal direction, the equation... The discrete form is:

[0084]

[0085] Step 5: Based on the monotonicity of the VOF function series, superimpose the corresponding HF function increment onto the VOF function of the interface unit to complete the correction of the VOF function field.

[0086] Specifically, the increment of the HF function Once determined, the VOF function field needs to be corrected by superimposing the corresponding HF function increment onto the VOF function of the interface unit based on the monotonicity of the VOF function column. For example... Figure 4 As shown, the specific handling is divided into the following two cases (assuming y- is the principal normal direction):

[0087] Step 51: If Then, using the HF function in the corresponding column, search in a monotonically decreasing direction to find the first C. m,j ∈[∈,1.-∈] unit, if Then the new VOF function End; if Then C' m,j =1., let Continue searching for the next C m,j For units ∈[0.,1.], perform the operation of superimposing the increments of the HF function until... ∈ is the cutoff value, and it is generally taken as ∈ = 10. -5～-11 .

[0088] Step 52: If Then, using the HF function in the corresponding column, search in a monotonically increasing direction to find the first C. m,j ∈[∈,1.-∈] unit, if Then the new VOF function End; if Then C' m,j =0., let Continue searching for the next C m,j For units ∈[0.,1.], perform the operation of superimposing the increments of the HF function until...

[0089] Step 6: Determine whether the smoothness requirement has been met.

[0090] Specifically, calculation If the convergence error requirement is met, the process ends; otherwise, return to step 2 to continue iterative calculation.

[0091] S4: Solve the discrete pressure Poisson equation;

[0092] Specifically, the pressure gradient discrete pressure Poisson equation is constructed using the S3-step virtual fluid method.

[0093]

[0094] S5: Correcting the velocity field using the pressure gradient;

[0095] The velocity field is corrected based on the pressure gradient constructed using the virtual fluid method in step S3, making the velocity field U... n+1 Satisfy the continuity condition:

[0096] S6: Determine if the calculation has ended.

[0097] Specifically, determine whether the calculation has ended. If not, return to step S1; if yes, the calculation ends.

[0098] Figure 5 To illustrate the application of this invention, a simulation diagram of propellant jet atomization at a rectangular fuel injection nozzle in a rocket engine is provided. Figure 6 The droplet size distribution at different times, obtained by simulating jet atomization using this invention, is presented. As shown in the figure, with time, atomization becomes more complete, and the droplet size distribution gradually approaches a power-law distribution of -4.5, which is basically consistent with experimental results. Table 1 below shows the geometric characteristics of the fuel injection port and the two-phase flow fluid characteristic parameters.

[0099]

[0100] Table 1

[0101] Compared with existing technologies, this invention addresses the problem of local curvature calculation based on the VOF function. It employs the height function method as a local curvature prediction method, deriving a curvature correction equation based on the principle of volume conservation in average curvature flow. Solving this equation yields an incrementally corrected VOF function from the HF function, significantly improving the accuracy and smoothness of local curvature calculation. This enables highly stable numerical simulation of the propellant injection atomization process. Furthermore, the method of this invention, based on the principle of volume conservation in average curvature flow, improves the smoothness of local curvature calculation while ensuring the conservation characteristics of the VOF function field.

[0102] Specific embodiments of the present invention have been described above. It should be understood that the present invention is not limited to the specific embodiments described above, and those skilled in the art can make various changes or modifications within the scope of the claims, which do not affect the essence of the present invention. Unless otherwise specified, the embodiments and features described in this application can be arbitrarily combined with each other.

Claims

1. A method for simulating the curvature-corrected volume function of propellant injection atomization, characterized in that, The method includes the following steps: The VOF method is used to solve the interface transport equations and update the VOF function field. The estimated velocity is calculated using an explicit time-progression method; Discretizing surface tension terms based on the virtual fluid method specifically includes: discretizing the surface tension model based on the virtual fluid method. The surface tension effect is simulated in the gradient operator by constructing pressure discontinuities. in, For local curvature, considering the discontinuous pressure gradient due to surface tension, the following is given: in, Virtual point pressure that includes both gradient discontinuities and pressure discontinuities. for: ； The calculation steps for the local curvature specifically include: Step 1: Loop through and check if a cell is a UI element; Step 2: Use the HF function method to complete the local curvature prediction, specifically including: according to and Calculate the interface normal vector ; according to Determine the main orientation of the interface; Calculate based on the obtained HF function ( The curvature, in discrete form, is: ; Step 3: Calculate the local mean curvature field; Step 4: Calculate the HF function correction increment based on VPMCF; Step 5: Based on the monotonicity of the VOF function series, superimpose the corresponding HF function increment onto the VOF function of the interface unit to complete the correction of the VOF function field. Specifically, this includes: like Then, starting from the corresponding column, the HF function searches in a monotonically decreasing direction to find the first one. Unit, if Then the new VOF function End; if ,but ,make Continue searching for the next one The unit performs the operation of superimposing the increments of the HF function until... ; To be the cutoff value, take ; like Then, starting from the corresponding column, the HF function searches in a monotonically increasing direction to find the first... Unit, if Then the new VOF function End; if ,but ,make Continue searching for the next one The unit performs the operation of superimposing the increments of the HF function until... ; Step 6: Determine if the smoothness requirement has been met; Solve the discrete pressure Poisson equation; Correcting the velocity field using pressure gradients; Determine whether the calculation has ended.

2. The curvature correction volume function method for simulating propellant injection atomization according to claim 1, characterized in that, In the step of solving the interface transport equation and updating the VOF function field using the VOF method, the interface transport equation is: ,in The step function representing the distribution of free surfaces within the grid has its integral within the grid as the VOF volume fraction. .

3. The curvature correction volume function method for simulating propellant injection atomization according to claim 1, characterized in that, In the step of calculating the estimated velocity using an explicit time-progression method, the estimated velocity... The calculation formula is: ,in, , and These are the convection term, diffusion term, and external force term, respectively. , ,in, For the mass flux of the interface unit, Represents the acceleration due to gravity. This represents the surface tension term.

4. The curvature-corrected volume function method for simulating propellant injection atomization according to claim 1, characterized in that, The steps for solving the discrete pressure Poisson equation specifically include: constructing the discrete pressure Poisson equation based on the pressure gradient using the virtual fluid method. .

5. The curvature correction volume function method for simulating propellant injection atomization according to claim 1, characterized in that, The step of correcting the velocity field using the pressure gradient specifically includes: correcting the velocity field according to the pressure gradient constructed by the virtual fluid method, so that the velocity field... Satisfy the continuity condition: .

Images