Non-coal mine slope stability evaluation method based on multi-working condition analysis

Abstract

The present invention discloses a non-coal mine slope stability evaluation method based on multi-operating condition analysis, comprising: obtaining a geotechnical dataset, a multi-operating condition dataset, and a monitoring dataset of the non-coal mine slope, clustering the geotechnical dataset to obtain parameter partitions, obtaining a patch model based on the parameter partitions and the multi-operating condition dataset using a patch algorithm, performing stability evaluation using the parameter partitions and the patch model using the limit equilibrium method, outputting a first stability result of the slope, screening out relatively dangerous areas based on the first stability result of the slope, performing stability evaluation on the relatively dangerous areas using the finite element method, outputting a second stability result of the slope, and using the second stability result as the final evaluation result. This method provides efficient and reliable evaluation results for non-coal mine slope disaster prevention and control by combining a double-patch model with a step-by-step verification process of the limit equilibrium method and the finite element method.

Description

Technical Field

[0001] The present invention belongs to the field of stability evaluation, and in particular relates to a non-coal mine slope stability evaluation method based on multi-working condition analysis. Background Art

[0002] With the mining and operation of non-coal mines, slope stability is a key factor in ensuring mine safety, personnel safety and surrounding environmental safety. With the continuous expansion of the scale of non-coal mining in my country and the gradual increase in mining depth, the geological conditions and stress conditions of the slopes are becoming more and more complex. Slope instability accidents occur frequently, causing huge losses to mining companies and society. Therefore, accurately evaluating the stability of non-coal mine slopes under multiple working conditions has important practical significance for mine planning, design, mining and safety management.

[0003] Traditional non-coal mine slope stability evaluation methods have certain limitations when faced with complex geotechnical data sets, multi-operating condition data sets and monitoring data sets. On the one hand, a single clustering method is often used to process geotechnical data, which fails to fully consider the parameter similarity and spatial proximity of data points, resulting in inaccurate parameter partitioning and an inability to accurately reflect the actual characteristics of the slope rock and soil. On the other hand, when considering the impact of multiple operating conditions on slope stability, there is a lack of effective methods for reasonable setting and comprehensive analysis of different operating conditions. At the same time, it is easy to ignore the key role of local special geological areas in slope stability, which affects the accuracy and reliability of the evaluation results.

[0004] In order to overcome the shortcomings of traditional methods and meet the actual needs of modern non-coal mine slope stability evaluation, a stability evaluation method is needed that can comprehensively process multi-source data, accurately divide parameter zones, fully consider the influence of multiple working conditions, and effectively identify local special geological zones. Against this background, a non-coal mine slope stability evaluation method based on multi-working condition analysis is proposed. The geotechnical data are processed by system clustering and fuzzy kernel clustering, the working condition patches are set using the patch algorithm, and the limit equilibrium method and the finite element method are combined to perform multi-stage stability evaluation to improve the accuracy and reliability of the evaluation results and provide a more scientific basis for the safety management of non-coal mine slopes. Summary of the Invention

[0005] The purpose of the present invention is to provide a non-coal mine slope stability evaluation method based on multi-operating condition analysis.

[0006] To achieve the above object, the present invention is implemented according to the following technical solutions:

[0007] The first aspect of the present invention provides a non-coal mine slope stability evaluation method based on multi-operating condition analysis, comprising:

[0008] Collecting non-coal mine slope data in a preset area and preprocessing the non-coal mine slope data; the non-coal mine slope data includes a geotechnical data set, a multi-working condition data set, and a monitoring data set;

[0009] Performing systematic clustering on the geotechnical dataset to obtain a plurality of first categories, dividing regions of the geotechnical dataset with similar parameters and spatial proximity into second categories based on the first categories by fuzzy kernel clustering, and using the second categories as parameter partitions;

[0010] Based on the parameter partitioning and the multi-operating condition data set, a first patch is set for different operating conditions using a patch algorithm; local special geological areas under different operating conditions are extracted according to the monitoring data set; a second patch is generated for the local special geological areas; and a patch model is obtained by combining the first patch and the second patch;

[0011] Using the parameter partition and the patch model to perform stability evaluation through a limit equilibrium method, and outputting a first stability result of the slope;

[0012] A relatively dangerous area is screened out according to the first stability result of the slope, stability evaluation is performed on the relatively dangerous area using the finite element method, a second stability result of the slope is output, and the second stability result is used as the final evaluation result.

[0013] As a further method, the pretreatment method comprises:

[0014] The density, cohesion, internal friction angle, elastic modulus and Poisson's ratio of different positions on non-coal mine slopes under saturated conditions are collected as geotechnical parameters; the fault zone distribution, homogeneous rock formations, special structural rock formations and occurrence of non-coal mine slopes are collected as geological layer data; the meteorological historical records, seismic network data, groundwater level information, mine excavation design and construction records of non-coal mine slope areas are collected; rainfall intensity, earthquake acceleration, groundwater depth and excavation load are obtained as multi-condition data; the slope displacement monitoring values, stress monitoring values ​​and pore water pressure monitoring values ​​are obtained as monitoring data; the geotechnical parameters, geological layer data, multi-condition data and monitoring data are formatted in a unified manner, time axis aligned, outliers eliminated and missing values ​​interpolated to obtain geotechnical data sets, multi-condition data sets and monitoring data sets.

[0015] As a further method, a method for performing systematic clustering on the geotechnical dataset to obtain a plurality of first categories includes:

[0016] A three-dimensional model is constructed based on the spatial location of the data points in the geotechnical dataset, geological layer data, and geotechnical parameters. The geotechnical parameters of the data points are used as clustering dimensions according to the three-dimensional model, and the Euclidean distance matrix is ​​calculated for the data points using the clustering dimensions.

[0017] Based on the sampling spatial position of the data points, the data points are connected into a triangular grid through triangulation, and the distance between the data points in the triangular grid is used as the spatial proximity weight. The spatial proximity weight is multiplied by the Euclidean distance of the Euclidean distance matrix to obtain the Euclidean distance matrix of the superimposed spatial proximity;

[0018] Hierarchical clustering is performed based on the Euclidean distance matrix of superimposed spatial proximity. The data points are taken as the initial class and the data points are gradually merged until the increment of the sum of squared deviations within the class after merging is minimized, and a clustering hierarchical tree is obtained.

[0019] The silhouette coefficient of the data points is calculated based on the clustering hierarchy tree, and the mean of the silhouette coefficient corresponding to different number of categories is obtained. The number of categories is traversed and the number of categories with the largest mean of the silhouette coefficient is selected as the final number of categories. The geotechnical data set is divided into multiple first categories according to the final number of categories.

[0020] The first type of data is imported into the geographic information system, and the dominant direction of the joint is determined according to the occurrence data corresponding to the geotechnical parameters. Kriging interpolation is performed on the geotechnical parameters based on the dominant direction of the joint. During the interpolation, the variation function is used to distinguish the range and base value of different joint dominant directions to obtain a heterogeneous parameter field.

[0021] As a further method, based on the first category, a method for dividing the regions of the geotechnical dataset with similar data point parameters and spatial proximity into a second category by fuzzy kernel clustering includes:

[0022] Based on the first type of geotechnical parameter vectors and spatial coordinate vectors of data points extracted in the heterogeneous parameter field, a multidimensional feature space containing parameter characteristics and spatial positions is obtained;

[0023] According to the multidimensional feature space, a kernel mapping matrix is ​​constructed by using the Gaussian kernel function. The initial fuzzy membership matrix and the initial cluster center matrix are randomly generated based on the kernel mapping matrix. The objective function is constructed according to the kernel mapping matrix, the initial fuzzy membership matrix and the initial cluster center matrix, which can be expressed as:

[0024]

[0025] Where J(U,V) is the objective function, U is the fuzzy membership matrix, V is the cluster center matrix, c is the number of cluster centers, n is the total number of data points, u ik The fuzzy membership of data point i to cluster center k, m is the fuzzy factor, p i Geomechanical parameter vector of data point i, l iThe spatial coordinate vector of data point i, h1 is the variance, which is the scale of kernel space to measure the similarity of geotechnical parameters. The smaller the value, the faster the parameter similarity decays as the parameter difference increases. Clustering is more sensitive to parameter differences. The larger the value, the data points with large parameter differences can still have high similarity. h2 is the variance, which is the scale of kernel space to measure the proximity of spatial coordinates. The smaller the value, the faster the spatial proximity decreases as the distance increases. Clustering emphasizes spatial compactness. The larger the value, the data points with farther distances can be classified into one category due to parameter similarity. w1 is the similarity weight coefficient of geotechnical parameters, w2 is the spatial proximity weight coefficient, w1+w2=1, d ik is the spatial Euclidean distance between data point i and cluster center k;

[0026] The Lagrange multiplier method is used to solve the minimization problem of the objective function. According to the update formula derived from the Lagrange multiplier method, the kernel mapping matrix, the current fuzzy membership matrix and the cluster center matrix are used to iterate the membership and cluster center through the update formula until the difference between the objective function values ​​of two adjacent iterations is less than 0.1%. The second category of fuzzy kernel clustering is obtained and the second category is used as the parameter partition.

[0027] As a further method, a method for setting a first patch for different working conditions by a patch algorithm based on the parameter partition and the multi-working condition data set includes:

[0028] Based on the unified spatial coordinates and data format of the multi-operating condition data set and parameter partition data, the multi-operating condition data set and parameter partition data are input into the geographic information system. The parameter partition data is used as one layer and the operating condition data is used as another layer in the geographic information system. Multiple layers are superimposed and displayed in sequence to obtain the spatial overlap range;

[0029] The cohesion variation under rainfall intensity conditions is obtained using the saturated-unsaturated seepage model based on the spatial overlap range. If the cohesion variation is greater than or equal to 5%, the corresponding parameter partition is marked as a sensitive area affected by rainfall conditions. The rock and soil elastic modulus variation under earthquake acceleration conditions is obtained using the dynamic response formula based on the spatial overlap range. If the rock and soil elastic modulus variation is greater than or equal to 10%, the corresponding parameter partition is marked as a sensitive area affected by earthquake conditions. The internal friction angle variation under excavation load conditions is obtained using the unloading theory based on the spatial overlap range. If the internal friction angle variation is greater than or equal to 8%, the corresponding parameter partition is marked as a sensitive area affected by excavation conditions.

[0030] A working condition-sensitive zoning mapping table is constructed based on the sensitive areas affected by different working conditions. The degree of working condition impact is obtained based on the change amplitude of cohesion, rock and soil elastic modulus, and internal friction angle. A change amplitude of less than 10% is considered mild, a change amplitude of 10% to 25% is considered moderate, and a change amplitude of more than 25% is considered severe.

[0031] Based on the working condition-sensitivity zoning mapping table and the working condition impact degree classification, the slightly affected area, moderately affected area, and severely affected area are obtained. The cohesion reduction coefficient of the severely affected area by rainfall is determined by the Mohr-Coulomb strength theory. The elastic modulus reduction coefficient of the severely affected area by earthquake is determined by the dynamic response theory. The internal friction angle reduction coefficient of the severely affected area by excavation is determined by the unloading theory. The cohesion reduction coefficient, elastic modulus reduction coefficient, and internal friction angle reduction coefficient are used as adjustment rules.

[0032] Based on the corresponding sensitive zones of the adjustment rules, the pore water pressure value is calculated based on the permeability characteristics and rainfall intensity under rainfall conditions. The magnitude and direction of the inertia force is calculated based on the seismic acceleration under earthquake conditions. The position of the free surface is extracted for excavation conditions. The pore water pressure value, the magnitude and direction of the inertia force, and the position of the free surface are used as boundary conditions to obtain a boundary condition correction list.

[0033] Parameter correction values ​​are generated based on the adjustment rules of the corresponding working conditions. The boundary parameters in the boundary condition correction list that match the spatial position of the parameter partition are embedded in the patch. The parameter correction values ​​and boundary conditions under each working condition are integrated into a data block, which is named the working condition-partition patch. The first patch set of the slope parameter partition is obtained. The first patch and the non-patch area are smoothed according to the first patch set to obtain the first patch model.

[0034] As a further method, the method of smoothing the first patch and the non-patch area according to the first patch set includes:

[0035] At the junction of the first patch and the non-patch area, the parameter gradual transition is performed through the transition area smoothing formula, which is expressed as:

[0036]

[0037] Where P(x,y,t) is the parameter value at time t of coordinate (x,y), I is the set of all working conditions, S is any subset of I, representing a certain working condition combination, P i0 is the initial parameter value of the action area of ​​working condition i, λ i is the time attenuation coefficient of working condition i, which is determined by least square fitting based on the monitoring data, t is the time range of the monitoring data, d i is the shortest distance from the coordinate (x, y) in the 3D model to the patch boundary of working condition i, d mi is the maximum distance threshold affected by the patch of working condition i.

[0038] As a further method, a method for extracting local special geological areas under different working conditions based on the monitoring data set and generating a second patch for the local special geological areas includes:

[0039] Based on the monitoring data set, a sliding window Z-score algorithm is used to detect mutation points. That is, the difference between the data at each time point and the mean in the window is calculated, and then divided by the standard deviation in the window to obtain the Z value. If the absolute value of the Z value is greater than 3, it is marked as an outlier.

[0040] Based on the three-dimensional slope model, outliers with a spatial distance of less than 5 meters were connected into regions to obtain local special geological zones. The overlap was calculated based on the occurrence time of different working conditions in the local special geological zones. Rainfall intensity conditions with an overlap greater than 80% were considered the main cause, and earthquake acceleration conditions or excavation load conditions with an extension of 15 minutes before and after the occurrence time period were considered the main cause.

[0041] Based on the local special geological area and the main inducement, a local finite element sub-model with a grid size of less than 1m was constructed. Bayesian inversion was performed based on the local finite element sub-model. The slope displacement monitoring values, stress monitoring values, and pore water pressure monitoring values ​​were used as observation data. The observed data were used to obtain the probability value of the inversion parameter through the likelihood function. The likelihood function formula is as follows:

[0042]

[0043] Where p(D|θ) is the likelihood function, that is, the probability of the observation data D appearing when the inversion parameter θ is given, T is the total number of time steps of the observation data, Δu t is the difference between the simulated and observed displacements at time step t, Δσ t is the difference between the simulated and observed stress values ​​at time step t, Δp t is the difference between the simulated and observed pore water pressure at time step t, Σ is the covariance matrix of the observation error, and det(Σ) is the determinant of the covariance matrix Σ;

[0044] The posterior distribution is calculated based on the inversion parameter probability value using the Markov chain Monte Carlo sampling algorithm. The mean of the posterior distribution is used as the inversion parameter value. The correction coefficient is obtained by dividing the inversion parameter value by the original parameter partition value. If the correction coefficient of the local special geological area is less than 0.7, the spatial coordinates, inversion parameter values ​​and main inducements of the local special geological area are used for packaging, and the packaged data block is named the second patch.

[0045] As a further method, a method for performing stability evaluation using the parameter partitioning and the patch model by a limit equilibrium method includes:

[0046] Based on parameter partitioning, the Swedish arc search method is used to extract the potential slip surface. The potential slip surface is divided into several soil strips. If the soil strip is in the area covered by the patch model, the geotechnical parameters corrected by the patch model are used. If the soil strip is not in the area covered by the patch model, the geotechnical parameters in the original parameter partition are used.

[0047] The sliding force and anti-sliding force of the soil strip are calculated according to the limit equilibrium method. The formulas for the sliding force and anti-sliding force are:

[0048] T i =W i sinα i

[0049] R i =c i l i +(W i cosα i -u i l i )tanφ i

[0050] The sliding force is T i , the anti-slip force is R i , W i is the gravity of the soil strip itself, sinα i is the sine of the angle between the bottom of the soil strip and the horizontal plane, c i is the cohesion of the soil mechanics parameters in the soil strip, l i is the length of the bottom of the soil strip, cosα i is the cosine value of the angle between the bottom of the soil strip and the horizontal plane, u i is the pore water pressure at the bottom, l i is the length of the bottom of the soil strip, tanφ i is the tangent of the internal friction angle of the soil mechanics parameter in the soil strip;

[0051] Based on the potential sliding surface, the sum of the anti-sliding forces of the soil strips is divided by the sum of the sliding forces to obtain the safety factor of the potential sliding surface. The potential sliding surfaces in each parameter partition of the non-coal mine slope are traversed and the safety factor is calculated. The smallest safety factor in the non-coal mine slope is taken as the first slope stability result.

[0052] As a further method, a method for performing stability evaluation on the relatively dangerous area using a finite element method and outputting a second stability result of the slope includes:

[0053] Based on the first stability results of the slope and national standards, the non-coal mine slopes were divided into unstable and stable areas. The safety factor under natural working conditions was used as the upper limit of the relative danger interval according to the stable area, and the lowest safety factor under multiple working conditions was used as the lower limit of the relative danger interval. The relative danger interval was used to divide the stable area and obtain the relative danger area.

[0054] Extracting geological layer data based on the geotechnical data set of the relatively dangerous area, obtaining the fault zone, special structural rock layer, and homogeneous rock layer in the relatively dangerous area according to the geological layer data, dividing the relatively dangerous area corresponding to the fault zone and special structural rock layer into tetrahedral grid units, and dividing the relatively dangerous area corresponding to the homogeneous rock layer into hexahedral grid units;

[0055] Based on the grid cells, geotechnical parameters are assigned through the finite element model, and boundary conditions are assigned to the grid cells based on the monitoring data set. The boundary conditions include slope displacement, stress, and pore water pressure. According to the finite element model, the cohesion and internal friction angle of the soil are continuously reduced through the finite element strength reduction method until the slope reaches the limit equilibrium state. The strength reduction coefficient of this limit equilibrium state is used as the second stability result of the slope, and the second stability result is used as the final evaluation result.

[0056] The second aspect of the present invention provides a non-coal mine slope stability evaluation system for multi-operating condition analysis, comprising:

[0057] Data preprocessing module: collects non-coal mine slope data in a preset area and preprocesses the non-coal mine slope data; the non-coal mine slope data includes geotechnical data sets, multi-working condition data sets and monitoring data sets;

[0058] A data clustering module is configured to systematically cluster the geotechnical dataset to obtain a plurality of first categories, and based on the first categories, divide the regions of the geotechnical dataset with similar parameters and spatial proximity into a second category by fuzzy kernel clustering, and use the second category as a parameter partition;

[0059] A patch algorithm module is configured to set a first patch for different working conditions using a patch algorithm based on the parameter partitions and the multi-working condition data set, extract local special geological areas under different working conditions according to the monitoring data set, generate a second patch for the local special geological areas, and obtain a patch model by combining the first patch and the second patch;

[0060] A first stability evaluation module is used to perform stability evaluation using the parameter partitioning and the patch model through a limit equilibrium method, and output a first stability result of the slope;

[0061] The second stability evaluation module is used to screen out relatively dangerous areas according to the first stability result of the slope, perform stability evaluation on the relatively dangerous areas by using the finite element method, output the second stability result of the slope, and use the second stability result as the final evaluation result.

[0062] Compared with the prior art, the embodiments of the present invention have at least the following advantages or beneficial effects:

[0063] (1) The present invention integrates geotechnical data, multi-operation data and monitoring data, and combines system clustering and fuzzy kernel clustering technology to achieve refined parameter partitioning of slope geological bodies, thereby improving data utilization efficiency and enhancing the accuracy of dangerous area identification.

[0064] (2) The present invention generates a first patch and a second patch through a patch algorithm to construct a patch model, which accurately depicts the impact of multiple working conditions and special geological conditions on the mechanical properties of the slope, reflects the effects of complex factors, and improves the accuracy and reliability of stability evaluation.

[0065] (3) The present invention preliminarily screens relatively dangerous areas through the limit equilibrium method, and then uses the finite element method to analyze the area in detail, performing rough screening first and then fine evaluation. The limit equilibrium method calculates and locates dangerous areas efficiently and quickly, and the finite element method accurately simulates complex mechanical behaviors. The combination of the two improves the evaluation efficiency and the evaluation accuracy of key areas, providing a scientific and efficient decision-making basis for slope management. BRIEF DESCRIPTION OF THE DRAWINGS

[0066] Figure 1 This is a flowchart of the steps of a non-coal mine slope stability evaluation method based on multi-working condition analysis in an embodiment of the present invention. DETAILED DESCRIPTION

[0067] The technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0068] Reference Figure 1 As shown, the present invention provides a non-coal mine slope stability evaluation method based on multi-working condition analysis, including:

[0069] Collecting non-coal mine slope data in a preset area and preprocessing the non-coal mine slope data; the non-coal mine slope data includes a geotechnical data set, a multi-working condition data set, and a monitoring data set;

[0070] In the actual evaluation, the saturation density of 21.5-24.0 g / cm was collected at 10 points on the open pit slope of an iron ore mine. 3, cohesion 320 ~ 670kPa, internal friction angle 31 ° ~ 37 °, elastic modulus 1.3 ~ 16GPa, Poisson's ratio 0.25 ~ 0.35 and other geotechnical parameters, its dip angle 28 ° ~ 38 °, at the same time, the fault zone, such as the F1 fault cohesion 300kPa, homogeneous rock layer is K1 sandstone, K2 shale, cave is filled with soft clay and other geological layer data, as well as the maximum daily rainfall 80mm / d, earthquake peak acceleration 0.05g, step excavation unloading stress 200kPa, etc. The data were collected from multiple working conditions and monitoring data such as hourly displacement of 2.1 to 5.8 mm, stress of 150 to 200 kPa, and pore water pressure of 100 to 170 kPa at each point. The data were then unified into the National 2000 Geodetic Coordinate System and CSV format, and the rainfall, earthquake, and excavation working conditions were aligned with the monitoring data time at 10-minute intervals. The missing stress values ​​of the points were supplemented by linear interpolation. Finally, a data set was formed that included standardized geotechnical parameters of 10 points, multi-working condition triggering time and parameters, and no missing monitoring data.

[0071] Performing systematic clustering on the geotechnical dataset to obtain a plurality of first categories, dividing regions of the geotechnical dataset with similar parameters and spatial proximity into second categories based on the first categories by fuzzy kernel clustering, and using the second categories as parameter partitions;

[0072] It should be explained that the purpose of using two clustering methods when generating parameter partitions is to take into account the parameter similarity and spatial proximity of the data, so as to more accurately characterize the heterogeneity of the geological body. The specific reason is that systematic clustering realizes preliminary rough classification. Systematic clustering takes geotechnical parameters as the core dimension and clusters data points with higher similarity into one category. This step can perform coarse-grained division of the geotechnical data set based on macroscopic parameter characteristics and preliminary spatial relationships. Based on the first category obtained by systematic clustering, fuzzy kernel clustering further quantifies the parameter similarity between data points through the Gaussian kernel function in the multidimensional feature space (parameter vector + spatial coordinate), and introduces spatial proximity weights to divide the areas with similar parameters and spatial proximity into the second category. This method makes up for the shortcomings of systematic clustering, which can not only retain parameter similarity but also strengthen the influence of spatial position on partitioning, so that the final parameter partitioning is more in line with the actual distribution of the geological body, providing a more refined and accurate parameter basis for subsequent stability evaluation.

[0073] In the actual assessment, a three-dimensional geological model was constructed based on the spatial coordinates and geotechnical parameters of the 10 monitoring points. Density, cohesion, internal friction angle, elastic modulus, and Poisson's ratio were used as clustering dimensions. The Euclidean distance matrix between points was calculated, and spatial proximity weights (e.g., P1 and P2, which are 15 meters apart, were assigned a weight of 0.8) were superimposed on the distance matrix. Hierarchical clustering was used to merge data points. When the incremental sum of squared deviations within a class was minimized and the mean silhouette coefficient was maximized, the number of classes (k = 3) was determined after traversal. The data were classified into sandstone classes, with P1, P2, P9, and P10 for sandstone, P6-P8 for shale, and P3-P5 for transitional layers, including the F1 fault. The three data classes were then imported into a GIS. The dominant joint direction was determined based on the sandstone dip angle of 35° and dip NE 30° in the occurrence data. Kriging interpolation was used to distinguish the ranges in different directions, such as 50 meters for sandstone and 30 meters for shale, and the sill value, to generate a heterogeneous parameter field.

[0074] In the constructed heterogeneous parameter field, the geotechnical parameter vectors are extracted for the data of 10 monitoring points P1-P10. These vectors contain information such as density, cohesion, internal friction angle, elastic modulus, etc. For example, the geotechnical parameter vector of point P1 is [23,670,35,15]. At the same time, the spatial coordinate vector of each point is extracted. This vector consists of longitude, latitude, and elevation. For example, the spatial coordinates of point P1 are [100,200,300]. This forms a multidimensional feature space that contains both parameter features and spatial positions. Then, the Gaussian kernel function is used to process this multidimensional feature space and construct a kernel mapping matrix. Based on this kernel mapping matrix, the initialized fuzzy membership matrix is ​​randomly generated, for example, it is initialized as And the initialized cluster center matrix, where the elements in the cluster center matrix are set according to the typical values ​​of different categories, such as V1 = [23,670,35,15] for sandstone, V2 = [22,400,33,8] for shale, and V3 = [21.5,320,31,1.3] for transition layer. Then, the geotechnical parameter similarity variance h1 = 50 is set to allow parameter differences within a certain range; h2 = 10 is the variance of the spatial coordinate proximity in the kernel space, emphasizing the compactness of the space; w1 = 0.6 is the geotechnical parameter similarity weight coefficient, and w2 = 0.4 is the spatial proximity coefficient. Weight coefficient, and satisfy w1+w2=1; fuzzy factor m=2, construct the objective function based on these parameters, and then use the Lagrange multiplier method to solve the minimization problem of the objective function. In the solution process, according to the derived update formula, the kernel mapping matrix, the current fuzzy membership matrix and the cluster center matrix are used to iteratively update the membership and cluster center. For example, after the first iteration, the cluster center V1 is adjusted to [22.4, 650, 34.1, 13.9]. After many such iterations, the iteration is stopped when the difference between the objective function values ​​of two adjacent iterations is less than 0.1%. For example, from the objective function value J of the nth iteration, the membership and cluster center are iteratively updated. n = 25.6 to J of the n+1th iteration n+1 =25.62, the iteration is stopped when the difference is much less than 0.1%, and the areas with similar parameters and spatial proximity are finally divided into the second type of parameter partitioning. For example, the sandstone type is further subdivided into two sub-areas, and the areas with closer parameters and spatial proximity among P1, P2, P9, and P10 are divided into one sub-area, while the areas near P10, which are slightly farther away in space but still have similar parameters, are adjusted and divided, thus obtaining the final parameter partitioning, which provides a more refined parameter distribution basis for the subsequent slope stability evaluation.

[0075] Based on the parameter partitioning and the multi-operating condition data set, a first patch is set for different operating conditions using a patch algorithm; local special geological areas under different operating conditions are extracted according to the monitoring data set; a second patch is generated for the local special geological areas; and a patch model is obtained by combining the first patch and the second patch;

[0076] In the actual assessment, for rainfall conditions on the slope, such as a rainstorm with an intensity of 30 mm / h, the seepage model calculation found that the cohesion change in the P3-P5 transition layer area exceeded 5%, and it was demarcated as a heavy rainfall affected area. Based on the Mohr-Coulomb theory, the cohesion was reduced from 650 kPa to 640 kPa, and the pore water pressure boundary condition was updated from 150 kPa to 170 kPa. For earthquake conditions on the slope, such as a magnitude 3.2 earthquake, the elastic modulus change was calculated near the F1 fault zone and found to be more than 10%. The elastic modulus was reduced from 6 GPa to 5 GPa, and the increase was 0. Add inertial force boundary conditions, and the inertial force boundary conditions are peak acceleration 0.05g; for the excavation condition of the slope, for example, the excavation depth of the slope step is 15m, the friction angle in the P1-P2 sandstone area changes by more than 8%, reducing from 38° to 35°, and set the free surface boundary conditions. The parameter corrections and boundary conditions of these heavily affected areas are encapsulated as "rainfall-transition layer patch", "earthquake-fault patch", and "excavation-sandstone patch". The first patch set of slope parameter partition is obtained. The first patch and non-patch area are smoothed according to the first patch set to obtain the first patch model.

[0077] For the rainstorm condition i1 of the first patch model, the initial cohesion P in the action area is i10 =650kPa, and the time decay coefficient λ is determined by fitting the monitoring data using the least squares method. i1 =0.05, the maximum distance threshold of patch influence d mi1 =10m; earthquake condition is recorded as i2, elastic modulus P i20 =6GPa,λ i2 =0.03,d mi2 =15m; excavation condition is recorded as i3, internal friction angle P i30 =38°, λ i3 =0.04,d mi3 =8m, taking a coordinate near point P3 as an example, time t is 10 days after the rainstorm, and the distance to the boundary of the rainstorm patch is d i1 =3m, distance to the earthquake patch boundary d i2 =8m, distance d to the excavation patch boundary i3 =12m, because d i3 >d mi3 , not affected by the excavation patch, the contribution of the subset S of different working condition combinations is calculated according to the transition zone smoothing formula. When S = i1, the contribution value is 1.25, when S = i2, the contribution value is 30.7, and when S = {i1, i2}, the contribution value is 389.7. All possible subsets are calculated and summed, and finally the parameter value at the coordinate (x, y) at time t is obtained, realizing a gradual transition of the parameters from the first patch, such as 640 kPa of cohesion in the rainstorm patch to 660 kPa in the non-patch area, avoiding sudden changes in the model boundary and ensuring the accuracy and rationality of the calculation.

[0078] During the actual assessment, monitoring data near point P5 showed a sudden change. For example, at a certain moment, the displacement Z value reached 11. Using a sliding window Z-score algorithm with a window size of 12, this was identified as a localized special geological zone, such as a karst cave. Examining the overlap between the excavation and subsequent rainfall periods revealed a time overlap greater than 80%, determining that the primary cause was excavation unloading combined with rainfall weakening. A local finite element sub-model was constructed, and the Markov Chain Monte Carlo sampling algorithm was used to invert the karst cave parameters. The original cohesion of 490 kPa was inverted to 320 kPa, with a correction factor of 320 / 490 ≈ 0.65 < 0.7. The parameters of this localized area were encapsulated as a second patch, accurately depicting the unique impact of special geological conditions on slope stability.

[0079] Using the parameter partition and the patch model to perform stability evaluation through a limit equilibrium method, and outputting a first stability result of the slope;

[0080] In the actual evaluation, the Swedish arc search method was used to extract a potential sliding surface in the transition layer area of ​​the slope parameter partition and divide it into 10 soil strips. Among them, soil strips 3, 5, and 7 are in the rainstorm patch coverage area. The cohesion is revised from the original 360kPa to 340kPa, and the pore water pressure is updated to 120kPa. The original parameter partition parameters are used for the remaining soil strips to calculate the sliding force and anti-sliding force of each soil strip. Taking soil strip 3 as an example, the gravity of the soil strip itself W3 = 10000kN, the angle between the bottom surface of the soil strip and the horizontal plane α3 = 45°, and the sliding force T3 = 1000 0kN*sin45°=7071kN, the cohesion of the soil strip is c3=400kPa, the length of the bottom surface of the soil strip is l3=10m, the pore water pressure at the bottom is u3=100kPa, the internal friction angle of the soil strip is φ3=33, then the anti-slip force R3=400×10+(7071-100×10)×0.6494=7942.5kN. For all soil strips, for example, the total sliding force on a potential slip surface is 5000kN, the total anti-slip force is 5800kN, and the safety factor is Continuing to traverse other potential sliding surfaces, we finally took the minimum safety factor of 1.16 as the first stability result of the slope. The standard stability coefficient of the slope under natural working conditions should be greater than or equal to 1.25, and it can be reduced to 1.23-1.18 under heavy rain conditions, indicating that this slope is unstable under heavy rain conditions, providing a preliminary basis for subsequent detailed evaluation.

[0081] A relatively dangerous area is screened out according to the first stability result of the slope, stability evaluation is performed on the relatively dangerous area using the finite element method, a second stability result of the slope is output, and the second stability result is used as the final evaluation result.

[0082] In the actual assessment, based on the first stability result of 1.197 of the iron ore slope, the slope was carefully divided according to the "Safety Regulations for Metal and Non-metal Mines" and other national standards, and unstable and stable areas were defined. In the stable area, the safety factor of 1.25 obtained through long-term monitoring and precise calculation under natural working conditions was used as the upper limit of the relative danger interval. The lowest safety factor of 1.15 obtained under multiple working conditions such as continuous heavy rainfall causing a significant increase in pore water pressure and frequent blasting vibration causing stress fluctuations was used as the lower limit of the relative danger interval. With the help of this relative danger interval, the stable area was divided, thereby accurately delineating the relative danger area. The geotechnical dataset of the domain was extracted to clearly identify geological horizon data within the region, such as the F1 fault, which is approximately 3-5 meters wide and filled with crushed rock debris with a cohesion as low as 380 kPa; special structural rock formations, such as sandstone layers with dense joints and an average joint spacing of less than 0.5 meters, which have poor rock integrity; and limestone layers with good integrity and uniform lithology. Based on these differences in geological characteristics, the relatively dangerous areas corresponding to the fault zone and special structural rock formations were divided into tetrahedral grid cells, and the relatively dangerous areas corresponding to the homogeneous rock formations were divided into hexahedral grid cells. Then, based on the divided grid cells, geotechnical parameters were assigned using a finite element model. For example, the cohesion of the fault zone is set to 320kPa and the internal friction angle is set to 25°, the cohesion of the homogeneous rock layer is set to 500kPa and the internal friction angle is set to 38°. According to the long-term on-site monitoring data set, the boundary conditions are assigned to the grid units, 20 monitoring points are arranged on the slope, and the monitoring data are integrated with the patch model to extract characteristic parameters such as displacement rate, stress change rate, and pore water pressure change rate. The finite element strength reduction method is used to continuously reduce the cohesion and internal friction angle of the soil, and simulate the mechanical response of the slope in the process of gradual weakening until the slope reaches the ultimate equilibrium state. For example, the displacement rate of a certain monitoring point is A rate of over 5 mm / d for three consecutive days was considered a signal of possible stability changes. During a monitoring session, it was discovered that the displacement rate at point P3 had increased suddenly. Analysis revealed that the pore water pressure had increased due to rainfall, so the pore water pressure boundary conditions in the area near P3 were updated to 150 kPa. The geotechnical parameters of the area were recalculated, and the updated patch model was applied to the finite element model. The strength reduction factor was recalculated, and the result dropped from 1.12 to 1.05. The strength reduction factor in this limit equilibrium state was used as the second stability result of the slope, indicating that the slope stability had further decreased. Based on this, reinforcement measures were initiated in a timely manner to ensure the safety of the slope.

[0083] In this embodiment, the preprocessing method includes:

[0084] The density, cohesion, internal friction angle, elastic modulus and Poisson's ratio at different locations of non-coal mine slopes under saturated conditions are collected as geotechnical parameters; the fault zone distribution, homogeneous rock layers, special structural rock layers, occurrence and lithological characteristics of non-coal mine slopes are collected as geological layer data; the meteorological historical records, seismic network data, groundwater level information, mine excavation design and construction records of non-coal mine slope areas are collected; rainfall intensity, earthquake acceleration, groundwater depth and excavation load are obtained as multi-condition data; the slope displacement monitoring values, stress monitoring values ​​and pore water pressure monitoring values ​​are obtained as monitoring data; the geotechnical parameters, geological layer data, multi-condition data and monitoring data are unified in format, time axis aligned, outliers eliminated and missing values ​​interpolated to obtain geotechnical data sets, multi-condition data sets and monitoring data sets.

[0085] In this embodiment, a method for performing systematic clustering on the geotechnical dataset to obtain a plurality of first categories includes:

[0086] A three-dimensional model is constructed based on the spatial location of the data points in the geotechnical dataset, geological layer data, and geotechnical parameters. The geotechnical parameters of the data points are used as clustering dimensions according to the three-dimensional model, and the Euclidean distance matrix is ​​calculated for the data points using the clustering dimensions.

[0087] Based on the sampling spatial position of the data points, the data points are connected into a triangular grid through triangulation, and the distance between the data points in the triangular grid is used as the spatial proximity weight. The spatial proximity weight is multiplied by the Euclidean distance of the Euclidean distance matrix to obtain the Euclidean distance matrix of the superimposed spatial proximity;

[0088] Hierarchical clustering is performed based on the Euclidean distance matrix of superimposed spatial proximity. The data points are taken as the initial class and the data points are gradually merged until the increment of the sum of squared deviations within the class after merging is minimized, and a clustering hierarchical tree is obtained.

[0089] The silhouette coefficient of the data points is calculated based on the clustering hierarchy tree, and the mean of the silhouette coefficient corresponding to different number of categories is obtained. The number of categories is traversed and the number of categories with the largest mean of the silhouette coefficient is selected as the final number of categories. The geotechnical data set is divided into multiple first categories according to the final number of categories.

[0090] The first type of data is imported into the geographic information system, and the dominant direction of the joint is determined according to the occurrence data corresponding to the geotechnical parameters. Kriging interpolation is performed on the geotechnical parameters based on the dominant direction of the joint. During the interpolation, the variation function is used to distinguish the range and base value of different joint dominant directions to obtain a heterogeneous parameter field.

[0091] In this embodiment, the method of dividing the regions of the geotechnical dataset with similar data point parameters and spatial proximity into the second category by fuzzy kernel clustering based on the first category includes:

[0092] Based on the first type of geotechnical parameter vectors and spatial coordinate vectors of data points extracted in the heterogeneous parameter field, a multidimensional feature space containing parameter characteristics and spatial positions is obtained;

[0093] According to the multidimensional feature space, a kernel mapping matrix is ​​constructed by using the Gaussian kernel function. The initial fuzzy membership matrix and the initial cluster center matrix are randomly generated based on the kernel mapping matrix. The objective function is constructed according to the kernel mapping matrix, the initial fuzzy membership matrix and the initial cluster center matrix, which can be expressed as:

[0094]

[0095] Where J(U,V) is the objective function, U is the fuzzy membership matrix, V is the cluster center matrix, c is the number of cluster centers, n is the total number of data points, u ik The fuzzy membership of data point i to cluster center k, m is the fuzzy factor, p i Geomechanical parameter vector of data point i, l i The spatial coordinate vector of data point i, h1 is the variance, which is the scale of kernel space to measure the similarity of geotechnical parameters. The smaller the value, the faster the parameter similarity decays as the parameter difference increases. Clustering is more sensitive to parameter differences. The larger the value, the data points with large parameter differences can still have high similarity. h2 is the variance, which is the scale of kernel space to measure the proximity of spatial coordinates. The smaller the value, the faster the spatial proximity decreases as the distance increases. Clustering emphasizes spatial compactness. The larger the value, the data points with farther distances can be classified into one category due to parameter similarity. w1 is the similarity weight coefficient of geotechnical parameters, w2 is the spatial proximity weight coefficient, w1+w2=1, d ik is the spatial Euclidean distance between data point i and cluster center k;

[0096] The Lagrange multiplier method is used to solve the minimization problem of the objective function. According to the update formula derived from the Lagrange multiplier method, the kernel mapping matrix, the current fuzzy membership matrix and the cluster center matrix are used to iterate the membership and cluster center through the update formula until the difference between the objective function values ​​of two adjacent iterations is less than 0.1%. The second category of fuzzy kernel clustering is obtained and the second category is used as the parameter partition.

[0097] In this embodiment, the method of setting the first patch for different working conditions by using a patch algorithm based on the parameter partition and the multi-working condition data set includes:

[0098] Based on the unified spatial coordinates and data format of the multi-operating condition data set and parameter partition data, the multi-operating condition data set and parameter partition data are input into the geographic information system. The parameter partition data is used as one layer and the operating condition data is used as another layer in the geographic information system. Multiple layers are superimposed and displayed in sequence to obtain the spatial overlap range;

[0099] The cohesion variation under rainfall intensity conditions is obtained using the saturated-unsaturated seepage model based on the spatial overlap range. If the cohesion variation is greater than or equal to 5%, the corresponding parameter partition is marked as a sensitive area affected by rainfall conditions. The rock and soil elastic modulus variation under earthquake acceleration conditions is obtained using the dynamic response formula based on the spatial overlap range. If the rock and soil elastic modulus variation is greater than or equal to 10%, the corresponding parameter partition is marked as a sensitive area affected by earthquake conditions. The internal friction angle variation under excavation load conditions is obtained using the unloading theory based on the spatial overlap range. If the internal friction angle variation is greater than or equal to 8%, the corresponding parameter partition is marked as a sensitive area affected by excavation conditions.

[0100] A working condition-sensitive zoning mapping table is constructed based on the sensitive areas affected by different working conditions. The degree of working condition impact is obtained based on the change amplitude of cohesion, rock and soil elastic modulus, and internal friction angle. A change amplitude of less than 10% is considered mild, a change amplitude of 10% to 25% is considered moderate, and a change amplitude of more than 25% is considered severe.

[0101] Based on the working condition-sensitivity zoning mapping table and the working condition impact degree classification, the slightly affected area, moderately affected area, and severely affected area are obtained. The cohesion reduction coefficient of the severely affected area by rainfall is determined by the Mohr-Coulomb strength theory. The elastic modulus reduction coefficient of the severely affected area by earthquake is determined by the dynamic response theory. The internal friction angle reduction coefficient of the severely affected area by excavation is determined by the unloading theory. The cohesion reduction coefficient, elastic modulus reduction coefficient, and internal friction angle reduction coefficient are used as adjustment rules.

[0102] Based on the corresponding sensitive zones of the adjustment rules, the pore water pressure value is calculated based on the permeability characteristics and rainfall intensity under rainfall conditions. The magnitude and direction of the inertia force is calculated based on the seismic acceleration under earthquake conditions. The position of the free surface is extracted for excavation conditions. The pore water pressure value, the magnitude and direction of the inertia force, and the position of the free surface are used as boundary conditions to obtain a boundary condition correction list.

[0103] Parameter correction values ​​are generated based on the adjustment rules of the corresponding working conditions. The boundary parameters in the boundary condition correction list that match the spatial position of the parameter partition are embedded in the patch. The parameter correction values ​​and boundary conditions under each working condition are integrated into a data block, which is named the working condition-partition patch. The first patch set of the slope parameter partition is obtained. The first patch and the non-patch area are smoothed according to the first patch set to obtain the first patch model.

[0104] In this embodiment, the method for smoothing the first patch and the non-patch area according to the first patch set includes:

[0105] At the junction of the first patch and the non-patch area, the parameter gradual transition is performed through the transition area smoothing formula, which is expressed as:

[0106]

[0107] Where P(x,y,t) is the parameter value at time t of coordinate (x,y), I is the set of all working conditions, S is any subset of I, representing a certain working condition combination, P i0 is the initial parameter value of the action area of ​​working condition i, λ i is the time attenuation coefficient of working condition i, which is determined by least square fitting based on the monitoring data, t is the time range of the monitoring data, d i is the shortest distance from the coordinate (x, y) in the 3D model to the patch boundary of working condition i, d mi is the maximum distance threshold affected by the patch of working condition i.

[0108] In this embodiment, a method for extracting local special geological areas under different working conditions based on the monitoring data set and generating a second patch for the local special geological areas includes:

[0109] Based on the monitoring data set, a sliding window Z-score algorithm is used to detect mutation points. That is, the difference between the data at each time point and the mean in the window is calculated, and then divided by the standard deviation in the window to obtain the Z value. If the absolute value of the Z value is greater than 3, it is marked as an outlier.

[0110] Based on the three-dimensional slope model, outliers with a spatial distance of less than 5 meters were connected into regions to obtain local special geological zones. The overlap was calculated based on the occurrence time of different working conditions in the local special geological zones. Rainfall intensity conditions with an overlap greater than 80% were considered the main cause, and earthquake acceleration conditions or excavation load conditions with an extension of 15 minutes before and after the occurrence time period were considered the main cause.

[0111] Based on the local special geological area and the main inducement, a local finite element sub-model with a grid size of less than 1m was constructed. Bayesian inversion was performed based on the local finite element sub-model. The slope displacement monitoring values, stress monitoring values, and pore water pressure monitoring values ​​were used as observation data. The observed data were used to obtain the probability value of the inversion parameter through the likelihood function. The likelihood function formula is as follows:

[0112]

[0113] Where p(D|θ) is the likelihood function, that is, the probability of the observation data D appearing when the inversion parameter θ is given, T is the total number of time steps of the observation data, Δu t is the difference between the simulated and observed displacements at time step t, Δσ t is the difference between the simulated and observed stress values ​​at time step t, Δp t is the difference between the simulated and observed pore water pressure at time step t, Σ is the covariance matrix of the observation error, and det(Σ) is the determinant of the covariance matrix Σ;

[0114] The posterior distribution is calculated based on the inversion parameter probability value using the Markov chain Monte Carlo sampling algorithm. The mean of the posterior distribution is used as the inversion parameter value. The correction coefficient is obtained by dividing the inversion parameter value by the original parameter partition value. If the correction coefficient of the local special geological area is less than 0.7, the spatial coordinates, inversion parameter values ​​and main inducements of the local special geological area are used for packaging, and the packaged data block is named the second patch.

[0115] In this embodiment, the method for evaluating stability using the parameter partitioning and the patch model through the limit equilibrium method includes:

[0116] Based on parameter partitioning, the Swedish arc search method is used to extract the potential slip surface. The potential slip surface is divided into several soil strips. If the soil strip is in the area covered by the patch model, the geotechnical parameters corrected by the patch model are used. If the soil strip is not in the area covered by the patch model, the geotechnical parameters in the original parameter partition are used.

[0117] The sliding force and anti-sliding force of the soil strip are calculated according to the limit equilibrium method. The formulas for the sliding force and anti-sliding force are:

[0118] T i =W i sinα i

[0119] R i =c i l i +(W i cosα i -u i l i )tanφ i

[0120] The sliding force is T i , the anti-slip force is R i , W i is the gravity of the soil strip itself, sinα i is the sine of the angle between the bottom of the soil strip and the horizontal plane, c i is the cohesion of the soil mechanics parameters in the soil strip, l i is the length of the bottom of the soil strip, cosα i is the cosine value of the angle between the bottom of the soil strip and the horizontal plane, u i is the pore water pressure at the bottom, l i is the length of the bottom of the soil strip, tanφ i is the tangent of the internal friction angle of the soil mechanics parameter in the soil strip;

[0121] Based on the potential sliding surface, the sum of the anti-sliding forces of the soil strips is divided by the sum of the sliding forces to obtain the safety factor of the potential sliding surface. The potential sliding surfaces in each parameter partition of the non-coal mine slope are traversed and the safety factor is calculated. The smallest safety factor in the non-coal mine slope is taken as the first slope stability result.

[0122] In this embodiment, the method of performing stability evaluation on the relatively dangerous area using the finite element method and outputting a second slope stability result includes:

[0123] Based on the first stability results of the slope and national standards, the non-coal mine slopes were divided into unstable and stable areas. The safety factor under natural working conditions was used as the upper limit of the relative danger interval according to the stable area, and the lowest safety factor under multiple working conditions was used as the lower limit of the relative danger interval. The relative danger interval was used to divide the stable area and obtain the relative danger area.

[0124] Extracting geological layer data based on the geotechnical data set of the relatively dangerous area, obtaining the fault zone, special structural rock layer, and homogeneous rock layer in the relatively dangerous area according to the geological layer data, dividing the relatively dangerous area corresponding to the fault zone and special structural rock layer into tetrahedral grid units, and dividing the relatively dangerous area corresponding to the homogeneous rock layer into hexahedral grid units;

[0125] Based on the grid cells, geotechnical parameters are assigned through the finite element model, and boundary conditions are assigned to the grid cells based on the monitoring data set. The boundary conditions include slope displacement, stress, and pore water pressure. According to the finite element model, the cohesion and internal friction angle of the soil are continuously reduced through the finite element strength reduction method until the slope reaches the limit equilibrium state. The strength reduction coefficient of this limit equilibrium state is used as the second stability result of the slope, and the second stability result is used as the final evaluation result.

[0126] The second aspect of the present invention further provides a non-coal mine slope stability evaluation system for multi-operating condition analysis, comprising:

[0127] Data preprocessing module: collects non-coal mine slope data in a preset area and preprocesses the non-coal mine slope data; the non-coal mine slope data includes geotechnical data sets, multi-working condition data sets and monitoring data sets;

[0128] A data clustering module is configured to systematically cluster the geotechnical dataset to obtain a plurality of first categories, and based on the first categories, divide the regions of the geotechnical dataset with similar parameters and spatial proximity into a second category by fuzzy kernel clustering, and use the second category as a parameter partition;

[0129] A patch algorithm module is configured to set a first patch for different working conditions using a patch algorithm based on the parameter partitions and the multi-working condition data set, extract local special geological areas under different working conditions according to the monitoring data set, generate a second patch for the local special geological areas, and obtain a patch model by combining the first patch and the second patch;

[0130] A first stability evaluation module is used to perform stability evaluation using the parameter partitioning and the patch model through a limit equilibrium method, and output a first stability result of the slope;

[0131] The second stability evaluation module is used to screen out relatively dangerous areas according to the first stability result of the slope, perform stability evaluation on the relatively dangerous areas by using the finite element method, output the second stability result of the slope, and use the second stability result as the final evaluation result.

[0132] The above content is merely an example and explanation of the structure of the present invention. Those skilled in the art may make various modifications or additions to the described specific embodiments or replace them in a similar manner. As long as they do not deviate from the structure of the invention or exceed the scope defined by the claims, they should all fall within the scope of protection of the present invention.

Claims

1. A non-coal mine slope stability evaluation method based on multi-working condition analysis is characterized by: The following steps are involved: Collecting non-coal mine slope data in a preset area and preprocessing the non-coal mine slope data; The non-coal mine slope data includes geotechnical data sets, multi-working condition data sets and monitoring data sets; Performing systematic clustering on the geotechnical dataset to obtain a plurality of first categories, dividing regions of the geotechnical dataset with similar parameters and spatial proximity into second categories based on the first categories by fuzzy kernel clustering, and using the second categories as parameter partitions; Based on the parameter partitioning and the multi-operating condition data set, a first patch is set for different operating conditions using a patch algorithm; local special geological areas under different operating conditions are extracted according to the monitoring data set; a second patch is generated for the local special geological areas; and a patch model is obtained by combining the first patch and the second patch; Using the parameter partition and the patch model to perform stability evaluation through a limit equilibrium method, and outputting a first stability result of the slope; Screening out relatively dangerous areas based on the first stability result of the slope, performing stability evaluation on the relatively dangerous areas using the finite element method, outputting a second stability result of the slope, and using the second stability result as the final evaluation result; The method for setting a first patch for different working conditions by a patch algorithm based on the parameter partition and the multi-working condition data set includes: Based on the unified spatial coordinates and data format of the multi-operating condition data set and parameter partition data, the multi-operating condition data set and parameter partition data are input into the geographic information system. The parameter partition data is used as one layer and the operating condition data is used as another layer in the geographic information system. Multiple layers are superimposed and displayed in sequence to obtain the spatial overlap range; According to the spatial overlap range, the cohesion variation under rainfall intensity conditions is obtained through the saturated-unsaturated seepage model. If the cohesion variation is greater than or equal to 5%, the corresponding parameter partition is marked as a sensitive area affected by rainfall conditions. According to the spatial overlap range, the rock and soil elastic modulus variation under earthquake acceleration conditions is obtained through the dynamic response formula. If the rock and soil elastic modulus variation is greater than or equal to 10%, the corresponding parameter partition is marked as a sensitive area affected by earthquake conditions. According to the spatial overlap range, the internal friction angle variation under excavation load conditions is obtained through the unloading theory. If the internal friction angle variation is greater than or equal to 8%, the corresponding parameter partition is marked as a sensitive area affected by excavation conditions. A working condition-sensitive zoning mapping table is constructed based on the sensitive areas affected by different working conditions. The degree of working condition impact is obtained based on the change amplitude of cohesion, rock and soil elastic modulus, and internal friction angle. A change amplitude of less than 10% is considered mild, a change amplitude of 10% to 25% is considered moderate, and a change amplitude of more than 25% is considered severe. Based on the working condition-sensitivity zoning mapping table and the working condition impact degree classification, the slightly affected area, moderately affected area, and severely affected area are obtained. The cohesion reduction coefficient of the severely affected area by rainfall is determined by the Mohr-Coulomb strength theory. The elastic modulus reduction coefficient of the severely affected area by earthquake is determined by the dynamic response theory. The internal friction angle reduction coefficient of the severely affected area by excavation is determined by the unloading theory. The cohesion reduction coefficient, elastic modulus reduction coefficient, and internal friction angle reduction coefficient are used as adjustment rules. Based on the corresponding sensitive zones of the adjustment rules, the pore water pressure value is calculated based on the permeability characteristics and rainfall intensity under rainfall conditions. The magnitude and direction of the inertia force is calculated based on the seismic acceleration under earthquake conditions. The position of the free surface is extracted for excavation conditions. The pore water pressure value, the magnitude and direction of the inertia force, and the position of the free surface are used as boundary conditions to obtain a boundary condition correction list. Generate parameter correction values ​​based on the adjustment rules of the corresponding working condition, embed boundary parameters in the boundary condition correction list that match the parameter partition spatial position into patches, integrate the parameter correction values ​​and boundary conditions under each working condition into a data block, name the data block working condition-partition patch, obtain a first patch set of slope parameter partitions, smooth the first patch and non-patch area according to the first patch set, and obtain a first patch model; A method for extracting local special geological areas under different working conditions according to the monitoring data set and generating a second patch for the local special geological areas includes: Based on the monitoring data set, a sliding window Z-score algorithm is used to detect mutation points. That is, the difference between the data at each time point and the mean in the window is calculated, and then divided by the standard deviation in the window to obtain the Z value. If the absolute value of the Z value is greater than 3, it is marked as an outlier. Based on the three-dimensional slope model, outliers with a spatial distance of less than 5 meters were connected into regions to obtain local special geological zones. The overlap was calculated based on the occurrence time of different working conditions in these local special geological zones. Rainfall intensity conditions with an overlap greater than 80% were considered the main cause, and earthquake acceleration conditions or excavation load conditions with an extension of 15 minutes before and after the occurrence time period were considered the main cause. Based on the local special geological area and the main inducement, a local finite element sub-model with a grid size of less than 1m was constructed. Bayesian inversion was performed based on the local finite element sub-model. The slope displacement monitoring values, stress monitoring values, and pore water pressure monitoring values ​​were used as observation data. The observed data were used to obtain the probability value of the inversion parameter through the likelihood function. The likelihood function formula is as follows: in is the likelihood function, that is, given the parameters to be inverted When the observation data The probability of occurrence, is the total number of time steps of observation data, For the time step The difference between the simulated and observed displacements is For the time step The difference between the simulated and observed stress values ​​at time , is in the time step The difference between the simulated and observed pore water pressure at time ∑ is the covariance matrix of the observation error, is the determinant of the covariance matrix Σ; The posterior distribution is calculated based on the inversion parameter probability value using the Markov chain Monte Carlo sampling algorithm. The mean of the posterior distribution is used as the inversion parameter value. The correction coefficient is obtained by dividing the inversion parameter value by the original parameter partition value. If the correction coefficient of the local special geological area is less than 0.7, the spatial coordinates, inversion parameter values ​​and main inducements of the local special geological area are used for packaging, and the packaged data block is named the second patch.

2. The non-coal mine slope stability evaluation method based on multi-operating condition analysis according to claim 1 is characterized in that: The pretreatment method comprises: The density, cohesion, internal friction angle, elastic modulus and Poisson's ratio of different positions on non-coal mine slopes under saturated conditions are collected as geotechnical parameters; the fault zone distribution, homogeneous rock formations, special structural rock formations and occurrence of non-coal mine slopes are collected as geological layer data; the meteorological historical records, seismic network data, groundwater level information, mine excavation design and construction records of non-coal mine slope areas are collected; rainfall intensity, earthquake acceleration, groundwater depth and excavation load are obtained as multi-condition data; the slope displacement monitoring values, stress monitoring values ​​and pore water pressure monitoring values ​​are obtained as monitoring data; the geotechnical parameters, geological layer data, multi-condition data and monitoring data are formatted in a unified manner, time axis aligned, outliers eliminated and missing values ​​interpolated to obtain geotechnical data sets, multi-condition data sets and monitoring data sets.

3. The non-coal mine slope stability evaluation method based on multi-operating condition analysis according to claim 1 is characterized in that: The method of performing systematic clustering on the geotechnical data set to obtain a plurality of first categories comprises: A three-dimensional model is constructed based on the spatial location of the data points in the geotechnical dataset, geological layer data, and geotechnical parameters. The geotechnical parameters of the data points are used as clustering dimensions according to the three-dimensional model, and the Euclidean distance matrix is ​​calculated for the data points using the clustering dimensions. Based on the sampling spatial position of the data points, the data points are connected into a triangular grid through triangulation, and the distance between the data points in the triangular grid is used as the spatial proximity weight. The spatial proximity weight is multiplied by the Euclidean distance of the Euclidean distance matrix to obtain the Euclidean distance matrix of the superimposed spatial proximity; Hierarchical clustering is performed based on the Euclidean distance matrix of superimposed spatial proximity. The data points are taken as the initial class and the data points are gradually merged until the increment of the sum of squared deviations within the class after merging is minimized, and a clustering hierarchical tree is obtained. The silhouette coefficient of the data points is calculated based on the clustering hierarchy tree, and the mean of the silhouette coefficient corresponding to different number of categories is obtained. The number of categories is traversed and the number of categories with the largest mean of the silhouette coefficient is selected as the final number of categories. The geotechnical data set is divided into multiple first categories according to the final number of categories. The first type of data is imported into the geographic information system, and the dominant direction of the joint is determined according to the occurrence data corresponding to the geotechnical parameters. Kriging interpolation is performed on the geotechnical parameters based on the dominant direction of the joint. During the interpolation, the variation function is used to distinguish the range and base value of different joint dominant directions to obtain a heterogeneous parameter field.

4. The non-coal mine slope stability evaluation method based on multi-operating condition analysis according to claim 1 is characterized in that: The method of dividing the areas of the geotechnical data set with similar data point parameters and spatial proximity into a second category by fuzzy kernel clustering based on the first category includes: Based on the first type of geotechnical parameter vectors and spatial coordinate vectors of data points extracted in the heterogeneous parameter field, a multidimensional feature space containing parameter characteristics and spatial positions is obtained; According to the multidimensional feature space, a kernel mapping matrix is ​​constructed by using the Gaussian kernel function. The initial fuzzy membership matrix and the initial cluster center matrix are randomly generated based on the kernel mapping matrix. The objective function is constructed according to the kernel mapping matrix, the initial fuzzy membership matrix and the initial cluster center matrix, which can be expressed as: in is the objective function, is the fuzzy membership matrix, is the cluster center matrix, is the number of cluster centers, is the total number of data points, Data Points Cluster centers The fuzzy membership of is the fuzzy factor, Data Points The geomechanical parameter vector of Data Points The space coordinate vector of is the variance, which is the scale of the kernel space to measure the similarity of geotechnical parameters. The smaller the value, the faster the parameter similarity decays as the parameter difference increases. Clustering is more sensitive to parameter differences. The larger the value, the data points with large parameter differences can still have high similarity. Variance is the kernel space's measure of spatial coordinate proximity. The smaller the value, the faster the spatial proximity decreases with increasing distance. Clustering emphasizes spatial compactness. The larger the value, the more distant data points are allowed to be grouped together due to parameter similarity. is the similarity weight coefficient of geotechnical parameters, is the spatial proximity weight coefficient, + , For data points and cluster centers The spatial Euclidean distance of The Lagrange multiplier method is used to solve the minimization problem of the objective function. According to the update formula derived from the Lagrange multiplier method, the kernel mapping matrix, the current fuzzy membership matrix and the cluster center matrix are used to iterate the membership and cluster center through the update formula until the difference between the objective function values ​​of two adjacent iterations is less than 0.1%. The second category of fuzzy kernel clustering is obtained and the second category is used as the parameter partition.

5. The non-coal mine slope stability evaluation method based on multi-operating condition analysis according to claim 1 is characterized in that: The method for smoothing the first patch and the non-patch area according to the first patch set includes: At the junction of the first patch and the non-patch area, the parameter gradual transition is performed through the transition area smoothing formula, which is expressed as: in For coordinates time The parameter value at is the set of all working conditions, for Any subset of represents a certain working condition combination, For working conditions The initial parameter value of the action area, For working conditions The time decay coefficient of is determined by least squares fitting based on the monitoring data. is the time range of the monitoring data, Coordinates in the 3D model To working conditions The shortest distance to the patch boundary, For working conditions The maximum distance threshold for a patch to be affected.

6. The non-coal mine slope stability evaluation method based on multi-operating condition analysis according to claim 1 is characterized in that: The method for performing stability evaluation by using the parameter partitioning and the patch model through the limit equilibrium method comprises: Based on parameter partitioning, the Swedish arc search method is used to extract the potential slip surface. The potential slip surface is divided into several soil strips. If the soil strip is in the area covered by the patch model, the geotechnical parameters corrected by the patch model are used. If the soil strip is not in the area covered by the patch model, the geotechnical parameters in the original parameter partition are used. The sliding force and anti-sliding force of the soil strip are calculated according to the limit equilibrium method. The formulas for the sliding force and anti-sliding force are: The sliding force is , the anti-slip force is , is the gravity of the soil strip itself, is the sine of the angle between the bottom of the soil strip and the horizontal plane, is the cohesion of the soil strip, is the length of the bottom surface of the soil strip, is the cosine of the angle between the bottom of the soil strip and the horizontal plane, is the pore water pressure at the bottom, is the length of the bottom surface of the soil strip, is the tangent of the internal friction angle of the soil mechanics parameter in the soil strip; Based on the potential sliding surface, the sum of the anti-sliding forces of the soil strips is divided by the sum of the sliding forces to obtain the safety factor of the potential sliding surface. The potential sliding surfaces in each parameter partition of the non-coal mine slope are traversed and the safety factor is calculated. The smallest safety factor in the non-coal mine slope is taken as the first slope stability result.

7. The non-coal mine slope stability evaluation method based on multi-operating condition analysis according to claim 1 is characterized in that: The method of performing stability evaluation on the relatively dangerous area by using the finite element method and outputting a second stability result of the slope includes: Based on the first stability results of the slope and national standards, the non-coal mine slopes were divided into unstable and stable areas. The safety factor under natural working conditions was used as the upper limit of the relative danger interval according to the stable area, and the lowest safety factor under multiple working conditions was used as the lower limit of the relative danger interval. The relative danger interval was used to divide the stable area and obtain the relative danger area. Extracting geological layer data based on the geotechnical data set of the relatively dangerous area, obtaining the fault zone, special structural rock layer, and homogeneous rock layer in the relatively dangerous area according to the geological layer data, dividing the relatively dangerous area corresponding to the fault zone and special structural rock layer into tetrahedral grid units, and dividing the relatively dangerous area corresponding to the homogeneous rock layer into hexahedral grid units; Based on the grid cells, geotechnical parameters are assigned through the finite element model, and boundary conditions are assigned to the grid cells based on the monitoring data set. The boundary conditions include slope displacement, stress, and pore water pressure. According to the finite element model, the cohesion and internal friction angle of the soil are continuously reduced through the finite element strength reduction method until the slope reaches the limit equilibrium state. The strength reduction coefficient of this limit equilibrium state is used as the second stability result of the slope, and the second stability result is used as the final evaluation result.

8. A non-coal mine slope stability evaluation system based on multi-operating condition analysis, for executing the non-coal mine slope stability evaluation method based on multi-operating condition analysis according to any one of claims 1 to 7, characterized in that: The system comprises: Data preprocessing module: collects non-coal mine slope data in a preset area and preprocesses the non-coal mine slope data; the non-coal mine slope data includes geotechnical data sets, multi-working condition data sets and monitoring data sets; A data clustering module is configured to systematically cluster the geotechnical dataset to obtain a plurality of first categories, and based on the first categories, divide the regions of the geotechnical dataset with similar parameters and spatial proximity into a second category by fuzzy kernel clustering, and use the second category as a parameter partition; A patch algorithm module is configured to set a first patch for different working conditions using a patch algorithm based on the parameter partitions and the multi-working condition data set, extract local special geological areas under different working conditions according to the monitoring data set, generate a second patch for the local special geological areas, and obtain a patch model by combining the first patch and the second patch; A first stability evaluation module is used to perform stability evaluation using the parameter partitioning and the patch model through a limit equilibrium method, and output a first stability result of the slope; The second stability evaluation module is used to screen out relatively dangerous areas according to the first stability result of the slope, perform stability evaluation on the relatively dangerous areas by using the finite element method, output the second stability result of the slope, and use the second stability result as the final evaluation result.

Images