CN111859529A

CN111859529A - Multi-grid disturbance domain updating acceleration method for aircraft streaming numerical simulation

Info

Publication number: CN111859529A
Application number: CN202010528125.7A
Authority: CN
Inventors: 蒋崇文; 胡姝瑶; 高振勋; 李椿萱
Original assignee: Beihang University
Current assignee: Beihang University
Priority date: 2020-06-11
Filing date: 2020-06-11
Publication date: 2020-10-30
Anticipated expiration: 2040-06-11
Also published as: CN111859529B

Abstract

The invention discloses a multiple grid disturbance domain updating and accelerating method for aircraft flow-around numerical simulation, which aims at the problem that a large amount of invalid calculation exists in the existing multiple grid method when aircraft flow-around is simulated, and solves a flow control equation in each iteration step in a convection dynamic calculation domain only containing unconverged disturbed units by adopting two dynamic calculation domains of convection and viscosity on a fine grid and a coarse grid, and considers the viscosity effect in a viscosity dynamic calculation domain only containing viscous effect leading units, thereby realizing the solution idea of only solving the unconverged disturbed units and only considering the viscosity effect in a local area. By effectively reducing the calculation amount of a single iteration step and the number of iteration steps required for reaching a convergence state, the method can obviously improve the calculation efficiency of simulating the aircraft flow by using a multi-grid method.

Description

Multi-grid disturbance domain updating acceleration method for aircraft streaming numerical simulation

Technical Field

The invention relates to the technical field of computational fluid mechanics, in particular to a multi-grid disturbance domain updating and accelerating method for aircraft bypass numerical simulation.

Background

With the rapid development of computer technology, computational fluid dynamics numerical simulation has become one of the important technical supports for the aerodynamic design of the modern aircraft. Although numerical simulation can accurately predict static and dynamic aerodynamic characteristics of an aircraft, the problem that numerical simulation calculation efficiency cannot meet the increasing practical application requirements still exists for the aerodynamic design of the aircraft. Therefore, the high-efficiency method suitable for simulating the aircraft bypass flow numerical value is of great significance and engineering application value for improving the aircraft design iteration efficiency and shortening the development period.

The multiple grid method is one of the mainstream acceleration techniques for the numerical simulation of aircraft bypass flow. The method can be used for iterative solution of the steady-state streaming flow field in aircraft static aerodynamic characteristic prediction and is also suitable for internal iterative solution in iterative propulsion of the unsteady streaming flow field in aircraft dynamic aerodynamic characteristic prediction. The multiple grid method is an accelerated convergence technology for solving a flow control equation circularly between a coarse grid and a fine grid, and the acceleration principle of the multiple grid method has two aspects: firstly, the coarse grid can advance the solution with larger iteration step length, so that the fine grid can be driven to reach the convergence state more quickly; second, the cyclic solution between the coarse and fine grids helps to attenuate high-frequency errors that hinder the convergence of the fine grid, so that the overall error of the fine grid can be attenuated more quickly.

At present, the multi-grid method still has a disadvantage when aiming at the aircraft flow-around numerical simulation, namely, each iteration step needs to update and solve all grid units in a preset calculation domain of an aircraft flow field grid. The global updating solution ignores the characteristics of the evolution of the aircraft flow-around flow field in the solution process, thereby causing invalid calculation in different degrees. For a subsonic aircraft, the mathematical property of the flow is elliptical, so that in order to avoid non-physical response caused by a limited boundary of numerical simulation, a preset calculation domain of an aircraft flow field grid needs to be simulated in a range which is tens of times of the reference length of the aircraft. However, as the aircraft bypass flow field gradually approaches to a steady state, the area which does not satisfy the flow control equation is correspondingly contracted, and the preset calculation domain of the aircraft flow field grid also tends to be redundant. For supersonic aircraft, the mathematical properties of the flow are hyperbolic, i.e. any point in the flow field is only affected by the flow in its upstream dependency domain. Therefore, in the iteration process, the area located at the downstream of the aircraft bypass flow field does not converge before the upstream flow field, and the upstream flow which has reached the steady state is not influenced by the updating of the unconverged downstream flow. Therefore, in the whole iterative solution process, the global updating solution of all grid units in the preset calculation domain of the aircraft flow field grid always generates a large amount of invalid calculations, so that the calculation efficiency of the aircraft bypass numerical simulation is obviously reduced. Therefore, if invalid calculation in the existing method can be eliminated, the calculation efficiency of the aircraft bypass numerical simulation is certainly and remarkably improved, and the problem that the numerical simulation calculation efficiency in the aircraft aerodynamic design cannot meet the increasing practical application requirements is effectively solved.

Disclosure of Invention

In view of this, the present invention provides a multiple-grid disturbance domain updating acceleration method for aircraft bypass numerical simulation, so as to solve the problem that a large amount of invalid calculations exist in the existing multiple-grid method when aircraft bypass is simulated.

The invention provides a method for updating and accelerating multiple grid disturbance domains for aircraft bypass numerical simulation, which is characterized by comprising the following steps of:

s1: reading data, including fine grids of an aircraft flow field, a preset calculation domain, boundary conditions and calculation settings;

s2: generating a coarse grid with a specified number of layers by thinning the fine grid;

s3: initializing flow fields in preset calculation domains on the coarse grid and the fine grid according to an incoming flow condition or a given flow field;

s4: establishing a convection dynamic calculation domain and a viscosity dynamic calculation domain according to a flow field initialization mode;

s5: dividing data to be stored for solving all layers of grids into two types; the first kind of data is inherent information of a unit, and a static data structure is adopted to store information of all grid units in the preset calculation domain, wherein the information comprises grid coordinates and flow field variables; the second kind of data is information related to solving and updating, and only the information of grid cells in the convection dynamic calculation domain is stored by adopting a dynamic data structure, wherein the information comprises a conservative updating quantity and a local iteration step length;

S6: assigning a value to the conservation quantity of the boundary virtual grid of the current layer grid according to the type of the boundary condition;

s7: dividing residual terms of a flow control equation into two types of non-sticky terms and sticky terms; wherein the inviscid terms comprise forcing functions of the convection flux and the coarse grid, the inviscid terms being computed in the convection dynamics computation domain; the viscosity term comprises a viscosity flux and turbulence model equation source term, the viscosity term being calculated in a viscosity dynamics computation domain;

s8: in the convection dynamic calculation domain, solving a conservative updating quantity and updating a flow field variable;

s9: judging whether the current grid layer is a fine grid or not; if yes, go to step S10; if not, go to step S21;

s10: judging whether the fine grid reaches a convergence condition; if not, go to step S11; if yes, go to step S23;

s11: judging whether all boundary units of the convection dynamic calculation domain are subjected to inviscid disturbance one by one; if yes, go to step S12; if not, returning to the step S11, and judging the next boundary unit of the convection dynamic calculation domain; after traversing all the boundary cells of the convection dynamic calculation domain, executing step S13;

s12: measuring the propagation direction of the inviscid disturbance, bringing the cells positioned in the propagation direction in the adjacent cells of the boundary cell into the convection dynamic calculation domain, returning to the step S11, and judging the next boundary cell in the convection dynamic calculation domain;

S13: judging whether all boundary units of the convection dynamic calculation domain meet four conditions of converged, in compressible flow, at the most upstream and no longer influenced by other units in the convection dynamic calculation domain; if yes, go to step S14; if not, returning to the step S13, and judging the next boundary unit of the convection dynamic calculation domain; after traversing all the boundary cells of the convection dynamic calculation domain, executing step S16;

s14: removing the boundary unit from the convection dynamic calculation domain, and judging whether the boundary unit exists in the viscosity dynamic calculation domain; if yes, go to step S15; if not, returning to the step S13, and judging the next boundary unit of the convection dynamic calculation domain;

s15: removing the boundary cell from the viscous dynamic calculation domain, returning to step S13, and determining the next boundary cell in the convective dynamic calculation domain;

s16: judging whether all boundary units of the viscosity dynamic calculation domain are dominated by viscosity effect one by one; if yes, go to step S17; if not, returning to the step S16, and judging the next boundary unit of the viscosity dynamic calculation domain; after traversing all boundary cells of the sticky dynamic computation domain, executing step S18;

S17: bringing the cell located in the convection dynamic calculation domain in the immediate cell of the boundary cell into the viscous dynamic calculation domain, returning to step S16, and determining the next boundary cell of the viscous dynamic calculation domain;

s18: updating the convection dynamic calculation domain and the viscosity dynamic calculation domain of all layers of the coarse grids according to the convection dynamic calculation domain and the viscosity dynamic calculation domain of the fine grids;

s19: reallocating storage spaces of second-class data in all the coarse grids and the fine grids according to the convection dynamic calculation domain;

s20: interpolating the conservative quantity to a coarser grid, returning to the step S6, and entering the calculation of the coarser grid;

s21: judging whether a coarser grid needs to be calculated; if yes, return to step S20; if not, go to step S22;

s22: interpolating the conservative quantity correction back to the finer grid, and judging whether the finer grid is a fine grid or not; if yes, returning to the step S6, and entering the calculation of the next iteration step; if not, returning to the step S6 to enter the calculation of the finer grid;

s23: and outputting the result.

In a possible implementation manner, in the method for accelerating updating of a disturbance domain of multiple grids provided by the present invention, step S3 is to initialize flow fields in preset computation domains on the coarse grid and the fine grid according to an incoming flow condition or according to a given flow field, and specifically includes:

Initializing according to incoming flow conditions, and assigning the conservative quantities of all grid units in a preset calculation domain on the coarse grid and the fine grid as incoming flow values;

and initializing according to a given flow field, assigning the conservation quantities of all grid units in a preset calculation domain on the fine grid to be given flow field values, and obtaining the conservation quantities of all grid units in the preset calculation domain on the coarse grid according to the conservation quantity interpolation of all grid units on the fine grid.

In a possible implementation manner, in the multiple-grid-disturbance-domain updating acceleration method provided by the present invention, step S4 is to establish a convection dynamic calculation domain and a viscosity dynamic calculation domain according to a flow field initialization manner, and specifically includes:

initializing according to incoming flow conditions, regarding a convection dynamic calculation domain, taking 1 layer of cells adjacent to a wall surface as initial cells by the coarsest layer of grids, taking L layers of adjacent cells of the wall surface as initial cells by the fine grids, wherein L satisfies the following conditions:

L＝2^l(1)

wherein l represents the number of coarse grid layers, is a positive integer, and is given by the calculation setting read in step S1; generating initial units of the coarse layer grids according to the initial units of the fine grids; for the viscous dynamic calculation domain, all grid layers take 1 layer of units adjacent to the wall surface as initial units;

According to initialization of a given flow field, for a convection dynamic calculation domain, an initial unit of a fine grid is an interfered unit with flow characteristics in the given flow field inconsistent with incoming flow conditions, and the following conditions are met:

||W-W_∞||/||ΔW⁽¹⁾||_max＞_a,c(2)

wherein W represents a conservative amount; w_∞A conservation quantity representing an incoming flow condition; | Δ W⁽¹⁾||_maxShowing the convection flow

The conservation quantity updating quantity of all units in the state calculation domain is the maximum value of the 1 st iteration step;_a,cthe newly added threshold value of convection is taken as 10^-4≤_a,c≤10^-6(ii) a For a sticky dynamic computing domain, the initial cells of the fine mesh are fine meshesAnd (3) a viscous effect dominant unit in the grid convection dynamic calculation domain satisfies the following conditions:

wherein Ψ represents a viscous effect measurement parameter as a ratio of viscous disturbance to non-viscous disturbance mass flow;

a viscosity effect measurement parameter representing the 1 st iteration; i, J and K respectively represent grid directions;

represents the spectral radius of the convection flux Jacobian matrix along the i direction,

representing the spectral radius of the viscous flux Jacobian matrix along the i direction;_a,vindicating the newly added threshold value of viscosity, take 10^-3≤_a,v≤10^-2(ii) a And the convection dynamic calculation domains of all the coarse grid layers are generated according to the convection dynamic calculation domain of the fine grid, and the viscosity dynamic calculation domains of all the coarse grid layers are generated according to the viscosity dynamic calculation domain of the fine grid.

In a possible implementation manner, in the multiple-grid-disturbance-domain updating acceleration method provided by the present invention, step S7 is to divide the residual terms of the flow control equation into two types, i.e., a non-sticky term and a sticky term; wherein the inviscid terms comprise forcing functions of the convection flux and the coarse grid, the inviscid terms being computed in the convection dynamics computation domain; the viscosity term comprises a viscosity flux and a turbulence model equation source term, and is calculated in a viscosity dynamic calculation domain, and specifically comprises the following steps:

The flow control equation is expressed as:

wherein t represents time, F_cDenotes the convective flux, F_vDenotes viscous flux, Q_FForcing function, Q, representing a coarse grid_TRepresenting the source term of the turbulence model equation, | Ω | representing the volume of the grid cell, N_fThe number of the unit surfaces is shown, and Delta S represents the area of the unit surfaces; the right-hand term of equation (5) is collectively referred to as the residual term of the flow control equation,

