Simulated fault stability test system
By using a simulated fault stability testing system, data collected from multiple types of sensors is used to identify the differences in physical properties of faults in coal-bearing reservoirs, calculate the stress migration coupling coefficient, simulate the stress distribution trajectory, determine the slip trend, and optimize the loading strategy. This solves the problem of large deviations in the fault slip limit value in existing technologies and enables accurate identification and control of fault slip risk.
Patent Information
- Application Number
- CN202511090915.0
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-08-05
- Publication Date
- 2025-11-25
AI Technical Summary
Existing simulation systems cannot accurately capture the asymmetric stress evolution process when faced with the heterogeneity of fault structures in coal-bearing reservoirs, resulting in large deviations in fault slip threshold values and making it difficult to accurately determine the risk of fault activation in actual engineering projects.
A simulated fault stability testing system is adopted, including a fault data acquisition module, a heterostructure identification module, a coupling factor calculation module, a stress evolution simulation module, and a critical evolution judgment module. Data is collected through multiple types of sensors to identify the differences in physical properties between the hanging wall and footwall of the fault, calculate the stress migration coupling coefficient, simulate the stress distribution trajectory, determine the slip trend, and optimize the loading strategy.
It enables accurate judgment and loading control of fault slip trends, enhances the forward-looking identification capability of fault slip risks, and significantly improves the accuracy of characterizing the stress evolution influence mechanism of fault systems in heterogeneous structures.
Smart Images

Figure CN121009689A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of fault testing technology, specifically to a simulated fault stability testing system. Background Technology
[0002] Geotechnical engineering and geological hazard research, as an important branch of geological engineering technology, has long focused on the stability of underground structures. Within this macro-technical framework, fault mechanics, due to its direct correlation with the mechanical continuity, fracturing characteristics, and slippage trends of the underground medium, is considered an important research object that cannot be ignored. Further focusing on the field of energy engineering, particularly in the green and safe development of coal resources, the fault structure within coal-bearing reservoirs, due to its dual role of "conduction / closure" in the occurrence, migration, and extraction of coalbed methane, has become one of the core bottlenecks in the efficiency and safety control of coalbed methane extraction.
[0003] Most current testing systems used for coal seam fault simulation make a simplified assumption in their model construction or algorithm settings: that the hanging wall and footwall of the fault are homogeneous structures with consistent mechanical properties. This assumption is barely applicable under stratigraphic conditions where lithological differences are not significant. However, in actual coal-bearing reservoir engineering, the hanging wall and footwall are often composed of combinations of different lithological layers. For example, the hanging wall may be a hard shale layer, while the footwall may be a porous coal-rock interlayer, with significant differences in elastic modulus, Poisson's ratio, creep characteristics, and fluid permeability. In such heterogeneous structures, stress does not diffuse proportionally to both sides, but is non-uniformly regulated by the interfacial properties during stress migration.
[0004] Current simulation systems are unable to accurately capture the asymmetric stress evolution process caused by lateral heterogeneity when dealing with this type of fault model. Especially under reloading or complex disturbance conditions, the fault slip limit values output by existing testing systems have significant deviations, making it difficult to accurately determine the risk of fault activation in actual engineering projects. Even with high-precision measurement modules, if the system software layer lacks the ability to model the heterogeneous stress migration mechanism, its evaluation results may still deviate from actual underground behavior. Summary of the Invention
[0005] To address the shortcomings of existing technologies, this invention provides a simulated fault stability testing system, which solves the problems mentioned in the background section.
[0006] To achieve the above objectives, the present invention provides the following technical solution: a simulated fault stability testing system, comprising a fault data acquisition module, a heterostructure identification module, a coupling factor calculation module, a stress evolution simulation module, a critical evolution judgment module, and a feedback adjustment evaluation module;
[0007] The fault data acquisition module collects the basic simulation parameters from the fault simulation device to form a basic dataset Dat, and establishes a spatial distribution mapping matrix Map for the physical property parameters of coal and rock in the hanging wall and footwall of the fault at different locations.
[0008] The heterostructure identification module identifies the structural boundaries of the hanging wall and footwall of the fault based on the basic dataset Dat and the spatial distribution mapping matrix Map, extracts heterogeneous difference parameters, and constructs the heterogeneous difference factor set Dif.
[0009] The coupling factor calculation module calculates the stress migration coupling coefficient Kst from the heterogeneous difference factor set Dif, and synchronously generates the coupling evolution dataset Seq in time series.
[0010] The stress evolution simulation module uses the coupled evolution dataset Seq as a basis to perform stress response simulation and calculate the stress distribution trajectory set Str in the fault region;
[0011] The critical evolution judgment module analyzes and obtains the slip trend index Sli based on the stress distribution trajectory set Str, and compares it with the preset structural migration sensitivity threshold Thr. Based on the comparison result, it determines to start the loading strategy optimization.
[0012] Preferably, the fault data acquisition module includes a simulation parameter acquisition unit;
[0013] The simulation parameter acquisition unit collects stress parameters, displacement parameters, temperature parameters, lithological parameters, and microstructure parameters during the simulation process using multiple types of sensors located at different positions in the fault simulation device.
[0014] The various types of sensors include stress-strain sensors, displacement sensors, thermal temperature sensors, and structural micro-probe sensors; simultaneously, these various types of sensors are synchronously connected to the system acquisition bus and perform synchronous sampling at preset time intervals;
[0015] After the data acquisition by the various types of sensors is completed, the raw sensor data is processed by a sliding window stability screening method to remove abnormal data. When the standard deviation of a continuous data segment exceeds a predetermined threshold, it is judged as abnormal and removed.
[0016] After removing outlier data, the linear unit conversion normalization method is used to eliminate the differences in the units of different parameters in the original sensor data. Then, the timestamps are linearly interpolated and aligned, and the data is packaged into structured frame data in chronological order to obtain the basic dataset Dat.
[0017] Preferably, the fault data acquisition module further includes a heterogeneous modeling processing unit;
[0018] The heterogeneous modeling processing unit extracts lithological parameters and microstructure parameters from the basic dataset Dat, and establishes a spatial attribute identification table by combining the spatial location index corresponding to each parameter point, dividing the parameter points into the hanging wall data group and the footwall data group.
[0019] Based on different lithological combinations and corresponding microstructural characteristics, the physical property estimation model built into the heterogeneous modeling unit is invoked to derive parameters and generate effectiveness indicators with equivalent mechanical significance. The effectiveness indicators include elastic modulus estimate, Poisson's ratio estimate, shear modulus estimate and density estimate, which are integrated into a representative physical property parameter set Est.
[0020] The lithological assemblage includes rocks, sandy mudstone, and shale; the microstructural features include fracture density and porosity.
[0021] After obtaining the representative physical property parameter set Est, the weighted inverse distance method built into the heterogeneous modeling processing unit performs spatial continuous modeling of the parameter points in the fault hanging wall and footwall data groups, realizes the smooth transition and internal interpolation of physical property indicators between different spatial regions, generates the function mapping model of each parameter point in the spatial structure, and outputs it as a multidimensional parameter matrix, which is marked as a spatial distribution mapping matrix Map.
[0022] Preferably, the heterostructure identification module includes a boundary identification unit;
[0023] The boundary identification unit is based on the spatial distribution mapping matrix Map. It extracts two parameters, the elastic modulus estimate and the shear modulus estimate, from the spatial distribution mapping matrix Map. It calculates the rate of change of parameters at continuous spatial locations using the first-order spatial difference method, which is used to identify regions of abrupt changes in physical properties in fault structures.
[0024] When there are several consecutive locations in space where the rate of change of the estimated elastic modulus or the estimated shear modulus simultaneously exceeds the preset abrupt change threshold and the direction of change is consistent, the several consecutive locations are determined to be candidate segments for structural abrupt change.
[0025] Meanwhile, the boundary identification unit calls the lithological combination in the basic dataset Dat to verify the stratigraphic composition of the spatial location of the candidate structural abrupt change section;
[0026] When the verification result shows that the lithological combination corresponding to the spatial location changes from a single coal-rock structure to a multi-component combination structure containing sandstone and shale, the spatial location is marked as the structural boundary determination point.
[0027] By integrating the spatial locations of all structural boundary determination points, a structural boundary index set Bidx is generated, which serves as the basis for dividing the hanging wall and footwall of the fault structure space.
[0028] Preferably, the heterostructure identification module includes a difference factor extraction unit;
[0029] The difference factor extraction unit uses the structural boundary index set Bidx as the basis for the spatial division of the fault structure, extracts the parameter content located in the hanging wall region and the footwall region of the fault respectively, extracts the elastic modulus estimate, Poisson's ratio estimate, shear modulus estimate and density estimate in the spatial distribution mapping matrix Map, calculates the difference between the average values of the hanging wall region and the footwall region of the fault, and generates an average value difference index; simultaneously calculates the ratio of the standard deviation of the hanging wall region and the footwall region of the fault, and generates a discrete difference index.
[0030] The average value difference indexes include the average value difference indexes of elastic modulus, Poisson's ratio, shear modulus, and density.
[0031] The dispersion difference indices include the dispersion difference index of elastic modulus, the dispersion difference index of Poisson's ratio, the dispersion difference index of shear modulus, and the dispersion difference index of density.
[0032] After integrating the mean difference index and the dispersion difference index, the heterogeneous difference factor set Dif is obtained.
[0033] Preferably, the coupling factor calculation module includes a path modeling unit;
[0034] The path modeling unit uses the input heterogeneous difference factor set Dif as the basic data structure, and presets a weight coefficient with physical meaning for each indicator in the heterogeneous difference factor set Dif. The specific value is set by the user, and they are uniformly combined into a weight set to obtain the parameter contribution score set Wset.
[0035] The parameter contribution score set Wset stores only the fixed weight coefficient of each difference index among the average difference index and the discrete difference index, which is used to represent the relative strength of the difference index in path coupling. It does not contain the index value itself and is a dimensionless value in the range of [0,1].
[0036] Then, the weight coefficients in the parameter contribution score set Wset are multiplied one-to-one with the corresponding difference indicators in the heterogeneous difference factor set Dif to obtain the indicator impedance value. After integration and processing, the weighted indicator result set Dres is formed.
[0037] Each item in the weighted index result set Dres represents the impedance strength influence corresponding to the difference index in the path space;
[0038] The weighted index result set Dres is then projected onto the spatial structure of the fault model. The path impedance superposition method is used to map all index impedance values onto the spatial connection path between the hanging wall region and the footwall region of the fault, forming a two-dimensional path influence matrix to obtain the stress migration path matrix Pmat.
[0039] Preferably, the coupling factor calculation module includes a coupling coefficient calculation unit;
[0040] The coupling coefficient calculation unit generates a stress migration coupling coefficient Kst that reflects the influence of the transverse heterogeneity of the structure on stress transmission based on the stress migration path matrix Pmat and the heterogeneity difference factor set Dif through a nonlinear mapping mechanism. The stress migration coupling coefficient Kst value obtained each time is numbered according to the time index and paired with the corresponding time point one by one, and summarized and recorded as a coupling evolution dataset Seq in time series format.
[0041] The stress migration coupling coefficient Kst is obtained by the following calculation formula:
[0042] Kst = log(1 + Davg) * (1 + Ztot);
[0043] In the formula, log represents the logarithmic function, Davg represents the normalized heterogeneity index, and Ztot represents the equivalent total impedance index.
[0044] Among them, the normalized heterogeneity index Davg is obtained by averaging the heterogeneity difference factor set Dif.
[0045] The equivalent total impedance index Ztot calculates the total impedance intensity of each coupling path in the fault structure by accumulating the index impedance values of each spatial connection path in the stress migration path matrix Pmat.
[0046] The nonlinear amplification mapping function is obtained by combining the normalized heterogeneous index Davg and the equivalent total impedance index Ztot.
[0047] Preferably, the stress evolution simulation module includes a stress loading scheduling unit;
[0048] The stress loading scheduling unit is based on the coupled evolution dataset Seq and combines the spatial node property estimation parameters contained in the spatial distribution mapping matrix Map to perform stress response simulation at each coupled evolution time step.
[0049] The stress response simulation process extracts the stress migration coupling coefficient Kst corresponding to each time step of the coupling evolution as the stress transfer modulation factor under the transverse heterostructure, and combines it with the four types of physical property estimation parameters at each spatial node position in the spatial distribution mapping matrix Map. The four types of physical property estimation parameters include elastic modulus estimation, Poisson's ratio estimation, shear modulus estimation and density estimation.
[0050] Based on the stress migration coupling coefficient Kst and four types of physical property estimation parameters, a coupling path network is constructed between spatially adjacent nodes. Stress transfer and aggregation calculations between nodes are performed to obtain the stress state update value of each node at the coupling evolution time step. The stress state update value includes principal stress and shear stress components.
[0051] The stress state update values of each node at each coupled evolution time step are continuously recorded and accumulated in the time dimension, and integrated into a stress change path structure with spatial location and time step as dual index dimensions, which is marked as stress distribution trajectory set Str.
[0052] Preferably, the critical evolution determination module includes a slip index extraction unit;
[0053] The slip index extraction unit, based on the slip trend index Sli, performs the following processing steps sequentially for all spatial nodes within the fault structure region:
[0054] Step 1: For each spatial node, extract the stress response record under the entire coupled evolution time series from the stress distribution trajectory set Str, and construct the principal stress change sequence Ser1 and the shear stress change sequence Ser2 respectively.
[0055] Step 2: Perform a first-order difference operation on the shear stress change sequence Ser2 to obtain its rate of change in the time direction, which is labeled as the shear stress growth rate sequence Seq1;
[0056] Step 3: At each coupled evolution time step, calculate the ratio of shear stress to principal stress at the corresponding time to form the shear principal stress ratio sequence Seq2;
[0057] Step 4: Perform sliding window abrupt change detection on the principal stress change sequence Ser1, extract the stress abrupt increase value within a local continuous time period, and obtain the principal stress abrupt increase factor sequence Seq3;
[0058] Step 5: Normalize the three trend factor sequences. Perform interval normalization on the shear stress growth rate sequence Seq1, the shear principal stress ratio sequence Seq2, and the principal stress surge factor sequence Seq3, respectively, and convert them to the [0,1] numerical interval. They are labeled as: Normalized Shear Stress Rate Factor Nor1, Normalized Shear Principal Stress Ratio Fluctuation Factor Nor2, and Normalized Principal Stress Surge Factor Nor3.
[0059] Step 6: Multiply the normalized shear stress rate factor Nor1, the normalized shear principal ratio fluctuation factor Nor2, and the normalized principal stress surge factor Nor3 by the weighting coefficients set by the user, perform a weighted average, and fuse them to obtain the final node-level index, which is marked as the index slip trend index Sli.
[0060] Preferably, the critical evolution judgment module includes a critical state judgment unit;
[0061] The critical state judgment unit compares the obtained slip trend index Sli with the preset structural migration sensitivity threshold Thr, and judges the current stress state based on the comparison result, outputting a warning state Sta, which is used to determine the start of loading strategy optimization.
[0062] The warning status Sta is obtained through the following determination method:
[0063] When the slippage trend index Sli ≥ the structural migration sensitivity threshold Thr, the node is marked as a slippage critical point, and the warning state of the node is marked Sta = 1;
[0064] When the slippage trend index Sli < the structural migration sensitivity threshold Thr, the node is marked as a non-critical point, and the warning state of the node is marked as Sta = 0.
[0065] Among them, the node status identifier Sta serves as the output for judging the local slip risk of each node;
[0066] Next, the nodes with node status identifier Sta=1 are counted. When the number exceeds the preset risk threshold, the loading strategy optimization is initiated, including adjusting the spatial distribution direction of the sliding critical point and adjusting the loading main axis direction to move it away from the ±15° range of the center vector of the sliding critical zone.
[0067] This invention provides a simulated fault stability testing system, which has the following advantages:
[0068] (1) The coupled model constructed based on the heterogeneous difference factor set Dif generates a stress migration coupling coefficient Kst that dynamically reflects the modulation effect of structural inconsistency on the stress transmission path. Using the coupled evolution dataset Seq driven by the stress migration coupling coefficient Kst, the system can obtain a complete stress distribution trajectory set Str during the simulated loading process, thus realistically reproducing the non-uniform stress concentration phenomenon formed in the fault region over time. Finally, by comparing the slip trend index Sli with the structural migration sensitivity threshold Thr derived from the stress migration coupling coefficient Kst, the system achieves active identification of the slip critical state and can trigger strategy-level loading path adjustment, completing the multi-stage coupled simulation of the entire process. This breaks the ideal assumption of consistent physical properties between the upper and lower walls of the fault in existing models, and has significant advantages in characterizing the stress evolution influence mechanism of transverse heterogeneous structures in faults, effectively improving the accuracy of fault slip trend judgment and the foresight of loading control.
[0069] (2) By assigning physically meaningful weight coefficients to various difference indices in the heterogeneous difference factor set Dif using path modeling units, a unified parameter contribution score set Wset is constructed. This set is then multiplied with the corresponding difference indices to generate a weighted index result set Dres, which characterizes the intensity of the difference's influence. This result set is then mapped to the fault spatial connection path using the path impedance superposition method, constructing an overall stress migration path matrix Pmat. This matrix realistically reflects the dynamic modulation effect of fault structural heterogeneity on the stress transmission process under different temporal loading conditions. This achieves continuous response modeling and adjustable parameter integration for structural differences on the stress transmission path, significantly improving the ability to capture the response of heterogeneous behavior within the fault system to the macroscopic mechanical evolution process. It provides precise driving parameters for subsequent slip criteria and loading strategy optimization.
[0070] (3) By constructing multiple characteristic trend sequences based on the stress distribution trajectory set Str, and generating a node-level slip trend index Sli through normalization fusion, and further comparing it with the structural migration sensitivity threshold Thr, an early warning state Sta is generated, which can realize spatially distributed critical slip determination. Through the fusion mechanism of trajectory evolution simulation and trend factor quantification calculation, synchronous modeling is achieved in three dimensions: mechanism identification, time response, and risk feedback, which can significantly improve the forward-looking identification capability and strategy response efficiency of fault system slip instability. Attached Figure Description
[0071] Figure 1 This is a schematic diagram of a simulated fault stability testing system according to the present invention;
[0072] Figure 2 This is a schematic diagram of the stress evolution trajectory at the nodes. Detailed Implementation
[0073] The technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without creative effort are within the scope of protection of the present invention.
[0074] Example 1
[0075] This invention provides a simulated fault stability testing system. Please refer to [link / reference]. Figure 1 It includes a fault data acquisition module, a heterostructure identification module, a coupling factor calculation module, a stress evolution simulation module, a critical evolution judgment module, and a feedback adjustment evaluation module;
[0076] The fault data acquisition module collects the basic simulation parameters from the fault simulation device to form a basic dataset Dat, and establishes a spatial distribution mapping matrix Map for the physical property parameters of coal and rock in the hanging wall and footwall of the fault at different locations.
[0077] The heterostructure identification module identifies the structural boundaries of the hanging wall and footwall of the fault based on the basic dataset Dat and the spatial distribution mapping matrix Map, extracts heterogeneous difference parameters, and constructs the heterogeneous difference factor set Dif.
[0078] The coupling factor calculation module calculates the stress migration coupling coefficient Kst from the heterogeneous difference factor set Dif, and synchronously generates the coupling evolution dataset Seq in time series.
[0079] The stress evolution simulation module uses the coupled evolution dataset Seq as a basis to perform stress response simulation and calculate the stress distribution trajectory set Str in the fault region;
[0080] The critical evolution judgment module analyzes and obtains the slip trend index Sli based on the stress distribution trajectory set Str, and compares it with the preset structural migration sensitivity threshold Thr. Based on the comparison result, it determines to start the loading strategy optimization.
[0081] In this embodiment, the collected basic dataset Dat and the established spatial distribution mapping matrix Map comprehensively cover the stress parameters, lithological parameters, and structural characteristics of the hanging wall and footwall at different spatial locations of the fault, effectively supporting subsequent identification of physical property inconsistencies and extraction of structural difference factors. Furthermore, through a coupling model constructed based on the heterogeneous difference factor set Dif, the resulting stress migration coupling coefficient Kst can dynamically reflect the modulation effect of structural inconsistencies on the stress transmission path. Using the coupling evolution dataset Seq driven by the stress migration coupling coefficient Kst, the system can obtain a complete stress distribution trajectory set Str during the simulated loading process, thereby realistically reproducing the non-uniform stress concentration phenomenon formed in the fault region over time. Finally, by comparing the slip trend index Sli with the structural migration sensitivity threshold Thr derived from the stress migration coupling coefficient Kst, the system achieves proactive identification of the slip critical state and can trigger strategy-level loading path adjustments, completing the multi-stage coupled simulation of the entire process. Breaking the ideal assumption of consistent physical properties between the hanging wall and footwall of a fault in existing models, this study achieves innovative breakthroughs in three core areas: structural identification, coupled modeling, and slip criteria. It has significant advantages, especially in characterizing the stress evolution mechanism of transverse heterogeneous structures in faults, and can effectively improve the accuracy of fault slip trend judgment and the foresight of loading control.
[0082] Example 2
[0083] Specifically: the fault data acquisition module includes a simulation parameter acquisition unit;
[0084] The simulation parameter acquisition unit collects stress parameters, displacement parameters, temperature parameters, lithological parameters, and microstructure parameters during the simulation process using multiple types of sensors located at different positions in the fault simulation device.
[0085] The various types of sensors include stress-strain sensors, displacement sensors, thermal temperature sensors, and structural micro-probe sensors; simultaneously, these various types of sensors are synchronously connected to the system acquisition bus and perform synchronous sampling at preset time intervals;
[0086] After the data acquisition by the various types of sensors is completed, the raw sensor data is processed by a sliding window stability screening method to remove abnormal data. The sliding window stability screening method evaluates the stability based on the change of standard deviation within a set time window. When the standard deviation of a continuous data segment exceeds a preset threshold, it is judged as abnormal and removed. This method can effectively remove data fluctuations caused by equipment micro-vibration and short-term interference, ensuring that the data reflects the true physical behavior.
[0087] After removing outlier data, the linear unit conversion normalization method is used to eliminate the differences in the units of different parameters in the original sensor data. Then, the timestamps are linearly interpolated and aligned, and the data is packaged into structured frame data in chronological order to obtain the basic dataset Dat.
[0088] The fault data acquisition module also includes a heterogeneous modeling processing unit;
[0089] The heterogeneous modeling processing unit extracts lithological parameters and microstructure parameters from the basic dataset Dat, and establishes a spatial attribute identification table by combining the spatial location index corresponding to each parameter point, dividing the parameter points into the hanging wall data group and the footwall data group.
[0090] Based on different lithological combinations and corresponding microstructural characteristics, the physical property estimation model built into the heterogeneous modeling unit is invoked to derive parameters and generate effectiveness indicators with equivalent mechanical significance. The effectiveness indicators include elastic modulus estimate, Poisson's ratio estimate, shear modulus estimate and density estimate, which are integrated into a representative physical property parameter set Est.
[0091] The lithological assemblage includes rocks, sandy mudstone, and shale; the microstructural features include fracture density and porosity.
[0092] After obtaining the representative physical property parameter set Est, the weighted inverse distance method built into the heterogeneous modeling processing unit performs spatial continuous modeling of the parameter points in the fault hanging wall and footwall data groups, realizes the smooth transition and internal interpolation of physical property indicators between different spatial regions, generates the function mapping model of each parameter point in the spatial structure, and outputs it as a multidimensional parameter matrix, which is marked as a spatial distribution mapping matrix Map.
[0093] In this embodiment, stress parameters, displacement parameters, temperature parameters, lithological parameters, and microstructure parameters are acquired simultaneously during the simulated loading process. A sliding window stability screening method is used to dynamically eliminate and optimize the original sensor data, ensuring the acquired basic dataset Dat possesses good stability and the ability to represent realistic physical responses. Based on this, a heterogeneous modeling unit performs spatial structural mapping of lithological and microstructure parameters. Combining hanging wall and footwall grouping and spatial location indexing, and utilizing a property estimation model and weighted inverse distance continuous interpolation, a representative set of physical property parameters Est with spatial consistency and equivalent mechanical significance is established. A spatial distribution mapping matrix Map, expressing the heterogeneous distribution characteristics of the fault, is then generated. This process not only significantly improves the construction accuracy of the spatial physical property model of the fault structure but also effectively solves the problem of the difficulty in coordinating the expression of lithological differences and microstructure parameters in traditional modeling methods. This makes the mechanical indices within the spatial range more continuous and complete, providing a stable and reliable input foundation for subsequent heterogeneous structure identification and coupled evolution calculations.
[0094] Example 3
[0095] Specifically: the heterostructure identification module includes a boundary identification unit;
[0096] The boundary identification unit is based on the spatial distribution mapping matrix Map. It extracts two parameters, the elastic modulus estimate and the shear modulus estimate, from the spatial distribution mapping matrix Map. It calculates the rate of change of parameters at continuous spatial locations using the first-order spatial difference method, which is used to identify regions of abrupt changes in physical properties in fault structures.
[0097] When there are several consecutive locations in space where the rate of change of the estimated elastic modulus or the estimated shear modulus simultaneously exceeds the preset abrupt change threshold and the direction of change is consistent, the several consecutive locations are determined to be candidate segments for structural abrupt change.
[0098] Meanwhile, the boundary identification unit calls the lithological combination in the basic dataset Dat to verify the stratigraphic composition of the spatial location of the candidate structural abrupt change section;
[0099] When the verification result shows that the lithological combination corresponding to the spatial location changes from a single coal-rock structure to a multi-component combination structure containing sandstone and shale, the spatial location is marked as the structural boundary determination point.
[0100] By integrating the spatial locations of all structural boundary determination points, a structural boundary index set Bidx is generated, which serves as the basis for dividing the hanging wall and footwall of the fault structure space.
[0101] The heterostructure identification module includes a differential factor extraction unit;
[0102] The difference factor extraction unit uses the structural boundary index set Bidx as the basis for the spatial division of the fault structure, extracts the parameter content located in the hanging wall region and the footwall region of the fault respectively, extracts the elastic modulus estimate, Poisson's ratio estimate, shear modulus estimate and density estimate in the spatial distribution mapping matrix Map, calculates the difference between the average values of the hanging wall region and the footwall region of the fault, and generates an average value difference index; simultaneously calculates the ratio of the standard deviation of the hanging wall region and the footwall region of the fault, and generates a discrete difference index.
[0103] The average value difference indexes include the average value difference indexes of elastic modulus, Poisson's ratio, shear modulus, and density.
[0104] The dispersion difference indices include the dispersion difference index of elastic modulus, the dispersion difference index of Poisson's ratio, the dispersion difference index of shear modulus, and the dispersion difference index of density.
[0105] After integrating the mean difference index and the dispersion difference index, the heterogeneous difference factor set Dif is obtained.
[0106] In this embodiment, a dual abrupt change identification mechanism based on the continuous variation behavior of physical property parameters and the combination of stratigraphic structure changes is employed. The boundary identification unit, based on the spatial distribution mapping matrix Map, performs first-order spatial difference on the estimated elastic modulus and shear modulus to accurately locate regions in space where the rate of parameter change continuously increases. This is verified by combining the compositional characteristics of lithological combinations in the basic dataset Dat, achieving cross-confirmation of abrupt change boundaries from two dimensions: structural behavior and compositional evolution. This ultimately forms the basic index set, structural boundary index set Bidx, for the spatial division of the hanging wall and footwall of the fault. Building upon this, the difference factor extraction unit further extracts four types of parameters—estimated elastic modulus, estimated Poisson's ratio, estimated shear modulus, and estimated density—within the hanging wall and footwall regions. By calculating the ratio of the difference in their average values to their standard deviations, a heterogeneous difference factor set Dif is formed, which can simultaneously measure the trend and dispersion of physical property migration. Compared to traditional methods that rely solely on a single parameter for qualitative analysis of fault heterogeneity, this approach achieves quantitative, index-based, and structured modeling representation of fault structural property inconsistencies, significantly enhancing the analytical accuracy and decision interpretation capabilities of subsequent coupled modeling and response analysis for lateral heterogeneity.
[0107] Example 4
[0108] Specifically: the coupling factor calculation module includes a path modeling unit;
[0109] The path modeling unit uses the input heterogeneous difference factor set Dif as the basic data structure, and presets a weight coefficient with physical meaning for each indicator in the heterogeneous difference factor set Dif. The specific value is set by the user, and they are uniformly combined into a weight set to obtain the parameter contribution score set Wset.
[0110] The parameter contribution score set Wset stores only the fixed weight coefficient of each difference index among the average difference index and the discrete difference index, which is used to represent the relative strength of the difference index in path coupling. It does not contain the index value itself and is a dimensionless value in the range of [0,1].
[0111] Then, the weight coefficients in the parameter contribution score set Wset are multiplied one-to-one with the corresponding difference indicators in the heterogeneous difference factor set Dif to obtain the indicator impedance value. After integration and processing, the weighted indicator result set Dres is formed.
[0112] Each item in the weighted index result set Dres represents the impedance strength influence corresponding to the difference index in the path space;
[0113] The weighted index result set Dres is then projected onto the spatial structure of the fault model. The path impedance superposition method is used to map all index impedance values onto the spatial connection path between the hanging wall region and the footwall region of the fault, and a two-dimensional path influence matrix is formed to obtain the stress migration path matrix Pmat.
[0114] The stress migration path matrix Pmat is a structure that describes the overall influence of the transverse heterogeneity of the fault structure on the stress transmission path, and is used to calculate the stress migration coupling coefficient Kst in subsequent calculations.
[0115] The coupling factor calculation module includes a coupling coefficient calculation unit;
[0116] The coupling coefficient calculation unit generates a stress migration coupling coefficient Kst that reflects the influence of the transverse heterogeneity of the structure on stress transmission based on the stress migration path matrix Pmat and the heterogeneity difference factor set Dif through a nonlinear mapping mechanism. The stress migration coupling coefficient Kst value obtained each time is numbered according to the time index and paired with the corresponding time point one by one, and summarized and recorded as a coupling evolution dataset Seq in time series format.
[0117] The stress migration coupling coefficient Kst is obtained by the following calculation formula:
[0118] Kst = log(1 + Davg) * (1 + Ztot);
[0119] In the formula, log represents the logarithmic function, Davg represents the normalized heterogeneity index, and Ztot represents the equivalent total impedance index.
[0120] Among them, the normalized heterogeneity index Davg is obtained by averaging the heterogeneity difference factor set Dif, forming a normalized dimensional index used to describe the overall level of lateral inconsistency in fault structure.
[0121] The equivalent total impedance index Ztot calculates the total impedance intensity of each coupling path in the fault structure by accumulating the index impedance values of each spatial connection path in the stress migration path matrix Pmat, and generates an overall path impedance characteristic index.
[0122] The nonlinear amplification mapping function is obtained by combining the normalized heterogeneous index Davg and the equivalent total impedance index Ztot.
[0123] In this embodiment, a path modeling unit assigns physically meaningful weight coefficients to various difference indices in the heterogeneous difference factor set Dif, constructing a unified parameter contribution score set Wset. This Wset is then multiplied with the corresponding difference indices to generate a weighted index result set Dres, which characterizes the intensity of the difference's influence. This Dres set is then mapped to the fault spatial connection path using a path impedance superposition method, constructing a holistic stress migration path matrix Pmat. Based on this, the system integrates path spatial impedance characteristics and the normalization level of the difference indices, calculating the stress migration coupling coefficient Kst through a nonlinear amplification mechanism. This forms a coupled evolution dataset Seq containing time-series structures, realistically reflecting the dynamic modulation effect of fault structural heterogeneity on the stress transmission process under different temporal loading conditions. Compared to traditional methods that treat heterogeneity as a static disturbance, this invention achieves continuous response modeling and adjustable parameter integration expression of structural differences on the stress transmission path. This significantly improves the ability to capture the response of heterogeneous behavior within the fault system to the macroscopic mechanical evolution process, providing precise driving parameters for subsequent slip criteria and loading strategy optimization.
[0124] Example 5
[0125] Please see Figure 1 and Figure 2 Specifically: the stress evolution simulation module includes a stress loading scheduling unit;
[0126] The stress loading scheduling unit is based on the coupled evolution dataset Seq and combines the spatial node property estimation parameters contained in the spatial distribution mapping matrix Map to perform stress response simulation at each coupled evolution time step.
[0127] The stress response simulation process extracts the stress migration coupling coefficient Kst corresponding to each time step of the coupling evolution as the stress transfer modulation factor under the transverse heterostructure, and combines it with the four types of physical property estimation parameters at each spatial node position in the spatial distribution mapping matrix Map. The four types of physical property estimation parameters include elastic modulus estimation, Poisson's ratio estimation, shear modulus estimation and density estimation.
[0128] Based on the stress migration coupling coefficient Kst and four types of physical property estimation parameters, a coupling path network is constructed between spatially adjacent nodes. Stress transfer and aggregation calculations between nodes are performed to obtain the stress state update value of each node at the coupling evolution time step. The stress state update value includes principal stress and shear stress components.
[0129] The stress state update values of each node at each coupled evolution time step are continuously recorded and accumulated in the time dimension, and integrated into a stress change path structure with spatial location and time step as dual index dimensions, which is marked as stress distribution trajectory set Str;
[0130] The stress distribution trajectory set Str represents the evolution path of stress at different nodes over time under the influence of transverse heterogeneous structures, and is used to reflect the non-uniform stress concentration trend of the fault during loading.
[0131] The critical evolution judgment module includes a slip index extraction unit;
[0132] The slip index extraction unit, based on the slip trend index Sli, performs the following processing steps sequentially for all spatial nodes within the fault structure region:
[0133] Step 1: For each spatial node, extract the stress response record under the entire coupled evolution time series from the stress distribution trajectory set Str, and construct the principal stress change sequence Ser1 and the shear stress change sequence Ser2 respectively.
[0134] Step 2: Perform a first-order difference operation on the shear stress change sequence Ser2 to obtain its rate of change in the time direction, which is labeled as the shear stress growth rate sequence Seq1, used to reflect the trend of accelerated growth of shear stress over time.
[0135] Step 3: At each coupling evolution time step, calculate the ratio of shear stress to principal stress at the corresponding time to form a shear principal stress ratio sequence Seq2. Perform maximum fluctuation amplitude calculation on this sequence to measure the degree of fluctuation in coupling strength.
[0136] Step 4: Perform sliding window abrupt change detection on the principal stress change sequence Ser1, extract the stress abrupt increase value within the local continuous time period, and obtain the principal stress abrupt increase factor sequence Seq3, which is used to describe the intensity and persistence of the local principal stress concentration area.
[0137] Step 5: Normalize the three trend factor sequences. Perform interval normalization on the shear stress growth rate sequence Seq1, the shear principal stress ratio sequence Seq2, and the principal stress surge factor sequence Seq3, respectively, and convert them to the [0,1] numerical interval. They are labeled as: Normalized Shear Stress Rate Factor Nor1, Normalized Shear Principal Stress Ratio Fluctuation Factor Nor2, and Normalized Principal Stress Surge Factor Nor3.
[0138] Step 6: Multiply the normalized shear stress rate factor Nor1, the normalized shear principal ratio fluctuation factor Nor2, and the normalized principal stress surge factor Nor3 by the weighting coefficients set by the user, perform a weighted average, and fuse them to obtain the final node-level index, which is marked as the index slip trend index Sli.
[0139] Table 1: Examples of three types of trend factor sequences:
[0140]
[0141]
[0142] The critical evolution judgment module includes a critical state judgment unit;
[0143] The critical state judgment unit compares the obtained slip trend index Sli with the preset structural migration sensitivity threshold Thr, and judges the current stress state based on the comparison result, outputting a warning state Sta, which is used to determine the start of loading strategy optimization.
[0144] The warning status Sta is obtained through the following determination method:
[0145] When the slippage trend index Sli ≥ the structural migration sensitivity threshold Thr, the node is marked as a slippage critical point, and the warning state of the node is marked Sta = 1;
[0146] When the slippage trend index Sli < the structural migration sensitivity threshold Thr, the node is marked as a non-critical point, and the warning state of the node is marked as Sta = 0.
[0147] Among them, the node status identifier Sta serves as the output for judging the local slip risk of each node;
[0148] Next, the nodes with node status identifier Sta=1 are counted. When the number exceeds the preset risk threshold, the loading strategy optimization is initiated, including adjusting the spatial distribution direction of the sliding critical point and adjusting the loading main axis direction to move it away from the ±15° range of the center vector of the sliding critical zone.
[0149] In this embodiment, the stress loading scheduling unit drives the loading process using the coupled evolution dataset Seq. Combining the physical property estimation parameters of each node in the spatial distribution mapping matrix Map, the current stress migration coupling coefficient Kst is used as a transmission modulation factor at each coupled evolution time step to construct a spatial coupled path network and perform stress transmission and aggregation operations between nodes. This generates a complete stress distribution trajectory set Str with both spatial and temporal indices, realistically reproducing the dynamic evolution of stress concentration trends in the fault region during loading. Based on this, the critical evolution judgment module, through the slip index extraction unit, constructs multiple characterizing trend sequences based on the stress distribution trajectory set Str, and generates a node-level slip trend index Sli through normalization fusion. This index is further compared with the structural migration sensitivity threshold Thr to generate a warning state Sta, enabling spatially distributed critical slip judgment. Compared to traditional schemes that cannot effectively describe the nonlinear stress change path and slip precursor behavior during fault loading, this method integrates trajectory evolution simulation and trend factor quantification to achieve simultaneous modeling in three dimensions: mechanism identification, time response, and risk feedback. This significantly improves the forward-looking identification capability and strategy response efficiency of fault system slip instability.
[0150] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and variations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.
Claims
1. A simulated fault stability testing system, characterized in that: It includes a fault data acquisition module, a heterostructure identification module, a coupling factor calculation module, a stress evolution simulation module, a critical evolution judgment module, and a feedback adjustment evaluation module; The fault data acquisition module collects the basic simulation parameters from the fault simulation device to form a basic dataset Dat, and establishes a spatial distribution mapping matrix Map for the physical property parameters of coal and rock in the hanging wall and footwall of the fault at different locations. The heterogeneous structure identification module identifies the structural boundaries of the hanging wall and footwall of the fault based on the basic dataset Dat and the spatial distribution mapping matrix Map, extracts the heterogeneous difference parameters, and constructs the heterogeneous difference factor set Dif. The coupling factor calculation module calculates the stress migration coupling coefficient Kst from the heterogeneous difference factor set Dif, and synchronously generates the coupling evolution dataset Seq in time series. The stress evolution simulation module uses the coupled evolution dataset Seq as a basis to perform stress response simulation and calculate the stress distribution trajectory set Str in the fault region; The critical evolution judgment module analyzes and obtains the slip trend index Sli based on the stress distribution trajectory set Str, and compares it with the preset structural migration sensitivity threshold Thr. Based on the comparison result, it determines to start the loading strategy optimization.
2. The simulated fault stability testing system according to claim 1, characterized in that: The fault data acquisition module includes a simulation parameter acquisition unit; The simulation parameter acquisition unit collects stress parameters, displacement parameters, temperature parameters, lithological parameters, and microstructure parameters during the simulation process using multiple types of sensors located at different positions in the fault simulation device. The various types of sensors include stress-strain sensors, displacement sensors, thermoelectric temperature sensors, and structural micro-probe sensors; Simultaneously, the various types of sensors are synchronously connected to the system's data acquisition bus and perform synchronous sampling at preset time intervals. After the data acquisition by the various types of sensors is completed, the raw sensor data is processed by a sliding window stability screening method to remove abnormal data. When the standard deviation of a continuous data segment exceeds a predetermined threshold, it is judged as abnormal and removed. After removing outlier data, the linear unit conversion normalization method is used to eliminate the differences in the units of different parameters in the original sensor data. Then, the timestamps are linearly interpolated and aligned, and the data is packaged into structured frame data in chronological order to obtain the basic dataset Dat.
3. The simulated fault stability testing system according to claim 2, characterized in that: The fault data acquisition module also includes a heterogeneous modeling processing unit; The heterogeneous modeling processing unit extracts lithological parameters and microstructure parameters from the basic dataset Dat, and establishes a spatial attribute identification table by combining the spatial location index corresponding to each parameter point, dividing the parameter points into the hanging wall data group and the footwall data group. Based on different lithological combinations and corresponding microstructural characteristics, the physical property estimation model built into the heterogeneous modeling unit is invoked to derive parameters and generate effectiveness indicators with equivalent mechanical significance. The effectiveness indicators include elastic modulus estimate, Poisson's ratio estimate, shear modulus estimate and density estimate, which are integrated into a representative physical property parameter set Est. The lithological assemblage includes rocks, sandy mudstone, and shale; the microstructural features include fracture density and porosity. After obtaining the representative physical property parameter set Est, the weighted inverse distance method built into the heterogeneous modeling processing unit performs spatial continuous modeling of the parameter points in the fault hanging wall and footwall data groups, realizes the smooth transition and internal interpolation of physical property indicators between different spatial regions, generates the function mapping model of each parameter point in the spatial structure, and outputs it as a multidimensional parameter matrix, which is marked as a spatial distribution mapping matrix Map.
4. The simulated fault stability testing system according to claim 3, characterized in that: The heterostructure identification module includes a boundary identification unit; The boundary identification unit is based on the spatial distribution mapping matrix Map. It extracts two parameters, the elastic modulus estimate and the shear modulus estimate, from the spatial distribution mapping matrix Map. It calculates the rate of change of parameters at continuous spatial locations using the first-order spatial difference method, which is used to identify regions of abrupt changes in physical properties in fault structures. When there are several consecutive locations in space where the rate of change of the estimated elastic modulus or the estimated shear modulus simultaneously exceeds the preset abrupt change threshold and the direction of change is consistent, the several consecutive locations are determined to be candidate segments for structural abrupt change. Meanwhile, the boundary identification unit calls the lithological combination in the basic dataset Dat to verify the stratigraphic composition of the spatial location of the candidate structural abrupt change section; When the verification result shows that the lithological combination corresponding to the spatial location changes from a single coal-rock structure to a multi-component combination structure containing sandstone and shale, the spatial location is marked as the structural boundary determination point. By integrating the spatial locations of all structural boundary determination points, a structural boundary index set Bidx is generated, which serves as the basis for dividing the hanging wall and footwall of the fault structure space.
5. The simulated fault stability testing system according to claim 4, characterized in that: The heterostructure identification module includes a differential factor extraction unit; The difference factor extraction unit uses the structural boundary index set Bidx as the basis for the spatial division of the fault structure, extracts the parameter content located in the hanging wall region and the footwall region of the fault respectively, extracts the elastic modulus estimate, Poisson's ratio estimate, shear modulus estimate and density estimate in the spatial distribution mapping matrix Map, calculates the difference between the average values of the hanging wall region and the footwall region of the fault, and generates an average value difference index; simultaneously calculates the ratio of the standard deviation of the hanging wall region and the footwall region of the fault, and generates a discrete difference index. The average value difference indexes include the average value difference indexes of elastic modulus, Poisson's ratio, shear modulus, and density. The dispersion difference indices include the dispersion difference index of elastic modulus, the dispersion difference index of Poisson's ratio, the dispersion difference index of shear modulus, and the dispersion difference index of density. After integrating the mean difference index and the dispersion difference index, the heterogeneous difference factor set Dif is obtained.
6. The simulated fault stability testing system according to claim 5, characterized in that: The coupling factor calculation module includes a path modeling unit; The path modeling unit uses the input heterogeneous difference factor set Dif as the basic data structure, and presets a weight coefficient with physical meaning for each indicator in the heterogeneous difference factor set Dif. The specific value is set by the user, and they are uniformly combined into a weight set to obtain the parameter contribution score set Wset. The parameter contribution score set Wset stores only the fixed weight coefficient of each difference index among the average difference index and the discrete difference index, which is used to represent the relative strength of the difference index in path coupling. It does not contain the index value itself and is a dimensionless value in the range of [0,1]. Then, the weight coefficients in the parameter contribution score set Wset are multiplied one-to-one with the corresponding difference indicators in the heterogeneous difference factor set Dif to obtain the indicator impedance value. After integration and processing, the weighted indicator result set Dres is formed. Each item in the weighted index result set Dres represents the impedance strength influence corresponding to the difference index in the path space; The weighted index result set Dres is then projected onto the spatial structure of the fault model. The path impedance superposition method is used to map all index impedance values onto the spatial connection path between the hanging wall region and the footwall region of the fault, forming a two-dimensional path influence matrix to obtain the stress migration path matrix Pmat.
7. The simulated fault stability testing system according to claim 6, characterized in that: The coupling factor calculation module includes a coupling coefficient calculation unit; The coupling coefficient calculation unit generates a stress migration coupling coefficient Kst that reflects the influence of the transverse heterogeneity of the structure on stress transmission based on the stress migration path matrix Pmat and the heterogeneity difference factor set Dif through a nonlinear mapping mechanism. The stress migration coupling coefficient Kst value obtained each time is numbered according to the time index and paired with the corresponding time point one by one, and summarized and recorded as a coupling evolution dataset Seq in time series format. The stress migration coupling coefficient Kst is obtained by the following calculation formula: Kst = log(1 + Davg) * (1 + Ztot); In the formula, log represents the logarithmic function, Davg represents the normalized heterogeneity index, and Ztot represents the equivalent total impedance index. Among them, the normalized heterogeneity index Davg is obtained by averaging the heterogeneity difference factor set Dif. The equivalent total impedance index Ztot calculates the total impedance intensity of each coupling path in the fault structure by accumulating the index impedance values of each spatial connection path in the stress migration path matrix Pmat. The nonlinear amplification mapping function is obtained by combining the normalized heterogeneous index Davg and the equivalent total impedance index Ztot.
8. The simulated fault stability testing system according to claim 7, characterized in that: The stress evolution simulation module includes a stress loading scheduling unit; The stress loading scheduling unit is based on the coupled evolution dataset Seq and combines the spatial node property estimation parameters contained in the spatial distribution mapping matrix Map to perform stress response simulation at each coupled evolution time step. The stress response simulation process extracts the stress migration coupling coefficient Kst corresponding to each time step of the coupling evolution as the stress transfer modulation factor under the transverse heterostructure, and combines it with the four types of physical property estimation parameters at each spatial node position in the spatial distribution mapping matrix Map. The four types of physical property estimation parameters include elastic modulus estimation, Poisson's ratio estimation, shear modulus estimation and density estimation. Based on the stress migration coupling coefficient Kst and four types of physical property estimation parameters, a coupling path network is constructed between spatially adjacent nodes. Stress transfer and aggregation calculations between nodes are performed to obtain the stress state update value of each node at the coupling evolution time step. The stress state update value includes principal stress and shear stress components. The stress state update values of each node at each coupled evolution time step are continuously recorded and accumulated in the time dimension, and integrated into a stress change path structure with spatial location and time step as dual index dimensions, which is marked as stress distribution trajectory set Str.
9. The simulated fault stability testing system according to claim 8, characterized in that: The critical evolution judgment module includes a slip index extraction unit; The slip index extraction unit, based on the slip trend index Sli, performs the following processing steps sequentially for all spatial nodes within the fault structure region: Step 1: For each spatial node, extract the stress response record under the entire coupled evolution time series from the stress distribution trajectory set Str, and construct the principal stress change sequence Ser1 and the shear stress change sequence Ser2 respectively. Step 2: Perform a first-order difference operation on the shear stress change sequence Ser2 to obtain its rate of change in the time direction, which is labeled as the shear stress growth rate sequence Seq1; Step 3: At each coupled evolution time step, calculate the ratio of shear stress to principal stress at the corresponding time to form the shear principal stress ratio sequence Seq2; Step 4: Perform sliding window abrupt change detection on the principal stress change sequence Ser1, extract the stress abrupt increase value within a local continuous time period, and obtain the principal stress abrupt increase factor sequence Seq3; Step 5: Normalize the three trend factor sequences. Perform interval normalization on the shear stress growth rate sequence Seq1, the shear principal stress ratio sequence Seq2, and the principal stress surge factor sequence Seq3, respectively, and convert them to the [0,1] numerical interval. They are labeled as: Normalized Shear Stress Rate Factor Nor1, Normalized Shear Principal Stress Ratio Fluctuation Factor Nor2, and Normalized Principal Stress Surge Factor Nor3. Step 6: Multiply the normalized shear stress rate factor Nor1, the normalized shear principal ratio fluctuation factor Nor2, and the normalized principal stress surge factor Nor3 by the weighting coefficients set by the user, perform a weighted average, and fuse them to obtain the final node-level index, which is marked as the index slip trend index Sli.
10. The simulated fault stability testing system according to claim 9, characterized in that: The critical evolution judgment module includes a critical state judgment unit; The critical state judgment unit compares the obtained slip trend index Sli with the preset structural migration sensitivity threshold Thr, and judges the current stress state based on the comparison result, outputting a warning state Sta, which is used to determine the start of loading strategy optimization. The warning status Sta is obtained through the following determination method: When the slippage trend index Sli ≥ the structural migration sensitivity threshold Thr, the node is marked as a slippage critical point, and the warning state of the node is marked Sta = 1; When the slippage trend index Sli < the structural migration sensitivity threshold Thr, the node is marked as a non-critical point, and the warning state of the node is marked as Sta = 0. Among them, the node status identifier Sta serves as the output for judging the local slip risk of each node; Next, the nodes with node status identifier Sta=1 are counted. When the number exceeds the preset risk threshold, the loading strategy optimization is initiated, including adjusting the spatial distribution direction of the sliding critical point and adjusting the loading main axis direction to move it away from the ±15° range of the center vector of the sliding critical zone.