Multi-component multiphase flow numerical simulation method based on prediction-correction lattice Boltzmann algorithm

Abstract

The invention discloses a multicomponent multiphase flow numerical simulation method based on a prediction-correction lattice Boltzmann algorithm, which comprises the following steps: solving a macroscopic equation based on an equilibrium state distribution function by adopting a prediction-correction step-by-step technology, and calculating an intermediate macroscopic variable through the equilibrium state distribution function in a prediction step; calculating to obtain an equilibrium state distribution function at an intermediate moment according to the intermediate macroscopic variable in the prediction step, and calculating to obtain a non-equilibrium item in combination with an original equilibrium state distribution function; substituting the macroscopic variable at the middle moment and the macroscopic variable at the previous time step into a discrete source item calculation formula, and respectively calculating an interface discrete source item and a flow field discrete force item; and correcting the macroscopic variable at the middle moment by using the unbalanced item and the discrete source item. Compared with a traditional lattice Boltzmann multiphase flow algorithm, the method has the advantages that the boundary processing format is simplified, the virtual memory occupation is reduced, the multiphase flow interface tracking precision is improved, and the method can be used for simulating and analyzing the large-density-ratio multicomponent multiphase flow problem.

Description

Technical Field

[0001] The invention relates to a multi-component multi-phase flow numerical simulation method, in particular to a multi-component multi-phase flow numerical simulation method based on a prediction-corrected lattice Boltzmann algorithm. Background Art

[0002] Multiphase fluid generally refers to a fluid system with two or more phases and a distinct interface. Multicomponent multiphase flow is a fluid mixture containing multiple types of substances and with a distinct interface. It involves the interaction between fluid dynamics and multiphases and is the most complex multiphase flow system. Multicomponent multiphase flow phenomena are widely present in life, nature, and many industrial fields such as energy, environment, and chemical industry, such as oil leaks in the ocean and the injection of carbon dioxide and other gases in oil production to increase the production rate. Due to the complexity of multicomponent multiphase flow, most of the existing research on "multiphase" fluid flow is about two-phase or three-phase fluid flow, and there are few research reports on multicomponent multiphase flow (equal to or more than four phases). In recent years, with the development and improvement of the theoretical system of the two-phase flow phase field-lattice Boltzmann model, the development of the phase field lattice Boltzmann multicomponent multiphase flow model has gradually been on the right track. Although the lattice Boltzmann method based on the phase field has achieved good results in the field of multiphase flow, these models still have major shortcomings. Among these problems, the most difficult one is the stability when dealing with multiphase flow problems with a density ratio greater than 1000. The simulation of multiphase flow problems with a large density ratio requires particularly high stability of the model. Because of this, the existing models can only simulate multiphase flow problems with a density ratio of about 10. Although there are many effective numerical models in the field of two-phase flow with a large density ratio, there is still a lack of numerical methods that can handle density ratios greater than 1000 in the field of multiphase flow numerical simulation. In addition, the existing phase field lattice Boltzmann multi-component multiphase flow model uses the distribution function as the evolution object, while the multi-component multiphase flow problem requires solving multiple physical quantities and phases. For each additional phase in the fluid system, a distribution function needs to be added to track the changes in the phase interface, which results in more virtual memory required in the calculation process, and the requirements for computing resources for solving complex multi-component multiphase flows are higher. Taking the D2Q9 velocity model as an example, in the case of extreme programming, 18 distribution functions need to be stored at each grid point for two-phase flow problems, 27 for three-phase flow, and so on. Therefore, the existing phase field lattice Boltzmann multi-component multi-phase flow model consumes a lot of virtual memory during calculation, which is not conducive to the numerical simulation analysis of actual engineering problems. On the other hand, after each time step is calculated, the distribution function inside the calculation area has been calculated, but there are still some unknown distribution functions on the boundary. Therefore, it is necessary to first deal with the evolution format of the unknown distribution function on the boundary to realize the boundary conditions. Designing a suitable distribution function boundary evolution format based on the physical problem to be solved is a critical and difficult task, and it will also seriously affect the stability of the model. Therefore, it is of great theoretical and practical significance to develop a lattice Boltzmann multi-component multi-phase flow model that can inherit the mesoscopic properties of the lattice Boltzmann method, optimize its shortcomings in virtual memory consumption and boundary condition realization, and stabilize the large density ratio multi-phase flow problem. Summary of the invention

[0003] Purpose of the invention: The purpose of the invention is to propose a multi-component multi-phase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm to realize the numerical simulation of multi-component multi-phase flow problems with large density ratios.

[0004] Technical solution: The present invention comprises the following steps:

[0005] The prediction-correction step-by-step technique is used to solve the macro equation based on the equilibrium distribution function, and the intermediate macro variables are calculated by the equilibrium distribution function in the prediction step;

[0006] The equilibrium distribution function at the intermediate moment is calculated based on the intermediate macro variables in the prediction step, and the non-equilibrium term is calculated by combining the original equilibrium distribution function;

[0007] The macro variables at the intermediate moment and the macro variables at the previous time step are introduced into the discrete source term calculation formula to calculate the interface discrete source term and the flow field discrete force term respectively;

[0008] The non-equilibrium terms and discrete source terms are used to correct the macro variables at intermediate moments.

[0009] The intermediate variable pressure p in the prediction step * , Sequence Parameters Speed * The solution process is:

[0010]

[0011] Where t is the current solution time, t-δ t is the previous moment, x is a grid point on the grid, δ t is the time interval between adjacent evolution steps, Among them, ρ i is the density of each fluid;

[0012] The equilibrium distribution function of the pressure is solved by using the density quantity ρ(x,t-δ t ) calculation, that is

[0013]

[0014] p ∨ =p(xe α δ t ,t-δ t )

[0015]

[0016] u ∨ =u(xe α δ t,t-δ t )

[0017] Where p, ρ and u are the pressure, density and velocity at the previous moment respectively; ω α is the weight coefficient in the grid direction; e α is the discrete velocity direction.

[0018] The interface equilibrium distribution function for solving the order parameter for:

[0019]

[0020] In the formula, φ i is the sequence parameter of phase i at the previous moment.

[0021] The velocity distribution function of the flow field equilibrium state is for:

[0022]

[0023] In the formula, ρ ∨ =ρ(xe α δ t ,t-δ t ) is the mesh density after migration.

[0024] The equilibrium distribution function at the intermediate moment and The calculation formulas are:

[0025]

[0026] where Γ α * The expression is:

[0027]

[0028] The non-equilibrium term The calculation formula is:

[0029]

[0030] In the formula, τ f is the dimensionless relaxation parameter of the flow field.

[0031] The interface discrete source term and the discrete force term of the flow field The calculation formula is:

[0032]

[0033] in,

[0034]

[0035]

[0036] Where n i is the unit gradient direction of the interface, W is the interface thickness, F s is the surface tension term.

[0037] The non-equilibrium term and the discrete source term are used to correct the macro variables at the intermediate moment, specifically:

[0038] The prediction-correction step-by-step technique is used to solve the macro equation based on the equilibrium distribution function restored by the Chapman-Enskog multiscale expansion analysis. In the correction step, the intermediate variable pressure p in the prediction step is * , Sequence Parameters and speed u * Further corrections.

[0039] The correction process is:

[0040] p=p *

[0041]

[0042] In the formula, τ g is the dimensionless relaxation parameter of the phase field; F b is the volume force.

[0043] Beneficial effects: The present invention has the following advantages:

[0044] (1) In the prediction step, the intermediate macro variables are calculated by the equilibrium distribution function, and in the correction step, the intermediate macro variables are corrected by using the non-equilibrium part and the source term to achieve stable and effective tracking of the interface and flow field of the multi-component multiphase flow;

