A method for detecting the vibration amplitude of a large-inertia flywheel in a high-speed rotating state

Abstract

This invention provides a method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation, belonging to the field of precision machining technology. The invention uses a thermal expansion finite element model to dynamically compensate the original displacement signal to separate the rapidly changing vibration components. It employs a kurtosis spectrum to automatically locate the optimal demodulation frequency band and extracts the envelope spectrum using Hilbert transform to construct a modulation amplitude matrix. A key phase sensor is used to resample the vibration signal to the angular domain, and Hilbert-Huang transform is combined to extract the instantaneous amplitude envelope in the angular domain. Sparse Bayesian learning is used to decompose the multidimensional vibration feature vector into unbalanced excitation components, bearing impact components, and oil film whirl components, and uncertainty is quantified. Finally, a two-layer Nash equilibrium closed-loop model is used to jointly optimize the bearing preload and thermal compensation correction, and an adjustment command is output. This method solves the problem of vibration amplitude extraction distortion in large-inertia flywheels under conditions of wide speed range and large temperature rise interference.

Description

Technical Field

[0001] This invention belongs to the field of precision machining technology, and specifically relates to a method for detecting the vibration amplitude of a high-inertia flywheel under high-speed rotation. Background Technology

[0002] Large-inertia flywheel energy storage devices are widely used in power grid peak shaving, rail transit braking energy recovery, and industrial uninterruptible power supplies. One of their core performance indicators is the vibration amplitude of the flywheel shaft system within a wide speed range of 0–10000 r / min. Traditional vibration detection methods typically combine piezoelectric accelerometers with time-domain peak statistics to perform single-point measurements of steady-state vibration under fixed speed conditions. For bearing fault feature extraction, engineering commonly employs spectral analysis combining fixed-bandwidth bandpass filtering and fast Fourier transform. For variable-speed, non-stationary conditions, some methods introduce short-time Fourier transform or wavelet transform for time-frequency analysis.

[0003] However, the above methods have significant limitations in practical applications of large-inertia flywheels. The flywheel bearing housing experiences a temperature rise of 60–90°C under high-speed rotation. The slow-varying displacement component caused by thermal expansion overlaps with the fast-varying vibration component in frequency range, making it impossible for traditional fixed-band filtering to effectively separate them. Simultaneously, the non-stationary acceleration and deceleration process over a wide speed range causes the frequencies of each order of vibration components to continuously drift with the rotational speed. Spectral analysis methods based on fixed-frequency coordinates produce frequency ambiguity effects, resulting in severely underestimated amplitude estimates. Furthermore, the vibration responses of three excitation mechanisms—unbalanced excitation, bearing impact, and oil film whirl—are superimposed in the frequency band, and existing methods lack a systematic approach to decouple and separate these three mechanisms.

[0004] In existing technologies, due to the combined effects of the aliasing of slow-varying thermal expansion components and fast-varying vibration components in the high-temperature-rise environment of large-inertia flywheel bearing housings, frequency ambiguity under non-stationary operating conditions over a wide speed range, and the aliasing of vibration responses from multiple excitation mechanisms, existing vibration detection methods cannot accurately extract vibration amplitudes under these complex interference conditions. In other words, existing technologies suffer from the technical problem of distorted vibration amplitude extraction for large-inertia flywheels under high-temperature-rise interference conditions over a wide speed range. Summary of the Invention

[0005] In view of this, the present invention provides a method for detecting the vibration amplitude of a large inertia flywheel under high-speed rotation. This method can solve the technical problem that existing vibration amplitude detection methods cannot effectively separate the slow-changing thermal expansion component and the fast-changing vibration component under wide speed range conditions where the large inertia flywheel rotates at high speed and the bearing housing has large temperature rise interference, resulting in distorted vibration amplitude extraction.

[0006] This invention is implemented as follows: This invention provides a method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation, comprising the following steps:

[0007] A piezoelectric sensor array is mounted on the flywheel bearing housing, and a high-precision thermocouple array is synchronously arranged at each measuring point on the flywheel bearing housing to collect the original displacement signal of the flywheel shaft in the range of 0 to 10000 r / min and the temperature field data of the flywheel bearing housing.

[0008] Based on the temperature field data of the flywheel bearing housing, the thermal expansion finite element model is called to calculate the static component of thermal deformation at each measuring point in real time, and the original displacement signal is dynamically thermally compensated and corrected. The slow thermal expansion component and the fast vibration component are separated, and the vibration displacement signal is output.

[0009] A multi-resolution bandpass filter bank is constructed for the vibration displacement signal. The kurtosis value is calculated for each sub-band signal and the kurtosis spectrum is plotted. The optimal band with the largest kurtosis value is automatically selected for bandpass filtering. The Hilbert transform is performed on the filtered signal to extract the instantaneous amplitude envelope. The fast Fourier transform is performed on the instantaneous amplitude envelope to obtain the envelope spectrum. The harmonic family amplitudes are searched in the envelope spectrum with the bearing fault characteristic frequencies as the fundamental frequencies and accumulated. The modulation amplitude matrix is ​​then output.

[0010] Using the key phase pulse signal acquired by the key phase sensor as a reference, the vibration displacement signal is resampled into an equal-angle domain signal. In the angle domain, a fast Fourier transform is performed on the equal-angle domain signal to eliminate the frequency ambiguity caused by the rotational speed change. Then, the result of the fast Fourier transform in the angle domain is combined with the Hilbert-Huang transform to extract the instantaneous amplitude envelope in the angle domain. The instantaneous amplitude envelope in the angle domain and the elements in the modulation amplitude matrix are arranged in order to form a multidimensional vibration feature vector.

[0011] The multidimensional vibration feature vector is sparsely decomposed under an overcomplete atom dictionary. A sparse prior is established using a Bayesian inference framework. The mean and variance of the posterior distribution of the sparse coefficients are iteratively updated using the expectation-maximization algorithm. Only atoms with significantly non-zero posterior variance are retained. The multidimensional vibration feature vector is decomposed into unbalanced excitation components, bearing impact components, and oil film whirl components. The amplitude envelopes of the unbalanced excitation components, bearing impact components, and oil film whirl components are extracted respectively. At the same time, the uncertainty quantification results of the amplitude envelopes of each component are output.

[0012] An upper-level optimization model is constructed with the goal of minimizing the amplitude of the unbalanced excitation component, and a lower-level optimization model is constructed with the goal of minimizing the bearing preload adjustment error. The upper-level optimization model and the lower-level optimization model are coupled through coupling terms to solve the Nash equilibrium and obtain the optimal bearing preload and the optimal thermal compensation correction amount. When any of the amplitude envelopes of the unbalanced excitation component, the bearing impact component, or the oil film whirl component exceeds a preset threshold, an adjustment command is output based on the optimal bearing preload and the optimal thermal compensation correction amount to correct the flywheel bearing preload state and the dynamic thermal compensation correction amount, thus completing the closed-loop detection of vibration amplitude.

[0013] Specifically, the dynamic thermal compensation correction refers to the real-time subtraction of the slow thermal expansion component in the original displacement signal by subtracting the static thermal deformation component of each measuring point calculated by the thermal expansion finite element model at the signal processing layer, so that the output vibration displacement signal retains only the fast dynamic vibration component.

[0014] Specifically, the thermal expansion finite element model is a calculation model that uses the geometric structure and material thermal expansion coefficient of the flywheel bearing housing as structural inputs and the temperature field data of the flywheel bearing housing as boundary conditions. It numerically solves the static components of thermal deformation generated at each measuring point under temperature rise conditions using the finite element method.

[0015] Specifically, the multi-resolution bandpass filter bank is a filtering structure composed of a group of bandpass filters with different center frequencies and bandwidths, used to decompose the vibration displacement signal into multiple sub-frequency band signals so as to calculate the kurtosis value for each frequency band.

