Numerical simulation method for streaming of rotationally symmetrical geometry

Abstract

The invention discloses a rotational symmetry geometry streaming numerical simulation method, and relates to the field of numerical simulation, and the method comprises the steps: generating a body-fitted structure grid containing an odd axis; virtual grid points are arranged on the outer side of the odd axis; giving space coordinates and flow field variables on the virtual grid points; calculating to obtain grid transformation derivatives and Jacobian determinants on the boundary and the internal grid points of the odd axis; a grid transformation derivative and a Jacobian determinant on the virtual grid point are given; giving a first derivative value of a flow field variable on a virtual grid point for the viscous streaming; discretizing a space convection term and a viscosity term of the control equation; updating the flow field variable of the flow around the rotationally symmetrical geometry; calculating to obtain a flow field variable on the odd axis; obtaining flow field data corresponding to the input data based on all calculation results; according to the method, the problem of adaptability of a high-precision flux reconstruction finite difference method and odd axis boundary processing is solved, numerical errors are reduced, and the flow field quality at the odd axis boundary during numerical simulation of the rotational symmetry geometry streaming problem is improved.

Description

Technical Field

[0001] The present invention relates to the field of numerical simulation, and in particular to a numerical simulation method for flow around a rotationally symmetric geometric body. Background Art

[0002] Singular axes are special boundaries introduced when dealing with rotationally symmetric geometric problems in computational fluid dynamics (CFD). In order to numerically simulate the flow around geometric bodies with rotational symmetry (such as the sharp spherical head or cylinder of a high-speed aircraft), it is usually chosen to generate an axisymmetric body-fitting structure grid around the axis of symmetry. In an axisymmetric body-fitting structure grid, the grid lines of the surface grid converge into a point at the axis of symmetry, which becomes a singular point. After forming a volume grid in space, these points constitute the singular axis. The singular axis appears as a straight line in physical space, but is mapped to a plane in computational space. Although the grid points on the singular axis theoretically rotate around the axis, the x, y, and z coordinate values ​​do not change. The corresponding geometric meaning is that the area of ​​the singular axis in physical space is zero. When performing coordinate transformation from the physical domain (x, y, and z) to the computational domain (ξ, η, ζ), the inverse of the determinant of the Jacobian matrix degenerates to zero at the singular axis (J -1 = 0). The mathematical singularity introduced by the coordinate system transformation itself (the grid derivative tends to zero or infinity) conflicts with the intrinsic smoothness requirement of the physical quantity in discrete numerical calculations, resulting in the singular axis problem, which requires special treatment in CFD calculations.

[0003] In the finite volume method, singular axes appear as degenerate boundary surfaces, resulting in zero flux integrals and eliminating the need for boundary values. However, in the finite difference method, grid points on the singular axes must be assigned values. If a grid contains a singular axis, traditional finite difference calculations typically employ a basic approach of simply setting the numerical flux of the grid points on the singular axis to zero. The flow field variable values ​​at the adjacent grid nodes on the inner side of the singular axis are then interpolated, summed circumferentially, and averaged to obtain the flow field variable values ​​at the grid points on the singular axis. While this approach seems reasonable when there's no angle of attack (AoA), when the mesh near the singular axis is sparse, the resulting calculation results often have large errors, with errors concentrating on the singular axis and fluctuations in the calculation. When there is an AoA, because the area of ​​the singular axis in physical space is zero, even though the actual flow passes through it, the zero flux singular axis treatment prevents variables on the axis from directly participating in the solution of the flow equations. This lacks the direct influence of the flow field on the singular axis, which can result in structural distortion or instability in the flow field near the singular axis. Therefore, traditional singular axis boundary treatment methods for finite-difference methods face challenges with numerical stability and low accuracy.

[0004] Recently, researchers have proposed using a hybrid algorithm based on traditional methods to solve the problem of singular axis calculation using the finite difference method. Specifically, a finite difference method based on reconstruction of conserved variables is used near the singular axis, and a finite difference method based on flux reconstruction is used in areas away from the singular axis. Reconstruction based on conserved variables avoids the problem of insufficient information transfer caused by the zero area of ​​the singular axis in flux calculation, and reflects the mutual influence between the variables on the singular axis and the surrounding flow field. Because the region needs to be divided for conservation variable reconstruction and flux reconstruction respectively, the hybrid algorithm is not only complex to implement, but also the region division is somewhat empirical. In addition, there is also the problem of smooth transition between different reconstruction methods. Summary of the Invention