[0045] (2) The equilibrium distribution function, non-equilibrium part and source term involved in the present invention can be directly calculated by macroscopic variables, and the direct evolution of macroscopic quantities such as order parameter and velocity is realized in the system of lattice Boltzmann method. The whole calculation process does not need to rely on the distribution function, and the direct evolution of macroscopic quantities is realized in the system of lattice Boltzmann.

[0046] (3) The present invention inherits the kinetic basis and mesoscopic characteristics of the lattice Boltzmann method, and improves the defects of the traditional lattice Boltzmann multiphase flow model, such as large virtual memory consumption, complex boundary condition implementation, and poor numerical calculation stability. It effectively improves the stability of numerical simulation calculations, greatly reduces the consumption of virtual memory in the calculation process, simplifies the distribution function processing steps on the boundary, improves the tracking accuracy of the multi-component multi-phase flow interface, and achieves the goal of stable and accurate simulation of large density ratio multi-component multi-phase flow problems. BRIEF DESCRIPTION OF THE DRAWINGS

[0047] Figure 1 is a flow chart of the present invention;

[0048] Figure 2 A schematic diagram of a simulation of a four-phase static droplet provided for an example of the present invention;

[0049] Figure 3 Comparison between simulated and theoretical values ​​of the distribution of the four order parameters of the horizontal central fluid;

[0050] Figure 4 A distribution diagram of non-physical measurement indicators of a fluid on the interface of three other fluids provided in an example of the present invention;

[0051] Figure 5 A schematic diagram of the simulation of Laplace's law provided for an example of the present invention;

[0052] Figure 6 A schematic diagram of the deformation of a three-phase single lens provided by an example of the present invention;

[0053] Figure 7 The distribution of various phases after the deformation of a single droplet lens stabilizes under different surface tension coefficients provided in the examples of the present invention;

[0054] Figure 8 A schematic diagram of the deformation of a four-phase flow double lens provided by an example of the present invention;

[0055] Fig. 9 The distribution of various phases after the double droplet lens is deformed and stabilized under different surface tension coefficients provided in the examples of the present invention. DETAILED DESCRIPTION

[0056] The present invention will be further described below in conjunction with the accompanying drawings.

[0057] The multi-component multiphase flow numerical simulation method based on the prediction-corrected lattice Boltzmann algorithm of this embodiment first calculates the macro variables at the intermediate moment according to the relationship between the equilibrium distribution function and the lattice tensor in the prediction step; further, the equilibrium distribution function at the intermediate moment is calculated according to the intermediate macro variables in the prediction step, and the non-equilibrium term is calculated by combining the original equilibrium distribution function; secondly, the macro variables at the intermediate moment and the macro variables at the previous time step are introduced into the discrete source term calculation formula to calculate the discrete source term; finally, the macro variables at the intermediate moment are corrected by using the non-equilibrium term and the discrete source term to achieve effective tracking of the interface and stable evolution of the flow field. Figure 1 As shown, the specific steps include:

[0058] S1. The prediction-correction step-by-step technique is used to solve the macro equation based on the equilibrium distribution function restored by the Chapman-Enskog multiscale expansion analysis. The intermediate macro variables are calculated by the equilibrium distribution function in the prediction step. The pressure p in the prediction step is * , Sequence Parameters and speed u * The solution process of intermediate variables is:

[0059]

[0060]

[0061] Where t is the current solution time, t-δ t is the previous moment, x is a grid point on the grid, δ t is the time interval between adjacent evolution steps, Among them, ρ i is the density of each fluid;

[0062] In the above formula, the density quantity ρ(x,t-δ) of the current grid is used in calculating the equilibrium distribution function of pressure t ),Right now

[0063]

[0064] p ∨ =p(xe α δ t ,t-δ t )

[0065]

[0066] u ˇ =u(xe α δ t ,t-δ t )

[0067] Where p, ρ and u are the pressure, density and velocity at the previous moment respectively; the weight coefficient ω in the grid direction is α and the discrete velocity direction e α The discrete model DdQq that depends on the lattice velocity, d is the spatial dimension, q is the number of discrete velocities, and this embodiment adopts the D2Q9 model, whose discrete velocities and weight coefficients are:

[0068]

[0069] Where c is the grid velocity, defined as the grid length δ x and time step δ t The ratio c = δ x / δ t , lattice sound speed

[0070] Solving the interface equilibrium distribution function of the intermediate order parameter And solve the equilibrium distribution function of the flow field at intermediate speeds They are:

[0071]

[0072] In the formula, φ i is the order parameter of phase i at the previous moment; ρ ˇ =ρ(xe α δ t ,t-δ t ) is the density after migration.

[0073] S2. The equilibrium distribution function at the intermediate moment is calculated based on the intermediate macro variables in the prediction step, and the non-equilibrium term is calculated by combining the original equilibrium distribution function.

[0074] Among them, the equilibrium distribution function at the intermediate time and The calculation formula is as follows:

[0075]

[0076] where Γ α * The expression is:

[0077]

[0078] Unbalanced Items The calculation formula is:

[0079]

[0080]

[0081] In the formula, τ f is the dimensionless relaxation parameter of the flow field, and its calculation is related to the fluid kinematic viscosity ν,

[0082] S3. Bring the macro variables at the intermediate moment and the macro variables at the previous time step into the discrete source term calculation formula to calculate the interface discrete source terms respectively. and the discrete force term of the flow field

[0083] in,

[0084]

[0085] Where n i is the unit gradient direction of the interface; W is the interface thickness; F s is the surface tension term, which is calculated based on the total free energy functional using the virtual work method:

[0086]

[0087] Where σ i,j represents the surface tension between fluid i and fluid j.

[0088] χ i is the indicator function, which is defined as

[0089]

[0090] Phase field gradient Calculated using second-order central difference:

[0091]

[0092] Laplace operator The calculation is:

[0093]

[0094] S4. Use non-equilibrium terms and discrete source terms to correct the macro variables at intermediate moments to achieve effective tracking of the interface and stable evolution of the flow field.

[0095] Specifically, the prediction-correction step-by-step technique is used to solve the macro equation based on the equilibrium distribution function restored by the Chapman-Enskog multiscale expansion analysis. In the correction step, the pressure p in the prediction step is * , Sequence Parameters and speed u * The intermediate variables are further corrected, and the correction process is:

[0096] p=p *

[0097]

[0098] In the formula, F b is the body force, and the dimensionless relaxation parameter of the phase field can be calculated based on Calculate; θ M is the mobility coefficient.

[0099] The application of the present invention to the problem of multi-component multi-phase flow with large density ratio is further described below with reference to specific examples.

[0100] Reference Example 1: Four-phase static droplet

[0101] The present invention relates to a multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm. First, the accuracy and precision of the model in tracking the interface between different fluids are verified. Figure 2 As shown in the figure, there are three stationary droplets of different types in the fluid. In a region with a length of L and a width of H, the fluid is filled, and three other droplets of different types but with a radius of R are infiltrated in it. The radius of the three droplets is R = 0.1L, and their respective positions are: droplet 1 (H / 4, H / 2), droplet 2 (3H / 4, H / 2), droplet 3 (5H / 4, H / 2). In the simulation, L = 300, H = 200, and periodic boundary conditions are applied to all boundaries. The mobility is set to θ m =0.01, that is, τ g =0.503. It should be noted that such a small mobility has a great impact on the stability of the standard lattice Boltzmann method, which will cause the simulation results to diverge. In order to verify the ability of the simplified LB model proposed in this invention to handle multi-component multiphase flow problems with large density ratios, the density ratio between the fluids is set to ρ1:ρ2:ρ3:ρ4=1000:500:300:1. The tension coefficient matrix between the fluids σ=[σ i,j ] 4×4 Set to

[0102]

