A tunnel surrounding rock global deformation prediction method of physical constraint sequential bayesian inversion

Abstract

The present application belongs to the field of deep tunnel and roadway engineering surrounding rock deformation monitoring and prediction, and discloses a kind of physical constraint sequential bayesian inversion's tunnel surrounding rock global deformation prediction method, solve the existing tunnel surrounding rock deformation monitoring data dispersion, difficult to realize global deformation deduction, lack of physical and mechanical support, low in computational efficiency, cannot quantitatively predict the problem of uncertainty.This method defines rock mass parameter vector according to working condition, establishes surrounding rock simulation grid model, divides observation and prediction node set and collects measured displacement data in each cycle;And based on numerical simulation, the physical constraint prior distribution of the real displacement field of surrounding rock is constructed, and the measured displacement likelihood distribution is established combined with monitoring noise;Through analytical edge simplification calculation, the rock mass parameter edge likelihood analytical expression is constructed, and the parameter posterior distribution is dynamically updated combined with sequential bayesian mechanism;The maximum a posteriori estimation of rock mass parameter is obtained by BFGS iterative optimization, and finally the global displacement prediction of surrounding rock is realized based on the optimal parameter.

Description

Technical Field

[0001] This invention belongs to the field of monitoring and prediction technology of surrounding rock deformation in deep tunnel and roadway engineering, specifically involving a physical constraint sequential Bayesian inversion method for predicting the global deformation of surrounding rock in tunnels. Background Technology

[0002] As my country's transportation tunnels, mine roadways, and underground chambers continue to develop towards deeper and longer scales, the early warning and prevention of surrounding rock deformation and instability disasters has become a key challenge restricting the safety of underground engineering construction. Currently, on-site monitoring of surrounding rock deformation mainly relies on convergence gauges, displacement gauges, strain gauges, or total stations deployed at limited locations on the tunnel walls. Whether using manual or automated monitoring methods to record displacement at measuring points, only discrete displacement data from a limited number of measuring points on a few monitoring sections can be obtained. However, surrounding rock deformation is a continuous three-dimensional spatial distribution. Numerous undetected areas exist between measuring points and outside the measuring section. The surrounding rock in these areas may experience local instability or large deformation before the measuring points, creating safety hazards. Engineering practice also shows that many tunnel (roadway) surrounding rock deformation and instability failures often occur in areas where deformation monitoring points are not deployed. Furthermore, even when the surrounding rock monitoring feedback at the monitored points approaches stability, the surrounding rock in areas without monitored points may still experience large deformation or even collapse and instability.

[0003] Due to cost and practical constraints, it is impossible to deploy an unlimited number of continuous monitoring points on the surrounding rock of tunnels. Therefore, accurately and efficiently sensing the deformation state of the surrounding rock of tunnels using deformation data observed from a limited number of monitoring points is of paramount importance for disaster prevention and control in deep-earth engineering. This requires breaking through the limitations of traditional discrete monitoring, inferring the potential deformation trends of unmonitored areas by exploring the spatiotemporal correlation patterns between data from a limited number of monitoring points, and integrating advanced algorithms and physical and mechanical models to construct a global deformation field prediction technology for the surrounding rock, from "points" to "surfaces" and even "volumes." This will provide a scientific basis for optimizing support design and formulating precise disaster prevention and control strategies, thereby improving the safety level of deep-earth engineering.

[0004] Existing technologies include convergent deformation prediction studies based on empirical models (such as Logistic, Sachs, and hyperbolic models), which estimate the final deformation by fitting model parameters to monitored data. In recent years, some studies have introduced Bayesian methods to dynamically update empirical model parameters, improving prediction accuracy to some extent. Furthermore, forward modeling of surrounding rock deformation based on numerical simulation is widely used in tunnel design. However, these existing technologies have certain limitations:

[0005] (1) Limited prediction range: Existing methods mostly focus on the deformation prediction of the locations where measuring points have been set up, and cannot perceive the deformation state of the surrounding rock in areas where no monitoring instruments have been set up outside the measuring points.

[0006] (2) Insufficient integration of physical mechanisms and monitoring data: The prediction method based on empirical models is completely divorced from the laws of surrounding rock mechanics, and the reliability of extrapolation in far-field areas without monitoring data is seriously insufficient; although numerical simulation has the physical basis of rock mechanics, it is mostly a one-time forward analysis that is divorced from the feedback of measured data during construction. The physical model and monitoring data are independent and have not been effectively integrated, making it difficult to dynamically calibrate the simulation model as the excavation progresses.

[0007] (3) Lack of quantitative evaluation of the reliability of the surrounding rock deformation prediction: Most existing methods only provide deterministic deformation prediction values ​​and fail to provide the credibility range of the prediction results, which is not conducive to engineers' scientific weighing in risk decision-making.

[0008] (4) Low efficiency of parameter calculation and update: Some methods that introduce the Bayesian update idea often use random sampling strategies such as MCMC for posterior calculation, which requires repeated numerical simulation, resulting in high calculation cost and long convergence time. This makes it difficult to meet the actual needs of tunnel construction for real-time dynamic and rapid updating of surrounding rock parameters.

[0009] In summary, there is an urgent need for a solution that combines a rigorous physical constraint mechanism with efficient parameter update capabilities, and can use measured data from limited measurement points to predict the deformation of the surrounding rock over the entire area and quantify the uncertainty of the prediction. Summary of the Invention

[0010] The technical problem to be solved by this invention is to provide a physical constraint sequential Bayesian inversion method for predicting the global deformation of tunnel surrounding rock, which solves the problems of existing tunnel surrounding rock deformation monitoring data being discrete, making it difficult to extrapolate the global deformation, and traditional prediction methods lacking physical and mechanical support, having low computational efficiency, and being unable to quantify prediction uncertainties.

[0011] The technical solution adopted by the present invention to solve the above-mentioned technical problems is as follows:

[0012] A method for predicting the global deformation of tunnel surrounding rock using physically constrained sequential Bayesian inversion includes the following steps:

[0013] S1. Define the rock mass parameter vector corresponding to each excavation cycle condition, and determine the initial rock mass parameter vector;

[0014] S2. Construct a numerical grid simulation model for the surrounding rock, discretize all grid nodes to form a complete set of nodes, and divide the complete set of nodes into a subset of observed nodes and a subset of predicted nodes; after each round of excavation and support procedures is completed, collect the measured displacement vectors of the measuring points corresponding to the observed nodes and define a global true displacement vector containing the displacement components of all grid nodes.

[0015] S3. Using the rock mass parameter vector of the current excavation cycle as input, mechanical forward modeling is performed using the surrounding rock numerical grid simulation model to obtain the simulation global displacement vector; using the simulation global displacement vector as the mean, a physical constraint Gaussian prior distribution of the global real displacement vector is constructed.

[0016] S4. Construct a binary observation matrix based on the spatial correspondence between the subset of observation nodes and the entire set of nodes. Use the binary observation matrix to extract the theoretical displacement components of the observation nodes from the true displacement vector of the entire domain. Combine the random measurement errors generated by field monitoring to establish the Gaussian likelihood distribution of the measured displacement vector of the measuring point relative to the true displacement vector of the entire domain.

[0017] S5. Use the maximum a posteriori estimate of the logarithmic vector of rock mass parameters obtained from the previous excavation cycle as the prior mean of the logarithmic vector of rock mass parameters for the current excavation cycle; perform analytical marginalization on the true displacement vector of the entire domain, and derive the marginal likelihood function containing only the logarithmic vector of rock mass parameters; combine the prior distribution and the marginal likelihood function, and update the posterior distribution of the logarithmic vector of rock mass parameters for the current excavation cycle according to Bayes' theorem;

[0018] S6. Construct a negative logarithmic objective function with the goal of maximizing the posterior probability of the parameters. Use the BFGS quasi-Newton method to iteratively optimize the objective function, solve for the maximum posterior estimate of the logarithmic vector of the rock mass parameters under the current excavation cycle, and restore the optimal rock mass parameter vector.

[0019] S7. Substitute the optimal rock mass parameter vector into the surrounding rock numerical grid simulation model to calculate the global displacement prediction field of the surrounding rock, which includes all grid positions of observation nodes and prediction nodes. Calculate the confidence interval of the displacement of each prediction node based on sensitivity propagation.

[0020] S8. Use the maximum a posteriori estimate of the logarithmic vector of the rock mass parameters under the current working condition as the prior benchmark for the next excavation cycle working condition, and repeat steps S5 to S7 until the full construction cycle iteration is completed, and output the prediction results of the global deformation of the surrounding rock.

[0021] This scheme utilizes a sequential Bayesian iterative mechanism to dynamically update and continuously correct rock mass parameters for each excavation cycle. It accumulates and inherits the surrounding rock geological information contained in monitoring data from each construction stage, accurately adapting to the dynamic parameter changes caused by stress release, surrounding rock deterioration, and uneven strata during tunnel excavation. Simultaneously, analytical edge-mapping simplifies high-dimensional inversion calculations, and quasi-Newton iterative optimization achieves efficient parameter solving. While ensuring the physical rationality of deformation prediction, it significantly improves the timeliness of surrounding rock parameter inversion and deformation prediction, ultimately outputting global deformation prediction results covering all surrounding rock grid nodes. This provides reliable data support for tunnel surrounding rock stability analysis, dynamic assessment of construction conditions, and optimized design of support structures.

[0022] Furthermore, in step S1, the rock mass parameter vector corresponding to each excavation cycle condition is defined, and the initial rock mass parameter vector is determined, specifically including:

[0023] Define the tunnel number Rock mass parameter vector under excavation-support cycle condition :

[0024] ;

[0025] in, The number of rock mass parameters to be inverted; , This refers to the number of tunnel excavation-support cycles. For the first The first excavation-support cycle condition The rock mass parameters to be inverted; This is the matrix transpose.

[0026] The rock mass parameters include one or more of the following: elastic modulus, cohesion, internal friction angle, compressive strength, integrity coefficient, and groundwater pressure.

[0027] Based on the preliminary engineering geological survey data, indoor rock mechanics tests, and in-situ mechanical testing results, combined with comparative experience from similar projects in the same area and the geological conditions of the surrounding rock revealed during actual excavation, the initial values ​​of the rock mass parameter vector before excavation were determined. :

[0028] ;

[0029] in, For the first The first excavation-support cycle condition The initial values ​​of the rock mass parameters to be inverted.

[0030] In this scheme, independent rock mass parameter vectors are constructed for each excavation-support cycle of the tunnel. Standardized vector mathematical forms are used to accurately represent the independent properties of the surrounding rock mechanical parameters at different construction stages. This approach can precisely adapt to the construction characteristics of step-by-step excavation, progressive disturbance, and dynamic stress adjustment in tunnels, and can finely depict the stage-by-stage changes in the mechanical properties of the surrounding rock throughout the entire tunneling process. This scheme explicitly defines multiple core surrounding rock mechanical parameters as inversion objects, covering key control indicators such as surrounding rock strength, deformation characteristics, rock mass integrity, and hydrogeological conditions. Inversion parameter combinations can be flexibly selected according to different tunnel geological conditions, adapting to various complex engineering scenarios such as soft rock with large deformation, fractured surrounding rock, and water-rich strata, effectively improving the adaptability and versatility of parameter inversion.

[0031] Furthermore, in step S2, a numerical grid simulation model of the surrounding rock is constructed, and all grid nodes are discretized to form a complete node set. The complete node set is then divided into an observation node subset and a prediction node subset, specifically including:

[0032] A numerical mesh simulation model of the surrounding rock is constructed, and the surrounding rock area within the excavation disturbance range is discretized into a finite difference mesh or a finite element mesh. The set of all mesh nodes in the model is denoted as . The total number of nodes is denoted as The total degrees of freedom of the displacement of nodes in the global grid are denoted as . ,in, For spatial dimensions, in two-dimensional analysis 3D analysis ;

[0033] Based on the spatial location of the on-site monitoring points, from the complete set of nodes The subset of observed nodes is obtained by matching. The remaining grid nodes constitute a subset of the predicted nodes. ;

[0034] Observation node subset The observation nodes are grid nodes that correspond one-to-one with the spatial locations of displacement monitoring instruments deployed on-site, including measuring points on the tunnel perimeter surface and measuring points in deep boreholes; assuming the total number of on-site monitoring measuring points is... Observation node subset The total degrees of freedom of the corresponding measuring point displacement is denoted as ;

[0035] Predicting node subset The predicted nodes are the entire set of nodes. All grid nodes except for the observation nodes correspond to the spatial locations of the surrounding rock where no monitoring instruments have been deployed on-site and deformation prediction is required; the prediction node subset The total degrees of freedom of the corresponding measuring point displacement is denoted as .

[0036] In this scheme, the surrounding rock calculation area is discretized into a grid to construct a standardized numerical simulation grid node system, achieving precise spatial matching between the numerical simulation model and the field monitoring system. Based on the location of field monitoring points, observation node subsets and prediction node subsets are divided, clarifying the physical definitions, spatial attributes, and degrees of freedom characteristics of both types of nodes, thus achieving precise division between the measured monitoring area and the unknown prediction area. The observation node subset can accurately connect to field monitoring data, providing effective observation information for data-constrained inversion; the prediction node subset corresponds to the surrounding rock area without monitoring equipment coverage, clarifying the core target area for global deformation extrapolation.

[0037] Furthermore, in step S2, after each round of excavation and support procedures is completed, the measured displacement vectors of the corresponding measuring points of the observation nodes are collected, and a global true displacement vector containing all grid node displacement components is defined, specifically including:

[0038] Complete the first After each excavation-support cycle, displacement monitoring instruments deployed on-site are used to collect data from a subset of observation nodes. The deformation data of the measuring points are sorted and combined according to the degrees of freedom of the measuring points to form... Dimensional observed node measured displacement vector ;

[0039] definition dimensional global true displacement vector , The total degrees of freedom for the displacement of nodes in the global grid. Includes the complete set of nodes All displacement components of all grid nodes, and the global true displacement vector. The displacement components of the corresponding subset of observed nodes are denoted as The true displacement vector of the entire domain The displacement components of the corresponding predicted node subset are denoted as .

[0040] In this scheme, deformation monitoring data are collected in strict accordance with the construction sequence of each excavation-support cycle to ensure that each set of measured displacement data corresponds one-to-one with the corresponding construction condition, avoiding data distortion caused by time sequence misalignment. At the same time, a global true displacement vector containing displacement information of all grid nodes is constructed, and the observed components and predicted components are decomposed into blocks to clearly distinguish between known measured deformation information and unknown deformation information to be determined.

[0041] Furthermore, in step S3, using the rock mass parameter vector of the current excavation cycle as input, mechanical forward modeling is performed using the surrounding rock numerical grid simulation model to obtain the simulated global displacement vector. Using this simulated global displacement vector as the mean, a physical constraint Gaussian prior distribution of the global true displacement vector is constructed, specifically including:

[0042] Read the first The construction condition parameters for each excavation-support cycle are used to vectorize the rock mass parameters under the current working condition. And the construction condition parameters are input into the surrounding rock numerical grid simulation model;

[0043] The numerical grid simulation model of the surrounding rock is used to perform forward simulation analysis until the model reaches mechanical equilibrium and converges.

[0044] Extract the complete set of nodes The simulation calculates the displacement components of all mesh nodes and compares them with the true displacement vector across the entire domain. Consistent degree-of-freedom ordering rules, assembled to obtain 3D simulation global displacement vector ;

[0045] Introducing the model error vector Establish the mapping relationship between actual displacement and simulated displacement: Among them, the model error vector The components are independent and follow a mean of 0 and a variance of . The Gaussian normal distribution is as follows: , for The identity matrix, The variance of the model error;

[0046] Determine the model error variance To simulate the global displacement vector As the mean vector, As the covariance matrix, it is constructed under a given working condition. and rock mass parameter vector Under the condition of global true displacement vector Physically constrained Gaussian prior distribution:

[0047] ;

[0048] in, The prior probability density function is the physical constraint. For a given working condition and rock mass parameter vector The simulation of the global displacement vector under the given conditions; It is the square of the Euclidean norm.

[0049] In this scheme, the parameters of each cycle of construction conditions and the rock mass parameters are used together as simulation inputs to ensure that the forward modeling process of numerical simulation fully matches the actual construction conditions on site, thereby improving the realism and reliability of the simulated displacement field. Through the model's mechanical equilibrium convergence judgment mechanism, the simulation results are guaranteed to meet the physical criteria of stress balance and deformation coordination in the surrounding rock, avoiding calculation deviations caused by non-converged simulation results. This scheme constructs a physically constrained Gaussian prior distribution using the displacement vector of the entire simulation domain as the mean, integrating the mechanical mechanism of the surrounding rock into the probabilistic inference process. It automatically suppresses unreasonable displacement fields that violate mechanical laws, effectively improving the rationality of displacement extrapolation in sparse and non-existent measurement point areas, and ensuring the reliability of the overall deformation prediction results from a physical perspective.

[0050] Furthermore, in step S4, a binary observation matrix is ​​constructed based on the spatial correspondence between the subset of observed nodes and the entire set of nodes. This binary observation matrix is ​​then used to extract the theoretical displacement components of the observed nodes from the true displacement vector across the entire domain. Combined with random measurement errors generated by on-site monitoring, a Gaussian likelihood distribution of the measured displacement vector at the measuring point relative to the true displacement vector across the entire domain is established. Specifically, this includes:

[0051] Constructing a binary observation matrix Binary observation matrix The dimension is ;in, The total degrees of freedom of displacement of the subset of observation nodes, The total degrees of freedom for the displacement of nodes in the global grid; binary observation matrix. The matrix elements consist only of 0s and 1s, and each row in the matrix corresponds to an observation degree of freedom. If the column index matches a subset of observation nodes... For internal grid nodes, the element at that position is set to 1, and the elements at other positions are set to 0.

[0052] Binary observation matrix Multiply by the global true displacement vector on the left The theoretical displacement components corresponding to the observed nodes are extracted. ;

[0053] Introducing a measurement error vector that follows a Gaussian distribution Establish the global true displacement vector With the measured displacement vector The observation equation between:

[0054] ;

[0055] Construct a given global true displacement vector Measured displacement vector under the condition Gaussian likelihood probability density function :

[0056] ;

[0057] in, This represents the variance of the measurement error of the monitoring instrument.

[0058] This scheme introduces random measurement errors from on-site monitoring and establishes standardized observation equations to accurately reproduce the actual engineering characteristics of noise and bias in on-site monitoring data, thus avoiding distortion of the inversion results. By constructing a Gaussian likelihood distribution, the matching degree between theoretical and measured displacements is quantitatively characterized, realizing the probabilistic expression of measured monitoring information. This effectively integrates real on-site engineering data into the Bayesian inversion system, ensuring that the parameter inversion and deformation prediction results closely match the actual deformation state of the surrounding rock on-site.

[0059] Furthermore, in step S5, the maximum a posteriori estimate of the logarithmic vector of rock mass parameters obtained from the previous excavation cycle is used as the prior mean of the logarithmic vector of rock mass parameters for the current excavation cycle. The true displacement vector across the entire domain is analytically marginalized, and a marginal likelihood function containing only the logarithmic vector of rock mass parameters is derived. Combining the prior distribution and the marginal likelihood function, the posterior distribution of the logarithmic vector of rock mass parameters for the current excavation cycle is updated according to Bayes' theorem. Specifically, this includes:

[0060] rock mass parameter vector The logarithmic vector of rock mass parameters is obtained by performing a natural logarithmic transformation on each component. :

[0061] ;

[0062] Determine the prior distribution of the logarithmic vector of the rock mass parameters under the current working condition:

[0063] when When the first excavation-support cycle begins, the initial logarithmic parameter vector is taken. ;

[0064] when At that time, the previous cycle of excavation conditions will be repeated. Maximum a posteriori estimation of the logarithmic vector of rock mass parameters As for the current working conditions Logarithmic vector of rock mass parameters The prior mean is used to obtain the current working condition. Rock mass parameters Prior probability distribution :

[0065] ;

[0066] Joint global true displacement vector Prior probability With the measured displacement vector likelihood probability For the global true displacement vector Perform marginal integration to construct a logarithmic vector with respect to rock mass parameters only. Measured displacement vector of the observation node Marginal likelihood function:

[0067] ;

[0068] displacement relationship Substitute into the observation equation The expression for error superposition is obtained as follows: ;

[0069] Based on the independent property of model error and measurement error, solve for the total error covariance matrix:

[0070] ;

[0071] in, It is an identity matrix with dimension 1. ;

[0072] By utilizing the marginalization analytical property of the Gaussian distribution and avoiding high-dimensional numerical integration calculations, an analytical expression for the marginal likelihood function is obtained:

[0073] ;

[0074] Combining the prior probability distribution and marginal likelihood function of the logarithmic vector of rock mass parameters, the unnormalized posterior distribution of the logarithmic vector of rock mass parameters under the current working condition is obtained according to Bayes' theorem:

[0075] ;

[0076] Substituting the analytical expressions for the marginal likelihood function and the prior probability distribution of the logarithmic vector of rock mass parameters into the unnormalized posterior distribution of the logarithmic vector of rock mass parameters under the current working condition, we obtain the complete analytical form of the unnormalized posterior distribution expression:

[0077] ;

[0078] in, Let be the total error covariance matrix.

[0079] This scheme employs a sequential prior update mechanism, using the optimal parameter estimate from the previous construction cycle as the prior mean for the current cycle. This enables the temporal accumulation and transfer of surrounding rock geological information and parameter cognition, accurately reflecting the gradual evolution characteristics of surrounding rock parameters during tunnel excavation. By analytically marginalizing the high-dimensional global displacement latent variables and replacing traditional high-dimensional numerical integration and random sampling with pure matrix analytical operations, the computational dimensionality and computational overhead are significantly reduced, improving the efficiency of inversion calculations.

[0080] Furthermore, in step S6, a negative logarithmic objective function is constructed with the goal of maximizing the posterior probability of the parameters. The BFGS quasi-Newton method is used to iteratively optimize the objective function, solving for the maximum posterior estimate of the logarithmic vector of the rock mass parameters under the current excavation cycle, and reconstructing the optimal rock mass parameter vector. Specifically, this includes:

[0081] Maximizing the posterior probability density is equivalent to minimizing the following objective function. :

[0082] ;

[0083] ;

[0084] in, To maximize the posterior probability density of the rock mass logarithmic parameter vector; The data fitting term is defined as the numerical simulation displacement of the surrounding rock deformation observation nodes. With the measured displacement vector The weighted sum of squared residuals between and:

[0085] ;

[0086] The prior regularization term is defined as the current operating condition. Logarithmic vector of rock mass parameters Compared to the previous working condition Maximum a posteriori estimation Weighted bias:

[0087] ;

[0088] in, The total error covariance matrix is... Its inverse matrix;

[0089] Set the initial conditions for the BFGS quasi-Newton method iteration: set the initial iteration point... Set as the previous construction cycle condition Maximum posterior estimate ,Right now ;when At that time, take the initial logarithmic parameter vector. Use it as the initial point for iteration; and initialize the approximate Hessian inverse matrix. , for The identity matrix; and setting the iteration convergence tolerance. ;

[0090] Perform BFGS quasi-Newton iterative optimization: Starting from the set initial point, use the BFGS quasi-Newton method to optimize the objective function. The solution is obtained through iterative steps. During the iteration, the approximate Hessian inverse matrix is ​​updated based on the gradient information, and the iteration step size is determined using a line search. The iteration terminates when the gradient of the objective function satisfies the convergence tolerance, and the maximum a posteriori estimate of the logarithmic vector of the rock mass parameters under the current working condition is output. ;

[0091] The optimal rock mass parameter vector is obtained by inverse exponential transform: .

[0092] In this scheme, a negative logarithmic objective function is constructed, which includes a data fitting term and a prior regularization term. The data fitting term ensures that the deformation results derived from the inversion parameters can closely match the field measured data, thus guaranteeing the engineering adaptability of the results. The prior regularization term constrains the stability of the parameter iteration update, avoids parameter mutations caused by disturbances in single-cycle monitoring data, effectively suppresses the inversion overfitting problem, and ensures the continuity of the temporal evolution of rock mass parameters.

[0093] Furthermore, in step S7, the optimal rock mass parameter vector is substituted into the surrounding rock numerical grid simulation model to calculate the global displacement prediction field of the surrounding rock, which includes all grid positions of the observation nodes and prediction nodes. Based on sensitivity propagation, the confidence interval of the displacement of each prediction node is calculated, specifically including:

[0094] The optimal rock mass parameter vector under the current working condition Input the surrounding rock numerical mesh simulation model, and substitute the first... Mechanical forward modeling was performed on the construction condition parameters corresponding to each excavation-support cycle. After the model reached mechanical convergence and equilibrium, the complete set of nodes was extracted. The displacement components of all grid nodes are assembled to obtain the global displacement prediction field of the surrounding rock. :

[0095] ;

[0096] The global displacement prediction field of the surrounding rock It also includes the displacement components of the subset of observation nodes. and predicted node subset displacement components ;

[0097] For any prediction node Construct the logarithmic vector of the nodal displacement with respect to the rock mass parameters. Sensitivity vector : ;

[0098] in, To simulate and predict displacement; For the first Logarithmic vector of rock mass parameters under excavation-support cycle conditions The One component;

[0099] Obtain the approximate Hessian inverse matrix after BFGS iteration convergence. It is used as the posterior covariance matrix of the logarithmic vector of rock mass parameters. The sensitivity vector is used to propagate parameter uncertainty into the displacement space, and the model error variance is superimposed. Calculate prediction nodes The total variance of the displacement:

[0100] ;

[0101] Based on the assumption that the predicted displacement values ​​follow a normal distribution, at a confidence level... Next, the node to be predicted true displacement The confidence interval is expressed as:

[0102] ;

[0103] in, The standard normal distribution is at a confidence level of 100%. The two-sided quantile values ​​below;

[0104] Traverse all predicted node subsets Within the nodes, obtain the displacement prediction values ​​and corresponding confidence intervals for all locations to be predicted.

[0105] This scheme constructs a sensitivity vector of displacement with respect to rock mass parameters to accurately quantify the impact of rock mass parameter fluctuations on surrounding rock deformation, thereby achieving a quantitative transfer of parameter uncertainty to the displacement field. By comprehensively superimposing the cognitive uncertainty of parameters and the uncertainty of model simulation, a total displacement variance calculation system is constructed to fully cover the error sources of the prediction results, providing comprehensive probabilistic data support for surrounding rock stability assessment.

[0106] Furthermore, the method also includes the following steps:

[0107] S9. Based on the prediction results of the global deformation of the surrounding rock, the probability distribution of the displacement of each node to be predicted is used, combined with the preset safety threshold for the deformation of the surrounding rock, to calculate the probability of deformation exceeding the limit for each node to be predicted; the probability of deformation exceeding the limit is compared with the preset graded early warning threshold, and when the probability of deformation exceeding the limit reaches the corresponding graded early warning threshold, the corresponding level of surrounding rock disaster early warning is triggered.

[0108] This solution enables simultaneous risk assessment of both monitored and unmonitored areas, addressing the challenge of providing early warnings for hidden or unmonitored areas of tunnel surrounding rock. By implementing tiered early warnings based on the probability of exceeding limits, it accurately distinguishes the risk levels of the surrounding rock, providing a basis for decision-making regarding on-site support optimization and risk classification management, thereby enhancing safety management during tunnel construction.

[0109] The beneficial effects of this invention are:

[0110] (1) A physical Gaussian prior distribution constrained by mechanical mechanism was constructed to improve the reliability and accuracy of global deformation prediction of surrounding rock under sparse measuring point conditions:

[0111] This invention uses the global displacement field of the surrounding rock output from numerical physics simulation as the mean center of the Gaussian prior distribution of the real displacement field, and establishes a prior probability density function with mechanical constraints. It automatically assigns extremely low prior probabilities to unreasonable displacement fields that severely deviate from the mechanical equilibrium and deformation coordination laws of the surrounding rock. Unlike traditional purely data-driven, empirical statistical prediction methods that lack physical constraints and rely on measurement point data, this invention forces the global displacement field to conform to the mechanical principles of the surrounding rock through physical priors. Even in the far-field region of the surrounding rock where monitoring points are sparse and no monitoring data is available, the deformation prediction results still possess rigorous physical rationality. This solves the problem of prediction distortion and large errors in the blank areas outside the measurement points in purely data-driven methods, improving the overall prediction accuracy and engineering reliability of the global deformation field of the tunnel surrounding rock.

[0112] (2) By utilizing the analytical marginalization property of Gaussian distribution, accurate dimensionality reduction of high-dimensional latent variables can be achieved, thus reducing the computational cost of Bayesian inversion:

[0113] This invention utilizes the analytical properties of Gaussian distribution linear transformation and marginalization to analytically marginalize the high-dimensional latent variable—the global true displacement field—from the joint distribution of prior and likelihood. It derives an accurate analytical expression for the marginal likelihood function specific to rock mechanics parameters. The Mahalanobis distance between simulated and measured displacements is quantitatively measured using the inverse of the total error covariance matrix, achieving precise quantification of parameter fitting accuracy. This technical solution completely abandons the cumbersome computational mode of traditional Bayesian inversion, which relies on high-dimensional numerical integration and large-scale random sampling. It simplifies the highly complex marginalization solution process into efficient matrix operations, avoiding the computational consumption of massive simulation calls and iterative sampling. Compared to traditional probabilistic inversion methods such as particle filtering and Monte Carlo sampling, it eliminates the need to maintain large-scale particle sets and repeatedly call numerical simulation models, improving overall computational efficiency by 1-2 orders of magnitude and meeting the rapid computational requirements of dynamic inversion in tunnel construction.

[0114] (3) Combining BFGS quasi-Newton optimization with sequential initial value iteration strategy, both parameter inversion accuracy and iteration convergence efficiency are considered:

[0115] This invention constructs a negative logarithmic objective function based on the maximum a posteriori probability criterion and employs the BFGS quasi-Newton method for deterministic iterative optimization. It implicitly updates the approximate Hessian matrix using iterative gradient information, achieving superlinear convergence with only first-order gradient information. It eliminates the need to solve and store complex second-order derivative matrices, making it suitable for low-dimensional rock mass parameter inversion scenarios. It exhibits strong iterative stability and high solution accuracy. Simultaneously, it constructs a sequential update mechanism for excavation cycles, using the maximum a posteriori estimate of rock mass parameters from the previous construction cycle as the initial point for the current iteration. This aligns with the engineering characteristics of tunnel surrounding rock parameters gradually changing with excavation, significantly shortening the iterative convergence path and effectively compressing the conventional iteration steps to 10–30. While ensuring the accuracy of dynamic calibration of rock mass parameters, it significantly improves the efficiency of single-condition inversion, providing core algorithmic support for the dynamic real-time prediction of surrounding rock deformation in deep-buried tunnels and roadways.

[0116] (4) An engineering closed-loop application system that deeply integrates physical mechanisms and monitoring data has been formed, promoting the upgrading of tunnel surrounding rock safety management and control model:

[0117] This invention, based on a sequential Bayesian recursive framework, automatically integrates newly added on-site monitoring data after each round of excavation and support procedures, dynamically iterates and updates the surrounding rock mechanical parameters, and simultaneously outputs a global grid node deformation prediction field for the surrounding rock. It also combines a sensitivity uncertainty quantification method to obtain the displacement confidence interval for key areas and achieves graded disaster early warning based on probability-based exceedance discrimination. This method organically integrates physical simulation with Bayesian updates, providing a "lightweight, automated, and real-time" global rock deformation perception solution for tunnel and roadway engineering. It helps to promote the transformation and upgrading of deep tunnel (roadway) surrounding rock deformation and disaster monitoring from "point-based monitoring and experience-based judgment" to "global perception and probabilistic decision-making." Attached Figure Description

[0118] Figure 1 The flowchart illustrates a method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion, as provided in this embodiment of the invention. Detailed Implementation

[0119] This invention aims to provide a physical constraint sequential Bayesian inversion method for predicting the global deformation of tunnel surrounding rock, which solves the problems of existing tunnel surrounding rock deformation monitoring data being discrete, making it difficult to extrapolate the global deformation, and traditional prediction methods lacking physical and mechanical support, having low computational efficiency, and being unable to quantify prediction uncertainties. Its core idea is as follows: Combining the temporal construction characteristics of tunnel excavation and support, the physical mechanism of surrounding rock numerical simulation is deeply coupled with the probabilistic inference of on-site monitoring. A Gaussian prior distribution with physical constraints is constructed based on the numerical mechanics forward modeling results to ensure that the global deformation field meets the conditions of surrounding rock mechanical equilibrium and deformation coordination. Based on the Gaussian analytical marginalization theory, the high-dimensional global displacement latent variables are reduced in dimensionality to derive the analytical expression of the marginal likelihood of rock mass parameters, avoiding the high-dimensional numerical integration and random sampling operations of traditional probabilistic inversion. Based on the maximum a posteriori probability criterion, a negative logarithmic objective function is constructed. The BFGS quasi-Newton method combined with a sequential iterative initial value strategy is used to quickly solve the optimal estimate of rock mass parameters for each excavation cycle, realizing the dynamic temporal update of surrounding rock mechanical parameters with the construction process. Furthermore, the dual uncertainty of parameters and model is quantified through the sensitivity propagation mechanism to output the global deformation prediction field of surrounding rock and the confidence interval of each node. Combined with the deformation exceeding probability, a graded risk warning is realized. Finally, a set of intelligent prediction system for global deformation of tunnel surrounding rock with physical rigor, computational efficiency and engineering practicality is constructed.

[0120] Example:

[0121] This embodiment provides a method for predicting the global deformation of tunnel surrounding rock using physically constrained sequential Bayesian inversion. (See also...) Figure 1 It includes the following implementation process:

[0122] S1. Define the rock mass parameter vector for each working condition and determine the initial parameter values:

[0123] In this step, rock mass parameter vectors corresponding to each excavation cycle are defined, and initial rock mass parameter vectors are determined. The purpose is to independently define rock mass parameter vectors for each condition according to the tunnel excavation-support sequence, achieving a refined characterization of the dynamic evolution of surrounding rock mechanical parameters with excavation disturbances, thus breaking away from the static assumption of fixed parameters throughout the traditional simulation process. Simultaneously, by comprehensively determining reasonable initial rock mass parameter vectors through multi-source geological information, reliable and realistic initial values ​​for subsequent sequential Bayesian inversion are provided, effectively avoiding iterative divergence and local optima problems caused by random selection of initial parameter values, and ensuring the stability and physical rationality of the parameter inversion throughout the entire process.

[0124] In one exemplary implementation, this step is carried out as follows:

[0125] Linear underground engineering projects such as tunnels and roadways often traverse different types of rock strata during excavation. Even within the same type of rock strata, the physical and mechanical properties of the surrounding rock and the values ​​of various rock mass parameters dynamically change with each excavation cycle. Furthermore, the deformation of the surrounding rock after each excavation cycle is closely related to the values ​​of these rock mass parameters. To accurately describe the evolution of rock mass parameters with each excavation step, this invention defines the... Rock mass parameter vectors used for numerical simulation calculation of surrounding rock deformation under a single cyclic excavation condition. :

[0126] (1.1);

[0127] in, The number of rock mass parameters to be inverted; Representing the tunnel number Step "excavation-support" cycle ( ), Represents the total number of construction cycles; Indicates the first The first step of the excavation cycle corresponds to the first step. Rock mass parameters.

[0128] Elements in the rock mass parameter vector The physical names and symbols are, but are not limited to, those shown in Table 1. The specific types and quantities of rock mass parameters are determined by the constitutive model selected in the physical simulation analysis, such as the Mohr-Coulomb constitutive model, the Drucker-Prager constitutive model, or the Kelvin viscoelastic constitutive model.

[0129] Table 1. Names of Rock Mass Parameter Vector Components

[0130]

[0131] Based on the preliminary engineering geological survey data, indoor rock mechanics tests, and in-situ mechanical testing results, combined with experience from similar projects in the same area and the surrounding rock geological conditions revealed during actual excavation, the initial values ​​of the rock mass parameter vectors corresponding to the various types of rock strata traversed by the tunnel (tunnel) before excavation were comprehensively determined. :

[0132] (1.2);

[0133] The initial parameter vector It will serve as the basis for determining the prior mean of parameters in the first excavation cycle in the subsequent sequential Bayesian inversion, and as the calculation benchmark for the dynamic iterative update of rock mass parameters throughout the tunnel construction cycle.

[0134] S2. Construct a numerical simulation mesh model, partition the nodes into subsets, and define the true displacement vector:

[0135] In this step, a numerical grid simulation model of the surrounding rock is constructed, and all grid nodes are discretized to form a complete node set. This complete node set is then divided into an observation node subset and a prediction node subset. After each round of excavation and support procedures, the measured displacement vectors of the corresponding measuring points of the observation nodes are collected, and a global true displacement vector containing the displacement components of all grid nodes is defined. The purpose is to establish a standardized numerical calculation system for the surrounding rock through grid discretization, accurately delineate the field-measured constraint area and the area to be predicted by distinguishing between the observation node subset and the prediction node subset, and achieve precise spatial mapping between monitoring measuring points and simulation grid nodes. By collecting measured displacements in each construction cycle and constructing a global true displacement vector, the global deformation information of the surrounding rock is structurally stored, achieving decoupling and unification of measured data, the true displacement field, and the displacement field to be predicted. This provides a standardized data and model foundation for subsequent physical prior construction, application of observation constraints, and global deformation deduction.

[0136] In one exemplary implementation, this step is carried out as follows:

[0137] (1) Construct a surrounding rock deformation simulation network and divide the node set. :

[0138] A three-dimensional simulation analysis model of tunnel surrounding rock deformation is established. The surrounding rock area within the excavation disturbance influence range is discretized into a finite difference mesh or a finite element mesh. The set of all mesh nodes in the model is denoted as . The total number of nodes is denoted as The total degrees of freedom of the mesh node displacements are denoted as . ,in, For spatial dimensions, in two-dimensional analysis 3D analysis .

[0139] Based on the surrounding rock deformation monitoring and measurement scheme at the tunnel construction site, surface displacement measuring points (for monitoring crown settlement, tunnel perimeter horizontal convergence, etc.) are arranged on the free face of the excavated section of the surrounding rock, while multiple displacement gauges are installed in boreholes inside the surrounding rock to obtain deep displacement data. A global grid node set is then established based on the spatial correspondence of the monitoring points. Divide into two non-overlapping subsets of nodes:

[0140] ① Observation node subset : Grid nodes that correspond one-to-one with the spatial locations of displacement monitoring elements deployed on-site, including measuring points on the tunnel perimeter surface and measuring points in deep boreholes; assuming the total number of on-site monitoring measuring points is . The total degrees of freedom of the displacement of the measuring point are denoted as .

[0141] ②Predicted node subset The total degrees of freedom of displacement for all grid nodes except the observation nodes, corresponding to the spatial locations of the surrounding rock where no monitoring equipment is deployed on site and deformation prediction needs to be carried out. .

[0142] The above two types of node subsets satisfy the set relation. The displacement degrees of freedom satisfy .

[0143] (2) Define the excavation cycle index working condition And collect the measured displacement vectors of the surrounding rock observation nodes. :

[0144] Define a single tunnel excavation-support construction cycle index as a working condition. ( , (Total number of tunnel construction cycles). Single working condition. This includes: the excavation advance of this cycle, the excavation method (such as the full-face method or the bench method), and the design support parameters.

[0145] Complete the first After the excavation cycle is completed, on-site monitoring and measurement of surrounding rock deformation are carried out, and displacement data at each observation node are collected and combined to form a dimension. Measured displacement vector of surrounding rock observation node .

[0146] Furthermore, define the first Global nodal displacement vector of surrounding rock after excavation cycle (vector dimension) ). It fully includes all degrees of freedom displacement components of the surrounding rock observation nodes and prediction nodes, serving as the core state variables to be solved and estimated in the subsequent physical simulation forward modeling and sequential Bayesian inference process.

[0147] S3. Constructing a Gaussian prior distribution of global displacement physical constraints based on numerical forward modeling:

[0148] In this step, the rock mass parameter vector of the current excavation cycle is used as input, and mechanical forward modeling is performed using the surrounding rock numerical grid simulation model to obtain the simulated global displacement vector. Using this simulated global displacement vector as the mean, a physically constrained Gaussian prior distribution of the true global displacement vector is constructed. The purpose is to construct a Gaussian prior distribution based on the simulated global displacement field obtained from the surrounding rock mechanics numerical forward modeling, forcibly integrating the physical and mechanical mechanisms of stress balance and deformation coordination in the surrounding rock into the probabilistic inference framework, ensuring that the probability distribution of the displacement field strictly conforms to the laws of geotechnical mechanics. By constraining the prior mean through numerical simulation results, non-physical interpretations in areas without monitoring or with sparse measuring points are suppressed, addressing the shortcomings of pure data-driven prediction methods in areas with scarce data and poor reliability, thus ensuring the rationality and accuracy of the global deformation prediction results from a physical perspective.

[0149] In one exemplary implementation, this step is carried out as follows:

[0150] To establish the probabilistic mapping relationship from rock mass parameters to the surrounding rock displacement field, in the first... One "excavation-support" cycle (working condition) After completion, the rock mass parameter vector to be inverted in this cycle will be used. As a condition variable, combined with the specific working condition The corresponding construction parameters (including excavation advance, excavation method and design support parameters, etc.) are used to carry out forward modeling analysis of physical simulation of surrounding rock deformation by calling finite difference or finite element numerical models.

[0151] After the numerical simulation system converges to mechanical equilibrium through iterative calculations, the output will be... Numerical solution of global nodal displacement vector of surrounding rock This serves as the physical basis for subsequent inference of the global displacement probability of the surrounding rock.

[0152] The displacement field of the surrounding rock obtained by physical simulation (numerical solution of the global nodal displacement vector of the surrounding rock) ) and the actual surrounding rock displacement field (actual global displacement vector of the surrounding rock) There are objective deviations between them; these deviations originate from factors such as grid discretization errors, assumptions introduced by the constitutive model, treatment of rock mass heterogeneity, and uncertainties in geostress boundary conditions. This embodiment introduces the dimension as follows: Model error vector to indicate and The difference between the two can be expressed as follows:

[0153] (3.1);

[0154] This embodiment uses the numerical solution of the global nodal displacement vector of the surrounding rock output by the physical simulation. Considered as the global displacement vector of the real surrounding rock The optimal estimate. Model error vector. The components are independent and follow a mean of 0 and a variance of . The Gaussian normal distribution is as follows: ,in, for The identity matrix.

[0155] Therefore, under a given working condition and rock mass parameter vector Under the condition of, the global nodal displacement vector of the real surrounding rock The physical constraint prior probability density function can be expressed as a Gaussian distribution centered on the simulation solution:

[0156] (3.2);

[0157] In the formula, The model error variance characterizes the logarithmic simulation results. The level of trust is determined by engineers based on statistical results from physical simulations and on-site measurement data, as well as their professional experience. As a preference, this embodiment presupposes... The value is taken as the square of 1% to 5% of the magnitude of the surrounding rock deformation. This is the square of the Euclidean norm, which is the sum of the squares of the deviations of the vector components. Its physical meaning is the total nodal displacement vector of the real surrounding rock. Numerical solution of global nodal displacement vectors of surrounding rock The total square distance between them.

[0158] The results obtained in this step The prior probability distribution of the true displacement field of the tunnel surrounding rock is characterized, which enforces the mechanical equilibrium and compatibility of the displacement field across the entire surrounding rock. Specifically, any displacement field that significantly deviates from the physical simulation results is assigned an extremely low prior probability density, and unreasonable displacement solutions are automatically excluded in subsequent Bayesian inference. Furthermore, this prior probability distribution function is the core foundation for constructing the marginal likelihood and completing the Bayesian iterative update of the rock mass parameters in subsequent step S5.

[0159] S4. Construct a binary observation matrix and establish the Gaussian likelihood distribution of the measured displacement:

[0160] In this step, a binary observation matrix is ​​constructed based on the spatial correspondence between the subset of observation nodes and the entire set of nodes. This binary observation matrix is ​​then used to extract the theoretical displacement components of the observation nodes from the global true displacement vector. Furthermore, a Gaussian likelihood distribution of the measured displacement vector at the measuring points relative to the global true displacement vector is established, taking into account random measurement errors generated during field monitoring. The aim is to achieve a precise linear projection of the high-dimensional global true displacement field onto the low-dimensional measured point space using the binary observation matrix, accurately extracting the theoretical displacement components of the observation nodes. By combining this with random noise from field monitoring to establish a Gaussian likelihood distribution, the probabilistic mapping relationship between the measured data and the actual deformation of the surrounding rock is accurately depicted. This introduces the measured deformation information into the Bayesian inference system in the form of probabilistic constraints, effectively constraining the model parameters and deformation field using the measured data.

[0161] In one exemplary implementation, this step is carried out as follows:

[0162] Completed working condition After the excavation cycle, deformation data of all observation nodes were collected using various displacement monitoring elements deployed on-site, resulting in a dimension of [missing information]. Measured displacement vector of surrounding rock observation node .

[0163] Due to random measurement errors caused by factors such as sensor equipment accuracy and manual reading during on-site monitoring, this embodiment introduces an observation matrix. Establish global nodal displacement vectors of real surrounding rock Measured displacement vector of surrounding rock observation node The relationships between them are as follows:

[0164] (4.1);

[0165] In the formula, For dimension The observation matrix is ​​used to obtain the displacement vector of all nodes in the surrounding rock. The displacement components at the observation nodes of surrounding rock deformation are extracted. The elements of this matrix consist only of 0s and 1s, and each row of the matrix corresponds to an observation degree of freedom. If the column index matches the subset of observation nodes... For internal grid nodes, the element at that position is set to 1, and the other positions are set to 0. The matrix is ​​uniquely determined by the field measurement point layout scheme and the geometric correspondence between the numerical grid coordinates. For dimension The measurement instrument error vector. Its elements are independent and follow a mean of 0 and a variance of . The Gaussian distribution, i.e. .variance The accuracy is determined based on the factory calibration of the displacement monitoring instrument and the statistical results of repeated on-site measurements.

[0166] Therefore, given the global node displacement vector of the real surrounding rock At that time, the measured displacement vector of the surrounding rock observation node The likelihood probability density function is:

[0167] (4.2);

[0168] The results obtained in this step The likelihood probability distribution of the measured displacement field of the tunnel surrounding rock is characterized, and is used to quantitatively measure the true global displacement field (the true global nodal displacement vector of the surrounding rock). ) and the actual measured data at the field measuring points (actual displacement vector of the surrounding rock observation node) The degree of matching between the observed node processing theory displacement and the measured value is: the smaller the deviation between the displacement calculated by the observation node processing theory and the measured value, the higher the corresponding likelihood probability value, indicating that the displacement field can better interpret the field monitoring data. This likelihood function provides the core foundation for the subsequent step S5 to construct the marginal likelihood and complete the sequential Bayesian iterative update of the rock mass parameters.

[0169] S5. Constructing marginal likelihoods and updating the posterior distribution of rock mass parameters based on analytical marginalization:

[0170] In this step, the maximum a posteriori estimate of the logarithmic vector of rock mass parameters obtained from the previous excavation cycle is used as the prior mean of the logarithmic vector of rock mass parameters for the current excavation cycle. The true displacement vector across the entire domain is analytically marginalized to derive a marginal likelihood function containing only the logarithmic vector of rock mass parameters. Combining the prior distribution and the marginal likelihood function, the posterior distribution of the logarithmic vector of rock mass parameters for the current excavation cycle is updated according to Bayes' theorem. The aim is to construct the prior parameters for the current excavation cycle using the optimal parameters from the previous cycle, achieving the temporal accumulation and transmission of surrounding rock geological information, and accurately reflecting the gradual evolution of surrounding rock parameters with excavation. By analytically marginalizing the high-dimensional global displacement latent variables, the high computational cost of traditional Bayesian inversion, which relies on high-dimensional numerical integration and random sampling, is avoided, resulting in a high-precision analytical expression for the marginal likelihood containing only rock mass parameters. Combining the prior and marginal likelihood to update the posterior distribution of parameters enables efficient and accurate correction of rock mass parameters by each round of construction monitoring information, completing the probabilistic inversion and dimensionality reduction solution for the low-dimensional parameter space.

[0171] In one exemplary implementation, this step is carried out as follows:

[0172] S51. Based on the previous cycle of excavation conditions Maximum a posteriori estimation of the logarithmic vector of rock mass parameters Construct the logarithmic vector of rock mass parameters for the current working condition t. and its prior probability distribution .

[0173] Since the physical and mechanical parameters of rock mass (such as elastic modulus and cohesion) are naturally non-negative, this embodiment addresses the need for Gaussian distribution operations in the negative value region of the real number space during subsequent Bayesian inference. The logarithmic vector of rock mass parameters is obtained by performing a natural logarithmic transformation on each component. :

[0174] (5.1-1);

[0175] Due to the logarithmic vector of rock mass parameters The components are independent random variables, and their prior probability density functions can be described using a multivariate Gaussian distribution model. In the Bayesian update method proposed in this embodiment, the previous cyclic excavation condition is... Maximum a posteriori estimation of the logarithmic vector of rock mass parameters As for the current working conditions Logarithmic vector of rock mass parameters The prior mean can be used to obtain the current working condition. Rock mass parameters Prior probability distribution :

[0176] (5.1-2);

[0177] The names and meanings of the symbols in the above formula are shown in Table 2:

[0178] Table 2 Formula Parameter Names and Meanings

[0179]

[0180] The prior probability of the true displacement field obtained by combining S5.2 with S3 Likelihood probability of the measured displacement field obtained from S4 Construct the measured displacement vector of the surrounding rock observation node Marginal likelihood function :

[0181] The geological conditions and mechanical parameters of the surrounding rock traversed by deep-earth engineering projects such as tunnels and roadways change dynamically as the excavation face advances. To achieve logarithmic vectorization of rock mass parameters... The posterior distribution is dynamically updated with the surrounding rock deformation observation data during construction, and needs to be combined with the true displacement field obtained in step S3. Prior probability The observed displacement field obtained in step S4 likelihood probability The latent variables (true displacement field) Marginalization integration is used to eliminate latent variables. Construct a logarithmic vector with respect to rock mass parameters only. Measured displacement vector of surrounding rock observation node Marginal likelihood function:

[0182] (5.2-1);

[0183] The displacement relationship in equation (3.1) of step S3. Substitute the observation equation into equation (4.1) in step S4. ,have to:

[0184] (5.2-2);

[0185] Among them, due to the physical simulation error vector (dimension) ) and measuring instrument error vector (dimension) The two are independent and both follow a zero-mean Gaussian distribution. Therefore, the total error vector after their superposition is... (dimension) The covariance matrix of ) (dimension) )for:

[0186] (5.2-3);

[0187] In the formula, It is an identity matrix with dimension 1. .

[0188] Based on the linear transformation and marginalization properties of the Gaussian distribution, the marginal likelihood function has an analytical form, namely: Equation (5.2-1) above requires high-dimensional numerical integration to complete the marginalization process, which has an accurate analytical expression. The working condition t and the logarithmic vector of the rock mass parameters... Under the condition, the measured displacement vector of the surrounding rock observation node The marginal likelihood probability function is expressed as:

[0189] (5.2-4);

[0190] In the formula, Covariance matrix The determinant of the matrix represents the scatter volume of the total error distribution in the observation space.

[0191] This embodiment uses the derived analytical expression of equation (5.2-4) to replace the integral operation of equation (5.2-1), avoiding the bottleneck of high-dimensional numerical integration calculation, greatly reducing the calculation cost, and providing strong support for the rapid sequential Bayesian update calculation of the logarithmic vector posterior probability of rock mass parameters for various working conditions.

[0192] S5.3 Based on the measured data of surrounding rock deformation under the current excavation cycle condition t, a logarithmic vector diagram of rock mass parameters is generated. Sequential Bayesian update of posterior probability:

[0193] remember This indicates the period up to the previous construction cycle (operating condition). This involves all surrounding rock deformation observation data acquired after the completion of the process. The data is then analyzed using the S5.1 method based on the previous... The current operating conditions can be obtained from the historical observation data. Rock mass parameters Prior probability .

[0194] Furthermore, when the latest observation data for the current operating condition t is obtained... Then, Bayes' theorem can be used to realize the current working condition. Logarithmic vector of rock mass parameters Posterior probability sequential Bayes update is as follows:

[0195] (5.3-1);

[0196] In the formula, The current operating conditions are based on historical observation data. Logarithmic vector of rock mass parameters The prior probability; For working conditions Measured displacement vector of surrounding rock observation node The marginal likelihood function is obtained using the analytical method derived in step S5.2. The normalization constant (also known as the evidence term) is defined in the parameter space. The m-fold multiple integral over the given surface is as follows:

[0197] (5.3-2);

[0198] Normalization constant The integral calculation involves traversing all possible rock mass parameter values, resulting in extremely high computational costs. Although it is a necessary term for the completeness of the probability distribution, due to A logarithmic vector of rock mass parameters Irrelevant positive real numbers are only used for probability normalization and do not affect parameter optimization calculations. Therefore, the operating conditions are determined in the subsequent step S6. The calculation of the maximum a posteriori (MAP) estimation of the logarithmic vector of rock mass parameters can be omitted, simplifying the problem to determining only the maximum value of the unnormalized posterior, i.e.:

[0199] (5.3-3);

[0200] Substituting the measured displacement margin likelihood function (Equation 5.2-4) from step S5.2 and the prior distribution of rock mass parameters (Equation 5.1-2) from step 5.1 into the above equation (5.3-3), we obtain the complete unnormalized posterior distribution expression:

[0201] (5.3-4);

[0202] S6. Construct a negative logarithmic objective function and use the BFGS algorithm to solve for the optimal rock mass parameter vector:

[0203] In this step, a negative logarithmic objective function is constructed with the goal of maximizing the posterior probability of the parameters. The BFGS quasi-Newton method is used to iteratively optimize the objective function, solving for the maximum a posteriori estimate of the logarithmic vector of the rock mass parameters under the current excavation cycle, and then reconstructing the optimal rock mass parameter vector. The purpose is to transform the maximum a posteriori probability solution into a problem of minimizing the negative logarithmic objective function. The data fitting term ensures that the parameter inversion results fit the measured deformation characteristics, and the prior regularization term ensures the smooth and continuous evolution of the parameters over time, suppressing the overfitting phenomenon. The BFGS quasi-Newton method achieves efficient superlinear convergence by relying only on the first-order gradient information, without requiring the second-order derivative matrix, resulting in low computational cost and strong robustness. Combined with the sequential iteration initial value strategy, the number of iteration steps is greatly reduced, and the optimal rock mass parameter vector for each working condition is quickly obtained, realizing the dynamic and accurate calibration of the surrounding rock mechanical parameters.

[0204] In one exemplary implementation, this step is carried out as follows:

[0205] S6.1 Objective function for optimizing structural rock mass parameters :

[0206] As shown in step S5.3, the posterior distribution of rock mass parameters is... Proportional to marginal likelihood With parameter prior The product of these. To facilitate the solution, we take the negative natural logarithm of equation (5.3-4) and ignore the product with respect to the parameter. Irrelevant constant terms (including normalization evidence) and the normalization factors of the two Gaussian distributions , (etc.), maximizing the posterior probability density is equivalent to minimizing the following objective function :

[0207] (6.1-1);

[0208] (6.1-2);

[0209] in, The data fitting term is defined as the numerical simulation displacement of the surrounding rock deformation observation nodes. With the measured displacement vector The weighted sum of squared residuals between and:

[0210] (6.1-3);

[0211] The prior regularization term is defined as the current operating condition. Logarithmic vector of rock mass parameters Compared to the previous working condition Maximum a posteriori estimation Weighted bias:

[0212] (6.1-4);

[0213] In equation (6.1-2), The total error covariance matrix is ​​defined in equation (5.2-3) of step S5.2. Its inverse matrix represents the corresponding fitting weights assigned to observation data with different levels of reliability; in equation (6.1-4) For the previous working condition The posterior accuracy matrix of rock mass parameters is used as a regularization effect to "anchor" the current rock mass parameter estimate to the vicinity of historical knowledge, thereby making the objective function... The goal is to find the optimal balance between "fitting current measured data" and "maintaining the coherence of historical data coherence".

[0214] S6.2 BFGS Iterative Solution for the Logarithmic Vector of Optimal Rock Mass Parameters under Current Working Condition t :

[0215] 1. Initialize parameter settings:

[0216] (1) Set the iteration step Since the rock mass parameters change gradually rather than abruptly between two adjacent construction cycles, the initial point parameters of the optimization iteration will be... Set as the previous construction cycle condition The maximum a posteriori estimate, i.e. This initial point is already at or very close to the optimal solution for the current working condition t, which can significantly reduce the number of subsequent iterations. For the first excavation cycle (t=1). Take the logarithmic mean of the initial rock mass parameters determined in step S1.

[0217] (2) Initialize the symmetric positive definite approximate Hessian inverse matrix ( for (identity matrix), in the absence of prior curvature information, an initial value strategy of unit search step size is adopted for each parameter direction.

[0218] (3) Set convergence tolerance for iterative calculation This example takes .

[0219] 2. Enter the main iteration loop:

[0220] For iteration steps Perform the following sub-steps (a) to (g) in sequence:

[0221] (a) Calculate the gradient vector of the objective function :

[0222] According to equations (6.1-2) to (6.2-4) and the chain rule, the objective function at the current point... gradient vector at The parsing expression is:

[0223] (6.2-1);

[0224] In the formula, The displacement matrix of the measuring point obtained from physical simulation calculation. Logarithmic vector of rock mass parameters The sensitivity matrix, with dimension Its (i,j)th element is The value is taken at that location.

[0225] The rock mass parameter vector to be inverted As a condition variable, in this working condition In the iterative calculation, the displacement vector is obtained by calling a finite difference or finite element numerical simulation model once.

[0226] (b) Convergence test:

[0227] Calculate the current iteration step Euclidean norm of gradient || If satisfied This indicates that the stationary point of the objective function has been reached, the algorithm has converged, the iteration terminates, and the optimal solution is output. .

[0228] (c) Calculate the BFGS search direction vector :

[0229] Using the approximate Hessian inverse matrix of the current iteration step For gradient vector Perform curvature correction to generate the descent search direction vector. :

[0230] (6.2-2);

[0231] when identity matrix At that time, in this search direction That is, the direction of gradient descent. As the iterations proceed, The BFGS correction formula (6.2-7) described later is used to update the objective function, accumulating the local curvature information of the objective function so that the search direction gradually approaches the Newton direction.

[0232] (d) Determine the search step size that satisfies the Armijo sufficient descent condition. :

[0233] Along the search direction In the process of performing iterative calculations to find extreme points, to ensure that each iteration brings a sufficiently significant decrease in function value, while avoiding oscillations caused by excessively large step sizes and inefficiencies caused by excessively small step sizes, and to ensure the stability and rapid convergence of the BFGS algorithm, this embodiment adopts a method that satisfies the Armijo sufficient descent condition. As the search step size.

[0234] Search Step Size The specific process is as follows:

[0235] From the initial trial step size Initially, the step size is continuously reduced by a shrinkage factor ρ∈(0,1) (e.g.: ... (until...) The Armijo descent condition described below must be satisfied:

[0236] (6.2-3);

[0237] in, To fully reduce the constant, this example takes... .

[0238] Each attempt at the search step size only requires recalculation. There is no need to re-call the numerical simulation, because the displacement field has already been obtained in the gradient calculation of the current iteration step, and a linear approximation calculation can be performed in the search direction.

[0239] (e) Update the logarithmic vector of rock mass parameters for the next iteration. :

[0240] Along the search direction A search step size satisfying the Armijo sufficient descent condition is adopted. Update the logarithmic vector of rock mass parameters for the next iteration. :

[0241] (6.2-4);

[0242] (f) At the updated parameter point Calculate the gradient vector of the objective function for the next iteration step. :

[0243] The logarithmic vector of rock mass parameters updated in the above sub-step (e) Calculate and update the gradient of the objective function for the next iteration step. .

[0244] (g) Update the approximate Hessian inverse matrix using BFGS correction. :

[0245] Constructing an approximate Hessian inverse matrix for BFGS correction update Two vectors and ,as follows:

[0246] (6.2-5);

[0247] (6.2-6);

[0248] (i) Judgment If so, skip this BFGS update and maintain This ensures the symmetric positive definiteness of the approximate Hessian inverse matrix and the numerical stability of the algorithm.

[0249] (ii) Judgment When updating, use the following BFGS correction formula.

[0250] (6.2-7);

[0251] In the formula, for The identity matrix is ​​[value missing]. This BFGS correction formula possesses symmetric positive definiteness, ensuring that the search direction is always the descent direction, and the calculation process only requires storing a [value missing] matrix of size [value missing]. matrix The rock mass parameter dimensions addressed in this embodiment In low-dimensional inversion scenarios, storage resources are consumed very little.

[0252] 3. Termination and Output:

[0253] Repeat steps (a) to (f) above for iterative calculation until the convergence criterion is met. When the loop terminates, output the final convergence point. For the current construction cycle (working condition) The maximum a posteriori parameter estimate, i.e., the logarithmic vector of rock mass parameters under the current working condition. optimal solution ,Will Expanding the logarithmically yields the operating condition. Optimal solution of rock mass parameter vector :

[0254] (6.2-8);

[0255] The iterative calculation method provided in this embodiment ensures that the total number of numerical simulation calls during the BFGS optimization process of each construction cycle is equal to the product of the number of iteration steps and the number of simulations required in a single iteration. Since the physical and mechanical properties of the surrounding rock in adjacent construction cycles within the same stratum are often quite similar, the excellent initial values ​​provided by the sequential Bayesian update in this embodiment are crucial. The solution is close to the optimal solution of the current loop, resulting in a typical convergence step count of only 10 to 30. Therefore, compared with existing techniques (such as the MCMC method) that require hundreds or thousands of simulation calls for random sampling, this method significantly reduces convergence time and saves computational resources while ensuring the accuracy of rock mass parameter estimation.

[0256] S7. Calculate the global deformation field based on optimal parameters and quantify the uncertainty of displacement prediction:

[0257] In this step, the optimal rock mass parameter vector is substituted into the surrounding rock numerical grid simulation model to calculate the global displacement prediction field of the surrounding rock, which includes all grid positions of the observation nodes and prediction nodes. The confidence interval of the displacement of each prediction node is then calculated based on sensitivity propagation. The aim is to re-perform numerical forward modeling using dynamically calibrated optimal rock mass parameters, generating a global deformation prediction field of the surrounding rock covering both the observed and unmonitored areas, thus overcoming the spatial limitations of traditional single-point monitoring. Through the sensitivity propagation mechanism, the uncertainty of rock mass parameters is transferred to the displacement space, and the inherent model error is superimposed to achieve dual uncertainty quantification. The confidence interval of the displacement of each prediction node is then calculated, overcoming the drawback of traditional single deterministic predictions that cannot evaluate the reliability of the results, and achieving simultaneous output of surrounding rock deformation prediction results and risk credibility.

[0258] In one exemplary implementation, this step is carried out as follows:

[0259] S7.1 Determine the predicted values ​​of nodal displacements across the entire surrounding rock area :

[0260] To obtain the optimal solution of the rock mass parameter vector for the t-th construction cycle (condition t). Then, substituting the finite difference or finite element numerical model described in step S3, a physical simulation forward modeling analysis is performed under the construction parameter conditions corresponding to working condition t to obtain the global node displacement vector of the tunnel surrounding rock. .

[0261] (7.1-1);

[0262] Displacement vector Characterizes the working condition The global displacement field of the surrounding rock includes the mesh node set of the surrounding rock model. All degrees of freedom displacement components. Among them, the subset of observation nodes. The displacement components are denoted as Predicting a subset of nodes The displacement components are denoted as ,Right now: .

[0263] S7.2 Reliable quantitative assessment of predicted displacement values ​​at any node in the surrounding rock:

[0264] Based on the global nodal displacement vector of the surrounding rock obtained in S7.1 above It can extract any point from the spatial prediction node subset of the surrounding rock. The displacement components are denoted as: The accuracy of this displacement prediction is affected by two factors: firstly, the error of the physical simulation model defined in step S3; and secondly, the residual parameter cognitive uncertainty in the posterior distribution of the rock mass parameters obtained in step S6. To quantitatively assess any node... displacement components To assess the reliability of the displacement prediction, this embodiment employs a sensitivity propagation mechanism to transfer the uncertainty in the parameter space to the displacement space, establishing a variance estimate of the displacement of the node to be predicted to evaluate the reliability of the displacement prediction. Specifically, this includes:

[0265] (1) For any node The displacement is defined by simulation prediction. Logarithmic vector of rock mass parameters The sensitivity vector composed of the partial derivatives of each component :

[0266] (7.2-1);

[0267] Wherein, the sensitivity vector The partial derivatives of each element value are determined using the central finite difference method: in the optimal estimation of the logarithmic vector of rock mass parameters. At each logarithmic parameter component Apply small positive and negative perturbations respectively. (This example takes) (and call the numerical simulation model to obtain the nodes under the corresponding working conditions) The displacement value is approximated by the following formula for each partial derivative component:

[0268] (7.2-2);

[0269] Therefore, any prediction node Displacement prediction value variance It can be obtained by superimposing the posterior covariance contribution term of the parameters and the model error variance, as shown in the following formula:

[0270] (7.2-3);

[0271] In the formula, The approximate inverse matrix after convergence of BFGS optimization in step S6. The final value; this matrix is ​​equivalent to the covariance matrix of the posterior distribution of the current operating condition parameters. ;

[0272] This represents the model error variance, specifically the physical simulation model error variance defined in step S3. It characterizes the numerical simulation results. The level of trust is determined by the engineer based on statistical results of physical simulation and on-site measurement data, as well as their professional experience. As a preferred approach, this embodiment presupposes... The value is taken as the square of 1% to 5% of the magnitude of the surrounding rock deformation.

[0273] Obtain Node displacement value variance Subsequently, since the predicted nodal displacement values ​​follow the predicted mean... If the distribution is a normal distribution centered at a certain confidence level, then at a certain confidence level... Next, the node to be predicted True displacement The confidence interval can be expressed as:

[0274] (7.2-4);

[0275] In the formula, The standard normal distribution is at a confidence level of 100%. The two-sided quantile values ​​below.

[0276] In engineering, a confidence level of 95% is typically used, that is: ,at this time Then the node The 95% confidence interval for the displacement prediction is [ , ].

[0277] S8. Multi-condition sequential iterative update outputs full-cycle global deformation prediction results:

[0278] In this step, the maximum a posteriori estimate of the logarithmic vector of the rock mass parameters under the current working condition is used as the prior benchmark for the next excavation cycle. Steps S5 to S7 are repeated until the entire construction cycle iteration is completed, and the prediction result of the global deformation of the surrounding rock is output. The purpose is to establish a Bayesian sequential iterative update mechanism for each excavation cycle, using the optimal parameters of a single cycle as the prior benchmark for the next cycle, continuously accumulating monitoring information and surrounding rock knowledge throughout the entire construction cycle, dynamically adapting to stress release, surrounding rock deterioration, and stratum changes during the tunneling process; through full-cycle iterative closed-loop updates of rock mass parameters and deformation fields, prediction deviations are continuously corrected, and finally, a high-precision prediction result of the global deformation of the surrounding rock that conforms to the evolution law of the entire tunnel construction is output, providing complete, dynamic, and reliable basic data for subsequent risk warning and construction optimization.

[0279] Based on the prediction results of the overall deformation of the surrounding rock, the probability distribution of the displacement of each node to be predicted is used, combined with the preset safety threshold for surrounding rock deformation, to calculate the probability of deformation exceeding the limit for each node. The probability of deformation exceeding the limit is compared with the preset graded early warning threshold. When the probability of deformation exceeding the limit reaches the corresponding graded early warning threshold, the corresponding level of surrounding rock disaster early warning is triggered. This guides the implementation of targeted prevention and control measures on site, such as denser support, optimized excavation progress, and evacuation of personnel.

[0280] This method fully considers the dual interference factors of uncertainties in rock mass parameter cognition and inherent errors in numerical models. It can not only achieve full-domain extrapolation and prediction of surrounding rock deformation in areas without monitoring points, but also quantify the reliability of deformation instability risk based on Bayesian posterior distribution and probability intervals. Compared with traditional deterministic threshold early warning methods that rely solely on single-point measured values, it can significantly reduce the probability of missed and false alarms in surrounding rock disaster early warning, and improve the reliability of geological risk management throughout the tunnel construction process.

[0281] Although embodiments of the present invention have been described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the present invention, and all such changes and alterations shall not depart from the protection scope of the present invention.

Claims

1. A method for predicting the global deformation of tunnel surrounding rock using physically constrained sequential Bayesian inversion, characterized in that, Includes the following steps: S1. Define the rock mass parameter vector corresponding to each excavation cycle condition, and determine the initial rock mass parameter vector; S2. Construct a numerical grid simulation model for the surrounding rock, discretize all grid nodes to form a complete set of nodes, and divide the complete set of nodes into a subset of observed nodes and a subset of predicted nodes; after each round of excavation and support procedures is completed, collect the measured displacement vectors of the measuring points corresponding to the observed nodes and define a global true displacement vector containing the displacement components of all grid nodes. S3. Using the rock mass parameter vector of the current excavation cycle as input, mechanical forward modeling is performed using the surrounding rock numerical grid simulation model to obtain the simulation global displacement vector; using the simulation global displacement vector as the mean, a physical constraint Gaussian prior distribution of the global real displacement vector is constructed. S4. Construct a binary observation matrix based on the spatial correspondence between the subset of observation nodes and the entire set of nodes. Use the binary observation matrix to extract the theoretical displacement components of the observation nodes from the true displacement vector of the entire domain. Combine the random measurement errors generated by field monitoring to establish the Gaussian likelihood distribution of the measured displacement vector of the measuring point relative to the true displacement vector of the entire domain. S5. Use the maximum a posteriori estimate of the logarithmic vector of rock mass parameters obtained from the previous excavation cycle as the prior mean of the logarithmic vector of rock mass parameters for the current excavation cycle; perform analytical marginalization on the true displacement vector across the entire domain, and derive the marginal likelihood function containing only the logarithmic vector of rock mass parameters; combine the prior distribution and the marginal likelihood function, and update the posterior distribution of the logarithmic vector of rock mass parameters for the current excavation cycle according to Bayes' theorem: rock mass parameter vector The logarithmic vector of rock mass parameters is obtained by performing a natural logarithmic transformation on each component. : ;in, This is the transpose of the matrix; For the first The first excavation-support cycle condition The rock mass parameters to be inverted; The number of rock mass parameters to be inverted; Determine the prior distribution of the logarithmic vector of the rock mass parameters under the current working condition: when When the first excavation-support cycle begins, the initial logarithmic parameter vector is taken. ; The initial values ​​for the rock mass parameter vector before excavation; when At that time, the previous cycle of excavation conditions will be repeated. Maximum a posteriori estimation of the logarithmic vector of rock mass parameters As for the current working conditions Logarithmic vector of rock mass parameters The prior mean is used to obtain the current working condition. Rock mass parameters Prior probability distribution : ； Joint global true displacement vector Prior probability With the measured displacement vector Likelihood probability For the global true displacement vector Perform marginal integration to construct a logarithmic vector with respect to rock mass parameters only. Measured displacement vector of the observation node Marginal likelihood function: ； displacement relationship Substitute into the observation equation The expression for error superposition is obtained as follows: ; in, The measurement error vector follows a Gaussian distribution; It is a binary observation matrix; To simulate the global displacement vector; This is the model error vector; Based on the independent property of model error and measurement error, solve for the total error covariance matrix: ； in, The total error covariance matrix; It is an identity matrix with dimension 1. ; The total degrees of freedom for the displacement of the measured points in the subset of observation nodes; The variance of the model error; To monitor the measurement error variance of the instrument; By utilizing the marginalization analytical property of the Gaussian distribution and avoiding high-dimensional numerical integration calculations, an analytical expression for the marginal likelihood function is obtained: ； Combining the prior probability distribution and marginal likelihood function of the logarithmic vector of rock mass parameters, the unnormalized posterior distribution of the logarithmic vector of rock mass parameters under the current working condition is obtained according to Bayes' theorem: ； Substituting the analytical expressions for the marginal likelihood function and the prior probability distribution of the logarithmic vector of rock mass parameters into the unnormalized posterior distribution of the logarithmic vector of rock mass parameters under the current working condition, we obtain the complete analytical form of the unnormalized posterior distribution expression: ； S6. Construct a negative logarithmic objective function with the goal of maximizing the posterior probability of the parameters. Use the BFGS quasi-Newton method to iteratively optimize the objective function, solve for the maximum posterior estimate of the logarithmic vector of the rock mass parameters under the current excavation cycle, and restore the optimal rock mass parameter vector. S7. Substitute the optimal rock mass parameter vector into the surrounding rock numerical grid simulation model to calculate the global displacement prediction field of the surrounding rock, which includes all grid positions of observation nodes and prediction nodes. Calculate the confidence interval of the displacement of each prediction node based on sensitivity propagation. S8. Use the maximum a posteriori estimate of the logarithmic vector of the rock mass parameters under the current working condition as the prior benchmark for the next excavation cycle working condition, and repeat steps S5 to S7 until the full construction cycle iteration is completed, and output the prediction results of the global deformation of the surrounding rock.

2. The method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion as described in claim 1, characterized in that, In step S1, the rock mass parameter vector corresponding to each excavation cycle condition is defined, and the initial rock mass parameter vector is determined, specifically including: Define the tunnel number Rock mass parameter vector under excavation-support cycle condition : ； in, The number of rock mass parameters to be inverted; , This refers to the number of tunnel excavation-support cycles. For the first The first excavation-support cycle condition The rock mass parameters to be inverted; This is the matrix transpose. The rock mass parameters include one or more of the following: elastic modulus, cohesion, internal friction angle, compressive strength, integrity coefficient, and groundwater pressure. Based on the preliminary engineering geological survey data, indoor rock mechanics tests, and in-situ mechanical testing results, combined with comparative experience from similar projects in the same area and the geological conditions of the surrounding rock revealed during actual excavation, the initial values ​​of the rock mass parameter vector before excavation were determined. : ； in, For the first The first excavation-support cycle condition The initial values ​​of the rock mass parameters to be inverted.

3. The method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion as described in claim 1, characterized in that, In step S2, a numerical grid simulation model of the surrounding rock is constructed, and all grid nodes are discretized to form a complete node set. The complete node set is then divided into an observation node subset and a prediction node subset, specifically including: A numerical mesh simulation model of the surrounding rock is constructed, and the surrounding rock area within the excavation disturbance range is discretized into a finite difference mesh or a finite element mesh. The set of all mesh nodes in the model is denoted as . The total number of nodes is denoted as The total degrees of freedom of the displacement of nodes in the global grid are denoted as . ,in, For spatial dimensions, in two-dimensional analysis 3D analysis ; Based on the spatial location of the on-site monitoring points, from the complete set of nodes The subset of observed nodes is obtained by matching. The remaining grid nodes constitute a subset of the predicted nodes. ; Observation node subset The observation nodes are grid nodes that correspond one-to-one with the spatial locations of displacement monitoring instruments deployed on-site, including measuring points on the tunnel perimeter surface and measuring points in deep boreholes; assuming the total number of on-site monitoring measuring points is... Observation node subset The total degrees of freedom of the corresponding measuring point displacement is denoted as ; Predicting node subset The predicted nodes are the entire set of nodes. All grid nodes except for the observation nodes correspond to the spatial locations of the surrounding rock where no monitoring instruments have been deployed on-site and deformation prediction is required; the prediction node subset The total degrees of freedom of the corresponding measuring point displacement is denoted as .

4. The method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion as described in claim 3, characterized in that, In step S2, after each round of excavation and support procedures is completed, the measured displacement vectors of the measuring points corresponding to the observation nodes are collected, and a global true displacement vector containing all grid node displacement components is defined, specifically including: Complete the first After each excavation-support cycle, displacement monitoring instruments deployed on-site are used to collect data from a subset of observation nodes. The deformation data of the measuring points are sorted and combined according to the degrees of freedom of the measuring points to form... Dimensional observed node measured displacement vector ; definition dimensional global true displacement vector , The total degrees of freedom for the displacement of nodes in the global grid. Includes the complete set of nodes All displacement components of all grid nodes, and the global true displacement vector. The displacement components of the corresponding subset of observed nodes are denoted as The true displacement vector of the entire domain The displacement components of the corresponding predicted node subset are denoted as .

5. The method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion as described in claim 4, characterized in that, In step S3, using the rock mass parameter vector of the current excavation cycle as input, mechanical forward modeling is performed using the surrounding rock numerical grid simulation model to obtain the simulation global displacement vector. Using this simulation global displacement vector as the mean, a physical constraint Gaussian prior distribution of the global true displacement vector is constructed, specifically including: Read the first The construction condition parameters for each excavation-support cycle are used to vectorize the rock mass parameters under the current working condition. And the construction condition parameters are input into the surrounding rock numerical grid simulation model; The numerical grid simulation model of the surrounding rock is used to perform forward simulation analysis until the model reaches mechanical equilibrium and converges. Extract the complete set of nodes The simulation calculates the displacement components of all mesh nodes and compares them with the true displacement vector across the entire domain. Consistent degree-of-freedom ordering rules, assembled to obtain 3D simulation global displacement vector ; Introducing the model error vector Establish the mapping relationship between real displacement and simulated displacement: Among them, the model error vector The components are independent and follow a mean of 0 and a variance of . The Gaussian normal distribution is as follows: , for The identity matrix, The variance of the model error; Determine the model error variance To simulate the global displacement vector As the mean vector, As the covariance matrix, it is constructed under a given working condition. and rock mass parameter vector Under the condition of global true displacement vector Physically constrained Gaussian prior distribution: ； in, The prior probability density function is the physical constraint. For a given working condition and rock mass parameter vector The simulation of the global displacement vector under the given conditions; It is the square of the Euclidean norm.

6. The method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion as described in claim 5, characterized in that, In step S4, a binary observation matrix is ​​constructed based on the spatial correspondence between the subset of observed nodes and the entire set of nodes. This binary observation matrix is ​​then used to extract the theoretical displacement components of the observed nodes from the true displacement vector across the entire domain. Combined with random measurement errors generated by on-site monitoring, a Gaussian likelihood distribution of the measured displacement vector at the measuring points relative to the true displacement vector across the entire domain is established. Specifically, this includes: Constructing a binary observation matrix Binary observation matrix The dimension is ;in, The total degrees of freedom of displacement of the subset of observation nodes. The total degrees of freedom for the displacement of nodes in the global grid; binary observation matrix. The matrix elements consist only of 0s and 1s, and each row in the matrix corresponds to an observation degree of freedom. If the column index matches a subset of observation nodes... For internal grid nodes, the element at that position is set to 1, and the elements at other positions are set to 0. Binary observation matrix Multiply by the global true displacement vector on the left The theoretical displacement components corresponding to the observed nodes are extracted. ; Introducing a measurement error vector that follows a Gaussian distribution Establish the global true displacement vector With the measured displacement vector The observation equation between: ； Construct a given global true displacement vector Measured displacement vector under the condition Gaussian likelihood probability density function : ； in, This represents the variance of the measurement error of the monitoring instrument.

7. The method for predicting the global deformation of tunnel surrounding rock using physically constrained sequential Bayesian inversion as described in claim 6, characterized in that, In step S6, a negative logarithmic objective function is constructed with the goal of maximizing the posterior probability of the parameters. The BFGS quasi-Newton method is used to iteratively optimize the objective function, solving for the maximum a posteriori estimate of the logarithmic vector of the rock mass parameters under the current excavation cycle, and reconstructing the optimal rock mass parameter vector. Specifically, this includes: Maximizing the posterior probability density is equivalent to minimizing the following objective function. : ； ； in, To maximize the posterior probability density of the rock mass logarithmic parameter vector; The data fitting term is defined as the numerical simulation displacement of the surrounding rock deformation observation nodes. With the measured displacement vector The weighted sum of squared residuals between and: ； The prior regularization term is defined as the current operating condition. Logarithmic vector of rock mass parameters Compared to the previous working condition Maximum a posteriori estimation Weighted bias: ； in, The total error covariance matrix is... Its inverse matrix; Set the initial conditions for the BFGS quasi-Newton method iteration: set the initial iteration point... Set as the previous construction cycle condition Maximum posterior estimate ,Right now ;when At that time, take the initial logarithmic parameter vector. This serves as the initial point for iteration; and the approximate Hessian inverse matrix is ​​initialized. , for The identity matrix; and setting the iteration convergence tolerance. ; Perform BFGS quasi-Newton iterative optimization: Starting from the set initial point, use the BFGS quasi-Newton method to optimize the objective function. The solution is obtained through iterative steps. During the iteration, the approximate Hessian inverse matrix is ​​updated based on the gradient information, and the iteration step size is determined by line search. The iteration terminates when the gradient of the objective function satisfies the convergence tolerance, and the maximum a posteriori estimate of the logarithmic vector of the rock mass parameters under the current working condition is output. ; The optimal rock mass parameter vector is obtained by inverse exponential transform: .

8. The method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion as described in claim 7, characterized in that, In step S7, the optimal rock mass parameter vector is substituted into the surrounding rock numerical grid simulation model to calculate the global displacement prediction field of the surrounding rock, which includes all grid positions of the observation nodes and prediction nodes. Based on sensitivity propagation, the confidence interval of the displacement of each prediction node is calculated, specifically including: The optimal rock mass parameter vector under the current working condition Input the surrounding rock numerical mesh simulation model, and substitute the first... Mechanical forward modeling was performed on the construction condition parameters corresponding to each excavation-support cycle. After the model reached mechanical convergence and equilibrium, the complete set of nodes was extracted. The displacement components of all grid nodes are assembled to obtain the global displacement prediction field of the surrounding rock. : ； The global displacement prediction field of the surrounding rock It also includes the displacement components of the subset of observation nodes. and predicted node subset displacement components ; For any prediction node Construct the logarithmic vector of the nodal displacement with respect to the rock mass parameters. Sensitivity vector : ; in, To simulate and predict displacement; For the first Logarithmic vector of rock mass parameters under excavation-support cycle conditions The One component; Obtain the approximate Hessian inverse matrix after BFGS iteration convergence. It is used as the posterior covariance matrix of the logarithmic vector of rock mass parameters. The sensitivity vector is used to propagate parameter uncertainty into the displacement space, and the model error variance is superimposed. Calculate prediction nodes The total variance of the displacement: ； Based on the assumption that the predicted displacement values ​​follow a normal distribution, at a confidence level... Next, the node to be predicted true displacement The confidence interval is expressed as: ； in, The standard normal distribution is at a confidence level of 100%. The two-sided quantile values ​​below; Traverse all predicted node subsets Within the nodes, obtain the displacement prediction values ​​and corresponding confidence intervals for all locations to be predicted.

9. A method for predicting the global deformation of tunnel surrounding rock using physical constraint sequential Bayesian inversion as described in any one of claims 1 to 8, characterized in that, The method also includes the following steps: S9. Based on the prediction results of the global deformation of the surrounding rock, the probability distribution of the displacement of each node to be predicted is used, combined with the preset safety threshold for the deformation of the surrounding rock, to calculate the probability of deformation exceeding the limit for each node to be predicted; the probability of deformation exceeding the limit is compared with the preset graded early warning threshold, and when the probability of deformation exceeding the limit reaches the corresponding graded early warning threshold, the corresponding level of surrounding rock disaster early warning is triggered.

Images