A method for quantitatively analyzing regional tectonic stress risks for tunnel engineering

Abstract

The present application relates to the technical field of geological engineering safety risk assessment, and particularly relates to a regional tectonic stress risk quantitative analysis method for tunnel engineering. The method comprises the following steps: collecting multi-source geomechanics parameters, simulating and inverting tectonic stress field and checking; building a multi-scale geomechanics model, coupling simulation of tectonic stress and gravity stress, and generating stress space basic layers; obtaining stress value and direction layers by using Kriging interpolation and inversion decomposition; calculating the angle between the maximum principal stress direction and the tunnel axis, and vector decomposition to obtain the vertical stress component; normalizing the factor, determining the weight to calculate the risk value, and finally generating a tectonic stress risk quantitative classification map and performing GIS visual rendering. The present application retains the directionality and spatial heterogeneity characteristics of the tectonic stress tensor, and builds a risk quantitative analysis mechanism for tunnel engineering, significantly improving the authenticity, resolution and engineering applicability of tectonic stress identification.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of geological engineering safety risk assessment, and particularly relates to a regional tectonic stress risk quantitative analysis method for tunnel engineering. BACKGROUND

[0002] A tectonic stress field refers to a geostress field that causes tectonic movement, generates tectonic strain and forms various tectonic traces, and is generally obtained through focal mechanism solution analysis, key tectonic site in-situ geostress testing (such as the hydraulic fracturing method, casing unclamping method and three-dimensional stress tensor ultrasonic measurement method) and numerical simulation inversion analysis. Traditional tectonic stress field analysis follows a "point-domain" progressive reconstruction method, in which discrete point data of geostress is obtained through in-situ testing means such as the hydraulic fracturing method and casing unclamping method, the regional stress field is reconstructed in combination with physical experiments or numerical inversion, and finally a tectonic stress factor layer is generated and used for regional stability evaluation. However, the existing method simplifies the three-dimensional geostress tensor into a scalar, such as the maximum principal stress intensity, which only reflects the stress size and ignores the relationship between the directionality and the engineering axis, resulting in systematic compression of the anisotropy characteristics of the tectonic stress, insufficient spatial heterogeneity of the evaluation results and risks such as underestimation of tunnel asymmetric damage and rock burst misjudgment. SUMMARY

[0003] Therefore, it is necessary to provide a regional tectonic stress risk quantitative analysis method for tunnel engineering to solve at least one of the above technical problems.

[0004] To achieve the above-mentioned purpose, a regional tectonic stress risk quantitative analysis method for tunnel engineering comprises the following steps:

[0005] Step S1: Collecting a multi-source geomechanical parameter data set of a target tunnel engineering region; simulating and inverting a regional tectonic stress field based on the multi-source geomechanical parameter data set and performing a measured check to obtain a tectonic stress tensor;

[0006] Step S2: Analyzing the spatial scale and geological unit of the target tunnel engineering region and constructing a regional multi-scale geomechanical simulation model;

[0007] Step S3: Numerically coupling and simulating the tectonic stress tensor and the gravity stress tensor based on the regional multi-scale geomechanical simulation model, and structurally processing the numerical coupling and simulation results to obtain a tectonic stress space basic layer;

[0008] Step S4: Kriging interpolation of the principal stress discrete points in the tectonic stress space basic layer to generate a continuous grid layer; inversion of the continuous grid layer and decomposition of the inversion results to obtain a stress value layer and a direction layer;

[0009] Step S5: Calculate the included angle between the maximum principal stress direction and the tunnel axis based on the direction layer and the preset tunnel route trend layer of the to-be-constructed tunnel, to obtain an included angle layer; perform vector decomposition on the stress value layer and the included angle layer, to obtain a vertical principal stress component layer;

[0010] Step S6: Perform factor standardization and normalization processing on the vertical principal stress component layer, to generate a geostress factor evaluation dataset; determine the weight coefficients of each factor based on the geostress factor evaluation dataset, and calculate a regional tectonic stress risk value, to obtain a tectonic stress risk quantitative classification layer; perform GIS visual rendering on the tectonic stress risk quantitative classification layer, to output a regional tectonic stress risk distribution map.

[0011] The present application enhances the physical real perception ability of the regional stress state by introducing the inversion and verification mechanism of the tectonic stress tensor, effectively improves the reliability and integrity of the tectonic stress field data; fuses the multi-scale geomechanical modeling system and the coordinate adaptive strategy, so that the model is closer to the spatial scale characteristics of the engineering scene, and the simulation accuracy and boundary adaptability are improved; a structured expression mode of the three-dimensional coupled stress field is constructed, the tensor-level fusion of the tectonic stress and the gravity stress is realized, and the limitation of simplifying the stress field to a scalar representation in the traditional method is broken through; the direction-retained tensor interpolation and spatial decomposition mechanism is adopted, the stress intensity information is retained, and the directional element is introduced to participate in the engineering axis coupling analysis, the spatial relationship difference between the tectonic stress and the tunnel trend is revealed, and the spatial resolution and engineering pertinence of the evaluation result are significantly improved; the vertical principal stress component is extracted through vector decomposition, a quantitative factor highly related to the engineering stability risk is constructed, and the evaluation index is more interpretable and responsive; combined with the normalization processing and the analytic hierarchy process, the multi-dimensional fusion and weight precise distribution of the risk factors are realized, the actual influence of the tectonic stress field on the potential damage mechanism of the tunnel engineering is effectively reflected; and finally, the tectonic stress risk distribution map output has high resolution and high recognition visual expression ability, provides scientific support for the tunnel routing optimization, high-risk section identification and construction risk control, and significantly improves the safety design and fine management level of the tunnel engineering in the geostress environment. BRIEF DESCRIPTION OF DRAWINGS