[0103] In order to quantitatively analyze this problem, Figure 3 The distribution of the order parameter φ4 of fluid 4 along the horizontal center line (y = H / 2) is plotted, and the initial distribution is given as a theoretical solution for comparison. It can be seen that the calculated results are in good agreement with the initial values. Furthermore, an indicator variable φ is introduced x =φ1φ2φ3φ4, as a measurement indicator to test the non-physical existence of one fluid at the interface of three other fluids. From a physical point of view, in an N-phase fluid system, except for the intersection point where N fluids coexist, φ xMust be zero throughout the region. Figure 4 The calculation area φ is shown in x As shown in the figure, φ x The maximum value is 10 -7 Such a small magnitude. It should be noted that the simulation accuracy of existing models for such problems is at 10 -3 It can be seen that the present invention has good accuracy in interface tracking and processing the boundary area between interfaces.

[0104] Reference Example 2: Laplace's Law

[0105] In a multi-component multi-phase flow system, the interaction between the surface tensions of the phases is more complex, so the model's ability to handle surface tension is required to be higher. The present invention further simulates the problem of three interlocking droplets immersed in another fluid to verify the model's ability to handle surface tension. Figure 5 As shown in Figure 1, the radius length ratio between the droplets is R1:R2:R3=1:2:3, where R1=R0. A grid of L×L=200×200 is used in the simulation. Periodic boundary conditions are applied to all boundaries. The density ratio and viscosity ratio are set to ρ1:ρ2:ρ3:=ρ4 1:3 and μ1:μ2:μ3:μ4=1:3:5:10, respectively. In addition, the other parameters are set to: θ m =0.05, μ1=0.01, W=5 and R0=30. Surface tension matrix σ=[σ i,j ] 4×4 Set to

[0106]

[0107] The centers of the three droplets are all located at (x c ,y c )=(L / 2,L / 2). In equilibrium, the droplet remains stationary and maintains its initial round shape. The pressure difference Δp between the inside and outside of the adjacent droplet interface i,j and surface tension coefficient σ i,j Satisfies Laplace's law

[0108] Δp i,j =σ i,j / R i,j

[0109] Among them, R i,j represents the radius of curvature between fluid i and fluid j.

[0110] Table 1 quantitatively compares the calculated values ​​of the pressure difference between the interface and the fluid phases at different initial radii, Δp i,j and the theoretical value R i,j σi,j , and the relative error between the two is listed. As can be seen from the table, the relative errors of all calculation results are within 4.4%, which well illustrates that the present invention can accurately express the surface tension between multi-component and multi-phase.

[0111] Table 1 Relative errors of the calculated results and analytical solutions of the internal and external pressure differences between the interfaces of the fluid phases at different initial radii

[0112]

[0113] Reference Example 3: Deformation of a single lens in three-phase flow

[0114] The problems discussed in the above two examples involve two intersecting interfaces, and there is no complex situation of the intersection of multiple phase interfaces. In this example, the present invention will simulate the problem of droplet lens deformation to test the model's ability to handle the problem of the intersection of multiple phase interfaces. When a droplet is located between the interfaces of other multiphase fluids, the droplet will deform under the action of surface tension and finally reach an equilibrium state. Its deformation is determined by the size of the surface tension between the phases. This example will simulate a three-phase flow problem, that is, a droplet is immersed in the interface of two other fluids. The schematic diagram of the problem is shown in Figure 6 In the experiment, when a droplet is placed on the interface of two other immiscible fluids, an area where the three fluids intersect will appear. In this area, the three-phase fluids will intersect at one point, namely the three-phase contact point. Due to the surface tension between the fluids, the immersed droplet will be deformed, and the final equilibrium state will be composed of two arcs. In the equilibrium state, static equilibrium is maintained near the three-phase contact point. At the contact point, the interface can be approximated as a straight line, and the angles between the interfaces are θ1, θ2 and θ3 respectively. In the equilibrium state, the distance d between the three-phase contact points is controlled by the surface tension coefficient and the interface contact angle, and its theoretical calculation expression is:

