A method and system for real-time high-precision acquisition of rock joint roughness coefficient based on a smooth joint model
Patent Information
- Application Number
- CN202611044446.3
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-07-14
- Publication Date
- 2026-10-16
AI Technical Summary
此类方法存在显著缺陷:一是时效性差,仅能获取某一静态时刻的JRC值,无法捕捉JRC随节理剪切磨损、粘结断裂、颗粒滑移等力学行为的动态演化过程;二是成本高昂,大型扫描设备的购置与维护成本高,且难以适配复杂工况下的原位测试
(1)本发明充分挖掘了SmoothJoint模型仅作用于预设节理面两侧颗粒接触的专属物理特性,将其作为节理接触的天然物理标签,结合颗粒分组归属校验与坐标范围校验构成双重校验机制,实现了节理接触的精准识别与筛选。相比现有仅依赖颗粒分组或几何边界的筛选方法,本发明有效避免了岩块内部颗粒接触被误判为节理接触的问题,显著降低了数据源中的冗余噪声,为后续形貌重构与JRC计算提供了高纯度的数据基础。
Smart Images

Figure CN122927005A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of numerical simulation technology in rock mechanics, specifically to a method and system for real-time high-precision acquisition of rock joint roughness coefficients based on the SmoothJoint model. Background Technology
[0002] Rock joints are ubiquitous discontinuities in rock masses, and their roughness coefficient (JRC) is a core parameter characterizing the geometry of joint surfaces and quantifying the shear strength, deformation characteristics, and permeability of jointed rock masses. Accurate acquisition of JRC is crucial for the safe design and risk management of rock mass engineering projects.
[0003] Traditional methods for obtaining joint morphology (JRC) primarily rely on offline detection. This involves acquiring point clouds of joint surface morphology using equipment such as 3D laser scanning, and then calculating the JRC value offline using specialized software. These methods have significant drawbacks: firstly, they lack timeliness, only providing JRC values at a static moment and failing to capture the dynamic evolution of JRC as joint shear wear, bond fracture, and particle slip occur; secondly, they are costly, with large scanning equipment incurring high purchase and maintenance costs, and are difficult to adapt to in-situ testing under complex conditions.
[0004] With the widespread application of Particle Flow Formula (PFC) in numerical simulation of rock mechanics, online JRC acquisition methods based on PFC have become a research hotspot. However, existing technologies still suffer from the following key technical bottlenecks, which limit the accuracy, real-time performance, and reliability of JRC acquisition, as follows: (1) High misjudgment rate of joint contact screening: Existing methods mostly screen joint contacts by particle grouping and geometric boundary range, without combining the dedicated mechanical model label of joint contact. This makes it easy to misjudge the particle contact inside the rock block as joint contact, resulting in a large amount of redundant noise in the data source of JRC calculation, which directly affects the subsequent morphology reconstruction and JRC calculation accuracy.
[0005] (2) Distortion of morphology construction benchmark: Existing technologies mostly fit the joint benchmark plane by discrete contact point cloud. The fitting process is easily affected by random outliers, and the fitted benchmark plane cannot be associated with the joint physical plane preset in PFC. This results in a systematic deviation between the macroscopic geometric benchmark of morphology reconstruction and the actual joint surface, affecting the accuracy of slope calculation results.
[0006] (3) It is difficult to balance real-time performance and accuracy: Existing methods adopt a fixed time step sampling strategy, which results in a large number of invalid calculations in large-scale particle models. This reduces the overall efficiency of PFC simulation and fails to focus on key mechanical stages where JRC changes significantly, such as the bonding fracture and shear slip stages. It is difficult to balance real-time monitoring and data accuracy.
[0007] The SmoothJoint (SJ) model is a contact model in PFC specifically designed to simulate the contact mechanics of rock joints. Its core characteristic is that it only acts on the particle contact on both sides of the preset joint surface, possessing a natural physical label for joint contact, which provides a technical basis for the accurate screening of joint contact.
[0008] However, there is currently no technology to fully exploit the unique characteristics of the SJ model and construct a real-time JRC acquisition method and system that balances accuracy, real-time performance, and robustness. Therefore, there is an urgent need for a real-time, high-precision acquisition method and system for rock joint roughness coefficients based on the SJ model to solve the aforementioned technical problems.
[0009] The information disclosed in this background section is intended only to enhance the understanding of the overall background of the invention and should not be construed as an admission or in any way implying that the information constitutes prior art known to those skilled in the art. Summary of the Invention
[0010] To address the aforementioned technical problems, this invention provides a method and system for real-time, high-precision acquisition of rock joint roughness coefficients based on the SmoothJoint model, thereby solving the problems mentioned in the background section. The method includes the following steps: A particle flow joint model containing a SmoothJoint model is constructed, and the SmoothJoint contact model is assigned to the particle contact on both sides of the joint plane to achieve the physical binding of SmoothJoint contact with the joint region. Collect and double-validate SmoothJoint contact data, filter out contacts with the model type SmoothJoint, perform double validation on the filtered contacts including particle grouping and coordinate range validation, remove abnormal contacts, and collect contact data that have passed double validation to form a structured dataset. To screen and purify effective contact points, the structured dataset is sequentially subjected to state filtering, statistical tests to remove outliers, dynamic mechanical threshold filtering, and clustering outlier removal to obtain a set of effective contact points. The joint morphology is adaptively digitally reconstructed by performing principal component analysis on the effective contact point set to automatically construct a local coordinate system, calculating the contact point density, and adaptively selecting an interpolation method based on the density threshold to generate an elevation matrix, followed by smoothing. To accurately calculate the rock joint roughness coefficient, the smoothed elevation matrix is extended by boundary, the root mean square of the slope is calculated using the second-order central difference method, and the rock joint roughness coefficient value is obtained by piecewise formula conversion based on localized calibration parameters. The three-level fault tolerance process handles abnormal data, monitors abnormal states in real time, and executes a three-level fault tolerance strategy based on the number of consecutive abnormal steps, including linear prediction, using historical valid values, and alarm degradation. Event-triggered sampling and real-time output: Real-time monitoring of mechanical event trigger thresholds; rock joint roughness coefficient calculation is initiated when any trigger threshold is met or the maximum sampling time step interval is reached; and data and visualization curves are output in real time. The calculation results are verified in multiple dimensions. Standard contour line verification, particle number convergence analysis and parameter sensitivity analysis are performed on the calculation results, and a verification report is output.
[0011] Preferably, the dual verification includes: Particle grouping and attribution verification: Verify whether two contact-related particles belong to the upper rock block group and the lower rock block group, respectively. Coordinate range verification verifies whether the coordinates of the contact point are within the preset geometric boundary of the model, and whether the accuracy of the coordinate values is consistent with the numerical calculation accuracy of the particle flow program. Furthermore, for contacts that pass dual verification, the coordinates of the contact point, normal contact force, tangential contact force, contact activation state, bond fracture state, and distance from the contact point to the joint plane are collected to form a structured dataset.
[0012] Preferably, the effective contact points for screening and purification specifically include: State filtering: Traverse the SmoothJoint contact dataset and retain only contacts that are active and have not broken. Grubbs test for outlier removal: Outliers were removed from the distance data from the contact point to the joint plane using a significance level of α=0.05; Dynamic mechanical threshold screening: Calculate the 5th percentile of the normal force of all SmoothJoint contacts at the current time step as the dynamic threshold, and only retain contacts with a normal force greater than 0 and not less than the threshold. DBSCAN clustering outlier removal: Set the neighborhood radius to twice the average particle radius and the minimum number of cluster samples to 3, retaining only the cluster core points and boundary points.
[0013] Preferably, in the adaptive digital reconstruction of joint morphology: After centering the effective contact point cloud, principal component analysis is performed. The centroid of the point cloud is taken as the origin, the first principal component is used as the joint orientation x-axis, the second principal component is used as the joint dip y-axis, and the third principal component is used as the joint normal z-axis. The coordinate system is orthogonalized and normalized to construct a local coordinate system. The approximate area of the joint surface is calculated based on the range of effective contact points on the x-axis and y-axis of the local coordinate system, and the contact point density is obtained by counting the number of contact points per unit area. Set the contact point density threshold to 2000 points / m 2 and 5000 points / m 2 When the density is ≥5000 points / m2 When using Kriging interpolation, with 2000 points / m 2 ≤Density<5000 points / m 2 Radial basis function interpolation is used when the density is <2000 points / m. 2 Nearest neighbor interpolation was used to generate the elevation matrix with a grid step size of 0.0005m. The elevation matrix is smoothed using a Gaussian filter with a 3×3 filter kernel and a standard deviation σ=1.
[0014] Preferably, in the accurate calculation of the rock joint roughness coefficient: Boundary mirror extension: The smoothed elevation matrix is mirrored and extended by the boundary with a width of 1 grid cell to eliminate the influence of boundary effects on slope calculation. The formula for calculating slope using the second-order central difference method is: x-direction slope ; y-direction slope ; in, grid_step =0.0005m is the grid step size; Total slope: ; Root mean square of slope: ; Where m and n are the number of rows and columns of the elevation matrix, respectively; Piecewise formula conversion of JRC: Based on the numerical simulation data fitting of the Barton standard profile, localized calibration parameters a and b are obtained. Piecewise formula is used to complete the JRC conversion, and the JRC value is limited to the range of 0~20.
[0015] Preferably, the three-level fault tolerance strategy is as follows: Level 1 fault tolerance: When the number of consecutive abnormal time steps is ≤3, extract the rock joint roughness coefficient data of the first 3 normal time steps to construct a linear fitting model JRC=k×step+b to predict the value of the current time step; Level 2 fault tolerance: When 3 < number of consecutive anomalies ≤ 7, the rock joint roughness coefficient value of the most recent effective time step is used; Level 3 fault tolerance: When the number of steps exceeds 7 during consecutive anomalies, the default value of 0 is output and an alarm is triggered. At the same time, the granular flow model status self-check is started. After the self-check passes, the calculation of rock joint roughness coefficient is restarted.
[0016] Preferably, in the event-triggered sampling: The mechanical event triggering thresholds include: shear displacement increment ≥ 0.00001m, and the percentage of the increase in the number of bond fractures to the total number of SmoothJoint contacts ≥ 5%; The maximum sampling time step interval is 5000 time steps; The method also includes automatically saving a snapshot of the particle flow model and intermediate data of rock joint roughness coefficients every 10,000 time steps for data recovery after simulation interruption.
[0017] A real-time, high-precision system for acquiring rock joint roughness coefficients based on the SmoothJoint model includes the following functional modules: Particle flow model building module: Used to integrate the particle flow program application programming interface, providing a visual configuration interface for SmoothJoint model parameters, and automatically completing the generation, grouping and allocation of rock particles and SmoothJoint models; Data acquisition and verification module: Built-in dual verification logic is used to automatically filter SmoothJoint contacts and complete the verification, and output a standardized SmoothJoint contact dataset; Contact point screening and purification module: encapsulates the Grubbs test and DBSCAN clustering algorithm to achieve automated screening and purification of effective contact points; Adaptive topography reconstruction module: Built-in principal component analysis coordinate system construction, multi-algorithm interpolation and Gaussian filtering functions, used to automatically output the smoothed elevation matrix; The precise calculation module for rock joint roughness coefficient is used to automate the execution of boundary extension, second-order difference slope calculation and piecewise formula conversion, and output the rock joint roughness coefficient value. Anomaly tolerance module: Used to monitor abnormal states in real time, perform fault tolerance processing according to hierarchical strategies, generate anomaly logs and trigger alarms; Intelligent sampling and output module: used to monitor mechanical event trigger thresholds, automatically perform sampling, output structured logs and dynamic curves, and manage model snapshots and intermediate data; Multi-dimensional verification module: used to automatically complete standard contour verification, particle number convergence analysis and parameter sensitivity analysis, and generate a visual verification report.
[0018] Preferably, the dual verification logic in the data acquisition and verification module includes: Particle grouping and attribution verification: Verify whether two contact-related particles belong to the upper rock block group and the lower rock block group, respectively; Coordinate range verification: Verify whether the coordinates of the contact point are within the preset geometric boundary of the model, and whether the accuracy of the coordinate values is consistent with the numerical calculation accuracy of the particle flow program; In addition, for contacts that pass dual verification, the coordinates of the contact point, normal contact force, tangential contact force, contact activation state, bond fracture state, and distance from the contact point to the joint plane are collected to form a structured dataset.
[0019] Preferably, the hierarchical strategy of the fault tolerance module is as follows: Level 1 fault tolerance: When the number of consecutive abnormal steps is ≤3, the current value is linearly estimated based on the data of the previous 3 normal time steps; Level 2 fault tolerance: When the number of consecutive anomalies is ≤7, the most recent valid value is used; Level 3 fault tolerance: When the number of steps exceeds 7 during consecutive anomalies, output the default value of 0 and trigger an alarm, while simultaneously starting the granular flow model state self-check.
[0020] The present invention provides a method and system for real-time high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model, which has the following beneficial effects: (1) This invention fully exploits the unique physical characteristic of the SmoothJoint model, which only applies to particle contact on both sides of a predefined joint surface. This characteristic is used as a natural physical label for joint contact. Combined with particle grouping and coordinate range verification, a dual verification mechanism is established, enabling accurate identification and screening of joint contact. Compared to existing screening methods that rely solely on particle grouping or geometric boundaries, this invention effectively avoids the problem of misidentifying particle contact within rock blocks as joint contact, significantly reduces redundant noise in the data source, and provides a high-purity data foundation for subsequent morphology reconstruction and JRC calculation.
[0021] (2) This invention uses principal component analysis to automatically construct a local coordinate system for the joint, with the centroid of the point cloud as the origin and the first, second, and third principal components as the strike axis, dip axis, and normal axis of the joint, respectively, thus avoiding subjective biases caused by manually setting the coordinate system. Simultaneously, based on the contact point density, Kriging interpolation, radial basis function interpolation, or nearest neighbor interpolation is adaptively selected to generate the elevation matrix, and Gaussian filtering is used for smoothing, achieving high-precision digital reconstruction of the joint morphology. This method significantly improves the geometric consistency of the elevation matrix and its matching degree with the actual joint surface, and significantly reduces the numerical error in slope calculation.
[0022] (3) This invention designs a three-level fault-tolerant processing mechanism to address potential data anomalies during numerical simulation. Based on the number of steps in case of continuous anomalies, three progressive fault-tolerant strategies are implemented: linear prediction, reuse of historical valid values, and alarm degradation. This effectively ensures the continuity of JRC monitoring and the stability of program operation. Even in the event of a temporary loss of data source or algorithm execution anomalies, the system can still output a reasonable JRC prediction value, avoiding simulation interruption and significantly reducing the failure rate during continuous operation.
[0023] (4) This invention adopts an event-triggered sampling strategy to replace the traditional fixed-step sampling. By monitoring two types of mechanical events in real time, namely the shear displacement increment and the bond fracture number increment, JRC calculation is started only at the critical stage when the joint morphology changes significantly. At the same time, a maximum sampling interval is set as a fallback mechanism, which avoids the impact of invalid calculation on simulation efficiency and ensures high-resolution capture of key mechanical processes, thus achieving a good balance between real-time performance and data accuracy.
[0024] (5) The present invention obtains localized calibration parameters based on numerical simulation data fitting of the Barton standard profile, uses a piecewise formula to complete the JRC conversion, and limits the JRC value to a reasonable range. This conversion method can adjust the parameters according to the indoor test data of different rock types, has strong adaptability, and controls the error between the calculated JRC value and the measured value within a small range, which significantly improves the accuracy compared with the general empirical formula.
[0025] (6) This invention provides a complete multi-dimensional verification system, including Barton standard contour verification, particle number convergence analysis, and parameter sensitivity analysis, ensuring the accuracy of JRC calculation results, model convergence, and robustness to parameter fluctuations. Meanwhile, the automatic saving mechanism for model snapshots and intermediate data enables lossless data recovery after simulation interruptions, providing reliable quality assurance for engineering applications. Attached Figure Description
[0026] Figure 1 This is a flowchart illustrating the overall technical framework of the method of the present invention; Figure 2 This is a functional module architecture diagram of the system of the present invention; Figure 3 This is a logical flowchart of the fault tolerance mechanism of the present invention. Detailed Implementation
[0027] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.
[0028] To verify and analyze the effectiveness of the foam RPB fluid dynamics prediction method based on the porous media resistance model provided by this invention, the following two embodiments are selected to verify the effectiveness of the method of this invention.
[0029] Example 1: Specific Implementation Steps of the Method of the Invention
[0030] S1. Construct a PFC granular flow joint model including the SJ model.
[0031] The core of this step is to achieve physical bonding between the SJ contact and the joint region, laying the foundation for accurate contact screening in the future. The specific operation is as follows: (1) SJ model parameter configuration: In the PFC program, first set the plane equation parameters of the joint plane (satisfying ax+by+cz+d=0), and assign a unique joint ID to the joint; then create an SJ contact model with parallel bonding, and configure the model mechanical parameters: normal stiffness 1e9N / m, tangential stiffness 5e8N / m, bond tensile strength 5e6Pa, bond cohesion 1e6Pa, internal friction angle 30° (which can be adjusted according to the indoor test data of the target rock type); (2) Formation and grouping of rock fragments: Rock fragments were formed above and below the joints, with a minimum radius of 0.001 m, a maximum radius of 0.0015 m, and a density of 2650 kg / m³. 3 The upper rock block particles are uniformly classified into the "upper_block" group, and the lower rock block particles are uniformly classified into the "lower_block" group, to ensure that the particles are distributed within a 0.02m range on both sides of the joint plane. (3) SJ model precise allocation: Traverse all particle contacts and only allocate the SJ model to contacts that simultaneously meet the following two conditions: First, the two particles in contact belong to the “upper_block” and “lower_block” groups (cross-block grouping); Second, the distance from the contact point to the joint plane is ≤ 1.2 times the average particle radius (1.2×(0.001+0.0015) / 2=0.0015m). (4) Model balancing: Enable PFC gravity field (gravitational acceleration 9.81 m / s²), run 10,000 time steps to achieve initial particle balancing and ensure the stability of the initial state of the model.
[0032] S2 and SJ contact data acquisition and dual verification
[0033] The core of this step is to obtain a highly reliable SJ contact data source, and the specific steps are as follows: SJ Contact Screening: All contacts are traversed through the PFCPythonAPI to screen out contacts with model type SJ and joint ID that matches the target joint, forming a preliminary contact set; Dual verification is performed: ① Particle grouping and attribution verification: Traverse each contact in the initial screening contact set and verify whether the two associated particles belong to the "upper_block" and "lower_block" groups respectively, and eliminate misassigned SJ contacts within the rock block; ② Coordinate range verification: Verify whether the contact point coordinates are within the preset geometric boundaries of the PFC model (-0.1m≤x≤0.1m, -0.1m≤y≤0.1m, -0.1m≤z≤0.1m), and whether the coordinate values retain 6 decimal places (consistent with the PFC numerical calculation precision), and eliminate abnormal contacts with distorted coordinate values; Data Acquisition: For contacts that pass dual verification, parameters such as actual contact point coordinates (average of the actual contact points of the two particles), normal contact force, tangential contact force, contact activation state, bond fracture state, and distance from the contact point to the joint plane are collected to form a structured SJ contact dataset. Simultaneously, the PFC simulation physical timestamp is recorded to ensure data time sequence consistency.
[0034] S3. Screening and purification of effective contact points
[0035] The core of this step is to remove noise and invalid data from the dataset and obtain a set of highly reliable and valid contact points. The specific steps are as follows: State filtering: Traverse the SJ contact dataset, retain only contacts that are active and have not broken, and remove failed contacts; Grubbs test for outlier removal: The distance data from the contact point to the joint plane is extracted, and the Grubbs test with a significance level of α=0.05 is used to identify and remove outliers to ensure the rationality of the geometric distribution of the contact points; Dynamic mechanical threshold screening: Calculate the 5th percentile of the normal force of all SJ contacts at the current time step as the dynamic threshold, and only retain contacts with normal force > 0 (excluding tensile contacts) and ≥ the threshold, while removing weak contacts; DBSCAN Clustering Outlier Removal: The contact point cloud that has been screened above is processed using the DBSCAN clustering algorithm. The neighborhood radius is set to twice the average particle radius (0.0025m) and the minimum number of cluster samples is 3. Only the core points and boundary points of the clusters are retained, and isolated outliers are removed to obtain the effective set of contact points.
[0036] S4. Adaptive digital reconstruction of joint morphology
[0037] The core of this step is to achieve high-precision digitization of joint morphology, and the specific operation is as follows: PCA automatically constructs a local coordinate system: the effective contact point cloud is centered (the centroid coordinates of the point cloud are subtracted); principal component analysis (PCA) is performed on the centered point cloud to extract the first, second, and third principal components; with the centroid of the point cloud as the origin, the first principal component is the joint strike x-axis, the second principal component is the joint dip y-axis, and the third principal component is the joint normal z-axis; The coordinate system is orthogonalized and normalized: x-axis = x-axis initial vector - (dot product of x-axis initial vector and z-axis) × z-axis, y-axis = cross product of z-axis and x-axis. The x-axis, y-axis and z-axis are normalized by dividing by their own magnitudes, and finally the joint local coordinate system is obtained. Contact point density calculation: Based on the range of effective contact points on the x and y axes of the local coordinate system, calculate the approximate area of the joint surface (area = (x_max - x_min) × (y_max - y_min)), count the number of contact points per unit area, and obtain the contact point density; Adaptive interpolation is used to generate the elevation matrix: Contact point density thresholds are set to 2000 points / m² and 5000 points / m². The interpolation method is selected according to the following rules: Density ≥ 5000 points / m²: Kriging interpolation is used to preserve the spatial correlation of joint morphology; 2000 points / m² ≤ density < 5000 points / m²: RBF interpolation is used to balance interpolation accuracy and computational efficiency; Density < 2000 points / m²: Nearest neighbor interpolation is used to ensure the robustness of the interpolation algorithm. The regular grid generated by interpolation has a step size of 0.0005m, resulting in the joint morphology elevation matrix. Elevation matrix smoothing: Gaussian filtering with a 3×3 filter kernel and a standard deviation σ=1 is used to smooth the elevation matrix, reducing the impact of numerical noise on subsequent slope calculations.
[0038] S5 and JRC precise calculations
[0039] The core of this step is to achieve accurate conversion from root mean square slope to JRC. The specific operation is as follows: Boundary mirror extension: The smoothed elevation matrix is mirrored and extended by the boundary with a width of 1 grid cell to eliminate the influence of boundary effects on slope calculation. Root mean square slope calculation: The slope is calculated using the second-order central difference method. x-direction slope ; y-direction slope ; in, grid_step =0.0005m is the grid step size; Total slope: ; Root mean square of slope: ; Where m and n are the number of rows and columns of the elevation matrix, respectively; Piecewise formula conversion of JRC: Based on the numerical simulation data of the Barton standard profile, the localized calibration parameters a=32.2 and b=32.47 were obtained by fitting. The piecewise formula was used to complete the JRC conversion: When the joint surface roughness reaches the upper limit, JRC=20. At the same time, the JRC value is limited to a reasonable range of 0 to 20.
[0040] S6. Three-level fault tolerance handling for abnormal data
[0041] The core of this step is to ensure the continuity of JRC monitoring. The specific fault tolerance strategy is as follows: Anomaly detection: Real-time monitoring of the number of valid contact points (<5 is an anomaly), interpolation algorithm execution status (error is an anomaly), and slope calculation results (NaN / Inf is an anomaly), and counting the number of steps when there are consecutive anomalies; Three-level fault tolerance execution: Level 1 fault tolerance (linear prediction, consecutive outliers ≤ 3): Extract the JRC data from the first 3 normal time steps, construct a linear fitting model JRC=k×step+b, predict the JRC value of the current time step, log it as "interpolation prediction", and the confidence level is "medium"; Level 2 fault tolerance (using historical values, 3 < consecutive outliers ≤ 7): Directly use the JRC value of the most recent valid time step, log it as "using historical valid values", and the confidence level is "low"; Level 3 fault tolerance (alarm degradation, consecutive anomalies > 7): Output JRC default value 0 and trigger pop-up / email alarm, and start PFC model status self-check. The self-check content includes particle overlap rate (>1% is abnormal), loading boundary displacement / force deviation (>10% is abnormal), and SJ model parameter validity. After the self-check passes, restart JRC calculation. Exception logging: Records exception time step, exception type, handling strategy, JRC output value, confidence level and other information in JSON format to enable exception tracing.
[0042] S7, Event-triggered sampling and real-time output
[0043] The core of this step is to balance real-time monitoring with computational efficiency, and the specific operation is as follows: Mechanical event monitoring: Real-time monitoring of two trigger thresholds: one is shear displacement increment ≥ 0.00001m (obtained through loading wall displacement sensor), and the other is the proportion of the increase in the number of bond fractures to the total number of SJ contacts ≥ 5% (obtained by statistical analysis of SJ contact bonding status). Intelligent sampling: JRC calculation is started when any trigger threshold is met, and the maximum sampling time step interval is set to 5000 to avoid monitoring interruption due to the failure of mechanical events to be triggered; Real-time output: Data log: JRC values, shear displacement, number of bond fractures, number of effective contact points, confidence level, and other data are written to a JSON structured log; Visualization curves: JRC-time step and JRC-shear displacement dynamic curves are plotted using Matplotlib to intuitively display the JRC evolution process; Data recovery: Every 10,000 hours, the PFC model snapshot (.sav) and JRC intermediate data (.npy) are automatically saved. After the simulation is interrupted, the most recent snapshot can be read to restore the model state and recalculate the JRC data for the interrupted period.
[0044] Multi-dimensional verification of S8 and JRC results
[0045] The core of this step is to ensure the accuracy and reliability of the JRC results. The specific verification content is as follows: Barton Standard Profile Validation: Construct Barton standard profile PFC models with JRC=0, 5, 10, 15, 20, and calculate JRC using this method. The absolute error and relative error between the calculated value and the standard value should be ≤1 and ≤5% respectively. Particle number convergence analysis: The particle number was adjusted to 1000, 5000, 10000, and 20000 in turn, and the JRC value was calculated for different particle sizes. The JRC change rate was ≤2% when the particle number was ≥10000, which was considered as model convergence. Parameter sensitivity analysis: By using the controlled variable method, the range of values of key parameters (normal stiffness, bond strength, internal friction angle) of the SJ model within ±10% was traversed, and the influence coefficient of parameter changes on JRC was calculated. The influence coefficient is required to be ≤2% to ensure the robustness of JRC results to parameter fluctuations. Verification report output: Integrate the above verification results and output a verification report on the accuracy, convergence, and robustness of the JRC results.
[0046] Example 2: Detailed Implementation of the Method of the Invention
[0047] (I) Experimental Preparation
[0048] Hardware and software environment: Hardware consists of an Intel Xeon Gold 6248R CPU, 128GB DDR4 memory, an NVIDIA RTX 3090 GPU, and a 1TB SSD; Software consists of PFC 6.0 granular flow program, Python 3.9, and Python dependencies include NumPy 1.23.5, SciPy 1.9.3, Scikit-learn 1.1.3, PyKrige 1.7.0, and Matplotlib 3.6.2. Model parameter configuration: joint plane z=0, SJ model normal stiffness 1e9N / m, tangential stiffness 5e8N / m, total number of particles 5000, normal stress 5MPa, shear rate 0.000005m / hour.
[0049] (II) Implementation Process
[0050] PFC Model Construction: Start the PFC model construction module, complete the SJ model configuration, rock particle generation and grouping, SJ model precise allocation, and run 10,000 time steps to achieve particle balance; Data Acquisition and Cleaning: The data acquisition and verification module automatically completes SJ contact screening and double verification, while the contact point screening and cleaning module performs multi-level screening to obtain a set of valid contact points; Topography Reconstruction and JRC Calculation: The adaptive topography reconstruction module constructs a local coordinate system using PCA with a contact point density of 6000 points / ㎡, automatically selects Kriging interpolation to generate an elevation matrix and completes smoothing; the JRC accurate calculation module calculates the root mean square slope as 0.12, and converts it to JRC=12.3 using a piecewise formula. Intelligent sampling and fault tolerance: The first sampling is triggered when the shear displacement increment reaches 0.00001m; when the simulation reaches 20000 time steps, there are only 3 effective contact points, and the first level of fault tolerance is activated, with linear prediction JRC=12.2; Results verification: The calculated value of the model with Barton standard JRC=10 is 9.8, with an absolute error of 0.2 and a relative error of 2%; when the number of particles increases to 10,000, the JRC change rate is 1.5%, which meets the convergence requirements; when the bond strength of the SJ model fluctuates by 10%, the JRC change rate is 1.8%, which meets the robustness standard.
[0051] (III) Case 1: Real-time monitoring of granite joints using JRC
[0052] For a shear numerical experiment on a deeply buried granite joint, the method of this invention was used for real-time JRC monitoring. Initially, JRC = 11.5. When the shear displacement reached 0.005 mm, the proportion of bonded fractures reached 5%, triggering JRC calculation. At this point, JRC dropped to 9.8, accurately reflecting the reduction in roughness of the joint surface due to shear wear. During the simulation, an interpolation algorithm failure occurred once. Secondary fault tolerance was activated, maintaining JRC = 9.6. After the anomaly was resolved, it automatically switched to the real-time calculated value. The entire process was uninterrupted. The final JRC evolution curve showed a correlation coefficient of 0.95 with the indoor test results, verifying the engineering applicability of the method.
[0053] Case 2: Localization and Calibration of JRC Formula for Sandstone Joints
[0054] For sandstone joints, five sets of numerical simulation data of Barton standard profiles were collected (RMS_slope were 0.02, 0.05, 0.1, 0.3, and 0.8, corresponding to standard JRC values of 2, 5, 8, 12, and 18). Through the calibration function of the multi-dimensional verification module, localized parameters a=30.5 and b=30.2 were obtained by fitting, with a goodness of fit R²=0.98. The JRC of sandstone joints calculated using these parameters had a root mean square error of 1.2 compared with the measured values from the laboratory test, which is 65.7% more accurate than the general formula.
[0055] While the specific embodiments of the present invention have been described above in conjunction with the accompanying drawings, this is not intended to limit the scope of protection of the present invention. Those skilled in the art should understand that various modifications or variations that can be made by those skilled in the art without creative effort based on the technical solutions of the present invention are still within the scope of protection of the present invention.
Claims
1. A method for real-time, high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model, characterized in that, Includes the following steps: A particle flow joint model containing a SmoothJoint model is constructed, and the SmoothJoint contact model is assigned to the particle contact on both sides of the joint plane to achieve the physical binding of SmoothJoint contact with the joint region. Collect and double-validate SmoothJoint contact data, filter out contacts with the model type SmoothJoint, perform double validation on the filtered contacts including particle grouping and coordinate range validation, remove abnormal contacts, and collect contact data that have passed double validation to form a structured dataset. To screen and purify effective contact points, the structured dataset is sequentially subjected to state filtering, statistical tests to remove outliers, dynamic mechanical threshold filtering, and clustering outlier removal to obtain a set of effective contact points. The joint morphology is adaptively digitally reconstructed by performing principal component analysis on the effective contact point set to automatically construct a local coordinate system, calculating the contact point density, and adaptively selecting an interpolation method based on the density threshold to generate an elevation matrix, followed by smoothing. To accurately calculate the rock joint roughness coefficient, the smoothed elevation matrix is extended by boundary, the root mean square of the slope is calculated using the second-order central difference method, and the rock joint roughness coefficient value is obtained by piecewise formula conversion based on localized calibration parameters. The three-level fault tolerance process handles abnormal data, monitors abnormal states in real time, and executes a three-level fault tolerance strategy based on the number of consecutive abnormal steps, including linear prediction, using historical valid values, and alarm degradation. Event-triggered sampling and real-time output: Real-time monitoring of mechanical event trigger thresholds; rock joint roughness coefficient calculation is initiated when any trigger threshold is met or the maximum sampling time step interval is reached, and data and visualization curves are output in real time. The calculation results are verified in multiple dimensions. Standard contour line verification, particle number convergence analysis and parameter sensitivity analysis are performed on the calculation results, and a verification report is output.
2. The method for real-time high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model according to claim 1, characterized in that, The dual verification includes: Particle grouping and attribution verification: Verify whether two contact-related particles belong to the upper rock block group and the lower rock block group, respectively. Coordinate range verification verifies whether the coordinates of the contact point are within the preset geometric boundary of the model, and whether the accuracy of the coordinate values is consistent with the numerical calculation accuracy of the particle flow program. Furthermore, for contacts that pass dual verification, the coordinates of the contact point, normal contact force, tangential contact force, contact activation state, bond fracture state, and distance from the contact point to the joint plane are collected to form a structured dataset.
3. The method for real-time high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model according to claim 1, characterized in that, The effective contact points for screening and purification specifically include: State filtering: Traverse the SmoothJoint contact dataset and retain only contacts that are active and have not broken. Grubbs test for outlier removal: Outliers were removed from the distance data from the contact point to the joint plane using a significance level of α=0.05; Dynamic mechanical threshold screening: Calculate the 5th percentile of the normal force of all SmoothJoint contacts at the current time step as the dynamic threshold, and only retain contacts with a normal force greater than 0 and not less than the threshold. DBSCAN clustering outlier removal: Set the neighborhood radius to twice the average particle radius and the minimum number of cluster samples to 3, retaining only the cluster core points and boundary points.
4. The method for real-time high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model according to claim 1, characterized in that, In the adaptive digital reconstruction of the joint morphology: After centering the effective contact point cloud, principal component analysis is performed. The centroid of the point cloud is taken as the origin, the first principal component is used as the joint orientation x-axis, the second principal component is used as the joint dip y-axis, and the third principal component is used as the joint normal z-axis. The coordinate system is orthogonalized and normalized to construct a local coordinate system. The approximate area of the joint surface is calculated based on the range of effective contact points on the x-axis and y-axis of the local coordinate system, and the contact point density is obtained by counting the number of contact points per unit area. Set the contact point density threshold to 2000 points / m 2 and 5000 points / m 2 When the density is ≥5000 points / m 2 When using Kriging interpolation, with 2000 points / m 2 ≤Density<5000 points / m 2 Radial basis function interpolation is used when the density is <2000 points / m. 2 Nearest neighbor interpolation was used to generate the elevation matrix with a grid step size of 0.0005m. The elevation matrix is smoothed using a Gaussian filter with a 3×3 filter kernel and a standard deviation σ=1.
5. The method for real-time high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model according to claim 1, characterized in that, In the precise calculation of the rock joint roughness coefficient: Boundary mirror extension: The smoothed elevation matrix is mirrored and extended by the boundary with a width of 1 grid cell to eliminate the influence of boundary effects on slope calculation. The formula for calculating slope using the second-order central difference method is: x-direction slope ; y-direction slope ; in, grid_step =0.0005m is the grid step size; Total slope: ; Root mean square of slope: ; Where m and n are the number of rows and columns of the elevation matrix, respectively; Piecewise formula conversion of JRC: Based on the numerical simulation data fitting of the Barton standard profile, localized calibration parameters a and b are obtained. Piecewise formula is used to complete the JRC conversion, and the JRC value is limited to the range of 0~20.
6. The method for real-time high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model according to claim 1, characterized in that, The three-level fault tolerance strategy is as follows: Level 1 fault tolerance: When the number of consecutive abnormal time steps is ≤3, extract the rock joint roughness coefficient data of the first 3 normal time steps to construct a linear fitting model JRC=k×step+b to predict the value of the current time step; Level 2 fault tolerance: When 3 < number of consecutive anomalies ≤ 7, the rock joint roughness coefficient value of the most recent effective time step is used; Level 3 fault tolerance: When the number of steps exceeds 7 during consecutive anomalies, the default value of 0 is output and an alarm is triggered. At the same time, the granular flow model status self-check is started. After the self-check passes, the calculation of rock joint roughness coefficient is restarted.
7. The method for real-time high-precision acquisition of rock joint roughness coefficient based on the SmoothJoint model according to claim 1, characterized in that, In the event-triggered sampling: The mechanical event triggering thresholds include: shear displacement increment ≥ 0.00001m, and the percentage of the increase in the number of bond fractures to the total number of SmoothJoint contacts ≥ 5%; The maximum sampling time step interval is 5000 time steps; The method also includes automatically saving a snapshot of the particle flow model and intermediate data of rock joint roughness coefficients every 10,000 time steps for data recovery after simulation interruption.
8. A real-time, high-precision system for acquiring rock joint roughness coefficients based on the SmoothJoint model, characterized in that, Includes the following functional modules: Particle flow model building module: Used to integrate the particle flow program application programming interface, providing a visual configuration interface for SmoothJoint model parameters, and automatically completing the generation, grouping and allocation of rock particles and SmoothJoint models; Data acquisition and verification module: Built-in dual verification logic is used to automatically filter SmoothJoint contacts and complete the verification, and output a standardized SmoothJoint contact dataset; Contact point screening and purification module: encapsulates the Grubbs test and DBSCAN clustering algorithm to achieve automated screening and purification of effective contact points; Adaptive topography reconstruction module: Built-in principal component analysis coordinate system construction, multi-algorithm interpolation and Gaussian filtering functions, used to automatically output the smoothed elevation matrix; The precise calculation module for rock joint roughness coefficient is used to automate the execution of boundary extension, second-order difference slope calculation and piecewise formula conversion, and output the rock joint roughness coefficient value. Anomaly tolerance module: Used to monitor abnormal states in real time, perform fault tolerance processing according to hierarchical strategies, generate anomaly logs and trigger alarms; Intelligent sampling and output module: used to monitor mechanical event trigger thresholds, automatically perform sampling, output structured logs and dynamic curves, and manage model snapshots and intermediate data; Multi-dimensional verification module: used to automatically complete standard contour verification, particle number convergence analysis and parameter sensitivity analysis, and generate a visual verification report.
9. The real-time high-precision acquisition system for rock joint roughness coefficient based on the SmoothJoint model according to claim 8, characterized in that, The dual verification logic in the data acquisition and verification module includes: Particle grouping and attribution verification: Verify whether two contact-related particles belong to the upper rock block group and the lower rock block group, respectively; Coordinate range verification: Verify whether the coordinates of the contact point are within the preset geometric boundary of the model, and whether the accuracy of the coordinate values is consistent with the numerical calculation accuracy of the particle flow program; In addition, for contacts that pass dual verification, the coordinates of the contact point, normal contact force, tangential contact force, contact activation state, bond fracture state, and distance from the contact point to the joint plane are collected to form a structured dataset.
10. The real-time high-precision acquisition system for rock joint roughness coefficient based on the SmoothJoint model according to claim 8, characterized in that, The hierarchical strategy of the fault tolerance module is as follows: Level 1 fault tolerance: When the number of consecutive abnormal steps is ≤3, the current value is linearly estimated based on the data of the previous 3 normal time steps; Level 2 fault tolerance: When the number of consecutive anomalies is ≤7, the most recent valid value is used; Level 3 fault tolerance: When the number of steps exceeds 7 during consecutive anomalies, output the default value of 0 and trigger an alarm, while simultaneously starting the granular flow model state self-check.