A method for monitoring the operating state of a sludge drying system

CN122797162APending Publication Date: 2026-09-22WUHAN BIQING ENVIRONMENTAL PROTECTION TECH CO LTD
View PDF 0 Cites 0 Cited by

Patent Information

Application Number
CN202611111475.7
Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
Filing Date
2026-07-24
Publication Date
2026-09-22

AI Technical Summary

Technical Problem

[0004]根据本发明的第一方面,提出一种污泥干化系统运行状态监测方法,用于解决现有方法在监测逻辑与状态评估机制上面临的问题,包括如下步骤:

Benefits of technology

[0007]本发明将污泥干化系统的过程变量划分为强时序性热工变量和弱时序性物料变量,符合工业过程变量时序特性差异明显的实际情况。通过构建滞后互相关积分核函数,利用变量与关键性能指标之间的关联程度对距离计算进行加权,克服了传统核函数中各维度权重均等、关键变量贡献不突出的缺陷;通过第一路投影系数调制第二路核函数,建立了热工状态对物料状态演化过程的融合约束关系。此外,本发明利用最大李雅普诺夫指数表征系统不稳定性,结合柯西施瓦茨散度提取融合偏离度指标,并在二维监测平面内构建非等轴高斯混合模型置信边界进行异常判定。由此能够识别系统内部非线性及混沌演化特征,放大微弱异常信号,提高污泥干化系统运行状态监测的灵敏度、准确性和稳定性。

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN122797162A_ABST
    Figure CN122797162A_ABST
Patent Text Reader

Abstract

The application provides a sludge drying system operation state monitoring method, including standardizing normal historical data sets, dividing process variables into a strong time sequence thermal variable subset and a weak time sequence material variable subset; constructing a lag cross-correlation integral kernel, weighting distance according to the correlation integral of variables and performance indicators, and performing first path kernel principal component analysis on thermal variables to obtain projection coefficients; modulating a second path kernel function with the projection coefficients, reducing the dimension of material variables, and constructing a time window phase space trajectory after double projection on online data; calculating the maximum Lyapunov index of the two trajectories to represent system instability, combining the Cauchy-Schwarz divergence of the projection coefficient norm sequence to obtain a fusion deviation degree; based on the fusion system instability representation quantity and the fusion deviation degree, a two-dimensional monitoring plane and a non-equiaxed Gaussian mixture model confidence boundary are constructed, and when the coordinate point is out of boundary, it is determined that the sludge drying system is abnormal, and thermal fluctuation abnormalities and material fusion mismatch abnormalities are distinguished.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This application belongs to the field of sludge drying, and in particular relates to a method for monitoring the operating status of a sludge drying system. Background Technology

[0002] With the continuous advancement of urbanization and industrialization, sludge treatment and disposal has become a crucial issue in environmental protection. Sludge drying is a key step in achieving sludge reduction, harmlessness, and resource utilization; its operational safety and stability directly impact the efficiency and cost of the entire treatment process. The sludge drying system is a typical nonlinear industrial process, its operating state influenced by various process parameters and external disturbances. The process variables within the system exhibit significant differences in their temporal characteristics, typically categorized into thermal variables such as temperature, pressure, and gas flow rate, which possess strong temporal evolution characteristics, and material variables such as sludge moisture content and feed / discharge rates, which are significantly affected by feed batches, moisture content fluctuations, and process loads, and do not exhibit sustained thermal memory characteristics within short time windows. These two types of variables differ in their evolution rates and response mechanisms and are coupled during continuous operation. Fluctuations in the thermal system's state affect the material drying efficiency and operating trajectory, while changes in material properties and load also negatively impact the energy balance of the thermal system. Therefore, precise monitoring of the operational status of the sludge drying system is an important requirement for preventing equipment failures and ensuring stable system operation.

[0003] Traditional monitoring technologies typically treat all process variables as a single high-dimensional set, failing to adequately consider the differences in characteristic attributes and time scales between strongly time-series thermal variables and weakly time-series material variables. They also struggle to identify the lagged correlations between different variables and key system performance indicators, resulting in core feature extraction results that fail to accurately reflect the underlying operational patterns of the system. Existing methods generally lack effective identification of interactions between subsystems, making it difficult to characterize the modulating influence of the thermal system state on the material system evolution process, and to measure the degree of imbalance and deviation in the fusion relationship between the thermal and material systems. Furthermore, the abnormal evolution of sludge drying systems is often accompanied by changes in internal nonlinear dynamic characteristics. Conventional monitoring indicators based on single static statistical features are insufficient to identify the deep-seated chaotic instability evolution patterns generated by external disturbances or internal degradation, leading to inaccurate anomaly judgment confidence boundary construction and a tendency to generate false alarms and missed alarms under varying operating conditions. This fails to meet the requirements of modern sludge drying processes for high-sensitivity and high-reliability state monitoring. Summary of the Invention

[0004] According to a first aspect of the present invention, a method for monitoring the operational status of a sludge drying system is proposed to address the problems faced by existing methods in terms of monitoring logic and status assessment mechanisms, comprising the following steps:

[0005] The original dataset of historical normal operation is obtained and standardized to obtain the modeling dataset. The process variables in the modeling dataset are divided into a subset of strong time-series thermal variables and a subset of weak time-series material variables. A lag cross-correlation integral kernel function is constructed, and the difference components of different dimensions are weighted when calculating the distance between sample points. The weights are determined by the absolute value integral of the lag cross-correlation function between the corresponding variable and the key performance index within a preset lag time interval. The first-way kernel principal component analysis is performed on the subset of strongly time-series thermal variables using the lag cross-correlation integral kernel function to obtain the first-way projection coefficients. A second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficient is constructed. The second kernel function is used to perform second kernel principal component analysis on the weakly temporal material variable subset to obtain the second projection coefficient. The two projection coefficients are calculated on the online real-time data and the two phase space trajectories are constructed within a preset time window. The maximum Lyapunov exponent of the two phase space trajectories within a preset time window is calculated as the instability representation of the thermal system and the fusion system. The Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients is calculated. The fusion deviation index is obtained by amplifying the thermal system instability representation after dimensionless transformation by the preset time scale parameter. A non-isoaxial Gaussian mixture model confidence boundary is established in the two-dimensional monitoring plane composed of the fusion system instability representation and the fusion deviation index. Anomalies are determined when the coordinate points fall outside the boundary.

[0006] According to a second aspect of the present invention, a sludge drying system operation status monitoring system is provided, comprising the following modules: The partitioning module is used to obtain the original dataset of historical normal operation and standardize it to obtain the modeling dataset. The process variables in the modeling dataset are divided into a subset of strong time-series thermal variables and a subset of weak time-series material variables. The module is used to construct a lag cross-correlation integral kernel function, which weights the difference components of different dimensions when calculating the distance between sample points. The weights are determined by the absolute value integral of the lag cross-correlation function between the corresponding variable and the key performance index within a preset lag time interval. The lag cross-correlation integral kernel function is used to perform a first-way kernel principal component analysis on the subset of strongly time-series thermal variables to obtain the first-way projection coefficients. The calculation module is used to construct a second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficient, use the second kernel function to perform second kernel principal component analysis on the weak temporal material variable subset to obtain the second projection coefficient, calculate the two projection coefficients on the online real-time data, and construct the two phase space trajectories within a preset time window; The determination module is used to calculate the maximum Lyapunov exponent of the two phase space trajectories within a preset time window as the instability representation of the thermal system and the fusion system. It calculates the Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients and amplifies the instability representation of the thermal system after it has been dimensionlessly transformed by the preset time scale parameter to obtain the fusion deviation index. A non-isoaxial Gaussian mixture model confidence boundary is established in the two-dimensional monitoring plane composed of the instability representation of the fusion system and the fusion deviation index. Anomalies are determined when the coordinate points fall outside the boundary.

[0007] This invention divides the process variables of a sludge drying system into strongly time-series thermal variables and weakly time-series material variables, reflecting the significant differences in the time-series characteristics of industrial process variables. By constructing a lagged cross-correlation integral kernel function, the distance calculation is weighted based on the correlation between variables and key performance indicators, overcoming the shortcomings of traditional kernel functions where all dimensions have equal weights and key variables do not contribute significantly. A second-path kernel function is modulated by the first-path projection coefficient, establishing a fusion constraint relationship between thermal state and material state evolution. Furthermore, this invention uses the maximum Lyapunov exponent to characterize system instability, combines Cauchy-Schwarz divergence to extract a fusion deviation index, and constructs a non-equiaxial Gaussian mixture model confidence boundary in a two-dimensional monitoring plane for anomaly detection. This enables the identification of nonlinear and chaotic evolution characteristics within the system, amplifies weak anomaly signals, and improves the sensitivity, accuracy, and stability of sludge drying system operation status monitoring. Attached Figure Description

[0008] Figure 1 A flowchart of the first embodiment; Figure 2 A schematic diagram of the hysteresis cross-correlation curve between thermal variables and moisture content at the drying outlet; Figure 3 This is a schematic diagram of the scatter points of samples and the distribution of non-isoaxial confidence ellipses in a two-dimensional monitoring plane. Detailed Implementation

[0009] The technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of this application, and not all embodiments. Based on the embodiments of this application, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this application.

[0010] In the first embodiment, the present invention proposes a method for monitoring the operating status of a sludge drying system, such as... Figure 1 This includes the following steps: S1. Obtain the original dataset of historical normal operation and standardize it to obtain the modeling dataset. Divide the process variables in the modeling dataset into a subset of strong time-series thermal variables and a subset of weak time-series material variables.

[0011] Continuous operational data of the sludge drying system under stable and fault-free conditions were exported from the distributed control system as the original historical normal operation dataset. The data format was a two-dimensional matrix containing timestamps and multiple sensor readings. This dataset was read, and each column was then standardized using Z-scores to ensure that the mean of each variable was 0 and the variance was 1, removing the influence of dimensions to obtain the modeling dataset. For variable partitioning, the augmented Dickey-Fuller test was not used as a direct criterion for rejecting the null hypothesis for strongly time-series thermal variables or weakly time-series material variables; the augmented Dickey-Fuller test was only used to determine whether there was a significant trend in the original sequence, and to assist in detrending or smoothing preprocessing when necessary. In actual classification, the decay rate of the autocorrelation function of each process variable's own time series was uniformly used to determine the decorrelation time scale, and the relationship between this decorrelation time scale and the time scale classification threshold was used to complete the variable partitioning.

[0012] Specifically, the autocorrelation function of each standardized process variable sequence under multiple positive time delays is calculated. The delay time corresponding to the first crossing of zero or the first time being less than or equal to a preset small threshold is found, and this delay time is used as the decorrelation time scale. Variables such as hot air temperature, carrier gas flow rate, and heating steam pressure with decorrelation time scales greater than the time scale classification threshold are classified into the strongly time-series thermal variable subset, while variables such as feed moisture content, sludge flow rate, and feed / discharge rate with decorrelation time scales less than or equal to the time scale classification threshold are classified into the weakly time-series material variable subset. This avoids mistaking the cross-correlation relationship between the original sequence and the first-order difference sequence for the autocorrelation coefficient, and also avoids confusing the stationarity test results with the long-term memory determination.

[0013] In an optional embodiment, the step of obtaining the historical normal operation raw dataset and standardizing it to obtain the modeling dataset, and dividing the process variables in the modeling dataset into a strongly time-series thermal variable subset and a weakly time-series material variable subset, includes: Calculate the time series autocorrelation function of each process variable in the modeling dataset; Determine the delay time corresponding to the first crossing of zero or reaching a preset small threshold for the autocorrelation function of each process variable, and use the delay time as the decorrelation time scale for the corresponding process variable; Set time-scale classification thresholds; Process variables with a decorrelated time scale greater than the time scale classification threshold are classified into the strongly time-series thermal variables subset, and process variables with a decorrelated time scale less than or equal to the time scale classification threshold are classified into the weakly time-series material variables subset.

[0014] A standardized N-dimensional process variable data matrix is ​​collected and analyzed at M time points during the continuous operation of the sludge drying system. For the j-th process variable, such as heat transfer oil temperature or sludge feed rate, the time series autocorrelation function is calculated. The preferred time delay step is 1 to 50 minutes, and the sampling interval is 1 minute. Subsequently, along the direction of increasing time delay, the minimum delay time that makes the autocorrelation function value first less than or equal to a preset small threshold (preferably 0.05 or zero) is searched. This minimum delay time is set as the decorrelation time scale of the variable, representing the duration of the memory effect within a single variable. If a process variable does not show a situation where the autocorrelation function first crosses zero or is less than or equal to the preset small threshold within the preset positive delay search range of 1 to 50 minutes, the upper limit of the search range is taken as the lower bound of the minimum decorrelation time scale of the variable. In actual classification, the decorrelation time scale of the variable can be recorded as a preset extension value greater than the upper limit of the search range, such as 50 minutes plus a sampling interval, or further extended to the preset maximum search upper limit for calculation. If a variable fails to cross the preset maximum search limit, it is considered to have long-term memory characteristics and is assigned to the side where the relevant time scale is greater than the time scale classification threshold. This ensures that slowly decaying variables can still be classified even if they do not cross the threshold, preventing the variable partitioning process from being interrupted.

[0015] Based on expert experience or statistical characteristics of data related to sludge drying processes, a global timescale classification threshold is set. The preferred range for this threshold is 10 to 15 minutes, such as 12 minutes in this example. The decorrelation timescales of each of the N process variables are compared with the classification threshold. For example, if the decorrelation timescale of the heat transfer oil inlet temperature is 18 minutes, which is greater than the 12-minute threshold, then the inlet temperature is considered to have long-term memory and is classified into the strongly time-series thermal variable subset. If the decorrelation timescale of the sludge feed screw speed is 5 minutes, which is less than the 12-minute threshold, then the sludge feed screw speed is classified into the weakly time-series material variable subset. The variable names mentioned above are merely examples; the actual classification is based on the decorrelation timescales calculated from the corresponding variables in historical normal operation data. When the same type of variable exhibits different autocorrelation decay characteristics under different operating conditions, it is classified according to the comparison results of its calculated decorrelation timescale and the timescale classification threshold. This achieves a hard threshold segmentation of the entire mixed process variable space into two independent subspaces with different characteristics. This variable space partitioning mechanism avoids the technical defect in traditional full-scale dimensionality reduction methods where slow-varying features are masked by fast-varying high-frequency noise.

[0016] S2, construct a lagged cross-correlation integral kernel function, and weight the difference components of different dimensions when calculating the distance between sample points. The weights are determined by the absolute value integral of the lagged cross-correlation function of the corresponding variable and the key performance index within a preset lag time interval. The lagged cross-correlation integral kernel function is used to perform a first-way kernel principal component analysis on the subset of strongly time-series thermal variables to obtain the first-way projection coefficients.

[0017] The moisture content of the dried sludge outlet was selected as a key performance indicator. The positive lag cross-correlation function between each variable in the strongly time-series thermal variable subset and the key performance indicator sequence was calculated. The positive lag was defined as the time difference between the change in the thermal variable occurring first and the response of the dried sludge outlet moisture content occurring later. The preset lag time interval was uniformly set to 5 to 30 minutes, covering the process delay interval from the transmission of front-end thermal parameter changes to the end-end moisture content response. Negative lag intervals were not used in the weight calculation. The absolute value of the lag cross-correlation function within this time interval was integrated. The obtained integral value was divided by the sum of the integral values ​​of all dimensions and normalized before being used as the weight of the data distance difference component for the corresponding variable dimension. The variation of the cross-correlation between the variable and the dried sludge moisture content index with lag time is as follows: Figure 2 As shown, a hysteresis cross-correlation integral kernel function is constructed based on the radial basis function kernel function framework. In this kernel function, the weights obtained earlier are used to perform a weighted summation of the squared differences of the corresponding dimensions of any two input sample points. This summation is then multiplied by -1 / 2 and divided by the square of the kernel width parameter before performing a natural exponential operation to obtain the kernel matrix element values. A kernel principal component analysis algorithm is used to perform nonlinear mapping and eigenvalue decomposition on a subset of strongly time-series thermal variables, extracting the top principal components with a cumulative contribution rate greater than 85%. The score vectors of each sample point in the directions of these principal components are the first-path projection coefficients.

[0018] In an optional embodiment, the construction of the lagged cross-correlation integral kernel function, when calculating the distance between sample points, weights the difference components of different dimensions. The weights are determined by the absolute value integral of the lagged cross-correlation function between the corresponding variable and the key performance index within a preset lag time interval. The lagged cross-correlation integral kernel function is then used to perform a first-way kernel principal component analysis on the subset of strongly time-series thermal variables to obtain the first-way projection coefficients, including: For each variable in the strongly time-series thermal variable subset, calculate the sequence of cross-correlation functions between it and the moisture content index of the sludge drying system within the preset lag time interval. Integrate the absolute values ​​of each cross-correlation function sequence within the preset lag time interval to obtain the initial weights of each variable; The initial weights of all strongly time-series thermal variables are normalized to obtain the final weights of each variable. When calculating the Gaussian kernel function, the squared difference between two sample points in each variable dimension is multiplied by the corresponding final weight, and the sum of the squared weighted differences in all dimensions is substituted into the exponential function to obtain the value of the lagged cross-correlation integral kernel function. The first-way kernel principal component analysis is performed on the subset of strongly time-series thermal variables using the hysteresis cross-correlation integral kernel function to obtain the first-way projection coefficients.

[0019] For a subset of strongly time-series thermal variables, such as features comprising P variables, the cross-correlation function between the i-th variable and the dry mud moisture content index is calculated. The preset lag time interval is preferably set to [5, 30] minutes, covering the transmission delay from changes in the front-end thermal parameters to the final moisture content response. This lag time interval is consistent with the previous embodiment, representing a positive causal interval where the thermal variables lead and the moisture content index lags. The absolute values ​​of the cross-correlation function sequence within this time interval are discretely summed and integrated to obtain the initial weight of the i-th thermal variable. This initial weight is then normalized by dividing it by the sum of the initial weights of all P strongly time-series thermal variables to obtain the final weight coefficient. This step calculates the degree of influence of each thermal variable on the final moisture content time lag; for example, the weight of the main steam pressure is 0.35, while the weight of the cooling water temperature is 0.05.

[0020] Before constructing the spatial distance metric, the equal-weighted Euclidean distance is replaced with a weighted distance when performing kernel principal component dimensionality reduction calculations. For any two strongly temporally sequential feature sample points in the input feature space, when calculating the squared distance between them, the difference between them in the i-th dimension is first extracted and squared. Then, the squared difference is multiplied by the final weight corresponding to that dimension to obtain a single-dimensional weighted component. The sum of the weighted components of all P dimensions is calculated as the total weighted squared distance between the two sample points. The total weighted squared distance is then negatively divided by twice the square of the kernel width parameter and substituted into the exponential function to calculate the kernel function value. The kernel width parameter of the Gaussian kernel function is preferably determined by the mean of the weighted distance between samples under normal operating conditions, and the preferred value range is 0.5 to 5.0, for example, 2.1. This hysteresis cross-correlation integral kernel function reduces the information interference of weaker dimensions and improves the feature extraction capability of spatiotemporal correlation in the drying process.

[0021] S3, construct a second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficient, use the second kernel function to perform second kernel principal component analysis on the weakly temporal material variable subset to obtain the second projection coefficient, calculate the two projection coefficients on the online real-time data and construct the two phase space trajectories within a preset time window.

[0022] The second kernel function is set as a Gaussian radial basis kernel function whose kernel width parameter is modulated by the temporal fluctuation characteristics of the first projection coefficients, instead of directly multiplying the time-varying modulation factor by the entire kernel function value. Specifically, the first k principal components of the first projection coefficients are extracted within a preset modulation time window, and the sum of squares or variance of the differences between adjacent time points of the first k principal components is calculated as the characteristic value of thermal state fluctuation. This characteristic value of thermal state fluctuation is substituted into a monotonically decreasing exponential decay function to obtain the uniform adjustment factor corresponding to the modulation time window. The preset basic kernel width parameter is multiplied by the uniform adjustment factor to obtain the uniform kernel width parameter used by the second kernel function within the modulation time window.

[0023] To ensure the online closed-loop application of the core principal component analysis (KPI) model, during the offline training phase, thermal state fluctuation characteristic values ​​under multiple modulation time windows are calculated based on historical normal operation data. Several preset modulation intervals are then defined according to the range of thermal state fluctuation characteristic values ​​or the unified adjustment factor. For each preset modulation interval, a second-way KPI sub-model is trained using the corresponding unified kernel width parameter. The corresponding training sample centering parameters, kernel matrix eigenvectors, eigenvalues, and projection calculation parameters are saved, forming a second-way KPI model library. During online monitoring, the thermal state fluctuation characteristic values ​​and the unified adjustment factor are first calculated based on the first-way projection coefficient of the current time window. Then, a second-way KPI sub-model matching the current modulation interval is selected to project the current weakly time-series material variable data. Therefore, during the online phase, the kernel width parameter of the already trained single KPI model is not arbitrarily changed. Instead, a model with a consistent kernel width is selected from the pre-trained and saved model library to complete the projection, ensuring consistency between the training kernel space and the online projection kernel space. When the unified adjustment factor calculated online is lower than the minimum modulation interval or higher than the maximum modulation interval in the model library, it is truncated to the nearest preset modulation interval, and the second-path kernel principal component analysis sub-model corresponding to the modulation interval is called; when it is near the boundary of two adjacent modulation intervals, the nearest modulation interval is selected according to the preset nearest neighbor interval selection rule to ensure that the online projection always corresponds to the trained and saved kernel space.

[0024] For online real-time data accessed via industrial communication protocols, standardized parameters are first applied for standardization. The data is then projected onto a first-path kernel principal component analysis model and a second-path kernel principal component analysis sub-model corresponding to the current modulation interval, respectively, to calculate the projection coefficients of both paths at the current moment. A preset time window length is set to 50 to 120 sampling periods, preferably 100 sampling periods. Within this preset time window, the L2 norms of the first and second-path projection coefficient vectors are calculated at each moment, forming a first-path one-dimensional time window projection coefficient norm sequence and a second-path one-dimensional time window projection coefficient norm sequence. Based on the Tukens embedding theorem, the phase space of the two one-dimensional time window projection coefficient norm sequences is reconstructed. The optimal delay time is determined using the mutual information method, and the minimum embedding dimension is determined using the spurious nearest neighbor method. Based on the determined delay time and embedding dimension, the one-dimensional norm sequence is converted into a high-dimensional phase space point set. These point sets are then connected in chronological order to construct the two-path phase space trajectories corresponding to the preset time window length. The Tukens embedding theorem states that for a dynamic system, a high-dimensional phase space that reflects the nonlinear dynamic characteristics of the original system can be reconstructed from the delayed coordinates of a single time series.

[0025] In an optional embodiment, the construction of a second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficients, and the use of the second kernel function to perform second-way kernel principal component analysis on the weakly temporal material variable subset to obtain the second projection coefficients, includes: Extract the first k kernel principal components from the first projection coefficients; Calculate the sum of squares or variance of the differences between the first k principal components at adjacent time points, and use it as the characteristic value of thermal state fluctuation; Substituting the thermal state fluctuation characteristic value into a monotonically decreasing exponential decay function, the unified adjustment factor corresponding to the preset modulation time window is obtained. Multiply the preset base kernel width parameter by the unified adjustment factor to obtain the unified kernel width parameter of the second kernel function within the preset modulation time window; The unified kernel width parameter is used to construct a second kernel function, so that any sample point pair within the same kernel matrix uses the same kernel width parameter, and the kernel width of the second kernel function is reduced when the thermal state time-series fluctuations are enhanced.

[0026] Obtain the strongly temporally sequential projection coefficient vector output in real time from the first-path kernel principal component analysis model, and extract the coefficients of the first k principal components. The number of principal components k is preferably determined by ensuring that the cumulative variance contribution rate is greater than 85%. In actual operating conditions, k is usually taken as an integer from 3 to 8, for example, k=5. At the current monitoring time t, the statistical measure of the sum of squares or variance of the differences between adjacent times of the first k kernel principal components within the preset modulation time window is used as the characteristic value of thermal state fluctuation, instead of directly using the norm of a single current-time projection vector as the overall scaling factor of the kernel function. This characteristic value represents the degree to which the current thermal subsystem deviates from the steady-state operating point and generates temporally sequential fluctuations through the magnitude of state changes in the reduced-dimensional space. Using a monotonically decreasing exponential decay function, the thermal state fluctuation characteristic value is multiplied by a negative decay coefficient to calculate the adjustment factor. The preferred range of the decay coefficient is set to 0.01 to 0.1, for example, 0.05. This coefficient is used to control the sensitivity of the subsequent kernel width to thermal state fluctuations and avoid numerical abrupt changes.

[0027] For constructing the Gaussian kernel function of the second-path weakly time-series material variable subset, a constant base kernel width parameter is pre-calculated using the median of the Euclidean distance between samples in the historical normal material training set, for example, set to 3.5. Within each modulation time window, the base kernel width parameter is multiplied by a uniform adjustment factor calculated for that time window to generate a modulated uniform kernel width parameter. This kernel width parameter is then substituted into the Gaussian kernel function expression of the second-path material data. Within the same kernel matrix, any pair of sample points uses the same uniform kernel width parameter, thus ensuring the kernel matrix maintains symmetry and positive semi-definiteness. Between different modulation intervals, kernel space matching is achieved by pre-training and saving the second-path kernel principal component analysis sub-model under the corresponding kernel width parameter. When fluctuations in the thermal state cause an increase in energy eigenvalues, the adjustment factor decays, resulting in a decrease in the Gaussian kernel width of the second-path material subset. At this point, the second kernel function improves the resolution of the data sample spacing difference within the material subspace, enabling it to identify minute feed blockages or flow interruptions caused by thermal instability and easily masked by process noise, thus improving the sensitivity of anomaly identification under all operating conditions.

[0028] In an optional embodiment, the step of calculating the two projection coefficients from online real-time data and constructing the two phase space trajectories within a preset time window includes: Calculate the norms of the first and second projection coefficient vectors within the preset time window to form two one-dimensional time window projection coefficient sequences. The optimal delay time is calculated using the mutual information method for the two one-dimensional time window projection coefficient sequences, and the minimum embedding dimension is calculated using the spurious nearest neighbor method. Based on the corresponding optimal delay time and minimum embedding dimension, the two one-dimensional time window projection coefficient sequences are reconstructed in phase space to obtain the first phase space trajectory matrix and the second phase space trajectory matrix.

[0029] An online monitoring sliding time window is set, with a length L preferably ranging from 60 to 120 consecutive sampling points. For example, if the data sampling interval is set to 1 minute, L is set to 100, corresponding to 100 minutes of historical information. As the real-time data stream progresses, the window is updated using a sliding step size of 1 sampling point. For each discrete moment within the current time window, the L2 norm of the multidimensional projection coefficient vectors of the first and second outputs is calculated, and these are sequentially arranged to form a one-dimensional first-path norm sequence and a second-path norm sequence, both of length L. This step aggregates the multidimensional, dispersed kernel principal component features into two one-dimensional sequences representing the energy evolution of the subsystem through Euclidean length compression. For the extracted one-dimensional sequences, the phase space time delay parameter is solved independently using the mutual information method. The mutual information is calculated as the delay step size increases, and the time step size corresponding to the first local minimum of the curve is extracted as the optimal delay time. This parameter is typically between 2 and 10, for example, 4 for the first path and 3 for the second path.

[0030] After obtaining the optimal delay time parameters, the spurious nearest neighbor method is used to determine the minimum embedding dimension required for high-dimensional unfolding of the sequence. Starting with dimension m=2, the test space is constructed by gradually increasing the dimension. The proportion of the increase in distance between two adjacent sample points is calculated as the reconstruction dimension increases. When the proportion of the increase in distance at a point exceeds a preset threshold, such as 15, it is determined to be a spurious nearest neighbor caused by low-dimensional projection folding. The dimension is continuously increased until the overall proportion of spurious nearest neighbors falls below a critical value of 5%, for example, decreasing to 1%. The corresponding m value at this point is selected as the minimum embedding dimension, with a preferred range generally between 3 and 8 dimensions; for example, 5 dimensions for the first path and 4 dimensions for the second path. The window length L, the optimal delay time τ, and the minimum embedding dimension m satisfy L-(m-1)τ greater than the preset minimum reconstruction point threshold. When the calculated τ or m results in insufficient reconstruction points, the window length L is increased first, or τ or m is adjusted within a preset range to ensure that the phase space trajectory matrix can be effectively constructed. The matrix transformation is performed based on the Packard-Takens delay coordinate reconstruction theorem: using the extracted optimal delay time and minimum embedding dimension as reconstruction parameters, a multidimensional row vector sequence is constructed by extracting elements of the original sequence at a specified step size, generating high-dimensional first-path phase space trajectory matrices and second-path phase space trajectory matrices respectively. This process expands the system's nonlinear dynamic attractor in the one-dimensional time-series signal into a multidimensional space.

[0031] S4. Calculate the maximum Lyapunov exponent of the two phase space trajectories within a preset time window as the instability representation of the thermal system and the fusion system. Calculate the Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients. Amplify the instability representation of the thermal system after dimensionless transformation by the preset time scale parameter to obtain the fusion deviation index. Establish a confidence boundary for a non-equiaxial Gaussian mixture model in a two-dimensional monitoring plane composed of the instability representation of the fusion system and the fusion deviation index. Anomalies are determined when the coordinate points fall outside the boundary.

[0032] For the reconstructed two phase space trajectories, the maximum Lyapunov exponents of the first and second phase space trajectories within a preset time window are calculated respectively. The maximum Lyapunov exponent of the first trajectory is used as the instability representation of the thermal system, and the maximum Lyapunov exponent of the second trajectory is used as the instability representation of the fused system. The L2 norms of the projection coefficients of the two trajectories at each moment within the preset time window are calculated to form two norm sequences of length L. The dot product of these two norm sequences is calculated, and the self-product of each norm sequence is also calculated. The dot product is divided by the square root of the product of the two self-products after adding a preset small normal number to obtain the normalized cosine value of the angle. The absolute value of this normalized cosine value is truncated with an upper limit of 1 and a lower limit of the preset small normal number, and then the negative natural logarithm is taken to obtain the Cauchy-Schwarz divergence. That is, the absolute cosine value under the negative logarithm is used as a unified definition, and the form of dividing the square of the inner product by the product of its own inner products is no longer used to avoid the double coefficient difference of the divergence index. The Cauchy-Schwarz divergence is multiplied by the dimensionless thermal system instability factor transformed by the natural exponential function, and the nonlinear amplification process is completed to obtain the fusion deviation index.

[0033] A two-dimensional monitoring plane is constructed using all instability indicators of the fusion system calculated under historical normal data as the x-axis and the corresponding fusion deviation index as the y-axis. The two-dimensional coordinate data is fitted to establish a non-isoaxial Gaussian mixture model, allowing each cluster to exhibit a non-isoaxial elliptical shape. The expectation-maximization algorithm is used to solve for the model parameters. The sample distribution and confidence ellipse boundary of the fusion system instability characteristics and fusion deviation are shown below. Figure 3As shown. The probability density values ​​of historical normal samples under the Gaussian mixture model are calculated, and a low-density contour line that can encompass a preset proportion of normal samples is selected as the confidence boundary. For example, when the confidence level is set to 99%, the 1st percentile of the probability density value of historical normal samples is selected as the density threshold, or equivalently, the 99th percentile of the negative log-likelihood value is selected as the anomaly threshold. If the probability density of the real-time two-dimensional coordinate point of the online data is lower than this low-density threshold, or the negative log-likelihood value is higher than the corresponding anomaly threshold, the real-time coordinate point is determined to fall outside the boundary, and an anomaly status signal is output. This ensures that the confidence boundary encompasses the vast majority of normal samples, rather than misjudging high-probability normal samples as anomalies. Furthermore, when the online coordinate point mainly shows an increase in the thermal system instability indicator while the fusion deviation index does not increase significantly, a thermal fluctuation anomaly warning is output; when the online coordinate point mainly shows an increase in the fusion deviation index, or when the fusion deviation index and the fusion system instability indicator increase synchronously, a material fusion mismatch anomaly warning is output. The anomaly type prompts are used to assist maintenance personnel in locating the source of anomalies and do not affect the anomaly determination results based on non-equiaxial confidence boundaries.

[0034] In an optional embodiment, the calculation of the Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients, and the amplification of the fusion deviation index using a dimensionless representation of thermal system instability after a preset time scale parameter, includes: Calculate the inner product of the norm sequence corresponding to the first projection coefficient and the norm sequence corresponding to the second projection coefficient, and calculate the inner product of the two norm sequences themselves respectively; The normalized cosine of the included angle, protected by numerical truncation, is calculated based on the inner product and its own inner product. The negative logarithm of the absolute value of the normalized cosine of the included angle is then taken to obtain the Cauchy-Schwarz divergence. The dimensionless thermal instability factor is obtained by multiplying the thermal system instability representation by a preset time scale parameter, and an exponential amplification factor is constructed using the dimensionless thermal instability factor. Multiplying the exponential amplification factor by the Cauchy-Schwarz divergence yields the fusion deviation index.

[0035] Within a preset time window of the current sliding process, for example, with a window data length of L=100, the norm sequence of the first projection coefficients is a one-dimensional vector A, and the norm sequence of the second projection coefficients is a one-dimensional vector B. The dot product of these two evolutionary sequences is calculated by multiplying and summing corresponding elements of vectors A and B, and this dot product is used as a cross-correlation term. The dot product of the first and second norm sequences themselves is calculated by multiplying and summing A and B, respectively. The numerator of this cross-correlation term is divided by the denominator obtained by taking the square root of the product of the two self-products after adding a preset small positive constant, to obtain the normalized cosine value of the angle in the vector space. This cosine value is in the range [-1, 1], representing the degree of similarity and synchronization of the energy change trajectories of the two subsystems in linear space. Extract the absolute value of the cosine and limit it to a range with an upper limit of 1 and a lower limit of a preset small positive constant, such as 10 to the power of -8, to avoid numerical singularities when taking the logarithm. Then, take the negative of the natural logarithm of the truncated absolute value to obtain the Cauchy-Schwarz divergence. If the divergence tends to 0, it indicates that the energy flow changes of the thermal and material subsystems are synchronized.

[0036] To highlight the disruptive characteristics of the two-subsystem fusion mechanism caused by thermal front-end disturbances under abnormal operating conditions, the maximum Lyapunov exponent is extracted from the phase space trajectory calculation results of the previous step as a representation of thermal system instability. This thermal system instability representation is multiplied by a preset time scale parameter to obtain a dimensionless thermal instability factor, which is then used as the base of the natural constant e to construct an exponential amplification factor. The exponential amplification factor can be expressed as exp(β·λ1·T0), where λ1 is the thermal system instability representation, T0 is the preset time scale parameter, and β is a dimensionless proportional adjustment coefficient. The preferred range of β is 1.0 to 2.5, for example, set to 1.5. Specifically, when λ1 is calculated per unit time, T0 takes the preset time scale corresponding to the time unit; when λ1 is calculated per sampling step, T0 takes the corresponding preset number of sampling steps, making λ1·T0 a dimensionless quantity. This mechanism ensures that when the thermal system enters the divergence or chaotic edge, i.e., when the maximum Lyapunov exponent is greater than 0, a nonlinear amplification gain can be generated. The generated exponential amplification factor is multiplied by the obtained Cauchy-Schwarz divergence to calculate the fusion deviation index. This index is used to amplify the characteristic decoupling divergence signal of the two subsystems under thermal instability. For material-end blockage, sudden load changes, or abnormal moisture content, the impact is first manifested in changes to the norm sequence of the second projection coefficient, the second phase space trajectory, and the instability representation of the fusion system. This is further reflected in the increased Cauchy-Schwarz divergence through a decrease in the synchronization relationship between the two norm sequences. When the thermal system exhibits instability trends simultaneously, the exponential amplification factor increases the sensitivity to this fusion mismatch. Therefore, this embodiment does not require material-end faults to necessarily propagate in reverse to thermal-end faults; instead, identification is achieved through the synchronous disruption of the two projection sequences and the joint offset of the two-dimensional monitoring plane.

[0037] In an optional embodiment, establishing a confidence boundary for the non-isoaxial Gaussian mixture model in a two-dimensional monitoring plane composed of the instability representation of the fused system and the fusion deviation index includes: Based on historical normal operation data, the instability representation and fusion deviation index of the fusion system are calculated to form a two-dimensional feature vector set; The Gaussian mixture model parameters of the two-dimensional feature vector set are estimated using the expectation-maximization algorithm, resulting in the mean vector of multiple Gaussian distribution components, the regularized off-diagonal covariance matrix, and the mixture weights. The probability density contour lines are calculated based on the Gaussian mixture model, and the probability density contour lines corresponding to the preset confidence threshold are selected as the non-isoaxial confidence boundaries of the normal operation area.

[0038] Data within the steady-state normal operating cycle is extracted, and a two-dimensional feature vector containing Q sample points is calculated through the aforementioned steps to obtain the standard model training set. Each feature vector contains two dimensional components: the first dimension is the maximum Lyapunov exponent obtained from the second phase space calculation, and the second dimension is the fusion deviation index. This two-dimensional feature vector set is fed into the expectation-maximization algorithm module to perform parameter fitting of the Gaussian mixture model. Before fitting, the optimal number of Gaussian components is preferably determined using the Bayesian information criterion, with the search range limited to 2 to 5, for example, 3 components. The EM algorithm is used to perform alternating iterations of E-step expectation estimation and M-step parameter maximization updates, with a maximum iteration limit of 500 times and a convergence tolerance of 10 to the power of -5, to estimate the mean vectors of the three independent Gaussian distribution components, the regularized off-diagonal covariance matrix, and the mixture weights. When the covariance matrix approaches singularity, a preset small regularization term is added to its diagonal to ensure the stability of the probability density calculation.

[0039] After estimating the model parameters, these components are fused to construct a global probability density function, which is then evaluated based on discrete grid points in the feature plane coordinate system, generating probability density contour lines representing the spatial data density. A confidence control threshold is set, preferably 95% or 99%. Specific density contours corresponding to this threshold condition are solved and extracted using numerical spatial integration, requiring that the region defined by these contours encompasses 99% of the normal state points in the training set. Specifically, historical normal sample probability density values ​​can be sorted from low to high, and the 1st percentile under a 99% confidence requirement is taken as the low density contour; alternatively, negative log-likelihood values ​​can be sorted from low to high, and the 99th percentile is taken as the anomaly detection threshold. These two methods are opposite but equivalent; the former determines anomalies based on probability density below the threshold, while the latter determines anomalies based on negative log-likelihood above the threshold. This boundary defined by the density contour is extracted, presenting a non-equiaxial confidence envelope on a two-dimensional plane. When deploying online monitoring, simply substitute the two-dimensional coordinates of the new moment into this distribution density function model. If the calculated probability density value is lower than the low density equal value threshold, or the negative log-likelihood value is higher than the corresponding anomaly judgment threshold, then the current coordinate point is determined to be outside the confidence boundary and an anomaly exceeding the limit alarm is triggered. This method realizes the monitoring and confidence protection of skewed and multimodal time-varying data.

[0040] In the second embodiment, the present invention also proposes a sludge drying system operation status monitoring system, comprising the following modules: The partitioning module is used to obtain the original dataset of historical normal operation and standardize it to obtain the modeling dataset. The process variables in the modeling dataset are divided into a subset of strong time-series thermal variables and a subset of weak time-series material variables. The module is used to construct a lag cross-correlation integral kernel function, which weights the difference components of different dimensions when calculating the distance between sample points. The weights are determined by the absolute value integral of the lag cross-correlation function between the corresponding variable and the key performance index within a preset lag time interval. The lag cross-correlation integral kernel function is used to perform a first-way kernel principal component analysis on the subset of strongly time-series thermal variables to obtain the first-way projection coefficients. The calculation module is used to construct a second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficient, use the second kernel function to perform second kernel principal component analysis on the weak temporal material variable subset to obtain the second projection coefficient, calculate the two projection coefficients on the online real-time data, and construct the two phase space trajectories within a preset time window; The determination module is used to calculate the maximum Lyapunov exponent of the two phase space trajectories within a preset time window as the instability representation of the thermal system and the fusion system. It calculates the Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients and amplifies the instability representation of the thermal system after it has been dimensionlessly transformed by the preset time scale parameter to obtain the fusion deviation index. A non-isoaxial Gaussian mixture model confidence boundary is established in the two-dimensional monitoring plane composed of the instability representation of the fusion system and the fusion deviation index. Anomalies are determined when the coordinate points fall outside the boundary.

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

Claims

1. A method for monitoring the operating status of a sludge drying system, characterized in that, Includes the following steps: The original dataset of historical normal operation is obtained and standardized to obtain the modeling dataset. The process variables in the modeling dataset are divided into a subset of strong time-series thermal variables and a subset of weak time-series material variables. A lag cross-correlation integral kernel function is constructed, and the difference components of different dimensions are weighted when calculating the distance between sample points. The weights are determined by the absolute value integral of the lag cross-correlation function between the corresponding variable and the key performance index within a preset lag time interval. The first-way kernel principal component analysis is performed on the subset of strongly time-series thermal variables using the lag cross-correlation integral kernel function to obtain the first-way projection coefficients. A second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficient is constructed. The second kernel function is used to perform second kernel principal component analysis on the weakly temporal material variable subset to obtain the second projection coefficient. The two projection coefficients are calculated on the online real-time data and the two phase space trajectories are constructed within a preset time window. The maximum Lyapunov exponent of the two phase space trajectories within a preset time window is calculated as the instability representation of the thermal system and the fusion system. The Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients is calculated. The fusion deviation index is obtained by amplifying the thermal system instability representation after dimensionless transformation by the preset time scale parameter. A non-isoaxial Gaussian mixture model confidence boundary is established in the two-dimensional monitoring plane composed of the fusion system instability representation and the fusion deviation index. Anomalies are determined when the coordinate points fall outside the boundary.

2. The method according to claim 1, characterized in that, The process involves obtaining and standardizing the original historical normal operation dataset to obtain the modeling dataset. The process variables in the modeling dataset are then divided into a strongly time-series thermal variable subset and a weakly time-series material variable subset, including: Calculate the time series autocorrelation function of each process variable in the modeling dataset; Determine the delay time corresponding to the first crossing of zero or reaching a preset small threshold for the autocorrelation function of each process variable, and use the delay time as the decorrelation time scale for the corresponding process variable; Set time-scale classification thresholds; Process variables with a decorrelated time scale greater than the time scale classification threshold are classified into the strongly time-series thermal variables subset, and process variables with a decorrelated time scale less than or equal to the time scale classification threshold are classified into the weakly time-series material variables subset.

3. The method according to claim 2, characterized in that, The construction of the lagged cross-correlation integral kernel function involves weighting the difference components of different dimensions when calculating the distance between sample points. The weights are determined by the absolute value integral of the lagged cross-correlation function between the corresponding variable and the key performance index over a preset lag time interval. The lagged cross-correlation integral kernel function is then used to perform a first-way kernel principal component analysis on the subset of strongly time-series thermal variables to obtain the first-way projection coefficients, including: For each variable in the strongly time-series thermal variable subset, calculate the sequence of cross-correlation functions between it and the moisture content index of the sludge drying system within the preset lag time interval. Integrate the absolute values ​​of each cross-correlation function sequence within the preset lag time interval to obtain the initial weights of each variable; The initial weights of all strongly time-series thermal variables are normalized to obtain the final weights of each variable. When calculating the Gaussian kernel function, the squared difference between two sample points in each variable dimension is multiplied by the corresponding final weight, and the sum of the squared weighted differences in all dimensions is substituted into the exponential function to obtain the value of the lagged cross-correlation integral kernel function. The first-way kernel principal component analysis is performed on the subset of strongly time-series thermal variables using the hysteresis cross-correlation integral kernel function to obtain the first-way projection coefficients.

4. The method according to claim 1, characterized in that, The construction of a second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficients, and the use of the second kernel function to perform second-way kernel principal component analysis on the weakly temporal material variable subset to obtain the second projection coefficients, includes: Extract the first k kernel principal components from the first projection coefficients; Calculate the sum of squares or variance of the differences between the first k principal components at adjacent time points, and use it as the characteristic value of thermal state fluctuation; Substituting the thermal state fluctuation characteristic value into a monotonically decreasing exponential decay function, the unified adjustment factor corresponding to the preset modulation time window is obtained. Multiply the preset base kernel width parameter by the unified adjustment factor to obtain the unified kernel width parameter of the second kernel function within the preset modulation time window; The unified kernel width parameter is used to construct a second kernel function, so that any sample point pair within the same kernel matrix uses the same kernel width parameter, and the kernel width of the second kernel function is reduced when the thermal state time-series fluctuations are enhanced.

5. The method according to claim 1 or 4, characterized in that, The step of calculating the projection coefficients of two channels from online real-time data and constructing the phase space trajectories of two channels within a preset time window includes: Calculate the norms of the first and second projection coefficient vectors within the preset time window to form two one-dimensional time window projection coefficient sequences. The optimal delay time is calculated using the mutual information method for the two one-dimensional time window projection coefficient sequences, and the minimum embedding dimension is calculated using the spurious nearest neighbor method. Based on the corresponding optimal delay time and minimum embedding dimension, the two one-dimensional time window projection coefficient sequences are reconstructed in phase space to obtain the first phase space trajectory matrix and the second phase space trajectory matrix.

6. The method according to claim 1, characterized in that, The calculation of the Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients, and the amplification of the thermal system instability representation after dimensionless transformation by a preset time scale parameter to obtain the fusion deviation index, includes: Calculate the inner product of the norm sequence corresponding to the first projection coefficient and the norm sequence corresponding to the second projection coefficient, and calculate the inner product of the two norm sequences themselves respectively; The normalized cosine of the included angle, protected by numerical truncation, is calculated based on the inner product and its own inner product. The negative logarithm of the absolute value of the normalized cosine of the included angle is then taken to obtain the Cauchy-Schwarz divergence. The dimensionless thermal instability factor is obtained by multiplying the thermal system instability representation by a preset time scale parameter, and an exponential amplification factor is constructed using the dimensionless thermal instability factor. Multiplying the exponential amplification factor by the Cauchy-Schwarz divergence yields the fusion deviation index.

7. The method according to claim 1, characterized in that, The establishment of a confidence boundary for a non-isoaxial Gaussian mixture model in a two-dimensional monitoring plane composed of the instability representation of the fused system and the fusion deviation index includes: Based on historical normal operation data, the instability representation and fusion deviation index of the fusion system are calculated to form a two-dimensional feature vector set; The Gaussian mixture model parameters of the two-dimensional feature vector set are estimated using the expectation-maximization algorithm to obtain the mean vector of multiple Gaussian distribution components, the regularized off-diagonal covariance matrix, and the mixture weights. The probability density contour lines are calculated based on the Gaussian mixture model, and the probability density contour lines corresponding to the preset confidence threshold are selected as the non-isoaxial confidence boundaries of the normal operation area.

8. A sludge drying system operation status monitoring system, characterized in that, Includes the following modules: The partitioning module is used to obtain the original dataset of historical normal operation and standardize it to obtain the modeling dataset. The process variables in the modeling dataset are divided into a subset of strong time-series thermal variables and a subset of weak time-series material variables. The module is used to construct a lag cross-correlation integral kernel function, which weights the difference components of different dimensions when calculating the distance between sample points. The weights are determined by the absolute value integral of the lag cross-correlation function between the corresponding variable and the key performance index within a preset lag time interval. The lag cross-correlation integral kernel function is used to perform a first-way kernel principal component analysis on the subset of strongly time-series thermal variables to obtain the first-way projection coefficients. The calculation module is used to construct a second kernel function modulated by the temporal fluctuation characteristics of the first projection coefficient, use the second kernel function to perform second kernel principal component analysis on the weak temporal material variable subset to obtain the second projection coefficient, calculate the two projection coefficients on the online real-time data, and construct the two phase space trajectories within a preset time window; The determination module is used to calculate the maximum Lyapunov exponent of the two phase space trajectories within a preset time window as the instability representation of the thermal system and the fusion system. It calculates the Cauchy-Schwarz divergence between the norm sequences corresponding to the two projection coefficients and amplifies the instability representation of the thermal system after it has been dimensionlessly transformed by the preset time scale parameter to obtain the fusion deviation index. A non-isoaxial Gaussian mixture model confidence boundary is established in the two-dimensional monitoring plane composed of the instability representation of the fusion system and the fusion deviation index. Anomalies are determined when the coordinate points fall outside the boundary.

9. The system according to claim 8, characterized in that, The process involves obtaining and standardizing the original historical normal operation dataset to obtain the modeling dataset. The process variables in the modeling dataset are then divided into a strongly time-series thermal variable subset and a weakly time-series material variable subset, including: Calculate the time series autocorrelation function of each process variable in the modeling dataset; Determine the delay time corresponding to the first crossing of zero or reaching a preset small threshold for the autocorrelation function of each process variable, and use the delay time as the decorrelation time scale for the corresponding process variable; Set time-scale classification thresholds; Process variables with a decorrelated time scale greater than the time scale classification threshold are classified into the strongly time-series thermal variables subset, and process variables with a decorrelated time scale less than or equal to the time scale classification threshold are classified into the weakly time-series material variables subset.

10. The system according to claim 8, characterized in that, The construction of the lagged cross-correlation integral kernel function involves weighting the difference components of different dimensions when calculating the distance between sample points. The weights are determined by the absolute value integral of the lagged cross-correlation function between the corresponding variable and the key performance index over a preset lag time interval. The lagged cross-correlation integral kernel function is then used to perform a first-way kernel principal component analysis on the subset of strongly time-series thermal variables to obtain the first-way projection coefficients, including: For each variable in the strongly time-series thermal variable subset, calculate the sequence of cross-correlation functions between it and the moisture content index of the sludge drying system within the preset lag time interval. Integrate the absolute values ​​of each cross-correlation function sequence within the preset lag time interval to obtain the initial weights of each variable; The initial weights of all strongly time-series thermal variables are normalized to obtain the final weights of each variable. When calculating the Gaussian kernel function, the squared difference between two sample points in each variable dimension is multiplied by the corresponding final weight, and the sum of the squared weighted differences in all dimensions is substituted into the exponential function to obtain the value of the lagged cross-correlation integral kernel function. The first-way kernel principal component analysis is performed on the subset of strongly time-series thermal variables using the hysteresis cross-correlation integral kernel function to obtain the first-way projection coefficients.