Method for predicting tunnel excavation water inflow through horizontal drilling

Abstract

The invention provides a method for predicting tunnel excavation water inflow through horizontal drilling, and relates to the technical field of tunnel construction water inflow prediction.The method comprises the steps that water inflow time sequence data and geological parameter data are obtained; calculating the permeability coefficient of the position where each monitoring point is located according to the water inflow time sequence data by using a Goodman formula based on rock integrity correction; decomposing the water inflow time sequence data to obtain a plurality of intrinsic mode functions and a residual term; forming a training data set; constructing a prediction model comprising a GRU neural network structure and a physical constraint mechanism, training the prediction model by using the training data set until the loss function converges, and obtaining a trained prediction model; and obtaining current monitoring data, inputting the current monitoring data into the trained prediction model, and predicting the water inflow and the water inflow risk level of the to-be-excavated section. According to the method, refined and high-precision prediction of the water inflow in the tunnel excavation process can be achieved, and potential safety hazards caused by coarse-grained prediction of a traditional method are avoided.

Description

Technical Field

[0001] The present invention relates to the technical field of water inflow prediction for tunnel construction, and in particular to a method for predicting water inflow during tunnel excavation through horizontal drilling. Background Art

[0002] Tunnel water inrush is one of the major geological hazards threatening the safety of tunnel projects and has a significant impact on construction safety. Sudden water inrush may not only cause casualties among construction workers and damage machinery and equipment, but also lead to a series of chain reactions, such as project delays and cost overruns. Existing methods for predicting tunnel water inrush mainly include empirical formulas, numerical simulations, and statistical analysis. These methods generally rely on vertical drilling and geophysical exploration to obtain geological parameters. However, the empirical formula method relies too much on empirical parameters and cannot accurately reflect the dynamic changes in hydrology under complex geological conditions. Although the numerical simulation method has a solid theoretical foundation, it is computationally intensive and difficult to obtain key parameters, limiting its practicality. The statistical analysis method, however, has low prediction accuracy and is difficult to achieve effective early warning of sudden water inrush events, thus exhibiting significant technical limitations.

[0003] Chinese invention patent 202111449519.4 discloses a method for predicting the amount of water inflow during tunnel excavation based on horizontal directional drilling exploration holes, attempting to solve some of the defects of traditional prediction methods. However, this method still has technical bottlenecks: it uses a formula method to predict the amount of water inflow in the tunnel in sections, and uses the distance between vertical boreholes as the basis for segmentation, but the vertical borehole spacing is usually 200-500 meters, or even larger. Considering that tunnel construction generally adopts a combination of drilling and blasting and TBM construction, and the advance of a single excavation cycle is about 2-3 meters, such a long prediction segment is bound to lead to the prediction result being too coarse and insufficient in accuracy, and it is impossible to timely warn of local water inflow risks, thereby posing potential hidden dangers to construction safety.

[0004] In summary, existing technologies suffer from three major technical deficiencies: First, traditional vertical borehole monitoring methods struggle to continuously and accurately reflect the changing characteristics of the tunnel's axial hydrogeological conditions; second, most prediction models rely on empirical formulas, lacking the support of real-time dynamic monitoring data and unable to adapt to complex changes in geological conditions; and third, prediction results are often too macroscopic and holistic, making it difficult to achieve high-precision localized segmental predictions, making them difficult to meet the practical requirements of tunnel construction safety. Therefore, there is an urgent need to develop a new tunnel water inflow prediction method that can achieve refined predictions. Summary of the Invention

[0005] In view of this, the present invention proposes a method for predicting water inflow during tunnel excavation by horizontal drilling, which realizes the refined and high-precision prediction of water inflow during tunnel excavation, avoids the safety hazards caused by coarse-grained prediction of traditional methods, provides forward-looking warning for water inflow risks during construction, and effectively improves the safety and efficiency of tunnel construction.

[0006] The technical solution of the present invention is achieved as follows:

[0007] The present invention provides a method for predicting water inflow during tunnel excavation by horizontal drilling, comprising:

[0008] S1. Horizontal directional drilling is arranged in sections along the tunnel axis and a monitoring system is installed to obtain time series data on water inflow and geological parameters at each monitoring point;

[0009] S2. Using the Goodman formula modified based on rock mass integrity, calculate the permeability coefficient at each monitoring point based on the water inflow time series data;

[0010] S3, using the complete set empirical mode decomposition method to decompose the water inflow time series data to obtain multiple intrinsic mode functions and a residual term;

[0011] S4. forming a training data set based on the intrinsic mode function, residual term, permeability coefficient and geological parameter data;

[0012] S5. Construct a prediction model that includes a GRU neural network structure and a physical constraint mechanism, and train the prediction model using a training data set until the loss function converges to obtain a trained prediction model;

[0013] S6. Obtain current monitoring data and input it into the trained prediction model to predict the water inflow and water inflow risk level of the section to be excavated.

[0014] Preferably, in step S1, arranging horizontal directional drilling in sections and installing a monitoring system includes:

[0015] Use geological radar to conduct advance geological surveys ahead of the tunnel to identify potential high-risk water inrush areas;

[0016] A monitoring point should be set up every 15-20 meters of drilling in ordinary areas, and every 5-10 meters of drilling in high-risk areas;

