CN102880589B - A kind of method that finite-difference modeling detecting mesh quality is flowed around complex configuration - Google Patents

Abstract

The invention discloses a kind of method that finite-difference modeling detecting mesh quality is flowed around complex configuration, it comprises: the generation of step one, computing grid, measuring and adjustation; Step 2, control differential equation discrete and solving; Step 3, postpositive disposal; It focuses on: generate complicated polylith docking structure grid to whole flowing zoning; The mesh transformations Jacobi calculated after the mathematical definition formula of mesh transformations Jacobi, symmetric form and symmetry conservation form adopt linear finite difference operator discrete is designated as V respectively ¹, V ²and V ³; According to described V ¹, V ²and V ³determine the Testing index of mesh quality; Testing index according to described mesh quality detects the computing grid generated, and determines that those do not reach the particular location of the grid cell of Testing index; Again detect after the described grid cell not reaching Testing index is adjusted, by this repeatably process until the grid cell after this adjustment all reaches Testing index.

Description

A kind of method that finite-difference modeling detecting mesh quality is flowed around complex configuration

Technical field

The present invention relates to the method for numerical simulation of a kind of complicated calculations region physical characteristics, particularly a kind of method detected around the finite-difference modeling Complex Flows of complex configuration zoning mess generation quality.

Background technology

Although wind tunnel test is still the important means predicting various space flight and aviation aerodynamic characteristics of vehicle up till now, numerical simulation technology plays more and more important effect in pneumatic design.In recent years, along with the continuous lifting of computing machine floating-point operation ability and the perfect gradually of numerical computation method, people more and more favor in high precision, high-resolution numerical computation method.The existing large quantifier elimination fact shows: the method for numerical simulation obtaining Numerical solution of partial defferential equatio based on traditional second order accuracy Finite Volume Method discrete differential equation is not well positioned to meet the demand of Practical Project problem to computational accuracy, especially the separated flow on a large scale of fighter plane generation when doing large attack angle flight, the Finite Volume Method of second order accuracy often can not provide gratifying numerical simulation result, needs to adopt high precision, high-resolution numerical method to simulate.The separated flow on a large scale that fighter plane produces when doing large attack angle flight is that one flows around complex configuration.Refer to that air or water or other fluid walk around the Complex Flows of the various true aviation aircraft complex configuration such as aircraft, guided missile or aircraft under water around complex configuration flowing.

When adopting higher order accuracy finite difference method to solve in polylith docking structure grid, how to ensure that the higher order accuracy characteristic of numerical simulation result seems particularly important.As everyone knows, the numerical solution of higher order accuracy be expected except needing to adopt the numerical discretization schemes of higher order accuracy, traditional second order accuracy numerical method is also different to the requirement of mess generation quality.It is important to note that the grid higher order accuracy numerical method generated can carry out calculating and higher order accuracy numerical method smoothly, can to obtain good numerical solution be two diverse concepts.The grid higher order accuracy numerical method generated can carry out smoothly calculating and first will meet, and just likely adopts higher order accuracy numerical method to obtain good numerical solution on this basis.The two is also not quite similar to the requirement of mess generation quality.Obviously, if want to utilize higher order accuracy numerical method to obtain good numerical solution, its to the requirement of mesh quality than being only that the requirement that higher order accuracy numerical method can calculate smoothly is high.So, how to judge generate grid can higher order accuracy numerical method discrete under obtain good numerical solution? the index that there is the mesh quality that certain can quantize can be good at the quality problems pointing out grid? problems rarely has domestic and foreign literature to mention.In fact the computing grid around complex configuration is normally made up of the grid cell of millions of even several hundred million, in a lot of situation the initial computing grid generated can not higher order accuracy numerical method discrete under obtain good numerical solution be only due to only a few (as several or tens) existence of inferior grid cell causes, and the inferior grid cell of these only a fews is not too easy to be identified from a large amount of computing grid unit.Mesh quality detection method disclosed by the invention can detect the quality problems of computing grid effectively, determine the particular location of quality grid cell not up to standard, by adjusting simply to ensure mesh quality, being conducive to simulation when complex configuration flows, adopting higher order accuracy numerical method to obtain good numerical solution smoothly.

Summary of the invention

The object of the invention is to solve and use the simulation of higher order accuracy finite difference method when complex configuration flow, how to ensure that mess generation quality is to obtain good or to meet the problem of numerical solution of Practical precision.

In order to achieve the above object, the technical solution used in the present invention is as follows:

The method that the finite-difference modeling detecting mesh quality is flowed around complex configuration, it comprises:

The generation of step one, computing grid, measuring and adjustation;

Step 2, control differential equation discrete and solving;

Step 3, postpositive disposal;

To it is characterized in that in described step one 1), complicated polylith docking structure grid is generated to whole flowing zoning, 2), the mesh transformations Jacobi that calculates after adopting linear finite difference operator discrete of the mathematical definition formula of mesh transformations Jacobi, symmetric form and symmetry conservation form is designated as V respectively ¹, V ²and V ³, 3), according to described V ¹, V ²and V ³determine the Testing index V of mesh quality _nand V _m, described V _nfor described V ¹, V ²and V ³minimum value; Described V _mfor V ¹², V ²³and V ¹³maximal value, wherein V ¹², V ²³and V ¹³by described V ¹, V ²and V ³by V ¹²=∣ V ¹?V ²∣/V ², V ²³=∣ V ²?V ³∣/V ³, V ¹³=∣ V ¹?V ³∣/V ³determine, 4), according to the Testing index of described mesh quality, the computing grid generated is detected, determine that those do not reach the particular location of the grid cell of Testing index, 5), the described grid cell not reaching Testing index is adjusted after again detect, by this repeatably process until the grid cell after this adjustment all reaches Testing index; In described step 2, the discrete control differential equation of higher order accuracy finite difference method is adopted to whole flowing zoning, and it is solved; In described step 3, logarithm value simulated flow pattern analyzes the physical significance of its correspondence.

The present invention has following technique effect:

1) because grid detection only have employed mesh transformations Jacobi variable, testing process facilitates feasible, and examination criteria is clearly easy to judge;

2) because grid detection can adopt the different form of calculation of mesh transformations Jacobi, the grid system that the form of calculation of often kind of mesh transformations Jacobi relates to is incomplete same, can effectively more comprehensively detect the quality of mess generation;

3) quality owing to adopting the grid detection technology in the present invention can generate than more comprehensive checking network lattice, by after grid detection, adopts the numerical method of higher order accuracy successfully can obtain the numerical solution of good higher order accuracy.

In addition, all there is many forms in three kinds of form of calculation of mesh transformations Jacobi of the present invention, and these many forms are discrete lower equivalent at linear finite difference operator.

Accompanying drawing explanation

Fig. 1 is the computing grid schematic diagram of the mesh transformations Jacobi of embodiment

Fig. 2 is the volume schematic diagram that the mathematical definition formula (7) of the mesh transformations Jacobi of embodiment calculates

Fig. 3 is the volume schematic diagram (volume fractiion relevant to ξ direction vector area) that the symmetric form (9) of the mesh transformations Jacobi of embodiment calculates

Fig. 4 is the volume schematic diagram (volume fractiion relevant to η direction vector area) that the symmetric form (9) of the mesh transformations Jacobi of embodiment calculates

Fig. 5 is the volume schematic diagram (volume fractiion relevant to ζ direction vector area) that the symmetric form (9) of the mesh transformations Jacobi of embodiment calculates

Fig. 6 is the volume schematic diagram (volume fractiion on the left side, ξ direction) that the conservation symmetric form (16) of the mesh transformations Jacobi of embodiment calculates

Fig. 7 is the volume schematic diagram (volume fractiion on the right of ξ direction) that the conservation symmetric form (16) of the mesh transformations Jacobi of embodiment calculates

Fig. 8 is the volume schematic diagram (volume fractiion on the left side, η direction) that the conservation symmetric form (16) of the mesh transformations Jacobi of embodiment calculates

Fig. 9 is the volume schematic diagram (volume fractiion on the right of η direction) that the conservation symmetric form (16) of the mesh transformations Jacobi of embodiment calculates

Figure 10 is the volume schematic diagram (volume fractiion on the left side, ζ direction) that the conservation symmetric form (16) of the mesh transformations Jacobi of embodiment calculates

Figure 11 is the volume schematic diagram (volume fractiion on the right of ζ direction) that the conservation symmetric form (16) of the mesh transformations Jacobi of embodiment calculates

Embodiment

Illustrate for the differential equation under following rectangular coordinate system by reference to the accompanying drawings:

$\frac{\∂Q}{\∂t}+\frac{\∂E}{\∂x}+\frac{\∂F}{\∂y}+\frac{\∂G}{\∂z}=0---\left(1\right)$

Wherein Q is the physical descriptor solved in complex configuration flowing, and E, F and G are the function about Q.Equation (1) carry out in polylith docking structure grid finite difference discrete time, under coordinates computed system need be transformed to, set up coordinates computed system (τ, ξ, η, ζ) and rectangular coordinate system (t, x, y, z) between one to one transformation relation be:

$\left\{\begin{array}{c}\mathrm{\τ}=t\\ \mathrm{\ξ}=\mathrm{\ξ}(t,x,y,z)\\ \mathrm{\η}=\mathrm{\η}(t,x,y,z)\\ \mathrm{\ζ}=\mathrm{\ζ}(t,x,y,z)\end{array}\right.---\left(2\right)$

Then the form of expression of equation (1) under coordinates computed system is:

$\frac{\∂\hat{Q}}{\∂\mathrm{\τ}}+\frac{\∂\hat{E}}{\∂\mathrm{\ξ}}+\frac{\∂\hat{F}}{\∂\mathrm{\η}}+\frac{\∂\hat{G}}{\∂\mathrm{\ζ}}=0---\left(3\right)$

Wherein:

$\begin{array}{c}\hat{Q}=J^{-1}Q\\ \hat{E}={\widehat{\mathrm{\ξ}}}_{t}Q+{\widehat{\mathrm{\ξ}}}_{x}E+{\widehat{\mathrm{\ξ}}}_{y}F+{\widehat{\mathrm{\ξ}}}_{z}G\\ \hat{F}={\widehat{\mathrm{\η}}}_{t}Q+{\widehat{\mathrm{\η}}}_{x}E+{\widehat{\mathrm{\η}}}_{y}F+{\widehat{\mathrm{\η}}}_{z}G\\ \hat{G}={\widehat{\mathrm{\ζ}}}_{t}Q+{\widehat{\mathrm{\ζ}}}_{x}E+{\widehat{\mathrm{\ζ}}}_{y}F+{\widehat{\mathrm{\ζ}}}_{z}G\end{array}---\left(4\right)$

In static grid, the mathematical definition formula of grid derivative is:

$\begin{array}{c}\left\{\begin{array}{c}{\widehat{\mathrm{\ξ}}}_x=J^{-1}{\mathrm{\ξ}}_x=y_{\mathrm{\η}}{z}_{\mathrm{\ζ}}-y_{\mathrm{\ζ}}{z}_{\mathrm{\η}}\\ {\widehat{\mathrm{\ξ}}}_y=J^{-1}{\mathrm{\ξ}}_y=x_{\mathrm{\ζ}}{z}_{\mathrm{\η}}-x_{\mathrm{\η}}{z}_{\mathrm{\ζ}}\\ {\widehat{\mathrm{\ξ}}}_z=J^{-1}{\mathrm{\ξ}}_z=x_{\mathrm{\η}}{y}_{\mathrm{\ζ}}-x_{\mathrm{\ζ}}{y}_{\mathrm{\η}}\end{array}\right.\\ \left\{\begin{array}{c}{\widehat{\mathrm{\η}}}_x=J^{-1}{\mathrm{\η}}_x=y_{\mathrm{\ζ}}{z}_{\mathrm{\ξ}}-y_{\mathrm{\ξ}}{z}_{\mathrm{\ζ}}\\ {\widehat{\mathrm{\η}}}_y=J^{-1}{\mathrm{\η}}_y=x_{\mathrm{\ξ}}{z}_{\mathrm{\ζ}}-x_{\mathrm{\ζ}}{z}_{\mathrm{\ξ}}\\ {\widehat{\mathrm{\η}}}_z=J^{-1}{\mathrm{\η}}_z=x_{\mathrm{\ζ}}{y}_{\mathrm{\ξ}}-x_{\mathrm{\ξ}}{y}_{\mathrm{\ζ}}\end{array}\right.\\ \left\{\begin{array}{c}{\widehat{\mathrm{\ζ}}}_x=J^{-1}{\mathrm{\ζ}}_x=y_{\mathrm{\ξ}}{z}_{\mathrm{\η}}-y_{\mathrm{\η}}{z}_{\mathrm{\ξ}}\\ {\widehat{\mathrm{\ζ}}}_y=J^{-1}{\mathrm{\ζ}}_y=x_{\mathrm{\η}}{z}_{\mathrm{\ξ}}-x_{\mathrm{\ξ}}{z}_{\mathrm{\η}}\\ {\widehat{\mathrm{\ζ}}}_z=J^{-1}{\mathrm{\ζ}}_z=x_{\mathrm{\ξ}}{y}_{\mathrm{\η}}-x_{\mathrm{\η}}{y}_{\mathrm{\ξ}}\end{array}\right.\end{array}---\left(5\right)$

Wherein subscript represents partial derivative, as x _ξdenotation coordination x is to the partial derivative in coordinates computed ξ direction.

In order to can strict meeting geometric conservation law, the geometrical property of computing grid can be reflected again exactly simultaneously, need the symmetry conservation form of calculation adopting grid derivative:

$\begin{array}{c}\left\{\begin{array}{c}{\widehat{\mathrm{\ξ}}}_x=\frac{1}{2}\[{\left({\mathrm{zy}}_{\mathrm{\η}}\right)}_{\mathrm{\ζ}}+{\left({\mathrm{yz}}_{\mathrm{\ζ}}\right)}_{\mathrm{\η}}-{\left({\mathrm{zy}}_{\mathrm{\ζ}}\right)}_{\mathrm{\η}}-{\left({\mathrm{yz}}_{\mathrm{\η}}\right)}_{\mathrm{\ζ}}\]\\ {\widehat{\mathrm{\ξ}}}_y=\frac{1}{2}\[{\left({\mathrm{xz}}_{\mathrm{\η}}\right)}_{\mathrm{\ζ}}+{\left({\mathrm{zx}}_{\mathrm{\ζ}}\right)}_{\mathrm{\η}}-{\left({\mathrm{xz}}_{\mathrm{\ζ}}\right)}_{\mathrm{\η}}-{\left({\mathrm{zx}}_{\mathrm{\η}}\right)}_{\mathrm{\ζ}}\]\\ {\widehat{\mathrm{\ξ}}}_z=\frac{1}{2}\[{\left({\mathrm{yx}}_{\mathrm{\η}}\right)}_{\mathrm{\ζ}}+{\left({\mathrm{xy}}_{\mathrm{\ζ}}\right)}_{\mathrm{\η}}-{\left({\mathrm{yx}}_{\mathrm{\ζ}}\right)}_{\mathrm{\η}}-{\left({\mathrm{xy}}_{\mathrm{\η}}\right)}_{\mathrm{\ζ}}\]\end{array}\right.\\ \left\{\begin{array}{c}{\widehat{\mathrm{\η}}}_x=\frac{1}{2}\[{\left({\mathrm{zy}}_{\mathrm{\ζ}}\right)}_{\mathrm{\ξ}}+{\left({\mathrm{yz}}_{\mathrm{\ξ}}\right)}_{\mathrm{\ζ}}-{\left({\mathrm{zy}}_{\mathrm{\ξ}}\right)}_{\mathrm{\ζ}}-{\left({\mathrm{yz}}_{\mathrm{\ζ}}\right)}_{\mathrm{\ξ}}\]\\ {\widehat{\mathrm{\η}}}_y=\frac{1}{2}\[{\left({\mathrm{xz}}_{\mathrm{\ζ}}\right)}_{\mathrm{\ξ}}+{\left({\mathrm{zx}}_{\mathrm{\ξ}}\right)}_{\mathrm{\ζ}}-{\left({\mathrm{xz}}_{\mathrm{\ξ}}\right)}_{\mathrm{\ζ}}-{\left({\mathrm{zx}}_{\mathrm{\ζ}}\right)}_{\mathrm{\ξ}}\]\\ {\widehat{\mathrm{\η}}}_z=\frac{1}{2}\[{\left({\mathrm{yx}}_{\mathrm{\ζ}}\right)}_{\mathrm{\ξ}}+{\left({\mathrm{xy}}_{\mathrm{\ξ}}\right)}_{\mathrm{\ζ}}-{\left({\mathrm{yx}}_{\mathrm{\ξ}}\right)}_{\mathrm{\ζ}}-{\left({\mathrm{xy}}_{\mathrm{\ζ}}\right)}_{\mathrm{\ξ}}\]\end{array}\right.\\ \left\{\begin{array}{c}{\widehat{\mathrm{\ζ}}}_x=\frac{1}{2}\[{\left({\mathrm{zy}}_{\mathrm{\ξ}}\right)}_{\mathrm{\η}}+{\left({\mathrm{yz}}_{\mathrm{\η}}\right)}_{\mathrm{\ξ}}-{\left({\mathrm{zy}}_{\mathrm{\η}}\right)}_{\mathrm{\ξ}}-{\left({\mathrm{yz}}_{\mathrm{\ξ}}\right)}_{\mathrm{\η}}\]\\ {\widehat{\mathrm{\ζ}}}_y=\frac{1}{2}\[{\left({\mathrm{xz}}_{\mathrm{\ξ}}\right)}_{\mathrm{\η}}+{\left({\mathrm{xz}}_{\mathrm{\η}}\right)}_{\mathrm{\ξ}}-{\left({\mathrm{xz}}_{\mathrm{\η}}\right)}_{\mathrm{\ξ}}-{\left({\mathrm{zx}}_{\mathrm{\ξ}}\right)}_{\mathrm{\η}}\]\\ {\widehat{\mathrm{\ζ}}}_z=\frac{1}{2}\[{\left({\mathrm{yx}}_{\mathrm{\ξ}}\right)}_{\mathrm{\η}}+{\left({\mathrm{xy}}_{\mathrm{\η}}\right)}_{\mathrm{\ξ}}-{\left({\mathrm{yx}}_{\mathrm{\η}}\right)}_{\mathrm{\ξ}}-{\left({\mathrm{xy}}_{\mathrm{\ξ}}\right)}_{\mathrm{\η}}\]\end{array}\right.\end{array}---\left(6\right)$

And the mathematical definition formula of mesh transformations Jacobi is:

$J^{-1}=\left|\begin{array}{ccc}{x}_{\mathrm{\ξ}}& x_{\mathrm{\η}}& x_{\mathrm{\ζ}}\\ y_{\mathrm{\ξ}}& y_{\mathrm{\η}}& y_{\mathrm{\ζ}}\\ z_{\mathrm{\ξ}}& z_{\mathrm{\η}}& z_{\mathrm{\ζ}}\end{array}\right|=x_{\mathrm{\ξ}}{y}_{\mathrm{\η}}{z}_{\mathrm{\ζ}}-x_{\mathrm{\ξ}}{y}_{\mathrm{\ζ}}{z}_{\mathrm{\η}}+x_{\mathrm{\η}}{y}_{\mathrm{\ζ}}{z}_{\mathrm{\ξ}}-x_{\mathrm{\η}}{y}_{\mathrm{\ξ}}{z}_{\mathrm{\ζ}}+x_{\mathrm{\ζ}}{y}_{\mathrm{\ξ}}{z}_{\mathrm{\η}}-x_{\mathrm{\ζ}}{y}_{\mathrm{\η}}{z}_{\mathrm{\ξ}}---\left(7\right)$

The geometry essence of the mathematical definition formula (7) of mesh transformations Jacobi is the volume of computing grid unit.Consider the geometric definition of volume, the volume of any unit is surrounded by its area, uniquely can be determined by its area.The mathematical definition formula (7) of mesh transformations Jacobi is expressed as the form of area:

$\begin{array}{c}{J}^{-1}=x_{\mathrm{\ξ}}{\widehat{\mathrm{\ξ}}}_x+x_{\mathrm{\η}}{\widehat{\mathrm{\η}}}_x+x_{\mathrm{\ζ}}{\widehat{\mathrm{\ζ}}}_x\stackrel{\mathrm{\Δ}}{=}{J}^{-1}\left(x\right)\\ J^{-1}=y_{\mathrm{\ξ}}{\widehat{\mathrm{\ξ}}}_y+y_{\mathrm{\η}}{\widehat{\mathrm{\η}}}_y+y_{\mathrm{\ζ}}{\widehat{\mathrm{\ζ}}}_y\stackrel{\mathrm{\Δ}}{=}{J}^{-1}\left(y\right)\\ J^{-1}=z_{\mathrm{\ξ}}{\widehat{\mathrm{\ξ}}}_z+z_{\mathrm{\η}}{\widehat{\mathrm{\η}}}_z+z_{\mathrm{\ζ}}{\widehat{\mathrm{\ζ}}}_z\stackrel{\mathrm{\Δ}}{=}{J}^{-1}\left(z\right)\end{array}---\left(8\right)$

Formula (8) does not have the geometric meaning of volume when carrying out finite difference and being discrete, formula (8) is rewritten as symmetric form:

$J^{-1}=\frac{1}{3}\left(x_{\mathrm{\ξ}}{\widehat{\mathrm{\ξ}}}_x+y_{\mathrm{\ξ}}{\widehat{\mathrm{\ξ}}}_y+z_{\mathrm{\ξ}}{\widehat{\mathrm{\ξ}}}_z+x_{\mathrm{\η}}{\widehat{\mathrm{\η}}}_x+y_{\mathrm{\η}}{\widehat{\mathrm{\η}}}_y+z_{\mathrm{\η}}{\widehat{\mathrm{\η}}}_z+x_{\mathrm{\ζ}}{\widehat{\mathrm{\ζ}}}_x+y_{\mathrm{\ζ}}{\widehat{\mathrm{\ζ}}}_y+z_{\mathrm{\ζ}}{\widehat{\mathrm{\ζ}}}_z\right)---\left(9\right)$

Formula (8) is rewritten as conservation form:

$\begin{array}{c}{J}^{-1}\left(x\right)={\left(x{\widehat{\mathrm{\ξ}}}_x\right)}_{\mathrm{\ξ}}+{\left(x{\widehat{\mathrm{\η}}}_x\right)}_{\mathrm{\η}}+{\left(x{\widehat{\mathrm{\ζ}}}_x\right)}_{\mathrm{\ζ}}-{\mathrm{xI}}_x\\ J^{-1}\left(y\right)={\left(y{\widehat{\mathrm{\ξ}}}_y\right)}_{\mathrm{\ξ}}+{\left(y{\widehat{\mathrm{\η}}}_y\right)}_{\mathrm{\η}}+{\left(y{\widehat{\mathrm{\ζ}}}_y\right)}_{\mathrm{\ζ}}-{\mathrm{yI}}_y\\ J^{-1}\left(z\right)={\left(z{\widehat{\mathrm{\ξ}}}_z\right)}_{\mathrm{\ξ}}+{\left(z{\widehat{\mathrm{\η}}}_z\right)}_{\mathrm{\η}}+{\left(z{\widehat{\mathrm{\ζ}}}_z\right)}_{\mathrm{\ζ}}-{\mathrm{zI}}_z\end{array}---\left(10\right)$

Wherein I _x, I _yand I _zfor geometry conservation law remainder:

$\begin{array}{c}{I}_x=\frac{\∂{\widehat{\mathrm{\ξ}}}_x}{\∂\mathrm{\ξ}}+\frac{\∂{\widehat{\mathrm{\η}}}_x}{\∂\mathrm{\η}}+\frac{\∂{\widehat{\mathrm{\ζ}}}_x}{\∂\mathrm{\ζ}}\\ I_y=\frac{\∂{\widehat{\mathrm{\ξ}}}_y}{\∂\mathrm{\ξ}}+\frac{\∂{\widehat{\mathrm{\η}}}_y}{\∂\mathrm{\η}}+\frac{\∂{\widehat{\mathrm{\ζ}}}_y}{\∂\mathrm{\ζ}}\\ I_z=\frac{\∂{\widehat{\mathrm{\ξ}}}_z}{\∂\mathrm{\ξ}}+\frac{\∂{\widehat{\mathrm{\η}}}_z}{\∂\mathrm{\η}}+\frac{\∂{\widehat{\mathrm{\ζ}}}_z}{\∂\mathrm{\ζ}}\end{array}---\left(11\right)$

When adopting finite difference method discrete, the internal layer difference operator in the symmetry conservation form (6) of grid derivative is designated as δ ³, outer difference operator is designated as δ ², the difference operator in formula (11) is designated as δ ¹, then the difference operator in the conservation form of calculation (10) of mesh transformations Jacobi and (15) and its symmetry conservation form of calculation (16) is δ ¹.In order to can strict meeting geometric conservation law, δ ¹, δ ²and δ ³linear finite difference operator need be.So-called linear finite difference operator, refers to that difference operator δ has following character:

For arbitrary constant a, b and variable φ, following formula is all had to set up:

And the commutative character of differentiate order, that is:

δ _ξ(δ _ηφ)＝δ _η(δ _ξφ)(13)

The described linear finite difference operator that strictly meets the demands of geometry conservation law meets following constraint condition:

${\mathrm{\δ}}_{\mathrm{\ξ}}^1={\mathrm{\δ}}_{\mathrm{\ξ}}^2={\mathrm{\δ}}_{\mathrm{\ξ}}^3,{\mathrm{\δ}}_{\mathrm{\η}}^1={\mathrm{\δ}}_{\mathrm{\η}}^2={\mathrm{\δ}}_{\mathrm{\η}}^3,{\mathrm{\δ}}_{\mathrm{\ζ}}^1={\mathrm{\δ}}_{\mathrm{\ζ}}^2={\mathrm{\δ}}_{\mathrm{\ζ}}^3---\left(14\right)$

Under the condition that geometry conservation law strictly meets: I _x=I _y=I _z=0, formula (10) is written as:

$\begin{array}{c}{J}^{-1}\left(x\right)={\left(x{\widehat{\mathrm{\ξ}}}_x\right)}_{\mathrm{\ξ}}+{\left(x{\widehat{\mathrm{\η}}}_x\right)}_{\mathrm{\η}}+{\left(x{\widehat{\mathrm{\ζ}}}_x\right)}_{\mathrm{\ζ}}\\ J^{-1}\left(y\right)={\left(y{\widehat{\mathrm{\ξ}}}_y\right)}_{\mathrm{\ξ}}+{\left(y{\widehat{\mathrm{\η}}}_y\right)}_{\mathrm{\η}}+{\left(y{\widehat{\mathrm{\ζ}}}_y\right)}_{\mathrm{\ζ}}\\ J^{-1}\left(z\right)={\left(z{\widehat{\mathrm{\ξ}}}_z\right)}_{\mathrm{\ξ}}+{\left(z{\widehat{\mathrm{\η}}}_z\right)}_{\mathrm{\η}}+{\left(z{\widehat{\mathrm{\ζ}}}_z\right)}_{\mathrm{\ζ}}\end{array}---\left(15\right)$

In like manner, formula (15) does not have the geometric meaning of volume when carrying out finite difference and being discrete, formula (15) is rewritten the symmetry conservation form of mesh transformations Jacobi:

$J^{-1}=\frac{1}{3}\[{\left(x{\widehat{\mathrm{\ξ}}}_x+y{\widehat{\mathrm{\ξ}}}_y+z{\widehat{\mathrm{\ξ}}}_z\right)}_{\mathrm{\ξ}}+{\left(x{\widehat{\mathrm{\η}}}_x+y{\widehat{\mathrm{\η}}}_y+z{\widehat{\mathrm{\η}}}_z\right)}_{\mathrm{\η}}+{\left(x{\widehat{\mathrm{\ζ}}}_x+y{\widehat{\mathrm{\ζ}}}_y+z{\widehat{\mathrm{\ζ}}}_z\right)}_{\mathrm{\ζ}}\]---\left(16\right)$

The geometry essence of mesh transformations Jacobi is the volume of grid cell.The form of calculation different to mesh transformations Jacobi carry out finite difference discrete time, the symmetry conservation form (16) of the mathematical definition formula (7) of mesh transformations Jacobi, the symmetric form (9) of mesh transformations Jacobi, mesh transformations Jacobi all can keep the geometrical property of its grid cell volume, and expression formula (8) and conservation form (15) all can not keep its due geometrical property, this needs formula (8) and formula (15) to write the reason of rewriting its symmetric form (9) and formula (16) respectively just.

When adopting linear finite difference operator discrete, might as well remember that the mesh transformations Jacobi calculated by the mathematical definition formula (7) of mesh transformations Jacobi is V ¹; Remember that the mesh transformations Jacobi calculated by the symmetric form (9) of mesh transformations Jacobi is V ²; Remember that the mesh transformations Jacobi calculated by the symmetry conservation form (16) of mesh transformations Jacobi is V ³.If calculate the grid schematic diagram of employing as shown in Figure 1, then, when adopting Using Second-Order Central Difference form discrete, calculated the mesh transformations Jacobi V of an O by the mathematical definition formula (7) of mesh transformations Jacobi ¹the schematic diagram of net point coordinate information used as shown in Figure 2; The mesh transformations Jacobi V of an O is calculated by the symmetric form (9) of mesh transformations Jacobi ²the schematic diagram of net point coordinate information used is respectively the partial volume addition of three shown in Fig. 3, Fig. 4 and Fig. 5 and obtains; The mesh transformations Jacobi V of an O is calculated by the symmetry conservation form (16) of mesh transformations Jacobi ³six partial volumes of schematic diagram as shown in Fig. 6, Fig. 7, Fig. 8, Fig. 9, Figure 10 and Figure 11 of net point coordinate information used are added and obtain.It should be noted that after adopting Using Second-Order Central Difference form discrete especially, the mesh transformations Jacobi V at some O place ¹be not equal to the volume shown in Fig. 2, but 3/4ths of volume shown in Fig. 2; V ²also three partial volume sums shown in Fig. 3, Fig. 4 and Fig. 5 are not exclusively equaled, but 1/8th of three partial volume sums shown in Fig. 3, Fig. 4 and Fig. 5; In like manner V ³also six partial volume sums shown in Fig. 6 to Figure 11 are not exclusively equaled, but 1/8th of six partial volume sums shown in Fig. 6 to Figure 11.

Easily find out, in Fig. 6 to Figure 11, all six partial volume sums are just for comprising the volume of the macrolattice unit of eight grid cells shown in Fig. 1.Due to V ¹, V ²and V ³mathematical definition is equivalent, the volume of the macrolattice unit that all can be made up of eight grid cells shown in representative graph 1 to a certain extent after its volume calculated is multiplied by octuple, wherein V ¹volume (the V carrying out the macrolattice unit shown in representative graph 1 with six of volume shown in Fig. 2 times is equivalent to after being multiplied by octuple ¹for 3/4ths of volume shown in Fig. 2), and V ²and V ³all can the volume of the macrolattice unit shown in direct representation Fig. 1 after being multiplied by octuple.

As can be seen from the different volumes schematic diagram shown in Fig. 2 to Figure 11, V ¹, V ²and V ³the net point coordinate information adopted when calculating volume is incomplete same, wherein V ¹the net point coordinate information adopted is minimum, V ²secondly, V ³the net point coordinate information adopted is maximum, V ³use the coordinate information of all net points forming macrolattice unit shown in Fig. 1.Therefore V ¹, V ²and V ³between difference just can embody the difference of relative position between net point coordinate.Illustrate: six times of volume shown in Fig. 2 (calculate volume V ¹octuple) with the macrolattice volume shown in Fig. 1 (calculate volume V ³octuple) between difference just in time reflect in Fig. 1 form macrolattice unit eight differences between grid cell and parallelepipedon, that is close all quadrilaterals of macrolattice unit and the difference of parallelogram in Fig. 1.In other words, if eight grid cells forming macrolattice unit shown in Fig. 1 are parallelepipedon, then six times of volume shown in Fig. 2 just equal with the volume of the unit of macrolattice shown in Fig. 1.If eight grid cells forming macrolattice unit shown in Fig. 1 are not all parallelepipedons or are not parallelepipedon entirely, then six times of volume shown in Fig. 2 generally unequal with the volume of the unit of macrolattice shown in Fig. 1.

In like manner, the volume V shown in Fig. 3 to Fig. 5 ²the net point coordinate information adopted compares V ¹many, but compare V ³few, its respectively with V ¹and V ³difference can also reflect respectively in Fig. 1 and surround some face of macrolattice unit and the difference of parallelogram.If all quadrilaterals surrounding macrolattice unit in Fig. 1 are parallelogram, namely calculate the mesh quality fine (optimal computing grid) adopted, then V ¹, V ²and V ³all equal.Given this, the present invention designs following mesh quality Testing index:

$\begin{array}{c}\begin{array}{cc}{V}_N=\min \left(V_{i,j,k}^1,V_{i,j,k}^2,V_{i,j,k}^3\right)& foralli,j,k\end{array}\\ \begin{array}{cc}{V}_M=\max \left(V_{i,j,k}^{12},V_{i,j,k}^{23},V_{i,j,k}^{13}\right)& foralli,j,k\end{array}\\ V^{12}=\frac{\left|V^1-V^2\right|}{V^2}\×100\%,V^{23}=\frac{\left|V^2-V^3\right|}{V^3}\×100\%,V^{13}=\frac{\left|V^1-V^3\right|}{V^3}\×100\%\end{array}---\left(17\right)$

I.e. V _nfor mesh transformations Jacobi three kinds of form of calculation V ¹, V ²and V ³minimum value on all computing grids; V _mfor V ¹², V ²³and V ¹³each mesh transformations Jacobi calculates the maximal value of relative error between component.It should be noted that, although the above-mentioned analysis for mesh transformations Jacobi is all the central difference schemes based on second order accuracy, finite difference scheme for higher order accuracy is applicable equally, and it can be applied in the finite difference method of higher order accuracy equally to the detection of mesh quality.

For arbitrary polylith docking structure grid, first to ensure its volume non-negative, i.e. V _n> 0.If want to obtain good higher order accuracy numerical solution when adopting the numerical method discrete calculation territory of higher order accuracy, then mesh quality also needs to meet another requirement: V _m≤ ε.In fact index V _mreflect the degree of grid cell distortion, be the attribute of adopted computing grid itself, better then its grid cell of mesh quality is more regular, V _mvalue will be less; Index ε reflects the degree that numerical computation method can obtain the grid cell distortion that better numerical simulation result allows, and be the embodiment of adopted numerical calculations ability, ε value is less shows that the requirement of numerical computation method to mess generation quality is higher.According to the experience of practical engineering application, for higher order accuracy (and second order accuracy) numerical computation method, we provide ε=10%.

By detect determine the particular location not reaching Testing index grid cell, again detect after the described grid cell not reaching Testing index is adjusted, by this repeatably process until adjustment after grid cell all reach Testing index; Effectively can ensure that mesh quality meets the requirement of higher order accuracy finite difference method calculating like this.Then when adopting that higher order accuracy finite difference method is discrete to be solved, successfully can obtain good higher order accuracy numerical solution.Logarithm value analog result carries out the ideal higher order accuracy numerical simulation that postpositive disposal can realize physical descriptor Q in opposing connection complex configuration computational fields.

Three kinds of form of calculation of mesh transformations Jacobi all have the different forms of expression, as by V ¹, V ²and V ³calculation expression can be rewritten as the form of vector respectively:

The vector form of the mathematical definition formula (7) of mesh transformations Jacobi:

$J^{-1}=\left({\stackrel{\→}{r}}_{\mathrm{\ξ}}\×{\stackrel{\→}{r}}_{\mathrm{\η}}\right)\·{\stackrel{\→}{r}}_{\mathrm{\ζ}}=\left({\stackrel{\→}{r}}_{\mathrm{\η}}\×{\stackrel{\→}{r}}_{\mathrm{\ζ}}\right)\·{\stackrel{\→}{r}}_{\mathrm{\ξ}}=\left({\stackrel{\→}{r}}_{\mathrm{\ζ}}\×{\stackrel{\→}{r}}_{\mathrm{\ξ}}\right)\·{\stackrel{\→}{r}}_{\mathrm{\η}}---\left(18\right)$

The vector form of the symmetric form (9) of mesh transformations Jacobi:

$J^{-1}=\frac{1}{3}\left({\stackrel{\→}{r}}_{\mathrm{\ξ}}\·{\stackrel{\→}{S}}^{\mathrm{\ξ}}+{\stackrel{\→}{r}}_{\mathrm{\η}}\·{\stackrel{\→}{S}}^{\mathrm{\η}}+{\stackrel{\→}{r}}_{\mathrm{\ζ}}\·{\stackrel{\→}{S}}^{\mathrm{\ζ}}\right)---\left(19\right)$

The vector form of the symmetry conservation form (16) of mesh transformations Jacobi:

$J^{-1}=\frac{1}{3}\[{\left(\stackrel{\→}{r}\·{\stackrel{\→}{S}}^{\mathrm{\ξ}}\right)}_{\mathrm{\ξ}}+{\left(\stackrel{\→}{r}\·{\stackrel{\→}{S}}^{\mathrm{\η}}\right)}_{\mathrm{\η}}+{\left(\stackrel{\→}{r}\·{\stackrel{\→}{S}}^{\mathrm{\ζ}}\right)}_{\mathrm{\ζ}}\]---\left(20\right)$

Wherein:

$\begin{array}{c}{\stackrel{\→}{S}}^{\mathrm{\ξ}}={\widehat{\mathrm{\ξ}}}_x\stackrel{\→}{i}+{\widehat{\mathrm{\ξ}}}_y\stackrel{\→}{j}+{\widehat{\mathrm{\ξ}}}_x\stackrel{\→}{k}\\ {\stackrel{\→}{S}}^{\mathrm{\η}}={\widehat{\mathrm{\η}}}_x\stackrel{\→}{i}+{\widehat{\mathrm{\η}}}_y\stackrel{\→}{j}+{\widehat{\mathrm{\η}}}_z\stackrel{\→}{k}\\ {\stackrel{\→}{S}}^{\mathrm{\ζ}}={\widehat{\mathrm{\ζ}}}_x\stackrel{\→}{i}+{\widehat{\mathrm{\ζ}}}_y\stackrel{\→}{j}+{\widehat{\mathrm{\ζ}}}_z\stackrel{\→}{k}\end{array}$

Obviously, formula (18), formula (19) and formula (20) are respectively the another kind of form of expression of symmetry conservation form (16) of the mathematical definition formula (7) of mesh transformations Jacobi, the symmetric form (9) of mesh transformations Jacobi and mesh transformations Jacobi, its essence is identical with formula (16) with formula (7), formula (9) respectively, these many forms are discrete lower equivalent at described linear finite difference operator, embody the uniqueness that mesh transformations Jacobi calculates; The many forms of mesh transformations Jacobi is convenient to again the application in different physical model on the other hand.

Claims (3)

1. the method that the finite-difference modeling detecting mesh quality is flowed around complex configuration, it comprises:

The generation of step one, computing grid, measuring and adjustation;

Step 2, control differential equation discrete and solving;

Step 3, postpositive disposal;

To it is characterized in that in described step one 1), complicated polylith docking structure grid is generated to whole flowing zoning, 2), the mesh transformations Jacobi that calculates after adopting linear finite difference operator discrete of the mathematical definition formula of mesh transformations Jacobi, symmetric form and symmetry conservation form is designated as V respectively ¹, V ²and V ³, 3), according to described V ¹, V ²and V ³determine the Testing index V of mesh quality _nand V _m, described V _nfor described V ¹, V ²and V ³minimum value; Described V _mfor V ¹², V ²³and V ¹³maximal value, wherein V ¹², V ²³and V ¹³by described V ¹, V ²and V ³by V ¹²=∣ V ¹?V ²∣/V ², V ²³=∣ V ²?V ³∣/V ³, V ¹³=∣ V ¹?V ³∣/V ³determine, 4), according to the Testing index of described mesh quality, the computing grid generated is detected, determine that those do not reach the particular location of the grid cell of Testing index, 5), the described grid cell not reaching Testing index is adjusted after again detect, by this repeatably process until the grid cell after this adjustment all reaches Testing index; In described step 2, the discrete control differential equation of higher order accuracy finite difference method is adopted to whole flowing zoning, and it is solved; In described step 3, logarithm value simulated flow pattern analyzes the physical significance of its correspondence.

2. the method that flows around complex configuration of the finite-difference modeling of detection mesh quality according to claim 1, when it is characterized in that adopting higher order accuracy finite difference method discrete, the Testing index of described mesh quality is for meet V simultaneously _n> 0 and V _m≤ 10%.

3. the method that flows around complex configuration of the finite-difference modeling of detection mesh quality according to claim 1 and 2, it is characterized in that three kinds of form of calculation of described mesh transformations Jacobi all have many forms, these many forms are discrete lower equivalent at described linear finite difference operator.