[0115]

[0116] In the above formula, A represents the area of ​​the droplet lens

[0117] It should be noted that the three-phase contact point needs to be formed under the following conditions:

[0118] Σ i =σ ij +σ ik -σ jk >0.

[0119] If the surface tension coefficient does not satisfy the above formula, the three fluids will separate and there will be no three-phase contact angle. In the simulation of this problem, a grid of L0×L0=301×301 is used. In the simulation, the relevant parameters are set as follows: μ1=μ2=μ3=0.025, θ m = 0.05, R0 = 32, W = 4. The deformation problem of the droplet lens is independent of the density ratio between the fluids, but in order to test the ability of the present invention to handle large density ratio problems, the density ratio of each phase is set to ρ1:ρ2:ρ3 = 1:1.4:0.0014.

[0120] In equilibrium, the ratio of the four surface tension coefficients (σ 12 :σ 13 :σ 23 )=(1:1:1),(σ 12 :σ 13 :σ 23 )=(1:1.4:0.8),(σ 12 :σ 13 :σ 23 )=(1:0.8:1.4) and (σ 12 :σ 13 :σ 23 )=(1:1:3) when the droplet lens deformation is shown in Figure 7 The comparison between the simulated and theoretical values ​​of the contact point spacing is shown in Table 2. It can be seen that the simulated results are consistent with the theoretical results, with a relative error of only about 1%, demonstrating the ability of this model to handle multiphase flows with large density ratios.

[0121] Table 2 Simulated and theoretical values ​​and relative errors of the distance between the three-phase contact points after the single droplet lens is stabilized

[0122]

[0123] Reference Example 4: Four-phase flow double lens deformation

[0124] This example will simulate a four-phase flow problem, that is, two different types of droplets infiltrated at the interface of two other fluids. The geometric parameters and contact angle definitions of this problem are as follows Figure 8 To ensure the regular deformation of the droplet lens and the measurement of related parameters, σ is set 13 =σ 14 =σ 23 =σ 24 Under this condition, the droplet lens is symmetrically distributed along the horizontal center and the vertical center when in equilibrium. When the system reaches equilibrium, the distance d between the left and right contact points and the distance h between the upper and lower contact points have the following theoretical solution:

[0125] d=R(sinθ1+cosθ2),

[0126] h=R(sinθ2-cosθ1),

[0127] The expression of the theoretical value of the contact angle is:

[0128]

[0129] R is the radius of curvature of the stabilized rear lens, and its calculation formula is:

[0130]

[0131] In the above formula, S p Represents the area of ​​fluid 1.

[0132] For the simulation of this problem, other relevant parameters are set as follows: M = 0.05, W = 4, μ1 = 0.01, the density ratio and viscosity ratio between the phases are ρ1:ρ2:ρ3:ρ4 = 1:3:50:1000 and μ1:μ2:μ3:μ4 = 1:3:50:1000 respectively. This example will simulate and analyze four different tension coefficients. The specific tension coefficient matrix σ = [σ i,j ] 4×4 The settings are

[0133] Calculation example 1:

[0134] Calculation example 2:

[0135] Calculation example 3:

[0136] Calculation example 4:

[0137] Fig. 9 The density cloud map under equilibrium state is shown in Figure 3. It can be seen that the entire interface shape is symmetrically distributed. Table 3 compares the simulated values ​​and theoretical values ​​of the horizontal distance and vertical distance of the contact point in different cases in detail, and gives the relative error between the simulated value and the theoretical value. It can be seen that the maximum relative error between the simulated value and the theoretical value is only 4.64%, which proves that this model can effectively simulate multi-component multi-phase flow.

[0138] Table 3 Simulated and theoretical values ​​and relative errors of the distance between the three-phase contact points after the double droplet lens is stabilized

[0139]

