A centrifugal air compressor stall and surge precursor identification and anti-surge control method

Abstract

The application belongs to the technical field of fuel cells, and particularly relates to a centrifugal air compressor stall surge precursor identification and anti-surge control method. The method preprocesses dynamic operation data of a fuel cell centrifugal air compressor, constructs an observation characteristic vector and a vortex structure dictionary based on a sliding short-time window, solves a sparse vortex structure state code and determines an uncertainty index, constructs a precursor risk index and a risk level in combination with a grouping and aggregation result and a change amount of adjacent short-time windows, generates an anti-surge control instruction under the condition of meeting a gas supply demand and an actuator constraint, and gradually recovers after the risk is removed. According to the precursor risk index and the uncertainty, the vortex structure dictionary is updated or frozen. The application can realize early and low false alarm identification of the centrifugal air compressor stall surge precursor, and takes into account minimum intervention of anti-surge control, gas supply stability and system operation efficiency.

Description

Technical Field

[0001] This application belongs to the field of fuel cell technology, and specifically relates to a method for identifying and preventing surge precursors in centrifugal air compressors. Background Technology

[0002] The net efficiency and dynamic response performance of a fuel cell system largely depend on the operating status of the air supply system. As the core component of the fuel cell air supply system, the centrifugal air compressor is responsible for continuously supplying air at a certain pressure and flow rate to the fuel cell stack. Its operational stability directly affects the output performance, energy utilization efficiency, and overall reliability of the fuel cell system.

[0003] In applications such as fuel cell vehicles and distributed power generation, air supply systems typically need to adapt to rapid load changes, environmental boundary fluctuations, and wide operating conditions. Due to factors such as pipeline volume effects, valve actuation lag, sensor delays, and control response speed, centrifugal air compressors are prone to entering high-risk areas near stall or surge during dynamic operation. Once stall or surge occurs, it not only leads to a decrease in compression efficiency but can also cause severe pressure fluctuations, sudden changes in bearing loads, and current surges in the drive motor, thereby affecting the lifespan and operational reliability of both the air supply system and the fuel cell system. Therefore, stall and surge identification and anti-surge control for centrifugal air compressors has always been a crucial technical challenge in the field of fuel cell air management.

[0004] Currently, engineering solutions for surge control of centrifugal air compressors typically use the relationship between the operating point and the surge boundary as the criterion, or use characteristic parameters such as pressure pulsation as the basis for instability triggering, and suppress instability by adjusting the bypass valve opening and torque command.

[0005] A common approach is a threshold-based determination method based on surge lines or surge margins. This type of method typically pre-sets a surge line based on the air compressor's pressure-flow ratio characteristic diagram and estimates the distance or margin between the current operating point and the surge line in real time. When the margin falls below a preset threshold, anti-surge actions such as rapid opening of a bypass valve are triggered to move the operating point away from the surge region. For example, patent document CN201610860217.9 discloses a surge control method and system. This type of method is relatively straightforward in its implementation, but it essentially still falls under boundary threshold-triggered control.

[0006] Another approach is based on instability determination methods using pressure pulsation amplitude, characteristic frequency, or frequency domain energy distribution. These methods typically utilize pressure sensors at the outlet, volute, or corresponding pipeline locations to extract parameters such as the pulsation amplitude, peak power spectrum, or low-frequency energy ratio of the pressure signal. When these characteristic parameters exceed thresholds, the air compressor is determined to have entered an unstable state, and anti-surge control is then implemented. For example, patent document CN202211478186.2 discloses a surge control method, system, electronic equipment, and storage medium for a centrifugal air compressor.

[0007] Another type of approach is predictive control based on fixed models. This type of method typically uses linearized or fixed-parameter models to model the air compressor, piping, and related actuators. Under the premise of satisfying constraints, it optimizes the calculation of control quantities such as bypass valve commands and torque commands to achieve anti-surge control. For example, patent document CN202510591655.9 discloses an anti-surge control method for compressor units.

[0008] However, existing technologies still have the following shortcomings in the wide operating conditions and high dynamic application scenarios of fuel cell centrifugal air compressors:

[0009] First, threshold determination methods based on surge lines or surge margins typically require the operating point to be close to the surge boundary before triggering control actions, which can easily lead to response lag during rapid operating condition transitions. To ensure safety, engineering practices often set thresholds conservatively, resulting in frequent opening of bypass valves, increased gas supply energy consumption, and reduced net output power and overall operating efficiency of the fuel cell system.

[0010] Secondly, while methods based on pressure pulsation amplitude or single frequency domain characteristics can reflect instability phenomena to some extent, they typically struggle to effectively distinguish between different instability mechanisms such as pre-stallization precursors, enhanced diffuser separation, and low-frequency system oscillations. Furthermore, they are susceptible to factors like sensor noise, zero-point drift, installation location differences, asynchronous sampling, and pipeline resonance, leading to false alarms or missed alarms. In addition, these methods lack a structured explanation of pre-stabilization precursors and have limited adaptability across different machine types and operating conditions.

[0011] Furthermore, predictive control methods based on fixed models are highly dependent on model accuracy and parameter consistency. In fuel cell air supply systems, variations in pipeline volume, valve characteristic drift, sensor delay changes, and system boundary fluctuations are all common. Deviations can easily occur between the actual situation and the preset model, thus affecting the prediction results and control effectiveness, and weakening the robustness and generalization ability of anti-surge control.

[0012] In addition, some existing data-driven instability identification methods tend to focus on black-box feature extraction or single statistical index judgment, lacking interpretable characterization of the instability precursor evolution process, and are quite sensitive to data quality, sensor consistency and operating condition distribution, making it difficult to balance early identification capability, false alarm rate control and control intervention costs.

[0013] Therefore, existing technologies lack a stall / surge identification and anti-surge control technology solution that can address a wide range of operating conditions for fuel cell centrifugal air compressors, while also considering interpretable characterization of instability precursors, early identification, low false alarm warnings, and minimal intervention control, in order to achieve efficient, stable, and reliable operation of the air supply system. Summary of the Invention

[0014] To address the problems of existing technologies based on surge margin, single pulsation threshold, or fixed models in the wide operating conditions and high dynamic scenarios of fuel cell centrifugal air compressors, such as delayed precursor identification, difficulty in distinguishing different instability mechanisms, susceptibility to noise and operating condition drift, high false alarm and false alarm rates, and overly conservative anti-surge control leading to increased energy consumption and decreased system efficiency, this application provides a method for identifying and controlling stall surge precursors in centrifugal air compressors. This method constructs a sliding short-time window feature from dynamic operating data, establishes a vortex structure dictionary constrained by instability mechanism groups, and solves the problem... The sparse vortex structure state code with temporal prior constraints is used to construct a precursor risk index by integrating aggregation intensity, window variation, and uncertainty. Based on the precursor risk index, minimum intervention continuous anti-surge control is implemented to meet the air supply demand constraints and actuator constraints. At the same time, combined with the small-step online update and freeze mechanism of the dictionary, interpretable, early, and low false alarm identification of stall surge precursors and hierarchical coordinated control of bypass valves, motors, and / or guide vanes are realized. This improves the operational stability, control robustness, net efficiency of fuel cell system, and overall reliability of air supply system under wide operating conditions.

[0015] This application provides a method for identifying and preventing surge precursors in a centrifugal air compressor, the method comprising:

[0016] Step 1: Acquire dynamic operating data during the operation of the fuel cell centrifugal air compressor, and preprocess the dynamic operating data to obtain input data;

[0017] Step 2: Based on the input data, construct the observation feature vector according to the sliding short time window, and establish a vortex structure dictionary constrained by the instability mechanism grouping, and divide the dictionary primitives in the vortex structure dictionary into multiple instability mechanism groups;

[0018] Step 3: For the current short window, based on the observation feature vector, vortex structure dictionary and the status code corresponding to the previous short window, solve the sparse vortex structure status code that satisfies the residual consistency, sparsity constraint and cross-window consistency constraint, and determine the uncertainty index according to the reconstruction residual of the current short window and / or the change in the status code of adjacent short windows.

[0019] Step 4: Aggregate the sparse vortex structure state codes according to the group to which the dictionary primitives belong, obtain the aggregation intensity of each group, and construct the precursor risk index based on the aggregation intensity of each group, the change in aggregation intensity of adjacent short time windows, and the uncertainty index. Determine the risk level based on the comparison result between the precursor risk index and the preset risk threshold.

[0020] Step 5: Based on the precursor risk index and its changing trend, under the conditions of meeting the gas supply demand constraints and actuator constraints, generate continuous control commands for the bypass valve opening, selectively generate continuous control commands for motor torque and / or guide vane angle, and when the precursor risk index decreases to the preset recovery condition, perform gradual recovery of the control commands to restore the air compressor to a high-efficiency operating state.

[0021] Step 6: Based on the precursor risk index and uncertainty, when the preset update conditions are met, the vortex structure dictionary is updated in small steps. When the preset freeze conditions are met, the update of the vortex structure dictionary is frozen and the current dictionary is maintained to participate in the solution of the sparse vortex structure state code.

[0022] In a preferred implementation, step 3 further includes:

[0023] Step 3.1: For the current short window, the observed feature vector corresponding to the current short window is solved by non-negative sparse coding with time-series prior constraints on the vortex structure dictionary. The goal is to minimize the feature reconstruction error, and simultaneously introduces sparsity constraints, cross-window amplitude continuity constraints with the previous short window state code, and suppression constraints on newly added activation components to obtain the sparse vortex structure state code of the current short window.

[0024] Step 3.2: Based on the weighted reconstruction residual of the current short-term window observation feature vector and the sparse vortex structure state code, and combined with the support set change and code value difference intensity between the current short-term window state code and the previous short-term window state code, an uncertainty index is constructed to characterize the credibility of the current identification result.

[0025] In a preferred implementation, step 4 further includes:

[0026] Step 4.1: Group the dictionary primitives in the vortex structure dictionary according to the instability mechanism to which they belong, and perform grouping summation or weighted aggregation of the activation intensity of the corresponding dictionary primitives in the sparse vortex structure state code of the current short window to obtain the aggregation intensity of the rotating stall related primitive group, the separation enhancement primitive group, the system-level critical coupling primitive group, and the background and interpretable disturbance primitive group, and based on the change of the aggregation intensity of each group between the current short window and the previous short window between adjacent short windows;

[0027] Step 4.2: Based on the aggregation intensity of each instability mechanism group, the change in aggregation intensity of each group between adjacent short time windows, and the uncertainty index, construct a precursor risk index that jointly characterizes the precursors of rotating stall, separation enhancement, and system-level critical coupling, and performs confidence correction on the background disturbance. Then, compare the precursor risk index with at least one preset risk threshold to determine the risk level corresponding to the current short time window.

[0028] In a preferred implementation, step 5 further includes:

[0029] Step 5.1: Determine the anti-surge control target based on the fuel cell air supply demand, and establish a set of control constraints including target cathode pressure and / or target air flow, excess air coefficient constraints, and amplitude and rate of change constraints of bypass valves, motors and / or guide vane actuators;

[0030] Step 5.2: Based on the precursor risk index and its adjacent short-term window change trends, generate continuous control commands for bypass valve opening that satisfy the control constraint set, and selectively coordinate the generation of continuous control commands for motor torque and / or guide vane angle according to the risk level, so as to achieve minimum intervention and surge prevention control according to the risk increment.

[0031] Step 5.3: When the precursor risk index decreases below the preset recovery threshold and the preset stability conditions are continuously met, perform a gradual recovery of the continuous control commands for the bypass valve opening and the continuous control commands for the motor torque and / or guide vane angle, so as to restore the air compressor to a high-efficiency operating state.

[0032] In the preferred implementation, further, in step 3.1, the sparse vortex structure state code of the current short-time window is obtained by solving the following optimization problem:

[0033] ;

[0034] Where t represents the time sequence index of the current sliding short window; t-1 represents the previous short window adjacent to the current short window; This represents the observed feature vector corresponding to the t-th short time window; This represents the dictionary matrix of the vortex structure corresponding to the t-th short-time window; This represents the sparse vortex structure state code to be solved in the t-th short time window; This represents the optimal sparse vortex structure state code corresponding to the t-th short time window; This indicates the sparse vortex structure state code corresponding to the previous short time window; Let W represent the objective function for the t-th short-time window; W represents the feature domain weighting matrix. This is a constraint term for residual consistency; Indicates the sparse constraint weight coefficient; For the L1 norm sparse terms of the status code; This represents the weighting coefficient of the cross-window consistency constraint; This is a cross-window consistency constraint term; ρ represents the newly added activation suppression weight coefficient. This represents the support set mask vector constructed from the status codes of the previous short window. This represents a mask vector that is complementary to the support set of the previous short window; Add a new activation inhibition term; Indicates a nonnegativity constraint; Represents the L2 norm; Denotes the norm 1;

[0035] In step 3.2, the uncertainty index is:

[0036] ;

[0037] in, This represents the uncertainty index corresponding to the t-th short-time window; This represents the weighted reconstruction residual energy of the current short-term window; This indicates the status code corresponding to the t-th short time window. Support set; This indicates the status code corresponding to the (t-1)th short time window. Support set; △ represents symmetric difference operation; This represents the amount of change in the support set of adjacent short time windows; Indicates the strength of the code value difference between adjacent short-time window status codes; , , The weighting coefficient is greater than 0.

[0038] In the preferred implementation, further, in step 4.1, the polymerization strength is:

[0039] ;

[0040] in, The polymerization intensity of the g-th mechanism group within the current short time window t; Let i be the set of dictionary primitive indices corresponding to the g-th mechanism group; i is the mechanism group. The i-th dictionary primitive index; is the sparse vortex structure state code corresponding to the i-th dictionary primitive within the t-th short-time window; t is the current short-time window number.

[0041] The change in polymerization intensity, constructed from the difference between the polymerization intensity corresponding to the current short-term window and the previous short-term window, can be expressed as:

[0042] ;

[0043] in, This represents the change in polymerization intensity for the g-th mechanism group between two adjacent short time windows; The change in aggregation intensity of the g-th mechanism group in the current short time window relative to the previous short time window is represented by its positive or negative sign, respectively, indicating the strengthening or weakening trend of the instability characteristic corresponding to the mechanism group. This represents the polymerization intensity corresponding to the previous short time window;

[0044] In step 4.2, the precursor risk index is:

[0045] ;

[0046] in, It is a monotonic mapping function used to map the fusion result to a predetermined risk range; This represents the growth trend of the polymerization intensity of the rotational stall-related primitives over time. This represents the growth trend of the aggregation intensity of the critically coupled primitive group at the system level over time. , , The weighting coefficients corresponding to the aggregation intensity of each instability mechanism are used to characterize the contribution of the current activity level of different mechanisms to the total risk; , The weighting coefficients corresponding to the growth trend items; This is the correction weight for uncertainty.

[0047] In the preferred implementation, further, in step 5.1, the set of control constraints corresponding to the t-th control time is:

[0048] ;

[0049] in, This represents the control input vector at the t-th control time. The control input vector includes at least the bypass valve opening. Optionally includes motor torque and / or guide vane angle ; This represents the actual cathode pressure at the t-th control moment; This represents the target cathode pressure at the t-th control moment; This indicates that the cathode pressure tracking error does not exceed the preset pressure error upper limit. ; This represents the actual airflow rate at the t-th control time. This represents the target airflow rate at the t-th control time. This indicates that the airflow tracking error does not exceed the preset flow error limit. ; This represents the actual excess air coefficient at control time t; and These represent the lower and upper limits of the allowable excess air coefficient at the t-th control time, respectively. and These represent the control input vectors respectively. The lower limit and upper limit of amplitude; This represents the change in the control input vector at adjacent control moments; This indicates the maximum allowable variation in the control input vector.

[0050] In the preferred implementation, further, in step 5.2, the bypass valve opening command... for:

[0051] ;

[0052] in, This indicates the bypass valve opening command at the t-th control moment; Indicates the basic opening degree; This indicates the precursory risk index corresponding to the current short-term window; This indicates the precursor risk index corresponding to the previous short time window; This represents the static response gain of the risk index. Gain in response to changes in the risk index trend; This represents a saturation function.

[0053] In a preferred implementation, step 2 further includes:

[0054] Step 2.1: Slice the input data into a sliding short time window and extract observation features that characterize the frequency band energy distribution, cross-channel coupling relationship and the consistency of disturbance propagation, and construct an observation feature vector;

[0055] Step 2.2: Based on the observed feature vectors, a dictionary learning method constrained by instability mechanism is used to establish a vortex structure dictionary. According to different instability mechanisms, the dictionary primitives are divided into rotating stall related primitive group, separation enhancement primitive group, system-level critical coupling primitive group, and background and interpretable disturbance primitive group.

[0056] In a preferred implementation, step 6 further includes:

[0057] Step 6.1: Based on the precursor risk index and uncertainty, determine whether the vortex structure dictionary corresponding to the current short time window is in update mode or frozen mode. Specifically, when the precursor risk index is in a preset low-risk range and the uncertainty is below a preset threshold, the vortex structure dictionary is determined to be in update mode; when the precursor risk index is in a preset high-risk range and / or the uncertainty is above a preset threshold, the vortex structure dictionary is determined to be in frozen mode, and the current dictionary is maintained to participate in the sparse vortex structure state code solution.

[0058] Step 6.2: In update mode, based on the observed feature vector and sparse vortex structure state code corresponding to the current short time window, perform incremental update of the vortex structure dictionary with step size limit and forgetting factor, and apply instability mechanism grouping constraint projection to the dictionary primitives after the update.

[0059] The beneficial effects of this application are:

[0060] First, the centrifugal air compressor stall surge precursor identification and anti-surge control method of this application constructs a sliding short-time window feature from the dynamic operating data of the fuel cell centrifugal air compressor, and introduces a vortex structure dictionary constrained by instability mechanism grouping and a sparse vortex structure state code solution mechanism with cross-window consistency. This achieves a structured and interpretable characterization of different instability mechanisms such as rotating stall precursors, diffuser / volute separation enhancement, and low-frequency coupling of the system. Compared with schemes based on surge margin, single pressure pulsation threshold, or fixed model prediction, it can identify instability precursors earlier and effectively reduce false alarms and missed alarms caused by noise, drift, and changes in operating conditions. Simultaneously, this invention constructs a precursor risk index based on aggregation strength, cross-window variation, and uncertainty, and... Based on this, continuous control commands for bypass valve opening, motor torque, and / or guide vane angle are generated under the conditions of satisfying gas supply demand constraints and actuator constraints. After the risk is eliminated, progressive recovery is performed, thereby realizing hierarchical early warning and minimum intervention to prevent surge control in the early stage, reducing frequent bypasses, increased compression power consumption, and decreased system net efficiency caused by conservative control. In addition, by updating the vortex structure dictionary in small steps when preset conditions are met, and freezing the update under high-risk or high-uncertainty conditions while maintaining the current dictionary to participate in state code solving, this method has online adaptive capability, mechanistic semantic stability, and long-term operational robustness, making it more suitable for efficient, stable, and reliable operation of fuel cell centrifugal air compressors under wide operating conditions and high dynamic scenarios.

[0061] Second, in the preferred implementation, step 3 of this application introduces sparsity constraints, cross-window amplitude continuity constraints, and newly added activation suppression constraints simultaneously into the non-negative sparse coding solution of the current short time window. This makes the solution of the state code of the current short time window no longer an independent instantaneous matching process, but a continuous state estimation process that takes into account reconstruction accuracy, temporal smoothness, and support set stability. Therefore, it can effectively suppress code value jumps, frequent support set switching, and pseudo-activation phenomena caused by measurement noise, local resonance, sampling disturbances, and short-term anomalies, and improve the solution stability and result consistency of sparse vortex structure state codes in a wide range of dynamic processes. At the same time, step 3.2 further integrates the weighted reconstruction residual, support set change, and code value difference intensity to construct an uncertainty index, which can quantitatively characterize the credibility of the current identification result and provide a credible basis for subsequent precursor risk index calculation, risk classification, and control decisions.

[0062] Third, in the preferred implementation, step 4 of this application aggregates the sparse vortex structure state codes according to instability mechanisms and further integrates the intensity changes and uncertainty indices of each group between adjacent short time windows to construct a precursor risk index. This elevates the local activation information originally scattered across multiple dictionary primitives into a mechanism-level risk representation oriented towards rotating stall precursors, separation enhancement precursors, and system-level critical coupling precursors. It can simultaneously reflect the strength of the current existence of various instability precursors, avoiding the one-sidedness caused by making judgments based solely on instantaneous amplitudes or single features. At the same time, by introducing background and interpretable disturbance grouping and uncertainty correction mechanisms, it can perform confidence correction on non-instability anomalies caused by factors such as valve action, speed regulation, and measurement disturbances, reducing the interference of background disturbances on risk judgment and improving the authenticity and discriminability of risk classification results.

[0063] Fourth, in the preferred implementation, step 5 of this application first establishes a set of control constraints including cathode pressure, air flow, excess air coefficient, and amplitude and rate of change limits for bypass valves, motors, and / or guide vane actuators. Then, it generates continuous control commands based on the precursor risk index and its changing trend, and performs gradual recovery after the risk is eliminated. This transforms the anti-surge control from the traditional single-threshold triggering and large-amplitude correction method into a continuous constraint optimization intervention process oriented towards the air supply target and the physical boundary of the actuator. This allows for smooth release of control action in the early stages of instability risk, avoiding air supply fluctuations, control shocks, and actuator fatigue caused by sudden opening and closing of bypass valves, rapid torque jumps, or drastic adjustments to the guide vane angle. At the same time, coordinated control is implemented as the risk level changes, and intervention is gradually withdrawn during the recovery phase. This helps to ensure that the air compressor is away from the risk of stall surge while also taking into account the air supply stability of the fuel cell system.

[0064] Fifth, in the preferred implementation, step 6 of this application uses the precursor risk index and uncertainty as the dictionary update gating conditions, allowing incremental correction of the vortex structure dictionary with step size restrictions and forgetting factors only within a short time window with low risk and high confidence. After the update, the instability mechanism grouping constraint projection is applied, so that the dictionary can gradually absorb the drift and slow change information under normal operating conditions without destroying the existing mechanism semantics. At the same time, the dictionary update is frozen in time under high risk or high uncertainty conditions to avoid writing unstable samples, abnormal disturbance samples or low confidence samples into the dictionary, which would cause representation pollution and semantic shift.

[0065] Sixth, in the preferred implementation, step 2 of this application performs sliding short-time window slicing on the input data and extracts observational features such as frequency band energy distribution, cross-channel coupling relationship, and disturbance propagation consistency from short-time local dynamics, so that the original operating data is transformed from discrete and mixed multi-source time-series signals into a unified feature expression that can characterize the unstable evolution structure. Furthermore, by establishing a vortex structure dictionary constrained by the unstable mechanism grouping, and pre-dividing the dictionary primitives into rotating stall related primitive groups, separation enhancement primitive groups, system-level critical coupling primitive groups, and background and interpretable disturbance primitive groups, it is beneficial to improve the separability between different unstable precursor modes and reduce the confusing effect of background disturbances and operating condition changes on feature expression. Attached Figure Description

[0066] Figure 1 This is a flowchart of the method for identifying and preventing surge precursors in a centrifugal air compressor, as described in an embodiment of the present invention. Detailed Implementation

[0067] To enable those skilled in the art to better understand the technical solutions of this application, the following will provide a more detailed description of this application in conjunction with the accompanying drawings and embodiments.

[0068] The directional terms such as above, below, left, right, front, and back used in this application are based on the positional relationships shown in the attached drawings. Different attached drawings may result in different positional relationships, therefore they should not be interpreted as limitations on the scope of protection.

[0069] In this application, the terms "installation," "connection," "interlocking," "linking," and "fixing," etc., should be interpreted broadly. For example, they can refer to a fixed connection, a detachable connection, an integral connection, a mechanical connection, an electrical connection, or a connection that allows communication between components. They can also refer to a direct connection or an indirect connection through an intermediate medium. They can refer to the internal connection of two components or the interaction between two components. For those skilled in the art, the specific meaning of the above terms in this application can be understood according to the specific circumstances.

[0070] This invention describes a method for identifying and preventing surge precursors in centrifugal air compressors. It aims to address the shortcomings of existing technologies, such as reliance on single thresholds or black-box features for precursor identification, difficulty in balancing early detection, interpretability, and engineering stability, and the generally conservative intervention, significant bypass losses, and insufficient long-term adaptability in surge prevention control. This invention provides an integrated technical solution for precursor identification and control in a wide range of operating conditions for fuel cell air supply systems. The method constructs a sparse vortex structure coding mechanism to map multi-source unsteady observation signals during air compressor operation into low-dimensional sparse vortex structure state codes with mechanistic semantics. During the state code solution process, residual consistency constraints, cross-window consistency constraints, and a new activation suppression mechanism are introduced to improve the ability to identify the continuous evolution of instability precursors and reduce the activation of false precursors caused by noise disturbances, structural resonance, measurement point differences, and operating condition drift. Based on this, a precursor risk index is constructed by combining the state aggregation intensity, window growth characteristics and uncertainty correction information of different mechanism groups. This enables early identification, graded warning and reliable quantification of stall surge precursors. Furthermore, risk-driven minimum intervention adaptive anti-surge control is implemented based on the precursor risk index. Under the premise of meeting gas supply constraints and operational safety requirements, the anti-surge margin is dynamically maintained, reducing unnecessary bypass actions and the resulting compression power consumption and net output loss. At the same time, the stability of dictionary mechanism semantics and risk boundaries is maintained by online small-step updates and constraint projection mechanisms, improving the applicability and robustness of the method under long-term operating conditions.

[0071] As per the instruction manual Figure 1 This invention provides a method for identifying and preventing stall surge precursors in centrifugal air compressors. This method is applied to the air supply control process of fuel cell systems, particularly to the identification, risk assessment, and prevention of stall surge precursors in centrifugal air compressors and their air supply circuits under wide operating conditions and strong dynamic conditions. It aims to achieve early detection, graded warning, and minimal intervention control of instability precursors while meeting air supply requirements and operational safety constraints. This reduces unnecessary bypassing triggered by false alarms, decreases compression power consumption and air supply disturbances, and improves the operational stability, control adaptability, and net efficiency of the air supply system.

[0072] The method includes:

[0073] Step 1: Acquire dynamic operating data during the operation of the fuel cell centrifugal air compressor, and preprocess the dynamic operating data to obtain input data.

[0074] Step 1 includes:

[0075] Step 1.1: Collect dynamic operating data of the air compressor during operation by using multi-source sensors and a high-frequency multi-channel data acquisition system installed in the air supply circuit of the fuel cell centrifugal air compressor.

[0076] Taking a centrifugal air compressor air supply circuit used in a vehicle fuel cell system as an example, the process of acquiring dynamic operating data in step 1 is explained. The air supply circuit includes a centrifugal air compressor, an intake pipe, an outlet pipe, a bypass branch, a drive motor, and corresponding controllers. The centrifugal air compressor is used to supply compressed air to the fuel cell stack, and the bypass branch is used to unload and adjust the air supply circuit when the risk of stall and surge approaches.

[0077] In this embodiment, to obtain dynamic operating data that characterizes the precursors of stall and surge in a centrifugal air compressor, multi-source sensors are arranged on the air compressor body and its air supply circuit. Specifically, an inlet pressure sensor and an inlet temperature sensor are installed at the air compressor inlet pipe to collect intake air state parameters; an outlet pressure sensor and an outlet temperature sensor are installed at the air compressor outlet pipe or volute outlet to collect compressed air state parameters; a flow sensor is installed in the intake or outlet passage to directly measure air flow; for implementations without independent flow sensors, the equivalent flow can also be calculated using differential pressure, temperature, speed, and a calibration model. Simultaneously, a speed sensor, a motor current sensor, and a bus voltage sensor are installed on the drive motor side to collect the dynamic operating status of the air compressor's electric drive system; a bypass valve opening sensor is installed on the bypass branch to reflect the operating status of the anti-surge actuator. To improve the ability to capture early unsteady disturbances, in some embodiments, dynamic pressure sensors or vibration sensors may be added to the volute wall, the area near the diffuser, or the outer wall of key pipe sections to collect information on local pressure pulsations or structural vibrations.

[0078] In this embodiment, all the aforementioned sensors are connected to a high-frequency multi-channel data acquisition system, which uses a unified clock source for synchronous sampling. For signals that change rapidly, such as outlet pressure, volute dynamic pressure, motor current, and speed, a higher sampling frequency is used for acquisition; for signals that change relatively slowly, such as temperature and valve position, a lower sampling frequency is used, and these signals are aligned with the high-speed signals in subsequent processing. The data output from each acquisition channel includes timestamp information for subsequent timing synchronization and window slicing processing.

[0079] In this embodiment, the data acquisition process spans the entire process of the centrifugal air compressor, including startup, speed-up, steady-state operation, sudden load increases, sudden load decreases, bypass valve operation, and shutdown. For example, when the fuel cell vehicle accelerates, decelerates, and idles, the control system records parameters such as inlet pressure, outlet pressure, inlet temperature, outlet temperature, air flow rate, air compressor speed, motor current, bus voltage, and bypass valve opening in real time, thereby forming raw dynamic operating data covering different operating conditions. The raw data obtained thus includes not only pressure, temperature, and flow rate information characterizing the air path state, but also speed, current, voltage, and valve position information characterizing the drive state and actuator state, which can relatively completely reflect the dynamic behavior of the centrifugal air compressor under complex operating conditions.

[0080] In a specific application scenario, when a fuel cell system rapidly switches from low load to high load, the air compressor speed increases rapidly, and the outlet pressure and airflow change synchronously. The bypass valve opening and motor current also dynamically adjust according to control commands. Through the aforementioned multi-source sensors and high-frequency data acquisition system, the transient changes of various physical quantities during this process can be recorded simultaneously, providing a raw data foundation for subsequent steps such as stall surge precursor feature extraction, state coding, and risk assessment.

[0081] Step 1.2: Perform time synchronization, resampling alignment, detrending, normalization, noise reduction, and outlier suppression on the collected dynamic running data to obtain the input data.

[0082] After obtaining the raw dynamic running data in step 1.1, the dynamic running data is preprocessed to eliminate time offset, dimensional differences, noise interference and abnormal sampling point effects between different acquisition channels, so as to obtain input data that can be used by subsequent recognition models or judgment modules.

[0083] In this embodiment, the raw dynamic operating data of each acquisition channel is first synchronized and resampled for alignment. Specifically, since the signals such as inlet pressure, outlet pressure, air flow, air compressor speed, motor current, bus voltage, and bypass valve opening originate from different sensors and acquisition channels, there may be inconsistencies in sampling times, communication delays, or buffer lags between different channels. Therefore, the timestamps attached to the data of each channel are aligned using the controller's master clock or the unified time base of the data acquisition system as the reference time axis. For discrete sampling points that do not fall on the unified time axis, linear interpolation or zero-order hold is used for resampling to form multi-channel synchronous time-series data with a unified sampling interval. Through this processing, the multi-source asynchronous data obtained in step 1.1 is transformed into a time-consistent synchronous data sequence.

[0084] After time synchronization is completed, outlier identification and suppression are performed on the synchronized channel data. In this embodiment, outliers mainly include spikes, jumps, and missing values ​​caused by sensor transient distortion, sampling packet loss, communication jitter, and electromagnetic interference. For example, when a vehicle accelerates rapidly or a bypass valve suddenly actuates, some sampling points may exhibit transient pulses that significantly deviate from the normal trend. To address this, sliding window threshold discrimination, median filtering discrimination, or methods based on physical upper and lower limits can be used to identify outliers. For identified outliers, neighborhood interpolation, local mean replacement, or filling with the previous valid value can be used for correction. After this processing, the interference of outlier sampling points on subsequent feature extraction and state determination can be reduced.

[0085] Furthermore, to mitigate the impact of slow drift and differences in the dimensions of different physical quantities on the joint analysis, the synchronized data after outlier processing undergoes detrending and normalization. In this embodiment, detrending is primarily used to eliminate low-frequency baseline shifts caused by changes in ambient temperature, sensor zero-point drift, and slow transitions in system operating conditions. For example, for signals such as outlet pressure, airflow, and motor current, moving average detrending or high-pass filtering can be used to remove slowly changing baseline components, preserving dynamic fluctuation components related to stall and surge precursors. Normalization is mainly used to map different physical quantities to a unified numerical range, facilitating subsequent multivariate feature construction and joint determination. In this embodiment, parameters such as inlet pressure, outlet pressure, flow rate, rotational speed, current, voltage, and valve position can be normalized by scaling to rated values ​​or using Z-score normalization, ensuring comparability of variables with different dimensions in subsequent inputs.

[0086] In some implementations, the synchronized data after detrending and normalization is also subjected to noise reduction. For channels containing high-frequency measurement noise, such as outlet pressure, volute dynamic pressure, and motor current, low-pass filtering, moving average filtering, or wavelet denoising can be used to suppress random noise components. The filtering parameters are set based on the principle of preserving the dynamic fluctuation information corresponding to stall surge precursors, avoiding the masking of early weak disturbance characteristics due to excessive smoothing. For example, for the outlet pressure signal, while preserving low-frequency fluctuations and characteristic frequency band energy components, high-frequency measurement noise that significantly exceeds the target analysis frequency band is filtered out, thereby improving signal quality.

[0087] After the above preprocessing is completed, the multi-channel synchronous data is segmented. In this embodiment, according to the preset analysis window length and sliding step size, a sliding window method is used to extract data segments of consistent length from the continuous running data to form continuous input samples. Each input sample consists of a multi-dimensional time series corresponding to multiple physical quantities within the same time window. For example, during a rapid load increase phase, inlet pressure, outlet pressure, flow rate, rotational speed, current, and bypass valve opening data within that phase can be simultaneously extracted to form the input sample for the corresponding time window, which is used for subsequent short-time window feature construction and vortex structure state encoding.

[0088] After time synchronization, resampling alignment, outlier suppression, detrending, normalization, noise reduction, and sample segmentation, input data characterizing the dynamic operating state of the fuel cell centrifugal air compressor is obtained. This input data retains key dynamic information about stall surge precursors while reducing the impact of noise, drift, and channel inconsistencies on the identification results, thus providing a data foundation for subsequent steps such as constructing observation feature vectors, solving sparse vortex structure state codes, and assessing precursor risks.

[0089] It should be noted that the sensor type, installation location, sampling frequency, synchronization method and preprocessing method mentioned above do not constitute a limitation on this application. Those skilled in the art can make equivalent substitutions or adaptive adjustments according to the actual structure and control requirements of the fuel cell air supply system, and all of these should fall within the protection scope of this application.

[0090] Step 2: Based on the input data, construct the observation feature vector according to the sliding short time window, and establish a vortex structure dictionary constrained by the instability mechanism grouping. Divide the dictionary primitives in the vortex structure dictionary into multiple instability mechanism groups.

[0091] Step 2 includes:

[0092] Step 2.1: Slice the input data into a sliding short time window and extract observation features that characterize the energy distribution of frequency bands, cross-channel coupling relationships and the consistency of disturbance propagation, and construct an observation feature vector.

[0093] In this embodiment, the input data in step 1 is represented as a multi-channel time series composed of parameters such as inlet pressure, outlet pressure, inlet temperature, outlet temperature, air flow rate or equivalent flow rate, air compressor speed, motor current, bus voltage, and bypass valve opening. For the multi-channel time series, sliding short-time window slices are first performed according to a preset window length and window sliding step size to obtain continuous short-time dynamic segments. The window length is comprehensively set based on the air compressor dynamic response frequency band, the dominant frequency range of pressure pulsation, and the control cycle, so that each short-time window can cover a local dynamic process while maintaining the ability to distinguish the details of the precursor evolution; the window sliding step size is smaller than the window length, so that adjacent short-time windows maintain partial overlap, thereby facilitating continuous tracking of the evolution process of instability precursors. For example, in a specific implementation scenario, the preprocessed multi-channel input data can be sequentially truncated according to a fixed duration to form dynamic data segments corresponding to the first short-time window, the second short-time window, and so on up to the t-th short-time window.

[0094] For each data segment within a short time window, observational features used to characterize precursors of instability are further extracted. In this embodiment, the observational features include at least frequency band energy distribution features, cross-channel coupling relationship features, and disturbance propagation consistency features.

[0095] The frequency band energy distribution characteristics are used to describe the energy distribution of each dynamic signal across different frequency bands. Specifically, time-frequency analysis is performed on dynamic signals such as outlet pressure, volute dynamic pressure, air flow, air compressor speed, and motor current within a short time window to extract frequency band energy information related to instability precursors. Time-frequency analysis can be implemented using short-time Fourier transform, wavelet decomposition, or discrete wavelet packet decomposition. By calculating parameters such as spectral energy, energy proportion, energy concentration near the main peak frequency, and the ratio of low-frequency to mid-high-frequency energy for each signal in the low-frequency, mid-frequency, and high-frequency bands, the differences in frequency domain distribution of rotating stall precursors, localized separation enhancement, and system-level low-frequency oscillations can be reflected. For operating conditions with a significant dynamic enhancement trend, the growth rate of the dominant frequency band energy over time can be further extracted to enhance the ability to characterize the evolution trend of precursors.

[0096] Cross-channel coupling characteristics are used to describe the degree of synchronous change and mutual coupling between different physical quantities. Since fuel cell centrifugal air compressors often exhibit changes in the coupling relationships between variables such as pressure, flow rate, speed, current, and valve position when nearing instability, this embodiment calculates the correlation, coherence, and phase relationship characteristics of different signal channels within the same short time window. Specifically, this includes the peak cross-correlation between outlet pressure and air flow rate, the amplitude squared coherence coefficient between outlet pressure and speed, the average phase difference between motor current and bypass valve opening, the phase fluctuation amplitude of different channels, and phase consistency indices. These characteristics reflect the degree of coordinated change across multiple channels when local disturbances intensify, and help distinguish between structural changes caused by actual flow instability and non-instability fluctuations caused by actuator action and boundary disturbances.

[0097] In implementations with multiple spatially distributed measuring points, disturbance propagation consistency characteristics are further extracted. Specifically, when multiple pressure measuring points are arranged circumferentially around the volute, in the area adjacent to the diffuser, or at different pipe sections, characteristic parameters reflecting the disturbance propagation path and spatial organization can be constructed based on the time delay, phase propagation direction, propagation velocity consistency, and correlation strength between adjacent measuring points. For example, the estimated time delay, phase lag angle, and consistency index along the circumferential propagation between adjacent measuring points can be calculated to describe the continuous spatial propagation characteristics of rotational disturbances; the spatial organization of disturbance distribution during the expansion of the local separation zone can also be characterized based on the degree of consistency in the enhancement of pressure pulsations at different measuring points. These characteristics help enhance the ability to identify precursors related to rotational stall.

[0098] In a specific application scenario, when a fuel cell system rapidly increases from a medium load to a high load, the air compressor outlet pressure, air flow rate, and motor current change synchronously within a short period. By extracting features from short-time window data within this period, a set of information that simultaneously includes frequency band energy changes, enhanced inter-channel coupling, and changes in propagation consistency can be obtained. This transforms the original multi-source time-series signal into structured features that can characterize the current local dynamic state.

[0099] To improve the stability and comparability of features under different operating conditions, this embodiment can also perform secondary standardization on the various features extracted above. Specifically, the frequency band energy features, cross-channel coupling features, and propagation consistency features can be normalized, limited, or mapped to rated operating conditions to reduce the impact of sensor range differences and operating condition drift on the feature value range. If necessary, preset feature screening rules can be used to remove feature items with excessively low variance, excessive redundancy, or weak engineering interpretation significance, retaining only feature components that are highly sensitive to stall surge precursors.

[0100] After the above processing, the frequency band energy distribution features, cross-channel coupling relationship features, and disturbance propagation consistency features extracted within each short time window are concatenated in a predetermined order to form the observation feature vector corresponding to that short time window. Taking the t-th short time window as an example, the t-th observation feature vector is composed of all the features within that short time window. The observation feature vectors corresponding to multiple consecutive short time windows are arranged in chronological order to form the input feature sequence used in subsequent steps for vortex structure dictionary construction, status code solving, and precursor risk assessment.

[0101] Step 2.2: Based on the observed feature vectors, a dictionary learning method constrained by instability mechanism is used to establish a vortex structure dictionary. According to different instability mechanisms, the dictionary primitives are divided into rotating stall related primitive group, separation enhancement primitive group, system-level critical coupling primitive group, and background and interpretable disturbance primitive group.

[0102] Specifically, in this embodiment, a dictionary matrix is ​​first constructed based on the observed feature vectors corresponding to multiple short time windows. Where M is the total number of dictionary primitives, and each column Let represent a dictionary primitive. Let each observed feature vector... If all elements have the same dimension, the dictionary matrix consists of multiple dictionary primitives. Each dictionary primitive corresponds to a typical structural pattern in the feature space, used to characterize the observed feature distribution of a fuel cell centrifugal air compressor under a certain type of local flow state, coupled state, or disturbance state. For any short time window, the observed feature vector... All of them can be approximated by a linear combination of a small number of dictionary primitives, thereby transforming the complex changes in the original feature space into activation combinations of several typical structural patterns.

[0103] In one specific embodiment, the dictionary initialization process can utilize observed feature vectors from historical operating data under different operating conditions as candidate samples. For example, representative observed feature vectors can be selected from the steady-state air supply stage of the air compressor, the rapid load switching stage, the bypass valve intervention stage, the near-instability stage, and the stage where significant pressure pulsation has occurred. An initial dictionary is formed by extracting cluster centers, drawing prototype samples, or selecting random representative samples.

[0104] In this embodiment, to ensure that the established dictionary not only completes the sparse representation of observed feature vectors but also has a clear ability to explain instability mechanisms, instability mechanism grouping constraints are imposed on the dictionary primitives. Specifically, the dictionary primitives are divided into a rotational stall-related primitive group, a separation enhancement primitive group, a system-level critical coupling primitive group, and a background and interpretable disturbance primitive group.

[0105] The rotating stall-related primitive set is used to characterize the typical feature patterns corresponding to the enhancement of rotating disturbances and precursors of rotating stall. In this embodiment, this primitive set mainly absorbs the observation features corresponding to enhanced propagation consistency among multiple measurement points, ordered phase changes, enhanced circumferential propagation characteristics, and enhanced coherence structure in specific frequency bands. For example, when the propagation consistency index of multiple pressure measurement points, the estimated time delay of adjacent measurement points, and the coherence coefficient of the dominant frequency band extracted in step 2.1 show obvious ordered changes within certain short time windows, the dictionary primitive that best matches this type of observation feature is classified into the rotating stall-related primitive set. Through this primitive set, the formation and evolution of precursors of rotating instability can be structurally characterized.

[0106] The separation enhancement primitive group is used to characterize the typical feature patterns corresponding to enhanced separation flow in the diffuser region, volute region, or local flow channels. In this embodiment, this primitive group mainly reflects features such as broadband energy increase, enhanced mid-to-high frequency pulsation, weakened local channel coherence, and changes in coupling relationships. For example, when the energy proportion of the outlet pressure or volute pressure signal in the mid-to-high frequency band increases, and the coupling relationship between airflow and local pressure measurement points degrades, the dictionary primitives corresponding to this type of observation feature can be classified into the separation enhancement primitive group. This primitive group is used to characterize the precursor states of local separation zone expansion or enhanced non-uniform flow.

[0107] The system-level critical coupling primitive set is used to characterize the system-level low-frequency coupling modes formed by the combined effects of the air compressor body, inlet and outlet pipe volumes, bypass branches, load changes, and control loops. In this embodiment, this primitive set mainly corresponds to observation modes characterized by enhanced multi-channel low-frequency synchronous fluctuations, improved global coherence, and increased systemic coupling characteristics. For example, during rapid load changes, if the outlet pressure, air flow rate, air compressor speed, and motor current show a synchronous increasing trend in the low-frequency range, then the dictionary primitives matching this characteristic can be classified into the system-level critical coupling primitive set. Through this primitive set, the instability risk caused by system-level resonance or critical coupling effects can be characterized, rather than being limited to a certain local flow structure.

[0108] The background and interpretable disturbance primitive group is used to characterize background changes and interpretable disturbances that do not directly originate from the instability mechanism but affect the observed feature vector. In this embodiment, this primitive group mainly absorbs observed feature patterns caused by factors such as rapid bypass valve action, step changes in drive motor speed, switching of electric drive control modes, changes in ambient temperature, and differences in sensor installation. For example, when the bypass valve opening undergoes a significant transition within a short time window, causing transient disturbances in pressure and flow signals, the dictionary primitive corresponding to this type of disturbance is classified into the background and interpretable disturbance primitive group. The purpose of setting up this primitive group is to separate feature changes caused by non-instability factors from true instability precursors, thereby reducing misjudgments in the subsequent identification process.

[0109] Step 3: For the current short window, based on the observation feature vector, vortex structure dictionary and the status code corresponding to the previous short window, solve for the sparse vortex structure status code that satisfies the residual consistency, sparsity constraint and cross-window consistency constraint, and determine the uncertainty index based on the reconstruction residual of the current short window and / or the change in the status code of the adjacent short windows.

[0110] Step 3 includes:

[0111] Step 3.1: For the current short window, the observed feature vector corresponding to the current short window is solved by non-negative sparse coding with time-series prior constraints on the vortex structure dictionary. The goal is to minimize the feature reconstruction error, and simultaneously introduces sparsity constraints, cross-window amplitude continuity constraints with the previous short window state code, and suppression constraints on newly added activation components to obtain the sparse vortex structure state code of the current short window.

[0112] Specifically, for the t-th short-time window, the observation feature vector corresponding to the short-time window is denoted as... The vortex structure dictionary established in step 2.2 is denoted as... Record the current short-time window status code to be solved as M is the total number of dictionary primitives, and the status code corresponding to the previous short window is recorded as... Among them, status codes Each component corresponds to the activation intensity of each dictionary primitive in the vortex structure dictionary within the current short time window, and is used to characterize the degree of participation and combination relationship of rotating stall related modes, separation enhancement modes, system-level critical coupling modes, and background and interpretable disturbance modes at the current moment.

[0113] In this embodiment, to enable the observation feature vector of the current short window to be reconstructed from a small number of dictionary primitives, while ensuring the continuity and stability of state evolution between adjacent short windows, a system is constructed for state codes. The constrained optimization model aims to minimize the error between the current observed feature vector and the dictionary reconstruction result, while simultaneously introducing sparsity constraints, cross-window amplitude continuity constraints, and newly added activation suppression constraints. The status code is obtained by solving the following optimization problem:

[0114]

[0115] Where t represents the time sequence index of the current sliding short window; t-1 represents the previous short window adjacent to the current short window; This represents the observed feature vector corresponding to the t-th short time window; This represents the dictionary matrix of the vortex structure corresponding to the t-th short-time window; This represents the sparse vortex structure state code to be solved in the t-th short time window; This represents the optimal sparse vortex structure state code corresponding to the t-th short time window; This indicates the sparse vortex structure state code corresponding to the previous short time window; Let represent the objective function of the t-th short window, which is used to comprehensively characterize the reconstruction error, sparsity, cross-window consistency, and degree of suppression of new activations; W represents the feature domain weighting matrix, which is used to apply different weights to different feature components to reduce the impact of high-noise features, low-confidence features, or features with inconsistent dimensions on the solution results. This is a constraint term for residual consistency, used to measure the error between the current observed feature vector and the dictionary reconstruction result. The smaller this term is, the higher the degree of interpretation of the observed features by the current status code. This represents the sparse constraint weight coefficient, a scalar parameter greater than zero, used to adjust the balance between the sparsity of the status code and the reconstruction accuracy; The L1 norm sparse terms of the status code are used to ensure that only a small number of dictionary primitives are activated within the current short time window, thereby improving the mechanistic interpretability of the status code. γ represents the cross-window consistency constraint weight coefficient, which is a scalar parameter greater than zero. It is used to adjust the continuity constraint strength between the current state code and the previous short window state code. When γ increases, the state code change between adjacent short windows is more strongly suppressed. When γ decreases, the state code is more sensitive to rapid evolution precursors. The cross-window consistency constraint term is used to constrain the difference between the current short-term window status code and the previous short-term window status code, reducing irregular jumps caused by transient noise, sampling jitter, or local resonance; ρ represents the newly added activation suppression weight coefficient, which is a scalar parameter greater than zero, used to adjust the penalty intensity of the dictionary primitive component "the previous short-term window was not activated and the current short-term window is suddenly activated". This represents the support set mask vector constructed from the previous short-window state codes, whose dimension is the same as that of the state codes. Similarly, in some implementations, if the i-th component of the previous short-time window status code satisfies Then let Otherwise ,in The preset activation threshold; This represents a mask vector that is complementary to the support set of the previous short window, used to select the positions of dictionary primitives that were not activated in the previous short window. This represents element-wise multiplication, used to convert the current status code. Extract the portion belonging to the "newly activated component" from the data; The newly added activation suppression term is used to penalize activation components that were not activated in the previous short window but newly appeared in the current short window, in order to enhance the stability of the support set and reduce the probability of false precursor activation. The non-negative constraint requires that each component of the state code is not less than zero in order to ensure that it has a clear physical meaning. In this embodiment, each component of the state code can be understood as the activation intensity, energy ratio or participation degree of the corresponding dictionary primitive. Therefore, the non-negative constraint is adopted to avoid the occurrence of inexplicable negative activation phenomena. Represents the L2 norm; It represents the first norm.

[0116] The residual consistency term is used to constrain the observed eigenvectors. With dictionary reconstruction vector The fitting error between the two is minimized. By minimizing the reconstruction residual between them, the solved state code can accurately represent the structural mode corresponding to the current short window. The sparsity constraint term is used to ensure that the current state code activates only a few dictionary primitives. Since the dynamic characteristics of a centrifugal air compressor are usually determined by only a finite number of dominant structural modes within any short window, introducing sparsity regularization into the state code can avoid a large number of irrelevant primitives from participating in the reconstruction at the same time, thereby improving the structural identification and physical interpretation capabilities of the state code. The cross-window amplitude continuity constraint term is used to limit the change amplitude between the current short window state code and the previous short window state code. The newly added activation suppression constraint term is used to impose additional penalties on dictionary primitives that were not activated in the previous short window but suddenly appear in the current short window.

[0117] In one specific embodiment, the above solution process can be summarized as follows: under non-negative constraints, the optimal state code for the current short time window is obtained by comprehensively minimizing the weighted reconstruction error, sparse regularization term, cross-window difference term, and newly added activation suppression term. To meet the real-time requirements of online identification and control applications, this embodiment preferably uses an iterative optimization method suitable for embedded implementation, such as the proximal iteration method, coordinate descent method, alternating direction multiplier method, or equivalent variations thereof.

[0118] Furthermore, in this embodiment, the status code of the previous short time window can be... This serves as the initial value for the current short-time window solution, constructing a warm-start process. By using the solution result of the previous short-time window as the initial value for the current iteration, the iteration convergence time can be shortened, and the continuity of state evolution between adjacent short-time windows can be enhanced. During the specific iteration process, operations such as residual gradient update, sparse threshold shrinkage, cross-window difference correction, mask shrinkage of newly added activation components, and non-negative projection can be performed sequentially until a preset stopping condition is met. The iteration stopping condition may include at least one of the following: the change in state code between two consecutive iterations is lower than a set threshold. The objective function decreases by less than a set threshold. The number of iterations has reached the upper limit. Or the solution time for a single short window reaches the computation time limit allowed by the controller. ,in , , , The settings are determined by the controller's computing power budget and real-time requirements.

[0119] In a specific application scenario, when the load on a fuel cell system rapidly increases from a moderate level, and the observed feature vector extracted in step 2.1 simultaneously exhibits an increase in energy in a certain frequency band of outlet pressure, enhanced consistency in propagation between adjacent measuring points, and changes in the coupling relationship between rotational speed and current, the solution model can select a small number of best-matched dictionary primitives from the rotating stall-related primitive group and the system-level critical coupling primitive group for activation, while suppressing background interference primitives inconsistent with the current observation mode, thereby forming a state code with a sparse structure. If the activation intensity of such primitives shows a continuous upward trend within adjacent short time windows, the state code can continuously reflect this precursor evolution process without frequent jumps due to local noise interference.

[0120] Through the above solution process, the sparse vortex structure state code corresponding to the t-th short-time window is obtained in this embodiment. Status code It can not only represent the projection results of the current observed feature vector on each instability mechanism dictionary primitive in a sparse form, but also maintain good support set stability and amplitude continuity between adjacent short time windows, thus providing a state basis for subsequent steps such as uncertainty index calculation, group aggregation intensity extraction and precursor risk assessment.

[0121] Step 3.2: Based on the weighted reconstruction residual of the current short-term window observation feature vector and the sparse vortex structure state code, and combined with the support set change and code value difference intensity between the current short-term window state code and the previous short-term window state code, an uncertainty index is constructed to characterize the credibility of the current identification result.

[0122] Specifically, in this embodiment, the observation feature vector corresponding to the t-th short time window is denoted as... The current short-time window sparse vortex structure state code obtained in step 3.1 is denoted as... The status code corresponding to the previous short window is recorded as The vortex structure dictionary established in step 2.2 is denoted as Based on the above variables, the current short-term window status code is output. Simultaneously, the uncertainty index corresponding to this status code is further calculated. .

[0123] Uncertainty index It is used to reflect the degree of matching between the current short-term window observation features and the dictionary structure pattern, the consistency of the state evolution of adjacent short-term windows, and whether there are abnormal fluctuations in the current coding results, thereby providing a reliable basis for subsequent precursor risk assessment and control decisions.

[0124] In this embodiment, the uncertainty index It is constructed based on at least the following three types of information.

[0125] The first type of information is the residual information from the feature reconstruction of the current short-time window. Specifically, it is the residual information from the current observed feature vector. With dictionary matrix and current status code The reconstructed feature vector Compare the two and calculate the weighted residual energy. , If the current observed features can be reconstructed well from the existing dictionary primitives and their corresponding state codes, it indicates that there is a good matching relationship between the current observation pattern and the mechanism pattern represented by the dictionary. In this case, the state code... A high level of physical interpretation capability leads to a high degree of reliability in the corresponding identification results. Conversely, if the weighted reconstruction residual is too large, it may indicate the presence of abnormal perturbation patterns not fully covered by the dictionary, strong measurement noise interference, abrupt changes in boundary conditions, or a deviation from the existing mechanism pattern within the current short-term window, thus reducing the interpretability of the current state code for the observed features. In a specific embodiment, the feature domain weighting matrix W from step 3.1 can be used to weight the observed feature residuals, thereby constructing the weighted residual energy of the current short-term window. The higher the weighted residual energy, the higher the uncertainty of the current identification result.

[0126] The second type of information is the change in the state code support set between adjacent short-term windows. Since the precursors to stall surge in centrifugal air compressors typically exhibit a certain degree of continuous evolution, the state code support sets that truly demonstrate mechanistic continuity within adjacent short-term windows generally have a high degree of overlap. Based on this, in this embodiment, further adjustments are made based on the current short-term window state code... Compared with the previous short window status code Construct their respective support sets and ,in , = , and Let represent the sets of dictionary primitive indices activated in the current short window and the previous short window, respectively, and calculate the change between the support sets. Specifically, based on a preset activation threshold, the dictionary primitive positions corresponding to components in the status code above the activation threshold are defined as the support set for that short-term window. Then, one or more of the following are calculated: the intersection-union ratio (IU), the symmetry difference, the number of newly added components in the support set, and the number of components that have disappeared from the support set, to construct an index of support set change. If the overlap of the support sets between the current and previous short-term windows is high, it indicates that the state evolution is relatively stable and the current recognition result has good continuity. If the support set changes significantly between adjacent short-term windows, it may mean that the current status code is affected by noise, short-term external disturbances, or insufficient dictionary matching, and its reliability is relatively reduced.

[0127] The third type of information is the code value difference strength between the current short-term window state code and the previous short-term window state code. Even if the support set of the current short-term window and the previous short-term window have not changed significantly, if the activation amplitude of the corresponding dictionary primitive changes too quickly or too much between adjacent short-term windows, it may indicate that the current state has strong fluctuations, or that the current solution is sensitive to local noise and boundary changes. Therefore, in this embodiment, the current state code is further calculated. Compared to the previous status code The difference norm between them serves as a code value differential strength index characterizing the degree of change in state code amplitude, such as... This type of information can be used to supplement the judgment of cases where the support set remains unchanged but the code value fluctuates drastically, thus reflecting the stability of the current recognition result more comprehensively.

[0128] In one specific embodiment, the weighted reconstruction residual information, support set variation information, and code value differential strength information are normalized and mapped to a unified dimension interval. Subsequently, the three types of information are fused according to preset weights to obtain the uncertainty index corresponding to the current short-term window. The fusion method can be implemented using linear weighted summation, piecewise mapping, table lookup mapping, or rule combination judgment. For example, the weighted residual energy can be first... The support set change index and the code value difference strength index are converted into dimensionless scores between 0 and 1, and then different weights are assigned according to their impact on false alarm risk to obtain the comprehensive uncertainty index. In this embodiment, when the reconstruction residual increases, the change in the support set increases, or the code value differential strength increases, the uncertainty index... Simultaneous increase; when all three types of information are at low levels, the uncertainty index... Maintain a low level. The uncertainty index is:

[0129]

[0130] in, This represents the uncertainty index corresponding to the t-th short-time window; This represents the weighted reconstruction residual energy of the current short-term window; This indicates the status code corresponding to the t-th short time window. Support set; This indicates the status code corresponding to the (t-1)th short time window. Support set; △ represents symmetric difference operation; This represents the amount of change in the support set of adjacent short time windows; Indicates the strength of the code value difference between adjacent short-time window status codes; , , A weighting coefficient greater than 0 is used to adjust the degree of influence of each item on the uncertainty index.

[0131] In another specific embodiment, to improve the engineering applicability of the uncertainty index, a stability identification interval threshold, a cautionary judgment interval threshold, or a graded threshold can also be set for it. Specifically, when the uncertainty index When the status code is below the stable identification threshold, it can be determined that the current status code has a high degree of credibility, allowing it to directly participate in subsequent precursor risk assessments; when the uncertainty index When the uncertainty index falls between the stable identification threshold and the cautious judgment threshold, the current identification result can be marked as moderately confident, and a moderate confidence correction can be introduced in subsequent risk outputs; when the uncertainty index... When the threshold for cautious judgment is exceeded, it indicates that the current status code may be affected by strong noise, dictionary mismatch, or boundary abrupt changes, triggering a cautious judgment mechanism. This cautious judgment mechanism may include one or more of the following: lowering the precursor risk score for the current short-term window; delaying the output of a high-level warning conclusion; requiring multiple consecutive short-term windows to meet the precursor enhancement conditions before confirming the risk level; or marking the current identification result as pending confirmation. This approach reduces the probability of false alarms under conditions of strong noise, actuator abruptness, or rapid operating condition switching.

[0132] Through the above methods, step 3.2 achieves the quantification of the reliability of the current short-time window identification result, enabling the system to identify whether the current judgment is stable, whether there is a large interpretation bias, and whether a conservative decision needs to be made while outputting the sparse vortex structure state code. Compared with the method of only outputting a single state code, this implementation can better adapt to the precursor identification needs of fuel cell centrifugal air compressors under wide operating conditions, high dynamics, and weak observability. The sparse vortex structure state code obtained in step 3.1 and the uncertainty index obtained in step 3.2 can be further used for subsequent instability precursor intensity assessment, risk classification and early warning, and minimum intervention and surge control decision-making.

[0133] Step 4: Aggregate the sparse vortex structure state codes according to the group to which the dictionary primitives belong, obtain the aggregation intensity of each group, and construct the precursor risk index based on the aggregation intensity of each group, the change in aggregation intensity of adjacent short time windows, and the uncertainty index. Determine the risk level based on the comparison result between the precursor risk index and the preset risk threshold.

[0134] Step 4 includes:

[0135] Step 4.1: Group the dictionary primitives in the vortex structure dictionary according to the instability mechanism to which they belong, and perform grouping summation or weighted aggregation of the activation intensity of the corresponding dictionary primitives in the sparse vortex structure state code of the current short window to obtain the aggregation intensity of the rotating stall related primitive group, the separation enhancement primitive group, the system-level critical coupling primitive group, and the background and interpretable disturbance primitive group, and based on the change of the aggregation intensity of each group between the current short window and the previous short window between adjacent short windows.

[0136] Let the set of dictionary primitive indices corresponding to the g-th mechanism group be . Then, the status code components belonging to the group can be aggregated to obtain the overall response strength of the mechanism group within the current short time window. In this embodiment, the aggregated strength of the g-th mechanism group... It can be obtained by summing or weighted summing the status code components corresponding to each dictionary primitive in the group, specifically as follows:

[0137]

[0138] By aggregating each group, the aggregation strength of the rotational stall-related primitives can be obtained separately. Separation enhances the polymerization strength of the unit group The aggregation strength of system-level critical coupling units And the polymerization intensity of background and explainable interference units .in, Used to characterize the response intensity of rotational instability precursors within the current short time window. Used to characterize the activity level corresponding to localized separation enhancement. Used to characterize system-level low-frequency coupling or critical coupling trends. It is used to characterize the intensity of background interference introduced by non-target factors such as valve action, speed step, electric drive switching, measurement disturbance or environmental change.

[0139] After obtaining the polymerization intensity of each mechanism group, this embodiment further extracts its variation characteristics between adjacent short time windows to characterize whether the corresponding precursor is in a state of continuous enhancement. Specifically, for the g-th mechanism group, its polymerization intensity variation can be constructed based on the difference between the polymerization intensity corresponding to the current short time window and the previous short time window, expressed as:

[0140]

[0141] in, This characterizes the strengthening or weakening trend of the g-th mechanism group in the current short-term window relative to the previous short-term window. When the value is greater than zero, it indicates that the response intensity of this mechanism grouping shows an upward trend; when When the value is less than zero, it indicates that the activity level of this mechanism group has decreased; when When the value is close to zero, it indicates that the mechanism grouping remains basically stable between adjacent short time windows.

[0142] In a specific application scenario, when a fuel cell system rapidly transitions from medium to high load, if multiple state code components corresponding to the rotating stall-related primitives continuously increase within several short time windows, then the resulting data after grouping and aggregation... It will gradually increase, and It remains positive over several consecutive short time windows; conversely, if only the aggregation intensity of the background and the explainable interfering units exists in a short time window... A brief increase followed by a rapid recovery within an adjacent short time window suggests that the change is more likely to correspond to a transient disturbance caused by bypass valve action or control switching, rather than a sustained increase in a true precursor to instability.

[0143] Step 4.2: Based on the aggregation intensity of each instability mechanism group, the change in aggregation intensity of each group between adjacent short time windows, and the uncertainty index, construct a precursor risk index that jointly characterizes the precursors of rotating stall, separation enhancement, and system-level critical coupling, and performs confidence correction on the background disturbance. Then, compare the precursor risk index with at least one preset risk threshold to determine the risk level corresponding to the current short time window.

[0144] Specifically, in this embodiment, the polymerization strength of the rotational stall-related primitives obtained in step 4.1 is first determined. Separation enhances the polymerization strength of the unit group The aggregation strength of system-level critical coupling units And the polymerization intensity of background and explainable interference units And combined with the corresponding change in polymerization strength and Each component risk quantity corresponding to a different instability mechanism is constructed.

[0145] In one specific embodiment, the risk component of rotational stall precursors can be determined by the aggregation strength of rotational stall-related primitives. and its change Jointly determined; the precursor risk component of separation enhancement can be determined by the polymerization strength of the separation enhancement unit set. and its change Jointly determined; the precursory risk component of system-level critical coupling can be determined by the aggregation strength of system-level critical coupling primitives. and its change The risk is determined jointly. In this embodiment, the risk quantities of each sub-item can be constructed using a linear weighting method. For example, for a certain mechanism group, its sub-item risk can be obtained by linearly combining the aggregation intensity and change amount of that group according to a preset weight; alternatively, piecewise mapping or monotonic nonlinear mapping can be used as needed to convert the aggregation intensity and change amount into corresponding risk scores. Preferably, before the aggregation intensity and change amount enter the fusion calculation, they are first normalized based on historical steady-state statistical intervals, rated operating condition reference values, or offline calibration samples to improve the comparability of risk quantities under different operating boundary conditions.

[0146] After obtaining the above-mentioned risk quantities, the precursors of rotating stall, enhanced separation, and system-level critical coupling are further combined and fused to construct a unified precursor risk index corresponding to the current short time window. In one specific embodiment, the precursor risk index can be constructed using the following linearly modified form:

[0147]

[0148] in, It is a monotonic mapping function used to map the fusion result to a predetermined risk range; This represents the growth trend of the polymerization intensity of the rotational stall-related primitives over time. This represents the growth trend of the aggregation intensity of the critically coupled primitive group at the system level over time. , , The weighting coefficients corresponding to the aggregation intensity of each instability mechanism are used to characterize the contribution of the current activity level of different mechanisms to the total risk; , The weighting coefficients corresponding to the growth trend items; This is the correction weight for uncertainty.

[0149] Using the above construction method, the precursor risk index This comprehensively reflects the response intensity of spin-stall precursors, enhanced separation precursors, and system-level critical coupling precursors within the current short time window. Furthermore, spin-stall precursors and system-level critical coupling precursors participate in risk correction through their growth trend terms to characterize whether the corresponding mechanisms are in a state of continuous enhancement. Simultaneously, an uncertainty index is introduced. The credibility correction of the fusion results suppresses the precursor risk index when the current short-window identification results are unstable, the residuals are large, or the cross-window consistency is poor, thereby reducing the probability of false alarms.

[0150] Based on the early warning risk index of the current short-term window Then, the precursor risk index is compared with at least one preset risk threshold to determine the risk level corresponding to the current short-term window. In this embodiment, multiple risk classification thresholds can be set, for example... , and and satisfy .according to The current threshold range can be used to classify the current operating status into different risk levels. For example, when... When the air compressor is in normal operation or without obvious warning signs, it is determined that the air compressor is in normal operation or without obvious warning signs. When, it is determined to be a precursor state; when When, it is determined to be a precursor enhancement state; when When this occurs, it is determined to be a state nearing instability or a high-risk state. The above risk level classification is only an example; in other implementations, two, three, or more risk levels can be set according to specific application requirements.

[0151] To avoid frequent risk level switching caused by fluctuations in the current short-term risk index around a threshold, some preferred embodiments may introduce a hysteresis comparison mechanism and a continuous window confirmation mechanism for risk level determination. The hysteresis comparison mechanism uses a higher threshold for entering a higher risk level than for exiting that risk level; the continuous window confirmation mechanism confirms entry into the corresponding risk level only when the current precursor risk index exceeds a certain risk threshold for several consecutive short-term windows. These mechanisms improve the stability of risk classification results and prevent misjudgments caused by single-window spikes or short-term disturbances.

[0152] After the above processing, step 4.2 outputs the precursor risk index and risk level corresponding to the current short-term window. In some implementations, additional information such as the dominant mechanism category corresponding to the current risk level, the risk growth trend, and whether it is necessary to enter the early warning confirmation state can also be output simultaneously, so that the controller can execute graded anti-surge control, bypass intervention optimization, or monitoring and alarm strategies accordingly.

[0153] After processing in steps 4.1 and 4.2, the short-time window level sparse eddy structure state code can be converted into a precursor risk index and risk level with clear mechanistic meaning, providing a unified decision-making basis for subsequent minimum intervention prevention and asthma control based on risk level.

[0154] Step 5: Based on the precursor risk index and its changing trend, under the conditions of meeting the gas supply demand constraints and actuator constraints, generate continuous control commands for the bypass valve opening, selectively generate continuous control commands for motor torque and / or guide vane angle, and when the precursor risk index decreases to the preset recovery condition, perform gradual recovery of the control commands to restore the air compressor to a high-efficiency operating state.

[0155] Step 5 includes:

[0156] Step 5.1: Determine the anti-surge control target based on the fuel cell air supply demand, and establish a set of control constraints including target cathode pressure and / or target air flow, excess air coefficient constraints, and amplitude and rate of change constraints of bypass valves, motors and / or guide vane actuators.

[0157] Specifically, in this embodiment, the air supply control target is first determined based on the current load demand of the fuel cell system. The control target is a hard constraint to meet the fuel cell's air supply demand, including at least the tracking requirements of the target cathode pressure and / or target air flow rate, and may further include constraints on the excess air coefficient. Simultaneously, the optimization direction is to reduce bypass losses and compression power consumption, thereby improving the operating efficiency of the air supply system while ensuring air supply safety and stable compressor operation.

[0158] In one specific embodiment, the target cathode pressure within the current control cycle is determined based on the current stack current request, power request, ambient pressure, ambient temperature, and stack operating status. and / or target airflow The target cathode pressure and target air flow rate can be determined based on the stack current request, power request, ambient pressure, ambient temperature, stack operating point, and system calibration mapping relationship. When the system adopts a pressure-driven air management mode, the target cathode pressure is the primary tracking factor; when the system adopts a flow-driven air management mode, the target air flow rate is the primary tracking factor; when the system simultaneously considers cathode pressure and flow rate distribution, both are incorporated into the control objectives.

[0159] In this embodiment, the control input includes at least the bypass valve opening degree. Optionally includes motor torque With guide vane angle The bypass valve opening is used to adjust the degree of opening of the compressed air return channel to release pressure buildup when it is close to instability; the motor torque is used to adjust the speed change trend of the air compressor, thereby affecting the air supply capacity and operating point position; the guide vane angle is used to adjust the intake aerodynamic matching state to improve the compressor's operating characteristics under different operating conditions.

[0160] After determining the control objective and control input, a set of control constraints is further established. This set of constraints includes at least actuator amplitude constraints and rate of change constraints. Let the control input vector be denoted as u, then the actuator amplitude constraints satisfy:

[0161]

[0162] The rate of change constraint of the actuator satisfies:

[0163]

[0164] in, and These represent the lower and upper limits of the control input, respectively; This indicates the change in control input at adjacent control moments; This indicates the maximum permissible variation in the control input. These constraints ensure that the control commands for the bypass valve, motor, and guide vane actuator remain within their actual executable range and prevent rapid abrupt changes in control inputs within a short period, thereby reducing the risks of actuator shock, gas supply fluctuations, and control instability.

[0165] In one specific embodiment, when the control input vector is composed of the bypass valve opening, motor torque, and guide vane angle, it can be: .at this time, The amplitude constraints can be respectively corresponding to the upper and lower limits of the bypass valve opening, the upper and lower limits of the motor torque, and the upper and lower limits of the guide vane angle. The rate-of-change constraints correspond to the limits on the changes in bypass valve opening, motor torque, and guide vane angle within adjacent control cycles, respectively. This approach allows for a unified constraint expression for coordinated control of multiple actuators.

[0166] Furthermore, the set of control constraints corresponding to the t-th control time is:

[0167]

[0168] in, This represents the control input vector at the t-th control time. The control input vector includes at least the bypass valve opening. Optionally includes motor torque and / or guide vane angle ; This represents the actual cathode pressure at the t-th control moment; This represents the target cathode pressure at the t-th control moment; This indicates that the cathode pressure tracking error does not exceed the preset pressure error upper limit. ; This represents the actual airflow rate at the t-th control time. This represents the target airflow rate at the t-th control time. This indicates that the airflow tracking error does not exceed the preset flow error limit. ; This represents the actual excess air coefficient at control time t; and These represent the lower and upper limits of the allowable excess air coefficient at the t-th control time, respectively. and These represent the control input vectors respectively. The lower and upper amplitude limits are used to define the physical reach of bypass valves, motors, and guide vane actuators; This represents the change in the control input vector at adjacent control moments; This indicates the maximum allowable variation in the control input vector.

[0169] Step 5.2: Based on the precursor risk index and its adjacent short-term window change trends, generate continuous control commands for bypass valve opening that satisfy the control constraint set, and selectively coordinate the generation of continuous control commands for motor torque and / or guide vane angle according to the risk level, so as to achieve minimum intervention and surge prevention control according to the risk increment.

[0170] To achieve minimum intervention and surge prevention control that involves "small-scale action in the early stage, rapid increase in margin in the critical stage, and smooth recovery after risk is eliminated," this embodiment generates continuous control commands for bypass valve opening that satisfy the set of control constraints based on the early risk index and its adjacent short-term window changes. Furthermore, it selectively coordinates and generates continuous control commands for motor torque and / or guide vane angle according to the risk level, thereby achieving minimum intervention and surge prevention control with increasing risk.

[0171] Specifically, in this embodiment, the bypass valve is first used as the main rapid intervention actuator in the anti-surge control. The bypass valve opening control is used to adjust the degree of opening of the compressed air return path, so as to quickly release pressure buildup, increase flow margin, and drive the operating point away from the stall surge region when the air compressor operating point is close to the instability boundary. To balance risk response capability and air supply efficiency, in this embodiment, the bypass valve opening command is based on the current short-term window precursor risk index. The changing trends between the window and its adjacent short time windows are continuously generated, rather than using a single threshold-triggered full-open or sudden-open control method.

[0172] In one specific embodiment, the bypass valve opening command for:

[0173]

[0174] in, This indicates the bypass valve opening command at the t-th control moment; Indicates the basic opening degree; This indicates the precursory risk index corresponding to the current short-term window; This indicates the precursor risk index corresponding to the previous short time window; This represents the static response gain of the risk index. Gain in response to changes in the risk index trend; The saturation function is used to ensure that the generated bypass valve opening command meets the actuator amplitude constraint and rate of change constraint established in step 5.1.

[0175] With the above construction, the bypass valve opening command... In the formula This item is used to adjust the bypass valve opening based on the current risk level; that is, when the current risk index increases, the bypass valve opening increases accordingly. This term is used to dynamically adjust the bypass valve opening based on the risk change trend. Specifically, when the risk index continues to rise between adjacent short time windows, the bypass valve opening growth rate is further increased to enhance the response capability to scenarios with rapidly rising risks. When the risk index decreases or tends to level off, this trend term decreases or becomes negative, so that the bypass valve opening can be smoothly recovered, avoiding excessive bypass losses after the risk is eliminated.

[0176] In some implementations, the basic opening The default bypass valve opening can be set under normal operating conditions, or it can be dynamically set based on the air supply target determined in step 5.1, the current load condition, or environmental boundary conditions. For low-risk conditions, A smaller value or a value close to the closure value can be selected to reduce unnecessary bypasses; for operating conditions that have already entered the precursor enhancement stage, It can also be appropriately increased to reserve a shorter execution path for subsequent rapid capacity increases.

[0177] In this embodiment, the saturation function This not only ensures that the bypass valve opening command does not exceed the upper and lower limits of the physical stroke, but also preferably satisfies the rate of change constraint in step 5.1. Specifically, when the bypass valve opening command... When the target bypass valve opening calculated by the formula exceeds the allowable range, it is truncated to the corresponding boundary; when the opening change between adjacent control moments exceeds the upper limit of the allowable rate of change, a limit correction is performed according to the upper limit of the rate of change. Through the above processing, it can be ensured that the bypass valve action can effectively respond to the increase in risk, and will not cause mechanical shock, control instability or large disturbances in gas supply due to sudden changes in commands.

[0178] In addition to the bypass valve opening command, in this embodiment, continuous motor torque control commands and / or continuous guide vane angle control commands can be selectively generated based on the risk level. Motor torque control is used to adjust the compressor speed variation trend, thereby affecting the compressor's air supply capacity and operating point position; guide vane angle control is used to adjust the intake aerodynamic matching state, thereby improving the compressor's aerodynamic stability under different operating conditions. Compared to the bypass valve, the adjustment of motor torque and guide vane angle is more biased towards operational state shaping and efficiency coordination. Therefore, in this embodiment, selective coordination is preferably performed according to the principle of "prioritizing efficiency without compromising efficiency."

[0179] Specifically, in one embodiment, when the precursor risk level is in the low-risk range or the initial stage of the precursor, small-amplitude motor torque regulation and / or guide vane angle fine-tuning are preferentially used to slightly correct the compressor operating point while minimizing bypass losses. For example, the tendency for the motor to accelerate further can be moderately reduced, or the guide vane angle can be slightly adjusted to keep the compressor operating point in a region more conducive to stable operation, thereby reserving a certain anti-surge margin at a small efficiency cost. Within this risk range, the bypass valve opening can still be set according to the bypass valve opening command. Formulas are generated, but their increments are typically kept small to reflect the principle of minimal intervention.

[0180] When step 4.2 determines that the current risk level has entered the precursor enhancement zone, this embodiment implements coordinated control of the bypass valve, motor, and guide vane actuator. At this time, the bypass valve opening degree is determined according to the bypass valve opening degree command. The formula continuously increases, while simultaneously applying correction commands to the motor torque to suppress further approach to the instability boundary, and / or applying compensation commands to the guide vane angle that are beneficial to aerodynamic stability. In this way, the bypass valve is used to quickly increase the margin, and the motor and guide vane adjustments are used to slow down the operating point from continuing to approach the high-risk region, thereby reducing the continuous efficiency loss caused by complete reliance on bypass.

[0181] When the current risk level enters the high-risk zone or approaches the instability zone, this embodiment prioritizes increasing the control weight of the bypass valve, making it the dominant intervention channel to achieve rapid escape from the critical zone. Simultaneously, adjustments to the motor torque and / or guide vane angle can be supplemented based on system capacity and response speed, but efficiency is no longer the primary objective; instead, the priority is to ensure the air compressor quickly moves away from the stall and surge boundary. In one specific embodiment, when When the high-risk threshold is exceeded, the bypass valve opening command can be adjusted. Gain in the formula and / or By increasing the bypass valve opening in stages, the valve opening becomes more sensitive to the risk index and the upward trend of risk, thereby enabling rapid margin control during high-risk phases.

[0182] In other embodiments, the continuous control commands for motor torque and guide vane angle can also be generated in a risk-driven manner. For example, the motor torque correction amount Δ can be constructed based on the current risk level, the trend of risk index changes, and the control objective in step 5.1. And guide vane angle correction Δ This allows the system to coordinate and adjust the target airflow, cathode pressure, and aerodynamic stability margin while satisfying their respective amplitude and rate of change constraints. For systems without guide vane actuators, only bypass valve opening commands and motor torque commands are output; for systems without independent motor torque coordination, only bypass valve opening commands and guide vane angle commands are output.

[0183] In a specific application scenario, when the load on the fuel cell system rapidly increases from a moderate level, and the precursor risk index calculated in step 4.2... When multiple short time windows rise consecutively but have not yet entered the high-risk zone, the bypass valve opening command can be used as a guide. The formula applies small, continuous opening commands to the bypass valve and simultaneously suppresses further rapid increases in motor torque to slow the rate at which the operating point approaches the instability boundary; if subsequently... If the risk level continues to rise and exceeds the high-risk threshold, the bypass valve opening will increase rapidly under the combined effect of the risk and trend factors, causing the air compressor operating point to quickly move away from the critical zone; as the subsequent risk index gradually declines, the bypass valve opening command will... The trend term in the formula decreases accordingly or even turns negative, thereby driving the bypass valve opening to gradually recover at a limited rate, avoiding continuous energy loss caused by the bypass channel maintaining an excessively large opening for a long time.

[0184] In this embodiment, by directly incorporating the precursor risk index and its adjacent short-term window variation trends into the continuous control law, gradual intervention at different risk stages can be achieved: small-amplitude, continuous, and low-disturbance control is implemented at the initial stage of precursor occurrence; control intensity is promptly increased during the enhancement and high-risk stages of precursor occurrence; and the intervention state is smoothly exited after the risk is resolved. Compared with traditional single-threshold triggered bypass anti-surge control, this method is beneficial in reducing compression work loss and gas supply disturbance caused by excessive conservative bypass.

[0185] Through the above processing, in this embodiment, step 5.2 outputs a continuous control command for the bypass valve opening that satisfies the control constraint set described in step 5.1, and selectively outputs continuous control commands for motor torque and / or guide vane angle according to the risk level, thereby forming a risk-driven minimum anti-surge control quantity for the fuel cell centrifugal air compressor, providing direct input for subsequent online control of the actuator.

[0186] Step 5.3: When the precursor risk index decreases below the preset recovery threshold and the preset stability conditions are continuously met, perform a gradual recovery of the continuous control commands for the bypass valve opening and the continuous control commands for the motor torque and / or guide vane angle, so as to restore the air compressor to a high-efficiency operating state.

[0187] Specifically, in this embodiment, the precursor risk index output in step 4.2 is first... Perform a recovery criterion assessment. When Decreased and stabilized below the recovery threshold Instead of immediately withdrawing all the asthma prevention intervention control measures applied in step 5.2, the condition is further required to continue to be met for a preset number of short-term windows to constitute a recovery confirmation condition. Let the number of continuous short-term windows be... Then when continuous Each short time window meets the requirements. When this is the case, it can be determined that the air compressor has moved away from the precursor enhancement or near-instability state and entered the recoverable stage.

[0188] In some implementations, the preset stabilization conditions include, in addition to the precursor risk index being below the recovery threshold and remaining thereafter... In addition to the short-term window, one or more of the following conditions may be further included: the precursor risk index no longer rises continuously between adjacent short-term windows; the uncertainty index obtained in step 3.2 The risk level falls below the preset stability threshold; the change in polymerization intensity of the rotating stall-related primitive group and the system-level critical coupling primitive group in step 4.1 no longer shows a continuous positive increase; or the current risk level has dropped to the normal operation state or the initial stage of the precursor state. By setting the above recovery criteria, the false recovery phenomenon caused by short-term fluctuations near the risk threshold can be reduced.

[0189] After confirming the entry into the recoverable phase, this embodiment performs a gradual recovery of the continuous control command for the bypass valve opening. Specifically, instead of directly restoring the bypass valve opening to the basic opening or closed state, the bypass valve opening is gradually reduced according to a preset recovery slope, allowing the air compressor operating point to smoothly return to the high-efficiency operating region. Let the current bypass valve opening command output in step 5.2 be... Restore the target opening degree to In a specific embodiment, the bypass valve opening recovery process can be implemented recursively as follows:

[0190]

[0191] in, The bypass valve opening and recovery slope; This is the time interval between adjacent control moments. The recovery process still satisfies the actuator amplitude and rate of change constraints established in step 5.1 to avoid gas supply fluctuations or system shocks caused by rapid shut-off of the bypass valve.

[0192] In one specific embodiment, let the current motor torque control command be: The corresponding high-efficiency reference torque is Assume the current control command for the guide vane angle is... The corresponding high-efficiency reference guide vane angle is During the recoverable phase, updates can be performed using the constrained asymptotic recycling method:

[0193] → , →

[0194] Here, "→" indicates that, under the premise of satisfying the amplitude and rate of change constraints in step 5.1, the corresponding high-efficiency reference value is gradually approached with a finite step size, rather than jumping to the target value all at once. In this way, new gas supply disturbances or mechanical shocks can be avoided during the recovery process.

[0195] Through the above processing, in this embodiment, when step 5.3 confirms that the current precursor risk index has decreased to below the preset recovery threshold and continues to meet the preset stability conditions, the continuous control command for the bypass valve opening and the continuous control command for the motor torque and / or guide vane angle are gradually recovered, so that the air compressor smoothly recovers from the anti-surge protection state to the high-efficiency operating state, thereby completing the closed-loop recovery process of risk-driven minimum interference anti-surge control.

[0196] Step 6: Based on the precursor risk index and uncertainty, when the preset update conditions are met, the vortex structure dictionary is updated in small steps. When the preset freeze conditions are met, the update of the vortex structure dictionary is frozen and the current dictionary is maintained to participate in the solution of sparse vortex structure state codes.

[0197] Step 6 includes:

[0198] Step 6.1: Based on the precursor risk index and uncertainty, determine whether the vortex structure dictionary corresponding to the current short time window is in update mode or frozen mode. Specifically, when the precursor risk index is in the preset low-risk range and the uncertainty is lower than the preset threshold, the vortex structure dictionary is determined to be in update mode; when the precursor risk index is in the preset high-risk range and / or the uncertainty is higher than the preset threshold, the vortex structure dictionary is determined to be in frozen mode, and the current dictionary is maintained to participate in the sparse vortex structure state code solution.

[0199] Specifically, in this embodiment, let the precursor risk index corresponding to the t-th short-term window be . The uncertainty index is The vortex structure dictionary currently used in solving the sparse vortex structure state codes in step 3.1 is: To determine whether online dictionary updates are permitted within the current short-term window, pre-defined dictionary operating mode rules are established. The operating modes include at least an update mode and a freeze mode. The update mode allows the dictionary to perform small-step adaptive corrections based on reliable samples within the current short-term window, adapting to slow system drift, changes in actuator characteristics, and environmental boundary variations. The freeze mode stops dictionary updates when the risk is high or the recognition result lacks sufficient reliability, maintaining only the current dictionary to continue participating in the status code solving of subsequent short-term windows, thereby preserving the semantic stability of the dictionary mechanism.

[0200] In one specific embodiment, a low-risk threshold for triggering dictionary updates is preset. and low uncertainty threshold and high-risk thresholds used to trigger dictionary freezing. and high uncertainty threshold Among them, the threshold satisfies , This creates a decision margin between update mode and freeze mode, avoiding frequent switching of the working mode around a single threshold.

[0201] When the current short-term window satisfies the condition that the precursor risk index is in a preset low-risk range and the uncertainty is below a preset threshold, the current vortex structure dictionary is determined to be in update mode. Specifically, in some implementations, when the condition is met... , In this case, the observed features corresponding to the current short-term window are considered to mainly reflect the steady-state operation, slight disturbance, or low-risk interpretable state. At this point, the current sample shows high consistency with the existing mechanism dictionary, and the identification results are highly reliable. Therefore, it is permissible to use the current short-term window samples to perform subsequent small-step updates to the vortex structure dictionary. In this way, the dictionary can gradually absorb slowly varying drift information of the system under normal and low-risk operating conditions, improving its long-term adaptability.

[0202] When the current short-term window satisfies the condition that the precursor risk index is in a preset high-risk range and / or the uncertainty is higher than a preset threshold, the current vortex structure dictionary is determined to be in frozen mode. Specifically, in some implementations, when the condition is met... and / or At this point, it is believed that the current system may be in the stage of enhanced instability precursors, approaching the critical instability zone, or, although the risk index has not yet increased significantly, the reliability of the current identification results has decreased. At this time, the currently observed features may contain abnormal perturbation patterns, samples with rapid boundary switching, samples with strong actuator actions, or high-noise samples. If the dictionary continues to be updated online at this stage, abnormal samples may be written into the dictionary, causing a shift in the mechanistic semantics of the dictionary primitives. Therefore, in freeze mode, the current dictionary is not updated; only the current dictionary is maintained. Continue to participate in the sparse vortex structure state code solution in step 3.1 to maintain the stability of the precursor identification model structure.

[0203] In this embodiment, when the current short-term window does not meet either the explicit update condition or the explicit freeze condition, the dictionary operating mode corresponding to the previous short-term window can be maintained unchanged. That is, a mode-maintaining interval is set between the update threshold and the freeze threshold, allowing the system to continue the original operating mode during intermediate transition states, thereby reducing the frequency of mode switching. In this way, the continuity and robustness of the online dictionary management process can be further enhanced.

[0204] Through the above processing, in this embodiment, step 6.1 outputs the dictionary working mode corresponding to the current short-term window. If the current short-term window corresponds to the update mode, then proceed to step 6.2 to perform dictionary update; if the current short-term window corresponds to the freeze mode, then maintain the current dictionary unchanged and continue to use it for solving the sparse vortex structure state codes of subsequent short-term windows. In this way, while ensuring the semantic stability of the dictionary mechanism, the adaptive ability and engineering robustness of the method in long-term operation can be improved.

[0205] Step 6.2: In update mode, based on the observed feature vector and sparse vortex structure state code corresponding to the current short time window, perform incremental update of the vortex structure dictionary with step size limit and forgetting factor, and apply instability mechanism grouping constraint projection to the dictionary primitives after the update.

[0206] Specifically, in this embodiment, if step 6.1 determines that the vortex structure dictionary corresponding to the t-th short-term window is in update mode, then the observed feature vector corresponding to the current short-term window is used. and the sparse vortex structure state code obtained from step 3.1 Based on the current dictionary Perform a restricted incremental update to obtain the updated dictionary. If step 6.1 determines that the vortex structure dictionary corresponding to the current short-term window is in frozen mode, then this step is not executed, and the current dictionary is maintained. constant.

[0207] In one specific embodiment, let the current vortex structure dictionary be: The corresponding sparse vortex structure state code is First, the reconstruction result of the current short-term window is calculated based on the current dictionary and the sparse vortex structure state code. and reconstructed residuals Then, based on the activation strength of each component in the status code, select the one that satisfies the condition. The dictionary primitives are the active primitives in this round, among which To preset the activation threshold, inactive dictionary primitives remain unchanged. For each active primitive... Construct the original update amount according to its corresponding status code amplitude and the reconstructed residual:

[0208]

[0209] in, The update step size is based on this. To prevent dictionary mutations caused by a single short-time-window sample, a further step size constraint is imposed on the original update amount. At that time, it was truncated as:

[0210]

[0211] Otherwise take ,in The maximum single-step update magnitude is preset. Then, a forgetting factor is introduced. By fusing historical dictionary information with current incremental information, candidate update primitives are obtained:

[0212]

[0213] in , A smaller value indicates a greater emphasis on the stability of the historical dictionary. The larger the value, the greater the emphasis on adaptability to current reliable samples.

[0214] In this embodiment, to maintain the semantic stability of the dictionary mechanism, the updated candidate primitives are... Instead of directly writing it into the dictionary, the primitive is projected back to the corresponding constraint set based on the instability mechanism group to which it belongs. Specifically, if the i-th dictionary primitive belongs to the rotating stall-related primitive group, it is projected onto the constraint set reflecting enhanced propagation consistency and enhanced cross-channel coherent structure; if it belongs to the diffuser / volute separation enhancement primitive group, it is projected onto the constraint set reflecting broadband energy rise and coherence degradation trends; if it belongs to the system-level critical coupling primitive group, it is projected onto the constraint set reflecting low-frequency global consistency; if it belongs to the background and interpretable disturbance primitive group, it is projected onto the constraint set aligned with bypass valve action, torque step, or speed regulation time. The projected primitive is denoted as:

[0215]

[0216] in, This represents the projection operator onto the constraint set of the mechanism group to which the i-th primitive belongs. After updating all active primitives, the vortex structure dictionary for the new time step is obtained:

[0217]

[0218] This is then used to solve the sparse vortex structure state code for subsequent short-term windows; if step 6.1 determines that the current short-term window is in frozen mode, the above update is not performed, and the current state code remains unchanged. The dictionary remains unchanged, and the encoding reasoning continues using only the current dictionary, thereby avoiding the contamination of the dictionary by high-risk or high-uncertainty samples.

[0219] By using steps 6.1 and 6.2, the vortex structure dictionary can be made to perform controlled small-step adaptive updates only under low-risk and low-uncertainty conditions, and remain frozen under high-risk or high-uncertainty conditions. This avoids critical samples from contaminating the dictionary, while improving the environmental adaptability, mechanism interpretability, and cross-condition robustness of the vortex structure dictionary in long-term operation.

[0220] The above descriptions are merely embodiments of this application, and common knowledge regarding specific structures and characteristics in the solutions is not described in detail here. It will be apparent to those skilled in the art that this application is not limited to the details of the above exemplary embodiments, and that this application can be implemented in other specific forms without departing from the spirit or essential characteristics of this application. Therefore, the embodiments should be considered exemplary and non-limiting in all respects, and the scope of this application is defined by the appended claims rather than the foregoing description. Therefore, it is intended that all variations falling within the meaning and scope of equivalents of the claims be included within this application. No reference numerals in the claims should be construed as limiting the scope of the claims.

Claims

1. A method for identifying and preventing surge precursors in a centrifugal air compressor, characterized in that, The method includes: Step 1: Acquire dynamic operating data during the operation of the fuel cell centrifugal air compressor, and preprocess the dynamic operating data to obtain input data; Step 2: Based on the input data, construct the observation feature vector according to the sliding short time window, and establish a vortex structure dictionary constrained by the instability mechanism grouping, and divide the dictionary primitives in the vortex structure dictionary into multiple instability mechanism groups; Step 3: For the current short window, based on the observation feature vector, vortex structure dictionary and the status code corresponding to the previous short window, solve the sparse vortex structure status code that satisfies the residual consistency, sparsity constraint and cross-window consistency constraint, and determine the uncertainty index according to the reconstruction residual of the current short window and / or the change in the status code of adjacent short windows. Step 4: Aggregate the sparse vortex structure state codes according to the group to which the dictionary primitives belong, obtain the aggregation intensity of each group, and construct the precursor risk index based on the aggregation intensity of each group, the change in aggregation intensity of adjacent short time windows, and the uncertainty index. Determine the risk level based on the comparison result between the precursor risk index and the preset risk threshold. Step 5: Based on the precursor risk index and its changing trend, under the conditions of meeting the gas supply demand constraints and actuator constraints, generate continuous control commands for the bypass valve opening, selectively generate continuous control commands for motor torque and / or guide vane angle, and when the precursor risk index decreases to the preset recovery condition, perform gradual recovery of the control commands to restore the air compressor to a high-efficiency operating state. Step 6: Based on the precursor risk index and uncertainty, when the preset update conditions are met, the vortex structure dictionary is updated in small steps. When the preset freeze conditions are met, the update of the vortex structure dictionary is frozen and the current dictionary is maintained to participate in the solution of the sparse vortex structure state code.

2. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 1, characterized in that, Step 3 includes: Step 3.1: For the current short window, the observed feature vector corresponding to the current short window is solved by non-negative sparse coding with time-series prior constraints on the vortex structure dictionary. The goal is to minimize the feature reconstruction error, and simultaneously introduces sparsity constraints, cross-window amplitude continuity constraints with the previous short window state code, and suppression constraints on newly added activation components to obtain the sparse vortex structure state code of the current short window. Step 3.2: Based on the weighted reconstruction residual of the current short-term window observation feature vector and the sparse vortex structure state code, and combined with the support set change and code value difference intensity between the current short-term window state code and the previous short-term window state code, an uncertainty index is constructed to characterize the credibility of the current identification result.

3. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 1, characterized in that, Step 4 includes: Step 4.1: Group the dictionary primitives in the vortex structure dictionary according to the instability mechanism to which they belong, and perform grouping summation or weighted aggregation of the activation intensity of the corresponding dictionary primitives in the sparse vortex structure state code of the current short window to obtain the aggregation intensity of the rotating stall related primitive group, the separation enhancement primitive group, the system-level critical coupling primitive group, and the background and interpretable disturbance primitive group, and based on the change of the aggregation intensity of each group between the current short window and the previous short window between adjacent short windows; Step 4.2: Based on the aggregation intensity of each instability mechanism group, the change in aggregation intensity of each group between adjacent short time windows, and the uncertainty index, construct a precursor risk index that jointly characterizes the precursors of rotating stall, separation enhancement, and system-level critical coupling, and performs confidence correction on the background disturbance. Then, compare the precursor risk index with at least one preset risk threshold to determine the risk level corresponding to the current short time window.

4. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 1, characterized in that, Step 5 includes: Step 5.1: Determine the anti-surge control target based on the fuel cell air supply demand, and establish a set of control constraints including target cathode pressure and / or target air flow, excess air coefficient constraints, and amplitude and rate of change constraints of bypass valves, motors and / or guide vane actuators; Step 5.2: Based on the precursor risk index and its adjacent short-term window change trends, generate continuous control commands for bypass valve opening that satisfy the control constraint set, and selectively coordinate the generation of continuous control commands for motor torque and / or guide vane angle according to the risk level, so as to achieve minimum intervention and surge prevention control according to the risk increment. Step 5.3: When the precursor risk index decreases below the preset recovery threshold and the preset stability conditions are continuously met, perform a gradual recovery of the continuous control commands for the bypass valve opening and the continuous control commands for the motor torque and / or guide vane angle, so as to restore the air compressor to a high-efficiency operating state.

5. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 2, characterized in that, In step 3.1, the sparse vortex structure state code for the current short-time window is obtained by solving the following optimization problem: ； Where t represents the time sequence index of the current sliding short window; t-1 represents the previous short window adjacent to the current short window; This represents the observed feature vector corresponding to the t-th short time window; This represents the dictionary matrix of the vortex structure corresponding to the t-th short-time window; This represents the sparse vortex structure state code to be solved in the t-th short time window; This represents the optimal sparse vortex structure state code corresponding to the t-th short time window; This indicates the sparse vortex structure state code corresponding to the previous short time window; Let W represent the objective function for the t-th short-time window; W represents the feature domain weighting matrix. This is a constraint term for residual consistency; Indicates the sparse constraint weight coefficient; For the L1 norm sparse terms of the status code; This represents the weighting coefficient of the cross-window consistency constraint; This is a cross-window consistency constraint term; ρ represents the newly added activation suppression weight coefficient. This represents the support set mask vector constructed from the status codes of the previous short window. This represents a mask vector that is complementary to the support set of the previous short window; Add a new activation inhibition term; Indicates a nonnegativity constraint; Represents the L2 norm; Denotes the norm 1; In step 3.2, the uncertainty index is: ； in, This represents the uncertainty index corresponding to the t-th short-time window; This represents the weighted reconstruction residual energy of the current short-term window; This indicates the status code corresponding to the t-th short window. Support set; This indicates the status code corresponding to the (t-1)th short time window. Support set; △ represents symmetric difference operation; This represents the amount of change in the support set of adjacent short time windows; Indicates the strength of the code value difference between adjacent short-time window status codes; , , The weighting coefficient is greater than 0.

6. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 3, characterized in that, In step 4.1, the polymerization strength is: ； in, The polymerization intensity of the g-th mechanism group within the current short time window t; Let i represent the set of dictionary primitive indices corresponding to the g-th mechanism group; i is the mechanism group. The i-th dictionary primitive index; is the sparse vortex structure state code corresponding to the i-th dictionary primitive within the t-th short-time window; t is the current short-time window number. The change in polymerization intensity, constructed from the difference between the polymerization intensity corresponding to the current short-term window and the previous short-term window, can be expressed as: ； in, This represents the change in polymerization intensity for the g-th mechanism group between two adjacent short time windows; The change in aggregation intensity of the g-th mechanism group in the current short time window relative to the previous short time window is represented by its positive or negative sign, respectively, indicating the strengthening or weakening trend of the instability characteristic corresponding to the mechanism group. This represents the polymerization intensity corresponding to the previous short time window; In step 4.2, the precursor risk index is: ； in, It is a monotonic mapping function used to map the fusion result to a predetermined risk range; This represents the growth trend of the polymerization intensity of the rotational stall-related primitives over time. This represents the growth trend of the aggregation intensity of the critically coupled primitive group at the system level over time. , , The weighting coefficients corresponding to the aggregation intensity of each instability mechanism are used to characterize the contribution of the current activity level of different mechanisms to the total risk; , The weighting coefficients corresponding to the growth trend items; This is the correction weight for uncertainty.

7. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 4, characterized in that, In step 5.1, the set of control constraints corresponding to the t-th control time is: ； in, This represents the control input vector at the t-th control time. The control input vector includes at least the bypass valve opening. Optionally includes motor torque and / or guide vane angle ; This represents the actual cathode pressure at the t-th control moment; This represents the target cathode pressure at the t-th control moment; This indicates that the cathode pressure tracking error does not exceed the preset pressure error upper limit. ; This represents the actual airflow rate at the t-th control time. This represents the target airflow rate at the t-th control time. This indicates that the airflow tracking error does not exceed the preset flow error limit. ; This represents the actual excess air coefficient at control time t; and These represent the lower and upper limits of the allowable excess air coefficient at the t-th control time, respectively. and These represent the control input vectors respectively. The lower limit and upper limit of amplitude; This represents the change in the control input vector at adjacent control moments; This indicates the maximum allowable variation in the control input vector.

8. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 4, characterized in that, In step 5.2, the bypass valve opening command for: ； in, This indicates the bypass valve opening command at the t-th control moment; Indicates the basic opening degree; This indicates the precursory risk index corresponding to the current short-term window; This indicates the precursor risk index corresponding to the previous short time window; This represents the static response gain of the risk index. Gain in response to changes in the risk index trend; This represents a saturation function.

9. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 1, characterized in that, Step 2 includes: Step 2.1: Slice the input data into a sliding short time window and extract observation features that characterize the frequency band energy distribution, cross-channel coupling relationship and the consistency of disturbance propagation, and construct an observation feature vector; Step 2.2: Based on the observed feature vectors, a dictionary learning method constrained by instability mechanism is used to establish a vortex structure dictionary. According to different instability mechanisms, the dictionary primitives are divided into rotating stall related primitive group, separation enhancement primitive group, system-level critical coupling primitive group, and background and interpretable disturbance primitive group.

10. The method for identifying and preventing surge precursors in centrifugal air compressors according to claim 1, characterized in that, Step 6 includes: Step 6.1: Based on the precursor risk index and uncertainty, determine whether the vortex structure dictionary corresponding to the current short time window is in update mode or frozen mode. Specifically, when the precursor risk index is in a preset low-risk range and the uncertainty is below a preset threshold, the vortex structure dictionary is determined to be in update mode; when the precursor risk index is in a preset high-risk range and / or the uncertainty is above a preset threshold, the vortex structure dictionary is determined to be in frozen mode, and the current dictionary is maintained to participate in the sparse vortex structure state code solution. Step 6.2: In update mode, based on the observed feature vector and sparse vortex structure state code corresponding to the current short time window, perform incremental update of the vortex structure dictionary with step size limit and forgetting factor, and apply instability mechanism grouping constraint projection to the dictionary primitives after the update.

Images