Gearbox structure stiffness detection method based on multi-body dynamics

Abstract

This invention discloses a method for detecting the stiffness of gearbox structural components based on multibody dynamics, belonging to the field of gearbox structural component stiffness detection technology. The method includes the following steps: establishing a multibody dynamics topological model of the gearbox, and condensing the gearbox finite element model into a reduced-order flexible body matrix; generating time-varying input boundary loads based on the operating load, and jointly solving the time-varying contact stiffness of the gears, the nonlinear support force of the bearings, and the flexible response of the gearbox within the same time step to obtain the dynamic contact force and support reaction force; transferring the load to the finite element mesh nodes through cross-domain mapping to complete the transient stress-strain calculation; and finally determining the stiffness compliance by combining the dynamic deformation and dynamic stiffness influence factor. This invention solves the problem of the inability to capture dynamic coupling effects caused by unidirectional series connection.

Description

Technical Field

[0001] This invention relates to the field of gearbox structural component stiffness testing technology, specifically a method for testing gearbox structural component stiffness based on multibody dynamics. Background Technology

[0002] As the capacity of single wind turbine units continues to increase, the dynamic loads borne by gearbox structural components are becoming increasingly complex. The gear meshing stiffness varies periodically with the meshing phase, and the bearing support stiffness varies nonlinearly with the eccentric load state. There is a strong coupling effect between these two and the flexibility of the gearbox structure, making it difficult to accurately characterize the true deformation state of the structural components through static analysis. At the same time, the contradiction between the engineering requirements for simulation of dynamic response under all working conditions and the limited computing resources is becoming increasingly prominent.

[0003] To address the aforementioned issues, existing technologies typically employ two approaches: one is based on the static equivalence approach, applying the transmission chain excitation force to the finite element model of the box under constant boundary conditions for stiffness verification; the other is to use commercial multibody dynamics software, such as ADAMS and Simpack, to extract the transmission chain excitation force, which is then manually converted and unidirectionally input into the finite element model as a quasi-static load for structural analysis.

[0004] However, the above-mentioned methods have the following shortcomings: the static equivalent method cannot capture the dynamic coupling effect between the time-varying contact stiffness of the gear and the nonlinear support of the bearing, resulting in a significant deviation between the verification results and the actual dynamic deformation; in the unidirectional serial process of multibody dynamics and finite element method, modal truncation often leads to distortion of local high-frequency contact characteristics, and the artificial conversion process of discrete dynamic load to continuous medium mesh has the problem of boundary condition distortion; in addition, the existing methods lack a quantitative evaluation mechanism that integrates dynamic deformation, vibration amplitude and fatigue damage, and the stiffness compliance judgment relies on engineering experience and lacks objective basis.

[0005] Therefore, this invention proposes a method for detecting the stiffness of gearbox structural components based on multibody dynamics. Summary of the Invention

[0006] The purpose of this invention is to provide a method for detecting the stiffness of gearbox structural components based on multibody dynamics, so as to solve the existing problems mentioned in the background art.

[0007] To achieve the above objectives, the present invention provides the following technical solution: a method for detecting the stiffness of gearbox structural components based on multibody dynamics, comprising the following steps: S1. Establish a multibody dynamics topology model, extract and matrix condense the boundary degrees of freedom at the connection between the box and the transmission chain in the box finite element model, generate a reduced-order flexible body matrix that retains the principal modal characteristics, and import the reduced-order flexible body matrix into the multibody dynamics topology model; S2. Configure load fingerprint parameters including impact coefficient and time-varying attenuation factor, and generate system input boundary load that varies with operating conditions based on the load fingerprint parameters; S3. Under the drive of the system input boundary load, the gear is discretized into multiple micro-element slices along the tooth width direction based on the slice theory and the time-varying contact stiffness of each micro-element is calculated in real time. The bearing nonlinear contact mechanics model is called, and the time-varying contact stiffness, the bearing nonlinear contact mechanics model, and the reduced-order flexible body matrix are solved in the same time step to obtain the dynamic contact force and support reaction force sequence in the full time domain. S4. Spatial interpolation is performed on the dynamic contact force and support reaction force sequence through the cross-domain mapping interface, and the interpolation results are mapped to the mesh nodes of the box finite element model to perform transient stress and strain calculation. S5. Based on the results of the transient stress-strain calculation, extract the dynamic deformation of the structural component, calculate the difference between the dynamic deformation and the preset design allowable deformation limit to obtain the remaining deformation space, and determine the stiffness compliance in conjunction with the dynamic stiffness influence factor evaluation system.

[0008] A further improvement of this invention is that the generation of the reduced-order flexible body matrix that preserves the principal modal characteristics includes: A finite element mesh is generated from the three-dimensional geometric model of the box to form a finite element model of the box, and the boundary degree-of-freedom nodes at each connection between the box and the transmission chain are identified and calibrated. Based on the boundary degree-of-freedom nodes, a matrix reduction method is used to perform matrix reduction. During the reduction process, the main modal characteristics of the box are preserved and the degrees of freedom corresponding to high-frequency redundant nodes are removed. The reduced-order flexible body matrix is ​​subjected to modal deviation verification. After confirming that the deviation between the main modal frequency and the full-order finite element model is within the preset allowable range, the reduced-order flexible body matrix is ​​output.

[0009] A further improvement of this invention is that the process of generating the system input boundary load that varies with operating conditions based on load fingerprint parameters includes: The system obtains the rotational speed time series and torque time series from the wind load parameters, performs time-domain alignment and filtering on the rotational speed time series and torque time series, configures the impact coefficient to characterize the amplification characteristics of transient load abrupt change amplitude and the time-varying attenuation factor to characterize the attenuation evolution characteristics of the load in the time domain, and modulates the reference load sequence based on the load fingerprint parameters composed of the impact coefficient and the time-varying attenuation factor to generate the system input boundary load.

[0010] A further improvement of this invention is that the equation solving within the same time step includes: S31. Simultaneously perform the calculation of gear time-varying contact stiffness based on slice theory and the calculation of bearing dynamic support force based on bearing nonlinear contact mechanics model within each time step. S32. The calculation results of the time-varying contact stiffness of the gear and the calculation results of the dynamic support force of the bearing are used as internal force inputs, and the dynamic equations are solved together with the reduced-order flexible body matrix within the current time step. S33. Run the joint dynamic equation solution process step by step until the entire simulation time domain is covered, and obtain the dynamic contact force sequence and support reaction force sequence.

[0011] A further improvement of this invention is that the calculation process for the time-varying contact stiffness of the gear includes: The tooth width, module, and pressure angle parameters of the gear pair are obtained, and the gear is divided into a preset number of micro-element slices at equal intervals along the tooth width direction. For each micro-element slice, the instantaneous contact stiffness of the micro-element slice is calculated based on the meshing phase corresponding to the current time step. The instantaneous contact stiffness of each micro-element slice is integrated and superimposed along the tooth width direction to obtain the time-varying contact stiffness of the gear pair at the current time step. The number of micro-element slices is determined according to the tooth width size and calculation accuracy requirements, and the instantaneous contact stiffness is updated in real time with the change of meshing phase.

[0012] A further improvement of this invention is that the specific construction process of the bearing nonlinear contact mechanics model includes: The number of rolling elements, contact angle, and clearance parameters of the bearing are obtained, and the relative displacement of the inner and outer rings of the bearing at the current time step is obtained. Based on the number of rolling elements, contact angle, clearance parameters, and relative displacement, the contact deformation of each rolling element is calculated. The contact force of each rolling element is calculated by applying a nonlinear stiffness-deformation relationship to the contact deformation of each rolling element. The contact forces of each rolling element are vector-superimposed to obtain the dynamic support resultant force and resultant torque of the bearing. When the contact deformation of the rolling element is negative, the corresponding contact force is set to zero.

[0013] A further improvement of this invention is that the step of performing spatial interpolation and mapping to the mesh nodes of the box finite element model through the cross-domain mapping interface includes: S41. Using the spatial coordinates of each force point in the multibody dynamics topology model and the spatial coordinates of each mesh node in the box finite element model as input, establish the mapping topology relationship between the force point and the mesh node and the corresponding spatial interpolation weight coefficient based on spatial proximity search, and output the mapping topology table. S42. Using the mapping topology table and the full-time dynamic contact force and support reaction force sequence as input, perform spatial interpolation calculation on the force vector of each time step according to the spatial interpolation weight coefficient, and distribute the interpolated force vector to the mesh node corresponding to the box finite element model, and output the distributed equivalent node force load sequence of each time step. S43. Using the distributed equivalent nodal force load sequence and the dynamic contact force and support reaction force sequence before mapping as input, perform numerical comparison and verification on the resultant force and resultant moment before and after mapping respectively. When the deviation exceeds the preset conservation tolerance, return to step S42 to adjust the spatial interpolation weight coefficient. When the deviation is within the preset conservation tolerance, output the distributed equivalent nodal force load sequence that has passed the conservation verification. S44. The distributed equivalent nodal force load sequence that has passed the conservation check is used as the load boundary condition input to the box finite element model, and transient stress and strain calculation is performed to output the transient stress field and strain field of the box structure in the whole time domain.

[0014] A further improvement of this invention is that the step of extracting the dynamic deformation of the structural component and calculating the remaining deformation space to determine stiffness compliance includes: Extract the maximum dynamic deformation of the structural component in the entire time domain from the results of the transient stress-strain calculation; obtain the preset design allowable deformation upper limit, calculate the difference between the maximum dynamic deformation and the design allowable deformation upper limit to obtain the remaining deformation space; determine whether the remaining deformation space meets the preset minimum remaining deformation space threshold, and determine that the stiffness is non-compliant when the remaining deformation space is less than the minimum remaining deformation space threshold.

[0015] A further improvement of this invention is that, after determining the stiffness compliance, it also includes a step of performing closed-loop correction on the multibody dynamics topology model, specifically including: Acquire measured data from the physical test bench, including fiber optic strain data and laser vibration data; The measured data and the calculation output of the simulation model under the same working conditions are time-domain aligned and feature-compared. A system parameter dimensionality reduction matrix is ​​constructed to map the parameters to be corrected to a low-dimensional space. Based on the dimensionality reduction matrix, a proxy mapping relationship model between the system parameters and the simulation output is established. Taking the deviation between the measured data and the calculated output as the target, the stiffness parameters and damping parameters in the simulation model are updated through inverse iteration. Finally, the corrected stiffness and damping parameters are substituted back into the reduced-order flexible body matrix and the bearing nonlinear contact mechanics model.

[0016] A further improvement of this invention is that the stiffness compliance determination process further includes: Transient stiffness variables are extracted from the full-time dynamic contact force sequence and support reaction force sequence. The vibration amplitude obtained from the transient stress-strain calculation and the fatigue life reduction coefficient of the bearing and gear are obtained. A dynamic stiffness influence factor evaluation matrix is ​​constructed. The transient stiffness variables, the vibration amplitude and the fatigue life reduction coefficient are used as input elements to perform fusion calculation to generate a quantitative ranking index. Based on the quantitative ranking index, the stiffness weakness of each structural region is ranked, and the comprehensive stiffness compliance judgment result is jointly output with the judgment result of the remaining deformation space.

[0017] Compared with the prior art, the beneficial effects of the present invention are: 1. This invention first solves the problem that the dynamic coupling effect cannot be captured in the traditional method because the three factors of housing flexibility, gear time-varying contact and bearing nonlinear contact mechanics model are independent and unidirectionally connected, by jointly solving the reduced-order flexible body matrix, the time-varying contact stiffness of gear based on slice theory and the nonlinear contact mechanics model of bearing within the same time step. This realizes real-time dynamic feedback between housing structure deformation and transmission chain contact characteristics, and improves the calculation reliability of dynamic response under all working conditions.

[0018] 2. By using a cross-domain mapping interface, the dynamic contact force and support reaction force in the discrete domain of multibody dynamics are spatially interpolated and then mapped to the continuous mesh nodes of the finite element method after applying a resultant force conservation check. This solves the problem of load boundary distortion introduced by the traditional method of manually applying concentrated force, and realizes the high-fidelity transfer of multibody dynamics solution results to finite element transient stress and strain calculation, providing a reliable load input basis for the accurate extraction of dynamic deformation of structural components. Attached Figure Description

[0019] Figure 1 This is a flowchart of the gearbox structural component stiffness detection method based on multibody dynamics according to the present invention. Detailed Implementation

[0020] The technical solution of the present invention will be described in detail below with reference to the accompanying drawings and specific embodiments. It should be understood that the embodiments of the present invention and the specific features in the embodiments are detailed descriptions of the technical solution of the present invention, rather than limitations thereof. In the absence of conflict, the embodiments of the present invention and the technical features in the embodiments can be combined with each other.

[0021] The term "and / or" is merely a description of the relationship between related objects, indicating that there can be three relationships. For example, A and / or B can represent three cases: A exists alone, A and B exist simultaneously, and B exists alone.

[0022] Example 1 This embodiment provides a method for stiffness testing of gearbox structural components based on multibody dynamics, targeting the dynamic stiffness verification of gearbox structural components in large-megawatt wind turbine generators. The implementation platform is a high-performance computing workstation equipped with an Intel Xeon W-2295 processor and 512GB ECC memory. Multibody dynamics simulation is performed in the MSC Adams solution environment, and finite element analysis is performed in the ANSYS Mechanical solution environment. Data transfer between the two is achieved through a self-developed cross-domain mapping interface program. This embodiment selects an 8MW wind turbine gearbox as the verification object, which includes a single-stage planetary gear train and two stages of parallel shaft gear pairs.

[0023] When gearbox power ratings reach 8MW and above, the time-varying nonlinear effects of gear meshing, the nonlinear stiffness of bearings varying with eccentric loads, and the dynamic coupling between the gearbox housing flexibility and the transmission chain intertwine to jointly determine the actual deformation state of the structural components. At this point, if a static equivalent method is used, these coupling effects are completely ignored; if a full-order finite element transient analysis is employed, the computational cost of millions of degrees of freedom under multiple operating conditions and long time domains becomes unacceptable. To address this issue, this implementation provides the following methodological steps, such as... Figure 1 As shown: Step S1: Establish a multibody dynamics topology model, extract and matrix condense the boundary degrees of freedom at the connection between the box and the transmission chain in the box finite element model, generate a reduced-order flexible body matrix that retains the principal modal characteristics, and import the reduced-order flexible body matrix into the multibody dynamics topology model.

[0024] First, the complete topology of the gearbox is built in a multibody dynamics solution environment, covering the sun gear, planet gears, ring gear, and planet carrier in the planetary gear train, each gear pair and drive shaft in the parallel stage, as well as all bearing support positions and gearbox structural components. The components are connected by revolute pairs, prismatic pairs, or coupling constraints to form a complete kinematic transmission path. Then, the finite element mesh of the gearbox's three-dimensional geometric model is generated, using ten-node tetrahedral elements (C3D10). After the mesh is generated, the boundary degree-of-freedom nodes at each connection interface between the gearbox and the transmission chain are identified and calibrated, including the inner surface nodes of each bearing housing hole, the mating surface nodes of the ring gear, and the bolt hole nodes of the base.

[0025] Subsequently, the modal synthesis condensation method was used to reduce the order of the full-order model. The system degrees of freedom were split into boundary degrees of freedom and internal degrees of freedom. The boundary degrees of freedom were fully preserved to maintain the accuracy of force transmission, while the internal degrees of freedom were compressed through modal truncation under fixed interface conditions. During the condensation process, the dominant modal characteristics of the box in the 0 to 2000 Hz frequency band were preserved, while redundant degrees of freedom corresponding to higher frequency bands were eliminated. After the order reduction was completed, the reduced-order flexible body matrix with several internal generalized coordinates added to the compressed boundary degree of freedom nodes was imported into the multibody dynamics topology model in MNF format to replace the rigid box components.

[0026] The reduced-order flexible body matrix output in step S1 provides the mathematical basis for the box's flexible participation in the dynamic response in the subsequent joint solution, while the external excitation source required for the joint solution is constructed in step S2.

[0027] Step S2: Configure load fingerprint parameters including impact coefficient and time-varying attenuation factor, and generate system input boundary load that varies with operating conditions based on load fingerprint parameters.

[0028] The rotational speed and torque time series from the wind load parameters are obtained. In this embodiment, a 600-second gust wind condition segment is extracted, with a sampling frequency of 100Hz. After performing time-domain alignment and low-pass filtering on the two sequences, the load fingerprint parameters are configured, including: Impact coefficient This value is used to characterize the magnitude amplification during transient load changes. In this embodiment, a value of 1.35 is preferred. This value is obtained through statistical analysis of the ratio of peak transient torque to rated torque in gusts in the service data of similar units. It can cover most typical gust events without being overly conservative. Time-varying decay factor The load fingerprint parameter is used to characterize the evolution of the load from its peak value to its descent in the time domain. In this embodiment, a value of 0.92 exponential decay rate per second is preferred. This value is determined through analysis of the torsional vibration mode damping ratio of the transmission system, reflecting the process by which the torsional vibration energy of the transmission chain decays to less than 5% of its initial value within approximately 12 seconds. The reference load sequence is modulated based on the load fingerprint parameters composed of the impact coefficient and the time-varying decay factor to generate the system input boundary load.

[0029] The time-varying load generated in step S2 provides an external excitation input with transient characteristics for the joint solution in step S3. If the characterization of the impact and decay process is lacking, the subsequent calculation will lack transient peak response information.

[0030] Step S3: Under the drive of the system input boundary load, the gear is discretized into multiple micro-element slices along the tooth width direction based on the slice theory and the time-varying contact stiffness of each micro-element is calculated in real time. The bearing nonlinear contact mechanics model is called, and the time-varying contact stiffness and the bearing nonlinear contact mechanics model are solved in the same time step with the reduced-order flexible body matrix to obtain the dynamic contact force and support reaction force sequence in the full time domain.

[0031] Within each time step, this embodiment can take... First, based on the slicing theory, the gear is discretized into multiple micro-element slices along the tooth width direction. Taking a planetary gear as an example, with a tooth width of 280mm, it is discretized into 28 micro-element slices, each with a width of 10mm. For each micro-element slice, its instantaneous contact stiffness is calculated based on the current meshing phase. The instantaneous contact stiffness of each micro-element is integrated and superimposed along the tooth width direction to obtain the time-varying contact stiffness of the gear pair at the current time step. At the same time, the nonlinear contact mechanics model of the bearing calculates the contact deformation of each rolling element based on the relative displacement of the inner and outer rings of the bearing at the current time step and obtains the dynamic support resultant force of the bearing. Then, the calculation results of the time-varying contact stiffness of the gear and the calculation results of the dynamic support force of the bearing are used as internal force inputs and solved together with the reduced-order flexible body matrix within the current time step for solving the joint dynamic equations. The process is run step by step until the entire simulation time domain is covered, and the dynamic contact force sequence and support reaction force sequence of the entire time domain are output.

[0032] By solving the problem simultaneously within the same time step, the changes in gear contact force are immediately fed back to the flexible deformation of the housing. The housing deformation then instantly changes the relative displacement of the bearings and the meshing state of the gears, forming a dynamic closed-loop coupling. However, the output of step S3 is the force vector at the discrete hinge points in the multibody dynamics model. These force application points do not correspond to the finite element mesh nodes. Therefore, the cross-domain mapping in step S4 is needed to convert them into distributed nodal forces acceptable to the finite element.

[0033] Step S4: Spatial interpolation of the dynamic contact force and support reaction force sequence is performed through the cross-domain mapping interface. The interpolation results are then mapped onto the mesh nodes of the box finite element model for transient stress and strain calculation.

[0034] Using the spatial coordinates of each force application point in the multibody dynamics topology model and the spatial coordinates of each mesh node in the box finite element model as input, a mapping topology relationship between the force application points and mesh nodes and the corresponding spatial interpolation weight coefficients are established based on spatial nearest neighbor search, and a mapping topology table is output. Based on the mapping topology table, spatial interpolation calculations are performed on the force vectors at each time step, and the interpolated force vectors are distributed to the corresponding mesh nodes of the box finite element model, outputting the distributed equivalent nodal force load sequence for each time step. Conservation checks are performed on the resultant force and resultant moment before and after mapping. In this embodiment, the resultant force deviation tolerance is 0.5%, and the resultant moment deviation tolerance is 1.0%. If the deviation exceeds the limit, the interpolation weight coefficients are adjusted until the check passes. The distributed equivalent nodal force load sequence that has passed the conservation check is used as the load boundary condition input to the box finite element model, and transient stress-strain calculations are performed, outputting the transient stress field and strain field of the box structure in the entire time domain.

[0035] Step S5: Extract the dynamic deformation of the structural component based on the results of transient stress-strain calculation, calculate the difference between the dynamic deformation and the preset design allowable deformation limit to obtain the remaining deformation space, and determine the stiffness compliance in conjunction with the dynamic stiffness influence factor evaluation system.

[0036] The maximum dynamic deformation of the structural components over the entire time domain is extracted from the transient stress-strain calculation results, with a focus on the deformation of the planetary carrier cantilever end face, the ellipticity deformation of the front housing bearing seat hole, and the normal deformation of the intermediate housing mating surface. A preset design allowable deformation upper limit is obtained (0.15mm for the planetary carrier cantilever end face, 0.05mm for the front housing bearing seat hole ellipticity, and 0.08mm for the normal deformation of the intermediate housing mating surface). The difference between the maximum dynamic deformation and the design allowable deformation upper limit is calculated to obtain the remaining deformation space. It is then determined whether the remaining deformation space meets a preset minimum remaining deformation space threshold. In this embodiment, 30% of the design allowable deformation upper limit can be used. If the remaining deformation space is less than this threshold, the stiffness is deemed non-compliant. Furthermore, transient stiffness variables are extracted from the full-time dynamic contact force sequence and support reaction force sequence. The vibration amplitude obtained from the transient stress-strain calculation and the fatigue life reduction coefficient of the bearing and gear are obtained. A dynamic stiffness influence factor evaluation matrix is ​​constructed. The above elements are integrated to generate a quantitative ranking index. Based on this index, the stiffness weakness of each structural region is ranked, and the comprehensive stiffness compliance judgment result is jointly output with the remaining deformation space judgment result.

[0037] This embodiment combines the reduced-order flexible body matrix, the time-varying contact stiffness of the gear, and the nonlinear contact mechanics model of the bearing within the same time step for joint solution, establishing a real-time dynamic feedback loop between the flexible deformation of the housing and the contact characteristics of the transmission chain. The cross-domain mapping interface replaces the manual application of concentrated force with automated interpolation and torque conservation verification, eliminating the hidden danger of load boundary distortion. The combined application of the residual deformation space criterion and the dynamic stiffness influence factor evaluation system provides a quantitative technical criterion for stiffness design that does not rely on subjective judgment based on engineering experience.

[0038] Example 2 This embodiment further refines and limits the key algorithms and parameters of each step in Embodiment 1, and is a preferred implementation. The 8MW wind turbine gearbox is still used as the object, and the computing platform is the same as in Embodiment 1.

[0039] Detailed implementation of generating the reduced-order flexible body matrix that preserves the principal modal characteristics in step S1: A finite element mesh is generated from the 3D geometric model of the housing to create a finite element model of the housing. Boundary degree-of-freedom nodes at each connection point between the housing and the transmission chain are identified and calibrated. Based on these boundary degree-of-freedom nodes, a modal synthesis reduction method is used for matrix reduction. During the reduction process, the principal modal characteristics of the housing are preserved, and the degrees of freedom corresponding to high-frequency redundant nodes are eliminated. This embodiment preserves the principal modes in the range of 0 to 2000 Hz. The cutoff frequency of 2000 Hz is selected based on the following: planetary-level meshing frequencies are distributed between 200 and 800 Hz, while parallel-stage high-speed-stage meshing frequencies can reach 1200 to 1500 Hz. The 2000 Hz cutoff frequency covers all meshing fundamental frequencies and major harmonics, ensuring that the dynamic response within the frequency band of concern in the engineering is not affected by order reduction.

[0040] Modal deviation verification is performed on the reduced-order flexible body matrix. After confirming that the deviation between the principal modal frequency and the full-order finite element model is within a preset allowable range, the reduced-order flexible body matrix is ​​output. In this embodiment, the allowable range is set to 3%—IEC61400-4 recommends a frequency accuracy of 5%, which is tightened to 3% in this embodiment to ensure that the cumulative error is controllable in long-time-domain transient solutions. If the deviation of a certain order frequency exceeds 3%, the number of boundary nodes in that region is appropriately increased, and the matrix is ​​re-condensed until all orders of interest meet the allowable requirements.

[0041] Detailed implementation of generating system input boundary loads based on load fingerprint parameters in step S2: The wind speed and torque time series from the wind load parameters are obtained, and time-domain alignment and filtering are performed on the two series. Time-domain alignment is based on the torque channel, and the speed channel is resampled after applying cross-correlation delay estimation. Filtering uses a fourth-order Butterworth low-pass filter with a cutoff frequency of 50Hz. The energy changes in the wind turbine gearbox input are concentrated within 50Hz, and this cutoff frequency retains all effective information while filtering out high-frequency noise.

[0042] Impact coefficient The calibration method is as follows: The ratio of the peak transient torque to the average steady-state torque under gust events is statistically analyzed from historical operating data of similar units, and the upper limit of the 95% confidence interval of this ratio is taken. This embodiment preferably uses... The reason for choosing 1.35 instead of a higher quantile is that 95% coverage already covers the vast majority of actual operating conditions, and further increasing it would lead to an overestimation of structural redundancy quality. When the operating environment deteriorates significantly, it can be adjusted upward to 1.45 to 1.55.

[0043] Based on the above parameters, the load modulation formula is as follows: in For a moment The system input boundary load (in Nm); This is the filtered baseline torque sequence (in Nm). The impact coefficient is dimensionless and is taken as 1.35. For the Heaviside step function, in The value is 1 if the condition is met and 0 otherwise, and is used to calibrate the first... The starting point of the shock event; For the first The time of occurrence of each impact event (in seconds), through the analysis of... Peak detection is applied, and a flag is triggered when the instantaneous torque exceeds 1.1 times the rated torque. The 1.1 times threshold defines the normal fluctuation range as within 10% above and below the rated value. Exceeding this value is considered the start of an impact. This represents the total number of impact events in the simulation time domain; The time-varying attenuation factor is dimensionless. This formula degenerates to a reference load in the absence of impact, but the load is instantaneously amplified to [value missing] upon impact. After exponential decay, the cumulative impact of continuous gusts can be reproduced by summing and superimposing multiple events.

[0044] Regarding fault tolerance: when an abnormally high number of events is detected, the trigger threshold is automatically increased to 1.15 to 1.2 times; when a certain event... Decay to below Remove it from the summation list to avoid floating-point precision issues.

[0045] Regarding the calculation process of the time-varying contact stiffness of the gear in step S3: Obtain the tooth width of the gear pair Modulus With pressure angle The gear is distributed at equal intervals along the tooth width direction. Each micro-element slice, slice width . Based on the tooth width dimensions and calculation accuracy requirements, the tooth width for the planetary star model in this embodiment is determined to be 280mm. (Slice width 10mm). The 10mm width was chosen because it is sufficient to capture the gradient of tooth surface bulging along the tooth width. Comparative verification showed that further refining to 5mm resulted in a contact force distribution deviation of less than 1% while doubling the computational load. 10mm represents a balance between accuracy and efficiency. For each micro-element slice... Based on the current engagement phase Calculate instantaneous contact stiffness Taking into account bending, Hertz contact, shear, and the stiffness components of the tooth root foundation, the stiffness is set to zero in the non-contact state to reflect the alternating single and double tooth meshing characteristics. The stiffness of each micro-element is integrally superimposed along the tooth width, and the calculation formula is expressed as: in The time-varying contact stiffness of the gear pair at the current time step (N / mm); For the first A microelement at time Instantaneous contact stiffness (N / mm / mm); The width of the slice is (mm). This integral can reflect the uneven distribution of load in the tooth width direction caused by profile modification, manufacturing errors or shaft deflection on a slice-by-slice basis, which is a key feature that constant average stiffness cannot characterize.

[0046] Regarding the construction process of the bearing nonlinear contact mechanics model in step S3: Obtain the number of rolling elements of the bearing Contact angle With radial clearance In this embodiment, the planetary gear bearing uses a double-row cylindrical roller bearing. (Separate column) , mm (C3 clearance group). Obtain the relative displacement vector between the inner and outer rings of the bearing at the current time step, and perform the following calculations for the first time step. A rolling element, Calculate its position angle radial contact deformation at the location ;when The contact force is calculated using a nonlinear stiffness-deformation relationship. The load-deformation index for the cylindrical roller is taken as 10 / 9. The contact stiffness coefficient is determined using the ISO / TS 16281 Palmgren formula combined with roller dimensions. The contact force is set to zero. The dynamic support resultant force and resultant moment are obtained by superimposing the contact force vectors of each rolling element. The setting of the unloading criterion allows the dynamic off-center load effect of the bearing caused by the flexible deformation of the housing to be naturally captured.

[0047] Detailed implementation of the cross-domain mapping interface in step S4: Step S41: For each point of force application, the kd-tree nearest neighbor algorithm is used to search for neighboring grid nodes. The search radius is three times the average side length of the grid in that region. This ratio ensures that enough nodes in the region participate in the interpolation without introducing far-field irrelevant contributions. The interpolation weights are inversely proportional to the distance. ,in From the point of application of the force to the first The distance (mm) between neighboring nodes is calculated, and the mapping topology table is output.

[0048] Step S42, based on the mapping topology table, the time... Force vector Assigned to neighboring nodes: Output the distributed equivalent nodal force load sequence.

[0049] Step S43: Verify that the resultant force deviation before and after mapping does not exceed 0.5% and the resultant torque deviation does not exceed 1.0%. The torque tolerance is slightly wider due to the lever arm amplification effect. When the limit is exceeded, the search radius of the force point that contributes the most to the deviation is expanded or a higher-order interpolation basis function is switched until it passes.

[0050] Step S44: The fully constrained mounting surface of the box base is superimposed with the verified load sequence, and the transient stress and strain calculation is performed using Newmark implicit integration with a time step of 0.0001s, outputting the full-time transient stress field and strain field.

[0051] Detailed implementation of dynamic deformation extraction and stiffness compliance determination in step S5: Extract the maximum dynamic deformation of each region of interest across the entire time domain. The upper limit of allowable deformation for each region and its source are as follows: planetary carrier cantilever end face. mm, obtained through sensitivity analysis of bearing life to planetary carrier deflection, shows that the life of bearing L10 is reduced by more than 15% when the limit is exceeded; front housing bearing seat ovality mm, when the limit is exceeded, the rolling element contact stress deviates from the design value by more than 20%; the normal direction of the intermediate box joint surface. mm, exceeding the limit, the gasket compression is insufficient to maintain the oil seal function. Calculate the remaining deformation space. The criterion is The selection of the 30% threshold needs to cover the combined effects of simulation simplification error (5% to 8%), material property dispersion (5% to 10%), and service load uncertainty (10% to 15%).

[0052] Further construct the dynamic stiffness influence factor evaluation matrix: extract transient stiffness variables, represented as The ratio of peak-to-valley difference to mean, vibration amplitude, peak value of nodal acceleration response, and fatigue life reduction factor for bearings and gears are calculated by inputting the dynamic load spectrum into the ISO281 bearing life model and the ISO 6336 gear strength model. The weights of each element are determined using the analytic hierarchy process (AHP): transient stiffness variable 0.4, vibration amplitude 0.35, and fatigue life reduction factor 0.25. After weighted normalization, a quantitative ranking index is generated, which, together with the remaining deformation space determination result, outputs a comprehensive stiffness compliance determination result. Simultaneously, the quantitative degree of non-compliant indicators is output for the design team to optimize accordingly.

[0053] Example 3 This embodiment provides alternative implementation methods for several steps S1 to S5, which constitute parallel optional variants with Embodiment 1 and Embodiment 2.

[0054] Regarding alternatives to the order reduction method in step S1: When the box geometry is relatively regular and the frequency range of interest is narrow, the Guyan static condensation method can be used instead of the modal synthesis condensation method. The Guyan method only retains the static constraint modes of the boundary degrees of freedom and does not retain the internal dynamic modes, resulting in faster order reduction and lower matrix dimensions, which can further accelerate the joint solution in step S3. However, its applicability is premised on the existence of at least a 20% frequency gap between the box's natural frequency and the meshing frequency; otherwise, the dynamic amplification effect will be significantly underestimated. Another alternative is a multi-level condensation method based on substructure partitioning, which is suitable for complex boxes consisting of multiple independent shell components connected by flange bolts: the box is divided into multiple substructures according to components, and each substructure is reduced in order and then assembled through interface coupling. When the design of a certain substructure changes, only that part needs to be updated and reassembled, reducing the amount of iterative repetitive calculations.

[0055] Regarding alternative solutions for the load modulation method in step S2: When the operating conditions are dominated by frequent periodic load fluctuations rather than occasional impact events, such as the tidal current cycle load of a tidal energy gearbox, the load fingerprint parameters can be replaced with a two-parameter periodic modulation model consisting of "amplitude modulation coefficient plus frequency modulation coefficient". The amplitude modulation coefficient characterizes the envelope of the load as the cycle changes, and the frequency modulation coefficient characterizes the frequency shift of the frequency component with the flow velocity. The output format of both solutions after step S2 is the same and does not affect subsequent steps. When historical operating data is lacking for calibrating the impact coefficient, the extreme operating condition time history of the design load conditions specified in Appendix IEC 61400-4 can be directly used as the benchmark. Both the impact coefficient and the attenuation factor are set to 1.0 to ensure consistency with the international standard system, which is suitable for compliance verification during the product certification stage.

[0056] Regarding alternative solutions for calculating gear contact stiffness in step S3: When the tooth surface profile is complex, including tooth-direction bulging, tooth-end thinning, and tooth profile edge trimming, a lookup table method can be used instead of analytical calculation. Finite element contact analysis is performed beforehand for different meshing positions and load levels, creating a two-dimensional stiffness lookup table indexed by meshing phase and load level. During runtime, each micro-element interpolates stiffness from the table, achieving higher accuracy than the analytical method, but requiring prior finite element calculations. Another alternative is an adaptive slicing strategy: When the tooth surface load is significantly biased towards one side of the tooth width, the resolution of equidistant slices is insufficient. Instead, the slice density is automatically adjusted based on the load gradient of each slice, refining areas with large gradients and widening areas with gentle gradients, improving the calculation accuracy of key areas without significantly increasing the total number of slices.

[0057] Regarding alternatives to the spatial interpolation method in step S4: When the mesh around the force application point is non-uniformly distributed, manifesting as a transition zone between a fine-mesh and a coarse-mesh region, inverse distance interpolation may result in unreasonable weight allocation. It can be replaced by an interpolation method based on radial basis functions (RBF), using thin-plate spline functions as basis functions to solve a system of linear equations to obtain the weights. RBF interpolation produces a more reasonable force distribution in non-uniform mesh regions, but its computational cost is higher than the inverse distance method.

[0058] Regarding alternatives to the evaluation system in step S5: When sufficient historical fault case data is available, a data-driven approach can be used to replace the analytic hierarchy process (AHP) in determining weights. A logistic regression classification model is established using fault occurrence as a label, with normalized regression coefficients serving as weights. This model is continuously optimized as fault data accumulates, avoiding the subjectivity of expert scoring. The minimum remaining deformation space threshold can also be set differently based on the importance level of structural regions: 35% for the planetary carrier region directly affecting transmission accuracy, 25% for the mating surface region primarily responsible for sealing, and 30% for the bearing seat bore region affecting bearing life, making the margin allocation more reasonable. Each alternative in this embodiment can independently or in combination replace the corresponding steps in Embodiments 1 and 2 without affecting the overall logical chain integrity of the method.

[0059] Example 4 This embodiment, after completing the stiffness compliance determination using the method in Embodiment 1, further performs closed-loop correction on the multibody dynamics topology model, constituting a complete implementation of the scheme described in claim 9. The application scenario is that a gearbox prototype has been installed on a test bench, and measured data can be obtained for model calibration.

[0060] After completing steps S1 to S5, the reliability of the stiffness determination results output by the simulation model depends on the estimation accuracy of the parameters in the initial modeling. Parameters such as the actual material damping of the housing, the actual clearance of the bearings after installation, and the actual contact stiffness of the gear pair deviate from their design nominal values. If calibration is not based on measured data, the accumulated deviations may lead to inaccurate determinations. Therefore, this embodiment introduces an experimental-simulation closed-loop correction mechanism.

[0061] The test bench was used to acquire measured data, including fiber optic strain data and laser vibration data. Fiber optic strain data was collected using fiber Bragg grating strain sensors mounted on key areas of the test chamber, including the circumferential direction of bearing housing holes, the root of the planetary carrier cantilever, and the mating flange. A total of 48 channels were used, with a sampling frequency of 5kHz and a measurement range of [missing information]. resolution Laser vibration data was collected using a scanning laser Doppler vibrometer (Polytec PSV-500) to assess the vibration velocity field of the outer wall surface of the enclosure. At least 200 scanning points were used, covering a frequency range of 0 to 3000 Hz, with a velocity resolution better than [previous value missing]. The locations of the two types of sensors correspond one-to-one with the observation points output by the simulation model.

[0062] The measured data and the simulation model's calculated output under the same operating conditions were time-domain aligned and their characteristics compared. Time-domain alignment applied the same timestamp alignment to the simulation output based on the key phase pulse of the test rotational speed signal. Characteristic comparison included: the time-domain waveform, frequency-domain amplitude spectrum, and amplitude and phase of key meshing orders at each strain measurement point; the effective value of vibration velocity and the spatial distribution of operational deformation at each vibration measurement point. By comparing and identifying the distribution characteristics of the deviation between simulation and measurement, the parameter domain requiring correction was determined.

[0063] A system parameter dimensionality reduction matrix is ​​constructed to map the parameters to be corrected to a lower-dimensional space. In this embodiment, sensitivity analysis is used to screen out 12 parameters that have a significant impact on the output, including the elastic modulus correction factor of each region of the housing, the material damping coefficient, the contact stiffness coefficient of each bearing, and the clearance correction amount. Then, principal component analysis is used to reduce the dimensionality. The significance of dimensionality reduction is to reduce the dimensionality of the inverse optimization search space, improve numerical stability, and avoid the risk of overfitting in high-dimensional space.

[0064] A surrogate mapping model between system parameters and simulation outputs is established based on a dimension-reduced matrix. This embodiment employs the Kriging surrogate model: several sample points are generated in the dimension-reduced parameter space through Latin hypercube sampling. For each sample point, complete steps S1 to S4 are performed, and the output is recorded. The Kriging model is trained using parameters as input and the response as output. The covariance function is selected as anisotropic Matérn 5 / 2, and the hyperparameters are determined through maximum likelihood estimation. Leave-one-out cross-validation is used to ensure the model accuracy is acceptable when the normalized root mean square error is below 5%.

[0065] The objective is to update the stiffness and damping parameters in the simulation model through inverse iterative updates, using the deviation between measured data and calculated output as the target. The problem is transformed into an optimization problem: the weighted least squares deviation between the Kriging predicted output and the measured data is used as the objective function, and the dimensionality-reduced low-dimensional coordinates are used as the design variables. Global optimization is performed within the physically feasible region. A differential evolution algorithm is used for the search. If the optimized objective function is still above the convergence threshold, sample points are added near the optimum, the Kriging model is updated, and optimization is repeated, with a maximum of 5 iterations.

[0066] Fault tolerance: When the objective function has multiple local minima, perform no less than 3 independent optimization runs from different initial populations and take the lowest result of the objective function; when the Kriging model's prediction variance exceeds 30% of the mean near the parameter boundary, supplement the sampling in that region to improve local accuracy and avoid the optimization being guided to a high uncertainty region.

[0067] The corrected stiffness and damping parameters are mapped back to the original 12-dimensional space through an inverse transformation of the dimension reduction matrix. These are then substituted back into the reduced-order flexible body matrix and the bearing nonlinear contact mechanics model. For the corrected values ​​of the housing material parameters, the material properties of the corresponding region in the finite element model are updated, and the dimension reduction is re-executed to generate a new flexible body matrix, replacing the old one. For the corrected values ​​of the bearing parameters, the corresponding parameters in the bearing model are directly updated. After the correction is complete, steps S3 to S5 are run again. At this point, the simulation output is the corrected result calibrated with measured data.

[0068] In some possible implementation scenarios, test bench testing includes not only steady-state conditions but also order-based acceleration and transient impact conditions. In this case, the measured data for closed-loop correction should cover multiple conditions, and the Kriging model training samples should also cover the corresponding dimensions of the condition variables to ensure that the corrected parameters are robust across the entire condition domain.

[0069] The closed-loop correction mechanism established in this embodiment upgrades the simulation model from open-loop prediction to closed-loop prediction after experimental calibration, eliminating the uncertainty of the stiffness determination conclusion due to the initial modeling parameter estimation deviation. The Kriging proxy model compresses the number of full model simulation calls to the level of sample points, and the parameter dimensionality reduction matrix ensures efficient convergence of inverse optimization in low-dimensional space. The synergy of these two mechanisms enables the gearbox model to complete parameter calibration from the design state to the physical state with reasonable computational resource consumption, providing a reliable model foundation for full-condition stiffness evaluation based on the calibrated model.

[0070] The threshold, weight, coefficient and other settings involved in the above embodiments can be set by default according to the present invention, or can be set by those skilled in the art based on experience and actual needs.

[0071] Those skilled in the art will understand that embodiments of the present invention can be provided as methods, systems, or computer program products. Therefore, the present invention can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, the present invention can take the form of a computer program product embodied on one or more computer-usable storage media (including, but not limited to, disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0072] This invention is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, generate instructions for implementing the flowchart illustrations and / or block diagrams. Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0073] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including instruction means, which are implemented in a process Figure 1 One or more processes and / or boxes Figure 1 The function specified in one or more boxes.

[0074] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0075] The embodiments of the present invention have been described above with reference to the accompanying drawings. However, the present invention is not limited to the specific embodiments described above. The specific embodiments described above are merely illustrative and not restrictive. Those skilled in the art can make many other forms under the guidance of the present invention without departing from the spirit and scope of the claims. All of these forms are within the protection scope of the present invention.

Claims

1. A method for detecting the stiffness of gearbox structural components based on multibody dynamics, characterized in that: Includes the following steps: S1. Establish a multibody dynamics topology model, extract and matrix condense the boundary degrees of freedom at the connection between the box and the transmission chain in the box finite element model, generate a reduced-order flexible body matrix that retains the principal modal characteristics, and import the reduced-order flexible body matrix into the multibody dynamics topology model; S2. Configure load fingerprint parameters including impact coefficient and time-varying attenuation factor, and generate system input boundary load that varies with operating conditions based on the load fingerprint parameters; S3. Under the drive of the system input boundary load, the gear is discretized into multiple micro-element slices along the tooth width direction based on the slice theory and the time-varying contact stiffness of each micro-element is calculated in real time. The bearing nonlinear contact mechanics model is called, and the time-varying contact stiffness, the bearing nonlinear contact mechanics model, and the reduced-order flexible body matrix are solved in the same time step to obtain the dynamic contact force and support reaction force sequence in the full time domain. S4. Spatial interpolation is performed on the dynamic contact force and support reaction force sequence through the cross-domain mapping interface, and the interpolation results are mapped to the mesh nodes of the box finite element model to perform transient stress and strain calculation. S5. Based on the results of the transient stress-strain calculation, extract the dynamic deformation of the structural component, calculate the difference between the dynamic deformation and the preset design allowable deformation limit to obtain the remaining deformation space, and determine the stiffness compliance in conjunction with the dynamic stiffness influence factor evaluation system.

2. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 1, characterized in that: The generation of the reduced-order flexible body matrix that preserves the principal modal characteristics includes: A finite element mesh is generated from the three-dimensional geometric model of the box to form a finite element model of the box, and the boundary degree-of-freedom nodes at each connection between the box and the transmission chain are identified and calibrated. Based on the boundary degree-of-freedom nodes, a matrix reduction method is used to perform matrix reduction. During the reduction process, the main modal characteristics of the box are preserved and the degrees of freedom corresponding to high-frequency redundant nodes are removed. The reduced-order flexible body matrix is ​​subjected to modal deviation verification. After confirming that the deviation between the main modal frequency and the full-order finite element model is within the preset allowable range, the reduced-order flexible body matrix is ​​output.

3. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 1, characterized in that: The process of generating time-varying system input boundary loads based on load fingerprint parameters includes: The system obtains the rotational speed time series and torque time series from the wind load parameters, performs time-domain alignment and filtering on the rotational speed time series and torque time series, configures the impact coefficient to characterize the amplification characteristics of transient load abrupt change amplitude and the time-varying attenuation factor to characterize the attenuation evolution characteristics of the load in the time domain, and modulates the reference load sequence based on the load fingerprint parameters composed of the impact coefficient and the time-varying attenuation factor to generate the system input boundary load.

4. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 1, characterized in that: Solving the equations within the same time step includes: S31. Simultaneously perform the calculation of gear time-varying contact stiffness based on slice theory and the calculation of bearing dynamic support force based on bearing nonlinear contact mechanics model within each time step. S32. The calculation results of the time-varying contact stiffness of the gear and the calculation results of the dynamic support force of the bearing are used as internal force inputs, and the dynamic equations are solved together with the reduced-order flexible body matrix within the current time step. S33. Run the joint dynamic equation solution process step by step until the entire simulation time domain is covered, and obtain the dynamic contact force sequence and support reaction force sequence.

5. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 4, characterized in that: The calculation process for the time-varying contact stiffness of the gear includes: The tooth width, module, and pressure angle parameters of the gear pair are obtained, and the gear is divided into a preset number of micro-element slices at equal intervals along the tooth width direction. For each micro-element slice, the instantaneous contact stiffness of the micro-element slice is calculated based on the meshing phase corresponding to the current time step. The instantaneous contact stiffness of each micro-element slice is integrated and superimposed along the tooth width direction to obtain the time-varying contact stiffness of the gear pair at the current time step. The number of micro-element slices is determined according to the tooth width size and calculation accuracy requirements, and the instantaneous contact stiffness is updated in real time with the change of meshing phase.

6. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 4, characterized in that: The specific construction process of the bearing nonlinear contact mechanics model includes: The number of rolling elements, contact angle, and clearance parameters of the bearing are obtained, and the relative displacement of the inner and outer rings of the bearing at the current time step is obtained. Based on the number of rolling elements, contact angle, clearance parameters, and relative displacement, the contact deformation of each rolling element is calculated. The contact force of each rolling element is calculated by applying a nonlinear stiffness-deformation relationship to the contact deformation of each rolling element. The contact forces of each rolling element are vector-superimposed to obtain the dynamic support resultant force and resultant torque of the bearing. When the contact deformation of the rolling element is negative, the corresponding contact force is set to zero.

7. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 1, characterized in that: The process of spatial interpolating and mapping the mesh nodes of the box finite element model through the cross-domain mapping interface includes: S41. Using the spatial coordinates of each force point in the multibody dynamics topology model and the spatial coordinates of each mesh node in the box finite element model as input, establish the mapping topology relationship between the force point and the mesh node and the corresponding spatial interpolation weight coefficient based on spatial proximity search, and output the mapping topology table. S42. Using the mapping topology table and the full-time dynamic contact force and support reaction force sequence as input, perform spatial interpolation calculation on the force vector of each time step according to the spatial interpolation weight coefficient, and distribute the interpolated force vector to the mesh node corresponding to the box finite element model, and output the distributed equivalent node force load sequence of each time step. S43. Using the distributed equivalent nodal force load sequence and the dynamic contact force and support reaction force sequence before mapping as input, perform numerical comparison and verification on the resultant force and resultant moment before and after mapping respectively. When the deviation exceeds the preset conservation tolerance, return to step S42 to adjust the spatial interpolation weight coefficient. When the deviation is within the preset conservation tolerance, output the distributed equivalent nodal force load sequence that has passed the conservation verification. S44. The distributed equivalent nodal force load sequence that has passed the conservation check is used as the load boundary condition input to the box finite element model, and transient stress and strain calculation is performed to output the transient stress field and strain field of the box structure in the whole time domain.

8. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 1, characterized in that: The process of extracting the dynamic deformation of structural components and calculating the remaining deformation space to determine stiffness compliance includes: Extract the maximum dynamic deformation of the structural component in the entire time domain from the results of the transient stress-strain calculation; obtain the preset design allowable deformation upper limit, calculate the difference between the maximum dynamic deformation and the design allowable deformation upper limit to obtain the remaining deformation space; determine whether the remaining deformation space meets the preset minimum remaining deformation space threshold, and determine that the stiffness is non-compliant when the remaining deformation space is less than the minimum remaining deformation space threshold.

9. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 1, characterized in that: After determining the stiffness compliance, the process also includes a closed-loop correction step for the multibody dynamics topology model, specifically including: Acquire measured data from the physical test bench, including fiber optic strain data and laser vibration data; The measured data and the calculation output of the simulation model under the same working conditions are time-domain aligned and feature-compared. A system parameter dimensionality reduction matrix is ​​constructed to map the parameters to be corrected to a low-dimensional space. Based on the dimensionality reduction matrix, a proxy mapping relationship model between the system parameters and the simulation output is established. Taking the deviation between the measured data and the calculated output as the target, the stiffness parameters and damping parameters in the simulation model are updated through inverse iteration. Finally, the corrected stiffness and damping parameters are substituted back into the reduced-order flexible body matrix and the bearing nonlinear contact mechanics model.

10. The method for detecting the stiffness of gearbox structural components based on multibody dynamics according to claim 1, characterized in that: The stiffness compliance determination process also includes: Transient stiffness variables are extracted from the full-time dynamic contact force sequence and support reaction force sequence. The vibration amplitude obtained from the transient stress-strain calculation and the fatigue life reduction coefficient of the bearing and gear are obtained. A dynamic stiffness influence factor evaluation matrix is ​​constructed. The transient stiffness variables, the vibration amplitude and the fatigue life reduction coefficient are used as input elements to perform fusion calculation to generate a quantitative ranking index. Based on the quantitative ranking index, the stiffness weakness of each structural region is ranked, and the comprehensive stiffness compliance judgment result is jointly output with the judgment result of the remaining deformation space.

Images