Method for evaluating fatigue life of aircraft structure repaired by laser considering interface effect

Abstract

This invention relates to the field of fault prediction technology, specifically to a fatigue life assessment method for laser-repaired aircraft structural components considering interface influence. The method includes the following steps: reading the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly; separating the node coordinate array of the repair layer from the node coordinate array of the matrix; generating a two-dimensional cohesive element coordinate tensor; reading the tangential and normal displacement components of the nodes within the two-dimensional cohesive element coordinate tensor; and calculating the constitutive characteristic tensor of the repair layer-matrix interface. This invention constrains the proportion of artificially dissipated energy by comparing the dissipated energy error array set, and uses the reduction of the system equilibrium traction force parameter to the numerical zero distribution boundary as the basis for life determination. This reduces the interference of artificial damping on fatigue life results, making the predicted debonding fatigue life of aircraft structural components closer to the actual interface damage evolution process.

Description

Technical Field

[0001] This invention relates to the field of fault prediction technology, and in particular to a method for assessing the fatigue life of laser-repaired aircraft structural components that takes into account the influence of interfaces. Background Technology

[0002] The fatigue life assessment method for laser-repaired aircraft structural components involves performance evaluation, interface bonding state identification, fatigue damage evolution analysis, and remaining life prediction of aircraft structural components after laser repair. It is usually aimed at evaluating the mechanical response of load-bearing structural components such as wings, skins, beams, frames, and joints under service loads, cyclic loads, thermal effects, and the coupling effect between the repair layer and the substrate interface.

[0003] Existing technologies in practical operation mainly focus on the overall performance evaluation, interface bonding state identification, fatigue damage evolution analysis, and remaining life prediction after laser repair. They typically assess the mechanical response of load-bearing structural components such as wings, skins, beams, frames, and joints under service loads, cyclic loads, thermal effects, and the coupling effect of the repair layer and substrate interface. However, the interface location is often simplified to an equivalent connection region or an overall material transition region, making it difficult to establish a continuous correspondence between the repair layer node coordinates, substrate node coordinates, interface displacement separation, cross-sectional contact stress, and local tangential stiffness changes. During cyclic load evaluation, if judgments are based solely on macroscopic stress amplitude, overall crack propagation trend, or averaged remaining life indicators, details such as local plastic slip accumulation, interface stiffness reduction, and fracture energy decay are easily smoothed out. For example, when slip concentration has occurred in a small area at the interface between the repair layer and the substrate, the overall structural response may still appear to meet safety requirements, leading to a lag in identifying the debonding initiation location. Therefore, improvements are needed. Summary of the Invention

[0004] The purpose of this invention is to overcome the shortcomings of existing technologies and to propose a fatigue life assessment method for laser-repaired aircraft structural components that takes into account the influence of interfaces.

[0005] To achieve the above objectives, the present invention adopts the following technical solution: a fatigue life assessment method for laser repair of aircraft structural components considering interface effects, comprising the following steps:

[0006] Read the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly, separate the node coordinate array of the repair layer and the node coordinate array of the matrix, and generate a two-dimensional cohesive element coordinate tensor; read the tangential displacement components and normal displacement components of the internal nodes of the two-dimensional cohesive element coordinate tensor, and calculate the constitutive feature tensor of the repair layer-matrix interface.

[0007] An external fatigue cyclic load spectrum is introduced, and the cumulative plastic slip parameters corresponding to the internal node variables of the constitutive feature tensor of the repair layer matrix interface are extracted to generate a cumulative distribution map of interface plastic slip. Based on the associated values ​​of the coordinate nodes of the corresponding regions of the cumulative distribution map of interface plastic slip, a new set of parameters is generated. The internal element array of the new set of parameters covering the constitutive feature tensor of the repair layer matrix interface is extracted and compiled to generate a fatigue cyclic degradation damage evolution matrix.

[0008] Extract the fatigue cycle degradation damage evolution matrix, perform incremental iterative calculation of the equation, record the Boolean flag variable whose diagonal element value drops below the extreme value limit during the iteration process, and calculate when the Boolean flag variable corresponds to the true state to generate an adaptive viscous regularized stiffness control matrix.

[0009] Based on the adaptive viscous regularized stiffness control matrix, a set of dissipated energy error comparison arrays is obtained; when the quotient of the set of dissipated energy error comparison arrays is lower than a set threshold, the artificial viscous damping factor value is decayed by calling the time variable formula to generate the debonding fatigue life prediction result of the aircraft structural component.

[0010] Preferably, the step of obtaining the constitutive feature tensor of the repair layer matrix interface is as follows:

[0011] Read the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly. According to the repair layer identifier, base identifier, and adjacency relationship at the interface, match the node number, spatial coordinates, and unit connection relationship item by item. Separate the repair layer node coordinate array and the base node coordinate array. Insert a dimension-two thickness zero node mesh element array at the same projection position at the interface. Write the initial interface penalty stiffness value and critical fracture energy value to each zero thickness node mesh element to obtain the dimension-two thickness zero node mesh element array.

[0012] Based on the aforementioned two-dimensional zero-thickness mesh element array, structural boundary conditions and initial static loads are applied. According to the node number index of the zero-thickness node mesh element, the characteristic vector of the mesh shape function node and the node displacement tensor at the interface are extracted. The node displacement tensor is mapped to the local coordinate direction of the zero-thickness node mesh element. The coordinate order, degree of freedom number, and displacement component correspondence of adjacent nodes are checked to generate a two-dimensional cohesive element coordinate tensor.

[0013] Read the tangential and normal displacement components of the internal nodes of the coordinate tensor of the two-dimensional cohesive unit, calculate the quotient of the cross-sectional contact stress component and the interface displacement separation component according to the local coordinate direction, update the constant term of the local tangential stiffness matrix of each internal node item by item according to the quotient, and write it into the interface constitutive feature position according to the node number order to obtain the interface constitutive feature tensor of the repair layer matrix.

[0014] Preferably, the steps for obtaining the new-order parameter set are as follows:

[0015] An external fatigue cyclic load spectrum is introduced. According to the load peak value, load valley value and cyclic loading order, the internal node variables of the constitutive feature tensor of the matrix interface of the repair layer are read one by one. The tangential slip increment of each internal node in adjacent load cycles is extracted. The slip increment of the same internal node is accumulated and the corresponding load cycle number is recorded. The load cycle number and cumulative plastic slip parameters are written into the hysteresis dissipation energy statistics list line by line. The hysteresis dissipation energy statistics list variable values ​​are mapped according to the coordinate node order of the boundary area to generate the interface plastic slip cumulative distribution map.

[0016] The associated values ​​of the coordinate nodes corresponding to the cumulative distribution map of plastic slip at the interface are read item by item. The change range of the cumulative plastic slip parameter of the same coordinate node under the variable of the number of continuous load cycles is statistically analyzed. The coordinate nodes whose change range reaches the preset fatigue damage triggering threshold are marked as decreasing control nodes. The dynamic damage decreasing control coefficient is calculated based on the associated values ​​of the decreasing control nodes. The dynamic damage decreasing control coefficient is multiplied by the initial interface penalty stiffness value and the critical fracture energy value to generate a new set of parameters.

[0017] Preferably, the step of obtaining the fatigue cyclic degradation damage evolution matrix is ​​as follows:

[0018] According to the node number, coordinate node position, and interface constitutive feature position of the internal element array of the interface constitutive feature tensor of the repair layer matrix, the corresponding updated penalty stiffness value and updated fracture energy value in the new order parameter set are extracted item by item. The updated penalty stiffness value is used to cover the constant term parameter of the local tangent stiffness matrix, and the updated fracture energy value is used to cover the interface damage energy position. The covered internal element array is then matrix-arranged according to the numerical variable of the load cycle number to generate the fatigue cycle degradation damage evolution matrix.

[0019] Preferably, the steps for obtaining the Boolean flag variable are as follows:

[0020] The fatigue cyclic degradation damage evolution matrix is ​​extracted and divided into continuous incremental step input items according to the numerical variable of the load cycle number, the node number of the interface region, and the position of the internal element array. For each continuous incremental step input item, the degradation stiffness component, damage evolution component, and node displacement component of the interface region node are read. The residual balance of the interface region node is calculated. The node displacement component is corrected according to the residual balance. When the residual balance reaches the preset balance convergence limit, the diagonal element of the local tangent stiffness matrix of the corresponding incremental step is recorded to obtain the diagonal element of the local tangent stiffness matrix.

[0021] Based on the diagonal elements of the local tangent stiffness matrix, the values ​​of the diagonal elements of the local tangent stiffness matrix are compared item by item with the set negative instability extreme value limit according to the increment step number, interface region node number, and local tangent stiffness matrix position index. If the value of the diagonal elements of the local tangent stiffness matrix drops below the negative instability extreme value limit, the corresponding interface region node, increment step number, and diagonal element value are written to the true state. If the value of the diagonal elements of the local tangent stiffness matrix does not drop below the negative instability extreme value limit, it is written to the false state, and a Boolean flag variable is obtained.

[0022] Preferably, the steps for obtaining the adaptive viscous regularized stiffness control matrix are as follows:

[0023] Filter the interface area nodes corresponding to the true state of the Boolean flag variable, read the corresponding incremental step number, specified base factor, and time-related variable, perform exponential multiplication operation on the specified base factor according to the time-related variable, retain the original diagonal element records for the interface area nodes corresponding to the false state of the Boolean flag variable, extract the exponential multiplication operation results corresponding to the true state interface area nodes, and obtain the value of the artificial viscous damping factor.

[0024] According to the node number of the true state interface region, the incremental step number, and the position of the diagonal constant term element of the local tangent stiffness matrix, the value of the artificial viscous damping factor is superimposed into the position of the diagonal constant term element of the local tangent stiffness matrix item by item. For the nodes of the false state interface region, the original value of the diagonal constant term element of the local tangent stiffness matrix is ​​retained. The superimposed local tangent stiffness matrix is ​​rearranged according to the node number of the interface region to generate an adaptive viscous regularized stiffness control matrix.

[0025] Preferably, the step of obtaining the dissipated energy error comparison array set is as follows:

[0026] The adaptive viscous regularized stiffness control matrix is ​​extracted. According to the position of the internal constant term element, the node number of the interface region, and the incremental step number, the artificial dissipation energy data and the interface damage dissipation energy data associated with the internal constant term are read one by one. The artificial dissipation energy data is used as the numerator and the interface damage dissipation energy data is used as the denominator. The division operation is performed node by node. The quotient values ​​are arranged according to the interface region node number to generate a dissipation energy error comparison array set.

[0027] Preferably, the steps for obtaining the predicted debonding fatigue life of the aircraft structural components are as follows:

[0028] Based on the dissipated energy error comparison array set, the quotient value corresponding to each interface region node is read item by item, the quotient value is compared with a set threshold, interface region nodes with quotient values ​​lower than the set threshold are filtered, the artificial viscous damping factor value corresponding to the interface region node with quotient value lower than the set threshold is read, the artificial viscous damping factor value is gradually decayed according to the time variable, and the decayed artificial viscous damping factor value is rewritten into the position of the diagonal constant term element of the local tangent stiffness matrix to generate the reconstructed control matrix.

[0029] Based on the reconstructed control matrix, the node traction force balance terms associated with the control matrix are read item by item according to the node number of the interface area, the connection relationship between adjacent coordinate nodes, and the position of the diagonal constant term element of the local tangent stiffness matrix. The difference of the node traction force balance terms of each adjacent coordinate node is checked, and the node traction force balance terms that satisfy the equilibrium convergence limit are determined as the system equilibrium traction force parameters, thus obtaining the system equilibrium traction force parameters.

[0030] Based on the connection relationship between adjacent coordinate nodes, the system balance traction force parameter values ​​are identified segment by segment until they drop to the range of the zero-value distribution boundary. The node number, increment step number, and load cycle number of the interface area covered by the zero-value distribution boundary are extracted. The total number of fatigue cycle load operations corresponding to the zero-value distribution boundary is used as the basis for determining the occurrence of debonding, and the debonding fatigue life prediction result of the aircraft structural component is generated.

[0031] Compared with the prior art, the advantages and positive effects of the present invention are as follows:

[0032] This invention extends fatigue life assessment from the overall structural response to the interface between the repair layer and the matrix by reading the geometric topology data node coordinates of the assembly and separating the repair layer node coordinate array and the matrix node coordinate array. The two-dimensional cohesive element coordinate tensor, tangential displacement components, normal displacement components, and the constitutive characteristic tensor of the repair layer-matrix interface form a continuous chain, allowing for the step-by-step expression of interface displacement separation, cross-sectional contact stress, and local tangential stiffness changes. The external fatigue cyclic load spectrum is combined with cumulative plastic slip parameters to form an interface plastic slip cumulative distribution map, transforming interface damage under cyclic loading from a single static response assessment into a regional distribution result that updates with the number of load cycles. The new parameter set covers the constitutive characteristic tensor of the repair layer-matrix interface. The system employs a partial element array and compiles a fatigue cyclic degradation damage evolution matrix, enabling the initial interface penalty stiffness and critical fracture energy values ​​to be reconstructed with the damage state, thus improving the accuracy of identifying stiffness degradation, interface debonding precursors, and fatigue failure locations. It performs extreme value boundary judgments on the diagonal elements of the local tangent stiffness matrix and introduces an artificial viscous damping factor when the Boolean flag variable is true, which can suppress iterative abrupt changes caused by interface instability, ensuring the computability of the debonding expansion stage. By constraining the proportion of artificial dissipation energy through a dissipation energy error comparison array set and using the system equilibrium traction force parameter reduced to the numerical zero distribution boundary as the basis for life determination, the interference of artificial damping on fatigue life results can be reduced, making the predicted fatigue life of aircraft structural components closer to the actual interface damage evolution process. Attached Figure Description

[0033] Figure 1 This is a schematic diagram of the steps of the present invention;

[0034] Figure 2 This is a cumulative distribution map of interfacial plastic slip;

[0035] Figure 3 This is a graph showing the decay of the viscous damping factor with the number of cycles. Detailed Implementation

[0036] 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.

[0037] Please see Figure 1-3 This invention provides a technical solution for fatigue life assessment of laser-repaired aircraft structural components, considering the influence of interfaces, comprising the following steps:

[0038] Read the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly, separate the node coordinate array of the repair layer from the node coordinate array of the matrix, and generate a two-dimensional cohesive element coordinate tensor; read the tangential displacement components and normal displacement components of the internal nodes of the two-dimensional cohesive element coordinate tensor, and calculate the constitutive characteristic tensor of the repair layer-matrix interface.

[0039] An external fatigue cyclic load spectrum is introduced, and the cumulative plastic slip parameters corresponding to the internal node variables of the constitutive feature tensor of the repair layer matrix interface are extracted to generate a cumulative distribution map of interface plastic slip. Based on the associated values ​​of the coordinate nodes of the corresponding region of the cumulative distribution map of interface plastic slip, a new set of parameters is generated. The internal element array of the constitutive feature tensor of the repair layer matrix interface covered by the new set of parameters is extracted and compiled to generate a fatigue cyclic degradation damage evolution matrix.

[0040] Extract the fatigue cycle degradation damage evolution matrix, perform incremental iterative calculation of the equation, record the Boolean flag variables whose diagonal element values ​​drop below the extreme value limit during the iteration process, and perform calculations when the Boolean flag variables correspond to the true state to generate an adaptive viscous regularized stiffness control matrix.

[0041] Based on the adaptive viscous regularized stiffness control matrix, the set of dissipation energy error comparison arrays is obtained; when the quotient of the set of dissipation energy error comparison arrays is lower than the set threshold, the artificial viscous damping factor value is decayed by calling the time variable formula, and the debonding fatigue life prediction result of the aircraft structural component is generated.

[0042] The steps for obtaining the constitutive feature tensor of the repair layer matrix interface are as follows:

[0043] Read the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly. According to the repair layer identifier, base identifier, and adjacency relationship at the interface, match the node number, spatial coordinates, and unit connection relationship item by item. Separate the repair layer node coordinate array and the base node coordinate array. Insert a dimension-two thickness zero node mesh element array at the same projection position at the interface. Write the initial interface penalty stiffness value and critical fracture energy value to each zero thickness node mesh element to obtain the dimension-two thickness zero node mesh element array.

[0044] Based on a two-dimensional zero-thickness mesh element array, structural boundary conditions and initial static loads are applied. According to the node number index of the zero-thickness node mesh element, the characteristic vector of the mesh shape function node and the node displacement tensor at the interface are extracted. The node displacement tensor is mapped to the local coordinate direction of the zero-thickness node mesh element. The coordinate order, degree of freedom number, and displacement component correspondence of adjacent nodes are checked to generate a two-dimensional cohesive element coordinate tensor.

[0045] Read the tangential and normal displacement components of the internal nodes of the two-dimensional cohesive element coordinate tensor, calculate the quotient of the cross-sectional contact stress component and the interface displacement separation component according to the local coordinate direction, update the constant parameters of the local tangential stiffness matrix of each internal node item by item according to the quotient, and write them into the interface constitutive feature position according to the node number order to obtain the interface constitutive feature tensor of the repair layer matrix.

[0046] Specifically, the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly are read. This coordinate data originates from the finite element preprocessing software and includes node numbers and three-dimensional spatial coordinates (X, Y, ...). The data list of Z) and unit connection relationships is first traversed. Based on the pre-defined component identifier of each node (e.g., "1" for repair layer, "0" for substrate), all node data are split into two independent storage structures: a repair layer node coordinate array and a substrate node coordinate array. Next, spatial position matching is performed on the nodes in the two arrays. Specifically, the projection coordinates of each repair layer node are calculated using the substrate surface normal as the projection direction. Nodes with Euclidean distances less than a preset tolerance (e.g., 0.01 mm) are searched in the substrate node coordinate array. These nodes (one repair layer node, one substrate node) are marked as interface node pairs. Two-dimensional four-node zero-thickness cohesive elements are batch-inserted between the original coordinate positions of all interface node pairs to form a complete interface mesh layer. Then, initial mechanical performance parameters are assigned to each newly generated zero-thickness node mesh element. The initial interface penalty stiffness value is set with reference to the material properties of adjacent solid elements; for example, its set value is the Young's modulus of the adjacent substrate material (e.g., 70). Dividing 50 times the GPa by the characteristic length of the adjacent solid element (e.g., 0.5 mm) yields 7.0 × 10⁻⁶. 6 The value of N / mm³ is written into the element properties, while the critical fracture energy value is based on the experimental test data of the material. For example, the critical fracture energy value of type I obtained by double cantilever beam test on the same material adhesive joint is 2.5 N / mm. This value is written into the corresponding element properties. After completing the parameter assignment of all zero-thickness elements, all element data are finally summarized to obtain a dimension two thickness zero-node mesh element array.

[0047] Based on a 2D thickness zero-node mesh element array, precise structural boundary conditions and loads are first applied to the entire aircraft structural component model. For example, all six degrees of freedom (three translational and three rotational) of all nodes at the bottom edge of the structural component are completely fixed, i.e., U1=U2=U3=UR1=UR2=UR3=0. Simultaneously, a uniformly distributed initial static load, such as a pressure of 10 MPa, is applied to the top surface of the repair layer. After loading, a linear static analysis is performed using a finite element solver to obtain the global displacement tensor of all nodes in the model. This tensor contains the displacement components of each node in the global coordinate system (X,Y,Z). Next, according to the node numbers recorded in the 2D thickness zero-node mesh element array, each zero-thickness element at the boundary is indexed one by one, and its corresponding global displacement is extracted. At the same time, a local coordinate system is established for each zero-thickness element, with its local normal axis... Perpendicular to the interface, local tangential axis Parallel to the interface, the angle of rotation of this local coordinate system relative to the global coordinate system (X, Y) The global displacement vector is then calculated using the coordinates of the two endpoint nodes of the element. Through a rotation matrix Mapping to a local coordinate system, its calculation method is as follows: ,in, It is the local normal and tangential relative displacement separation of the upper and lower surfaces of the interface. and These are the global displacement vectors and rotation matrices of the nodes on the upper and lower surfaces of the zero-thickness element, respectively. The definition of During the mapping process, the degree of freedom numbers and displacement components of adjacent node pairs for each element are matched and checked. After confirming that the data correspondence is correct, the node numbers, local coordinates, and local normal displacement components of all zero-thickness elements are then mapped. and tangential displacement components Integrate and store the data to generate a two-dimensional cohesive unit coordinate tensor.

