CN111859530B

CN111859530B - Iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation

Info

Publication number: CN111859530B
Application number: CN202010528415.1A
Authority: CN
Inventors: 蒋崇文; 胡姝瑶; 高振勋; 李椿萱
Original assignee: Beihang University
Current assignee: Beihang University
Priority date: 2020-06-11
Filing date: 2020-06-11
Publication date: 2022-04-22
Anticipated expiration: 2040-06-11
Also published as: CN111859530A

Abstract

The invention discloses an iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation, which adopts three types of dynamic calculation domains of unsteady, convection and viscosity to provide a solution idea for solving an unconverged unsteady unit and considering a viscosity effect only in a local area aiming at the problem of a large amount of invalid calculations in the existing aircraft dynamic aerodynamic characteristic numerical simulation method. The unsteady dynamic calculation domain only comprises grid units with the dominant unsteady effect at each physical moment, the convective dynamic calculation domain only comprises grid units which are disturbed in each iteration step and are not converged in solution, the viscous dynamic calculation domain only comprises grid units with the dominant viscous effect in the convective dynamic calculation domain, and the three types of dynamic calculation domains are updated in each iteration step in a self-adaptive manner. By avoiding invalid calculation in the unsteady numerical simulation and reducing the iteration steps required by the internal iteration to reach the convergence state, the method can obviously improve the calculation efficiency of the aircraft dynamic aerodynamic characteristic numerical simulation.

Description

Iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation

Technical Field

The invention relates to the technical field of computational fluid mechanics, in particular to an iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation.

Background

Dynamic aerodynamic properties are key performance parameters for aircraft design. Due to the rapid development of computer hardware performance, unsteady flow numerical simulation is widely applied to dynamic aerodynamic characteristic prediction of various aircrafts. However, the flow problem is complex, which results in large calculation amount, low solving efficiency of the numerical simulation method, and the like, the numerical simulation of the dynamic aerodynamic characteristics of the aircraft still needs to consume a large amount of calculation resources and an ultra-long solving period, and the problem that the numerical simulation calculation efficiency is difficult to meet the increasing practical application requirements of the aerodynamic design of the aircraft exists. Therefore, it has been a hot research point of computational fluid dynamics to propose an efficient numerical simulation technique suitable for predicting the dynamic aerodynamic characteristics of an aircraft.

Currently, the iterative propulsion method is one of the mainstream methods for aircraft dynamic aerodynamic characteristic numerical simulation, such as a two-time step method based on second-order backward difference. The method introduces internal iteration in the time-marching solution of the aircraft dynamic aerodynamic characteristic numerical simulation, so that the numerical simulation precision of the aircraft dynamic aerodynamic characteristic can be ensured, and the mature aircraft static aerodynamic characteristic numerical simulation method can be popularized to the numerical simulation of the dynamic aerodynamic characteristic.

At present, an acceleration technology suitable for the dynamic aerodynamic characteristic numerical simulation of the iterative propulsion aircraft is mainly realized by three ways of implicit iteration format, spatial parallelization and reasonable grid generation. The implicit iteration format can break the stability condition limit of the iteration step length, and the step number required by the inner iteration is reduced by greatly increasing the iteration step length. The spatial parallelization is to divide the aircraft flow field grid into a plurality of grid blocks, namely to decompose the computation task into a plurality of subtasks for parallel processing, so that the computation amount of serial processing can be reduced. The rationalization of grid generation is to reasonably control the density of the aircraft flow field grid, so that the calculation amount can be reduced by reducing the grid amount while the spatial precision is ensured.

The technology realizes acceleration of the aircraft dynamic aerodynamic characteristic numerical simulation from the aspects of discrete formats, solving algorithms, grids and the like, but global updating and solving are carried out on all grid units in the static calculation domain of the aircraft flow field grid at each physical moment and in each iteration step. This updating method has two disadvantages: on one hand, the static calculation domain of the aircraft flow field grid is crucial to the precision and efficiency of dynamic aerodynamic characteristic numerical simulation, but the selection of the static calculation domain completely depends on experience, the solution fails due to the insufficient selected calculation domain, and the selected calculation domain is too large and wastes calculation resources; on the other hand, in the static calculation domain of the aircraft flow field grid, the numerical simulation can not utilize the characteristic that the required solution range in the solution continuously changes along with the evolution and the solution convergence, so that a large amount of invalid calculation can be caused, and even more numerical errors are introduced. Therefore, if invalid calculation caused by adopting static calculation domain simulation of an aircraft flow field grid in the conventional iterative propulsion aircraft dynamic aerodynamic characteristic numerical simulation method can be eliminated, the calculation amount required by aircraft dynamic aerodynamic characteristic numerical simulation is obviously reduced, the calculation efficiency of aircraft dynamic aerodynamic characteristic numerical simulation is improved, and the problem that the numerical simulation calculation efficiency cannot meet the increasing practical application requirement of aircraft aerodynamic design is effectively solved.

Disclosure of Invention

In view of the above, the present invention provides an iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation, so as to solve the problem of a large amount of invalid calculations existing in the existing aircraft dynamic aerodynamic characteristic numerical simulation method.

The invention provides an iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation, which comprises the following steps:

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

s2: initializing the flow field in the preset calculation domain according to the incoming flow condition or according to a given flow field;

s3: establishing an unsteady dynamic calculation domain, a convection dynamic calculation domain and a viscosity dynamic calculation domain according to a flow field initialization mode;

s4: dividing data required to be stored for solving 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 of the current moment and the last two moments; 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;

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

s6: dividing residual terms of a flow control equation into two types of non-sticky terms and sticky terms; wherein the sticky-free terms comprise source terms of convection flux and iterative marching methods, the sticky-free 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;

s7: in the convection dynamic calculation domain, solving a conservative quantity updating quantity and updating a flow field variable at the current moment;

s8: judging whether all boundary units of the convection dynamic calculation domain are subjected to inviscid disturbance one by one; if yes, go to step S9; if not, returning to the step S8, 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 S10;

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

s10: judging whether the unsteady dynamic calculation domain is increased or not; if yes, executing step S11 and step S12; if not, go to step S12;

s11: reallocating the storage space of the second class of data according to the unsteady dynamic computation domain;

s12: 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 S13; if not, returning to the step S12, 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 S15;

s13: 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 S14; if not, returning to the step S12, and judging the next boundary unit of the convection dynamic calculation domain;

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

s15: judging whether all boundary units of the viscosity dynamic calculation domain are dominated by viscosity effect one by one; if yes, go to step S16; if not, returning to the step S15, 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 S17;

s16: 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 S15, and determining the next boundary cell of the viscous dynamic calculation domain;

s17: judging whether the inner iteration reaches a convergence condition; if yes, go to step S18; if not, returning to the step S5, entering the calculation of the next iteration step of the internal iteration, and repeatedly executing the steps S5-S17;

s18: measuring the influence of the unsteady effect by using a physical time derivative, and judging whether all boundary units of the unsteady dynamic calculation domain are no longer dominated by the unsteady effect one by one; if yes, go to step S19; if not, returning to the step S18, and judging the next boundary unit of the unsteady dynamic calculation domain; after traversing all the boundary cells of the unsteady dynamic computation domain, executing step S20;

s19: removing the boundary cell from the unsteady dynamic computation domain, returning to step S18, and determining the next boundary cell in the unsteady dynamic computation domain;

s20: judging whether the unsteady dynamic calculation domain is reduced or not; if yes, executing step S21 and step S22; if not, go to step S22;

s21: reallocating the storage space of the second class of data according to the unsteady dynamic computation domain;

s22: updating flow field variables at the last two moments, resetting the convection dynamic calculation domain to the range of the unsteady dynamic calculation domain, and resetting the viscous dynamic calculation domain to the region where the viscous effect plays a dominant role in the convection dynamic calculation domain;

s23: judging whether the specified final state physical time is reached; if yes, go to step S24; if not, returning to the step S5, entering the calculation of the next physical time, and repeatedly executing the steps S5-S23;

s24: and outputting the result.

In a possible implementation manner, in the iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation provided by the present invention, step S2 is to initialize a flow field in the preset calculation domain 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 the preset calculation domain as incoming flow values;

and according to the initialization of the given flow field, assigning the conservation quantities of all grid units in the preset calculation domain as the given flow field value.

In a possible implementation manner, in the iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation provided by the present invention, step S3 is to establish an unsteady dynamic calculation domain, a convective dynamic calculation domain, and a viscous dynamic calculation domain according to a flow field initialization manner, and specifically includes:

initializing according to the incoming flow conditions, taking a plurality of layers of adjacent units of the wall surface as initial units in the unsteady dynamic calculation domain and the convection dynamic calculation domain, and taking 1 layer of units adjacent to the wall surface as initial units in the viscous dynamic calculation domain;

