An aircraft life control method and system based on digital twinning

Abstract

The present application relates to the technical field of digital twinning, in particular to a kind of airplane life control method and system based on digital twinning.The real-time boundary condition is constructed by mapping the motion data collected by multiple source sensors to virtual environment in the present application, which can accurately calculate the geometric nonlinear stiffness change caused by large deformation, combine the instantaneous stiffness characteristics with external excitation frequency to realize the accurate identification of resonance risk and the superposition of dynamic stress amplification effect, convert complex stress state into quantifiable damage potential index using strain energy density criterion containing stress triaxiality correction, define dynamic safety envelope according to material fatigue tolerance threshold to realize the real-time visualization of structure safety boundary, perform spatial intersection operation on asymmetric safety area to generate accurate executable flight constraint range, predict flight state based on driving intention and actively cut off out-of-bound maneuvering signal, which avoids dangerous maneuvering that leads to high strain energy density accumulation in flight process.

Description

Technical Field

[0001] This invention relates to the field of digital twin technology, and in particular to an aircraft life control method and system based on digital twins. Background Technology

[0002] The field of digital twin technology involves constructing a digital model consistent with the real object in virtual space by collecting and mapping the status of physical entities through multi-source data and dynamic modeling, and realizing operation monitoring, status assessment, behavior prediction and decision support throughout its entire life cycle. This includes modeling the entity's geometric structure, physical mechanism, collecting operating condition parameters, associating historical operating data, and establishing a virtual-real synchronous update mechanism.

[0003] Among them, aircraft life control methods refer to a class of technical solutions adopted for the management and control of the service life of aircraft structural components and key systems during long-term service. Under the circumstances of constantly changing flight load environmental conditions and usage frequency, the usable life of the aircraft is determined and controlled. Typically, life calculation benchmarks based on structural parts are established during the aircraft structural design stage. During use, parameters such as the number of takeoffs and landings, flight duration, speed, altitude, and load changes for each flight are recorded. The structural life consumption is accumulated item by item according to a pre-determined fatigue accumulation rule. At the same time, the remaining life is manually corrected by combining periodic inspection records and maintenance and replacement records, thereby completing the control and management of aircraft life.

[0004] Existing aircraft life control methods mainly rely on static load benchmarks and general fatigue accumulation rules determined during the design phase, ignoring the individualized structural evolution characteristics of a single aircraft in complex and variable service environments. Simply accumulating takeoffs and landings or flight durations is insufficient to capture the instantaneous stress redistribution caused by geometric nonlinear deformation due to high-maneuver flight. Pre-set theoretical fatigue models cannot reflect the complex mechanical response generated by the coupling of aerodynamic loads and structural stiffness in real time, leading to the continuous accumulation of potential structural damage without monitoring. The reliance on periodic manual inspections and maintenance records has significant time lags and cannot promptly identify resonance hazards caused by component stiffness degradation or specific frequency coupling. This passive life management mode cannot constrain pilot control behavior in real time during flight and cannot prevent the input of ultimate loads that lead to premature structural failure, ultimately resulting in deviations in the assessment of the remaining life of the airframe and insufficient flight safety assurance capabilities. Summary of the Invention

[0005] The purpose of this invention is to overcome the shortcomings of existing technologies and propose an aircraft life control method and system based on digital twins.

[0006] To achieve the above objectives, the present invention adopts the following technical solution: an aircraft lifespan control method based on digital twins, comprising the following steps: S1: Collect motion data of the left and right wings and key load-bearing components of the aircraft, and map them into the digital twin environment to construct a real-time boundary condition parameter set for the digital twin; S2: Based on the digital twin real-time boundary condition parameter set, construct the finite element model of the left wing and the finite element model of the right wing of the aircraft, perform nonlinear finite element tangential stiffness calculation, and generate instantaneous stiffness characteristic data of the left and right wings. S3: Using the instantaneous stiffness characteristic data of the left and right wings, perform eigenvalue solving operations on the finite element models of the left and right wings respectively, and determine the value of the resonant additional stress by combining the digital twin real-time boundary condition parameter set. S4: Set virtual flight speed and virtual overload values ​​covering the entire flight envelope, and superimpose them onto the finite element models of the left and right wings. Extract the first and second safety envelope regions to obtain the asymmetric structure safety envelope dataset. S5: Determine the overlapping dynamic safety flight envelope based on the first safety envelope region and the second safety envelope region in the asymmetric structure safety envelope dataset. If the expected flight state coordinates are outside the overlapping dynamic safety flight envelope, generate a flight safety boundary control command.

[0007] As a further aspect of the present invention, the digital twin real-time boundary condition parameter set includes real-time strain values ​​of the left wing structure, real-time strain values ​​of the right wing structure, acceleration values ​​of the key load-bearing component structure motion state, engine speed values, flight speed values, and aerodynamic control surface deflection angle values. The instantaneous stiffness characteristic data of the left and right wings includes the sum of the geometric stiffness increment value of the left wing and the initial stiffness value of the left wing, and the sum of the geometric stiffness increment value of the right wing and the initial stiffness value of the right wing. The resonance additional stress value is specifically the additional vibration stress amplitude generated under resonance state calculated based on the retrieved dynamic amplification factor. The asymmetric structure safety envelope dataset includes the first safety envelope region of the left wing structure and the second safety envelope region of the right wing structure. The flight safety boundary control command includes truncation of electrical signals exceeding the overlapping dynamic safety flight envelope and the aerodynamic control surface deflection angle value after forced adjustment.

[0008] As a further aspect of the present invention, the step of obtaining the digital twin real-time boundary condition parameter set specifically includes: S111: Collect the real-time strain values ​​of the left wing structure output by the sensor array arranged on the surface of the left wing of the aircraft, collect the real-time strain values ​​of the right wing structure output by the sensor array arranged on the surface of the right wing of the aircraft, and simultaneously collect the acceleration values ​​of the key load-bearing component structure motion state output by the sensor array on the key load-bearing component, and generate multi-channel structural motion state sensing data. S112: Detect the current engine speed value through the flight control bus, obtain the real-time flight speed value and aerodynamic control surface deflection angle value, and integrate the multi-channel structural motion state sensing data to obtain a multi-source heterogeneous state dataset of the whole aircraft. S113: Map the real-time strain values ​​of the left wing structure, the real-time strain values ​​of the right wing structure, the acceleration values ​​of the key load-bearing component structure, the engine speed values, the flight speed values, and the aerodynamic control surface deflection angle values ​​in the full-aircraft multi-source heterogeneous state dataset as input variables, establish the boundary constraints required to drive the virtual model to run, and obtain the digital twin real-time boundary condition parameter set.

[0009] As a further aspect of the present invention, the step of obtaining the instantaneous stiffness characteristic data of the left and right wings specifically includes: S211: Based on the digital twin real-time boundary condition parameter set, construct the left wing finite element model and the right wing finite element model in the virtual space. Calculate the inertial load of the key load-bearing component based on the numerical acceleration of the structural motion state of the key load-bearing component. Calculate the aerodynamic distributed load based on the numerical deflection angle of the aerodynamic control surface. Apply the inertial load and aerodynamic distributed load of the key load-bearing component to the left wing finite element model and the right wing finite element model respectively, and calculate the internal stress distribution of the model under the load to generate the current model stress distribution value. S212: Based on the real-time strain values ​​of the left wing structure and the right wing structure, calculate the geometric deformation of the left wing and the right wing respectively. Combined with the stress distribution in the current model stress distribution values, calculate the geometric stiffness increment of the left wing caused by the geometric deformation of the left wing and the geometric stiffness increment of the right wing caused by the geometric deformation of the right wing, and obtain the nonlinear geometric stiffness increment dataset. S213: Calculate the sum of the geometric stiffness increment value of the left wing and the initial stiffness value of the left wing in the nonlinear geometric stiffness increment data, and calculate the sum of the geometric stiffness increment value of the right wing and the initial stiffness value of the right wing to obtain the instantaneous stiffness characteristic data of the left and right wings respectively.

[0010] As a further aspect of the present invention, the step of obtaining the value of the resonant additional stress is specifically as follows: S311: Using the instantaneous stiffness characteristic data of the left and right wings, update and eigenvalue solving operations are performed on the finite element models of the left and right wings respectively, extracting the natural frequency values ​​of the left and right wings, calculating the excitation frequency value based on the engine speed value in the digital twin real-time boundary condition parameter set, calculating the aerodynamic frequency value based on the flight speed value, and generating a set of wing modes and external excitation frequencies; S312: Calculate the frequency differences between the natural frequency values ​​of the left and right wings and the excitation frequency and aerodynamic frequency values ​​in the set of wing modes and external excitation frequencies, respectively, and filter the frequency differences that fall within the preset resonance sensitive bandwidth range threshold to obtain sensitive resonance frequency deviation index data. S313: Based on the sensitive resonance frequency deviation index data, retrieve the corresponding dynamic amplification factor, and combine the dynamic amplification factor to proportionally amplify the stress distribution in the current model stress distribution value to generate the resonance additional stress value.

[0011] As a further aspect of the present invention, the steps for obtaining the asymmetric structure safety envelope dataset are as follows: S411: Set virtual flight speed values ​​and virtual overload values ​​covering the entire flight envelope, apply the virtual flight speed values ​​and virtual overload values ​​to the updated left wing finite element model and right wing finite element model respectively, superimpose the resonant additional stress values ​​to perform multi-physics load composite, and generate full envelope virtual loading response state data. S412: Based on the full envelope virtual loading response state data, calculate the strain energy density values ​​of the key parts of the left wing and the key parts of the right wing to obtain the local strain energy density distribution sets of the left and right wings. S413: Compare the local strain energy density distribution sets of the left and right wings with the material fatigue damage tolerance threshold respectively, screen out safe points where the strain energy density value is lower than the material fatigue damage tolerance threshold, construct the first safety envelope region of the left wing structure and the second safety envelope region of the right wing structure, and generate an asymmetric structure safety envelope dataset.

[0012] As a further aspect of the present invention, the step of obtaining the flight safety boundary control command specifically includes: S511: Perform a spatial Boolean intersection operation on the first safety envelope region and the second safety envelope region in the asymmetric structure safety envelope dataset, retain the part of the two regions that completely overlap in the velocity and overload coordinate system, and generate an overlapping dynamic safety flight envelope. S512: Collect the pilot's input control commands and aerodynamic derivatives to calculate the expected flight state coordinates caused by the control commands, determine whether the expected flight state coordinates are within the geometric range defined by the overlapping dynamic safety flight envelope, identify potential dangerous states that exceed the geometric range, and obtain the flight state boundary judgment result. S513: For the flight state exceeding the boundary judgment result in the dangerous state, perform the electrical signal cut-off operation for the excess part, adjust the aerodynamic control surface deflection angle value to the allowable range, and generate flight safety boundary control command.

[0013] A digital twin-based aircraft lifecycle control system, the system comprising: The boundary acquisition module collects motion data of the aircraft's left wing, right wing, and key load-bearing components, and maps it to the digital twin environment to construct a real-time boundary condition parameter set for the digital twin. The stiffness calculation module constructs finite element models of the left and right wings of the aircraft based on the digital twin real-time boundary condition parameter set, performs nonlinear finite element tangential stiffness calculation, and generates instantaneous stiffness characteristic data of the left and right wings. The stress assessment module uses the instantaneous stiffness characteristic data of the left and right wings to perform eigenvalue solving calculations on the finite element models of the left and right wings respectively, and determines the value of the resonant additional stress by combining the real-time boundary condition parameter set of the digital twin. The envelope construction module sets virtual flight speed and virtual overload values ​​that cover the entire flight envelope and superimposes them onto the finite element models of the left and right wings. It then extracts the first and second safety envelope regions to obtain the asymmetric structure safety envelope dataset. The boundary control module determines the overlapping dynamic safety flight envelope based on the first and second safety envelope regions in the asymmetric structure safety envelope dataset. If the expected flight state coordinates are outside the overlapping dynamic safety flight envelope, a flight safety boundary control command is generated.

[0014] Compared with the prior art, the advantages and positive effects of the present invention are as follows: In this invention, by mapping motion data collected by multi-source sensors to a virtual environment to construct real-time boundary conditions, the geometric nonlinear stiffness changes caused by large deformations can be accurately calculated. By combining instantaneous stiffness characteristics with external excitation frequency, the accurate identification of resonance risk and the superposition of dynamic stress amplification effects are achieved. By using a strain energy density criterion that includes stress triaxiality correction, complex stress states are transformed into quantifiable damage potential indicators. A dynamic safety envelope is defined based on the material fatigue tolerance threshold, realizing real-time visualization of the structural safety boundary. Spatial intersection operations are performed on asymmetric safety regions to generate precise executable flight constraint ranges. Flight states are predicted based on pilot intent and overstepping control signals are actively intercepted. Dangerous maneuvers that lead to high strain energy density accumulation are avoided in real time during flight, preventing irreversible structural damage from the source. While ensuring the structural integrity throughout the entire life cycle, the maximum utilization of aircraft service efficiency is achieved. Attached Figure Description

[0015] Figure 1 This is a schematic diagram of the workflow of the present invention; Figure 2 This is a flowchart illustrating the construction process of the real-time boundary condition parameter set for the digital twin of this invention. Figure 3 This is a flowchart illustrating the generation process of instantaneous stiffness characteristic data for the left and right wings in this invention. Figure 4 This is a flowchart illustrating the determination of the numerical value of the resonant additional stress in this invention. Figure 5 This is a flowchart illustrating the generation process of the asymmetric structure safety envelope dataset for this invention. Figure 6 This is a flowchart illustrating the process of generating flight safety boundary control commands for this invention. Detailed Implementation

[0016] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0017] Please see Figure 1 This invention provides a technical solution, a method for aircraft life control based on digital twins, comprising the following steps: S1: Collect motion data of the left and right wings and key load-bearing components of the aircraft, and map them into the digital twin environment to construct a real-time boundary condition parameter set for the digital twin; S2: Based on the digital twin real-time boundary condition parameter set, construct the finite element model of the left wing and the right wing of the aircraft, perform nonlinear finite element tangential stiffness calculation, and generate instantaneous stiffness characteristic data of the left and right wings. S3: Using the instantaneous stiffness characteristic data of the left and right wings, perform eigenvalue solving calculations on the finite element models of the left and right wings respectively, and combine the digital twin real-time boundary condition parameter set to determine the value of the resonance additional stress. S4: Set virtual flight speed and virtual overload values ​​covering the entire flight envelope, and superimpose them onto the finite element models of the left and right wings. Extract the first and second safety envelope regions to obtain the asymmetric structure safety envelope dataset. S5: Determine the overlapping dynamic safety flight envelope based on the first and second safety envelope regions in the asymmetric structure safety envelope dataset. If the expected flight state coordinates are outside the overlapping dynamic safety flight envelope, generate flight safety boundary control commands.

[0018] The digital twin real-time boundary condition parameter set includes real-time strain values ​​of the left wing structure, real-time strain values ​​of the right wing structure, acceleration values ​​of the structural motion state of key load-bearing components, engine speed values, flight speed values, and aerodynamic control surface deflection angle values. The instantaneous stiffness characteristic data of the left and right wings include the sum of the geometric stiffness increment value of the left wing and the initial stiffness value of the left wing, and the sum of the geometric stiffness increment value of the right wing and the initial stiffness value of the right wing. The resonance additional stress value is specifically the additional vibration stress amplitude generated under resonance state calculated based on the retrieved dynamic amplification factor. The asymmetric structure safety envelope dataset includes the first safety envelope region of the left wing structure and the second safety envelope region of the right wing structure. The flight safety boundary control commands include electrical signals that truncate the portion exceeding the overlapping dynamic safety flight envelope and the aerodynamic control surface deflection angle values ​​after forced adjustment.

[0019] Please see Figure 2 The specific steps for obtaining the real-time boundary condition parameter set of a digital twin are as follows: S111: Collect the real-time strain values ​​of the left wing structure output by the sensor array arranged on the surface of the left wing of the aircraft, collect the real-time strain values ​​of the right wing structure output by the sensor array arranged on the surface of the right wing of the aircraft, and simultaneously collect the acceleration values ​​of the key load-bearing component structure motion state output by the sensor array on the key load-bearing component, and generate multi-channel structural motion state sensing data. The fiber Bragg grating sensor array deployed on the surface skin and internal beam structure of the left wing, along with a corresponding array of the same type deployed on the right wing, initiates the data acquisition process. The data acquisition controller drives the fiber Bragg grating demodulator to emit laser pulses at a sampling frequency of 2000 Hz and receive reflected light signals, demodulating the center position of the reflected wavelength in real time. Based on the photoelastic effect, the demodulation converts the acquired wavelength drift into micro-strain values ​​using an internally stored wavelength-strain sensitivity coefficient. This sensitivity coefficient is determined by performing uniaxial tensile calibration experiments on the same batch of sensors in a constant-temperature laboratory environment, recording wavelength changes under different stress levels, and using the slope of a linear regression equation fitted using the least squares method. This process generates real-time strain values ​​for both the left and right wing structures. Simultaneously, a microelectromechanical system (MEMS) accelerometer array installed on a key load-bearing component at the fuselage-wing connection is activated, synchronously acquiring vibration signals along the three axes of the aircraft's axisymmetric system. The original acceleration signal was denoised using a low-pass filter. The cutoff frequency was set based on the Nyquist sampling theorem and the frequency distribution characteristics of the main structural modes of the airframe. Noise spectrum energy distribution data from the aircraft's historical flight test database were statistically analyzed, and a frequency point (500 Hz) covering 98% of the effective signal energy was selected as the cutoff frequency to filter out high-frequency environmental noise interference, thereby extracting the acceleration values ​​of the key load-bearing components. Subsequently, the synchronous controller, using a precise time protocol and the GPS timing signal as a unified clock source, strictly aligned the strain data of the left and right wings with the acceleration data of the key load-bearing components on the time axis, generating multi-channel structural motion state sensing data containing spatial distribution and temporal information. Taking a cruise monitoring scenario as an example, the wavelength drift collected by the sensor at the root of the left wing was 1.2 nanometers, and the sensitivity coefficient determined through calibration experiments was 0.001 nanometers per microstrain. The data acquisition controller performs the following calculation: 1.2 / 0.001 = 1200, yielding a real-time strain value of 1200 microstrain for the left wing structure. Similarly, a drift of 1.15 nanometers is measured at the corresponding position on the right wing, and the calculation is: 1.15 / 0.001 = 1150, yielding a real-time strain value of 1150 microstrain for the right wing structure. The output voltage signal from the accelerometer of the critical load-bearing component in the vertical direction is converted to 2.5 times the gravitational acceleration, meaning the acceleration value of the critical load-bearing component structure in motion is 2.5 times the gravitational acceleration.

[0020] S112: Detects the current engine speed value through the flight control bus, obtains the real-time flight speed value and aerodynamic control surface deflection angle value, and integrates the multi-channel structural motion state sensing data to obtain the full-aircraft multi-source heterogeneous state dataset. The aircraft's onboard bus interface controller reads flight status data frames broadcast by the flight control computer in real time via the MIL-STD-1553B data bus. The bus parser first extracts the high-pressure and low-pressure rotor speed signals of the engine from specific bits in the data frame, and obtains the current engine speed value through weighted averaging. Simultaneously, it analyzes the vacuum speed and Mach number output by the atmospheric data computer to obtain the real-time flight speed value; and extracts the actual deflection angles of the ailerons, flaps, and rudder from the feedback signals of the flight control actuation system to obtain the aerodynamic control surface deflection angle values. Next, a multi-source heterogeneous data fusion operation is performed, associating the above bus data with multi-channel structural motion state sensor data. Since the bus data is typically refreshed at a low frequency (e.g., 50 Hz), while the sensor data is refreshed at 2000 Hz, a weighted time synchronization interpolation algorithm based on a Gaussian kernel function is used to upsample the bus data to ensure accurate alignment in the time dimension. This algorithm introduces time proximity and data reliability as weighting factors to calculate the precise value at the target time. The calculation formula is as follows: In the formula, The time obtained by interpolation calculation The flight speed value; The number of sampling points within the time window for interpolation is set to 4 based on real-time requirements and computational resource limitations. For the first time window Data values ​​of one original bus sampling point; For the first The timestamp corresponding to each sampling point; This is the Gaussian kernel bandwidth parameter, which is set based on the average refresh interval statistics of the bus data (e.g., set to 0.01 seconds) and is used to control the impact of time distance on the weights. For the first The signal quality confidence weight of each data point is set according to the integrity of the check bit of that data point. If the check passes, it is set to 1.0; otherwise, it is set to 0.1. The trend compensation coefficient is set to 0.05 to compensate for the smoothing lag caused by interpolation. This represents the gradient value of local data changes. Under the aforementioned cruise state, the bus parser detected the flight speed at timestamp 20 milliseconds. It was detected at Mach 0.78 in 40 milliseconds. The target interpolation time is Mach 0.782. 21 milliseconds. Set bandwidth parameters. For 10 milliseconds, confidence level All are 1, gradient The result is 0.0001. Substituting this into the above formula, the numerator is used to calculate the weighted sum, the denominator is used to calculate the weight normalization factor, and a trend compensation term is added. Calculation is then performed: The real-time flight speed at 21 milliseconds was obtained as Mach 0.7803, thus unifying the temporal and spatial references of the data and integrating it to obtain a multi-source heterogeneous state dataset of the whole aircraft containing multi-dimensional information such as aerodynamics, motion, and structural response.

[0021] S113: Map the real-time strain values ​​of the left wing structure, the real-time strain values ​​of the right wing structure, the acceleration values ​​of the key load-bearing components, the engine speed values, the flight speed values, and the aerodynamic control surface deflection angle values ​​in the multi-source heterogeneous state dataset of the whole aircraft as input variables, establish the boundary constraints required to drive the virtual model to run, and obtain the digital twin real-time boundary condition parameter set. The multi-source heterogeneous state dataset of the entire aircraft is imported into the digital twin engine interface. First, the real-time strain values ​​of the left and right wing structures, the acceleration values ​​of key load-bearing components, engine speeds, flight speeds, and aerodynamic control surface deflection angles in the dataset are normalized, converting them into dimensionless input variables between 0 and 1. The normalization reference values ​​are selected from the design limit values ​​of various parameters in the aircraft design manual. Then, according to the input specifications of the finite element simulation software, the mapping relationship between these variables and the boundary conditions of the virtual model is established. Specifically, flight speed and aerodynamic control surface deflection angles are mapped to flow field boundary conditions, acceleration values ​​are mapped to inertial force boundary conditions, and strain values ​​are mapped to structural displacement boundary conditions. The script generator generates simulation configuration files containing these constraint definitions by writing automated scripts, establishing the boundary constraints required to drive the virtual model, and finally obtaining a formatted set of real-time boundary condition parameters for the digital twin. Taking flight speed as an example, its design limit value is Mach 2.0. The calculation 0.7803 / 2.0 = 0.39015 maps the normalized values ​​to virtual wind tunnel inlet velocity boundary conditions. For the 5-degree aerodynamic control surface deflection angle, it is converted into a rotation constraint for the wing trailing edge mesh nodes, with the rotation angle set to 5 degrees. For the 2.5 gravitational acceleration value, it is converted into an inertial load vector applied vertically downwards to the model's center of mass.

[0022] Please see Figure 3 The specific steps for obtaining the instantaneous stiffness characteristic data of the left and right wings are as follows: S211: Based on the real-time boundary condition parameter set of digital twin, construct the finite element models of the left wing and the right wing in virtual space. Calculate the inertial load of the key load-bearing components based on the numerical acceleration of the structural motion state of the key load-bearing components, and calculate the aerodynamic distributed load based on the numerical deflection angle of the aerodynamic control surfaces. Apply the inertial load and aerodynamic distributed load of the key load-bearing components to the finite element models of the left and right wings respectively, and calculate the internal stress distribution of the models under the load, generating the current model stress distribution values. Based on the real-time boundary condition parameter set of the digital twin, a pre-built high-fidelity wing structure model is called by the finite element solver in the virtual simulation platform. The solver loads the finite element models of the left and right wings respectively, both of which are modeled using a hybrid of hexahedral solid elements and shell elements. The load calculation logic reads the acceleration values ​​of the structural motion state of the key load-bearing components. According to Newton's second law, the inertial load of the key load-bearing components is calculated by multiplying the mass matrix of each node of the model with the acceleration values, and then applied as a volume load to the model. At the same time, the aerodynamic load calculation logic reads the deflection angle values ​​of the aerodynamic control surfaces. By calling a fast aerodynamic calculation program based on potential flow theory, the pressure distribution data of the wing surface is obtained, the aerodynamic distributed load is calculated, and it is mapped to the surface nodes of the wing skin elements. After the load is applied, the nonlinear statics solution module is started to calculate the mechanical response of the model under the combined action of inertial and aerodynamic loads, solves the stress tensor of each element node, and generates the current model stress distribution values ​​containing full-field von Mises stress data. In the numerical example, the total mass of the model is known to be 2000 kg, and the acceleration of the key load-bearing component is 25 m / s² (approximately 2.5 g). The load calculation logic executes the calculation: 2000 × 25 = 50000, yielding an inertial load of 50000 Newtons for the key load-bearing component. Simultaneously, based on a 5-degree control surface deflection, the aerodynamic calculation program integrates to determine the total lift generated by the aerodynamic distributed loads acting on the wing surface as 45000 Newtons. After applying these loads, the solver calculates the maximum local stress at the lower edge of the main spars at the wing root to be 80 MPa.

[0023] S212: Based on the real-time strain values ​​of the left wing structure and the right wing structure, calculate the geometric deformation of the left wing and the right wing respectively. Combined with the stress distribution in the current model stress distribution values, calculate the geometric stiffness increment of the left wing caused by the geometric deformation of the left wing and the geometric stiffness increment of the right wing caused by the geometric deformation of the right wing, and obtain the nonlinear geometric stiffness increment dataset. First, the real-time strain values ​​of the left and right wing structures were read. Using a pre-constructed strain-displacement transformation matrix, the geometric deformation of the left and right wings was calculated. This transformation matrix was determined through full-scale ground static tests on the wing structure, recording strain and displacement data under hundreds of loading conditions, and then using least squares regression analysis. Next, the stiffness increment calculation logic, combined with the current model stress distribution values, introduced geometric nonlinear theory. This logic quantifies the strengthening or softening effect of in-plane stress on the lateral stiffness of the structure by calculating the initial stress matrix (geometric stiffness matrix). Specifically, the membrane force component in the stress distribution was extracted and integrated into the element geometric stiffness matrix. The additional stiffness term caused by deformation was calculated, thus yielding the geometric stiffness increment values ​​of the left and right wings caused by the geometric deformation of the left and right wings, respectively. For example, based on the inversion of 1200 microstrain, the geometric deformation (wingtip deflection) of the left wing was determined to be 150 mm. The stiffness increment calculation logic, combined with the tensile stress distribution of 80 MPa, calculates based on geometric stiffness theory that the "tensioning effect" generated by this tensile stress increases the wing stiffness. Assuming the baseline unit stiffness increment coefficient is 19.44 (this coefficient is derived from the material's elastic modulus and cross-sectional geometric properties), the calculation logic executes: 80 × 19.44 = 1555 (approximate simplified calculation), resulting in a geometric stiffness increment of 1555 Newtons per meter for the left wing.

[0024] S213: Calculate the sum of the geometric stiffness increment value of the left wing and the initial stiffness value of the left wing in the nonlinear geometric stiffness increment data, and calculate the sum of the geometric stiffness increment value of the right wing and the initial stiffness value of the right wing to obtain the instantaneous stiffness characteristic data of the left and right wings respectively. The linear elastic stiffness of the wing under unloaded conditions was retrieved from the material property database, specifically the initial stiffness values ​​of the left and right wings. These initial stiffness values ​​were baseline values ​​determined through modal analysis based on the wing's CAD design model. Next, matrix addition was performed to sum the incremental geometric stiffness values ​​of the left wing extracted from the nonlinear geometric stiffness increment dataset with the initial stiffness values ​​of the left wing, yielding the corrected instantaneous stiffness characteristics of the left wing. Similarly, the incremental geometric stiffness values ​​of the right wing were summed with the initial stiffness values ​​of the right wing to obtain the instantaneous stiffness characteristics of the right wing. This process updated the model from a linear stiffness model to a nonlinear real-time stiffness model. The initial stiffness value of the left wing was set to 50,000 N / m. The calculation 50,000 + 1555 = 51,555 yielded an instantaneous stiffness characteristic of 51,555 N / m for the left wing. This result indicates that the actual stiffness of the wing increased by 3% due to the geometric nonlinear effects caused by aerodynamic loads.

[0025] Please see Figure 4 The specific steps for obtaining the value of the resonant additional stress are as follows: S311: Using the instantaneous stiffness characteristic data of the left and right wings, update and eigenvalue solving operations are performed on the finite element models of the left and right wings respectively. The natural frequency values ​​of the left and right wings are extracted. The excitation frequency value is calculated based on the engine speed value in the real-time boundary condition parameter set of the digital twin. The aerodynamic frequency value is calculated based on the flight speed value. The wing mode and external excitation frequency set are generated. The instantaneous stiffness characteristics of the left and right wings are assigned to the stiffness matrix of the finite element model, and eigenvalue solving (using the Lanzos algorithm) is performed to extract the first 10 modal parameters, obtaining the natural frequencies of the left and right wings respectively. Simultaneously, the frequency calculation logic, based on the engine speed values ​​in the real-time boundary condition parameter set of the digital twin, calculates the excitation frequency by dividing the speed by 60 and multiplying by the number of engine blades or the rotor order. Furthermore, based on the flight speed and wing characteristic chord length, the vortex shedding frequency is calculated according to the Strouhal number principle, generating the aerodynamic frequency. The frequency calculation logic summarizes these frequency data to generate a set of wing modes and external excitation frequencies that include structural properties and external excitation sources. For example, analyzing the updated left wing model after engine solving, the natural frequency value of the left wing in its first-order bending mode is obtained as 12.5 Hz. At the same time, the frequency calculation logic processes the engine speed of 8500 rpm, and for the first-order rotor imbalance, performs the calculation: 8500 / 60 = 141.6, obtaining the high-pressure rotor excitation frequency. For the low-pressure rotor (assuming a speed of 900 rpm), the calculation is: 900 / 60=15, yielding an excitation frequency of 15 Hz for the low-pressure rotor. Based on the flight speed of Mach 0.7803 (approximately 265 m / s), a characteristic chord length of 2 meters, and a Strouhal number of 0.1, the calculation is: 0.1×265 / 2=13.25, rounded down, yielding an aerodynamic vibration frequency of 13 Hz.

[0026] S312: The frequency differences between the natural frequency values ​​of the left and right wings and the excitation frequency and the aerodynamic frequency values ​​in the computer wing mode and external excitation frequency set are selected and the frequency differences falling within the preset resonance sensitive bandwidth range threshold to obtain sensitive resonance frequency deviation index data. A frequency difference matrix is ​​constructed by traversing the set of wing modes and external excitation frequencies. Specifically, the absolute difference between the natural frequency and excitation frequency of the left wing, and the absolute difference between the natural frequency and aerodynamic frequency of the left wing, are calculated separately. The same difference calculation is then performed on the right wing. Subsequently, a preset resonance-sensitive bandwidth threshold is introduced. This threshold is determined by statistically analyzing the standard deviation of the dispersion of each modal frequency during a full-aircraft ground resonance test of the same aircraft model, and taking three times the standard deviation as the bandwidth range (e.g., set to 1.0 Hz). All calculated frequency differences are compared with this threshold. The comparison results fall into two possibilities: first, the frequency difference is less than or equal to the threshold, indicating that the excitation frequency falls within the structure's resonance band, posing a resonance risk; second, the frequency difference is greater than the threshold, indicating that the excitation frequency is far from the structure's natural frequency and is within a safe range. The indices and deviation values ​​of the dangerous frequency pairs belonging to the first possibility are recorded, obtaining the sensitive resonance frequency deviation index data. In the example, the natural frequency of the left wing is 12.5 Hz, and the aerodynamic frequency is 13 Hz. Calculate the absolute difference: |12.5-13.0|=0.5, resulting in a difference of 0.5 Hz. The preset resonance sensitivity bandwidth threshold is known to be 1.0 Hz. The comparison logic performs a judgment: 0.5 is less than 1.0, this value belongs to the first possibility, and a resonance risk is determined.

[0027] S313: Based on the sensitive resonance frequency deviation index data, retrieve the corresponding dynamic amplification factor, combine the dynamic amplification factor to proportionally amplify the stress distribution in the current model stress distribution value, and generate the resonance additional stress value. Based on the sensitive resonance frequency deviation index data, the corresponding damping ratio is determined by retrieving structural dynamics manuals or pre-determined damping characteristic curves. This damping ratio is usually obtained through logarithmic decay experiments on material specimens. The calculation logic uses the single-degree-of-freedom forced vibration theory formula, inputting the frequency ratio (the ratio of the excitation frequency to the natural frequency) and the damping ratio, to calculate the dynamic amplification factor. For conditions with resonance risk, the stress distribution in the current model stress distribution value is proportionally amplified by combining the dynamic amplification factor, i.e., the static stress is multiplied by the dynamic amplification factor to estimate the true stress level including the resonance effect, generating the resonance additional stress value. The excitation frequency is 13 Hz and the natural frequency is 12.5 Hz. The calculation logic yields a frequency ratio of 1.04. The experimentally determined structural damping ratio is set to 0.02. Substituting into the amplification factor formula, the dynamic amplification factor is calculated to be 2.0. The calculation logic applies this factor to the current model stress distribution value (80 MPa) and performs the calculation: 80 × 2.0 = 160, resulting in a resonance additional stress value of 160 MPa.

[0028] Please see Figure 5The specific steps for obtaining the asymmetric structure safety envelope dataset are as follows: S411: Set virtual flight speed and virtual overload values ​​covering the entire flight envelope, apply the virtual flight speed and virtual overload values ​​to the updated left wing finite element model and right wing finite element model respectively, superimpose the resonance additional stress value to perform multi-physics load composite, and generate full envelope virtual loading response state data. The process of setting virtual flight speed and virtual overload values ​​that cover the entire flight envelope is as follows: Obtain the design limit flight speed boundary data and the design limit maneuver overload boundary data of the aircraft; The velocity discretization sampling interval is determined based on the design limit flight speed boundary data, and the overload discretization sampling interval is determined based on the design limit maneuver overload boundary data; Set the speed sampling step size and the overload sampling step size; By using the velocity sampling step size value, equal-interval interpolation is performed on the velocity discretization sampling interval to generate a discretized virtual flight velocity numerical sequence; By using the overload sampling step size value, equal-interval interpolation processing is performed on the overload discretized sampling interval to generate a discretized virtual overload numerical sequence; First, the design limit flight speed boundary data (Mach 0.2 to 0.9) and design limit maneuver overload boundary data (gravitational acceleration from -1 to +9) of the aircraft are acquired. The environment generator determines the velocity discretization sampling interval based on the design limit flight speed boundary data and sets a velocity sampling step size (e.g., Mach 0.05). Using this step size, it performs equal-interval interpolation on the interval to generate a discretized virtual flight speed numerical sequence. Similarly, using an overload sampling step size (e.g., 0.5 gravitational acceleration), it performs equal-interval interpolation on the overload discretization sampling interval to generate a discretized virtual overload numerical sequence. Subsequently, the environment generator combines these two sequences to generate thousands of operating points covering the entire flight envelope. For each operating point, the simulation engine applies it as input to the updated left and right wing finite element models, calculates the basic response, superimposes resonant additional stress values, performs multiphysics load composite calculations, and finally generates virtual loading response state data covering the entire flight envelope. For example, a virtual flight speed of Mach 0.8 and a virtual overload of 5 gravitational acceleration can be selected as an operating point. The simulation engine calculates the base stress under this condition to be 160 MPa. Adding the additional stress (160 MPa) under the current resonance state, the calculation is: 160 + 160 = 320, resulting in a total stress response of 320 MPa under this condition.

[0029] S412: Based on the virtual loading response state data of the full envelope, calculate the strain energy density values ​​of the key parts of the left wing and the key parts of the right wing to obtain the local strain energy density distribution sets of the left and right wings. Based on the virtual loading response state data of the full envelope, the stress and strain tensors of key parts (such as the wing root and pylon connections) of the left and right wing models at each discrete working point are extracted. Tensor inner product calculations are performed, and to more accurately characterize the structural damage potential under complex stress states, a correction factor considering stress triaxiality is introduced to calculate the corrected strain energy density per unit volume. The improved strain energy density algorithm formula used in the calculation process is as follows: In the formula, The calculated corrected strain energy density is expressed in megajoules per cubic meter. For the extracted von Mises equivalent stress; This corresponds to the equivalent plastic strain; This is the triaxiality sensitivity coefficient of the material, which is determined based on the plane strain fracture toughness test of the material (for example, it is set to 1.5 for high-strength aluminum alloys). A higher value indicates that the material is more sensitive to hydrostatic pressure; The hydrostatic pressure (mean normal stress) is obtained from the trace of the stress tensor. This is the yield strength constant of the material (e.g., 500 MPa). The damage nonlinearity index describes the rate at which damage worsens with stress state. It is determined by fitting fatigue crack propagation rate test data and is typically set to 2.0. This formula is used to calculate the strain energy density values ​​of critical components on the left and right wings, respectively. In specific calculations, for a given operating point, the equivalent stress at the critical components is... The equivalent strain is 320 MPa. The hydrostatic pressure is 0.004. Set the yield strength to 100 MPa. It is 500 MPa. It is 1.5. The value is 2.0. Perform the calculation: The strain energy density at this point was found to be 1.0816 megajoules per cubic meter, a value that comprehensively reflects the combined effect of tension and volumetric deformation on lifespan. The calculation results for all operating points were summarized to obtain the local strain energy density distribution sets for the left and right wings.

[0030] S413: Compare the local strain energy density distribution sets of the left and right wings with the material fatigue damage tolerance threshold respectively, screen out safe points where the strain energy density value is lower than the material fatigue damage tolerance threshold, construct the first safety envelope region of the left wing structure and the second safety envelope region of the right wing structure, and generate an asymmetric structure safety envelope dataset. The material fatigue damage tolerance threshold for wing skin materials (such as 7050 aluminum alloy) is retrieved from a materials science database. This threshold is calculated by converting the critical stress intensity factor of the material, determined through fracture toughness experiments, into the critical strain energy density value using the crack propagation rate formula. The logic unit compares each value in the local strain energy density distribution of the left and right wings with this threshold. The comparison results have two possible classifications: the first is a "low damage risk" safe point, i.e., a condition where the strain energy density value is strictly below the material fatigue damage tolerance threshold; the second is a "high damage risk" dangerous point, i.e., a condition where the strain energy density value is higher than or equal to the threshold. The logic unit filters out the condition points belonging to the first category, connects all the safe points of the left wing in a velocity-overload coordinate system to construct the first safety envelope region of the left wing structure; similarly, the second safety envelope region of the right wing structure is constructed, generating an asymmetric structure safety envelope dataset. For example, the set material fatigue damage tolerance threshold is 1.2 megajoules per cubic meter. The logic unit compares 1.0816 with 1.2 and determines that 1.0816 is less than 1.2, classifying the operating condition as a Class I condition and marking it as a safe point. If the calculated value for another high overload condition is 1.35 megajoules per cubic meter, then it determines that 1.35 is greater than 1.2, classifying it as a Class II condition and marking it as a dangerous point.

[0031] Please see Figure 6 The specific steps for obtaining flight safety boundary control commands are as follows: S511: Perform a spatial Boolean intersection operation on the first and second safety envelope regions in the asymmetric structure safety envelope dataset, retaining the parts of the two regions that completely overlap in the velocity and overload coordinate system, and generating an overlapping dynamic safety flight envelope. Extract the first safety envelope region (left wing safety domain) and the second safety envelope region (right wing safety domain) from the asymmetric structure safety envelope dataset. These two regions appear as two polygonal regions in the velocity-overload plane. Perform a spatial Boolean intersection operation on these two geometric regions, retaining only points that belong to both the first and second regions. The result generates a new, narrower closed region, namely the overlapping dynamic safety flight envelope. A numerical example illustrates that the first safety envelope region covers overloads from -1g to 6.5g, and the second safety envelope region covers from -1g to 5.0g. Performing the intersection operation and comparing the upper limits of the two intervals: min(6.5, 5.0) = 5.0, we find that the upper limit of the overlapping dynamic safety flight envelope in the overload dimension is 5.0g.

[0032] S512: Collects the pilot's input control commands and aerodynamic derivatives to calculate the expected flight state coordinates caused by the control commands, determines whether the expected flight state coordinates are within the geometric range defined by the overlapping dynamic safety flight envelope, identifies potential dangerous states that exceed the geometric range, and obtains the flight state boundary judgment result. The system acquires pilot input commands via a joystick sensor and, combined with current flight parameters, performs predictive calculations using a stored aerodynamic derivative matrix. This matrix contains parameters such as lift slope and moment coefficient under different flight conditions, obtained through fitting and analysis of extensive wind tunnel test data. The monitor derives the expected flight state coordinates that the command will cause within a short time (e.g., 0.5 seconds). Subsequently, the monitor uses a ray casting or polygon inclusion test algorithm to determine whether these expected flight state coordinates are within the geometrically defined overlapping dynamic safety flight envelope. The judgment result is divided into two possibilities: the first is a "normal" state, where the expected coordinates are within the geometric envelope; the second is a "crossing the boundary" dangerous state, where the expected coordinates are outside the geometric envelope. If the coordinates fall outside the range, it is identified that the command will cause the aircraft to enter a potentially dangerous state, and a flight state crossing the boundary judgment result for the dangerous state is output. For example, the pilot requests a 6.0g overload. The monitor calculates the expected flight state coordinates as (Mach 0.8, 6.0g). The monitor compares this coordinate with the overlapping dynamic safety flight envelope (upper limit 5.0g). The judgment is: 6.0 is greater than 5.0, the coordinate point belongs to the second possibility, and is determined to be outside the geometric range of the envelope, falling into the "out of bounds" category.

[0033] S513: For flight status exceeding the boundary judgment result in a dangerous situation, perform electrical signal cut-off operation for the excess part, adjust the aerodynamic control surface deflection angle value to the allowable range, and generate flight safety boundary control command; When a flight status exceeding the safety envelope is received, the actuator activates the boundary protection logic. The actuator performs an electrical signal cutoff operation for the portion exceeding the safety envelope, effectively limiting the pilot's command signals within the fly-by-wire system. Based on the boundary values ​​of the overlapping dynamic safety flight envelope, the actuator calculates the corresponding maximum permissible control surface deflection angle and forcibly adjusts the aerodynamic control surface deflection angle to within the permissible range. This process generates a modified flight safety boundary control command. This control command serves as a hard constraint input to the fly-by-wire system, physically limiting the deflection amplitude of the control surfaces. This effectively prevents the aircraft from entering the dangerous region of excessive strain energy density caused by high-G maneuverability, thus avoiding irreversible plastic deformation or resonance damage accumulation, ensuring the integrity of the airframe structure and fatigue life safety during escort operations. Faced with a 6.0g command, the actuator detects a safety boundary of 5.0g. The actuator performs a cutoff operation. Assuming 6.0g corresponds to a command voltage of 6 volts and 5.0g corresponds to a command voltage of 5 volts, the actuator forcibly clamps the output voltage to 5 volts. If 6.0g corresponds to an elevator deflection of 15 degrees, and 5.0g corresponds to 12 degrees, the actuator calculates a correction: 12 + 0 = 12, forcibly adjusting the aerodynamic control surface deflection angle to 12 degrees. This correction ensures that the actual flight overload is strictly limited to the safe upper limit of 5.0g, successfully preventing excessive maneuvering that could cause structural damage.

[0034] A digital twin-based aircraft lifecycle control system, comprising: The boundary acquisition module collects motion data of the aircraft's left wing, right wing, and key load-bearing components, and maps it to the digital twin environment to construct a real-time boundary condition parameter set for the digital twin. The stiffness calculation module constructs finite element models of the left and right wings of the aircraft based on the real-time boundary condition parameter set of the digital twin, performs nonlinear finite element tangential stiffness calculation, and generates instantaneous stiffness characteristic data of the left and right wings. The stress assessment module uses the instantaneous stiffness characteristic data of the left and right wings to perform eigenvalue solving calculations on the finite element models of the left and right wings respectively, and combines the real-time boundary condition parameter set of the digital twin to determine the value of the resonant additional stress. The envelope construction module sets virtual flight speed and virtual overload values ​​that cover the entire flight envelope and superimposes them onto the finite element models of the left and right wings. It then extracts the first and second safety envelope regions to obtain the asymmetric structure safety envelope dataset. The boundary control module determines the overlapping dynamic safety flight envelope based on the first and second safety envelope regions in the asymmetric structure safety envelope dataset. If the expected flight state coordinates are outside the overlapping dynamic safety flight envelope, a flight safety boundary control command is generated.

[0035] The above are merely preferred embodiments of the present invention and are not intended to limit the present invention in any other way. Any person skilled in the art may make changes or modifications to the above-disclosed technical content to create equivalent embodiments that can be applied to other fields. However, any simple modifications, equivalent changes, and modifications made to the above embodiments based on the technical essence of the present invention without departing from the scope of the present invention shall still fall within the protection scope of the present invention.

Claims

1. A method for aircraft lifespan control based on digital twins, characterized in that, Includes the following steps: S1: Collect motion data of the left and right wings and key load-bearing components of the aircraft, and map them into the digital twin environment to construct a real-time boundary condition parameter set for the digital twin; S2: Based on the digital twin real-time boundary condition parameter set, construct the finite element model of the left wing and the finite element model of the right wing of the aircraft, perform nonlinear finite element tangential stiffness calculation, and generate instantaneous stiffness characteristic data of the left and right wings. The specific steps for obtaining the instantaneous stiffness characteristic data of the left and right wings are as follows: S211: Based on the digital twin real-time boundary condition parameter set, construct the left wing finite element model and the right wing finite element model in the virtual space. Calculate the inertial load of the key load-bearing component based on the numerical acceleration of the structural motion state of the key load-bearing component. Calculate the aerodynamic distributed load based on the numerical deflection angle of the aerodynamic control surface. Apply the inertial load and aerodynamic distributed load of the key load-bearing component to the left wing finite element model and the right wing finite element model respectively, and calculate the internal stress distribution of the model under the load to generate the current model stress distribution value. S212: Based on the real-time strain values ​​of the left wing structure and the right wing structure, calculate the geometric deformation of the left wing and the right wing respectively. Combined with the stress distribution in the current model stress distribution values, calculate the geometric stiffness increment of the left wing caused by the geometric deformation of the left wing and the geometric stiffness increment of the right wing caused by the geometric deformation of the right wing, and obtain the nonlinear geometric stiffness increment dataset. S213: Calculate the sum of the geometric stiffness increment value of the left wing and the initial stiffness value of the left wing in the nonlinear geometric stiffness increment data, and calculate the sum of the geometric stiffness increment value of the right wing and the initial stiffness value of the right wing to obtain the instantaneous stiffness characteristic data of the left and right wings respectively. S3: Using the instantaneous stiffness characteristic data of the left and right wings, perform eigenvalue solving operations on the finite element models of the left and right wings respectively, and determine the value of the resonant additional stress by combining the digital twin real-time boundary condition parameter set. The specific steps for obtaining the value of the resonant additional stress are as follows: S311: Using the instantaneous stiffness characteristic data of the left and right wings, update and eigenvalue solving operations are performed on the finite element models of the left and right wings respectively, extracting the natural frequency values ​​of the left and right wings, calculating the excitation frequency value based on the engine speed value in the digital twin real-time boundary condition parameter set, calculating the aerodynamic frequency value based on the flight speed value, and generating a set of wing modes and external excitation frequencies; S312: Calculate the frequency differences between the natural frequency values ​​of the left and right wings and the excitation frequency and aerodynamic frequency values ​​in the set of wing modes and external excitation frequencies, respectively, and filter the frequency differences that fall within the preset resonance sensitive bandwidth range threshold to obtain sensitive resonance frequency deviation index data. S313: Based on the sensitive resonance frequency deviation index data, retrieve the corresponding dynamic amplification factor, and combine the dynamic amplification factor to proportionally amplify the stress distribution in the current model stress distribution value to generate the resonance additional stress value. S4: Set virtual flight speed and virtual overload values ​​covering the entire flight envelope, and superimpose them onto the finite element models of the left and right wings. Extract the first and second safety envelope regions to obtain the asymmetric structure safety envelope dataset. The specific steps for obtaining the asymmetric structure safety envelope dataset are as follows: S411: Set virtual flight speed values ​​and virtual overload values ​​covering the entire flight envelope, apply the virtual flight speed values ​​and virtual overload values ​​to the updated left wing finite element model and right wing finite element model respectively, superimpose the resonant additional stress values ​​to perform multi-physics load composite, and generate full envelope virtual loading response state data. S412: Based on the full envelope virtual loading response state data, calculate the strain energy density values ​​of the key parts of the left wing and the key parts of the right wing to obtain the local strain energy density distribution sets of the left and right wings. S413: Compare the local strain energy density distribution sets of the left and right wings with the material fatigue damage tolerance threshold respectively, screen out safe points where the strain energy density value is lower than the material fatigue damage tolerance threshold, construct the first safety envelope region of the left wing structure and the second safety envelope region of the right wing structure, and generate an asymmetric structure safety envelope dataset. S5: Determine the overlapping dynamic safety flight envelope based on the first safety envelope region and the second safety envelope region in the asymmetric structure safety envelope dataset. If the expected flight state coordinates are outside the overlapping dynamic safety flight envelope, generate a flight safety boundary control command.

2. The aircraft lifespan control method based on digital twins according to claim 1, characterized in that, The digital twin real-time boundary condition parameter set includes real-time strain values ​​of the left wing structure, real-time strain values ​​of the right wing structure, acceleration values ​​of the key load-bearing component structure motion state, engine speed values, flight speed values, and aerodynamic control surface deflection angle values. The instantaneous stiffness characteristic data of the left and right wings includes the sum of the geometric stiffness increment value of the left wing and the initial stiffness value of the left wing, and the sum of the geometric stiffness increment value of the right wing and the initial stiffness value of the right wing. The resonance additional stress value is specifically the additional vibration stress amplitude generated under resonance state calculated based on the retrieved dynamic amplification factor. The asymmetric structure safety envelope dataset includes the first safety envelope region of the left wing structure and the second safety envelope region of the right wing structure. The flight safety boundary control command includes truncation of electrical signals exceeding the overlapping dynamic safety flight envelope and the aerodynamic control surface deflection angle value after forced adjustment.

3. The aircraft lifespan control method based on digital twins according to claim 1, characterized in that, The specific steps for obtaining the digital twin real-time boundary condition parameter set are as follows: S111: Collect the real-time strain values ​​of the left wing structure output by the sensor array arranged on the surface of the left wing of the aircraft, collect the real-time strain values ​​of the right wing structure output by the sensor array arranged on the surface of the right wing of the aircraft, and simultaneously collect the acceleration values ​​of the key load-bearing component structure motion state output by the sensor array on the key load-bearing component, and generate multi-channel structural motion state sensing data. S112: Detect the current engine speed value through the flight control bus, obtain the real-time flight speed value and aerodynamic control surface deflection angle value, and integrate the multi-channel structural motion state sensing data to obtain a multi-source heterogeneous state dataset of the whole aircraft. S113: Map the real-time strain values ​​of the left wing structure, the real-time strain values ​​of the right wing structure, the acceleration values ​​of the key load-bearing component structure, the engine speed values, the flight speed values, and the aerodynamic control surface deflection angle values ​​in the full-aircraft multi-source heterogeneous state dataset as input variables, establish the boundary constraints required to drive the virtual model to run, and obtain the digital twin real-time boundary condition parameter set.

4. The aircraft lifespan control method based on digital twins according to claim 1, characterized in that, The specific steps for obtaining the flight safety boundary control command are as follows: S511: Perform a spatial Boolean intersection operation on the first safety envelope region and the second safety envelope region in the asymmetric structure safety envelope dataset, retain the part of the two regions that completely overlap in the velocity and overload coordinate system, and generate an overlapping dynamic safety flight envelope. S512: Collect the pilot's input control commands and aerodynamic derivatives to calculate the expected flight state coordinates caused by the control commands, determine whether the expected flight state coordinates are within the geometric range defined by the overlapping dynamic safety flight envelope, identify potential dangerous states that exceed the geometric range, and obtain the flight state boundary judgment result. S513: For the flight state exceeding the boundary judgment result in the dangerous state, perform the electrical signal cut-off operation for the excess part, adjust the aerodynamic control surface deflection angle value to the allowable range, and generate flight safety boundary control command.

5. The aircraft life control method based on digital twin according to claim 3, characterized in that, For real-time flight speed values, the formula is used: ； In the formula, The value represents the flight velocity at time t obtained through interpolation, and M is the number of sampling points within the time window used for interpolation. This represents the data value of the k-th raw bus sampling point within the time window. The timestamp corresponding to the kth sampling point. For Gaussian kernel bandwidth parameters, The signal quality confidence weight for the k-th data point is set based on the integrity of the data point's checksum bit; if the checksum passes, it is set to 1.0, otherwise it is set to 0.

1. This is the trend compensation coefficient. This represents the gradient value of local data changes.

6. The aircraft life control method based on digital twin according to claim 1, characterized in that, The process of setting virtual flight speed and virtual overload values ​​that cover the entire flight envelope is as follows: Obtain the design limit flight speed boundary data and the design limit maneuver overload boundary data of the aircraft; The velocity discretization sampling interval is determined based on the design limit flight speed boundary data, and the overload discretization sampling interval is determined based on the design limit maneuver overload boundary data; Set the speed sampling step size and overload sampling step size; By using the velocity sampling step size value, equal-interval interpolation is performed on the velocity discretization sampling interval to generate a discretized virtual flight velocity numerical sequence; By using the overload sampling step size value, equal-interval interpolation is performed on the overload discretized sampling interval to generate a discretized virtual overload numerical sequence.

7. An aircraft lifespan control system based on digital twins, characterized in that, The aircraft life control method based on digital twins according to any one of claims 1-6, the system comprising: The boundary acquisition module collects motion data of the aircraft's left wing, right wing, and key load-bearing components, and maps it to the digital twin environment to construct a real-time boundary condition parameter set for the digital twin. The stiffness calculation module constructs finite element models of the left and right wings of the aircraft based on the digital twin real-time boundary condition parameter set, performs nonlinear finite element tangential stiffness calculation, and generates instantaneous stiffness characteristic data of the left and right wings. The stress assessment module uses the instantaneous stiffness characteristic data of the left and right wings to perform eigenvalue solving calculations on the finite element models of the left and right wings respectively, and determines the value of the resonant additional stress by combining the real-time boundary condition parameter set of the digital twin. The envelope construction module sets virtual flight speed and virtual overload values ​​that cover the entire flight envelope and superimposes them onto the finite element models of the left and right wings. It then extracts the first and second safety envelope regions to obtain the asymmetric structure safety envelope dataset. The boundary control module determines the overlapping dynamic safety flight envelope based on the first and second safety envelope regions in the asymmetric structure safety envelope dataset. If the expected flight state coordinates are outside the overlapping dynamic safety flight envelope, a flight safety boundary control command is generated.

Images