[0017] High-precision flow sensors and data transmission devices are deployed at each monitoring point.

[0018] Preferably, the geological parameter data include: groundwater level, borehole diameter, water-bearing body length, rock mass quality index RQD, fault fracture zone width, karst development degree and environmental influencing factors.

[0019] Preferably, step S3 includes:

[0020] S31, adding white noise of different amplitudes to the water inflow time series data to generate multiple sets of data;

[0021] S32, performing empirical mode decomposition on each set of noise-added data;

[0022] S33. Average the decomposition results of all groups to obtain the intrinsic mode function and residual term.

[0023] Preferably, the amplitude range of the white noise in step S31 is 0.1-0.3 times the standard deviation of the original data, and the number of noise groups is 50-100 groups.

[0024] Preferably, the Goodman formula corrected based on rock mass integrity in step S2 is:

[0025]

[0026] Where K is the permeability coefficient; Q is the time series data of water inflow; L is the length of the horizontal borehole through the water-bearing body; H is the head difference, which represents the vertical height difference between the groundwater level and the tunnel drainage surface; d is the cross-sectional diameter of the horizontal borehole; α(RQD) is the correction factor based on rock mass integrity, which is calculated based on the rock mass quality index RQD, as shown in the following formula:

[0027]

[0028] Where C1 is the basic correction value, and its value range is 0.15-0.25; C2 is the gain coefficient, and its value range is 0.75-0.85; k is the curve steepness coefficient, which controls the sensitivity of the correction coefficient to changes in RQD, and its value range is 0.04-0.06; R0 is the inflection point parameter, which indicates the critical value of the correction coefficient change, and its value range is 45-55.

[0029] Preferably, the prediction model includes:

[0030] Feature preprocessing layer, which combines the intrinsic mode function and residual term to form time series features, and combines the permeability coefficient and geological parameters to form static features;

[0031] The feature encoding layer extracts temporal features and static features, and concatenates them in the feature dimension to form a fused feature vector;

[0032] The temporal modeling layer includes multiple stacked bidirectional GRU units to capture temporal dependencies. The output of each GRU unit is weighted by the attention mechanism to obtain the context feature representation.

[0033] The prediction layer maps the context feature representation to the water inflow prediction value through a fully connected network;

[0034] At the risk assessment layer, the ratio of the predicted water inflow volume to the safety threshold is calculated, and the water inflow risk level is divided according to the size of the ratio.

[0035] Preferably, the loss function is:

[0036] L total =L MSE +λ1(N)·L physical +λ2(N)·L temporal

[0037]

[0038] Where, L total is the total loss function; L MSE is the mean square error loss; L physical is the physical constraint loss; L temporal is the time continuity loss; m is the number of samples; y i is the actual water inflow of the i-th sample; is the predicted water inflow of the i-th sample; A tunnel is the tunnel cross-sectional area; K i is the permeability coefficient of the monitoring point corresponding to the i-th sample; is the hydraulic gradient of the monitoring point corresponding to the i-th sample, indicating the change in hydraulic head per unit distance; Δt i represents the time interval between two samples i and i+1; λ1(N) and λ2(N) are weight coefficients that change with the accumulated data volume N.

[0039] Preferably, step S4 includes:

[0040] S41. Taking the intrinsic mode function and the residual term as time series features, and the permeability coefficient and geological parameter data as static features;

[0041] S42. For each monitoring point, select the time series features and corresponding static features of T consecutive time points to form an input feature vector; and use the water inflow at time point T+1 as the output prediction target;

[0042] S43, using a sliding window method with a window length of M and a sliding step size of 1, sequentially construct multiple input-output sample pairs;

[0043] S44. Randomly divide the sample pairs into a training set, a validation set, and a test set to form a training data set.

[0044] Preferably, the determination of the water inrush risk level in step S6 includes:

[0045] Calculate the ratio of the predicted water inflow to the safety threshold to obtain the risk index P;

[0046] The risk level is divided according to the risk index P: P < 0.3 is low risk, 0.3 ≤ P < 0.6 is medium risk, 0.6 ≤ P < 0.8 is high risk, and P ≥ 0.8 is extremely high risk.

[0047] The present invention has the following beneficial effects compared to the prior art:

[0048] (1) The present invention achieves high-precision prediction of tunnel excavation water inflow by arranging horizontal directional boreholes in sections along the tunnel axis and installing a monitoring system. It combines the complete set empirical mode decomposition method to process the time series data of water inflow, and adopts a prediction model that integrates the GRU neural network and the physical constraint mechanism. Compared with traditional vertical drilling and empirical formula prediction methods, the present invention can capture the continuous changes in the hydrogeological conditions along the tunnel axis, provide more refined prediction results, and provide a forward-looking warning for water inflow risks during tunnel construction.

[0049] (2) The present invention adopts a segmented horizontal directional drilling strategy, setting a monitoring point every 15-20 meters in ordinary areas and every 5-10 meters in high-risk areas, breaking through the 200-500 meter spacing limit of traditional vertical drilling. This dense layout makes the distribution of monitoring points more reasonable, accurately reflecting the spatial variation characteristics of hydrogeological conditions along the tunnel, and the monitoring granularity is more closely matched with the single-cycle footage of tunnel construction, effectively improving the prediction model's ability to identify local small-scale water inrush risks;

[0050] (3) The present invention introduces the Goodman formula based on rock mass integrity correction to calculate the permeability coefficient, and uses the α(RQD) correction coefficient to correct the permeability characteristics under different rock mass quality grades. This correction coefficient is constructed using a Sigmoid function and can accurately describe the nonlinear characteristics of the permeability coefficient changing with rock mass integrity, effectively solving the problem of poor adaptability of traditional empirical formulas to rock masses with different integrity levels.

[0051] (4) The present invention uses the complete set empirical mode decomposition method to process the water inflow time series data. By adding white noise with an amplitude of 0.1-0.3 times the standard deviation to the original data and generating 50-100 groups of data, each group of data is subjected to empirical mode decomposition and averaged, effectively overcoming the modal aliasing and endpoint effect problems existing in the traditional EMD method;

[0052] (5) The prediction model constructed by the present invention adopts a hierarchical network structure. The feature preprocessing layer extracts temporal features and static features respectively, and combines the feature encoding layer and the temporal modeling layer to achieve effective fusion of multi-source heterogeneous data. In particular, the temporal modeling layer adopts a bidirectional GRU structure with an attention mechanism, which can simultaneously capture the forward and backward dependencies of the water flow time series data and adaptively assign different weights to features at different time points, effectively suppressing noise interference while retaining long-term memory. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments or the description of the prior art. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative work.

[0054] Figure 1 is a flow chart of the method of the present invention;

[0055] Figure 2 This is a technical implementation diagram of the present invention. DETAILED DESCRIPTION

[0056] The following will be combined with the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the embodiments described are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of the present invention.

[0057] like Figure 1 As shown, the present invention provides a method for predicting water inflow during tunnel excavation by horizontal drilling, comprising:

[0058] S1. Horizontal directional drilling is arranged in sections along the tunnel axis and a monitoring system is installed to obtain time series data on water inflow and geological parameters at each monitoring point;

[0059] S2. Using the Goodman formula modified based on rock mass integrity, calculate the permeability coefficient at each monitoring point based on the water inflow time series data;

[0060] S3, using the complete set empirical mode decomposition method to decompose the water inflow time series data to obtain multiple intrinsic mode functions and a residual term;

[0061] S4. forming a training data set based on the intrinsic mode function, residual term, permeability coefficient and geological parameter data;

[0062] S5. Construct a prediction model that includes a GRU neural network structure and a physical constraint mechanism, and train the prediction model using a training data set until the loss function converges to obtain a trained prediction model;

[0063] S6. Obtain current monitoring data and input it into the trained prediction model to predict the water inflow and water inflow risk level of the section to be excavated.

[0064] like Figure 2 As shown in the figure, the technical implementation process of the present invention is as follows: first, horizontal directional drilling is arranged in sections along the tunnel axis and a monitoring system is installed. Monitoring points are set up every 15-20 meters in ordinary areas and every 5-10 meters in high-risk areas to obtain time series data of water inflow and geological parameter data; secondly, the Goodman formula based on rock integrity correction is used to calculate the permeability coefficient, and the correction coefficient α(RQD) constructed by the Sigmoid function is introduced to fine-tune the rock masses with different integrity; then, the complete set empirical mode decomposition method is used to process the time series data of water inflow, and multiple groups of data are generated by adding white noise and averaging them to overcome the traditional EM The modal aliasing problem of D is solved to obtain stable intrinsic mode functions and residual terms; then, the intrinsic mode functions and residual terms are used as time series features, and the permeability coefficient and geological parameters are used as static features to construct a training data set; then, a prediction model integrating the GRU neural network structure and the physical constraint mechanism is designed, and a composite loss function including mean square error loss, physical constraint loss based on Darcy's law, and time continuity loss is used for training; finally, the current monitoring data is input into the trained model to predict the water inflow of the section to be excavated, and the risk level is divided according to the ratio of the prediction result to the safety threshold, providing accurate early warning for tunnel construction.

[0065] Specifically, in one embodiment of the present invention, in step S1, arranging horizontal directional drilling in sections and installing a monitoring system includes:

[0066] Use geological radar to conduct advance geological surveys ahead of the tunnel to identify potential high-risk water inrush areas;

[0067] A monitoring point should be set up every 15-20 meters of drilling in ordinary areas, and every 5-10 meters of drilling in high-risk areas;

[0068] High-precision flow sensors and data transmission devices are deployed at each monitoring point.

[0069] In this embodiment, geological radar is first used to conduct advance surveys of the geology ahead of the tunnel to identify potential high-risk water inflow areas, such as fault fracture zones and karst development areas. By emitting electromagnetic waves and receiving reflected signals, geological radar can effectively detect geological structural anomalies ahead of the tunnel, providing guidance for subsequent drilling layout. After obtaining the geological radar detection results, horizontal directional boreholes are laid out along the tunnel axis. The boreholes have a certain slope to allow for natural groundwater outflow, facilitating water inflow monitoring. A differentiated drilling strategy is employed: in standard areas, a monitoring point is set every 15-20 meters of drilling progress; in high-risk areas (such as fault fracture zones and karst development areas identified by geological radar), a monitoring point is set every 5-10 meters of drilling progress, forming a rational monitoring network layout. For hard rock formations, drilling can be carried out simultaneously with TBM (Transportation Machine) excavation machinery to improve construction efficiency. At each monitoring point, a high-precision flow sensor and data transmission device are deployed. The flow sensor accurately measures the amount of water inflow from the borehole, and the data transmission device transmits the monitoring data in real time to the data processing center. The monitoring system uses wireless transmission, and each monitoring point is equipped with an RFID tag for easy positioning and identification.

[0070] Specifically, the depth and flow information of each monitoring point is recorded simultaneously, and geological parameter data is collected, including:

[0071] The groundwater level is obtained by a water level gauge in the borehole, which represents the vertical height difference between the groundwater level and the tunnel drainage surface.

[0072] The diameter of the borehole is the cross-sectional diameter of the horizontal directional borehole, which is determined by the specifications of the drill bit used during drilling.

[0073] The length of the water body is the length of the horizontal borehole through the aquifer, which is obtained through geological records during the drilling process and electrical logging technology.

[0074] The rock quality index (RQD) is obtained through core sampling and analysis of drilled rocks. It is an indicator of the integrity of the rock mass. The calculation formula is the percentage of the cumulative length of the core sections greater than 10 cm to the total drilling length.

[0075] The width of the fault fracture zone is determined by geological radar detection and the fault fracture zone conditions encountered during drilling, and the width of the fault fracture zone is recorded.

[0076] The degree of karst development can be divided into three levels: weak, moderate and strong through comprehensive evaluation such as core analysis, geological radar detection and hydrogeological survey during drilling.

[0077] Environmental influencing factors, including rainfall, atmospheric pressure changes, surrounding construction activities and other external factors that may affect the amount of water inflow, are obtained through meteorological stations and environmental monitoring equipment.

[0078] At the same time, the time series data of water inflow is collected, that is, the numerical sequence of water inflow measured at different time points at each monitoring point.

[0079] All data are marked with timestamps to establish a unified spatiotemporal database.

[0080] Specifically, in one embodiment of the present invention, step S2 calculates the permeability coefficient of each section based on the water inflow of each section of the borehole, and uses the improved Goodman formula for derivation:

[0081]

[0082] Where K is the permeability coefficient; Q is the time series data of water inflow; L is the length of the horizontal borehole through the water-bearing body; H is the head difference, which represents the vertical height difference between the groundwater level and the tunnel drainage surface; d is the cross-sectional diameter of the horizontal borehole; α(RQD) is the correction factor based on rock mass integrity, which is calculated based on the rock mass quality index RQD, as shown in the following formula:

[0083]

[0084] Where C1 is the basic correction value, and its value range is 0.15-0.25; C2 is the gain coefficient, and its value range is 0.75-0.85; k is the curve steepness coefficient, which controls the sensitivity of the correction coefficient to changes in RQD, and its value range is 0.04-0.06; R0 is the inflection point parameter, which indicates the critical value of the correction coefficient change, and its value range is 45-55.

[0085] This embodiment introduces a correction factor α(RQD) constructed using a sigmoid function to accurately describe the nonlinear characteristics of the permeability coefficient as it varies with rock mass integrity. A low RQD value indicates high rock fragmentation, with the correction factor approaching C1+C2 and a high permeability coefficient. A high RQD value indicates high rock mass integrity, with the correction factor approaching C1 and a low permeability coefficient. This correction effectively addresses the adaptability issues of the traditional Goodman formula for rock masses of varying integrity.

[0086] In the actual calculation process, parameters such as the water inflow time series data Q, water-bearing body length L, head difference H, borehole diameter d, and rock mass quality index RQD are first extracted from the monitoring data obtained by S1. The correction factor α(RQD) is then calculated based on the rock mass quality index RQD. Finally, the modified Goodman formula is substituted into the calculated permeability coefficient K at each monitoring point.

[0087] Specifically, in one embodiment of the present invention, step S3 includes:

[0088] S31. Add white noise of different amplitudes to the water inflow time series data to generate multiple sets of data.

[0089] Specifically, white noise of varying amplitudes is added to the original water inflow time series data X(t) to generate multiple sets of noise superposition data. The amplitude of the white noise ranges from 0.1 to 0.3 times the standard deviation of the original data, and the number of noise groups is 50 to 100. The formula is as follows:

[0090] X j (t) = X(t) + ε j ·σ·N j (t),j=1,2,...,M

[0091] Where, X j (t) is the data after adding noise to the jth group; X(t) is the original water inflow time series data; ε j is the amplitude coefficient, ranging from 0.1 to 0.3; σ is the standard deviation of the original data; N j (t) is the jth standard white noise sequence (mean is 0, standard deviation is 1); M is the number of noise groups, ranging from 50 to 100 groups.

[0092] S32. Perform empirical mode decomposition on each set of noise-added data.

[0093] Specifically, for each set of data X after adding noise j (t) Perform traditional empirical mode decomposition and decompose it into a series of intrinsic mode functions (IMFs) and a residual term.

[0094] The basic idea of empirical mode decomposition is to decompose a complex time series signal into a series of intrinsic mode functions (IMFs) with different characteristic scales. Each IMF represents an oscillation mode in the original signal. The decomposition process is implemented through an iterative "screening" process, which consists of the following steps:

[0095] 1. Identify X j All local maxima and minima in (t);

[0096] 2. Fit the upper envelope through the maximum point and the lower envelope through the minimum point;

[0097] 3. Calculate the average value of the upper and lower envelopes, denoted as m(t);

[0098] 4. Calculate X j The difference between (t) and m(t) gives h(t) = X j (t)-m(t);

[0099] 5. Determine whether h(t) satisfies the two IMF conditions: (a) the difference between the number of maximum and minimum points does not exceed 1; (b) the local mean is zero;

[0100] 6. If h(t) does not satisfy the IMF condition, treat h(t) as a new signal and repeat steps 1-5 until the IMF condition is satisfied.

[0101] 7. Define h(t) that meets the conditions as the first IMF component, denoted as IMF1;

[0102] 8. Calculate the residual r(t) = X j (t)-IMF1;

[0103] 9. Use r(t) as the new input signal and repeat the above steps to obtain IMF2, IMF3, ..., until the residual function becomes a monotonic function or its amplitude is less than a preset threshold.

[0104] For each set of data X j (t), after the above decomposition, we get:

[0105]

[0106] Where, is the i-th eigenmode function obtained by decomposing the j-th group of data; r j (t) is the residual term of the jth group of data; N j is the number of intrinsic mode functions obtained by decomposing the j-th group of data.

[0107] S33. Average the decomposition results of all groups to obtain the intrinsic mode function and residual term.

[0108] Specifically:

[0109]

[0110] Where, IMF i (t) is the final i-th eigenmode function; r(t) is the final residual term; N is the number of eigenmode functions, which ranges from 5 to 7.

[0111] Specifically, in one embodiment of the present invention, step S4 includes:

[0112] S41, take the intrinsic mode function and residual term as time series features, and take the permeability coefficient and geological parameter data as static features; the time series features usually include 5-7 intrinsic mode functions (IMF1, IMF2, ..., IMF n) and a residual term (r), each of which is a time-varying sequence. These features reflect the variation of water inflow at different time scales. The high-frequency IMF reflects short-term fluctuations, the low-frequency IMF reflects medium-term changes, and the residual term reflects long-term trends. Static features include the permeability coefficient K and geological parameter data (groundwater level, borehole diameter, water-bearing body length, rock mass quality index (RQD), fault fracture zone width, karst development level, and environmental factors). These features are relatively stable at each monitoring point and reflect the spatial geohydrological characteristics.

[0113] S42. For each monitoring point, select the time series features and corresponding static features of T consecutive time points to form an input feature vector; take the water inflow at time point T+1 as the output prediction target; specifically expressed as:

[0114] Input feature vector: X i =[IMF1(t i-T+1 :t i ),IMF2(t i-T+1 :t i ),…,IMF n (t i-T+1 :t i ),r(t i-T+1 :t i ),K i ,G i ];

[0115] Output prediction target: y i =Q(t i+1 ).

[0116] Among them, the IMF j (t i-T+1 :t i ) represents the jth eigenmode function at t i-T+1 to t i The value at a time point; r(t i-T+1 :t i ) indicates the residual term at t i-T+1 to t i The value at the time point; K i represents the permeability coefficient of the i-th monitoring point; G i represents the geological parameter data of the i-th monitoring point; Q(t i+1 ) represents t i+1 The actual water inflow at a point in time.

[0117] S43. Use a sliding window method with a window length of M and a sliding step of 1 to sequentially construct multiple input-output sample pairs. The sliding window process is as follows:

[0118] First, the window is placed at the beginning of the time series to construct the first sample pair. Then, the window is slid forward by one unit to construct the second sample pair. This process is repeated until the window reaches the end of the time series. This method fully utilizes time series data and generates a large number of training samples from limited monitoring data. For example, if a monitoring point has data for 1500 time points, T = 24, and M = 1440, then 1500 - 24 - 1 + 1 = 1476 sample pairs can be constructed. Each sample pair contains the input features of 24 time points and the output target of 1 time point.

[0119] S44. Randomly divide the sample pairs into a training set, a validation set, and a test set to form a training data set.

[0120] The division ratio is: training set 70%, validation set 15%, and test set 15%.

[0121] Specifically, in one embodiment of the present invention, the prediction model includes:

[0122] Feature preprocessing layer, the intrinsic mode function and the residual term are combined to form the time series feature, and the permeability coefficient and geological parameters are combined to form the static feature; specifically, the multiple intrinsic mode functions (IMF1, IMF2, ..., IMF n ) and a residual term (r) are combined to form a time series feature matrix with a dimension of T×(n+1), where T is the length of the time window and n+1 is the number of features. The calculated permeability coefficient K and the obtained geological parameter data (groundwater level, borehole diameter, water-bearing body length, rock quality index RQD, fault fracture zone width, karst development degree, and environmental influencing factors) are combined to form a static feature vector. All features are standardized to eliminate dimensionality effects.

[0123] The feature encoding layer extracts temporal features and static features, and concatenates them in the feature dimension to form a fused feature vector.

[0124] The feature encoding layer encodes different types of features independently and then fuses them into a unified representation:

[0125] For temporal features, a one-dimensional convolutional network (Conv1D) is used for feature extraction with a convolution kernel size of 3 and a step size of 1 to extract local temporal patterns. For static features, a fully connected network is used for feature mapping to capture the nonlinear relationship between features. The processed temporal features and static features are concatenated in the feature dimension to form a fused feature vector, which is used as the input of the temporal modeling layer.

[0126] The temporal modeling layer includes multiple stacked bidirectional GRU units to capture temporal dependencies. The output of each GRU unit is weighted by the attention mechanism to obtain contextual feature representation.

[0127] Specifically, the time series modeling layer contains 3-4 stacked bidirectional GRU units, each of which contains 64-128 hidden neurons. The bidirectional GRU can consider both past and future information, enhancing its ability to capture time series patterns. The calculation process of each GRU unit is as follows:

[0128] z t =σ(W z ·[h t-1 ,x t ]+b z )

[0129] r t =σ(W r ·[h t-1 ,x t ]+b r )

[0130]

[0131] Among them, z t is the update gate, which controls the degree of retention of the previous state information; r t To reset the gate, control the influence of the previous state on the current candidate state; is the candidate hidden state; h t is the current hidden state; x t is the input feature, including hydrological parameters, geological parameters, etc.; W z 、W r , W is the network weight matrix; b z 、b r , b is the bias term.

[0132] The output of each GRU unit is weighted by the attention mechanism to highlight the impact of important time points:

[0133] Attention weight calculation: α t =softmax(v T tanh(W a ·h t +b a ));

[0134] Weighted context vector: c = ∑α t ·h t ;

[0135] Finally, a contextual feature representation containing rich temporal information is obtained.

[0136] The prediction layer maps the context feature representation to the water inflow prediction value through a fully connected network.

[0137] The prediction layer maps the context feature representation to the water inflow prediction value:

[0138] It contains 2-3 fully connected layers, and the number of hidden layer neurons decreases successively (such as 128→64→32); the activation function uses ReLU to avoid the gradient vanishing problem; the latter layer is the output layer of a single neuron, which represents the predicted water inflow value; the output layer does not use an activation function and directly outputs the linear combination result.

[0139] At the risk assessment layer, the ratio of the predicted water inflow volume to the safety threshold is calculated, and the water inflow risk level is divided according to the size of the ratio.

[0140] The risk assessment layer calculates the ratio of the predicted water inflow value to the safety threshold and divides the risk level accordingly:

[0141] The safety threshold is determined based on tunnel construction specifications and site conditions; the risk index P = predicted water inflow / safety threshold; the risk level classification is: P < 0.3 is low risk, 0.3 ≤ P < 0.6 is medium risk, 0.6 ≤ P < 0.8 is high risk, and P ≥ 0.8 is extremely high risk.

[0142] After the prediction model is built, it is trained using the training data set. The loss function during training introduces a physical constraint mechanism so that the model prediction results conform to both the data rules and the laws of physics. The loss function is:

[0143] L total =L MSE +λ1(N)·L physical +λ2(N)·L temporal

[0144]

[0145] Where, L total is the total loss function; L MSE is the mean square error loss; L physical is the physical constraint loss; L temporal is the time continuity loss; m is the number of samples; y i is the actual water inflow of the i-th sample; is the predicted water inflow of the i-th sample; A tunnel is the tunnel cross-sectional area; K i is the permeability coefficient of the monitoring point corresponding to the i-th sample; is the hydraulic gradient of the monitoring point corresponding to the i-th sample, indicating the change in hydraulic head per unit distance; Δt i represents the time interval between two samples i and i+1; λ1(N) and λ2(N) are weight coefficients that change with the accumulated data volume N. The weight coefficients are calculated as follows:

[0146] λ i (N) = λ i0·e -γN

[0147] Among them, λ i0 is the initial weight; γ is the decay coefficient, and N is the amount of accumulated data. When the amount of data is small, the physical constraint weight is large. As the data accumulates, the constraint weight gradually decreases, enhancing the model's ability to fit the measured data.

[0148] The specific training process of the model is as follows:

[0149] 1. Initialization: Randomly initialize model parameters and set hyperparameters such as learning rate, batch size, and number of training rounds;

[0150] 2. Forward propagation: Input data is fed into the model and the prediction results are obtained through calculations at each layer;

[0151] 3. Loss calculation: Calculate the loss value under the current parameters according to the above composite loss function;

[0152] 4. Backpropagation: Calculate the gradient of the loss function with respect to each parameter and use the optimizer to update the parameters; the optimizer uses the Adam algorithm, with the initial learning rate set to 0.001; a learning rate decay strategy is used, with the learning rate decaying to 0.9 times the original value every 10 training rounds;

[0153] 5. Early stopping strategy: Monitor model performance on the validation set and stop training early if the validation set loss does not decrease for 10 consecutive rounds;

[0154] 6. Model saving: Save the model parameters with the best performance on the validation set.

[0155] During training, the overall loss function and its three components are monitored to ensure that the model learns both the statistical patterns in the data and the laws of hydrophysics. Training ends when the loss function converges to a stable value or reaches the preset maximum number of training rounds, resulting in a trained prediction model.

[0156] Specifically, in one embodiment of the present invention, step S6 includes:

[0157] First, the latest water inflow time series data and geological parameter data are collected in real time from the horizontal directional drilling monitoring system. Then, the complete set empirical mode decomposition (CEMD) of the water inflow time series data is performed to obtain the intrinsic mode function and residual term. Then, combined with the calculated permeability coefficient, these features are organized into an input vector according to the training format. The input vector is then fed into the trained GRU neural network prediction model to obtain the predicted water inflow value at a future time point. The ratio of the predicted water inflow to the safety threshold is calculated to obtain the risk index P, and the water inflow risk level is divided according to the threshold (P < 0.3 is low risk, 0.3 ≤ P < 0.6 is medium risk, 0.6 ≤ P < 0.8 is high risk, and P ≥ 0.8 is extremely high risk).

[0158] Specifically, in another embodiment of the present invention, after the prediction model is trained, the ensorRT optimization technology is used to compress and accelerate the trained GRU prediction model to achieve efficient deployment on edge computing devices. The details are as follows:

[0159] The model's FP32 (32-bit floating point) parameters are converted to INT8 (8-bit integer) or FP16 (16-bit floating point) format, achieving improved performance while reducing parameter precision and reducing model size by 50%-75%. A calibration dataset is used during the quantization process to ensure that accuracy loss is controlled within 2%. Through vertical and horizontal layer fusion technology, multiple adjacent network layers (such as convolutional layers, batch normalization layers, and activation layers) are merged into a single computing unit, reducing the storage of intermediate results and the number of memory accesses, thereby reducing computational latency. The optimal CUDA kernel implementation is automatically selected based on the deployment hardware characteristics (such as GPU model and memory bandwidth) to maximize hardware resource utilization.

[0160] The lightweight model is deployed on the edge computing device. The model deployment adopts a microservice architecture, which includes the following components: data acquisition service, responsible for obtaining real-time data from the monitoring system; feature processing service, responsible for data preprocessing and feature extraction; inference service, loading the optimized model to perform predictive calculations; result distribution service, transmitting the prediction results to the construction control system.

[0161] This embodiment also designs an incremental learning strategy to enable the model to continuously learn from newly collected data and continuously improve prediction accuracy. The specific implementation is as follows:

[0162] Incremental learning is triggered by the following two conditions:

[0163] Periodic triggering: Incremental learning is automatically triggered every time an excavation cycle (usually 2-3 meters) is completed;

[0164] Error trigger: Incremental learning is triggered when the prediction error exceeds a preset threshold ε (ε=15%) for N consecutive times (N=5).

[0165] A sliding time window mechanism is used to select data for incremental learning. The window size W is defined as the monitoring records of the most recent M time points (usually M = 1440, corresponding to 1 day of minute-level data). Each time incremental learning is triggered, the window slides forward, incorporating newly collected data into the training set while removing the oldest data. The sliding window dataset is represented as:

[0166] D window =(X t ,y t )|t∈[t current-M+1 ,t current ]

[0167] Where, X t represents the input feature vector at time t; y t represents the actual water inflow at time t; t current Indicates the current time point.

[0168] Incremental learning updates model parameters through the gradient descent method, and its mathematical expression is:

[0169]

[0170] Where W t+1 is the updated network weight parameter set; W t is the current network weight parameter set; η t is the dynamic learning rate, which decays over time; is the gradient of the loss function calculated based on the window data. Dynamic learning rate η t Using the time decay strategy, the calculation formula is:

[0171]

[0172] Where η0 is the initial learning rate, which is set to 0.001; δ is the decay coefficient, which is set to 0.0005; and t is the cumulative number of incremental learning cycles.

[0173] Specifically, the loss function L window Consistent with the composite loss function used in the model training phase, it also includes three parts: mean square error loss, physical constraint loss, and time continuity loss.

[0174] To prevent catastrophic forgetting (i.e., the model forgets previously learned knowledge) during incremental learning, the present invention employs the following protection mechanisms:

[0175] Experience replay strategy: Maintain an experience pool with a capacity of K (usually K = 5000) to store historical important samples. During each incremental learning, randomly extract P (usually P = 1000) samples from the experience pool and mix them with new data for training.

[0176] Parameter regularization: Add a regularization term to the incremental learning loss function to limit the deviation of the new parameters from the original parameters:

[0177] L regularized =L window +λ reg ||W t+1 -W t || 2

[0178] Among them, λ reg is the regularization coefficient, usually set to 0.01.

[0179] Knowledge distillation: The original model prediction results are used as "soft labels" to assist in training, and the loss function adds distillation terms:

[0180] L distill =α·L window +(1-α)·KL(f(X;W t ),f(X;W t+1 ))

[0181] Among them, α is the balance coefficient, which is set to 0.7; KL represents the KL divergence, which measures the difference between two distributions; f(X; W) represents the predicted output of the model with parameter W for input X.

[0182] In summary, the present invention collects water inflow time-series data and geological parameter data through horizontal drilling. It then uses CEEMD decomposition to extract multi-scale time-series features. Combined with the physically calculated permeability coefficient, these features are input into a designed prediction model to predict water inflow and assess risk levels. Edge computing achieves millisecond-level response, and an incremental learning mechanism enables the model to continuously adapt to changing geological conditions. This effectively addresses the issues inherent in traditional methods, such as the disconnect between data and physical models, low prediction accuracy, and poor adaptability.

[0183] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A method for predicting water inflow during tunnel excavation by horizontal drilling, characterized in that: include: S1. Horizontal directional drilling is arranged in sections along the tunnel axis and a monitoring system is installed to obtain time series data on water inflow and geological parameters at each monitoring point; S2. Using the Goodman formula modified based on rock mass integrity, calculate the permeability coefficient at each monitoring point based on the water inflow time series data; S3, using the complete set empirical mode decomposition method to decompose the water inflow time series data to obtain multiple intrinsic mode functions and a residual term; S4. forming a training data set based on the intrinsic mode function, residual term, permeability coefficient and geological parameter data; S5. Construct a prediction model that includes a GRU neural network structure and a physical constraint mechanism, and train the prediction model using a training data set until the loss function converges to obtain a trained prediction model; S6. Obtain current monitoring data and input it into the trained prediction model to predict the water inflow and water inflow risk level of the section to be excavated.

2. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 1, characterized in that: In step S1, the step of laying out horizontal directional drilling holes in sections and installing a monitoring system includes: Use geological radar to conduct advance geological surveys ahead of the tunnel to identify potential high-risk water inrush areas; A monitoring point is set up every 15-20 meters of drilling in ordinary areas, and a monitoring point is set up every 5-10 meters of drilling in high-risk areas; High-precision flow sensors and data transmission devices are deployed at each monitoring point.

3. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 1, characterized in that: Geological parameter data include: groundwater level, borehole diameter, water-bearing body length, rock mass quality index RQD, fault fracture zone width, karst development degree and environmental influencing factors.

4. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 1, characterized in that: Step S3 includes: S31, adding white noise of different amplitudes to the water inflow time series data to generate multiple sets of data; S32, performing empirical mode decomposition on each set of noise-added data; S33. Average the decomposition results of all groups to obtain the intrinsic mode function and residual term.

5. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 4, characterized in that: In step S31 , the amplitude range of the white noise is 0.1-0.3 times the standard deviation of the original data, and the number of noise groups is 50-100 groups.

6. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 3, characterized in that: The Goodman formula corrected based on rock mass integrity in step S2 is: Where K is the permeability coefficient; Q is the time series data of water inflow; L is the length of the horizontal borehole through the water-bearing body; H is the head difference, which represents the vertical height difference between the groundwater level and the tunnel drainage surface; d is the cross-sectional diameter of the horizontal borehole; α(RQD) is the correction factor based on rock mass integrity, which is calculated based on the rock mass quality index RQD, as shown in the following formula: Where C1 is the basic correction value, and its value range is 0.15-0.25; C2 is the gain coefficient, and its value range is 0.75-0.85; k is the curve steepness coefficient, which controls the sensitivity of the correction coefficient to changes in RQD, and its value range is 0.04-0.06; R0 is the inflection point parameter, which indicates the critical value of the correction coefficient change, and its value range is 45-55.

7. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 1, characterized in that: Predictive models include: Feature preprocessing layer, which combines the intrinsic mode function and residual term to form time series features, and combines the permeability coefficient and geological parameters to form static features; The feature encoding layer extracts temporal features and static features, and concatenates them in the feature dimension to form a fused feature vector; The temporal modeling layer includes multiple stacked bidirectional GRU units to capture temporal dependencies. The output of each GRU unit is weighted by the attention mechanism to obtain the context feature representation. The prediction layer maps the context feature representation to the water inflow prediction value through a fully connected network; At the risk assessment layer, the ratio of the predicted water inflow volume to the safety threshold is calculated, and the water inflow risk level is divided according to the size of the ratio.

8. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 7, characterized in that: The loss function is: L total =L MSE +λ1(N)·L physical +λ2(N)·L temporal Where, L total is the total loss function; L MSE is the mean square error loss; L physical is the physical constraint loss; L temporal is the time continuity loss; m is the number of samples; y i is the actual water inflow of the i-th sample; is the predicted water inflow of the i-th sample; A tunnel is the tunnel cross-sectional area; K i is the permeability coefficient of the monitoring point corresponding to the i-th sample; is the hydraulic gradient of the monitoring point corresponding to the i-th sample, indicating the change in hydraulic head per unit distance; Δt i represents the time interval between two samples i and i+1; λ1(N) and λ2(N) are weight coefficients that change with the accumulated data volume N.

9. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 1, characterized in that: Step S4 includes: S41. Taking the intrinsic mode function and the residual term as time series features, and the permeability coefficient and geological parameter data as static features; S42. For each monitoring point, select the time series features and corresponding static features of T consecutive time points to form an input feature vector; and use the water inflow at time point T+1 as the output prediction target; S43, using a sliding window method with a window length of M and a sliding step size of 1, sequentially construct multiple input-output sample pairs; S44. Randomly divide the sample pairs into a training set, a validation set, and a test set to form a training data set.

10. The method for predicting water inflow during tunnel excavation by horizontal drilling according to claim 1, characterized in that: The determination of the water inrush risk level in step S6 includes: Calculate the ratio of the predicted water inflow to the safety threshold to obtain the risk index P; The risk level is divided according to the risk index P: P < 0.3 is low risk, 0.3 ≤ P < 0.6 is medium risk, 0.6 ≤ P < 0.8 is high risk, and P ≥ 0.8 is extremely high risk.

Images