according to initialization of a given flow field, an initial unit of a convection dynamic calculation domain 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 (1)

wherein W represents a conservative amount; w_∞A conservation quantity representing an incoming flow condition; | Δ W⁽¹⁾||_maxRepresenting the maximum value of the 1 st iteration step of the conservative updating quantity of all units in the convection dynamic calculation domain at the 1 st physical moment; epsilon_a,cThe newly added threshold value of convection is taken as 10^-4≤ε_a,c≤10^-6(ii) a The initial unit of the unsteady dynamic computation domain is consistent with the initial unit of the convection dynamic computation domain; the initial unit of the viscosity dynamic calculation domain is a viscosity effect dominant unit in the convection dynamic calculation domain, and 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 step at the current physical moment; 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; epsilon_a,vIndicating the newly added threshold value of viscosity, take 10^-3≤ε_a,v≤10^-2。

In a possible implementation manner, in the iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation provided by the invention, step S6 is to divide the residual terms of the flow control equation into two types, namely a non-sticky term and a sticky term; wherein the sticky-free terms comprise source terms of convection flux and iterative marching methods, the sticky-free 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 τ represents a virtual time; Δ t represents a physical time step; f_cDenotes the convective flux, F_vDenotes viscous flux, Q_TEquation of a model representing turbulenceThe source of the item is,

is a source term of the iterative boosting method; | Ω | represents 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; superscripts (n-2), (n-1) and (n) respectively represent the last moment, the last moment and the current moment; the right-hand term of equation (4) 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 iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation provided by the present invention, step S8 is to determine, one by one, whether all boundary cells of the convective dynamic computation domain have been subjected to inviscid disturbance, 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 (5)

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 S9, measuring the propagation direction without viscous disturbance, and incorporating the cells located in the propagation direction in the cells immediately adjacent to the boundary cell into the convection dynamic calculation domain and the unsteady 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 (6)

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

In a possible implementation manner, in the above iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation provided by the present invention, step S12 is to determine, one by one, whether all boundary cells of the convective dynamic computation domain satisfy four conditions of converged, in compressible flow, in the most upstream, and no longer affected by other cells in the convective dynamic computation domain; if yes, go to step S13; if not, returning to the step S12, and judging the next boundary unit of the convection dynamic calculation domain; the method specifically comprises the following steps:

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

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

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

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

s123: 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, execute the stepS124; if not, returning to the step S121, and judging a next boundary unit of the convection dynamic calculation domain;

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

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

||ΔW||+||Δ(ΔW)||＜ε_d (9)

wherein, | Δ (Δ W) | | represents a module value of the boundary unit conservative update quantity affected by the adjacent unit of the boundary unit in the convection dynamic calculation domain; Δ 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 S13; if not, returning to step S121, and determining a next boundary unit of the convection dynamic calculation domain.

In a possible implementation manner, in the iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation provided by the present invention, step S15 is to determine, one by one, whether all boundary units of the viscous dynamic calculation domain are dominated by a viscous effect, which 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)_minAnd the minimum value of the maximum value of the single-step conservation quantity updating quantity from the 1 st iteration step of the internal iteration in the current moment to the current iteration step is represented.

In a possible implementation manner, in the iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation provided by the present invention, in step S18, it is determined one by one whether all boundary units of the unsteady dynamic computation domain are no longer dominated by the unsteady effect, specifically including:

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

according to the iterative propulsion disturbance domain updating method for the aircraft dynamic aerodynamic characteristic simulation, residual estimation is only used for calculating the non-sticky item of the residual in the convection dynamic calculation domain, and the sticky item of the residual is only calculated in the sticky dynamic calculation domain; therefore, the iterative propulsion disturbance domain updating method provided by the invention can efficiently avoid invalid calculation caused by global updating of all grid units in the preset calculation domain of the aircraft flow field grid in the existing method by adopting three dynamic calculation domains of unsteady, convective and viscous, thereby effectively improving the calculation efficiency of the aircraft dynamic aerodynamic characteristic numerical simulation and solving the problem that the numerical simulation calculation efficiency is difficult to meet the increasing practical application requirements of the aircraft aerodynamic design. 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 grid cells in a preset calculation domain in a static data structure; for the second type of data, only the information of the grid cells in the convection dynamic calculation domain is stored in a dynamic data structure, and the storage space of the second type of data is adjusted in each iteration step and each physical moment in time according to the unsteady dynamic domain, so that the memory requirement of the aircraft dynamic aerodynamic characteristic numerical simulation can be effectively reduced by the iterative propulsion disturbance domain updating method provided by the invention. According to the disturbance propagation characteristics, only adding a disturbed unit and a viscous effect leading unit into a convection dynamic calculation domain and a viscous dynamic calculation domain respectively; removing the converged unit from the solution update in time on the premise of not influencing the convergence rate by using the solution convergence characteristics; removing the constant unit from the solution update in time by measuring the influence of the unsteady effect; therefore, the iterative propulsion disturbance domain updating method provided by the invention can always solve only the unconverged unsteady units, and only considers the viscosity in the viscous effect leading unit, thereby effectively avoiding invalid calculation caused by global updating of all grid units in the preset calculation domain of the aircraft flow field grid in the existing method, and solving the problem that the numerical simulation calculation efficiency is difficult to meet the ever-increasing practical application requirement of the aircraft aerodynamic design.