[0005] In view of the shortcomings of the existing finite difference processing method for singular axis boundaries in the numerical simulation of flow around rotationally symmetric geometric bodies, the present invention proposes a finite difference processing method for singular axis boundaries based on virtual grid points and adopts unified high-precision flux reconstruction in the entire computational domain. It solves the compatibility problem between the high-precision flux reconstruction finite difference method and the singular axis boundary processing, reduces numerical errors, and improves the flow field quality at the singular axis boundary when numerically simulating the flow around rotationally symmetric geometric bodies.

[0006] To achieve the above-mentioned object, the present invention provides a method for numerical simulation of flow around a rotationally symmetric geometric body, the method comprising:

[0007] Step 1: Generate a body-fitting mesh containing singular axes required for numerical simulation of flow around a rotationally symmetric geometry using mesh generation software;

[0008] Step 2: Set virtual grid points outside the singular axis of the body-fitted structure grid;

[0009] Step 3: Based on the spatial coordinate values ​​of the internal grid points of the body-fitted structure grid, the spatial coordinate values ​​of the virtual grid points are given;

[0010] Step 4: Based on the spatial coordinate values ​​of the virtual grid points, the singular axis boundary grid points, and the internal grid points, the grid transformation derivatives and Jacobian determinants at the singular axis boundary grid points and the internal grid points are calculated and used to discretize the spatial convection derivative terms and spatial viscous derivative terms of the governing equations.

[0011] Step 5: Based on the grid transformation derivatives and Jacobian determinant at the internal grid points of the body-fitting structure grid, the grid transformation derivatives and Jacobian determinant at the given virtual grid points are used to discretize the spatial convection derivative terms and spatial viscous derivative terms of the governing equations;

[0012] Step 6: Based on the flow field variable values ​​at the internal grid points of the body-fitting structure grid, the flow field variable values ​​at the virtual grid points are given;

[0013] Step 7: For the viscous flow problem around a rotationally symmetric body, the first-order derivatives of the flow field variables at the boundary grid points of the singular axis and the internal grid points are calculated based on the flow field variable values ​​at the boundary grid points of the singular axis and the internal grid points. The first-order derivatives of the flow field variables at the virtual grid points are given and used to discretize the spatial viscous derivative terms of the governing equations.

[0014] Step 8: Discretize the spatial convection derivative and spatial viscous derivative of the governing equation based on the grid transformation derivatives, Jacobian determinant, flow field variable values, and first-order derivative values ​​of flow field variables at the virtual grid points, singular axis boundary grid points, and internal grid points;

[0015] Step 9: Based on the spatial convection derivatives and spatial viscous derivatives of the discretized governing equations, update the flow field variables of the flow around the rotationally symmetric geometry through time-marching iterations.

[0016] Step 10: Based on the updated flow field variable values ​​of the flow around the rotationally symmetric geometric body, the flow field variable values ​​on the singular axis of the rotationally symmetric geometric body are calculated;

[0017] Step 11: Determine whether the total number of calculation time steps or the calculation convergence condition has been reached. If not, return to step 6. If so, end the simulation and obtain the flow field data of the rotationally symmetric geometry corresponding to the input data based on all calculation results.

[0018] The principle of this method is as follows: by setting virtual grid points outside the singular axis boundary of the body-fitting structure grid of the rotationally symmetric geometry, and reasonably setting the grid coordinates, grid transformation derivatives, flow field variables and their first-order derivatives on the virtual grid points, a unified high-precision flux reconstruction difference method is used to discretize the spatial derivative terms of the control equation in the entire computational domain including the singular axis, thereby solving the compatibility problem between the high-precision flux reconstruction finite difference method and the singular axis boundary processing of the body-fitting structure grid.

[0019] Preferably, the rotationally symmetrical geometric body is a spherical head or a cylindrical body of an aircraft or spacecraft.

[0020] Preferably, the step 1 specifically includes:

[0021] Based on the mesh generation software, a body-fitting structure grid with singular axes is generated for the numerical simulation of the flow around the rotationally symmetric geometric body. The coordinates of the body-fitting structure grid in the physical space are (x, y, z), and the singular axes are assumed to fall on the x-axis. The coordinates of the body-fitting structure grid in the computational space are (ξ, η, ζ), and the ξ direction corresponds to the normal grid line and is represented by the index i, i=I maxrepresents the object surface, the η direction corresponds to the flow grid line and is represented by the index j, j = 1 represents the singular axis, the ζ direction corresponds to the circumferential grid line and is represented by the index k; the number of grid points on the normal grid line is I max ; The number of grid points on the flow grid line is J max ; The number of grid points on the circumferential grid line is K max , K max is an odd value, and the grid points on the circumferential grid line are evenly distributed.

[0022] Preferably, the step 2 specifically includes:

[0023] The grid points on the singular axis are expanded circumferentially in the computational space to form a plane j = 1. For each grid point on the plane, several layers of virtual grid points are set in the negative direction of the index j (i.e., the direction in which the index j decreases).

[0024] Preferably, the spatial coordinate values ​​and flow field variable values ​​on the virtual grid points are obtained as follows:

[0025] The spatial coordinate values ​​and flow field variable values ​​on the virtual grid points are obtained as follows:

[0026] The spatial coordinate values ​​and flow field variable values ​​at the current virtual grid point are set based on the spatial coordinate values ​​and flow field variable values ​​at the internal grid points of the body-fitting structure grid. The specific corresponding relationships are as follows:

[0027] ;

[0028] Among them, the flow field variables 、 、 、 and are density, velocity component in the x direction, velocity component in the y direction, velocity component in the z direction, and pressure respectively; jg is the index of the virtual grid point in the η direction, jt and kt are the indexes of the internal grid points corresponding to the virtual grid point in the η direction and the ζ direction respectively, is the number of virtual grid points.

[0029] Preferably, the method for calculating the grid transformation derivative from the physical space to the computational space in step 4 is:

[0030] ;

[0031] Where J is the Jacobian determinant of the coordinate transformation, 、 、 、 、 、 、 、 and They are the calculation space coordinates ( , , ) is the derivative of the mesh transformation in the physical space (x, y, z) direction, is the partial derivative operator, 、 、 、 、 、 、 、 and They are the physical space coordinates (x, y, z) in the computational space ( , , ) direction of the grid transformation derivative;

[0032] Among them, at the boundary grid points of the singular axis, let:

[0033] ;

[0034] in, , , , , , .

[0035] Preferably, the step 5 is specifically to determine the grid transformation derivative and Jacobian determinant at the virtual grid point according to the following formula:

[0036] ;

[0037] in, 、 、 、 、 、 、 、 and are the grid transformation derivatives 、 、 、 、 、 、 、 and The value at the virtual grid point (i, jg, k), 、 、 、 、 、 、 、 and are the grid transformation derivatives 、 、 、 、 、 、 、 and The value at the internal grid point (i, jt, k), is the value of the Jacobian determinant at the virtual grid point (i, jg, k), is the value of the Jacobian determinant at the interior grid point (i, jt, k).

[0038] Preferably, in step 7, if the flow in the viscous flow problem around a rotationally symmetric geometric body is an axisymmetric flow, the first-order derivative value of the flow field variable at the virtual grid point is obtained according to the following formula:

[0039] ;

[0040] in, and are the values ​​of the first-order derivative of the velocity component u in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component u in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component u in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component v in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component v in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component v in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component w in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component w in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component w in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of temperature t in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first derivative of temperature t in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of temperature t in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), respectively;

[0041] In step 7, if the flow in the viscous flow problem around the rotationally symmetric geometric body is non-axisymmetric, the first-order derivative value of the flow field variable at the virtual grid point is calculated by the unilateral difference method based on the flow field variable values ​​at the virtual grid point and the internal grid point.

[0042] Preferably, step 8 includes: based on the grid transformation derivatives, Jacobian determinants, flow field variable values ​​and first-order derivative values ​​of flow field variables on the virtual grid points, the singular axis boundary grid points and the internal grid points, at the singular axis boundary and its adjacent areas, a flux reconstruction difference method consistent with the method used for the internal grid points is used to discretize the spatial convection derivative terms of the control equation, and a central difference method consistent with the method used for the internal grid points is used to discretize the spatial viscosity derivative terms of the control equation to obtain the discretized spatial derivative terms of the control equation.

[0043] Preferably, the step 10 specifically includes:

[0044] Based on the updated flow field variable values ​​of the flow around the rotationally symmetric geometric body, the flow field variable values ​​on the singular axis are interpolated along the j direction for the index position (i, k); the interpolated flow field variable values ​​are replaced by the average of the circumferential numbers as the flow field variable values ​​on the singular axis of the rotationally symmetric geometric body.

[0045] One or more technical solutions provided by the present invention have at least the following technical effects or advantages:

[0046] The present invention introduces virtual grid points outside the singular axis boundary of the body-fitting structure grid of a rotationally symmetric geometric body, reasonably sets the grid coordinates, grid transformation derivatives, flow field variables and their first-order derivatives on the virtual grid points, and discretizes the spatial terms at the singular axis boundary and its adjacent areas using the same flux reconstruction method as the internal grid points. The influence and contribution of the flow field variables on the other side of the singular axis to the local flow field are taken into account, and the singular axis boundary processing of the body-fitting structure grid based on the high-precision flux reconstruction finite difference method is realized, thereby reducing numerical errors and improving the flow field quality at the singular axis boundary when numerically simulating the flow around a rotationally symmetric geometric body. BRIEF DESCRIPTION OF THE DRAWINGS

[0047] The drawings described herein are used to provide a further understanding of the embodiments of the present invention, constitute a part of the present invention, and do not constitute a limitation of the embodiments of the present invention;

[0048] Figure 1 This is a flow chart of a numerical simulation method for flow around a rotationally symmetric geometric body;

[0049] Figure 2 Schematic diagram of coordinate axes and indices of the supersonic spherical head with singular axis body-fitting structure grid;

[0050] Figure 3 Schematic diagram of the wall pressure distribution in the calculation results of the inviscid flow around a supersonic spherical head in implementation case 1;

[0051] Figure 4 Schematic diagram of the pressure distribution on the symmetric surface in the calculation results of the inviscid flow around a supersonic spherical head in implementation case 1;

[0052] Figure 5 Schematic diagram of the wall pressure distribution in the calculation results of the supersonic viscous flow around a spherical head in implementation case 2;

[0053] Figure 6 Schematic diagram of the pressure distribution on the symmetric surface in the calculation results of the supersonic viscous flow around a spherical head in implementation case 2. DETAILED DESCRIPTION

[0054] In order to more clearly understand the above-mentioned objects, features and advantages of the present invention, the present invention is further described in detail below in conjunction with the accompanying drawings and specific embodiments. It should be noted that the embodiments of the present invention and the features therein can be combined with each other without conflict.

[0055] In the following description, many specific details are set forth to facilitate a full understanding of the present invention. However, the present invention may also be implemented in other ways different from those described herein. Therefore, the scope of protection of the present invention is not limited to the specific embodiments disclosed below.

[0056] Embodiment 1;

[0057] Please refer to Figure 1 , Figure 1 Schematic diagram of a flow chart of a method for numerically simulating flow around a rotationally symmetric geometric body. Embodiment 1 of the present invention provides a method for numerically simulating flow around a rotationally symmetric geometric body, the method comprising:

[0058] Step 1: Generate a body-fitting mesh containing singular axes required for numerical simulation of flow around a rotationally symmetric geometry using mesh generation software;

[0059] Step 2: Set virtual grid points outside the singular axis of the body-fitted structure grid;

[0060] Step 3: Based on the spatial coordinate values ​​of the internal grid points of the body-fitted structure grid, the spatial coordinate values ​​of the virtual grid points are given;

[0061] Step 4: Based on the spatial coordinate values ​​of the virtual grid points, the singular axis boundary grid points, and the internal grid points, the grid transformation derivatives and Jacobian determinants at the singular axis boundary grid points and the internal grid points are calculated and used to discretize the spatial convection derivative terms and spatial viscous derivative terms of the governing equations.

[0062] Step 5: Based on the grid transformation derivatives and Jacobian determinant at the internal grid points of the body-fitting structure grid, the grid transformation derivatives and Jacobian determinant at the given virtual grid points are used to discretize the spatial convection derivative terms and spatial viscous derivative terms of the governing equations;

[0063] Step 6: Based on the flow field variable values ​​at the internal grid points of the body-fitting structure grid, the flow field variable values ​​at the virtual grid points are given;

[0064] Step 7: For the viscous flow problem around a rotationally symmetric body, the first-order derivatives of the flow field variables at the boundary grid points of the singular axis and the internal grid points are calculated based on the flow field variable values ​​at the boundary grid points of the singular axis and the internal grid points. The first-order derivatives of the flow field variables at the virtual grid points are given and used to discretize the spatial viscous derivative terms of the governing equations.

[0065] Step 8: Discretize the spatial convection derivative and spatial viscous derivative of the governing equation based on the grid transformation derivatives, Jacobian determinant, flow field variable values, and first-order derivative values ​​of flow field variables at the virtual grid points, singular axis boundary grid points, and internal grid points;

[0066] Step 9: Based on the spatial convection derivatives and spatial viscous derivatives of the discretized governing equations, update the flow field variables of the flow around the rotationally symmetric geometry through time-marching iterations.

[0067] Step 10: Based on the updated flow field variable values ​​of the flow around the rotationally symmetric geometric body, the flow field variable values ​​on the singular axis of the rotationally symmetric geometric body are calculated;

