Tunnel construction data processing method based on real-time monitoring

Abstract

This invention relates to the field of electronic digital data processing technology and discloses a method for processing tunnel construction data based on real-time monitoring. The method involves acquiring real-time data from displacement, stress, and micro-vibration sensors; constructing a spatial graph adjacency matrix based on coordinates and rock strata topology; and plotting the dynamic features of nodes in the time series graph. The data is then input into a physically constrained spatiotemporal graph convolutional network. A displacement-stress transmission constraint term based on the elasticity equation is introduced as a weight attenuation penalty factor in the graph convolutional propagation layer to restrict cross-node information transmission in accordance with the laws of surrounding rock mechanics. The dynamic features of the time dimension are extracted through a time-gated cyclic unit, and a tensor data stream is output. This invention restricts weight updates within physical boundaries, eliminates non-physical interference features from spatial aggregation, avoids feature confusion during multi-source data fusion, and overcomes the problem of distorted precursor feature extraction due to a lack of physical constraints. Furthermore, improvements are made to the coupling processing of multi-source heterogeneous data under geological topology, taking into account classification numbers.

Description

Technical Field

[0001] This invention relates to the field of electronic digital data processing technology, and discloses a method for processing tunnel construction data based on real-time monitoring. Background Technology

[0002] During tunnel construction, displacement sensors, stress sensors, and micro-vibration sensors are required to acquire data on the surrounding rock condition. Existing data processing solutions, when handling multi-source heterogeneous monitoring data, typically construct a spatial graph structure using the geographic coordinates of each sensor as nodes. A graph neural network is then used to aggregate data from adjacent nodes to extract spatial features, followed by a recurrent neural network to extract dynamic features over time. In this conventional approach, information transfer between nodes in the graph convolutional propagation layer relies entirely on data-driven topological distance weights. During training, the network automatically updates the connection weights between nodes based on the backpropagation gradient of the loss function, achieving the fusion processing of data from different sensors.

[0003] When processing tunnel monitoring data under complex geological conditions, the aforementioned existing technical solutions lack physical boundary constraints in updating the connection weights between adjacent nodes in the graph convolution propagation layer. Rock deformation and stress transmission are constrained by the equations of elasticity and have specific physical transmission paths. However, conventional graph neural networks, when updating weights based on data, introduce spatially close but mechanically unrelated sensor data into the propagation path. This leads to the introduction of non-physically related interference features when aggregating multi-source heterogeneous data, causing feature confusion and making it impossible to extract true precursor features of surrounding rock deformation in complex fault areas. Summary of the Invention

[0004] The purpose of this invention is to provide a method for processing tunnel construction data based on real-time monitoring, which can effectively solve the problems mentioned in the background art.

[0005] To achieve the above objectives, the technical solution adopted by the present invention is as follows: The data processing method for tunnel construction based on real-time monitoring includes: acquiring real-time monitoring data from displacement sensors, stress sensors, and micro-vibration sensors in the tunnel construction environment; Based on the geographic coordinates of each sensor and the topological connections of the rock strata, an adjacency matrix of the spatial graph structure is constructed, and the time series sampling points are used as the dynamic features of the graph nodes. The real-time monitoring data is input into a pre-constructed physical constraint spatiotemporal graph convolutional network. A displacement-stress transmission constraint term based on the elasticity equation is introduced as a weight attenuation penalty factor in the graph convolution propagation layer of the physical constraint spatiotemporal graph convolutional network to restrict the cross-node information transmission path to conform to the surrounding rock mechanics transmission law. Dynamic features of the time dimension are extracted by the time-gated recurrent unit embedded in the physical constraint spatiotemporal graph convolutional network. Output a unified tensor data stream that integrates heterogeneous data from multiple sources and eliminates spatiotemporal misalignment errors.

[0006] Preferably, the step of constructing the adjacency matrix of the spatial graph structure based on the geographic coordinates of each sensor and the topological connection relationship of the rock strata includes: using the displacement sensor, stress sensor and micro-vibration sensor as graph nodes of the spatial graph structure, respectively. Calculate the Euclidean distance between adjacent graph nodes, and use the reciprocal of the Euclidean distance as the initial spatial connectivity weight; Obtain a ground-penetrating radar scan profile of the tunnel excavation face and identify the locations of rock strata interfaces and fault fracture zones in the ground-penetrating radar scan profile. When the line connecting two graph nodes penetrates the rock strata interface or fault fracture zone, the initial spatial connection weight is multiplied by a preset geological weakening coefficient to generate the adjacency matrix.

[0007] Preferably, the dynamic feature of using time series sampling points as graph nodes includes: standardizing the real-time monitoring data of the displacement sensor, stress sensor and micro-vibration sensor respectively to generate standard time series of each sensor; The standard time series is truncated using a fixed time window, and the standard time series of each sensor within the same time window are concatenated into a high-dimensional feature vector; The high-dimensional feature vector is assigned to the corresponding graph node as the dynamic feature of the graph node at the current time step, and the high-dimensional feature vector is updated using a sliding window in the next time step.

[0008] Preferably, the step of introducing a displacement-stress transmission constraint term based on the elasticity equation as a weight attenuation penalty factor into the graph convolution propagation layer of the physically constrained spatiotemporal graph convolutional network includes: extracting the learnable weight matrix between adjacent graph nodes in the graph convolution propagation layer; Establish an expression for the relationship between displacement and stress partial derivatives based on the three-dimensional elasticity equilibrium differential equation; Calculate the sum of squared residuals between the updated graph node feature distribution and the expression relating displacement and stress partial derivatives; The sum of squared residuals is added as the weight decay penalty factor to the total loss function of the physically constrained spatiotemporal graph convolutional network.

[0009] Preferably, the step of extracting dynamic change features of the time dimension through the time-gated recurrent unit embedded in the physically constrained spatiotemporal graph convolutional network includes: inputting the graph node feature sequence containing spatial topological associations output by the graph convolutional propagation layer into the time-gated recurrent unit; A causal convolution kernel is introduced into the update gate and reset gate of the time-gated loop unit, and the causal convolution kernel is used to expand the local temporal receptive field of the graph node feature sequence. The update gate controls the retention ratio of the hidden state of the historical time step, and the reset gate filters out the historical time step interference information that is irrelevant to the current surrounding rock deformation state, outputting a time hidden state vector containing spatiotemporal coupling characteristics.

[0010] Preferably, the output of a unified tensor data stream that integrates multi-source heterogeneous data and eliminates spatiotemporal misalignment errors includes: rearranging the temporal hidden state vector output by the time-gated loop unit according to the graph node number to generate a two-dimensional spatiotemporal feature matrix. Tensor expansion is performed on the two-dimensional spatiotemporal feature matrix to map the multi-source sensor channel dimension, time step dimension, and spatial node dimension into a three-dimensional unified tensor. Singular value decomposition is performed on the three-dimensional unified tensor along the sensor channel dimension, suboptimal components smaller than a preset singular value threshold are truncated, and the dimension-reduced unified tensor data stream is output.

[0011] Preferably, after generating the adjacency matrix, the process further includes a dynamic correction process for the adjacency matrix: real-time extraction of the location coordinates of the surrounding rock micro-seismic events monitored by the micro-vibration sensor; A spatial attenuation function is constructed centered on the location coordinates of the microseismic events in the surrounding rock. The spatial connectivity weights of graph nodes within the influence radius of microseismic events in the adjacency matrix are enhanced and superimposed using the spatial attenuation function. When a preset number of microseismic events in the surrounding rock are triggered consecutively in the same area, a virtual stress transmission node is added to the adjacency matrix, and a bidirectional connection edge is established between the virtual stress transmission node and the affected graph node.

[0012] Preferably, when the residual sum of squares is superimposed as the weight decay penalty factor into the total loss function of the physically constrained spatiotemporal graph convolutional network, a dynamic adaptive coefficient adjustment method is adopted: the variance ratio of the data-driven loss term in the current training batch to the residual sum of squares during gradient backpropagation is calculated. When the variance ratio is greater than the preset imbalance threshold, a dual optimization equation for the data-driven loss term and the weight decay penalty factor is constructed using the Lagrange multiplier method. By solving the dual optimization equation, the weight decay penalty factor in the current iteration weight coefficient of the total loss function is reduced.

[0013] Preferably, after filtering out historical time step interference information that is irrelevant to the current surrounding rock deformation state through the reset gate, a time series anomaly focusing mechanism is also included: calculating the cosine similarity matrix of the graph node feature sequence between the current time step and the historical time step; The cosine similarity matrix is ​​input into a multilayer perceptron, and the attention score distribution representing the probability of state change at each time step is output. The attention score distribution is used to weight the time hidden state vector output by the update gate element by element, thereby amplifying the characteristic amplitude in the time hidden state vector corresponding to the sudden release of surrounding rock stress.

[0014] Preferably, after reducing the weight decay penalty factor in the current iteration weight coefficient of the total loss function, an early stopping mechanism based on the physical boundary violation rate is also included: after each weight coefficient update, the learnable weight matrix is ​​extracted and the equivalent elastic modulus transfer value between adjacent graph nodes is calculated. Determine whether the equivalent elastic modulus transfer value exceeds the preset extreme value range of the physical material elastic modulus, and count the proportion of the number of node pairs that exceed the extreme value range of the physical material elastic modulus to the total number of node pairs as the physical boundary violation rate. When the physical boundary violation rate shows an upward trend for a continuous preset window period, the model training is forcibly stopped early and the current network weight parameters are locked.

[0015] Compared with the prior art, the beneficial effects of the present invention are as follows: This invention introduces a displacement-stress transmission constraint term based on the equation of elasticity as a weight attenuation penalty factor into the graph convolution propagation layer of a physically constrained spatiotemporal graph convolutional network. This ensures that the cross-node information transmission path between graph nodes conforms to the transmission law of surrounding rock mechanics. This technique restricts the data-driven weight update process within physical boundaries, eliminates non-physical interference features generated during spatial aggregation, avoids feature confusion during the fusion of multi-source heterogeneous data, and overcomes the problem of distortion in the extraction of precursor features of surrounding rock deformation caused by the lack of physical constraints in conventional graph neural networks.

[0016] This invention adjusts the initial spatial connectivity weights of the adjacency matrix using geological weakening coefficients based on the locations of rock strata interfaces and fault fracture zones obtained from ground-penetrating radar (GPR) scan profiles, ensuring that the graph structure's topological connectivity closely matches the actual physical distribution of rock strata. During network training, a dual optimization equation is constructed between the data-driven loss term and the weight decay penalty factor to dynamically adjust the penalty factor's weight coefficients, preventing network non-convergence caused by physical constraints. Combined with an early stopping mechanism based on physical boundary violation rates, the equivalent elastic modulus transfer values ​​between adjacent graph nodes are monitored during training to ensure whether they exceed the physical material's extreme range, guaranteeing the continuity of the surrounding rock mechanics transfer in the final output tensor data stream. Attached Figure Description

[0017] Figure 1 This is a flowchart illustrating the overall process of the tunnel construction data processing method based on real-time monitoring according to the present invention. Figure 2 This is a flowchart illustrating the construction and dynamic correction process of the spatial graph structure adjacency matrix in this invention. Figure 3 This is a flowchart of the dynamic feature extraction process for graph nodes in this invention; Figure 4 This is a flowchart illustrating the construction of physical constraint terms, dynamic adaptive adjustment, and early stopping mechanism of the present invention. Figure 5 This is a flowchart of the time-gated loop unit feature extraction and anomaly focusing process of the present invention; Figure 6 This is a flowchart illustrating the unified tensor data stream generation and output process of this invention. Detailed Implementation

[0018] The technical solution of the present invention will be further described in detail below with reference to specific embodiments.

[0019] Please refer to the attached document. Figure 1 This embodiment provides displacement sensors, stress sensors, and micro-vibration sensors deployed along the tunnel excavation outline and the depth direction within the surrounding rock within the tunnel construction environment. Each sensor corresponds to a unique geographic coordinate, which is calibrated using a three-dimensional construction coordinate system established with the tunnel excavation portal as the origin. During the continuous tunnel excavation operation, real-time monitoring data collected by each sensor is acquired. The real-time monitoring data collected by the displacement sensors are the axial, radial, and vertical deformation displacements of the surrounding rock monitoring points in the three-dimensional construction coordinate system. The real-time monitoring data collected by the stress sensors are the radial, tangential, and axial stress values ​​at the surrounding rock monitoring points. The real-time monitoring data collected by the micro-vibration sensors are the three-dimensional vibration acceleration time-domain sequence generated by micro-fracture events in the surrounding rock. The real-time monitoring data of each sensor are synchronized and aligned according to a unified time reference to eliminate time deviations in the sampling trigger times of different sensors, providing a time-consistent basic data for subsequent multi-source data fusion processing.

[0020] Please refer to the attached document. Figure 3Using each sensor as an independent graph node, a spatial graph structure is constructed to characterize the spatial topological relationships of the surrounding rock. The node set of the spatial graph structure includes nodes corresponding to all deployed displacement sensors, stress sensors, and micro-vibration sensors. The connection edges between nodes are constructed based on the geographic coordinates of the sensors and the topological connection relationship of the rock strata, generating a corresponding adjacency matrix. The dimension of the adjacency matrix is ​​consistent with the total number of nodes, and the element values ​​in the adjacency matrix correspond to the spatial connection weight between two nodes. The spatial connection weight characterizes the strength of the surrounding rock mechanical transmission correlation between two corresponding sensor monitoring points. While constructing the adjacency matrix, time-series sampling points are used as dynamic features of the graph nodes. Specifically, the real-time monitoring data after synchronization and alignment of each sensor is divided into sequences according to the time step. The time-series sampling point data collected by each sensor within the same time step is mapped to the feature vector of the corresponding graph node, so that the features of the graph node are dynamically updated as the time step progresses, realizing the coupling of spatial topology and time-series features.

[0021] Real-time monitoring data, synchronized and aligned, is input into a pre-constructed physically constrained spatiotemporal graph convolutional network (PCN). This network comprises cascaded graph convolutional propagation layers, time-gated recurrent units, and tensor output layers. During feature propagation in the graph convolutional propagation layers, a displacement-stress transmission constraint term based on the elasticity equation is introduced as a weight decay penalty factor. This penalty factor acts on the learnable weight matrix between adjacent nodes in the graph convolutional propagation layers. During the backpropagation gradient update process, a penalty is applied to weight update directions that deviate from the mechanical transmission laws of the surrounding rock, restricting cross-node information transmission paths to conform to the elasticity transmission laws of the surrounding rock. Specifically, the deformation displacement and stress transmission of the surrounding rock follow a three-dimensional elasticity equilibrium differential equation. This equation defines the intrinsic physical relationship between the displacement field and stress field within the surrounding rock. The displacement-stress transmission constraint term established based on this physical relationship provides a clear physical boundary for information transmission between nodes during the graph convolutional propagation process. This avoids introducing spatially close but mechanically unrelated node features during data-driven weight updates, eliminating the influence of non-physically related interference information on the feature aggregation process.

[0022] The graph convolutional propagation layer outputs a sequence of graph node features containing spatial topological relationships. This sequence is input into a time-gated recurrent unit embedded within a physically constrained spatiotemporal graph convolutional network. The time-gated recurrent unit extracts the dynamic change features of the monitoring data over time. Through its internal gating structure, the time-gated recurrent unit selectively retains and filters feature information from historical time steps. The update gate controls the proportion of hidden states from historical time steps transmitted to the current time step, while the reset gate filters redundant information from historical time steps that is unrelated to the current rock deformation state. This allows the network to focus on the temporal evolution of rock deformation and extract time-dimensional features that characterize the dynamic changes in the rock state.

[0023] The hidden state vector, containing spatiotemporal coupling characteristics, output by the time-gated loop unit is input into the tensor output layer for unified processing. Specifically, the hidden state vector is rearranged according to the node numbering order of the graph, generating a two-dimensional spatiotemporal feature matrix containing spatial node dimensions and time step dimensions. Tensor expansion is performed on the two-dimensional spatiotemporal feature matrix, mapping the channel dimensions, time step dimensions, and spatial node dimensions corresponding to the multi-source sensors to a unified three-dimensional tensor space, generating a unified three-dimensional tensor. Singular value decomposition is performed on the unified three-dimensional tensor along the sensor channel dimensions, truncating suboptimal components smaller than a preset singular value threshold to eliminate spatiotemporal misalignment errors generated during the fusion of multi-source heterogeneous data, and outputting a dimension-reduced unified tensor data stream. The unified tensor data stream contains the spatial topological correlation and temporal dynamic evolution characteristics of the multi-source sensor monitoring data, and all features conform to the physical laws of surrounding rock mechanics transmission, which can be directly used for surrounding rock stability analysis and disaster early warning during tunnel construction.

[0024] Table 1. Mapping Table of Sensor Node Types and Corresponding Feature Dimensions

[0025] Specifically, Table 1 clarifies the feature dimensions and physical meanings corresponding to different types of sensor nodes, providing a unified dimension mapping standard for assigning values ​​to the dynamic features of graph nodes, and avoiding the problem of dimension mismatch in the feature stitching process of multi-source heterogeneous data. Here, T is the number of sampling points within a fixed time window, which corresponds to the sampling time step of the sensor.

