Wing boundary simulation method, terminal device and storage medium for optimizing isoparametric units

Abstract

The present invention discloses a wing boundary simulation method, terminal device, and storage medium for optimizing isoparametric units. After obtaining wing curved boundary information, the method determines the curved unit shape based on this information. Based on the curved unit shape, the Jacobian matrix of the wing physical unit and the isoparametric unit is constructed to obtain surface integrals and volume integrals. This method eliminates the need to define the distribution of parameter points, significantly improving the calculation accuracy of surface integrals and volume integrals.

Description

Technical Field

[0001] The present invention relates to aircraft wing simulation technology, in particular to a wing boundary simulation method for optimizing isoparametric units, terminal equipment and storage medium. Background Art

[0002] The airfoil, commonly known as the wing profile or blade cross-section, is one of the core factors affecting the overall performance of an aircraft. The airfoil is a fundamental element in the design of aerodynamic surfaces such as aircraft wings, tail planes, helicopter rotors, propellers, and wind turbine blades. It directly influences the aerodynamic performance of an aircraft, including key aerodynamic characteristics such as lift, drag, and pitching moment. The aerodynamic characteristics of the airfoil also directly affect the maneuverability and stability of an aircraft. A well-designed airfoil can provide sufficient control torque, ensuring stability and maneuverability under various flight conditions, thereby improving flight safety.

[0003] In addition, airfoils are closely related to energy conservation, emission reduction, and consumption reduction of aircraft. Optimizing airfoil design is one of the effective measures to reduce drag, energy consumption, and fuel consumption. Airfoil design determines the lift-to-drag ratio of an aircraft [Li Yunpeng, Han Yongzhi. Research Progress on Design of Lift-Enhancing Devices Based on Laminar Flow Wings [J]. Advances in Aeronautical Engineering, 2021, 12(4): 1-11.], that is, the ratio of lift to drag. The greater the lift-to-drag ratio, the higher the flight efficiency, which can provide a longer range under limited fuel conditions while reducing environmental pollution (Allison, E., Kroo, I., Sturdza, P., Suzuki, Y., & Martins-Rivas, H. (2010). Aircraft conceptual design with natural laminar flow. In Proceedings of the 27th Congress of the International Council of the Aeronautical Sciences (pp. 428–436). Nice, France: ICAS.). Excellent airfoil design can achieve a high lift-to-drag ratio, effectively reduce fuel consumption, lower carbon emissions, increase range and endurance, which is crucial to improving the flight efficiency and economy of aircraft.

[0004] In order to economically and efficiently obtain excellent airfoils with strong stability and good drag reduction effects, computational fluid dynamics methods are often used in engineering applications to achieve the design, simulation, and optimization of airfoils. Computational fluid dynamics is a numerical simulation method that discretizes the physical space into a computational space composed of grids, and then uses time and space numerical formats to simulate the complex flow problems of the wing during flight. The process of discretizing the physical space is called grid generation. Commonly used grids in engineering applications include structured grids, unstructured grids, and Cartesian grids. Structured grids have high generation efficiency in simple geometric shapes, but grid generation becomes increasingly difficult in complex situations. Cartesian grids are simple to generate and can provide high computational accuracy within a regular computational domain, but complex geometric shapes are prone to mismatch between the grid and the boundary, requiring special processing methods such as the immersed boundary method, which increases the complexity of implementation. Unstructured grids have strong automatic generation capabilities and can flexibly adapt to arbitrarily complex geometric shapes. They have been widely used in practical engineering problems in recent decades. Numerical schemes for unstructured grids include finite volume (FV), spectral volume (SV), discontinuous Galerkin (DG), PNPM, and correction procedure via reconstruction (CPR). Finite volume schemes are relatively simple to derive and implement, and can be easily extended to high-order accuracy. In recent years, the introduction of the compact finite volume (CFV) scheme and mature shock wave capture techniques have continuously increased the influence of finite volume schemes.

[0005] In order for unstructured grid finite volume numerical simulation methods to truly reflect the flight state of an airfoil, two aspects require special attention. First, the mesh edges of unstructured grids generated in engineering applications are all line segments, while airfoils generally have curved boundaries. Replacing curved boundaries with line segments results in a slight difference between the numerical shape of the wing and the actual shape, which may lead to inaccurate numerical simulations. Therefore, it is necessary to express the airfoil with a high-order curved boundary. Second, to improve computational accuracy and solution efficiency, high-precision numerical formats are often used for solutions in engineering applications. However, high-precision formats are often used for calculations within the calculation region. Numerical formats on the curved boundaries of airfoils often have difficulty obtaining high-order accuracy. The values ​​at the airfoil boundaries directly affect the calculation results of the entire calculation region. In other words, numerical simulations of airfoils with curved boundaries cannot achieve high-order accuracy, which will inevitably affect the overall accuracy.

[0006] Finite volume schemes discretize conservation law equations in integral form. Therefore, achieving high-order accuracy relies on matching the scheme's accuracy with high-order numerical integration, such as for the surface integral of fluxes and the volume integral of source terms. High-order numerical integration is relatively straightforward for straight-edge elements, and accurate results can be obtained by following the integration points and weights described in the literature (Ollivier-Gooch C, Nejat A, Michalak K. Obtaining and verifying high-order unstructured finite volume solutions to the Euler equations [J]. AIAA Journal. 2009, 47:2105–2120.). However, for elements with curved edges, particularly those with curved airfoils, the numerical integration process presents numerous challenges. Currently, two common approaches exist. The first is the isoparametric element method, commonly used in the finite element field. This method achieves equivalent computations on isoparametric elements by constructing a mapping between physical space and isoparametric space. However, this method requires defining the distribution of parameter points, solving and storing the Lagrange interpolation coefficients and Jacobian, so the process is somewhat cumbersome, especially in high-order isoparametric elements, which will further increase the computational complexity. In 2006, Krivodonova and Berger proposed a new curved element processing method for the DG format. It is based on straight-edge element calculations and uses the normal on the surface at the integration point, which significantly reduces the complexity of curved element processing, but the overall accuracy cannot be guaranteed. Li and Nishikawa used this method in their work and achieved high-order accuracy in the DG / FV and EB3 (Edge-Based Third-Oirder) formats, but the accuracy of the calculation results is not as good as the isoparametric method, and the accuracy can be further improved. Summary of the Invention

[0007] The technical problem to be solved by the present invention is to provide a wing boundary simulation method, terminal device and storage medium for optimizing isoparametric units in view of the shortcomings of the existing technology, so that the unstructured grid finite volume numerical simulation method can truly reflect the flight state of the wing.

[0008] To solve the above technical problems, the technical solution adopted by the present invention is: a wing boundary simulation method for optimizing isoparametric units, comprising the following steps:

[0009] S1. Obtain wing curve boundary information;

[0010] S2. Determine the shape of the curved element according to the curved boundary information. If it is a triangular mesh, construct a 10-point cubic element shape function of the triangular mesh; if it is a quadrilateral mesh, determine the isoparametric element order used by the quadrilateral curved element of the airfoil. If the order is three, construct a 10-point cubic element shape function of the quadrilateral mesh; if the order is five, construct a 17.5-order element shape function of the quadrilateral mesh;

[0011] S3, construct the Jacobian matrix of the wing physical unit and the isoparametric unit;

[0012] S4. Calculate the surface integral and / or volume integral on the isoparametric element using the Jacobian matrix.

[0013] In CFD methods, obtaining airfoil curve information is necessary for subsequent flow simulation. The key to flow simulation is performing volume and surface integrals on the flow governing equations. After obtaining the wing's curved boundary information, the present invention determines the curved unit shape based on this information. Based on the curved unit shape, the Jacobian matrix of the wing's physical unit and isoparametric unit is constructed to obtain the surface and volume integrals. This method eliminates the need to define the distribution of parameter points, significantly improving the accuracy of the surface and volume integral calculations.

[0014] In step S2, the shape function expression of the ten-point cubic element of the triangular mesh is:

[0015] ;

[0016] in, is the reference point coordinate of the standard isoparametric unit, is the shape function of the i-th parameter point of the isoparametric element.

[0017] In step S2, the shape function expression of the ten-point cubic element of the quadrilateral mesh is:

[0018] ;

[0019] in, The reference point coordinates of the standard isoparametric unit.

[0020] In step S2, the shape function expression of the 17.5-order element of the quadrilateral mesh is:

[0021] ;

[0022] in, The reference point coordinates of the standard isoparametric unit.

[0023] In step S3, the expression of the Jacobian matrix is:

[0024] ;

[0025] Where J is the Jacobian matrix, is the reference point coordinate of the standard isoparametric unit, The reference point coordinates are in actual physical units.

[0026] The surface integral expression is:

[0027] ;

[0028] in, and They refer to the two endpoints of the corresponding edge of the isoparametric unit in the isoparametric space. The coordinates of the direction, , represents the number of surface integrals, is the surface integral coefficient, f is the convection or viscous flux in the space outside the wing, which is related to physical quantities such as speed and temperature. is the flux of the kth Gaussian integration point (related to physical quantities such as speed and temperature), and is the Jacobian at the Gaussian integration point.

[0029] The volume integral expression is:

[0030] ;

[0031] in, Indicates the number of volume points, is the volume fraction coefficient, is the source function, is the corresponding reference point coordinate of the physical unit, n is the number of reference points of the isoparametric unit, is the shape function, is the coordinate of the qth Gaussian integral point in the isoparametric space, is the Jacobian matrix at the qth Gaussian integration point.

[0032] As an inventive concept, the present invention also provides a terminal device, including a memory, a processor, and a computer program stored in the memory; the processor executes the computer program to implement the steps of the above method.

[0033] As an inventive concept, the present invention also provides a computer-readable storage medium having a computer program / instruction stored thereon; the computer program / instruction implements the steps of the above method when executed by a processor.

[0034] As an inventive concept, the present invention also provides a computer program product, comprising a computer program / instruction; when the computer program / instruction is executed by a processor, the steps of the above method are implemented.

[0035] Compared with existing technologies, the present invention has the following advantages: after obtaining wing curved boundary information, the present invention determines the curved unit shape based on this information, constructs the Jacobian matrix of the wing physical unit and the isoparametric unit based on the curved unit shape, and then obtains the surface integral and volume integral. This invention does not require the definition of the distribution of parameter points, greatly improving the calculation accuracy of the surface integral and volume integral. BRIEF DESCRIPTION OF THE DRAWINGS

[0036] Figure 1 It is a quadrilateral mesh NACA0012 airfoil;

[0037] Figure 2 It is a schematic diagram of the standard cubic triangle physical unit and the isoparametric unit;

[0038] Figure 3 To simplify the schematic diagram of cubic triangle physical unit and isoparametric unit;

[0039] Figure 4 It is a simplified schematic diagram of a cubic quadrilateral isoparametric unit;

[0040] Figure 5 Possible reference point arrangements for standard cubic quadrilateral isoparametric elements;

[0041] Figure 6 It is a schematic diagram of a standard cubic quadrilateral isoparametric unit;

[0042] Figure 7 It is a schematic diagram of a 17-point quintic isoparametric unit in a quadrilateral;

[0043] Figure 8 This is a flow chart of a method according to an embodiment of the present invention;

[0044] Figure 9 The sparsest triangular and quadrilateral meshes for inviscid flow around a circular cylinder; (a) Schematic diagram of the overall mesh for flow around a cylinder, (b) an enlarged diagram of the cylinder boundary;

[0045] Figure 10 The errors obtained by different curved element discretization methods on quadrilateral meshes under third-order and fourth-order formats; (a) overall error of all elements, (b) boundary element error;

[0046] Figure 11 Entropy error statistics of different discretization methods on triangular meshes; (a) overall error of all cells, (b) boundary cell error;

[0047] Figure 12 The pressure coefficient distribution of the cylinder surface obtained by different curved unit discretization methods; (a) the pressure coefficient distribution of the upper surface of the cylinder, (b) the pressure coefficient distribution of the local magnification of the upper surface of the cylinder; Figure 12 The horizontal axis is the horizontal coordinate of the wing boundary point;

[0048] Figure 13 Integral of the drag coefficient obtained by different methods on the quadrilateral grid C D ; (a) sparse grid, (b) medium grid, (c) dense grid;

[0049] Figure 14 The position relationship of the integration points based on the straight-edge and curved-edge unit discretization methods; (a) The position relationship of the integration points based on the straight-edge unit discretization method, (b) The position of the integration points based on the curved-edge unit discretization method;

[0050] Figure 15 Integral of the drag coefficient obtained by different methods on the triangular mesh C D ; (a) sparse grid, (b) medium grid, (c) dense grid. DETAILED DESCRIPTION

[0051] To make the objectives, technical solutions, and advantages of the embodiments of the present invention more clear, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts shall fall within the scope of protection of the present invention.

[0052] Example 1

[0053] The airfoil curved boundary expression requires the use of a curve function to represent the boundary of the airfoil's geometric model, which involves the mapping from the geometric model to the computational model.

[0054] Mapping the geometric model to the computational model requires three steps: generating a mesh, building a coordinate system, and calculating coordinates. Figure 1 Taking the NACA0012 airfoil shown in the figure as an example, the three steps are explained.

[0055] 1. Generate a mesh

[0056] like Figure 1 As shown in the figure, firstly, the geometric model boundary of the NACA0012 airfoil is filled with points for discretization of the geometric boundary. These points are called grid points, and adjacent grid points are connected to form grid edges. After the grid edges on the geometric boundary are formed, the mesh generation software can be used to generate the mesh of the external calculation area of ​​the NACA0012 airfoil, as shown in the figure. Figure 1As shown, a quadrilateral mesh is generated, and a triangular mesh can also be generated. Although the mesh edges of unstructured meshes are composed of line segments, the airfoil's curved boundary can be accurately expressed using the high-order curved boundary representation method (Gao H, Wang Z, Liu Y. A study of curved boundary representations for 2D high-order Euler solvers [J]. Journal of Scientific Computing. 2010, 44: 323–336.). This is the basis for constructing a high-order precision format for airfoil de-elementing and performing high-order numerical integration.

[0057] 2. Constructing the coordinate system

[0058] You need to construct the XOY plane. Set the leading edge of the NACA0012 airfoil to point O. Set the x-axis along the chord line of the NACA0012 airfoil, the line connecting the leading and trailing edges, with the positive x-axis direction. Set the y-axis perpendicular to the x-axis through point O, with the positive y-axis direction pointing upward.

[0059] 3. Calculate coordinates

[0060] The construction of the numerical format is related to the grid point coordinates. For each grid point in the calculation model, the grid coordinates of the grid point in the XOY plane are calculated according to the grid scale. .

[0061] After completing the above three steps, the high-order precision numerical format of the airfoil curved element can be constructed.

[0062] The standard isoparametric element discretization method is introduced below.

[0063] The isoparametric element method (Isoparametric Element Method) is an effective and widely used technique. It simplifies the computation and analysis process by mapping complex curved elements to standard isoparametric elements. The key to the Isoparametric Element Method lies in the construction of the isoparametric mapping.

[0064] For two-dimensional problems, a standard reference unit is usually selected, such as an equilateral triangle, an isosceles right triangle, or a unit square, to map the curved unit onto the reference unit. The mapping function can be defined by an interpolation function. Assume that the reference point coordinates of the standard isoparametric unit are , the reference point coordinates of the actual physical unit are The isoparametric mapping can be expressed as:

[0065] ;

[0066] Where, is the shape function, n is the number of reference points of the isoparametric unit, is the corresponding reference point coordinate of the physical unit.

[0067] The mapping between airfoil curved elements and isoparametric elements is introduced below.

[0068] Wang and Liu discussed quadratic and cubic triangular isoparametric elements for the spectral volume format, and considering that in most cases a triangular element will contain at most one curved edge, they further simplified the standard quadratic and cubic triangular isoparametric elements.

[0069] Since the embodiments of the present invention mainly explore the spatial discretization and numerical integration on curved elements in third-order and fourth-order accuracy unstructured finite volume schemes, quadratic isoparametric elements are no longer considered. Figure 2 and Figure 3 Schematics of standard and simplified cubic triangle physical units and isoparametric units are shown respectively.

[0070] As can be seen, a standard cubic triangle isoparametric element requires 10 reference points, in which case a unique shape function can be obtained at each point. While the simplified isoparametric element omits the reference points on the two straight edges, there are theoretically an infinite number of possible shape functions at each point, introducing significant uncertainty into the calculation process.

[0071] In the embodiment of the present invention, the default method for triangular units is Figure 2 The standard case is shown in , and in this arrangement the shape function at each point is:

[0072] ;

[0073] Furthermore, the mapping relationship between the song unit and the isoparametric unit can be obtained according to the isoparametric mapping.

[0074] Currently, triangular isoparametric elements are a more common choice, while quadrilateral isoparametric elements are less discussed. Li used simplified quadrilateral isoparametric elements in his work, such as Figure 4 The standard quadrilateral isoparametric elements are shown, but not discussed.

[0075] The following introduces the spatial discretization method based on isoparametric elements.

[0076] Based on the shape functions of the airfoil triangle and quadrilateral isoparametric elements, the Jacobian between the physical element and the isoparametric element of the airfoil can be obtained:

[0077] ;

[0078] Flux f and source terms s Taking as an example, the forms of surface integral and volume integral on isoparametric units are given.

[0079] Surface integral over isoparametric elements

[0080] The surface integral on the curved element can be split into two parts: the integral on the straight edge and the integral on the curved edge:

[0081] ;

[0082] The integral on the straight edge is easy to calculate, and the following mainly discusses the integral on the curved edge. For the curved infinitesimal element:

[0083] ;

[0084] Combining the relationship between the tangent vector and the normal vector, we can know the outward normal on the curved edge:

[0085] ;

[0086] Therefore, the flux integral over the curved edge can be expressed as:

[0087] ;

[0088] At this point, the integral on the physical unit is converted to the isoparametric unit, where and Refers to the two endpoints of the corresponding edge of the isoparametric unit ( ). Gaussian numerical integration scheme for combined conservation law equations:

[0089] ;

[0090] in, is the conserved variable in the control volume The average value on and are the convection and viscous fluxes along the normal direction outside the unit surface. According to the discretization scheme of the high-order finite volume format, in the control volume unit i superior, The semi-discrete form of can be expressed as:

[0091] ;

[0092] Represents the control volume element i The number of unit surfaces, and represent the number of surface integral and volume integral points respectively, and and are the coefficients of surface integral and volume integral. The integral of the flux on the curved edge can finally be expressed as:

[0093] ;

[0094] In the formula and is the Jacobian at the Gaussian integration point, which can be obtained according to Calculated.

[0095] Volume integral over isoparametric elements:

[0096] For volume integrals, there is the following relationship between isoparametric units and physical units:

[0097] ;

[0098] Therefore, the volume integral of the source function s(x,y) can be expressed as:

[0099] ;

[0100] Similarly, we can get:

[0101] ;

[0102] The following describes an optimization isoparametric unit according to an embodiment of the present invention.

[0103] The triangular isoparametric elements of the airfoil are still mapped according to the above method. For the quadrilateral isoparametric elements, if they are arranged according to the standard cubic quadrilateral isoparametric elements, that is, according to Figure 2 The most intuitive way to arrange the reference points of the quadrilateral elements is as follows: Figure 5 shown.

[0104] Figure 5 In the figure, reference point 10 is the geometric center of the quadrilateral. Unfortunately, based on this reference point arrangement, the shape function cannot be obtained. Therefore, the embodiment of the present invention adopts Figure 6 The reference point arrangement shown in .

[0105] At this time, reference point 10 is located at the point where the vertical center line is divided into three equal parts near the bottom. Based on this reference point arrangement, the number of shapes obtained is:

[0106] ;

[0107] On this basis, the degree of the quadrilateral isoparametric element can be further improved, that is, on the basis of the quartic shape function, two quintic terms are added to construct a more accurate 17 point quintic quadrilateral isoparametric element. The shape function used is:

[0108] ;

[0109] The parameter points are arranged as follows Figure 7 As shown. We found that this reference point arrangement is more symmetrical, and point 17 inside the quadrilateral is exactly the geometric center of the quadrilateral:

[0110] ;

[0111] In the embodiment of the present invention, the airfoil triangle element adopts the standard 10-point cubic isoparametric element; the airfoil quadrilateral element adopts the standard 10-point cubic and 17-point quintic isoparametric elements respectively.

[0112] The high-order discretization method of optimizing the airfoil curved boundary of the isoparametric unit adopts the method of optimizing the isoparametric unit to construct the airfoil triangle unit and quadrilateral unit isoparametric unit, and then adopts the spatial discretization method to perform surface integral and volume integral. The flow chart is as follows Figure 8 shown.

[0113] 1. First, read the airfoil curve boundary information after high-order representation;

[0114] 2. Determine the shape of the curved element. If it is a triangular mesh, turn 3; if it is a quadrilateral mesh, turn 4;

[0115] 3. Construct the ten-point cubic element shape function of the triangular mesh and go to 7;

[0116] 4. Determine the order of the isoparametric element used for the quadrilateral curved element of the airfoil. The third order is converted to 5, and the fifth order is converted to 6;

[0117] 5. Construct the 10-point cubic element shape function of the quadrilateral mesh, go to 7;

[0118] 6. Construct the 17.5-order element shape function of the quadrilateral mesh, go to 7;

[0119] 7. Construct the Jacobian of the airfoil physical unit and the isoparametric unit;

[0120] 8. Solve the surface integral or volume integral of isoparametric elements;

[0121] 9. End.

[0122] The high-order curved boundary expression method of airfoil is used to test the calculation results of the curved unit discretization method based on quintic polynomial curve.

[0123] Inviscid flow around a circular cylinder example:

[0124] This example is a typical isentropic outflow problem. In this section, we consider the incoming flow Mach number as The calculation result is the same as that of 0.38. The calculation domain far field is 30 times the cylinder diameter, and the cylinder diameter D= 1. This section uses six sets of triangular and quadrilateral meshes with different densities to verify different methods.

[0125] The incoming flow Mach number is Ma=0.1, and the calculation domain far field is 30 times the cylinder diameter. In this paper, the cylinder diameter is taken as D = 1. This part uses seven sets of triangular and quadrilateral meshes with different densities for verification. The number of triangular mesh elements ranges from 2400 to 17496, and the number of quadrilateral mesh elements ranges from 2500 to 12100. Figure 9 The sparsest triangle and quad mesh is shown in .

[0126] Figure 10 The entropy error of the embodiment of the present invention and the curved unit discretization method and the straight-edge unit discretization method in the third-order format is shown in the figure. The third-order format - cubic isoparametric unit and the third-order format - quintic isoparametric unit are respectively the third-order cubic and quintic isoparametric unit curved unit discretization methods of the embodiment of the present invention, the fourth-order format - cubic isoparametric unit and the fourth-order format - quintic isoparametric unit are respectively the fourth-order cubic and quintic isoparametric unit curved unit discretization methods of the embodiment of the present invention, the third-order format - simplified curved-edge unit and the fourth-order format - simplified curved-edge unit are respectively the third-order and fourth-order curved unit discretization methods, the third-order format - straight-edge unit and the fourth-order format - straight-edge unit are respectively the straight-edge unit discretization methods. It can be seen that for the third-order format, the results obtained by using the two curved unit discretization methods on all units and wall units are very close, and there is no obvious loss of accuracy as the grid is continuously encrypted. Even when the discretization is completely based on straight-edge units, the obtained results can still be maintained near the third-order accuracy.

[0127] When the accuracy is increased to the fourth order, the fourth-order format experiences a faster rate of error reduction as the mesh becomes denser. Using the embodiment of the present invention and the two curved element discretization schemes, the calculation results consistently maintain accuracy above the fourth order for all elements and wall elements. However, combining the specific data in Tables 1 and 2, it can be seen that using straight-edge element discretization results in accuracy loss as the mesh becomes denser. In particular, between the last two sets of dense meshes, the accuracy of the wall elements decreases to 2.5 order, far less than the embodiment of the present invention and the fourth-order accuracy. This demonstrates the necessity of using the curved element discretization scheme for high-precision numerical formats.

[0128] Table 1 Discretization errors of different discretization methods on the densest quadrilateral grid

[0129]

[0130] Table 2 Calculation accuracy of different discretization methods between the last two sets of dense grids

[0131]

[0132] On this basis, Figure 11 The entropy error obtained by different methods on triangular meshes is further demonstrated in the figure. The results show that for the third-order format, the results obtained using the two curved element discretization methods and the straight-edge element discretization method consistently maintain near third-order accuracy for all element types. Furthermore, for wall elements, the computational accuracy of the three methods gradually approaches third-order as the mesh becomes denser.

[0133] However, in the fourth-order format, the two curved element discretization methods yielded similar results. Combined with the specific data presented in Table 3, it can be seen that the entropy error accuracy for all elements exceeded the fourth order, and for wall elements, the calculation accuracy between the two densest meshes exceeded 3.7 orders. However, based on the straight-edge element discretization, as the mesh density increases, the entropy error calculation results for wall elements degrade significantly, reaching only 1.065 orders between the two densest meshes. Therefore, the correctness of the embodiment of the present invention is preliminarily verified.

[0134] Table 3 Discretization error and calculation accuracy of different methods on the densest grid on triangular mesh

[0135]

[0136] Pressure coefficient distribution:

[0137] For the inviscid flow around a circular cylinder discussed in this section, the pressure coefficient has an analytical form. The exact and numerical solutions for the pressure coefficient distribution on the cylindrical surface are analyzed below.

[0138] Since the numerical example used in the test is the flow Mach number , it can be regarded as an incompressible flow. Considering that the pressure coefficient is defined as:

[0139] ;

[0140] According to Bernoulli's equation, the pressure distribution on the cylindrical surface can be obtained as:

[0141] ;

[0142] Will Substituting the Bernoulli equation and simplifying it, we can get:

[0143] ;

[0144] Substituting the above formula into the definition of pressure coefficient given in the definition of pressure coefficient, we can obtain:

[0145] ;

[0146] Under the conditions of third-order and fourth-order formats, the distribution results of the cylindrical surface pressure coefficient obtained by three different discretization methods on a 30×120 sparse grid are as follows Figure 12 As shown, the exact solution of the pressure coefficient is the formula The exact solution for the calculated pressure coefficient.

[0147] from Figure 12 The results show that on the sparse grid, the results obtained by the three discretization schemes in the third-order format are close, with some deviation from the exact value of the pressure coefficient; however, in the fourth-order format, the calculation accuracy is further improved, and the pressure coefficient curves obtained by the three methods are almost completely consistent with the exact curve. Figure 12 The zoomed-in view in the right figure further reveals that, compared with the embodiment of the present invention and the two curved element discretization methods, the discretization based entirely on straight-edge elements in the fourth-order format results in a significant error in the cylindrical surface pressure coefficient. This is consistent with the pattern observed in the accuracy test.

[0148] Example of flow around a viscous cylinder:

[0149] This part is about the subsonic viscous flow around a circular cylinder. The given Mach number and Reynolds number of the incoming flow are , .

[0150] Considering that for finite volume schemes, third-order accuracy schemes can only achieve second-order accuracy for viscosity problems, this section only considers fourth-order schemes. The following two sections will present the results on quadrilateral and triangular meshes.

[0151] Quad mesh:

[0152] Figure 13 The pressure coefficient is shown in C D The curve of changes with the number of iteration steps, where the cubic isoparametric unit and the quintic isoparametric unit are respectively the cubic and quintic optimized isoparametric unit curved unit discretization methods of the embodiments of the present invention, the simplified curved edge unit is the curved unit discretization method, the straight edge unit refers to the straight edge unit discretization method, and the surface integral / volume integral results are the standard Gaussian surface integral and volume integral results based on the cylindrical analytical expression.

[0153] The drag coefficient curves further reveal that the results of the curved element discretization method using the optimized isoparametric element are nearly identical to those of the standard discretization method, while the other two methods exhibit significant deviations, particularly the straight-edge element discretization method based on the curved boundary normal. Preliminary analysis shows that for both the isoparametric element method and the standard discretization method, when calculating the pressure coefficient curve and the drag coefficient integral curve on the wall boundary, the sampling points are located at the integration points of the reconstructed polynomial curve and the analytical cylindrical curve; whereas the sampling points for the other two methods are all located at the integration points of the straight-edge elements.

[0154] The two integral points are used C Gand S G To express, combined Figure 14 It can be seen that the pressure coefficient can be expressed as:

[0155] ;

[0156] The drag coefficient is

[0157] ;

[0158] In the formula and are the lengths of the straight and curved boundaries of the element, respectively. Therefore, even if the discretization scheme of straight-edge elements based on the curved boundary normal can guarantee the flow results of corresponding accuracy, the calculation of pressure and drag coefficients involves the actual integration point location and integration along the boundary surface. In this case, the integration points and the numerical integration process are all along the straight boundary, resulting in significant differences in the calculation results compared with the discretization method based on curved-edge elements.

[0159] Triangular mesh:

[0160] Figure 15 From the curve showing the change of the integral of the cylindrical surface resistance coefficient with the number of iteration steps, it can be seen that the pressure coefficient results obtained by the optimized isoparametric unit method in the embodiment of the present invention are close to those obtained by the standard discretization method based on the cylindrical analytical curve, while there is a certain deviation between the results obtained by the two discretization methods based on straight-edge units and the curved unit discretization method.

[0161] In summary, the numerical performance of different methods was further compared on quadrilateral and triangular meshes using a subsonic viscous flow example around a circular cylinder. The results show that the use of a quintic polynomial curve combined with the optimized isoparametric element discretization scheme according to the present invention yields results very close to those of the standard finite volume discretization scheme based on an analytical representation of the cylinder.

[0162] Example 2

[0163] Embodiment 2 of the present invention provides a measurement system corresponding to the above-mentioned embodiment 1. The measurement system may be a processing device for a client, such as a mobile phone, a laptop, a tablet computer, a desktop computer, etc., to execute the method of the above-mentioned embodiment.

[0164] The measurement system of this embodiment includes a memory, a processor, and a computer program stored in the memory; the processor executes the computer program in the memory to implement the steps of the method in the above-mentioned embodiment 1.

[0165] In some implementations, the memory may be a high-speed random access memory (RAM), and may also include a non-volatile memory, such as at least one disk storage.

[0166] In other implementations, the processor may be a central processing unit (CPU), a digital signal processor (DSP), or other general-purpose processors, which are not limited herein.

[0167] Although the preferred embodiments of the present application have been described, those skilled in the art may make additional changes and modifications to these embodiments once they have learned the basic creative concept. Therefore, the appended claims are intended to be interpreted as including the preferred embodiments and all changes and modifications that fall within the scope of the present application.

[0168] Obviously, those skilled in the art may make various changes and modifications to this application without departing from the spirit and scope of this application. Thus, if these modifications and variations of this application fall within the scope of the claims of this application and their equivalents, this application is intended to include these modifications and variations.

Claims

1. A wing boundary simulation method for optimizing isoparametric units, characterized in that: The following steps are involved: S1. Obtain wing curve boundary information; S2. Determine the shape of the curved element according to the curved boundary information. If it is a triangular mesh, construct a 10-point cubic element shape function of the triangular mesh; if it is a quadrilateral mesh, determine the isoparametric element order used by the quadrilateral curved element of the airfoil. If the order is three, construct a 10-point cubic element shape function of the quadrilateral mesh; if the order is five, construct a 17.5-order element shape function of the quadrilateral mesh; S3, construct the Jacobian matrix of the wing physical unit and the isoparametric unit; S4. Calculate the surface integral and / or volume integral on the isoparametric element using the Jacobian matrix.

2. The wing boundary simulation method of optimizing isoparametric units according to claim 1, characterized in that: In step S2, the shape function expression of the ten-point cubic element of the triangular mesh is: ； in, is the reference point coordinate of the standard isoparametric unit, is the shape function of the i-th parameter point of the isoparametric element.

3. The wing boundary simulation method of optimizing isoparametric units according to claim 1, characterized in that: In step S2, the shape function expression of the ten-point cubic element of the quadrilateral mesh is: ； in, The reference point coordinates of the standard isoparametric unit.

4. The wing boundary simulation method of optimizing isoparametric units according to claim 1, characterized in that: In step S2, the shape function expression of the 17.5-order element of the quadrilateral mesh is: ； in, The reference point coordinates of the standard isoparametric unit.

5. The wing boundary simulation method of optimizing isoparametric units according to claim 1, characterized in that: In step S3, the expression of the Jacobian matrix is: ； Where J is the Jacobian matrix, is the reference point coordinate of the standard isoparametric unit, The reference point coordinates are in actual physical units.

6. The wing boundary simulation method of optimizing isoparametric units according to claim 1, characterized in that: The surface integral expression is: ； in, and They refer to the two endpoints of the corresponding edge of the isoparametric unit in the isoparametric space. The coordinates of the direction, , represents the number of surface integrals, is the surface integral coefficient, f is the convection or viscous flux outside the wing, is the flux of the kth Gaussian integration point, and is the Jacobian at the Gaussian integration point.

7. The wing boundary simulation method of optimizing isoparametric units according to claim 1, characterized in that: The volume integral expression is: ； in, Indicates the number of volume points, is the volume fraction coefficient, is the source function, , is the coordinate of the corresponding reference point of the physical unit, n is the number of reference points of the isoparametric unit, is the shape function, is the coordinate of the qth Gaussian integral point in the isoparametric space, is the Jacobian matrix at the qth Gaussian integration point.

8. A terminal device comprising a memory, a processor, and a computer program stored in the memory; characterized in that: The processor executes the computer program to implement the steps of the method according to any one of claims 1 to 7.

9. A computer-readable storage medium having a computer program / instruction stored thereon; characterized in that: When the computer program / instructions are executed by a processor, the steps of the method according to any one of claims 1 to 7 are implemented.

10. A computer program product comprising a computer program / instructions; characterized in that When the computer program / instructions are executed by a processor, the steps of the method according to any one of claims 1 to 7 are implemented.

Images