Beidou intelligent early warning analysis method and system applied to basin-level reservoir dam group

Abstract

The application provides a Beidou intelligent early warning analysis method and system applied to a basin-level reservoir dam group, time sequence records of dam body deformation and corresponding time sequence records of hydrological environment collected by a Beidou satellite positioning terminal of a basin-level reservoir dam group monitoring system are acquired, self-encoding processing is performed on associated time sequence data sets in a dam body deformation response mode, a deformation response characteristic element set representing a dam body structure state is generated, the deformation response characteristic element set is input into a pre-constructed topological relationship graph of the basin-level reservoir dam group, abnormal state propagation calculation is performed through state dependence relationships between nodes in the graph, a cascading risk transmission path set is obtained, risk evolution trend deduction is performed on reservoir dam nodes in each path based on this, a risk evolution index matrix is generated, and a precursor risk identification sequence for key reservoir dam nodes is generated. Through the application, the structuralization and response operability of early warning information can be improved, and the scientificity and practicality of basin-level reservoir dam group early warning analysis can be comprehensively improved.

Description

Technical Field

[0001] This invention relates to the field of data processing technology, and in particular to a BeiDou intelligent early warning analysis method and system applied to a basin-level reservoir and dam group. Background Technology

[0002] With the expansion of water conservancy projects, safety monitoring and early warning analysis of basin-level reservoir-dam groups have become crucial means to ensure the safe operation of water conservancy facilities. Currently, dam deformation monitoring and hydrological environment data collection based on BeiDou satellite positioning technology have been applied in practice, providing a data foundation for reservoir-dam safety assessment. However, existing early warning analysis methods still have room for improvement in practical applications. For example, the coupling relationship between dam data and hydrological data is not sufficiently captured, making it difficult to accurately extract key features reflecting the dam's structural state. At the basin-level reservoir-dam group level, the systematic analysis of the state dependencies and risk transmission paths between reservoir-dam nodes is insufficient, affecting the overall effectiveness of basin-level reservoir-dam group early warning analysis. Summary of the Invention

[0003] In view of this, the present invention provides a BeiDou intelligent early warning analysis method and system for basin-level reservoir and dam groups.

[0004] On one hand, embodiments of the present invention provide a BeiDou intelligent early warning analysis method applied to a basin-level reservoir-dam group. The method includes: acquiring time-series records of dam deformation and corresponding hydrological environment collected by a basin-level reservoir-dam group monitoring system via a BeiDou satellite positioning terminal, wherein the dam deformation time-series records and the hydrological environment time-series records are associated and bound through a unified timestamp to obtain an associated time-series dataset; performing self-encoding processing on the associated time-series dataset to generate a set of deformation response feature primitives characterizing the dam structure state by capturing the coupling relationship between deformation and the hydrological environment; inputting the set of deformation response feature primitives into a pre-constructed topological relationship map of the basin-level reservoir-dam group, performing anomaly state propagation calculation through the state dependency relationship between map nodes to obtain a set of chain risk transmission paths containing risk propagation intensity and path probability; performing risk evolution trend deduction on reservoir-dam nodes in each path based on the set of chain risk transmission paths to generate a risk evolution index matrix containing node risk accumulation rate and critical threshold spacing; and generating a precursor risk identification sequence for key reservoir-dam nodes based on the risk evolution index matrix.

[0005] On the other hand, embodiments of the present invention provide an early warning system, including a memory and a processor. The memory stores a computer program that can run on the processor, and the processor executes the program to implement the steps in the above-described method.

[0006] This invention establishes a coupled data foundation of structural response and environment-driven processes by acquiring and binding time-series records of dam deformation and hydrological environment, thereby improving the correlation and analytical effectiveness of the original data. Through self-encoding of the dam deformation response patterns in the associated time-series dataset, the coupling relationship between deformation and the hydrological environment is captured, and a set of deformation response feature primitives characterizing the dam's structural state is generated. This represents a fundamental improvement from high-dimensional time-series data to structured state features, enhancing the expressive power of complex dam response patterns. By inputting the feature primitive set into the topological relationship map of the watershed reservoir-dam group and performing anomaly state propagation calculations, a network-level risk analysis framework for inter-node state dependencies is constructed, overcoming the limitations of static assessment of isolated nodes and revealing... The study of risk transmission patterns along hydraulic and spatial paths enhances the systematicness and comprehensiveness of risk assessment. By performing risk evolution trend deduction on reservoir and dam nodes along the path and generating a risk evolution index matrix containing node risk accumulation rates and critical threshold intervals, quantitative analysis from the current state to dynamic trends is achieved, overcoming the lack of foresight in static threshold comparisons and enhancing the dynamism and predictability of risk assessment. By generating a sequence of precursor risk identifiers containing risk warning levels and spatiotemporal triggering conditions based on the index matrix, multidimensional risk characteristics are integrated and key nodes are focused on, avoiding the problem of fragmented risk information making it difficult to guide decision-making. This improves the structuring of warning information and the operability of response, thereby comprehensively enhancing the scientificity and practicality of basin-level reservoir and dam group early warning analysis. Attached Figure Description

[0007] Figure 1 This is a schematic diagram illustrating the implementation process of a BeiDou intelligent early warning analysis method applied to a watershed-level reservoir and dam group, as provided in an embodiment of the present invention.

[0008] Figure 2 This is a schematic diagram of the hardware entity of an early warning system provided in an embodiment of the present invention. Detailed Implementation

[0009] This invention provides a BeiDou intelligent early warning analysis method applied to a basin-level reservoir and dam group. This method can be executed by the processor of an early warning system. The early warning system can refer to devices with data processing capabilities, such as servers, laptops, tablets, and desktop computers.

[0010] Figure 1 This is a schematic diagram illustrating the implementation process of a BeiDou intelligent early warning analysis method applied to a basin-level reservoir and dam group, as provided in an embodiment of the present invention. Figure 1 As shown, the method includes:

[0011] Step S100: Obtain the dam deformation time series records and corresponding hydrological environment time series records collected by the basin-level reservoir and dam group monitoring system through the Beidou satellite positioning terminal. The dam deformation time series records and hydrological environment time series records are associated and bound through a unified timestamp to obtain an associated time series dataset.

[0012] Dam deformation time-series records are data sets that record the deformation of the dam body in chronological order over a period of time. They reflect changes in the dam's structural state at different times, such as settlement, tilting, and displacement. Hydrological environment time-series records are data that record hydrological environmental factors within a watershed in chronological order over the same time period, such as changes in water level and flow velocity. For example, BeiDou satellite positioning terminals are installed at key locations in each reservoir and dam, continuously collecting dam deformation data. Simultaneously, hydrological monitoring equipment, such as water level gauges and current meters, is deployed at appropriate locations within the watershed to collect hydrological environmental data such as water level and flow velocity. This collected data is transmitted to the server of the watershed-level reservoir and dam group monitoring system, where a precise timestamp is added to each data point. Then, through a database management system, dam deformation data and hydrological environmental data with the same timestamp are associated and bound to form a correlated time-series dataset.

[0013] Step S200: Perform autoencoding processing on the associated time series dataset to generate a set of deformation response feature primitives that characterize the dam structure state by capturing the coupling relationship between deformation and hydrological environment.

[0014] The coupling relationship between deformation and hydrological environment reflects the mutual influence and interaction between the two. For example, a rise in water level may cause the dam body to bear greater pressure, thereby causing an increase in dam settlement.

[0015] In one implementation, step S200 may include the following steps S210-S260:

[0016] Step S210: Divide the associated time series dataset into time windows. Divide the continuous time series data into blocks with time continuity according to a preset time interval. Each block contains dam deformation record segments and hydrological environment record segments within the corresponding window. The number of blocks is determined according to the total duration of the dataset and the window length.

[0017] Time window partitioning is the operation of dividing a continuous, correlated time-series dataset into segments according to a set time length. The preset time interval is a fixed, pre-defined time length used to determine the size of each time window.

[0018] For example, the preset time interval is first determined based on the analysis requirements and data characteristics. For instance, if the goal is to analyze the dam's response to short-term hydrological changes, a shorter preset time interval can be set; if the focus is on long-term trends, a longer preset time interval can be set. Then, the number of blocks is calculated based on the total dataset duration and the preset time interval. Specifically, the total dataset duration is divided by the preset time interval. For example, if the total dataset duration is a long period and the preset time interval is a short period, the number of blocks can be obtained through division. Afterward, the associated time-series dataset is segmented sequentially according to the preset time interval, dividing the continuous data into multiple temporally continuous blocks. The dam deformation record segment and the hydrological environment record segment in each block correspond temporally.

[0019] Step S220: Extract deformation features and hydrological features from the segmented data. Deformation features include the rate of change of settlement and the amplitude of horizontal displacement fluctuation. Hydrological features include the amplitude of water level rise and fall and the gradient of flow velocity change. The feature extraction process retains the timestamp correlation of the original data.

[0020] The rate of change in settlement refers to the speed at which the dam body's settlement changes over a certain period of time, reflecting the dam's stability in the vertical direction. The amplitude of horizontal displacement fluctuation refers to the difference between the maximum and minimum values ​​of the dam body's horizontal displacement, reflecting the range of horizontal displacement changes and helping to determine if the dam body exhibits horizontal instability. Hydrological characteristics are indicators extracted from hydrological environmental data that reflect the characteristics of the hydrological environment. The amplitude of water level rise and fall refers to the height difference between the rise and fall of water level over a certain period of time; positive values ​​indicate rising water levels, and negative values ​​indicate falling water levels, reflecting changes in water volume within the basin. The velocity change gradient refers to the degree of change in flow velocity within a certain time window, usually expressed as the standard deviation of flow velocity, reflecting the stability of water flow velocity.

[0021] During feature extraction, for the extraction of the settlement change rate, it is first necessary to obtain the initial and final values ​​of the dam settlement in the data blocks, then calculate the difference between them, and finally divide the difference by the corresponding time length of the data block to obtain the settlement change rate. For the extraction of the horizontal displacement fluctuation amplitude, it is necessary to find the maximum and minimum values ​​of the dam's horizontal displacement in the data blocks, and then calculate their difference. For the extraction of the water level rise and fall amplitude, the initial and final values ​​of the water level in the data blocks are obtained, and the difference between them is calculated. For the extraction of the velocity change gradient, the standard deviation of the velocity in the data blocks is first calculated; the larger the standard deviation, the greater the velocity change gradient. Throughout the entire feature extraction process, special attention should be paid to preserving the timestamp correlation of the original data.

[0022] Step S230: Model the mutual correlation between the extracted deformation features and hydrological features, calculate the feature mutual information value under different time lags, select the lag with the largest mutual information value as the optimal response delay, and generate a spatiotemporal correlation matrix containing response intensity and delay duration.

[0023] Interrelation modeling is the process of establishing the correlation between deformation characteristics and hydrological characteristics. Its purpose is to identify the intrinsic links between these characteristics through data analysis, and to understand how changes in the hydrological environment affect dam deformation. A time lag refers to the delay in time between one variable and another when analyzing the correlation between two variables. The optimal response delay is the time lag that maximizes the mutual information value between the deformation and hydrological characteristics, representing the best response time of the dam deformation to changes in the hydrological environment.

[0024] In one implementation, step S230 may include the following steps S231-S236:

[0025] Step S231: Perform stationarity tests on the extracted deformation and hydrological features, verify the stationarity of the sequence using the unit root test method, and perform differencing on non-stationary sequences until the stationarity requirements are met.

[0026] For example, the extracted deformation and hydrological features are first treated as time series data and subjected to unit root tests. For each feature sequence, a series of tests are performed according to the rules of the unit root test method. If the test results indicate that the sequence has a unit root, i.e., the sequence is non-stationary, then the sequence needs to be differencing. The differencing process involves calculating the difference between adjacent data points to obtain a new sequence. Then, the new sequence is subjected to a unit root test again to determine whether it is stationary. If the new sequence is still non-stationary, the differencing process continues until the sequence meets the stationarity requirement.

[0027] Step S232: Obtain the preset time lag search range. The lag value ranges from the initial value to the maximum lag value. The maximum lag value is determined based on the sampling frequency of the feature sequence and the hydrological response period. The lag interval is consistent with the sampling period.

[0028] For example, the maximum lag value is first determined based on the sampling frequency of the feature sequence and the hydrological response period. For instance, if the sampling frequency of the feature sequence is high, indicating denser data collection, the maximum lag value can be appropriately increased to cover more possible time lag scenarios. If the hydrological response period is long, the maximum lag value also needs to be set larger to ensure that the complete response process of dam deformation to changes in the hydrological environment can be captured. After determining the initial value and the maximum lag value, the search range for the time lag is obtained. Then, according to the lag interval consistent with the sampling period, different time lag values ​​are sequentially selected within the search range for subsequent mutual information value calculation. For example, if the sampling period is a fixed time length, the lag interval is also set to that time length. Starting from the initial value, one lag interval is added each time to obtain a series of different time lag values.

[0029] Step S233: Calculate the mutual information value between deformation features and hydrological features for each hysteresis. The mutual information value reflects the correlation strength between features. The larger the value, the stronger the correlation. The calculation process takes into account the estimation error of the feature probability density function.

[0030] For example, for each selected time lag, the deformation feature sequence and the hydrological feature sequence are aligned according to that lag. That is, the hydrological feature sequence is shifted backward relative to the deformation feature sequence by the time length corresponding to that lag, and then the mutual information value is calculated for the two shifted sequences. When calculating the mutual information value, it is first necessary to estimate the probability density functions of the deformation and hydrological features. Some commonly used probability density function estimation methods, such as kernel density estimation, can be used. During the estimation process, the characteristics and distribution of the data should be fully considered to minimize estimation errors. Then, based on the estimated probability density functions, the mutual information value is calculated according to the rules for calculating mutual information values. During the calculation process, estimation errors should be appropriately handled; for example, some correction methods can be used to reduce the impact of errors on the mutual information value.

[0031] Step S234: Record the mutual information value corresponding to each hysteresis quantity, and plot the mutual information-hysteresis quantity relationship curve. The hysteresis quantity corresponding to the peak value of the curve is the optimal response delay, and the peak mutual information value is the response intensity. Each feature pair corresponds to a set of optimal delay and intensity.

[0032] The peak value of the curve refers to the maximum point in the mutual information-hysteresis relationship curve. The hysteresis corresponding to this peak value is the optimal response delay, indicating that under this time lag, the correlation between deformation characteristics and hydrological characteristics is the strongest, and the response of dam deformation to changes in the hydrological environment is most significant. The peak mutual information value is the response intensity, reflecting the degree of correlation between deformation characteristics and hydrological characteristics. Each feature pair (i.e., each combination of deformation characteristics and hydrological characteristics) corresponds to a set of optimal delays and intensities. This is because the correlation between different feature pairs may be different, and the optimal response time and response intensity will also vary.

[0033] For example, the recorded hysteresis and corresponding mutual information values ​​are input into a plotting tool to draw a mutual information-hysteresis relationship curve. Then, the peak point of the curve is identified by analysis. After determining the peak point, the hysteresis and mutual information values ​​corresponding to that peak point are recorded, which are respectively used as the optimal response delay and response strength of that feature pair.

[0034] Step S235: Based on the optimal response delay and response intensity of all feature pairs, construct a spatiotemporal correlation matrix. The matrix row index is the deformation feature type, the column index is the hydrological feature type, the matrix element is the response intensity of the corresponding feature pair, and the element position corresponds to the optimal delay duration.

[0035] For example, firstly, based on the optimal response delay and response intensity of each feature pair calculated previously, the size and element values ​​of the matrix are determined. The number of rows in the matrix equals the number of deformation feature types, and the number of columns equals the number of hydrological feature types. For each feature pair, its response intensity is used as the element value at the corresponding position in the matrix, with the row and column indices corresponding to the deformation feature type and the hydrological feature type, respectively. Simultaneously, the specific position of the element in the matrix is ​​determined based on the optimal response delay, reflecting the optimal response time of the feature pair. For instance, for the feature pair of settlement change rate and water level fluctuation amplitude, its response intensity is used as the element value in the matrix corresponding to the row of settlement change rate and the column of water level fluctuation amplitude, and the specific position of the element in the matrix is ​​determined based on the optimal response delay of the feature pair.

[0036] Step S236: Sparsify the spatiotemporal correlation matrix by retaining matrix elements with response strength greater than a preset threshold and setting the remaining elements to zero to generate a sparse spatiotemporal correlation matrix.

[0037] For example, a preset threshold is first determined. The preset threshold can be adjusted according to the needs of the analysis and the characteristics of the data. For example, if more correlation information is desired, the preset threshold can be set lower; if only feature pairs with very strong correlations are considered, the preset threshold can be set higher. Then, each element in the spatiotemporal correlation matrix is ​​iterated, and elements with response strength greater than the preset threshold are retained, while elements with response strength less than the preset threshold are set to zero.

[0038] Step S240: Based on the spatiotemporal correlation matrix, perform multi-scale feature separation on the block data, decompose it into feature components of different vibration modes, extract short-term deformation features from the instantaneous response components, extract long-term deformation features from the trend influence components, and match the feature dimensions with the dimensions of the block data.

[0039] For example, firstly, based on the information in the spatiotemporal correlation matrix, the correlation relationships and response patterns between different features in the segmented data are determined. Then, a suitable multi-scale analysis method, such as wavelet analysis, is used to decompose the segmented data. Wavelet analysis can decompose the segmented data into components of different frequencies and scales, corresponding to different vibration modes. By analyzing the decomposed components, they are divided into instantaneous response components and trend influence components. For the instantaneous response components, short-term deformation characteristics are extracted. For example, by calculating and processing the data in the components, characteristic indicators reflecting the instantaneous deformation of the dam body can be obtained. For the trend influence components, long-term deformation characteristics are extracted. For example, by smoothing the data in the components and performing trend analysis, characteristic indicators reflecting the long-term deformation trend of the dam body can be obtained.

[0040] Step S250: Combine short-term deformation features with long-term hydrological features to form joint features. Input the joint features into an autoencoder processing unit. Compress the joint features into a fixed-dimensional latent space vector through a multi-layer nonlinear mapping network. Reconstruct the original features through a decoding network. When the reconstruction error is less than a preset convergence threshold, save the current latent space vector as a candidate feature primitive.

[0041] For example, the extracted short-term deformation features and long-term hydrological features are first combined to form joint features. Then, the joint features are input into a multi-layer nonlinear mapping network of the autoencoder processing unit. The multi-layer nonlinear mapping network consists of multiple nonlinear layers, each performing a nonlinear transformation on the input data, progressively compressing the data into a fixed-dimensional latent space to obtain a latent space vector. Next, the latent space vector is input into a decoding network, which reconstructs the original joint features from the latent space vector through a series of inverse transformation operations. The difference between the reconstructed features and the original features is calculated to obtain the reconstruction error. The reconstruction error is compared with a preset convergence threshold. If the reconstruction error is less than the preset convergence threshold, the autoencoder processing unit is considered to have converged, and the current latent space vector is saved as a candidate feature primitive. If the reconstruction error is greater than the preset convergence threshold, the parameters of the autoencoder processing unit, such as the network weights and biases, need to be adjusted, and the compression and reconstruction operations are performed again until the reconstruction error is less than the preset convergence threshold.

[0042] Step S260: Perform density clustering screening on all candidate feature primitives, divide feature clusters using a clustering algorithm based on spatial distribution density, calculate the spatial distance between the center vectors of each cluster, merge adjacent clusters with a distance less than the clustering threshold, extract the center vector of each cluster as the final deformation response feature primitive, and combine all center vectors to obtain the deformation response feature primitive set.

[0043] A feature cluster is a set of candidate feature primitives that are close to each other in space and have similar characteristics. The cluster center vector is the vector corresponding to the center position of each feature cluster, representing the main feature of that cluster. The clustering threshold is a pre-set standard value used to determine whether to merge adjacent feature clusters. If the spatial distance between the center vectors of two clusters is less than the clustering threshold, the two feature clusters are considered adjacent and similar, and they need to be merged into a larger feature cluster. The final deformation response feature primitives are vectors that represent the deformation response characteristics of the dam body after screening and merging.

[0044] For example, firstly, all candidate feature primitives are taken as input and clustered using a spatial distribution density-based clustering algorithm. This algorithm divides them into different feature clusters based on the spatial distribution density of the feature primitives. Then, the center vector of each feature cluster is calculated, resulting in cluster center vectors. Next, the spatial distance between the cluster center vectors is calculated, and adjacent clusters with a distance less than a clustering threshold are merged. The merging process involves combining candidate feature primitives from adjacent clusters into a new feature cluster and recalculating the center vector of the new feature cluster. Finally, the center vectors of each cluster are extracted as the final deformation response feature primitives, and all center vectors are combined to obtain the deformation response feature primitive set.

[0045] In another implementation, step S200 may include the following steps S201-S206:

[0046] Step S201: Divide the associated time series dataset into segments at equal intervals along the time axis to obtain multiple continuous data segments. Each data segment contains complete dam deformation records and hydrological environment records for the corresponding time period. There is no time overlap between data segments and the data segments cover the entire dataset.

[0047] For example, the time intervals are first determined based on the analytical needs and the characteristics of the data. For instance, the time interval can be set to a suitable value based on the periodic characteristics of dam deformation and hydrological environment changes. Then, the associated time-series dataset is segmented sequentially according to the determined time intervals, dividing the continuous data into multiple continuous data segments. The dam deformation records and hydrological environment records in each data segment are complete and correspond in time.

[0048] Step S202: Extract multi-dimensional feature sequences from each data segment. Deformation dimension features include vertical settlement fluctuation sequences and horizontal displacement change sequences. Hydrological dimension features include water level fluctuation sequences and flow velocity gradient sequences. The sampling frequency of each sequence is consistent with the original record.

[0049] In one implementation, step S202 may include the following steps S2021-S2026:

[0050] Step S2021: Extract vertical settlement data from the dam deformation record, arrange them in chronological order to form the original settlement sequence, calculate the settlement difference between adjacent sampling points to obtain the vertical settlement fluctuation sequence, which reflects the change in settlement rate per unit time.

[0051] For example, vertical settlement data is first extracted from the dam deformation record. This data can be collected by vertical settlement monitoring equipment installed on the dam. Then, the extracted vertical settlement data is arranged in chronological order to form an original settlement sequence. Next, the original settlement sequence is processed to calculate the settlement difference between adjacent sampling points. For example, for the first and second sampling points in the original settlement sequence, their settlement difference is calculated, and this difference is used as the first data point in the vertical settlement fluctuation sequence. Similarly, the settlement difference between all adjacent sampling points in the original settlement sequence is calculated to obtain the vertical settlement fluctuation sequence. By analyzing the vertical settlement fluctuation sequence, the change in the settlement rate of the dam per unit time can be understood, and the vertical stability of the dam can be determined.

[0052] Step S2022: Extract the horizontal displacement data from the dam deformation record, decompose it into lateral and longitudinal displacement components, calculate the fluctuation amplitude of each component (i.e., the difference between the maximum and minimum values), arrange them in chronological order to form a horizontal displacement change sequence, reflecting the changing trend of displacement direction and range.

[0053] For example, firstly, horizontal displacement data is extracted from the dam deformation record. This data can be collected by horizontal displacement monitoring equipment installed on the dam. Then, the horizontal displacement data is decomposed into lateral and longitudinal displacement components. A coordinate transformation method can be used to decompose the horizontal displacement data into components in the lateral and longitudinal directions within a pre-defined coordinate system. Next, for each displacement component, its maximum and minimum values ​​over a certain period of time are calculated, and their difference is calculated to obtain the fluctuation amplitude of that component. Finally, the fluctuation amplitudes of each component are arranged in chronological order to form a horizontal displacement change sequence. By analyzing the horizontal displacement change sequence, the displacement trend of the dam in the horizontal direction can be understood, and the stability of the dam in the horizontal direction can be determined.

[0054] Step S2023: Extract water level monitoring data from hydrological environmental records, calculate the water level difference between adjacent time points, positive values ​​indicate rising water levels and negative values ​​indicate falling water levels, arrange them in chronological order to form a water level fluctuation sequence, and keep the sequence length consistent with the deformation record segment.

[0055] For example, water level monitoring data is first extracted from hydrological environmental records. This data can be collected by water level monitoring equipment installed within the watershed. Then, the water level monitoring data is processed to calculate the water level difference between adjacent time points. For instance, for the first and second time points in the water level monitoring data, their water level difference is calculated, and this difference is used as the first data point in the water level fluctuation sequence. Similarly, the water level difference between all adjacent time points in the water level monitoring data is calculated to obtain the water level fluctuation sequence. During the calculation process, it is essential to ensure that the length of the water level fluctuation sequence is consistent with the length of the deformation record segment, so that water level changes and dam deformation correspond in time. By analyzing the water level fluctuation sequence, the trend of water level changes and the potential impact of water level changes on dam deformation can be understood.

[0056] Step S2024: Extract flow velocity monitoring data from hydrological environmental records, calculate the standard deviation of flow velocity within the window. The standard deviation sequence reflects the intensity of flow velocity fluctuations, and is arranged in chronological order to form a flow velocity gradient sequence. The window size is dynamically adjusted according to the frequency of flow velocity changes.

[0057] For example, flow velocity monitoring data is first extracted from hydrological environmental records. This data can be collected by flow velocity monitoring equipment installed within the watershed. Then, the window size is dynamically adjusted based on the frequency of flow velocity changes. For instance, if frequent flow velocity changes are observed, the window size is set smaller; if slow changes occur, the window size is set larger. Next, the standard deviation of the flow velocity is calculated within each time window. For the flow velocity data within each time window, the standard deviation is calculated using the standard deviation calculation method to obtain the flow velocity standard deviation for that time window. Finally, the flow velocity standard deviations within different time windows are arranged in chronological order to form a flow velocity gradient sequence. By analyzing the flow velocity gradient sequence, the intensity of flow velocity fluctuations can be understood, as well as the potential impact of flow velocity changes on dam deformation.

[0058] Step S2025: Timestamp alignment is performed on the vertical settlement fluctuation sequence, horizontal displacement change sequence, water level fluctuation sequence, and velocity gradient sequence. Based on the unified timestamp, the sampling points of each sequence are matched one by one to generate a multi-dimensional feature matrix.

[0059] For example, firstly, the timestamp information for each feature sequence is obtained. This timestamp information records the collection time of each sampling point. Then, based on a unified timestamp, the sampling points of the vertical settlement fluctuation sequence, horizontal displacement change sequence, water level fluctuation sequence, and velocity gradient sequence are matched one by one. For each timestamp, the sampling points corresponding to the timestamp in all feature sequences are found, and these sampling points are combined together to form a row of a matrix. All timestamps are processed sequentially to finally generate a multi-dimensional feature matrix.

[0060] Step S2026: Standardize the multi-dimensional feature matrix and map the values ​​of each feature sequence to a preset interval to make each feature sequence comparable.

[0061] For example, a preset interval is first determined. The preset interval can be adjusted according to the needs of the analysis and the characteristics of the data. For example, the preset interval can be set as a numerical range, such as [0,1] or [-1,1]. Then, for each feature sequence in the multi-dimensional feature matrix, a standardization method is used to map its values ​​to the preset interval. The standardization method can be something like min-max standardization or z-score standardization.

[0062] Step S203: Dynamically fuse the multi-dimensional feature sequences, assign fusion weights based on the mutual information values ​​between features, assign higher weights to features with higher mutual information values, and generate a fused feature sequence containing deformation-hydrological coupling relationship. The sequence length is the same as the data segment length.

[0063] For example, the mutual information value between features in the multi-dimensional feature sequence is first calculated. A mutual information calculation method can be used to calculate the mutual information value for each pair of features. Then, fusion weights are assigned based on the mutual information values. Features with high mutual information values ​​are assigned higher weights, and features with low mutual information values ​​are assigned lower weights. Next, the multi-dimensional feature sequences are weighted and summed according to the fusion weights. For each time point, the value of each feature sequence at that time point is multiplied by its corresponding fusion weight, and the results are summed to obtain the value of the fused feature sequence at that time point. This process is repeated for all time points to finally generate a fused feature sequence containing the deformation-hydrological coupling relationship. For example, for a vertical settlement fluctuation sequence, a horizontal displacement change sequence, a water level fluctuation sequence, and a velocity gradient sequence, their mutual information values ​​are calculated, fusion weights are assigned based on the mutual information values, and then the four sequences are weighted and summed according to the fusion weights to obtain a fused feature sequence.

[0064] Step S204: Perform trend separation processing on the fused feature sequence to separate the short-term dynamic component and the long-term trend component. The short-term dynamic component reflects the instantaneous deformation response characteristics, and the long-term trend component reflects the cumulative impact characteristics of the hydrological environment.

[0065] For example, trend separation methods such as moving averages or wavelet decomposition can be used. For moving averages, a suitable moving window size is first selected. The choice of moving window size can be adjusted according to the analysis requirements and data characteristics. For example, if short-term dynamic changes are desired, a smaller moving window size can be chosen; if long-term trend changes are the focus, a larger moving window size can be chosen. Then, a moving average is calculated on the fused feature sequence to obtain a smoothed sequence representing the long-term trend component of the fused feature sequence. Next, the long-term trend component is subtracted from the original fused feature sequence to obtain the short-term dynamic component. For wavelet decomposition, the fused feature sequence is decomposed using wavelets to obtain wavelet coefficients at different scales. By selecting appropriate wavelet coefficients, the short-term dynamic component and long-term trend component are reconstructed. For example, for a fused feature sequence containing information on dam deformation and hydrological environment, the moving average method can be used to separate it into short-term dynamic components and long-term trend components.

[0066] Step S205: Input the separated feature components into the encoding processing module, and convert them into fixed-dimensional feature vectors through a multilayer perceptron network. The vector dimension is determined according to the number of data segments and the feature complexity.

[0067] The encoding module, composed of a multilayer perceptron network, encodes and transforms the feature components. The multilayer perceptron network consists of an input layer, hidden layers, and an output layer. The input layer receives the separated feature components as input data, the hidden layers perform nonlinear transformations on the input data, and the output layer converts the transformed result into a feature vector of fixed dimensions. The vector dimension is the number of dimensions of the feature vector, determined by the number of data segments and the feature complexity. The number of data segments reflects the size of the data, while feature complexity reflects the complexity of the features. If the number of data segments is large or the feature complexity is high, a larger vector dimension may be needed to retain more feature information; conversely, if the number of data segments is small or the feature complexity is low, a smaller vector dimension can be used to reduce computation.

[0068] For example, the separated short-term dynamic component and long-term trend component are first input into the input layer of a multilayer perceptron network. The input layer passes the input data to the hidden layers, where neurons perform non-linear transformations on the input data. Each neuron has an activation function, such as the sigmoid function or the ReLU function, which non-linearly maps the input data, converting it into a new output value. There can be multiple hidden layers, each further processing and transforming the input data. Finally, the output layer converts the output of the hidden layers into a fixed-dimensional feature vector. The vector dimension is adjusted according to the number of data segments and the feature complexity.

[0069] Step S206: Calculate the similarity between the feature vector and the preset standard feature template. When the similarity exceeds the matching threshold, it is determined to be a valid feature primitive. Filter all valid feature primitives and arrange them in chronological order to obtain a set of deformation response feature primitives. The size of the set is positively correlated with the number of data segments.

[0070] The preset standard feature template is a pre-defined set of standard features used to measure the similarity of feature vectors. Similarity is an indicator that measures the degree of similarity between a feature vector and the preset standard feature template, and can be calculated using methods such as cosine similarity and Euclidean distance. The matching threshold is a pre-defined standard value used to determine whether a feature vector is a valid feature primitive. When the similarity between a feature vector and the preset standard feature template exceeds the matching threshold, the feature vector is considered a valid feature primitive. The deformation response feature primitive set is a set obtained by arranging all valid feature primitives in chronological order, representing the dam's response characteristics to changes in the hydrological environment at different times. The size of the set is positively correlated with the number of data segments, meaning that the more data segments there are, the more valid feature primitives can be obtained, and the larger the size of the deformation response feature primitive set will be.

[0071] Step S300: Input the set of deformation response feature primitives into the pre-constructed topological relationship map of the watershed reservoir-dam group, and perform anomaly state propagation calculation through the state dependency relationship between map nodes to obtain a set of chain risk transmission paths containing risk propagation intensity and path probability.

[0072] A watershed reservoir-dam group topology map is a graph used to describe the topological structure and interrelationships among reservoirs and dams within a watershed. Nodes represent reservoir / dam entities, and edges represent the hydraulic connections between them. The state dependencies between nodes reflect how a change in the state of one dam affects the states of other dams. Anomaly propagation calculus simulates the propagation of anomalies within the reservoir-dam group based on these state dependencies. Risk propagation intensity refers to the degree to which an anomaly affects other dams during propagation, while path probability refers to the likelihood of an anomaly propagating along a specific path. The cascading risk transmission path set is a set of all possible risk propagation paths obtained through anomaly propagation calculus; each path contains both risk propagation intensity and path probability information.

[0073] In one implementation, step S300 may include the following steps S310-S360:

[0074] Step S310: Call the pre-constructed topological relationship map of the watershed reservoir-dam group. In the map, nodes represent reservoir-dam entities, edges represent hydraulic connection relationships between reservoirs and dams, node attributes include dam structural parameters and bound deformation response feature primitives, and edge attributes include hydraulic conduction coefficients and connection directions.

[0075] For example, a pre-constructed topological map of the watershed reservoir-dam group can be obtained by calling a relevant graph management system or database. This map can be stored in the database in the form of a graph data structure and loaded into memory through query and retrieval operations.

[0076] Step S320: Assign the set of deformation response feature primitives to the corresponding reservoir nodes by matching spatial coordinates. Calculate the spatial distance between the grid coordinates of the feature primitives and the geographic coordinates of the reservoir nodes. The node with the smallest distance is the matching node. Generate a primitive-node association table, which contains the primitive index and the corresponding node identifier.

[0077] Spatial distance is the distance between the grid coordinates of a feature primitive and the geographic coordinates of a reservoir / dam node. It can be calculated using distance calculation methods, such as the Euclidean distance method. A matching node is the reservoir / dam node with the smallest spatial distance to the feature primitive. The primitive-node association table is a table used to record the matching relationship between feature primitives and reservoir / dam nodes. It contains primitive indexes and corresponding node identifiers. This table allows for quick lookup of the reservoir / dam node corresponding to each feature primitive.

[0078] For example, firstly, the grid coordinates of each feature element in the deformation response feature element set and the geographic coordinates of each reservoir-dam node in the watershed reservoir-dam group topology map are obtained. Then, for each feature element, its spatial distance to all reservoir-dam nodes is calculated. The Euclidean distance method can be used to calculate the distance between the feature element and the reservoir-dam node based on the coordinate values. Next, the reservoir-dam node with the smallest spatial distance to the feature element is identified and designated as the matching node. Finally, the index of the feature element and the identifier of the matching node are recorded in the element-node association table. For example, for a set containing multiple deformation response feature elements, the matching node for each feature element is found by calculating the spatial distance between each feature element and all reservoir-dam nodes, and an element-node association table is generated.

[0079] Step S330: Initialize the state parameters of the map nodes. The state parameters include deformation anomaly, hydrological sensitivity, and propagation probability. The deformation anomaly is calculated based on the deviation between the feature primitive and the standard primitive. The propagation probability is initialized to a value that is positively correlated with the deformation anomaly.

[0080] The standard primitives are a pre-defined set of characteristic primitives representing the normal state. By comparing the differences between the characteristic primitives and the standard primitives, the deformation anomaly degree can be obtained. Hydrological sensitivity refers to the degree to which a reservoir / dam node is sensitive to changes in the hydrological environment, reflecting the extent to which the state of the reservoir / dam node is affected by the hydrological environment. Propagation probability refers to the likelihood that a reservoir / dam node will propagate its abnormal state to other reservoir / dam nodes. It is initialized to a value positively correlated with the deformation anomaly degree; that is, the higher the deformation anomaly degree, the greater the propagation probability.

[0081] For example, firstly, the deformation response feature primitives bound to each reservoir-dam node are obtained. Then, each feature primitive is compared with a standard primitive, and the deviation value between them is calculated. The deformation anomaly degree is calculated based on the deviation value, which can be achieved by calculating the distance or similarity between the feature primitive and the standard primitive. Next, the hydrological sensitivity of each reservoir-dam node is determined based on experience or historical data. Finally, the propagation probability is initialized to a value positively correlated with the deformation anomaly degree. For example, for a watershed reservoir-dam group topology map containing multiple reservoir-dam nodes, for each reservoir-dam node, the deviation value between its bound deformation response feature primitive and the standard primitive is calculated to obtain the deformation anomaly degree. The hydrological sensitivity of each reservoir-dam node is determined empirically, and the propagation probability is initialized to a value positively correlated with the deformation anomaly degree.

[0082] Step S340: Identify initial abnormal nodes. When the abnormality of a node deformation exceeds a preset safety threshold, it is marked as an initial abnormal node. Its propagation probability is set to the initial propagation probability. The propagation probability of other nodes is initialized to zero. Record the list of initial abnormal nodes and their corresponding state parameters.

[0083] For example, first, each reservoir-dam node in the topological relationship map of the watershed reservoir-dam group is traversed to obtain its deformation anomaly degree. Then, the deformation anomaly degree of each node is compared with a preset safety threshold. If the deformation anomaly degree of a node exceeds the preset safety threshold, the node is marked as an initial anomalous node, and its propagation probability is set to the initial propagation probability. For nodes whose deformation anomaly degree does not exceed the preset safety threshold, their propagation probability is initialized to zero. Finally, the list of initial anomalous nodes and their corresponding state parameters are recorded, including deformation anomaly degree, hydrological sensitivity, and propagation probability. For example, for a topological relationship map of a watershed reservoir-dam group containing multiple reservoir-dam nodes, each node is traversed, and it is determined whether its deformation anomaly degree exceeds the preset safety threshold. Nodes exceeding the threshold are marked as initial anomalous nodes, and their propagation probability is set to the initial propagation probability. Nodes whose deformation anomaly degree does not exceed the threshold are initialized to zero.

[0084] Step S350: Starting from the initial abnormal node, perform multiple rounds of abnormal state propagation simulation. In each round, the node attempts to activate the adjacent node according to the propagation probability. After successful activation, the state parameters of the adjacent node are updated. The update magnitude is positively correlated with the edge propagation coefficient. Record the node sequence and activation time of each propagation path.

[0085] The propagation probability is the probability that each node will transmit an abnormal state to its neighboring nodes, which has been initialized in previous steps based on factors such as the degree of deformation anomaly of the nodes. Neighboring nodes refer to other nodes in the topological relationship graph of the watershed reservoir-dam group that are connected to the current node through edges; these edges represent the hydraulic connections between the reservoirs and dams.

[0086] During each round of propagation, for each initial anomalous node or a node activated in the previous round, attempts are made to activate neighboring nodes based on their propagation probability. The specific activation operation can be understood as a probabilistic event; for example, a random number can be generated and compared with the node's propagation probability. If the random number is less than the propagation probability, activation is considered successful. Once a neighboring node is successfully activated, its state parameters need to be updated. State parameters include deformation anomaly degree, hydrological sensitivity, and propagation probability. The update magnitude is positively correlated with the edge conduction coefficient, an attribute of the edges in the topological relationship graph of the watershed reservoir-dam group, reflecting the ability of water flow to conduct between reservoirs and dams. The larger the edge conduction coefficient, the closer the hydraulic connection between the two reservoirs and dams, and the greater the propagation impact of the anomalous state; therefore, the larger the update magnitude of the state parameters of neighboring nodes.

[0087] Step S360: Calculate the cumulative propagation intensity and path probability of each propagation path. The propagation intensity is the product of the transmission coefficients of each side of the path, and the path probability is the product of the activation probabilities of each node. Filter the paths whose propagation intensity is greater than the intensity threshold and whose path probability is greater than the probability threshold, and sort them in ascending order of propagation steps to obtain the set of chain risk transmission paths.

[0088] Cumulative propagation intensity is an indicator that measures the impact of an abnormal state propagation along a propagation path. It is obtained by multiplying the conduction coefficients of each edge of the path. The edge conduction coefficient reflects the water flow conduction capacity between the reservoir and dam; the larger the product of the conduction coefficients, the stronger the propagation impact of the abnormal state along this path. Path probability refers to the likelihood of an abnormal state propagating along a certain path, obtained by multiplying the activation probabilities of each node along the path. The node activation probability is the probability that each node will be activated during the propagation simulation of the abnormal state; the larger the product of the activation probabilities, the higher the likelihood of the abnormal state propagating along this path. Intensity threshold and probability threshold are two pre-set standard values ​​used to screen paths with high propagation impact and probability. Only when the cumulative propagation intensity of a propagation path is greater than the intensity threshold and the path probability is greater than the probability threshold is the path considered a meaningful and risk transmission path that requires attention. Sort by propagation steps in ascending order means sorting the screened paths according to the number of nodes contained in the path minus one. The fewer the propagation steps, the shorter the distance and the faster the abnormal state propagates. For paths with the same number of propagation steps, they are then sorted in descending order of cumulative propagation intensity, with priority given to paths with higher propagation intensity. The final set of chain risk transmission paths includes those with significant propagation impact, high probability of propagation, and relatively fast propagation speed.

[0089] In one implementation, step S360 may include the following steps S361-S366:

[0090] Step S361: Traverse all propagation paths and extract the edge set and node set contained in each propagation path. The edge set contains the hydraulic conduction coefficient of each segment of the path, and the node set contains the activation probability of each node. The path identifier is associated with the propagation round.

[0091] Path identifiers are information used to uniquely identify each propagation path and are associated with the propagation round, thus clearly indicating in which round each path was formed. By associating path identifiers with propagation rounds, the propagation process and temporal patterns of abnormal states can be better analyzed.

[0092] For example, by traversing the database or data structure that records the propagation path, for each path, the edge and node information contained therein can be extracted from the topological relationship map of the watershed reservoir-dam group. The hydraulic conduction coefficient of the edge and the activation probability of the node can be collected into the edge set and the node set respectively, and the path identifier and propagation round association information can be recorded.

[0093] Step S362: Calculate the cumulative propagation intensity of the path. From the starting node to the ending node of the path, multiply the hydraulic conduction coefficients of each side in turn. The product is the cumulative propagation intensity of the path, which reflects the risk transmission capacity of the path.

[0094] The cumulative propagation strength of a propagation path is a process of quantitatively assessing the risk transmission capacity of each path. Starting from the initial node of the path, each edge is visited sequentially along the path, and the hydraulic conductivity coefficients of each edge are multiplied together. The final product is the cumulative propagation strength of the path. The greater the cumulative propagation strength, the closer the hydraulic connection between the reservoir and dam along the path, and the greater the propagation impact of the anomalous state along the path; that is, the stronger the risk transmission capacity of the path. For example, for each propagation path with an extracted edge set, starting from the first edge of the edge set, its hydraulic conductivity coefficient is used as the initial value, and then it is multiplied by the hydraulic conductivity coefficients of subsequent edges in the edge set until all edges have been traversed. The final result is the cumulative propagation strength of the path.

[0095] Step S363: Calculate the path probability of the path. From the starting node to the ending node of the path, multiply the activation probabilities of each node in turn. The product is the path probability of the path, which reflects the likelihood of the path occurring.

[0096] Calculating the path probability is a process of quantifying the likelihood of each propagation path occurring. Starting from the initial node of the path, each node is visited sequentially along the path, and the activation probabilities of each node are multiplied together. The final product is the path probability. A higher path probability indicates a higher probability that an abnormal state is propagating along this path. For example, if the activation probabilities of all nodes on a path are relatively high, their product will also be relatively high, indicating that an abnormal state is very likely to propagate along this path.

[0097] For example, for each propagation path that has extracted a set of nodes, starting from the first node in the set, its activation probability is used as the initial value, and then multiplied by the activation probabilities of subsequent nodes in the set in turn, until all nodes have been traversed. The final result is the path probability of that path.

[0098] Step S364: Obtain the preset intensity screening threshold and probability screening threshold. The intensity screening threshold is determined based on the statistical analysis of the transmission intensity of historical risk events in the watershed, and the probability screening threshold is set based on the path reliability requirements.

[0099] The probability screening threshold is set based on path reliability requirements. Path reliability requirements are determined according to actual risk assessment and management needs. For example, if the tolerance for risk is low and it is desired to focus only on paths with a very high probability of occurrence, then the probability screening threshold can be set higher; if the tolerance for risk is high, the probability screening threshold can be set lower.

[0100] For example, by querying relevant historical data records or risk assessment models, one can obtain the intensity screening threshold determined based on the statistical analysis of the transmission intensity of historical risk events in the watershed, as well as the probability screening threshold set according to the path reliability requirements.

[0101] Step S365: Perform double screening on the propagation paths, retaining paths whose cumulative propagation intensity is greater than the intensity screening threshold and whose path probability is greater than the probability screening threshold, and deleting paths that do not meet either threshold condition.

[0102] For each propagation path, its cumulative propagation intensity is compared with an intensity screening threshold, and its path probability is compared with a probability screening threshold. If both conditions are met, it indicates that the path has both sufficient risk transmission capacity and a high probability of occurrence, and it is retained; if neither condition is met, it indicates that the path has a low risk level or a low probability of occurrence, and it is removed from the propagation path set.

[0103] For example, all propagation paths are traversed. For each path, its cumulative propagation strength and path probability are checked to see if they meet the filtering criteria. If the criteria are met, the path is marked as retained; otherwise, it is marked as deleted. Finally, all paths marked as deleted are removed from the propagation path set, resulting in the filtered propagation path set.

[0104] Step S366: Sort the filtered paths in ascending order by the number of propagation steps (i.e., the number of nodes contained in the path minus one), and sort the paths with the same number of propagation steps in descending order by cumulative propagation intensity. Extract the sorted path set as the chain risk transmission path set.

[0105] Sort the filtered paths to better organize and analyze these high-risk propagation paths. The propagation step count is the number of nodes in the path minus one, reflecting the number of intermediate nodes the abnormal state passes through from the starting node to the ending node. The fewer the propagation steps, the shorter the distance the abnormal state propagates and the faster it travels.

[0106] For example, a sorting algorithm can be used to sort the filtered set of paths. First, the paths are sorted according to the number of propagation steps. If the number of propagation steps is the same, a second sort is performed based on the cumulative propagation strength. After sorting, the sorted set of paths is extracted as the set of cascading risk transmission paths.

[0107] As another implementation of step S300, the following steps S301-S306 may be included:

[0108] Step S301: Load the topological relationship map data of the watershed reservoir-dam group. The map is stored in a graph structure, including a node set and an edge set. The node set stores the spatial location and deformation feature primitives of the reservoir-dam group, and the edge set stores the hydraulic connection relationship and conduction parameters between the nodes.

[0109] The node set stores relevant information about each reservoir and dam within the watershed, including their spatial location and associated deformation feature primitives. The spatial location information of the reservoirs and dams can be represented using geographic coordinates, which helps determine their relative positions. The deformation feature primitives, generated in previous steps, represent the deformation response characteristics of the reservoirs and dams under different hydrological environments. Storing them in the nodes provides a basis for subsequent anomaly analysis. The edge set stores the hydraulic connectivity relationships and conduction parameters between nodes. Hydraulic connectivity relationships describe which reservoirs and dams are hydraulically connected, i.e., water can flow between them. Conduction parameters are attributes of the edges, reflecting the ability of water flow to be conducted between reservoirs and dams; for example, they can be represented by edge conduction coefficients.

[0110] Step S302: Based on the spatial coordinates of the reservoir node, bind the feature primitives in the deformation response feature primitive set to the corresponding node. The matching process calculates the Euclidean distance between the feature primitive and the node. Primitives with a distance less than the matching threshold are assigned to the node, and a primitive assignment table is generated.

[0111] The matching threshold is a pre-defined standard value used to determine whether a feature primitive should be assigned to a reservoir node. If the Euclidean distance between a feature primitive and a reservoir node is less than the matching threshold, it indicates that the feature primitive and the node are spatially close, suggesting a correlation, and the feature primitive is assigned to that node. The primitive assignment table records the matching relationships between feature primitives and reservoir nodes, containing the identifier of the feature primitive and the corresponding reservoir node identifier. Using the primitive assignment table, the reservoir node corresponding to each feature primitive can be quickly queried, providing a foundation for subsequent abnormal state propagation calculations.

[0112] For example, first, obtain the spatial coordinates of each feature element in the deformation response feature element set and the spatial coordinates of each reservoir-dam node in the topological relationship map of the watershed reservoir-dam group. Then, for each feature element, calculate its Euclidean distance to all reservoir-dam nodes. Reservoir-dam nodes with Euclidean distances less than a matching threshold are selected as candidate matching nodes. If multiple candidate matching nodes exist, the node with the closest distance is selected as the final matching node. Record the identifiers of the feature elements and the final matching node in the element allocation table.

[0113] Step S303: Obtain the preset initial value of the node state vector. The state vector includes three dimensions: deformation deviation, environmental correlation, and propagation capability. The deformation deviation is calculated based on the similarity between the feature primitive and the standard primitive. The propagation capability is positively correlated with the deformation deviation.

[0114] The preset initial values ​​of the node state vector are a set of values ​​set for the initial state of each reservoir-dam node in the topological relationship map of the watershed reservoir-dam group. The state vector includes three dimensions: deformation deviation, environmental correlation, and propagation capability.

[0115] Deformation deviation is an indicator that measures the degree of deviation between the deformation of a reservoir / dam node and the standard condition. It is calculated based on the similarity between feature elements and standard elements. Standard elements are a pre-defined set of feature elements representing the normal state. By comparing the similarity between feature elements and standard elements, the deformation deviation can be obtained. The lower the similarity, the greater the difference between the deformation of the reservoir / dam node and the standard condition, and the higher the deformation deviation.

[0116] Environmental correlation reflects the degree of connection between a reservoir / dam node and its surrounding hydrological environment, and can be determined based on factors such as the reservoir / dam's geographical location and hydrological monitoring data. Propagation capacity refers to the ability of a reservoir / dam node to propagate abnormal conditions to other reservoir / dam nodes, and is positively correlated with deformation deviation. That is, the higher the deformation deviation, the more severe the abnormal condition of the reservoir / dam node, and the stronger its ability to propagate the abnormal condition, resulting in a higher propagation capacity value.

[0117] Step S304: Determine the initial source of abnormal propagation. When the deformation deviation of a node exceeds the design safety threshold, it is marked as an abnormal propagation source. Set its propagation capability value to the baseline propagation probability. The propagation capability values ​​of other nodes are initialized to zero. Record the node identifier and state vector of the abnormal propagation source.

[0118] The initial source of anomaly propagation refers to the reservoir / dam node in the watershed's reservoir-dam group whose deformation deviation exceeds the design safety threshold. The design safety threshold is a pre-set standard value used to determine whether a reservoir / dam node is in an abnormal state. When the deformation deviation of a node exceeds the design safety threshold, it indicates that the node's deformation exceeds the normal range and an abnormal state may exist; therefore, it is marked as an anomaly propagation source. The baseline propagation probability is the initial likelihood that the anomaly propagation source will propagate the abnormal state to other reservoir / dam nodes. For nodes marked as anomaly propagation sources, their propagation capability value is set to the baseline propagation probability to initiate the propagation process of the abnormal state. Other nodes refer to reservoir / dam nodes whose deformation deviation does not exceed the design safety threshold. In the initial stage, they are considered to be in a normal state and will not propagate the abnormal state; therefore, their propagation capability value is initialized to zero.

[0119] Recording the node identifier and state vector of the anomaly propagation source is essential for accurately tracking and analyzing the propagation process of the anomaly in subsequent anomaly propagation calculations. The node identifier uniquely identifies each reservoir node, and the state vector contains information such as the node's deformation deviation, environmental correlation, and propagation capability value.

[0120] Step S305: Simulate the propagation of abnormal states based on the independent cascade model. The propagation source node attempts to activate adjacent nodes according to the propagation capability value. After successful activation, the state vector of the adjacent node is updated. The deformation deviation increases the deviation contribution of the propagation source. The environmental correlation is adjusted according to the edge attribute. The propagation capability value is set to the current activation probability.

[0121] The process of a source node attempting to activate its neighboring nodes based on its propagation capability value can be understood as a probabilistic event. For example, a random number can be generated and compared to the source node's propagation capability value; if the random number is less than the propagation capability value, activation is considered successful. Once a neighboring node is successfully activated, its state vector needs to be updated. The state vector contains three dimensions: deformation deviation, environmental correlation, and propagation capability value. Deformation deviation increases the deviation contribution of the propagation source; that is, the deformation deviation of a neighboring node increases due to the abnormal propagation of the source node, and the magnitude of this increase is determined by the source node's deformation deviation. Environmental correlation is adjusted based on edge attributes, which are the relevant properties of the edges connecting the source node and its neighboring nodes, such as edge conduction coefficients. These edge attributes can adjust the degree of correlation between the neighboring node and its surrounding hydrological environment. The propagation capability value is set as the current activation probability, calculated based on the source node's propagation capability value and edge conduction coefficients, reflecting the likelihood that the neighboring node will continue to propagate the abnormal state after activation.

[0122] In one implementation, step S305 may include the following steps S3051-S3056:

[0123] Step S3051: Obtain preset propagation simulation parameters, including maximum propagation rounds, propagation time step, and probability decay factor. The maximum propagation rounds are dynamically determined based on the number of reservoirs and dams in the basin. The time step is consistent with the sampling period of BeiDou monitoring data. The decay factor decreases as the propagation rounds increase.

[0124] The preset propagation simulation parameters are crucial for simulating the propagation of abnormal states. These parameters include the maximum number of propagation rounds, the propagation time step, and the probability decay factor. The maximum number of propagation rounds refers to the maximum number of rounds allowed during the abnormal state propagation simulation, dynamically determined based on the number of reservoirs and dams in the basin. The more reservoirs and dams in the basin, the higher the probability and complexity of abnormal state propagation, thus increasing the maximum number of propagation rounds accordingly. The propagation time step refers to the time interval between each propagation simulation, consistent with the sampling period of the BeiDou monitoring data.

[0125] The probability decay factor is a parameter used to control the change in the probability of an abnormal state propagation with each propagation round; it decreases as the number of propagation rounds increases. As the number of propagation rounds increases, the propagation influence of the abnormal state gradually weakens. This decay effect can be simulated using the probability decay factor. For example, in the first round of propagation, the propagation probability may be relatively high, but as the number of propagation rounds increases, the propagation probability will gradually decrease because the abnormal state is affected by various factors during propagation, causing it to gradually spread and weaken.

[0126] Step S3052: During the first round of propagation, the initial abnormal propagation source node is added to the set of active nodes, its state vector is recorded, the deformation deviation is the deviation between the current monitoring value and the safety value, the environmental correlation is the hydrological coupling coefficient in the primitive, and the propagation capability value is the initial baseline probability.

[0127] The first round of propagation is the initial stage of the anomaly propagation simulation. The initial source nodes for anomaly propagation are determined in the previous steps, and these nodes are added to the active node set. The active node set is a collection used to store nodes that are currently propagating an anomaly.

[0128] The state vector of the initial anomaly propagation source node is recorded. The state vector includes three dimensions: deformation deviation, environmental correlation, and propagation capability. Deformation deviation is the deviation between the current monitoring value and the safety value. The current monitoring value is the actual deformation value of the reservoir node obtained through BeiDou monitoring data and other means. The safety value is a pre-set deformation value representing the normal state of the reservoir node. The deviation between the two reflects the degree of anomaly of the reservoir node.

[0129] The environmental correlation degree is the hydrological coupling coefficient in the primitive. The hydrological coupling coefficient is an attribute in the deformation response characteristic primitive, which reflects the degree of correlation between the reservoir / dam node and the surrounding hydrological environment.

[0130] The propagation capability value is the initial baseline probability, which is the propagation capability value set for the abnormal propagation source node in the previous steps. It represents the probability that the abnormal propagation source node will propagate the abnormal state in the first round of propagation.

[0131] Step S3053: In each round of propagation, traverse the nodes in the active node set, extract all outgoing neighbor nodes (i.e., downstream nodes along the water flow direction), calculate the activation probability of the current node for each neighbor node, and the activation probability is the product of the current node's propagation capability value and the edge propagation coefficient, multiplied by the attenuation factor.

[0132] In each round of propagation, each node in the active node set needs to be processed. The nodes in the active node set are those currently propagating the abnormal state. For each active node, all its outgoing neighbor nodes need to be extracted. Outgoing neighbor nodes refer to nodes located downstream of the node along the water flow direction in the topological relationship graph of the watershed reservoir-dam group. These nodes are the target nodes that the abnormal state may propagate to.

[0133] The activation probability of the current node to each of its neighbors is calculated. This activation probability represents the likelihood that the current node will propagate its abnormal state to its neighbors. The activation probability is calculated by multiplying the current node's propagation capability value by the conduction coefficient of the edges connecting them, and then multiplying by a decay factor. The current node's propagation capability value reflects its ability to propagate abnormal states, the edge conduction coefficient reflects the strength of the hydraulic connection between the two nodes, and the decay factor controls how the propagation probability decreases with each propagation round. For example, if an active node has a large propagation capability value, a large edge conduction coefficient connecting it to its neighbors, and a large decay factor at shorter propagation rounds, then the activation probability will be relatively high.

[0134] Step S3054: Generate a random value for each neighbor node. When the random value is less than the activation probability, the neighbor node is successfully activated and added to the new active node set. Update its state vector. The deformation deviation is the weighted sum of the deformation deviation of the current node and its original deviation. The environmental correlation is the hydrological sensitivity in the edge attribute. The propagation capability value is the current activation probability.

[0135] Generating random values ​​for each neighboring node is to simulate the probabilistic nature of abnormal state propagation. The generated random values ​​are generally in the range [0,1]. The generated random values ​​are compared with the calculated activation probability; if the random value is less than the activation probability, it means that the neighboring node has been successfully activated.

[0136] Once a neighboring node is successfully activated, it is added to the newly active node set. The newly active node set is a collection used to store newly activated nodes in this round, and these nodes will participate in the next round of abnormal state propagation.

[0137] Updating the state vector of neighboring nodes involves updating three dimensions: deformation deviation, environmental correlation, and propagation capability. Deformation deviation is a weighted sum of the current node's deformation deviation and its original deviation. The weighted sum is calculated by combining the current node's deformation deviation with the original deviation of its neighboring nodes according to set weights, thus comprehensively considering the impact of the current node's abnormal propagation and the initial abnormality of the neighboring nodes. Environmental correlation is the hydrological sensitivity attribute in the edge properties. Hydrological sensitivity is an attribute of the edge connecting the current node and its neighboring nodes, reflecting the degree of correlation between the neighboring node and its surrounding hydrological environment.

[0138] Step S3055: The propagation rounds are incremented, the active node set is updated to a new active node set, and the propagation process is repeated until the maximum propagation rounds are reached or the new active node set is empty. Record the node activation status and state vector change trajectory of each propagation round.

[0139] Increasing the number of propagation rounds indicates that the simulation of abnormal state propagation has entered the next round. The set of active nodes is updated to the new set of active nodes, meaning that the newly activated nodes in this round become the nodes for propagating abnormal states in the next round.

[0140] The repeated propagation process refers to executing steps S3053-S3054 again, which involves traversing the nodes in the active node set, calculating the activation probability, attempting to activate neighboring nodes, and updating the state vector of the activated nodes. This process repeats until the maximum propagation rounds are reached or the new active node set is empty. Reaching the maximum propagation rounds indicates that the abnormal state propagation simulation has run for a sufficient number of rounds, and the simulation process ends; an empty new active node set indicates that no new nodes have been activated in the current round, and the abnormal state propagation has stopped.

[0141] Step S3056: After the propagation is completed, organize all propagation paths. Each propagation path includes the sequence of activated nodes, the activation time of each node, and state parameters. Calculate the cumulative propagation strength and path probability of the propagation path. The calculation method is similar to the previous steps.

[0142] The calculation of the cumulative propagation strength and path probability of a propagation path is similar to the previous steps. The cumulative propagation strength is the result of multiplying the propagation coefficients of each edge from the starting node to the ending node, reflecting the risk propagation capability of the path. The path probability is the result of multiplying the activation probabilities of each node from the starting node to the ending node, reflecting the likelihood of the path occurring.

[0143] Step S306: During the propagation simulation, record the node sequence, activation time, and state vector changes of each propagation path, calculate the cumulative propagation intensity and path probability of the path, and filter the paths that meet the threshold conditions by grouping and sorting them according to the number of propagation steps to obtain a set of chain risk transmission paths.

[0144] The methods for calculating the cumulative propagation strength and path probability of a path are consistent with those described in the previous steps. The cumulative propagation strength is obtained by multiplying the propagation coefficients of each edge on the path, reflecting the risk propagation capability of the path; the path probability is obtained by multiplying the activation probabilities of each node on the path, reflecting the likelihood of the path occurring.

[0145] Filtering paths that meet the threshold conditions means comparing the calculated cumulative propagation intensity and path probability with the preset intensity filtering threshold and probability filtering threshold, and only retaining paths whose cumulative propagation intensity is greater than the intensity filtering threshold and whose path probability is greater than the probability filtering threshold.

[0146] Grouping and sorting by propagation steps involves grouping the filtered paths according to the number of propagation steps (i.e., the number of nodes in the path minus one), and then sorting them in descending order of cumulative propagation intensity within each group. Paths with fewer propagation steps indicate that the abnormal state propagates quickly, while paths with higher cumulative propagation intensity indicate stronger risk transmission capabilities.

[0147] The final set of chain risk transmission paths includes those with high risk transmission capacity, high probability of occurrence, and fast propagation speed. These paths are of great significance for risk early warning and management of river basin reservoir and dam groups.

[0148] Step S400: Based on the set of chain risk transmission paths, perform risk evolution trend deduction on the reservoir nodes in each path, and generate a risk evolution index matrix that includes the node risk accumulation rate and critical threshold spacing.

[0149] Risk evolution trend projection is a process of predicting and analyzing the changing trends of risk at reservoir / dam nodes along a set of cascading risk transmission paths over time. This set of cascading risk transmission paths includes a series of possible risk propagation paths, and each reservoir / dam node along each path may be affected by abnormal conditions, causing its risk status to change over time. The node risk accumulation rate refers to the speed at which the risk at a reservoir / dam node accumulates per unit time, reflecting the rate of risk growth. The critical threshold gap refers to the distance between the current risk value of a reservoir / dam node and a safety critical threshold. The safety critical threshold is a pre-set standard value used to determine whether a reservoir / dam node is in a dangerous state; the smaller the critical threshold gap, the closer the node is to a dangerous state.

[0150] In one implementation, step S400 may include the following steps S410-S460:

[0151] Step S410: Decompose the set of chain risk transmission paths, split each path into a sequence of reservoir and dam nodes connected in sequence, record the node connection order, the corresponding path propagation strength and path probability value, and generate a path-node association table. Each row in the table contains the path identifier, node number and node name.

[0152] Corresponding path propagation intensity refers to the cumulative propagation intensity of each path, reflecting the risk transmission capacity of that path. Path probability value refers to the probability of each path occurring, reflecting the likelihood of that path occurring. Recording path propagation intensity and path probability values ​​can provide important reference information for subsequent risk analysis.

[0153] The path-node association table is used to record the relationships between paths and nodes. Each row in the table contains a path identifier, a node sequence number, and a node name. The path identifier uniquely identifies each path, the node sequence number indicates the node's order within the path, and the node name identifies the specific reservoir / dam node. Using the path-node association table, you can quickly query which path each node belongs to and its position within that path.

[0154] Step S420: For each reservoir node in the path-node association table, retrieve the deformation monitoring record of the corresponding node from the historical monitoring database. The record contains the settlement, horizontal displacement and tilt angle data within the preset time period. Sort the data in ascending order by timestamp to obtain the node's historical deformation sequence. The sequence length matches the historical time period.

[0155] Deformation monitoring records include data on settlement, horizontal displacement, and tilt angle over a predetermined time period. Settlement reflects the vertical subsidence of the dam, horizontal displacement reflects its horizontal movement, and tilt angle reflects its degree of tilt. This data helps analyze the structural stability and deformation trends of the dam.

[0156] Sort by timestamp in ascending order, this sorts the acquired deformation monitoring records according to the chronological order of their acquisition time, resulting in a historical deformation sequence for each node. The sequence length matches the historical duration, meaning the number of data points in the sequence corresponds to the number of monitoring sessions within the preset historical duration.

[0157] In one implementation, step S420 may include the following steps S421-S426:

[0158] Step S421: Query the deformation monitoring records of the node from the BeiDou monitoring historical database using the node's unique identifier, filter valid records with timestamps within a preset historical period, and exclude record entries with a data missing rate exceeding the allowable threshold.

[0159] A node unique identifier is information used to uniquely identify each reservoir / dam node. This identifier allows for accurate retrieval of deformation monitoring records for the corresponding node from the BeiDou monitoring historical database. The BeiDou monitoring historical database is a dedicated database storing reservoir / dam deformation monitoring data collected by BeiDou satellite positioning terminals.

[0160] For example, using the node's unique identifier as a query condition, deformation monitoring records for the corresponding node are retrieved from the BeiDou monitoring historical database. Records within a preset historical time period are filtered based on timestamps, the data missing rate for each record is calculated, and records with a data missing rate exceeding the allowable threshold are excluded.

[0161] Step S422: Perform time standardization processing on the selected valid records, and use linear interpolation to fill in the missing data points with uneven sampling time intervals, so that the sampling timestamps of settlement, horizontal displacement and tilt angle are completely consistent, and obtain a time-aligned monitoring dataset.

[0162] Time standardization of the selected valid records is necessary to ensure temporal consistency among different types of deformation monitoring data (settlement, horizontal displacement, and tilt angle). In actual monitoring, the sampling times of different monitoring devices may vary, leading to inconsistent timestamps, which can complicate subsequent analysis.

[0163] Linear interpolation is a method used to fill in missing data points. When the sampling time interval is uneven, data gaps may occur. Linear interpolation can estimate the value of the missing data point based on adjacent valid data points. For example, if there is a missing data point between two known data points, the value of the missing data point can be estimated using a linear relationship between the values ​​of the two known data points and the time interval.

[0164] Step S423: Extract key deformation parameters from the time-aligned monitoring dataset, including vertical settlement, horizontal displacement and tilt angle. Each parameter is arranged in ascending order by timestamp to form a single-parameter sequence, and the parameter units are consistent with the original records.

[0165] Extracting key deformation parameters from the time-aligned monitoring dataset aims to further focus on parameters that significantly impact the stability of the reservoir / dam structure. Vertical settlement reflects the dam's vertical subsidence, horizontal displacement reflects its horizontal movement, and tilt angle reflects the degree of tilt. Each parameter is arranged in ascending order by timestamp to form a single-parameter sequence; that is, vertical settlement, horizontal displacement, and tilt angle are arranged chronologically, forming three independent sequences. The units of the parameters are consistent with the original records to ensure data accuracy and comparability.

[0166] Step S424: Perform outlier detection and correction on each single-parameter sequence. Use an outlier identification method based on interquartile range to mark outlier data points. Replace outlier points with a weighted average of neighboring normal data points to preserve the temporal continuity and trend characteristics of the sequence.

[0167] The interquartile range (IVR) is the difference between the 75th percentile and the 25th percentile in a data sequence. By calculating the IVR, a range of normal data can be determined. If a data point exceeds this range, it is marked as an outlier. Outliers are replaced using a weighted average of their neighboring normal data points. This means that a reasonable value is estimated based on the values ​​of the normal data points surrounding the outlier and used to replace it. This preserves the temporal continuity and trend of the sequence, preventing outliers from disrupting the overall trend of the data sequence.

[0168] Step S425: Combine the corrected single-parameter sequences in timestamp order to obtain a node historical deformation sequence containing multiple parameters. The sequence length is the number of valid sampling points within the historical period, and each sampling point contains the measured values ​​of multiple deformation parameters.

[0169] The modified single-parameter sequences are combined in timestamp order to merge the three single-parameter sequences of vertical settlement, horizontal displacement, and tilt angle into a single sequence containing multiple parameters. This allows for a more comprehensive consideration of information from multiple deformation parameters, providing a more complete picture of the reservoir's deformation. The length of the node's historical deformation sequence is equal to the number of valid sampling points within the historical period; that is, the sequence contains all valid sampling data within the preset historical period. Each sampling point contains measurements of multiple deformation parameters; for example, at a certain time point, the sampling point contains measurements of vertical settlement, horizontal displacement, and tilt angle at that moment.

[0170] Step S426: Perform smoothing and denoising processing on the historical deformation sequence of the nodes. The size of the filtering window is determined according to the sampling frequency of the monitoring data to ensure that the smoothed sequence retains the original trend characteristics and filters out short-term random fluctuations.

[0171] The size of the filtering window is determined based on the sampling frequency of the monitoring data. A higher sampling frequency indicates a shorter data acquisition time interval, allowing for a relatively smaller filtering window; conversely, a lower sampling frequency allows for a relatively larger filtering window. The size of the filtering window affects the smoothing effect. A suitable filtering window size can effectively filter out short-term random fluctuations while preserving the original trend characteristics. For example, the filtering window size is determined based on the sampling frequency of the monitoring data, and an appropriate filtering method (such as moving average filtering) is used to smooth and denoise the historical deformation sequence of the nodes. This ensures that the smoothed sequence retains the trend characteristics of the original sequence, such as upward or downward trends, while filtering out short-term random fluctuations.

[0172] Step S430: Correlate the node’s historical deformation sequence with the propagation intensity of the corresponding path in time, calculate the node’s activation time in the path based on the number of propagation steps, and use the activation time as a benchmark to extract deformation data of a preset duration before and after as an analysis window. The window length is dynamically adjusted according to the number of propagation steps in the path.

[0173] The activation time of a node in the path is calculated based on the number of propagation steps. The number of propagation steps refers to the number of nodes that the abnormal state passes through from the initial abnormality source node to the current node minus one. Based on the number of propagation steps and the propagation time step (set in the previous steps), the specific time when the node is activated can be calculated.

[0174] The analysis window uses deformation data for a preset time period before and after the activation time as a baseline. The preset time period is a pre-defined time range. By extracting deformation data within this time period, the deformation changes of nodes under the influence of abnormal state propagation can be focused on. The window length is dynamically adjusted according to the number of propagation steps. The more propagation steps, the farther the abnormal state propagates and the longer the time, and the longer the analysis window can be increased accordingly to more comprehensively observe the deformation change trend of nodes.

[0175] Step S440: Construct a risk trend prediction model, taking the historical deformation sequence of nodes and the path propagation intensity within the analysis window as input, capturing the temporal dependence of the sequence through a recurrent neural network, predicting the node deformation sequence within a preset time period in the future, and outputting the predicted values ​​containing each monitoring dimension.

[0176] The analysis window uses the node's historical deformation sequence and path propagation intensity as input. The historical deformation sequence contains the node's past deformation under the influence of anomalous state propagation, while the path propagation intensity reflects the impact of the anomalous state propagation along the path. Using these two pieces of information as input allows for a comprehensive consideration of the node's historical deformation and the impact of anomalous state propagation. The system predicts the node's deformation sequence within a preset timeframe. This preset timeframe is a pre-defined period, and the risk trend prediction model can predict the node's deformation within this future timeframe. The output includes predicted values ​​for various monitoring dimensions, such as vertical settlement, horizontal displacement, and tilt angle.

[0177] Step S450: Calculate the node risk accumulation rate from the predicted deformation sequence, calculate the deformation difference between adjacent predicted points in chronological order, sum all positive differences and divide by the prediction duration to obtain the average risk growth rate per unit time. A positive rate value indicates that the risk continues to accumulate.

[0178] The predicted deformation sequence is a node deformation sequence within a preset future timeframe obtained through a risk trend prediction model. It includes predicted values ​​for multiple monitoring dimensions, such as vertical settlement, horizontal displacement, and tilt angle. The deformation difference between adjacent predicted points is calculated chronologically; that is, for each monitoring dimension, the difference between the values ​​of two adjacent predicted points is calculated sequentially. A positive difference indicates increasing deformation, while a negative difference indicates decreasing deformation. All positive differences are summed up; only positive differences are considered because they indicate increasing risk, while negative differences indicate decreasing risk. The summation of all positive differences is then divided by the predicted timeframe to obtain the average risk growth rate per unit time. The predicted timeframe is the preset future timeframe set in the previous steps. A positive rate value indicates continuous risk accumulation, meaning that the risk at the node is constantly increasing over time.

[0179] In one implementation, step S450 may include the following steps S451-S456:

[0180] Step S451: Extract the predicted values ​​of each monitoring dimension from the deformation prediction sequence, and arrange them in chronological order to form a vertical settlement prediction subsequence, a horizontal displacement prediction subsequence, and a tilt angle prediction subsequence. The length of each subsequence is consistent with the number of prediction cycles.

[0181] Extracting predicted values ​​for each monitoring dimension from the deformation prediction sequence is to process the predicted data for different monitoring dimensions separately. The deformation prediction sequence is a sequence of node deformations within a preset future time period obtained through a risk trend prediction model, containing predicted values ​​for multiple monitoring dimensions such as vertical settlement, horizontal displacement, and tilt angle. Arranged chronologically, it forms vertical settlement prediction subsequences, horizontal displacement prediction subsequences, and tilt angle prediction subsequences; that is, the predicted values ​​for each monitoring dimension are arranged in chronological order to form independent subsequences. The length of each subsequence is consistent with the number of prediction periods, which refers to the number of sampling points within the preset future time period; that is, each subsequence contains the same number of predicted data points.

[0182] Step S452: Perform time difference processing on each predicted subsequence, calculate the deformation increment of adjacent periods, i.e., the predicted value of the next period minus the predicted value of the previous period, to obtain the vertical settlement increment sequence, the horizontal displacement increment sequence, and the tilt angle increment sequence.

[0183] Temporal differencing is performed on each predicted subsequence to calculate the deformation change of each monitoring dimension between adjacent periods. The vertical settlement increment sequence, horizontal displacement increment sequence, and tilt angle increment sequence record the changes in vertical settlement, horizontal displacement, and tilt angle between adjacent periods, respectively. Positive increment values ​​indicate an increase in deformation for that monitoring dimension, while negative increment values ​​indicate a decrease in deformation.

[0184] Step S453: Identify positive increment values ​​in each increment sequence, which indicate risk growth, and ignore negative increment values, which indicate risk reduction. Sum the positive increment values ​​of each dimension to obtain the total increment in the vertical direction, the total increment in the horizontal direction, and the total increment in the tilt direction.

[0185] Identifying positive increments in each increment sequence is to focus on situations where risk is increasing. A positive increment indicates increasing deformation in that monitoring dimension, meaning the node's risk is growing; a negative increment indicates decreasing deformation, meaning the node's risk is decreasing, and negative increments are ignored when calculating the risk accumulation rate. The positive increments for each dimension are summed separately: the positive increments in the vertical settlement increment sequence, the horizontal displacement increment sequence, and the tilt angle increment sequence are added together to obtain the total increments in the vertical, horizontal, and tilt directions, respectively. These three total increments reflect the risk accumulation at the node in the vertical, horizontal, and tilt directions over a predetermined future time period.

[0186] Step S454: The total increment in each direction is weighted and summed. The weights are determined based on the degree of influence of each dimension on the safety of the dam body. Vertical settlement has the highest weight, followed by horizontal displacement, and tilt angle has the lowest weight, thus obtaining the total increment of comprehensive risk.

[0187] The weighted summation of the total increments in each direction is used to comprehensively consider the impact of different monitoring dimensions on dam safety. Different monitoring dimensions (vertical settlement, horizontal displacement, and tilt angle) have varying degrees of impact on dam safety. By assigning different weights to each dimension, the overall risk of the node can be assessed more accurately. Vertical settlement has the highest weight because an increase in vertical settlement may lead to foundation instability, significantly impacting dam safety. Horizontal displacement is next, as changes in horizontal displacement may affect the overall structural stability of the dam. Tilting angle has the lowest weight; although changes in tilt angle also affect dam safety, the impact is relatively minor. The total increment of overall risk is a value obtained by weighted summation of the total increments in each direction, reflecting the cumulative overall risk of the node across multiple monitoring dimensions.

[0188] Step S455: Calculate the node risk accumulation rate. The average risk growth rate per unit time is obtained by dividing the total risk increment by the number of prediction periods. The rate unit matches the deformation unit and the period unit.

[0189] Calculating the node risk accumulation rate is the final step in quantitatively assessing the overall risk growth of a node. The total increase in overall risk, calculated in previous steps, reflects the overall risk accumulation of a node across multiple monitoring dimensions.

[0190] The prediction period number is the number of sampling points within a predetermined future timeframe. By dividing the total increase in comprehensive risk by the prediction period number, the average risk growth rate per unit time can be obtained. The rate unit matches the deformation unit and the period unit.

[0191] Step S456: Mark the directional characteristics of the node risk accumulation rate. When the vertical settlement increment accounts for more than a preset proportion of the total increment, it is marked as vertically dominant. When the horizontal displacement increment accounts for more than a preset proportion, it is marked as horizontally dominant. The marking is used as an auxiliary reference for judging the warning level.

[0192] The percentage of vertical settlement increment refers to the proportion of the total vertical increment in the total increase of comprehensive risk, while the percentage of horizontal displacement increment refers to the proportion of the total horizontal increment in the total increase of comprehensive risk. The preset percentage is a pre-defined standard value. When the percentage of vertical settlement increment exceeds the preset percentage, it indicates that the risk increase at the node is mainly caused by vertical settlement, and it is marked as vertically dominant. When the percentage of horizontal displacement increment exceeds the preset percentage, it indicates that the risk increase at the node is mainly caused by horizontal displacement, and it is marked as horizontally dominant.

[0193] These markers can serve as an auxiliary reference for judging the warning level. Different dominant directions may correspond to different levels of risk and scope of impact. When judging the warning level, these markers can be combined to more accurately assess the risk status of nodes.

[0194] Step S460: Retrieve the design safety deformation threshold of the reservoir node, calculate the difference between the current deformation value (i.e., the end value of the analysis window) and the safety threshold as the critical threshold spacing. The smaller the spacing value, the closer the node is to the dangerous state. Arrange the risk accumulation rate of all nodes and the critical threshold spacing according to the path and node number to generate a two-dimensional risk evolution index matrix.

[0195] The current deformation value is the endpoint value of the analysis window, which is the node deformation value at the last time point within the analysis window captured in the previous steps. The critical threshold spacing is obtained by calculating the difference between the current deformation value and the safety threshold. This spacing value reflects the distance between the current deformation state of the reservoir node and its critical safety state. The smaller the spacing value, the closer the node is to a dangerous state, requiring more attention.

[0196] A two-dimensional risk evolution index matrix is ​​generated by arranging the risk accumulation rate and critical threshold spacing of all nodes according to their paths and node numbers. Path information indicates the specific route of anomalous state propagation, while node numbers specify the position of each node within that path. This arrangement clearly demonstrates the risk status of each reservoir / dam node under different propagation paths. Rows in the matrix can correspond to different nodes, and columns can represent the risk accumulation rate and critical threshold spacing, respectively. This matrix structure facilitates a systematic analysis and comparison of the risk evolution of the entire watershed's reservoir / dam group.

[0197] For example, firstly, the design safety deformation threshold for each reservoir / dam node is retrieved from the relevant design documents or database. Then, for each node, the endpoint value of its analysis window is obtained as the current deformation value, and the difference between this value and the safety threshold is calculated to obtain the critical threshold interval. Next, the risk accumulation rate and critical threshold interval for each node, calculated previously, are organized according to the path and node number. A database management system or data processing software can be used to sort the data and generate the matrix. For example, using the database's sorting function, the data can be sorted in ascending order according to the path and node number, and then the sorted data can be output in matrix form to form a two-dimensional risk evolution index matrix.

[0198] Step S500: Based on the risk evolution index matrix, generate a sequence of precursor risk identifiers for key reservoir and dam nodes.

[0199] The precursor risk identification sequence is an important set of information used for early warning of potential risks to reservoirs and dams. The risk evolution index matrix contains key information such as the risk accumulation rate and critical threshold spacing of each reservoir and dam node. Through further analysis and processing of this information, key reservoir and dam nodes can be identified, and corresponding precursor risk identifications can be generated for them.

[0200] In one implementation, step S500 may include the following steps S510-S560:

[0201] Step S510: Calculate the comprehensive risk index for each reservoir node in the risk evolution index matrix. The comprehensive risk index is obtained by dynamic weighted fusion of the node risk accumulation rate and the critical threshold distance. The rate weight and the distance weight are assigned based on the correlation of the impact factors in the basin's historical risk event database. The correlation is derived by inverting the node damage degree of the risk event and environmental parameters. The index range is transformed to a preset range through nonlinear mapping. The mapping process retains the monotonicity of the risk trend.

[0202] The node risk accumulation rate reflects the speed at which node risk increases, while the critical threshold distance reflects the distance between the node's current state and its dangerous state. By combining these two factors through dynamic weighted fusion, the degree of node risk can be measured more accurately.

[0203] The allocation of rate weights and spacing weights is based on the correlations of influencing factors in the historical risk event database of the watershed. This database records various past risk events, including the nodes where the events occurred, the extent of damage to those nodes, and the environmental parameters at the time. Through analysis and inverse derivation of this data, the influence of node risk accumulation rates and critical threshold spacing on risk events can be determined, thus identifying appropriate weights. For example, if the node risk accumulation rate plays a more critical role in the occurrence of a historical risk event, then the rate weight will be relatively larger.

[0204] For example, relevant data is first extracted from the historical risk event database of the watershed, and rate weights and spacing weights are determined through data analysis and inversion derivation. Then, for each reservoir-dam node in the risk evolution index matrix, the node risk accumulation rate and critical threshold spacing are weighted and fused according to the determined weights to calculate the comprehensive risk index. Finally, a nonlinear mapping method is used to transform the comprehensive risk index to a preset interval. For example, some common nonlinear functions, such as the sigmoid function, can be used to map the comprehensive risk index to the interval [0,1].

[0205] Step S520: The comprehensive risk index of all reservoir and dam nodes is screened using a density-based clustering algorithm. Risk feature clusters are divided according to the distribution density of index values ​​in the feature space. The average deviation between the index value of the center of each cluster and the index values ​​of all nodes in the cluster is calculated. Nodes with deviations exceeding the average deviation within the cluster are extracted as key reservoir and dam nodes. The screening process iteratively adjusts the cluster radius to make the risk characteristics of nodes within the cluster consistent. The initial value of the cluster radius is dynamically set based on the spatial distribution density of the nodes.

[0206] For example, the comprehensive risk index of all reservoir / dam nodes is first input into a density-based clustering algorithm, and the initial cluster radius is dynamically set according to the spatial distribution density of the nodes. Then, clustering is performed to obtain clusters with different risk characteristics. For each cluster, the average deviation between the cluster center index value and the index values ​​of the nodes within the cluster is calculated. Nodes with deviations exceeding the average deviation within the cluster are selected as key reservoir / dam nodes. During the selection process, the cluster radius is iteratively adjusted until the requirement of consistent risk characteristics among nodes within the cluster is met.

[0207] Step S530: Determine the warning level classification threshold through the multi-inflection point detection algorithm of the risk index distribution curve, dynamically divide the comprehensive risk index sequence into intervals, calculate the first and second derivatives of the curve after smoothing the index sequence, identify the inflection points of the curve by the sign change of the second derivative, divide multiple continuous intervals with the inflection points as the boundary, the interval with the highest index value corresponds to the highest warning level, and the interval with the lowest index value corresponds to the lowest warning level, and the number of intervals matches the number of inflection points.

[0208] The risk index distribution curve is formed by arranging the comprehensive risk indices of all reservoir and dam nodes in sequence. The multi-inflection point detection algorithm is a method used to identify multiple inflection points on this curve. An inflection point is a point where the slope of the curve changes, representing a turning point in the trend of the comprehensive risk index. Dynamic interval division of the comprehensive risk index sequence involves dividing the sequence into multiple continuous intervals based on the inflection points of the curve. Before division, the index sequence needs to be smoothed to eliminate noise and fluctuations in the data, making the curve smoother and facilitating accurate calculation of the first and second derivatives. The first derivative reflects the slope of the curve, and the second derivative reflects the rate of change of the slope. Inflection points can be identified by the change in the sign of the second derivative; when the sign of the second derivative changes from positive to negative or from negative to positive, the corresponding point is the inflection point of the curve.

[0209] Multiple consecutive intervals are divided by inflection points, each corresponding to a different warning level. The interval with the highest index value corresponds to the highest warning level, indicating that the risk situation of the nodes within that interval is the most severe; the interval with the lowest index value corresponds to the lowest warning level, indicating that the risk situation of the nodes within that interval is relatively mild. The number of intervals matches the number of inflection points; that is, the number of intervals corresponds to the number of inflection points.

[0210] Step S540: Calculate the early warning trigger time for each key reservoir / dam node. The trigger time is derived by the ratio of the critical threshold interval to the risk accumulation rate. During the derivation process, the risk accumulation rate is smoothed by a sliding window to eliminate the interference of short-term fluctuations on the trigger time. When the risk accumulation rate is zero, the trigger time is set to the maximum duration of the time window of the risk evolution index matrix, and the time unit is consistent with the monitoring data sampling period.

[0211] For example, for each critical reservoir node, its critical threshold spacing and risk accumulation rate are obtained. The risk accumulation rate is smoothed using a sliding window smoothing method. Then, the ratio of the critical threshold spacing to the smoothed risk accumulation rate is calculated to obtain the warning trigger time. If the risk accumulation rate is zero, the trigger time is set to the maximum duration of the time window of the risk evolution index matrix.

[0212] Step S550: Map the relative spatial coordinates of key reservoir nodes to the absolute coordinate system of the watershed through the spatial transformation rules of the topological relationship map. The coordinate transformation is based on the spatial adjacency matrix between nodes in the map. During the transformation process, the transformation error is corrected by verifying the coordinates of adjacent nodes. Spatiotemporal triggering conditions are generated in combination with the triggering time. The spatiotemporal triggering conditions include the correlation between the absolute coordinates and the triggering time.

[0213] During coordinate transformation, transformation errors are corrected by verifying the coordinates of adjacent nodes. Since there may be inherent errors during coordinate transformation, verifying the results using the coordinate information of adjacent nodes can correct the transformation and improve its accuracy. For example, if the transformed coordinates of a node do not conform to the spatial adjacency relationship in the topological graph, then the transformed coordinates of that node need to be adjusted.

[0214] For example, firstly, spatial transformation rules and the spatial adjacency matrix between nodes are obtained from the topological relationship graph. For each key reservoir / dam node, its relative spatial coordinates are transformed according to the spatial transformation rules to obtain its coordinates in the absolute coordinate system of the watershed. During the transformation process, the coordinates of adjacent nodes are used for verification and correction. Then, the transformed absolute coordinates are combined with the previously calculated warning trigger time to generate spatiotemporal triggering conditions.

[0215] Step S560: Arrange key reservoir nodes in descending order of comprehensive risk index, and associate them with their corresponding early warning level, unique node identifier, spatiotemporal triggering conditions and generation timestamp in sequence to form a precursor risk identifier sequence containing core risk characteristics. Each entry in the sequence contains four dimensions: early warning level, node code, spatiotemporal triggering conditions and generation timestamp. The dimensions are logically linked through the comprehensive risk index to ensure that the sequence can be directly used for early warning response scheduling.

[0216] For example, the key reservoir / dam nodes are first sorted in descending order based on their comprehensive risk index. Then, for each node, its corresponding warning level, unique node identifier, spatiotemporal triggering conditions, and the timestamp of its generation are associated to form an entry containing four dimensions of information. All the node entries are then combined in the ordered sequence to form a precursor risk identifier sequence. This sequence can be stored in a database for the early warning response scheduling system to query and use at any time.

[0217] Figure 2 This is a schematic diagram of a hardware entity of an early warning system provided in an embodiment of the present invention, such as... Figure 2 As shown, the hardware entity of the early warning system 1000 includes a processor 1001 and a memory 1002, wherein the memory 1002 stores a computer program that can run on the processor 1001, and the processor 1001 executes the program to implement the steps in the method of any of the above embodiments.

Claims

1. A BeiDou intelligent early warning analysis method applied to a basin-level reservoir and dam group, characterized in that, The method includes: The basin-level reservoir and dam group monitoring system acquires dam deformation time-series records and corresponding hydrological environment time-series records collected by the Beidou satellite positioning terminal. The dam deformation time-series records and the hydrological environment time-series records are associated and bound together by a unified timestamp to obtain an associated time-series dataset. The associated time-series dataset undergoes autoencoding processing of dam deformation response patterns. By capturing the coupling relationship between deformation and the hydrological environment, a set of deformation response feature primitives characterizing the dam's structural state is generated. Specifically, this includes: dividing the associated time-series dataset into time-continuous blocks according to a preset time interval. Each block contains a dam deformation record segment and a hydrological environment record segment within the corresponding window. The number of blocks is determined based on the total duration of the dataset and the window length. Deformation features and hydrological features are extracted from the block data. Deformation features include the rate of change of settlement and the amplitude of horizontal displacement fluctuations. Hydrological features include the amplitude of water level rise and fall and the gradient of flow velocity changes. The feature extraction process retains the timestamp correlation of the original data. The extracted deformation features and hydrological features are modeled to correlate with each other. The mutual information value of the features under different time lags is calculated. The lag with the largest mutual information value is selected as the optimal response delay, generating a spatiotemporal model containing response intensity and delay duration. The spatiotemporal correlation matrix is ​​used to perform multi-scale feature separation on the segmented data, decomposing it into feature components of different vibration modes. Short-term deformation features are extracted from the instantaneous response components, and long-term deformation features are extracted from the trend influence components. The feature dimensions are matched with the dimensions of the segmented data. The short-term deformation features and the long-term deformation features are input into an autoencoder processing unit. The joint features formed by combining the short-term deformation features and the long-term hydrological features are compressed into a fixed-dimensional latent space vector through a multi-layer nonlinear mapping network. The original features are then reconstructed through a decoding network. When the reconstruction error is less than a preset convergence threshold, the current latent space vector is saved as a candidate feature primitive. Density clustering is used to screen all candidate feature primitives. A clustering algorithm based on spatial distribution density is used to divide the feature clusters. The spatial distance between the center vectors of each cluster is calculated. Adjacent clusters with a distance less than the clustering threshold are merged. The center vectors of each cluster are extracted as the final deformation response feature primitives. All center vectors are combined to obtain the set of deformation response feature primitives. The set of deformation response feature primitives is input into a pre-constructed topological relationship map of the watershed reservoir-dam group. Anomaly propagation calculation is performed through the state dependency relationship between the nodes in the map to obtain a set of chain risk transmission paths that include risk propagation intensity and path probability. In the map, nodes represent reservoir-dam entities, edges represent hydraulic connection relationships between reservoirs and dams, node attributes include dam structural parameters and bound deformation response feature primitives, and edge attributes include hydraulic conduction coefficients and connection directions. Based on the set of chain risk transmission paths, risk evolution trend is extrapolated for reservoir nodes in each path, generating a risk evolution index matrix that includes node risk accumulation rate and critical threshold spacing. Based on the aforementioned risk evolution index matrix, a sequence of precursor risk identifiers is generated for key reservoir and dam nodes.

2. The method according to claim 1, characterized in that, The step of performing autoencoding processing on the associated time-series dataset to generate a set of deformation response feature primitives characterizing the dam's structural state by capturing the coupling relationship between deformation and the hydrological environment includes: The associated time series dataset is divided into segments at equal intervals along the time axis to obtain multiple continuous data segments. Each data segment contains complete dam deformation records and hydrological environment records within the corresponding time period. There is no time overlap between the data segments and the data segments cover the entire dataset. Multi-dimensional feature sequences were extracted from each data segment. Deformation dimension features included vertical settlement fluctuation sequences and horizontal displacement change sequences, while hydrological dimension features included water level fluctuation sequences and flow velocity gradient sequences. The sampling frequency of each sequence was consistent with the original record. Multi-dimensional feature sequences are dynamically fused, and fusion weights are assigned according to the mutual information values ​​between features. Features with higher mutual information values ​​are given higher weights, generating a fused feature sequence that includes deformation-hydrological coupling relationship. The sequence length is the same as the data segment length. The fused feature sequence is subjected to trend separation processing to separate short-term dynamic components and long-term trend components. The short-term dynamic components reflect instantaneous deformation response characteristics, while the long-term trend components reflect the cumulative impact characteristics of the hydrological environment. The separated feature components are input into the encoding processing module, and converted into fixed-dimensional feature vectors through a multilayer perceptron network. The vector dimension is determined according to the number of data segments and the feature complexity. The similarity between the feature vector and the preset standard feature template is calculated. When the similarity exceeds the matching threshold, it is determined to be a valid feature primitive. All valid feature primitives are filtered and arranged in chronological order to obtain the set of deformation response feature primitives. The size of the set is positively correlated with the number of data segments.

3. The method according to claim 2, characterized in that, The process extracts multi-dimensional feature sequences from each data segment. Deformation dimension features include vertical settlement fluctuation sequences and horizontal displacement change sequences; hydrological dimension features include water level fluctuation sequences and velocity gradient sequences. The sampling frequency of each sequence remains consistent with the original record, including: Vertical settlement data is extracted from the dam deformation record, arranged in chronological order to form the original settlement sequence, and the settlement difference between adjacent sampling points is calculated to obtain the vertical settlement fluctuation sequence. Horizontal displacement data are extracted from the dam deformation record, decomposed into lateral and longitudinal displacement components, and the fluctuation amplitude of each component, i.e. the difference between the maximum and minimum values, is calculated and arranged in chronological order to form a horizontal displacement change sequence. Water level monitoring data is extracted from hydrological environmental records, and the water level difference between adjacent time points is calculated. Positive values ​​indicate rising water levels, and negative values ​​indicate falling water levels. The data are arranged in chronological order to form a water level fluctuation sequence. Extract flow velocity monitoring data from hydrological environmental records, calculate the standard deviation of flow velocity within the window, the standard deviation sequence reflects the intensity of flow velocity fluctuations, arrange them in chronological order to form a flow velocity gradient sequence, and dynamically adjust the window size according to the frequency of flow velocity changes. The vertical settlement fluctuation sequence, horizontal displacement change sequence, water level fluctuation sequence and flow velocity gradient sequence are time-stamp aligned, and the sampling points of each sequence are matched one by one based on the unified timestamp to generate a multi-dimensional feature matrix. The multi-dimensional feature matrix is ​​standardized to map the values ​​of each feature sequence to a preset interval.

4. The method according to claim 1, characterized in that, The process involves modeling the correlation between extracted deformation features and hydrological features, calculating the feature mutual information value under different time lags, selecting the lag with the largest mutual information value as the optimal response delay, and generating a spatiotemporal correlation matrix containing response intensity and delay duration, including: The extracted deformation and hydrological features are tested for stationarity, and non-stationary sequences are differentially processed until the stationarity requirements are met. Obtain the preset time lag search range. The lag value ranges from the initial value to the maximum lag value. The maximum lag value is determined based on the sampling frequency of the feature sequence and the hydrological response cycle. The lag interval is consistent with the sampling cycle. The mutual information value between deformation features and hydrological features is calculated for each hysteresis. The mutual information value reflects the correlation strength between features. The larger the value, the stronger the correlation. The estimation error of the feature probability density function is taken into account in the calculation process. Record the mutual information value corresponding to each lag, plot the mutual information-lag relationship curve, the lag corresponding to the peak of the curve is the optimal response delay, the peak mutual information value is the response intensity, and each feature pair corresponds to a set of optimal delay and intensity. Based on the optimal response delay and response intensity of all feature pairs, a spatiotemporal correlation matrix is ​​constructed. The matrix row index represents the deformation feature type, the column index represents the hydrological feature type, the matrix element represents the response intensity of the corresponding feature pair, and the element position corresponds to the optimal delay duration. The spatiotemporal correlation matrix is ​​sparsified by retaining matrix elements with response strength greater than a preset threshold and setting the remaining elements to zero to generate a sparse spatiotemporal correlation matrix.

5. The method according to claim 1, characterized in that, The process involves inputting the set of deformation response feature primitives into a pre-constructed topological map of the watershed reservoir-dam group, and performing anomaly propagation calculations based on the state dependencies between map nodes to obtain a set of cascading risk transmission paths that includes risk propagation intensity and path probability, including: Call the pre-built topological relationship map of the watershed reservoir-dam group; The set of deformation response feature primitives is assigned to the corresponding reservoir nodes by spatial coordinate matching. The spatial distance between the grid coordinates of the feature primitives and the geographical coordinates of the reservoir nodes is calculated. The node with the smallest distance is the matching node. A primitive-node association table is generated, which contains the primitive index and the corresponding node identifier. Initialize the state parameters of the map nodes. The state parameters include deformation anomaly degree, hydrological sensitivity, and propagation probability. The deformation anomaly degree is calculated based on the deviation value between the feature primitive and the standard primitive. The propagation probability is initialized to a value that is positively correlated with the deformation anomaly degree. Identify initial abnormal nodes. When the abnormality of a node deformation exceeds a preset safety threshold, mark it as an initial abnormal node, set its propagation probability as the initial propagation probability, initialize the propagation probability of other nodes to zero, and record the list of initial abnormal nodes and their corresponding state parameters. Starting from the initial abnormal node, multiple rounds of abnormal propagation simulation are performed. In each round, the node attempts to activate the adjacent node according to the propagation probability. After successful activation, the state parameters of the adjacent node are updated. The update magnitude is positively correlated with the edge propagation coefficient. The node sequence and activation time of each propagation path are recorded. Calculate the cumulative propagation intensity and path probability for each propagation path. The propagation intensity is the product of the transmission coefficients of each side of the path, and the path probability is the product of the activation probabilities of each node. Paths with propagation intensity greater than the intensity threshold and path probability greater than the probability threshold are selected and sorted in ascending order by the number of propagation steps to obtain the set of chain risk transmission paths.

6. The method according to claim 1, characterized in that, The process involves inputting the set of deformation response feature primitives into a pre-constructed topological map of the watershed reservoir-dam group, and performing anomaly propagation calculations based on the state dependencies between map nodes to obtain a set of cascading risk transmission paths that includes risk propagation intensity and path probability, including: Load the topological relationship map data of the watershed reservoir-dam group. The map is stored in a graph structure, including node sets and edge sets. The node sets store the spatial location and deformation feature primitives of the reservoir-dam group, and the edge sets store the hydraulic connection relationship and conduction parameters between the nodes. Based on the spatial coordinates of the reservoir nodes, the feature primitives in the set of deformation response feature primitives are bound to the corresponding nodes. The matching process calculates the Euclidean distance between the feature primitives and the nodes. Primitives with a distance less than the matching threshold are assigned to the node, and a primitive allocation table is generated. Obtain the preset initial value of the node state vector. The state vector includes three dimensions: deformation deviation, environmental correlation, and propagation capability. The deformation deviation is calculated based on the similarity between the feature primitive and the standard primitive. The propagation capability is positively correlated with the deformation deviation. Determine the initial source of anomaly propagation. When the deformation deviation of a node exceeds the design safety threshold, mark it as an anomaly propagation source and set its propagation capability value as the baseline propagation probability. Initialize the propagation capability values ​​of other nodes to zero and record the node identifier and state vector of the anomaly propagation source. Anomaly propagation simulation is performed based on an independent cascade model. The propagation source node attempts to activate adjacent nodes according to the propagation capability value. After successful activation, the state vector of the adjacent nodes is updated. During the propagation simulation, the node sequence, activation time, and state vector changes of each propagation path are recorded. The cumulative propagation intensity and path probability of the path are calculated. Paths that meet the threshold conditions are selected and sorted by the number of propagation steps to obtain the set of chain risk transmission paths.

7. The method according to claim 6, characterized in that, The simulation of abnormal state propagation based on the independent cascade model involves the source node attempting to activate adjacent nodes according to its propagation capability value. Upon successful activation, the state vector of the adjacent nodes is updated. The deformation deviation increases the deviation contribution of the propagation source, and the environmental correlation is adjusted according to the edge attributes. The propagation capability value is set as the current activation probability, including: Obtain preset propagation simulation parameters, including maximum propagation rounds, propagation time step, and probability decay factor. The maximum propagation rounds are dynamically determined based on the number of reservoirs and dams in the basin. The time step is consistent with the sampling period of BeiDou monitoring data. The decay factor decreases as the propagation rounds increase. During the first round of propagation, the initial abnormal propagation source node is added to the set of active nodes, its state vector is recorded, the deformation deviation is the deviation between the current monitoring value and the safe value, the environmental correlation is the hydrological coupling coefficient in the primitive, and the propagation capability value is the initial baseline probability. During each round of propagation, the nodes in the active node set are traversed, and all outgoing neighbor nodes, i.e. downstream nodes along the water flow direction, are extracted. The activation probability of the current node for each neighbor node is calculated. For each neighboring node, a random value is generated. When the random value is less than the activation probability, the neighboring node is successfully activated and added to the new active node set. Its state vector is updated. The deformation deviation is the weighted sum of the deformation deviation of the current node and its original deviation. The environmental correlation is the hydrological sensitivity in the edge attribute. The propagation ability value is the current activation probability. The propagation rounds increase, the set of active nodes is updated to a new set of active nodes, and the propagation process is repeated until the maximum propagation rounds are reached or the new set of active nodes is empty. The node activation status and state vector change trajectory of each round of propagation are recorded. After the propagation ends, all propagation paths are organized. Each propagation path includes the sequence of activated nodes, the activation time of each node, and state parameters. The cumulative propagation strength and path probability of the propagation path are calculated.

8. The method according to claim 5, characterized in that, The calculation of the cumulative propagation strength and path probability for each propagation path, where propagation strength is the product of the transmission coefficients of each edge of the path and path probability is the product of the activation probabilities of each node, filters out paths with propagation strength greater than a strength threshold and path probability greater than a probability threshold, and sorts them in ascending order of propagation steps to obtain the set of chain risk transmission paths, including: Traverse all propagation paths, extract the edge set and node set contained in each propagation path. The edge set contains the hydraulic conduction coefficient of each segment of the path, and the node set contains the activation probability of each node. The path identifier is associated with the propagation round. The cumulative propagation intensity of the path is calculated by multiplying the hydraulic conductivity coefficients of each side from the starting node to the ending node. The product is the cumulative propagation intensity of the path, which reflects the risk transmission capacity of the path. To calculate the path probability, start from the starting node and end node of the path, multiply the activation probabilities of each node in turn. The product is the path probability of the path, which reflects the likelihood of the path occurring. Obtain preset intensity screening thresholds and probability screening thresholds. The intensity screening thresholds are determined based on the statistical analysis of the transmission intensity of historical risk events in the watershed, while the probability screening thresholds are set based on path reliability requirements. The propagation paths are double-screened, retaining paths whose cumulative propagation intensity is greater than the intensity screening threshold and whose path probability is greater than the probability screening threshold, and deleting paths that do not meet either threshold condition. The filtered paths are sorted in ascending order by the number of propagation steps (i.e., the number of nodes contained in the path minus one). Paths with the same number of propagation steps are sorted in descending order by cumulative propagation intensity. The sorted set of paths is extracted as the set of chain risk transmission paths.

9. An early warning system, comprising a memory and a processor, wherein the memory stores a computer program executable on the processor, characterized in that, When the processor executes the program, it implements the steps of the method according to any one of claims 1 to 8.

Images