Drawings

Fig. 1 is a flowchart of an iterative perturbation domain updating method provided in embodiment 1 of the present invention;

fig. 2 is an evolution process diagram of a flow field and two types of steady and convective dynamic calculation domains for solving a transonic speed NACA0012 airfoil oscillation problem by using the iterative propulsion disturbance domain updating method provided in embodiment 1 of the present invention;

fig. 3 is an unsteady and convective dynamic calculation domain grid quantity change curve and convergence curve for solving the transonic speed NACA0012 airfoil oscillation problem by using the iterative propulsion disturbance domain updating method provided in embodiment 1 of the present invention;

fig. 4 is an evolution process diagram of a flow field and two types of steady and convective dynamic calculation domains for solving the pitch problem of the supersonic winged missile by using the iterative propulsion disturbance domain updating method provided by the embodiment 1 of the invention;

fig. 5 is an unsteady and convective dynamic calculation domain grid quantity change curve and convergence curve for solving the pitch problem of the supersonic winged missile by using the iterative propulsion disturbance domain updating 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 iterative propulsion solving process for the aircraft dynamic aerodynamic characteristic simulation 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; (5) the unsteady effect dominant region gradually changes along with the evolution of the flow field. In order to fully utilize the 5-point characteristics and avoid invalid calculation in the conventional numerical simulation method, the invention provides an iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation, which is characterized in that only an unconverged unsteady unit is solved, viscosity is considered only in a viscous effect leading unit, the calculation amount of a single iteration step in internal iteration and the number of iteration steps required for reaching a convergence state are effectively reduced, and the efficiency of aircraft dynamic aerodynamic characteristic numerical simulation is improved.

s2: initializing a flow field in a preset calculation domain according to an incoming flow condition or a given flow field;

s4: dividing data required to be stored for solving 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 a preset calculation domain, wherein the information comprises grid coordinates and flow field variables of the current moment and the last two moments; 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;

s6: dividing residual terms of a flow control equation into two types of non-sticky terms and sticky terms; wherein the inviscid terms comprise source terms of convection flux and an iterative propulsion method, and the inviscid terms are 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;

s7: in the convection dynamic calculation domain, solving the conservative updating quantity and updating the flow field variable at the current moment;

s8: 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 S9; if not, returning to the step S8, 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 S10;

s9: measuring the propagation direction without viscous disturbance, respectively bringing the cells positioned in the propagation direction in the adjacent cells of the boundary cell into a convection dynamic calculation domain and an unsteady dynamic calculation domain, returning to the step S8, and judging the next boundary cell of the convection dynamic calculation domain;

s12: 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 S13; if not, returning to the step S12, 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 S15;

s13: 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 S14; if not, returning to the step S12, and judging the next boundary unit of the dynamic calculation domain of the convection current;

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

s15: 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 S16; if not, returning to the step S15 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 S17;

s16: 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 S15, and determining the next boundary cell in the viscous dynamic calculation domain;

s18: measuring the influence of the unsteady effect by using a physical time derivative, and judging whether all boundary units of the unsteady dynamic calculation domain are no longer dominated by the unsteady effect one by one; if yes, go to step S19; if not, returning to the step S18 to judge the next boundary unit of the unsteady dynamic calculation domain; after traversing all the boundary cells of the unsteady dynamic computation domain, executing step S20;

s22: updating flow field variables of the last two moments, resetting the convection dynamic calculation domain to be the range of the unsteady dynamic calculation domain, and resetting the viscous dynamic calculation domain to be the region which plays a dominant role in viscous effect in the flow dynamic calculation domain;

s24: and outputting the result.

The following describes a specific implementation of the above iterative perturbation domain updating method provided by the present invention in detail through a specific embodiment.

Example 1: as shown in fig. 1, the iterative perturbation domain updating method 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 grids of the aircraft flow field, a preset calculation domain, boundary conditions and calculation settings.

The second step is that: 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 initializing according to incoming flow conditions, assigning the conservative quantities of all grid units in a preset calculation domain as incoming flow values; when the method is initialized according to a given flow field, the conservative quantities of all grid units in a preset calculation domain are assigned to be given flow field values.

The third step: three dynamic calculation domains of unsteady, convection and viscosity are established.

(1) When the method is initialized according to the incoming flow conditions, three dynamic calculation domains of unsteady state, convection and viscosity are all established according to the wall surface boundary. Disturbance carrying flow field change information is generated from the surface of an object, and a unit dominated by viscosity is also close to the surface of the object, so that two types of dynamic calculation domains of unsteady and convection take a plurality of layers of adjacent units of a wall surface as initial units, for example, 5-10 layers of adjacent grid units of the wall surface can be taken, and the dynamic calculation domain of viscosity only takes 1 layer of units close to the wall surface as initial units.

(2) When the method is initialized according to a given flow field, three types of dynamic calculation domains of unsteady state, convection and viscosity are established according to the flow characteristics of the given flow field. The initial units of the unsteady and convective dynamic calculation domains are disturbed units with flow characteristics in a given flow field inconsistent with incoming flow conditions, and the initial unit of the viscous dynamic calculation domain is a viscous effect leading unit in the convective dynamic calculation domain.

The initial units of the convection dynamic calculation domain 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 disturbed units meet the following conditions:

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

wherein, W_∞A conservation quantity representing an incoming flow condition; | Δ W⁽¹⁾||_maxRepresenting the maximum value of the 1 st iteration step of the conservative updating quantity of all units in the convection dynamic calculation domain at the 1 st physical moment; epsilon_a,cThe newly added threshold value of convection is taken as 10^-4≤ε_a,c≤10^-6。

The initial unit of the unsteady dynamic computation domain is consistent with the initial unit of the convection dynamic computation domain.

The initial unit of the viscous dynamic calculation domain is the viscous effect dominant unit in the convective dynamic calculation 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:

The fourth step: and allocating a storage space.

Dividing data required to be stored for solving 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 a preset calculation domain, wherein the information comprises grid coordinates and flow field variables of the current moment and the last two moments; 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 one row in three types of dynamic computation domains, including: a. a single linked list-records the I label range of the line unit in the current iteration step in the updating of the three 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 three 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 three types of dynamic calculation domains are updated; c. a real two-dimensional array storing the second type of data of the row of cells in the convective dynamic domain, corresponding to the range of the I-label 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 conventional method, all data of all grid cells in a static preset calculation domain are stored in a static data structure. In contrast, the iterative perturbation domain updating method provided by the invention divides the data to be stored into two types, namely the inherent information of all units and the 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 stream dynamic computation domain is stored in a dynamic data structure, and the ninth step and the fourteenth step can adjust the storage space of the second kind of data according to the unsteady dynamic computation domain in time at each iteration step and each physical moment. Therefore, compared with the existing method, the iterative propulsion disturbance domain updating method provided by the invention can effectively reduce the memory requirement of the aircraft dynamic aerodynamic characteristic numerical simulation.

The fifth step: and (5) processing boundary conditions.

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

And a sixth step: and estimating a residual error.

In order to effectively avoid invalid calculation in the conventional unsteady numerical simulation method and reduce the calculation amount of a single iteration step, the method divides the residual error term of the flow control equation into an inviscid term and a viscous term, solves the inviscid term of the residual error term only in the convection dynamic calculation domain, and solves the viscous term of the residual error term in the viscous dynamic calculation domain.

The flow control equation is expressed as:

wherein τ represents a virtual time; Δ t represents a physical time step; f_cDenotes the convective flux, F_vDenotes viscous flux, Q_TThe source terms of the turbulence model equation are represented,

is a source term of the iterative boosting method; | Ω | represents 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; superscripts (n-2), (n-1) and (n) respectively represent the last moment, the last moment and the current moment; the right-hand term of equation (4) 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 sticky term.

The seventh step: and (4) time integration.

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