[0026] In this embodiment, by constructing a spatial graph structure that couples spatial topology and temporal dynamic features, physical constraint terms based on the equation of elasticity are introduced during the graph convolution propagation process. Temporal dynamic features are extracted by combining time-gated cyclic units, and finally a unified tensor data stream that fuses multi-source heterogeneous data is output. This realizes the spatiotemporal coupling processing of tunnel construction monitoring data, makes the information transmission between nodes conform to the transmission law of surrounding rock mechanics, and eliminates the influence of non-physical interference features on the data fusion process.

[0027] In one optional embodiment, when constructing the adjacency matrix of the spatial graph structure, displacement sensors, stress sensors, and micro-vibration sensors are used as graph nodes in the spatial graph structure. Each graph node corresponds to a unique three-dimensional geographic coordinate system, which is calibrated in the tunnel construction coordinate system. The origin of the coordinate system is the design center point of the tunnel excavation portal. The X-axis points along the tunnel axis towards the excavation direction, the Y-axis points along the tunnel's horizontal radial direction towards the right side of the tunnel, and the Z-axis points vertically upwards. The Euclidean distance between adjacent graph nodes is calculated based on the three-dimensional geographic coordinates corresponding to the two graph nodes. The reciprocal of the Euclidean distance is used as the initial spatial connection weight between the two graph nodes. The initial spatial connection weight is negatively correlated with the spatial distance between the two nodes, indicating that the closer the monitoring points are, the stronger the mechanical transmission correlation.

[0028] A ground-penetrating radar (GPR) profile of the tunnel excavation face is obtained. This profile is acquired point-by-point by the GPR equipment along the contour line of the tunnel excavation face. The horizontal axis of the profile corresponds to the horizontal position of the tunnel excavation face, and the vertical axis corresponds to the depth within the surrounding rock. The grayscale values ​​of the profile correspond to the differences in dielectric constant of the surrounding rock medium. The GPR profile is preprocessed, including grayscale normalization, adaptive filtering for noise reduction, and background signal removal to eliminate interference from environmental noise and system errors during the scanning process. An edge detection algorithm is used to extract the continuous contour lines of the rock strata interfaces from the preprocessed GPR profile. These interfaces correspond to grayscale jumps in dielectric constant. A continuous sequence of rock strata interface coordinates is obtained through non-maximum suppression and dual threshold filtering. Simultaneously, a threshold segmentation algorithm is used to identify the location of fault fracture zones from the preprocessed GPR profile. Fault fracture zones correspond to low-grayscale regions with significantly lower dielectric constants than intact surrounding rock. A region growing algorithm is used to extract the spatial contour and coordinate range of the fault fracture zones.

[0029] When the connection between two graph nodes penetrates the identified rock strata interface or fault fracture zone, the initial spatial connection weights corresponding to the two graph nodes are multiplied by a preset geological weakening coefficient to generate the final adjacency matrix. The geological weakening coefficient is less than 1, representing the blocking effect of the rock strata interface or fault fracture zone on the transmission of stress and deformation in the surrounding rock. The more rock strata interfaces or fault fracture zones that are penetrated, the smaller the geological weakening coefficient, and the lower the spatial connection weight between the corresponding nodes.

[0030] Table 2. Values ​​of Geological Weakening Coefficient for Different Geological Structure Types

[0031] Specifically, Table 2 clarifies the range of values ​​for the geological weakening coefficient and the weight correction rules corresponding to different geological structures, providing a quantitative basis for the construction of the adjacency matrix that conforms to real geological conditions. This enables the topological connection relationship of the spatial graph structure to accurately reflect the mechanical transmission path inside the surrounding rock, avoiding distortion of node connection weight settings caused by geological structure obstruction.

[0032] Real-time monitoring data from displacement sensors, stress sensors, and micro-vibration sensors were standardized using the Z-score standardization method. Based on the mean and standard deviation of historical monitoring data for each sensor, the real-time monitoring data were mapped to a standard distribution interval with a mean of 0 and a variance of 1, generating standard time series for each sensor. Standardization eliminated dimensional and numerical differences between monitoring data from different types of sensors, providing a unified numerical basis for feature fusion of multi-source data. The standard time series were truncated within a fixed time window, containing a preset number of continuous time step sampling points. The length of the time window was set according to the evolution rate of surrounding rock deformation, ensuring that the truncated time series could fully characterize the short-term dynamic changes in the surrounding rock state. The standard time series from each sensor within the same time window were concatenated according to a preset dimensional order to generate high-dimensional feature vectors for the corresponding sensors. The dimensions of the high-dimensional feature vectors corresponded to the number of sampling points within the time window and the dimension of the physical quantities monitored by the sensors. The generated high-dimensional feature vector is assigned to the corresponding graph node as the dynamic feature of the graph node at the current time step. In the next time step, a sliding window is used to re-truncate the standard time series. The sliding window has a sliding step size of a single time step. The high-dimensional feature vector corresponding to the graph node is updated through the sliding window to realize the real-time update of the dynamic features of the graph node.

[0033] After generating the adjacency matrix, it is dynamically corrected based on real-time monitoring data from micro-vibration sensors. Micro-seismic events in the surrounding rock detected by the sensors are extracted in real-time, and their location is processed to obtain their three-dimensional coordinates. The location of the micro-seismic events is calculated based on the time difference of arrival of vibration signals from the multi-channel micro-vibration sensors. A spatial attenuation function is constructed centered on the location coordinates of the micro-seismic events. The value of the spatial attenuation function monotonically decreases with increasing distance from the center of the micro-seismic event, characterizing the attenuation of the influence of the micro-seismic events on the mechanical transmission properties of the surrounding rock with increasing distance. The expression for the spatial attenuation function is:

[0034] in Let the spatial attenuation function be taken at a distance d from the center of the microseismic event. The baseline influence coefficient for the center of microseismic events. denoted as the surrounding rock attenuation coefficient, and d is the Euclidean distance between the graph node and the microseismic event location coordinates.

[0035] The spatial connectivity weights of graph nodes within the influence radius of microseismic events in the adjacency matrix are enhanced by using a spatial attenuation function. The enhanced spatial connectivity weights are the sum of the original weights and the corresponding values ​​of the spatial attenuation function, thus increasing the node connectivity weights within the microseismic event influence area and reflecting the strengthening effect of rock fissure propagation caused by the microseismic event on the mechanical transmission path. When a preset number of rock microseismic events are triggered consecutively within the same spatial area, virtual stress transmission nodes are added to the adjacency matrix. The coordinates of the virtual stress transmission nodes are the average of the location coordinates of all microseismic events within the area. Bidirectional connections are established between the virtual stress transmission nodes and all affected graph nodes within the area. The weights of the bidirectional connections are set based on the reciprocal of the Euclidean distance between the corresponding graph node and the virtual node. By adding virtual nodes, the mechanical transmission path in the rock fissure development area is supplemented, enabling the adjacency matrix to respond in real time to the dynamic changes in the internal structure of the surrounding rock.

[0036] In this embodiment, please refer to the appendix. Figure 2 The initial spatial connectivity weights of the adjacency matrix are corrected by identifying rock strata interfaces and fault fracture zones through ground-penetrating radar cross-sections, ensuring that the topological connections of the spatial graph structure closely match the actual geological distribution of the surrounding rock. Standardization and sliding window methods are used to construct dynamic features of graph nodes, achieving unified feature mapping of multi-source heterogeneous data. Dynamic correction of the adjacency matrix through microseismic events enables the spatial graph structure to respond in real-time to dynamic changes in the internal structure of the surrounding rock, improving the accuracy of the graph structure's representation of the mechanical transmission path of the surrounding rock.

[0037] In another embodiment, the graph convolutional propagation layer of the physically constrained spatiotemporal graph convolutional network adopts a Chebyshev graph convolutional structure. The input of the graph convolutional propagation layer is the dynamic feature matrix and adjacency matrix of the graph nodes, and the output is the graph node feature matrix after spatial topology aggregation. During the feature propagation process of the graph convolutional propagation layer, a displacement-stress transmission constraint term based on the elasticity equation is introduced as a weight attenuation penalty factor. Specifically, a learnable weight matrix between adjacent graph nodes is extracted in the graph convolutional propagation layer. The dimension of the learnable weight matrix is ​​the same as that of the adjacency matrix, and the element values ​​in the learnable weight matrix correspond to the information transfer weight between two adjacent nodes. During the training process of the network, iterative updates are performed through backpropagation gradients.

[0038] Establish an expression for the relationship between displacement and stress partial derivatives based on the three-dimensional elasticity equilibrium differential equations. The three-dimensional elasticity equilibrium differential equations are as follows:

[0039] in These represent the normal stresses of the surrounding rock medium in the x, y, and z directions, respectively. These are the shear stress components of the surrounding rock medium. These are the volume force components of the surrounding rock medium in the x, y, and z directions, respectively. x, y, and z correspond to the three coordinate axes of the tunnel construction coordinate system.

[0040] Based on the generalized Hooke's law, the expression for the partial derivative relationship between displacement and stress is established:

[0041] in These are the stress tensor components. This is the coordinate system direction index, with values ​​1, 2, and 3 corresponding to the x, y, and z directions, respectively. and Let Lame constant be the value of the surrounding rock medium. Let Kronecker function be used when hour ,when hour , Let be the component of the displacement tensor in the i-direction. Let k be the coordinate component in the k-direction.

[0042] The expression for updating node features in the graph convolutional propagation layer is:

[0043] in For the first The node feature matrix output by the layered graph convolutional propagation layer. For the first The input node feature matrix of the layer, It is a non-linear activation function. Let be the learnable weight matrix corresponding to the k-th order Chebyshev polynomial. For k-order normalized Laplace matrix The Chebyshev polynomial, where K is the order of the Chebyshev polynomial. For the normalized graph Laplace matrix, L is the graph Laplace matrix. Let D be the degree matrix of the graph, A be the adjacency matrix, and I be the identity matrix. Let L be the largest eigenvalue of the graph Laplacian matrix L.

[0044] The sum of squared residuals between the updated graph node feature distribution and the expression relating displacement and stress partial derivatives after the learnable weight matrix is ​​updated is given by:

[0045] in Let N be the sum of squared residuals corresponding to the displacement-stress transmission constraint term, and N be the total number of nodes in the graph. The predicted stress tensor components are calculated based on the displacement features output by graph convolution at the i-th node. The measured stress tensor components at the i-th node are collected by the stress sensor. The predicted stress tensor components are calculated from the node displacement characteristics output by graph convolution based on the displacement-stress partial derivative relationship.

[0046] The sum of squared residuals is added as a weight decay penalty factor to the total loss function of the physically constrained spatiotemporal graph convolutional network. The expression for the total loss function is:

[0047] in Let be the total loss function of the network. The data-driven loss term uses mean squared error loss and is calculated based on the predicted surrounding rock deformation value output by the network and the measured monitoring value. The weight coefficient of the weight decay penalty factor. This is the weight decay penalty factor based on the sum of squared residuals.

[0048] Table 3. Range of extreme values ​​of elastic modulus corresponding to common surrounding rock types.

[0049] Specifically, Table 3 clarifies the extreme range of elastic modulus and the range of Lamé constant values ​​for different lithologies of surrounding rocks, providing a basis for the calculation of displacement-stress transmission constraint terms and physical parameters of the surrounding rocks. It also provides a physical extreme value judgment standard for the subsequent calculation of physical boundary violation rate, ensuring that the calculation of physical constraint terms conforms to the actual material properties of the surrounding rocks.

[0050] When the sum of squared residuals is added as a weight decay penalty factor to the total loss function, a dynamic adaptive coefficient adjustment method is used to adjust the weight coefficients. Make real-time adjustments. Specifically, calculate the data-driven loss term in the current training batch. With the sum of squared residuals The variance ratio during gradient backpropagation is calculated based on the variance of the gradient values ​​of the two loss terms across all training samples in the current batch. When the variance ratio exceeds a preset imbalance threshold, a dual optimization equation for the data-driven loss term and the weight decay penalty factor is constructed using the Lagrange multiplier method. The corresponding Lagrange function is:

[0051] in Let Lagrangian function be the function corresponding to the dual optimization equation. For Lagrange multipliers, The gradient variance of the data-driven loss term. denoted as the gradient variance of the weight decay penalty factor, and r is the preset gradient variance balance ratio.

[0052] By solving the dual optimization equation, the optimal weighting coefficients that minimize the Lagrange function are obtained. This reduces the weight decay penalty factor in the current iteration weight coefficient of the total loss function, so that the gradient update magnitudes of the two loss terms reach a balanced state.

[0053] After reducing the weight coefficients of the current iteration by the weight decay penalty factor, an early stopping mechanism based on the physical boundary violation rate is implemented. Specifically, after each network weight coefficient update, the learnable weight matrix in the graph convolutional propagation layer is extracted. Based on the learnable weight matrix and the displacement-stress characteristics of adjacent nodes, the equivalent elastic modulus transfer value between adjacent graph nodes is calculated. The equivalent elastic modulus transfer value is calculated based on the mapping relationship between the displacement-stress relationship and the learnable weight matrix, representing the equivalent elastic modulus of the surrounding rock medium between two adjacent nodes. It is determined whether the calculated equivalent elastic modulus transfer value exceeds the preset extreme range of the physical material elastic modulus. The extreme range of the elastic modulus is obtained from Table 3 based on the surrounding rock lithology type of the corresponding region. The proportion of node pairs exceeding the extreme range of the physical material elastic modulus is counted as the physical boundary violation rate, expressed as:

[0054] in For physical boundary violation rate, This represents the number of node pairs whose equivalent elastic modulus transfer value exceeds the extreme range of the physical material's elastic modulus. This represents the total number of all node pairs in the spatial graph structure that have connecting edges.

[0055] When the physical boundary violation rate shows a continuous upward trend within a preset window period, the model training is forcibly stopped early to lock the current network weight parameters and prevent the network from overfitting the training data and violating the physical boundary conditions.

[0056] In this embodiment, please refer to the appendix. Figure 4 Based on the three-dimensional elasticity equilibrium differential equation and the generalized Hooke's law, a displacement-stress transmission constraint term was constructed and superimposed on the total loss function of the network as a weight decay penalty factor, thereby realizing the physical constraint on the information transmission process of graph convolution. The gradient update magnitude of the data-driven loss term and the physical constraint term was balanced by a dynamic adaptive coefficient adjustment method, avoiding the non-convergence problem in the network training process. The network weight parameters that meet the physical boundary conditions were locked by an early stopping mechanism based on the physical boundary violation rate, ensuring the physical rationality of the network output features.

[0057] In yet another alternative embodiment, please refer to the appendix. Figure 5 The graph convolutional propagation layer of the physically constrained spatiotemporal graph convolutional network outputs a sequence of graph node features containing spatial topological relationships. This sequence is then input to an embedded time-gated recurrent unit (TLU), which extracts the dynamic change features of the monitoring data over time. The graph node feature sequence has the following dimensions: [number of time steps, number of nodes, feature dimension]. The number of time steps corresponds to the number of sampling points within the sliding window, the number of nodes corresponds to the total number of graph nodes in the spatial graph structure, and the feature dimension corresponds to the dimension of the dynamic features of the graph nodes.

[0058] A causal convolution kernel is introduced into the update and reset gates of the time-gated loop unit. This kernel is a one-dimensional causal convolution structure, and its output at the current time step depends only on the input features of the current time step and all previous time steps, without introducing information from future time steps. This aligns with the temporal causality required in real-time monitoring of tunnel construction and avoids interference from future data in the extraction of current state features. The causal convolution kernel is then used to perform convolution operations on the input graph node feature sequence, expanding the local temporal receptive field of the feature sequence and enabling the gating structure to capture the temporal evolution features of the surrounding rock state over a longer time range.

[0059] Table 4. Causal Convolution Kernel Parameters and Corresponding Temporal Receptive Field Range

[0060] Specifically, Table 4 clarifies the temporal sensing field range and applicable scenarios corresponding to different combinations of causal convolution kernel parameters, providing a quantitative basis for setting the causal convolution kernel parameters of the time-gated cyclic unit, enabling the gating structure to adjust the temporal sensing field range according to different monitoring needs and accurately capture the dynamic changes in the surrounding rock state at the corresponding time scale.

[0061] The formula for calculating the update gate of the time-gated loop unit is:

[0062] in The update gate output for the current time step t. It is the sigmoid activation function. To update the learnable weight matrix of the gate, This is a one-dimensional causal convolution operation. The graph node feature matrix is ​​input at the current time step t. For the previous time step The hidden state vector, To update the bias term of the gate.

[0063] The formula for calculating the reset gate of the time-gated loop unit is:

[0064] in The reset gate output for the current time step t. To reset the learnable weight matrix of the gate, The offset term for resetting the door is used; the physical meaning of the other parameters is the same as that of the door update calculation formula.

[0065] The formulas for calculating the candidate hidden state and the hidden state update of the time-gated loop unit are as follows:

[0066] in Let be the candidate hidden state vector at the current time step t. This is the hidden state vector output at the current time step t. The hyperbolic tangent activation function is used. Let be the learnable weight matrix of the candidate hidden states. The bias term for the candidate hidden state. For element-wise multiplication, the physical meaning of the remaining parameters is consistent with the calculation formula of the update gate.

[0067] The update gate controls the retention ratio of hidden states from historical time steps. The closer the update gate's output value is to 1, the more important the feature information of the current time step, and the higher the proportion of candidate hidden states retained for the current time step. Conversely, the closer the update gate's output value is to 0, the more important the hidden states of historical time steps, and the higher the proportion of hidden states retained. The reset gate filters out interfering information from historical time steps unrelated to the current surrounding rock deformation state. The closer the reset gate's output value is to 0, the more irrelevant the corresponding historical time step information is to the current state, and the higher the proportion filtered out. The closer the reset gate's output value is to 1, the higher the correlation between the corresponding historical time step information and the current state, and the higher the proportion retained. After processing by the gating structure, the time-gated cyclic unit outputs a time-hidden state vector containing spatiotemporal coupling characteristics.

[0068] After filtering out historical time-step interference information irrelevant to the current surrounding rock deformation state through a reset gate, a temporal anomaly focusing mechanism is executed. Specifically, the cosine similarity matrix of the graph node feature sequences between the current time step and all historical time steps is calculated. The formula for calculating the cosine similarity is:

[0069] in Let t be the current time step and the historical time step. Cosine similarity between them This is the feature matrix of the graph nodes at the current time step t. For historical time steps The feature matrix of graph nodes, Let f be the Frobenius norm of the matrix.

[0070] The elements in the cosine similarity matrix represent the similarity between the current time step and the historical time steps in terms of the surrounding rock state characteristics. A lower cosine similarity value indicates a sudden change in the surrounding rock state relative to the historical time step, corresponding to a higher probability of sudden stress release or accelerated deformation. The generated cosine similarity matrix is ​​input into a multilayer perceptron (MLP), which contains two cascaded fully connected layers. The first fully connected layer uses the ReLU activation function, and the second fully connected layer uses the sigmoid activation function. The MLP outputs an attention score distribution representing the probability of state abrupt change at each time step. The output attention score distribution is used to element-wise weight the time-hidden state vector output by the update gate. The weighted hidden state vector is the element-wise product of the original hidden state vector and the attention score distribution. This weighting operation amplifies the feature amplitude corresponding to sudden stress release in the time-hidden state vector, enabling the network to focus on the abrupt changes in the surrounding rock state.

[0071] The time-gated cyclic unit outputs the hidden state vectors according to the graph node numbers, generating a two-dimensional spatiotemporal feature matrix. The rows of the two-dimensional spatiotemporal feature matrix correspond to the graph node numbers, the columns correspond to the time steps, and the elements in the matrix are the hidden state feature values ​​of the corresponding nodes at the corresponding time steps. Tensor expansion is performed on the two-dimensional spatiotemporal feature matrix, mapping the multi-source sensor channel dimension, time step dimension, and spatial node dimension to a three-dimensional unified tensor. The three dimensions of the three-dimensional unified tensor are the spatial node dimension, the time step dimension, and the sensor channel dimension, where the sensor channel dimension corresponds to the three data channels of the displacement sensor, stress sensor, and micro-vibration sensor. Singular value decomposition is performed on the three-dimensional unified tensor along the sensor channel dimension. The formula for singular value decomposition is:

[0072] Where M is a two-dimensional matrix obtained by expanding the three-dimensional unified tensor along the sensor channel dimension, and U is a left singular matrix. It is a diagonal singular value matrix, where the elements on the diagonal are singular values ​​arranged in descending order. It is the transpose of the right singular matrix.

[0073] Please refer to the attached document. Figure 6The process truncates components corresponding to second-best singular values ​​in the diagonal singular value matrix that are less than a preset singular value threshold, retaining only the components corresponding to the top K largest singular values, where K is the preset number of principal components. This truncation operation eliminates redundant information and spatiotemporal misalignment errors generated during the fusion of multi-source heterogeneous data, outputting a unified tensor data stream with reduced dimensionality. This unified tensor data stream contains the spatial topological correlation features and temporal dynamic evolution features of the multi-source sensor monitoring data, and all features conform to the physical laws of rock mechanics transmission, making it directly applicable to subsequent rock stability analysis and disaster early warning.

[0074] In this embodiment, by introducing a causal convolution kernel into the gating structure of the time-gated loop unit, the receptive field of the time dimension is expanded, ensuring the causality of temporal feature extraction; the feature amplitude corresponding to the sudden change in the surrounding rock state is amplified through the temporal anomaly focusing mechanism, improving the ability to capture the precursor information of surrounding rock instability; and the spatiotemporal misalignment error in the process of multi-source heterogeneous data fusion is eliminated through tensor expansion and singular value decomposition, outputting a tensor data stream with unified dimensions, and realizing the efficient fusion of multi-source monitoring data.

Claims

1. A method for processing tunnel construction data based on real-time monitoring, characterized in that, include: To acquire real-time monitoring data from displacement sensors, stress sensors, and micro-vibration sensors in the tunnel construction environment; Based on the geographic coordinates of each sensor and the topological connections of the rock strata, an adjacency matrix of the spatial graph structure is constructed, and the time series sampling points are used as the dynamic features of the graph nodes. The real-time monitoring data is input into a pre-constructed physical constraint spatiotemporal graph convolutional network. A displacement-stress transmission constraint term based on the elasticity equation is introduced as a weight attenuation penalty factor in the graph convolution propagation layer of the physical constraint spatiotemporal graph convolutional network to restrict the cross-node information transmission path to conform to the surrounding rock mechanics transmission law. Dynamic features of the time dimension are extracted by the time-gated recurrent unit embedded in the physical constraint spatiotemporal graph convolutional network. Output a unified tensor data stream that integrates heterogeneous data from multiple sources and eliminates spatiotemporal misalignment errors.

2. The tunnel construction data processing method based on real-time monitoring according to claim 1, characterized in that, The process of constructing an adjacency matrix for a spatial graph structure based on the geographic coordinates of each sensor and the topological connection relationship of the rock strata includes: using the displacement sensor, stress sensor, and micro-vibration sensor as graph nodes of the spatial graph structure, respectively. Calculate the Euclidean distance between adjacent graph nodes, and use the reciprocal of the Euclidean distance as the initial spatial connectivity weight; Obtain a ground-penetrating radar scan profile of the tunnel excavation face and identify the locations of rock strata interfaces and fault fracture zones in the ground-penetrating radar scan profile. When the line connecting two graph nodes penetrates the rock strata interface or fault fracture zone, the initial spatial connection weight is multiplied by a preset geological weakening coefficient to generate the adjacency matrix.

3. The tunnel construction data processing method based on real-time monitoring according to claim 1, characterized in that, The dynamic feature of using time series sampling points as graph nodes includes: standardizing the real-time monitoring data of the displacement sensor, stress sensor and micro-vibration sensor respectively to generate standard time series of each sensor; The standard time series is truncated using a fixed time window, and the standard time series of each sensor within the same time window are concatenated into a high-dimensional feature vector; The high-dimensional feature vector is assigned to the corresponding graph node as the dynamic feature of the graph node at the current time step, and the high-dimensional feature vector is updated using a sliding window in the next time step.

4. The tunnel construction data processing method based on real-time monitoring according to claim 1, characterized in that, The method of introducing a displacement-stress transmission constraint term based on the elasticity equation as a weight attenuation penalty factor in the graph convolution propagation layer of the physically constrained spatiotemporal graph convolutional network includes: extracting the learnable weight matrix between adjacent graph nodes in the graph convolution propagation layer. Establish an expression for the relationship between displacement and stress partial derivatives based on the three-dimensional elasticity equilibrium differential equation; Calculate the sum of squared residuals between the updated graph node feature distribution and the expression relating displacement and stress partial derivatives; The sum of squared residuals is added as the weight decay penalty factor to the total loss function of the physically constrained spatiotemporal graph convolutional network.

5. The tunnel construction data processing method based on real-time monitoring according to claim 1, characterized in that, The step of extracting dynamic change features in the time dimension through the time-gated recurrent unit embedded in the physically constrained spatiotemporal graph convolutional network includes: inputting the graph node feature sequence containing spatial topological associations output by the graph convolutional propagation layer into the time-gated recurrent unit; A causal convolution kernel is introduced into the update gate and reset gate of the time-gated loop unit, and the causal convolution kernel is used to expand the local temporal receptive field of the graph node feature sequence. The update gate controls the retention ratio of the hidden state of the historical time step, and the reset gate filters out the historical time step interference information that is irrelevant to the current surrounding rock deformation state, outputting a time hidden state vector containing spatiotemporal coupling characteristics.

6. The tunnel construction data processing method based on real-time monitoring according to claim 1, characterized in that, The output is a unified tensor data stream that integrates multi-source heterogeneous data and eliminates spatiotemporal misalignment errors, including: rearranging the time hidden state vector output by the time-gated loop unit according to the graph node number to generate a two-dimensional spatiotemporal feature matrix; Tensor expansion is performed on the two-dimensional spatiotemporal feature matrix to map the multi-source sensor channel dimension, time step dimension, and spatial node dimension into a three-dimensional unified tensor. Singular value decomposition is performed on the three-dimensional unified tensor along the sensor channel dimension, suboptimal components smaller than a preset singular value threshold are truncated, and the dimension-reduced unified tensor data stream is output.

7. The tunnel construction data processing method based on real-time monitoring according to claim 2, characterized in that, After generating the adjacency matrix, a dynamic correction process for the adjacency matrix is ​​also included: real-time extraction of the location coordinates of surrounding rock micro-seismic events monitored by micro-vibration sensors; A spatial attenuation function is constructed centered on the location coordinates of the microseismic events in the surrounding rock. The spatial connectivity weights of graph nodes within the influence radius of microseismic events in the adjacency matrix are enhanced and superimposed using the spatial attenuation function. When a preset number of microseismic events in the surrounding rock are triggered consecutively in the same area, a virtual stress transmission node is added to the adjacency matrix, and a bidirectional connection edge is established between the virtual stress transmission node and the affected graph node.

8. The tunnel construction data processing method based on real-time monitoring according to claim 4, characterized in that, When the residual sum of squares is superimposed as the weight decay penalty factor into the total loss function of the physically constrained spatiotemporal graph convolutional network, a dynamic adaptive coefficient adjustment method is adopted: the variance ratio of the data-driven loss term in the current training batch to the residual sum of squares during gradient backpropagation is calculated. When the variance ratio is greater than the preset imbalance threshold, a dual optimization equation for the data-driven loss term and the weight decay penalty factor is constructed using the Lagrange multiplier method. By solving the dual optimization equation, the weight decay penalty factor in the current iteration weight coefficient of the total loss function is reduced.

9. The tunnel construction data processing method based on real-time monitoring according to claim 5, characterized in that, After filtering out historical time step interference information that is irrelevant to the current surrounding rock deformation state through the reset gate, a time anomaly focusing mechanism is also included: calculating the cosine similarity matrix of the graph node feature sequence between the current time step and the historical time step; The cosine similarity matrix is ​​input into a multilayer perceptron, and the attention score distribution representing the probability of state change at each time step is output. The attention score distribution is used to weight the time hidden state vector output by the update gate element by element, thereby amplifying the characteristic amplitude in the time hidden state vector corresponding to the sudden release of surrounding rock stress.

10. The tunnel construction data processing method based on real-time monitoring according to claim 8, characterized in that, After reducing the weight decay penalty factor in the current iteration weight coefficient of the total loss function, an early stopping mechanism based on the physical boundary violation rate is also included: after each weight coefficient update, the learnable weight matrix is ​​extracted and the equivalent elastic modulus transfer value between adjacent graph nodes is calculated. Determine whether the equivalent elastic modulus transfer value exceeds the preset extreme value range of the physical material elastic modulus, and count the proportion of the number of node pairs that exceed the extreme value range of the physical material elastic modulus to the total number of node pairs as the physical boundary violation rate. When the physical boundary violation rate shows an upward trend for a continuous preset window period, the model training is forcibly stopped early and the current network weight parameters are locked.

Images