is a non-sticky item, and is provided with a sticky part,

is a viscosity term; solving the non-sticky terms of the residual terms in the flow dynamics computation domain and solving the sticky terms of the residual terms in the sticky dynamics computation domain.

In a possible implementation manner, in the method for accelerating updating of a multiple grid disturbance domain provided by the present invention, step S11 is to determine whether all boundary cells of the convection dynamic calculation domain have been subjected to inviscid disturbance one by one, and specifically includes:

judging whether all boundary units of the convection dynamic calculation domain meet the following requirements one by one:

||ΔW||/||ΔW⁽¹⁾||_max＞_a,c(6)

wherein, | Δ W | | | represents a modulus of a constant update quantity of a boundary unit of the convection dynamic calculation domain in the current iteration step;

step S12, measuring the propagation direction without viscous disturbance, and incorporating the cell located in the propagation direction in the cells immediately adjacent to the boundary cell into the convection dynamic calculation domain, which specifically includes:

The inviscid disturbance is transmitted at the sound velocity relative to the flow, and the transmission speed of the inviscid disturbance is positive along the transmission direction of the inviscid disturbance; let q denote the unit direction vector, then the propagation of the inviscid perturbation along the q direction is expressed as:

u·q+a＞0 (7)

wherein u represents a flow velocity vector, and a represents a sound velocity; and q is a unit vector of a grid point of the boundary unit of the convection dynamic calculation domain pointing to a certain grid point of the boundary unit, and if q satisfies the formula (7), an adjacent unit sharing the grid point with the boundary unit of the convection dynamic calculation domain is brought into the convection dynamic calculation domain.

In a possible implementation manner, in the aforementioned multiple-grid-disturbance-domain updating acceleration method provided by the present invention, step S13 is to determine, one by one, whether all boundary cells of the convection dynamic computation domain satisfy four conditions that have converged, are located in a compressible flow, are located at the most upstream, and are no longer affected by other cells in the convection dynamic computation domain; if yes, go to step S14; if not, returning to the step S13, and judging the next boundary unit of the convection dynamic calculation domain; the method specifically comprises the following steps:

s131: judging whether all boundary units of the convection dynamic calculation domain meet the following requirements one by one:

||ΔW||/||ΔW⁽¹⁾||_max＜_d(8)

Wherein,_dindicating a deletion threshold, take_d≤10^-7(ii) a If yes, go to step S132; if not, returning to the step S131, and judging a next boundary unit of the convection dynamic calculation domain;

s132: judging whether the Mach number of the boundary unit is larger than 0.3; if yes, go to step S133; if not, returning to the step S131, and judging a next boundary unit of the convection dynamic calculation domain;

s133: judging whether the boundary unit satisfies the following conditions:

wherein, q is a unit vector of grid center of the convection dynamic calculation domain boundary unit pointing to the adjacent cell center of the boundary unit, theta_dRepresenting the tolerance angle of the upstream unit, and taking theta more than or equal to 5 degrees_dLess than or equal to 10 degrees; if yes, go to step S134; if not, returning to the step S131, and judging a next boundary unit of the convection dynamic calculation domain;

s134: for supersonic non-viscous flow, step S14 is directly executed without further judgment;

for subsonic and viscous flows, it is determined whether the boundary cell satisfies:

||ΔW||+||Δ(ΔW)||＜_d(10)

wherein, | Δ (Δ W) | | represents a module value of the boundary unit conservative update quantity affected by the adjacent unit in the convection dynamic calculation domain of the boundary unit; Δ t denotes the iteration step, C_CFLNumber of CFLs, Δ R, representing a time advance format ⁱRepresenting the influence of the residual error items of the adjacent units on the boundary units along the i direction in the flow dynamic calculation domain;

for subsonic inviscid cells, the effect Δ R in the i direction of the residual terms of the adjacent cells to the boundary cell in the convection dynamics computation domainⁱExpressed as:

wherein, Δ F_cRepresenting the amount of change in convective flux; the subscript i +1 represents the immediate vicinity of the boundary cell of the convective dynamics computation domain along the positive i direction, the subscript i-1 represents the immediate vicinity of the boundary cell of the convective dynamics computation domain along the negative i direction, Δ S_i+1/2The area of a boundary unit of the convection dynamic calculation domain along a positive i direction unit surface is represented; delta S_i-1/2The area of a boundary unit of the convection dynamic calculation domain along a negative i direction unit surface is represented;

for sticky cells, the effect Δ R in the i direction of the residual term of the immediate cell on the boundary cell in the flow dynamics computation domainⁱExpressed as:

if yes, go to step S14; if not, the process returns to step S131 to determine the next boundary unit of the convection dynamic calculation domain.

In a possible implementation manner, in the multiple-grid-disturbance-domain updating acceleration method provided by the present invention, step S16 is to determine whether all boundary units of the sticky dynamic calculation domain are dominated by a sticky effect one by one, and specifically includes:

Judging whether all boundary units of the viscosity dynamic calculation domain meet the following conditions one by one:

wherein phi represents the scaling factor of the newly added viscosity threshold, | | Δ W | | magnetism_maxRepresenting the maximum value of the conservative update quantity of all units in the convection dynamic calculation domain in the current iteration step (| | Δ W | | count_max)_minRepresents the minimum value of the maximum value of the conservation constant updating quantity of the single step from the 1 st iteration step to the current iteration step.

According to the method for updating and accelerating the multiple grid disturbance domains for simulating the aircraft flow-around numerical value, the residual estimation only calculates the inviscid terms of the residual in the convection dynamic calculation domain of each layer of grid, only calculates the viscosity terms of the residual in the viscosity dynamic calculation domain of each layer of grid, and the time integration is only executed in the convection dynamic calculation domain of each layer of grid. Dividing the data to be stored into two types of inherent information of all units and information related to solving and updating; for the first type of data, storing information of all units in a preset calculation domain in a static data structure; for the second type of data, only the information of the grid units in the convection dynamic calculation domain of each grid layer is stored in a dynamic data structure, and the storage space of the second type of data is adjusted in each iteration step according to the convection dynamic calculation domain of each grid layer, so that compared with the existing multi-grid method, the method can effectively reduce the storage space requirement of the aircraft bypass numerical simulation. Aiming at the fine grid of the aircraft flow field, only adding the disturbed unit into a convection dynamic calculation domain according to disturbance propagation characteristics; adding only the grid cells dominated by the viscosity effect into the viscosity dynamic calculation domain of the fine grid by measuring the influence of the viscosity effect; removing the converged unit from the two dynamic calculation domains of convection and viscosity by using the convergence solving characteristic on the premise of not influencing the convergence rate; updating two dynamic calculation domains of convection and viscosity of all coarse grids of the aircraft flow field according to the dynamic calculation domain of the fine grid; therefore, the method can always solve all grid layers of the aircraft flow field only for the disturbed units which are not converged in the solution, and only considers the viscosity in the viscous effect leading unit, thereby effectively avoiding the invalid calculation of the existing multi-grid method in the aircraft flow-around numerical simulation and solving the problem that the numerical simulation calculation efficiency is difficult to meet the increasing practical application requirement of the aircraft aerodynamic design.

