CN102831304B - A kind of method using the finite-difference modeling of geometry conservation to flow around complex configuration - Google Patents

Abstract

The invention discloses a kind of method using the finite-difference modeling of geometry conservation to flow around complex configuration, it comprises: the generation of step one, computing grid; Step 2, control differential equation discrete and solving; Step 3, postpositive disposal; Mainly need to meet in described step 2: the symmetry conservation form of calculation 1, adopting grid derivative, make the calculating of grid derivative in the geometry essence ensureing embody again while geometry conservation law strictly meets its grid cell vector area, 2, the symmetry conservation form of calculation of mesh transformations Jacobi is adopted, geometry that the volume making discrete mesh transformations Jacobi can embody grid cell is surrounded by its vector area essence, 3, the space difference operator of discrete differential equation differentiate is designated as δ ¹, the outer difference operator of conservation grid derivative is designated as δ ², the internal layer difference operator of conservation grid derivative is designated as δ ³, then three described difference operators are linear difference operator, and they are equal respectively on each coordinates computed direction.The present invention greatly strengthen higher order accuracy finite difference method in the computational stability in complex configuration flowing and adaptive faculty.

Description

A kind of method using the finite-difference modeling of geometry conservation to flow 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 finite-difference modeling method of the higher order accuracy around complex configuration flowing.

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 configuration such as aircraft, guided missile or aircraft under water around complex configuration flowing.

In existing numerical discretization schemes, the Finite Volume Method of second order accuracy obtains fairly large application because of its good computational stability; And ancient finite difference method due to computational stability bad, the computational stability in the structured grid especially docked at complicated polylith is poor, and limits its large-scale application in Practical Project problem; The excessive large-scale application being also unfavorable for Practical Project problem of calculated amount of Finite Element Method.If but consider the implementation procedure of high precision in three dimensions, high resolving power numerical computation method, finite difference method reaches high precision, high resolving power owing to can adopt the mode by dimension discreet derivative, and Finite Volume Method and Finite Element Method just can reach high precision because needs adopt multidimensional to reconstruct, its calculated amount magnitude larger than the finite difference method of same precision is even more.Under existence conditions, consider the restrictive condition of computational resource, become the best-of-breed technology approach solved around complex configuration flow numerical simulation in Practical Project problem based on the high order accurate numerical method that finite difference is discrete.

The stability calculated and adaptive faculty refer to the dependency degree of method for numerical simulation to grid enable quality.The computational stability of Finite Volume Method is good, its main cause is that Finite Volume Method is when computing grid derivative, grid surface sum volume is all calculate based on the geometric meaning of computing grid, can strict meeting geometric conservation law (its vector area of the arbitrary closed cell cube of meaning and space and should be zero, then show as the physical characteristics that governing equation freely flows conservation in the calculation), the geometrical property of computing grid can be reflected again preferably simultaneously.And the poor major reason of computational stability in the structured grid that finite difference method docks at complicated polylith is just that its geometry conservation law is not easy to meet, especially its geometry conservation law of finite difference method of higher order accuracy is more difficult meets.Although the people such as Deng little Gang are at " Journal of ComputationalPhysics " Volume 230.Issue 4, 20February 2011, the conservation grid derivative method that Pages 1100-1115 proposes can strict meeting geometric conservation law, strengthen High Resolution Finite Difference method to the adaptive faculty of complicated polylith docking structure grid, but still there is following problem: 1, conservation grid derivative method only discuss a kind of conservation grid derivative of form, and do not possess the geometrical property of the computing grid vector area that grid derivative should possess, 2, the form of calculation of grid derivative fails in fact uniquely to determine, the form of calculation of different grid derivatives is caused also to be difficult to assessment to the impact of result of calculation, 3, in conservation grid derivative method, fully strong proof is not provided to choosing of conservation grid derivative internal layer difference operator in grid derivative calculations process, there is ambiguity, easily misunderstand, cause computation process can not by correct enforcement, 4, the form of calculation of conservation grid derivative method to mesh transformations Jacobi is not mentioned, can there is the mesh transformations Jacobi of various ways in actual computation process, and it is also difficult to pass judgment on the impact of result of calculation.Above problem causes the simulation of current High Resolution Finite Difference method still computational stability in urgent need to be improved and adaptive faculty when complex configuration flows.

Summary of the invention

The object of the invention is to solve the computational stability and adaptive faculty problem of using the simulation of High Resolution Finite Difference method when complex configuration flows.

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

Use the method that the finite-difference modeling of geometry conservation is flowed around complex configuration, it comprises:

The generation of step one, computing grid;

Step 2, control differential equation discrete and solving;

Step 3, postpositive disposal;

It is characterized in that in described step one, complicated polylith docking structure grid is generated to whole flowing zoning; In described step 21, adopt the symmetry conservation form of calculation of grid derivative, make the calculating of grid derivative in the geometry essence ensureing embody again while geometry conservation law strictly meets its vector area, 2, the symmetry conservation form of calculation of mesh transformations Jacobi is adopted, make discrete after the mesh transformations Jacobi volume that can the embody grid cell geometry essence of being surrounded by its vector area, 3, the space difference operator of discrete differential equation differentiate is designated as δ ¹, the outer difference operator of conservation grid derivative is designated as δ ², the internal layer difference operator of conservation grid derivative is designated as δ ³, then three described difference operators are linear difference operator, and they are equal respectively on each coordinates computed direction; 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 the calculating of grid derivative have employed its symmetry conservation form, so can ensure that the strict of geometry conservation law meets;

2), after adopting linear finite difference operator discrete in polylith docking structure grid, grid derivative can be expressed as the form of the vector area of grid cell or the linear combination of grid cell vector area;

3), after adopting linear finite difference operator discrete in polylith docking structure grid, mesh transformations Jacobi can be expressed as the linear combination of grid cell volume or grid cell volume;

4), after adopting linear finite difference operator discrete in polylith docking structure grid, mesh transformations Jacobi is surrounded by its vector area and is also drawn by its vector areal calculation;

5) because grid derivative have employed the form of differentiating operator, so can directly adopt high precision, high-resolution finite difference operator to obtain high precision, high-resolution numerical simulation result by the differentiating operator in dimension discrete grid block derivative;

6) because mesh transformations Jacobi have employed the form of differentiating operator, so can directly adopt high precision, high-resolution finite difference operator to obtain high precision, high-resolution numerical simulation result by the differentiating operator in dimension discrete grid block conversion Jacobi;

7) when less calculated amount, high precision, high resolving power finite difference method is substantially increased in the stability of numerical simulation and adaptive faculty in complex configuration flowing.

In addition, there is many forms in the symmetry conservation form of calculation of grid derivative of the present invention, and these many forms are discrete lower equivalent at described linear difference operator.

There is many forms in the symmetry conservation form of calculation of mesh transformations Jacobi of the present invention, these many forms are discrete lower equivalent at described linear difference operator.

Accompanying drawing explanation

Fig. 1 is the computing grid schematic diagram of the grid derivative of embodiment

Fig. 2 is the schematic diagram (Using Second-Order Central Difference form) of the vector area that the symmetry conservation grid derivative (8) of embodiment calculates

Fig. 3 is the schematic diagram (Using Second-Order Central Difference form) of the vector area that the original mesh derivative (5) of embodiment calculates

Fig. 4 is the schematic diagram (Using Second-Order Central Difference form) of the vector area that class Finite Volume Method calculates

Fig. 5 is the schematic diagram (before single order difference form) of the vector area that the symmetry conservation grid derivative (8) of embodiment calculates

Fig. 6 is the schematic diagram (after single order difference form) of the vector area that the symmetry conservation grid derivative (8) of embodiment calculates

Fig. 7 be the vector area that the symmetry conservation grid derivative (8) of embodiment calculates schematic diagram ( direction: difference form before single order; direction: Using Second-Order Central Difference form)

Fig. 8 be the vector area that the symmetry conservation grid derivative (8) of embodiment calculates schematic diagram ( direction: Using Second-Order Central Difference form; direction: difference form after single order)

Fig. 9 is the volume schematic diagram that 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:

$\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---\left(4\right)$

$\hat{G}={\widehat{\mathrm{\ζ}}}_{t}Q+{\widehat{\mathrm{\ζ}}}_{x}E+{\widehat{\mathrm{\ζ}}}_{y}F+{\widehat{\mathrm{\ζ}}}_{z}G$

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

$\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(5\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.$

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

And the mathematical definition formula of mesh transformations Jacobi is:

J ^-1＝x _ξy _ηz _ζ-x _ξy _ζz _η+x _ηy _ζz _ξ-x _ηy _ξz _ζ+x _ζy _ξz _η-x _ζy _ηz _ξ(6)

The I that strictly meets the demands of geometry conservation law _x=I _y=I _z=0, wherein:

$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{\ζ}}---\left(7\right)$

$I_z=\frac{\∂{\widehat{\mathrm{\ξ}}}_z}{\∂\mathrm{\ξ}}+\frac{\∂{\widehat{\mathrm{\η}}}_z}{\∂\mathrm{\η}}+\frac{\∂{\widehat{\mathrm{\ζ}}}_z}{\∂\mathrm{\ζ}}$

In order to can strict meeting geometric conservation law, can reflect again the geometrical property of computing grid exactly, the present invention provides the form of calculation of the symmetry conservation of following grid derivative simultaneously:

$\left\{\begin{array}{c}{\widehat{\mathrm{\ξ}}}_x=\frac{1}{2}[{\left({\mathrm{zy}}_{\mathrm{\η}}\right)}_{\mathrm{\ζ}}+{\left(y{z}_{\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(y{z}_{\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(y{z}_{\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(8\right)$

As can be seen from the mathematical definition formula (5) of grid derivative, discrete grid leads the vector area that geometry of numbers essence is grid cell, and each geometric sense in formula (5) is real is the component of grid cell vector area under rectangular coordinate system.The form of calculation (8) of the grid derivative provided in the present invention can also keep the geometrical property of its grid cell vector area when adopting suitable linear finite difference operator discrete.

The geometry essence of the mathematical definition formula (6) of mesh transformations Jacobi is the volume of grid cell.Consider the geometric definition of volume, the volume of any unit is surrounded by its vector area, also uniquely can be determined by its vector area.Given this, the symmetry conservation form of calculation of mesh transformations Jacobi that provides of the present invention:

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

The geometry essence of mesh transformations Jacobi is the volume of grid cell.When adopting linear finite difference operator discrete to the symmetry conservation form of calculation of mesh transformations Jacobi, if adopt the grid derivative of symmetry conservation simultaneously, then the volume that formula (9) calculates is just the whole volume in Fig. 9.

In the present invention, the finite difference operator δ of indication is only for linear difference operator, and linear difference operator δ has following characteristic:

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

Further derivation can prove: linear difference operator also meets the commutative character of differentiate order, that is:

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

The outer difference operator of the grid derivative in symmetry conservation form (8) is designated as δ ², internal layer difference operator is designated as δ ³, the differential derivative in formula (8) is expressed as the form of difference operator:

$\left\{\begin{array}{c}{\widehat{\mathrm{\ξ}}}_x=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\η}}^{3}y\right)+{\mathrm{\δ}}_{\mathrm{\η}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}z\right)-{\mathrm{\δ}}_{\mathrm{\η}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}y\right)-{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\η}}^{3}z\right)]\\ {\widehat{\mathrm{\ξ}}}_y=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\η}}^{3}z\right)+{\mathrm{\δ}}_{\mathrm{\η}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}x\right)-{\mathrm{\δ}}_{\mathrm{\η}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}z\right)-{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\η}}^{3}x\right)]\\ {\widehat{\mathrm{\ξ}}}_z=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\η}}^{3}x\right)+{\mathrm{\δ}}_{\mathrm{\η}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}y\right)-{\mathrm{\δ}}_{\mathrm{\η}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}x\right)-{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\η}}^{3}y\right)]\end{array}\right.$

$\left\{\begin{array}{c}{\widehat{\mathrm{\η}}}_x=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}y\right)+{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}z\right)-{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}y\right)-{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}z\right)]\\ {\widehat{\mathrm{\η}}}_y=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}z\right)+{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}x\right)-{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}z\right)-{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}x\right)]\\ {\widehat{\mathrm{\η}}}_z=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}x\right)+{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}y\right)-{\mathrm{\δ}}_{\mathrm{\ζ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}x\right)-{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ζ}}^{3}y\right)]\end{array}\right.$

$\left\{\begin{array}{c}{\widehat{\mathrm{\ζ}}}_x=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\η}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}y\right)+{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\η}}^{3}z\right)-{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\η}}^{3}y\right)-{\mathrm{\δ}}_{\mathrm{\η}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}z\right)]\\ {\widehat{\mathrm{\ζ}}}_y=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\η}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}z\right)+{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(z{\mathrm{\δ}}_{\mathrm{\η}}^{3}x\right)-{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\η}}^{3}z\right)-{\mathrm{\δ}}_{\mathrm{\η}}^2\left(z{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}x\right)]\\ {\widehat{\mathrm{\ζ}}}_z=\frac{1}{2}[{\mathrm{\δ}}_{\mathrm{\η}}^2\left(y{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}x\right)+{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(x{\mathrm{\δ}}_{\mathrm{\η}}^{3}y\right)-{\mathrm{\δ}}_{\mathrm{\ξ}}^2\left(y{\mathrm{\δ}}_{\mathrm{\η}}^{3}x\right)-{\mathrm{\δ}}_{\mathrm{\η}}^2\left(x{\mathrm{\δ}}_{\mathrm{\ξ}}^{3}y\right)]\end{array}\right.---\left(12\right)$

In the present invention, when difference operator meets following condition:

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

Then formula (12) always can be expressed as the form of grid cell vector area after linear finite difference operator is discrete.Illustrate for the computing grid schematic diagram shown in Fig. 1: if δ ²and δ ³be Using Second-Order Central Difference at coordinates computed system ξ and η direction, then the vector area that calculates of formula (12) as shown in Figure 2 (the vector area in Fig. 2 to Fig. 8 all refers to the area of dash area).Adopt identical difference operator, then the vector area that calculates of the mathematical definition formula (5) of grid derivative as shown in Figure 3.And the vector area that class Finite Volume Method provides as shown in Figure 4.As difference operator δ ²and δ ³not during difference, as δ ²and δ ³when coordinates computed system ξ and η direction are single order forward difference form, then the vector area that calculates of formula (12) as shown in Figure 5; In like manner Fig. 6 is then δ ²and δ ³the vector area of single order backward difference form calculating is at coordinates computed system ξ and η direction.Notice in formula (13) and do not require that the difference operator in different coordinates computed systems direction wants the same, when the difference operator of two coordinate directions is different, as δ ²and δ ³be single order form forward in coordinates computed system ξ direction, and be Using Second-Order Central Difference form in coordinates computed system η direction, as shown in Figure 7, Fig. 8 then represents the vector area in another situation to its vector area: δ ²and δ ³be Using Second-Order Central Difference form in coordinates computed system ξ direction, and be single order backward difference form in coordinates computed system η direction.

When carrying out finite difference to equation (3) and being discrete, the space difference operator in equation (3) is designated as δ ¹, then have and freely flow conservation (geometry conservation law) requirement:

$I_x={\mathrm{\δ}}_{\mathrm{\ξ}}^1\left({\widehat{\mathrm{\ξ}}}_x\right)+{\mathrm{\δ}}_{\mathrm{\η}}^1\left({\widehat{\mathrm{\η}}}_x\right)+{\mathrm{\δ}}_{\mathrm{\ζ}}^1\left({\widehat{\mathrm{\ζ}}}_x\right)$

$I_y={\mathrm{\δ}}_{\mathrm{\ξ}}^1\left({\widehat{\mathrm{\ξ}}}_y\right)+{\mathrm{\δ}}_{\mathrm{\η}}^1\left({\widehat{\mathrm{\η}}}_y\right)+{\mathrm{\δ}}_{\mathrm{\ζ}}^1\left({\widehat{\mathrm{\ζ}}}_y\right)---\left(14\right)$

$I_z={\mathrm{\δ}}_{\mathrm{\ξ}}^1\left({\widehat{\mathrm{\ξ}}}_z\right)+{\mathrm{\δ}}_{\mathrm{\η}}^1\left({\widehat{\mathrm{\η}}}_z\right)+{\mathrm{\δ}}_{\mathrm{\ζ}}^1\left({\widehat{\mathrm{\ζ}}}_z\right)$

If the grid derivative in formula (14) adopts such as formula the conservation form of calculation in (8), then according to the inherent characteristic (11) of linear finite difference operator, the adequate condition easily providing strict meeting geometric conservation law is:

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

Convolution (13) and formula (15), the linear finite difference operator that the present invention provides should meet:

${\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(16\right)$

Can find out that the symmetry conservation form of grid derivative is the vector area of grid cell after discrete by Fig. 2 and Fig. 9, mesh transformations Jacobi is surrounded by its vector area completely and is also drawn by its vector areal calculation, grid derivative and mesh transformations Jacobi geometry essence is separately fully reflected, and is conducive to the numerical evaluation stability and the adaptive faculty that strengthen finite difference method.When less calculated amount, substantially increase high precision, high resolving power finite difference method in the stability of numerical simulation and adaptive faculty in complex configuration flowing.Logarithm value analog result is carried out postpositive disposal and can be realized the good numerical simulation of physical descriptor Q in calculating basin.

The symmetry conservation form of calculation (8) of grid derivative can also simply be written as:

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

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

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

And the vector calculus form of its correspondence:

${\stackrel{\→}{S}}^{\mathrm{\ξ}}={\widehat{\mathrm{\ξ}}}_x\stackrel{\→}{i}+{\widehat{\mathrm{\ξ}}}_y\stackrel{\→}{j}+{\widehat{\mathrm{\ξ}}}_x\stackrel{\→}{k}=\frac{1}{2}[{({\stackrel{\→}{r}}_{\mathrm{\η}}\×\stackrel{\→}{r})}_{\mathrm{\ζ}}-{({\stackrel{\→}{r}}_{\mathrm{\ζ}}\×\stackrel{\→}{r})}_{\mathrm{\η}}]$

${\stackrel{\→}{S}}^{\mathrm{\η}}={\widehat{\mathrm{\η}}}_x\stackrel{\→}{i}+{\widehat{\mathrm{\η}}}_y\stackrel{\→}{j}+{\widehat{\mathrm{\η}}}_z\stackrel{\→}{k}=\frac{1}{2}[{({\stackrel{\→}{r}}_{\mathrm{\ζ}}\×\stackrel{\→}{r})}_{\mathrm{\ξ}}-{({\stackrel{\→}{r}}_{\mathrm{\ξ}}\×\stackrel{\→}{r})}_{\mathrm{\ζ}}]---\left(18\right)$

${\stackrel{\→}{S}}^{\mathrm{\ζ}}={\widehat{\mathrm{\ζ}}}_x\stackrel{\→}{i}+{\widehat{\mathrm{\ζ}}}_y\stackrel{\→}{j}+{\widehat{\mathrm{\ζ}}}_z\stackrel{\→}{k}=\frac{1}{2}[{({\stackrel{\→}{r}}_{\mathrm{\ξ}}\×\stackrel{\→}{r})}_{\mathrm{\η}}-{({\stackrel{\→}{r}}_{\mathrm{\η}}\×\stackrel{\→}{r})}_{\mathrm{\ξ}}]$

Obvious formula (17) and (18) to be all equal to the discrete form (12) of formula (8) afterwards the above-mentioned linear difference operator of employing is discrete.In like manner, formula (19) is equal to after adopting above-mentioned linear difference operator discrete with formula (9) mesh transformations Jacobi.

Therefore, in the present invention, can say that the symmetry conservation form of calculation of grid derivative is in fact well-determined, the symmetry conservation form of mesh transformations Jacobi is in fact also well-determined, so just substantially increases stability and the adaptive faculty of numerical simulation calculation; There is the application that many forms is conducive in different physical model in the symmetry conservation form of calculation of grid derivative and the symmetry conservation form of calculation of mesh transformations Jacobi.

Uniquely determining the form of calculation of grid derivative and mesh transformations Jacobi, after eliminating the uncertainty of grid derivative in computation process and mesh transformations Jacobi, logarithm value analog result carries out the ideal numerical simulation that postpositive disposal can realize calculating physical descriptor Q in basin.

Claims (3)

1. use the method that the finite-difference modeling of geometry conservation is flowed around complex configuration, it comprises:

The generation of step one, computing grid;

Step 2, control differential equation discrete and solving;

Step 3, postpositive disposal;

It is characterized in that in described step one, complicated polylith docking structure grid is generated to whole flowing zoning; In described step 21, adopt the symmetry conservation form of calculation of grid derivative, make the calculating of grid derivative in the geometry essence ensureing embody again while geometry conservation law strictly meets its vector area, 2, the symmetry conservation form of calculation of mesh transformations Jacobi is adopted, make discrete after the mesh transformations Jacobi volume that can the embody grid cell geometry essence of being surrounded by its vector area, 3, the space difference operator of discrete differential equation differentiate is designated as δ ¹, the outer difference operator of conservation grid derivative is designated as δ ², the internal layer difference operator of conservation grid derivative is designated as δ ³, then three described difference operators are linear difference operator, and they are equal respectively on each coordinates computed direction; 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 utilization geometry conservation according to claim 1, it is characterized in that the symmetry conservation form of calculation of described grid derivative exists many forms, these many forms are discrete lower equivalent at described linear difference operator.

3. the method that flows around complex configuration of the finite-difference modeling of utilization geometry conservation according to claim 1, it is characterized in that the symmetry conservation form of calculation of described mesh transformations Jacobi exists many forms, these many forms are discrete lower equivalent at described linear difference operator.