In the existing method, the most time-consuming residual estimation and time integration steps are required to be performed on all grid units in a static preset calculation domain, and the total solution time is about 99%. In contrast, in the above iterative perturbation domain updating method provided by the present invention, the residual estimation in the sixth step only calculates the non-sticky term of the residual in the convection dynamic calculation domain, and only calculates the sticky term of the residual in the sticky dynamic calculation domain; the time integration of the seventh step is also performed only in the convection dynamics calculation domain. Therefore, the iterative propulsion disturbance domain updating method provided by the invention can effectively improve the calculation efficiency of the aircraft dynamic aerodynamic characteristic numerical simulation.

Eighth step: and increasing two types of dynamic calculation domains of unsteady and convective.

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 (5)

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 (5), 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 and the unsteady dynamic calculation domain, returning to the eighth 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 (5), returning to the eighth step, and judging the next boundary unit of the convection dynamic calculation domain; and after traversing all boundary units of the convection dynamic calculation domain, executing the ninth step.

u·q+a＞0 (6)

wherein u represents a flow velocity vector, and a represents a sound velocity; and q is a unit vector of a lattice point of the boundary unit of the convection dynamic calculation domain pointing to a certain lattice point of the boundary unit, and if q satisfies the formula (6), an adjacent unit sharing the lattice point with the boundary unit of the convection dynamic calculation domain is brought into the convection dynamic calculation domain and the unsteady dynamic calculation domain. The unsteady dynamic computation domains should always contain the convection dynamic computation domain.

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

Judging whether the unsteady dynamic calculation domain is increased or not; if yes, reallocating the storage space of the second type of data according to the unsteady dynamic calculation domain, and then executing the tenth step; if not, executing the tenth step.

The tenth step: and reducing two types of dynamic calculation domains of convection and viscosity.

In order to reduce the calculation amount of a single iteration step as much as possible without influencing the convergence rate, the convection 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 (7)

wherein epsilon_dIndicating the deletion threshold, take ε_d≤10^-7(ii) a If the boundary cell satisfies the formula (7), 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 (8), 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 (9)

wherein, | Δ (Δ W) | | represents a module value of the boundary unit conservative update quantity affected by the adjacent unit of the boundary unit in the convection dynamic calculation domain; Δ 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; no viscosity to subsonic velocityA unit for calculating influence delta R of residual error term of boundary unit along i direction next to the unit in flow dynamic calculation 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 unit satisfies the formula (9) and the formula (10), the boundary unit simultaneously satisfies four conditions of being converged, being located in the compressible flow, being located at the most upstream and no longer being influenced by other units in the convection dynamic calculation domain, removing the boundary unit from the convection dynamic calculation domain, if the boundary unit simultaneously exists in the viscous dynamic calculation domain, removing the boundary unit from the viscous dynamic calculation domain to ensure that the convection dynamic calculation domain always contains the viscous dynamic calculation domain, returning to the tenth step, and judging the next boundary unit of the convection dynamic calculation domain; if the boundary cell does not satisfy the formula (9) and the formula (10), the tenth step is returned, and the next boundary cell in the dynamic calculation domain of the convection current is judged. And executing the eleventh step after traversing all boundary units of the convection dynamic calculation domain.

The eleventh 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 formulas (13) and (14), the unit positioned in the flow dynamic calculation domain in the adjacent unit of the boundary unit is taken into the viscous dynamic calculation domain, the tenth step 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 formulas (13) and (14), returning to the tenth 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 twelfth step is executed. 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 twelfth step: judging whether the inner iteration reaches a convergence condition; if yes, executing the thirteenth step; if not, returning to the fifth step, entering the calculation of the next iteration step of the internal iteration, and repeatedly executing the fifth step to the twelfth step.

The thirteenth step: and reducing the unsteady dynamic calculation domain.

The influence of the unsteady effect is measured by utilizing the physical time derivative, whether all boundary units of the unsteady dynamic calculation domain are not dominated by the unsteady effect is judged one by one, namely whether all boundary units of the unsteady dynamic calculation domain meet the following requirements is judged one by one:

if the boundary unit of the unsteady dynamic calculation domain satisfies the formula (15), removing the boundary unit from the unsteady dynamic calculation domain, returning to the thirteen step, and judging the next boundary unit of the unsteady dynamic calculation domain; if the boundary unit of the unsteady dynamic calculation domain does not satisfy the formula (15), returning to the thirteenth step, and judging the next boundary unit of the unsteady dynamic calculation domain; and executing the fourteenth step after traversing all boundary units of the unsteady dynamic computation domain.

In the eighth step and the tenth step, only a disturbed unit and a viscous effect leading unit are respectively added into a convection dynamic calculation domain and a viscous dynamic calculation domain according to disturbance propagation characteristics; in the tenth step, the converged units are removed from the solution update in time on the premise of not influencing the convergence rate by using the solution convergence characteristics; in the tenth step, the influence of the unsteady effect is measured, and the unsteady unit is removed from the solution updating in time; therefore, the iterative propulsion disturbance domain updating method provided by the invention can always solve only the non-convergent non-stationary unit, only the viscosity effect is considered in the viscous effect leading unit, and the invalid calculation in the existing non-stationary numerical simulation method can be effectively avoided.

The fourteenth step is that: and reallocating the storage space of the second type data.

Judging whether the unsteady dynamic calculation domain is reduced or not; if yes, reallocating the storage space of the second type of data according to the unsteady dynamic calculation domain, and then executing the fifteenth step; if not, executing the fifteenth step.

The fifteenth step: resetting two dynamic calculation domains of convection and viscosity.

Updating the flow field variables of the last two moments; and resetting the convection dynamic calculation domain to the range of the unsteady dynamic calculation domain, and resetting the viscosity dynamic calculation domain to the viscosity effect dominant region in the convection dynamic calculation domain.

Sixteenth, step: judging whether the specified final state physical time is reached; if yes, executing the seventeenth step; if not, returning to the fifth step, entering the calculation of the next physical moment, and repeatedly executing the fifth step to the sixteenth step.

Seventeenth step: and outputting the result.

Simulation is performed below for the NACA0012 airfoil pitch oscillation problem at mach number 0.755 and the winged missile pitch motion problem at mach number 1.58, respectively.

(1) A simulation was performed of the NACA0012 airfoil pitch oscillation problem at mach number 0.755.

Fig. 2 illustrates the unsteady dynamic computation domain and the minimum convection dynamic computation domain in the internal iteration corresponding to the above iterative perturbation domain updating method provided by the present invention at different physical moments. In FIG. 2,. eta._uRepresenting the ratio of the grid quantities, η, of the unsteady dynamic computational domain to the predetermined computational domain_c,minAnd representing the ratio of the grid quantity of the inner iteration minimum convection dynamic calculation domain and the preset calculation domain. The flow field is initialized according to the incoming flow condition, the disturbance is generated from the wall surface boundary, and 10 layers of adjacent units of the wall surface are taken as initial units of the unsteady dynamic calculation domain. At the 1 st physical moment, the unsteady dynamic computation domain extends from the wall surface to the whole preset computation domain and keeps in the range all the time. In the inner iteration of all the time, the convection dynamic calculation domain is reset according to the unsteady 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 area and a downstream area thereof are left.

Fig. 3 shows the unsteady and convective dynamic computation domain grid volume change curves of the iterative perturbation domain updating method provided by the present invention at the last 6 physical moments, and compares the unsteady and convective dynamic computation domain grid volume change curves with the convergence curve of the existing global updating method. The grid quantity variation curve of the convection dynamic calculation domain shown in FIG. 3 shows that the convection dynamic calculation domain is monotonically reduced, eta, in one iteration_cFrom 1.0 down to 0.44, i.e. at the end of the inner iteration, only 44% of the pre-set computing field units participate in the solution. In fig. 2, a comparison between the convergence curves of the existing global update method and the above iterative-boosting perturbation-domain update method provided by the present invention shows that the present inventionThe provided iteration propulsion disturbance domain updating method not only can not generate adverse effect on the convergence rate of the unsteady numerical simulation, but also can reduce the number of solving steps required by the inner iteration to reach the convergence state. Compared with the existing method, the iterative boosting perturbation domain updating method provided by the invention can save 26.5% of time.

(2) A simulation was performed of the winged missile pitch problem at mach 1.58.

