Treatment method for curved boundaries in numerical simulation using the lattice Boltzmann method

Abstract

The present invention proposes a method for processing curved boundaries in numerical simulations using the lattice Boltzmann method. This method comprehensively combines the Mei format, Chen format, and Guo format based on the actual numerical simulation calculation requirements, and can be used to develop a computational platform for high-precision and high-efficiency numerical simulations of complex flow fields, for practical application fields such as high-performance airfoil design and new-generation engine design, promoting the development and innovation of this discipline. The present invention is more accurate than the three existing curved boundary processing formats. Since each step of the calculation is more accurate, the convergence speed is also faster, achieving high precision and high efficiency in the calculation.

Description

Technical Field

[0001] The present invention relates to the technical field of fluid mechanics, and specifically to a method for processing the physical boundary of relevant numerical simulations in computational fluid dynamics, which is used for numerical simulations such as airfoil, cylinder flow around, and internal flow in engines. Background Art

[0002] The lattice Boltzmann method is a mesoscopic computational fluid dynamics algorithm that solves the distribution function based on the collision and migration theory of fluid molecules at the grid points, so as to obtain the macroscopic physical quantities at the grid points, such as velocity, density, pressure, etc. In actual numerical simulations, for the curved boundary grid points (such as Figure 1 shown) commonly found in complex shapes such as airfoils, cylinder flow around, and internal flow in engines, there is no way to participate in the evolution process (collision and migration) of the entire flow field. Therefore, interpolation formats are usually used to manually assign values to these curved boundary points (curved boundary processing). There are three relatively popular existing processing formats, namely the Mei format, the Chen format, and the Guo format. Among them, the Mei format and the Guo format divide the distribution function into an equilibrium state and a non-equilibrium state, and then interpolate to solve the distribution function of the physical surface points by creating virtual points. Its disadvantages are as follows: First, it does not conform to physical reality because the positions of the virtual points are actually inside the object, and there is no migration and collision of fluid molecules; second, the interpolation formats of the physical quantities (velocity and density) of the virtual points are relatively simple and do not contain enough information. The Chen format directly solves the distribution function of the physical surface points through linear interpolation of fluid points. The interpolation format is also relatively simple, does not contain enough information near the physical surface points, and does not process the distribution function according to the non-equilibrium state and the equilibrium state, so it is not accurate enough at the physical level. The above three curved boundary formats perform excellently in dealing with general numerical simulation problems, but there is still great room for optimization or development in terms of pursuing high-precision and high-efficiency numerical simulations, especially in high-performance airfoil design, new-generation engine design, etc. In order to meet the actual application requirements of increasingly demanding airfoil flow around numerical simulations and internal flow numerical simulations in the engine cavity, and to improve the numerical simulation accuracy and efficiency for specific engineering application backgrounds, it is necessary to develop High precision and high efficiency Efficiency computational methods. Summary of the Invention

[0003] To solve the problems existing in the prior art, meet the requirements of high-precision and high-efficiency numerical simulation while conforming to the physical reality, the present invention comprehensively combines the above three surface treatment formats based on the actual numerical simulation calculation requirements, and proposes a new surface boundary optimization processing method applicable to LBM numerical simulation, which can be used to develop a calculation platform for high-precision and high-efficiency numerical simulation of complex flow fields. Based on this method, corresponding numerical simulation methods for airfoil flow around and internal flow in the engine cavity are further proposed, and are applied to practical application fields such as high-performance airfoil design and new-generation engine design, improving the accuracy and efficiency of numerical simulation of airfoil flow around and internal flow in the engine cavity, and promoting the development and innovation of this discipline.

[0004] The technical solution of the present invention is as follows:

[0005] A method for processing the curved surface boundary in the numerical simulation of the lattice Boltzmann method, and the distribution function of point f along the direction on the physical surface boundary is solved through the following process:

[0006] Step 1: Divide the distribution function of point f along the direction into the equilibrium distribution function and the non-equilibrium distribution function :

[0007]

[0008] Among them, the equilibrium distribution function is obtained through the calculation formula of the equilibrium distribution function ;

[0009] Step 2: Solve the non-equilibrium distribution function of point f by the following steps:

[0010] Step 2.1: Use linear interpolation to synthesize the information of points ff and w

[0011]

[0012] where is the embedding depth;

[0013] Step 2.2: Use linear extrapolation to synthesize the information of points fff and ff

[0014]

[0015] Step 2.3: Average the results of linear interpolation and extrapolation to obtain:

[0016]

[0017] Among them: The information of w is calculated by the following process :

[0018] The distribution function of point w is divided into two parts: non-equilibrium state and equilibrium state

[0019]

[0020] It is determined by the no-slip condition on the object surface , and we get ; among them, the equilibrium state distribution function of point w is obtained through the calculation formula of the equilibrium state distribution function ;

[0021] At the same time, linear interpolation inside and outside is used to solve the distribution function of point w

[0022]

[0023]

[0024] Thus, we get , where the distribution function of point b along direction is migrated from the distribution function of point f along direction after collision migration. The distribution function of point f along direction is migrated from the distribution function of point ff along direction after collision migration. The distribution function of point ff along direction is migrated from the distribution function of point fff along direction after collision migration;

[0025] Thus, we get:

[0026]

[0027] Furthermore, we get

[0028] .

[0029] Furthermore, in step 1, when calculating the equilibrium state distribution function , the density and velocity of point f are interpolated by the points around point f:

[0030] .

[0031] Furthermore, in step 3, the information of point ff and f is used to interpolate the velocity and density of point w, and then the equilibrium state distribution function of point w is calculated through the equilibrium state distribution function .

[0032] Beneficial effects

[0033] The method for handling curved boundaries in the numerical simulation of the lattice Boltzmann method proposed by the present invention can meet the problems of handling curved boundaries encountered in the numerical simulation of complex flow fields and complex shapes by LBM. At the same time, it is more accurate than the above three existing curved boundary handling formats. Since each calculation step is more accurate, the convergence speed is also faster, achieving high precision and high efficiency in numerical simulation. In practical applications such as airfoil design and engine design, good application effects can be obtained.

[0034] Additional aspects and advantages of the present invention will be given in part in the following description, become apparent in part from the following description, or be learned through the practice of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS

[0035] The above and / or additional aspects and advantages of the present invention will become apparent and be readily understood from the description of the embodiments in conjunction with the following drawings, in which:

[0036] Figure 1 The curved solid surface boundary in the rectangular grid computational domain is shown;

[0037] Figure 2 The basic structure of the curved solid surface boundary handling format for the lattice Boltzmann method is described;

[0038] Figure 3 The streamline diagram of the steady state (Re = 40) of the flow around a cylinder is described;

[0039] Figure 4 The streamline diagram of the unsteady periodic state (Re = 3900) of the flow around a cylinder is described;

[0040] Figure 5 The numerical simulation results of the NACA0012 airfoil (Re = 3000, ) are described; (1) The pressure coefficient on the upper surface, (2) The streamline diagram and the pressure contour map;

[0041] Figure 6 The numerical simulation results of the NACA0012 airfoil (Re = 100000, ) are described; (1) The lift coefficient, (2) The drag coefficient, (3) The streamline diagram, (4) The pressure contour map;

[0042] Figure 7 The numerical simulation results of the SD7003 airfoil (Re = 60000, ) and the U-shaped cavity (Re = 1000) are described; (1) The SD7003 airfoil, Re = 60000, The streamline diagram, (2) The streamline diagram of the flow in the U-shaped square cavity, Re = 1000. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0043] In view of the problem of handling the curved surface boundary in the numerical simulation of complex shapes based on the Lattice Boltzmann method (LBM) algorithm, the present invention proposes a new curved surface boundary handling method, which can not only handle the external flow fields such as airfoil flow around in practical engineering calculations, but also be applicable to the numerical simulation of the internal flow field in a cavity with a curved physical boundary.

[0044] The following takes the numerical simulation of cylinder flow around, the numerical simulation of the turbulent flow field of NACA0012 airfoil, the numerical simulation of the turbulent flow field of SD7003 airfoil, and the numerical simulation of the internal flow field in a U-shaped cavity as examples to illustrate the new curved surface boundary handling method.

[0045] Example 1, numerical simulation of cylinder flow around:

[0046] When performing the numerical simulation of cylinder flow around, after performing conventional CFD operations such as establishing a three-dimensional model of the cylinder and mesh generation, the distribution function of the points on the cylinder surface boundary is solved through the following process:

[0047] As Figure 2 shown, taking the distribution function of point f along the direction as an example for solution, points fff, ff, and f are fluid points. After each evolution, the information of points fff and ff has been updated. However, the distribution function of point f along the direction is unknown because it migrates from the virtual point b, and point b is a virtual point without a distribution function. Therefore, we can only construct the distribution function of point f along the direction through the information of points fff, ff, and w.

[0048] First, the distribution function of point f along the direction is divided into the equilibrium distribution function and the non-equilibrium distribution function for processing.

[0049] (1)

[0050] Among them, the equilibrium distribution function can be obtained through the calculation formula of the equilibrium distribution function . Therefore, we only need to know the density and velocity of point f. In this method, the density and velocity of point f can be interpolated by the points around point f

[0051] (points marked with black triangles)

[0052] In this way, we can obtain the equilibrium distribution function , Next, we will solve the non - equilibrium distribution function of point f .

[0053] First, we use linear interpolation to synthesize the information of points ff and w

[0054] (2)

[0055] where is the embedding depth, and then we use linear extrapolation to synthesize the information of points fff and ff

[0056] (3)

[0057] Then, we average the results of linear interpolation and extrapolation to obtain

[0058] (4)

[0059] After each evolution, in formula (4), except for this term, everything else is known. Therefore, next we need to solve .

[0060] First, we know that the distribution function of point w can be divided into two parts: non - equilibrium state and equilibrium state

[0061] (5)

[0062] or

[0063] (6)

[0064] At the same time, we know that due to the no - slip condition restriction on the solid surface

[0065] (7)

[0066] Therefore, formula (6) becomes

[0067] (8)

[0068] Since the equilibrium distribution function of point w can be obtained through the equilibrium distribution function calculation formula, we only need to use the information of points ff and f to interpolate the velocity and density of point w. If we consider a stationary solid surface, the surface velocity is 0. If the surface is not stationary, we can determine the velocity of point w according to the specific situation. Then, we use linear interpolation and extrapolation again to solve the distribution function of point w

[0069] (9)

[0070] (10)

[0071] Combining formula (9) and (10)

[0072] (11)

[0073] where the distribution function of point b along the direction migrates from the distribution function of point f along the direction after collision migration. Similarly, the distribution function of point f along the direction migrates from the distribution function of point ff along the direction after collision migration. The distribution function of point ff along the direction migrates from the distribution function of point fff along the direction after collision migration. Substitute formula (11) into formula (8)

[0074] (12)

[0075] Substitute formula (12) into formula (4) and we get

[0076] (13)

[0077] Up to this point, the distribution function of point f along the direction has been solved. The non-equilibrium distribution function of this point is constructed by comprehensively using the information of points fff, ff, f and w, and the equilibrium distribution function is constructed by using the comprehensive information of the red-circle punctuation.

[0078] Compared with the Chen format, the present invention treats the non-equilibrium and equilibrium distribution functions differently when constructing the distribution function. And the new format simultaneously uses linear interpolation and extrapolation to comprehensively consider the information and influence of surrounding points. Compared with the Mei format and the Guo format, more reasonable surrounding points are selected as the interpolation information sources of the target point velocity and density when constructing the equilibrium state, and the velocities and densities of virtual points (inconsistent with physical reality) are not forcibly constructed. The distribution function of point b along the direction is just a concept for data exchange (assignment). At the same time, when calculating the non-equilibrium distribution function, the Mei format and the Guo format only use the information of point f as the interpolation information source, while the present format selects internal and external linear interpolation, which is more reasonable and simultaneously considers the comprehensive information of multiple surrounding points.

[0079]

[0080]

[0081] Table 1. Comparison of the numerical calculation results of the flow around a cylinder using the present invention with other data (Re = 40)

[0082]

[0083] Table 1 presents the comparison of the flow data of the steady flow around a cylinder with other results when the new curved surface boundary treatment format proposed by the present invention is used at Re = 40. As can be seen from the table, the calculation accuracy of the new curved surface boundary treatment format is better than that of the other three popular curved surface boundary treatment formats. Since the new format is different from the other formats and is divided into two parts: pre - treatment and post - construction, the time required for each iteration step is longer than that of the other formats. However, due to the high accuracy of the new format itself, it is closer to the theoretical solution after each iteration step, so it can accelerate convergence, and thus the overall time consumed is basically the same as that of the other three formats.

[0084] Figure 3 The streamline diagram of the steady flow around a cylinder at Re = 40 is shown, and two stable attached vortices can be clearly observed.

[0085] Table 2. Comparison of the Hopf flow bifurcation points of the flow around a cylinder using the new format with other results

[0086]

[0087] Table 2 presents the comparison of the predictions of the Hopf flow bifurcation points of the flow around a cylinder by different literatures and three popular formats with the calculation results of the new format. It can be seen that the calculation results of the new format are in good agreement with those of other literatures, indicating that the calculation accuracy of the new format is relatively high.

[0088] Table 3. Comparison of the calculation results of the flow around a cylinder using the new format with other literatures (Re = 3900).

[0089]

[0090] Figure 4 And Table 3 shows the calculation results of the unsteady periodic flow state of the flow around a cylinder using the new format. Figure 4 (a) Describes the curve of the lift coefficient varying with time, , , , , represents a complete cycle, corresponding to each time node Figure 4 (b), (c), (d), (e), (f) successively show the instantaneous streamline diagrams of the flow around a cylinder at Re = 3900.

[0091] Example 2, numerical simulation of the turbulent flow field of the NACA0012 airfoil, including two states: Re = 3000, angle of attack = 0° and Re = 100000, angle of attack = 0°;

[0092] Example 3, numerical simulation of the turbulent flow field of the SD7003 airfoil, with the state of Re = 60000, angle of attack = 14°;

[0093] Example 4, numerical simulation of the internal flow field of the U-shaped cavity, with the state of Re = 1000;

[0094] In the above Examples 2, 3, and 4, the distribution function is solved for the points on the physical surface boundaries of the NACA0012 airfoil, SD7003 airfoil, and U-shaped cavity using the same method as in Example 1, and the numerical simulation results obtained are as Figures 5 to 7 shown.

[0095] Figure 5 The numerical simulation of the turbulent flow field of the NACA0012 airfoil is presented, showing the upper surface pressure coefficient, streamline diagram, and pressure contour map at Re = 3000, angle of attack = 0°. Figure 6 The numerical simulation of the turbulent flow field of the NACA0012 airfoil is described, showing the lift and drag coefficients, streamline diagram, and pressure contour map at Re = 100000, angle of attack = 0°. Figure 7 The numerical simulation of the turbulent flow field of the SD7003 airfoil is presented, showing the streamline diagram at Re = 60000, angle of attack = 14°, and the streamline diagram of the internal flow field of the U-shaped cavity (Re = 1000) is also given. It can be seen that the curved surface boundary format proposed in the present invention is applicable to the numerical simulation of the external flow field around the airfoil and the internal flow field of the U-shaped cavity.

[0096] Although the embodiments of the present invention have been shown and described above, it can be understood that the above embodiments are exemplary and should not be construed as limiting the present invention. Those of ordinary skill in the art can make changes, modifications, substitutions, and variations to the above embodiments within the scope of the present invention without departing from the principles and spirit of the present invention.

Claims

1. A numerical simulation method for an airfoil spoiler flow field, the numerical simulation including establishing an airfoil model, dividing grids, and processing curved surface boundaries, characterized in that: During the numerical simulation process, the following treatment method is used for the airfoil surface boundary: For the point f on the object surface boundary, the distribution function of point f along the α direction is solved through the following process Step 1: The distribution function of point f along the α direction is divided into an equilibrium distribution function and a non-equilibrium distribution function Among them, the equilibrium distribution function is obtained through the calculation formula of the equilibrium distribution function Function(ρ,u); Step 2: The following steps are used to solve the non-equilibrium distribution function of point f Step 2.1: Use linear interpolation to synthesize the information of points ff and w wherein is the embedding depth; x f is the abscissa of point f, x ff is the abscissa of point ff, x w is the abscissa of point w, x b is the abscissa of point b; Step 2.2: Use linear extrapolation to synthesize the information of points fff and ff Step 2.3: Average the results of linear interpolation and extrapolation to obtain: Among them: The information of w is calculated by the following process The distribution function of point w is divided into two parts: non-equilibrium state and equilibrium state Determined by the no-slip condition on the object surface Obtained Among them, the equilibrium distribution function of point w Obtained through the calculation formula of the equilibrium distribution function Function(ρ,u); At the same time, use linear interpolation and extrapolation to solve the distribution function of point w Thus, wherein the distribution function of point b along the -α direction is migrated from the distribution function of point f along the -α direction after collision migration, the distribution function of point f along the -α direction is migrated from the distribution function of point ff along the -α direction after collision migration, and the distribution function of point ff along the -α direction is migrated from the distribution function of point fff along the -α direction after collision migration; Thus, obtain: Furthermore, obtain 2. The method according to claim 1, wherein: In step 1, when calculating the equilibrium distribution function the density and velocity of point f are interpolated from the points around point f:

3. The method according to claim 1, wherein: In step 3, the information of points ff and f is used to interpolate the velocity and density of point w, and then the equilibrium distribution function of point w is calculated through the equilibrium distribution function 4. A numerical simulation method for the flow field in an engine cavity, wherein the numerical simulation includes establishing a cavity model, meshing, and curved surface boundary treatment, and is characterized in that: During the numerical simulation process, the following treatment method is used for the cavity surface boundary: For the point f on the object surface boundary, the distribution function of point f along the α direction is solved through the following process Step 1: The distribution function of point f along the α direction is divided into an equilibrium state distribution function and a non-equilibrium state distribution function Among them, the equilibrium distribution function is obtained by calculating the formula of the equilibrium distribution function Function(ρ,u); Step 2: The following steps are used to solve the non-equilibrium distribution function of point f Step 2.1: Use linear interpolation to synthesize the information of points ff and w Among them is the embedding depth; x f is the abscissa of point f, x ff is the abscissa of point ff, x w is the abscissa of point w, x b is the abscissa of point b; Step 2.2: Use linear extrapolation to synthesize the information of points fff and ff Step 2.3: Average the results of linear interpolation and extrapolation to obtain: Among them: The information of w is calculated by the following process The distribution function of point w is divided into two parts: non-equilibrium state and equilibrium state Determined by the no-slip condition on the object surface Obtained Among them, the equilibrium distribution function of point w Obtained through the calculation formula of the equilibrium distribution function Function(ρ,u); At the same time, use linear interpolation and extrapolation to solve the distribution function of point w Thus, Among them, the distribution function of point b along the -α direction is migrated from the distribution function of point f along the -α direction after collision migration. The distribution function of point f along the -α direction is migrated from the distribution function of point ff along the -α direction after collision migration. The distribution function of point ff along the -α direction is migrated from the distribution function of point fff along the -α direction after collision migration; Thus, obtain: Furthermore, obtain 5. The method according to claim 4, characterized in that: In step 1, the equilibrium distribution function is calculated When, the density and velocity of point f are interpolated using the points around point f:

6. The method according to claim 4, wherein: In step 3, the information of points ff and f is used to interpolate the velocity and density of point w, and then the equilibrium distribution function of point w is calculated through the equilibrium distribution function.

Images