[0140] The present invention first obtains the intermediate macroscopic quantity through the equilibrium distribution function in the prediction step, and then makes corrections in the correction step to obtain a stable and accurate numerical value. The entire calculation process only involves the equilibrium distribution function and the non-equilibrium part, wherein the equilibrium distribution function can be directly calculated from the macroscopic quantity according to the formula, and the non-equilibrium part is also calculated by two equilibrium distribution functions adjacent in time and space steps. Therefore, the present invention realizes the direct evolution of macroscopic quantities in the lattice Boltzmann system. During the calculation process, only the macroscopic quantities need to be stored, which directly eliminates the large part of the virtual memory occupied by the distribution function, greatly improving the efficiency of the computer. Macroscopic physical boundary conditions can be directly realized through macroscopic quantities, which greatly simplifies the boundary condition processing format. The present invention combines the mesoscopic characteristics of the lattice Boltzmann method and has high numerical stability, providing a stable and reliable numerical simulation tool for the study of multi-component multi-phase flow problems with large density ratios.

Claims

1. A multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm, characterized in that: The following steps are involved: The prediction-correction step-by-step technique is used to solve the macro equation based on the equilibrium distribution function, and the intermediate macro variables are calculated by the equilibrium distribution function in the prediction step; The equilibrium distribution function at the intermediate moment is calculated based on the intermediate macro variables in the prediction step, and the non-equilibrium term is calculated by combining the original equilibrium distribution function; The macro variables at the intermediate moment and the macro variables at the previous time step are introduced into the discrete source term calculation formula to calculate the interface discrete source term and the flow field discrete force term respectively; The non-equilibrium terms and discrete source terms are used to correct the macro variables at intermediate moments.

2. The method for numerical simulation of multi-component multi-phase flow based on the predicted-corrected lattice Boltzmann algorithm according to claim 1, characterized in that: The intermediate variable pressure p in the prediction step * , Sequence Parameters Speed * The solution process is: Where t is the current solution time, t-δ t is the previous moment, x is a grid point on the grid, δ t is the time interval between adjacent evolution steps, Among them, ρ i is the density of each fluid.

3. The multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 2 is characterized in that: The equilibrium distribution function of the pressure is solved by using the density quantity ρ(x,t-δ t ) calculation, that is p ∨ =p(xe α d t ,t-d t ) you ∨ =u(xe α d t ,t-d t ) Where p, ρ and u are the pressure, density and velocity at the previous moment respectively; ω α is the weight coefficient in the grid direction; e α is the discrete velocity direction.

4. The multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 2 is characterized in that: Solving the interface equilibrium distribution function for order parameters for: In the formula, φ i is the sequence parameter of phase i at the previous moment.

5. The multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 2, characterized in that: Solve the velocity equilibrium distribution function of the flow field for: In the formula, ρ ∨ =ρ(xe α δ t ,t-δ t ) is the mesh density after migration.

6. The multi-component multi-phase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 2, characterized in that: The equilibrium distribution function at the intermediate moment and The calculation formulas are: in The expression is:

7. The multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 6, characterized in that: The non-equilibrium term The calculation formula is: In the formula, τ f is the dimensionless relaxation parameter of the flow field.

8. The multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 7, characterized in that: The interface discrete source term and the discrete force term of the flow field The calculation formula is: in, Where n i is the unit gradient direction of the interface, W is the interface thickness, F s is the surface tension term.

9. The multi-component multiphase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 8, characterized in that: The non-equilibrium term and the discrete source term are used to correct the macro variables at the intermediate moment, specifically: The prediction-correction step-by-step technique is used to solve the macro equation based on the equilibrium distribution function restored by the Chapman-Enskog multiscale expansion analysis. In the correction step, the intermediate variable pressure p in the prediction step is * , Sequence Parameters and speed u * Further corrections.

10. The multi-component multi-phase flow numerical simulation method based on the predicted-corrected lattice Boltzmann algorithm according to claim 9, characterized in that: The correction process is: p=p * In the formula, τ g is the dimensionless relaxation parameter of the phase field; F b is the volume force.

Images