Fig. 4 shows the unsteady dynamic computation domain and the minimum convection dynamic computation domain in the internal iteration, which correspond to the above iterative perturbation domain updating method provided by the present invention at different physical moments. For supersonic flow, the unsteady dynamic calculation domain is not always kept in the whole preset calculation domain range as the subsonic velocity condition, but shows periodic change along with the pitching motion of the missile. When the 1 st physical moment begins, the unsteady dynamic computation domain is established according to the wall surface, and only 10 layers of units adjacent to the wall surface boundary are taken as initial units; and at the end of the 1 st physical moment, the eta of the unsteady dynamic calculation domain_uHas reached 0.999. In the first half period, the unsteady dynamic calculation domain keeps a larger range all the time; then, the unsteady dynamic computation domain repeats the reduction and the increase every half cycle. When the missile returns to the initial state, the unsteady dynamic computation domain begins to shrink. The reduction process of the unsteady dynamic calculation domain is firstly rapid and then slow, and the reduction is slower as the approach angle is closer to the maximum value of the attack angle after a small half period. In contrast, the increase from the minimum range to the maximum range is very rapid and is only done in one moment.

Fig. 5 shows the unsteady and convective dynamic computation domain grid volume change curves of the iterative perturbation domain updating method provided by the present invention at the last several physical moments, and the change curves are compared with the convergence curve of the existing global updating method. The comparison result shows that the iterative propulsion disturbance domain updating method provided by the invention is the same as the existing global updating method, and the solution convergence is realized in each internal iteration step. For the supersonic three-dimensional problem, the acceleration effect of the invention is more obvious; when the iterative update disturbance domain updating method is finished solving, the existing method still has 16 moments to be calculated. By virtue of the computation amount saving and the iteration step number saving obtained by reducing the unsteady and convective dynamic computation domains, the iterative boost perturbation domain updating method provided by the invention can obtain 38.9% of time saving compared with the existing method.

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. An iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation is characterized by comprising the following steps:

s24: and outputting the result.

2. The iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation according to claim 1, wherein step S2, based on the incoming flow condition or based on the given flow field, initializes the flow field in the preset calculation domain, specifically including:

3. The iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation according to claim 1, wherein step S3 is to establish an unsteady dynamic calculation domain, a convective dynamic calculation domain and a viscous dynamic calculation domain according to a flow field initialization mode, and specifically includes:

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

wherein W represents a conservative amount; w_∞A conservation quantity representing an incoming flow condition; | Δ W⁽¹⁾||_maxRepresenting the maximum value of the 1 st iteration step of the conservative updating quantity of all units in the convection dynamic calculation domain at the 1 st physical moment; epsilon_a,cThe newly added threshold value of convection is taken as 10^-6≤ε_a,c≤10^-4(ii) a The initial unit of the unsteady dynamic computation domain is consistent with the initial unit of the convection dynamic computation domain; the initial unit of the viscosity dynamic calculation domain is a viscosity effect dominant unit in the convection dynamic calculation domain, and satisfies the following conditions:

4. The iterative propulsion disturbance domain updating method for the dynamic aerodynamic characteristic simulation of the aircraft according to claim 3, wherein in step S6, the residual terms of the flow control equation are divided into two types, namely a non-viscous term and a viscous term; wherein the sticky-free terms comprise source terms of convection flux and iterative marching methods, the sticky-free 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 iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation according to claim 4, wherein step S8 is performed to determine whether all boundary cells of the convective dynamic computation domain have been subjected to inviscid disturbance one by one, and specifically includes:

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

u·q+a＞0 (6)

6. The iterative-propulsion-disturbance-domain updating method for aircraft dynamic aerodynamic characteristic simulation according to claim 5, wherein step S12 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; if yes, go to step S13; if not, returning to the step S12, and judging the next boundary unit of the convection dynamic calculation domain; the method specifically comprises the following steps:

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

s123: 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 S124; if not, returning to the step S121, and judging a next boundary unit of the convection dynamic calculation domain;

||ΔW||+||Δ(ΔW)||＜ε_d (9)

wherein, Δ F_cRepresenting the amount of change in convective flux; the subscript i +1 indicates the convection flowThe immediate vicinity of the boundary cell of the state computation domain along the positive i direction, the subscript i-1 denotes the immediate vicinity of the boundary cell of the convection 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 iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation according to claim 6, wherein step S15 is performed to determine whether all boundary cells of the viscous dynamic calculation domain are dominated by viscous effect one by one, and specifically includes:

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 1 st iteration step of the internal iteration in the current moment to the maximum value of the single-step conservative quantity updating quantity in the current iteration stepIs measured.

8. The iterative propulsion disturbance domain updating method for aircraft dynamic aerodynamic characteristic simulation according to claim 7, wherein in step S18, it is determined one by one whether the unsteady effect is no longer dominant for all boundary elements of the unsteady dynamic computation domain, specifically including: