A method for evaluating dynamic stability of a hypersonic vehicle

Abstract

The application discloses a kind of hypersonic vehicle dynamic stability evaluation methods, scheme includes: for the aircraft to be evaluated, the numerical simulation working condition including center of gravity position and reference length is constructed;Under the condition that oscillation excitation is not applied, steady-state flow field calculation is carried out to initialize transition model;The coupling solution framework of unsteady RANS solver and intermittent factor type transition model is established, and the aerodynamic moment time history under pitch oscillation is obtained by the synchronous iteration of each time step;Real-time extraction of transition position is carried out in oscillation period, and its periodic evolution characteristics are analyzed;Adaptive identification method combining multi-harmonic fitting and energy criterion is used to extract dynamic derivative;Finally, based on the coupling relationship between transition position and center of gravity position, the dynamic stability risk level is comprehensively judged, and the risk map is output.The application can realize the time-domain coupling evaluation of transition evolution and dynamic response, and improve the accuracy and engineering applicability of dynamic stability prediction under hypersonic condition.

Description

Technical Field

[0001] This application relates to the field of dynamic stability assessment technology for hypersonic vehicles, and in particular to a method for assessing the dynamic stability of hypersonic vehicles. Background Technology

[0002] In the field of dynamic stability assessment for hypersonic vehicles, engineers sometimes rely on the dynamic derivative as a key parameter to characterize aerodynamic damping characteristics for the purpose of accurately predicting the attitude stability of the vehicle during reentry. For example, they obtain the sensitivity coefficient of pitch moment to angle of attack and its rate of change through forced oscillation wind tunnel tests or numerical simulations. In existing technologies, the test team usually provides oscillation data of the model at a preset frequency / amplitude, or the numerical simulation team calculates the aerodynamic moment time history based on the RANS / URANS framework, and then uses Fourier decomposition or least squares regression to extract the dynamic derivative, thus providing a quantitative basis for stability criteria. However, under typical reentry conditions at high Reynolds numbers and high Mach numbers, the transition of the boundary layer from laminar to turbulent flow is often unavoidable. This transition significantly changes wall friction, heat / momentum transport mechanisms, and the intensity of shock-boundary layer interaction, thereby affecting the phase relationship between aerodynamic forces and moments, as well as the equivalent damping. However, in existing evaluation schemes, the common practice is to simply assume that the boundary layer is "fully turbulent" or "fully laminar," or to use only conventional turbulence models for calculation. As a result, the transition effect is equated to noise or nonlinear error in the dynamic derivative identification stage. This makes it difficult to explain the sudden dynamic instability phenomena that occur within certain height / Reynolds number windows, such as the situation where dynamic instability is significantly amplified when the transition region coincides with the center of gravity. Summary of the Invention

[0003] This specification provides an embodiment of a method for evaluating the dynamic stability of a hypersonic vehicle to solve at least one of the technical problems mentioned above.

[0004] To solve the above-mentioned technical problems, the embodiments in this specification are implemented as follows: According to an embodiment of the present invention, a method for evaluating the dynamic stability of a hypersonic vehicle is provided, comprising: For the hypersonic vehicle to be evaluated, a numerical simulation condition is constructed. The numerical simulation condition includes input vehicle geometry and reference parameters, external flow conditions, wall boundary conditions, and pitch oscillation excitation parameters. The vehicle geometry and reference parameters include the vehicle's center of gravity position and reference length. Based on the numerical simulation conditions, steady-state or quasi-steady-state flow field calculations are performed without applying pitch oscillation excitation defined by the pitch oscillation excitation parameters to obtain the initial flow field, and the relevant variables of the transition model are initialized to match the initial flow field. A coupled solution framework is established between the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor transition model. The flow field variables and transition state variables are updated synchronously in each time step and iterated until convergence, thereby obtaining the aerodynamic torque time history of the aircraft under the pitch oscillation excitation. Within the oscillation period, the transition position is extracted in real time along a specified path on the aircraft wall to obtain a sequence of transition positions changing over time. Periodic analysis is then performed on this sequence to obtain migration characteristic parameters characterizing the periodic evolution of the transition position. The migration characteristic parameters include the periodic average transition position, migration amplitude, and migration phase. The obtained aerodynamic moment time history is preprocessed to remove the initial transient period, and the dynamic derivative or equivalent damping index of the aircraft is extracted by an adaptive identification method that combines multi-harmonic fitting, energy criteria, or a combination thereof. Based on the periodic average transition position and the aircraft's center of gravity position, a transition-center of gravity coupling index is constructed. Combined with the dynamic derivative or equivalent damping index, the dynamic stability risk level of the aircraft under the current operating conditions is comprehensively determined, and a risk window or risk map is output.

[0005] In some alternative implementations, the pitch oscillation excitation parameters include the initial angle of attack. Oscillation amplitude oscillation frequency and initial phase The variation law of the aircraft's pitch angle or equivalent angle of attack over time is as follows: ;in, Define the angular frequency; define the reduced frequency. ,in, For reference length, The incoming flow velocity; The aircraft geometry and reference parameters include reference length. Reference area, center of gravity position And the definition of the coordinate system.

[0006] In some optional implementations, obtaining the migration characteristic parameters characterizing the periodic evolution of the transition position includes: Define the transition point along the specified wall path. ,in, For the distribution of intermittent factors along the wall path, The threshold value is a preset value used to determine the transition. For the transition position Harmonic fitting is performed, expressed as:

[0007] in, The oscillation angular frequency is used to obtain the periodic average transition position. Migration amplitude and migration phase ; in, , Calculated using orthogonal projection:

[0008] The oscillation period is [the period of time].

[0009] In some optional implementations, the transition determination threshold is... The transition location is determined by one or more combinations of the intermittent factor threshold, the heat flux change threshold along the friction, and the friction change threshold along the friction; and the transition location is extracted along multiple azimuth paths or multiple meridians to form the spatial distribution characteristics of the transition migration.

[0010] In some optional implementations, the adaptive identification method employing multi-harmonic fitting, energy criteria, or a combination thereof to extract the dynamic derivative or equivalent damping index of the aircraft includes: The pitch moment coefficient is obtained by using multi-harmonic fitting. Represented as:

[0011] in, For the fitting order; when When, the static stability derivative is obtained. and the combined dynamic derivative ,in, For pitch angular velocity damping derivative, The derivative of the rate of change of angle of attack, This represents the oscillation amplitude. And / or, using the energy criterion, define a period. Equivalent damping work within ,in, For pitch angular velocity, when It is determined to be dynamically stable when It is determined to be dynamically unstable.

[0012] In some optional implementations, the preprocessing of the obtained aerodynamic moment time history further includes: calculating the proportion of fundamental frequency components, the proportion of harmonic energy, or the regression residual to form an identification quality index; and adaptively selecting, based on the identification quality index, to use multi-harmonic fitting or energy criteria, or to use both for cross-validation.

[0013] In some optional implementations, the construction transition-centroid coupling index includes: Define the periodic average transition position and the center of gravity position. normalized distance

[0014] in, For reference length; according to With dynamic derivative The overall risk level output includes: when Less than or equal to the preset risk assessment threshold and At that time, it was determined to be a high-risk, dynamically unstable window; when Greater than the preset threshold for risk assessment and At that time, it was determined to be a low-risk, dynamically stable window; Other situations are classified as medium risk.

[0015] In some optional implementations, the output risk window or risk map includes associating the risk level with corresponding height, Reynolds number, wall temperature, oscillation amplitude, and frequency parameters to form a risk window or risk map.

[0016] In some optional implementations, the transition model is an intermittent factor transition model, which is solved synchronously with the unsteady Reynolds-averaged Navier-Stokes solver at each time step. It is used to update the transition state in real time during the unsteady process, obtain the wall intermittent factor distribution by solving the transport equation containing the intermittent factor, and then determine the transition location and its evolution over time.

[0017] In some optional implementations, the establishment of a coupled solution framework between the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor transition model, which synchronously updates the flow field variables and transition state variables at each time step and iterates until convergence, includes: completing closed-loop iterations of flow field advancement, transition variable update, turbulent viscosity correction, and boundary layer state update at each time step until the convergence criteria for each step are met.

[0018] One embodiment of this specification can achieve at least the following beneficial effects: The technical solution of this application first constructs a numerical simulation condition for the hypersonic vehicle to be evaluated, including the center of gravity position and reference length, to establish complete input conditions and spatial references for subsequent dynamic stability analysis. Then, steady-state flow field calculations are performed under conditions without oscillatory excitation to initialize the transition model, ensuring the physical continuity of the evolution of the transition state in the subsequent unsteady solution process. Furthermore, a coupled solution framework is established between the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor-type transition model, synchronously updating the flow field variables and transition state quantities at each time step, thereby obtaining the aerodynamic moment time history under pitch oscillation excitation. This coupling method differs from the existing approach of using transition as a fixed position or steady-state result input for dynamic calculation, enabling consistent coupling between transition evolution and aerodynamic response in the time domain. Based on this, the transition position is extracted in real time along a specified path on the wall within the oscillation period and subjected to periodic analysis to obtain migration characteristic parameters such as the period-averaged transition position, migration amplitude, and migration phase. This transforms the previously invisible transition dynamic process into a quantitative indicator that can be used for engineering characterization. Then, an adaptive identification method combining multi-harmonic fitting and energy criteria is used to extract dynamic derivatives or equivalent damping indices. This effectively addresses the nonlinearity and higher-order harmonic enhancement problems caused by transition migration, improving the repeatability and robustness of derivative identification. Finally, a transition-centroid coupling index is constructed based on the period-averaged transition position and the center of gravity position. Combined with the dynamic derivative, the dynamic stability risk level is comprehensively determined, and a risk window or risk map is output. This determination mechanism incorporates the spatial relationship between the transition position and the center of gravity as a key criterion into the evaluation process, expanding the dynamic stability conclusion from a single derivative output to a risk expression with clear engineering orientation, facilitating the rapid identification of high-risk working condition windows. Attached Figure Description

[0019] To more clearly illustrate the technical solutions in the embodiments or prior art of this specification, the drawings used in the description of the embodiments or prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments recorded in this application. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0020] Figure 1 This is a flowchart of a method for evaluating the dynamic stability of a hypersonic vehicle provided by the present invention; Figure 2 This paper presents a schematic diagram illustrating the geometric parameters of the blunt cone model used as an embodiment in the technical solution of this application. Figure 3 Given Figure 2 The diagram shows the computational mesh of the blunt cone model in the figure, which is used to illustrate the mesh topology and spatial discretization method used in the numerical simulation. Figure 4 It displays data on the change of angle of attack over time as measured in flight tests; Figure 5 The paper demonstrates the distribution of wall intermittency factor obtained by numerical calculation at three different altitudes of 20km, 24km, and 28km in the simulation experiment. Detailed Implementation

[0021] To make the objectives, technical solutions, and advantages of one or more embodiments of this specification clearer, the technical solutions of one or more embodiments of this specification will be clearly and completely described below in conjunction with specific embodiments and corresponding drawings. Obviously, the described embodiments are only a part of the embodiments of this specification, and not all of them. Based on the embodiments in this specification, all other embodiments obtained by those skilled in the art without creative effort are within the protection scope of one or more embodiments of this specification.

[0022] It should be understood that although the terms first, second, third, etc., may be used in this application to describe various information, this information should not be limited to these terms. These terms are only used to distinguish information of the same type from one another.

[0023] To address the technical problems in the prior art mentioned in the background section, this application establishes an evaluation framework that synchronously couples the boundary layer transition model with the flow solution at each time step during the unsteady numerical simulation of forced pitch oscillation of a hypersonic vehicle. This framework tracks and extracts the periodic evolution characteristics of the transition position over time in real time, combining this with dynamic derivatives identified by multi-harmonic fitting or energy criteria based on aerodynamic moment time histories. This constructs a quantitative coupling index between the transition position and the vehicle's center of gravity, ultimately outputting a risk window or risk map that encompasses dynamic stability assessment and risk level. This enables integrated analysis from transition dynamic evolution to quantitative stability evaluation.

[0024] The technical solution of this application will be described in detail below based on the accompanying drawings, such as... Figure 1 As shown, Figure 1 A flowchart illustrating a method for evaluating the dynamic stability of a hypersonic vehicle, as provided in this application, is shown below. Figure 1 As shown, the method may include: Step 102: For the hypersonic vehicle to be evaluated, construct a numerical simulation condition. The numerical simulation condition includes inputting the vehicle's geometry and reference parameters, external flow conditions, wall boundary conditions, and pitch oscillation excitation parameters. The vehicle's geometry and reference parameters include the vehicle's center of gravity position and reference length.

[0025] In the embodiments of this specification, constructing a numerical simulation case for the hypersonic vehicle to be evaluated can be understood as establishing a complete set of input conditions and computational environment for subsequent dynamic stability analysis. This step aims to transform the actual physical parameters of the vehicle, the flight environment, and the preset disturbance excitations into mathematical expressions that the numerical solver can recognize and process. Specifically, the inputs involved in constructing the numerical simulation case can first include the vehicle's geometry and reference quantities. The vehicle geometry is used to generate the computational mesh or define the object surface boundary, while the reference quantities can include basic parameters for dimensionless processing, such as reference length, reference area, and coordinate system definition. The geometry and reference quantities can explicitly include the vehicle's center of gravity position, because the center of gravity is a key reference point for subsequent evaluation of transition effects and calculation of pitching moment. Its relative relationship with the aerodynamic load distribution directly determines the magnitude and sign of the moment. In practical applications, constructing the numerical simulation case also requires inputting external flow conditions. These conditions can include parameters such as Mach number, Reynolds number per unit length, incoming flow temperature, pressure, density, and incoming flow turbulence intensity, which together describe the hypersonic flow field environment in which the vehicle operates. Furthermore, wall boundary conditions are also an essential component of the operational condition construction. These can be isothermal wall temperatures, adiabatic walls, or given heat flux distributions, used to simulate the heat exchange process on the aircraft surface under different thermal protection conditions. Finally, the input of pitch oscillation excitation parameters is to define the motion modes used for dynamic stability assessment. These parameters can include at least the initial angle of attack, oscillation amplitude, oscillation frequency, and oscillation duration. They specify how the aircraft performs periodic pitch motion around a certain average attitude in the numerical simulation, such as a sinusoidally varying angle of attack input, thus providing an excitation source for subsequent extraction of unsteady aerodynamic responses and identification of dynamic derivatives.

[0026] Step 104: Based on the numerical simulation conditions, perform steady-state or quasi-steady-state flow field calculations without applying pitch oscillation excitation defined by the pitch oscillation excitation parameters to obtain the initial flow field, and initialize the relevant variables of the transition model to match the initial flow field.

[0027] In the embodiments described in this specification, after constructing the numerical simulation conditions, the first step is to solve the steady-state or quasi-steady-state flow field without applying any pitch oscillation excitation. At this time, the aircraft maintains a fixed attitude, and the solver iteratively calculates based on the given inflow conditions and wall boundary conditions until convergence, thereby obtaining an initial flow field distribution that matches the current flight state. Based on this, the relevant variables in the transition model need to be initialized. For example, in the intermittent factor type transition model, the transport variables such as the intermittent factor γ, fluctuating kinetic energy k, and its specific dissipation rate ω need to be assigned initial values ​​based on the local flow characteristics of the initial flow field or obtained by solving the steady-state transition transport equation, so that they achieve a state of complete matching with the initial flow field. Through the above two-stage preprocessing process, it can be ensured that when forced oscillation excitation is applied and unsteady coupling solution is started, the evolution of the transition state does not start from zero or make arbitrary assumptions, but is based on a steady state with physical continuity. This avoids the introduction of additional numerical transient disturbances due to initial field mismatch and improves the computational stability and result reliability of the entire dynamic stability assessment process.

[0028] Step 106: Establish a coupled solution framework of unsteady Reynolds-averaged Navier-Stokes solver and intermittent factor transition model, synchronously update flow field variables and transition state variables in each time step, and iterate until convergence, thereby obtaining the aerodynamic torque time history of the aircraft under the pitch oscillation excitation.

[0029] In the embodiments of this specification, the coupled solution framework of the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor transition model can refer to a numerical computation architecture established to achieve synchronous solution of the flow control equations and transition transport equations in the time domain. Specifically, in unsteady calculations under pitch oscillation excitation, the updates of flow field variables and transition state variables need to be completed simultaneously in each physical time step, forming a tight coupling relationship rather than a sequential unidirectional transfer. In the actual calculation process, the compressible Navier-Stokes equations are first solved based on the current flow field state to obtain initially updated flow field variables such as density, velocity, and pressure. Then, based on the updated flow field information, the transport equations contained in the intermittent factor transition model are solved synchronously. On this basis, the updated intermittent factor is used... For effective viscosity coefficient A correction is made, and this corrected viscosity coefficient is fed back into the turbulence term in the Navier-Stokes equations, thus affecting the flow field solution in the next iteration. This process is repeated within the same time step until both the flow field variables and transition state variables meet the preset convergence criteria before proceeding to the next physical time step. This tightly coupled solution method, which realizes the real-time evolution of the transition state in each time step, can accurately capture the dynamic migration of the transition position caused by the periodic change in angle of attack during forced pitch oscillations. Ultimately, it obtains the time-varying history data of the aerodynamic torque of the aircraft under a given oscillatory excitation, providing a fundamental input for subsequent dynamic stability analysis.

[0030] Step 108: During the oscillation period, the transition position is extracted in real time along the specified path on the aircraft wall to obtain a sequence of transition positions changing over time. Periodic analysis is performed on the sequence to obtain migration characteristic parameters characterizing the periodic evolution of the transition position. The migration characteristic parameters include the periodic average transition position, migration amplitude, and migration phase.

[0031] In the embodiments of this specification, the real-time extraction and periodic analysis of the transition position can refer to a series of processing methods used in the unsteady calculation of forced pitch oscillations to transform the originally continuous transition evolution process within the boundary layer into discrete, quantifiable characteristic parameters. Specifically, at each physical moment or every several time steps, the solver needs to determine the transition position at that moment along a predetermined path on the aircraft wall based on the current intermittent factor distribution. For example, the intermittent factor along the path can be used to determine the transition position when it first reaches a predetermined threshold. The flow direction coordinate (e.g., a value of 0.5) is defined as the instantaneous transition position at that moment. By repeating this extraction operation throughout the entire oscillation period, a discrete sequence showing the continuous time-varying transition position can be obtained. Based on this, periodic mathematical analysis can be performed on the sequence to extract migration characteristic parameters with engineering significance. Typically, harmonic fitting methods can be used to approximate the time-varying law of the transition position as a fundamental frequency oscillation, i.e. It is decomposed into a superposition of sine and cosine components with the same frequency as the pitch motion. Through this fitting process, the periodic average transition position can be further calculated. Migration amplitude and migration phase These are the three core migration feature parameters; among them, the periodic average transition position. It can characterize the average axial position and migration amplitude of the transition region within a complete oscillation cycle. It can quantify the degree to which the transition position oscillates around this average value, and the shift phase. It can reflect the lag or lead of the transition position response relative to the pitch motion excitation. Through the above extraction and analysis process, the transition dynamic process originally hidden inside the flow field is transformed into a set of concise numerical indicators, which can provide direct quantitative basis for explaining the variation law of dynamic derivatives from the perspective of physical mechanism and constructing the relationship between transition and center of gravity position.

[0032] Step 110: Preprocess the obtained aerodynamic moment time history, remove the initial transient period, and use an adaptive identification method that combines multi-harmonic fitting, energy criteria, or a combination thereof to extract the dynamic derivative or equivalent damping index of the aircraft.

[0033] In the embodiments of this specification, since non-physical transient oscillations may exist in the initial stage of unsteady calculations due to the transition of the initial flow field, the original torque time history can be preprocessed first to remove data from the first few cycles to eliminate the influence of initial transients and ensure that subsequent analysis is based on a stable periodic response. Based on this, a multi-harmonic fitting method can be used to expand the torque time history into a Fourier series form, and the static stability derivative and combined dynamic derivative can be extracted by fitting the fundamental frequency and higher-order harmonic components. Simultaneously, an energy criterion can be used to directly determine the sign and magnitude of the equivalent damping by calculating the integral of the pitching moment with respect to the pitching angular velocity over a complete cycle. In practical applications, based on quality indicators such as the proportion of the fundamental frequency component or the regression residuals obtained in the preprocessing stage, multi-harmonic fitting, the energy criterion, or a combination of both can be adaptively selected for cross-validation to improve the robustness and reliability of derivative identification under nonlinear responses caused by transition migration.

[0034] Step 112: Based on the periodic average transition position and the aircraft's center of gravity position, construct a transition-center of gravity coupling index, and combine it with the dynamic derivative or equivalent damping index to comprehensively determine the dynamic stability risk level of the aircraft under the current operating conditions, and output a risk window or risk map.

[0035] This step quantitatively correlates the transition dynamic features extracted in the previous steps with the aircraft's center of gravity position, and based on this, makes a final determination of the dynamic stability under the current flight conditions. Specifically, it first needs to be based on the obtained periodic average transition position. relative to the center of gravity of the aircraft We can construct a normalized coupling index to characterize the spatial proximity of the two components, which can be defined as the ratio of their absolute difference to the reference length. Then, it is combined with the combined dynamic derivatives obtained from the previous steps. Or equivalent damping index, which is comprehensively output according to preset risk assessment rules, for example when Small enough and When the display system is in a positive damping state, it can be determined that the current operating condition has a high risk of dynamic instability. Conversely, when the transition zone is far from the center of gravity and the damping is negative, it is determined to be a low-risk stable state. By correlating the above risk levels with corresponding operating parameters such as flight altitude, Reynolds number, and oscillation frequency, engineering output forms such as risk windows or risk maps can be generated. This transforms complex flow transition effects and attitude stability problems into risk distribution maps that designers can directly view and use, thus providing intuitive quantitative basis for flight envelope planning, center of gravity position optimization, and transition control measures decisions.

[0036] The technical solution of this application first constructs a numerical simulation condition for the hypersonic vehicle to be evaluated, including the center of gravity position and reference length, to establish complete input conditions and spatial references for subsequent dynamic stability analysis. Then, steady-state flow field calculations are performed under conditions without oscillatory excitation to initialize the transition model, ensuring the physical continuity of the evolution of the transition state in the subsequent unsteady solution process. Furthermore, a coupled solution framework is established between the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor-type transition model, synchronously updating the flow field variables and transition state quantities at each time step, thereby obtaining the aerodynamic moment time history under pitch oscillation excitation. This coupling method differs from the existing approach of using transition as a fixed position or steady-state result input for dynamic calculation, enabling consistent coupling between transition evolution and aerodynamic response in the time domain. Based on this, the transition position is extracted in real time along a specified path on the wall within the oscillation period and subjected to periodic analysis to obtain migration characteristic parameters such as the period-averaged transition position, migration amplitude, and migration phase. This transforms the previously invisible transition dynamic process into a quantitative indicator that can be used for engineering characterization. Then, an adaptive identification method combining multi-harmonic fitting and energy criteria is used to extract dynamic derivatives or equivalent damping indices. This effectively addresses the nonlinearity and higher-order harmonic enhancement problems caused by transition migration, improving the repeatability and robustness of derivative identification. Finally, a transition-centroid coupling index is constructed based on the period-averaged transition position and the center of gravity position. Combined with the dynamic derivative, the dynamic stability risk level is comprehensively determined, and a risk window or risk map is output. This determination mechanism incorporates the spatial relationship between the transition position and the center of gravity as a key criterion into the evaluation process, expanding the dynamic stability conclusion from a single derivative output to a risk expression with clear engineering orientation, facilitating the rapid identification of high-risk working condition windows.

[0037] Based on the technical solutions described above, this application also provides some more specific technical solutions, which are described below.

[0038] In an optional embodiment, the pitch oscillation excitation parameters include the initial angle of attack. Oscillation amplitude oscillation frequency and initial phase The variation law of the aircraft's pitch angle or equivalent angle of attack over time is as follows: ;in, Define the angular frequency; define the reduced frequency. ,in, For reference length, The incoming flow velocity; The aircraft geometry and reference parameters include reference length. Reference area, center of gravity position And the definition of the coordinate system.

[0039] In the technical solution of this application, the accurate definition of pitch oscillation excitation parameters is the foundation for implementing forced oscillation numerical simulation and ensuring the repeatability of subsequent dynamic derivative identification. Specifically, these parameters may include at least the initial angle of attack. Oscillation amplitude oscillation frequency and initial phase These factors collectively define how an aircraft performs periodic pitch motions around a given average attitude in numerical simulations. During implementation, the change in the aircraft's pitch angle or equivalent angle of attack over time can be expressed as... ,in The angular frequency, a sinusoidal motion input, can simulate the attitude disturbances that an aircraft may encounter in actual flight, and provides a clear kinematic reference for subsequently extracting static and dynamic derivatives from the aerodynamic moment time history. To facilitate comparative analysis between different flight speeds and reference scales, this embodiment can also introduce a reduced frequency. As a dimensionless frequency parameter, it is defined as follows: In the formula The reference length (which can be the characteristic chord length or the configuration reference length) is used. Let be the incoming flow velocity. Introducing a reduced frequency allows the effect of the oscillation frequency to be decoupled from the incoming flow velocity, thus revealing more clearly the relative contributions of viscous effects and added mass effects in the unsteady aerodynamic response as a function of frequency.

[0040] On the other hand, the aircraft geometry and reference inputs constitute the spatial reference and dimensionless basis for the entire numerical simulation. These inputs can include at least the reference length. Reference area, center of gravity position And the definition of the coordinate system. Among them, the reference length... The reference area is used to dimensionlessly process the calculated aerodynamic and moment coefficients, making the results comparable across different scales or operating conditions. The technical solution in this embodiment takes the center of gravity position into account. The center of gravity serves as a reference point for subsequent calculations of the pitching moment and is also a key element in constructing the transition-center-of-gravity coupling index and evaluating dynamic stability. Essentially, the magnitude and sign of the moment depend on the lever arm of the aerodynamic load distribution relative to the center of gravity. The migration of the transition zone causes dynamic changes in this lever arm. Therefore, in this embodiment, the center of gravity position is explicitly defined as a fixed input during the condition construction phase. The coordinate system definition (e.g., using an absolute coordinate system or a body coordinate system) ensures that the application of the motion boundary conditions and the interpretation of the aerodynamic force output have a consistent spatial reference, thereby guaranteeing complete physical self-consistency throughout the entire evaluation process from input to output.

[0041] In optional embodiments, obtaining migration characteristic parameters characterizing the periodic evolution of transition positions may include: Define the transition point along the specified wall path. ,in, For the distribution of intermittent factors along the wall path, The threshold value is a preset value used to determine the transition. For the transition position Harmonic fitting is performed, expressed as:

[0042] in, The oscillation angular frequency is used to obtain the periodic average transition position. Migration amplitude and migration phase ; in, , Calculated using orthogonal projection:

[0043] The oscillation period is [the period of time].

[0044] In the technical solution of this application, obtaining the migration characteristic parameters characterizing the periodic evolution of the transition position is an intermediate step in realizing the quantitative characterization of the transition dynamic effect. This step first relies on extracting the time-varying sequence of the transition position from the unsteady coupling solution results. The basic definition is as follows: along a specified path on the aircraft wall (e.g., a meridian or a given azimuth profile), at each physical moment... Intermittent factor The distribution along the route reaches the preset threshold for the first time. The smallest The coordinates are determined as the turning point at that moment. ,Right now Among them, the intermittent factor These are field variables obtained by directly solving the transition transport equation. Their values ​​change continuously from 0 to 1, representing the probability or state of the flow being in fully laminar, intermittent transition, and fully turbulent flow, respectively. The preset transition threshold can be selected based on the configuration characteristics and verification data in specific implementations. For example, its value can be 0.5 or a preset value within the range of 0.4 to 0.6. Its physical significance lies in defining the dominant transition point between laminar and turbulent flow. In this embodiment, this definition method can transform the continuous intermittent factor distribution field into a discrete position signal that varies with time, thereby enabling the forward and backward migration behavior of the transition region within the oscillation period to be quantitatively captured and subsequently analyzed.

[0045] When the turning point is obtained Based on this, the technical solution of this embodiment further performs periodic harmonic analysis to extract migration characteristic parameters with engineering significance. Since the forced pitch oscillation motion has a clear periodicity, the response at the transition position is usually dominated by the fundamental frequency component; therefore, it can be approximately expressed as a first-order Fourier series, i.e. ,in The coefficients in this expression represent the angular frequency of the pitch oscillation. and These represent the amplitudes of the components in phase with and orthogonal to the motion in the transition position response, respectively. They can be calculated from the original time history data using the orthogonal projection formula, i.e. , In the formula The oscillation period is given. Based on the above fitting results, the three key migration characteristic parameters mentioned above can be further derived, namely, the average transition position of the period. It reflects the average axial position of the transition zone within one oscillation cycle; migration amplitude It quantifies the degree to which the transition position oscillates around its average value; the migration phase It characterizes the phase lag or lead relationship between the transition position response and the pitch motion excitation. These three parameters together constitute a concise and complete engineering description of the dynamic evolution process of the transition, which can compress the originally complex unsteady intermittent factor field into several scalars with clear physical meaning.

[0046] In an optional embodiment, the transition determination threshold is... The transition location is determined by one or more combinations of the intermittent factor threshold, the heat flux change threshold along the friction, and the friction change threshold along the friction; and the transition location is extracted along multiple azimuth paths or multiple meridians to form the spatial distribution characteristics of the transition migration.

[0047] In this embodiment, the determination of the transition position can adopt a definition method adapted to different flow characteristics. Specifically, the transition judgment threshold... The selection of the threshold value can be determined based on the actual flow characteristics and the output type of the numerical solution, using one or more combinations of the interval factor threshold, the heat flux friction change threshold, or the friction friction change threshold. Among these, the interval factor threshold is the interval factor directly output from the transition model. The distribution along the wall path is used to determine the transition location; for example, it can be set... The threshold value is 0.5 or a preset value within the range of 0.4 to 0.6. When the friction loss intermittent factor first reaches or exceeds this threshold, the location is determined to be the transition point at the current moment. The friction loss abrupt change threshold is determined using the wall heat flow. The transition location is determined by the gradient change characteristics along the flow direction. For example, the location where the heat flux gradient reaches an extreme value or abruptly changes can be used as a marker of transition. Similarly, the friction coefficient abrupt change threshold can be determined by the friction coefficient gradient characteristics along the flow direction. Furthermore, the above three determination methods can be combined according to the specific configuration and incoming flow conditions. For example, the transition location can be finally determined when both the intermittent factor threshold and the heat flux abrupt change threshold are satisfied, thus enhancing the robustness and physical consistency of the criteria. Moreover, to comprehensively capture the influence of three-dimensional effects on the transition process, the extraction of the transition location can be carried out simultaneously along multiple azimuth paths or multiple meridians along the aircraft wall, rather than being limited to a single meridian or profile. This allows for obtaining the distribution of the transition location at different circumferential angles and its changes over time, thus forming the spatial distribution characteristics of the transition migration. These spatial distribution characteristics reflect the three-dimensional morphology of the transition region on the aircraft surface as it evolves with oscillating motion.

[0048] In optional embodiments, the adaptive identification method using multi-harmonic fitting, energy criteria, or a combination thereof to extract the dynamic derivative or equivalent damping index of the aircraft may include: The pitch moment coefficient is obtained by using multi-harmonic fitting. Represented as:

[0049] in, For the fitting order; when When, the static stability derivative is obtained. and the combined dynamic derivative ,in, For pitch angular velocity damping derivative, The derivative of the rate of change of angle of attack, This represents the oscillation amplitude. And / or, using the energy criterion, define a period. Equivalent damping work within ,in, For pitch angular velocity, when It is determined to be dynamically stable when It is determined to be dynamically unstable.

[0050] In this embodiment, to accurately extract dynamic stability characteristics from the aerodynamic moment time history obtained from forced oscillation numerical simulation, a multi-harmonic fitting method is employed. This method can adapt to nonlinear aerodynamic responses caused by complex flow phenomena such as boundary layer transition and shock wave-boundary layer interference. Specifically, the pitching moment coefficient... Represented as time The function is approximated using Fourier series form:

[0051] in, The angular frequency of the pitch oscillation. This is the periodic average value of the torque coefficient. , These are the sine and cosine coefficients of the fundamental frequency component, respectively. , These are the coefficients of the higher-order harmonic components. The fitting order used. When only the fundamental frequency component is considered (i.e., When this occurs, the method degenerates into a traditional linearized identification method, which can then be determined by the fundamental frequency sinusoidal coefficients. With oscillation amplitude The statically stable derivative is obtained directly from the ratio. Simultaneously, the combined dynamic derivative used to characterize the pitch damping properties is extracted. ,in For pitch angular velocity damping derivative, This is the derivative of the rate of change of angle of attack. In practical applications, an appropriate fitting order can be selected based on the degree of nonlinearity of the torque response. This allows for more accurate capture of higher-order harmonic contributions, thereby improving the reliability and repeatability of dynamic derivative identification.

[0052] On the other hand, to avoid errors that may be introduced by the linear assumption and to directly assess dynamic stability from the perspective of energy dissipation, this embodiment also introduces an energy criterion as a supplementary or alternative method. This criterion is calculated by taking a complete oscillation period. Pitching moment coefficient pitch angular velocity The integral of the equivalent damping work is defined. ,Right now From a physical perspective, this integral value is proportional to the net work done by the aerodynamic pitching moment on the aircraft during a single cycle, and its sign directly reflects the trend of increase or decrease in system energy. This indicates that aerodynamic damping consumes vibrational energy, and the system exhibits dynamic stability; if This indicates that the system absorbs or dissipates zero energy from the airflow, corresponding to a dynamically unstable or critical state. This energy criterion does not depend on the harmonic structure of the torque response and is applicable to strongly nonlinear, multi-frequency coupling situations caused by transition migration. It can be cross-validated with the results obtained from multi-harmonic fitting methods.

[0053] In an optional embodiment, the preprocessing of the obtained aerodynamic moment time history may further include: calculating the proportion of fundamental frequency components, the proportion of harmonic energy, or the regression residual to form an identification quality index; and adaptively selecting, based on the identification quality index, to use multi-harmonic fitting or energy criteria, or to use both for cross-validation.

[0054] In this embodiment, after eliminating the initial transient period, spectral analysis and statistical evaluation can be performed on the remaining periodic torque data. Specifically, this can include calculating quantitative indicators such as the proportion of the fundamental frequency component, the proportion of harmonic energy, and the regression residual. The proportion of the fundamental frequency component reflects the relative magnitude of the energy of the component with the same frequency as the pitch oscillation in the aerodynamic torque response. When this proportion is close to 1, it indicates that the response is dominated by linear components. The proportion of harmonic energy measures the weight of higher-order harmonic components in the total response energy. An increase in its value often indicates a transition or nonlinear enhancement caused by unsteady shock motion. The regression residual can characterize the magnitude of the unexplained residual error when fitting the torque time history using the fundamental frequency linear model, and can directly reflect the applicability of the linear assumption.

[0055] Based on the aforementioned identification quality indicators, this embodiment further introduces an adaptive selection mechanism to employ the most suitable dynamic derivative extraction method under different flow response characteristics. Specifically, when the fundamental frequency component accounts for a high proportion and the regression residual is small, it indicates that the aerodynamic torque response is close to linear. In this case, the fundamental frequency fitting method can be directly used to efficiently extract static and dynamic derivatives. When the proportion of harmonic energy increases significantly or the regression residual exceeds a preset threshold, it indicates the presence of non-negligible nonlinear components in the response. In this case, a multi-harmonic fitting method is needed to more accurately capture the contribution of higher-order harmonics, or an energy criterion independent of harmonic structure can be directly used for equivalent damping determination. When the identification quality is in the critical or uncertain range, both multi-harmonic fitting and energy criterion methods can be used simultaneously for cross-validation. By analyzing the consistency of the results obtained from the two methods, the reliability of the dynamic stability determination can be further improved.

[0056] In optional embodiments, the construction of the transition-centroid coupling index may include: Define the periodic average transition position and the center of gravity position. normalized distance

[0057] in, For reference length; according to With dynamic derivative The overall risk level output includes: when Less than or equal to the preset risk assessment threshold and At that time, it was determined to be a high-risk, dynamically unstable window; when Greater than the preset threshold for risk assessment and At that time, it was determined to be a low-risk, dynamically stable window; Other situations are classified as medium risk.

[0058] In this embodiment, the step of constructing a transition-center-of-gravity coupling index and outputting a risk level based on it aims to transform the flow characteristic quantity of the transition position into a stability criterion with clear engineering orientation. The construction of this index relies on the periodic average transition position extracted from unsteady calculations. This value can be obtained by determining the transition position within the oscillation period in the preceding steps. Harmonic fitting or periodic statistics are performed to obtain the average axial position of the transition region during the dynamic process. Furthermore, this average transition position is correlated with the aircraft's center of gravity position. Establish quantitative relationships and define normalized distance. ,in As a reference length (such as fuselage length or reference chord length), this dimensionless parameter characterizes the average deviation of the transition region from the center of gravity. Its physical significance lies in quantifying the spatial relationship between the distribution of turbulent loads caused by transition and the center of gravity. Based on this, combined with the combined dynamic derivative also obtained from the previous steps... (Used to characterize equivalent pitch damping capacity), a comprehensive risk level judgment logic can be constructed. Specifically, when the average position of the transition zone is sufficiently close to the center of gravity (i.e., Furthermore, the dynamic derivative display system is in an unstable state. When a window is identified as a high-risk, dynamically unstable window, this determination logic corresponds to, for example, ... Figure 5 The revealed physical mechanism, namely Figure 5 As shown, Figure 5The paper presents the wall intermittency factor distributions obtained through numerical calculations at three different altitudes (20km, 24km, and 28km) for the proposed technical solution in simulation experiments. Comparison of these distributions reveals that at 24km, the transition point is precisely near the aircraft's center of gravity. At this altitude, the dynamic derivative reaches its maximum value among the three altitudes, indicating the worst dynamic stability of the aircraft. Specifically, when the transition zone coincides with the center of gravity, the additional aerodynamic loads generated by the turbulent boundary layer create an unbalanced distribution before and after the center of gravity, making it difficult for the positive and negative contributions of the pitching moment to cancel each other out. This amplifies the sensitivity of the aerodynamic moment to angle-of-attack disturbances, leading to increased dynamic instability. Conversely, when the transition zone is far from the center of gravity (…),… And the dynamic derivative display system is stable. When a certain condition is met, it is determined to be a low-risk dynamic stability window. All other combinations falling between these two conditions are uniformly classified as medium-risk. This determination rule includes a preset threshold. The risk assessment threshold can be set according to the specific aircraft configuration, mission profile, and engineering experience. For example, the risk assessment threshold for a specific configuration can be summarized by analyzing the influence of the relative relationship between the transition position and the center of gravity at different altitudes on the dynamic derivative. Through this coupling index and classification rule, the technical solution of this application can ultimately condense the complex dynamic evolution process of transition into a risk level that can be directly used for engineering screening, enabling designers to quickly identify high-risk operating condition windows induced by transition-center of gravity coupling.

[0059] In an optional embodiment, the output risk window or risk map includes associating the risk level with the corresponding height, Reynolds number, wall temperature, oscillation amplitude, and frequency parameters to form a risk window or risk map.

[0060] In this embodiment, the dynamic stability determination result obtained in the preceding steps is explicitly correlated with the specific flight condition parameters that produce the result, thereby forming an intuitive, searchable, and easily applied risk expression. Specifically, this step binds the risk level determined based on the transition-center-of-gravity coupling index and dynamic derivative with the corresponding flight condition parameters that induce the risk level. These flight condition parameters may include at least flight altitude, Reynolds number, wall temperature conditions (such as isothermal or adiabatic wall temperature), and pitch oscillation amplitude. and oscillation frequency (or reduce frequency) ), where the reduction frequency is defined as , For reference length, The incoming flow velocity is used as the reference point. In this way, for each calculated or analyzed operating point, the assessment method not only determines whether the point is high-risk, medium-risk, or low-risk, but also maps this determination to a specific combination of parameters (e.g., a certain oscillation frequency at a specific altitude and Reynolds number). By repeatedly executing the aforementioned assessment process for a series of continuously changing operating conditions (e.g., changing altitude or Reynolds number along the reentry trajectory, or scanning the oscillation frequency at the same altitude), the distribution pattern of risk levels with these key parameters can be obtained. Organizing these distribution patterns in graphical form creates risk windows or risk maps. For example, a two-dimensional contour map can be used to mark the distribution range of high-risk areas (windows) on a parameter plane with altitude (or Reynolds number) as the horizontal axis and oscillation frequency as the vertical axis. Alternatively, a multi-dimensional table or graph can be used to show the changing trend of risk levels with multiple parameter combinations such as altitude, Reynolds number, wall temperature, oscillation amplitude, and frequency.

[0061] In an optional embodiment, the transition model is an intermittent factor transition model, which is solved synchronously with the unsteady Reynolds-averaged Navier-Stokes solver at each time step. It is used to update the transition state in real time during the unsteady process, obtain the wall intermittent factor distribution by solving the transport equation containing the intermittent factor, and then determine the transition location and its evolution over time.

[0062] The transition model used in this embodiment can specifically be an intermittent factor-type transition model, which obtains the intermittent factor characterizing the intermittent state of the flow between laminar and turbulent flow by solving additional transport equations. This allows for a quantitative description of the initiation and development of transitions. As a preferred implementation of this type of model, the Fu-Wang transition model can be adopted. This model, within the unsteady Reynolds-averaged Navier-Stokes (URANS) framework, couples the solution to the pulsating kinetic energy... Specific dissipation rate and intermittent factors The three transport equations have the following governing equations:

[0063]

[0064]

[0065] Among them, the intermittent factor The introduction of this allows for a fine characterization of the flow state, that is, when Time indicates fully laminar flow. The time interval represents fully turbulent flow, while the intermediate value corresponds to the transition region. By solving the above set of equations, the distribution of the intermittent factor at each location and time within the computational domain can be obtained, thus laying the foundation for accurately capturing the dynamic process of transition.

[0066] In this embodiment, the intermittent factor transition model and the URANS solver are not simply used in series, but rather they are solved and interact with each other closely and synchronously within each time step. Specifically, during the iteration of each physical time step, the flow field progression solution is first completed to obtain updated basic flow field variables such as velocity, pressure, and density. Subsequently, based on the updated flow field information, the intermittent factor is solved synchronously. pulsating kinetic energy and specific dissipation rate The transport equations are then used to update the transition-related variables. Next, the updated intermittent factors are used. Based on the effective viscosity coefficient Modeling formulas (such as) The turbulent viscosity is corrected, and this corrected viscosity coefficient will again affect the turbulent terms in the flow field solution. Finally, it is determined whether the flow field variables and transition variables in this time step meet the convergence criteria. If not, the above solution process is repeated until convergence, thus forming a closed-loop iteration of "flow field advancement - transition variable update - turbulent viscosity correction - boundary layer state update" within one time step. Through this coupling method that realizes real-time updates of the transition state in each time step, the evaluation method in this embodiment can accurately reflect the dynamic migration of the transition position caused by the periodic changes of parameters such as angle of attack and pressure gradient during forced pitch oscillation, providing a physical basis for subsequently extracting the time-varying characteristics of the transition position and analyzing its impact on dynamic stability.

[0067] In an optional embodiment, the establishment of a coupled solution framework between an unsteady Reynolds-averaged Navier-Stokes solver and an intermittent factor transition model, which synchronously updates flow field variables and transition state variables at each time step and iterates until convergence, may include: completing closed-loop iterations of flow field advancement, transition variable update, turbulent viscosity correction, and boundary layer state update at each time step until the convergence criteria for each step are met.

[0068] In this embodiment, the established coupled solution framework of the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor transition model aims to achieve tightly coupled solution of the flow control equations and transition transport equations within each physical time step. This coupled solution framework does not treat the transition model as a post-processing step or merely perform unidirectional data transfer; instead, it requires that flow field variables and transition state variables be synchronously updated and mutually influence each other during the iteration process of each time step. Specifically, at the beginning of each time step, the compressible Navier-Stokes equations are solved based on the current flow field variables (such as density, velocity, pressure, etc.) to obtain an initial update of the flow field. Subsequently, based on the updated flow field information, the transport equations included in the intermittent factor transition model are solved synchronously to obtain the current step values ​​of the transition-related variables. Based on this, the updated intermittent factor is used... The turbulent viscosity coefficient is corrected based on the effective viscosity coefficient. The modeling relationship is used to weight and fuse the contributions of laminar and turbulent flow. This corrected effective viscosity coefficient is fed back into the turbulence term in the Navier-Stokes equations, thus affecting the flow field solution in the next iteration. This process is repeated within the same time step until the flow field variables and transition state variables all meet the preset convergence criteria, thereby forming a complete closed-loop iteration.

[0069] This closed-loop iterative process can be further refined into four sequentially executed and coupled sub-steps. Within a single time step, the flow field advancement, transition variable update, turbulent viscosity correction, and boundary layer state update are sequentially performed until convergence. First, in the flow field advancement step, based on the effective viscosity coefficient of the previous time step or iteration, the unsteady Reynolds-averaged Navier-Stokes equations are solved to obtain an approximate solution for the flow field variables in the current iteration. Second, in the transition variable update step, the intermittent factor is solved using the latest obtained flow field variables. pulsating kinetic energy and specific dissipation rate The transport equations allow the transition state to respond in real time to changes in the flow field. For example, when an increase in the angle of attack leads to an enhancement of the adverse pressure gradient, the intermittent factor distribution can be adjusted accordingly to reflect the shift in the transition position. Then, in the turbulent viscosity correction step, based on the updated intermittent factor... Through effective viscosity coefficient The turbulent viscosity is recalculated using a modeled formula. This corrected viscosity coefficient directly reflects the intermittent characteristics of the flow at the current moment, i.e., the spatial distribution ratio of laminar and turbulent flow. Finally, in the boundary layer state update step, the corrected turbulent viscosity and other turbulent flows are substituted into the boundary layer-related calculations or directly used as the source term for the next flow field propulsion, thereby completing the update of momentum and energy transport characteristics within the boundary layer. If the iteration results of the above four sub-steps have not yet reached the convergence criterion, the closed-loop process is repeated based on the updated value of the current step until the flow field variables and transition state variables remain stable within a given tolerance before proceeding to the solution of the next physical time step. Through this tightly coupled solution method that realizes the real-time evolution of the transition state in each time step, the evaluation method in this application can accurately capture the dynamic migration of the transition position caused by the periodic change of the angle of attack during forced pitch oscillation and its coupling effect with shock wave-boundary layer disturbance.

[0070] To further illustrate the effectiveness of this technical solution, a detailed explanation is provided below with reference to a specific embodiment.

[0071] This embodiment is based on flight test data of a blunt cone with a 5° semi-cone angle obtained from NASA. It studies a blunt cone with a 5° semi-cone angle configuration, and its geometric model parameters are as follows: Figure 2 As shown, that is, the total length Head radius The position of the centroid is Place. Figure 3 The computational grid used for numerical simulation is shown. Flight test data indicates that, as Figure 4 As shown, the aircraft vibrates most violently at an altitude of approximately 24 km, indicating poor dynamic stability under this condition. This technical solution is applied to evaluate the blunt cone model. First, a numerical simulation condition is constructed according to step 102, with input including... Figure 2 The center of gravity position shown ( ) and reference length ( Determine the geometric parameters of the flow and set the incoming flow conditions and wall boundary conditions. Select... Figure 4 The vibration amplitude and frequency shown are (0.67° and 6°). Dynamic stability calculations were performed at altitudes of 20km, 24km, and 28km. Subsequently, steps 104 to 110 were followed for steady-state initialization, unsteady coupled solution, transition feature extraction, and dynamic derivative identification. Through the unsteady-transition coupled solution of this technical solution, dynamic stability evaluation results at different altitudes can be obtained. Figure 5 The distribution of wall intermittent factors, obtained through numerical calculations, is presented at three different altitudes: 20km, 24km, and 28km. From... Figure 5 It can be clearly seen that at altitudes of 20km and 28km, there are transition zones (intermittent factors). The transition zone from 0 to 1 is located at the distance from the center of gravity of the aircraft. The turning point is relatively far away. At an altitude of 24km, the turning point is located near the center of gravity of the aircraft.

[0072] Further analysis of the calculated dynamic derivatives revealed that at an altitude of 24 km, the combined dynamic derivatives... The highest value among the three altitudes indicates that the aircraft is most unstable at this altitude, which is consistent with... Figure 4 The flight test observations shown are highly consistent. When the transition occurs near the center of gravity, the additional aerodynamic loads generated by the turbulent boundary layer form a highly unbalanced distribution before and after the center of gravity, making it difficult for the positive and negative contributions of the pitching moment to cancel each other out. This amplifies the change in aerodynamic moment on the angle-of-attack disturbance, leading to a significant increase in dynamic instability. When the transition position moves further upstream of the center of gravity, the lever arm of the turbulent loads shortens, and the additional aerodynamic moment partially cancels out, thus mitigating instability to some extent. Similarly, when the transition is downstream of the center of gravity, the turbulent loads are concentrated behind the center of gravity, which also disrupts the moment balance. However, its resultant moment and lever arm are relatively small, so the destabilizing effect is weaker compared to the case where the transition is at the center of gravity.

[0073] As can be seen from this embodiment, the technical solution of this application extracts the average periodic transition position through step 108. Combined with the transition-centroid coupling index constructed in step 112 (such as...) This method can effectively identify the 24km high-risk dynamic instability window. This verifies the effectiveness of this technical solution in integrating and coupling transition dynamic evolution with stability analysis, providing more accurate quantitative basis for engineering design than traditional methods.

[0074] It should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and are not intended to limit them. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A method for evaluating the dynamic stability of a hypersonic vehicle, characterized in that, include: For the hypersonic vehicle to be evaluated, a numerical simulation condition is constructed. The numerical simulation condition includes input vehicle geometry and reference parameters, external flow conditions, wall boundary conditions, and pitch oscillation excitation parameters. The vehicle geometry and reference parameters include the vehicle's center of gravity position and reference length. Based on the numerical simulation conditions, steady-state or quasi-steady-state flow field calculations are performed without applying pitch oscillation excitation defined by the pitch oscillation excitation parameters to obtain the initial flow field, and the relevant variables of the transition model are initialized to match the initial flow field. A coupled solution framework is established between the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor transition model. The flow field variables and transition state variables are updated synchronously in each time step and iterated until convergence, thereby obtaining the aerodynamic torque time history of the aircraft under the pitch oscillation excitation. Within the oscillation period, the transition position is extracted in real time along a specified path on the aircraft wall to obtain a sequence of transition positions changing over time. Periodic analysis is then performed on this sequence to obtain migration characteristic parameters characterizing the periodic evolution of the transition position. The migration characteristic parameters include the periodic average transition position, migration amplitude, and migration phase. The obtained aerodynamic moment time history is preprocessed to remove the initial transient period, and the dynamic derivative or equivalent damping index of the aircraft is extracted by an adaptive identification method that combines multi-harmonic fitting, energy criteria, or a combination thereof. Based on the periodic average transition position and the aircraft's center of gravity position, a transition-center of gravity coupling index is constructed. Combined with the dynamic derivative or equivalent damping index, the dynamic stability risk level of the aircraft under the current operating conditions is comprehensively determined, and a risk window or risk map is output.

2. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 1, characterized in that, The pitch oscillation excitation parameters include the initial angle of attack. Oscillation amplitude oscillation frequency and initial phase The variation law of the aircraft's pitch angle or equivalent angle of attack over time is as follows: ;in, Define the angular frequency; define the reduced frequency. ,in, For reference length, The incoming flow velocity; The aircraft geometry and reference parameters include reference length. Reference area, center of gravity position And the definition of the coordinate system.

3. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 1, characterized in that, The migration feature parameters used to characterize the periodic evolution of the transition position include: Define the transition point along the specified wall path. ,in, For the distribution of intermittent factors along the wall path, The threshold value is a preset value used to determine the transition. For the transition position Harmonic fitting is performed, expressed as: in, The oscillation angular frequency is used to obtain the periodic average transition position. Migration amplitude and migration phase ; in, , Calculated using orthogonal projection: The oscillation period is [the period of time].

4. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 3, characterized in that, The transition determination threshold The transition location is determined by one or more combinations of the intermittent factor threshold, the heat flux change threshold along the friction, and the friction change threshold along the friction; and the transition location is extracted along multiple azimuth paths or multiple meridians to form the spatial distribution characteristics of the transition migration.

5. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 1, characterized in that, The adaptive identification method using multi-harmonic fitting, energy criteria, or a combination thereof to extract the dynamic derivative or equivalent damping index of the aircraft includes: The pitching moment coefficient is obtained by using multi-harmonic fitting. Represented as: in, For the fitting order; when When, the static stability derivative is obtained. and the combined dynamic derivative ,in, For pitch angular velocity damping derivative, The derivative of the rate of change of angle of attack, This represents the oscillation amplitude. And / or, using the energy criterion, define a period. Equivalent damping work within ,in, For pitch angular velocity, when It is determined to be dynamically stable when It is determined to be dynamically unstable.

6. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 5, characterized in that, The preprocessing of the obtained aerodynamic moment time history further includes: calculating the proportion of fundamental frequency components, the proportion of harmonic energy, or the regression residual to form an identification quality index; and adaptively selecting, based on the identification quality index, to use multi-harmonic fitting or energy criteria, or to use both for cross-validation.

7. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 1, characterized in that, The construction transition-centroid coupling index includes: Define the periodic average transition position and the center of gravity position. normalized distance in, For reference length; according to With dynamic derivative The overall risk level output includes: when Less than or equal to the preset risk assessment threshold and At that time, it was determined to be a high-risk, dynamically unstable window; when Greater than the preset threshold for risk assessment and At that time, it was determined to be a low-risk, dynamically stable window; Other situations are classified as medium risk.

8. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 7, characterized in that, The output risk window or risk map includes associating the risk level with the corresponding height, Reynolds number, wall temperature, oscillation amplitude, and frequency parameters to form a risk window or risk map.

9. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 1, characterized in that, The transition model is an intermittent factor type transition model, which is solved synchronously with the unsteady Reynolds-averaged Navier-Stokes solver at each time step. It is used to update the transition state in real time during the unsteady process. By solving the transport equation containing the intermittent factor, the distribution of the wall intermittent factor is obtained, thereby determining the transition location and its evolution over time.

10. The method for evaluating the dynamic stability of a hypersonic vehicle according to claim 1, characterized in that, The establishment of a coupled solution framework between the unsteady Reynolds-averaged Navier-Stokes solver and the intermittent factor transition model, which synchronously updates the flow field variables and transition state variables at each time step and iterates until convergence, includes: completing closed-loop iterations of flow field advancement, transition variable update, turbulent viscosity correction, and boundary layer state update at each time step until the convergence criteria for each step are met.

Images