CN104268322A - Boundary processing technology of WENO difference method - Google Patents

Abstract

The invention provides a boundary processing technology when the flow field calculation is performed by a WENO difference method. The consistent high accuracy can be achieved when the flow field calculation is performed and the complicated boundary can be processed. The boundary processing technology mainly comprises step 1, grid division, wherein grids are divided into an internal area and a boundary area, the boundary area only comprises a layer of grids, and accordingly the complicated calculation of the boundary area DG is minimized; step 2, space dispersion, wherein the space dispersion is performed on the internal area by the WENO finite difference method; step 3, boundary processing, wherein the space dispersion is performed on the boundary area by a DG method; step 4, coupling processing, wherein the coupling is performed on the boundary processing portion and the internal area and the derivative information is provided by the boundary area through DG polynomial approximation; step 5, time dispersion, wherein the time dispersion is displayed by a TVD-Runge-Kutta method; step 6, post-processing, wherein the visual stimulation is performed on a flow field result through commercial software such as Tecplot or Paraview.

Description

A kind of boundary treatment technology of WENO difference method

Technical field

The present invention relates to Fluid Mechanics Computation numerical method field, particularly relate to a kind of boundary treatment technology solving the High Resolution Finite Difference method of Hyperbolic Conservation equation.

Background technology

The numerical simulation of aircraft three-dimensional Complex Flows and relevant multi-objective optimization question are the forward position hot issues in current Fluid Mechanics Computation, are also the application problems of an Engineering Oriented actual demand simultaneously.But at current computer size with under the condition solving ability, the numerical method of current main flow can't meet the needs of this practical implementation problem in counting yield, one of key addressed this problem improves the efficiency of fluid diagnosis.

There is the coupling process of the multiple advantage in order to make full use of popular high order accurate numerical method at present, such as WENO limiter Discontinuous Finite Element Method, the SV method of high-order, the DG/FV method of coupling, DG and the WENO method of multizone coupling, wherein DG and the WENO method of multizone coupling has efficiently, high precision and the advantage of easy process complex boundary, but most effective in current all higher-order methods is still finite difference method, but the method generally all needs to calculate on structured grid, and can not complex boundary be processed, compactness is not had during process border yet, thus lack a kind of can help finite difference method compact process complex boundary higher order values method.

Summary of the invention

In order to overcome the deficiencies in the prior art, we have proposed a kind of a kind of boundary treatment technology of WNEO finite difference method as process WENO difference method that can process complex boundary, the method has developed DG and the WENO-FD method of multizone coupling, make use of the advantage of finite difference method and the advantage of DG more fully, almost whole region all utilizes WENO-FD method to calculate, thus efficiency is higher, only use DG method at one deck grid place on border, boundary treatment utilizes the process of DG method, thus the process complex boundary that can compact.

Technical scheme of the present invention is:

Solved problem is carried out structured grid subdivision, then institute domain is divided into two kinds of unit areas, one is internal element region, and another kind is boundary element region, and boundary element region only has one deck computing grid.

Utilize traditional WENO-FD method to carry out spatial spreading to inner unit area, utilize traditional DG method to carry out spatial spreading to unit area, border, and process border.

The multizone coupling scheme be coupled in DG and WENO method are utilized to realize the coupling processing of WENO-FD method and DG method.Need the numerical flux constructing internal element zone boundary place in coupling process, the building method of WENO numerical flux uses special HWENO construction process to realize.

DG method in boundary element place calculates the n rank approximation by polynomi-als result that this unit is separated, and can obtain the value at this unit center point place and be no more than the derivative value on n rank by this approximation by polynomi-als result.

After spatial spreading terminates, the TVD Runge-Kutta method of high-order is utilized to solve the numerical result obtaining subsequent time.

The result of the visual softwares such as Tecplot to checking examples such as two-dimentional scalar Hyperbolic Conservation equation and Two-dimensional Euler Equations groups is utilized to carry out visualization processing.

The invention has the beneficial effects as follows:

The present invention takes full advantage of the efficient simple advantage of High Order WENO finite difference method, combine DG method and can process complex boundary and the feature with locality and compactness, develop multizone coupling DG and WENO method, propose a kind of WNEO-FD method based on DG process border, the method only uses DG method at one deck grid place, border, process extensive problem be that efficiency is almost equivalent to WENO-FD method, and there is the ability of process complex boundary, and no longer need to construct multiple virtual net lattice point when processing complex boundary, template size when reducing WENO-FD boundary treatment thus process time there is compactness and locality.

Accompanying drawing explanation

Fig. 1 is one dimension example stress and strain model schematic diagram

Fig. 2 is internal element place tradition WENO template schematic diagram

Fig. 3 is internal element boundary HWENO template schematic diagram

Fig. 4 is one dimension scalar Burgers equation example schematic diagram

Fig. 5 is one dimension precision test form schematic diagram

Fig. 6 is that two-dimentional double Mach reflection example structured grid divides schematic diagram

Fig. 7 is double Mach reflection example description figure

Fig. 8 is two-dimentional double Mach reflection example density isoline schematic diagram

Fig. 9 is algorithm realization process flow diagram flow chart

Embodiment

First we illustrate our algorithm main process for one dimension scalar Burgers equation, and the Region dividing as Fig. 1 first step is two parts, and wherein WENO-FD part in the left side is internal element region, and the right one deck DG grid is one deck boundary element region.

The specific implementation process of interior zone WENO-FD is as follows:

First we do not consider thinking position, suppose f'(u) > 0, so at unit I _ithe finite difference scheme form of conservation is

$\frac{d\left(u_i\left(t\right)\right)}{\mathrm{dt}}+\frac{1}{\mathrm{\Δx}}({\hat{f}}_{i+\frac{1}{2}}-{\hat{f}}_{i-\frac{1}{2}})$

Wherein u _i(t)=u (x _i, t) be unit I _ithe point value of midpoint, for numerical flux can be expressed as

$\hat{f}=f(u_{i-r},...,u_{i+s})$

If Fig. 2 is except the flux of two, internal element border, most of internal element flux we select the template S=(S of Fig. 2 ₁, S ₂, S ₃), flux $\hat{f}=f\left(u_{i+\frac{1}{2}}\right),$ Wherein $u_{i+\frac{1}{2}}=w_1{u}_{i+\frac{1}{2}}^1+w_2{u}_{i+\frac{1}{2}}^2+w_3{u}_{i+\frac{1}{2}}^3,$ Then $u_{i+\frac{1}{2}}^s,s=\mathrm{1,2,3}$ Be obtained by Lagrange's interpolation by the value that each little template comprises unit, the little template of the WENO-FD for 5rd is drawn by 3rd Lagrange's interpolation, and final interpolation result is:

$\left\{\begin{array}{c}{u}_{i+\frac{1}{2}}^1=-\frac{1}{6}\left(u_{i-1}\right)+\frac{5}{6}\left(u_i\right)+\frac{1}{3}\left(u_{i+1}\right)\\ u_{i+\frac{1}{2}}^2=\frac{1}{3}\left(u_i\right)+\frac{5}{6}\left(u_{i+1}\right)-\frac{1}{6}\left(u_{i+2}\right)\\ u_{i+\frac{1}{2}}^3=\frac{11}{6}\left(u_{i+1}\right)-\frac{7}{6}\left(u_{i+2}\right)+\frac{1}{3}\left(u_{i+3}\right)\end{array}\right.$

Nonlinear weight w ₁, w ₂, w ₃computation process in meet w _r=d _r+ O (Δ x ³):

$w_r=\frac{{\mathrm{\α}}_r}{\underset{s=1}{\overset{3}{\mathrm{\Σ}}}{\mathrm{\α}}_s},{\mathrm{\α}}_r=\frac{d_r}{{(\mathrm{\ϵ}+{\mathrm{\β}}_r)}^2},r=\mathrm{1,2},3$

Wherein β _rbe referred to as " smoothing factor ", select smoothing factor here:

${\mathrm{\β}}_r=\underset{l=1}{\overset{2}{\mathrm{\Σ}}}{\∫}_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}\mathrm{\Δ}{x}^{2l-1}{\left(\frac{{\∂}^l{p}_r\left(x\right)}{{\∂}_{l}x}\right)}^2\mathrm{dx}.$

The computation process of nonlinear weight needs to calculate linear power: $u_{i+\frac{1}{2}}=d_1{u}_{i+\frac{1}{2}}^1+d_2{u}_{i+\frac{1}{2}}^2+d_3{u}_{i+\frac{1}{2}}^3=u\left(x_{i+\frac{1}{2}}\right)+O\left({\mathrm{\Δx}}^5\right),$ Such linear coefficient $d_1=\frac{3}{10},d_2=\frac{3}{5},d_3=\frac{1}{10}.$

As Fig. 3, if we will calculate two elementary boundary flux on internal element border, less template is used in our method, virtual net lattice point in figure we no longer need, but adopt special Hermite interpolation to obtain the interpolation result of each little template that we want in these two relatively little templates the result of Hermite interpolation is:

$\left\{\begin{array}{c}{u}_{N-\frac{3}{2}}^1=-\frac{1}{6}\left(u_{N-3}\right)+\frac{5}{6}\left(u_{N-2}\right)+\frac{1}{3}\left(u_{N-1}\right)\\ u_{N-\frac{3}{2}}^2=\frac{1}{3}\left(u_{N-2}\right)+\frac{5}{6}\left(u_{N-1}\right)-\frac{1}{6}\left(u_N\right)\\ u_{N-\frac{3}{2}}^3=\frac{13}{6}\left(u_{N-1}\right)-\frac{7}{6}\left(u_N\right)+\frac{2}{3}\mathrm{\Δx}\left(u_N^{\′}\right)\end{array}\right.$

With

$\left\{\begin{array}{c}{u}_{N-\frac{1}{2}}^1=-\frac{1}{6}\left(u_{N-2}\right)+\frac{5}{6}\left(u_{N-1}\right)+\frac{1}{3}\left(u_N\right)\\ u_{N-\frac{1}{2}}^2=\frac{1}{6}\left(u_{N-1}\right)+\frac{5}{6}\left(u_N\right)-\frac{1}{3}\left(u_N\right)\\ u_{N-\frac{1}{2}}^3=\left(u_N\right)-\frac{1}{2}\mathrm{\Δx}\left(u_N^{\′}\right)+\frac{1}{12}\mathrm{\Δ}{x}^2\left(u_N^{\′\′}\right)\end{array}\right.$

And then obtained nonlinear weight and the smoothing factor of our needs by same mode, thus complete solution can be obtained by solving equation.

According to the step in technical scheme, we also need to be obtained as Fig. 3 boundary element I by DG method ₁andI _nthe value of place's unit midpoint and higher derivative value, such as I _nthe approximation by polynomi-als solution that unit obtains is:

u(x)＝u ₀p ₀(x)+u ₁p ₁(x)+u ₂p ₂(x)

Coefficient u ₀, u ₁, u ₂approach by DG method the degree of freedom obtained.Higher derivative information also draws easily via polynomial derivation

u′(x)＝u ₀p′ ₀(x)+u ₁p′ ₁(x)+u ₂p′ ₂(x)

u″(x)＝u ₀p″ ₀(x)+u ₁p″ ₁(x)+u ₂p″ ₂(x)

According to the process in technical scheme, we need the border of two kinds of methods to be coupled to together, DG and HWENO coupling process is as follows:

$\left\{\begin{array}{c}\frac{d}{\mathrm{dt}}{u}_1-\frac{1}{\mathrm{\Δ}{x}_1}({\tilde{f}}_{\frac{1}{2}}-{\tilde{f}}_{i+\frac{1}{2}})=0,\\ \frac{d}{\mathrm{dt}}{u}_i-\frac{1}{\mathrm{\Δ}{x}_i}({\hat{f}}_{i-\frac{1}{2}}-{\hat{f}}_{i+\frac{1}{2}})=0,i=2,...,N-1,\\ \frac{d}{\mathrm{dt}}{u}_N-\frac{1}{\mathrm{\Δ}{x}_N}({\tilde{f}}_{N-\frac{1}{2}}-{\tilde{f}}_{N+\frac{1}{2}})=0,\end{array}\right.$

Wherein representative reconstructs by common WENO the flux obtained, represent the flux that DG method obtains, that the mode reconstructed by above-mentioned HWENO is obtained.

In order to verify above-mentioned algorithm, as Fig. 4, we calculate intense shock wave example for the Burgers equation of one dimension scalar, and the one dimension example of our linear has carried out the precision test of algorithm as Fig. 3.

In order to verify above-mentioned algorithm, we verify for two-dimentional double Mach reflection example, and the mesh generation of double Mach reflection problem is as Fig. 6, and the meshing wherein represents internal element grid, and bottom one deck grid represents border unit grid.

As Fig. 7, suppose initial time, the normal shock wave of Mach number Ma=10 is to treadmill exercise.60 ° of angles are become between shock wave with flat board.By the evolution of the WENO finite difference method computational flow on the discontinuous Galerkin process border of foregoing description.

As Fig. 3 and Fig. 7, we are according to the result of calculation of the WENO finite difference method on the discontinuous Galerkin process border described in us, Visualization Demo has been carried out by the Flow Field Calculation result of Tecplot visual software to us, the change of flow field density in double Mach reflection problem is demonstrated in Fig. 7, computing grid is 320*80, initial Mach 2 ship Ma=10, to demonstrate from ρ=1.3965 to ρ=22.682 30 rod density isoline altogether in figure.

The content described in detail is not had all to belong to the general knowledge that the art personnel know in the present invention.

Claims (7)

1. A kind of boundary treatment technology of 1.WENO difference method, is characterized in that the performing step on described method structured grid is as follows:

   First solved problem is carried out structured grid subdivision by A, then institute domain is divided into two kinds of unit areas, and one is internal element region, and another kind is boundary element region, and boundary element region only has one deck computing grid.

   B utilizes traditional WENO-FD method to carry out spatial spreading to inner unit area, utilizes traditional DG method to carry out spatial spreading to unit area, border, and processes border.

   C utilizes the multizone coupling scheme be coupled in DG and WENO method to realize the coupling processing of WENO-FD method and DG method.Need the numerical flux constructing internal element zone boundary place in coupling process, the building method of WENO numerical flux uses special HWENO construction process to realize.

   D boundary element place DG method calculates the n rank approximation by polynomi-als result that this unit is separated, and can obtain the value at this unit center point place and be no more than the derivative value on n rank by this approximation by polynomi-als result.

   E, after spatial spreading terminates, utilizes the TVD Runge-Kutta method of high-order to solve the numerical result obtaining subsequent time.

   F utilizes the visual softwares such as Tecplot or Paraview to carry out the Visual Implementation to the flow field result utilizing said method to obtain, and comprises the density in flow field, speed, pressure etc.
2. 2. a kind of boundary treatment technology of WENO difference method according to claim 1, when it is characterized in that we carry out mesh generation in step, by the boundary element region that one deck stress and strain model of boundary is independent, when carrying out spatial spreading in step B, boundary treatment process is born by DG computation process completely, and no longer need to construct multiple virtual grid unit, thus ensure compactness.
3. 3. a kind of boundary treatment technology of WENO difference method according to claim 1, when it is characterized in that us two kinds of methods be coupled in step C, the structure of coupled interface place WENO numerical flux, we utilize the thought of HWENO, utilize Hermite interpolation to construct WENO-FD internal element border flux f ^wENO.
4. 4. according to a kind of boundary treatment technology of the WENO difference method described in claim 1, it is characterized in that the information needed when carrying out Hermite interpolation in we step C, a part is calculated the unit center place of gained value by internal element WENO-FD method provides, and the value of the boundary element center that another part is obtained by DG method approximation by polynomi-als result in step D and higher derivative value provide.
5. 5., according to a kind of boundary treatment technology of the WENO difference method described in claim 1, it is characterized in that we need the DG numerical flux f required for tectonic boundary cells D G method at WENO-FD method internal element boundary ^dG(u ^-, u ⁺), this numerical flux f ^dG(u ^-, u ⁺) in u ^-then need to be provided by WENO-FD method, we use the u that the mode of Lagrange's interpolation provides here ^-corresponding point estimated value, the information that wherein Lagrange's interpolation needs is provided by internal element completely.
6. 6. according to a kind of boundary treatment technology of the WENO difference method described in claim 1, it is characterized in that we use the TVD Runge-Kutta method of display to carry out time discrete in step e, in order to keep form to be stable in departure process, time step is selected to have certain restriction, Δ x≤CFL* Δ x, CFL=min (WENOCFL, DGCFL).
7. 7., according to a kind of boundary treatment technology of the WENO difference method described in claim 1, it is characterized in that in step F, we utilize the result of the visual softwares such as Tecplot to checking examples such as two-dimentional scalar Hyperbolic Conservation equation and Two-dimensional Euler Equations groups to carry out visualization processing.