Compact high-order reconstruction method based on multiple discrete moments on unstructured grid

Abstract

The invention relates to a compact high-order reconstruction method based on multiple discrete moments on an unstructured grid. The compact high-order reconstruction method comprises the steps of obtaining a grid topological relation; calculating and storing the volume integral point of each grid unit and the area division point of each surface; constructing a linear polynomial based on nodes on the boundary surface and the centroids of adjacent grid units, and calculating and storing a least square coefficient; solving and storing the body center value of each grid unit; solving a linear least square overdetermined equation of each surface to obtain an average value and a first derivative value of calculation variables of each boundary surface; obtaining first-order and second-order derivative values of the calculation variables according to the Gaussian divergence theorem; calculating a variable value of each boundary surface according to the Taylor series; mapping each grid unit to a standard grid on a local coordinate system to realize quadratic polynomial reconstruction of point values; calculating flux; updating the volume integral average value and the point value; and completing time stepping according to a third-order Runge-Kutta method. According to the method, the third-order precision can be achieved based on the compact template on the unstructured grid.

Description

Technical Field

[0001] The present invention relates to the technical field of computational fluid dynamics, and particularly to a compact high-order reconstruction method based on multi-discrete moments on unstructured grids. Background Art

[0002] In the field of computational fluid dynamics (CFD), high-order numerical methods (referring to numerical algorithms of order three or higher, such as high-order finite volume methods, discontinuous Galerkin methods, etc.) have always been the focus of research due to their significant advantages in numerical accuracy and computational efficiency. Although second-order methods are widely used in commercial software and engineering applications, their inherent numerical dissipation problems make it difficult to accurately simulate flow problems involving complex physical phenomena such as turbulence, eddies, and vortices. In contrast, high-order methods can provide higher numerical accuracy and efficiency because their truncation errors decrease faster than second-order methods.

[0003] However, the implementation of high-order methods on unstructured grids is not easy. The main difficulty lies in that the irregularity of the grid topology increases the complexity of constructing templates and calculating gradients or high-order derivatives. Different from structured grids, the node connection topological relationships of unstructured grids are uneven and difficult to predict, and complex data structures and interpolation techniques are required to handle the spatial distribution of adjacent cells, which increases the computational cost and reduces the computational efficiency. In addition, due to the fact that grid cells with poor quality (judged by orthogonality, aspect ratio, skewness, etc.) will reduce the accuracy and lead to numerical instability, the consistency and stability of high-order methods on unstructured grids have high requirements for grid quality.

[0004] Currently, common high-order numerical methods include: high-order finite volume (FV) methods, discontinuous Galerkin (DG) methods, spectral volume (SV) methods, spectral difference (SD) methods, etc.

[0005] Due to its inherent conservation and adaptability to unstructured grids, the finite volume method has become the main method for most production-level computational codes. High-order accuracy in the finite volume method is achieved through polynomial approximation to perform high-order reconstruction of physical variables within each cell. However, since only one value (i.e., the volume integral average, VIA) is defined for each cell in the finite volume method, high-order reconstruction in the finite volume method needs to be achieved by expanding the template. The problems brought by large templates are mainly as follows:

[0006] (1) Cache misses: Since it is necessary to access the solution data stored in distant memory locations, the cache miss rate increases;

[0007] (2) Reduced parallel computing efficiency: A large amount of data exchange is required between parallel partitions, reducing the efficiency of parallel computing;

[0008] (3) Large memory occupation: The high-order polynomial reconstruction process requires a large amount of memory, especially for the storage of the coefficient matrix.

[0009] (4) Matrix condition number sensitivity: In the least squares formula, the condition number of the matrix depends on the template size and is sensitive to the grid quality, which may lead to matrix ill-conditioning or even singularity.

[0010] In the framework of the traditional finite volume method, the most widely used high-order method is the k-exact least squares method. For a desired convergence order, this method assumes that the variable values within each cell can be represented by a polynomial of the corresponding order. Usually, a k-th degree polynomial can achieve (k + 1)-order spatial reconstruction. This method constructs a system of equations based on a template with sufficient degrees of freedom (DOF), and then uses the least squares method to solve the overdetermined system of equations to determine the coefficients of the reconstructed polynomial. Through this method, third-order, fourth-order, fifth-order, and even higher-order accuracies can be achieved theoretically. However, the wide application of the k-exact least squares method is restricted by two key issues: template selection and the solution of the overdetermined system of equations. The template must contain more cells than the number of coefficients to be solved, which introduces problems related to the selection of the surrounding grid, the storage of topological relationships, boundary conditions, and parallel implementation strategies. In addition, solving the overdetermined system of equations also raises problems in terms of computational efficiency and numerical robustness. Although the later proposed compact least squares finite volume method (CLSFV) constructs the reconstructed polynomial and its spatial derivatives by using face-adjacent cells, this method generates a large sparse linear system. This system must be solved implicitly using iterative methods or in combination with implicit time integration methods, which is computationally too expensive for unsteady flow problems.

[0011] As an alternative to the traditional finite volume method, in recent years, high-order methods that utilize additional degrees of freedom within compact templates have received extensive attention. Among them, the most well-known method is the discontinuous Galerkin (DG) method. The DG method adopts a discontinuous approximation within each cell and weakly enforces conservation, and its numerical trace plays a key role in determining the consistency, stability, and accuracy of the scheme. In addition to the DG method, relatively more research has been conducted on the spectral volume (SV) method and the spectral difference (SD) method. The SV method introduces the concept of "spectral volume", divides each spectral volume into control volumes, and uses the volume integral average (VIAs) of these sub-cells. The SD method places solution points and flux points within each cell, usually aligned with the quadrature points. Although the basic principles of the DG, SV, and SD schemes are similar, they differ in how to select and update the degrees of freedom. Huynh further generalizes these methods through the flux reconstruction (FR) framework, from which the DG, SV, and SD methods can be derived. These methods have the following problems in practical flow problems with complex structures:

[0012] (1) Low computational efficiency: The computational cost per unit is relatively high, leading to an increase in the overall computational cost and storage requirements;

[0013] (2) Strict stability conditions: The CFL condition for explicit time integration is strict, which limits the time step;

[0014] (3) Lack of robustness: There is a lack of reliable limiting techniques near discontinuities to suppress spurious oscillations.

[0015] To achieve an optimal balance among numerical accuracy, computational cost, and numerical robustness, a series of new numerical schemes applicable to unstructured grids based on the Multi-moment Finite Volume Method have been developed. Different from traditional compact stencil high-order methods, this method introduces the Point Values (PVs) at cell vertices as additional degrees of freedom on the basis of the traditional finite volume method. Both the Volume Integral Averages (VIAs) and the Point Values (PVs) are regarded as computational variables and updated synchronously in time. The VIA, as a conserved quantity, is updated by applying the finite volume method in the integral form of the governing equations; while the PVs are calculated by the point value numerical scheme derived from the differential form. This method is called the VPM method (Volume-integrated average and Point value based Multi-moment), which improves the numerical efficiency and robustness while ensuring strict conservation by sharing the PVs of adjacent cells.

[0016] The original VPM method was proposed in 2014, which established a high-order multi-moment format in the local coordinate system. This format uses the VIAs and PVs of the target element and its directly adjacent elements as degrees of freedom. Although this format achieves third-order accuracy on most mesh types, on tetrahedral meshes, due to insufficient degrees of freedom provided by the compact face adjacent template, quadratic polynomial reconstruction cannot be achieved. To solve this problem, the VPM-Ex method (VPM method based on extended template reconstruction) was proposed. This method directly constructs a quadratic polynomial by using an extended template (including adjacent elements sharing vertices) in the global coordinate system, and updates the PV moments using the derivative of this polynomial, thereby improving the numerical accuracy on mixed unstructured meshes. To extend the convergence order beyond third order, the VPM-CLS method (VPM method based on constrained least squares) was further proposed. This method can theoretically achieve arbitrary high-order accuracy by making full use of the degrees of freedom on the compact template. The VPM-CLS method constructs a reconstruction polynomial using the constrained least squares method in the global coordinate system, eliminating the need for mesh transformation in the multi-moment format. On the compact template containing all directly adjacent elements, the VPM-CLS scheme can achieve quartic polynomial reconstruction, thus reaching fifth-order accuracy. This significantly improves the accuracy in terms of numerical error and convergence behavior, while not increasing the algorithm complexity.

[0017] Although the existing VPM methods have made significant progress, there is still room for improvement in terms of numerical accuracy, algorithm simplicity, and computational efficiency. Further improvements are needed to optimize the balance of various performance aspects to make it more suitable for complex flows. Summary of the Invention

[0018] The object of the present invention is to provide a compact high-order reconstruction method based on multi-discrete moments on unstructured meshes, which can achieve third-order accuracy based on a compact template on unstructured meshes and solve compressible flow problems including shocks, discontinuities, etc.

[0019] To achieve the above object, the present invention provides a compact high-order reconstruction method based on multi-discrete moments on unstructured meshes, including the following steps:

[0020] Obtain the mesh topology relationship and define the initial volume integral average value and point value of each mesh element;

[0021] Calculate and store the volume integral points of each mesh element and the area integral points of each face;

[0022] Construct a linear polynomial based on the nodes on the boundary surface and the centroids of adjacent mesh elements, and calculate and store the least square coefficients on each mesh face;

[0023] Based on the volume integral average values and point values of each stored grid cell, solve and store the centroid values of each grid cell according to the conversion relationship between the centroid value, the volume integral average value, and the point value.

[0024] Solve the linear least squares overdetermined equations for each face to obtain the average value and the first derivative value of the computational variables on each boundary face.

