Wing boundary simulation method for optimizing isoparametric unit, terminal equipment and storage medium

Abstract

The invention discloses a wing boundary simulation method for optimizing an isoparametric unit, a terminal device and a storage medium, and the method comprises the steps: obtaining wing curved boundary information, judging the shape of a curved unit according to the wing curved boundary information, constructing a Jacobian matrix of a physical unit and an isoparametric unit of a wing according to the shape of the curved unit, and obtaining an area component and a volume component; according to the method, distribution of parameter points does not need to be defined, and the calculation precision of the area component and the volume component is greatly improved.

Description

Technical Field

[0001] The present invention relates to the airfoil simulation technology of an aircraft, in particular to a wing boundary simulation method, a terminal device and a storage medium for optimizing isoparametric elements. Background Art

[0002] The airfoil, commonly known as the wing section or blade section, is one of the core factors affecting the comprehensive performance of an aircraft. The airfoil is the basic element for the aerodynamic surface shape design of aircraft wings, tails, helicopter rotors, propellers, and wind turbine blades. It directly affects the aerodynamic performance of the aircraft, including key aerodynamic characteristics such as lift, drag, and pitching moment. The aerodynamic characteristics of the airfoil can also directly affect the maneuverability and stability of the aircraft. A good airfoil design can provide sufficient control moments to ensure the stability and maneuverability of the aircraft under various flight conditions, thereby improving flight safety.

[0003] In addition, the airfoil is also closely related to the energy conservation, emission reduction, and consumption reduction of the aircraft. Optimizing the airfoil design is one of the effective measures to reduce drag and energy consumption and reduce fuel consumption. The airfoil design determines the lift-to-drag ratio of the aircraft [Li Yunpeng, Han Yongzhi. Research progress on the design of lift augmentation devices based on laminar wings [J]. Advances in Aeronautical Science and Engineering, 2021, 12(4): 1-11.], that is, the ratio of lift to drag. The higher the lift-to-drag ratio, the higher the flight efficiency, and it 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.). An excellent airfoil design can achieve a high lift-to-drag ratio, effectively reduce fuel consumption, reduce carbon emissions, increase the range and endurance time, which is crucial for improving the flight efficiency and economy of the aircraft.

[0004] In order to obtain excellent airfoils with strong stability and good drag reduction effects economically and efficiently, computational fluid dynamics methods are usually adopted 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 during the flight of the wing. The process of discretizing the physical space is called mesh generation. The commonly used meshes in engineering applications include structured meshes, unstructured meshes, and Cartesian meshes. Structured meshes have a high generation efficiency in simple geometric shapes, but it becomes increasingly difficult to generate meshes in complex situations. The generation of Cartesian meshes is simple and can provide high computational accuracy in regular computational domains, but in complex geometric shapes, the situation of mesh and boundary mismatch is likely to occur, and special treatment methods such as the immersed boundary method are required, which increases the complexity of implementation. Unstructured meshes have strong automatic generation capabilities and can flexibly adapt to any complex geometric shape, and have been widely used in practical engineering problems in recent decades. The numerical formats on unstructured meshes include Finite Volume (FV), Spectral Volume (SV), Discontinuous Galerkin (DG), PNPM, and Correction Procedure via Reconstruction (CPR), etc. Among them, the derivation and implementation process of the finite volume format is relatively simple and is easy to be extended to high-order accuracy. In recent years, the emergence of the Compact Finite Volume (CFV) format and the mature shock capturing technology has continuously enhanced the influence of the finite volume format.

[0005] In order to enable the unstructured mesh finite volume numerical simulation method to truly reflect the flight state of the wing, two aspects need to be particularly noted. One is that for the unstructured meshes generated in engineering applications, the mesh edges are all line segments, while the airfoil generally has a curved boundary. Using line segments to replace the curved boundary makes the numerical shape of the wing slightly different from the actual shape, which may lead to inaccurate numerical simulation results. Therefore, high-order curved boundary representation of the airfoil is required. The other is that in order to improve the computational accuracy and solution efficiency, high-precision numerical formats are often used for solution in engineering applications. However, high-precision formats are often used for the calculation inside the computational domain, and it is often difficult to obtain high-order accuracy for the numerical format on the curved boundary of the airfoil. And the values at the airfoil boundary directly affect the calculation results of the entire computational domain, that is, the numerical simulation of the curved boundary of the airfoil cannot reach high-order accuracy, and the overall accuracy will surely be affected.

