Coal pillar bearing performance evaluation method based on cyclic loading and unloading data

Abstract

The application discloses a coal column bearing performance evaluation method based on cyclic loading and unloading data, and relates to the technical field of coal column bearing performance evaluation. The specific implementation scheme is as follows: stress response sequence and strain response sequence of the coal column under the condition of cyclic loading and unloading are obtained; unloading stiffness, residual strain and hysteresis energy consumption of each cycle are extracted and respectively converted into stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics; the respective weights are determined according to the data dispersion degree in multiple cycles to obtain single-cycle damage values; the single-cycle damage values of multiple cycles are subjected to product type cumulative calculation to obtain the current cumulative damage degree; the residual bearing stress is calculated based on the current cumulative damage degree and the stress response sequence; and the safety level is output according to the residual bearing stress. Through the application, multiple types of damage response characteristics can be comprehensively reflected, the feature weights can be objectively determined, the laboratory results can be effectively transferred to the field working conditions, and the accuracy of the coal column bearing performance evaluation can be improved.

Description

Technical Field

[0001] This application relates to the field of coal pillar bearing capacity assessment technology, and more specifically, to a method for assessing coal pillar bearing capacity based on cyclic loading and unloading data. Background Technology

[0002] In deep coal mining, coal pillars such as working face pillars, section protection pillars, and side pillars of roadways are subjected to periodic mining disturbances or cyclic loading and unloading for a long time. The damage evolution of these pillars under cyclic loading has significant cumulative and irreversible characteristics. The degradation of their bearing capacity is directly related to the safety and stability of the surrounding rock and roadway system. Therefore, quantitatively evaluating the bearing capacity of coal pillars under cyclic loading and unloading conditions is a technical problem that urgently needs to be solved in engineering practice.

[0003] Currently, the performance evaluation of coal pillars under cyclic loading mainly involves the following technical approaches: One approach is based on a single monitoring quantity, such as monitoring only the stress or displacement changes of the coal pillar and issuing an early warning when the monitored quantity exceeds a preset threshold. While this method is simple to implement, a single monitoring quantity cannot comprehensively reflect the evolution of damage within the coal pillar. During cyclic loading and unloading, the damage to the coal pillar is not only manifested as a decrease in stress level, but also accompanied by changes in various response characteristics such as stiffness degradation, accumulation of permanent deformation, and energy dissipation. Relying solely on a single monitoring quantity for evaluation easily overlooks crucial information about damage evolution, leading to a significant deviation between the evaluation results and the actual state of the coal pillar.

[0004] The second method is the analogy evaluation based on laboratory test data. This type of method usually conducts cyclic loading and unloading tests on standard coal pillar samples in the laboratory to obtain the damage evolution law of the coal pillar samples, and then directly applies the test results to the evaluation of coal pillars in the field. However, there are significant differences between laboratory standard coal pillar samples and field coal pillars in terms of geometric dimensions, confining pressure conditions, moisture content, and degree of fracture development. These differences will directly affect the mechanical response and damage evolution path of the coal pillar. If the laboratory test results are transferred to field conditions without correction, the accuracy and reliability of the evaluation results are difficult to guarantee.

[0005] Another type of evaluation method is based on empirical formulas or simplified mechanical models. These methods typically describe the degradation of the coal pillar's bearing capacity as a simple function of the number of cycles, such as linear reduction models or models with fixed reduction coefficients. However, the damage evolution of a coal pillar under cyclic loading and unloading exhibits significant nonlinear characteristics. The contribution of different cycles to the damage is not the same, and the cumulative effect of damage weakens the bearing capacity of subsequent cycles. Using simplified models with linear or fixed coefficients makes it difficult to accurately characterize the true degradation behavior of a coal pillar under cyclic loading.

[0006] Furthermore, in existing technologies, the weights of each indicator are typically determined manually when constructing damage assessment indicators. Manual weighting relies on the experience and judgment of the assessors, which is highly subjective. Different assessors may provide different weight allocation schemes, leading to insufficient repeatability and objectivity of the assessment results.

[0007] Therefore, in the existing technology, under cyclic loading and unloading conditions, there is a lack of a quantitative assessment method for coal pillar damage assessment that can comprehensively reflect the characteristics of multiple damage responses, objectively determine the characteristic weights, and effectively transfer laboratory results to field working conditions. As a result, the accuracy, objectivity, and engineering practicality of the assessment results are difficult to meet actual needs. Summary of the Invention

[0008] In view of this, the purpose of this application is to provide a method for evaluating the bearing capacity of coal pillars based on cyclic loading and unloading data. This method can simultaneously extract three types of damage features from the cyclic loading and unloading response: unloading stiffness, residual strain, and hysteretic energy dissipation. These features are then objectively weighted and fused into a single-cycle damage value. Finally, the residual bearing stress and safety level are output after scenario correction, thereby achieving a quantitative evaluation of the bearing capacity of coal pillars.

[0009] According to a first aspect of this application, a method for evaluating the bearing capacity of a coal pillar based on cyclic loading and unloading data is provided, comprising: Obtain the stress response sequence and strain response sequence of the coal pillar under cyclic loading and unloading conditions; Multiple cycles are identified from the stress response sequence. In each cycle, the unloading stiffness, residual strain, and hysteresis energy are extracted based on the stress response sequence combined with the strain response sequence. Based on the predetermined critical residual strain and critical hysteresis energy dissipation, the unloading stiffness, residual strain and hysteresis energy dissipation of each cycle are converted into stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics, respectively. Based on the degree of data dispersion of stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics in multiple cycles, the weights corresponding to stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics are determined. By using the weights corresponding to stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics, the stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics of each cycle are weighted and summed to obtain the single-cycle damage value of each cycle. The current cumulative damage level is obtained by performing a product-type cumulative calculation on the single-cycle damage values ​​of multiple cycles; The residual bearing stress is calculated based on the current cumulative damage, the initial bearing stress determined by the stress response sequence, the damage sensitivity coefficient, and the scenario correction coefficient. The safety level is output based on the relationship between the residual bearing stress and the required minimum bearing stress.

[0010] In some embodiments, the stress response sequence is smoothed to obtain a smoothed stress sequence; Based on the trend of the smooth stress sequence over time, distinguish between the loading and unloading segments; The point where the smooth stress sequence changes from rising to falling and reaches a local maximum value is taken as the peak point of each cycle; The starting point of each cycle is the point before the peak point of each cycle and whose stress is a first preset proportion of the peak stress of the previous cycle. The ending point of each cycle is the point after the peak point of each cycle and whose stress is a second preset proportion of the peak stress of the current cycle. In this way, multiple complete cycles are divided from the stress response sequence.

[0011] In some embodiments, the unloading stiffness, residual strain, and hysteretic energy dissipation for each cycle are extracted based on the stress response sequence combined with the strain response sequence, including: Linear fitting is performed on the stress and strain data of the unloading section in each cycle. The slope of the fitted line is taken as the unloading stiffness of that cycle, and the intercept of the fitted line on the strain axis is taken as the residual strain of that cycle. The area enclosed between the stress-strain curve of the loading section and the stress-strain curve of the unloading section in each cycle is taken as the hysteresis energy dissipation of that cycle.

[0012] In some embodiments, the unloading stiffness, residual strain, and hysteretic energy dissipation of each cycle are converted into stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics, respectively, including: Calculate the difference between the initial cyclic unloading stiffness and the unloading stiffness of each cycle, and use the proportion of this difference to the initial cyclic unloading stiffness as the stiffness degradation characteristic of that cycle. Calculate the difference between the residual strain and the initial residual strain in each cycle, and use the proportion of this difference to the difference between the critical residual strain and the initial residual strain as the residual deformation characteristic of that cycle. Calculate the difference between the hysteresis energy consumption of each cycle and the initial hysteresis energy consumption, and use the proportion of this difference to the difference between the critical hysteresis energy consumption and the initial hysteresis energy consumption as the energy dissipation characteristic of that cycle; the initial cyclic unloading stiffness, initial residual strain and initial hysteresis energy consumption are specifically the cyclic unloading stiffness, residual strain and hysteresis energy consumption under the loading and unloading conditions of the first cycle.

[0013] In some embodiments, the critical residual strain and critical hysteresis energy dissipation are predetermined by: Cyclic loading and unloading tests were conducted on coal pillar samples from the same seam to obtain the peak stress, residual strain and hysteretic energy consumption of the coal pillar samples under each cycle. Rounds in which the peak stress drops below the preset peak stress reduction ratio of the maximum peak stress during the test and the residual strain increment exceeds the preset residual strain increment are judged as unstable rounds. The average residual strain of the instability cycle and several cycles before and after it is taken as the critical residual strain, and the average hysteresis energy consumption of the instability cycle and several cycles before and after it is taken as the critical hysteresis energy consumption.

[0014] In some embodiments, the weights corresponding to each feature are determined based on the data dispersion of the stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics of the multi-cycle stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics, including: The stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics of multiple cycles are normalized to obtain the normalized characteristic sequences corresponding to each characteristic. Based on the normalized feature sequence corresponding to each feature, calculate the information entropy corresponding to each feature; The information entropy corresponding to each feature is converted into the weight corresponding to each feature, so that the feature with greater data dispersion receives greater weight; By using the weights corresponding to each feature, the stiffness degradation feature, residual deformation feature, and energy dissipation feature of each cycle are weighted and summed to obtain the single-cycle damage value of each cycle.

[0015] In some embodiments, the single-cycle damage values ​​from multiple cycles are calculated using a product-based cumulative calculation to obtain the current cumulative damage level, including: For each cycle of damage, calculate the damage survival value. The damage survival value is specifically represented as a constant 1 minus the single-cycle damage value of that cycle. Multiply the damage survival values ​​of each cycle up to the current round to obtain the current survival product; The difference between the constant 1 and the current survival product is used as the current cumulative damage level; When obtaining the single-cycle damage value for each cycle, if the single-cycle damage value obtained by weighted summation is not greater than the preset lower limit, the single-cycle damage value for that cycle is set to the preset lower limit; if the single-cycle damage value obtained by weighted summation is greater than the preset upper limit, the single-cycle damage value for that cycle is set to the preset upper limit.

[0016] In some embodiments, the residual bearing stress is calculated based on the current cumulative damage, initial bearing stress, damage sensitivity coefficient, and scenario correction coefficient, including: Multiply the current cumulative damage level by the damage sensitivity coefficient to obtain the damage attenuation index; The attenuation factor is calculated using the natural constant as the base and the negative of the damage attenuation index as the exponent. Multiplying the initial bearing stress by the attenuation factor yields the damage-reduced bearing stress; The residual bearing stress is obtained by multiplying the damage-reduced bearing stress by the scenario correction factor.

[0017] In some embodiments, the scene correction coefficient is predetermined in the following manner: The power function of the ratio of the width-to-height ratio of the coal pillar to the standard width-to-height ratio is used as the width-to-height ratio correction factor, and the exponent of the power function is the width-to-height ratio correction parameter. The power function of the ratio of the equivalent confining pressure on the coal pillar to the standard equivalent confining pressure is used as the confining pressure correction factor, and the exponent of the power function is the confining pressure correction parameter. Substituting the difference between the moisture content of the coal pillar and the standard moisture content into the exponential decay function yields the moisture content correction factor, which includes the moisture content correction parameter. Substituting the difference between the fracture density and the standard fracture density of the coal pillar into the exponential decay function yields the fracture density correction factor, which includes the fracture density correction parameter. Multiply the aspect ratio correction factor, confining pressure correction factor, moisture content correction factor, and fracture density correction factor together to obtain the scene correction coefficient.

[0018] In some embodiments, the deviation between the field-measured bearing stress based on historical samples and the bearing stress after correction by the laboratory calibration bearing stress using the scenario correction factor is minimized, and the aspect ratio correction parameter, confining pressure correction parameter, moisture content correction parameter, and crack density correction parameter are determined by least squares fitting. The damage sensitivity coefficient was determined by regression fitting. The goal of regression fitting was to minimize the deviation between the bearing stress at which the coal column sample became unstable and the bearing stress calculated based on the damage sensitivity coefficient, the initial bearing stress, and the cumulative damage. The initial bearing stress is determined in the following way: The average value of the peak stress from multiple consecutive cycles before the start of the monitoring phase is taken as the initial bearing stress.

[0019] One embodiment of the above application has the following advantages or beneficial effects: This application uses unloading stiffness to reflect the recovery capacity of the coal pillar during the unloading process, residual strain to reflect the permanent deformation of the coal pillar that cannot be recovered after the cycle ends, and hysteretic energy dissipation to reflect the irreversible energy consumed by the coal pillar during the cycle due to microcrack propagation, particle slippage, and structural reorganization. These three features characterize the damage state of the coal pillar from three dimensions: stiffness degradation, permanent deformation, and energy dissipation, avoiding the problem of incomplete damage evolution information from a single monitoring quantity. Furthermore, the weights are determined based on the data dispersion of each feature across multiple cycles; features with greater dispersion are more sensitive to damage evolution and thus receive greater weights, eliminating the subjectivity of subjective weighting. The product-based cumulative calculation method reflects how damage in one cycle weakens the bearing capacity of subsequent cycles, better reflecting the actual degradation behavior of the coal pillar under cyclic loading. A scenario correction coefficient is used to eliminate the difference between laboratory conditions and field conditions, enabling the assessment results to be effectively transferred to the field coal pillar.

[0020] Other effects of the above-mentioned alternative methods will be described below in conjunction with specific embodiments. Attached Figure Description

[0021] The accompanying drawings are provided for a better understanding of this solution and do not constitute a limitation of this application. Wherein: Figure 1 This is a schematic diagram of a process for evaluating the bearing capacity of a coal pillar based on cyclic loading and unloading data, provided in an embodiment of this application. Figure 2 This is a line graph showing the evolution of unloading stiffness with the number of cycles provided in the embodiments of this application; Figure 3 This is a line graph showing the evolution of residual strain with the number of cycles provided in the embodiments of this application; Figure 4 This is a line graph showing the evolution of damage characteristic values ​​with the number of cycles provided in the embodiments of this application; Figure 5 This is a comparison chart of cumulative damage analysis provided in the embodiments of this application; Figure 6 This is a security level analysis attribution diagram provided in the embodiments of this application. Detailed Implementation

[0022] In the following description, specific details such as particular system architectures and techniques are set forth for illustrative purposes and not for limitation, in order to provide a thorough understanding of the embodiments of this application. However, those skilled in the art will understand that this application may also be implemented in other embodiments without these specific details. In other instances, detailed descriptions of well-known systems, apparatuses, circuits, and methods have been omitted so as not to obscure the description of this application with unnecessary detail.

[0023] See Figure 1This is a flowchart illustrating a method for evaluating the bearing capacity of a coal pillar based on cyclic loading and unloading data, provided in an embodiment of this application. Figure 1 The execution subject of the method shown can be a combination of software and / or hardware, specifically, it can be one or more of various types of terminals, hardware systems, cloud computing, etc.

[0024] Figure 1 The method for evaluating the bearing capacity of a coal pillar based on cyclic loading and unloading data includes steps S101 to S108, as detailed below: Step S101: Obtain the stress response sequence and strain response sequence of the coal pillar under cyclic loading and unloading conditions.

[0025] In this embodiment, a sensor system is deployed at the mine site for the target coal pillar that needs to be evaluated. The sensor system includes stress sensors, displacement sensors, strain sensors, and a data acquisition terminal. The stress sensors are used to acquire stress changes inside or on the surface of the coal pillar, the displacement sensors are used to acquire the convergence of the coal pillar surface or the proximity of the roof and floor, the strain sensors are used to acquire local deformation of the coal pillar, and the data acquisition terminal is used to synchronously record the output signals of each sensor.

[0026] After sensor deployment, zero-point calibration and no-load correction are performed to ensure that the initial output value of all sensors is zero or close to zero before formal monitoring, thereby reducing drift error. The preferred sensor placement and quantity are as follows: stress sensors are preferentially placed at the contact zone between the coal pillar and the roof and floor, and on the middle surface of the coal pillar's sidewalls, spaced evenly along the coal pillar's direction, typically 2 to 4 sensors. The spacing between adjacent sensors is 5 to 10 meters. The roof and floor contact zone is the main load transfer interface between the coal pillar and the surrounding rock, where stress changes are most sensitive to cyclic pressure. Measurement points in the middle of the sidewalls are used to reflect the overall load level of the coal pillar. Averaging multiple measurement points reduces the impact of coal heterogeneity and local stress concentration on single-point measurement results. Displacement sensors are deployed on both sides of the coal pillar along the coal seam normal, with one set of convergence measurement lines, typically 2 to 4 sets, to measure the roof and floor convergence and the sidewall surface convergence, respectively, to capture the coal pillar deformation response caused by cyclic loading and unloading. Strain sensors are arranged in pairs near the stress sensor along the direction perpendicular to the principal stress. The number of pairs is generally 4 to 8. The pair arrangement can distinguish between shallow deformation on the surface of the coal pillar and overall compression deformation.

[0027] Stress gauge sensors are embedded inside the coal pillar at a depth of no less than 0.5 meters, avoiding surface fracture zones and areas with obvious crack development. This ensures that the measuring points reflect the effective stress inside the coal pillar rather than local fluctuations caused by surface loosening. Through this deployment scheme, data from each measuring point are cross-checked spatially and recorded synchronously in time, improving the representativeness of the monitoring data for the overall damage state of the coal pillar and reducing the interference of local anomalies on subsequent feature extraction and damage assessment.

[0028] In the laboratory, standard cylindrical coal pillar samples are prepared from raw coal of the same stratum, type, and structural section as the target coal pillar. The diameter-to-height ratio of the coal pillar sample meets the requirements of conventional rock mechanics testing, such as a diameter of 50 mm and a height of 100 mm, or uses a similar size according to the specifications of the testing equipment, but all calibration coal pillar samples maintain consistent dimensions. After sample preparation, the coal pillar samples are numbered, and the moisture content, end face flatness, and fracture distribution are uniformly controlled. Subsequently, cyclic loading and unloading tests are carried out on a triaxial or uniaxial cyclic loading testing machine. The loading rate, unloading rate, and number of cycles are designed to simulate the loading characteristics of cyclic pressure in the field as closely as possible, ensuring the transferability of the calibration results.

[0029] Let the acquired original stress sequence and strain sequence be respectively and ,in Number the sampling points. The unit is megapascal (MPa). Dimensionless, sampling frequency is fixed. The unit is Hertz. The value of the sampling frequency is determined based on the capability of the field data acquisition terminal and the rate of change of the cycle. It is generally taken as 50 Hertz to 200 Hertz to ensure that a sufficient number of data points can be collected in each cycle for subsequent fitting and integration calculations.

[0030] Identifying multiple cycles from a stress response sequence includes: smoothing the stress response sequence to obtain a smoothed stress sequence; distinguishing between loading and unloading segments based on the time-varying trend of the smoothed stress sequence; taking the point where the smoothed stress sequence changes from rising to falling and reaches a local maximum as the peak point of each cycle; taking the point before the peak point of each cycle and whose stress is a first preset proportion of the peak stress of the previous cycle as the starting point of the current cycle; and taking the point after the peak point of each cycle and whose stress is a second preset proportion of the peak stress of the current cycle as the ending point of the current cycle. In this way, multiple complete cycles are divided from the stress response sequence.

[0031] In this embodiment, to eliminate high-frequency noise, mechanical vibration, and occasional sensor pulse interference, the original stress sequence and strain sequence are smoothed using a five-point moving average. The smoothed stress sequence and strain sequence are represented as follows:

[0032]

[0033] in, The stress value is the smoothed value, in megapascals (MPa). The strain value is the smoothed strain value, dimensionless. Number the sampling points. and The first Stress and strain data at each sampling point This represents the window offset, with a window size of 5, meaning it takes the average of the current sampling point and the two sampling points before and after it. The five-point moving average was chosen because it is simple to implement, computationally inexpensive, suitable for online processing, and effectively preserves the main morphological characteristics of the cyclic curve.

[0034] Setting the window size to 5 is a balance between smoothing effect and shape preservation: if the window is too small, the smoothing effect will be insufficient, and if the window is too large, the smoothing will be over-smoothed, which will weaken the peak and inflection point information of the cycle curve.

[0035] After obtaining the smoothed stress sequence, the loading and unloading segments are distinguished based on the trend of smoothed stress variation over time. Since the collected data consists of equally spaced discrete sampling points, the continuous derivative cannot be directly calculated. Therefore, in actual calculations, the first-order difference between adjacent sampling points is used. Approximate number 1 Stress change rate at each sampling point.

[0036] The sampling interval is fixed, and the sign of the difference value is consistent with the sign of the derivative. A difference value greater than 0 indicates that the stress is in the rising stage, and a difference value less than 0 indicates that the stress is in the falling stage. To avoid misjudging the stage transition due to residual noise causing the sign to jitter at individual sampling points, the current loading or unloading stage is only confirmed when the signs of the difference values ​​at three or more consecutive sampling points are consistent.

[0037] Based on this, the loading and unloading segments are distinguished according to the derivative criterion: when the derivative of the smooth stress with respect to time satisfies... When this stage is determined to be the loading stage, and the conditions are met... At this point, the stage is determined to be the unloading stage. The point where the stress derivative changes from positive to negative and the stress reaches a local maximum is determined to be the [number missing]th [stage missing]. The peak point of the cycle, the peak stress is denoted as... The unit is megapascal (MPa), and the peak strain is denoted as . , dimensionless.

[0038] To avoid interference from noise in the low-stress zone on the identification of the cycle boundary, this embodiment defines the start and end points of each cycle as points corresponding to 5% of the peak stress. The specific determination formula is as follows:

[0039]

[0040] in, For the first The starting time of the cycle, in seconds. For the first The end time of the cycle, in seconds. For the first Peak stress in a cycle, measured in megapascals (MPA). For the first Peak stress in a cycle, measured in megapascals (MPA). This is the index for loading and unloading loops.

[0041] The reason for choosing 5% as the start and end thresholds is as follows: If the threshold is set to 0%, requiring the stress to completely return to zero before determining the end of the cycle, noise fluctuations in the low-stress zone will cause the cycle boundary to be repeatedly triggered, making it impossible to stably separate each cycle. If the threshold is set too high, such as 20%, damage information in the low-stress segment below 20% of the peak stress will be excluded from the cycle, resulting in the omission of this part of the damage. The 5% threshold strikes a balance between avoiding noise interference in the low-stress zone and retaining damage information in the low-stress segment, thereby stably separating each cycle.

[0042] Step S102: Identify multiple cycles from the stress response sequence. In the multiple cycles, extract the unloading stiffness, residual strain and hysteresis energy of each cycle based on the stress response sequence combined with the strain response sequence.

[0043] Extract the unloading stiffness, residual strain, and hysteretic energy dissipation for each cycle, including: Linear fitting is performed on the stress and strain data of the unloading section in each cycle. The slope of the fitted line is taken as the unloading stiffness of that cycle, and the intercept of the fitted line on the strain axis is taken as the residual strain of that cycle. The area enclosed between the stress-strain curve of the loading section and the stress-strain curve of the unloading section in each cycle is taken as the hysteresis energy dissipation of that cycle.

[0044] In this embodiment, for the first The unloading section data from the cycle was linearly fitted using the least squares method. The stress-strain relationship of the unloading section is expressed as:

[0045] in, The stress in the unloading section is expressed in megapascals (MPa). The strain of the unloading segment is dimensionless. For the first The unloading stiffness of the wheel cycle, measured in megapascals (MPa). The residual strain corresponding to this cycle is dimensionless.

[0046] The least squares method was used to perform linear fitting on the stress and strain data of the unloading section in each cycle to obtain the fitted straight line.

[0047] The formula for calculating unloading stiffness is:

[0048] in, For the first The unloading stiffness of the wheel cycle, measured in megapascals (MPa). and The first Sampling numbers of the start and end points of the unloading section in the cycle. For the first The strain values ​​at each sampling point are dimensionless. For the first The stress value at each sampling point is expressed in megapascals (MPa). For the first The average strain of the unloading section during the wheel cycle is dimensionless. For the first The average stress during the cyclic unloading section is expressed in megapascals (MPa). The least squares method was chosen because it provides stable fitting results, requires minimal computation, and is suitable for online processing. For details, please refer to [link to relevant documentation]. Figure 2 ,according to Figure 2 It can be seen that the unloading stiffness continuously decreases from the initial value of about 4510MPa to about 1500MPa, showing a clear monotonically decreasing trend. The degradation process of the coal pillar structure stiffness was quantified, indicating that as the number of cycles increases, the number of internal cracks in the coal pillar increases, the structural integrity decreases, and the elastic recovery ability continues to weaken. Moreover, the decrease in unloading stiffness is not linear, but shows a nonlinear decay trend of fast at first and slow at the end.

[0049] The basis for the approximate linearity of the stress-strain relationship in the unloading section is that the coal body primarily undergoes elastic recovery during unloading, with minimal additional plastic deformation. The stress-strain relationship in the unloading section approximates the behavior of an elastic linear line, thus a single straight line can be used to fit the unloading section data as a whole. The slope of the fitted straight line represents the unloading stiffness of that cycle, reflecting the elastic recovery capacity of the coal pillar under its current damage state. As the number of cycles increases, internal cracks develop and expand, leading to a continuous decrease in unloading stiffness. This characteristic can directly quantify the degradation process and rate of structural stiffness.

[0050] The further reason for using overall linear fitting instead of just taking the slope of the line connecting the start and end points of the unloading segment is that least squares fitting utilizes all the sampling point information of the unloading segment, which has an averaging effect on the random errors of individual measurement points, and the slope obtained by fitting is more stable than the slope of the line connecting the endpoints.

[0051] The residual strain is obtained from the intercept of the linear fit of the unloading section, and the calculation formula is as follows:

[0052] in, For the first The residual strain of the wheel cycle is dimensionless. The strain is the average value of the unloading section, dimensionless. This represents the average stress in the unloading section, in megapascals (MPA). For the first The unloading stiffness of the cycle is expressed in megapascals (MPa). Residual strain is the amount of permanent deformation that a coal pillar cannot recover when the stress returns to zero. The residual strain is not directly taken as the strain value at the unloading endpoint, but rather extrapolated from the intercept of the fitted straight line on the strain axis during the unloading segment. This is because, in field monitoring, each unloading cycle often begins before the stress has completely returned to zero, leaving some unreleased elastic deformation at the unloading endpoint. Directly using the strain at the unloading endpoint would mix this elastic deformation into the residual strain, leading to an overestimation of permanent deformation. By extrapolating the strain corresponding to zero stress using the fitted straight line, unrecoverable permanent deformation can be accurately separated from recoverable elastic deformation, ensuring that the residual strain only reflects the degree of plastic failure. For details, please refer to [link to relevant documentation]. Figure 3 ,according to Figure 3 Analysis shows that the residual strain increased from an initial value of about 0.506‰ to about 6.996‰, reflecting the irreversible accumulation process of permanent deformation of the coal pillar; the larger the residual strain, the worse the fracture closure and recovery ability, and the more serious the accumulation of permanent deformation.

[0053] The hysteresis energy dissipation for each cycle is obtained as follows: Discretely integrate the stress-strain curves of the loading and unloading sections for each cycle to obtain the area under the curve of the loading section and the area under the curve of the unloading section. The difference between the area under the curve of the loading section and the area under the curve of the unloading section is taken as the hysteresis energy dissipation for that cycle.

[0054] The continuous integral definition of hysteresis energy dissipation is:

[0055] in, For the first The hysteresis energy consumption of the wheel cycle, measured in megajoules per cubic meter. For the first Stress function of the wheel-loaded segment, in megapascals (MPa). For the first Stress function of the wheel unloading section, in megapascals (MPa). The strain at the starting point of this cycle is dimensionless. For peak strain, dimensionless. The residual strain is dimensionless.

[0056] Since the actual collected data consists of discrete sampling points, the trapezoidal formula is used to perform discrete integration on the loading and unloading segments respectively. The calculation formula is as follows:

[0057] in, For the first The hysteresis energy consumption of the wheel cycle, measured in megajoules per cubic meter. and These are the th elements within the current integration segment. The and the first The stress value at each sampling point is expressed in megapascals (MPa). and These are the th elements within the current integration segment. The and the first The strain values ​​at each sampling point are dimensionless. and These are the sampling numbers of the start and end points of the current integration segment. When integrating the loading segment, the start and end points of the loading segment are used; when integrating the unloading segment, the start and end points of the unloading segment are used. The trapezoidal integral results of the loading segment and the unloading segment are calculated using the above formula. The integration result of the loading segment is then subtracted from the integration result of the unloading segment to obtain the [number of segments]. The hysteresis energy consumption of the cycle is consistent with the definition of continuous integral. The trapezoidal rule is chosen because the data consists of discrete sampling points, and the trapezoidal rule achieves a good balance between computational accuracy and computational complexity. It is also simple to implement and suitable for online processing in the field.

[0058] The three types of characteristics each have a clear physical meaning. Unloading stiffness reflects the recovery ability of the coal pillar during the unloading process. If the unloading stiffness continues to decrease, it indicates an increase in internal cracks and a significant degradation in structural stiffness. Residual strain reflects the degree of plastic failure of the coal pillar. The larger the residual strain, the worse the crack closure recovery ability and the more severe the accumulation of permanent deformation. The mechanical meaning of hysteresis energy dissipation is that the area under the curve of the loading segment corresponds to the total work density done by the external force on the coal pillar, and the area under the curve of the unloading segment corresponds to the elastic energy density released by the coal pillar during the unloading process. The difference between the two is the energy that has not been recovered in this cycle. This part of the energy is mainly consumed by irreversible processes such as microcrack propagation, particle slip, and structural reorganization. Therefore, the larger the hysteresis energy dissipation, the more severe the damage evolution in this cycle. The three types of characteristics characterize the damage state of the coal pillar from three dimensions: stiffness degradation, permanent deformation, and energy dissipation, avoiding the problem that a single monitoring quantity cannot fully reflect the damage evolution information.

[0059] Step S103: Based on the predetermined critical residual strain and critical hysteresis energy dissipation, the unloading stiffness, residual strain and hysteresis energy dissipation of each cycle are converted into stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics, respectively.

[0060] The unloading stiffness, residual strain, and hysteretic energy dissipation of each cycle are converted into stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics, respectively, including: Calculate the difference between the initial cyclic unloading stiffness and the unloading stiffness of each cycle, and use the proportion of this difference to the initial cyclic unloading stiffness as the stiffness degradation characteristic of that cycle. Calculate the difference between the residual strain and the initial residual strain in each cycle, and use the proportion of this difference to the difference between the critical residual strain and the initial residual strain as the residual deformation characteristic of that cycle. Calculate the difference between the hysteresis energy consumption of each cycle and the initial hysteresis energy consumption, and use the proportion of this difference to the difference between the critical hysteresis energy consumption and the initial hysteresis energy consumption as the energy dissipation characteristic of that cycle; the initial cyclic unloading stiffness, initial residual strain and initial hysteresis energy consumption are specifically the cyclic unloading stiffness, residual strain and hysteresis energy consumption under the loading and unloading conditions of the first cycle.

[0061] In this embodiment, the calculation formulas for the three types of dimensionless damage characteristics are as follows:

[0062]

[0063]

[0064] in, For the first The stiffness degradation characteristics of a wheel cycle are dimensionless. For the first The residual deformation characteristics of the wheel cycle are dimensionless. For the first The energy dissipation characteristics of a cycle are dimensionless. For the first The unloading stiffness of the wheel cycle, measured in megapascals (MPa). The initial cyclic unloading stiffness is expressed in megapascals (MPA). For the first The residual strain of the wheel cycle is dimensionless. The initial residual strain is dimensionless. The critical residual strain is dimensionless. For the first The hysteresis energy consumption of the wheel cycle, measured in megajoules per cubic meter. This represents the initial hysteresis energy consumption, expressed in megajoules per cubic meter. This refers to critical hysteresis energy consumption, measured in megajoules per cubic meter. Please refer to [link / reference] for details. Figure 4 ,according to Figure 4 It can be seen that all three characteristics show a monotonically increasing trend, reflecting the irreversible accumulation process of damage. When the coal pillar is in its initial state, all three characteristics are 0, and when it reaches the critical state of instability, the residual deformation characteristic and energy dissipation characteristic are close to 1.

[0065] The purpose of dimensionless transformation is to eliminate the dimensional differences among unloading stiffness, residual strain, and hysteretic energy dissipation, making the three types of features numerically comparable, thus enabling weighted fusion in subsequent steps. The reason for selecting a relative comparison benchmark for each feature is as follows: the stiffness degradation feature uses the initial cyclic unloading stiffness as the benchmark, reflecting the relative attenuation of the current unloading stiffness relative to the initial state; the residual deformation feature uses the difference between the critical residual strain and the initial residual strain as the benchmark, reflecting the degree to which the current residual strain approaches the instability state; and the energy dissipation feature uses the difference between the critical hysteretic energy dissipation and the initial hysteretic energy dissipation as the benchmark, reflecting the degree to which the current hysteretic energy dissipation approaches the instability state. According to the above structure, all three characteristics are 0 when the coal pillar is in its initial state. When the instability critical state determined by the laboratory is reached, the residual deformation characteristic and energy dissipation characteristic are close to 1. Therefore, the three characteristics share a unified 0 to 1 scale. The single-cycle damage value obtained by subsequent weighted summation naturally falls within the bounded interval. The weight calculation and damage value fusion are carried out under the same metric scale, avoiding the artificial amplification of a certain characteristic in the fusion result due to the large difference in the absolute magnitude of each characteristic. The ratio form construction also makes the gain differences between similar sensors cancel each other out in the numerator and denominator, thus making it insensitive to the differences in the calibration accuracy of the monitoring instruments.

[0066] The critical residual strain and critical hysteresis energy consumption are predetermined by the following method: cyclic loading and unloading tests are carried out on coal pillar samples from the same seam to obtain the peak stress, residual strain and hysteresis energy consumption of the coal pillar samples in each cycle. The cycle in which the peak stress drops below the preset peak value reduction ratio of the maximum peak stress during the test and the residual strain increment exceeds the preset residual strain increment is determined as the unstable cycle. The average residual strain of the unstable cycle and several cycles before and after it is taken as the critical residual strain, and the average hysteresis energy consumption of the unstable cycle and several cycles before and after it is taken as the critical hysteresis energy consumption.

[0067] In this embodiment, during the laboratory cyclic loading and unloading test, loading was continued until the coal pillar sample experienced a significant peak and then became unstable. The cycle in which the peak bearing capacity of the coal pillar sample decreased significantly and was accompanied by a sudden increase in residual strain was defined as the instability cycle. The following conditions are used to determine instability:

[0068] and

[0069] in, The peak stress at the instability stage is expressed in megapascals (MPa). The maximum peak stress observed during the test is expressed in megapascals (MPa). The residual strain of the instability cycle is dimensionless. The residual strain of the previous cycle before the instability cycle is dimensionless.

[0070] The reason for setting the peak stress reduction ratio to 90% is that when the peak stress drops below 90% of the maximum peak stress, it indicates that the coal pillar sample has experienced a significant decrease in bearing capacity. This criterion is consistent with the judgment standard for coal pillar instability in engineering practice. The reason for setting the residual strain increment to 0.002 is that a residual strain increment of 0.002 corresponds to the identifiable threshold for a sudden increase in residual strain, which can effectively distinguish between residual strain fluctuations during normal cycles and sudden increases in residual strain before instability. Instability judgment adopts a joint judgment method that requires both conditions to be met simultaneously. If only the peak stress reduction is used as the criterion, the coal pillar sample may experience fluctuations such as a single-cycle peak drop followed by partial recovery due to unloading disturbances or local damage, which can easily be misjudged as instability. If only the sudden increase in residual strain is used as the criterion, the residual strain in the pre-peak stage also accumulates naturally with the increase in the number of cycles, making it difficult to clearly distinguish from the sudden increase before instability. When both conditions are met simultaneously, it means that the decrease in peak stress and the sudden increase in residual strain occur in the same cycle. That is, the coal pillar sample loses its bearing capacity at the same time and the irreversible deformation intensifies. Based on this, the unstable cycle can reduce the probability of misjudgment and omission, and make the calibration results of the critical threshold more reliable.

[0071] After determining the instability cycle, the average residual strain of the instability cycle and the two cycles before and after it (a total of three cycles) is taken as the critical residual strain, and the average hysteresis energy dissipation of the instability cycle and the two cycles before and after it (a total of three cycles) is taken as the critical hysteresis energy dissipation. The calculation formula is as follows:

[0072]

[0073] in, The critical residual strain is dimensionless. This refers to critical hysteresis energy consumption, measured in megajoules per cubic meter. , , The residual strains are given by the first two cycles, the previous cycle, and the cycle after instability, respectively. They are dimensionless. , , These are the hysteresis energy consumption of the two rounds before the instability, the previous round, and the instability round, respectively, in megajoules per cubic meter.

[0074] The reason for taking the average of three rounds is to avoid a single outlier causing a shift in the threshold, making the critical parameter more stable and reliable. Critical cumulative damage degree The damage was determined using laboratory cyclic test data following this procedure. First, the data from each cycle of the laboratory coal pillar sample underwent the same evaluation process as that used for field coal pillars. This involved linearly fitting the unloading section of each cycle to obtain the unloading stiffness and residual strain, integrating to obtain the hysteresis energy dissipation, constructing dimensionless damage characteristics, weighting them using the entropy weighting method, and then summing them to obtain the single-cycle damage value. Finally, a product-type accumulation was performed to obtain the cumulative damage sequence of the coal pillar sample across each cycle. Then, the instability cycle was determined. The cumulative damage value at the point of instability is used as the critical cumulative damage value. The rationale for using the value at the instability cycle is that this value corresponds to the actual damage level of the coal pillar sample at the moment of transition from cyclic damage accumulation to macroscopic instability. Exceeding this level, the coal pillar sample enters the rapid post-peak instability stage; therefore, this value can be directly used as the upper limit of damage that the coal pillar can reach. To reduce the influence of the dispersion of a single coal pillar sample, the average critical cumulative damage value of multiple coal pillar samples from the same seam is further used as the final critical cumulative damage value for predicting the remaining cycle life.

[0075] Step S104: Based on the data dispersion of stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics in multiple cycles, determine the weights corresponding to stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics.

[0076] The stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics of multiple cycles are normalized to obtain the normalized feature sequence corresponding to each feature. Based on the normalized feature sequence corresponding to each feature, the information entropy corresponding to each feature is calculated. The information entropy corresponding to each feature is converted into the weight corresponding to each feature, so that the feature with greater data dispersion receives a larger weight. The stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics of each cycle are weighted and summed using the weights corresponding to each feature to obtain the single-cycle damage value of each cycle.

[0077] In this embodiment, the three types of features are first normalized, and the normalization formula is as follows:

[0078] in, For the first Cycle 1 Normalized values ​​of each feature, dimensionless. For the first Cycle 1 The value of each feature is dimensionless. For the first The minimum value of a feature across all cycles, dimensionless. For the first The maximum value of a feature across all cycles, dimensionless. We choose 1, 2, and 3, which correspond to stiffness degradation, residual deformation, and energy dissipation characteristics, respectively. The reason for choosing linear normalization is that linear mapping ensures the normalized value is non-negative, meeting the requirement of non-negativity for subsequent probability distribution construction. It also does not change the relative proportions within the feature sequence, does not introduce additional distribution assumptions, has low computational cost, and is suitable for online processing. If a certain feature takes a constant value in all loops, i.e. If the normalization denominator is zero, the constant feature remains unchanged throughout the cycle and lacks the ability to distinguish damage evolution. In this case, the weight corresponding to that feature is directly set to 0 and no longer participates in the calculation of information entropy. The remaining features are still weighted according to the information entropy method. After this processing, the normalization process is defined for any feature sequence.

[0079] Then, a probability distribution is constructed, and the probability formula is:

[0080] in, For the first Cycle 1 The probability values ​​of each feature are dimensionless. For the first Cycle 1 Normalized values ​​of each feature, dimensionless. It is a very small positive number, used to prevent the logarithm from being undefined when some terms are zero. This represents the total number of cycles, i.e., the number of monitoring cycles that have been completed up to the present. Pick The reason is that this value is much smaller than the typical order of magnitude of the normalized eigenvalues, providing only a small positive shift to the zero-value probability term to eliminate the obstacle of logarithmic calculation, and basically not changing the relative proportions between the non-zero terms. Taking excessively large values ​​will artificially inflate the probability of zero terms and compress the differences between the probabilities of non-zero terms, thus distorting the information entropy estimation. The value should be minimized as much as possible while ensuring the stability of numerical calculations.

[0081] Next, calculate the first... The information entropy of a feature is given by the formula:

[0082] in, For the first The information entropy of a feature is dimensionless. The total number of cycles, For the first Cycle 1 The probability values ​​of each feature are dimensionless. This is a normalization factor used to constrain the range of information entropy between 0 and 1. The smaller the information entropy, the greater the dispersion of the feature and the more sensitive it is to damage evolution.

[0083] Finally, the complementary amount of information entropy is converted into weights, and the weight formula is:

[0084] in, For the first The weights of each feature are dimensionless. For the first The information entropy of each feature is dimensionless. The weights satisfy... Constraints.

[0085] The reason for choosing the entropy weight method is that it allocates weights entirely based on data variability, objectively reflecting which features are more discriminative and avoiding subjective errors caused by manual weighting. If equal weights are assigned manually, for example, by setting the weight of each of the three features to 1 / 3, then when the dispersion of a certain feature is significantly greater than that of the other features, the discriminative role of that feature cannot be highlighted, leading to a decrease in the sensitivity of damage assessment.

[0086] Step S105: Using the weights corresponding to the stiffness degradation feature, the residual deformation feature, and the energy dissipation feature, the stiffness degradation feature, residual deformation feature, and energy dissipation feature of each cycle are weighted and summed to obtain the single-cycle damage value of each cycle.

[0087] By using the weights corresponding to each feature, the stiffness degradation feature, residual deformation feature, and energy dissipation feature of each cycle are weighted and summed to obtain the single-cycle damage value for each cycle. The calculation formula is as follows:

[0088] in, For the first Single-cycle damage value of the cycle, dimensionless. , , These represent the weights corresponding to stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics, respectively, and are dimensionless. , , The first The stiffness degradation, residual deformation, and energy dissipation characteristics of the wheel cycle are dimensionless. The single-cycle damage value ranges from 0 to 1, with larger values ​​indicating more severe damage caused by the wheel cycle.

[0089] Step S106: Perform a product-type cumulative calculation on the single-cycle damage values ​​of multiple cycles to obtain the current cumulative damage degree.

[0090] For each cycle, calculate the damage survival value based on the single-cycle damage value. The damage survival value is 1 minus the single-cycle damage value of that cycle. Multiply the damage survival values ​​of all cycles up to the current cycle to obtain the current survival product. Subtract the current survival product from 1 as the current cumulative damage level.

[0091] The current cumulative damage level is obtained by performing a product-type cumulative calculation on the single-cycle damage values ​​from multiple cycles, including: For each cycle of damage, calculate the damage survival value. The damage survival value is specifically represented as a constant 1 minus the single-cycle damage value of that cycle. Multiply the damage survival values ​​of each cycle up to the current round to obtain the current survival product; The difference between the constant 1 and the current survival product is taken as the current cumulative damage level; in this embodiment, up to the [number]th ... The formula for calculating the cumulative damage degree of a cycle is:

[0092] in, As of the date The cumulative damage of the cycle, dimensionless. For the first Single-cycle damage value of the cycle, dimensionless. As of the date The product of the survival values ​​of each round of damage is the current survival product.

[0093] The reason for choosing product-based accumulation over nonlinear superposition is that damage in one cycle weakens the load-bearing foundation of subsequent cycles. Subsequent damage accumulates on top of the reduced load-bearing capacity caused by previous damage. Product-based accumulation better reflects the nonlinear characteristics of damage evolution under cyclic loading. If linear superposition were used, i.e. When the damage value is large, the cumulative damage degree may exceed 1, losing its physical meaning and failing to reflect the nonlinear acceleration effect of damage accumulation. Please refer to [link / reference needed] for details. Figure 5 , Figure 5 The product-based accumulation method (solid blue line) ensures that the cumulative damage value remains within the [0,1] interval, possessing a clear physical meaning: D=0 represents no damage, and D=1 represents complete damage. In contrast, the linear superposition method (dashed orange line) may result in a cumulative damage value exceeding 1 after multiple cycles, losing its physical meaning. For example, when the single-cycle damage value is 0.4 for three consecutive cycles, linear superposition yields 1.2, exceeding a reasonable range; product-based accumulation yields 0.784, remaining within a reasonable range.

[0094] When obtaining the single-cycle damage value for each cycle, if the single-cycle damage value obtained by weighted summation is not greater than the preset lower limit, the single-cycle damage value for that cycle is set to the preset lower limit; if the single-cycle damage value obtained by weighted summation is greater than the preset upper limit, the single-cycle damage value for that cycle is set to the preset upper limit.

[0095] In this embodiment, the preset lower limit is set to 0, and the preset upper limit is set to 1. When the weighted summation result is less than 0 due to numerical error, the single-cycle damage value of that cycle is set to 0; when the weighted summation result is greater than 1 due to numerical error, the single-cycle damage value of that cycle is set to 1. The reason for the truncation is to prevent numerical error from causing the damage value to exceed the range of physical meaning, ensuring that the single-cycle damage value always remains within a reasonable range of 0 to 1.

[0096] Step S107: Calculate the residual bearing stress based on the current cumulative damage degree, the initial bearing stress determined by the stress response sequence, the damage sensitivity coefficient, and the scenario correction coefficient.

[0097] The scenario correction factor is determined by four types of differences between field and laboratory conditions: aspect ratio of the coal pillar, equivalent confining pressure, moisture content, and fracture density. Among them, aspect ratio and equivalent confining pressure have an amplification effect on bearing capacity, while moisture content and fracture density have an attenuation effect on bearing capacity.

[0098] In this embodiment, the complete formula for the scene correction coefficient is:

[0099] in, The scene correction coefficient is dimensionless. This refers to the width of the coal pillar on site, in meters. This refers to the height of the coal pillar at the site, in meters. To calibrate the standard aspect ratio of coal pillar samples in the laboratory, dimensionless. The equivalent confining pressure on site is expressed in megapascals (MPa). The standard equivalent confining pressure is the confining pressure used in laboratory calibration tests, measured in megapascals (MPa). The moisture content of the coal seam at the site is expressed as a percentage. Standard moisture content, expressed as a percentage. This represents the density of fractures on site, expressed in units of fractures per meter. Standard fracture density, in units of fractures per meter. This is a dimensionless parameter for aspect ratio correction. The confining pressure correction parameter is dimensionless. This is a moisture content correction parameter, and its value is calculated for every 1 percentage point difference in moisture content. This is the correction parameter for fracture density, expressed in meters per fracture.

[0100] The aspect ratio of the coal pillar to the standard aspect ratio is used as a power function as the aspect ratio correction factor, and the exponent of this power function is the aspect ratio correction parameter. The confining pressure of the coal pillar to the standard equivalent confining pressure is used as a power function as the confining pressure correction factor, and the exponent of this power function is the confining pressure correction parameter. The difference between the moisture content of the coal pillar and the standard moisture content is substituted into the exponential decay function to obtain the moisture content correction factor, which includes the moisture content correction parameter. The difference between the fracture density of the coal pillar and the standard fracture density is substituted into the exponential decay function to obtain the fracture density correction factor, which includes the fracture density correction parameter. The aspect ratio correction factor, confining pressure correction factor, moisture content correction factor, and fracture density correction factor are multiplied together to obtain the scene correction coefficient.

[0101] The reason for choosing the combination of power functions and exponential decay functions is that aspect ratio and confining pressure have a scale amplification effect on bearing capacity. Using power functions can reflect this continuously amplified effect with scale changes. Moisture content and fracture density mainly exhibit a decay effect on bearing capacity. Using exponential decay functions can reflect the accelerated decline in bearing capacity as moisture content or fracture density increases. Power functions and exponential functions respectively adapt to these two different influence mechanisms. The aspect ratio correction factor and confining pressure correction factor are constructed using the ratio of field values ​​to standard values ​​because the ratio can eliminate the influence of absolute size magnitude, retaining only the degree of deviation of field conditions from laboratory conditions. The moisture content correction factor and fracture density correction factor are constructed by substituting the difference between field values ​​and standard values ​​into an exponential decay form because the weakening of bearing capacity by moisture content and fracture density is reflected in the extent exceeding the standard state; the greater the deviation, the more obvious the decay. All four correction factors are 1, which corresponds to the situation where the field coal pillar and the laboratory calibration conditions are completely consistent. In this case, the scenario correction coefficient automatically degenerates to 1, indicating that no correction is needed. Therefore, the construction takes 1 as the no-correction benchmark. When the field conditions are enhanced relative to the laboratory conditions, the correction coefficient is greater than 1, and when they are weakened, the correction coefficient is less than 1. The correction direction is intuitive and clear.

[0102] The deviation between the field-measured bearing stress based on historical samples and the bearing stress after correction by the laboratory calibration bearing stress using the scenario correction factor is minimized. The aspect ratio correction parameter, confining pressure correction parameter, moisture content correction parameter, and crack density correction parameter are determined by least squares fitting.

[0103] In this embodiment, the correction parameters are determined by least-squares fitting of historical samples from the same mining area, and the objective function for fitting is:

[0104] in, For the first The measured stable bearing stress of a single field sample, in megapascals (MPa). To correspond to the laboratory calibration bearing stress, the unit is megapascals (MPa). The correction coefficients are for corresponding scenarios and are dimensionless. To determine the sample size, , , , These are the correction parameters to be fitted. The fitting objective is to minimize the deviation between the field-measured bearing stress and the laboratory-calibrated bearing stress after scenario correction. The historical samples are constructed as follows: coal pillars from the same mining area that have completed their service cycle or have complete monitoring records during mining are selected as samples. The aspect ratio, equivalent confining pressure, moisture content, and fracture density of each sample are recorded. The bearing stress measured in the stable bearing stage of the sample is used as the field-measured bearing stress, and the corresponding laboratory-calibrated bearing stress is determined according to the laboratory cyclic test results of coal pillar samples from the same seam. The initial values ​​for iteration in the fitting calculation can be calculated as follows: , , , The iteration convergence condition is set to the parameter change in two consecutive iterations being less than 1 / 2. Or the relative decrease in the objective function value is less than Based on mining engineering experience, aspect ratio correction parameters. The typical value range is 0.2 to 0.5, for the confining pressure correction parameter. The typical value range is 0.3 to 0.7, for the moisture content correction parameter. The typical value range is 0.1 to 0.3, for the fracture density correction parameter. The typical value range is 0.05 to 0.2. When the aspect ratio correction parameter and confining pressure correction parameter are less than 1, the bearing capacity increases with increasing aspect ratio or confining pressure, but the rate of increase gradually slows down. The larger the moisture content correction parameter and fracture density correction parameter, the more significant the marginal weakening of the bearing capacity beyond the standard state. (Number of calibration samples) The number of samples should not be less than three. If the number of samples is too small, the fitted correction parameters will be unstable. The samples should cover different coal pillar widths, burial depths, and water-bearing states to enhance the representativeness of the fitted parameters under different working conditions. After the correction parameters are calibrated, they will be fixed for subsequent on-site coal pillar assessments in the mining area. If the geological conditions in the mining area change significantly, the calibration should be carried out again.

[0105] If on-site equivalent confining pressure is not feasible with actual measurement, it can be estimated based on burial depth and lateral pressure coefficient, i.e. ,in This represents the average unit weight of the overlying rock strata, expressed in kilonewtons per cubic meter. Burial depth, in meters. The lateral pressure coefficient is dimensionless. The unit conversion factor is used to make the estimated result in megapascals. When actual measurement conditions are available, the equivalent confining pressure in the field is preferably obtained by measuring using the hydraulic fracturing method or the casing stress relief method. The horizontal stress component perpendicular to the coal pillar axis at the measuring point is taken as the equivalent confining pressure in the field. When actual measurement conditions are not available, it is determined according to the above estimation formula. The basis of the estimation formula is that the vertical component of the initial stress field of the deep coal pillar is generally close to the self-weight stress field, and the vertical stress is approximately... The lateral constraint of the coal pillar mainly comes from the horizontal ground stress, which is calculated using the lateral pressure coefficient. Therefore, the converted horizontal in-situ stress is used as the equivalent confining pressure. The lateral pressure coefficient is determined based on measured in-situ stress data of the mining area. When measured data is lacking, a value range of 0.8 to 1.2 is selected in conjunction with the burial depth. The upper limit is used in mining areas with strong tectonic stress, and the lower limit is used in mining areas with weak tectonic stress. The greater the confining pressure, the higher the post-peak strength of the coal body and the stronger the residual bearing capacity. Therefore, the confining pressure correction factor has a gain effect on bearing capacity.

[0106] The damage sensitivity coefficient is determined by regression fitting. The goal of regression fitting is to minimize the deviation between the bearing stress at which the coal column sample becomes unstable and the bearing stress calculated based on the damage sensitivity coefficient, the initial bearing stress, and the cumulative damage degree.

[0107] In this embodiment, the objective function for the regression fitting of the damage sensitivity coefficient is:

[0108] in, The damage sensitivity coefficient is dimensionless. For the first Peak stress at instability of a coal pillar sample, in megapascals (MPa). The initial bearing stress of the coal pillar sample is expressed in megapascals (MPa). This represents the cumulative damage at the time of instability, and is dimensionless. The number of calibration samples is determined. The regression objective is to minimize the deviation between the bearing stress at the instability of the coal pillar sample and the bearing stress calculated based on the damage sensitivity coefficient, initial bearing stress, and cumulative damage. A globally unified damage sensitivity coefficient is determined using this formula. The physical meaning of the damage sensitivity coefficient is the relative rate of bearing capacity reduction caused by unit cumulative damage. For the exponential decay model, the larger the damage sensitivity coefficient, the more significant the reduction in bearing capacity under the same cumulative damage, and the more sensitive the coal pillar is to damage. To ensure that the residual bearing stress decreases monotonically with increasing cumulative damage, the damage sensitivity coefficient should be greater than 0. If the regression result is negative, it indicates that the model does not match the calibration data, and the sample data needs to be reviewed and recalibrated. Since the field coal pillar and laboratory coal pillar samples use the same type of exponential decay model to calculate the bearing capacity reduction, the damage sensitivity coefficient determined by the laboratory regression can be directly transferred to the calculation of the residual bearing stress in the field. The consistency of the model structure ensures the transferability of the calibration parameters. When regressing multiple sets of coal pillar samples, a single damage sensitivity coefficient is fitted with all coal pillar sample data to ensure a consistent evaluation caliber within the mining area.

[0109] The initial bearing stress is determined by taking the average of the peak stresses from multiple consecutive cycles before the start of monitoring as the initial bearing stress.

[0110] In this embodiment, the formula for calculating the initial bearing stress is:

[0111] in, The initial bearing stress is expressed in megapascals (MPa). , , These represent the peak stresses of the first three cycles during the initial monitoring phase, expressed in megapascals (MPa). The average of the first three cycles is used to reduce the impact of abnormal fluctuations in a single cycle on the initial bearing capacity, allowing the initial bearing capacity stress to more stably reflect the bearing capacity of the coal pillar at the beginning of monitoring.

[0112] The residual bearing stress is calculated based on the current cumulative damage, the initial bearing stress determined by the stress response sequence, the damage sensitivity coefficient, and the scenario correction coefficient.

[0113] Multiply the current cumulative damage degree by the damage sensitivity coefficient to obtain the damage attenuation index. Calculate the attenuation factor with the natural constant as the base and the negative of the damage attenuation index as the exponent. Multiply the initial bearing stress by the attenuation factor to obtain the damage-reduced bearing stress. Multiply the damage-reduced bearing stress by the scenario correction coefficient to obtain the residual bearing stress.

[0114] In this embodiment, the formula for calculating the residual bearing stress is:

[0115] in, For the first Residual bearing stress after cycle, in megapascals (MPa). The initial bearing stress is expressed in megapascals (MPa). The damage sensitivity coefficient is dimensionless. The current cumulative damage is dimensionless. This is a scene correction coefficient, dimensionless.

[0116] The reason for choosing the exponential decay form is that the bearing capacity decreases exponentially with the degree of damage, which is consistent with the nonlinear characteristics of the bearing capacity degradation of coal pillars under cyclic loading. When the cumulative damage is small, the reduction in bearing capacity is small. As the cumulative damage increases, the reduction in bearing capacity accelerates, which is consistent with the actual behavior of accelerated degradation of coal pillars due to damage accumulation under cyclic loading.

[0117] Considering on-site monitoring errors, coal body heterogeneity, and parameter fitting errors, the residual bearing stress can also be estimated over intervals. The error propagation formula is:

[0118] in, This is an error estimate of the residual bearing stress, in megapascals (MPa). The input parameters that affect the residual bearing stress include , , , wait, For the corresponding parameter error, For the number of parameters, The confidence coefficient is 1.96, corresponding to a 95% confidence level. The reason for using the error propagation formula is that when the errors of each input parameter are independent and relatively small in magnitude, the error of the residual bearing stress can be approximated by combining the components of each parameter error under the first-order Taylor expansion. The first-order partial derivative in the formula is the error amplification factor for the corresponding parameter. For example, residual bearing stress on The partial derivatives are ,right The partial derivatives are ,right The partial derivatives are ,right The partial derivatives are The partial derivatives can be analytically obtained; when an analytical expression cannot be written, numerical difference approximation is used. The errors of each parameter are estimated according to their respective sources of uncertainty: the error of the initial bearing stress is determined by a combination of the calibration accuracy of the stress sensor and monitoring fluctuations; the errors of the damage sensitivity coefficient and the scenario correction coefficient are taken from the parameter uncertainties output after regression fitting convergence; and the error of the current cumulative damage is estimated by weighting and fusing the fitting residual of the unloading section and the integral error of hysteresis energy consumption. A confidence coefficient of 1.96 corresponds to a 95% confidence level, indicating that the actual residual bearing stress falls within the given interval with approximately a 95% probability. When determining the safety level based on this interval, the lower limit of the interval is compared with the required minimum bearing stress. The determination result is conservative, preserving a safety margin for engineering decisions. Therefore, the residual bearing stress interval is obtained as follows:

[0119]

[0120] in, This represents the lower limit of residual bearing stress, expressed in megapascals (MPa). This represents the upper limit of residual bearing stress, expressed in megapascals (MPa). The output range, rather than a single value, more accurately reflects the bearing uncertainty under complex geological conditions.

[0121] Step S108: Output the safety level based on the relationship between the residual bearing stress and the required minimum bearing stress.

[0122] In this embodiment, the required minimum bearing stress is set according to the mine design requirements or support needs. The unit is megapascal (MPa). The third-level criterion is: when the lower limit of the residual bearing stress meets... When the coal pillar is determined to be in a safe bearing zone, and the following conditions are met... At that time, the coal pillar is determined to be in the warning zone, when the following conditions are met. At that time, the coal pillar was determined to be in a high-risk area for instability.

[0123] The safe bearing capacity zone indicates that the current bearing capacity of the coal pillar meets the engineering requirements, and mining operations can proceed normally. The warning zone indicates that the bearing capacity of the coal pillar is approaching the required minimum bearing stress, and monitoring and support reinforcement are needed. The high-risk instability zone indicates that the bearing capacity of the coal pillar is below the required minimum bearing stress, and mining and tunneling must be restricted or emergency measures must be taken.

[0124] The remaining number of cycles can be estimated based on the current cumulative damage. The formula for the remaining number of cycles is:

[0125] in, To estimate the remaining tolerable number of cycles, The critical cumulative damage degree is dimensionless. The current cumulative damage is dimensionless. The average damage increment over the most recent 3 cycles, dimensionless. This indicates rounding up. This result can be used to predict the bearing capacity of the coal pillar before the next cycle of pressure, thus achieving early warning. The specific calculation method for the average damage increment is as follows: [The rest of the text appears to be incomplete and requires further context.] The damage increment of a cycle is defined as the difference in cumulative damage between two adjacent cycles, i.e. Average damage increment in the last 3 cycles The reason for using the average of the most recent three cycles instead of all historical cycles is that the damage evolution rate changes with the damage accumulation stage. In the later stages of cyclic loading, the damage increment per cycle gradually increases. If the average increment of all historical cycles were used, the smaller increments in the early stages would dilute the evolution rate of the current stage, leading to an overestimation of the remaining lifetime. The average of the most recent three cycles reflects the actual damage evolution rate of the current stage, and the prediction result can be tightened in a timely manner as damage accelerates, resulting in higher safety. If the calculated remaining cycle count is not an integer, it is rounded up to ensure that the lifetime remaining calculated based on the complete cycle count is not underestimated, and the judgment result is biased towards safety. When the average damage increment is zero or negative, it indicates that damage no longer accumulates, and the remaining cycle lifetime is not limited by the current damage rate.

[0126] If the cumulative damage calculation method is replaced with linear superposition, that is... If the cumulative damage value exceeds 1 after multiple cycles, it loses its physical meaning. For example, when the single-cycle damage values ​​for three consecutive cycles are 0.4, 0.4, and 0.4 respectively, the cumulative damage value obtained by linear superposition is 1.2, which exceeds the reasonable range of damage values. However, the cumulative damage value obtained by the product-type accumulation in this embodiment is... The values ​​should remain within a reasonable range of 0 to 1. Furthermore, linear superposition cannot reflect the nonlinear accelerating effect of damage accumulation, i.e., the weakening effect of previous damage on the foundation bearing capacity of subsequent damage. This embodiment avoids the above problems by using a product-type accumulation method.

[0127] If the weights of the three features are all fixed at 1 / 3, the discriminative role of a feature cannot be highlighted when the dispersion of a certain feature is significantly greater than that of the other features. For example, when the dispersion of the stiffness degradation feature in multiple cycles is much greater than that of the residual deformation feature and the energy dissipation feature, the stiffness degradation feature has the strongest ability to distinguish damage evolution. However, equal weighting cannot reflect this difference, resulting in a decrease in the sensitivity of damage assessment. In this embodiment, the entropy weight method is used to automatically allocate weights according to the dispersion of the data, so that the feature with a greater dispersion receives a greater weight, thereby improving the sensitivity of damage assessment.

[0128] If the laboratory calibration results are directly applied to the field coal pillar without correction, the differences between the laboratory coal pillar sample and the field coal pillar in terms of aspect ratio, equivalent confining pressure, moisture content, and fracture density will be ignored. For example, the aspect ratio of the field coal pillar is usually much larger than that of the laboratory standard coal pillar sample. Without correction, the bearing stress calibrated in the laboratory will not accurately reflect the actual bearing capacity of the field coal pillar. In this embodiment, the laboratory calibration results are corrected by a scenario correction coefficient, so that the evaluation results can be effectively transferred to the field working conditions.

[0129] The application focuses on coal pillars in deep coal mine working faces. Coal pillar width... The coal pillar is 20 meters high. It is 4 meters deep and buried at a depth of 4 meters. The depth is 800 meters, and the average unit weight of the overlying strata is... The lateral pressure coefficient is 25 kN / m³. If it is 0.8, then the on-site equivalent confining pressure is... Megapascals.

[0130] On-site coal moisture content The in-situ fracture density is 3.5%. Five strips per meter. Laboratory calibration of the standard width-to-height ratio of coal pillar samples. The standard equivalent confining pressure is 2. 10 MPa, standard moisture content 2%, standard fracture density The parameters are 2 per meter. Corrected parameters are obtained through least-squares fitting of historical samples from the same mining area. , , , Substituting into the scene correction coefficient formula, the result is obtained. Calculated It is approximately 0.79.

[0131] Critical residual strain was obtained through laboratory calibration. The critical hysteresis energy consumption is 0.008. The damage sensitivity coefficient is 0.35 megajoules per cubic meter. The critical cumulative damage is 2.5. The value is 0.85. During the initial stage of on-site monitoring, the peak stresses in the first three cycles were 18.5 MPa, 18.2 MPa, and 18.6 MPa, respectively, and the initial bearing stress was... Megapascals.