[0068] Step 11: Determine whether the total number of calculation time steps or the calculation convergence condition has been reached. If not, return to step 6. If so, end the simulation and obtain the flow field data of the rotationally symmetric geometry corresponding to the input data based on all calculation results.

[0069] The following is a detailed description of this method:

[0070] This paper proposes a singular axis boundary finite difference method for numerical simulation of flow around rotationally symmetric geometric bodies based on virtual grid points and unified high-precision flux reconstruction in the entire computational domain. The method is implemented in the following 11 steps:

[0071] Step 1: Based on mesh generation software such as POINTWISE, generate the body-fitting mesh including singular axes required for numerical simulation of flow around rotationally symmetric geometric bodies. Specific requirements are as follows: Figure 2 As shown in the example of the spherical head flow grid, the coordinates in the physical space are (x, y, z), and the singular axis is assumed to fall on the x-axis; the coordinates in the computational space are (ξ, η, ζ), and ξ is assumed to correspond to the normal grid line and is represented by the index i, where i=I max represents the object surface, η corresponds to the flow grid line and is represented by the index j, where j = 1 represents the singular axis, and the ζ direction corresponds to the circumferential grid line and is represented by the index k. The number of grid points on the normal grid line is I max , i.e. i=1, …, I max ; The number of grid points on the flow grid line is J max , that is, j=1,…,J max ; The number of grid points on the circumferential grid line is K max , that is, k=1, …, K max Here, the grid points on the circumferential grid line are required to be evenly distributed, and the number of circumferential grid points K is required to be max Expected to be an odd value.

[0072] Step 2: Set virtual grid points outside the singularity axis of the body-fitting mesh of the rotationally symmetric geometry. The grid points on the singularity axis are expanded circumferentially within the computational space to form a plane with j = 1. For each grid point on the plane, set several layers of virtual grid points in the negative direction of the index j. The number of virtual grid point layers is determined by the accuracy of the finite difference method. For the fifth-order finite difference WENO method, the number of virtual grid point layers is set to Nghost = 3.

[0073] Step 3: Given the spatial coordinate values ​​of the virtual grid points outside the singular axis of the rotationally symmetric geometric body, which are used to calculate the grid transformation derivatives, the spatial coordinate values ​​of the current virtual grid points are set based on the grid coordinate values ​​of the grid points inside the body-fitting structure grid. The specific corresponding relationship is as follows:

[0074] ; (1)

[0075] Where jg is the index of the virtual grid point in the η direction, jt and kt are the indexes of the internal grid points corresponding to the virtual grid point in the η direction and the ζ direction, respectively. is the number of virtual grid points.

[0076] Step 4: Calculate the grid transformation derivatives and Jacobian at the grid points along the singularity axis and the interior of the rotationally symmetric geometry. These are used to discretize the spatial convection and viscous derivatives of the governing equations. Using a consistent high-precision central difference method (such as fourth- or sixth-order central difference) based on the virtual grid points and the spatial coordinates of the grid points along the singularity axis and the interior, the grid transformation derivatives and Jacobian at the grid points along the singularity axis and the interior are calculated.

[0077] The derivative of the grid transformation from physical space to computational space is calculated as follows:

[0078] ; (2)

[0079] Where J is the Jacobian determinant of the coordinate transformation, 、 、 、 、 、 、 、 and is the corresponding mesh transformation derivative.

[0080] It should be noted that at the grid points on the boundary of the singular axis, the grid transformation derivative has the following expression: , , , , , and In actual numerical simulations, in order to avoid the situation of calculating the inverse of zero, at the boundary grid points of the singular axis, let:

[0081] . (3)

[0082] Step 5: Given the grid transformation derivatives and Jacobian determinants at the virtual grid points outside the singular axis of the rotationally symmetric geometry, the spatial convection derivatives and viscous derivatives of the governing equations are discretized. Based on the grid transformation derivatives and Jacobian determinants at the internal grid points, the grid transformation derivatives and Jacobian determinants at the virtual grid points are determined according to the following correspondence:

[0083] ; (4)

[0084] Step 6: Given the flow field variable value at the virtual grid point outside the singular axis of the rotationally symmetric geometric body, the flow field variable value at the current virtual grid point is set based on the flow field variable value at the internal grid point. The specific corresponding relationship is as follows:

[0085] ; (5)

[0086] Among them, the flow field variables 、 、 、 and They are density, x-direction velocity component, y-direction velocity component, z-direction velocity component and pressure.

[0087] Step 7: For viscous flow problems around rotationally symmetric geometries, the first-order derivatives of the flow field variables at the singular axis boundary grid points, internal grid points, and virtual grid points are calculated based on the flow field variable values. The first-order derivatives of the flow field variables at the singular axis boundary grid points and internal grid points are also given, and the first-order derivatives of the flow field variables at the virtual grid points are used to discretize the spatial viscous derivative terms of the governing equations. For viscous flow problems around geometries, the first-order derivatives of the flow field variables are required to construct the viscous flux in the viscous terms of the governing equations. If the flow is axisymmetric, the first-order derivatives of the flow field variables at the virtual grid points are given according to the following relationship:

[0088] ; (6)

[0089] If the flow is non-axisymmetric, the above relationship does not hold. The first-order derivatives of the flow field variables at the virtual grid points can be calculated using high-precision single-sided difference methods (such as fourth-order or sixth-order single-sided difference methods) based on the flow field variable values ​​at the virtual grid points and internal grid points.

[0090] Step 8: Discretize the spatial convection derivative term and the viscous derivative term of the control equation based on the grid transformation derivatives, the Jacobian determinant, the flow field variable values and the first-order derivative values of the flow field variables on the virtual grid points, the singular axis boundary grid points and the internal grid points. At the singular axis boundary and its adjacent area, the convection derivative term is discretized by using a high-precision flux reconstruction difference method (such as a fifth-order precision finite difference WENO method) consistent with the method used for the internal grid points, and the viscous derivative term is discretized by using a high-precision central difference method (such as a fourth-order or sixth-order central difference method) consistent with the method used for the internal grid points, to obtain the spatial derivative term of the discretized control equation. The fifth-order precision finite difference WENO method can refer to the literature, and the embodiments of the present application do not perform corresponding elaboration:

[0091] 1. G.S. Jiang, C.-W. Shu, Efficient implementation of weighted ENOschemes, J. Comput. Phys. 126 (1996) 202–228.

[0092] 2. R. Borges, M. Carmona, B. Costa, W.S. Don, An improved weightedessentially non-oscillatory scheme for hyperbolic conservation laws, J.Comput. Phys. 227 (2008) 3101–3211.

[0093] Step 9: Update the flow field variable values of the flow around the rotationally symmetric geometric body by time marching iteration based on the discretized spatial convection derivative term and the viscous derivative term of the control equation.

[0094] Step 10: Calculate the flow field variable values on the singular axis of the rotationally symmetric geometric body by interpolation and then averaging based on the updated flow field variable values.

[0095] Specifically, first, for a certain (i, k) index position, the flow field variable on the singular axis (j=1) is obtained by interpolation along the j direction, and different k positions may obtain different on-axis flow field variables; second, the obtained value is replaced by averaging along the circumferential direction (all points of k values are traversed) to obtain the on-axis flow field variable value. Taking the variable u as an example, there are two ways of interpolation: (1) zero gradient condition: For inviscid flow, a second-order difference is used for discretization, and for viscous flow, a third-order difference is used for discretization. (2) Equal rate of change, that is, the second-order derivative of the flow direction of the second point outside the singular axis is zero, and the specific implementation is u i =2ui+1 - u i+2 .

[0096] Step 11: judging whether the total time step or the convergence condition is reached, if not, returning to execute step 6, if yes, ending the simulation, and obtaining the flow field data of the rotational symmetric geometric body around the flow based on all the calculation results.

[0097] The application solves the adaptability problem of the high-precision flux reconstruction finite difference method and the singular axis boundary treatment of the body-fitted structure grid by setting virtual grid points outside the singular axis boundary of the rotational symmetric geometric body body-fitted structure grid, and reasonably setting the grid coordinates, grid transformation derivatives, flow field variables and their first-order derivatives on the virtual grid points, and using a unified high-precision flux reconstruction difference method to disperse the spatial derivative terms of the control equation in the whole calculation domain including the singular axis.

[0098] Example two;

[0099] On the basis of example one, example two is introduced and explained by a specific implementation case.

[0100] Implementation case 1 is supersonic inviscid sphere head flow around;

[0101] Supersonic inviscid sphere head flow around calculation condition: the Mach number of the incoming flow is M=3.0, and the attack angle is zero. The surface of the sphere head is a slip boundary, and the supersonic inlet and outlet conditions are adopted. The number of circumferential grid points is 31, which is uniformly distributed; the number of normal grid points is 101; and the number of flow direction grid points is 61. The control equation is the Euler equation under the general curvilinear coordinate system, the spatial dispersion of the convection term adopts the five-order precision finite difference WENO method based on flux reconstruction, and the time advancement adopts the three-order precision TVD Runge-Kutta method. Figure 3-Figure 4 is the calculation result, including the wall surface pressure distribution and the symmetric plane pressure distribution.

[0102] Implementation case 2 is supersonic viscous sphere head flow around;

[0103] Supersonic viscous sphere head flow around calculation condition: laminar flow, attack angle is zero, incoming flow Mach number is M=3.0, Reynolds number based on the sphere head radius is Re=1×10 5 . The surface of the sphere head adopts non-slip condition, and the surface pressure is obtained by the zero pressure gradient along the surface normal. Supersonic inlet and outlet conditions are adopted. The number of circumferential grid points is 31, which is uniformly distributed; the number of normal grid points is 101; and the number of flow direction grid points is 61. The control equation is the Navier-Stokes equation under the general curvilinear coordinate system, the spatial dispersion of the convection term adopts the five-order precision finite difference WENO method based on flux reconstruction, the spatial dispersion of the viscous term adopts the six-order precision central difference method, and the time advancement adopts the three-order precision TVD Runge-Kutta method. Figure 5-Figure 6 is the calculation result, including the wall pressure distribution and the symmetric surface pressure distribution.

[0104] Although the preferred embodiments of the present invention 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 invention.

[0105] Obviously, those skilled in the art may make various changes and modifications to the present invention without departing from the spirit and scope of the present invention. Thus, if such changes and modifications fall within the scope of the claims and their equivalents, the present invention is intended to include such changes and modifications.

Claims

1. A numerical simulation method for flow around a rotationally symmetric geometric body, characterized in that: The method comprises: Step 1: Generate a body-fitting mesh containing singular axes required for numerical simulation of flow around a rotationally symmetric geometry using mesh generation software; Step 2: Set virtual grid points outside the singular axis of the body-fitted structure grid; Step 3: Based on the spatial coordinate values ​​of the internal grid points of the body-fitted structure grid, the spatial coordinate values ​​of the virtual grid points are given; Step 4: Based on the spatial coordinate values ​​of the virtual grid points, the singular axis boundary grid points, and the internal grid points, the grid transformation derivatives and Jacobian determinants at the singular axis boundary grid points and the internal grid points are calculated and used to discretize the spatial convection derivative terms and spatial viscous derivative terms of the governing equations. Step 5: Based on the grid transformation derivatives and Jacobian determinant at the internal grid points of the body-fitting structure grid, the grid transformation derivatives and Jacobian determinant at the given virtual grid points are used to discretize the spatial convection derivative terms and spatial viscous derivative terms of the governing equations; Step 6: Based on the flow field variable values ​​at the internal grid points of the body-fitting structure grid, the flow field variable values ​​at the virtual grid points are given; Step 7: For the viscous flow problem around a rotationally symmetric body, the first-order derivatives of the flow field variables at the boundary grid points of the singular axis and the internal grid points are calculated based on the flow field variable values ​​at the boundary grid points of the singular axis and the internal grid points. The first-order derivatives of the flow field variables at the virtual grid points are given and used to discretize the spatial viscous derivative terms of the governing equations. Step 8: Discretize the spatial convection derivative and spatial viscous derivative of the governing equation based on the grid transformation derivatives, Jacobian determinant, flow field variable values, and first-order derivative values ​​of flow field variables at the virtual grid points, singular axis boundary grid points, and internal grid points; Step 9: Based on the spatial convection derivatives and spatial viscous derivatives of the discretized governing equations, update the flow field variables of the flow around the rotationally symmetric geometry through time-marching iterations. Step 10: Based on the updated flow field variable values ​​of the flow around the rotationally symmetric geometric body, the flow field variable values ​​on the singular axis of the rotationally symmetric geometric body are calculated; Step 11: Determine whether the total number of calculation time steps or the calculation convergence condition has been reached. If not, return to step 6. If so, end the simulation and obtain the flow field data of the rotationally symmetric geometry corresponding to the input data based on all calculation results.

2. The numerical simulation method for flow around a rotationally symmetric geometric body according to claim 1, characterized in that: The rotationally symmetric geometric body is the spherical head and cylindrical body of an aircraft or spacecraft.

3. The numerical simulation method for flow around a rotationally symmetric geometric body according to claim 1, characterized in that: The step 1 specifically includes: Based on the mesh generation software, a body-fitting structure grid with singular axes is generated for the numerical simulation of the flow around the rotationally symmetric geometric body. The coordinates of the body-fitting structure grid in the physical space are (x, y, z), and the singular axes are assumed to fall on the x-axis. The coordinates of the body-fitting structure grid in the computational space are (ξ, η, ζ), and the ξ direction corresponds to the normal grid line and is represented by the index i, i=I max represents the object surface, the η direction corresponds to the flow grid line and is represented by the index j, j = 1 represents the singular axis, the ζ direction corresponds to the circumferential grid line and is represented by the index k; the number of grid points on the normal grid line is I max ; The number of grid points on the flow grid line is J max ; The number of grid points on the circumferential grid line is K max , K max is an odd value, and the grid points on the circumferential grid line are evenly distributed.

4. The numerical simulation method for flow around a rotationally symmetric geometric body according to claim 3, characterized in that: The step 2 specifically includes: The grid points on the singular axis are expanded circumferentially in the computational space to form a plane with j=1. For each grid point on the plane, several layers of virtual grid points are set in the negative direction of the index j.

5. The numerical simulation method for flow around a rotationally symmetric geometric body according to claim 3, characterized in that: The spatial coordinate values ​​and flow field variable values ​​on the virtual grid points are obtained as follows: The spatial coordinate values ​​and flow field variable values ​​at the current virtual grid point are set based on the spatial coordinate values ​​and flow field variable values ​​at the internal grid points of the body-fitting structure grid. The specific corresponding relationships are as follows: ； Among them, the flow field variables and They are density, x-direction velocity component, y-direction velocity component, z-direction velocity component and pressure; jg is the index of the virtual grid point in the η direction, jt and kt are the indices of the internal grid points corresponding to the virtual grid points in the η direction and the ζ direction, is the number of virtual grid points.

6. The method for numerical simulation of flow around a rotationally symmetric geometric body according to claim 3, characterized in that: The method for calculating the grid transformation derivative from the physical space to the computational space in step 4 is: ； Where J is the Jacobian determinant of the coordinate transformation, Calculate the spatial coordinates The derivatives of the mesh transformation in the physical space (x, y, z) directions, is the partial derivative operator, They are the physical space coordinates (x, y, z) in the computational space ( ) direction of the grid transformation derivative; Among them, at the boundary grid points of the singular axis, let: ； in, , .

7. The method for numerical simulation of flow around a rotationally symmetric geometric body according to claim 5, characterized in that: The step 5 specifically determines the grid transformation derivative and Jacobian determinant at the virtual grid point according to the following formula: ； in, 、 、 、 、 、 、 、 and are the grid transformation derivatives and The value at the virtual grid point (i, jg, k), 、 、 、 、 、 、 、 and are the grid transformation derivatives and The value at the internal grid point (i, jt, k), is the Jacobian determinant at the virtual grid point The value on is the value of the Jacobian determinant at the interior grid point (i, jt, k).

8. The method for numerical simulation of flow around a rotationally symmetric geometric body according to claim 5, characterized in that: In step 7, if the flow in the viscous flow problem around a rotationally symmetric geometry is axisymmetric, the first-order derivative of the flow field variable at the virtual grid point is obtained according to the following formula: ； in, and are the values ​​of the first-order derivative of the velocity component u in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component u in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component u in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component v in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component v in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component v in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component w in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component w in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of the velocity component w in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of temperature t in the x direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first derivative of temperature t in the y direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), and are the values ​​of the first-order derivative of temperature t in the z direction at the virtual grid point (i, jg, k) and the internal grid point (i, jt, k), respectively; In step 7, if the flow in the viscous flow problem around the rotationally symmetric geometric body is non-axisymmetric, the first-order derivative value of the flow field variable at the virtual grid point is calculated by the unilateral difference method based on the flow field variable values ​​at the virtual grid point and the internal grid point.

9. The method for numerical simulation of flow around a rotationally symmetric geometric body according to claim 1, characterized in that: The step 8 includes: based on the grid transformation derivatives, Jacobian determinants, flow field variable values ​​and first-order derivative values ​​of the flow field variables at the virtual grid points, the singular axis boundary grid points and the internal grid points, the spatial convection derivative terms and viscous derivative terms of the discretized control equation are discretized. At the singular axis boundary and its adjacent areas, the convection derivative terms of the control equation are discretized using a flux reconstruction difference method consistent with the method used for the internal grid points, and the viscous derivative terms of the control equation are discretized using a central difference method consistent with the method used for the internal grid points to obtain the discretized spatial derivative terms of the control equation.

10. The method for numerical simulation of flow around a rotationally symmetric geometric body according to claim 3, characterized in that: The step 10 specifically includes: Based on the updated flow field variable values ​​of the flow around the rotationally symmetric geometric body, the flow field variable values ​​on the singular axis are interpolated along the j direction for the index position (i, K); the interpolated flow field variables are replaced along the circumferential direction and averaged as the flow field variable values ​​on the singular axis of the rotationally symmetric geometric body.

Images