Drawings

Fig. 1 is a flowchart of a multiple-grid-disturbance-domain update acceleration method according to embodiment 1 of the present invention;

fig. 2 is a schematic diagram of multiple grids generating a coarse grid from a fine grid according to embodiment 1 of the present invention;

FIG. 3 is a diagram illustrating a multi-trellis loop in example 1 of the present invention;

fig. 4 is an evolution process diagram for solving a transonic speed NACA0012 airfoil non-viscous flow field and a convection dynamic calculation domain by using the multiple grid disturbance domain updating acceleration method provided in embodiment 1 of the present invention;

fig. 5a is a graph illustrating a grid quantity change curve of a dynamic convection calculation domain for transonic speed NACA0012 airfoil non-viscous flow by using a multiple grid disturbance domain updating acceleration method provided in embodiment 1 of the present invention;

fig. 5b is a graph illustrating a convergence curve of transonic NACA0012 airfoil non-viscous flow solved by using the multiple-grid perturbation domain updating acceleration method provided in embodiment 1 of the present invention;

fig. 6 is an evolution process diagram of a flow field, a convection field and a viscosity dynamic calculation field for solving transonic speed RAE2822 airfoil turbulent flow by using the multiple-grid disturbance domain updating acceleration method provided in embodiment 1 of the present invention;

fig. 7a is a graph for solving a convective and viscous dynamic computation domain grid quantity change curve of transonic speed RAE2822 airfoil turbulent flow by using the multiple grid disturbance domain updating acceleration method provided in embodiment 1 of the present invention;

Fig. 7b is a graph for solving the convergence curve of the transonic speed RAE2822 airfoil turbulent streaming by using the multiple grid disturbance domain updating acceleration method provided in embodiment 1 of the present invention.

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 illustrative and are not intended to limit the present invention.

The solving process for simulating the aircraft streaming by adopting the multi-grid method has the following characteristics: (1) disturbance of the flow field begins where the governing equation cannot be satisfied; (2) the disturbance carrying the discontinuous information is gradually transmitted to the surrounding flow field at a limited speed along with time, and the area which is not disturbed still keeps the initial state; (3) in compressible flow, the convergence of the upstream flow field is not later than that of the downstream flow field; (4) whether or not detached, the viscous effect only dominates over a limited area. In order to fully utilize the 4-point characteristics and avoid invalid calculation in the conventional multiple grid method, the invention provides a multiple grid disturbance domain updating acceleration method for simulating the aircraft flow-around numerical value.

The invention provides a method for updating and accelerating multiple grid disturbance domains for aircraft streaming numerical simulation, which comprises the following steps:

s2: generating a coarse grid with a specified number of layers through the sparse fine grid;

s5: dividing data to be stored for solving all layers of grids into two types; the first kind of data is inherent information of units, and a static data structure is adopted to store information of all grid units in a preset calculation domain, wherein the information comprises grid coordinates and flow field variables; the second kind of data is information related to solving and updating, and only stores information of grid cells in a convection dynamic calculation domain by adopting a dynamic data structure, wherein the information comprises a conservative updating quantity and a local iteration step length;

s7: dividing residual terms of a flow control equation into two types of non-sticky terms and sticky terms; wherein the inviscid term comprises a forcing function of the convection flux and the coarse grid, and the inviscid term is calculated in a convection dynamic calculation domain; the viscosity term comprises viscosity flux and turbulence model equation source terms, and is calculated in a viscosity dynamics calculation domain;

S8: in the convection dynamic calculation domain, solving the conservative updating quantity and updating the flow field variable;

s11: judging whether all boundary units of the opposite flow dynamic calculation domain are subjected to inviscid disturbance one by one; if yes, go to step S12; if not, returning to the step S11, and judging the next boundary unit of the dynamic calculation domain of the convection current; after traversing all the boundary cells of the convection dynamic calculation domain, executing step S13;

s12: measuring the propagation direction without viscous disturbance, bringing the cell positioned in the propagation direction in the adjacent cell of the boundary cell into a convection dynamic calculation domain, returning to the step S11, and judging the next boundary cell of the convection dynamic calculation domain;

s13: judging whether all boundary units of the convection dynamic calculation domain meet four conditions of converged, positioned in the compressible flow, positioned at the most upstream and not influenced by other units in the convection dynamic calculation domain; if yes, go to step S14; if not, returning to the step S13, and judging the next boundary unit of the dynamic calculation domain of the convection current; after traversing all the boundary cells of the convection dynamic calculation domain, executing step S16;

S14: removing the boundary unit from the convection dynamic calculation domain, and judging whether the boundary unit exists in the viscosity dynamic calculation domain; if yes, go to step S15; if not, returning to the step S13, and judging the next boundary unit of the dynamic calculation domain of the convection current;

s15: removing the boundary cell from the viscous dynamic computation domain, returning to step S13, and determining the next boundary cell in the dynamic computation domain;

s16: judging whether all boundary units of the viscous dynamic calculation domain are dominated by a viscous effect one by one; if yes, go to step S17; if not, returning to the step S16 to judge the next boundary unit of the viscosity dynamic calculation domain; after traversing all boundary cells of the sticky dynamic computation domain, executing step S18;

s17: bringing the cell located in the pair flow dynamic calculation domain in the immediate neighboring cell of the boundary cell into the viscous dynamic calculation domain, returning to step S16, and determining the next boundary cell in the viscous dynamic calculation domain;

s23: and outputting the result.

The following describes a specific implementation of the above-mentioned multiple-grid-perturbation-domain update acceleration method provided by the present invention in detail through a specific embodiment.

Example 1: as shown in fig. 1, the method for accelerating the update of the multiple grid perturbation domains provided by the present invention specifically includes the following steps:

the first step is as follows: and (6) reading in data.

The read-in data comprises fine grids of the aircraft flow field, a preset calculation domain, boundary conditions and calculation settings.

The second step is that: a coarse mesh is generated.

And generating a coarse grid with a specified number of layers through the sparse fine grid. For the structural mesh, a common coarse mesh generation method is to sparsely generate the fine mesh by one time in all mesh directions, and fig. 2 illustrates the coarse mesh generation method by taking one dimension as an example. In fig. 2, h represents a fine mesh, and 2h and 4h represent a layer 1 coarse mesh and a layer 2 coarse mesh, respectively. That is, each cell of the layer 1 coarse mesh contains 2 fine mesh cells, and each cell in the layer 2 coarse mesh contains 2 layer 1 coarse mesh cells.

The third step: and initializing a flow field.

The initialization of the flow field is divided into two modes of initialization according to the inflow condition and initialization according to the given flow field. When the initialization is carried out according to the incoming flow conditions, the conservative quantities of all grid units in the preset calculation domains on the coarse grid and the fine grid are assigned as incoming flow values. When the method is initialized according to a given flow field, the conservation quantities of all grid units in a preset calculation domain on the fine grid are assigned to a given flow field value, and the conservation quantities of all grid units in the preset calculation domain on the coarse grid are obtained according to the conservation quantity interpolation of all grid units on the fine grid. Taking one dimension as an example, the conservation quantity of the coarse grid unit satisfies:

wherein, W_2hRepresenting the conservation quantity, W, of the grid cells of the coarse grid_hRepresenting the conservation quantity of the grid cells of the fine grid; (| omega |)_2hRepresenting the volume of the coarse grid cells, | Ω -_hRepresenting a volume of the fine-mesh grid cells; 2I-1 denotes the 2I-1 st cell along the grid direction, and 2I denotes the 2I-th cell along the grid direction.

The fourth step: and establishing a convection dynamic calculation domain and a viscosity dynamic calculation domain.

(1) When the method is initialized according to the incoming flow condition, the convection dynamic calculation domain and the viscosity dynamic calculation domain are established according to the wall surface boundary. The disturbance carrying the flow field change information is generated from the surface of the object, and the unit dominated by viscosity is also close to the surface of the object, so that a plurality of layers of adjacent units of the wall surface are taken as initial units in the convection dynamic calculation domain. In order to avoid the influence of the increase speed of the convection dynamic calculation domain on the convergence rate of the multiple-grid method, the coarsest grid takes 1 layer of cells adjacent to the wall surface as initial cells, and for the fine grid under the structural grid framework, the convection dynamic calculation domain of the fine grid needs to take L layers of adjacent cells of the wall surface as initial cells, wherein L meets the following requirements:

L＝2^l(2)

Wherein l represents the number of coarse grid layers, is a positive integer, and is given by the calculation setting read in step S1; generating initial units of the coarse layer grids according to the initial units of the fine grids; for the sticky dynamic calculation domain, all grid layers take 1 layer of units adjacent to the wall surface as initial units so as to avoid invalid calculation caused by too fast growth of the sticky dynamic calculation domain.

(2) When initialized according to a given flow field, a convection dynamic calculation domain and a viscosity dynamic calculation domain are established according to flow characteristics of the given flow field. For the convection dynamic calculation domain, the initial units of the fine grid are disturbed units in a given flow field, the flow characteristics of the disturbed units are not consistent with the conditions of the incoming flow, and the conditions are as follows:

||W-W_∞||/||ΔW⁽¹⁾||_max＞_a,c(3)

wherein, W_∞A conservation quantity representing an incoming flow condition; | Δ W⁽¹⁾||_maxRepresenting all sheets of a convective dynamic computing Domain

The updating amount of the conservation constant of the element is at the maximum value of the 1 st iteration step;_a,cthe newly added threshold value of convection is taken as 10^-4≤_a,c≤10^-6。

For the sticky dynamic computation domain, the initial cell of the fine mesh is the sticky effect dominant cell in the fine mesh convection dynamic computation domain. The viscous effect metric Ψ is defined to represent the ratio of viscous disturbance to non-viscous disturbance mass flow. The viscous effect measurement parameter of the viscous effect dominant cell should be at a larger magnitude, so the cells in the viscous dynamics computation domain should satisfy:

representing viscous flux JacobiThe spectral radius of an matrix along the i direction;_a,vindicating the newly added threshold value of viscosity, take 10^-3≤_a,v≤10^-2。

And the convection dynamic calculation domains of all the coarse grid layers are generated according to the convection dynamic calculation domain of the fine grid, and the viscosity dynamic calculation domains of all the coarse grid layers are generated according to the viscosity dynamic calculation domain of the fine grid.

The fifth step: and allocating a storage space.

Dividing data to be stored for solving all layers of grids into two types; the first kind of data is inherent information of units, and a static data structure is adopted to store information of all grid units in a preset calculation domain, wherein the information comprises grid coordinates and flow field variables; the second kind of data is information related to solving and updating, and a dynamic data structure is adopted to only store information of grid cells in the flow dynamic calculation domain, wherein the information comprises a conservative updating quantity and a local iteration step length.

And for the first type of data, storing the information of all grid cells in a preset calculation domain by adopting a static data structure. For example, let I, J, K denote the grid direction, I _max,J_max,K_maxRespectively representing the number of cells in each grid direction, for a block I_max×J_max×K_maxIf there are 10 quantities to be recorded per cell, a real four-dimensional array A can be used₁(10,I_max,J_max,K_max) And storing.

For the second type of data, only the information of the grid cells in the flow dynamic calculation domain is stored by adopting a dynamic data structure. For example, let I, J, K denote the grid direction, and the cells with the same J, K labels are collectively called a row, and a special data type is defined for storing information of a row in the dynamic computation domains of stream and viscosity, and the special data type includes: a. a single linked list-records the I label range of the line unit of the current iteration step in the updating of two types of dynamic calculation domains; b. an integral two-dimensional array, wherein the I label range of the line unit in the last iteration step is stored when the two types of dynamic calculation domains are updated, and the I label range of the line unit in the current iteration step is stored after the two types of dynamic calculation domains are updated; c. a real two-dimensional arrayAnd storing the second type data of the unit of the row in the dynamic calculation domain of the stream, which corresponds to the I label range of the integer two-dimensional array. For a block I_max×J_max×K_maxIf there are 10 quantities to record per cell, a special data type two-dimensional array A is used ₂(J_max,K_max)，A₂The medium real number type two-dimensional array is B (10, I)₁:I₂)，I₁、I₂Respectively representing the minimum and maximum values of the I indices of the elements in the rows (J, K), 1 < I₁,I₂＜I_max。

In the existing multi-grid method, all data of all grid cells in all grid layers are stored in a static data structure. In contrast, the multiple-grid disturbance domain updating acceleration method provided by the invention divides the data to be stored into two types of inherent information of all units and information related to solving and updating; for the first type of data, storing information of all grid cells in a preset calculation domain in a static data structure; for the second kind of data, only the information of the grid cells in the convection dynamic calculation domain of each grid layer is stored in a dynamic data structure, and the fifteenth step adjusts the storage space of the second kind of data according to the convection dynamic calculation domain of each grid layer in each iteration step. Therefore, compared with the existing multi-grid method, the method can effectively reduce the storage space requirement.

And a sixth step: and (5) processing boundary conditions.

And assigning a value to the conservation quantity of the boundary virtual grid of the current layer grid according to the type of the boundary condition.

The seventh step: and estimating a residual error.

The flow control equation is expressed as:

wherein t represents time, F _cDenotes the convective flux, F_vDenotes viscous flux, Q_FForcing function, Q, representing a coarse grid_TRepresenting the source term of the turbulence model equation, N_fPresentation unitThe number of faces, Δ S, represents the area of the unit face; the right-hand term of equation (6) is collectively referred to as the residual term of the flow control equation,

is a non-sticky item, and is provided with a sticky part,

Eighth step: and (4) time integration.

In the convection dynamic calculation domain, calculating the left-end term of the formula (6) by adopting a time-marching format; the time integration method used is consistent with the existing method.

In the existing multi-grid method, the two steps of residual estimation and time integration which are most time-consuming need to be executed on all grid units of each layer of grid in a static preset calculation domain, and the total solving time is about 99%. In contrast, in the multiple-grid disturbance domain updating acceleration method provided by the invention, the residual estimation in the seventh step only calculates the non-sticky item of the residual in the convection dynamic calculation domain of each layer of grid, and only calculates the sticky item of the residual in the sticky dynamic calculation domain of each layer of grid; the time integration of the eighth step is also performed only in the convective dynamics computation domain of the respective layer meshes. Therefore, the invention can simultaneously reduce the calculation amount of each layer of grid and obviously improve the calculation efficiency of the multi-grid method.

The ninth step: judging whether the current grid layer is a fine grid or not; if yes, executing the tenth step; if not, executing the seventeenth step.

The tenth step: judging whether the fine grid reaches a convergence condition; if not, executing the eleventh step; if yes, executing the nineteenth step.

The eleventh step: the convection dynamics calculation domain is increased.

Judging whether all boundary units of the pair flow dynamic calculation domain are subjected to inviscid disturbance one by one, namely judging whether all boundary units of the pair flow dynamic calculation domain meet the following conditions one by one:

||ΔW||/||ΔW⁽¹⁾||_max＞_a,c(7)

wherein, | Δ W | | | represents a modulus of a constant update quantity of a boundary unit of the convection dynamic calculation domain in the current iteration step; if the boundary unit of the dynamic convection calculation domain satisfies the formula (7), measuring the propagation direction without viscous disturbance, bringing the unit positioned in the propagation direction in the adjacent unit of the boundary unit into the dynamic convection calculation domain, returning to the tenth step, and judging the next boundary unit of the dynamic convection calculation domain; if the boundary unit of the convection dynamic calculation domain does not satisfy the formula (7), returning to the tenth step, and judging the next boundary unit of the convection dynamic calculation domain; and executing the twelfth step after traversing all boundary units of the convection dynamic calculation domain.

u·q+a＞0 (8)

wherein u represents a flow velocity vector, and a represents a sound velocity; and q is a unit vector of a grid point of the boundary unit of the convection dynamic calculation domain pointing to a certain grid point of the boundary unit, and if q satisfies the formula (8), an adjacent unit sharing the grid point with the boundary unit of the convection dynamic calculation domain is brought into the convection dynamic calculation domain.

The twelfth step: and reducing the convection dynamic calculation domain and the viscosity dynamic calculation domain.

In order to reduce the calculation amount of a single iteration step as much as possible without influencing the convergence rate, the boundary unit of the flow dynamic calculation domain can be reduced. And (3) judging the following conditions for all boundary units of the convection dynamic calculation domain one by one to realize the reduction of the convection dynamic calculation domain:

(1) judging whether convergence is achieved on all boundary units of the pair flow dynamic calculation domain one by one, namely judging whether the convergence is met on all boundary units of the pair flow dynamic calculation domain one by one:

||ΔW||/||ΔW⁽¹⁾||_max＜_d(9)

wherein,_dindicating a deletion threshold, take_d≤10^-7(ii) a If it isIf the boundary cell satisfies the formula (9), continuing to judge (2); if not, judging a next boundary unit of the dynamic calculation domain of the pair of flows;

(2) Judging whether the boundary unit is positioned in the compressible flow, namely judging whether the Mach number of the boundary unit is more than 0.3; if yes, continuing to judge (3); if not, returning to the judgment (1) and judging the next boundary unit of the dynamic calculation domain of the pair of flows;

(3) judging whether the boundary unit is positioned at the most upstream, namely judging whether the boundary unit satisfies the following conditions:

wherein, q is a unit vector of grid center of the convection dynamic calculation domain boundary unit pointing to the adjacent cell center of the boundary unit, theta_dRepresenting the tolerance angle of the upstream unit, and taking theta more than or equal to 5 degrees_dLess than or equal to 10 degrees; if the boundary cell satisfies the formula (10), continuing to judge (4); if not, returning to the judgment (1) and judging the next boundary unit of the dynamic calculation domain of the pair of flows;

(4) judging whether the boundary unit is not influenced by other units in the convection dynamic calculation domain any more, and if the boundary unit does not flow at the supersonic speed, judging is not needed and the boundary unit can automatically meet the requirement; for subsonic and viscous flows, it is determined whether the boundary cell satisfies:

||ΔW||+||Δ(ΔW)||＜_d(11)

wherein, | Δ (Δ W) | | represents a module value of the boundary unit conservative update quantity affected by the adjacent unit in the convection dynamic calculation domain of the boundary unit; Δ t denotes the iteration step, C_CFLNumber of CFLs, Δ R, representing a time advance format ⁱRepresenting the influence of the residual error items of the adjacent units on the boundary units along the i direction in the flow dynamic calculation domain; for subsonic inviscid cells, the effect Δ R in the i direction of the residual terms of the adjacent cells to the boundary cell in the convection dynamics computation domainⁱExpressed as:

wherein, Δ F_cRepresenting the amount of change in convective flux; the subscript i +1 represents the immediate vicinity of the boundary cell of the convective dynamics computation domain along the positive i direction, the subscript i-1 represents the immediate vicinity of the boundary cell of the convective dynamics computation domain along the negative i direction, Δ S_i+1/2The area of a boundary unit of the convection dynamic calculation domain along a positive i direction unit surface is represented; delta S_i-1/2The area of a boundary unit of the convection dynamic calculation domain along a negative i direction unit surface is represented; for sticky cells, the effect Δ R in the i direction of the residual term of the immediate cell on the boundary cell in the flow dynamics computation domainⁱExpressed as:

if the boundary cell satisfies the formula (11) and the formula (12), the boundary cell simultaneously satisfies four conditions of converged, in compressible flow, at the most upstream and no longer influenced by other cells in the convection dynamic calculation domain, removing the boundary cell from the convection dynamic calculation domain, if the boundary cell simultaneously exists in the viscous dynamic calculation domain, removing the boundary cell from the viscous dynamic calculation domain, returning to the twelfth step, and judging the next boundary cell in the convection dynamic calculation domain; if the boundary cell does not satisfy the formula (11) and the formula (12), the process returns to the twelfth step, and the next boundary cell in the dynamic calculation domain of the convection current is determined. And executing the thirteenth step after traversing all boundary units of the convection dynamic calculation domain.

The thirteenth step: the sticky dynamic computation domain is increased.

Judging whether all boundary units of the viscosity dynamic calculation domain are dominated by the viscosity effect one by one, namely judging whether all boundary units of the viscosity dynamic calculation domain meet the following requirements one by one:

wherein phi represents the scaling factor of the newly added viscosity threshold, | | Δ W | | magnetism_maxRepresenting the maximum value of the conservative update quantity of all units in the convection dynamic calculation domain in the current iteration step (| | Δ W | | count_max)_minRepresenting the minimum value of the maximum value of the conservation-constant updating quantity in the single step from the 1 st iteration step to the current iteration step; if the boundary unit of the viscous dynamic calculation domain satisfies the formula (15) and the formula (16), the unit located in the flow dynamic calculation domain in the unit immediately adjacent to the boundary unit is brought into the viscous dynamic calculation domain, the step three is returned, and the next boundary unit of the viscous dynamic calculation domain is judged; if the boundary unit of the viscous dynamic calculation domain does not satisfy the formula (15) and the formula (16), returning to the thirteenth step, and judging the next boundary unit of the viscous dynamic calculation domain; after traversing all the boundary cells of the sticky dynamic computation domain, the fourteenth step is performed. The viscous term of the flow control equation is elliptical, so that all the adjacent units of the units which are disturbed by the viscosity can be added into the viscosity dynamic calculation domain without increasing the convection dynamic calculation domain.

The fourteenth step is that: and updating the coarse grid dynamic calculation domain.

And updating the convection dynamic calculation domains of all layers of the coarse grids according to the convection dynamic calculation domains of the fine grids, and updating the viscosity dynamic calculation domains of all layers of the coarse grids according to the viscosity dynamic calculation domains of the fine grids. For the structural grid, the corresponding relationship between the dynamic computation domains of the coarse grid and the fine grid is shown in fig. 2, and if a certain fine grid cell is included in the dynamic computation domain of the fine grid, the coarse grid cell corresponding to the cell is also added into the dynamic computation domain of the coarse grid.

In the tenth step, aiming at the fine grid, only adding disturbed units into a convection dynamic calculation domain according to disturbance propagation characteristics; in the tenth step, only the grid cells dominated by the viscosity effect are added into the viscosity dynamic calculation domain of the fine grid by measuring the influence of the viscosity effect; in the twelfth step, the converged units are removed from the two dynamic calculation domains of convection and viscosity in time by using the convergence solving characteristic on the premise of not influencing the convergence rate; in the fourteenth step, according to the dynamic calculation domain of the fine grid, updating two types of dynamic calculation domains of convection and stickiness of all the coarse grids. Therefore, the method can always solve only the disturbed units which are not converged in the solution in all the grid layers, and only considers the viscosity in the viscous effect leading unit, thereby effectively avoiding the invalid calculation in the solution of the existing multi-grid method.

The fifteenth step: and reallocating the storage space of the second type data.

And reallocating the storage space of the second type data in all the coarse grids and the fine grids according to the convection dynamic calculation domain.

Sixteenth, step: and (4) interpolating the conservative quantity to a coarser grid, returning to the sixth step, and entering the calculation of the coarser grid.

The interpolation mode of the conservative quantity from the fine grid to the coarse grid is the same as the interpolation mode initialized according to the given flow field in the third step.

Seventeenth step: judging whether a coarser grid needs to be calculated; if the coarser grid needs to be calculated, returning to the sixteenth step; if the coarser grid does not need to be calculated, executing the eighteenth step;

the common cyclic mode of multiple grids is shown in fig. 3, the left subgraph represents "V" type cycle, the right subgraph represents "W" type cycle, h represents fine grid in fig. 3, 2h represents layer 1 coarse grid, and 4h represents layer 2 coarse grid. Taking a V-shaped cycle as an example, if the current calculation layer is a fine grid, selecting a coarser grid to be calculated in the seventeenth step; if the current calculation layer is the coarsest grid, the seventeenth step selects to not need to calculate a coarser grid.

And eighteenth step: interpolating the conservative quantity correction back to the finer grid, and judging whether the finer grid is a fine grid or not; if the finer grid is the fine grid, returning to the sixth step and entering the calculation of the next iteration step; and if the finer grid is not the fine grid, returning to the sixth step and entering the calculation of the finer grid.

The nineteenth step: and outputting the result.

The NACA0012 airfoil non-turbulent problem at Mach 0.755 and the RAE2822 airfoil turbulent flow problem at Mach 0.729 were simulated separately below.

(1) The NACA0012 airfoil no-stick-flow problem at mach 0.755 was simulated.

Fig. 4 illustrates the evolution process of the convective dynamic computation domain of the above-mentioned multiple-grid-disturbance-domain update acceleration method provided by the present invention. In FIG. 4,. eta._cRepresenting the ratio of the grid quantities of the convection dynamic calculation Domain and the Preset calculation Domain, N_maxRepresenting the total number of iteration steps required to reach the converged state. The flow field is initialized with incoming flow conditions and disturbances are generated from the wall boundaries. In the initial stage of solving, the convection current dynamic calculation domain is gradually expanded to the whole preset calculation domain from the wall surface. In the middle and later stages of the solution, the convection dynamic calculation domain gradually shrinks from a far field to a wall surface and from upstream to downstream along with the convergence of the solution, and finally only a near-wall region and a downstream region thereof are left. As can be seen from the grid quantity change curve of the convection dynamic calculation domain of the multiple grid disturbance domain update acceleration method shown in fig. 5a, the multiple grid disturbance domain update acceleration method provided by the present invention is applicable to multiple grids with different grid layer numbers, and the evolution processes of the convection dynamic calculation domain under different grid layer numbers are basically the same.

Fig. 5b is a convergence curve comparing the above-mentioned multiple-grid perturbation domain updating acceleration method provided by the present invention with the existing multiple-grid method of global updating. The comparison result shows that for the multiple grid method with different layers, compared with the existing method, the multiple grid disturbance domain updating acceleration method provided by the invention can obviously improve the numerical simulation efficiency by reducing the calculation amount of a single iteration step and the total iteration step number. When the method adopts 1-layer coarse grids, the calculation time can be saved by 52.1 percent compared with the conventional multi-grid method; when the invention adopts 2-layer coarse grids, 45.6 percent of calculation time can be saved compared with the existing multi-grid method.

(2) A transonic RAE2822 airfoil turbulent streaming problem with a Mach number of 0.729 is simulated.

Fig. 6 shows a dynamic computation domain evolution process of the multiple-grid-disturbance-domain updating acceleration method provided by the invention. Flow ofThe fields are initialized with incoming flow conditions, and thus, the flow dynamics computation domain and the viscosity dynamics computation domain are still built from the walls. The rate of increase of the flow dynamics computation domain is much greater than the rate of increase of the viscosity dynamics computation domain. From 1% N in FIG. 6_maxIt can be seen that when the convection dynamics calculation domain has covered the whole preset calculation domain, the viscosity dynamics calculation domain remains near the wall surface, and only includes the near-wall region and its downstream region. As can be seen from the change curves of the grid volume of the convection dynamic calculation domain and the viscosity dynamic calculation domain of the multiple grid disturbance domain update acceleration method shown in fig. 7a, the multiple grid disturbance domain update acceleration method provided by the present invention is suitable for viscous flow, and both the convection dynamic calculation domain and the viscosity dynamic calculation domain follow the evolution law of first increasing and then shrinking.

Fig. 7b is a convergence curve comparing the above-mentioned multiple-grid-perturbation-domain updating acceleration method provided by the present invention with the existing multiple-grid method of global updating. The comparison result shows that the multiple-grid disturbance domain updating acceleration method provided by the invention can obtain a convergence solution as the existing multiple-grid method. By virtue of the reduction of two dynamic calculation domains of convection and viscosity and the consideration of the viscosity only in local parts, the method for accelerating the updating of the multiple grid disturbance domains provided by the invention has the advantages that the calculation efficiency is remarkably improved, and the calculation time can be saved by 40.7%.

It will be apparent to those skilled in the art that various changes and modifications may be made in the present invention without departing from the spirit and scope of the invention. Thus, if such modifications and variations of the present invention fall within the scope of the claims of the present invention and their equivalents, the present invention is also intended to include such modifications and variations.

Claims

1. A multiple grid disturbance domain updating and accelerating method for aircraft bypass flow numerical simulation is characterized by comprising the following steps:

s23: and outputting the result.

2. The method for accelerating updating of multiple grid disturbance domains according to claim 1, wherein step S3, initializing flow fields in preset computation domains on the coarse grid and the fine grid according to an incoming flow condition or according to a given flow field, specifically comprises:

3. The method for accelerating updating of multiple grid disturbance domains according to claim 1, wherein step S4 is to establish a convection dynamic calculation domain and a viscosity dynamic calculation domain according to a flow field initialization method, and specifically includes:

L＝2^l(1)

||W-W_∞||/||ΔW⁽¹⁾||_max＞_a,c(2)

wherein W represents a conservative amount; w_∞A conservation quantity representing an incoming flow condition; | Δ W⁽¹⁾||_maxRepresenting the maximum value of the conservative updating quantity of all units in the convection dynamic calculation domain in the 1 st iteration step;_a,cthe newly added threshold value of convection is taken as 10^-4≤_a,c≤10^-6(ii) a For the viscous dynamic calculation domain, the initial unit of the fine grid is a viscous effect dominant unit in the fine grid convection dynamic calculation domain, and the following conditions are met:

4. The multiple mesh perturbation domain updating acceleration method according to claim 3 characterized in that step S7, the residual terms of the flow control equation are divided into two types, non-sticky terms and sticky terms; wherein the inviscid terms comprise forcing functions of the convection flux and the coarse grid, the inviscid terms being computed in the convection dynamics computation domain; the viscosity term comprises a viscosity flux and a turbulence model equation source term, and is calculated in a viscosity dynamic calculation domain, and specifically comprises the following steps:

The flow control equation is expressed as:

is a non-sticky item, and is provided with a sticky part,

5. The method as claimed in claim 4, wherein the step S11 of determining whether all boundary cells of the convective dynamic computation domain have been disturbed without sticking one by one includes:

||ΔW||/||ΔW⁽¹⁾||_max＞_a,c(6)

u·q+a＞0 (7)

6. The multiple-grid-disturbance-domain updating acceleration method according to claim 5, wherein step S13 is implemented to determine whether four conditions of converged, in compressible flow, in most upstream, and no longer influenced by other units in the convective dynamic computation domain are simultaneously satisfied for all boundary cells of the convective dynamic computation domain one by one; if yes, go to step S14; if not, returning to the step S13, and judging the next boundary unit of the convection dynamic calculation domain; the method specifically comprises the following steps:

||ΔW||/||ΔW⁽¹⁾||_max＜_d(8)

wherein,_dindicating a deletion threshold, take _d≤10^-7(ii) a If yes, go to step S132; if not, returning to the step S131, and judging a next boundary unit of the convection dynamic calculation domain;

s133: judging whether the boundary unit satisfies the following conditions:

||ΔW||+||Δ(ΔW)||＜_d(10)

wherein, | Δ (Δ W) | | represents a module value of the boundary unit conservative update quantity affected by the adjacent unit in the convection dynamic calculation domain of the boundary unit; Δ t denotes the iteration step, C_CFLNumber of CFLs, Δ R, representing a time advance formatⁱRepresenting the influence of the residual error items of the adjacent units on the boundary units along the i direction in the flow dynamic calculation domain;

wherein, Δ F_cRepresenting the amount of change in convective flux; the subscript i +1 represents the immediate vicinity of the boundary cell of the convective dynamics computation domain along the positive i direction, the subscript i-1 represents the immediate vicinity of the boundary cell of the convective dynamics computation domain along the negative i direction, Δ S_i+1/2The area of a boundary unit of the convection dynamic calculation domain along a positive i direction unit surface is represented; delta S_i-_1/2The area of a boundary unit of the convection dynamic calculation domain along a negative i direction unit surface is represented;

7. The method for accelerating updating of multiple grid disturbance domains according to claim 6, wherein step S16, determining whether all boundary cells of the sticky dynamic computation domain are dominated by a sticky effect one by one, specifically comprises: