CN111046615B - Riemann solver shock wave instability suppression method and system - Google Patents

Riemann solver shock wave instability suppression method and system Download PDF

Info

Publication number
CN111046615B
CN111046615B CN201911374043.5A CN201911374043A CN111046615B CN 111046615 B CN111046615 B CN 111046615B CN 201911374043 A CN201911374043 A CN 201911374043A CN 111046615 B CN111046615 B CN 111046615B
Authority
CN
China
Prior art keywords
pressure
function
shock wave
solver
riemann
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201911374043.5A
Other languages
Chinese (zh)
Other versions
CN111046615A (en
Inventor
谢文佳
李桦
田正雨
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
National University of Defense Technology
Original Assignee
National University of Defense Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by National University of Defense Technology filed Critical National University of Defense Technology
Priority to CN201911374043.5A priority Critical patent/CN111046615B/en
Publication of CN111046615A publication Critical patent/CN111046615A/en
Application granted granted Critical
Publication of CN111046615B publication Critical patent/CN111046615B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Abstract

The invention discloses a method and a system for inhibiting shock wave instability of a Riemann solver, and a method and a system for simulating shock waves, wherein the method for inhibiting shock waves comprises the following steps: constructing a universal pressure dissipation function through the product of the pressure weight function and the pressure dissipation function; the pressure weight function value increases with increasing pressure gradient and the pressure weight function value varies between intervals [0,1 ]; establishing a relation between the general pressure dissipation function and a result obtained by a general Riemannian solver solving a convection term in a Venenstorx equation so as to change the dissipation performance of the Riemannian solver in a shock wave region, and activating the pressure dissipation function through a pressure weight function when the Riemannian solver is close to the shock wave region; when the shock wave is far away from the shock wave area, the pressure dissipation function is closed through the pressure weight function, and the shock wave instability phenomenon is restrained. The problems of inaccurate simulation, low calculation reliability and universality and the like in the prior art are solved, the change of the dispersity of the solver in different areas is realized, the influence on the numerical simulation precision is reduced, and the calculation reliability and universality are improved.

Description

Riemann solver shock wave instability suppression method and system
Technical Field
The invention relates to the technical field of shock wave suppression, in particular to a shock wave instability suppression method and system for a Riemann solver and a shock wave simulation method and system.
Background
With the development of computer technology and the progress of numerical methods, computational fluid dynamics has been developed. Because the ground experiment and the flight experiment are long in time consumption and high in cost, the numerical simulation increasingly becomes a main means for designing the aerodynamic layout of the aircraft. In the case of an aircraft flying in the near space, due to the high cruising speed and the high flight mach number, the high-speed airflow is often severely compressed to form shock waves around the aircraft. The occurrence of shock waves can seriously affect the aerodynamic performance of the aircraft, thereby affecting the layout design and optimization of the aircraft. Therefore, how to accurately simulate the shock wave is one of the key contents for numerically calculating the aerodynamic characteristics of the near space vehicle.
The peripheral flow field of the high-altitude aircraft is a typical compressible flow field, and one of the numerical simulation methods commonly used in engineering at present adopts a control equation of finite volume discrete fluid flow, namely a Navier-Stokes equation. The equation is a typical hyperbolic parabolic equation, and can be divided into a time term, a convection term and a viscosity term. For time terms and viscosity terms, efficient and accurate calculation methods exist, and for flow terms, the existing solution methods are numerous and have large difference in solution effect. Shock waves are one of typical phenomena of compressible flow fields, and the properties of shock waves are closely related to convection terms of the Navier-Stokes equation. Therefore, the difficulty in solving the equation numerically discretely is how to accurately and stably calculate the convection term of the equation.
At present, an approximate Riemann solver is commonly used in engineering calculation to disperse convection terms of a Navier-Stokes equation, and the commonly used Riemann solvers can be generally divided into two categories, namely a dissipative Riemann solver and a low-dissipation Riemann solver according to the difference of dissipation properties of the Riemann solvers. The difference between the two is that the dissipative Riemann solver has large numerical viscosity, and the accuracy of solving the convection term is low, so that the shock wave thickness of numerical simulation is large; the low-dissipation Riemann solver has low numerical viscosity, and when the convection term is calculated by using the Riemann solver, unstable phenomena such as oscillation or shock wave shape deformation are often generated in the calculated shock wave, so that the accuracy and the calculation efficiency of numerical simulation are seriously influenced. Thus, engineering calculations require a Riemann solver that is less dissipative and shock stable. Therefore, on the premise of keeping low dissipation, how to inhibit the shock wave instability of the Riemann solver and improve the calculation stability of the Riemann solver are of great significance for realizing efficient numerical simulation of the compressible flow field, and meanwhile, the improvement of the numerical simulation stability has great application value for engineering practical problems such as the optimization design of the aerodynamic layout of the aircraft and the like.
Currently, there are numerous improved methods for the problem of Riemann solver shock instability. The most popular method is to mix a dissipative riemann resolver with a dissipative riemann resolver by virtue of the characteristic that the dissipative riemann resolver can stably calculate a shock wave. Namely, a dissipative Riemann solver is adopted during the calculation of the shock wave, and the stable calculation of the shock wave is realized by means of the numerical viscosity of the dissipative Riemann solver; and a low-dissipation Riemann solver is adopted to calculate the flow field in a flow field area far away from the shock wave, so that the calculation precision is improved. To achieve a mix of these two riemann solvers, some shockwave probe functions are often employed to achieve the application of the two riemann solvers in different regions. However, this has the disadvantages that: currently, there is still a lack of an effective shock wave detection method, and therefore, this hybrid method may introduce additional numerical stickiness in some cases, thereby affecting the accuracy of numerical simulation. In addition, the method for restraining the shock wave instability of the Riemann solver also faces the numerical inconsistency caused by mixing two different Riemann solvers, and the reliability of calculation is possibly influenced. It is worth mentioning that this type of shock wave instability suppression method lacks versatility and can only be applied to a specific riemann solver. Therefore, this type of method for suppressing the shock instability of the riemann solver is not a good solution.
As described above, the riemann resolver, particularly the low-dissipation riemann resolver, generates a shock instability phenomenon when calculating a shock, and a general method capable of effectively suppressing the shock instability without losing the calculation accuracy is still lacking at present.
Disclosure of Invention
The invention provides a method and a system for inhibiting shock wave instability of a Riemann solver, and a shock wave simulation method and a system, which are used for overcoming the defects of lower precision and reliability, poorer universality and the like in the prior art.
In order to achieve the above object, the present invention provides a method for suppressing shock wave instability of a riemann resolver, comprising the steps of:
constructing a pressure weight function based on pressure gradient change according to pressure gradient around a grid interface obtained by initializing a compressible flow field through a finite volume method and reconstructing a processed flow field variable by adopting a limiter function; the pressure weight function value increases with increasing pressure gradient and the pressure weight function value varies between intervals [0,1 ];
constructing a universal pressure dissipation function through the product of the pressure weight function and the pressure dissipation function;
establishing a relation between a general pressure dissipation function and a convection term solving result in a Venus-Tox equation obtained by solving flow field variables inside the grid unit through a general Riemann solver so as to change the dissipative property of the Riemann solver in a shock wave area and update a numerical flux function value on a grid interface, and activating the pressure dissipation function through a pressure weight function when the pressure dissipation function is close to the shock wave area; the pressure dissipation function is closed by the pressure weight function when away from the shock region.
In order to achieve the above object, the present invention further provides a system for suppressing the shock wave instability of the riemann resolver, which includes a memory and a processor, wherein the memory stores a riemann resolver shock wave instability suppression program, and the processor executes the steps of the riemann resolver shock wave instability suppression method when the riemann resolver shock wave instability suppression program is run.
In order to achieve the above object, the present invention further provides a shock wave simulation method based on a riemann solver, including the following steps:
initializing the inflow conditions of the compressible flow field as variables by adopting a finite volume method;
reconstructing initialized variables in the flow field by adopting a limiter function to obtain flow field variable distribution in each grid unit;
calculating a flow field in the grid unit by adopting a Riemann solver to obtain a numerical flux function value on a grid interface, namely a discrete value of a flow term of a Naviserstokes equation;
updating the numerical flux function value on the grid interface by the Riemann solver shock wave instability suppression method;
judging whether the iteration time reaches an iteration condition or not, and returning to the flow field variable reconstruction step when the iteration time does not reach the iteration condition; and outputting a flow field result when the iteration time is reached when the iteration condition is reached.
In order to achieve the above object, the present invention further provides a shock wave simulation system based on a riemann resolver, which includes a memory and a processor, wherein the memory stores a shock wave simulation program based on the riemann resolver, and the processor executes the steps of the shock wave simulation method based on the riemann resolver when the processor runs the shock wave simulation program based on the riemann resolver.
When a compressible flow field around a high-altitude aircraft is simulated, firstly, initializing the compressed flow field by taking an incoming flow condition as a variable through a finite element volume method, and then reconstructing the initialized variable in the flow field by adopting a limiter function to obtain the variable distribution of the flow field inside each grid unit; calculating a flow field in the grid unit by adopting a Riemann solver to obtain a numerical flux function value on a grid interface, namely a pressure value, constructing a pressure weight function according to the pressure gradient change of the adjacent grid interface in the flow field, and judging a shock wave region; and then, a universal pressure dissipation function is constructed according to the pressure weight function and the pressure dissipation function, the dissipation property of the original Riemann solver near the shock wave is changed, the shock wave instability of the original solver is inhibited, and the shock wave instability phenomenon is further inhibited. In order to control the pressure dissipation function to only act near the shock wave without influencing the solving precision of other areas of the flow field, the invention defines the pressure weight function to control the size of the pressure dissipation function in different areas of the flow field. In the vicinity of the shock wave, the pressure weight function is unit 1, and at the moment, the pressure dissipation function acts to inhibit the shock wave from being unstable; in other areas far away from the shock wave, the value of the pressure dissipation function is zero, and the pressure dissipation function disappears; and the specific action of the pressure dissipation function is determined by the pressure gradient in the region close to the shock wave to different extents.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the embodiments or the prior art will be briefly described below, it is obvious that the drawings in the following description are only some embodiments of the present invention, and for those skilled in the art, other drawings can be obtained according to the structures shown in the drawings without creative efforts.
Fig. 1 is a basic flow diagram of a finite volume method in a method for suppressing shock wave instability of a riemann solver according to an embodiment of the present invention;
FIG. 2 is a flow chart of an improved flux calculation in the first embodiment;
FIG. 3 is a schematic diagram of the calculation of the pressure weight function f;
FIG. 4 is a density cloud of a blunt streaming flow with the original low dissipation Riemann solver calculating the incoming stream Mach number 20;
FIG. 5 is a density flow field result calculated by the Riemann solver shock wave instability suppression method provided by the invention;
FIG. 6 is a density cloud plot of the original low-dissipation Riemann solver-computed double Mach number reflection problem;
fig. 7 is a density flow field result calculated by using the method for suppressing the shock wave instability of the riemann resolver provided by the invention.
The implementation, functional features and advantages of the present invention will be further described with reference to the accompanying drawings.
Detailed Description
The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
It should be noted that all directional indicators (such as upper, lower, left, right, front and rear … …) in the embodiment of the present invention are only used to explain the relative position relationship between the components, the movement situation, etc. in a specific posture (as shown in the drawing), and if the specific posture is changed, the directional indicator is changed accordingly.
In addition, the descriptions related to "first", "second", etc. in the present invention are only for descriptive purposes and are not to be construed as indicating or implying relative importance or implicitly indicating the number of technical features indicated. Thus, a feature defined as "first" or "second" may explicitly or implicitly include at least one such feature. In the description of the present invention, "a plurality" means at least two, e.g., two, three, etc., unless explicitly specified otherwise.
In the present invention, unless otherwise expressly stated or limited, the terms "connected," "secured," and the like are to be construed broadly, and for example, "secured" may be a fixed connection, a removable connection, or an integral part; the connection can be mechanical connection, electrical connection, physical connection or wireless communication connection; they may be directly connected or indirectly connected through intervening media, or they may be connected internally or in any other suitable relationship, unless expressly stated otherwise. The specific meanings of the above terms in the present invention can be understood by those skilled in the art according to specific situations.
In addition, the technical solutions in the embodiments of the present invention may be combined with each other, but it must be based on the realization of the technical solutions by those skilled in the art, and when the technical solutions are contradictory to each other or cannot be realized, such a combination of the technical solutions should be considered to be absent, and is not within the protection scope claimed by the present invention.
Example one
As shown in fig. 1 to 4, an embodiment of the present invention provides a method for suppressing shock wave instability of a riemann resolver, including:
the discrete NaViStokes equation using the finite volume method is as follows:
Figure BDA0002340439250000061
in the formula (1), U represents a conservation variable vector in a flow field, omega represents the volume of a computational grid, N is the number of grids, Delta S represents the length of a grid interface, and F c And F v The convection flux vector and the viscous flux vector are respectively expressed, are functions of the conservation variable U, can be directly obtained from state variables in the flow field, and the flow term in equation (1) is solved by using a riemann solver, which can be generally expressed as:
Figure BDA0002340439250000062
a typical riemann solver can be represented in the form of equation (2). In the formula (2), U L And U R Representing the conservation variable vector, function phi (U), on the left and right sides of the computational grid interface, respectively L ,U R ) For the dissipation factor matrix is U L And U R Determines the dissipative properties of the Riemann solver, the dissipative properties of different Riemann solvers being different in function phi (U) L ,U R ) By different definitions of (c). In order to inhibit the instability of the function when solving the laser, the invention provides a general improved method, namely, the F is calculated c An additional pressure dissipation function is added, which can be specifically expressed as:
Figure BDA0002340439250000063
wherein F p Representing a pressure dissipation function, and f represents a pressure weight function, which is unity near the shock and zero in other regions; the function f varies between zero and one in the shock transition region and increases with increasing pressure gradient.
Constructing a pressure weight function based on pressure gradient change according to the pressure gradient; the pressure weight function value increases with increasing pressure gradient and the pressure weight function value varies between intervals [0,1 ]. When the pressure gradient takes the maximum value, the pressure weight function value tends to 1; when the pressure gradient is minimized, the pressure weight function value approaches 0.
The step of obtaining the pressure weight function value on the determined grid interface i comprises the following steps:
acquiring a plurality of adjacent grid interfaces k adjacent to the grid interface i;
for the grid interface i and each adjacent grid interface k, the self left-right pressure gradient is respectively obtained (see the following formula)
Figure BDA0002340439250000071
) Right and left pressure ladderDegree (see in the following formula
Figure BDA0002340439250000072
) Taking the m-th power of the minimum value of the pressure gradient as a calculation result, wherein m is more than or equal to 3;
taking the minimum value of the calculation result in the grid interface i and each adjacent grid interface k;
and subtracting the minimum value of the calculation result from the value 1 to obtain the pressure weight function value on the grid interface i.
Preferably, the expression of the pressure weight function is:
Figure BDA0002340439250000073
Figure BDA0002340439250000074
wherein the variable p L And p R Respectively representing the pressure on the left side and the right side of a grid interface k, wherein i is any one natural number from 1 to 7, and the rest 6 grid interfaces are adjacent to the grid interface i.
Referring to FIG. 3, the function f near the shock wave is large due to the large pressure gradient variation k The value of (b) is close to 0, and the value of the pressure weight function f is close to 1 at the moment; in the smooth region away from the shock wave, the function f is not much changed due to the pressure gradient k The value of (a) is close to 1, and the value of the pressure weight function f is close to 0.
Preferably, the pressure dissipation function F p The specific expression of (a) is then:
Figure BDA0002340439250000075
in equation (7), Δ p is the pressure difference across the grid interface, i.e., Δ p ═ p R -p L ,R 2 Is the second eigenvector of the right eigenvector of the jacobian matrix of euler's equations,
Figure BDA0002340439250000076
is the speed of sound calculated using the mean amount of Roe,
Figure BDA0002340439250000077
is the Mach number calculated from the mean amount of Roe, and S L And S R Then represent wave velocities, which are calculated as follows:
S L =min(u L -c L ,u R -c R ) S R =max(u L +c L ,u R +c R ) (8)
in the formula (8), u L 、u R Representing the velocities of the left and right sides of the grid interface, c L 、c R And sound velocities of the left side and the right side of the grid interface are represented, and the variables can be directly obtained by calculating state variables in the calculation grid in the flow field of the last time step. The problem of numerical instability of the shock is caused by the propagation of disturbance errors generated in the vicinity of the shock downstream and laterally of the shock. Pressure dissipation function F p The method has the advantages that the propagation of pressure disturbance errors near the shock wave is restrained through the newly constructed pressure dissipation terms near the shock wave, and therefore the purpose of restraining the shock wave instability is achieved. The pressure dissipation function meets the requirement of conservation of constancy, and can not generate negative influence on the numerical calculation of viscous flow in a smooth area while inhibiting the problem of unstable shock waves. From equation (3), it can be seen that the function F is unit 1 in the vicinity of the shock wave, when the convection term is from the original function F c Changing to a modified function
Figure BDA0002340439250000081
So that the pressure dissipation function F p Plays a role in inhibiting possible shock wave instability; and in the region far away from the shock wave, the function f is zero, and the correction function of the convection term
Figure BDA0002340439250000082
Is changed into F c At this time, the pressure dissipation function F p And thus, does not affect the computational accuracy of the original method.
The technical scheme of the invention has the following technical effects:
(1) according to the invention, the shock wave instability problem of the original Riemann solver is inhibited by adopting the additional pressure dissipation function near the shock wave, the method is not limited by the original Riemann solver, has universality, and is suitable for various different types of Riemann solvers;
(2) the improved method of the pressure dissipation function adopted by the invention can not lose the solving precision of the original Riemann solver;
(3) the improved method provided by the invention can be conveniently implanted into the existing algorithm and is used for inhibiting the shock wave instability problem of the existing Riemann solver.
Example two
Fig. 1 shows a basic flow chart for solving the navistokes equation in the finite volume format of the riemann solver, and a specific embodiment is as follows. The method comprises the following steps:
initializing the inflow conditions of the compressible flow field as variables by adopting a finite volume method;
reconstructing the initialized variables in the flow field by adopting a limiter function to obtain the flow field variable distribution in each grid unit;
calculating a flow field in the grid unit by adopting a Riemann solver to obtain a numerical flux function value on a grid interface, namely a discrete value of a convection term of a NavieStokes equation;
updating a numerical flux function value on a grid interface by a Riemann solver shock wave instability suppression method of any embodiment I;
judging whether the iteration time reaches an iteration condition or not, and returning to the flow field variable reconstruction step when the iteration time does not reach the iteration condition; and outputting a flow field result when the iteration time is reached when the iteration condition is reached.
Firstly, a finite volume method (known technology) is used to perform initialization processing on a flow field for some compressible flow problem needing to be calculated, that is, the flow field property variables such as flow speed, pressure, temperature and the like in the flow field are set as the numerical values such as the speed, pressure, temperature and the like of incoming flow, and then initialization is performed. Then, a limiter function is used to carry out the variable in the flow fieldReconstructing (one of steps of a computational fluid dynamics finite volume method, which is used for redistributing variables in a flow field) to obtain flow field variable distribution inside each grid unit (a grid unit refers to a certain grid, such as a cuboid, a grid interface refers to a boundary surface of the grid, such as the surface of the cuboid, and parameters such as pressure gradient of the grid interface can be obtained under the condition that the flow field variable distribution inside the grid unit is known). Then, a numerical flux function F on a grid interface is calculated by adopting a Riemann solver c I.e. the dispersion of the equation to the flow term. At this time, a pressure dissipation function for suppressing the shock instability is calculated from the flow field variables around the mesh interface. And after the numerical flux of the grid interface is obtained through calculation, iterative solution is carried out by adopting a time format to obtain the flow field information of the next time step, so that the updating of the flow field variable is completed. And finally, judging whether the numerical calculation is finished or not according to the specific flow problem of the calculation, and outputting a flow field result.
Figure 2 shows a flow diagram of an improved numerical flux function calculation. Firstly, according to flow field information of the left and right sides of a grid interface, a Riemann solver is adopted to calculate a numerical flux function F at the interface c (ii) a Then, a pressure dissipation function F is calculated according to flow field information around the grid interface p And a pressure weight function f (according to the characteristic that the pressure gradient near the shock wave changes greatly, the pressure weight function is designed, the pressure weight function is a function of the pressure gradient near the shock wave, the value near the shock wave is approximate to 1, and the value far away from the shock wave is approximate to 0); finally, the corrected numerical flux is calculated
Figure BDA0002340439250000091
And finally, updating the flow field variable by adopting a time iterative algorithm to obtain a flow field result after time advancing, and outputting a calculation result of the flow field. Embedding the improved flux calculation method of fig. 2 into the calculation framework of the original finite volume method (as shown in fig. 1) results in a method for suppressing shock instability of the riemann solver. The invention realizes the problem of restraining the shock wave of the original Riemann solver from being unstable by adding an additional pressure dissipation term near the shock wave.
Figure 4 shows a density cloud of blunt streaming flow at mach number 20 of the incoming stream calculated by the original low dissipation riemann solver. Taking a typical low-dissipation HLLC Riemann solver as an example, it can be seen that a shock wave profile calculated by the original HLLC Riemann solver is distorted near a stagnation point line, and a shock wave abnormal phenomenon occurs in a density flow field. Fig. 5 shows the density flow field result calculated by the improved method for suppressing shock instability, which is provided by the present invention, and it can be seen that after the improved scheme provided by the present invention is adopted, the shock profile obtained by calculation is accurate, and no obvious shock abnormal phenomenon occurs. The results show that the method for inhibiting shock instability provided by the invention is effective.
Figure 6 shows a density cloud of the original low-dissipation double mach number reflection problem calculated by the riemann solver. Still taking a typical low-dissipation HLLC Riemannian solver as an example, the phenomenon of an inclined Mach rod appearing on a shock wave profile calculated by the original HLLC Riemannian solver can be seen, and an abnormity appears in a density flow field. Figure 7 shows the density flow field results calculated using the improved method of suppressing shock instability proposed by the present invention. It can be seen that after the improved scheme provided by the invention is adopted, the shock wave shape obtained by calculation is accurate, and no obvious shock wave abnormal phenomenon occurs. The results show that the method for inhibiting shock instability provided by the invention is effective.
EXAMPLE III
Based on the first embodiment, the invention provides a system for suppressing the shock instability of a riemann resolver, which comprises a memory and a processor, wherein the memory stores a riemann resolver shock instability suppression program, and the processor executes the steps of any embodiment of the method when the riemann resolver shock instability suppression program is operated.
Based on the second embodiment, the invention provides a shock wave simulation system based on a riemann solver, which comprises a memory and a processor, wherein the memory stores a shock wave simulation program based on the riemann solver, and the processor executes the steps of the method of any second embodiment when running the shock wave simulation program based on the riemann solver.
The above description is only a preferred embodiment of the present invention, and is not intended to limit the scope of the present invention, and all modifications and equivalents of the present invention, which are made by the contents of the present specification and the accompanying drawings, or directly/indirectly applied to other related technical fields, are included in the scope of the present invention.

Claims (9)

1. A method for suppressing shock wave instability of a Riemann solver is characterized by comprising the following steps:
constructing a pressure weight function based on pressure gradient change according to pressure gradient around a grid interface obtained by initializing a compressible flow field through a finite volume method and reconstructing a processed flow field variable by adopting a limiter function; the pressure weight function value increases with increasing pressure gradient and the pressure weight function value varies between intervals [0,1 ];
constructing a universal pressure dissipation function through the product of the pressure weight function and the pressure dissipation function;
establishing a relation between a general pressure dissipation function and a convection term solving result in a Venus-Tox equation obtained by solving flow field variables inside grid units through a general Riemann solver so as to change the dissipation performance of the Riemann solver in a shock wave region and update a numerical flux function value on a grid interface, and activating the pressure dissipation function through the pressure weight function when the pressure dissipation function is close to the shock wave region; the pressure dissipation function is closed by the pressure weight function when away from the shock region.
2. The Riemann solver shock instability mitigation method of claim 1, the step of obtaining the pressure weight function values for a determined grid interface i comprising:
acquiring a plurality of adjacent grid interfaces k adjacent to the grid interface i;
respectively acquiring the pressure gradient of the grid interface i and each adjacent grid interface k thereof, and taking the m-th power of the minimum value of the pressure gradient as a calculation result, wherein m is more than or equal to 3;
taking the minimum value of the calculation result in the grid interface i and each adjacent grid interface k;
and subtracting the minimum value of the calculation result from the value 1 to obtain the pressure weight function value on the grid interface i.
3. The Riemann solver shock instability suppression method as claimed in claim 2, wherein the pressure weight function is expressed as:
Figure FDA0002340439240000011
Figure FDA0002340439240000012
wherein the variable p L And p R Respectively representing the pressure on the left side and the pressure on the right side of a grid interface k, wherein i is any one natural number from 1 to 7, and the rest 6 grid interfaces are all adjacent to the grid interface i.
4. The Riemann solver shock instability suppression method as claimed in any one of claims 1 to 3, wherein the pressure dissipation function F p The specific expression of (A) is as follows:
Figure FDA0002340439240000021
Δ p is the pressure difference across the grid interface, R 2 Is the second eigenvector of the right eigenvector of the jacobian matrix of euler's equations,
Figure FDA0002340439240000022
is the speed of sound calculated using the mean amount of Roe,
Figure FDA0002340439240000023
is the Mach number calculated from the mean amount of Roe, and S L And S R Then represent wave velocities, which are asThe following way is calculated:
S L =min(u L -c L ,u R -c R )S R =max(u L +c L ,u R +c R )
in the formula u L 、u R Representing the velocities of the left and right sides of the grid interface, c L 、c R Representing the sound velocities on the left and right sides of the grid interface.
5. The Riemann solver shock instability mitigation method of claim 4, wherein the generic pressure dissipation function is:
Figure FDA0002340439240000024
wherein F is a pressure weight function, F p As a function of pressure dissipation.
6. The Riemann solver shock instability suppression method of claim 5, wherein the step of relating the generic pressure dissipation function to the result of the generic Riemann solver solving the convection term in the Venenstorx equation to change the dissipation of the Riemann solver in the shock region comprises:
the discrete NaviStokes equation by the finite volume method is as follows:
Figure FDA0002340439240000025
in the formula (1), U represents a conservation variable vector in a flow field, omega represents the volume of a computational grid, N is the number of grids, Delta S represents the length of a grid interface, and F c And F v The flow flux vector and the viscous flux vector are respectively expressed and are functions of a conservation variable U, the functions are directly obtained by state variables in a flow field, and a Riemann solver is adopted to solve the flow term in the equation (1), and the expression is as follows:
Figure FDA0002340439240000026
in formula (2), U L And U R Representing the conservation variable vector, function phi (U), on the left and right sides of the computational grid interface, respectively L ,U R ) For the dissipation coefficient matrix is U L And U R A function of (a); in the calculation of F c An additional pressure dissipation function is added, which can be specifically expressed as:
Figure FDA0002340439240000031
wherein F p Representing a pressure dissipation function; f represents a pressure weight function.
7. A Riemann solver shock instability suppression system comprising a memory having a Riemann solver shock instability suppression program stored therein and a processor that performs the steps of the method of any one of claims 1 to 6 when running the Riemann solver shock instability suppression program.
8. A shock wave simulation method based on a Riemann solver is characterized by comprising the following steps:
initializing the incoming flow condition of the compressible flow field as a variable by adopting a finite volume method;
reconstructing initialized variables in the flow field by adopting a limiter function to obtain flow field variable distribution in each grid unit;
calculating a flow field in the grid unit by adopting a Riemann solver to obtain a numerical flux function value on a grid interface, namely a discrete value of a flow term of a Naviserstokes equation;
updating the numerical flux function value on the grid interface by the Riemann solver shock wave instability suppression method according to any one of claims 1-6;
judging whether the iteration time reaches an iteration condition or not, and returning to the flow field variable reconstruction step when the iteration time does not reach the iteration condition; and outputting a flow field result when the iteration time is reached when the iteration condition is reached.
9. A riemann solver-based shock wave simulation system comprising a memory storing a riemann solver-based shock wave simulation program and a processor that, when running the riemann solver-based shock wave simulation program, performs the steps of the method of claim 8.
CN201911374043.5A 2019-12-27 2019-12-27 Riemann solver shock wave instability suppression method and system Active CN111046615B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911374043.5A CN111046615B (en) 2019-12-27 2019-12-27 Riemann solver shock wave instability suppression method and system

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911374043.5A CN111046615B (en) 2019-12-27 2019-12-27 Riemann solver shock wave instability suppression method and system

Publications (2)

Publication Number Publication Date
CN111046615A CN111046615A (en) 2020-04-21
CN111046615B true CN111046615B (en) 2022-09-23

Family

ID=70240445

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911374043.5A Active CN111046615B (en) 2019-12-27 2019-12-27 Riemann solver shock wave instability suppression method and system

Country Status (1)

Country Link
CN (1) CN111046615B (en)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113486453B (en) * 2021-09-07 2021-11-19 中国空气动力研究与发展中心低速空气动力研究所 Numerical simulation method of compressible fluid
CN113987691B (en) * 2021-11-25 2022-08-23 中国人民解放军国防科技大学 High-precision hybrid calculation method, device, equipment and storage medium for shock wave instability

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102890751A (en) * 2012-09-18 2013-01-23 天津空中代码工程应用软件开发有限公司 Numerical method for solving two-dimensional Riemannian problem and simulating subsonic non-viscous stream
CN103823916B (en) * 2013-10-23 2016-09-14 沈智军 A kind of arbitary Lagrangian-Eularian based on multidimensional Riemann Solution

Also Published As

Publication number Publication date
CN111046615A (en) 2020-04-21

Similar Documents

Publication Publication Date Title
CN108563843B (en) Method for updating disturbance area of steady compressible flow
Aftosmis et al. Behavior of linear reconstruction techniques on unstructured meshes
CN111046615B (en) Riemann solver shock wave instability suppression method and system
CN115238397B (en) Method and device for calculating thermal environment of hypersonic aircraft and computer equipment
JP2008165804A (en) Flow simulation calculating method and system
CN113850008B (en) Self-adaptive grid disturbance domain updating acceleration method for aircraft aerodynamic characteristic prediction
CN105122248A (en) Method for modelling a part, in particular a turbine blade
CN116245049B (en) Node type non-structural grid boundary correction method, device, equipment and medium
CN111079228A (en) Pneumatic shape optimization method based on flow field prediction
Swanson et al. An efficient solver for the RANS equations and a one-equation turbulence model
Peng et al. Numerical investigation of the effects of structural geometric and material nonlinearities on limit-cycle oscillation of a cropped delta wing
CN115774968A (en) Curved surface mesh generation method and system based on recursive decomposition and computer equipment
US20130191087A1 (en) Method of calculating dynamic pressure at the level of an aircraft surface
Su Accurate and robust adaptive mesh refinement for aerodynamic simulation with multi‐block structured curvilinear mesh
Mei et al. Implicit numerical simulation of transonic flow through turbine cascades on unstructured grids
Park et al. Validation of an output-adaptive, tetrahedral cut-cell method for sonic boom prediction
Gou et al. A high‐order element based adaptive mesh refinement strategy for three‐dimensional unstructured grid
Galbraith et al. Full potential revisited: A medium fidelty aerodynamic analysis tool
Roy et al. A finite volume method for viscous incompressible flows using a consistent flux reconstruction scheme
JP6163897B2 (en) Numerical calculation program, numerical calculation method, and information processing apparatus
JP3587827B2 (en) Airfoil performance estimation method and airfoil performance estimation program
Diskin et al. Evaluation of multigrid solutions for turbulent flows
WO2013078912A1 (en) Numerical method for simulating wake field of ship propeller
Park et al. Output-adaptive tetrahedral cut-cell validation for sonic boom prediction
Pandya et al. Assessment of Preconditioner for a USM3D Hierarchical Adaptive Nonlinear Method (HANIM)

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant