Fan modal identification method based on harmonic comb positioning and symplectic geometry harmonic elimination

Abstract

The application discloses a fan modal identification method based on harmonic comb positioning and symplectic geometry harmonic elimination, and belongs to the technical field of structural health monitoring of wind generating sets, and comprises the following steps: collecting multi-direction vibration responses, constructing fusion spectrum and background normalization; positioning the fundamental frequency and harmonic frequency through harmonic comb scoring; performing symplectic matrix decomposition on the responses to obtain time domain components, determining harmonic dominant components according to the consistency of the main frequency and the harmonic frequency and the predominance of the harmonic energy, and eliminating the harmonic dominant components, and reconstructing the responses after the elimination of the harmonics; constructing a two-dimensional parameter domain of system order-covariance time lag, performing SSI-COV on the responses after the elimination of the harmonics, recursively subdividing a high divergence area according to the pole divergence degree of adjacent parameter units, fusing stable poles of effective parameter combinations, and then extracting physical modal families through DBSCAN clustering. The application does not need to preset harmonic frequencies or target frequency bands, can suppress harmonic interference, and improves the adaptability and reliability of modal identification under variable speed conditions.

Description

Technical Field

[0001] This invention relates to the fields of wind turbine generator structural health monitoring, operational modal analysis, vibration signal processing, and wind turbine tower dynamic characteristic identification. Specifically, it relates to an automatic wind turbine structural modal identification method based on harmonic comb scoring, symplectic geometric harmonic elimination, pole divergence feedback, SSI parameter domain recursive subdivision, and DBSCAN clustering. Background Technology

[0002] With the rapid development of offshore wind power, a large number of offshore wind turbines have entered long-term service. Due to complex environmental loads such as wind, waves, and currents, as well as the limitations of offshore operation and maintenance conditions, the long-term safety status of the towers has become a major concern. Structural health monitoring can characterize the overall dynamic characteristics of the structure by analyzing measured vibration responses and identifying parameters such as modal frequencies, damping ratios, and mode shapes, providing a basis for condition assessment, numerical model correction, and operation and maintenance. Since external excitations are often difficult to measure directly and completely during actual operation, operational modal analysis can identify modal parameters based solely on the structural response under unknown input conditions, and is therefore widely used for monitoring the dynamic characteristics of offshore wind turbines. Among existing operational modal analysis methods, Covariance-Driven Stochastic Subspace Identification (SSI-COV) is widely used to directly extract wind turbine structural modal parameters from measured vibration responses. SSI-COV typically assumes that the unmeasured input can be approximated as a broadband random excitation. However, the wind turbine operating response also contains deterministic or quasi-deterministic narrowband components related to rotor rotation and blade passage. The operating harmonics manifest in the power spectrum as sharp, narrow-band peaks, peaks that drift with rotational speed, or peak clusters with a certain width, and can generate pole branches that repeat with the system order in the stability diagram. When the operating harmonic frequency is close to the structural mode frequency, it is difficult to distinguish between harmonic-induced poles and physical mode poles based solely on pole stability, which can easily lead to mode misidentification.

[0003] Existing methods typically address this issue by using fixed band-stop filtering, manually specifying harmonic frequencies, directly deleting PSD peaks, or setting a target frequency band in the SSI results. However, these methods have the following shortcomings: Fixed band-stop filters cannot adapt to harmonic drift and bandwidth expansion caused by variable fan speeds; PSD peak detection alone is insufficient to distinguish between the modal peaks of the wind turbine tower structure and the harmonic peaks of wind turbine operation; Setting a target frequency band at the SSI output would weaken the effectiveness of the method in automatically identifying structural modes.

[0004] Therefore, a method is needed that can adaptively identify operating harmonics directly from multi-directional measured responses, selectively remove harmonic components at the response component level, and adaptively explore the system order-covariance time-delay parameter domain based on SSI pole results feedback. This method should avoid pre-specifying operating harmonic frequencies, target mode bands, and fixed parameter sampling sequences. Instead, it should recursively refine local parameter regions with significant discrepancies in pole results, and automatically fuse stable poles after pole results convergence. Then, it should employ a density-based spatial clustering of applications with noise (DBSCAN) algorithm for clustering, thus balancing the reliability of mode identification with computational specificity. Summary of the Invention

[0005] The purpose of this invention is to provide a wind turbine structural mode identification method based on harmonic comb localization and symplectic geometric harmonic elimination, aiming to solve the problem of rotational harmonics (1P, 3P, etc.) being misidentified as tower structural modes by SSI-COV during wind turbine operation. Through methods such as harmonic comb scoring for automatic localization of harmonic families, symplectic geometric harmonic elimination for automatic elimination of harmonic components, recursive subdivision of the SSI parameter domain with pole bifurcation feedback, and DBSCAN clustering, this method can automatically and reliably identify the true structural modes of the tower.

[0006] To achieve the above objectives, this invention provides a wind turbine mode identification method based on harmonic comb localization and symplectic geometric harmonic rejection, comprising the following steps: S1. Under normal operating conditions of the wind turbine generator set, vibration sensors are installed on the tower, nacelle, foundation or support structure to collect vibration responses in all directions, synchronize and preprocess the vibration responses in all directions, record the actual sampling frequency, and obtain multi-channel response signals. S2. Based on the power spectral density of the vibration response in each direction, construct a multi-directional fusion spectrum to retain significant narrowband peaks in any direction, and construct a background normalization spectrum based on the slowly varying spectral background of the multi-directional fusion spectrum. S3. Within the candidate fundamental frequency range, construct a harmonic comb score based on the degree of peak enhancement in the background normalized spectrum according to the integer harmonic positions of the candidate fundamental frequency, and determine the operating fundamental frequency; using the integer harmonics of the operating fundamental frequency as the nominal harmonic frequencies, perform local peak correction and effective harmonic order screening in the neighborhood of each nominal harmonic frequency to obtain the set of operating harmonic frequencies and the influence band corresponding to each operating harmonic; S4. Perform symplectic matrix decomposition on each channel of the multi-channel response signal, and reconstruct multiple time-domain response components by anti-angle averaging. S5. Calculate the dominant frequency, spectral energy within the operating harmonic influence band, and spectral energy outside the operating harmonic influence band for each time-domain response component; when the dominant frequency of a time-domain response component satisfies the frequency consistency condition with the set of operating harmonic frequencies, and the spectral energy within its operating harmonic influence band satisfies the energy dominance condition, the time-domain response component is determined to be the harmonic dominant component. S6. Remove the dominant harmonic components and superimpose the time-domain response components that were not determined to be dominant harmonic components to obtain the harmonic-removed responses in each direction. S7. Construct a two-dimensional parameter domain consisting of system order and covariance time delay. Perform covariance-driven random subspace identification on the response after harmonic removal at each parameter combination in the two-dimensional parameter domain. Adaptively adjust the parameters according to the degree of divergence of the pole set obtained by adjacent parameter combinations, and fuse the stable poles corresponding to the effective parameter combinations to form a global stable pole set. S8. The density-based noise spatial clustering algorithm (DBSCAN) is used to cluster the globally stable pole set. The pole clusters that meet the requirements of neighborhood density and minimum number of samples are determined as physical mode families. The modal information of the wind turbine support structure is output according to each physical mode family.

[0007] Preferably, S1 specifically includes: S11. The vibration response in each direction includes at least two orthogonal vibration responses, which are referred to as the X-direction response and the Y-direction response; the vibration responses in each direction are synchronized in time so that the sampling points of each channel correspond to each other in time. S12. Identify and correct abnormal sampling points that exceed the preset amplitude threshold; S13. Remove, interpolate, or segment missing sampling points; S14. When it is necessary to reduce the sampling frequency, the sampling frequency after downsampling the vibration response in each direction meets the Nyquist sampling requirements of the target analysis frequency band, and a multi-channel response signal is obtained.

[0008] Preferably, S2 specifically includes: S21. Calculate the power spectral density of the vibration response in the two orthogonal directions, and denot them as follows: and ; S22. Construct a multidirectional fusion spectrum for retaining significant narrowband peaks in any direction according to the following formula, wherein the multidirectional fusion spectrum is a bidirectional maximum power spectrum: S23. The bidirectional maximum power spectrum is processed using moving median, moving quantile, morphological filtering, or smoothing fitting to obtain the slowly varying spectral background. S24. Construct the background normalized spectrum according to the following formula: ; S25, among which, The slowly varying spectral background, To prevent positive numbers from being divided by zero.

[0009] Preferably, S3 specifically includes: S31, Regarding candidate fundamental frequencies According to the upper limit of the analysis frequency Determine the available harmonic order ; S32. Construct a harmonic comb score that normalizes the available harmonic orders according to the following formula: ; S33, among which For the first Weights of first harmonics; selection of factors that make... The candidate fundamental frequency that yields the maximum value is used as the operating fundamental frequency. .

[0010] S34, with For the first The nominal harmonic frequency, in The local maximum value of the background normalized spectrum is searched to obtain the corrected harmonic frequency. ;in, For spectral resolution, This is a multiple of the resolution. S35. Determine whether the modified harmonic frequency forms a local spectral peak, and calculate its integer harmonic deviation relative to the corresponding nominal harmonic frequency. S36. Only retain the harmonic orders that simultaneously satisfy the local peak condition and the integer multiple frequency deviation condition to form the set of operating harmonic frequencies; S37. Based on the spectral resolution, local peak width, or a combination of both, adaptively determine the influence band corresponding to each operating harmonic.

[0011] Preferably, S4 includes: S41. Based on delay embedding, each single-channel vibration response is constructed as a trajectory matrix; S42. Construct a symplectic correlation matrix based on the trajectory matrix; S43. Perform complex Schur decomposition or equivalent symmetric orthogonal decomposition on the symmetric correlation matrix to obtain a symmetric orthogonal basis, and project the symmetric orthogonal basis onto the trajectory matrix to obtain multiple projection matrices; S44. Perform anti-angle averaging on each projection matrix to obtain multiple time-domain response components that satisfy the original response reconstructibility relationship.

[0012] Preferably, S5 specifically includes: S51, Extract the first Direction The main frequency of each time domain response component The set of the main frequency and the operating harmonic frequencies is calculated according to the following formula. Distance between: S52, when When the time-domain response component satisfies the frequency consistency condition, it is determined that the time-domain response component satisfies the frequency consistency condition. This is a multiple of the frequency resolution. S53. The sum of the spectral energy within each operating harmonic influence band is denoted as... The spectral energy within the remaining analysis frequency band is denoted as ;when When the time-domain response component satisfies the energy dominance condition, it is determined that the time-domain response component satisfies the energy dominance condition. The energy discrimination coefficient is not less than zero; S54. The time-domain response component that simultaneously satisfies the frequency consistency condition and the energy dominance condition is determined as the harmonic dominant component.

[0013] Preferably, S6 includes: S61, Read the sequence of... Direction Each time-domain response component ; S62, if Simultaneously satisfying the frequency consistency condition and the energy dominance condition, Set to zero or multiply by an attenuation coefficient less than 1; otherwise, retain. S63. Superimpose all the processed time-domain response components to obtain the first... The response after directional harmonic removal.

[0014] Preferably, S7 specifically includes: S71. Determine the range of system order values ​​and the range of covariance time delay values ​​to form a two-dimensional parameter domain; S72. Divide the two-dimensional parameter domain into multiple initial parameter units according to a preset initial division scale, and select at least one representative parameter combination that satisfies the integer constraint from each initial parameter unit. S73. Perform covariance-driven random subspace identification (SSI-COV) on the response after harmonic removal using each representative parameter combination to obtain the candidate pole set corresponding to each parameter unit. S74. Match the candidate poles corresponding to adjacent parameter elements of the shared boundary or vertex, and calculate the pole divergence based on at least two of the following: frequency difference, damping ratio difference, mode shape difference, and proportion of unmatched poles. S75. The parameter unit with a pole divergence degree greater than the preset divergence threshold is determined as the parameter unit to be subdivided. The parameter unit to be subdivided is subdivided along the system order direction, the covariance time delay direction or both directions, and a representative parameter combination is reselected in the obtained sub-parameter units. S76. Repeat S73 to S75 until the maximum pole divergence of the parameter unit is lower than the preset divergence threshold, the change of the physical modality family cluster center obtained by the clustering algorithm in two consecutive rounds is lower than the preset convergence threshold, the proportion of newly added stable poles in the current round is lower than the preset proportion threshold, or at least one of the following is reached: the preset maximum number of recursion layers or the preset maximum number of recognitions. S77. Collect all candidate poles corresponding to all effective parameter combinations during the recursive subdivision process, screen and merge stable poles based on frequency stability, damping ratio stability and mode shape consistency, form a global stable pole set, and generate a stability graph.

[0015] Preferably, S8 specifically includes: S81. The modal frequencies and damping ratios of the global stable poles are standardized to form a feature vector for the clustering algorithm, or the extended distance is constructed by combining modal similarity. S82. Execute the clustering algorithm according to the preset neighborhood radius and minimum number of samples, and determine the poles that are not assigned to any pole cluster as noise or pseudo-modes, and determine the remaining pole clusters as physical mode families. S83. Using the mean, weighted mean, median or robust estimation of the poles within the cluster, obtain the representative modal frequencies and damping ratios of each physical mode family, and finally output the structural modal parameters of the wind turbine tower. The structural modal parameters include at least the modal frequencies, damping ratios and modal directions determined by the mode shapes of the responses in each direction.

[0016] Therefore, compared with the prior art, the wind turbine structural modal identification method based on harmonic comb localization, symplectic matrix decomposition, pole bifurcation feedback and DBSCAN clustering in this invention has the following beneficial effects: (1) By fusing the power spectrum information of vibration responses in two orthogonal directions and using background normalization to highlight local narrowband peaks, the interference of broadband random response on the identification of operating harmonics is reduced. On this basis, the candidate fundamental frequency and its integer harmonic positions are matched as a whole by using harmonic comb scoring. Combined with local peak correction and effective harmonic order screening, the operating harmonics under different operating conditions can be adaptively identified without pre-specifying the operating harmonic frequency, reducing the possibility of misjudging the structural modal peak as the operating harmonic by single peak detection.

[0017] (2) The original vibration response is decomposed into multiple time-domain components with different spectral characteristics by symplectic matrix decomposition. Based on the frequency consistency between the dominant frequency of each time-domain component and the identified operating harmonics, as well as the degree of energy dominance in the harmonic neighborhood, the dominant harmonic component is jointly determined. Only the components that simultaneously meet the conditions of frequency consistency and harmonic energy dominance are eliminated. The remaining components are used to reconstruct the response after harmonic elimination, thereby selectively weakening the operating harmonic components and reducing the impact on the modal information of the adjacent structure and the broadband random response components.

[0018] (3) By constructing a two-dimensional parameter domain consisting of system order and covariance time delay, and recursively subdividing local parameter units with large pole divergence according to the degree of divergence between pole sets obtained from adjacent parameter units, the SSI-COV calculation is concentrated in the region where the modality recognition result is more sensitive to parameter changes; the parameter region with consistent pole results is no longer repeatedly encrypted, and the subdivision stops after the pole results meet the preset convergence conditions, thereby reducing the dependence of the modality recognition result on a single system order or covariance time delay setting, and reducing the repeated calculation of low information gain parameter regions.

[0019] (4) By integrating all effective parameter combinations obtained in the recursive subdivision process, and using the density-based noise application space clustering algorithm DBSCAN to automatically cluster the stable poles, it is possible to distinguish between high-density physical modal pole clusters and discrete noise poles without pre-setting the number of modal categories. The modal frequency, damping ratio and mode shape of the wind turbine support structure are automatically obtained from each physical modal pole cluster, so that the stable poles induced by the operating harmonics are effectively suppressed, and the automation, stability and reliability of structural modal identification under the wind turbine operating state are improved.

[0020] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description

[0021] Figure 1 This is a flowchart illustrating the wind turbine mode identification method based on harmonic comb localization and symplectic geometric harmonic elimination of the present invention. Figure 2 These are comparison images of PSDs before and after harmonic removal in an embodiment of the present invention: (a) is the original image, and (b) is the image after harmonic removal. Figure 3 These are comparison diagrams of stable pole identification in embodiments of the present invention: (a) is the result diagram using the conventional method, and (b) is the result diagram using the method of the present invention. Figure 4 The following are the modal parameter evolution trend diagrams of the embodiments of the present invention: (a) is the time domain variation diagram of the identified modal frequencies, and (b) is the time domain variation diagram of the identified damping ratios. Detailed Implementation

[0022] The technical solution of the present invention will be further described below with reference to the accompanying drawings and embodiments.

[0023] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0024] Example like Figure 1 As shown, this embodiment provides a wind turbine mode identification method based on harmonic comb localization and symplectic geometric harmonic rejection, including the following steps: S1. Under normal operating conditions of the wind turbine generator set, vibration sensors are installed on the tower, nacelle, foundation or support structure to collect vibration responses in all directions, synchronize and preprocess the vibration responses in all directions, record the actual sampling frequency, and obtain multi-channel response signals. Specifically, it includes: S11. The vibration response in each direction includes at least two orthogonal vibration responses, which are called the X-direction response and the Y-direction response. The vibration responses in each direction are synchronized in time so that the sampling points of each channel correspond to each other in time. S12. Identify and correct abnormal sampling points that exceed the preset amplitude threshold; S13. Remove, interpolate, or segment missing sampling points; S14. When it is necessary to reduce the sampling frequency, the sampling frequency after downsampling the vibration response in each direction meets the Nyquist sampling requirements of the target analysis frequency band, and a multi-channel response signal is obtained.

[0025] S2. Based on the power spectral density of the vibration response in each direction, construct a multi-directional fusion spectrum to retain significant narrowband peaks in any direction, and construct a background normalization spectrum based on the slowly varying spectral background of the multi-directional fusion spectrum. Specifically, it includes: S21. Calculate the power spectral density of the vibration response in the two orthogonal directions, and denot them as follows: and ; S22. Construct a multidirectional fusion spectrum to retain significant narrowband peaks in any direction according to the following formula. The multidirectional fusion spectrum is a bidirectional maximum power spectrum: ; S23. The bidirectional maximum power spectrum is processed by moving median, moving quantile, morphological filtering or smoothing fitting to obtain the slowly varying spectral background. S24. Construct the background normalized spectrum according to the following formula: ; S25, among which, For the slowly varying spectral background, To prevent positive numbers from being divided by zero.

[0026] S3. Within the candidate fundamental frequency range, construct a harmonic comb score based on the enhancement degree of the spectral peaks in the background normalized spectrum at the integer harmonic positions of the candidate fundamental frequency to determine the operating fundamental frequency; using the integer harmonics of the operating fundamental frequency as the nominal harmonic frequencies, perform local peak correction and effective harmonic order screening in the neighborhood of each nominal harmonic frequency to obtain the set of operating harmonic frequencies and the influence band corresponding to each operating harmonic. The core concept of S3 lies in the fact that traditional methods typically search for individual peaks in the power spectrum in isolation, making it difficult to distinguish between structural mode peaks and wind turbine rotation harmonic peaks. This invention utilizes the physical characteristic that wind turbine rotation harmonics have a strict octave relationship, and constructs a weighted harmonic comb scoring function to perform overall matching and comprehensive scoring of the enhancement degree of spectral peaks at all octave positions of the candidate fundamental frequency.

[0027] This "harmonic comb" localization method is equivalent to using a comb with adjustable frequency spacing to match the signal-to-noise ratio spectrum. The spacing between the comb teeth corresponds to the candidate fundamental frequency, and the enhancement degree of the spectral peaks at each tooth location is weighted and summed. When the comb spacing matches the true harmonic fundamental frequency, all harmonic peaks are simultaneously "combed," and the scoring function reaches its maximum value. This method does not search for a single isolated PSD peak, but rather matches a family of harmonics with harmonic relationships, thus reliably distinguishing rotational harmonics from structural modes and avoiding the misjudgment problem of traditional single-peak detection methods.

[0028] Specifically, it includes: S31, Regarding candidate fundamental frequencies According to the upper limit of the analysis frequency Determine the available harmonic order ; S32. Construct a harmonic comb score that normalizes the available harmonic orders according to the following formula: ; S33, among which For the first Weights of first harmonics; selection of factors that make... The candidate fundamental frequency that yields the maximum value is used as the operating fundamental frequency. .

[0029] S34, with For the first The nominal harmonic frequency, in The local maxima of the background normalized spectrum are searched to obtain the corrected harmonic frequencies. ;in, For spectral resolution, This is a multiple of the resolution. S35. Determine whether the corrected harmonic frequency forms a local spectral peak, and calculate its integer harmonic deviation relative to the corresponding nominal harmonic frequency. S36. Only retain the harmonic orders that simultaneously satisfy the local peak condition and the integer multiple frequency deviation condition to form the set of operating harmonic frequencies; S37. Based on the spectral resolution, local peak width, or a combination of both, adaptively determine the influence band corresponding to each operating harmonic.

[0030] S4. Perform symplectic matrix decomposition on each channel of the multi-channel response signal, and reconstruct multiple time-domain response components by anti-angle averaging. The specific process is as follows: S41. Based on delay embedding, each single-channel vibration response is constructed as a trajectory matrix; S42. Construct the symplectic correlation matrix based on the trajectory matrix; S43. Perform complex Schur decomposition or equivalent symmetric orthogonal decomposition on the symmetric correlation matrix to obtain a symmetric orthogonal basis, and project the symmetric orthogonal basis onto the trajectory matrix to obtain multiple projection matrices. S44. Perform anti-angle averaging on each projection matrix to obtain multiple time-domain response components that satisfy the original response reconstructibility relationship.

[0031] S5. Calculate the dominant frequency, spectral energy within the operating harmonic influence band, and spectral energy outside the operating harmonic influence band for each time-domain response component. When the dominant frequency of a time-domain response component satisfies the frequency consistency condition with the set of operating harmonic frequencies, and the spectral energy within its operating harmonic influence band satisfies the energy dominance condition, the time-domain response component is determined to be the harmonic dominant component. Specifically, it includes: S51, Extract the first Direction The main frequency of each time domain response component And calculate the set of the main frequency and the operating harmonic frequencies according to the following formula. Distance between: ; S52, when When the time-domain response components satisfy the frequency consistency condition, This is a multiple of the frequency resolution. S53. The sum of the spectral energy within each operating harmonic influence band is denoted as... The spectral energy within the remaining analysis frequency band is denoted as ;when When the time-domain response component satisfies the energy dominance condition, The energy discrimination coefficient is not less than zero; S54. The time-domain response component that simultaneously satisfies the frequency uniformity condition and the energy dominance condition is determined as the harmonic dominant component.

[0032] S6. Remove the dominant harmonic components and superimpose the time-domain response components that were not identified as dominant harmonic components to obtain the harmonic-removed responses in each direction. Specifically, it includes: S61, Read the sequence of... Direction Each time-domain response component ; S62, if Simultaneously satisfying both the frequency uniformity condition and the energy dominance condition, Set to zero or multiply by an attenuation coefficient less than 1; otherwise, retain. S63. Superimpose all the processed time-domain response components to obtain the first... The response after directional harmonic removal.

[0033] Figure 2 This is a comparison image of the PSD before and after harmonic removal in an embodiment of the present invention, wherein: Figure 2 (a) in the figure shows the power spectral density diagram corresponding to the original vibration response. As can be seen from the figure, the original response exhibits obvious narrow-band peaks at the 1P, 2P, and 3P harmonic frequencies, with some operating harmonic peaks being particularly prominent and distributed within the analysis frequency band along with the broadband vibration components of the structure. When the operating harmonic frequency is close to the tower structure's modal frequency, it is difficult to accurately distinguish between the operating harmonics and the structural modes based solely on local peaks in the power spectral density, which can easily interfere with subsequent modal identification.

[0034] Figure 2 Figure (b) shows the power spectral density diagram of the response after harmonic removal, obtained by harmonic comb localization and symplectic matrix decomposition of the original vibration response according to the present invention. As can be seen from the figure, the narrowband spectral peaks at the operating harmonic frequencies such as 1P, 2P, and 3P have been significantly weakened, while the main spectral features near the actual structural modes of the tower are still preserved. This indicates that the present invention can specifically identify and remove dominant harmonic components while reducing the impact on structural modal components and broadband random response components.

[0035] pass Figure 2The before-and-after comparison shows that the amplitude of the spectral peak corresponding to the running harmonics is significantly reduced after processing by the present invention, while the frequency band of the structural mode does not shift significantly. This proves that the present invention can retain the effective modal information in the structural dynamic response while suppressing the interference of running harmonics, providing more reliable input data for subsequent modal parameter identification.

[0036] S7. Construct a two-dimensional parameter domain consisting of the system order and covariance time delay. Perform covariance-driven random subspace identification on the response after harmonic removal at each parameter combination in the two-dimensional parameter domain. Adaptively adjust the parameters according to the degree of divergence of the pole set obtained by adjacent parameter combinations, and fuse the stable poles corresponding to the effective parameter combinations to form a global stable pole set. Specifically, it includes: S71. Determine the range of system order values ​​and the range of covariance time delay values ​​to form a two-dimensional parameter domain; S72. Divide the two-dimensional parameter domain into multiple initial parameter units according to the preset initial division scale, and select at least one representative parameter combination that satisfies the integer constraint from each initial parameter unit. S73. Using each representative parameter combination, perform covariance-driven random subspace identification on the response after harmonic removal to obtain the candidate pole set corresponding to each parameter unit. S74. Match the candidate poles corresponding to adjacent parameter elements of the shared boundary or vertex, and calculate the pole divergence based on at least two of the following: frequency difference, damping ratio difference, mode shape difference, and proportion of unmatched poles. S75. The parameter unit with a pole divergence greater than the preset divergence threshold is identified as the parameter unit to be subdivided. The parameter unit to be subdivided is subdivided along the system order direction, the covariance time delay direction, or both directions. Representative parameter combinations are then selected from the sub-parameter units. S76. Repeat S73 to S75 until the maximum pole divergence of the parameter unit is lower than the preset divergence threshold, the change of the physical modality family cluster center obtained by DBSCAN for two consecutive rounds is lower than the preset convergence threshold, the proportion of newly added stable poles in the current round is lower than the preset proportion threshold, or at least one of the following is reached: the preset maximum number of recursive layers or the preset maximum number of recognitions. S77. Collect all candidate poles corresponding to all effective parameter combinations during the recursive subdivision process, screen and merge stable poles based on frequency stability, damping ratio stability and mode shape consistency, form a global stable pole set, and generate a stability graph.

[0037] Figure 3 This is a comparison chart of stable pole identification in an embodiment of the present invention, wherein: Figure 3(a) in the figure shows the stability plot obtained by directly inputting the original response into the SSI-COV program. As can be seen from the figure, there are numerous scattered unstable poles at the 1P, 2P, and other harmonic frequencies. These harmonic-induced pseudo-poles form a dense distribution on the stability plot, severely interfering with the identification of the true structural modes and making it difficult for engineers to directly determine which poles represent the true structural modes of the tower from the original stability plot.

[0038] Figure 3 Figure (b) shows the stability diagram obtained after processing according to this invention. As can be seen from the figure, a large number of unstable poles at the harmonic locations have been effectively eliminated, and the number of interfering poles in the stability diagram has been significantly reduced. At the same time, the vertical column of dark-colored stable poles corresponding to the actual structural modes of the tower (located around 0.3Hz in the figure) is completely preserved without attenuation or shift due to harmonic removal. Figure 3 The left-right comparison clearly demonstrates that while eliminating false harmonic poles, the present invention does not damage the true structural modal information of the tower, and the clarity and engineering interpretability of the stability diagram are significantly improved.

[0039] S8. The density-based noise spatial clustering algorithm (DBSCAN) is used to cluster the globally stable pole set. The pole clusters that meet the requirements of neighborhood density and minimum number of samples are determined as physical mode families. The modal information of the wind turbine support structure is output according to each physical mode family.

[0040] Specifically, it includes: S81. Standardize the modal frequencies and damping ratios of the global stable poles to form feature vectors for clustering algorithms, or further combine them with modal similarity to construct extended distances; S82. Perform DBSCAN based on the preset neighborhood radius and minimum number of samples. Determine the poles that are not assigned to any pole cluster as noise or pseudo-modes. Determine the remaining pole clusters as physical mode families. S83. Using the mean, weighted mean, median or robust estimation of the poles within the cluster, obtain the representative modal frequencies and damping ratios of each physical mode family, and finally output the structural modal parameters of the wind turbine tower. The structural modal parameters include at least the modal frequencies, damping ratios and modal directions determined by the mode shapes of the responses in each direction.

[0041] Figure 4 The modal parameter evolution trend diagram of this invention shows the stability change of the structural modal parameters identified by S8 over time during a wind turbine operation monitoring cycle of approximately two consecutive days, wherein: Figure 4Figure (a) shows the variation of the identified modal frequencies over time. As can be seen from the figure, the modal frequencies of each structural order can be continuously identified during the continuous monitoring period. The identification results are mainly concentrated within their respective narrow frequency ranges, without large-scale jumps or obvious misidentification caused by operating harmonic frequency drift. The modal frequencies obtained in different time windows show good consistency, indicating that the present invention can stably track the modal frequencies of the wind turbine support structure under changes in wind turbine operating status and external excitation.

[0042] Figure 4 Figure (b) shows the variation of the identified damping ratio over time. As can be seen from the figure, the damping ratio of each mode fluctuates to a certain extent within the continuous time window, and its variation amplitude is greater than that of the mode frequency. This is because the damping ratio is more sensitive to environmental excitation, wind turbine operating conditions, response noise, and pole dispersion. Nevertheless, the damping ratios of each mode remain within a relatively concentrated and reasonable range, without a large number of abnormal discrete values, indicating that this invention can effectively reduce the influence of operating harmonic-induced poles on the damping ratio identification results.

[0043] pass Figure 4 The results of modal frequency and damping ratio variations show that the present invention can continuously obtain stable structural modal parameters under continuous operation and monitoring conditions. The modal frequencies exhibit good temporal continuity and concentration, and the damping ratio variations remain within an interpretable range, demonstrating that the present invention has good parameter identification stability and long-term monitoring applicability, providing a reliable basis for tracking the dynamic characteristics and evaluating the operational status of wind turbine support structures.

[0044] Therefore, this invention employs the aforementioned wind turbine modal identification method based on harmonic comb localization and symplectic geometric harmonic removal, which can automatically locate 1P, 3P, and other rotating harmonics and adaptively determine their bandwidth to adapt to wind turbine variable speed operating conditions. By decomposing the vibration response into multiple time-domain components through symplectic matrix decomposition, and identifying and removing dominant harmonic components based on the consistency between the dominant frequency of the component and the frequency of the operating harmonics, as well as the degree of energy dominance in the harmonic neighborhood, the invention weakens the interference of operating harmonics while preserving the structural modal-related response. Furthermore, a parameter domain recursive subdivision method based on pole divergence feedback is used to adaptively explore the system order and covariance time delay parameter domains, fusing stable poles obtained under different effective parameter combinations, and using DBSCAN to cluster the stable poles, thereby suppressing the interference of harmonic-induced poles and discrete noise poles. Thus, this invention can automatically obtain the modal frequencies, damping ratios, and mode shapes of the wind turbine support structure without pre-setting a target identification frequency band, improving the automation, stability, reliability, and engineering applicability of structural modal identification under wind turbine operating conditions.

[0045] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit them. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the technical solutions of the present invention, and these modifications or equivalent substitutions cannot cause the modified technical solutions to deviate from the spirit and scope of the technical solutions of the present invention.

Claims

1. A wind turbine mode identification method based on harmonic comb localization and symplectic geometric harmonic rejection, characterized in that, Includes the following steps: S1. Under normal operating conditions of the wind turbine generator set, vibration sensors are installed on the tower, nacelle, foundation or support structure to collect vibration responses in all directions, synchronize and preprocess the vibration responses in all directions, record the actual sampling frequency, and obtain multi-channel response signals. S2. Based on the power spectral density of the vibration response in each direction, construct a multi-directional fusion spectrum to retain significant narrowband peaks in any direction, and construct a background normalization spectrum based on the slowly varying spectral background of the multi-directional fusion spectrum. S3. Within the candidate fundamental frequency range, construct a harmonic comb score based on the degree of peak enhancement in the background normalized spectrum according to the integer harmonic positions of the candidate fundamental frequency, and determine the operating fundamental frequency. Using the integer multiples of the operating fundamental frequency as nominal harmonic frequencies, local peak correction and effective harmonic order screening are performed in the neighborhood of each nominal harmonic frequency to obtain the set of operating harmonic frequencies and the influence band corresponding to each operating harmonic. S4. Perform symplectic matrix decomposition on each channel of the multi-channel response signal, and reconstruct multiple time-domain response components by anti-angle averaging. S5. Calculate the dominant frequency, spectral energy within the operating harmonic influence band, and spectral energy outside the operating harmonic influence band for each time domain response component. When the dominant frequency of a time-domain response component satisfies the frequency consistency condition with the set of operating harmonic frequencies, and the spectral energy within its operating harmonic influence band satisfies the energy dominance condition, the time-domain response component is determined to be a harmonic dominant component. S6. Remove the dominant harmonic components and superimpose the time-domain response components that were not determined to be dominant harmonic components to obtain the harmonic-removed responses in each direction. S7. Construct a two-dimensional parameter domain consisting of system order and covariance time delay. Perform covariance-driven random subspace identification on the response after harmonic removal at each parameter combination in the two-dimensional parameter domain. Adaptively adjust the parameters according to the degree of divergence of the pole set obtained by adjacent parameter combinations, and fuse the stable poles corresponding to the effective parameter combinations to form a global stable pole set. S8. A density-based noise spatial clustering algorithm is used to cluster the globally stable pole set. The pole clusters that meet the requirements of neighborhood density and minimum number of samples are determined as physical mode families. The modal information of the wind turbine support structure is output according to each physical mode family.

2. The method according to claim 1, characterized in that, S1 specifically includes: S11. The vibration response in each direction includes at least two orthogonal vibration responses, which are referred to as the X-direction response and the Y-direction response; the vibration responses in each direction are synchronized in time so that the sampling points of each channel correspond to each other in time. S12. Identify and correct abnormal sampling points that exceed the preset amplitude threshold; S13. Remove, interpolate, or segment missing sampling points; S14. When it is necessary to reduce the sampling frequency, the sampling frequency after downsampling the vibration response in each direction meets the Nyquist sampling requirements of the target analysis frequency band, and a multi-channel response signal is obtained.

3. The method according to claim 2, characterized in that, S2 specifically includes: S21. Calculate the power spectral density of the vibration response in the two orthogonal directions, and denot them as follows: and ; S22. Construct a multidirectional fusion spectrum for retaining significant narrowband peaks in any direction according to the following formula, wherein the multidirectional fusion spectrum is a bidirectional maximum power spectrum: S23. The bidirectional maximum power spectrum is processed using moving median, moving quantile, morphological filtering, or smoothing fitting to obtain the slowly varying spectral background. S24. Construct the background normalized spectrum according to the following formula: ； S25, among which, The slowly varying spectral background, To prevent positive numbers from being divided by zero.

4. The method according to claim 3, characterized in that, S3 specifically includes: S31, Regarding candidate fundamental frequencies According to the upper limit of the analysis frequency Determine the available harmonic order ; S32. Construct a harmonic comb score that normalizes the available harmonic orders according to the following formula: ； S33, among which... For the first Weights of first harmonics; selection of weights that make the first harmonic harmonics... The candidate fundamental frequency that yields the maximum value is used as the operating fundamental frequency. ; S34, with For the first The nominal harmonic frequency, in The local maximum value of the background normalized spectrum is searched to obtain the corrected harmonic frequency. ;in, For spectral resolution, This is a multiple of the resolution. S35. Determine whether the modified harmonic frequency forms a local spectral peak, and calculate its integer harmonic deviation relative to the corresponding nominal harmonic frequency; S36. Only retain the harmonic orders that simultaneously satisfy the local peak condition and the integer multiple frequency deviation condition to form the set of operating harmonic frequencies; S37. Based on the spectral resolution, local peak width, or a combination of both, adaptively determine the influence band corresponding to each operating harmonic.

5. The method according to claim 4, characterized in that, S4 includes: S41. Based on delay embedding, each single-channel vibration response is constructed as a trajectory matrix; S42. Construct a symplectic correlation matrix based on the trajectory matrix; S43. Perform complex Schur decomposition or equivalent symmetric orthogonal decomposition on the symmetric correlation matrix to obtain a symmetric orthogonal basis, and project the symmetric orthogonal basis onto the trajectory matrix to obtain multiple projection matrices; S44. Perform anti-angle averaging on each projection matrix to obtain multiple time-domain response components that satisfy the original response reconstructibility relationship.

6. The method according to claim 5, characterized in that, S5 specifically includes: S51, Extract the first Direction The main frequency of each time domain response component The set of the main frequency and the operating harmonic frequencies is calculated according to the following formula. Distance between: S52, when When the time-domain response component is satisfied, it is determined that the frequency consistency condition is met. This is a multiple of the frequency resolution. S53. The sum of the spectral energy within each operating harmonic influence band is denoted as... The spectral energy within the remaining analysis frequency band is denoted as ;when When the time-domain response component satisfies the energy dominance condition, it is determined that the time-domain response component satisfies the energy dominance condition. The energy discrimination coefficient is not less than zero; S54. The time-domain response component that simultaneously satisfies the frequency consistency condition and the energy dominance condition is determined as the harmonic dominant component.

7. The method according to claim 5, characterized in that, S6 includes: S61, Read the sequence of... Direction Each time-domain response component ; S62, if Simultaneously satisfying the frequency consistency condition and the energy dominance condition, Set to zero or multiply by an attenuation coefficient less than 1; otherwise, retain. S63. Superimpose all the processed time-domain response components to obtain the first... The response after directional harmonic removal.

8. The method according to claim 7, characterized in that, Specifically, S7 includes: S71. Determine the range of system order values ​​and the range of covariance time delay values ​​to form a two-dimensional parameter domain; S72. Divide the two-dimensional parameter domain into multiple initial parameter units according to the preset initial division scale, and select at least one representative parameter combination that satisfies the integer constraint from each initial parameter unit. S73. Using each representative parameter combination, perform covariance-driven random subspace identification on the response after harmonic removal to obtain the candidate pole set corresponding to each parameter unit. S74. Match the candidate poles corresponding to adjacent parameter elements of the shared boundary or vertex, and calculate the pole divergence based on at least two of the following: frequency difference, damping ratio difference, mode shape difference, and proportion of unmatched poles. S75. The parameter unit with a pole divergence degree greater than the preset divergence threshold is determined as the parameter unit to be subdivided. The parameter unit to be subdivided is subdivided along the system order direction, the covariance time delay direction or both directions, and a representative parameter combination is reselected in the obtained sub-parameter units. S76. Repeat S73 to S75 until the maximum pole divergence of the parameter unit is lower than the preset divergence threshold, the change of the physical modality family cluster center obtained by the clustering algorithm in two consecutive rounds is lower than the preset convergence threshold, the proportion of newly added stable poles in the current round is lower than the preset proportion threshold, or at least one of the following is reached: the preset maximum number of recursion layers or the preset maximum number of recognitions. S77. Collect all candidate poles corresponding to all effective parameter combinations during the recursive subdivision process, screen and merge stable poles based on frequency stability, damping ratio stability and mode shape consistency, form a global stable pole set, and generate a stability graph.

9. The method according to claim 1, characterized in that, S8 specifically includes: S81. The modal frequencies and damping ratios of the global stable poles are standardized to form a feature vector for the clustering algorithm, or the extended distance is constructed by combining modal similarity. S82. Execute the clustering algorithm according to the preset neighborhood radius and minimum number of samples, and determine the poles that are not assigned to any pole cluster as noise or pseudo-modes, and determine the remaining pole clusters as physical mode families. S83. Using the mean, weighted mean, median or robust estimation of the poles within the cluster, obtain the representative modal frequencies and damping ratios of each physical mode family, and finally output the structural modal parameters of the wind turbine tower. The structural modal parameters include at least the modal frequencies, damping ratios and modal directions determined by the mode shapes of the responses in each direction.

Images