[0016] Specifically, the kurtosis spectrum is a two-dimensional spectrum plotted with the sub-band center frequency and sub-band bandwidth as two-dimensional coordinate axes and the kurtosis value of each sub-band signal as the vertical height. The optimal frequency band is the sub-band with the largest kurtosis value.

[0017] Specifically, the modulation amplitude matrix is ​​a matrix composed of the characteristic frequencies of each bearing fault and the sum of their harmonic amplitudes, with rows corresponding to the bearing fault type and columns corresponding to the amplitudes of each order of harmonics.

[0018] Specifically, the key phase sensor is a sensor installed on the flywheel shaft that outputs a key phase pulse signal every time the flywheel shaft rotates through a fixed angle. The key phase pulse signal is used to mark the rotation angle position of the flywheel shaft so as to achieve uniform resampling of the vibration displacement signal to the equal angle domain signal.

[0019] Specifically, the equal-angle domain signal is a signal sequence obtained by uniformly resampling the vibration displacement signal with the flywheel shaft rotation angle as the horizontal axis, so that the frequency of each order vibration component in the angle domain remains fixed under variable speed conditions, thereby eliminating the frequency ambiguity effect caused by speed fluctuation.

[0020] Specifically, the Hilbert-Huang transform is a method that first performs empirical mode decomposition on the signal to obtain several intrinsic mode functions, and then performs Hilbert transform on each intrinsic mode function to obtain the instantaneous amplitude and instantaneous frequency. It is suitable for the instantaneous feature extraction of non-stationary vibration signals.

[0021] The overcomplete atomic dictionary is specifically a dictionary matrix composed of structural mode vectors and modulation waveforms, with the number of atoms exceeding the dimension of the signal. The structural mode vectors represent the natural vibration modes of each order of the flywheel shaft system, and the modulation waveforms represent the periodic modulation components of bearing impact and oil film vortex.

[0022] Specifically, the sparse prior is a prior probability distribution assumption applied to each atomic coefficient within the Bayesian inference framework. Each atomic coefficient follows a zero-mean Laplace distribution, such that only a few atomic coefficients in the decomposition result are significantly non-zero.

[0023] Specifically, the expectation-maximization algorithm is an optimization algorithm that iteratively solves the mean and variance of the posterior distribution of sparse coefficients by alternately executing the expectation step and the maximization step. The expectation step calculates the expectation of the posterior distribution of sparse coefficients, and the maximization step updates the hyperparameters to maximize the log-marginal likelihood.

[0024] Specifically, the uncertainty quantification result is the variance of the posterior distribution of each atomic coefficient after the expectation-maximization algorithm converges. It is used to evaluate the reliability of the extraction results of the amplitude envelope of the unbalanced excitation component, the amplitude envelope of the bearing impact component, and the amplitude envelope of the oil film whirl component.

[0025] The objective function of the upper-level optimization model is: The input includes the amplitude envelope of the unbalanced excitation component. With flywheel shaft speed The output is the target value of the upper layer optimization. The constraint condition is that the envelope of the unbalanced excitation component amplitude does not exceed a preset threshold.

[0026] Wherein, the objective function of the lower-level optimization model is Inputs include bearing preload. Temperature of flywheel bearing housing The constraint condition is the bearing preload. Within the allowable adjustment range and with the flywheel bearing housing temperature rise not exceeding 90℃, the upper-level optimization model and the lower-level optimization model are coupled through a coupling term. The interaction is used to solve the joint configuration of the optimal bearing preload and the optimal thermal compensation correction using Nash equilibrium.

[0027] This invention organically integrates five steps: real-time dynamic thermal compensation of thermal expansion finite element model, adaptive resonance demodulation of kurtosis spectrum, joint envelope extraction of angle domain resampling and Hilbert-Huang transform, separation of multi-excitation components by sparse Bayesian learning, and closed-loop control of double-layer Nash equilibrium. This constructs a complete vibration amplitude detection framework for large inertia flywheels operating under wide speed range and large temperature rise conditions.

[0028] This invention eliminates slowly varying interference caused by thermal expansion at the signal source, ensuring that the signals processed in subsequent signal processing stages only contain real dynamic vibration components; kurtosis spectrum automatically locates the frequency band where impact energy is most concentrated, eliminating reliance on manual experience; angular domain resampling transforms frequency drift caused by variable rotation speed into a fixed order, eliminating frequency ambiguity at its source; sparse Bayesian learning utilizes the sparsity differences in the atomic dictionary of three excitation mechanisms—unbalanced excitation, bearing impact, and oil film eddy—to achieve physical decoupling and provides a quantification of the uncertainty of the amplitude envelope of each component, ensuring the reliability of the extraction results; a two-layer Nash equilibrium model coordinates the joint optimization of bearing preload and thermal compensation correction in a closed-loop manner, enabling the system to continuously maintain the optimal detection state under dynamic operating conditions.

[0029] In summary, this invention solves the technical problem of distortion in vibration amplitude extraction of large inertia flywheels under conditions of large temperature rise interference over a wide speed range, as mentioned in the background art. Attached Figure Description

[0030] Figure 1 This is a flowchart of the method of the present invention.

[0031] Figure 2 This is a schematic diagram of a large inertia flywheel device.

[0032] Figure 3 This is a graph showing the distribution of the temperature field of the flywheel bearing housing as a function of rotational speed.

[0033] Figure 4 This is the kurtosis spectrum of the vibration displacement signal.

[0034] Figure 5 This is a graph showing the change of the instantaneous amplitude envelope in the angle domain with the rotation angle.

[0035] Figure 6 The graph shows the evolution of the amplitude envelope of each component of the sparse Bayesian decomposition as a function of the acceleration process. Detailed Implementation

[0036] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below.

[0037] like Figure 1 The diagram shows a flowchart of a method for detecting the vibration amplitude of a high-inertia flywheel under high-speed rotation, provided by this invention. The method includes the following steps:

[0038] S01. Install the piezoelectric sensor array on the flywheel bearing housing, and synchronously arrange the high-precision thermocouple array at each measuring point of the flywheel bearing housing to collect the original displacement signal of the flywheel shaft in the range of 0 to 10000 r / min and the temperature field data of the flywheel bearing housing.

[0039] S02. Based on the temperature field data of the flywheel bearing housing, the thermal expansion finite element model is called to calculate the static component of thermal deformation at each measuring point in real time, and the original displacement signal is dynamically thermally compensated and corrected. The slow thermal expansion component and the fast vibration component are separated, and the vibration displacement signal is output.

[0040] S03. Construct a multi-resolution bandpass filter bank for the vibration displacement signal, calculate the kurtosis value for each sub-band signal and draw the kurtosis spectrum, automatically select the optimal band with the largest kurtosis value for bandpass filtering, perform Hilbert transform on the filtered signal to extract the instantaneous amplitude envelope, perform fast Fourier transform on the instantaneous amplitude envelope to obtain the envelope spectrum, search for harmonic family amplitudes in the envelope spectrum with the bearing fault characteristic frequency as the fundamental frequency and accumulate them, and output the modulation amplitude matrix.

[0041] S04. Using the key phase pulse signal acquired by the key phase sensor as a reference, the vibration displacement signal is resampled into an equal-angle domain signal. In the angle domain, a fast Fourier transform is performed on the equal-angle domain signal to eliminate the frequency ambiguity caused by the rotational speed change. Then, the result of the fast Fourier transform in the angle domain is combined with the Hilbert-Huang transform to extract the instantaneous amplitude envelope in the angle domain. The instantaneous amplitude envelope in the angle domain and each element in the modulation amplitude matrix are arranged in order to form a multidimensional vibration feature vector.

[0042] S05. The multidimensional vibration feature vector is sparsely decomposed under an overcomplete atom dictionary. A sparse prior is established using a Bayesian inference framework. The mean and variance of the posterior distribution of the sparse coefficients are iteratively updated using the expectation-maximization algorithm. Only atoms with significantly non-zero posterior variance are retained. The multidimensional vibration feature vector is decomposed into unbalanced excitation components, bearing impact components, and oil film whirl components. The amplitude envelopes of the unbalanced excitation components, bearing impact components, and oil film whirl components are extracted respectively. At the same time, the uncertainty quantification results of the amplitude envelopes of each component are output.

[0043] S06. An upper-level optimization model is constructed with the goal of minimizing the amplitude of the unbalanced excitation component, and a lower-level optimization model is constructed with the goal of minimizing the bearing preload adjustment error. The upper-level optimization model and the lower-level optimization model are coupled through coupling terms to solve the Nash equilibrium and obtain the optimal bearing preload and the optimal thermal compensation correction amount. When any one of the amplitude envelopes of the unbalanced excitation component, the bearing impact component, or the oil film whirl component exceeds the preset threshold, an adjustment command is output based on the optimal bearing preload and the optimal thermal compensation correction amount to correct the flywheel bearing preload state and the dynamic thermal compensation correction amount, thus completing the closed-loop detection of vibration amplitude.

[0044] Among them, dynamic thermal compensation correction refers to the real-time subtraction of the slow thermal expansion component in the original displacement signal by subtracting the static thermal deformation component of each measuring point calculated by the thermal expansion finite element model in the signal processing layer, so that the output vibration displacement signal retains only the fast dynamic vibration component.

[0045] Among them, the thermal expansion finite element model refers to a calculation model that uses the geometric structure and material thermal expansion coefficient of the flywheel bearing housing as structural input and the temperature field data of the flywheel bearing housing as boundary conditions, and numerically solves the static components of thermal deformation generated at each measuring point under temperature rise conditions using the finite element method.

[0046] Among them, the multi-resolution bandpass filter bank refers to a filtering structure composed of a group of bandpass filters with different center frequencies and bandwidths, which is used to decompose the vibration displacement signal into multiple sub-frequency band signals so that the kurtosis value can be calculated for different frequency bands.

[0047] The kurtosis spectrum refers to a two-dimensional spectrum plotted with the center frequency and bandwidth of the sub-band as the two-dimensional coordinate axes and the kurtosis value of each sub-band signal as the vertical height. It is used to automatically locate the optimal frequency band in the vibration displacement signal where the impact component is most concentrated. The optimal frequency band is the sub-band with the largest kurtosis value.

[0048] The modulation amplitude matrix is ​​a matrix composed of the sum of the characteristic frequencies of each bearing fault and their harmonic family amplitudes. The rows correspond to the bearing fault types and the columns correspond to the amplitudes of each order of harmonics. It is used to quantitatively characterize the vibration amplitude level corresponding to each bearing fault type.

[0049] Among them, the key phase sensor refers to a sensor installed on the flywheel shaft that outputs a key phase pulse signal every time the flywheel shaft rotates through a fixed angle. The key phase pulse signal is used to mark the rotation angle position of the flywheel shaft so as to achieve uniform resampling of the vibration displacement signal to the equal angle domain signal.

[0050] Among them, the equal-angle domain signal refers to the signal sequence obtained by uniformly resampling the vibration displacement signal with the rotation angle of the flywheel shaft as the horizontal axis, so that the frequency of each order vibration component in the angle domain remains fixed under variable speed conditions, thereby eliminating the frequency ambiguity effect caused by speed fluctuation.

[0051] Among them, the instantaneous amplitude envelope in the angle domain refers to the instantaneous amplitude envelope extracted by performing a fast Fourier transform on the equiangular domain signal in the angle domain and then combining it with the Hilbert-Huang transform, reflecting the instantaneous amplitude changes of each order of vibration of the flywheel shaft during the variable speed process.

[0052] The Hilbert-Huang transform is a method that first performs empirical mode decomposition on the signal to obtain several intrinsic mode functions, and then performs Hilbert transform on each intrinsic mode function to obtain the instantaneous amplitude and instantaneous frequency. It is suitable for the instantaneous feature extraction of non-stationary vibration signals.

[0053] Among them, the multidimensional vibration feature vector refers to the vector formed by arranging each element in the modulation amplitude matrix and the instantaneous amplitude envelope value in the angle domain in sequence, which comprehensively characterizes the vibration amplitude characteristics of the flywheel shaft under the current operating state.

[0054] Among them, the overcomplete atomic dictionary refers to a dictionary matrix composed of structural mode vectors and modulation waveforms, with the number of atoms exceeding the dimension of the signal. The structural mode vectors represent the natural vibration modes of each order of the flywheel shaft system, and the modulation waveforms represent the periodic modulation components of bearing impact and oil film vortex.

[0055] Sparse prior refers to the prior probability distribution assumptions applied to each atomic coefficient under the Bayesian inference framework. Each atomic coefficient follows a zero-mean Laplace distribution, so that only a few atomic coefficients in the decomposition result are significantly non-zero, thereby realizing the sparse representation of the multidimensional vibrational feature vector.

[0056] Among them, the expectation maximization algorithm is an optimization algorithm that iteratively solves the mean and variance of the posterior distribution of sparse coefficients by alternately executing the expectation step (calculating the expectation of the posterior distribution of sparse coefficients) and the maximization step (updating hyperparameters to maximize the log marginal likelihood).

[0057] Among them, the uncertainty quantification result refers to the variance of the posterior distribution of each atomic coefficient after the expectation-maximization algorithm converges. It is used to evaluate the reliability of the extraction results of the amplitude envelope of the unbalanced excitation component, the amplitude envelope of the bearing impact component, and the amplitude envelope of the oil film whirl component.

[0058] The objective function of the upper-level optimization model is expressed as follows: The objective function of the upper-level optimization model is used to minimize the dimensionless index of vibration amplitude caused by unbalanced excitation. The input includes the amplitude envelope of the unbalanced excitation component. (Unit: μm) and flywheel shaft speed (Unit: r / min), the output is the upper-level optimization objective value. ;in , , These are dimensionless weighting coefficients. The reference value for the amplitude envelope of the unbalanced excitation component (unit: μm). This is the reference value for rotational speed (unit: r / min). The baseline value of the upper-level objective function, coupling terms For bearing preload (Unit: N) and flywheel bearing housing temperature Dimensionless results of the coupling function (unit: °C) The reference value for the coupling function is given; the constraint condition is that the envelope of the unbalanced excitation component amplitude does not exceed a preset threshold.

[0059] The objective function of the lower-level optimization model is expressed as follows: The objective function of the lower-level optimization model is used to minimize the degree to which the bearing preload deviates from the temperature-adaptive optimal preload. The input includes the bearing preload. (Unit: N) and flywheel bearing housing temperature (Unit: °C), the output is the target value for the lower-level optimization. ;in To maintain the temperature of the flywheel bearing housing Temperature-adaptive optimal preload under given conditions (unit: N). This is the reference value for bearing preload (unit: N). The baseline value for the lower-level objective function. These are dimensionless weighting coefficients, and the coupling terms are... Shared with the upper-level optimization model; the constraint condition is the bearing preload. Within the allowable adjustment range and with the flywheel bearing housing temperature rise not exceeding 90℃; the upper-level optimization model and the lower-level optimization model are coupled through coupling terms. The interaction is used to solve the joint configuration of the optimal bearing preload and the optimal thermal compensation correction using Nash equilibrium.

[0060] The combined introduction of kurtosis spectrum adaptive resonance demodulation and sparse Bayesian learning vibration signal compression reconstruction enables the detection method to automatically locate the optimal demodulation frequency band without relying on human experience. At the same time, it separates the multidimensional vibration feature vector into unbalanced excitation components, bearing impact components, and oil film whirl components, and provides quantitative results of the uncertainty of the amplitude envelope of each component. This significantly improves the accuracy and robustness of vibration amplitude extraction of large inertia flywheels in a wide speed range and under the interference of large temperature rise in flywheel bearing housing. It also effectively suppresses amplitude deviation under non-stationary acceleration and deceleration conditions.

[0061] The high-inertia flywheel device comprises a flywheel housing, flywheel shaft, bearings, flywheel bearing housing, protective cover, input flange, and output flange. The flywheel housing serves as the overall base of the device, weighing approximately 1.4 tons. Its overall dimensions (including the flywheel housing base) are 1700mm × 900mm × 1906mm (length × width × height), providing the mounting foundation for the flywheel shaft system. The flywheel shaft is a single forging, with a center height of (1300±0.5) mm (including the flywheel housing base). The designed moment of inertia of the flywheel shaft is... The maximum design speed is 10,000 r / min. The flywheel shaft and flywheel disc constitute the inertia flywheel assembly, and the total weight of the inertia flywheel assembly (including the flywheel housing base) is approximately 4 tons. One set of bearings is installed on each end of the flywheel shaft, and the two sets of bearings are respectively installed in the corresponding flywheel bearing housings. The flywheel bearing housings are fixed on the flywheel housing, and the temperature rise range of the flywheel bearing housings is 60-90℃. A protective cover is installed on the top of the flywheel housing to ensure normal oil return of the flywheel shaft system when the flywheel shaft rotates at high speed. The input flange is located at one end of the flywheel shaft, and the output flange is located at the other end of the flywheel shaft. The input flange and the output flange are used to connect the input and output ends of the engine transmission system, respectively, to realize torque transmission.

[0062] The specific implementation of step S01 is as follows: At least three piezoelectric sensor measuring points are arranged axially and radially on both sides of the flywheel bearing housing of the large inertia flywheel device, forming a piezoelectric sensor array. This array collects the original displacement signal of the flywheel shaft within the speed range of 0–10000 r / min, with a sampling frequency of not less than 50 kHz to meet the Nyquist theorem's sampling requirements for high-frequency vibration components. Simultaneously, high-precision thermocouples are synchronously installed near each piezoelectric sensor measuring point, forming a thermocouple array. This array collects real-time temperature field data from each measuring point on the flywheel bearing housing. The thermocouple accuracy class is not less than 0.5, and the temperature sampling period is not greater than 1 second, ensuring that the temperature field data accurately reflects the thermal distribution state of the bearing housing within the 60–90℃ temperature rise range. The original displacement signal and temperature field data are stored in a data buffer after being timestamped by a synchronous acquisition card for subsequent processing.

[0063] The specific implementation of step S02 is as follows: using the three-dimensional geometric model of the flywheel bearing housing and the linear expansion coefficient of the material (the reference value for steel is approximately...) Using the structural parameters and the real-time temperature field data acquired by the thermocouple array as boundary conditions, the finite element method is employed to numerically solve the static thermal deformation components of each measuring point on the bearing housing under the current temperature rise, forming a thermal expansion finite element model. This model outputs the static thermal deformation displacement values ​​of each measuring point within each temperature sampling cycle. The dynamic thermal compensation correction is implemented by subtracting the static thermal deformation components of each measuring point from the original displacement signal of the corresponding measuring point in real time using a sampling-point interpolation method. Since the time constant of the thermal expansion component (on the order of tens of seconds) is much larger than the vibration period (on the order of milliseconds), the two are naturally separable on a time scale. The residual signal after subtraction is the vibration displacement signal, which retains only the fast-changing components of dynamic vibration, completely separating the slow-changing interference of thermal expansion.

[0064] The specific implementation of step S03 is as follows: A multi-resolution bandpass filter bank is constructed for the vibration displacement signal, and a fast spectral kurtosis algorithm is used to decompose the vibration displacement signal into several sub-band signals within the two-dimensional frequency band space formed by the center frequency and bandwidth. The center frequency coverage range is 0– ( (sampling frequency), bandwidth according to ( To decompose the signal into layers (reference values ​​range from 1 to 6), the signal is refined layer by layer. The kurtosis value is calculated for each sub-band signal. The kurtosis value is defined as the ratio of the fourth-order cumulant to the square of the variance; a larger value indicates a more concentrated impact component within the band. The center frequency, bandwidth, and corresponding kurtosis value of each sub-band are plotted as a kurtosis spectrum. The sub-band with the largest kurtosis value is automatically selected as the optimal band. The kurtosis threshold reference value is 6; signals below this value are considered to have insignificant impact characteristics. A Hilbert transform is performed on the filtered signal within the optimal band to extract the magnitude of the analytic signal as the instantaneous amplitude envelope. A Fast Fourier Transform is then performed on the instantaneous amplitude envelope to obtain the envelope spectrum. In the envelope spectrum, the bearing inner ring fault characteristic frequency, outer ring fault characteristic frequency, rolling element fault characteristic frequency, and cage fault characteristic frequency are used as the fundamental frequencies. Harmonic family amplitudes are searched and accumulated in higher-order directions to form a modulation amplitude matrix. The matrix rows correspond to the four types of bearing faults, and the columns correspond to the 1st to 5th order harmonic amplitudes.

[0065] The specific implementation of step S04 is as follows: A key phase sensor is installed on the flywheel shaft, outputting a fixed number of key phase pulses per revolution (reference value is 1024 pulses per revolution), with each pulse precisely corresponding to a fixed angular increment in the flywheel shaft rotation. Using the key phase pulse signal as the resampling trigger reference, the vibration displacement signal is resampled at equal angular intervals, transforming the time axis into an angular axis to obtain an equal angular domain signal. Since the horizontal axis of the equal angular domain signal is the rotation angle rather than time, the vibration components of each order do not drift with changes in rotational speed on the angular domain frequency coordinate (unit: order), thereby eliminating the frequency ambiguity effect caused by rotational speed fluctuations in the time-domain spectrum analysis. A fast Fourier transform is performed on the equal angular domain signal to obtain the order spectrum. Subsequently, empirical mode decomposition (EMD) is performed on the equiangular domain signal, decomposing it into several intrinsic mode functions (EMFs). A Hilbert transform is applied to each EMF to extract its instantaneous amplitude. These instantaneous amplitudes are then synthesized according to energy weights to obtain the angular domain instantaneous amplitude envelope, reflecting the dynamic changes in the instantaneous amplitude of each order of vibration of the flywheel shaft during variable speed. The sampled values ​​of the angular domain instantaneous amplitude envelope are arranged sequentially with all elements of the modulation amplitude matrix to form a multidimensional vibration feature vector. The vector dimension is equal to the sum of the number of sampling points in the angular domain instantaneous amplitude envelope and the total number of elements in the modulation amplitude matrix.

[0066] The specific implementation of step S05 is as follows: The overcomplete atomic dictionary consists of two types of atoms. The first type is the structural mode vector, which obtains the natural vibration modes of each order through modal analysis of the large inertia flywheel shaft system, with the first 10 modes taken as reference orders. The second type is the modulation waveform, which constructs periodic modulation atoms according to the typical modulation frequency and modulation depth parameters of three excitation mechanisms: unbalanced excitation, bearing impact, and oil film vortex. The total number of atoms in the dictionary matrix is ​​much larger than the dimension of the multidimensional vibration feature vector, thus constituting an overcomplete dictionary. Under the Bayesian inference framework, a zero-mean Laplace prior distribution is applied to each atomic coefficient to reflect the sparsity assumption of each component in the dictionary space. The expectation-maximization algorithm is used for iterative solution: the expectation step estimates the posterior mean and posterior variance of each atomic coefficient using the current hyperparameters; the maximization step updates the hyperparameters by maximizing the log marginal likelihood; the two steps are iterated alternately until the posterior mean converges, and the convergence criterion is that the relative change of the posterior mean between two adjacent iterations is less than 1. After convergence, atoms with significantly non-zero posterior variance correspond to real vibration components. The multidimensional vibration feature vector is reconstructed into unbalanced excitation components, bearing impact components, and oil film whirl components according to atom type, and the amplitude envelopes of the three types of components are extracted respectively. The posterior variance of each atom coefficient is directly output as the uncertainty quantification result of the corresponding component amplitude envelope. The upper limit of uncertainty is 0.05μm. If it exceeds this value, it indicates that the confidence of the component extraction result is insufficient.

[0067] The specific implementation of step S06 is as follows: the upper-level optimization model aims to minimize the amplitude of the unbalanced excitation component, and the objective function is... Envelope of unbalanced excitation component amplitude Flywheel shaft speed Coupling terms Weighted composition, weight coefficient , , The reference values ​​are 0.5, 0.3, and 0.2, and the constraints are as follows: The value should not exceed a preset threshold (reference value is 50μm). The lower-level optimization model aims to minimize the degree to which the bearing preload deviates from the temperature-adaptive optimal preload, with the objective function... It is composed of a weighted average of preload deviation and coupling terms, with weighting coefficients... The reference value is 0.4, and the constraints are that the bearing preload is within the allowable range and the flywheel bearing housing temperature rise does not exceed 90℃. The two-layer optimization model uses shared coupling terms. Interaction, coupling function This system reflects the combined effect of bearing preload and bearing housing temperature on bearing stiffness. A Nash equilibrium solution method is used to iterate through two layers of the model. Each layer solves for its own optimal solution under the condition that the decision variables of the other layer are fixed. This process is repeated alternately until the objective functions of both layers no longer decrease, yielding the optimal bearing preload and optimal thermal compensation correction. When the amplitude envelope of any component exceeds a preset threshold, the system outputs an adjustment command based on the Nash equilibrium result, driving the bearing preload actuator and the thermal compensation correction module to adjust to their optimal states, thus completing the closed-loop detection of vibration amplitude.

[0068] The key technologies of this invention will be described below.

[0069] The first key technical approach is real-time dynamic thermal compensation using a thermal expansion finite element model. Traditional methods for handling thermal drift typically rely on low-pass filters to cut off low-frequency components. However, under the temperature rise conditions of 60–90°C for large-inertia flywheels, the slow-varying displacement components caused by thermal expansion overlap with the low-frequency vibration components in frequency. Low-pass filtering simultaneously weakens the true low-frequency vibration information. This invention directly calculates and subtracts the static component of thermal deformation using physical modeling, preserving the complete fast-varying vibration components and ensuring signal quality for subsequent analysis from the signal source.

[0070] The second key technical approach is the joint envelope extraction of angle-domain resampling and Hilbert-Huang transform. Traditional time-domain spectral analysis suffers from spectral broadening and amplitude underestimation due to frequency drift under variable speed conditions. While short-time Fourier transform can partially alleviate this, it is limited by the time-frequency resolution of the window function. Angle-domain resampling transforms the time axis into an angle axis, making the angle-domain frequencies of each order component independent of the rotational speed, fundamentally eliminating frequency ambiguity. The Hilbert-Huang transform adaptively decomposes non-stationary signals through empirical mode decomposition, without relying on fixed basis functions, thus adapting to the extraction of nonlinear and non-stationary vibration features.

[0071] The third key technical approach is sparse Bayesian learning for multi-excitation component separation and uncertainty quantification. Traditional blind source separation methods typically assume statistical independence of source signals, but the three components—unbalanced excitation, bearing impact, and oil film vortex—are not strictly statistically independent in terms of their physical mechanisms. Sparse Bayesian learning utilizes the sparsity differences of the three types of excitations in the atomic dictionary space for separation. The physical prior is embedded in the dictionary structure, and the separation accuracy does not depend on the statistical independence assumption. The posterior variance after the expectation-maximization algorithm converges directly provides the uncertainty quantification, enabling the separation results to have statistical reliability assessment capabilities.

[0072] The synergistic effect of the three key technical approaches is reflected in the following aspects: dynamic thermal compensation provides a clean vibration signal for angle domain analysis, eliminating the contamination of the order spectrum by thermal drift; angle domain resampling provides a fixed-frequency feature vector for sparse Bayesian decomposition, making the construction of the atomic dictionary deterministic; and the component amplitude envelope and uncertainty quantification of the output of sparse Bayesian decomposition provide high-quality input for two-layer Nash equilibrium optimization, ensuring the reliability of closed-loop control commands. These three approaches progress sequentially in the causal chain, and the accuracy improvement at any stage can be inherited and amplified by subsequent stages, ultimately achieving a synergistic gain in overall detection accuracy.

[0073] It should be noted that this invention also solves the following technical problem: In the vibration detection of large inertia flywheels, when there is an early, weak fault in the bearing, the fault impact amplitude is much lower than the unbalanced excitation amplitude, and the two overlap in the frequency band. Traditional single-band filtering methods cannot separate the bearing impact component from the strong background vibration, leading to missed detection of early faults. This invention automatically locates the optimal frequency band where the impact energy is most concentrated through kurtosis spectrum, and then independently extracts the bearing impact component from the multidimensional vibration feature vector through sparse Bayesian learning. This achieves effective separation of weak impact components in a strong background, solving the technical problem that the impact characteristics of early bearing faults are submerged and cannot be detected.

[0074] Furthermore, this invention addresses the technical problem that during the acceleration and deceleration of the flywheel shaft over a wide speed range, there is a coupling effect between the bearing preload and the bearing housing temperature rise, and that adjusting the bearing preload or thermal compensation parameters alone cannot simultaneously achieve the optimal vibration suppression state of the system. This invention constructs a two-layer Nash equilibrium optimization model with coupled terms, incorporating the optimization of bearing preload adjustment and thermal compensation correction into a unified game framework. Under the condition that both objective functions reach Nash equilibrium, a joint optimal adjustment command is output, overcoming the defect of suboptimal adjustment commands caused by neglecting coupling relationships in single-layer independent optimization.

[0075] Specifically, the principle of this invention is as follows: The fundamental reason why this invention can solve the above-mentioned technical problems lies in the fact that its technical solution establishes corresponding physical models and algorithm mechanisms for the three independent causes of vibration amplitude distortion, and integrates the three into the same detection framework through closed-loop feedback. First, the thermal expansion finite element model calculates the static component of thermal deformation at each measuring point in real time using temperature field data as boundary conditions. This static component is completely independent of dynamic vibration in terms of physical mechanism, so direct subtraction can eliminate thermal expansion interference without losing any vibration information. The logical rationality of this operation is based on the finite element method's ability to accurately describe thermal deformation numerically. Second, the angle domain resampling guided by the key phase sensor transforms the time axis into the angle axis, so that the vibration order always corresponds to a fixed angle domain frequency component during the variable speed process, thereby completely eliminating frequency ambiguity at the Fourier transform level. The physical basis of this transformation is the deterministic correspondence between the order characteristics of rotating machinery vibration and the angle. Furthermore, sparse Bayesian learning decouples components by applying sparse priors to an overcomplete atom dictionary and utilizing the physical facts of different atoms occupying different spaces in the dictionary space through three types of excitation mechanisms. The posterior variance after the expectation-maximization algorithm converges directly provides a quantification of uncertainty, logically guaranteeing the statistical reliability of the decomposition results. These three mechanisms are sequentially progressive in the causal chain, jointly ensuring the accuracy and robustness of vibration amplitude extraction. Therefore, the technical solution of this invention is logically self-consistent and complete.

[0076] The following provides a specific embodiment 1 of the present invention, and the specific implementation of each step in this embodiment 1 is described in detail below.

[0077] The specific implementation of step S01 is as follows: A piezoelectric sensor array is installed at each measuring point on the flywheel bearing housing, and a high-precision thermocouple array is simultaneously arranged to collect the original displacement signal of the flywheel shaft within the speed range of 0–10000 r / min and the temperature field data of the flywheel bearing housing. The original displacement signal output by the piezoelectric sensor array is denoted as... The temperatures at each measuring point output by the thermocouple array are denoted as follows: ,in Sampling time, This is the measurement point number.

[0078] The specific implementation of step S02 is as follows: based on the flywheel bearing housing temperature field data The thermal expansion finite element model is invoked to calculate the static components of thermal deformation at each measuring point in real time. The thermal expansion finite element model uses the geometry of the flywheel bearing housing and the thermal expansion coefficient of the material as structural inputs. As the boundary conditions, the first... static component of thermal deformation at measuring point The dynamic thermal compensation correction formula is expressed as follows:

[0079] ;

[0080] In the formula, For the first Vibration displacement signal output from the measuring point (unit: μm); For the first Original displacement signal of the measuring point (unit: μm); The first result obtained from the thermal expansion finite element model Static component of thermal deformation at measuring point (unit: μm); The displacement reference value (unit: μm) has an empirical value of 1 μm. For the first The residual error term for thermal compensation at the measuring point (unit: μm) reflects the fitting deviation between the finite element model and the measured temperature field, and its range is usually ±0.5μm.

[0081] The specific implementation of step S03 is as follows: [The text abruptly ends here, likely due to an incomplete sentence or a formatting error.] Construct a multi-resolution bandpass filter bank to decompose the signal into Sub-band signal. The center frequency of each sub-band is (Unit: Hz), bandwidth is (Unit: Hz), the corresponding sub-band signal is denoted as (Unit: μm), where Calculate the kurtosis value for each sub-band signal. kurtosis spectrum and Two-dimensional coordinate axes, with The vertical height is plotted. The formula for calculating the kurtosis value is as follows:

[0082] ;

[0083] In the formula, For the first The dimensionless kurtosis value of the sub-band signal; This represents the number of sampling points for the sub-band signal. For the first Sub-band signal Amplitude of each sampling point (unit: μm); For the first Mean value of sub-band signal (unit: μm); For the first Standard deviation of sub-band signal (unit: μm); This is the baseline value for kurtosis, defaulting to 3, corresponding to the Gaussian kurtosis; all items in parentheses are dimensionless ratios. Overall dimensionless. Automatically selects the sub-band with the largest kurtosis value. For the optimal frequency band, Extracting the instantaneous amplitude envelope using Hilbert transform The formula is expressed as follows:

[0084] ;

[0085] In the formula, for The Hilbert transform result (unit: μm), i.e. ; This represents the instantaneous amplitude envelope (unit: μm); the meanings of other variables are the same as before. Then... The envelope spectrum is obtained by performing a fast Fourier transform, and the characteristic frequencies of each bearing fault are then analyzed within the envelope spectrum. (Unit: Hz) The search is performed to find and accumulate the amplitudes of each harmonic order at the fundamental frequency. Type 1 Fault The amplitude of the first harmonic wave is denoted as (Unit: μm), forming the modulation amplitude matrix The formula is expressed as follows:

[0086] ;

[0087] In the formula, This is a dimensionless modulation amplitude matrix, with rows corresponding to bearing fault types and columns corresponding to each order of harmonics. This represents the total number of bearing failure types. To search for harmonic orders; The amplitude reference value (unit: μm) has an empirical value of 1 μm.

[0088] The specific implementation of step S04 is as follows: using the bond phase pulse signal acquired by the bond phase sensor as a reference, ... Resampling to equal-angle domain signal (Unit: μm), where This represents the flywheel shaft rotation angle (unit: rad). Perform a Fast Fourier Transform in the angle domain to eliminate frequency ambiguity, and then combine it with a Hilbert-Huang Transform to extract the instantaneous amplitude envelope in the angle domain. Hilbert-Huang transform first... Empirical mode decomposition is performed, and the formula is expressed as follows:

[0089] ;

[0090] In the formula, For the first One intrinsic mode function (unit: μm); Residual components of empirical mode decomposition (unit: μm); This represents the total number of intrinsic mode functions. The instantaneous amplitude envelope reference value in the angle domain is 1 μm, with an empirical value of 1 μm. The instantaneous amplitude is then obtained by performing a Hilbert transform on each eigenmode function, as expressed in the following formula:

[0091] ;

[0092] In the formula, for Hilbert transform results (unit: μm); For the first The instantaneous amplitude (unit: μm) of each eigenmode function. Take the dominant eigenmode function... As (Unit: μm). [The following appears to be a separate, unrelated sentence:] Will... The sampled values ​​at each angle and The elements in the vector are arranged in order to form a multidimensional vibration characteristic vector. The formula is expressed as follows:

[0093] ;

[0094] In the formula, This is a multidimensional vibration feature vector with dimension . Each element is dimensionless; For the first Angle sampling points (unit: rad) ); Total number of sampling points in the angle domain; superscript This indicates transpose.

[0095] The specific implementation of step S05 is as follows: The multidimensional vibration feature vector... In a complete atomic dictionary Sparse decomposition is then performed. It is composed of the structural mode vector and the modulation waveform, and the formula is expressed as follows:

[0096] ;

[0097] In the formula, For an overcomplete atomic dictionary matrix, the dimension is... ; For the first Atom vectors ( ), dimension is Each element is dimensionless and consists of structural mode vectors or modulation waveforms. The total number of atoms satisfies The sparse decomposition model formula is expressed as follows:

[0098] ;

[0099] In the formula, Let be a sparse coefficient vector with dimension O(n). The elements are dimensionless; To decompose the residual vector, the dimension is The elements are dimensionless; , and All elements are dimensionless, and both sides of the equation have the same dimensions. Within the Bayesian inference framework, for each atomic coefficient... Apply a sparse prior with a zero-mean Laplace distribution, i.e. ,in This is the regularization parameter (dimensionless), with an empirical value range of 0.01 to 0.1; express The absolute value. The expectation step (calculation) is performed alternately using the expectation-maximization algorithm. Posterior distribution expectation The maximization step (updating hyperparameters to maximize the log-marginal likelihood) iteratively updates the mean vector of the posterior distribution of sparse coefficients. (Dimension is) ) and posterior covariance matrix (Dimension is) ),in for The The diagonal element represents the th element. The posterior variance of each atomic coefficient; only retaining Significantly non-zero atoms will Decomposed into unbalanced excitation components Bearing impact component With oil film eddy component (all dimensions are) (Each element is dimensionless), and the corresponding amplitude envelopes are extracted respectively. (Unit: μm) (Unit: μm) (Unit: μm), the uncertainty quantification result is the posterior variance of the corresponding atom. .

[0100] The specific implementation of step S06 is as follows: An upper-level optimization model is constructed with the objective of minimizing the amplitude of the unbalanced excitation components. The objective function formula is expressed as follows:

[0101] ;

[0102] In the formula, The value of the upper-level objective function (unit: μm); The baseline value for the upper-level objective function (unit: μm); The amplitude envelope of the unbalanced excitation component (unit: μm); Reference value for the amplitude envelope of the unbalanced excitation component (unit: μm); The flywheel shaft speed (unit: r / min); This is the reference value for rotational speed (unit: r / min); , , These are dimensionless weighting coefficients, with empirical values ​​of 0.5, 0.3, and 0.2, respectively. Bearing preload (unit: N); Temperature of the flywheel bearing housing (unit: °C); This is a coupling function of bearing preload and flywheel bearing housing temperature (unit: N·℃). Reference value for the coupling function (unit: N·℃); This is the reference value for bearing preload (unit: N); The reference temperature value for the flywheel bearing housing (unit: °C) is 20 °C (empirical value); the constraint condition is... ,in A preset threshold (unit: μm) is set for the amplitude envelope of the unbalanced excitation component. A lower-level optimization model is constructed with the objective of minimizing the bearing preload adjustment error. The objective function formula is as follows:

[0103] ;

[0104] In the formula, The value of the lower-level objective function (unit: N); The baseline value for the lower-level objective function (unit: N); To maintain the temperature of the flywheel bearing housing Temperature-adaptive optimal preload under given conditions (unit: N); The weight is a dimensionless weighting coefficient with an empirical value of 0.2; the constraint condition is... Within the allowable adjustment range and with the flywheel bearing housing temperature rise not exceeding 90℃, the upper-level optimization model and the lower-level optimization model are coupled through a coupling term. Interactions, using Nash equilibrium to solve for the optimal bearing preload (Unit: N) and optimal thermal compensation correction amount (Unit: μm) joint configuration; when , or If any one of them exceeds the preset threshold, according to and Output adjustment commands to correct the flywheel bearing preload and dynamic thermal compensation correction amount, and complete the closed-loop detection of vibration amplitude.

[0105] To better understand and implement this invention, the following is a specific application scenario of the invention, Example 2: To verify the effectiveness of the invention, technicians set up a test environment, using a high-inertia flywheel energy storage device in actual operation as the object. The flywheel housing of this device weighs approximately 1.4t, and its overall dimensions are 1700mm × 900mm × 1906mm. The designed rotational inertia of the flywheel shaft is... The maximum design speed is 10,000 r / min, the total weight of the inertial flywheel assembly is approximately 4 tons, and the center height of the flywheel shaft is 1300 mm. A piezoelectric sensor array is arranged with four measuring points (two axial and two radial) on each side of the flywheel bearing housing. A thermocouple array is synchronously installed near each measuring point, with a sampling frequency set to 100 kHz and a thermocouple sampling period of 0.5 s. A key phase sensor is installed at one end of the flywheel shaft, outputting 1024 key phase pulses per revolution. The test condition is the entire acceleration process of the flywheel shaft from 0 r / min to 10,000 r / min over an acceleration time of 120 s, to fully cover non-stationary conditions over a wide speed range.

[0106] During the signal acquisition phase, the raw displacement signal and temperature field data are stored in the buffer after being synchronized and timestamped by the acquisition card. The flywheel bearing housing temperature rises from an initial 25℃ to a maximum of 82℃ during acceleration, a temperature increase of 57℃. Figure 3 The figure shows the temperature field distribution curve of the flywheel bearing housing as a function of rotational speed, intuitively reflecting the nonlinear characteristics of the temperature rise process. The thermal expansion finite element model uses the linear expansion coefficient of the flywheel bearing housing steel. Using material parameters, the static component of thermal deformation at each measuring point was numerically solved in real time. The peak value of the static component of thermal deformation at the radial measuring point of the left bearing housing was approximately 38 μm, and the corresponding value on the right was approximately 35 μm. After dynamic thermal compensation correction, the thermal drift component in the vibration displacement signal was completely separated, and the peak-to-peak value of the output vibration displacement signal was significantly reduced. The comparison of the peak-to-peak values ​​of vibration displacement at each measuring point before and after thermal compensation is shown in Table 1.

[0107] Table 1 Comparison of peak-to-peak vibration displacement values ​​at various measuring points before and after thermal compensation

[0108]

[0109] In the optimal frequency band localization stage, a multi-resolution bandpass filter bank was constructed for the vibration displacement signal, with a decomposition layer of 6, generating a total of 63 sub-frequency bands. Kujicic spectrum calculation results show that within the sub-frequency band with a center frequency of approximately 12800Hz and a bandwidth of approximately 3200Hz, the kurtosis value reaches 23.7, far exceeding that of other frequency bands, such as... Figure 4 The diagram shows the two-dimensional distribution of the kurtosis spectrum. The horizontal axis represents the center frequency of the sub-band, and the vertical axis represents the bandwidth of the sub-band. The color intensity indicates the magnitude of the kurtosis value, and the location of the optimal frequency band is clearly identifiable. A Hilbert transform is performed on the filtered signal of the optimal frequency band, and a Fast Fourier Transform is performed on the instantaneous amplitude envelope to obtain the envelope spectrum. The amplitudes of each harmonic family at the characteristic frequency of the bearing outer ring fault (approximately 87.4 Hz) are searched within the envelope spectrum. The amplitudes of the 1st to 5th harmonics are accumulated to obtain an outer ring fault modulation amplitude of 12.3 μm and an inner ring fault modulation amplitude of 4.1 μm. The modulation amplitude matrix is ​​shown in Table 2.

[0110] Table 2 Modulation amplitude matrix for various bearing fault types (unit: μm)

[0111]

[0112] In the angle domain analysis stage, using the 1024 pulse / revolution signal output by the key phase sensor as a reference, the vibration displacement signal was resampled at equal angle intervals to convert the entire time domain signal of the 120s acceleration process into an equal angle domain signal. After performing a fast Fourier transform on the equal angle domain signal, the spectral lines of each order component were clear, with a first harmonic (unbalanced excitation principal order) amplitude of 47.6 μm and a second harmonic amplitude of 8.3 μm, eliminating the frequency ambiguity effect. Subsequently, empirical mode decomposition was performed on the equal angle domain signal, yielding eight intrinsic mode functions. After performing a Hilbert transform on each intrinsic mode function, the instantaneous amplitude envelope in the angle domain was synthesized according to energy weights, as shown below. Figure 5 The curve shown is the change curve of the instantaneous amplitude envelope in the angular domain with the rotation angle, reflecting the instantaneous evolution law of the vibration amplitude of the flywheel shaft during the acceleration process.

[0113] In the sparse Bayesian decomposition stage, the overcomplete atom dictionary consists of a 10th-order structural mode vector and modulation waveform atoms corresponding to three types of excitation mechanisms, with a total of 152 atoms. The multidimensional vibrational feature vector has a dimension of 64, and the dictionary overcompleteness ratio is 2.375. The expectation-maximization algorithm converges after 37 iterations, and the relative change in the posterior mean is reduced to... The convergence criterion is met. The decomposition results show that the mean amplitude envelope of the unbalanced excitation component is 45.2 μm, with a corresponding uncertainty of 0.031 μm; the mean amplitude envelope of the bearing impact component is 13.1 μm, with a corresponding uncertainty of 0.044 μm; and the mean amplitude envelope of the oil film whirl component is 6.7 μm, with a corresponding uncertainty of 0.028 μm. The uncertainties of all three types of components are lower than the upper reference limit of 0.05 μm, and the confidence level of the extraction results meets the requirements.

[0114] During the closed-loop control phase, the mean value of the bearing impact component amplitude envelope (13.1 μm) exceeded the preset threshold of 10 μm, triggering the system adjustment process. The upper-level and lower-level optimization models are coupled by a term; when the current flywheel bearing housing temperature is 78℃, the temperature-adaptive optimal preload is determined. The initial preload was 4200N, while the current actual preload was 3600N, with a deviation of 600N. After 11 iterations of Nash equilibrium, the system converged, outputting an optimal bearing preload of 4180N and an optimal thermal compensation correction value of +2.3μm. Based on the output command, the system drove the preload actuator to adjust the bearing preload to 4180N and updated the thermal compensation correction value. After adjustment, the mean value of the bearing impact component amplitude envelope decreased to 7.8μm, which is lower than the preset threshold. The closed-loop detection was completed. The comparison of each component before and after adjustment is shown in Table 3.

[0115] Table 3 Comparison of amplitude envelopes of each component before and after closed-loop control (unit: μm)

[0116]

[0117] like Figure 3 As shown, the temperature field of the flywheel bearing housing changes with rotational speed, clearly demonstrating the nonlinear characteristics of temperature rise, which confirms the necessity of using a finite element model for real-time dynamic thermal compensation. During the full testing of the embodiments, compared with the traditional fixed-band filtering plus time-domain spectrum analysis scheme, the present invention demonstrates the following improvements at the principle level: The traditional scheme sends thermal drift and vibration signals into the filter together, and the interference of the slowly varying thermal expansion component on the low-frequency vibration amplitude cannot be eliminated. However, the present invention directly separates the two types of components through physical modeling, ensuring the physical purity of the vibration signal from the signal source. The traditional scheme leads to underestimation of amplitude due to spectral line broadening under variable speed conditions. The present invention converts the influence of speed change into a fixed order that strictly corresponds to the angle through angle domain resampling, and the amplitude estimation accuracy is not affected by speed fluctuations. The traditional scheme cannot physically decouple the three types of components: unbalanced excitation, bearing impact, and oil film vortex. Threshold alarms under aliasing conditions have the risk of false alarms or missed alarms. The present invention achieves independent extraction of the three types of components through sparse Bayesian learning. Each type of component has an independent amplitude envelope and uncertainty quantification result, and the physical basis of the alarm logic is more sufficient.

[0118] It should be noted that the variables involved in this invention are explained in detail in Tables 4 and 5.

[0119] Table 4. Variable Explanation Table (Part 1)

[0120]

[0121] Table 5. Variable Explanation Table (Part Two)

[0122]

[0123] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for detecting the amplitude of vibration of a high-inertia flywheel in a high-speed rotating state, characterized by, Includes the following steps: A piezoelectric sensor array is mounted on the flywheel bearing housing, and a high-precision thermocouple array is synchronously arranged at each measuring point on the flywheel bearing housing to collect the original displacement signal of the flywheel shaft in the range of 0 to 10000 r / min and the temperature field data of the flywheel bearing housing. Based on the temperature field data of the flywheel bearing housing, the thermal expansion finite element model is called to calculate the static component of thermal deformation at each measuring point in real time, and the original displacement signal is dynamically thermally compensated and corrected. The slow thermal expansion component and the fast vibration component are separated, and the vibration displacement signal is output. A multi-resolution bandpass filter bank is constructed for the vibration displacement signal. The kurtosis value is calculated for each sub-band signal and the kurtosis spectrum is plotted. The optimal band with the largest kurtosis value is automatically selected for bandpass filtering. The Hilbert transform is performed on the filtered signal to extract the instantaneous amplitude envelope. The fast Fourier transform is performed on the instantaneous amplitude envelope to obtain the envelope spectrum. The harmonic family amplitudes are searched in the envelope spectrum with the bearing fault characteristic frequencies as the fundamental frequencies and accumulated. The modulation amplitude matrix is ​​then output. Using the key phase pulse signal acquired by the key phase sensor as a reference, the vibration displacement signal is resampled into an equal-angle domain signal. In the angle domain, a fast Fourier transform is performed on the equal-angle domain signal to eliminate the frequency ambiguity caused by the rotational speed change. Then, the result of the fast Fourier transform in the angle domain is combined with the Hilbert-Huang transform to extract the instantaneous amplitude envelope in the angle domain. The instantaneous amplitude envelope in the angle domain and the elements in the modulation amplitude matrix are arranged in order to form a multidimensional vibration feature vector. The multidimensional vibration feature vector is sparsely decomposed under an overcomplete atom dictionary. A sparse prior is established using a Bayesian inference framework. The mean and variance of the posterior distribution of the sparse coefficients are iteratively updated using the expectation-maximization algorithm. Only atoms with significantly non-zero posterior variance are retained. The multidimensional vibration feature vector is decomposed into unbalanced excitation components, bearing impact components, and oil film whirl components. The amplitude envelopes of the unbalanced excitation components, bearing impact components, and oil film whirl components are extracted respectively. At the same time, the uncertainty quantification results of the amplitude envelopes of each component are output. An upper-level optimization model is constructed with the goal of minimizing the amplitude of the unbalanced excitation component, and a lower-level optimization model is constructed with the goal of minimizing the bearing preload adjustment error. The upper-level optimization model and the lower-level optimization model are coupled through coupling terms to solve the Nash equilibrium and obtain the optimal bearing preload and the optimal thermal compensation correction amount. When any of the amplitude envelopes of the unbalanced excitation component, the bearing impact component, or the oil film whirl component exceeds a preset threshold, an adjustment command is output based on the optimal bearing preload and the optimal thermal compensation correction amount to correct the flywheel bearing preload state and the dynamic thermal compensation correction amount, thus completing the closed-loop detection of vibration amplitude.

2. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 1, characterized in that, The dynamic thermal compensation correction specifically refers to the real-time subtraction of the slow thermal expansion component in the original displacement signal by calculating the static thermal deformation component of each measuring point based on the thermal expansion finite element model at the signal processing layer, so that the output vibration displacement signal retains only the dynamic vibration fast component.

3. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 2, characterized in that, The thermal expansion finite element model is specifically a calculation model that uses the geometric structure and material thermal expansion coefficient of the flywheel bearing housing as structural inputs and the temperature field data of the flywheel bearing housing as boundary conditions. It numerically solves the static components of thermal deformation generated at each measuring point under temperature rise conditions using the finite element method.

4. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 3, characterized in that, The multi-resolution bandpass filter bank is specifically a filtering structure composed of a group of bandpass filters with different center frequencies and bandwidths, used to decompose the vibration displacement signal into multiple sub-band signals so that kurtosis values ​​can be calculated for different frequency bands.

5. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 4, characterized in that, The kurtosis spectrum is specifically a two-dimensional spectrum plotted with the sub-band center frequency and sub-band bandwidth as two-dimensional coordinate axes and the kurtosis value of each sub-band signal as the vertical height. The optimal frequency band is the sub-band with the largest kurtosis value.

6. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 5, characterized in that, The modulation amplitude matrix is ​​specifically a matrix composed of the characteristic frequencies of each bearing fault and the sum of their harmonic amplitudes. The rows correspond to the bearing fault types, and the columns correspond to the amplitudes of each order of harmonics.

7. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 6, characterized in that, The key phase sensor is specifically a sensor installed on the flywheel shaft that outputs a key phase pulse signal every time the flywheel shaft rotates through a fixed angle. The key phase pulse signal is used to mark the rotation angle position of the flywheel shaft so as to achieve uniform resampling of the vibration displacement signal to the equal angle domain signal.

8. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 7, characterized in that, The angular domain signal is specifically a signal sequence obtained by uniformly resampling the vibration displacement signal with the flywheel shaft rotation angle as the horizontal axis. This ensures that the frequency of each order vibration component in the angular domain remains fixed under variable speed conditions, thereby eliminating the frequency ambiguity effect caused by speed fluctuations.

9. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 8, characterized in that, The Hilbert-Huang transform is a method that first performs empirical mode decomposition on the signal to obtain several intrinsic mode functions, and then performs Hilbert transform on each intrinsic mode function to obtain the instantaneous amplitude and instantaneous frequency. It is suitable for the instantaneous feature extraction of non-stationary vibration signals.

10. The method for detecting the vibration amplitude of a large-inertia flywheel under high-speed rotation as described in claim 9, characterized in that, The overcomplete atomic dictionary is specifically a dictionary matrix composed of structural mode vectors and modulation waveforms, with the number of atoms exceeding the dimension of the signal. The structural mode vectors represent the natural vibration modes of each order of the flywheel shaft system, and the modulation waveforms represent the periodic modulation components of bearing impact and oil film vortex.

Images