Method and system for monitoring internal stresses of a dam

Abstract

The present application belongs to the technical field of stress monitoring, and particularly relates to a dam internal stress monitoring method and system, comprising: acquiring full-matrix capture data, real-time temperature data and real-time stress data; constructing an initial non-uniform anisotropic sound velocity field according to a finite element model and an acoustic-elastic constitutive relation; performing iterative optimization, calculating group velocity propagation path and transmission-reception delay time by ray tracing, and obtaining a sound pressure amplitude graph by coherent superposition; extracting an instantaneous phase to generate a normalized coherence coefficient graph, updating the sound velocity field based on a coherence coefficient maximization rule, and generating an ultrasonic image after convergence; extracting a gray level co-occurrence matrix feature to construct a high-dimensional texture feature vector, dynamically updating an abnormality judgment threshold, and outputting a warning through an isolation forest abnormality detection model. The present application realizes adaptive follow-up adjustment of the abnormality judgment threshold, eliminates the interference of environmental factors on the warning accuracy, and enhances the warning reliability under complex working conditions.

Description

Technical Field

[0001] This invention relates to the field of stress monitoring technology. More specifically, this invention relates to a method and system for monitoring stress inside dams. Background Technology

[0002] Reservoir dams are crucial water conservancy infrastructures, and their safe operation directly impacts the safety of life and property and public safety in downstream areas. With increasing service life, dam structures are subjected to the combined effects of water erosion, reservoir water pressure, geological settlement, temperature changes, and dynamic and static loads, which can gradually lead to hidden damage such as microcracks, voids, and cavitation within the dam body. Such damage is difficult to detect through visual inspections in its early stages, and if not identified and warned of in a timely manner, it may further develop into serious safety accidents such as leakage, local instability, or even dam failure. In recent years, IoT monitoring technology and non-destructive testing technology have been widely applied to reservoir dam safety monitoring. IoT systems can continuously collect multi-source status data such as dam body temperature, stress, and strain, while ultrasonic testing, due to its strong penetrating power and high sensitivity to internal defects, is used to assess the structural integrity of the dam body. By combining multi-source sensor data with ultrasonic testing data, a data foundation can be provided for dam operation status analysis and early warning.

[0003] However, existing technologies still have shortcomings in high-precision imaging and intelligent alarm systems. As large-volume load-bearing structures, dams experience continuous changes in their internal temperature field, stress field, and material state. The propagation speed and path of sound waves within the dam body vary with spatial location and direction, exhibiting significant non-uniformity and anisotropy. Traditional ultrasonic imaging methods typically assume uniform and isotropic sound velocity in the medium, making it difficult to fully consider the impact of actual stress distribution on sound wave propagation paths, speeds, and transmission / reception delays. This can easily lead to image defocusing, increased artifacts, and inaccurate defect location. Furthermore, conventional anomaly detection methods often rely on manually set static thresholds for alarm judgment, which are ill-suited to normal data shifts caused by seasonal temperature differences, load variations, and stress fluctuations. In complex service environments, fixed thresholds can easily lead to false alarms or missed alarms, and are also insufficient for effectively extracting structural degradation information hidden in high-dimensional ultrasonic texture features. Therefore, how to combine real-time environmental and stress state adjustments to correct ultrasonic imaging models and establish an anomaly early warning mechanism that adapts to service conditions has become a pressing technical problem to be solved in the field of reservoir dam safety monitoring. Summary of the Invention

[0004] To address the technical problems in existing dam ultrasonic imaging technologies, such as insufficient consideration of non-uniform and anisotropic changes in sound velocity caused by temperature and stress, leading to image defocusing, defect location errors, and false alarms or missed alarms due to difficulty in adapting alarm thresholds to changes in service conditions, this invention provides a method and system for monitoring internal stress of dams.

[0005] In a first aspect, the present invention provides a method for monitoring stress inside a dam, comprising: S1: acquiring full-matrix capture data of the dam surface, collecting real-time temperature data and real-time stress data from embedded sensors within the dam via the Internet of Things, and constructing an initial non-uniform anisotropic sound velocity field based on the real-time temperature data, real-time stress data, finite element model, and acoustoelastic constitutive relation; S2: performing iterative optimization processing, using ray tracing to calculate the group velocity propagation path and transmission / reception delay time from each array element to the pixel, coherently superimposing the full-matrix capture data based on the transmission / reception delay time to obtain a sound pressure amplitude map, and extracting the instantaneous phase to calculate the inverse of the dispersion to generate a normalized coherence coefficient map. If the change in the coherence coefficient map is not lower than the convergence threshold, the direction-correlated sound velocity perturbation is determined based on the coherence coefficient maximization rule to update the non-uniform anisotropic sound velocity field and the iteration is repeated. If the change in the coherence coefficient map is lower than the convergence threshold, the sound pressure amplitude map and the coherence coefficient map are fused to generate an ultrasonic image. S3: Wavelet transform is performed on the ultrasonic image to extract gray-level co-occurrence matrix features and concatenate them into a high-dimensional texture feature vector. The mapping relationship between real-time temperature data and real-time stress data and the upper limit of the benchmark anomaly score is established and the anomaly judgment threshold is dynamically updated. The high-dimensional texture feature vector is input into the isolated forest anomaly detection model. When the anomaly score exceeds the updated anomaly judgment threshold, the dam early warning information is output.

[0006] By employing the above technical solutions, this invention achieves deep alignment between the ultrasonic imaging model and the actual physical state of the dam by constructing an initial non-uniform anisotropic sound velocity field based on real-time temperature and stress correction and performing iterative optimization. This effectively eliminates image artifacts caused by sound wave propagation path deviations and significantly improves the positioning accuracy of concealed damage. Simultaneously, by extracting high-dimensional texture features and establishing an anomaly scoring mapping mechanism based on environmental parameters, adaptive adjustment of the anomaly judgment threshold is achieved. This eliminates the interference of seasonal temperature differences and load fluctuations on the accuracy of early warning, enhancing the reliability of the system's early warning under complex operating conditions.

[0007] Preferably, the step of constructing an initial non-uniform anisotropic sound velocity field based on real-time temperature data, real-time stress data, a finite element model, and an acoustoelastic constitutive relation includes: dividing the finite element model into a three-dimensional mesh; calculating the stress tensor of each mesh node based on the real-time stress data; determining the temperature distribution and temperature correction of each mesh node based on the real-time temperature data; obtaining the inherent acoustoelastic coefficient tensor of the dam material; combining the initial reference sound velocity under zero stress and reference temperature conditions, and based on the acoustoelastic constitutive relation, calculating the sound velocity offset of each mesh node under the three-dimensional propagation direction parameters characterized by the propagation azimuth angle and propagation pitch angle through the inherent acoustoelastic coefficient tensor, the stress tensor of each mesh node, and the temperature correction of each mesh node; and superimposing these to obtain the initial non-uniform anisotropic sound velocity field.

[0008] By adopting the above technical solution, this invention, starting from the laws of physical evolution, determines the coupled influence of stress tensor and temperature correction on sound propagation characteristics, providing an initial model that closely approximates the actual physical field for subsequent iterative optimization, which can accelerate the convergence speed of sound velocity field reconstruction and improve the accuracy of calculation results.

[0009] Preferably, the step of calculating the group velocity propagation path and transmit / receive delay time from each array element to the pixel using the ray tracing method includes: discretizing the imaging space into a grid, obtaining the coordinates of each array element and the coordinates of the pixel, using each array element as the starting point, solving the equation of motion in the current non-uniform anisotropic sound velocity field using the fast travel method to obtain the wavefront travel time field; combining the deviation relationship between phase velocity and group velocity, determining the group velocity direction based on the three-dimensional propagation direction parameters and the current non-uniform anisotropic sound velocity field, and tracing backward from the pixel to the corresponding array element along the group velocity direction to obtain the group velocity propagation path; performing a line integral on the slowness along the group velocity propagation path, calculating the propagation delay time when each array element is used as a transmitter and receiver respectively, and adding them together to obtain the transmit / receive delay time.

[0010] Preferably, the step of using the fast-moving method to solve the anisotropic equation in the current non-uniform anisotropic sound velocity field to obtain the wavefront travel time field, and tracing backward from the pixel to the corresponding array element along the group velocity direction, includes: using the fast-moving method to solve the anisotropic equation in the current non-uniform anisotropic sound velocity field, generating the arrival time of each pixel and constructing the wavefront travel time field matrix; using the fourth-order Runge-Kutta method to perform step-by-step backward tracing from the target pixel along the opposite direction of the group velocity vector until returning to the position of the corresponding array element, generating a discretized set of group velocity propagation path points.

[0011] Preferably, the step of extracting the instantaneous phase, calculating the reciprocal of the discreteness, and generating a normalized coherence coefficient map includes: selecting a spatial neighborhood centered on the pixel, obtaining the delayed received signals corresponding to all array channels within the spatial neighborhood; performing a Hilbert transform on the delayed received signals to construct an analytical signal, and extracting the instantaneous phase angle of each array channel; calculating the circular square difference of the instantaneous phase angles of all array channels as the phase discreteness, taking the reciprocal of the sum of the phase discreteness and a preset minimum constant, and mapping it to the interval between 0 and 1 through extreme value normalization to generate the normalized coherence coefficient map.

[0012] Preferably, the step of determining the direction-dependent sound velocity perturbation based on the coherence coefficient maximization rule to update the non-uniform anisotropic sound velocity field includes: using a functional chain calculus derivation from the coherence coefficient to the instantaneous phase, the transmit / receive delay time, and even the global ray path, calculating the functional gradient of the coherence coefficient with respect to the direction-dependent sound velocity parameters at each spatial grid node; normalizing the functional gradient, and multiplying the normalized update direction by the update step size to perform gradient ascent, obtaining the sound velocity perturbation at the corresponding spatial grid node and the corresponding three-dimensional propagation direction parameters; and superimposing the sound velocity perturbation onto the corresponding direction-dependent sound velocity parameters in the current non-uniform anisotropic sound velocity field to complete the model update.

[0013] By adopting the above technical solution, the present invention can directly convert the quality of imaging into correction instructions for the sound velocity parameters of each grid point in space, enabling the sound velocity field model to spontaneously compensate for the local inhomogeneity of the material based on the measured data, and realizing closed-loop optimization of ultrasonic detection data and physical field model.

[0014] Preferably, the step of extracting gray-level co-occurrence matrix features from the ultrasonic image using wavelet transform and concatenating them into a high-dimensional texture feature vector includes: decomposing the ultrasonic image using multi-scale redundant wavelet transform to obtain low-frequency sub-band images and high-frequency sub-band images at different scales; calculating the gray-level co-occurrence matrix for each low-frequency sub-band image and high-frequency sub-band image in four directions: 0 degrees, 45 degrees, 90 degrees, and 135 degrees; extracting energy features, contrast features, and correlation features from each gray-level co-occurrence matrix; and concatenating the features extracted at all scales and in all directions in a preset order to form the high-dimensional texture feature vector.

[0015] Preferably, the step of establishing the mapping relationship between real-time temperature data and real-time stress data to the upper limit of the benchmark anomaly score and dynamically updating the anomaly judgment threshold includes: collecting multiple sets of temperature data, multiple sets of stress data, and the corresponding benchmark anomaly scores output by the isolated forest anomaly detection model under historical normal conditions; determining the upper limit of the benchmark anomaly score based on the maximum value of the benchmark anomaly score within a preset time window; training a regression prediction model using a support vector regression algorithm, with the multiple sets of temperature data and multiple sets of stress data as input and the upper limit of the benchmark anomaly score as output; inputting the real-time collected temperature data and stress data into the regression prediction model, outputting the upper limit of the benchmark anomaly score under the current environmental condition, and adding a preset fault tolerance margin as the updated anomaly judgment threshold.

[0016] By adopting the above technical solution, the present invention effectively decouples signal fluctuations under normal service conditions from structural anomaly signals. By dynamically outputting the upper limit of the benchmark anomaly score under the current temperature and pressure environment and adding a fault tolerance margin, false alarms caused by changes in the normal environment are eliminated, ensuring that the early warning system can maintain stable detection sensitivity under different seasons and load levels.

[0017] Preferably, after the iterative optimization process meets the convergence condition, the method further includes a full-field stress field inversion step inside the dam: based on the temperature distribution data of each spatial grid node, the sound velocity offset caused by temperature is subtracted from the final converged sound velocity field, and the direction-related sound velocity deviation caused by the stress change inside the dam is separated; the direction-related sound velocity deviation is substituted into the acoustoelastic constitutive relation, and the actual principal stress magnitude and direction at each spatial grid node are calculated by inverse solution using the nonlinear least squares method, and the actual principal stress of each grid node is mapped to three-dimensional space to generate a full-field stress distribution cloud map inside the dam.

[0018] By adopting the above technical solution, the present invention realizes the extension from qualitative damage imaging to quantitative stress inversion. The generated three-dimensional full-field stress distribution cloud map can intuitively reveal the stress concentration area and structural safety evolution trend inside the dam, providing data support for the overall stability assessment of the dam.

[0019] Secondly, the present invention provides a dam internal stress monitoring system, including a processor and a memory, wherein the memory stores computer program instructions, and when the computer program instructions are executed by the processor, the above-mentioned dam internal stress monitoring method is implemented.

[0020] By adopting the above technical solution, the above-mentioned method for monitoring internal stress of dams is generated into a computer program and stored in a memory so that it can be loaded and executed by a processor. Based on the memory and processor, a terminal device can be made for convenient use.

[0021] The technical solution of the present invention has the following beneficial technical effects: This invention combines real-time temperature and stress data acquired via the Internet of Things (IoT) and utilizes a finite element model and acoustoelastic constitutive relations to construct an initial non-uniform anisotropic sound velocity field, enabling the ultrasonic propagation model to closely match the actual service conditions of the dam. By employing ray tracing to calculate the group velocity propagation path and transmission / reception delay time from each array element to the pixel, and iteratively optimizing and correcting the non-uniform anisotropic sound velocity field based on the coherence coefficient maximization rule, the generated sound pressure amplitude map and the normalized coherence coefficient map are fused to generate the final ultrasonic image. This effectively reduces imaging defocusing and artifacts caused by temperature changes and stress field distribution fluctuations.

[0022] Furthermore, this invention extracts high-dimensional texture feature vectors through multi-scale redundant wavelet transform and gray-level co-occurrence matrix, and establishes a mapping relationship between real-time temperature data and real-time stress data and the upper limit of the benchmark anomaly score to achieve dynamic updating of the anomaly judgment threshold. Then, it combines the isolated forest anomaly detection model to perform anomaly judgment, reducing the negative impact of environmental factor fluctuations on alarm results and improving the accuracy and stability of reservoir safety monitoring and early warning. Attached Figure Description

[0023] Figure 1 This is a flowchart of the method for monitoring internal stress of a dam in this invention; Figure 2 This is a schematic diagram of the mean curve of the global coherence coefficient during optimization iteration; Figure 3 It is a bar chart showing the numerical distribution of contrast features under different wavelet decomposition scales and different gray-level co-occurrence matrix directions; Figure 4 This is a schematic diagram showing the comparison curves between the dynamically updated anomaly detection threshold and the static baseline threshold within a 24-hour monitoring period. Detailed Implementation

[0024] The technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments of the present invention. Obviously, the described embodiments are some embodiments of the present invention, but not all embodiments.

[0025] This invention discloses a method for monitoring internal stress in dams, referring to... Figure 1 This includes steps S1-S3: S1: Construct an initial non-uniform sound velocity field driven by multi-source data.

[0026] Full matrix capture data of the dam surface is obtained, and real-time temperature and stress data of embedded sensors in the dam body are collected through the Internet of Things. Based on the real-time temperature data, real-time stress data, finite element model and acoustoelastic constitutive relation, an initial non-uniform anisotropic sound velocity field is constructed.

[0027] In the specific implementation process, the transmission and reception control board of the ultrasonic phased array instrument controls the operation of each array element, enabling each element to transmit ultrasonic signals individually in a preset order. During each transmission, all array elements synchronously receive the echo signals, thereby generating full-matrix acquisition data. The acquired full-matrix acquisition data is organized according to three dimensions: transmitting array element number, receiving array element number, and time sampling point, forming a three-dimensional data array. Simultaneously, real-time temperature data is collected using fiber optic temperature sensors deployed inside the dam body, and real-time stress data is collected using vibrating wire stress gauges. The acquired real-time temperature and stress data are then uploaded to a cloud server for centralized processing via a message queue telemetry transmission protocol.

[0028] After receiving the aforementioned sensor data, the cloud server calls the application programming interface (API) of the Abaqus finite element analysis software to read the dam's three-dimensional structural mesh and imports the received real-time stress data into the dam's three-dimensional structural mesh for interpolation calculation of the stress tensor at the mesh nodes. Simultaneously, real-time temperature data is mapped to the corresponding mesh nodes to obtain the temperature distribution data for each mesh node. Based on this, the system establishes the relationship equation between principal stress and sound velocity changes based on acoustoelastic theory and solves for the phase velocity and group velocity in the three-dimensional anisotropic medium using the Christoffel equation. According to the difference between the temperature distribution data of each mesh node and the preset reference temperature, the temperature-velocity correction is calculated using a pre-calibrated temperature-velocity correction coefficient. Subsequently, a meshed initial non-uniform anisotropic sound velocity field model is constructed using a spatial interpolation algorithm based on non-uniform rational B-splines. This allows the sound velocity at any three-dimensional coordinate point in the space to be expressed as a spherical harmonic function expansion with the spatial propagation pitch angle and azimuth angle as independent variables. This spherical harmonic function is parameterized and stored for subsequent use.

[0029] In an optional embodiment, an initial non-uniform anisotropic sound velocity field is constructed based on real-time temperature data, real-time stress data, a finite element model, and an acoustoelastic constitutive relation. Specifically, the following processing logic is included: First, the finite element model is divided into a three-dimensional mesh. The stress tensor of each mesh node is calculated based on the real-time stress data, and the temperature distribution and temperature correction of each mesh node are determined based on the real-time temperature data. Second, the inherent acoustoelastic coefficient tensor of the dam material is obtained. Finally, combining the initial reference sound velocity under zero stress and reference temperature conditions, and based on the acoustoelastic constitutive relation, the sound velocity offset of each mesh node under the three-dimensional propagation direction parameters characterized by the propagation azimuth and pitch angles is calculated using the inherent acoustoelastic coefficient tensor, the stress tensor of each mesh node, and the temperature correction of each mesh node. The various sound velocity offsets are then superimposed to obtain the initial non-uniform anisotropic sound velocity field.

[0030] In the specific operation of dividing the finite element model into a three-dimensional mesh, eight-node hexahedral elements are used to discretize the dam finite element model. To balance computational accuracy and processing efficiency, the mesh side length is set within the range of 0.1m to 0.5m, with 0.2m being preferred in practice. The selection is based on the following: the dam is a large-volume concrete structure of hundreds of meters in diameter, and the gradient of its internal macroscopic temperature and stress fields is relatively gentle in space. Furthermore, the actual physical spacing of the sensors embedded within the dam is typically on the order of tens of centimeters to several meters. If the mesh side length is greater than 0.5m, the mesh nodes are too sparse, making accurate interpolation and capturing subtle technical state fluctuations between sensors impossible, thus reducing the accuracy of stress tensor representation. If the mesh side length is less than 0.1m, the number of global mesh nodes increases exponentially, easily leading to convergence difficulties in the global finite element matrix of the dam during interpolation calculations or causing computational overload. Therefore, limiting the grid side length to 0.1m to 0.5m (e.g., 0.2m) enables the grid density spatial resolution to be highly matched with the spatial deployment density of embedded sensors and the gradient of macroscopic stress field, achieving a closed-loop balance between computational accuracy and convergence efficiency in practical engineering applications.

[0031] Based on real-time stress data collected by sensors embedded within the dam, and through finite element mesh interpolation, stress tensor data at each mesh node is obtained. This stress data measurement range covers -5MPa to 5MPa with an accuracy of 0.01MPa. Then, the six independent stress tensor components of each node are calculated using a principal stress interpolation algorithm. When obtaining the inherent acoustoelastic coefficient tensor of the dam concrete material, taking a dam material with a density of 2400kg / m³ and a Poisson's ratio of 0.2 as an example, the corresponding second- and third-order acoustoelastic coefficient tensor components are obtained through pre-conducted indoor triaxial compression ultrasonic testing or material acoustoelastic calibration tests. The reference P-wave velocity of the dam under zero stress and reference temperature conditions is extracted, with a reference P-wave velocity of 4000m / s as an example. The aforementioned inherent acoustoelastic coefficient tensor and reference longitudinal wave velocity were obtained by conducting ultrasonic propagation velocity tests on concrete samples with the same mix proportion, age, and similar water content as the dam under different stress states, and by accurately fitting and calibrating the stress change, sound velocity change, and zero-stress sound velocity obtained from the tests.

[0032] Based on the acoustoelastic constitutive relation, the obtained intrinsic acoustoelastic coefficient tensor and the real-time principal stress tensor of each grid node are substituted into the Christoffel equation for solution. Simultaneously, based on the actual difference between the real-time temperature data of each grid node and the reference temperature, the corresponding temperature-velocity correction is calculated using a pre-calibrated temperature-velocity correction coefficient. This temperature-velocity correction coefficient is obtained by conducting ultrasonic propagation velocity tests on the same dam material sample under different constant temperature conditions, and then performing linear regression or piecewise fitting calibration on the observed temperature change and sound velocity change.

[0033] When determining the three-dimensional propagation direction, the azimuth and elevation angles of the ultrasonic wave are set as parameters for the three-dimensional propagation direction. The azimuth and elevation angles are discretized and combined according to a preset angular step size. Based on this, the phase velocity change caused by stress variation under each three-dimensional propagation direction parameter is calculated. This phase velocity change, along with the aforementioned temperature-based sound velocity correction, is used as the sound velocity offset. Taking a grid node subjected to a compressive stress of 2.5 MPa and the sound wave propagation direction parallel to the principal compressive stress direction as an example, the system calculates the corresponding stress-based sound velocity offset and combines it with the temperature-based sound velocity correction for that specific grid node. Both are then superimposed onto the initial reference sound velocity under the reference state. By linearly superimposing the initial reference sound velocity under zero stress and reference temperature conditions with the comprehensive sound velocity offset calculated for the target grid node under each three-dimensional propagation direction parameter, the system ultimately generates a non-uniform anisotropic sound velocity field parameter array containing three-dimensional spatial coordinates and three-dimensional propagation direction parameter indices.

[0034] S2: Iterative optimization and imaging based on maximizing coherence coefficient.

[0035] Iterative optimization is performed, and the group velocity propagation path and transmit / receive delay time from each array element to the pixel are calculated using the ray tracing method. Based on the transmit / receive delay time, the full matrix captured data is coherently superimposed to obtain the sound pressure amplitude map. The instantaneous phase is extracted and the inverse of the dispersion is calculated to generate a normalized coherence coefficient map. If the change in the coherence coefficient map is not lower than the convergence threshold, the direction-dependent sound velocity perturbation is determined based on the coherence coefficient maximization rule to update the non-uniform anisotropic sound velocity field and the iteration is repeated. If the change in the coherence coefficient map is lower than the convergence threshold, the sound pressure amplitude map and the coherence coefficient map are fused to generate an ultrasonic image.

[0036] In the specific iterative optimization process, the system first establishes a two-dimensional imaging pixel grid in the initial non-uniform anisotropic sound velocity field. Based on this grid, the ray tracing method is used to calculate the group velocity propagation path of the ultrasonic wave from each transmitting element to the target pixel and back to the corresponding receiving element, as well as the corresponding two-way transmission and reception delay time. After obtaining the transmission and reception delay time, the system extracts the echo data at the corresponding time from the full matrix acquisition data, and uses a full focusing algorithm to perform delay compensation and coherent superposition on the data of each transmission and reception channel to obtain the sound pressure amplitude of each pixel, thereby generating a sound pressure amplitude map reflecting the internal structure of the dam.

[0037] Simultaneously, signal analysis is performed on the full matrix captured data to extract the instantaneous phase. For each pixel, the dispersion of the instantaneous phase of all transmit and receive channels is calculated, and the dispersion is converted into a normalized coherence coefficient map. Subsequently, the difference between the coherence coefficient maps of two adjacent iterations is calculated using the normalized norm to evaluate the convergence state.

[0038] If the convergence condition is not met, the direction-dependent sound velocity perturbation is calculated based on the coherence coefficient maximization rule, superimposed on the current sound velocity field to complete the update and trigger the next round of calculation; if the convergence condition is met, the iteration stops, and the converged sound pressure amplitude map is multiplied pixel by pixel with the coherence coefficient map to generate an ultrasonic image used to enhance the display of the boundary of the internal structure of the dam.

[0039] In an optional embodiment, the group velocity propagation path and transmit / receive delay time from each array element to the pixel are calculated using ray tracing. Specifically, the following steps are included: First, the imaging space is discretized into a grid to obtain the coordinates of each array element and the pixel. Starting from each array element, the wavefront travel time field is solved using the equation of process function in the current non-uniform anisotropic sound velocity field. Then, combining the deviation relationship between phase velocity and group velocity, the group velocity direction is determined based on the three-dimensional propagation direction parameters and the current non-uniform anisotropic sound velocity field. The propagation path is then traced backward from the pixel to the corresponding array element along the group velocity direction. Finally, the slowness along the group velocity propagation path is linearly integrated to calculate the propagation delay time for each array element as both a transmitter and receiver, and these delay times are summed to obtain the transmit / receive delay time.

[0040] In the specific operation of the above-mentioned ray tracing, when discretizing the imaging space into a grid, the imaging field of view is constructed with the detection area as the center, and the grid is divided according to the set spatial resolution to obtain the coordinates of each array element and the coordinates of each pixel. The system reads the array element configuration of the ultrasonic phased array instrument. Taking a 64-element linear probe with a center frequency of 1MHz and an array element spacing of 1.5mm as an example, the position coordinates of each array element in its three-dimensional coordinate system are obtained. The spatial coordinates of the first array element are set as (0,0,0), and the spatial coordinates of the sixty-fourth array element are set as (94.5,0,0). At the same time, the coordinates of each pixel are obtained. When solving the wavefront travel time field, the i-th array element, where i=1,2,…,64 is the starting point of the point source, is used. The fast travel method based on the second-order upwind difference scheme is used to solve the anisotropic equation in the current non-uniform anisotropic sound velocity field. This equation is used to describe the travel time relationship when the sound wave propagation speed changes with spatial position and propagation direction. The arrival time required for ultrasonic waves to propagate to each pixel in the imaging grid is calculated by using each transmitting element as a wave source, and the arrival times of each pixel are arranged according to their spatial positions to generate the corresponding wavefront travel time field matrix.

[0041] In this embodiment, the technical basis for explicitly limiting the use of a second-order upwind difference scheme in the fast travel method is as follows: Dams are complex media with large volumes and non-uniform anisotropic sound velocity inside. During the propagation of sound waves, there are significant local gradient changes due to spatial position variations. If a traditional first-order upwind difference scheme is used for spatial derivative discretization, severe numerical dissipation is introduced during the calculation, causing smooth truncation of the anisotropic wavefront travel time field at the spatial grid boundaries. This results in the accumulation of millisecond-level calculation errors in the two-way transmission and reception delay time, leading to edge defocusing and artifacts in the final ultrasonic imaging. However, the second-order upwind difference scheme can utilize wider spatial neighborhood node information for higher-order difference approximation, improving the numerical calculation accuracy of the travel time field to the microsecond level without increasing computational complexity, thus meeting the high-resolution pixel imaging requirements of 1mm×1mm.

[0042] After acquiring the wavefront travel time field, the system calculates the deviation angle between the phase velocity and group velocity in the current propagation direction caused by anisotropy. This deviation angle varies from 0 to 5 degrees. Based on this deviation angle, the wavefront normal gradient is corrected to obtain the group velocity vector direction. A fourth-order Runge-Kutta method is then used to perform step-wise backward tracking from the target pixel along the opposite direction of the group velocity vector. The tracking step size is set to 0.1 mm until the tracking path returns to the location of the transmitting or receiving element, thus generating a discretized set of group velocity propagation path points. For each step-size segment on this path point set, the corresponding slowness (i.e., the reciprocal of the local group velocity at the current position) is calculated using a path line integral, where the slowness value ranges from 2.0 × 10⁻⁶. -4 s / m to 2.6×10 -4 s / m. The forward propagation delay time from the transmitting element m to the pixel P is obtained by line integration. and the echo delay time of pixel P scattered to receiving element n .Will and The summation yields a transmit / receive delay time matrix containing all transmit / receive array element combinations and all pixels. To reduce the initial storage pressure, the transmit / receive delay time matrix can be calculated and cached in batches according to pixel blocks or channel blocks.

[0043] In this embodiment, the specific scheme for batch computation and caching is as follows: Global pixels are divided into multiple independent imaging pixel blocks, and effective channel combinations are divided into multiple channel blocks. During GPU parallel acceleration computation, only the delayed received signal corresponding to one current channel block and the group velocity propagation path of one current pixel block are read into high-speed video memory at a time. The line integral of the path slowness is executed in parallel using a single-instruction multi-threaded architecture. The two-way transmit / receive delay time of the current local block is calculated and cached in video memory for full-focus coherent superposition. Once the sound pressure amplitude of the local block is accumulated, the corresponding delay time cache is immediately released, and the calculation of the next pixel block and channel block is moved in. Through this dual-block batch iteration mechanism, the instantaneous video memory usage in a single iteration is effectively controlled, ensuring the feasibility of real-time computation during massive data throughput at both the hardware and software levels.

[0044] In an optional embodiment, the instantaneous phase is extracted, and the reciprocal of the discreteness is calculated to generate a normalized coherence coefficient map. The specific steps include: First, a spatial neighborhood is selected centered on the pixel, and the delayed received signals corresponding to all array channels within the spatial neighborhood are obtained. Then, a Hilbert transform is performed on the delayed received signals to construct an analytical signal, and the instantaneous phase angle of each array channel is extracted. Next, the circular square difference of the instantaneous phase angles of all array channels is calculated as the phase discreteness. The reciprocal of the sum of the phase discreteness and a preset minimum constant is taken, and then normalized to the interval between 0 and 1 through extreme value normalization to generate a normalized coherence coefficient map.

[0045] In the specific implementation details of the above feature extraction, a spatial neighborhood is selected centered on the currently calculated pixel, and a circular region with a radius set to half the wavelength of the ultrasonic wave is used as the spatial neighborhood window. At a sound speed of 4000 m / s and a frequency of 1 MHz, this half wavelength is 2 mm. Based on the aforementioned calculated transmit / receive delay time, delay correction is performed on the acquired full matrix acquisition data, and a time-domain signal segment with a preset time width is extracted from the echo signal of each transmit / receive channel, centered on the sampling point corresponding to the transmit / receive delay time. This preset time width is set to a sampling length covering 3 to 5 center frequency cycles, and a Hanning window is used to window the extracted segment to reduce the truncation effect. The system extracts the time-domain delayed received signal segments corresponding to all N transmit-receive channel combinations within this spatial neighborhood. For a 64-element system, the effective channel combination N is 4096. The above-mentioned windowed delayed received signal segments are then processed. Perform Hilbert transform to obtain orthogonal signal components As the imaginary part, it is combined with the original signal as the real part to construct an analytic signal in complex form. The analytical amplitude of the analytical signal is extracted, and the instantaneous phase angle of each array channel at the corresponding time point is calculated using the arctangent function. Its value range is distributed in - To + between.

[0046] To measure the phase consistency between channels, the extracted 4096 instantaneous phase angles were projected onto a unit circle on the complex plane, and the results were calculated using the vector summation formula. Calculate the length of the average phase vector .in, This represents the length of the average phase vector. This represents the total number of valid transmit-receive channel combinations. The imaginary unit, Indicates the first The instantaneous phase angle of each array channel. Further using the formula... The circular square error of this set of phases is calculated as the phase dispersion. The value range is between 0 and 1. Subsequently, the aforementioned phase dispersion is calculated. With the preset minimum constant The reciprocal of the sum of the values ​​is used to perform Min-Max extremum normalization across the global pixel image. This involves subtracting the global minimum from the reciprocal of the dispersion of all pixels, then dividing by the difference between the global maximum and the global minimum, to accurately map the values ​​to a standardized range of 0 to 1. This process is repeated for all pixels to generate a two-dimensional coherence coefficient map representing the spatial consistency of the scattered signals from the dam's internal structure.

[0047] In an optional embodiment, the direction-dependent sound velocity perturbation is determined based on the coherence coefficient maximization rule to update the non-uniform anisotropic sound velocity field. This mainly includes the following functional processing steps: Using a functional chain calculus derivation from the coherence coefficient to the instantaneous phase, transmit / receive delay time, and finally the global ray path, the functional gradient of the coherence coefficient with respect to the direction-dependent sound velocity parameters at each spatial grid node is calculated. This functional gradient is normalized, and the normalized update direction is multiplied by the update step size to perform gradient ascent, obtaining the sound velocity perturbation for the corresponding spatial grid node and the corresponding three-dimensional propagation direction parameters. Finally, the sound velocity perturbation is superimposed onto the corresponding direction-dependent sound velocity parameters in the current non-uniform anisotropic sound velocity field, completing the model update.

[0048] In the specific calculation process of model update, the global mean coherence coefficient is used as the optimization objective to construct the coherence coefficient objective function in each iteration cycle. To determine the update direction of the sound velocity field, the adjoint state method is used to perform functional chain calculus derivation to calculate the sensitivity of the objective function to direction-dependent sound velocity parameters. Specifically, according to the dependencies of coherence coefficient, instantaneous phase, transmit / receive delay time, and ray path slowness, the physical dependency functional chain of the objective function on coherence coefficient, instantaneous phase, transmit / receive delay time, and path slowness is calculated sequentially. The partial derivatives are propagated forward step by step using the multivariate calculus chain rule to obtain the degree of influence of sound velocity changes on the objective function under each spatial grid node and the corresponding three-dimensional propagation direction parameters. Then, the sensitivity information is backpropagated using the adjoint travel time field to form the multidimensional functional gradient of the coherence coefficient objective function with respect to the direction-dependent sound velocity parameters. This adjoint travel time field is composed of the back sensitivity distribution constructed based on the coherence coefficient objective function, which specifically represents the influence weight generated when the coherence coefficient change of each pixel is backpropagated along the ray path to the sound velocity field grid node. After normalizing the gradient of the multidimensional functional, the dimensionless update direction for subsequent sound velocity field update is obtained.

[0049] Following the optimization principle of maximizing the coherence coefficient to approximate the true sound field distribution, the learning rate is set as the update step size. For example, the initial update step size is set to 10 m / s. The reason for limiting the initial update step size to 10 m / s is that the reference longitudinal wave velocity of dam concrete material under zero stress and reference temperature is usually about 4000 m / s. Choosing 10 m / s as the initial update step size makes it about 0.25% of the reference sound velocity. This can provide sufficient gradient update sensitivity in the early stage of iteration to accelerate model convergence, while effectively avoiding numerical oscillations caused by excessive step size and the optimization process crossing the global optimum.

[0050] Based on this, specific step size adjustment algorithms are selected according to actual computing resources and sound field characteristics: when system memory resources are sufficient and a smooth sound velocity field with high-precision convergence is required, the L-BFGS optimization algorithm is adopted, and its storage history step number parameter for approximating the inverse Hessian matrix is ​​set to 5 to 10 steps; when the sound velocity distribution has many local extrema and the algorithm needs to have strong stability to quickly jump out of local optima, the Adam optimization algorithm is adopted, and the exponential decay rate of the first moment estimation is set. The exponential decay rate estimated by the second moment is 0.9. The value is 0.999. Using the selected optimization algorithm, the step size of the dimensionless update direction is adjusted and iteratively updated. Gradient ascent is performed to solve for the direction-dependent sound velocity perturbation for each grid node and the corresponding three-dimensional propagation direction parameters. To prevent divergence, the amplitude of a single direction-dependent sound velocity perturbation generated in each iteration is strictly limited to the range of -20 m / s to +20 m / s. This is based on the fact that strictly controlling the single fluctuation within 0.5% of the reference sound velocity ensures that when the sound velocity field model is iteratively corrected, the local sound velocity gradient between spatially adjacent grid nodes will not undergo non-physical drastic changes. This ensures that when subsequently solving the anisotropic equation using the fast travel method and performing step-by-step ray backward tracking using the fourth-order Runge-Kutta method, the wavefront travel time field matrix remains highly continuous, and the group velocity vector direction remains smooth and stable. This eliminates ray crossing or severe defocusing artifacts caused by local distortion of the sound velocity field from both physical and geometric perspectives, ensuring the numerical robustness of the imaging closed-loop optimization.

[0051] The system linearly superimposes the generated direction-dependent sound velocity perturbation onto the corresponding direction-dependent sound velocity parameters in the non-uniform anisotropic sound velocity field of the previous iteration. For example, if the original sound velocity at a grid point under a certain three-dimensional propagation direction parameter is 4020 m / s, and the calculated sound velocity perturbation is +5 m / s, then the updated sound velocity at that grid point under that three-dimensional propagation direction parameter is 4025 m / s. After completing a single iteration of model update, the new parameters are reintroduced into the next forward wave field propagation calculation. To improve the stability of convergence judgment, the system sets a minimum number of iterations and a maximum number of iterations. The difference between the coherence coefficient diagrams of two adjacent iterations is calculated using the normalized Frobenius norm, and this difference is compared with a preset convergence threshold. Comparison is performed. The preset convergence threshold is limited to... The basis for this is that the pixel values ​​of the coherence coefficient map in this invention have all been mapped to the standardized interval [0,1]. When the normalized Frobenius norm is used to measure the overall difference of its matrix, The threshold ensures that the spatial pixel-level structural consistency fluctuation is less than one ten-thousandth, at which point the boundary focusing and defect defocus compensation of the ultrasound image have reached substantial stability in a physical sense; if this threshold is greater than... If the optimization iteration terminates prematurely, the correction of the direction-dependent sound velocity field will be insufficient, and artifacts will still remain in the imaging; if this threshold is less than Due to the inherent random white noise in the sensing system, the matrix difference between two images is extremely difficult to overcome this minute step, easily causing the algorithm to fall into an invalid infinite loop and exhaust computational resources. Therefore, [the algorithm] is selected... As a convergence threshold, it can effectively balance imaging accuracy and algorithm convergence efficiency.

[0052] When the number of iterations reaches the minimum number of iterations, if the difference is lower than the convergence threshold for multiple consecutive iteration cycles, the convergence condition is determined to be met and the iteration stops; if the maximum number of iterations is reached but the convergence condition is still not met, the iteration result with the highest coherence coefficient is output as the convergence approximation result.

[0053] In an optional embodiment, after the iterative optimization process meets the convergence condition, the method further includes a full-field stress field inversion step inside the dam. Based on the temperature distribution data of each spatial grid node, the sound velocity offset caused by temperature is subtracted from the converged final sound velocity field, separating the direction-dependent sound velocity deviation caused by stress changes inside the dam. The direction-dependent sound velocity deviation is substituted into the acoustoelastic constitutive relation, and the actual principal stress magnitude and direction at each spatial grid node are calculated by inverse nonlinear least squares method. The actual principal stress of each grid node is then mapped to three-dimensional space to generate a full-field stress distribution cloud map inside the dam.

[0054] At this stage, the system acquires the final sound velocity values ​​of each spatial grid node in each direction after correction. Since the change in sound velocity is dominated by both temperature and stress effects, the system, based on the pre-acquired temperature distribution data of each grid node, subtracts the sound velocity offset caused by temperature, accurately separating the direction-dependent sound velocity deviation caused purely by the stress changes inside the dam. This deviation is then substituted as a known variable into the aforementioned acoustoelastic constitutive equations for inverse constraint solution, thereby obtaining the actual principal stress tensor at each pixel grid node inside the dam and completing the three-dimensional spatial mapping, ultimately achieving the visualization output of the full-field stress distribution cloud map inside the dam.

[0055] S3: Anomaly detection and early warning based on texture features and dynamic thresholds.

[0056] Wavelet transform is performed on ultrasonic images to extract gray-level co-occurrence matrix features, which are then concatenated into a high-dimensional texture feature vector. A mapping relationship is established between real-time temperature data and real-time stress data and the upper limit of the benchmark anomaly score, and the anomaly judgment threshold is dynamically updated. The high-dimensional texture feature vector is input into the isolated forest anomaly detection model, and dam early warning information is output when the anomaly score exceeds the updated anomaly judgment threshold.

[0057] In its implementation, the system first performs a multi-scale two-dimensional discrete wavelet transform on the generated ultrasonic image and uses the Daubechies wavelet basis function to decompose the image into a low-frequency approximate subband and multiple high-frequency detail subbands. Then, at each decomposition scale, it calculates the gray-level co-occurrence matrices in four directions: 0 degrees, 45 degrees, 90 degrees, and 135 degrees, and extracts three texture features: energy features, contrast features, and correlation features. The feature scalars of all scales and directions are then tiled, concatenated, and stitched together to form a one-dimensional high-dimensional texture feature vector. Simultaneously, a nonlinear regression prediction model is established using a support vector regression algorithm. This model is pre-trained using a large amount of real-time temperature and stress data under historical normal conditions as input and the corresponding upper limit of the isolated forest benchmark anomaly score as the output label. During real-time detection, the currently acquired real-time temperature and stress data are input into this trained support vector regression prediction model to calculate and output the upper limit of the benchmark anomaly score under the current environmental condition. Based on this upper limit, the anomaly judgment threshold is updated.

[0058] During the training phase of the anomaly detection model, the system constructs a training sample set based on high-dimensional texture feature vectors extracted from historical normal operating conditions, thereby establishing an isolated forest anomaly detection model. This isolated forest anomaly detection model consists of multiple random isolated trees. During training, a portion of samples and features are randomly extracted from each isolated tree, and the sample space is recursively divided by randomly selecting features and random splitting points until the samples are isolated or the preset tree depth is reached, thus forming an unsupervised anomaly detection model.

[0059] In this embodiment, the key training parameters of the isolated forest anomaly detection model are specifically limited as follows: the number of random isolated trees is set to 100 to 200, preferably 150; the subsampling sample size of each random isolated tree is set to 256 samples; and the maximum tree depth is set to 8 layers. The above structural parameters are limited because: 150 random isolated trees can statistically ensure efficient convergence of anomaly scores, avoiding occasional score fluctuations caused by random segmentation of a single tree; limiting the subsampling sample size to 256 ensures that the spatial partitioning of a single tree maintains local sparsity, thereby effectively combating the hidden noise interference that may be included in the historical normal operating data of the dam, allowing the true structural damage feature vectors to be quickly isolated as a very small number of extremely outlier paths, thus improving the anomaly detection rate of high-dimensional texture features.

[0060] After entering the real-time detection phase, the one-dimensional high-dimensional texture feature vector corresponding to the current ultrasonic image is input into the trained isolated forest anomaly detection model. The number of node layers traversed by the current sample from the root node to the leaf node in each isolated tree is counted, and the average path length of all isolated trees is calculated to obtain the average path length. The system then normalizes this average path length using a path length normalization coefficient corresponding to the number of training samples, and converts it into anomaly score using an exponential function. The shorter the average path length, the easier it is for the sample to be isolated by the random segmentation process, and the higher its normalized anomaly score and the greater the degree of anomaly. When the normalized anomaly score is greater than the updated anomaly judgment threshold, the system triggers an alarm process and pushes dam structure damage warning information containing the anomaly score, detection time, and warning type to the reservoir safety management terminal via a short message point-to-point protocol service gateway. If combined with the aforementioned full-field stress distribution cloud map inside the dam, the pushed warning information can further include stress over-limit warnings and stress gradient mutation warnings.

[0061] In an optional embodiment, wavelet transform is performed on the ultrasonic image to extract gray-level co-occurrence matrix features, which are then concatenated into a high-dimensional texture feature vector. Specifically, the feature extraction process includes the following steps: Multi-scale redundant wavelet transform is used to decompose the ultrasonic image to obtain low-frequency and high-frequency sub-band images at different scales. Then, gray-level co-occurrence matrices are calculated for each low-frequency and high-frequency sub-band image in four directions: 0 degrees, 45 degrees, 90 degrees, and 135 degrees. Finally, energy features, contrast features, and correlation features are extracted from each gray-level co-occurrence matrix, and the features extracted at all scales and in all directions are concatenated in a preset order to form a high-dimensional texture feature vector.

[0062] In the implementation details of the above feature extraction, the fused ultrasonic 2D image generated after iterative convergence is subjected to a continuous three-level multi-scale wavelet decomposition operation using the Db4 translation invariant redundant wavelet basis function. During this non-downsampling process, each level generates one low-frequency approximate sub-band image representing the overall image structure and three high-frequency detail sub-band images representing edge details in different directions. Specifically, the high-frequency detail sub-band images include horizontal, vertical, and diagonal sub-bands. The system extracts one lowest-frequency approximate image and nine multi-scale high-frequency detail images. For these 10 sub-band images with unchanged dimensions, their gray levels are uniformly divided into 64 levels, specifically representing discrete gray levels from 0 to 63. The system sets the spatial displacement distance of co-occurring pixels to 2 pixels and calculates the joint probability of adjacent pixel gray-level pairs along four spatial directions: 0 degrees, 45 degrees, 90 degrees, and 135 degrees. This generates four 64x64 gray-level co-occurrence matrices for each sub-band image, resulting in a total of 40 gray-level co-occurrence matrices.

[0063] After obtaining all gray-level co-occurrence matrices, the system calculates three core texture parameters by normalizing the matrix elements. First, an energy feature representing the uniformity of gray-level distribution in the image is extracted. This feature is obtained by summing the squares of each element in the gray-level co-occurrence matrix. Second, a contrast feature representing the drasticness of local gray-level changes is extracted. This feature is obtained by weighting and summing the squares of the gray-level differences between each gray-level combination according to their occurrence probabilities. Finally, a correlation feature representing the degree of linear correlation between gray levels is extracted. This feature is obtained by calculating the product of the deviations of paired gray levels relative to their respective means, weighting and summing them according to their occurrence probabilities, and then normalizing by the standard deviation. It is mainly used to characterize the consistency and directionality of gray-level changes between adjacent pixels. The system flattens the three types of texture features extracted from 10 sub-band images in four calculation directions and seamlessly stitches them together in a fixed order: first sub-band, then direction, then feature type. This generates a high-dimensional texture feature vector of 120 dimensions, which is used as input feature data to accurately represent the structural anomaly state of ultrasonic scattering inside the dam.

[0064] In an optional embodiment, a mapping relationship is established between real-time temperature data and real-time stress data and a baseline anomaly score upper limit, and the anomaly determination threshold is dynamically updated. This includes the following prediction and update steps: Collect multiple sets of temperature data, multiple sets of stress data, and the corresponding baseline anomaly scores output by the isolated forest anomaly detection model under historical normal conditions. Determine the baseline anomaly score upper limit based on the maximum baseline anomaly score within a preset time window. Use a support vector regression algorithm, taking multiple sets of temperature data and multiple sets of stress data as input and the baseline anomaly score upper limit as output, to train a regression prediction model. Input the real-time collected temperature data and stress data into the regression prediction model, output the baseline anomaly score upper limit under the current environmental condition, and add a preset fault tolerance margin as the updated anomaly determination threshold.

[0065] During the initialization and training phases of the dynamic threshold model, the system continuously extracts multiple sets of environmental variables from the historical database of the dam's Internet of Things (IoT) for periods without structural anomaly labels over at least one full year. Specifically, taking 8760 hours of continuous data as an example, the system extracts average temperature data from multiple points on the dam body, including feedback from a distributed fiber optic temperature measurement system, at a sampling frequency of once per hour. The typical fluctuation range is between -15℃ and 40℃. Simultaneously, it extracts comprehensive stress data measured in the corresponding area, with a typical fluctuation range between -3MPa and 3MPa. The system inputs the ultrasonic high-dimensional texture feature vectors from the same period into a pre-trained isolated forest anomaly detection model to calculate anomaly scores, collecting benchmark anomaly scores corresponding to these 8760 sets of normal operating conditions. These values ​​are mainly distributed in the range of 0.35 to 0.55. The system sets a sliding time window of 24 hours and takes the maximum anomaly score within this window as the corresponding label, thus constructing a standardized dataset that correlates the two-dimensional independent variables of temperature and stress with the upper limit of the anomaly score as a one-dimensional dependent variable. Using this normal-state dataset, a radial basis function was selected as the kernel function, with a penalty parameter of 10 and a kernel parameter of 0.01. Support vector regression was then used to train the model, thereby establishing a nonlinear regression prediction model for the impact of environmental load on the baseline for isolated forest assessment. It should be noted that this regression prediction model was trained solely on historical temperature and stress data and corresponding anomaly scores under normal operating conditions, specifically designed to characterize the upper limit of anomaly scores as environmental conditions change under normal service conditions. This model is not involved in ultrasonic image generation or isolated forest anomaly score calculation; it is only used for the dynamic calibration of the end-point alarm threshold.

[0066] In this embodiment, the penalty parameter of 10 and the kernel parameter of 0.01 were determined through joint optimization using 5-fold cross-validation and grid search algorithms on a historical normal state dataset. The specific selection criteria are as follows: The penalty parameter is used to balance the training error and structural complexity of the regression model. Setting it to 10 can tolerate the small fitting residuals caused by the inherent random measurement white noise of IoT temperature and stress sensors, while also applying sufficient penalty to outliers that deviate from the general trend of normal seasonal operating conditions, effectively avoiding model overfitting. The kernel parameter determines the spatial influence range of the radial basis function. Setting it to 0.01 ensures that the influence radius of the support vector is highly matched with the gradual change period of the dam's annual temperature and pressure load on a physical scale. This ensures that the regression prediction model can still accurately output the upper limit of the baseline outlier score under the current environmental conditions when facing future extreme alternating temperature difference conditions, and has generalization calibration capabilities.

[0067] During real-time online monitoring of reservoir safety, the IoT gateway synchronously collects real-time dam temperature and compressive stress data every minute. Assuming the measured temperature is 22.5℃ and the compressive stress is 1.2MPa, after preprocessing, the system inputs this set of environmental parameter vectors into a pre-trained support vector regression prediction model to calculate and output the upper limit of the baseline anomaly score under the current environmental conditions, assuming the model's predicted value is 0.52. Considering the dam safety level and the inherent random white noise effect of the sensing system, a fixed constant proportionality coefficient is set as a preset fault tolerance margin. In practice, an absolute margin constant of 0.05 can be selected, or 10% of the baseline value can be used as the relative fluctuation range. Adding the preset fault tolerance margin of 0.05 to the predicted upper limit of the baseline anomaly score (0.52) yields an updated anomaly judgment threshold of 0.57. This updated threshold is used to determine the real-time anomaly score output by the isolated forest, thereby effectively reducing the probability of false alarms triggered by seasonal environmental changes in safety warnings. In this embodiment, the technical basis for limiting the absolute margin constant to 0.05 is that, under normal historical operating conditions, the statistical fluctuation variance of the isolated forest benchmark anomaly score caused by measurement white noise from the IoT sensing system and uncontrollable minor environmental disturbances typically ranges from 0.01 to 0.02. Introducing an absolute margin of 0.05 is equivalent to establishing a dynamic confidence protection band of approximately three standard deviations for the peak score under normal operating conditions. This value can perfectly absorb the artificially high scores caused by hardware noise and seasonal environmental fluctuations, thereby eliminating false alarms; at the same time, it will not severely dull the sensitivity of the algorithm due to an excessively large margin setting, such as greater than 0.10, ensuring that when the dam structure suffers early, hidden damage, the small increment of the anomaly score can promptly break through this threshold, guaranteeing the detection sensitivity of the early warning system.

[0068] To verify the actual effectiveness of the aforementioned monitoring system and methods, the experiment was conducted using a long-term online safety monitoring system for a large reservoir's concrete gravity dam. The experimental environment realistically covered alternating temperatures from -15℃ to 40℃ and alternating stress states from -3MPa to 3MPa. The experiment divided the verification system configuration into four comparative models: a baseline group, ablation group one, ablation group two, and a complete scheme group. The baseline group adopted a traditional uniform isotropic sound velocity field combined with single-scale gray-level features and a static anomaly threshold detection algorithm. Ablation group one, based on the baseline group, input a non-uniform anisotropic sound velocity field jointly corrected by temperature and stress, and a coherence coefficient functional gradient update model. Ablation group two, based on ablation group one, added a multi-scale redundant wavelet transform and a high-dimensional feature extraction algorithm for gray-level co-occurrence matrix. The complete scheme group included all the aforementioned innovative features and added an anomaly determination threshold update model jointly driven by support vector regression and multi-source environmental data.

[0069] The experiment collected 8760 hours of continuous ultrasonic full-matrix capture data and synchronous environmental load monitoring data from the dam. A total of 525,600 valid detection samples were provided, based on one valid detection record per minute. During threshold regression model training, the minute-level data was aggregated hourly to form 8760 sets of environmental state samples. The test set included 240 real structural anomaly samples, such as microcracks and internal pores, from indoor experiments and historical data calibration. The test platform was equipped with a dual-channel high-performance graphics processor, with a uniform ultrasonic center frequency of 1MHz, an imaging spatial grid resolution of 1mm, and an iterative update step size of 10m / s. The performance of each core module was comprehensively evaluated by calculating the peak signal-to-noise ratio of the images, the spatial localization error of microdefects, the accurate detection rate of anomaly samples, and the number of false alarms without faults throughout the year when the four models processed this dataset.

[0070] Test results showed that the baseline group had a peak signal-to-noise ratio (PSNR) of only 22.4 dB, a micro-defect localization error of 4.5 mm, an anomaly detection rate of 71.3%, and a high number of false alarms (145) throughout the year. Ablation group one improved its PSNR to 28.1 dB, reduced the localization error to 2.2 mm, increased the detection rate to 84.6%, and reduced the number of false alarms to 121. Ablation group two maintained a PSNR of 28.3 dB, a localization error of 2.1 mm, increased the detection rate to 93.2%, and reduced the number of false alarms to 108. The complete solution group, employing all steps, demonstrated the best performance, with an ultrasound image PSNR as high as 32.5 dB, a significant reduction in micro-defect spatial localization error to 0.8 mm, an anomaly detection rate of 98.9%, and a substantial reduction in the number of false alarms throughout the year to 4.

[0071] The experimental data demonstrates that utilizing a non-uniform anisotropic sound velocity field and applying ray tracing and functional gradient models effectively corrects sound velocity drift and travel time errors under temperature and stress, improving image signal-to-noise ratio and enhancing positioning accuracy. The multi-scale gray-level co-occurrence matrix feature extraction scheme fully exploits edge scattering texture features across different spatial dimensions, enhancing the model's ability to identify and detect early-stage minute defects. Finally, applying support vector regression to fuse environmental stress and temperature data for dynamic adjustment of the anomaly detection threshold effectively mitigates background noise interference caused by seasonal temperature differences and load fluctuations, resulting in a sharp decrease in the system's false alarm rate. This strongly proves that the entire scheme possesses excellent monitoring accuracy and overall stability under actual service conditions.

[0072] Figure 2This diagram illustrates the optimization of the global coherence coefficient mean curve. The circular scatter points connected by dashed lines represent the main distribution trajectory. Observing the evolution of the circular scatter points' positions in the coordinate system reveals a significant upward trend in the trajectory as the number of iterations on the horizontal axis increases, before entering a state of gentle local fluctuations after forty iterations. This phenomenon corresponds to the iterative optimization process described in the embodiment, proving that the coherence coefficient maximization rule can drive the non-uniform anisotropic sound velocity field model to continuously converge towards the actual physical state of the dam.

[0073] Figure 3 This is a bar chart showing the numerical distribution of contrast features at different wavelet decomposition scales and different gray-level co-occurrence matrix directions. Above each specific directional scale, three bar symbols are arranged side-by-side. In each group of side-by-side bar symbols, the first, second, and third bars, arranged from left to right, correspond to the first, second, and third decomposition scales of the multi-scale redundant wavelet transform, respectively. Observation reveals a significant step-like difference in the height of the bars at different arrangement positions, and the distribution height of the bars varies at different directional scales. This demonstrates that high-dimensional texture feature vectors can extract subtle differences in the internal structure of the dam from multiple spatial dimensions and frequency scales.

[0074] Figure 4 This is a schematic diagram comparing the changes in the dynamically updated anomaly detection threshold and the static baseline threshold over a 24-hour monitoring period. The thick solid line, which remains horizontally extended, represents the static baseline threshold. The square scatter points connected by the dotted line constitute the trajectory of the dynamically updated anomaly detection threshold. Observing the distribution of the square scatter points reveals that the line connecting the points exhibits a significant fluctuating trend over time, and is located above or below the static baseline threshold for most of the time periods. This corresponds to the physical characteristics described in the specific implementation method, which utilizes a support vector regression model to respond to environmental load changes in real time and dynamically correct the threshold, demonstrating that the system can effectively eliminate interference caused by environmental fluctuations.

[0075] This invention also discloses a dam internal stress monitoring system, including a processor and a memory. The memory stores computer program instructions, which, when executed by the processor, implement the dam internal stress monitoring method according to the invention.

[0076] The system also includes other components well known to those skilled in the art, such as communication buses and communication interfaces, the settings and functions of which are known in the art and will not be described in detail here.

[0077] It should be noted that those skilled in the art can make various modifications and improvements without departing from the inventive concept, and these all fall within the scope of protection of this invention. Therefore, the scope of protection of this patent should be determined by the appended claims.

Claims

1. A method for monitoring internal stress of a dam, characterized in that, include: S1: Acquire full-matrix capture data of the dam surface, collect real-time temperature and stress data from embedded sensors within the dam body via the Internet of Things, and construct an initial non-uniform anisotropic sound velocity field based on the real-time temperature data, real-time stress data, finite element model, and acoustoelastic constitutive relation; specifically including: dividing the finite element model into a three-dimensional mesh, calculating the stress tensor of each mesh node based on the real-time stress data, and determining the temperature distribution and temperature correction of each mesh node based on the real-time temperature data; acquiring the inherent acoustoelastic coefficient tensor of the dam material; combining the initial reference sound velocity under zero stress and reference temperature conditions, and based on the acoustoelastic constitutive relation, calculating the sound velocity offset of each mesh node under the three-dimensional propagation direction parameters characterized by the propagation azimuth and pitch angles through the inherent acoustoelastic coefficient tensor, the stress tensor of each mesh node, and the temperature correction of each mesh node, and superimposing them to obtain the initial non-uniform anisotropic sound velocity field; S2: Execute The process involves iterative optimization, using ray tracing to calculate the group velocity propagation path and transmission / reception delay time from each array element to the pixel. Based on the transmission / reception delay time, the captured data of the entire matrix is ​​coherently superimposed to obtain the sound pressure amplitude map. The instantaneous phase is extracted, and the inverse of the discreteness is calculated to generate a normalized coherence coefficient map. If the change in the coherence coefficient map is not lower than the convergence threshold, the direction-dependent sound velocity perturbation is determined based on the coherence coefficient maximization rule to update the non-uniform anisotropic sound velocity field and the iteration is repeated. If the change in the coherence coefficient map is lower than the convergence threshold, the sound pressure amplitude map and the coherence coefficient map are fused to generate an ultrasonic image. S3: Wavelet transform is performed on the ultrasonic image to extract gray-level co-occurrence matrix features and concatenate them into a high-dimensional texture feature vector. A mapping relationship is established between real-time temperature data and real-time stress data and the upper limit of the benchmark anomaly score, and the anomaly judgment threshold is dynamically updated. The high-dimensional texture feature vector is input into the isolated forest anomaly detection model. When the anomaly score exceeds the updated anomaly judgment threshold, dam early warning information is output.

2. The method for monitoring internal stress of a dam according to claim 1, characterized in that, The method of calculating the group velocity propagation path and transmission / reception delay time from each array element to the pixel using ray tracing includes: discretizing the imaging space into a grid, obtaining the coordinates of each array element and the coordinates of the pixel, using each array element as the starting point, solving the equation of motion in the current non-uniform anisotropic sound velocity field using the fast travel method to obtain the wavefront travel time field; combining the deviation relationship between phase velocity and group velocity, determining the group velocity direction based on the three-dimensional propagation direction parameters and the current non-uniform anisotropic sound velocity field, and tracing backward from the pixel to the corresponding array element along the group velocity direction to obtain the group velocity propagation path; performing line integration on the slowness along the group velocity propagation path, calculating the propagation delay time when each array element is used as a transmitter and receiver respectively, and adding them together to obtain the transmission / reception delay time.

3. The method for monitoring internal stress of a dam according to claim 2, characterized in that, The method of solving the anisotropic equation in the current non-uniform anisotropic sound velocity field using the fast travel method to obtain the wavefront travel time field, and tracing backward from the pixel to the corresponding array element along the group velocity direction, includes: solving the anisotropic equation in the current non-uniform anisotropic sound velocity field using the fast travel method to generate the arrival time of each pixel and construct the wavefront travel time field matrix; and using the fourth-order Runge-Kutta method to perform step-by-step backward tracing from the target pixel along the opposite direction of the group velocity vector until returning to the position of the corresponding array element, generating a discretized set of group velocity propagation path points.

4. The method for monitoring internal stress of a dam according to claim 1, characterized in that, The step of extracting the instantaneous phase, calculating the reciprocal of the discreteness, and generating a normalized coherence coefficient map includes: selecting a spatial neighborhood centered on the pixel, obtaining the delayed received signals corresponding to all array channels within the spatial neighborhood; performing a Hilbert transform on the delayed received signals to construct an analytical signal, and extracting the instantaneous phase angle of each array channel; calculating the circular square difference of the instantaneous phase angles of all array channels as the phase discreteness, taking the reciprocal of the sum of the phase discreteness and a preset minimum constant, and mapping it to the interval from 0 to 1 through extreme value normalization to generate the normalized coherence coefficient map.

5. The method for monitoring internal stress of a dam according to claim 1, characterized in that, The method of determining the direction-dependent sound velocity perturbation based on the coherence coefficient maximization rule to update the non-uniform anisotropic sound velocity field includes: using a functional chain calculus derivation from the coherence coefficient to the instantaneous phase, the transmit / receive delay time, and even the global ray path, to calculate the functional gradient of the coherence coefficient with respect to the direction-dependent sound velocity parameters at each spatial grid node; normalizing the functional gradient and multiplying the normalized update direction by the update step size to perform gradient ascent, thereby obtaining the sound velocity perturbation at the corresponding spatial grid node and the corresponding three-dimensional propagation direction parameters; and superimposing the sound velocity perturbation onto the corresponding direction-dependent sound velocity parameters in the current non-uniform anisotropic sound velocity field to complete the model update.

6. The method for monitoring internal stress of a dam according to claim 1, characterized in that, The step of extracting gray-level co-occurrence matrix features from the ultrasonic image using wavelet transform and concatenating them into a high-dimensional texture feature vector includes: decomposing the ultrasonic image using multi-scale redundant wavelet transform to obtain low-frequency sub-band images and high-frequency sub-band images at different scales; calculating the gray-level co-occurrence matrix for each low-frequency sub-band image and high-frequency sub-band image in four directions: 0 degrees, 45 degrees, 90 degrees, and 135 degrees; extracting energy features, contrast features, and correlation features from each gray-level co-occurrence matrix; and concatenating the features extracted at all scales and in all directions in a preset order to form the high-dimensional texture feature vector.

7. The method for monitoring internal stress of a dam according to claim 1, characterized in that, The process of establishing a mapping relationship between real-time temperature data, real-time stress data, and the upper limit of the baseline anomaly score, and dynamically updating the anomaly judgment threshold, includes: collecting multiple sets of temperature data, multiple sets of stress data, and the corresponding baseline anomaly scores output by the isolated forest anomaly detection model under historical normal conditions; determining the upper limit of the baseline anomaly score based on the maximum value of the baseline anomaly score within a preset time window; training a regression prediction model using a support vector regression algorithm, with the multiple sets of temperature data and multiple sets of stress data as input and the upper limit of the baseline anomaly score as output; inputting the real-time collected temperature data and stress data into the regression prediction model, outputting the upper limit of the baseline anomaly score under the current environmental condition, and adding a preset fault tolerance margin as the updated anomaly judgment threshold.

8. The method for monitoring internal stress of a dam according to claim 1, characterized in that, After the iterative optimization process meets the convergence condition, the method further includes a full-field stress field inversion step inside the dam: based on the temperature distribution data of each spatial grid node, the sound velocity offset caused by temperature is subtracted from the final converged sound velocity field, and the direction-related sound velocity deviation caused by the stress change inside the dam is separated; the direction-related sound velocity deviation is substituted into the acoustoelastic constitutive relation, and the actual principal stress magnitude and direction at each spatial grid node are calculated by inverse solution using the nonlinear least squares method, and the actual principal stress of each grid node is mapped to three-dimensional space to generate a full-field stress distribution cloud map inside the dam.

9. A dam internal stress monitoring system, characterized in that, include: A processor and a memory, the memory storing computer program instructions that, when executed by the processor, implement the method for monitoring internal stress of a dam according to any one of claims 1-8.

Images