[0048] Read the tangential displacement components of the internal nodes of the coordinate tensor of the two-dimensional cohesive element. With normal displacement components To determine the interface constitutive relationship, a calculation model based on the traction-separation rule is adopted here. First, according to the local coordinate direction, the initial interface penalty stiffness value set in the previous step is used. and (Corresponding to normal and tangential stiffness respectively), calculate the initial cross-sectional contact stress components for each node, and the calculation expression is as follows: and ,in For normal contact stress, For tangential contact stress, this process occurs, for example, under initial static load when the interface is in the linear elastic stage and no damage occurs. The calculated quotient of the contact stress component and the interface displacement separation component is... and The result is the initial stiffness value. and Next, based on this quotient, the local tangent stiffness matrix of each internal node is... The constant term parameters are updated by assigning values. This local tangent stiffness matrix is ​​a 2x2 matrix that describes the relationship between the local stress increment and the displacement separation increment, and its form is: In the initial linear elastic phase, the off-diagonal coupling term and If the value is set to 0, the diagonal items are updated directly using the calculated quotient. and After the update is complete, the local tangent stiffness matrix of each node will be... According to the node numbering order, the interface constitutive feature positions are written into a predefined data structure. This data structure associates the node number with the corresponding constitutive matrix, thereby organizing the discrete node constitutive information into an ordered whole, and obtaining the interface constitutive feature tensor of the repair layer matrix.

[0049] The steps to obtain the new-order parameter set are as follows:

[0050] An external fatigue cyclic load spectrum is introduced. According to the load peak value, load valley value and cyclic loading order, the internal node variables of the constitutive feature tensor of the repair layer matrix interface are read one by one. The tangential slip increment of each internal node in adjacent load cycles is extracted. The slip increment of the same internal node is accumulated and the corresponding load cycle number is recorded. The load cycle number and cumulative plastic slip parameters are written into the hysteresis dissipation energy statistics list line by line. The hysteresis dissipation energy statistics list variable values ​​are mapped according to the coordinate node order of the boundary area to generate the cumulative distribution map of interface plastic slip.

[0051] The associated values ​​of the coordinate nodes corresponding to the cumulative distribution map of plastic slip at the interface are read item by item. The change range of the cumulative plastic slip parameter of the same coordinate node under the variable of the number of continuous load cycles is statistically analyzed. The coordinate nodes whose change range reaches the preset fatigue damage triggering threshold are marked as decreasing control nodes. The dynamic damage decreasing control coefficient is calculated based on the associated values ​​of the decreasing control nodes. The dynamic damage decreasing control coefficient is multiplied by the initial interface penalty stiffness value and the critical fracture energy value to generate a new set of parameters.

[0052] Specifically, an external fatigue cyclic load spectrum with constant amplitude is introduced. For example, a sinusoidal load waveform with a stress ratio of 0.1 and a maximum stress of 180 MPa is set and discretized into a time series containing 100 increment steps for each load cycle. Then, a cyclic iterative process is initiated. Starting from the first load cycle, a finite element solution is performed for each increment step. Based on the node numbers stored in the constitutive characteristic tensor of the repair layer matrix interface, the displacement and stress variables of each internal node of the interface are read sequentially. At the end of the Nth load cycle, the maximum and minimum tangential displacement separation components experienced by each internal node in that cycle are extracted, and the absolute value of the difference between the two is calculated as the tangential slip increment of the current cycle. Next, this increment is added to the plastic slip parameter accumulated over all previous cycles for that node, with the update rule being: ,in It is the cumulative plastic slip parameter at the end of the Nth cycle. It is the accumulated value from the previous loop, the initial value. Set the value to 0, and simultaneously record the current load cycle count variable N. Then, combine the load cycle count variable N of each node with its corresponding cumulative plastic slip parameter. As a data pair, a dynamically updated list of hysteresis dissipation energy statistics is written line by line. Each line of this list contains the node number, the number of cycles, and the cumulative slip value. At a preset cycle interval, such as every 1000 cycles, the latest cumulative plastic slip parameter value of each node in the list of hysteresis dissipation energy statistics is mapped onto a two-dimensional grid according to the physical spatial coordinates of the nodes in the boundary area. The X and Y axes of the grid correspond to the spatial coordinates of the interface nodes, and the Z axis or color scale represents the value of the cumulative plastic slip parameter, generating a cumulative distribution map of interface plastic slip.

[0053] Read the values ​​associated with the coordinate nodes of each interface region in the cumulative distribution map of plastic slip, i.e., the cumulative plastic slip parameters. Then, for each coordinate node, its load is monitored in consecutive load cycles. The value is calculated and compared with a preset fatigue damage triggering threshold, which is set based on the fracture mechanics properties of the material, for example, the critical slip displacement when the interface material fails completely under pure tangential loading. 15%, if the experiment measured If the cumulative plastic slip parameter at a certain coordinate node is 0.4 mm, then the fatigue damage triggering threshold is set to 0.06 mm. When the value first exceeds 0.06 mm, the node is marked as a decreasing control node. Then, based on the associated values ​​of all nodes marked as decreasing control nodes, i.e., their current cumulative plastic slip parameters... Calculate a dynamic damage reduction control coefficient D. Here, an energy-based linear damage evolution model is used, and the formula for calculating this coefficient is as follows: ,in, This is the dynamic damage reduction control coefficient, with a value ranging from 0 to 1. It is the cumulative plastic slip parameter of the current node. This is the slip value at which damage begins to accumulate, i.e., the fatigue damage trigger threshold of 0.06 mm. This is the cumulative slip displacement when the material completely fails. For example, based on experimental data, it is set to 0.35 mm. For instance, if the current cumulative slip of a decreasing control node is 0.15 mm, then its dynamic damage decreasing control coefficient D is calculated as (0.15 - 0.06) / (0.35 - 0.06) ≈ 0.31. Finally, the dynamic damage decreasing control coefficient D is multiplied by the initial interface penalty stiffness value and critical fracture energy value of the node. The calculation formula is: updated penalty stiffness value = (1 - D) × initial interface penalty stiffness value, updated fracture energy value = (1 - D) × critical fracture energy value. The node numbers of all decreasing control nodes, the updated penalty stiffness values, and the updated fracture energy values ​​are collected to generate a new set of parameters.

[0054] The steps to obtain the fatigue cycle degradation damage evolution matrix are as follows:

[0055] Based on the node number, coordinate node position, and interface constitutive feature position of the internal element array of the constitutive feature tensor of the repair layer matrix interface, the corresponding updated penalty stiffness value and updated fracture energy value are extracted from the new order parameter set item by item. The updated penalty stiffness value is used to cover the constant term parameter of the local tangent stiffness matrix, and the updated fracture energy value is used to cover the interface damage energy position. The covered internal element array is then matrix-arranged according to the numerical variable of the load cycle number to generate the fatigue cycle degradation damage evolution matrix.

[0056] Specifically, according to the organization structure of the internal element array of the constitutive feature tensor of the repair layer matrix interface, that is, based on the node number, coordinate node position, and index of the interface constitutive feature position, the newly generated parameter set in the previous step is first traversed. For each entry in the set, the node number, updated penalty stiffness value, and updated fracture energy value are extracted item by item. Then, the extracted node number is used to locate the corresponding internal element array in the constitutive feature tensor of the repair layer matrix interface. The extracted updated penalty stiffness value, for example, for tangential stiffness, is calculated using the updated value (1 - D) × Overwrite the corresponding constant term parameter in the local tangent stiffness matrix of that node. Similarly, normal stiffness The corresponding overwrite is also performed. Then, the extracted updated fracture energy value overwrites the value of the interface damage energy location recorded in the node's attributes. For example, the original critical fracture energy is overwritten. Updated to (1 - D) × For nodes that do not appear in the new order parameter set, all their constitutive parameters remain unchanged. After completing the parameter overwrite operation for all damaged nodes, the modified constitutive feature tensor of the repair layer matrix interface, which reflects the material degradation state under the current load cycle, is stored as a whole data block and associated with the current load cycle number variable. This process is repeated at each preset damage assessment cycle interval (e.g., every 1000 cycles). The degraded constitutive feature tensor generated at each assessment point is stacked and arranged in the order of the load cycle number to form a time-series data structure, generating a fatigue cycle degradation damage evolution matrix.

[0057] The steps to obtain the Boolean flag variable are as follows:

[0058] The fatigue cyclic degradation damage evolution matrix is ​​extracted. According to the numerical variable of the load cycle number, the node number of the interface region, and the position of the internal element array, the fatigue cyclic degradation damage evolution matrix is ​​split into continuous incremental step input items. For each continuous incremental step input item, the degradation stiffness component, damage evolution component, and node displacement component of the interface region node are read. The residual balance of the interface region node is calculated. The node displacement component is corrected according to the residual balance. When the residual balance reaches the preset balance convergence limit, the diagonal element of the local tangent stiffness matrix of the corresponding incremental step is recorded to obtain the diagonal element of the local tangent stiffness matrix.

[0059] Based on the diagonal elements of the local tangent stiffness matrix, the values ​​of the diagonal elements of the local tangent stiffness matrix are compared with the set negative instability extreme value limit according to the increment step number, interface region node number, and local tangent stiffness matrix position index. If the value of the diagonal elements of the local tangent stiffness matrix drops below the negative instability extreme value limit, the corresponding interface region node, increment step number, and diagonal element value are written to the true state. If the value of the diagonal elements of the local tangent stiffness matrix does not drop below the negative instability extreme value limit, it is written to the false state, and a Boolean flag variable is obtained.

[0060] Specifically, the fatigue cycle degradation damage evolution matrix is ​​extracted. This matrix is ​​essentially a three-dimensional array, with dimensions representing the number of load cycles, interface node numbers, and constitutive parameters, respectively. First, the matrix is ​​decomposed along the time dimension in ascending order of the number of load cycles. For each load cycle, it is further decomposed into a series of discrete incremental step inputs. Each input corresponds to the structural state at the beginning of that incremental step, including the degradation stiffness components and damage evolution components of all interface nodes. For each incremental step, a nonlinear iterative solution process using the Newton-Raphson method is initiated. In the i-th iteration, the current displacement components of each interface region node are read. By combining the degraded stiffness component and the damage evolution component, the resistance inside the node is calculated through integration. Then, calculate the residual balance of the nodes in the interface region, i.e., the unbalanced force vector. ,in It is the external load vector applied in the current increment step, and then, based on the residual balance... and the local tangent stiffness matrix of the current node Solving for displacement correction And update the node displacement to Repeat this iterative process until the L2 norm of the residual balance is reached. Less than a preset equilibrium convergence limit, which is set as the L2 norm of the external load vector. 0.01%, for example, if If the value is 1000N, then the convergence limit is 0.1N. When this condition is met, the incremental step is considered convergent, and the two diagonal elements of the local tangent stiffness matrix of each interface region node in the convergent state are recorded. and This yields the diagonal elements of the local tangent stiffness matrix.

[0061] Based on the diagonal elements of the local tangent stiffness matrix, this data is a set containing the diagonal stiffness values ​​of all incremental steps and all interface nodes. A three-level index traversal is performed according to the incremental step number, interface region node number, and the position of the local tangent stiffness matrix (i.e., the normal or tangential diagonal element). Each diagonal element value is compared with a preset negative instability extreme value limit. This limit is set to identify numerical instability early, before computational divergence occurs due to stiffness becoming negative during the material softening stage. Its value is determined based on the initial undamaged interface penalty stiffness value; for example, if the initial normal interface penalty stiffness value is 7.0 × 10⁻⁶. 6If the negative instability extreme limit is N / mm³, then it can be set to -0.1% of its absolute value, i.e., -7000 N / mm³. This is an empirical value determined through numerical experiments that can ensure computational stability without excessively suppressing the actual softening behavior of the material. During the traversal, if the normal stiffness value of a diagonal element of a local tangent stiffness matrix, such as the node with node number 1051, is in the 25th increment step of the 5000th iteration, then... When the value drops to -7500 N / mm³, since -7500 is less than -7000, the node 1051 of the interface region, the increment step number 25, and its specific diagonal element value of -7500 N / mm³ are packaged together, and a Boolean flag is set to true for recording. Conversely, if the diagonal element value of another node is -6500 N / mm³, and its value has not dropped below the negative value instability extreme limit, then only its corresponding Boolean flag is set to false. After completing the judgment of the diagonal elements of all nodes under all increment steps, the Boolean state records of all nodes are summarized to obtain the Boolean flag variable.

[0062] The steps for obtaining the adaptive viscous regularized stiffness control matrix are as follows:

[0063] Filter the interface area nodes corresponding to the true state of the Boolean flag variable, read the corresponding incremental step number, specified base factor, and time-related variable, perform exponential multiplication on the specified base factor according to the time-related variable, retain the original diagonal element records for the interface area nodes corresponding to the false state of the Boolean flag variable, extract the exponential multiplication result corresponding to the true state interface area node, and obtain the value of the artificial viscous damping factor.

[0064] According to the node number of the true state interface region, the incremental step number, and the position of the diagonal constant term element of the local tangent stiffness matrix, the artificial viscous damping factor value is superimposed into the position of the diagonal constant term element of the local tangent stiffness matrix item by item. For the nodes of the false state interface region, the original value of the diagonal constant term element of the local tangent stiffness matrix is ​​retained. The superimposed local tangent stiffness matrix is ​​rearranged according to the node number of the interface region to generate the adaptive viscous regularized stiffness control matrix.

[0065] Specifically, the process involves filtering all records corresponding to true states in the Boolean flag variable dataset, extracting the associated interface area node numbers and increment step numbers, and for each filtered true state node, reading a pre-defined cardinality factor in the corresponding increment step. A time-related variable The specified base factor is a physical viscosity coefficient, the value of which is set according to the material type and the expected numerical stability requirements. For example, for aluminum alloy repair interfaces, it can be set empirically to [value missing]. N·s / mm³, this value needs to be small enough to avoid having a significant impact on the physical results, while also effectively suppressing numerical oscillations. The time-related variable is directly taken as the time step of the current increment step. ,For example The exponential multiplication operation 's' calculates a viscosity regularization coefficient, i.e., the numerical value of the artificial viscosity damping factor. Its calculation method follows the basic definition of viscous damping, and is... Based on the example values, it can be calculated that N / mm³, this is the calculated value The value represents the viscous damping factor for the true state node. For the interface region nodes corresponding to the Boolean flag variable in the false state, their physical behavior is stable and no artificial damping is required; therefore, their artificial viscous damping factor is set to 0. Finally, the calculation results for all true state interface region nodes are... The 0 values ​​corresponding to the pseudo-state nodes are summarized to obtain the artificial viscous damping factor value.

[0066] Following the numbering order of all interface area nodes, each node is processed one by one. First, its status is queried in the Boolean flag variable according to its number. If it is true, the value of its artificial viscous damping factor in the corresponding increment step is read. (For example (N / mm³), and simultaneously extract the current local tangent stiffness matrix of the node. Then, the artificial viscous damping factor value was... The elements are sequentially added to the diagonal constant term elements of the matrix. Specifically, this involves adding the original diagonal elements... and Updated to and For example, if the normal stiffness of a node Due to the material softening to -7500 N / mm³, its new value becomes after the addition of an artificial viscous damping factor. The positive increment of N / mm³ pulls the negative stiffness value back above the instability extreme limit, thus stabilizing the calculation process. The off-diagonal elements of the matrix remain unchanged. For interface region nodes with a false query result, no modification is made to their local tangent stiffness matrix, and the original values ​​of their diagonal constant terms are retained. After completing the conditional modification of the local tangent stiffness matrix of all interface region nodes, these modified or unchanged local tangent stiffness matrices are reorganized and arranged according to the original numbering order of the interface region nodes to form an updated stiffness matrix set containing the adaptive viscosity regularization term, generating the adaptive viscosity regularization stiffness control matrix.

[0067] The steps for obtaining the dissipated energy error comparison array set are as follows:

[0068] Extract the adaptive viscous regularized stiffness control matrix. According to the position of the internal constant term element, the node number of the interface region, and the incremental step number, read the artificial dissipation energy data and the interface failure dissipation energy data associated with the internal constant term item by item. Take the artificial dissipation energy data as the numerator and the interface failure dissipation energy data as the denominator. Perform division operation on each node and arrange the quotient values ​​according to the interface region node number to generate a dissipation energy error comparison array set.

[0069] Specifically, the adaptive viscous regularized stiffness control matrix is ​​extracted, and each interface node in each increment step is traversed according to the interface region node number and increment step sequence number. At each node, two key energy values ​​are calculated, the first being the artificial dissipation energy data. It is the energy consumed in a single increment step due to the introduction of artificial viscous damping, and its calculation formula is: ,in This is the value of the artificial viscous damping factor at that node. It is the interface displacement separation increment within this increment step. It is the duration of the increment step. and These are the start and end times of the incremental step, followed by the energy dissipation data due to interface destruction. It represents the energy consumed physically due to the evolution of material damage, and is calculated as the area under the traction force-displacement separation curve within the current increment step, which can be approximated as... ,in and These are the interface traction forces at the start and end of the incremental step, respectively. Then, the calculated artificial energy dissipation data... As a numerator, the data on energy dissipation due to interface failure will be used. As the denominator, a division operation is performed on each node with artificial viscous damping (i.e., the node whose Boolean flag variable is true) to obtain the energy ratio. For nodes without artificial adhesion, this ratio is recorded as 0. Finally, the energy ratios of all interface region nodes are calculated. Arrange them according to their node numbers to form a one-dimensional array, and generate a set of energy dissipation error comparison arrays.

[0070] The steps for obtaining the prediction results of debonding fatigue life of aircraft structural components are as follows:

[0071] Based on the dissipated energy error comparison array set, the quotient value corresponding to each interface region node is read item by item. The quotient value is compared with a set threshold, and interface region nodes with quotient values ​​lower than the set threshold are filtered out. The artificial viscous damping factor value corresponding to the interface region node with a quotient value lower than the set threshold is read. The artificial viscous damping factor value is gradually decayed according to the time variable. The decayed artificial viscous damping factor value is rewritten into the position of the diagonal constant term element of the local tangent stiffness matrix to generate the reconstructed control matrix.

[0072] Based on the reconstructed control matrix, the node traction force balance terms associated with the control matrix are read item by item according to the node number of the interface area, the connection relationship between adjacent coordinate nodes, and the position of the diagonal constant term element of the local tangent stiffness matrix. The difference of the node traction force balance terms of each adjacent coordinate node is checked, and the node traction force balance terms that meet the equilibrium convergence limit are determined as the system equilibrium traction force parameters, thus obtaining the system equilibrium traction force parameters.

[0073] Based on the connection relationship between adjacent coordinate nodes, the values ​​of the system balance traction force parameters are identified segment by segment and reduced to the range of the zero-distribution boundary. The node number, incremental step number, and load cycle number of the interface area covered by the zero-distribution boundary are extracted. The total number of fatigue cycle load operations corresponding to the zero-distribution boundary is used as the basis for determining the occurrence of debonding, and the debonding fatigue life prediction result of the aircraft structural component is generated.

[0074] Specifically, based on the dissipation energy error comparison array set, the quotient value, i.e. the energy ratio, is read item by item for each interface region node. This quotient is then compared to a preset threshold. This threshold is designed to control the proportion of artificial viscous energy in the total dissipated energy to ensure the physical accuracy of the calculation. Its value is typically set based on numerical accuracy requirements and experience; for example, it could be set to 5%, or 0.05. During the comparison, all interface nodes with quotients lower than the preset threshold are selected. For example, if a node's calculated quotient is 0.03, which is lower than 0.05, then that node is selected. For these selected nodes, the artificial viscous damping factor value used in the current increment step is read. And based on the time variable, i.e., the increment step, a decay operation is performed incrementally. The purpose of decay is to gradually remove the artificial damping after the system recovers stability, so as to avoid its continued impact on subsequent calculations. The decay formula can adopt linear decay, for example... The attenuation coefficient It is a small positive number, such as 0.1, which means that the attenuation is 10% at each step. The new artificial viscous damping factor value is obtained after attenuation calculation. This is done by rewriting the diagonal constant element of the local tangent stiffness matrix corresponding to that node, i.e., updating... and For nodes with quotient values ​​higher than or equal to a set threshold, the artificial viscous damping factor value is kept unchanged. After conditionally updating the artificial viscous damping factor values ​​of all nodes, an adjusted stiffness matrix set is obtained, and the reconstructed control matrix is ​​generated.

[0075] Based on the reconstructed control matrix, following the node numbering order of the interface region and combining the predefined connection relationships between adjacent coordinate nodes (e.g., on a straight line, node i is adjacent to nodes i-1 and i+1), the node traction force balance term associated with each node in the reconstructed control matrix is ​​read item by item. This balance term is actually the interface traction force T borne by the node at the end of the current increment step. This traction force is calculated based on the final converged displacement and damage state through the interface constitutive relation. Then, for each node, the difference between the node traction force balance terms of its adjacent coordinate nodes is checked. Specifically, the traction force of node i is calculated. traction force with adjacent node i-1 The difference between This difference is then compared to a preset equilibrium convergence limit, which is used to determine whether the spatial distribution of the traction force is smooth, indicating that stress concentration has been released. Numerically, this limit can be set to 0.1% of the nominal strength of the interface material; for example, if the material strength is 50 MPa, the limit is 0.05 MPa. If the difference... If the traction force of the two adjacent nodes is less than the equilibrium convergence limit, it is considered that the traction force of the two adjacent nodes has reached spatial equilibrium. The traction force values ​​of these nodes that meet the equilibrium conditions are determined as the system equilibrium traction force parameters. These parameters together constitute a stress distribution diagram describing the interface reaching a stable or failed state under fatigue load, thus obtaining the system equilibrium traction force parameters.

[0076] Based on the connection relationship between adjacent coordinate nodes, the system equilibrium traction force parameters obtained in the previous process are scanned and identified segment by segment. The system equilibrium traction force parameters are a series of traction force values ​​distributed along the interface. The segment-by-segment identification process involves checking the changing trend of the traction force values ​​along the interface path. The purpose is to find the boundary where the traction force value transitions from a non-zero value region to a region with a value close to zero, i.e., the zero-value distribution boundary. Here, "zero value" does not mean absolute zero, but rather that its value is less than a preset judgment threshold. This threshold is usually set to be a tiny amount much smaller than the material strength, such as 0.01% of the material strength, to distinguish between the load-bearing area and the complete failure area. Once such a zero-value distribution boundary is identified, the following steps are taken: Take the node numbers of all nodes within the interface region covered by the boundary, as well as the current incremental step number and the corresponding load cycle number. Here, "coverage" refers to the continuous node region along the crack propagation direction where the traction force is continuously below the zero threshold, starting from the node where the traction force first drops below the zero threshold. Then, record the total number of fatigue cycle load runs corresponding to the first appearance of this zero-value distribution boundary. This total number of cycles is used as the criterion for determining the occurrence of debonding, because it indicates that a macroscopically defined crack region no longer transmits load has appeared on the interface. Use this recorded total number of cycles as the final prediction result to generate the debonding fatigue life prediction result for aircraft structural components.

[0077] 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 fatigue life assessment method for laser-repaired aircraft structural components considering interface effects, characterized in that, Includes the following steps: Read the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly, separate the repair layer node coordinate array and the matrix node coordinate array, and generate a two-dimensional cohesive element coordinate tensor. The tangential and normal displacement components of the internal nodes of the coordinate tensor of the two-dimensional cohesive unit are read, and the constitutive characteristic tensor of the repair layer matrix interface is calculated. An external fatigue cyclic load spectrum is introduced, and the cumulative plastic slip parameters corresponding to the internal node variables of the constitutive feature tensor of the repair layer matrix interface are extracted to generate a cumulative distribution map of interface plastic slip. Based on the associated values ​​of the coordinate nodes of the corresponding regions of the cumulative distribution map of interface plastic slip, a new set of parameters is generated. The internal element array of the new set of parameters covering the constitutive feature tensor of the repair layer matrix interface is extracted and compiled to generate a fatigue cyclic degradation damage evolution matrix. Extract the fatigue cycle degradation damage evolution matrix, perform incremental iterative calculation of the equation, record the Boolean flag variable whose diagonal element value drops below the extreme value limit during the iteration process, and calculate when the Boolean flag variable corresponds to the true state to generate an adaptive viscous regularized stiffness control matrix. Based on the adaptive viscous regularized stiffness control matrix, obtain the dissipated energy error comparison array set; If the quotient of the dissipated energy error comparison array set is lower than a set threshold, the artificial viscous damping factor value is decayed by calling the time variable formula to generate the debonding fatigue life prediction result of the aircraft structural component.

2. The fatigue life assessment method for laser-repaired aircraft structural components considering interface effects according to claim 1, characterized in that, The steps for obtaining the constitutive feature tensor of the repair layer matrix interface are as follows: Read the geometric topology data node coordinates of the laser-repaired aircraft structural component assembly. According to the repair layer identifier, base identifier, and adjacency relationship at the interface, match the node number, spatial coordinates, and unit connection relationship item by item. Separate the repair layer node coordinate array and the base node coordinate array. Insert a dimension-two thickness zero node mesh element array at the same projection position at the interface. Write the initial interface penalty stiffness value and critical fracture energy value to each zero thickness node mesh element to obtain the dimension-two thickness zero node mesh element array. Based on the aforementioned two-dimensional zero-thickness mesh element array, structural boundary conditions and initial static loads are applied. According to the node number index of the zero-thickness node mesh element, the characteristic vector of the mesh shape function node and the node displacement tensor at the interface are extracted. The node displacement tensor is mapped to the local coordinate direction of the zero-thickness node mesh element. The coordinate order, degree of freedom number, and displacement component correspondence of adjacent nodes are checked to generate a two-dimensional cohesive element coordinate tensor. Read the tangential and normal displacement components of the internal nodes of the coordinate tensor of the two-dimensional cohesive unit, calculate the quotient of the cross-sectional contact stress component and the interface displacement separation component according to the local coordinate direction, update the constant term of the local tangential stiffness matrix of each internal node item by item according to the quotient, and write it into the interface constitutive feature position according to the node number order to obtain the interface constitutive feature tensor of the repair layer matrix.

3. The fatigue life assessment method for laser-repaired aircraft structural components considering interface effects according to claim 1, characterized in that, The steps for obtaining the new-order parameter set are as follows: An external fatigue cyclic load spectrum is introduced. According to the load peak value, load valley value and cyclic loading order, the internal node variables of the constitutive feature tensor of the matrix interface of the repair layer are read one by one. The tangential slip increment of each internal node in adjacent load cycles is extracted. The slip increment of the same internal node is accumulated and the corresponding load cycle number is recorded. The load cycle number and cumulative plastic slip parameters are written into the hysteresis dissipation energy statistics list line by line. The hysteresis dissipation energy statistics list variable values ​​are mapped according to the coordinate node order of the boundary area to generate the interface plastic slip cumulative distribution map. The associated values ​​of the coordinate nodes corresponding to the cumulative distribution map of plastic slip at the interface are read item by item. The change range of the cumulative plastic slip parameter of the same coordinate node under the variable of the number of continuous load cycles is statistically analyzed. The coordinate nodes whose change range reaches the preset fatigue damage triggering threshold are marked as decreasing control nodes. The dynamic damage decreasing control coefficient is calculated based on the associated values ​​of the decreasing control nodes. The dynamic damage decreasing control coefficient is multiplied by the initial interface penalty stiffness value and the critical fracture energy value to generate a new set of parameters.

4. The fatigue life assessment method for laser-repaired aircraft structural components considering interface effects according to claim 1, characterized in that, The steps for obtaining the fatigue cycle degradation damage evolution matrix are as follows: According to the node number, coordinate node position, and interface constitutive feature position of the internal element array of the interface constitutive feature tensor of the repair layer matrix, the corresponding updated penalty stiffness value and updated fracture energy value in the new order parameter set are extracted item by item. The updated penalty stiffness value is used to cover the constant term parameter of the local tangent stiffness matrix, and the updated fracture energy value is used to cover the interface damage energy position. The covered internal element array is then matrix-arranged according to the numerical variable of the load cycle number to generate the fatigue cycle degradation damage evolution matrix.

5. The fatigue life assessment method for laser-repaired aircraft structural components considering interface effects according to claim 1, characterized in that, The steps for obtaining the Boolean flag variable are as follows: The fatigue cyclic degradation damage evolution matrix is ​​extracted and divided into continuous incremental step input items according to the numerical variable of the load cycle number, the node number of the interface region, and the position of the internal element array. For each continuous incremental step input item, the degradation stiffness component, damage evolution component, and node displacement component of the interface region node are read. The residual balance of the interface region node is calculated. The node displacement component is corrected according to the residual balance. When the residual balance reaches the preset balance convergence limit, the diagonal element of the local tangent stiffness matrix of the corresponding incremental step is recorded to obtain the diagonal element of the local tangent stiffness matrix. Based on the diagonal elements of the local tangent stiffness matrix, the values ​​of the diagonal elements of the local tangent stiffness matrix are compared item by item with the set negative instability extreme value limit according to the increment step number, interface region node number, and local tangent stiffness matrix position index. If the value of the diagonal elements of the local tangent stiffness matrix drops below the negative instability extreme value limit, the corresponding interface region node, increment step number, and diagonal element value are written to the true state. If the value of the diagonal elements of the local tangent stiffness matrix does not drop below the negative instability extreme value limit, it is written to the false state, and a Boolean flag variable is obtained.

6. The fatigue life assessment method for laser-repaired aircraft structural components considering interface effects according to claim 1, characterized in that, The steps for obtaining the adaptive viscous regularized stiffness control matrix are as follows: Filter the interface area nodes corresponding to the true state of the Boolean flag variable, read the corresponding incremental step number, specified base factor, and time-related variable, perform exponential multiplication operation on the specified base factor according to the time-related variable, retain the original diagonal element records for the interface area nodes corresponding to the false state of the Boolean flag variable, extract the exponential multiplication operation results corresponding to the true state interface area nodes, and obtain the value of the artificial viscous damping factor. According to the node number of the true state interface region, the incremental step number, and the position of the diagonal constant term element of the local tangent stiffness matrix, the value of the artificial viscous damping factor is superimposed into the position of the diagonal constant term element of the local tangent stiffness matrix item by item. For the nodes of the false state interface region, the original value of the diagonal constant term element of the local tangent stiffness matrix is ​​retained. The superimposed local tangent stiffness matrix is ​​rearranged according to the node number of the interface region to generate an adaptive viscous regularized stiffness control matrix.

7. The fatigue life assessment method for laser-repaired aircraft structural components considering interface effects according to claim 1, characterized in that, The steps for obtaining the dissipated energy error comparison array set are as follows: The adaptive viscous regularized stiffness control matrix is ​​extracted. According to the position of the internal constant term element, the node number of the interface region, and the incremental step number, the artificial dissipation energy data and the interface damage dissipation energy data associated with the internal constant term are read one by one. The artificial dissipation energy data is used as the numerator and the interface damage dissipation energy data is used as the denominator. The division operation is performed node by node. The quotient values ​​are arranged according to the interface region node number to generate a dissipation energy error comparison array set.

8. The fatigue life assessment method for laser-repaired aircraft structural components considering interface effects according to claim 1, characterized in that, The steps for obtaining the predicted debonding fatigue life of the aircraft structural components are as follows: Based on the dissipated energy error comparison array set, the quotient value corresponding to each interface region node is read item by item, the quotient value is compared with a set threshold, interface region nodes with quotient values ​​lower than the set threshold are filtered, the artificial viscous damping factor value corresponding to the interface region node with quotient value lower than the set threshold is read, the artificial viscous damping factor value is gradually decayed according to the time variable, and the decayed artificial viscous damping factor value is rewritten into the position of the diagonal constant term element of the local tangent stiffness matrix to generate the reconstructed control matrix. Based on the reconstructed control matrix, the node traction force balance terms associated with the control matrix are read item by item according to the node number of the interface area, the connection relationship between adjacent coordinate nodes, and the position of the diagonal constant term element of the local tangent stiffness matrix. The difference of the node traction force balance terms of each adjacent coordinate node is checked, and the node traction force balance terms that satisfy the equilibrium convergence limit are determined as the system equilibrium traction force parameters, thus obtaining the system equilibrium traction force parameters. Based on the connection relationship between adjacent coordinate nodes, the system balance traction force parameter values ​​are identified segment by segment until they drop to the range of the zero-value distribution boundary. The node number, increment step number, and load cycle number of the interface area covered by the zero-value distribution boundary are extracted. The total number of fatigue cycle load operations corresponding to the zero-value distribution boundary is used as the basis for determining the occurrence of debonding, and the debonding fatigue life prediction result of the aircraft structural component is generated.

Images