[0006] The object discretized by the finite volume method is the conservation law equation in integral form. Therefore, one of the keys to achieving high-order accuracy lies in the high-order numerical integration that matches the format accuracy, such as the flux area integral and the source term volume integral. On straight-sided elements, the implementation of high-order numerical integration is relatively simple, and accurate results can be obtained according to the integration points and weight coefficients 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, in the case of elements with curved sides, especially in the problem of the curved boundary of an airfoil, the numerical integration process faces many challenges. There are currently two common treatment methods. The first is the isoparametric element method commonly used in the finite element field, which realizes equivalent calculations on isoparametric elements by constructing a mapping relationship between the physical space and the isoparametric space. However, this method requires defining the distribution of parameter points, solving and storing Lagrangian interpolation coefficients and Jacobians, so the process is slightly cumbersome, especially in high-order isoparametric elements, which will further increase the computational complexity. Krivodonova and Berger proposed a new curved element treatment method for the DG format in 2006. It is based on the calculation of straight-sided elements and uses the normal direction on the surface at the integration points, significantly reducing the complexity of curved element treatment, but the overall accuracy cannot be guaranteed. Li and Nishikawa followed this method in their work and achieved high-order accuracy in the DG / FV and EB3 (Edge-Based Third-Order) formats, but the accuracy of the calculation results is not as good as that of 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, a terminal device, and a storage medium that optimize isoparametric elements in view of the deficiencies of the prior art, 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 problem, the technical solution adopted by the present invention is: a wing boundary simulation method that optimizes isoparametric elements, including the following steps: S1. Obtain the wing curved boundary information; S2. Determine the shape of the curved element based on the curved boundary information. If it is a triangular mesh, construct the ten-node cubic element shape functions for the triangular mesh. If it is a quadrilateral mesh, determine the order of the isoparametric element used for the airfoil quadrilateral curved element. If the order is three, construct the ten-node cubic element shape functions for the quadrilateral mesh. If the order is five, construct the seventeen-node quintic element shape functions for the quadrilateral mesh. S3. Construct the Jacobian matrix of the wing physical element and the isoparametric element. S4. Use the Jacobian matrix to calculate the surface integral and / or volume integral on the isoparametric element.

[0009] In the CFD method, obtaining the curve information of the airfoil is for subsequent flow simulation. The key to flow simulation is to perform volume integration and surface integration on the flow control equations. After obtaining the wing curved boundary information in the present invention, the shape of the curved element is judged according to the wing curved boundary information, and the Jacobian matrix of the wing physical element and the isoparametric element is constructed according to the shape of the curved element, and then the surface integral and volume integral are obtained. The present invention does not need to define the distribution of parameter points, which greatly improves the calculation accuracy of the surface integral and volume integral.

[0010] In step S2, the expression of the ten-node cubic element shape functions for the triangular mesh is: ; where is the reference point coordinates of the standard isoparametric element, is the shape function of the i-th parameter point of the isoparametric element.

[0011] In step S2, the expression of the ten-node cubic element shape functions for the quadrilateral mesh is: ; where is the reference point coordinates of the standard isoparametric element.

[0012] In step S2, the expression of the seventeen-node quintic element shape functions for the quadrilateral mesh is: ; where is the reference point coordinates of the standard isoparametric element.

[0013] In step S3, the expression of the Jacobian matrix is: ; where J is the Jacobian matrix, is the reference point coordinates of the standard isoparametric element, is the reference point coordinates of the actual physical element.

[0014] The expression of the surface integral is: ; Among them, and respectively refer to the coordinates of the two endpoints of the corresponding side of the isoparametric element in the isoparametric space in the direction, , represents the number of surface integrals, is the surface integral coefficient, f is the convective or viscous flux in the external space of the wing, which is related to physical quantities such as velocity and temperature, is the flux at the k-th Gauss integration point (related to physical quantities such as velocity and temperature), and is the Jacobian at the Gauss integration point.

[0015] The volume integral expression is: ; Among them, represents the number of volume integration points, is the volume integral coefficient, is the source term function, is the coordinate of the corresponding reference point of the physical element, n is the number of reference points of the isoparametric element, is the shape function, is the coordinate of the q-th Gauss integration point in the isoparametric space, is the Jacobian matrix at the q-th Gauss integration point.

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

[0017] As an inventive concept, the present invention also provides a computer-readable storage medium, on which a computer program / instructions are stored; when the computer program / instructions are executed by a processor, the steps of the above method are implemented.

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

[0019] Compared with the prior art, the beneficial effects of the present invention are as follows: After obtaining the wing curved boundary information, the present invention judges the shape of the curved element according to the wing curved boundary information, constructs the Jacobian matrix of the wing physical element and the isoparametric element according to the shape of the curved element, and then obtains the surface integral and the volume integral. The present invention does not need to define the distribution of parameter points, greatly improving the calculation accuracy of the surface integral and the volume integral. Description of the Drawings

[0020] Figure 1 is a quadrilateral grid NACA0012 airfoil; Figure 2 is a schematic diagram of a standard cubic triangular physical element and an isoparametric element; Figure 3 is a schematic diagram of a simplified cubic triangular physical element and an isoparametric element; Figure 4 is a schematic diagram of a simplified cubic quadrilateral isoparametric element; Figure 5 is a possible arrangement of reference points for a standard cubic quadrilateral isoparametric element; Figure 6 is a schematic diagram of a standard cubic quadrilateral isoparametric element; Figure 7 is a schematic diagram of a 17 - point quintic isoparametric element for a quadrilateral; Figure 8 is a flowchart of the method according to an embodiment of the present invention; Figure 9 is the sparsest triangular and quadrilateral grids for inviscid cylinder flow; (a) overall schematic of the cylinder - flow grid, (b) magnified schematic of the cylinder boundary; Figure 10 is the error obtained by different curved - element discretization methods in third - order and fourth - order formats on a quadrilateral grid; (a) overall error of all elements, (b) error of boundary elements; Figure 11 is the entropy - error statistics of different discretization methods on a triangular grid; (a) overall error of all elements, (b) error of boundary elements; Figure 12 is the distribution of the cylinder - surface pressure coefficient obtained by different curved - element discretization methods; (a) distribution of the pressure coefficient on the upper surface of the cylinder, (b) locally magnified distribution of the pressure coefficient on the upper surface of the cylinder; Figure 12 The horizontal axis of Figure 13 is the drag - coefficient integral obtained by different methods on a quadrilateral grid C D ; (a) sparse grid, (b) medium grid, (c) dense grid; Figure 14 is the relationship of the integration - point positions based on straight - edge and curved - edge element discretization methods; (a) relationship of the integration - point positions based on straight - edge element discretization method, (b) positions of the integration points based on curved - edge element discretization method; Figure 15 is the drag - coefficient integral obtained by different methods on a triangular grid C D ; (a) sparse grid, (b) medium grid, (c) dense grid. Detailed implementation manners

[0021] In order to make the purpose, technical solution and advantages of the embodiments of the present invention clearer, the technical solution in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are 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 creative work are within the scope of protection of the present invention.

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

[0023] 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.

[0024] 1. Generate a mesh 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 grid generation software can be used to generate grids for the external calculation area of the NACA0012 airfoil, such as Figure 1 As shown in the figure, a quadrilateral mesh is generated, and a triangular mesh can also be generated. Although the mesh edges of the unstructured mesh are composed of line segments, the curved boundary of the airfoil can be accurately expressed through the high-order curved boundary representation method of the airfoil (Gao H, Wang Z, Liu Y. A study of curved boundary representationsfor 2D high-orderEuler solvers [J]. Journal of Scientific Computing. 2010,44: 323–336.), which is the basis for constructing a high-order precision format for airfoil removal and performing high-order numerical integration.

[0025] 2. Constructing the coordinate system The XOY plane needs to be constructed. Set the leading edge point of the NACA0012 airfoil as point O. Set the x-axis along the chord line of the NACA0012 airfoil, that is, the straight line connecting the leading edge and the trailing edge, and the direction from the leading edge to the trailing edge is the positive direction of the x-axis. Through point O, the direction perpendicular to the x-axis is set to the y-axis, and the upward direction is the positive direction of the y-axis.

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

[0027] After completing the above three steps, a high-order accurate numerical format for the airfoil curved element can be constructed.

[0028] The following introduces the standard isoparametric element discretization method.

[0029] The Isoparametric Element Method is an effective and widely used technique. This method simplifies the calculation 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.

[0030] For two-dimensional problems, a standard reference element is usually selected, such as an equilateral triangle, an isosceles right triangle, or a unit square, etc., and the curved element is mapped onto this reference element. The mapping function can be defined by interpolation functions. Assuming the reference point coordinates of the standard isoparametric element are , and the reference point coordinates of the actual physical element are . The isoparametric mapping can be expressed as: ; where is the shape function, n is the number of reference points of the isoparametric element, are the corresponding reference point coordinates of the physical element.

[0031] The following introduces the mapping between the airfoil curved element and the isoparametric element.

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

[0033] Since the embodiments of the present invention mainly explore the spatial discretization and numerical integration on curved elements in the third- and fourth-order accurate unstructured finite volume format, quadratic isoparametric elements are no longer considered. Figure 2 And Figure 3 respectively show the schematic diagrams of the standard and simplified cubic triangular physical elements and isoparametric elements.

[0034] It can be seen that the standard cubic triangle isoparametric unit requires a total of 10 reference points. In this case, a unique shape function can be obtained at each point. Although the simplified isoparametric unit omits the reference points on the two straight edges, there are theoretically an infinite number of possible shape functions at each point, which introduces a lot of uncertainty into the calculation process.

[0035] In the embodiment of the present invention, the default method for the triangle unit is Figure 2 The standard case is shown in , and in this arrangement the shape function at each point is: ; Furthermore, the mapping relationship between the song unit and the isoparametric unit can be obtained according to the isoparametric mapping.

[0036] Currently, triangular isoparametric elements are a 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.

[0037] The following introduces the spatial discretization method based on isoparametric units.

[0038] Based on the shape functions of the airfoil triangle and quadrilateral isoparametric units, the Jacobian between the physical unit and the isoparametric unit of the airfoil can be obtained: ; Flux f and source term s Taking the example, the surface integral and volume integral forms on the isoparametric unit are given.

[0039] Surface integral over isoparametric elements 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: ; 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: ; Combining the relationship between the tangent vector and the normal vector, we can know the external normal on the curved edge: ; Therefore, the flux integral on the curved edge can be expressed as: ; At this point, the integral on the physical unit is transformed into the isoparametric unit, where and They refer to the two endpoints of the corresponding edges of the isoparametric unit ( ). Gaussian numerical integration scheme for combined conservation law equations: ; in, is the conserved variable in the control volume The average value on and are the convection and viscous fluxes along the normal direction out of the element surface, respectively. According to the discretization scheme of the high-order finite volume format, in the control volume element i superior, The semi-discrete form of can be expressed as: ; Represents the control volume unit i The number of unit surfaces, and denote the number of surface integral and volume integral points respectively, and and are the surface integral and volume integral coefficients. Then the integral of the flux on the curved edge can be finally expressed as: ; In the formula and is the Jacobian at the Gaussian integration point, which can be obtained according to Calculated.

[0040] Volume integral over isoparametric cells: For volume integral, there is the following relationship between isoparametric units and physical units: ; Therefore, the volume integral of the source function s(x,y) can be expressed as: ; Similarly, we can get: ; The following describes an optimization isoparametric unit according to an embodiment of the present invention.

[0041] 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.

[0042] Figure 5 In the example, 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 .

[0043] At this time, the reference point No. 10 is located at the trisection point near the lower edge of the vertical midline. Based on this reference point layout, the shape number obtained is: ; On this basis, the order of the quadrilateral isoparametric element can be further improved, that is, on the basis of the quartic shape function, 2 quintic terms are added to construct a more accurate 17-node quintic quadrilateral isoparametric element. The form of the shape function adopted is: ; The layout of the parameter points is as Figure 7 shown. We find that this reference point layout has better symmetry, and the No. 17 point located inside the quadrilateral is exactly the geometric center of the quadrilateral: ; In the embodiments of the present invention, the airfoil triangular element adopts a standard 10-node cubic isoparametric element; the airfoil quadrilateral element will respectively adopt a standard 10-node cubic and a 17-node quintic isoparametric element.

[0044] The high-order discretization method for optimizing the curved boundary of the isoparametric element constructs the isoparametric elements of the airfoil triangular element and the quadrilateral element by using the method of optimizing the isoparametric element, and then uses the spatial discretization method for area integration and volume integration. Its flow chart is as Figure 8 shown.

[0045] 1. First, read the information of the airfoil curved boundary after high-order representation; 2. Judge the shape of the curved element. If it is a triangular mesh, go to 3; if it is a quadrilateral mesh, go to 4; 3. Construct the shape function of the ten-node cubic element of the triangular mesh, and go to 7; 4. Judge the order of the isoparametric element adopted by the airfoil quadrilateral curved element. If it is the third order, go to 5; if it is the fifth order, go to 6; 5. Construct the shape function of the ten-node cubic element of the quadrilateral mesh, and go to 7; 6. Construct the shape function of the seventeen-node quintic element of the quadrilateral mesh, and go to 7; 7. Construct the Jacobian of the airfoil physical element and the isoparametric element; 8. Solve the area integral or volume integral of the isoparametric element; 9. End.

[0046] Adopt the airfoil high-order curved boundary expression method to test the calculation results of the curved element discretization method based on the quintic polynomial curve.

[0047] Inviscid cylinder flow around example: This example is a typical isentropic external flow problem. In this part, we respectively consider the incoming flow Mach number to be Calculation results at 0.38. The far field of the computational domain is taken as 30 times the cylinder diameter, and the cylinder diameter D = 1. In this part, 6 sets of triangular and quadrilateral meshes with different densities are used to verify different methods.

[0048] The incoming flow Mach number is Ma = 0.1. The far field of the computational domain is taken as 30 times the cylinder diameter. In this paper, the cylinder diameter D = 1. In this part, 7 sets of triangular and quadrilateral meshes with different densities are used 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 coarsest triangular and quadrilateral meshes are shown in

[0049] Figure 10 shows the entropy error situations of the curved element discretization method and the straight - edge element discretization method in the third - order format using the embodiments of the present invention. The third - order format - cubic isoparametric element and the third - order format - quintic isoparametric element are respectively the curved element discretization methods of the third - order cubic and quintic isoparametric elements of the embodiments of the present invention. The fourth - order format - cubic isoparametric element and the fourth - order format - quintic isoparametric element are respectively the curved element discretization methods of the fourth - order cubic and quintic isoparametric elements of the embodiments of the present invention. The third - order format - simplified curved - edge element and the fourth - order format - simplified curved - edge element are respectively the curved element discretization methods of the third - order and fourth - order. The third - order format - straight - edge element and the fourth - order format - straight - edge element are respectively the straight - edge element discretization methods. It can be seen that for the third - order format, the results obtained by using the two curved element discretization methods on all elements and wall elements are very close, and no obvious loss of accuracy occurs as the mesh is continuously refined. Even when discretization is completely based on straight - edge elements, the obtained results can still remain near the third - order accuracy.

[0050] When the accuracy is improved to the fourth - order, since the error of the fourth - order format decreases faster as the mesh is refined, at this time, using the embodiments of the present invention and the two curved element discretization schemes, the calculation results always maintain an accuracy above the fourth - order on all elements and wall elements. However, combining the specific data in Table 1 and Table 2, it can be seen that when using straight - edge element discretization, accuracy loss occurs as the mesh is refined. Especially between the last two sets of dense meshes, the solution accuracy on the wall elements will be reduced to 2.5 - order, which is much lower than the embodiments of the present invention and the fourth - order accuracy. Therefore, it is proved that for high - precision numerical formats, the necessity of using the curved element discretization scheme.

[0051] Table 1 Discretization errors of different discretization methods on the densest quadrilateral mesh

[0052] Table 2 Computational accuracy between the last two sets of dense meshes for different discretization methods

[0053] On this basis, Figure 11 Figure 11 further shows the entropy error situations obtained by different methods on a triangular grid. From the results, for the third-order scheme, the results obtained by using the two curved element discretization methods and the straight-edge element discretization method can always maintain around the third-order accuracy in all elements, and on the wall elements, as the grid is refined, the calculation accuracies of the three methods gradually approach the third order.

[0054] However, in the fourth-order scheme, the results obtained by the two curved element discretization methods are close. Combining with the specific data given in Table 3, it can be seen that the entropy error accuracy on all elements exceeds the fourth order, and on the wall elements, the calculation accuracy between the two densest grids exceeds 3.7 orders; but based on the straight-edge element discretization, as the grid is refined, the calculation results of the entropy error on the wall elements are severely lost, and it is only 1.065 orders between the two densest grids. Therefore, the correctness of the embodiment of the present invention is preliminarily proved.

[0055] Table 3 Discretization errors and calculation accuracies of different methods on the densest grid of triangular grid

[0056] Pressure coefficient distribution: For the inviscid flow around a circular cylinder discussed in this part, the pressure coefficient has an analytical form. The exact solution and numerical solution of the pressure coefficient distribution on the surface of the circular cylinder are analyzed below.

[0057] Since the incoming flow Mach number of the numerical example used in the test , it can be regarded as an incompressible flow. Considering the definition of the pressure coefficient as: ; Through the Bernoulli equation, the pressure distribution on the surface of the circular cylinder can be obtained as: ; Substitute into the Bernoulli equation and simplify to obtain: ; Substitute the above formula into the pressure coefficient definition given in the definition of the pressure coefficient to obtain: ; Under the conditions of the third-order and fourth-order schemes, the results of the pressure coefficient distribution on the surface of the circular cylinder obtained by three different discretization methods on a sparse grid of 30×120 are as Figure 12 shown, where the exact solution of the pressure coefficient is the exact solution of the pressure coefficient calculated by the formula .

[0058] From Figure 12As can be seen from the results, on the sparse grid, the results obtained by the three discretization schemes in the third-order format are close, with a certain 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 almost completely coincide with the exact curve. But from Figure 12 the further enlarged view shown in the right figure, we further find that compared with the embodiment of the present invention and the two curved element discretization methods, there are obvious errors in the pressure coefficient on the cylindrical surface obtained by completely discretizing based on straight-edge elements in the fourth-order format. This is consistent with the law presented in the accuracy test.

[0059] Viscous circular cylinder flow example: This part is about the subsonic viscous circular cylinder flow problem, and the given incoming flow Mach number and Reynolds number are respectively , .

[0060] Considering that for the finite volume format, the third-order accuracy format can only achieve second-order accuracy in calculating viscous problems, only the fourth-order format is considered in this part. The calculation results on quadrilateral and triangular meshes will be given in two parts in turn below.

[0061] Quadrilateral mesh: Figure 13 shows the variation curve of the pressure coefficient C D with the number of iteration steps. Among them, the cubic isoparametric element and the quintic isoparametric element are the cubic and quintic optimized isoparametric element curved element discretization methods of the embodiment of the present invention, the simplified curved edge element is the curved element discretization method, the straight edge element refers to the straight edge element discretization method, and the surface integral / volume integral result is the standard Gaussian surface integral and volume integral result based on the cylindrical analytical expression.

[0062] From the drag coefficient curve, it is further found that the curved element discretization of the optimized isoparametric element is almost the same as the result of the standard discretization method, and there are obvious deviations in the other two methods, especially the straight edge element discretization method based on the normal direction of the curved boundary. Preliminary analysis shows that for the isoparametric element method and the standard discretization method, when statistically analyzing the pressure coefficient curve and the drag coefficient integral curve on the wall boundary, the sampling points are located at the integral points on the reconstructed polynomial curve and the analytical cylindrical curve; while the sampling points of the other two methods are all located at the integral points of the straight edge element.

[0063] Represent the two integral points with C G and S G respectively. Combining Figure 14 shows that the pressure coefficient can be expressed as: ; The drag coefficient is ; In the formula and are the lengths of the straight and curved boundaries of the unit, respectively. Therefore, even though the straight-edge element discretization scheme based on the normal of the curved boundary can ensure the flow results with corresponding accuracy, the calculation of the pressure and drag coefficients involves the actual integration point positions and the integration along the boundary surface. At this time, both the integration points and the numerical integration process are along the straight boundary, resulting in obvious differences between the calculation results and the discretization method based on curved elements.

[0064] Triangular mesh: Figure 15 From the curve of the change of the integral of the drag coefficient on the cylinder surface shown in

[0065] with the number of iteration steps, it can be seen that the pressure coefficient results obtained by the optimized isoparametric element method of the embodiment of the present invention are close to those of the standard discretization method based on the analytical curve of the cylinder, while there are certain deviations between the results obtained by the two straight-edge element-based discretization methods and the curved element discretization method.

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

[0067] 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 on the memory to implement the steps of the method of Embodiment 1 above.

[0068] In some implementations, the memory can be a high-speed random access memory (RAM: Random Access Memory), and may also include a non-volatile memory, such as at least one disk memory.

[0069] In other implementations, the processor can be a general-purpose processor of various types such as a central processing unit (CPU), a digital signal processor (DSP), etc., which is not limited here.

[0070] Although the preferred embodiments of the present application have been described, additional changes and modifications can be made to these embodiments by those skilled in the art once they learn the basic creative concept. Therefore, the appended claims are intended to be construed to include the preferred embodiments as well as all changes and modifications that fall within the scope of the present application.

[0071] Obviously, those skilled in the art can make various changes and modifications to the present application without departing from the spirit and scope of the present application. Thus, if these modifications and variations of the present application fall within the scope of the claims of the present application and their equivalent technologies, the present application is also intended to include these modifications and variations.

Claims

1. A wing boundary simulation method for optimizing isoparametric elements, characterized in that, It includes the following steps: S1. Obtain the wing curved boundary information; S2. Judge the shape of the curved element according to the curved boundary information. If it is a triangular mesh, construct a ten-node cubic element shape function for the triangular mesh. If it is a quadrilateral mesh, judge the isoparametric element order adopted by the airfoil quadrilateral curved element. If the order is three, construct a ten-node cubic element shape function for the quadrilateral mesh. If the order is five, construct a seventeen-node quintic element shape function for the quadrilateral mesh; S3. Construct the Jacobian matrix of the wing physical element and the isoparametric element; S4. Use the Jacobian matrix to calculate the surface integral and / or volume integral on the isoparametric element.

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

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

4. The wing boundary simulation method for optimizing isoparametric elements according to claim 1, characterized in that In step S2, the expression of the seventeen-node quintic element shape function of the quadrilateral mesh is: ； Among them, are the reference point coordinates of the standard isoparametric element.

5. The wing boundary simulation method for optimizing isoparametric elements 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 element, is the reference point coordinate of the actual physical element.

6. The wing boundary simulation method for optimizing isoparametric elements according to claim 1, characterized in that The surface integral expression is: ； Among them, and respectively refer to the coordinates of the two endpoints of the corresponding side of the isoparametric element in the direction in the isoparametric space, , represents the number of surface integrals, is the surface integral coefficient, f is the convective or viscous flux in the external space of the wing, is the flux at the k-th Gauss integration point, and is the Jacobian at the Gauss integration point.

7. The method for simulating the wing boundary of an optimized isoparametric element according to claim 1, characterized in that, The volume integral expression is: ； Among them, represents the number of volume integration points, is the volume integration coefficient, is the source term function, , is the coordinate of the corresponding reference point of the physical element, and n is the number of reference points of the isoparametric element, is the shape function, is the coordinate of the q-th Gauss integration point in the isoparametric space, is the Jacobian matrix at the q-th Gauss integration point.

8. A terminal device, comprising a memory, a processor, and a computer program stored on 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 above.

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

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

Images