[0132] After several rounds of iterations, the current cumulative damage level is calculated. The value is 0.42. Substituting this into the residual bearing stress formula, we obtain... Calculated Approximately 5.1 MPa. If the required minimum bearing stress... If the residual bearing stress is 8 MPa, and the lower limit of the interval estimation is 4.8 MPa and the upper limit is 5.6 MPa, then the condition is satisfied. The coal pillar was determined to be in a high-risk area of ​​instability, requiring restricted mining or emergency measures.

[0133] This embodiment's method starts from on-site cyclic loading and unloading data, and through feature extraction, objective weighted fusion, product-type accumulation, scene correction, and residual bearing stress calculation, ultimately outputs a clear safety level determination result. It is feasible in engineering practice; please refer to [link / reference needed]. Figure 6 This specifically explains that the three-level judgment system can provide differentiated engineering decision-making suggestions based on the actual bearing state of the coal pillar, making engineering management more refined and scientific.

[0134] The above embodiments can be implemented, in whole or in part, by software, hardware, firmware, or any other combination thereof. When implemented in software, the above embodiments can be implemented, in whole or in part, as a computer program product. Those skilled in the art will recognize that the units and algorithm steps of the various examples described in conjunction with the embodiments disclosed herein can be implemented by electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are implemented in hardware or software depends on the specific application and design constraints of the technical solution.

[0135] The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; they may be located in one place or distributed across multiple network units. Some or all of the units can be selected to achieve the purpose of this embodiment, depending on actual needs.

[0136] The above description is merely a specific embodiment of this application, but the scope of protection of this application is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in this application should be included within the scope of protection of this application.

Claims

1. A method for evaluating the bearing capacity of coal pillars based on cyclic loading and unloading data, characterized in that, include: Obtain the stress response sequence and strain response sequence of the coal pillar under cyclic loading and unloading conditions; Multiple cycles are identified from the stress response sequence. In each cycle, the unloading stiffness, residual strain, and hysteresis energy are extracted based on the stress response sequence combined with the strain response sequence. Based on the predetermined critical residual strain and critical hysteresis energy dissipation, the unloading stiffness, residual strain and hysteresis energy dissipation of each cycle are converted into stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics, respectively. Based on the degree of data dispersion of stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics in multiple cycles, the weights corresponding to stiffness degradation characteristics, residual deformation characteristics and energy dissipation characteristics are determined. By using the weights corresponding to stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics, the stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics of each cycle are weighted and summed to obtain the single-cycle damage value of each cycle. The current cumulative damage level is obtained by performing a product-type cumulative calculation on the single-cycle damage values ​​of multiple cycles; The residual bearing stress is calculated based on the current cumulative damage, the initial bearing stress determined by the stress response sequence, the damage sensitivity coefficient, and the scenario correction coefficient. The safety level is output based on the relationship between the residual bearing stress and the required minimum bearing stress.

2. The method according to claim 1, characterized in that, The stress response sequence is smoothed to obtain a smoothed stress sequence. Based on the trend of the smooth stress sequence over time, distinguish between the loading and unloading segments; The point where the smooth stress sequence changes from rising to falling and reaches a local maximum value is taken as the peak point of each cycle; The starting point of each cycle is the point before the peak point of each cycle and whose stress is a first preset proportion of the peak stress of the previous cycle. The ending point of each cycle is the point after the peak point of each cycle and whose stress is a second preset proportion of the peak stress of the current cycle. In this way, multiple complete cycles are divided from the stress response sequence.

3. The method according to claim 2, characterized in that, Based on the stress response sequence combined with the strain response sequence, the unloading stiffness, residual strain, and hysteretic energy dissipation for each cycle are extracted, including: Linear fitting is performed on the stress and strain data of the unloading section in each cycle. The slope of the fitted line is taken as the unloading stiffness of that cycle, and the intercept of the fitted line on the strain axis is taken as the residual strain of that cycle. The area enclosed between the stress-strain curve of the loading section and the stress-strain curve of the unloading section in each cycle is taken as the hysteresis energy dissipation of that cycle.

4. The method according to claim 1, characterized in that, The unloading stiffness, residual strain, and hysteretic energy dissipation of each cycle are converted into stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics, respectively, including: Calculate the difference between the initial cyclic unloading stiffness and the unloading stiffness of each cycle, and use the proportion of this difference to the initial cyclic unloading stiffness as the stiffness degradation characteristic of that cycle. Calculate the difference between the residual strain and the initial residual strain in each cycle, and use the proportion of this difference to the difference between the critical residual strain and the initial residual strain as the residual deformation characteristic of that cycle. Calculate the difference between the hysteresis energy consumption of each cycle and the initial hysteresis energy consumption, and use the proportion of this difference to the difference between the critical hysteresis energy consumption and the initial hysteresis energy consumption as the energy dissipation characteristic of that cycle; the initial cyclic unloading stiffness, initial residual strain and initial hysteresis energy consumption are specifically the cyclic unloading stiffness, residual strain and hysteresis energy consumption under the loading and unloading conditions of the first cycle.

5. The method according to claim 4, characterized in that, The critical residual strain and critical hysteresis energy dissipation are predetermined by the following method: Cyclic loading and unloading tests were conducted on coal pillar samples from the same seam to obtain the peak stress, residual strain and hysteretic energy consumption of the coal pillar samples under each cycle. Rounds in which the peak stress drops below the preset peak stress reduction ratio of the maximum peak stress during the test and the residual strain increment exceeds the preset residual strain increment are judged as unstable rounds. The average residual strain of the instability cycle and several cycles before and after it is taken as the critical residual strain, and the average hysteresis energy consumption of the instability cycle and several cycles before and after it is taken as the critical hysteresis energy consumption.

6. The method according to claim 1, characterized in that, Based on the data dispersion of stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics across multiple cycles, the weights corresponding to each characteristic are determined, including: The stiffness degradation characteristics, residual deformation characteristics, and energy dissipation characteristics of multiple cycles are normalized to obtain the normalized characteristic sequences corresponding to each characteristic. Based on the normalized feature sequence corresponding to each feature, calculate the information entropy corresponding to each feature; The information entropy corresponding to each feature is converted into the weight corresponding to each feature, so that the feature with greater data dispersion receives greater weight; By using the weights corresponding to each feature, the stiffness degradation feature, residual deformation feature, and energy dissipation feature of each cycle are weighted and summed to obtain the single-cycle damage value of each cycle.

7. The method according to claim 6, characterized in that, The current cumulative damage level is obtained by performing a product-type cumulative calculation on the single-cycle damage values ​​from multiple cycles, including: For each cycle of damage, calculate the damage survival value. The damage survival value is specifically represented as a constant 1 minus the single-cycle damage value of that cycle. Multiply the damage survival values ​​of each cycle up to the current round to obtain the current survival product; The difference between the constant 1 and the current survival product is used as the current cumulative damage level; When obtaining the single-cycle damage value for each cycle, if the single-cycle damage value obtained by weighted summation is not greater than the preset lower limit, the single-cycle damage value for that cycle is set to the preset lower limit; if the single-cycle damage value obtained by weighted summation is greater than the preset upper limit, the single-cycle damage value for that cycle is set to the preset upper limit.

8. The method according to claim 1, characterized in that, Based on the current cumulative damage, initial bearing stress, damage sensitivity coefficient, and scenario correction coefficient, the residual bearing stress is calculated, including: Multiply the current cumulative damage level by the damage sensitivity coefficient to obtain the damage attenuation index; The attenuation factor is calculated using the natural constant as the base and the negative of the damage attenuation index as the exponent. Multiplying the initial bearing stress by the attenuation factor yields the damage-reduced bearing stress; The residual bearing stress is obtained by multiplying the damage-reduced bearing stress by the scenario correction factor.

9. The method according to claim 8, characterized in that, The scene correction coefficients are predetermined in the following way: The power function of the ratio of the width-to-height ratio of the coal pillar to the standard width-to-height ratio is used as the width-to-height ratio correction factor, and the exponent of the power function is the width-to-height ratio correction parameter. The power function of the ratio of the equivalent confining pressure on the coal pillar to the standard equivalent confining pressure is used as the confining pressure correction factor, and the exponent of the power function is the confining pressure correction parameter. Substituting the difference between the moisture content of the coal pillar and the standard moisture content into the exponential decay function yields the moisture content correction factor, which includes the moisture content correction parameter. Substituting the difference between the fracture density and the standard fracture density of the coal pillar into the exponential decay function yields the fracture density correction factor, which includes the fracture density correction parameter. Multiply the aspect ratio correction factor, confining pressure correction factor, moisture content correction factor, and fracture density correction factor together to obtain the scene correction coefficient.

10. The method according to claim 9, characterized in that, The deviation between the field-measured bearing stress based on historical samples and the bearing stress after correction by the laboratory calibration bearing stress using the scenario correction factor is minimized. The aspect ratio correction parameter, confining pressure correction parameter, moisture content correction parameter and crack density correction parameter are determined by least squares fitting. The damage sensitivity coefficient was determined by regression fitting. The goal of regression fitting was to minimize the deviation between the bearing stress at which the coal column sample became unstable and the bearing stress calculated based on the damage sensitivity coefficient, the initial bearing stress, and the cumulative damage. The initial bearing stress is determined in the following way: The average value of the peak stress from multiple consecutive cycles before the start of the monitoring phase is taken as the initial bearing stress.

Images