Separation column blockage early warning method and system based on time series data

Abstract

The application provides a separation column blockage early warning method and system based on time series data, and relates to the technical field of data processing.The method comprises the following steps: based on a feature vector sequence, the dynamic mutual information correlation degree of the column axial pressure difference change rate and the outlet volume flow fluctuation coefficient in a sliding time window is calculated, the dynamic mutual information correlation degree is compared with a preset dynamic threshold value, and an initial abnormal feature is obtained; based on the time window determined according to the initial abnormal feature, the corresponding shell vibration signal is extracted and frequency domain transformation is performed to obtain an energy spectrum density evolution trend; the energy spectrum density evolution trend is matched with the initial abnormal feature, and a local accumulation state and a spatial distribution probability are obtained.The application realizes the identification and early warning of the local blockage state of the separation column, and improves the stable operation of the separation process.

Description

Technical Field

[0001] This invention relates to the field of data processing technology, and in particular to a method and system for early warning of column blockage based on time-series data. Background Technology

[0002] In the operation monitoring of separation columns, such as fixed-bed adsorption towers in chemical processes, real-time monitoring of the internal fluid state is crucial for maintaining process stability. Most existing monitoring solutions focus on threshold monitoring of macroscopic process parameters such as inlet and outlet pressure difference and flow rate. However, in actual industrial scenarios, blockage inside the separation column often exhibits progressive and non-uniform characteristics. Early local accumulation may not immediately cause the overall pressure difference to exceed the preset safety limit.

[0003] Taking a molecular sieve adsorption tower used by a petrochemical company as an example, when the raw gas carries trace amounts of heavy component impurities, these impurities tend to slowly deposit in specific areas at the bottom of the tower, forming localized narrowing of the flow channels. In this initial stage, the overall axial pressure difference of the tower changes relatively little and remains within the range of normal fluctuations. Traditional alarm systems based on a single pressure difference threshold are mostly unable to trigger warnings. At the same time, this localized change in flow resistance will excite specific micro-vibration signals in the shell, whose frequency and energy distribution are slightly different from those under normal operating conditions. Existing data processing methods mostly analyze pressure data and vibration data separately or use simple linear statistical methods to find the correlation between the two, which may be difficult to effectively capture the complex nonlinear dynamic coupling characteristics between the rate of change of pressure difference and flow fluctuation within a short-term sliding window. Summary of the Invention

[0004] This invention provides a method and system for early warning of blockage in a separation column based on time-series data, which enables the identification and early warning of local blockage in the separation column and improves the stable operation of the separation process.

[0005] To solve the above-mentioned technical problems, the technical solution of the present invention is as follows:

[0006] Firstly, a method for early warning of column blockage based on time-series data, the method comprising:

[0007] Step 1: Process the multidimensional time-series monitoring data to obtain a sequence of feature vectors;

[0008] Step 2: Based on the feature vector sequence, calculate the dynamic mutual information correlation degree between the rate of change of axial pressure difference of the column and the fluctuation coefficient of outlet volume flow rate within the sliding time window, and compare the dynamic mutual information correlation degree with the preset dynamic threshold to obtain the initial abnormal features.

[0009] Step 3: Based on the time window determined by the initial anomaly features, extract the corresponding shell vibration signal and perform frequency domain transformation to obtain the energy spectral density evolution trend; match the energy spectral density evolution trend with the initial anomaly features to obtain the local accumulation state and spatial distribution probability.

[0010] Step 4: Input the local accumulation state and spatial distribution probability into the pre-trained discretized state-space model to obtain the predicted hydraulic state value. Compare the predicted hydraulic state value with the actual monitoring value to generate a state estimation error vector. Map the state estimation error vector to a resistance potential energy deviation sequence to establish a dynamic equilibrium equation, calculate the correction gain matrix and generate a state-space correction vector. Use the state-space correction vector to adaptively update the discretized state-space model to obtain the regional flow resistance distortion correction parameters.

[0011] Step 5: Use the regional flow resistance distortion correction parameter to calibrate the local accumulation state, obtain the calibrated state information, and perform nonlinear time-series extrapolation based on the calibrated state information to obtain the pressure drop increase trajectory; compare the pressure drop increase trajectory with the preset safe operation threshold sequence to obtain the graded blockage early warning command.

[0012] Secondly, a separation column blockage early warning system based on time-series data includes:

[0013] The module is used to process multidimensional time-series monitoring data to obtain a sequence of feature vectors;

[0014] The identification module is used to calculate the dynamic mutual information correlation degree between the rate of change of axial pressure difference of the column and the fluctuation coefficient of outlet volume flow rate within the sliding time window based on the feature vector sequence, and compare the dynamic mutual information correlation degree with the preset dynamic threshold to obtain the initial abnormal features.

[0015] The discrimination module is used to extract the corresponding shell vibration signal based on the time window determined by the initial anomaly features and perform frequency domain transformation to obtain the energy spectral density evolution trend; the energy spectral density evolution trend is matched with the initial anomaly features to obtain the local accumulation state and spatial distribution probability;

[0016] The generation module is used to input the local accumulation state and spatial distribution probability into the pre-trained discretized state-space model to obtain the predicted hydraulic state value. The predicted hydraulic state value is compared with the actual monitoring value to generate a state estimation error vector. The state estimation error vector is mapped to the drag potential energy deviation sequence to establish a dynamic equilibrium equation, calculate the correction gain matrix and generate a state-space correction vector. The state-space correction vector is used to adaptively update the discretized state-space model to obtain the regional flow resistance distortion correction parameters.

[0017] The early warning module is used to calibrate the local accumulation state using regional flow resistance distortion correction parameters, obtain calibrated state information, and perform nonlinear time-series extrapolation calculation based on the calibrated state information to obtain the pressure drop increase trajectory; the pressure drop increase trajectory is compared with the preset safe operation threshold sequence to obtain a graded blockage early warning command.

[0018] The above-described solution of the present invention has at least the following beneficial effects:

[0019] By acquiring multi-dimensional time-series data in real time and performing synchronous alignment and noise filtering, the constructed standard feature vector sequence can characterize the instantaneous hydraulic state of the separation column, avoiding misjudgment and omission due to data deviation. Based on the joint analysis of dynamic mutual information correlation degree and frequency domain energy spectral density, the initial abnormal characteristics of flow field distribution instability are captured, the state and spatial distribution probability of local impurity accumulation are clarified, and early identification of anomalies is achieved. Through adaptive updating of the discretized state space model and flow resistance distortion correction, the impurity accumulation state is calibrated in real time. By predicting the pressure drop rise trajectory through nonlinear time-series extrapolation and generating graded early warning instructions, the risk of blockage can be predicted in advance, and targeted disposal work can be carried out to reduce operation and maintenance costs and production losses. Attached Figure Description

[0020] Figure 1 This is a flowchart illustrating the separation column blockage early warning method based on time-series data provided in an embodiment of the present invention.

[0021] Figure 2 This is a schematic diagram of a separation column blockage early warning system based on time-series data provided in an embodiment of the present invention. Detailed Implementation

[0022] Exemplary embodiments of the present disclosure will now be described in more detail with reference to the accompanying drawings. While exemplary embodiments of the present disclosure are shown in the drawings, it should be understood that the present disclosure may be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that this disclosure will be thorough and complete, and will fully convey the scope of the disclosure to those skilled in the art.

[0023] like Figure 1 As shown, an embodiment of the present invention proposes a separation column blockage early warning method based on time-series data, the method comprising the following steps:

[0024] Step 1: Process the multidimensional time-series monitoring data to obtain a sequence of feature vectors;

[0025] Step 2: Based on the feature vector sequence, calculate the dynamic mutual information correlation degree between the rate of change of axial pressure difference of the column and the fluctuation coefficient of outlet volume flow rate within the sliding time window, and compare the dynamic mutual information correlation degree with the preset dynamic threshold to obtain the initial abnormal features.

[0026] Step 3: Based on the time window determined by the initial anomaly features, extract the corresponding shell vibration signal and perform frequency domain transformation to obtain the energy spectral density evolution trend; match the energy spectral density evolution trend with the initial anomaly features to obtain the local accumulation state and spatial distribution probability.

[0027] Step 4: Input the local accumulation state and spatial distribution probability into the pre-trained discretized state-space model to obtain the predicted hydraulic state value. Compare the predicted hydraulic state value with the actual monitoring value to generate a state estimation error vector. Map the state estimation error vector to a resistance potential energy deviation sequence to establish a dynamic equilibrium equation, calculate the correction gain matrix and generate a state-space correction vector. Use the state-space correction vector to adaptively update the discretized state-space model to obtain the regional flow resistance distortion correction parameters.

[0028] Step 5: Use the regional flow resistance distortion correction parameter to calibrate the local accumulation state, obtain the calibrated state information, and perform nonlinear time-series extrapolation based on the calibrated state information to obtain the pressure drop increase trajectory; compare the pressure drop increase trajectory with the preset safe operation threshold sequence to obtain the graded blockage early warning command.

[0029] In this embodiment of the invention, by acquiring multi-dimensional time-series data in real time and performing synchronous alignment and noise filtering, the constructed standard feature vector sequence can characterize the instantaneous hydraulic state of the separation column, avoiding misjudgment and missed judgment caused by data deviation; based on the joint analysis of dynamic mutual information correlation degree and frequency domain energy spectral density, the initial abnormal characteristics of flow field distribution instability are captured, the state and spatial distribution probability of local impurity accumulation are clarified, and early identification of anomalies is achieved; through the adaptive update of the discretized state space model and flow resistance distortion correction, the impurity accumulation state is calibrated in real time, and the pressure drop rise trajectory is predicted by nonlinear time-series extrapolation and graded early warning instructions are generated, which can predict the blockage risk in advance, carry out targeted disposal work, and reduce operation and maintenance costs and production losses.

[0030] In a preferred embodiment of the present invention, multi-dimensional time-series monitoring data under the operating conditions of the separation column are collected in real time. The multi-dimensional time-series monitoring data includes at least inlet static pressure, outlet volumetric flow rate, column axial pressure difference, medium temperature, and shell vibration acceleration signals, and may include:

[0031] In this embodiment of the invention, the collected multi-dimensional time-series monitoring data includes at least the inlet static pressure, outlet volumetric flow rate, column axial pressure difference, medium temperature, and shell vibration acceleration signal. These five types of parameters can comprehensively reflect the instantaneous hydraulic state, medium characteristics, and equipment vibration of the separation column, covering key characterizing parameters of flow field distortion caused by impurity accumulation. This overcomes the limitations of single-point parameter monitoring in existing technologies. Sensors that meet industrial-grade accuracy requirements, have fast response speeds, and are compatible with the operating medium of the separation column are selected. Following the principles of data acquisition, non-interference with the flow field, and ease of maintenance, the sensors are fixedly deployed. The inlet static pressure sensor is installed on the vertical section of the separation column inlet pipe, 1.5-2.0m away from the inlet flange, avoiding interference from pipe bends, valves, etc., to ensure accurate data acquisition. The inlet static pressure data accurately reflects the inlet operating conditions; the outlet volumetric flow sensor is installed on the outlet pipe of the separation column, with a distance of no less than 3 times the pipe diameter from the outlet flange to ensure stable fluid flow and improve flow acquisition accuracy; the column axial differential pressure sensor adopts a dual-point deployment, installed at the top and bottom of the separation column respectively, with the line connecting the two points parallel to the column axis, to collect the axial pressure difference value of the column and reflect the internal flow resistance change; the medium temperature sensor is inserted into the internal flow channel of the separation column to a depth of no less than 50mm, ensuring full contact with the medium and avoiding temperature measurement deviation caused by the sensor probe being exposed to the air; the shell vibration acceleration sensor is attached to the middle and bottom of the separation column shell, and the corresponding position of the vortex cavity, with a total of 2-3 measurement points deployed to ensure that vibration signals from different areas can be captured.

[0032] A multi-channel data acquisition card is selected, supporting at least 8 analog signal inputs and simultaneously receiving analog signals from five types of sensors. It converts analog signals to digital signals with a sampling bit depth of ≥16 bits. Industrial Ethernet (Profinet protocol) is used for data transmission with a transmission rate of ≥100Mbps to ensure real-time, latency-free data transmission. A backup wireless transmission module (4G / 5G) is also provided to prevent data loss due to wired transmission interruptions. The data caching unit is equipped with an industrial-grade cache server with a cache capacity of at least 1TB for temporary storage of real-time acquired data, preventing data overflow due to processing delays. Cache data retention is at least 72 hours for easy data traceability. A unified sampling frequency of 10-50Hz is set based on the operating cycle of the separation column, with the shell vibration acceleration signal sampling frequency... The sampling frequency is set to 50Hz, and the sampling frequency of other parameters is set to 10Hz. This ensures the timeliness of the data while avoiding data redundancy caused by excessively high sampling frequencies. A unified global timestamp is configured for all sensors and data acquisition modules, and GPS time synchronization or industrial clock synchronization protocols are used to ensure that the data collected by different sensors are synchronized in the time dimension, eliminating time phase deviation caused by sensor response delay. According to the rated operating parameters of the separation column, the acquisition range of each sensor is set. For example, the inlet static pressure range is set to 0-1.6MPa, the outlet volumetric flow rate range is set to 0-50m³ / h, the column axial pressure difference range is set to 0-0.5MPa, the medium temperature range is set to 0-100℃, and the shell vibration acceleration range is set to 0-10g. This ensures that the acquired data is within the range and avoids measurement errors caused by signal saturation.

[0033] During data acquisition, the transmitted data undergoes preliminary verification in real time to remove invalid data. The verification rules are as follows: Range verification: The system checks whether the acquired data is within the preset range. Data exceeding the range is marked as invalid and automatically removed. Sudden change verification: The maximum allowable rate of change for each parameter is set. If the sudden change in data exceeds the allowable range at any given moment and the duration is less than 3 sampling periods, it is determined to be invalid data caused by sensor malfunction or interference and is removed. Zero value verification: For parameters such as outlet volumetric flow rate and column axial pressure difference, if a continuous zero value occurs under normal operating conditions of the separation column and the duration exceeds 5 sampling periods, it is determined to be invalid data, removed, and recorded. An anomaly is detected, prompting a check of the sensor and equipment operating status. Valid multi-dimensional time-series data, after initial verification, is stored in the data cache unit and backend database in the format of timestamp-parameter type-sensor number-data value. A partitioned storage strategy is adopted, storing data by time (one partition per hour) to facilitate rapid data retrieval and retrieval. Stored data is also backed up periodically to prevent data loss. The acquisition system provides real-time feedback on the operating status of each sensor and data acquisition status. If anomalies occur, such as sensor malfunction, data transmission interruption, or an excessively high percentage of invalid data (over 10%), an immediate alert signal is issued to prompt timely troubleshooting and ensure continuous and stable data acquisition.

[0034] In a preferred embodiment of the present invention, step 1 above, which processes the multidimensional time-series monitoring data to obtain a feature vector sequence, may include:

[0035] In this embodiment of the invention, step 110 involves cleaning the real-time acquired multidimensional time-series monitoring data to remove abnormal jump points and noise interference, resulting in preprocessed monitoring data. Specifically, this includes cleaning the data of five types of parameters—inlet static pressure, outlet volumetric flow rate, column axial pressure difference, medium temperature, and shell vibration acceleration—separately based on the pre-verified multidimensional time-series monitoring data to ensure the validity of each type of parameter. For the removal of abnormal jump points, an adjacent data comparison method is used. Using the data at the current sampling time as a benchmark, the difference between the current data and the data at the previous and next sampling times is calculated, and a reasonable jump threshold is set. This threshold is determined based on the normal operating fluctuation range of each type of parameter. For example, the jump threshold for column axial pressure difference is set to 5% of its rated value, and the jump threshold for outlet volumetric flow rate is set to 3% of its rated value. If the difference between the current data and the preceding and following data both exceed the jump threshold of the corresponding parameter, and the duration of the abnormal data does not exceed 3 sampling periods, it is determined to be an abnormal jump point. The data is then removed, and the average of the two valid data before and after is used to replace the abnormal point to ensure the continuity of the data sequence. If the duration of the abnormal data exceeds 3 sampling periods, it is determined to be an anomaly caused by sensor failure. The abnormal data segment is removed, and the fault period is recorded. For noise interference elimination, corresponding processing methods are adopted for the noise characteristics of different parameters. The high-frequency noise of the shell vibration acceleration signal is processed using the moving average method. Five consecutive sampling points are selected as the sliding window, and the average value of the five data within the window is calculated. This average value is used to replace the data of the middle sampling point in the window, and the noise elimination of the entire sequence is completed by sliding sequentially. The data of inlet static pressure, outlet volumetric flow rate, column axial pressure difference, and medium temperature are smoothed. The current data is replaced by calculating the weighted average of the data at the current time and the data at the previous two sampling times (the current data has a weight of 0.6, and the previous two data each have a weight of 0.2), thereby reducing the impact of random noise on the data. After the above processing, preprocessed monitoring data with no anomalies and low noise is obtained.

[0036] Step 111: The preprocessed monitoring data is extracted using a sliding time window of a preset length to obtain data segments within the continuous time window. Specifically, this includes: determining the preset length and sliding step size of the sliding time window. The window length is set based on the operating characteristics of the separation column and the sampling frequency. Considering the early characteristics of local accumulation in the separation column have a certain time duration, and avoiding excessively long windows leading to feature lag or excessively short windows leading to indistinct features, the window length is set to 10-30 seconds. The number of sampling points corresponding to different parameters is as follows: 500-1500 sampling points for shell vibration acceleration signal (sampling frequency 50Hz); and for inlet static pressure, outlet volumetric flow rate, column axial pressure difference, and medium temperature (sampling frequency 10Hz). (Hz) corresponds to 100-300 sampling points. The sliding step size is set to 1 / 5 of the window length. For example, when the window length is 10 seconds, the sliding step size is 2 seconds to ensure that there is a certain amount of data overlap between adjacent windows, avoid missing key features, and reduce data redundancy. When truncating, start from the first sampling point of the preprocessed monitoring data and truncate the first data segment according to the set window length. Then slide backward according to the sliding step size to truncate the next data segment. Repeat this process until all preprocessed data is truncated. Finally, a series of continuous and ordered time window data segments are obtained. Each data segment corresponds to a fixed time interval and contains complete data of all five types of monitoring parameters within that time interval.

[0037] Step 112: Extract time-domain statistical features from the data segments within each time window. These features include at least the mean, variance, peak value, rate of change, and fluctuation coefficient. Specifically, for each type of monitoring parameter data segment within each time window, extract the aforementioned five types of time-domain statistical features to ensure that each type of parameter has a corresponding feature value within each time window. The extraction method is as follows: the mean reflects the average operating level of a certain type of parameter within the time window, and is calculated as the sum of all sampled data for that type of parameter within the window, divided by the sampled value of that type of parameter within the window. The number of sampling points, for example, if there are 100 sampling data points for the outlet volumetric flow rate within a certain window, summing these 100 data points and dividing by 100, yields the mean of the outlet volumetric flow rate within that window. Variance reflects the degree of dispersion of a certain type of parameter data within that time window; the greater the dispersion, the more drastic the parameter fluctuation. For example, if there are 100 sampling data points for the axial pressure difference of a column within a certain window, first calculate the mean of that window, then subtract the mean from each sampling data point to obtain the deviation of each data point, square each deviation, sum them up, and then divide by 99 (100 minus 1) to obtain the variance of the axial pressure difference of the column within that window. The peak value reflects the maximum fluctuation amplitude of a certain type of parameter within the time window. It is calculated by subtracting the minimum sampled data value from the maximum sampled data value of the parameter within the window. For example, if the maximum sampled value of shell vibration acceleration within a certain window is 5g and the minimum sampled value is 1g, then the peak value of shell vibration acceleration within that window is 4g. The rate of change reflects the overall trend of a certain type of parameter within the time window. It is calculated by subtracting the first sampled data value from the last sampled data value of the parameter within the window, and then dividing by the time length of the window (the window time length is equal to the number of sampling points in the window divided by the time length of the window). The sampling frequency of the parameters), for example, if the first sampled value of the inlet static pressure in a certain window is 0.8 MPa and the last sampled value is 0.9 MPa, and the window length is 10 seconds, then the rate of change of the inlet static pressure in this window is 0.01 MPa / second; the fluctuation coefficient reflects the relative fluctuation degree of a certain type of parameter in this time window, and is calculated by dividing the standard deviation of a certain type of parameter in this window by the mean. According to the above method, the five types of time-domain statistical features of all five types of monitoring parameters in each time window are extracted, and each time window corresponds to 25 feature values ​​(5 types of parameters × 5 types of features).

[0038] Step 113: Normalize the time-domain statistical features extracted from each time window and combine them in chronological order to form a feature vector sequence. Specifically, this includes: uniformly normalizing the same type of time-domain statistical features extracted from all time windows. The purpose is to eliminate dimensional differences between different parameters and features, avoiding feature weight imbalances caused by different dimensions, which could affect subsequent analysis results. The normalization process uses the min-max normalization method; for example, the minimum mean of the outlet volumetric flow rate within all time windows is 10 m³ / s. 3 / h, with a maximum value of 40m3 / h, the average outlet volumetric flow rate within a certain window is 25m³ / h. 3 / h, then the normalized value of the mean is (25-10)÷(40-10)=0.5; after completing the normalization of all features, the 25 normalized feature values ​​corresponding to each time window are arranged in a fixed order, namely the mean, variance, peak value, rate of change, and fluctuation coefficient of the inlet static pressure, followed by the five types of features: outlet volumetric flow rate, column axial pressure difference, medium temperature, and shell vibration acceleration, forming a 25-dimensional feature vector. Then, according to the order of the time windows, the feature vectors corresponding to all time windows are combined in sequence to form a continuous feature vector sequence.

[0039] By using a sliding time window to extract data segments, the system captures the operational characteristics of the separation columns within different time intervals, adapting to the gradual and non-uniform nature of local accumulation. By extracting and normalizing multiple types of time-domain statistical features, the system comprehensively captures the operational level, fluctuation degree, and trend of various monitoring parameters, thereby improving the accuracy and timeliness of early warnings.

[0040] In a preferred embodiment of the present invention, step 2 above, based on the feature vector sequence, calculates the dynamic mutual information correlation between the rate of change of axial pressure difference of the column and the fluctuation coefficient of outlet volumetric flow rate within the sliding time window, and compares the dynamic mutual information correlation with a preset dynamic threshold to obtain initial abnormal features, which may include:

[0041] In this embodiment of the invention, step 220 involves extracting the feature components corresponding to the rate of change of the axial pressure difference of the column from the feature vector sequence to form a time series of the rate of change of pressure difference, and simultaneously extracting the feature components corresponding to the fluctuation coefficient of the outlet volumetric flow rate to form a time series of the flow rate fluctuation coefficient. Specifically, this includes: combining the feature vector sequence generated in step 113, determining the fixed arrangement order of each feature vector, where the mean, variance, peak value, rate of change, and fluctuation coefficient of the inlet static pressure are, in order, five types of time-domain statistical features of the outlet volumetric flow rate, the axial pressure difference of the column, the medium temperature, and the shell vibration acceleration. Each feature vector has 25 dimensions, where the rate of change of the axial pressure difference of the column corresponds to the 14th feature component in the feature vector (the axial pressure difference of the column is the third type of parameter, and the rate of change is the fourth feature of this type of parameter, 5×2+4=1). 4) The 10th feature component in the feature vector corresponding to the outlet volumetric flow rate fluctuation coefficient (outlet volumetric flow rate is the second type of parameter, and the fluctuation coefficient is the fifth feature of this type of parameter, 5×1+5=10); according to the order of the time window, the corresponding 14th feature component in each feature vector is extracted one by one, and all the extracted components are arranged in time order to form a pressure difference change rate time series. This series can completely reflect the changing trend of the axial pressure difference of the column in different time intervals; at the same time, the corresponding 10th feature component in each feature vector is extracted one by one and arranged in time order to form a flow rate fluctuation coefficient time series. This series can reflect the relative fluctuation degree of the outlet volumetric flow rate in different time intervals. The length of the two time series is consistent with the feature vector sequence, and the time dimension is completely synchronized.

[0042] Step 221 involves segmenting the time series of differential pressure change rate and the time series of flow fluctuation coefficient using a sliding time window of preset length, resulting in multiple consecutive time window segments corresponding to the differential pressure change rate and flow fluctuation coefficient. Specifically, the sliding time window used is adapted to the sliding time window in step 111, balancing time resolution and feature integrity to avoid omissions or redundancy of coupled features due to unreasonable window settings. The length of the sliding time window is set, taking into account the duration of early features of local accumulation in the separation column and the time interval of the feature vector sequence. Each feature vector corresponds to a sliding time window in step 111, with the time interval being the sliding step size. The sliding time window length is set to 5-10 feature vectors, corresponding to an actual time length of 10-60 seconds, ensuring that the dynamic coupling relationship between the differential pressure change rate and the flow fluctuation coefficient can be captured. The sliding step size is set to half the window length; for example, when the window length is 8 feature vectors, the sliding step size is 4 feature vectors. This ensures a certain amount of data overlap between adjacent windows, avoiding omission of key coupled features, while also reducing data redundancy and computational pressure.

[0043] During segmented processing, the process starts simultaneously from the beginning of the time series of differential pressure change rate and the time series of flow fluctuation coefficient. According to the set window length, the first window segment of each series is extracted, and the two window segments correspond to the same time interval. Then, according to the set sliding step size, the window slides backward simultaneously to extract the next window segment. This process is repeated until all segmented processing of the two time series is completed, resulting in multiple continuous and synchronous window segments of differential pressure change rate and flow fluctuation coefficient. Each time window corresponds to a set of two synchronous window segments.

[0044] Step 222: For each time window, calculate the dynamic mutual information correlation degree between the pressure difference change rate window segment and the flow fluctuation coefficient window segment to obtain the dynamic mutual information correlation degree sequence. Specifically, the dynamic mutual information correlation degree is used to quantify the degree of nonlinear dynamic coupling between two time series segments. The larger the correlation degree value, the tighter the coupling relationship between the pressure difference change rate and the flow fluctuation coefficient, and the more likely there is abnormal flow field distortion, i.e., flow resistance change caused by local accumulation. The calculation formula is as follows: ,in This represents the dynamic mutual information correlation between a window segment representing the rate of change of pressure difference and a window segment representing the flow fluctuation coefficient within a certain time window. Represents a window segment of the rate of change of pressure difference The information entropy reflects the uncertainty of the differential pressure change rate data within the window; Represents a window segment of the flow fluctuation coefficient The information entropy reflects the uncertainty of the flow fluctuation coefficient data within this window; Represents a window segment of the rate of change of pressure difference With flow fluctuation coefficient window segment The joint information entropy reflects the shared uncertainty of the two sequence segments.

[0045] The specific calculation method for information entropy is as follows: for the discretized sequence segment, ,in For sequence fragments Information entropy For sequence fragments The number of value categories after discretization For sequence fragments The value of the middle is The probability is calculated by dividing the number of times the value occurs by the total number of data points in the window segment. The specific calculation process involves first dividing the pressure difference change rate window segment corresponding to each time window... and flow fluctuation coefficient window segment Discretize the data into two segments, mapping each segment to 5-10 discrete intervals, counting the frequency of data occurrences within each interval, and calculating the probability of each value. and and the probability of taking the joint value. Then calculate according to the above formula respectively. , and Finally, the three entropy values ​​are substituted into the dynamic mutual information correlation formula to calculate the dynamic mutual information correlation of the time window. Following the above method, the dynamic mutual information correlation of all time windows is calculated in sequence. All correlation values ​​are arranged in the order of the time windows to obtain the dynamic mutual information correlation sequence. This sequence can reflect the changing trend of the coupling relationship between the pressure difference change rate and the flow fluctuation coefficient in different time intervals.

[0046] Step 223: Compare each value in the dynamic mutual information correlation sequence with a preset dynamic threshold. If the correlation exceeds the dynamic threshold, mark the time window as a potential anomaly window. Specifically, the preset dynamic threshold is not a fixed value, but is determined by combining the statistical results of dynamic mutual information correlation under normal operating conditions of the separation column. This ensures that the threshold can adapt to different operating conditions and reduce false alarms and false negatives. First, collect multi-dimensional time-series monitoring data under normal operating conditions of the separation column (no local accumulation, stable flow field). Following the methods in steps 1 and 220 to 222, calculate the dynamic mutual information correlation sequence under normal operating conditions, and statistically analyze the mean and standard deviation of the sequence. The dynamic threshold is set as follows: The threshold is the mean correlation during normal operation plus twice the standard deviation. This threshold can cover the maximum fluctuation range during normal operation and can effectively identify abnormal correlations that exceed the normal range.

[0047] During the comparison process, each correlation value in the dynamic mutual information correlation sequence is read one by one and compared with a preset dynamic threshold. If the correlation value of a certain time window is greater than the dynamic threshold, it indicates that the coupling relationship between the pressure difference change rate and the flow fluctuation coefficient within that time window exceeds the normal range, and there may be an anomaly in the flow field (flow resistance distortion caused by local accumulation). In this case, the time window is marked as a potential abnormal window, and the time interval, correlation value, and difference between the correlation and the threshold are recorded. If the correlation value is less than or equal to the dynamic threshold, the time window is determined to be a normal window and is not marked.

[0048] Step 224 involves merging and filtering all potential abnormal windows, eliminating isolated short-term anomalies, and extracting continuous abnormal time intervals and their corresponding correlation deviations to obtain initial abnormal features. Specifically, this includes: first, organizing all marked potential abnormal windows and arranging them in chronological order; then determining whether adjacent potential abnormal windows are continuous, i.e., whether the end time of the previous window and the start time of the next window are connected without any time interval; if two adjacent potential abnormal windows are continuous, they are merged into one continuous abnormal time interval; if there are normal windows between two adjacent potential abnormal windows, and the number of normal windows between them does not exceed two, they are also merged into one continuous abnormal time interval to ensure that no continuous abnormal features are missed.

[0049] Subsequently, isolated short-term anomalies are eliminated. Screening criteria are set: if a potential anomaly window is not merged with any other potential anomaly windows, and the actual duration of that window is less than 30 seconds (i.e., the number of windows is less than 5), combined with the window length setting, it is determined to be an isolated short-term anomaly. Such anomalies are mostly caused by external interference or instantaneous sensor fluctuations, not by local accumulation, and are therefore eliminated. After screening, for each merged continuous anomaly time interval, corresponding initial anomaly features are extracted. These initial anomaly features include at least the anomaly start time, anomaly end time, anomaly peak intensity, and anomaly duration. The anomaly start time is the start time of the first window in the continuous anomaly time interval, and the anomaly end time is the end time of the last window in that interval. The anomaly peak intensity is the maximum value among all dynamic mutual information correlation values ​​within that interval. The anomaly duration is the anomaly end time minus the anomaly start time. Simultaneously, the correlation deviation is calculated, i.e., the anomaly peak intensity within that interval minus a preset dynamic threshold, used to characterize the severity of the anomaly. The initial anomaly features corresponding to all continuous anomaly time intervals are then compiled and summarized to obtain complete initial anomaly features.

[0050] By calculating the correlation degree of dynamic mutual information, the complex nonlinear dynamic coupling characteristics between two types of parameters are captured. The setting of dynamic thresholds is combined with the statistical results of normal operating conditions to adapt to different operating scenarios and reduce false alarms and missed alarms. By merging and filtering potential abnormal windows, interfering isolated short-term anomalies are eliminated to ensure the authenticity of the initial abnormal characteristics and improve the early warning of separation column blockage.

[0051] In a preferred embodiment of the present invention, step 3 above, based on the time window determined by the initial abnormal features, extracts the corresponding shell vibration signal, performs frequency domain transformation to obtain the energy spectral density evolution trend, and matches it with the initial abnormal features to determine the local accumulation state and spatial distribution probability of impurities in the filling layer inside the separation column or on the wall of the swirling cavity, may include:

[0052] In this embodiment of the invention, step 330 involves extracting the shell vibration acceleration signal within the corresponding time period from the multidimensional time-series monitoring data based on the abnormal start time and abnormal end time included in the initial abnormal features, thereby obtaining a vibration signal segment. Specifically, this includes: using the initial abnormal features obtained in step 224 as the core basis, first determining the abnormal start time and abnormal end time corresponding to each initial abnormal feature. This time information is completely synchronized with the global timestamp of the multidimensional time-series monitoring data collected in step 1, ensuring time dimension matching. From the multidimensional time-series monitoring data collected in real-time and preliminarily verified in step 1, parameters such as shell vibration acceleration signals are selected. Using the abnormal start time of each initial abnormal feature as the starting node and the abnormal end time as the ending node, all vibration signals are analyzed one by one. The shell vibration acceleration signal at the measuring points is extracted. A total of 2-3 measuring points are deployed for the shell vibration acceleration signal, located in the middle, bottom and corresponding positions of the vortex cavity of the separation column shell. When extracting, it is necessary to ensure that the signal of each measuring point is complete and without loss during the abnormal period. Blank data, duplicate data and signal breakpoints within the extraction range are removed. At the same time, the original sampling frequency (50Hz) and timestamp information of the signal are retained, so that the duration of the vibration signal segment after extraction at each measuring point is completely consistent with the abnormal duration of the initial abnormal feature. Moreover, the vibration signal segment of each measuring point is completely synchronized with the signal segments of other measuring points on the time axis. A vibration signal segment corresponding to each initial abnormal feature and containing the synchronized data of all vibration measuring points is obtained. This segment is used as the only input data for the frequency domain transformation in step 331.

[0053] Step 331 involves performing a frequency domain transformation on the vibration signal segment, converting the time-domain signal into a spectral distribution, and calculating the energy spectral density in different frequency intervals to obtain an energy spectral density evolution sequence. Specifically, this includes: using the vibration signal segment obtained in step 330 as the processing object, preprocessing the vibration signal segment at each measuring point to avoid spectral leakage during the frequency domain transformation. The preprocessing uses a Hanning window, multiplying the Hanning window with the vibration signal segment point by point to smooth the attenuation at both ends of the signal segment and reduce spectral distortion caused by signal abrupt changes. After the windowing process is completed, a fast Fourier transform is performed on the windowed vibration signal segment at each measuring point to convert the time-domain vibration signal that originally varied with time into frequency-domain spectral data that varied with frequency, realizing the conversion from the time-acceleration dimension to the frequency-amplitude dimension, and obtaining the complex spectral data corresponding to the vibration signal segment at each measuring point.

[0054] The energy spectral density at each frequency point is calculated based on the obtained complex spectrum data. Energy spectral density characterizes the energy of the vibration signal within a unit frequency range and is a core indicator reflecting vibration characteristics. Its calculation formula is as follows: ,in Represents frequency Energy spectral density at that location Represents the frequency after Fast Fourier Transform The corresponding complex spectrum amplitude, The sampling frequency of the shell vibration acceleration signal is indicated (fixed at 50Hz). This represents the total number of sampling points within the vibration signal segment (the total number of sampling points equals the sampling frequency multiplied by the duration of the anomaly). After calculation, the frequency range is reasonably divided. Combining the vibration frequency characteristics of the separation column during normal operation and local accumulation, the frequency is divided into three intervals: low frequency interval (0-100Hz), mid frequency interval (100-500Hz), and high frequency interval (500-1000Hz). The micro-amplitude vibrations caused by local accumulation are mostly concentrated in the mid-frequency and high-frequency intervals. According to the time sequence of the vibration signal segment, the sum of the energy spectral density of all frequency points in each frequency interval is calculated time by time to obtain the energy spectral density value corresponding to each time and each frequency interval. The energy spectral density values ​​of all times are arranged in chronological order to form the energy spectral density evolution sequence corresponding to each measuring point. This sequence fully reflects the change trend of vibration energy in each frequency interval during the abnormal period and serves as the input data for step 332 to extract the evolution trend of vibration characteristics.

[0055] Step 332 involves extracting the dominant frequency migration trajectory, high-frequency energy ratio changes, and specific frequency band energy abrupt change points from the energy spectral density evolution sequence to form the vibration characteristic evolution trend. Specifically, this includes: using the energy spectral density evolution sequence of each measurement point obtained in step 331 as the processing object, extracting features from the sequence of each measurement point one by one to ensure that the extracted vibration features accurately reflect the vibration changes caused by local accumulation. The specific extraction process is as follows: First, extract the dominant frequency migration trajectory by traversing all frequency points at each moment in the energy spectral density evolution sequence and selecting the frequency with the largest energy spectral density value at that moment as the dominant frequency. The dominant frequency is the frequency at which energy is most concentrated in the vibration signal. Its change is directly related to the location and degree of local accumulation. Arranging the dominant frequencies of each moment at the measuring point in chronological order forms the dominant frequency migration trajectory of the measuring point, which fully reflects the change law of the dominant frequency with time during the abnormal period. Secondly, the change of high-frequency energy proportion is extracted. For each moment in the energy spectral density evolution sequence, the sum of energy spectral density of the entire frequency range (0-1000Hz) is first calculated, and then the sum of energy spectral density of the high-frequency interval (500-1000Hz) is calculated. The sum of energy spectral density of the high-frequency interval is divided by the sum of energy spectral density of the entire frequency range. The sum of the energy spectral density of the surrounding area is used to obtain the high-frequency energy ratio at that moment. The higher the high-frequency energy ratio, the more obvious the micro-amplitude vibration caused by local accumulation. The high-frequency energy ratios at each moment are arranged in chronological order to obtain the high-frequency energy ratio change curve of that measurement point. Third, specific frequency band energy mutation points are extracted, and a preset energy mutation threshold is set. This threshold is determined based on the fluctuation range of the high-frequency range energy spectral density during normal operation of the separation column. It is usually set to 3 times the average energy spectral density of the high-frequency range during normal operation. The energy spectral density evolution sequence is traversed. If the energy spectral density of a certain frequency band (mainly the mid-frequency and high-frequency ranges) exceeds a certain threshold, the energy mutation threshold is determined. If the increase in amplitude within a sampling period exceeds the preset energy mutation threshold, and the upward trend continues for at least two sampling periods, then the frequency range of that frequency band and the corresponding time are marked as specific frequency band energy mutation points. These mutation points are important features that indicate the beginning or aggravation of local accumulation. After extraction, the main frequency migration trajectory, high-frequency energy ratio change curve, and specific frequency band energy mutation points of each measurement point are integrated. Redundant features between different measurement points are eliminated, and vibration features with consistency and representativeness are retained to form a complete vibration feature evolution trend. This trend serves as the core input for spatiotemporal correlation matching in step 333.

[0056] Step 333 involves performing spatiotemporal correlation matching between the vibration feature evolution trend and the abnormal peak intensity and duration in the initial anomaly features. Based on the alignment of the energy mutation point in a specific frequency band with the abnormal peak intensity on the time axis, and the correlation between the rising trend of the high-frequency energy proportion and the anomaly duration, it is determined whether there is local accumulation of solid particles within the time window, thus obtaining the local accumulation state. Specifically, this includes: using the vibration feature evolution trend obtained in step 332 and the initial anomaly features obtained in step 224 as common inputs, performing spatiotemporal correlation matching on a unified time axis to ensure that the vibration features and hydraulic anomaly features are consistent. The specific matching and judgment process is as follows: Time-dimensional alignment matching is performed. The occurrence time of the energy mutation point in a specific frequency band in the evolution trend of vibration characteristics is compared one by one with the occurrence time of the abnormal peak intensity in the initial abnormal characteristics. Since both are based on the same global timestamp, if the occurrence time of the energy mutation point in a specific frequency band and the occurrence time of the abnormal peak intensity completely coincide, or the time deviation does not exceed one sampling period (0.02 seconds for a 50Hz sampling frequency), it is determined that the two are highly coupled in time, indicating that the vibration anomaly and the hydraulic parameter anomaly are caused by the same reason (i.e., local accumulation).

[0057] Correlation matching of trends is performed by comparing the change curve of high-frequency energy proportion in the vibration feature evolution trend with the duration of the anomaly in the initial anomaly features. The change law of high-frequency energy proportion with the duration of the anomaly is observed. If the high-frequency energy proportion gradually increases with the increase of the duration of the anomaly, and the increase is positively correlated with the magnitude of the anomaly peak intensity (the greater the peak intensity, the greater the increase in the high-frequency energy proportion), then the two trends are determined to be highly correlated, further verifying that the vibration anomaly is caused by local accumulation. Based on the above two matching results, the local accumulation state is determined. If both time alignment matching and trend correlation matching are satisfied, then it is determined that there is local accumulation of solid particles within the time window; if only one matching is satisfied, or neither is satisfied, then it is determined that there is no local accumulation within the time window (the anomaly is caused by external interference or sensor fluctuations).

[0058] Based on the degree of matching, the accumulation location and degree of accumulation are further distinguished. If the energy mutation point of a specific frequency band is concentrated in the mid-frequency range and the main frequency migration trajectory is biased towards the low frequency, the accumulation location is determined to be the filling layer inside the separation column. If the energy mutation point of a specific frequency band is concentrated in the high-frequency range and the main frequency migration trajectory is biased towards the high frequency, the accumulation location is determined to be the wall of the vortex cavity. Based on the increase in the proportion of high-frequency energy, the degree of accumulation is divided (an increase of 10%-30% is slight accumulation, 30%-60% is moderate accumulation, and more than 60% is severe accumulation). The local accumulation state containing three core information items, namely the presence or absence of accumulation, accumulation location, and accumulation degree, is obtained. This state serves as the input basis for solving the spatial distribution probability in step 334.

[0059] Step 334: Based on the local accumulation state, and according to the spatial directivity of the dominant frequency migration trajectory in the vibration characteristic evolution trend and the phase difference of vibration signals from different measuring points, calculate the spatial distribution probability of the accumulation area in the radial and axial directions of the column to obtain the spatial distribution probability. Specifically, this includes: taking the local accumulation state obtained in step 333 (confirming the existence of local accumulation) as a premise, combining the vibration characteristic evolution trend (dominant frequency migration trajectory) obtained in step 332 and the multi-measuring-point vibration signal segments (phase information) obtained in step 330, calculating the spatial distribution probability. The specific process is as follows: based on the spatial directivity of the dominant frequency migration trajectory, initially locate the axial direction of the accumulation area. Based on the location and structural characteristics of the separation column, the main frequency migration patterns caused by accumulation at different axial positions are different. Specifically, low-frequency main frequency (0-100Hz) migration corresponds to accumulation at the bottom of the column, mid-frequency main frequency (100-500Hz) migration corresponds to accumulation in the middle of the column, and high-frequency main frequency (500-1000Hz) migration corresponds to accumulation at the top of the column or in the vortex cavity region. Based on the duration and amplitude of the main frequency migration, the axial range of the accumulation area can be preliminarily determined. The longer the duration and the greater the amplitude of the main frequency migration, the wider the axial accumulation range. The initial accumulation probability values ​​for each axial region (top, middle, bottom, and vortex cavity) can be preliminarily determined.

[0060] Based on the phase difference of vibration signals from different measuring points, the radial position of the accumulation area is located and the axial probability is corrected. The phase information of vibration signal segments from different measuring points (middle, bottom, and vortex cavity) obtained in step 330 is extracted. At the same time, the phase difference of vibration signals from each measuring point is compared. Measuring points with phase leading indicate that they are closer to the accumulation area (vibration signals generated in the accumulation area will be transmitted to nearby measuring points first), while measuring points with phase lagging are far from the accumulation area. According to the degree of phase leading, the initial value of accumulation probability is assigned to each radial region (central region and wall region). The greater the phase leading, the higher the accumulation probability of the region near the corresponding measuring point.

[0061] The initial values ​​of axial and radial accumulation probabilities are normalized so that the sum of probabilities of all regions is 1, eliminating the influence of dimensional differences. At the same time, the probability values ​​are corrected by combining the accumulation location and accumulation degree obtained in step 333. For example, when the accumulation location is determined to be the filling layer (middle), the accumulation probability of the middle of the column is increased and the probability of other regions is decreased. When it is determined to be heavy accumulation, the probability distribution range of the corresponding region is expanded to obtain the accumulation probability of each region including the radial (central region, wall region) and axial (top, middle, bottom, vortex cavity) regions of the column, that is, the complete spatial distribution probability.

[0062] Targeted extraction of shell vibration signals and frequency domain analysis were performed, achieving deep integration of vibration signals and hydraulic anomaly characteristics. By extracting the evolution trend of energy spectral density and matching spatiotemporal correlation, the micro-amplitude vibration characteristics in the early stage of local accumulation were captured. Combined with the directivity of the main frequency migration and the phase difference of multiple measurement points, the radial and axial spatial distribution probability of the accumulation area was calculated, improving the early warning of separation column blockage.

[0063] In a preferred embodiment of the present invention, step 4 above involves inputting the local accumulation state and spatial distribution probability into a pre-trained discretized state-space model to obtain a predicted hydraulic state value. The predicted hydraulic state value is then compared with the actual monitored value to generate a state estimation error vector. This state estimation error vector is mapped to a drag potential energy deviation sequence to establish a dynamic equilibrium equation. A correction gain matrix is ​​calculated, and a state-space correction vector is generated. The state-space correction vector is then used to adaptively update the discretized state-space model to obtain regional flow resistance distortion correction parameters, which may include:

[0064] In this embodiment of the invention, step 440 involves constructing an input vector for a discretized state-space model based on the local accumulation state and spatial distribution probability, combined with the current operating parameters. The input vector is then substituted into the pre-trained discretized state-space model to perform forward recursive calculations, yielding the predicted hydraulic state value at the current moment. Specifically, this includes completing the construction and pre-training of the discretized state-space model to ensure that the model reflects the correlation between the hydraulic state of the separation column and the local accumulation state and spatial distribution probability. The specific construction and training process is as follows: the discretized state-space model is constructed based on the hydraulic dynamic characteristics of the separation column, consisting of two parts: a state equation and an observation equation. The state equation describes the dynamic relationship between the change in flow resistance inside the separation column and the hydraulic state (axial pressure difference of the column, outlet volumetric flow rate), adapting to the flow resistance distortion scenario caused by local accumulation. The state equation describes the changing law of the hydraulic state inside the separation column, and its discretized form is... ,in, for The state vector at time t, containing The two core hydraulic state parameters are the axial pressure difference of the column and the outlet volumetric flow rate. for The state vector at any given time; for The state transition matrix at time step 1, with a dimension of 2×2, is used to represent... The hydraulic state at any time The influence of the hydraulic state at time +1 is related to the flow resistance characteristics of the separation column. for The input matrix at time step has a dimension of 2×m (where m is the number of input parameters) and is used to characterize the influence of the input parameters on the hydraulic state. for The input vector at time step 1 is consistent with the input vector constructed subsequently. The process noise vector follows a normal distribution with a mean of 0 and a pre-defined variance, used to characterize random disturbances (such as small fluctuations in medium viscosity) that the model cannot capture. The observation equation describes the relationship between the observed and actual hydraulic states, and its discretized form is... ,in, for The observation vector at time t, namely the measured values ​​of the axial pressure difference of the column and the measured values ​​of the outlet volumetric flow rate after preprocessing in step 1; for The observation matrix at time step has a dimension of 2×2 and is used to represent the mapping relationship between the state vector and the observation vector. Its initial value is set as the identity matrix. The observed noise vector follows a normal distribution with a mean of 0 and a variance of a preset constant, and is used to characterize the sensor measurement error.

[0065] The model training utilizes historical monitoring data from different operating conditions of the separation column to ensure that the trained model can adapt to the hydraulic characteristics of normal operating conditions and different degrees of local accumulation. The specific training process is as follows: Training data is collected, including multi-dimensional time-series monitoring data for four operating conditions: normal operation, slight accumulation, moderate accumulation, and severe accumulation. The monitoring data duration for each condition is no less than 72 hours, covering different medium temperatures, inlet static pressures, and other operating conditions. All data undergoes preprocessing in step 1 (data cleaning, truncation, and feature extraction), and the local accumulation state, spatial distribution probability, and measured hydraulic state values ​​for the corresponding operating conditions are labeled. Model parameters and the state transition matrix are then initialized. Input matrix Observation matrix The initial value, the process noise vector and observation noise vector The variance is initially set based on the rated operating parameters of the separation column and engineering experience. Iterative optimization training uses local accumulation states, spatial distribution probabilities, and operating condition parameters from the training data as input, and measured hydraulic state values ​​as labels. These are substituted into the discretized state-space model, and the sum of squared errors between the predicted and measured hydraulic state values ​​is minimized using the least squares method. The model is then iteratively updated. , , The parameter values ​​are calculated until the sum of squared errors is less than the preset training convergence threshold (usually set to 0.001). For model validation, test data (covering four operating conditions) that were not used in training are selected and substituted into the trained model. The error between the predicted hydraulic state value and the measured value is calculated. If the mean error is less than the preset allowable error (axial pressure difference error of the column ≤ 0.01 MPa, outlet volumetric flow rate error ≤ 0.5 m³ / h), the model training is complete. If the error does not meet the requirements, the process returns to step three to continue iterative optimization until the requirements are met, resulting in a pre-trained discretized state-space model.

[0066] After the model pre-training is completed, the input vector is constructed and forward recursive calculation is performed. The local accumulation state obtained in step 333 and the spatial distribution probability obtained in step 334 are used as the core inputs. Combined with the current operating condition parameters (the measured values ​​of inlet static pressure and medium temperature after preprocessing in step 1), the input vector of the discretized state-space model is constructed. In this process, the degree of accumulation in the local accumulation state is quantified as a value between 0 and 1 (0.25 for slight accumulation, 0.5 for moderate accumulation, and 0.75 for severe accumulation). The accumulation location is assigned a value according to different regions (0.4 for the filling layer and 0.6 for the vortex cavity wall). The spatial distribution probability is directly adopted from the radial and axial probability values ​​of each region obtained in step 334. The operating condition parameters are the measured normalized values ​​at the current moment. These parameters are arranged in a fixed order to form an input vector with a dimension of 1×m. The constructed input vector Substituting the pre-trained discretized state-space model, forward recursive calculation is performed according to the state equation and observation equation, i.e., through... State vector at time step Input vector and state transition matrix Input matrix Calculations yielded The state vector at time +1 Then through the observation matrix The predicted hydraulic state value at the current moment is obtained by mapping. This predicted hydraulic state value includes at least the predicted axial pressure difference of the column and the predicted outlet volumetric flow rate, which serves as the comparison benchmark for step 441.

[0067] Step 441 involves comparing the predicted hydraulic state value with the measured values ​​of the column axial pressure difference and outlet volumetric flow rate at the same moment in the preprocessed actual monitoring data, calculating the deviation in each dimension and combining them to obtain the state estimation error vector. Specifically, this includes: using the predicted hydraulic state value obtained in step 440 as the processing object, simultaneously retrieving the actual monitoring data after data cleaning and preprocessing in step 1, and extracting the measured values ​​of the column axial pressure difference and outlet volumetric flow rate at the same moment as the predicted hydraulic state value (based on global timestamp alignment), ensuring that the two values ​​are completely synchronized in the time dimension to avoid inaccurate calculations due to time deviations; then performing a dimension-by-dimensional comparison calculation, with the first dimension comparing the column axial pressure difference... The deviation is calculated by subtracting the predicted axial pressure difference of the column from the measured value of the column axial pressure difference at the same time. The second dimension is the comparison of the outlet volumetric flow rate, which is obtained by subtracting the predicted outlet volumetric flow rate from the measured value of the outlet volumetric flow rate at the same time. The calculated deviations of the column axial pressure difference and the outlet volumetric flow rate are combined in sequence to form a 2-dimensional state estimation error vector. Each element of this vector corresponds to the prediction deviation of a hydraulic parameter. A positive deviation value indicates that the measured value is greater than the predicted value, and a negative deviation value indicates that the measured value is less than the predicted value. The larger the absolute value of the deviation, the lower the prediction accuracy of the discretized state space model. This state estimation error vector is used as the input data for step 442.

[0068] Step 442 involves converting the pressure difference deviation and flow rate deviation in the state estimation error vector into changes in drag potential energy distributed along the column axis, forming a drag potential energy deviation sequence, and constructing a dynamic equilibrium equation based on this sequence. Specifically, this includes using the state estimation error vector obtained in step 441 as the processing object, first converting the two deviation components in the vector into changes in drag potential energy distributed along the column axis. The conversion logic is based on the hydraulic dynamics principle of the separation column, where the change in drag potential energy is directly related to the pressure difference deviation and flow rate deviation. The specific conversion method is as follows: the drag potential energy change corresponding to the pressure difference deviation along the column axis... The change in force potential energy is calculated by multiplying the axial pressure difference deviation of the column by the cross-sectional area of ​​the separation column, and then by the axial length of the column. This yields the change in resistance potential energy along the column's axial direction, which is equal to the axial pressure difference deviation of the column × the cross-sectional area of ​​the separation column × the axial length of the column. The change in resistance potential energy corresponding to the outlet volumetric flow rate deviation is calculated by multiplying the square of the outlet volumetric flow rate deviation by the average flow resistance coefficient along the column's axial direction, and then by the axial length of the column. This yields the corresponding change in resistance potential energy, which is equal to the outlet volumetric flow rate deviation × the outlet volumetric flow rate deviation × the average flow resistance coefficient × the axial length of the column.

[0069] Following a chronological order, the changes in drag potential energy at the current moment and several previous moments (usually the first five moments) are arranged sequentially to form a drag potential energy deviation sequence. This sequence fully reflects the recent trend of drag potential energy deviation and can demonstrate the impact of flow resistance distortion caused by local accumulation on drag potential energy. Based on the obtained drag potential energy deviation sequence, a dynamic equilibrium equation is constructed. This equation characterizes the balance between drag potential energy deviation and model parameter deviation, ensuring that the solved correction gain matrix can correct model deviations. The form of the dynamic equilibrium equation is: ,in, for The drag potential energy deviation sequence vector at time t has a dimension of 1×n (n is the sequence length, usually taken as 6). for The coefficient matrix at time step has a dimension of n×2 and is used to characterize the relationship between the model parameter deviation and the drag potential energy deviation. for The model parameter deviation vector at time step, corresponding to the state transition matrix. and input matrix Deviation; The equilibrium error vector follows a normal distribution with a mean of 0 and is used to characterize the random disturbances in the equation. This dynamic equilibrium equation serves as the basis for solving the correction gain matrix in step 443.

[0070] Step 443: Solve the dynamic equilibrium equation to obtain the correction gain matrix. Use the correction gain matrix to perform a weighted transformation on the state estimation error vector to obtain the state-space correction vector. Specifically, this includes: using the dynamic equilibrium equation obtained in step 442 as the core, solving the equation using the least squares method to find the final correction gain matrix that minimizes the model parameter deviation, thereby achieving model correction. The specific solution process is as follows: First, organize the dynamic equilibrium equation, and combine the drag potential energy deviation sequence vector and the model parameter deviation vector at multiple time points (k-n+1 to k) into a matrix form to obtain... ,in, It is a matrix composed of the drag potential energy deviation sequence at multiple time points, with a dimension of n×n; It is a matrix consisting of coefficient matrices at multiple time points, with a dimension of n×2; It is a matrix consisting of the model parameter bias vectors at multiple time points, with a dimension of 2×n; It is a matrix consisting of balance error vectors at multiple time points, with dimensions n×n.

[0071] Solving the above matrix equation using the least squares method yields the coefficient matrix. The final solution is used as the correction gain matrix. The dimension of the correction gain matrix is ​​2×2, and its calculation formula is as follows: ,in, Coefficient matrix The transpose of the matrix, for The inverse matrix of the formula is used to find the correction gain that minimizes the drag potential energy deviation through matrix operations, ensuring that the subsequent correction vector can correct the model deviation. After the correction gain matrix is ​​solved, the state estimation error vector obtained in step 441 is weighted and transformed using this matrix. Specifically, the state space correction vector is obtained by multiplying the correction gain matrix by the state estimation error vector. That is, the state space correction vector is equal to the correction gain matrix × the state estimation error vector. The dimension of the state space correction vector is the same as that of the state estimation error vector (2-dimensional). Each element corresponds to the correction amount of a model parameter. The magnitude and direction of the correction amount are jointly determined by the correction gain matrix and the state estimation error vector. It is used to update the parameters of the discretized state space model. The state space correction vector is used as the input data for step 444.

[0072] Step 444 involves adaptively updating the internal state parameters of the discretized state-space model using the state-space correction vector to obtain updated model parameters, which are then used as parameters for correcting regional flow resistance distortion. Specifically, this includes: using the state-space correction vector obtained in step 443 as the processing object, adaptively updating the core internal state parameters of the pre-trained discretized state-space model to ensure the model can adapt to flow resistance distortion caused by local accumulation. The specific update process is as follows: determining the internal state parameters of the model that need to be updated, mainly including the state transition matrix. and input matrix These two matrices directly determine the model's accuracy in predicting the hydraulic state, and are also the parameters most affected by local flow resistance distortion; observation matrix Since it is mainly related to the sensor's observation characteristics, it will not be updated for the time being, and the parameter values ​​after the initial training will be retained.

[0073] Perform parameter update calculations using the current state transition matrix. Adding the first dimension of the state-space correction vector, we obtain the updated state transition matrix. ,Right now (in (The first dimension of the state-space correction vector); using the input matrix at the current time step. Adding the second dimension of the state-space correction vector, we obtain the updated input matrix. ,Right now (in The second dimension of the state-space correction vector is the correction value; after the update is completed, the updated model parameters ( , To verify the prediction, substitute the values ​​into the discretized state-space model, recalculate the predicted hydraulic state values, and compare them with the measured values ​​at the same time. If the prediction error is less than the error before the update and meets the preset allowable error, the parameter update is valid. If the error does not decrease, adjust the weights of the correction gain matrix and update again until the requirements are met. Then, update the state transition matrix. and input matrix As a parameter for regional flow resistance distortion correction.

[0074] By constructing an input vector that includes local accumulation states, spatial distribution probabilities, and operating conditions, the relevance and accuracy of model predictions are improved. Through error analysis, solving dynamic equilibrium equations, and adaptive model updates, the prediction bias caused by flow resistance distortion is effectively corrected, and the obtained regional flow resistance distortion correction parameters can quantify the impact of flow resistance distortion.

[0075] In a preferred embodiment of the present invention, step 5 above, which involves calibrating the local accumulation state using regional flow resistance distortion correction parameters to obtain calibrated state information, and performing nonlinear time-series extrapolation based on the calibrated state information to obtain a pressure drop increase trajectory; comparing the pressure drop increase trajectory with a preset safe operating threshold sequence to obtain a graded blockage warning instruction, may include:

[0076] In this embodiment of the invention, step 550 involves correcting multiple local characteristic parameters included in the local accumulation state based on the regional flow resistance distortion correction parameters, eliminating the deviation introduced by local flow resistance distortion, and obtaining calibrated state information characterizing the overall blockage process of the separation column. Specifically, this includes the regional flow resistance distortion correction parameters obtained in step 444 and the local accumulation state obtained in step 333. The regional flow resistance distortion correction parameters are the updated state transition matrix and input matrix. These two matrices quantify the degree and distribution characteristics of flow resistance distortion caused by local accumulation and are the core basis for correcting the deviation of the local accumulation state. The local accumulation state includes four core local characteristic parameters: presence or absence of accumulation, accumulation location, accumulation degree, and accumulation range. Although these parameters have been determined in step 3, they are affected by local flow resistance distortion and have certain deviations. They need to be corrected one by one in this step to ensure complete consistency with the actual operating state of the separation column.

[0077] The specific correction operations are as follows: each correction strictly relies on the regional flow resistance distortion correction parameter, and the corrected data serves as the basis for the correction. For the correction of the accumulation degree, firstly, the flow resistance distortion quantified value is extracted from the regional flow resistance distortion correction parameter. This value is the average of the diagonal elements of the updated state transition matrix, reflecting the overall degree of local flow resistance distortion. The larger the value, the more severe the flow resistance distortion. During correction, first multiply the accumulation degree value in the original local accumulation state (quantized to values ​​between 0 and 1, with 0.25 for slight accumulation, 0.5 for moderate accumulation, and 0.75 for severe accumulation) by this flow resistance distortion quantified value to obtain the preliminary correction value; then, based on the rated flow resistance of the separation column and the current actual flow... The difference in flow resistance is used to determine the correction compensation amount. The compensation amount is calculated as (current actual flow resistance - rated flow resistance) ÷ rated flow resistance, with a range of 0.01-0.1. The larger the flow resistance deviation, the larger the compensation amount. Finally, the initial correction value is added to the correction compensation amount to obtain the calibrated accumulation level. The calibrated value remains between 0 and 1. If it exceeds the range, the boundary value is taken (0 for less than 0, 1 for greater than 1). For example, if the original accumulation level is 0.5 (moderate accumulation), the distorted flow resistance value is 1.12, the current actual flow resistance is 0.32 MPa, and the rated flow resistance is 0.28 MPa, the correction compensation amount is (0.32 - 0.28) ÷ 0.28 ≈ 0.14. If the compensation exceeds 0.1, it is set to 0.1. The calibrated accumulation degree is 0.5 × 1.12 + 0.1 = 0.66, which is judged as moderate to heavy accumulation. For the correction of the accumulation location, the updated value of the input matrix is ​​extracted from the regional flow resistance distortion correction parameters. The parameters in different columns of the input matrix correspond to the flow resistance influence coefficients of different regions of the separation column (filling layer, vortex cavity wall, column top, and column bottom). The larger the parameter update amplitude, the more obvious the flow resistance distortion in that region, meaning the accumulation is more likely to be concentrated in that region. During correction, the update amplitude of the corresponding parameters in each region of the input matrix is ​​first calculated (updated parameter value - original parameter value), and the region with the largest update amplitude is found. If the original accumulation location... If the original accumulation position is consistent with the region, the original position remains unchanged; if the original accumulation position is inconsistent with the region, the original accumulation position is fine-tuned towards the region with the largest update amplitude. The fine-tuning amplitude is determined according to the difference in update amplitude. The larger the difference in update amplitude, the larger the fine-tuning amplitude, ensuring that the calibrated accumulation position completely corresponds to the actual flow resistance distortion region. For example, if the original accumulation position is determined to be the filling layer, the update amplitude of the corresponding parameter of the vortex cavity wall in the input matrix is ​​0.23, and the update amplitude of the corresponding parameter of the filling layer is 0.08. The update amplitude of the vortex cavity wall is larger, indicating that the flow resistance distortion is mainly concentrated on the vortex cavity wall. At this time, the original accumulation position is fine-tuned to the junction of the vortex cavity wall and the filling layer to match the actual accumulation situation.

[0078] The accumulation range is corrected, encompassing both radial and axial ranges. The radial range is divided into the central region and the wall region, while the axial range includes the top, middle, and bottom of the column, as well as the vortex cavity region. Correction relies on the updated state transition matrix in the region flow resistance distortion correction parameters. For the radial range, the updated value of the corresponding radial parameter in the state transition matrix is ​​multiplied by the original range size of each radial region to obtain the calibrated radial range size. If the updated value is greater than 1, it indicates that flow resistance distortion has caused an expansion of the radial accumulation range; if it is less than 1, it indicates that the radial accumulation range is smaller than the original judgment value. For the axial range, the same method is used: the updated value of the corresponding axial parameter in the state transition matrix is ​​multiplied by the original range size of each axial region to obtain the calibrated axial range size. For example, if the original accumulated radial wall region range is 50-80mm, and the updated radial parameter value in the state transition matrix is ​​1.08, the calibrated radial wall region range is (50×1.08)-(80×1.08)=54-86.4mm, with a slight increase in the accumulated radial range; the original axial middle region range is 1000-1500mm, and the updated value of the axial parameter in the state transition matrix is... The parameter update value is 0.97. The calibrated axial mid-area range is (1000×0.97)-(1500×0.97)=970-1455mm, and the axial accumulation range has slightly decreased. A final verification of the presence or absence of accumulation is performed based on the calibration results of the above three parameters. The accumulation degree judgment threshold is set at 0.2 (i.e., the minimum quantification value for slight accumulation). If the calibrated accumulation degree is greater than 0.2, and the calibrated accumulation range is greater than the preset minimum accumulation range (radial not less than 10mm, axial not less than 50mm), then... The original determination of accumulation is maintained; if the degree of accumulation after calibration is less than 0.2, or the range of accumulation after calibration is less than the preset minimum accumulation range, it is corrected to no accumulation, and it is determined to be a false anomaly caused by external interference and instantaneous fluctuations of the sensor in the early stage, so as to ensure that the determination of whether there is accumulation is accurate. After the above four local characteristic parameters are corrected and verified, they are integrated to form the calibration status information. This information fully includes the presence or absence of accumulation after calibration, accumulation location, accumulation degree, and accumulation range, which can characterize the true state of local accumulation of impurities inside the separation column, as well as the development trend of the overall blockage process.

[0079] Step 551: Based on the calibrated state information, determine the starting point and dynamic evolution step size of the nonlinear time-series extrapolation. Calculate the pressure drop value of the separation column at each of the future continuous sampling times to obtain the pressure drop prediction value at each time point. All pressure drop prediction values ​​constitute a complete pressure drop increase trajectory. Specifically, this includes: the calibrated state information obtained in step 550. The core purpose of the extrapolation operation is to extrapolate the upward trajectory of the axial pressure difference of the separation column over a future period based on the current and historical pressure drop change patterns, combined with the calibrated local accumulation state. This intuitively reflects the development speed of the blockage process and provides a basis for threshold comparison and early warning. The specific operation is divided into four steps, each closely connected, and incorporates the complete implementation of the fusion surface patch stitching algorithm. The first step is to determine the starting point and dynamic evolution step size of the nonlinear time-series extrapolation. The extrapolation starting point is strictly set to the current sampling time. The base pressure drop value at this time uses the measured axial pressure value of the column after data cleaning, anomaly removal, and smoothing in step 1. To ensure the accuracy and reliability of the starting point data, free from noise interference and abnormal deviations, the system retrieves the calibrated accumulation level at the current moment. Based on the accumulation level, the dynamic evolution step size is determined. The evolution step size is negatively correlated with the accumulation level; the more severe the accumulation, the faster the blockage evolves, and the shorter the evolution step size. Specifically, when the calibrated accumulation level is 0.2-0.4 (slight accumulation), the dynamic evolution step size is set to 5 sampling periods. Combined with a sampling frequency of 50Hz, the corresponding time interval is 0. The dynamic evolution step size is set to 0.1 seconds; when the accumulation level after calibration is 0.4-0.7 (moderate accumulation), the dynamic evolution step size is set to 3 sampling periods, corresponding to a time interval of 0.06 seconds; when the accumulation level after calibration is 0.7-1.0 (severe accumulation), the dynamic evolution step size is set to 1 sampling period, corresponding to a time interval of 0.02 seconds. This setting ensures both the timeliness of extrapolation, enabling timely capture of rapid congestion changes during severe accumulation, and the accuracy of calculation, avoiding trajectory distortion caused by excessively large step sizes.

[0080] The second step involves preprocessing the data for the fusion surface patch stitching algorithm. To ensure extrapolation accuracy, the basic data is preprocessed first, extracting historical pressure drop data associated with the calibrated state information. Specifically, the measured axial pressure values ​​of the cylinder at the current moment and the previous 30 sampling moments are selected, totaling 31 data points. This number is chosen to cover sufficient historical variation patterns while avoiding computational redundancy caused by excessive data. During preprocessing, each of the 31 historical pressure drop data points is checked, and abnormal fluctuation data is removed, i.e., the difference between the data and two adjacent data exceeds a preset fluctuation threshold, which is set at 0.01 MPa. If abnormal fluctuation data exists, the average of the two adjacent valid data points is used to replace the abnormal data, ensuring the continuity and reliability of historical data. Subsequently, the 31 processed historical pressure drop data points are arranged sequentially in chronological order, and the sampling time corresponding to each data point is labeled to form a complete historical pressure drop time series dataset, which serves as the basic data source for the fusion surface patch stitching algorithm.

[0081] The third step involves the core implementation of the surface patch stitching algorithm (surface patch segmentation, fitting, and fusion stitching). Surface patch segmentation divides the preprocessed historical voltage drop time series dataset into segments according to time sequence. The voltage drop data at every 10 sampling times is considered as a data segment. The last remaining data point is merged with the last data segment, resulting in three consecutive surface patch data segments: sampling times 1-10, sampling times 11-20, and sampling times 21-31. Each surface patch corresponds to a continuous time interval, reflecting the voltage drop trend within the corresponding time period, such as slow increase, rapid increase, or stable fluctuation. During segmentation, it is ensured that there is one data point overlap between adjacent surface patches. That is, the 10th data point belongs to both the first and second surface patches, and the 20th data point belongs to both the second and third surface patches. This lays the foundation for subsequent surface patch fusion stitching and avoids discontinuities at the stitching points.

[0082] The surface patch fitting process performs smooth surface fitting on each data segment. The core purpose of the fitting is to accurately reflect the actual change pattern of pressure drop within each time period, while ensuring the smoothness of the surface patch and avoiding abrupt changes or distortions during the fitting process. During fitting, a point-by-point smoothing fitting method is adopted. For each data point within each surface patch, the slope of the fitted surface is adjusted by combining the values ​​of the two data points before and after it, so that the fitted surface can pass through the vicinity of each data point, and the slope of the fitted curve between adjacent data points transitions smoothly without obvious abrupt changes. After the fitting is completed, each data segment corresponds to an independent smooth surface patch, and each surface patch can accurately reflect the change trend of pressure drop within the corresponding time period. For example, the first surface patch corresponds to the slow increase phase of pressure drop in the early historical period, the second surface patch corresponds to the accelerated increase phase of pressure drop, and the third surface patch corresponds to the recent rapid increase phase of pressure drop.

[0083] The surface patch fusion and stitching process merges three independently fitted surface patches to form a complete and smooth pressure drop variation surface. The core of the stitching is to eliminate the connection traces between adjacent surface patches, ensuring that the stitched surface can completely and continuously reflect the historical pressure drop variation pattern and can smoothly extend to future moments. During stitching, the focus is on processing the overlapping data points of adjacent surface patches. The fitted values ​​of the overlapping data points on the two adjacent surface patches are calculated, and the average of the two fitted values ​​is taken as the final fitted value of the data point. At the same time, the slope at the junction of adjacent surface patches is adjusted so that the slope at the end of the previous surface patch is consistent with the slope at the beginning of the next surface patch. A gradient transition method is used to gradually adjust the slope difference until there is no obvious inflection point at the junction, ensuring that the stitched surface is smooth, continuous, without discontinuities or abrupt changes. After the fusion and stitching is completed, a complete historical pressure drop variation surface is obtained.

[0084] The fourth step involves nonlinear temporal extrapolation and point-by-point calculation. Based on the complete pressure drop change surface after fusion and splicing, the estimated pressure drop value for future continuous sampling moments is calculated point by point according to the set dynamic evolution step size. During the calculation, the change trend of the surface is strictly followed. Starting from the surface position corresponding to the current moment, the calculation is gradually extended to future moments according to the dynamic evolution step size. For each extension step size, the corresponding value on the surface at that moment is read as the estimated pressure drop value for that moment. At the same time, combined with the accumulation degree in the calibrated state information, each estimated value is fine-tuned. The fine-tuning rule is that for every 0.1 increase in the accumulation degree after calibration, the increase in the estimated pressure drop value increases by 5%, ensuring that the estimated value conforms to the actual congestion evolution law. For example, if the accumulation degree after calibration is 0.6 (moderate to severe), the estimated value of the surface extension at a certain moment is 0.32 MPa, and the fine-tuned value is 0.32 × (1 + 0.6 × 5%) = 0.32 × 1.03 = 0.3296 MPa.

[0085] The duration of point-by-point extrapolation is determined based on the degree of accumulation after calibration. For slight accumulation, the estimated pressure drop for the next 24 hours is extrapolated; for moderate accumulation, the estimated pressure drop for the next 12 hours is extrapolated; and for severe accumulation, the estimated pressure drop for the next 6 hours is extrapolated. This ensures that the extrapolation duration can meet the warning preparation time while avoiding unnecessary computational redundancy. After point-by-point extrapolation is completed, the estimated pressure drop for all future sampling times is arranged in chronological order, and the future time corresponding to each estimated value is marked to form a complete trajectory of pressure drop increase.

[0086] Step 552 involves comparing the pressure drop increase trajectory sequentially with each threshold level in a preset safe operating threshold sequence, recording the time point at which the pressure drop increase trajectory first crosses each threshold level and the magnitude of exceeding the threshold, thus forming a threshold comparison result. Specifically, this includes the pressure drop increase trajectory obtained in step 551, and simultaneously retrieving the preset safe operating threshold sequence. This threshold sequence is based on the rated operating parameters of the separator column, industrial safety operating standards, and long-term engineering practice experience, comprehensively covering the entire range of the separator column from normal operation to severe blockage, and is divided into three distinct threshold levels, starting from low... The warning thresholds, from highest to lowest, are: mild warning threshold, moderate warning threshold, and severe warning threshold. Each threshold level corresponds to a fixed axial pressure difference value for the column. The specific setting standards are as follows: the mild warning threshold is set at 50% of the rated axial pressure difference of the separation column. For example, if the rated axial pressure difference of the separation column is 0.3 MPa, the mild warning threshold is 0.15 MPa; the moderate warning threshold is set at 100% of the rated axial pressure difference, i.e., 0.3 MPa; and the severe warning threshold is set at 150% of the rated axial pressure difference, i.e., 0.45 MPa. This setting meets the requirements for safe industrial operation and can distinguish the warning levels of different degrees of blockage.

[0087] The threshold comparison process is carried out strictly in chronological order to ensure consistency between the pressure drop increase trajectory and the safe operation threshold sequence in the time dimension. The specific operation steps are as follows: sort the pressure drop increase trajectory in chronological order and mark the future time corresponding to each pressure drop estimate; at the same time, sort the preset safe operation threshold sequence from low to high level, clarify the specific value of each threshold level, and ensure that the comparison order is consistent and there are no omissions or confusion; starting from the starting time of the pressure drop increase trajectory (the current time), compare the pressure drop estimate at this time with each threshold level in the safe operation threshold sequence moment by moment. The comparison order is first the mild warning threshold, then the moderate warning threshold, and finally the severe warning threshold. The comparison of all three levels is completed at each moment to ensure that no threshold crossing is missed.

[0088] For each threshold level, the key point to record is the time when the estimated pressure drop first exceeds that threshold level, i.e., the first crossing point. The exact moment of this recording should be precisely marked. Simultaneously, the difference between the estimated pressure drop at that moment and the corresponding threshold level is calculated as the magnitude of exceeding the threshold. The magnitude is calculated by subtracting the threshold level value from the estimated pressure drop. The larger the difference, the more severe the congestion, and the magnitude of the exceedance is positive. If no exceedance is achieved, a negative value is recorded (not recorded). For example, the pressure drop increase trajectory is within the next 10 hours (25...). At 10 hours, 25 minutes, and 30 seconds, the estimated pressure drop is 0.152 MPa, exceeding the mild warning threshold of 0.15 MPa for the first time. The first crossing time is recorded as 10 hours, 25 minutes, and 30 seconds, with an exceedance of 0.152 - 0.15 = 0.002 MPa. At 18 hours, 10 minutes, and 15 seconds, the estimated pressure drop is 0.305 MPa, exceeding the moderate warning threshold of 0.3 MPa for the first time. The first crossing time is recorded as 18 hours, 10 minutes, and 15 seconds, with an exceedance of 0.005 MPa.

[0089] If the pressure drop increase trajectory does not cross a certain threshold level, for example, only crosses the mild or moderate thresholds but does not reach the severe threshold, then the failure to cross that threshold level is clearly recorded, and the reason for not crossing is noted, such as the threshold not being reached within the estimated time. If the pressure drop increase trajectory directly crosses a certain threshold level, for example, jumping directly from 0.14 MPa below the mild warning threshold to 0.31 MPa above the moderate warning threshold, then the first crossing time of the crossed threshold level is recorded as the moment of the jump to that threshold. At the same time, the difference between the estimated pressure drop at that moment and the crossed threshold level is calculated as the exceedance magnitude. If the pressure drop increase trajectory exceeds the same threshold level at multiple times, only the first crossing time and the corresponding exceedance magnitude are recorded. Subsequent exceedances are not recorded again to avoid redundancy. All recorded information, including the first crossing time, exceedance magnitude, and failure to cross for each threshold level, is organized and summarized from low to high threshold level to form a complete threshold comparison result.

[0090] Step 553: Based on the threshold comparison results, identify the highest warning level reached by the pressure drop increase trajectory, and the urgency of reaching that level, to obtain a graded congestion warning instruction including level identifier and time information; specifically, this includes: the threshold comparison results obtained in step 552. The core purpose is to identify the highest warning level of the current congestion process based on the comparison results, determine the urgency of reaching that level, and generate a graded congestion warning instruction that can be directly used for industrial site safety management. The specific operation is as follows: Based on the crossing situation of each threshold level in the threshold comparison results, determine the highest warning level reached by the pressure drop increase trajectory. The identification rule is clear and unique. Specifically, if the threshold comparison results only record the first... If the first crossing information does not record crossing information at the moderate or severe warning thresholds (i.e., it does not cross the moderate or severe thresholds), the highest warning level is determined to be a mild warning. If the threshold comparison results record the first crossing information at the mild or moderate warning thresholds but not at the severe warning threshold, the highest warning level is determined to be a moderate warning. If the threshold comparison results record the first crossing information at all three threshold levels (mild, moderate, and severe), the highest warning level is determined to be a severe warning. If the threshold comparison results do not record crossing information at any threshold level (i.e., the pressure drop trajectory is always below the mild warning threshold), it is determined to be no warning, no warning command needs to be generated, only the current congestion process is recorded as being within the normal range, and regular monitoring is required.

[0091] The urgency of intervention is determined by the time difference between the current moment and the first crossing time of the highest warning level, combined with the exceedance magnitude corresponding to the highest warning level. The specific criteria are as follows: If the highest warning level is a mild warning, and the time difference between the current moment and the first crossing time of the mild warning is greater than 24 hours, it is considered a relaxed timeframe, indicating a slow blockage process, and monitoring the separation column's operating status every 6 hours is recommended; if the time difference is between 12 and 24 hours, it is considered a relatively urgent timeframe, and monitoring every 3 hours is recommended to prepare for intervention; if the time difference is less than 12 hours, it is considered an urgent timeframe, and monitoring every hour is recommended to prepare for mild intervention measures, such as adjusting the medium flow rate or increasing the backwashing frequency; if the highest warning level is a moderate warning, when... If the time difference between the current moment and the first crossing time of the moderate warning is greater than 12 hours, it is considered that the time is relatively lenient, and monitoring should be carried out every 2 hours; if the time difference is between 6 and 12 hours, it is considered that the time is urgent, and monitoring should be carried out every 30 minutes, and moderate intervention measures should be implemented immediately, such as adjusting the medium temperature and initiating local backwashing; if the time difference is less than 6 hours, it is considered that the time is extremely urgent, and real-time monitoring should be carried out, and preparations for shutdown should be made; if the highest warning level is a severe warning, if the time difference between the current moment and the first crossing time of the severe warning is greater than 6 hours, it is considered that the time is urgent, and real-time monitoring should be implemented, and severe intervention measures should be implemented immediately, such as shutdown for cleaning and replacement of the filling medium; if the time difference is less than 6 hours, it is considered that the time is extremely urgent, and immediate shutdown should be carried out to avoid equipment failure and process fluctuations caused by separation column blockage.

[0092] Meanwhile, if the exceedance magnitude corresponding to the highest warning level exceeds 10% of the threshold level, the time urgency level will be increased by one level. For example, the exceedance magnitude of a moderate warning is 0.03 MPa (the moderate warning threshold is 0.3 MPa, and 10% is 0.03 MPa). In this case, the time urgency level will be increased by one level. If the time difference is 12-24 hours, it will be increased to time urgency. The highest warning level, time urgency, the first crossing time point of each level, the exceedance magnitude, and the corresponding intervention suggestions will be integrated to generate a graded congestion warning instruction. The instruction specifically includes the following core information: warning level identifier (mild warning / moderate warning / severe warning / no warning), the first crossing time point of the highest warning level, the time urgency level, the exceedance magnitude of each level, and targeted intervention suggestions. The intervention suggestions are set in combination with the warning level and time urgency to ensure that the suggestions are operable.

[0093] By incorporating a fusion surface patch stitching algorithm into nonlinear time-series extrapolation, the inference accuracy of the pressure drop rise trajectory is improved. Through multi-threshold hierarchical comparison and time urgency judgment, a graded blockage warning command is generated to achieve accurate calibration of the separation column blockage and improve the safety and stability of the separation column operation.

[0094] like Figure 2As shown, embodiments of the present invention also provide a separation column blockage early warning system based on time-series data, including:

[0095] The module is used to process multidimensional time-series monitoring data to obtain a sequence of feature vectors;

[0096] The identification module is used to calculate the dynamic mutual information correlation degree between the rate of change of axial pressure difference of the column and the fluctuation coefficient of outlet volume flow rate within the sliding time window based on the feature vector sequence, and compare the dynamic mutual information correlation degree with the preset dynamic threshold to obtain the initial abnormal features.

[0097] The discrimination module is used to extract the corresponding shell vibration signal based on the time window determined by the initial anomaly features and perform frequency domain transformation to obtain the energy spectral density evolution trend; the energy spectral density evolution trend is matched with the initial anomaly features to obtain the local accumulation state and spatial distribution probability;

[0098] The generation module is used to input the local accumulation state and spatial distribution probability into the pre-trained discretized state-space model to obtain the predicted hydraulic state value. The predicted hydraulic state value is compared with the actual monitoring value to generate a state estimation error vector. The state estimation error vector is mapped to the drag potential energy deviation sequence to establish a dynamic equilibrium equation, calculate the correction gain matrix and generate a state-space correction vector. The state-space correction vector is used to adaptively update the discretized state-space model to obtain the regional flow resistance distortion correction parameters.

[0099] The early warning module is used to calibrate the local accumulation state using regional flow resistance distortion correction parameters, obtain calibrated state information, and perform nonlinear time-series extrapolation calculation based on the calibrated state information to obtain the pressure drop increase trajectory; the pressure drop increase trajectory is compared with the preset safe operation threshold sequence to obtain a graded blockage early warning command.

[0100] It should be noted that this system is a system corresponding to the above method. All implementation methods in the above method embodiments are applicable to this embodiment and can achieve the same technical effect.

[0101] The above description represents the preferred embodiments of the present invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.

Claims

1. A method for early warning of column blockage based on time-series data, characterized in that, The method includes: Step 1: Real-time acquisition of multi-dimensional time-series monitoring data under the operating conditions of the separation column. This data includes inlet static pressure, outlet volumetric flow rate, column axial pressure difference, medium temperature, and shell vibration acceleration signals. The multi-dimensional time-series monitoring data is processed to obtain a feature vector sequence. This process includes data cleaning of the real-time acquired multi-dimensional time-series monitoring data to remove abnormal jumps and noise interference, resulting in pre-processed monitoring data. A sliding time window of preset length is used to truncate the pre-processed monitoring data, obtaining data segments within continuous time windows. Time-domain statistical features are extracted from each data segment within a time window, including at least the mean, variance, peak value, rate of change, and fluctuation coefficient. The extracted time-domain statistical features from each time window are normalized and combined in chronological order to form a feature vector sequence. Step 2: Based on the feature vector sequence, calculate the dynamic mutual information correlation degree between the rate of change of axial pressure difference of the column and the fluctuation coefficient of outlet volume flow rate within the sliding time window, and compare the dynamic mutual information correlation degree with the preset dynamic threshold to obtain the initial abnormal features. Step 3: Based on the time window determined by the initial anomaly features, extract the corresponding shell vibration signal and perform frequency domain transformation to obtain the energy spectral density evolution trend; match the energy spectral density evolution trend with the initial anomaly features to obtain the local accumulation state and spatial distribution probability. Step 4: Input the local accumulation state and spatial distribution probability into the pre-trained discretized state-space model to obtain the predicted hydraulic state value. Compare the predicted hydraulic state value with the actual monitoring value to generate a state estimation error vector. Map the state estimation error vector to a resistance potential energy deviation sequence to establish a dynamic equilibrium equation, calculate the correction gain matrix and generate a state-space correction vector. Use the state-space correction vector to adaptively update the discretized state-space model to obtain the regional flow resistance distortion correction parameters. Step 5: Use the regional flow resistance distortion correction parameter to calibrate the local accumulation state, obtain the calibrated state information, and perform nonlinear time-series extrapolation based on the calibrated state information to obtain the pressure drop increase trajectory; compare the pressure drop increase trajectory with the preset safe operation threshold sequence to obtain the graded blockage early warning command.

2. The separation column blockage early warning method based on time-series data according to claim 1, characterized in that, Based on the feature vector sequence, the dynamic mutual information correlation degree between the rate of change of axial pressure difference in the column and the fluctuation coefficient of outlet volumetric flow rate within the sliding time window is calculated. The dynamic mutual information correlation degree is compared with a preset dynamic threshold to obtain initial abnormal features, including: The feature components corresponding to the rate of change of axial pressure difference of the column are extracted from the feature vector sequence to form the time series of pressure difference change rate. At the same time, the feature components corresponding to the outlet volume flow rate fluctuation coefficient are extracted to form the time series of flow rate fluctuation coefficient. The time series of differential pressure change rate and the time series of flow fluctuation coefficient are segmented by using a sliding time window of preset length to obtain multiple consecutive time window segments of differential pressure change rate and flow fluctuation coefficient. For each time window, the dynamic mutual information correlation degree between the pressure difference change rate window segment and the flow fluctuation coefficient window segment is calculated to obtain the dynamic mutual information correlation degree sequence; Each value in the dynamic mutual information correlation sequence is compared with a preset dynamic threshold. If the correlation exceeds the dynamic threshold, the time window is marked as a potential abnormal window. All potential anomaly windows are merged and filtered to remove isolated short-term anomalies. The time intervals of continuous anomalies and their corresponding correlation deviations are extracted to obtain the initial anomaly features.

3. The separation column blockage early warning method based on time-series data according to claim 2, characterized in that, The initial abnormal features include at least the abnormal start time, abnormal end time, abnormal peak intensity, and abnormal duration.

4. The separation column blockage early warning method based on time-series data according to claim 3, characterized in that, Determining the local accumulation state and spatial distribution probability of impurities in the packing layer or swirling chamber wall inside the separation column, including: Based on the abnormal start time and abnormal end time included in the initial abnormal characteristics, the shell vibration acceleration signal within the corresponding time period is extracted from the multi-dimensional time-series monitoring data to obtain vibration signal segments; The vibration signal segment is transformed in the frequency domain to convert the time domain signal into a spectral distribution, and the energy spectral density in different frequency intervals is calculated to obtain the energy spectral density evolution sequence. The evolution trend of vibration characteristics is formed by extracting the main frequency migration trajectory, the change in the proportion of high-frequency energy and the energy mutation points in specific frequency bands from the energy spectral density evolution sequence. The evolution trend of vibration characteristics is matched with the abnormal peak intensity and duration of abnormality in the initial abnormal characteristics in a spatiotemporal correlation. Based on the alignment of the energy mutation point of a specific frequency band with the abnormal peak intensity on the time axis and the correlation between the rising trend of high frequency energy proportion and the duration of abnormality, it is determined whether there is local accumulation of solid particles in the time window, and the local accumulation state is obtained. Based on the local accumulation state, and according to the spatial directivity of the dominant frequency migration trajectory in the vibration characteristic evolution trend and the phase difference of vibration signals at different measuring points, the spatial distribution probability of the accumulation area in the radial and axial directions of the column is calculated, and the spatial distribution probability is obtained.

5. The separation column blockage early warning method based on time-series data according to claim 4, characterized in that, Step 4 includes: Based on the local accumulation state and spatial distribution probability, combined with the operating parameters at the current moment, the input vector of the discretized state space model is constructed. The input vector is then substituted into the pre-trained discretized state space model to perform forward recursive calculation, thereby obtaining the predicted hydraulic state value at the current moment. The predicted hydraulic state value is compared dimension by dimension with the measured values ​​of the column axial pressure difference and the outlet volumetric flow rate at the same time in the pre-processed actual monitoring data. The deviation of each dimension is calculated and combined to obtain the state estimation error vector. The pressure difference deviation and flow rate deviation in the state estimation error vector are converted into the change in drag potential energy distributed along the axis of the column, forming a drag potential energy deviation sequence, and a dynamic equilibrium equation is constructed based on the drag potential energy deviation sequence. Solve the dynamic equilibrium equation to obtain the correction gain matrix. Use the correction gain matrix to perform a weighted transformation on the state estimation error vector to obtain the state space correction vector. The internal state parameters of the discretized state-space model are adaptively updated using the state-space correction vector to obtain the updated model parameters, which are then used as parameters for regional flow resistance distortion correction.

6. The separation column blockage early warning method based on time-series data according to claim 5, characterized in that, The predicted hydraulic state values ​​include at least the predicted axial pressure difference of the column and the predicted outlet volumetric flow rate.

7. The separation column blockage early warning method based on time-series data according to claim 6, characterized in that, The local accumulation state is calibrated using the regional flow resistance distortion correction parameter to obtain the calibrated state information. Based on the calibrated state information, a nonlinear time-series extrapolation operation is performed to obtain the pressure drop increase trajectory. The pressure drop increase trajectory is compared with a preset safe operating threshold sequence to obtain graded congestion early warning instructions, including: Based on the regional flow resistance distortion correction parameters, multiple local characteristic parameters included in the local accumulation state are corrected respectively to eliminate the deviation introduced by local flow resistance distortion and obtain calibrated state information characterizing the overall blockage process of the separation column. Based on the calibrated state information, the starting point and dynamic evolution step size of the nonlinear time series extrapolation are determined, and the pressure drop value of the separation column at each future continuous sampling time is calculated point by point to obtain the pressure drop prediction value at each time. All the pressure drop prediction values ​​constitute a complete pressure drop increase trajectory. The pressure drop increase trajectory is compared one by one with each threshold level in the preset safe operation threshold sequence according to the time sequence. The time point when the pressure drop increase trajectory first crosses each threshold level and the magnitude of exceeding the threshold are recorded to form the threshold comparison result. Based on the threshold comparison results, the highest warning level reached by the pressure drop increase trajectory is identified, as well as the urgency of reaching that level, and a graded congestion warning instruction including level identifier and time information is obtained.

8. A separation column blockage early warning system based on time-series data, wherein the system implements the method as described in any one of claims 1 to 7, characterized in that, include: The module is used to process multidimensional time-series monitoring data to obtain a sequence of feature vectors; The identification module is used to calculate the dynamic mutual information correlation degree between the rate of change of axial pressure difference of the column and the fluctuation coefficient of outlet volume flow rate within the sliding time window based on the feature vector sequence, and compare the dynamic mutual information correlation degree with the preset dynamic threshold to obtain the initial abnormal features. The discrimination module is used to extract the corresponding shell vibration signal based on the time window determined by the initial anomaly features and perform frequency domain transformation to obtain the energy spectral density evolution trend; the energy spectral density evolution trend is matched with the initial anomaly features to obtain the local accumulation state and spatial distribution probability; The generation module is used to input the local accumulation state and spatial distribution probability into the pre-trained discretized state-space model to obtain the predicted hydraulic state value. The predicted hydraulic state value is compared with the actual monitoring value to generate a state estimation error vector. The state estimation error vector is mapped to the drag potential energy deviation sequence to establish a dynamic equilibrium equation, calculate the correction gain matrix and generate a state-space correction vector. The state-space correction vector is used to adaptively update the discretized state-space model to obtain the regional flow resistance distortion correction parameters. The early warning module is used to calibrate the local accumulation state using regional flow resistance distortion correction parameters, obtain calibrated state information, and perform nonlinear time-series extrapolation calculation based on the calibrated state information to obtain the pressure drop increase trajectory; the pressure drop increase trajectory is compared with the preset safe operation threshold sequence to obtain a graded blockage early warning command.

Images