[0012] Other characteristics, objects and advantages of the present application will become more apparent from the following detailed description of non-restrictive embodiments made with reference to the accompanying drawings:

[0013] Figure 1 A step flowchart of a regional tectonic stress risk quantitative analysis method for tunnel engineering according to the present application;

[0014] Figure 2 A tectonic stress field analysis technical flowchart in the embodiment of the present application;

[0015] Figure 3 FIG. 1 is a schematic diagram of a GIS in an embodiment of the present application;

[0016] Figure 4 FIG. 4 is a schematic diagram of an asymmetric damage process in an embodiment of the present application;

[0017] Figure 5 FIG. 5 is a schematic diagram of a method for evaluating crustal stability considering the direction of tectonic stress field in an embodiment of the present application. DETAILED DESCRIPTION

[0018] The technical method of the present application will be described clearly and completely below in conjunction with the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor fall within the scope of protection of the present application.

[0019] In addition, the accompanying drawings are only schematic illustrations of the present application, and are not necessarily drawn to scale. The same reference numerals in the drawings denote the same or similar parts, and thus repeated descriptions thereof will be omitted. Some block diagrams shown in the drawings are functional entities, and do not necessarily correspond to physically or logically independent entities. The functional entities can be implemented in the form of software, or in one or more hardware modules or integrated circuits, or in different network and / or processor methods and / or microcontroller methods.

[0020] It should be understood that although the terms "first", "second", etc. can be used herein to describe various elements, these elements should not be limited by these terms. These terms are only used to distinguish one element from another. For example, without departing from the scope of the exemplary embodiments, a first element can be called a second element, and similarly a second element can be called a first element. The term "and / or" as used herein includes any and all combinations of one or more of the associated listed items.

[0021] To achieve the above-mentioned purpose, please refer to Figures 1 to 3 The present application provides a method for quantitative analysis of regional tectonic stress risk for tunnel engineering, which comprises the following steps:

[0022] Step S1: Collecting a multi-source geomechanical parameter data set of the target tunnel engineering region; simulating and inverting the regional tectonic stress field based on the multi-source geomechanical parameter data set, and performing measurement verification to obtain a tectonic stress tensor;

[0023] Step S2: Analyzing the spatial scale and geological unit of the target tunnel engineering region, and constructing a regional multi-scale geomechanical simulation model; Step S3: Calculating the regional tectonic stress risk based on the regional multi-scale geomechanical simulation model, and obtaining a tectonic stress risk map of the target tunnel engineering region.

[0024] Step S3: Numerically coupling simulation of the tectonic stress tensor and the gravity stress tensor based on the regional multi-scale geomechanical simulation model, and structuring the numerical coupling simulation results to obtain a tectonic stress space basic layer;

[0025] Step S4: Kriging interpolation of the principal stress discrete points in the tectonic stress space basic layer to generate a continuous grid layer; inversion of the continuous grid layer and decomposition of the inversion results to obtain a stress value layer and a direction layer;

[0026] Step S5: Calculating the angle between the maximum principal stress direction and the tunnel axis for each grid cell based on the direction layer and a preset tunnel route direction layer to obtain an angle layer; vector decomposition of the stress value layer and the angle layer to obtain a vertical principal stress component layer;

[0027] Step S6: Factor standardization and normalization processing of the vertical principal stress component layer to generate a geostress factor evaluation dataset; determining the weight coefficients of each factor based on the geostress factor evaluation dataset and calculating the regional tectonic stress risk value to obtain a tectonic stress risk quantitative classification layer; GIS visual rendering of the tectonic stress risk quantitative classification layer to output a regional tectonic stress risk distribution map.

[0028] Preferably, step S1 comprises the following steps:

[0029] Step S11: Collecting geological structure deformation data, lithology parameters, water pressure fracturing data, casing uncoring data, and GPS monitoring data of the target tunnel engineering region to obtain a multi-source geomechanical parameter dataset;

[0030] Step S12: Physical modeling of the regional tectonic stress evolution process based on the multi-source geomechanical parameter dataset and stress inversion to obtain a preliminary simulation result of the regional tectonic stress field;

[0031] Step S13: Verification of the preliminary simulation result of the tectonic stress field and residual error evaluation to obtain a tectonic stress tensor fitting residual result;

[0032] Step S14: Evaluating the reliability of the regional tectonic stress tensor based on the tectonic stress tensor fitting residual result to obtain a tectonic stress tensor.

[0033] In the embodiment of the present application, through systematic collection of tectonic deformation data, lithology parameters, hydraulic fracturing data, casing release data and GPS monitoring data in the target tunnel engineering area, a multi-source geomechanical parameter data set is constructed, wherein the tectonic deformation data needs to cover fault distribution, fold axial direction and geological structure combination relationship, the lithology parameters include Young's modulus, shear modulus, Poisson's ratio and density and other mechanical properties, the hydraulic fracturing data and casing release data must comply with the test procedures in the Specification for Field Test of Geomechanics, and directional drilling layout and hierarchical pressure increase scheme are adopted to ensure the accuracy of the calculation of the principal stress direction, and the GPS monitoring data should include horizontal and vertical displacement rate information for nearly three years, and spatial registration of the observation points is completed; based on the above multi-source geomechanical parameter data set, a physical model of the tectonic stress evolution process is constructed, a stress zoning method is used to divide the three-layer structure of the upper, middle and lower crust, the zoning boundary conditions are set according to the block boundary and the main fault distribution, the regional tectonic stress field is solved by the static equilibrium equation set and the stress boundary conditions, and the preliminary simulation result is obtained, which is output as a stress field data set in the form of a three-dimensional tensor; the preliminary simulation result is compared with the measured stress tensor values obtained by hydraulic fracturing and casing release point by point, the fitting error between the stress tensors of each measuring point is calculated using the Frobenius norm, and the 0.3 times normalized residual is used as the fitting accuracy threshold, and the tectonic stress tensor fitting residual result is output; the stress tensor reliability is evaluated according to the residual result distribution, for the area with fitting residual less than the threshold, the simulation tensor is retained as the effective result, and for the area with residual exceeding the limit, the tensor is corrected using local weighted interpolation, and finally the tectonic stress tensor with complete structure and controllable accuracy is obtained.

[0034] Preferably, step S2 comprises the following steps:

[0035] Step S21: identifying the spatial scale of the block distribution, latitude and longitude span and proposed route range of the target tunnel engineering area, to obtain a regional scale determination result;

[0036] Step S22: setting a model coordinate system based on the regional scale determination result, to obtain a regional mechanics model coordinate framework;

[0037] Step S23: analyzing the geological profile data, fault structure information and upper, middle and lower crust layering information of the target tunnel engineering area based on the regional mechanics model coordinate framework, and structurally modeling the region, to obtain a spatial division model;

[0038] Step S24: assigning corresponding physical properties to different geological units based on the spatial division model and tectonic stress tensor, and setting velocity constraint boundary conditions, to obtain a regional multi-scale geomechanical simulation model.

[0039] In the embodiment of the present application, the number of plots, boundary distribution, latitude and longitude span and projection length of the tunnel line in the target tunnel engineering area are comprehensively calculated by importing the plot boundary vector data of the target tunnel engineering area and the start and end coordinates of the proposed route into the ArcGIS platform, wherein the latitude and longitude span needs to be accurate to 0.001 If the latitude and longitude span is greater than 1.5 or the engineering route spans multiple first-order structural units, the area is identified as a large-scale engineering area, and the area scale determination result is output for coordinate system selection; the model coordinate system is set according to the area scale determination result, the spherical shell coordinate system is set if the determination result is a large-scale area, the spatial projection conversion is performed with the WGS-84 earth ellipsoid parameters, the plane rectangular coordinate system is set if the determination result is a small-scale area, and the unified Gauss projection zone coordinates are adopted, and finally the area mechanics model coordinate framework for three-dimensional tectonic mechanics analysis is obtained; after the coordinate framework is determined, the upper crust (0-15 km), the middle crust (15-30 km) and the lower crust (30-45 km) are divided into three layers by integrating the existing geological digital data, the layered boundary is set based on the seismic tomography and Moho depth data, the tectonic contact surface is constrained by the fault dislocation line and the fold axis surface, the geological unit boundary surface and volume unit are constructed by a three-dimensional modeling tool such as GOCAD, and a complete spatial division model is generated; based on the spatial division model and the tectonic stress tensor data, the mechanical parameters of each geological unit are assigned, wherein the Young's modulus is limited between 10-90 GPa, the shear modulus is limited between 5-45 GPa, the Poisson's ratio is set to 0.20-0.30, the rock density is set to 2.5-2.9 , and the model boundary velocity constraint is set according to the regional plate motion monitoring data, wherein the horizontal velocity of the north-south boundary is set to ±25 mm / a, and the vertical boundary is set to 0 displacement condition, and finally the regional multi-scale geomechanical simulation model satisfying the continuity of structure, coordinate and physical property is established.

[0040] Preferably, step S3 comprises the following steps:

[0041] Step S31: tensor synthesis is performed on the tectonic stress tensor and the gravity stress tensor, and physical coupling of the multi-source stress field is performed to obtain a coupled stress tensor field;

[0042] Step S32: numerical simulation is performed based on the regional multi-scale geomechanical simulation model and the coupled stress tensor field, stress propagation and deformation response are calculated, and a shallow coupled crustal stress numerical simulation result is obtained;

[0043] Step S33: residual calculation is performed on the shallow coupled crustal stress numerical simulation result, stress tensor data is extracted based on the residual calculation result, spatial coding and layer decomposition are performed according to the three-dimensional grid structure, and a tectonic stress spatial basic layer is obtained.

[0044] In the embodiment of the present application, the tectonic stress tensor and the gravity stress tensor are taken as inputs, the point-by-point summation processing of the two tensor fields is performed by using the tensor superposition principle, the gravity stress tensor is calculated according to , wherein is the rock density, is the gravity acceleration (9.81 ), is the depth, and the density parameter is valued between 2.5-2.9 according to different geological units. The calculation result maintains the spatial structure corresponding to the three-dimensional regional grid one by one. The Frobenius inner product form is used for tensor synthesis. The stress tensor of each spatial unit is synthesized into six components , , , , , to obtain the coupled stress tensor field. The regional multi-scale geomechanical simulation model is taken as the calculation domain, the coupled stress tensor field is taken as the initial stress field, the explicit difference method is applied for static simulation, the distribution response and deformation influence of the coupled stress in the shallow region (depth 0-1000 m) are calculated, the vertical layering structure is adopted in the simulation process, the stress conduction solution is calculated for each layer node, and multiple groups of indexes including the principal stress, shear stress and volumetric strain are output to form the numerical simulation result of the shallow coupled geostress; the simulation output result is compared with the original tectonic stress tensor field point by point, the residual error is calculated by using the Frobenius norm formula, the tensor data with the residual error less than 0.3 is taken as the effective value, the three-dimensional space grid (horizontal 500 m × 500 m, vertical 100 m) is divided at equal intervals, the tensor is spatially encoded, the field is marked , , and the direction angle, and the multi-channel layer structure is generated for each spatial unit by using the raster grid point mapping tool in ArcGIS, and finally the tectonic stress spatial basic layer is output.

[0045] Preferably, step S4 comprises the following steps:

[0046] Step S41: vector extraction is performed on the principal stress discrete points in the tectonic stress spatial basic layer to obtain a principal stress discrete point data set;

[0047] Step S42: spatial interpolation is performed on the principal stress size in the principal stress discrete point data set by using the Kriging method to obtain a continuous raster layer of the principal stress value;

[0048] Step S43: the main stress direction vector information in the main stress discrete point data set is directionally interpolated by using the weighted Kriging interpolation method, and a continuous grid layer of the main stress direction is obtained.

[0049] Step S44: the continuous grid layer of the main stress value and the continuous grid layer of the main stress direction are inversely synthesized to obtain a main stress tensor field layer.

[0050] Step S45: the tensor field is directionally decomposed based on the tensor principal axis direction information of each grid unit in the main stress tensor field layer, and a tensor modulus length component is extracted to obtain a stress value layer and a direction layer.

[0051] In the embodiment of the application, the main stress discrete points in the tectonic stress space basic layer are processed by using the vector field extraction tool of the ArcGIS platform, the three-dimensional direction vector components of the maximum principal stress in the tensor field (σ1) , , , ) and the main stress size value are called, the main stress discrete point data set is output according to the spatial coordinates of each point, the data set includes the tensor principal value and the corresponding direction information; based on the stress size field in the main stress discrete point data set, the ordinary Kriging method in the Geostatistical Analyst module is selected for spatial interpolation calculation, the search neighborhood radius is set to 1500 m, the minimum number of neighboring points is set to 10, the interpolation range covers the tunnel line control range plus the boundary buffer area of 1000 m, the continuous grid layer of the main stress value is output, and the grid unit size is set to 100 m*100 m; the main stress direction information is interpolated, each main stress direction vector is first normalized to a unit vector, then the three direction components (x, y, z) , , ) are respectively interpolated by using the weighted Kriging interpolation, the weight size is set according to the tensor modulus length, the multivariate joint interpolation tool of the Geostatistical Analyst is used to generate the continuous grid layers of the three direction components, then the direction field is restored to the complete main stress direction vector in each grid unit by using the three-dimensional vector synthesis method, and the continuous grid layer of the main stress direction is formed; the continuous grid layer of the main stress value and the continuous grid layer of the main stress direction are inversely synthesized point by point, and the value and the direction vector three components are restored in each grid. The spatial distribution of the tensor yields the principal stress tensor field layer, which is saved as multi-band floating-point raster data containing the magnitude and direction of the principal stresses. Using the principal stress tensor field layer as input, the principal axis direction information of each raster cell is extracted, the tensor vector magnitude is calculated, and the inverse cosine solution of the principal axis direction is performed to extract the maximum principal stress value as the stress value layer. At the same time, the principal axis direction angles (tilt and azimuth) are extracted to construct the direction layer, ultimately forming the stress value layer and direction layer for geostress direction analysis and engineering axis coupling calculation.

[0052] Preferably, step S5 includes the following steps:

[0053] In step S5, the angle between the direction of the maximum principal stress and the tunnel axis is calculated for each grid cell based on the direction layer and the preset proposed tunnel alignment layer.

[0054] Coordinate vectors are extracted from the direction layer, and the extracted results are standardized into principal stress direction angle layers within the range of [0°~360°].

[0055] The preset proposed tunnel route alignment layer is vectorized, and the route direction value of each grid position is extracted according to the tunnel axis direction to obtain the route alignment angle layer.

[0056] Based on the principal stress direction angle layer and the line direction angle layer, the direction difference is calculated for each grid cell to obtain the angle layer between the maximum principal stress direction and the line direction.

[0057] In this embodiment of the invention, firstly, coordinate vectors are extracted from the direction layer, and the extraction method is based on the three-dimensional components of the principal stress directions ( , Construct a horizontal projection vector in a two-dimensional plane, and use the Field Calculator tool in ArcGIS to calculate the azimuth angle of the direction vector relative to true north in each grid cell. Then, the Raster Calculator module was used to convert negative values ​​into standard orientation angles between 0° and 360°, generating a principal stress orientation angle layer. Next, the preset proposed tunnel alignment layer was vectorized. First, the alignment layer was segmented at 100m intervals. The starting and ending coordinates of each segment were extracted, and its orientation vector on the horizontal plane was calculated. Then, the orientation angle of each segment relative to true north was calculated. Use the Join Field tool to Values ​​are assigned to each spatially overlapping raster cell to generate a route orientation angle layer, with output values ​​limited to the range of 0° to 180°. Finally, the direction difference between the two angle layers is calculated using the MapAlgebra expression in ArcGIS, with the formula as follows: , if greater than 90° is replaced by 180° to ensure that the included angle is always within the range of 0°~90°, and the final output included angle layer is used as the basic data for analyzing the geometric relationship between the ground stress direction and the tunnel trend.

[0058] Preferably, the vector decomposition of the maximum principal stress value layer and the included angle layer in step S5 includes:

[0059] The included angle layer is subjected to a sine value conversion process to obtain a direction included angle sine layer;

[0060] The stress value layer is subjected to a standardization process to convert the original principal stress value into a vector module length form to obtain a standard principal stress module length layer;

[0061] Based on the standard principal stress module length layer and the direction included angle sine layer, the principal stress vector is directionally decomposed to extract a stress vector component layer perpendicular to the tunnel line trend, and the stress vector component layer is recorded as a vertical principal stress component layer.

[0062] In the embodiment of the present application, first, the included angle layer is subjected to a sine value conversion process, the Raster Calculator tool is called in the ArcGIS platform, the angle value of each grid cell in the included angle layer is converted from to radian value , and the expression is applied to obtain the included angle sine value layer, the value range of which is limited between 0 and 1, and is used to represent the orthogonal component proportion between the principal stress direction and the tunnel trend; secondly, the maximum principal stress value layer is subjected to a standardization process, the original stress value range is set to to (wherein , ), the linear normalization formula is used to standardize the original principal stress value to a dimensionless principal stress module length layer between 0 and 1, and the Reclassify tool is used to ensure that all abnormal values or null values are removed or filled to eliminate the influence of boundary interference on subsequent decomposition calculation; finally, the standard principal stress module length layer and the included angle sine value layer are subjected to a grid-by-grid multiplication operation, the expression is used to obtain the vertical principal stress component layer.

[0063] Preferably, the factor standardization and normalization process of the vertical principal stress component layer in step S6 includes:

[0064] The vertical principal stress component layer is subjected to a spatial reclassification process to divide the continuous values into a set of level intervals to obtain a vertical principal stress level layer;

[0065] obtain a lithology type layer of a target tunnel engineering area, tunnel buried depth data, a regional fault distribution layer and regional terrain data;

[0066] extract engineering response coefficients of different lithology categories in the lithology type layer, and perform encoding conversion to obtain a lithology response factor layer;

[0067] extract tunnel axis depth information in the regional terrain data, and perform spatial interpolation in combination with the tunnel buried depth data to obtain a buried depth factor layer;

[0068] perform buffer zone analysis on the regional fault distribution layer to identify a fault adjacent area risk zone, and obtain a tectonic activity zone factor layer;

[0069] perform discrete classification processing on the included angle layer, divide control levels according to the discrete classification processing result, and obtain an included angle control factor layer;

[0070] perform linear standardization processing on each layer value based on the vertical principal stress level layer, the lithology response factor layer, the buried depth factor layer, the tectonic activity zone factor layer and the included angle control factor layer to obtain a geostress factor evaluation data set.

[0071] In the embodiment of the present application, first, the Reclassify tool of ArcGIS is used to perform spatial reclassification processing on the vertical principal stress component layer, and the numerical range thereof is divided into five interval ranges, such as 0-0.2, 0.2-0.4, 0.4-0.6, 0.6-0.8 and 0.8-1.0, according to the threshold value of engineering experience, and the interval codes 1 to 5 are assigned, and a vertical principal stress level layer is output, reflecting the relative strength of the principal controlling lateral stress level in different regions; secondly, the lithology type layer, the tunnel buried depth data, the regional fault distribution layer and the regional terrain data of the target tunnel engineering region are obtained, wherein the lithology type layer is derived from the 1:50000 geological map digitization data, the Field Calculator is called to assign the coding to the lithology field, for example, hard rock (such as granite and diabase) is assigned as 1, medium hard rock (such as limestone and sandstone) is assigned as 2, and soft rock (such as shale and mudstone) is assigned as 3, and a lithology response level table is established in combination with the existing rock elastic modulus and shear modulus test results to generate a lithology response factor layer; the tunnel buried depth data is spatially interpolated to generate a terrain continuous surface by using the Topo to Raster method, and the difference between the minimum elevation of the tunnel axis on the terrain surface and the design buried depth data is extracted by using the Extract by Mask tool in ArcGIS Spatial Analyst, and a buried depth factor layer is output, with the value range set to 50 m to 800 m; the regional fault distribution layer is subjected to buffer zone analysis, the first-order fault buffer width is set to 1000 m and the second-order fault buffer width is set to 500 m according to the fault activity level, and a tectonic adjacent area is extracted by superposition to generate a tectonic activity belt factor layer; the included angle layer obtained in the previous step is discretely graded according to the range of 0°-90°, and five included angle levels are divided according to 0-15°, 15-30°, 30-45°, 45-60° and 60-90°, and the control levels 1 to 5 are assigned, and an included angle control factor layer is output; finally, the value ranges of the above five types of factor layers are uniformly converted to dimensionless data of 0 to 1 by using the linear standardization method, and the standardization calculation is performed on each layer by using the Raster Calculator, and finally the geostress factor evaluation data set composed of the vertical principal stress level layer, the lithology response factor layer, the buried depth factor layer, the tectonic activity belt factor layer and the included angle control factor layer is formed. =0°–90° range, and five included angle levels are divided according to 0-15°, 15-30°, 30-45°, 45-60° and 60-90°, and the control levels 1 to 5 are assigned, and an included angle control factor layer is output; finally, the value ranges of the above five types of factor layers are uniformly converted to dimensionless data of 0 to 1 by using the linear standardization method, and the standardization calculation is performed on each layer by using the Raster Calculator, and finally the geostress factor evaluation data set composed of the vertical principal stress level layer, the lithology response factor layer, the buried depth factor layer, the tectonic activity belt factor layer and the included angle control factor layer is formed.

[0072] Preferably, the step S6 of determining the weight coefficients of each factor based on the geostress factor evaluation data set and calculating the regional tectonic stress risk value comprises:

[0073] The weight influence factor scoring is performed on each factor layer in the geostress factor evaluation data set, a pair comparison matrix is constructed, and an analytic hierarchy process input matrix is obtained;

[0074] Eigenvalues and eigenvectors of each risk factor are calculated based on the input matrix of the analytic hierarchy process, consistency ratio is judged whether to meet the set threshold, and a weight coefficient set is obtained;

[0075] The regional tectonic stress risk index raster layer is constructed by weighted superposition calculation of the tectonic stress factor evaluation data set and the weight coefficient set;

[0076] The tectonic stress risk quantitative grading layer is obtained by reclassifying the tectonic stress risk index raster layer according to the preset risk level division standard.

[0077] In the embodiment of the present application, first, the weight influence factor scoring of each factor layer in the tectonic stress factor evaluation data set is performed, the expert scoring weight is set for five layers, i.e., the vertical principal stress level layer, the lithology response factor layer, the burial depth factor layer, the tectonic activity zone factor layer and the angle control factor layer, according to the influence degree of engineering mechanics, the scoring value represents the relative importance between factors in the pair ratio system of 1 to 9, a 5x5 pair comparison matrix is constructed, the diagonal elements in the matrix are 1, and the reciprocal elements are the reciprocal, which is used as the input matrix of the analytic hierarchy process; secondly, the maximum eigenvalue of the matrix is calculated by the eigenvalue method , and the corresponding eigenvector is extracted , wherein each represents the normalized weight coefficient of the corresponding factor, and the consistency ratio is used for consistency test, , RI is taken from the Saaty consistency random index table, when CR 0.1, the matrix consistency is determined to pass, and the weight coefficient set is output; then, each layer in the tectonic stress factor evaluation data set is weighted and superimposed pixel by pixel according to the corresponding weight by using the Raster Calculator tool of ArcGIS, the regional tectonic stress risk index raster layer is generated, and the layer value range is limited to 0 to 1; finally, the tectonic stress risk index raster layer is reclassified, the threshold interval is divided according to the preset five-level risk level standard, the grading interval is set as [0-0.2), [0.2-0.4), [0.4-0.6), [0.6-0.8), [0.8-1.0] according to experience, and the level identification 1 to 5 is given, which respectively corresponds to the low risk area, the lower risk area, the medium risk area, the higher risk area and the high risk area, and the tectonic stress risk quantitative grading layer is output.

[0078] Preferably, the GIS visualization rendering of the tectonic stress risk quantitative grading layer in step S6 comprises:

[0079] The color mapping setting is performed on the tectonic stress risk quantitative grading layer, the hierarchical color scale is configured for different risk levels, and the rendered tectonic stress risk layer is obtained.

[0080] The legend generation and map finishing are performed on the rendered tectonic stress risk layer to obtain a tectonic stress risk visualization map draft.

[0081] The high-risk area extraction is performed on the tectonic stress risk distribution map, the typical stress mutation zone and the tectonic active section are marked, and the regional tectonic stress risk distribution map is obtained.

[0082] In the embodiment of the present application, firstly, the tectonic stress risk quantitative grading layer is opened in the ArcGIS platform, the Symbology module is called to perform color mapping setting on the layer, the classification rendering is performed according to the risk level field assigned in the grading layer, a five-level color scale scheme is set, wherein the low-risk area is assigned a color of light green, the lower-risk area is assigned a color of yellow, the medium-risk area is assigned a color of orange, the higher-risk area is assigned a color of red, and the high-risk area is assigned a color of dark red, so as to ensure that the color contrast is clear and continuous, the rendered tectonic stress risk layer is output, and the grading field is bound in the layer attribute to ensure the dynamic linkage of the map; secondly, the Layout View function is used to perform map finishing on the rendered layer, the compass, the scale, the latitude and longitude network, the legend box and the title are added in the drawing, the legend is composed of the hierarchical color band and the corresponding level text description, the legend position is set to be fixed at the lower right corner of the drawing in the Legend Properties, the map border line width is set to 0.5 pt, and the map title is uniformly named as the regional tectonic stress risk visualization map, and the tectonic stress risk visualization map draft is output; finally, the space analysis processing is performed on the map, the high-risk raster unit with the grading value of 4 and 5 is extracted by calling the Raster Calculator, the continuous area is identified by using the Region Group tool, the boundary value variation range of the area is calculated by using the Zonal Statistics, the boundary paragraph with the gradient change rate greater than 0.4 is extracted as the typical stress mutation zone, the fault layer loaded in step S21 is combined, the spatial overlap area of the high-risk area and the fault active zone is identified by using the Intersect tool, and the text marking is performed, and finally the regional tectonic stress risk distribution map containing the risk grading, the typical mutation zone and the tectonic active section annotation information is generated.

[0083] Therefore, from any point of view, the embodiments should be regarded as exemplary and non-limiting, the scope of the present application is not limited by the above description, and therefore all changes falling within the meaning and scope of the equivalent elements of the application file are intended to be included in the present application.

[0084] The foregoing is considered as illustrative only of the principles of the application. Numerous modifications and changes will readily occur to those skilled in the art, and it is intended to embrace all such modifications and changes that fall within the scope of the application. Accordingly, the application is not to be restricted in scope to the specific embodiments disclosed herein but is to be accorded the full scope that the principles and novel features request appropriately granted.

Claims

1. A method for quantitatively analyzing regional tectonic stress risks for tunnel engineering, characterized in that, The method comprises the following steps: Step S1: collecting a multi-source geomechanical parameter data set of a target tunnel engineering area; Based on the multi-source geomechanical parameter data set, the regional tectonic stress field is simulated and inverted, and the measured calibration is performed to obtain the tectonic stress tensor; Step S2: analyzing the spatial scale and geological unit of the target tunnel engineering area, and constructing a regional multi-scale geomechanical simulation model; Step S3: based on the regional multi-scale geomechanical simulation model, the tectonic stress tensor and the gravity stress tensor are numerically coupled and simulated, and the numerical coupling simulation results are structured to obtain a tectonic stress space basic layer; Step S4: Kriging interpolation is performed on the principal stress discrete points in the tectonic stress space basic layer to generate a continuous grid layer; The continuous grid layer is inverted and the inversion results are decomposed to obtain a stress value layer and a direction layer; Step S5: based on the direction layer and a preset tunnel route direction layer, the angle between the maximum principal stress direction and the tunnel axis of each grid unit is calculated to obtain an angle layer; the stress value layer and the angle layer are vector decomposed to obtain a vertical principal stress component layer; Step S6: the vertical principal stress component layer is factor standardized and normalized to generate a geostress factor evaluation data set; Based on the geostress factor evaluation data set, the weight coefficients of each factor are determined, and the regional tectonic stress risk value is calculated to obtain a tectonic stress risk quantitative classification layer; The tectonic stress risk quantitative classification layer is rendered by GIS visualization to output a regional tectonic stress risk distribution map.

2. The method of quantitative analysis of regional tectonic stress risks for tunnel engineering according to claim 1, characterized in that, Step S1 comprises the following steps: Step S11: collecting geological structure deformation data, lithology parameters, water pressure fracturing data, casing core release data and GPS monitoring data of the target tunnel engineering area to obtain a multi-source geomechanical parameter data set; Step S12: based on the multi-source geomechanical parameter data set, a physical model of the regional tectonic stress evolution process is established, and stress inversion is performed to obtain a preliminary simulation result of the regional tectonic stress field; Step S13: the preliminary simulation result of the tectonic stress field is verified, and residual error is evaluated to obtain a fitting residual result of the tectonic stress tensor; Step S14: based on the fitting residual result of the tectonic stress tensor, the reliability of the regional tectonic stress tensor is evaluated to obtain the tectonic stress tensor.

3. The method of quantitative analysis of regional tectonic stress risks for tunneling according to claim 1, characterized in that, Step S2 comprises the following steps: Step S21: identifying the spatial scale of the block distribution, latitude and longitude span and the range of the proposed route of the target tunnel engineering area to obtain a regional scale determination result; Step S22: setting a model coordinate system based on the regional scale determination result to obtain a regional mechanical model coordinate framework; Step S23: based on the regional mechanical model coordinate framework, the geological profile data, fault structure information and upper, middle and lower crust layering information of the target tunnel engineering area are analyzed, and a structure model of the region is established to obtain a spatial division model; Step S24: based on the spatial division model and the tectonic stress tensor, the corresponding physical properties are given to different geological units, and the velocity constraint boundary conditions are set to obtain a regional multi-scale geomechanical simulation model.

4. The method for quantitatively analyzing regional tectonic stress risks for tunnel engineering according to claim 1, wherein, Step S3 comprises the following steps: Step S31: tensor synthesis is performed on the tectonic stress tensor and the gravity stress tensor, and physical coupling of the multi-source stress field is performed to obtain a coupled stress tensor field; Step S32: numerical simulation is performed based on the regional multi-scale geomechanical simulation model and the coupled stress tensor field, stress propagation and deformation response are calculated, and a shallow coupled crustal stress numerical simulation result is obtained; Step S33: residual calculation is performed on the shallow coupled crustal stress numerical simulation result, stress tensor data is extracted based on the residual calculation result, spatial coding and layer decomposition are performed according to a three-dimensional grid structure, and a tectonic stress space basic layer is obtained.

5. The method for quantitatively analyzing regional tectonic stress risks for tunnel engineering according to claim 1, characterized in that, Step S4 includes the following steps: Step S41: vector extraction is performed on the principal stress discrete points in the tectonic stress space basic layer to obtain a principal stress discrete point data set; Step S42: spatial interpolation is performed on the principal stress size in the principal stress discrete point data set by using the Kriging method to obtain a continuous grid layer of the principal stress value; Step S43: directional interpolation processing is performed on the principal stress direction vector information in the principal stress discrete point data set by using the weighted Kriging interpolation method to obtain a continuous grid layer of the principal stress direction; Step S44: the continuous grid layer of the principal stress value and the continuous grid layer of the principal stress direction are inversely synthesized to obtain a principal stress tensor field layer; Step S45: the tensor field is directionally decomposed based on the tensor principal axis direction information of each grid cell in the principal stress tensor field layer, the tensor module length component is extracted, and a stress value layer and a direction layer are obtained.

6. The method of quantitative analysis of regional tectonic stress risks for tunneling according to claim 1, characterized in that, The angle between the maximum principal stress direction and the tunnel axis is calculated for each grid cell based on the direction layer and the preset tunnel route direction layer in step S5, which includes: Coordinate vector extraction is performed on the direction layer, and the extraction result is standardized to a principal stress direction angle layer within [0°~360°]; Vectorization processing is performed on the preset tunnel route direction layer, the route direction value of each grid position is extracted according to the tunnel axis direction, and a route direction angle layer is obtained; Direction difference calculation is performed on each grid cell based on the principal stress direction angle layer and the route direction angle layer to obtain a maximum principal stress direction and route direction angle layer.

7. The method of quantitative analysis of regional tectonic stress risks for tunneling according to claim 1, characterized in that, The maximum principal stress value layer and the angle layer are vector decomposed in step S5, which includes: Sinusoidal value conversion processing is performed on the angle layer to obtain a direction angle sine layer; Standardization processing is performed on the stress value layer to convert the original principal stress value to a vector module length form to obtain a standard principal stress module length layer; Directional decomposition is performed on the principal stress vector based on the standard principal stress module length layer and the direction angle sine layer, a stress vector component layer perpendicular to the tunnel route direction is extracted, and the stress vector component layer is recorded as a vertical principal stress component layer.

8. The method of quantitative analysis of regional tectonic stress risks for tunneling according to claim 1, characterized in that, Factor standardization and normalization processing is performed on the vertical principal stress component layer in step S6, which includes: Spatial reclassification processing is performed on the vertical principal stress component layer to divide the continuous value into a set of level intervals to obtain a vertical principal stress level layer; Lithology type layer, tunnel burial depth data, regional fault distribution layer and regional terrain data of a target tunnel engineering area are obtained; The engineering response coefficients of different lithology categories in the lithology type layer are extracted and coded to obtain a lithology response factor layer; The depth information of the tunnel axis in the regional terrain data is extracted, and spatial interpolation is performed in combination with the tunnel burial depth data to obtain a burial depth factor layer; The buffer zone analysis is performed on the regional fault distribution layer to identify the fault adjacent zone risk belt to obtain a tectonic activity belt factor layer; The discrete classification processing is performed on the included angle layer, and the control level is divided according to the discrete classification processing result to obtain an included angle control factor layer; Based on the vertical principal stress level layer, the lithology response factor layer, the burial depth factor layer, the tectonic activity belt factor layer and the included angle control factor layer, the linear standardization processing is performed on each layer value to obtain a geostress factor evaluation data set.

9. The method of quantitative analysis of regional tectonic stress risks for tunneling according to claim 1, characterized in that, In step S6, the weight coefficients of each factor are determined based on the geostress factor evaluation data set, and the regional tectonic stress risk value is calculated, which includes: The weight influence factor scoring is performed on each factor layer in the geostress factor evaluation data set to build a pair-wise comparison matrix to obtain an analytic hierarchy process input matrix; Based on the analytic hierarchy process input matrix, the eigenvalues and eigenvectors of each risk factor are calculated, and it is judged whether the consistency ratio meets the set threshold to obtain a weight coefficient set; The weighted superposition calculation is performed on the geostress factor evaluation data set and the weight coefficient set to build a regional tectonic stress risk index raster layer; The reclassification processing is performed on the tectonic stress risk index raster layer, and the interval is divided according to the preset risk level division standard to obtain a tectonic stress risk quantitative classification layer.

10. The method of quantitative analysis of regional tectonic stress risks for tunneling according to claim 1, characterized in that, In step S6, the GIS visualization rendering is performed on the tectonic stress risk quantitative classification layer, which includes: The color mapping setting is performed on the tectonic stress risk quantitative classification layer, and the hierarchical color scale is configured for different risk levels to obtain a rendered tectonic stress risk layer; The legend generation and map finishing are performed on the rendered tectonic stress risk layer to obtain a tectonic stress risk visualization map sketch; The high-risk area extraction is performed on the tectonic stress risk distribution map, and the typical stress mutation zone and the tectonic active section are marked to obtain a regional tectonic stress risk distribution map.

Images