[0025] According to the Gauss divergence theorem, obtain the first and second derivative values of the computational variables within each grid cell.

[0026] Construct the quadratic reconstruction polynomial for each grid cell according to the Taylor series, and calculate the variable values on each boundary face.

[0027] Map each grid cell to the standard grid in the local coordinate system to achieve the quadratic polynomial reconstruction of the point value.

[0028] Calculate the flux according to the upwind method or the Riemann solver.

[0029] Update the volume integral average value of each grid cell by the finite volume method according to the integral form of the governing equation, and update the point value of each grid cell according to the differential form of the governing equation.

[0030] Complete the time stepping according to the third-order Runge-Kutta method.

[0031] Optionally, the calculation and storage of the least squares coefficients on each grid face are performed before the time loop.

[0032] Optionally, when constructing the quadratic reconstruction polynomial for each grid cell according to the Taylor series and calculating the variable values on each boundary face, a limiter is applied to the flow problems that need to be restricted.

[0033] Optionally, the limiter is not activated in the smooth region.

[0034] Optionally, the upwind method is used to solve the flux of the convection equation, the Riemann solver is used to solve the flux of the Euler and Navier-Stokes equations, and the isentropic correction is adopted.

[0035] Optionally, the grid topology relationship includes the self-index of each grid cell, the face index it contains, the node index it contains, the adjacent grid cell index, and the adjacent node index.

[0036] Based on the same inventive concept, the present invention also provides a readable storage medium, on which a computer program is stored, and when the computer program is executed, it can implement the compact high-order reconstruction method based on multi-discrete moments on the unstructured grid as described above.

[0037] The compact high-order reconstruction method based on multi-discrete moments on unstructured grids proposed by the present invention is different from the existing finite volume high-order methods in that:

[0038] 1) The method for obtaining the second-order derivative required for quadratic reconstruction is different. Instead of using the quadratic least squares method based on a large template, after obtaining the variable values and first-order derivative values of the boundary faces through the linear least squares method on each boundary face, the first-order and second-order derivatives of the cell are directly obtained through the Gauss divergence theorem;

[0039] 2) Only face-adjacent cells are required as templates, and the templates are very compact;

[0040] 3) The first-order derivative of the boundary face can be used not only for reconstruction but also for solving the flux of the Laplace term.

[0041] Therefore, the compact high-order reconstruction method based on multi-discrete moments on unstructured grids proposed by the present invention has at least the following beneficial effects compared with the prior art:

[0042] 1) Compact templates, only including face-adjacent grids, avoiding a series of problems brought by large templates;

[0043] 2) It is not necessary to directly calculate the second-order derivative in the global coordinates through quadratic least squares, greatly improving the calculation efficiency and robustness;

[0044] 3) It is not necessary to calculate and store the coefficient matrix for quadratic polynomial calculation, saving the calculation cost;

[0045] 4) Calculating the flux of the Laplace term through the first-order derivative of the boundary face further simplifies the calculation;

[0046] 5) Good stability and robustness;

[0047] 6) Combined with a compressible limiter, it can well capture structures such as shock waves and discontinuities. Description of the Drawings

[0048] Those of ordinary skill in the art should understand that the provided drawings are used to better understand the present invention and do not constitute any limitation to the scope of the present invention. Among them:

[0049] Figure 1 is the flowchart of the compact high-order reconstruction method based on multi-discrete moments on unstructured grids provided by an embodiment of the present invention;

[0050] Figure 2 is the schematic diagram of the grid type provided by an embodiment of the present invention;

[0051] Figure 3 is the schematic diagram of the template for face linear reconstruction of a two-dimensional triangular grid and a three-dimensional tetrahedral grid provided by an embodiment of the present invention;

[0052] Figure 4 A comparison diagram of the results of the VPM-FR method provided by an embodiment of the present invention and the existing finite volume high-order method on a quadrilateral;

[0053] Figure 5 A comparison diagram of the results of the VPM-FR method provided by an embodiment of the present invention and the existing finite volume high-order method on a triangle;

[0054] Figure 6 A comparison diagram of the results of the VPM-FR method provided by an embodiment of the present invention and the existing finite volume high-order method on a distorted grid;

[0055] Figure 7 The double Mach numerical calculation results with a grid size of 1 / 240 provided in one embodiment of the present invention;

[0056] Figure 8 The double Mach numerical calculation results with a grid size of 1 / 480 provided in one embodiment of the present invention;

[0057] Figure 9 An instantaneous pseudo-schlieren view of the calculation result of the VPM-FR method provided by an embodiment of the present invention when calculating a three-dimensional compressible practical problem;

[0058] Figure 10 A Q-standard diagram of the calculation results of a three-dimensional compressible practical problem calculated using the VPM-FR method provided by an embodiment of the present invention; DETAILED DESCRIPTION

[0059] In order to make the purpose, advantages and features of the present invention clearer, the present invention is further described in detail below in conjunction with the accompanying drawings and specific embodiments. It should be noted that the accompanying drawings are in a very simplified form and use non-precise proportions, which are only used to conveniently and clearly assist in explaining the purpose of the embodiments of the present invention. In order to make the purpose, features and advantages of the present invention more obvious and easy to understand, please refer to the accompanying drawings. It should be noted that the structure, proportion, size, etc. illustrated in the drawings of this specification are only used to match the content disclosed in the specification for people familiar with this technology to understand and read, and are not used to limit the limiting conditions for the implementation of the present invention. Any modification of the structure, change in the proportional relationship or adjustment of the size, under the same or similar conditions as the effects that can be produced by the present invention and the purposes that can be achieved, should still fall within the scope of the technical content disclosed by the present invention.

[0060] As used in this disclosure, the singular forms "a", "an", and "the" include plural referents unless the context clearly dictates otherwise. As used in this disclosure, the term "or" is generally used in a sense including "and / or" unless the context clearly dictates otherwise. As used in this disclosure,

[0061] Please refer to Figure 1 , this embodiment provides a compact high-order reconstruction method based on multi-discrete moments on unstructured grids, including the following steps:

[0062] S1. Obtain the grid topological relationship, and define the initial volume integral average value and point value of each grid cell;

[0063] S2. Calculate and store the volume integral points of each grid cell and the area integral points of each face;

[0064] S3. Construct a linear polynomial based on the nodes on the boundary surface and the centroid of the adjacent grid cells, and calculate and store the least square coefficients on each grid face;

[0065] S4. Based on the stored initial volume integral average value and point value of each grid cell, solve and store the centroid value of each grid cell according to the conversion relationship between the centroid value and the volume integral average value and the point value;

[0066] S5. Solve the linear least square overdetermined equations of each face to obtain the average value and the first derivative value of the calculated variables on each boundary surface;

[0067] S6. According to the Gauss divergence theorem, obtain the first and second derivative values of the calculated variables within each grid cell;

[0068] S7. According to the Taylor series, construct a quadratic reconstruction polynomial for each grid cell, and calculate the variable values on each boundary surface;

[0069] S8. Map each grid cell to the standard grid in the local coordinate system to realize the quadratic polynomial reconstruction of the point value;

[0070] S9. Calculate the flux according to the upwind method or the Riemann solver;

[0071] S10. According to the integral form of the governing equation, update the volume integral average value of each grid cell by the finite volume method, and update the point value of each grid cell according to the differential form of the governing equation;

[0072] S11. Complete the time stepping according to the third-order Runge-Kutta method.

[0073] Specifically, first execute S1 to obtain the grid topological relationship, and define the initial volume integral average value (VIA) and point value (PV s ). Among them, the grid topological relationship includes the self-index of each grid cell, the included face index, the included node index, the adjacent grid cell index, and the adjacent node index, etc.

[0074] The present invention is applicable to a variety of mesh types, including two-dimensional triangles, quadrilaterals, and three-dimensional tetrahedrons, hexahedrons, triangular prisms, and pyramidal meshes, as Figure 2 shown.

[0075] Figure 2 In, Ω i (i = 1,..., I) represents the control volume, where I represents the total number of meshes. The control volume Ω i has a boundary surface represented by Γ ij (j = 1,..., J), and the points of the control volume Ω i are represented by θ ik (k = 1,..., K). The centroid of the control volume Ω i is represented by θ ic , and the centroid of the surface Γ ij is represented by θ ij .

[0076] Using φ to represent the calculation variable, the definitions of the initial volume integral average value and the point value are as follows:

[0077]

[0078] φ ik (t) ≡ φ(x ik , y ik , z ik , t), k = 1, 2,..., K.

[0079] Then proceed to execute S2 to calculate and store the volume integral points of each mesh cell and the area integral points of each surface. This is because in the principle of the finite volume method, the volume integral average value needs to be calculated through volume averaging; and to obtain the reconstruction expression, the surface average value needs to be calculated through surface averaging.

[0080] Next, execute S3 to construct a linear polynomial based on the nodes on the boundary surface and the centroids of adjacent mesh cells, and calculate and store the least squares coefficients on each mesh surface. Figure 3 Taking two-dimensional triangular meshes and three-dimensional tetrahedral meshes as examples, the templates for surface linear reconstruction are shown.

[0081] The linear interpolation polynomial is:

[0082]

[0083] Among them, φ represents the calculation variable, and φ ij represents the value of the variable at the j-th surface of cell i, and (x ij , y ij , z ij ) is the centroid coordinate of the surface.

[0084] The system of equations formed by substituting into the template can be written in matrix form:

[0085] M·A = B

[0086] Where M is the linear least squares matrix for each face, A is the coefficient vector, and B is the vector composed of the calculated variables at each point;

[0087] It should be noted that this step is carried out before the time loop, that is, the linear least squares matrix M for each face only needs to be calculated once and saved for subsequent use, while the vector B composed of the calculated variables and the coefficient vector A need to be defined or solved within the time loop (S5).

[0088] Next, execute S4. Based on the volume integral average values and point values of each grid cell that have been stored, according to the conversion relationship between the centroid value and the volume integral average value and the point value, solve and store the centroid values of each grid cell. The conversion relationships for different grid cell types provided in this embodiment are shown in Table 1 below:

[0089]

[0090] Table 1

[0091] Where, φ ic represents the value of the variable at the centroid of cell i, represents the volume integral average value of cell i, and φ in represents the point values of cell i.

[0092] Next, execute S5 to solve the linear least squares overdetermined equations for each face, and obtain the average values and first derivative values of the calculated variables on each boundary face.

[0093] In this embodiment, the overdetermined system of equations obtained after substituting each value that makes up the template is solved by the following formula:

[0094] A = (M T M) -1 M T B

[0095] Then execute S6. According to the Gauss divergence theorem, obtain the first and second derivative values of the calculated variables within each grid cell. The specific form is:

[0096]

[0097] Where, |Ω i | represents the volume of cell i, |Γ ij | represents the area of the jth face of cell i, and n ij represents the outer normal vector of the face.

[0098] Next, execute S7. According to the Taylor series, construct the quadratic reconstruction polynomials for each grid cell, calculate the variable values at each interface, and apply a limiter to the flow problems that need to be restricted.

[0099] The quadratic reconstruction polynomial has the form:

[0100] Φ i (φ|; x) = c 000 + c 100 (x - x ic ) + c 010 (y - y ic ) + c 001 (z - z ic ) + c 200 (x - x ic ) 2 + c 020 (y - y ic ) 2 + c 002 (z - Z ic ) 2 + c 110 (x - x ic )(y - y ic ) + c 011 (y - y ic (z - z ic ) + c 101 (x - x ic (z - z ic ),

[0101] The high-order limiter adopted in this embodiment combines the MOGV limiter and the limiting function proposed by Nishikawa. A limiting function is defined for each cell, and the calculation method is:

[0102]

[0103] Among them, represents the slope limiter at the α-th Gauss integration point on the surface Γ ij , and the calculation method is:

[0104]

[0105] Among them, φ max and φ min represent the maximum and minimum calculated variable values in the reconstruction template of the current grid, and φ ij,α represents the calculated variable value at the surface Gauss point using the unconstrained calculation.

[0106] In order to maintain the high-order convergence rate, it is necessary to not activate the limiter in the smooth region and further modify the limiter:

[0107]

[0108] In this embodiment, while executing S7, S8 can also be executed to map each grid cell to a standard grid in the local coordinate system, realizing the quadratic polynomial reconstruction of the point values.

[0109] In this embodiment, the reconstruction function within any type of grid can be expressed as:

[0110]

[0111] where represents the local coordinates, and represent the first-order and second-order derivatives at the cell center, {ψ} represents the basis function, and {A} and {B} are parameters adjusted according to different grid types.

[0112] Taking the triangular grid as an example, A1 = 1, B1 = 0, A2 = B2 = 0, then the reconstruction expression is:

[0113]

[0114] where the basis function is:

[0115]

[0116] Then execute S9 to calculate the flux according to the upwind method or the Riemann solver.

[0117] The governing equation for the convection problem is:

[0118]

[0119] where φ represents the calculated variable, u represents the velocity vector, and the flux uφ at each face is calculated.

[0120] The governing equation for the compressible problem is the Navier - Stokes (NS) equation:

[0121]

[0122] where U represents the vector composed of calculated variables, and in this invention, conserved variables are adopted. F represents the inviscid flux, G represents the viscous flux, and the specific expressions are:

[0123] U = [ρ, ρu, ρv, ρw, E] T ,

[0124]

[0125]

[0126]

[0127]

[0128]

[0129]

[0130] Among them, ρ is the density, (u, v, w) are the velocity components in the (x, y, z) directions, E is the total energy, p is the pressure, T is the temperature, κ is the thermal conductivity, and the calculation method is as follows:

[0131]

[0132] The calculation method of the viscous stress tensor is as follows:

[0133]

[0134] Among them, μ is the dynamic viscosity coefficient, γ is the specific heat ratio, Pr is the Prandtl number, R is the gas constant, and δ αβ is the Kronecker function. Taking γ = 1.4, the total energy is:

[0135]

[0136] In the actual application process of the present invention, the upwind method is adopted when solving the flux of the convection equation, the Roe Riemann solver is adopted when solving the flux of the Euler and NS equations, and the isentropic correction is adopted.

[0137] The formula for calculating the inviscid flux of the volume integral average value is:

[0138]

[0139] Among them, L and R are the matrices composed of the left and right eigenvectors, Λ is the diagonal matrix composed of the eigenvalues, and A = RΛL is the Jacobian matrix. are the left and right values of the surface, that is, the values calculated by the reconstruction of the elements on both sides of the surface.

[0140] The formula for calculating the viscous flux of the volume integral average value is:

[0141]

[0142] Among them, the calculation formulas for each quantity are:

[0143]

[0144] Among them, the values of the calculation variables and the first-order derivatives on the boundary surface are obtained by step S5.

[0145] The formula for calculating the inviscid flux of the point value is:

[0146]

[0147] Among them, A, B, and C are Jacobian matrices;

[0148]

[0149] The Jacobian matrix is calculated after Roe-averaging based on the average value of all volume integrals around the point value. The Roe-averaging process is as follows:

[0150]

[0151] Then the non-viscous flux calculation formula for the point value can be written as:

[0152]

[0153] The left and right values in the formula are obtained from the reconstructed values of the surrounding cells according to the following weighting formula:

[0154]

[0155] Among them, M represents the Jacobian matrix for converting from the local coordinate system to the global coordinate system, and ω ikl represents the weight of the l-th cell around the point value. The calculation method is:

[0156]

[0157] Among them, represents the vector pointing from the centroid of the control volume Ω i to the point θ ik and n τ± represents the unit vectors in the x, y, and z directions.

[0158] The viscous flux calculation formula for the point value is:

[0159]

[0160] Among them, the first-order and second-order derivatives are obtained by quadratic least squares based on the average value of the volume integral around the point value and other point values.

[0161] Then execute S10. According to the control equation in integral form, update the average value of the volume integral of each grid cell by the finite volume method, and update the point value of each grid cell according to the control equation in differential form.

[0162] The update method for the integral form adopted by the average value of the volume integral is:

[0163]

[0164] The update method for the differential form adopted by the point values is as follows:

[0165]

[0166] After the update, the following restrictions are imposed on the point values:

[0167]

[0168] Among them, represents the point values around the volume integral average value.

[0169] Finally, execute S11 to complete the time stepping according to the third-order Runge-Kutta method.

[0170] The third-order Runge-Kutta formula is:

[0171]

[0172] Among them, L(U) is the sum of the fluxes calculated in S10.

[0173] In this embodiment, S1 - S3 are before the time loop, S4 - S10 are within the time loop, and S4 - S8 are essentially a compact high-order reconstruction process based on multi-discrete moments on an unstructured grid, while S9 - S10 are the control equation solving processes.

[0174] The compact high-order reconstruction method based on multi-discrete moments on an unstructured grid proposed by the present invention is different from the existing finite volume high-order methods in that:

[0175] 1) The method for obtaining the second-order derivative required for quadratic reconstruction is different. Instead of using the quadratic least squares method based on a large template, first, through the linear least squares method on each boundary surface, after obtaining the variable values and first-order derivative values on the boundary surface, the first and second derivatives of the unit are directly obtained through the Gauss divergence theorem;

[0176] 2) Only face-adjacent cells are required as templates, and the templates are very compact;

[0177] 3) The first derivative on the boundary surface can be used not only for reconstruction but also for solving the flux of the Laplace term.

[0178] Therefore, the compact high-order reconstruction method based on multi-discrete moments on an unstructured grid proposed by the present invention has the following beneficial effects compared with the prior art:

[0179] 1) Compact templates, only including face-adjacent grids, avoiding a series of problems brought by large templates;

[0180] 2) It is not necessary to directly calculate the second-order derivative in the global coordinates through quadratic least squares, greatly improving the calculation efficiency and robustness;

[0181] 3) There is no need to calculate and store the coefficient matrix for quadratic polynomial calculations, saving computational costs;

[0182] 4) Calculate the flux of the Laplace term through the first-order derivative of the boundary surface, further simplifying the calculation;

[0183] 5) Good stability and robustness;

[0184] 6) Combined with a compressible limiter, it can well capture structures such as shock waves and discontinuities.

[0185] Next, some typical numerical example calculation results are presented to demonstrate the advantages of the compact high-order reconstruction method based on multi-discrete moments on unstructured grids (hereinafter simply referred to as the VPM-FR method) proposed in the present invention.

[0186] 1) Sin wave convection test: Demonstrate the numerical error, convergence order of the convection equation solved by the present method, and its characteristics on distorted grids.

[0187]

[0188] Table 2

[0189] Table 2 shows the numerical errors, convergence orders, and computational times of the VPM-FR method and the previous version of VPM (VPM-Ex method) on quadrilateral and triangular grids. The computational time is obtained by averaging after calculating 100 time steps on a single-core 13th Gen Intel(R) Core(TM) i7-13700H, 2.40 GHz CPU. It can be seen that both can achieve third-order accuracy, but the numerical error of the VPM-FR method is only 1 / 6 - 1 / 3 of that of the VPM-Ex method, and the computational time of the VPM-FR method is only 1 / 3 of that of the VPM-Ex method.

[0190] Figures 4 - 5 The results comparison between VPM-FR and existing high-order finite volume methods is shown. Among them, VPM-FR-MOGN represents the VPM-FR method combined with a compressible limiter, FV-MLP represents the finite volume method combined with the MLP limiter, and FV-ROUND represents the finite volume method combined with the ROUND limiter. It can be seen that the VPM-FR method has the smallest numerical error and the highest convergence order. Although it is affected after adding the limiter, its performance is still better than other methods.

[0191] Figure 6 The results comparison of various methods on distorted grids is shown. The abscissa represents the grid distortion degree. The smaller the abscissa, the higher the grid distortion degree. It can be seen that the VPM-FR method not only has the smallest numerical error, but also has the smallest absolute value of the slope as the abscissa decreases, indicating that it is less affected by the grid distortion degree.

[0192] 2) Double Mach test: Demonstrate the VPM-FR method provided by the present invention to capture the characteristics of shock waves and complex flows in Euler equation calculations.

[0193] Figures 7 - 8 The double Mach numerical calculation results with grid sizes of 1 / 240 and 1 / 480 are shown. The main advantages include the ability to resolve complex curled vortex structures, generate compact shock wave profiles without numerical oscillations, and exhibit extremely low dissipation when capturing shear-driven instabilities. The results show that the VPM-FR method performs well in balancing high-resolution discontinuity capture with fidelity of small-scale flow characteristics, and is particularly suitable for complex flow problems involving shock-vortex interactions.

[0194] 3) Transonic flow through a pit test: Demonstrate the accuracy of the VPM-FR method provided by the present invention in calculating three-dimensional compressible practical problems. Figures 9 - 10 The instantaneous pseudo-schlieren view and Q-standard diagram of the calculation results are shown respectively, and the three vortices S1, S2, and S3 can be clearly seen.

[0195] Based on the same inventive concept, an embodiment of the present invention further proposes a readable storage medium having a computer program stored thereon, which, when executed, can implement the compact high-order reconstruction method based on multiple discrete moments on the unstructured grid as described above.

[0196] A readable storage medium can be a tangible device that can hold and store instructions used by an instruction execution device. For example, it can be, but is not limited to, an electrical storage device, a magnetic storage device, an optical storage device, an electromagnetic storage device, a semiconductor storage device, or any suitable combination of the above. More specific examples (non-exhaustive list) of the readable storage medium include: portable computer disks, hard disks, random access memory (RAM), read-only memory (ROM), erasable programmable read-only memory (EPROM or flash memory), static random access memory (SRAM), portable compact disk read-only memory (CD-ROM), digital versatile disks (DVD), memory sticks, floppy disks, mechanically encoded devices such as punch cards or raised structures in grooves having instructions stored thereon, and any suitable combination of the above. The computer programs described herein can be downloaded from the readable storage medium to various computing / processing devices, or downloaded to an external computer or external storage device via a network, such as the Internet, a local area network, a wide area network, and / or a wireless network. The network can include copper transmission cables, optical fiber transmission, wireless transmission, routers, firewalls, switches, gateway computers, and / or edge servers. A network adapter card or network interface in each computing / processing device receives the computer program from the network and forwards the computer program for storage in the readable storage medium in each computing / processing device. The computer program for performing the operations of the present invention can be assembly instructions, instruction set architecture (ISA) instructions, machine instructions, machine-related instructions, microcode, firmware instructions, state setting data, or source code or object code written in any combination of one or more programming languages, including object-oriented programming languages - such as Smalltalk, C++, etc., and conventional procedural programming languages - such as the "C" language or similar programming languages. The computer program can be executed entirely on the user's computer, partially on the user's computer, executed as a stand-alone software package, partially on the user's computer and partially on a remote computer, or entirely on a remote computer or server. In the case of a remote computer, the remote computer can be connected to the user's computer through any type of network, including a local area network (LAN) or a wide area network (WAN), or can be connected to an external computer (e.g., via an Internet service provider through the Internet). In some embodiments, by using the state information of the computer program to personalize an electronic circuit, such as a programmable logic circuit, a field programmable gate array (FPGA), or a programmable logic array (PLA), the electronic circuit can execute computer-readable program instructions to implement various aspects of the present invention.

[0197] Aspects of the present invention are described herein with reference to the flowcharts and / or block diagrams of methods, systems, and computer program products according to embodiments of the present invention. It should be understood that each block of the flowcharts and / or block diagrams, and the combinations of blocks in the flowcharts and / or block diagrams, can be implemented by computer programs. These computer programs can be provided to a processor of a general-purpose computer, a special-purpose computer, or other programmable data processing device, thereby producing a machine such that when these programs are executed by the processor of the computer or other programmable data processing device, a device is produced that implements the functions / actions specified in one or more blocks of the flowchart and / or block diagram. These computer programs can also be stored in a readable storage medium, and these computer programs cause the computer, programmable data processing device, and / or other devices to work in a specific manner. Thus, the readable storage medium storing the computer program includes a manufactured article that includes instructions for implementing various aspects of the functions / actions specified in one or more blocks of the flowchart and / or block diagram.

[0198] The computer programs can also be loaded onto a computer, other programmable data processing device, or other device, such that a series of operation steps are executed on the computer, other programmable data processing device, or other device to produce a computer-implemented process, so that the computer programs executed on the computer, other programmable data processing device, or other device implement the functions / actions specified in one or more blocks of the flowchart and / or block diagram.

[0199] Since the readable storage medium provided by the present invention and the compact high-order reconstruction method based on multi-discrete moments on an unstructured grid described above belong to the same inventive concept, the readable storage medium provided by the present invention has all the advantages of the compact high-order reconstruction method based on multi-discrete moments on an unstructured grid described above. Therefore, the beneficial effects of the readable storage medium provided by the present invention will not be elaborated one by one here.

[0200] In summary, the embodiments of the present invention provide a compact high-order reconstruction method based on multi-discrete moments on an unstructured grid. On the basis of the existing VPM method, the reconstruction method is improved. Compared with the prior art, it has at least the following beneficial effects:

[0201] 1) Compact template, only including face-adjacent grids, avoiding a series of problems caused by large templates;

[0202] 2) It is not necessary to directly calculate the second-order derivative in the global coordinates through quadratic least squares, greatly improving the calculation efficiency and robustness;

[0203] 3) It is not necessary to calculate and store the coefficient matrix for quadratic polynomial calculation, saving the calculation cost;

[0204] 4) The flux of the Laplacian term is calculated through the first derivative of the boundary surface, further simplifying the calculation;

[0205] 5) Good stability and robustness;

[0206] 6) Combined with a compressible limiter, it can well capture structures such as shock waves and discontinuities.

[0207] In addition, it should also be recognized that although the present invention has been disclosed above in preferred embodiments, the above embodiments are not intended to limit the present invention. For any person skilled in the art, without departing from the scope of the technical solution of the present invention, many possible variations and modifications can be made to the technical solution of the present invention by using the technical content disclosed above, or modified into equivalent embodiments with equivalent changes. Therefore, any simple modification, equivalent change and modification made to the above embodiments according to the technical essence of the present invention without departing from the content of the technical solution of the present invention still fall within the scope of protection of the technical solution of the present invention.

Claims

1. A compact high-order reconstruction method based on multi-discrete moments on unstructured grids, characterized in that, It includes the following steps: Obtain the grid topological relationship, and define the initial volume integral average value and point value of each grid cell; Calculate and store the volume integral points of each grid cell and the area integral points of each face; Construct a linear polynomial based on the nodes on the boundary face and the centroid of the adjacent grid cell, and calculate and store the least-square coefficients on each grid face; Based on the stored volume integral average value and point value of each grid cell, solve and store the centroid value of each grid cell according to the conversion relationship between the centroid value and the volume integral average value and the point value; Solve the linear least-square overdetermined equation of each face to obtain the average value and the first-order derivative value of the calculated variable on each boundary face; According to the Gauss divergence theorem, obtain the first-order and second-order derivative values of the calculated variable in each grid cell; Construct a quadratic reconstruction polynomial of each grid cell according to the Taylor series, and calculate the variable value of each boundary face; Map each grid cell to the standard grid in the local coordinate system to realize the quadratic polynomial reconstruction of the point value; Calculate the flux according to the upwind method or the Riemann solver; According to the integral form of the governing equation, update the volume integral average value of each grid cell by the finite volume method, and update the point value of each grid cell according to the differential form of the governing equation; Complete the time stepping according to the third-order Runge-Kutta method.

2. The compact high-order reconstruction method based on multi-discrete moments on unstructured grids according to claim 1, characterized in that, The calculation and storage of the least-square coefficients on each grid face are performed before the time loop.

3. The compact high-order reconstruction method based on multi-discrete moments on unstructured grids according to claim 1, characterized in that, When constructing the quadratic reconstruction polynomial of each grid cell according to the Taylor series and calculating the variable value of each boundary face, a limiter is applied to the flow problem that needs to be restricted.

4. The compact high-order reconstruction method based on multi-discrete moments on unstructured grids according to claim 3, characterized in that, The limiter is not activated in the smooth region.

5. The compact high-order reconstruction method based on multi-discrete moments on unstructured grids according to claim 1, characterized in that, The upwind method is used to solve the flux of the convection equation, the Riemann solver is used to solve the flux of the Euler and Navier-Stokes equations, and isentropic correction is adopted.

6. The compact high-order reconstruction method based on multi-discrete moments on unstructured grids according to claim 1, characterized in that, The grid topological relationship includes the self-index of each grid cell, the included face index, the included node index, the adjacent grid cell index, and the adjacent node index.

7. A readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed, it can implement the compact high-order reconstruction method based on multi-discrete moments on the unstructured grid according to any one of claims 1-6.

Images