Study on Shear Slip Characteristics of Sandstone Plane Joints under Normal Dynamic Load Disturbance

Rock joints are susceptible to slip instability due to dynamic load disturbances such as blasting, earthquakes, and fracturation. A series of direct shear tests under the dynamic load were conducted on sandstone plane joints using the RDS-200xl. The work investigated the effects of normal static loads and normal dynamic-load frequencies and amplitudes on plane joints. Besides, the following items were proposed, that is, the peak-to-valley response rate, shear velocity vibration dominant frequency, shear-stress reduction coefficient, and discrete element numerical simulation method for plane-joint direct shear tests. The results were as follows: (1) The normal dynamic load frequency played a role in attenuating the shear stress amplitude with a threshold value of 0.5 Hz. (2) The shear velocity of the plane joint was completely controlled by the high normal dynamic load frequency. Their vibrational dominant frequencies were identical. (3) The amplitude of shear stress increased, and the median stress decreased with the increased normal dynamic load amplitude. The reduction-coefficient equation for sandstone plane joints was proposed to evaluate the shear stress under the normal dynamic load disturbance. (4) The shear-stress hysteresis phenomenon existed in the plane joints under the normal dynamic load, which required excessive shear displacements to reach peak shear strength. The peak shear displacement increased with the increased normal static load. Numerical simulations and indoor tests showed that high- and low-shear-velocity regions were the main reason for shear-stress hysteresis. The findings are conducive to revealing the shear destabilization mechanism of rock joints under dynamic load disturbance.Rock is a complex geological body composed of joints and rock masses, and shear damage along the joints is one of the main damage modes of rocks [1]. The shear behavior of rock joints is an important basis for the project design and safety assessment in practical engineering. The force form of rock engineering is a combination of dynamic and static loads due to blasting, explosion, or seismic stimulation. Compared with static-load shear damage, the shear behavior of rocks under the dynamic load is more complicated. Research findings indicate that fault slip exhibits an inherent instability, characterized by a concurrent misalignment and a decline in the stress. When an earthquake occurs, normal and shear stresses change around the fault [2-4]. Therefore, it is of great significance to study the shear behavior of rock joints under the dynamic–static load combination.Tests consider real stresses in rock engineering, and research generally focuses on dynamic loading conditions in rocks’ normal and tangential directions. Guo et al. [5] investigated the fatigue damage and irreversible deformation of salt rocks under the uniaxial cyclic loading. The fatigue life of salt rocks is mainly affected by their structure and normal stress amplitude. Liu et al. [6] carried out axial cyclic loading tests on sandstone specimens. The increased frequency of the dynamic load significantly affects axial strain and the damage pattern of sandstone. Kilgore et al. [7, 8] employed a double direct-shear apparatus to investigate the frictional response of dry and bare granite surfaces under varying normal loads, including gradual increase, gradual decrease, and pulsating loading conditions. Alterations in the normal load can enhance or diminish the frictional resistance. Tao et al. [9] used the DJZ-500 shear apparatus to study the shear characteristics of saw-cut granite fractures infilled with quartz under different normal load amplitudes, cycles, and sliding velocities. Dang et al. [10] analyzed the effect of shear velocity, frequency, and amplitude of dynamic normal stress on the dynamic friction characteristics of dredged soil.These studies primarily focus on the shear velocity and cyclic shear associated with the tangential dynamic loading. Belem et al. [11] conducted cyclic shear tests on specimens with different shapes and proposed two models for degrading the surface roughness of rock joints. The degradation characteristics of structural surfaces are different under cyclic shear and primary shear. Dang et al. [12-14] investigated the flat structural surfaces under CNL conditions with direct-shear behavior. A new shear-resistant equation is proposed, indicating that the minimum coefficient of the dynamic friction increases with the increased slip velocity. Liu et al. [15] systematically investigated the strength, deformation characteristics, and macro- and fine-scale cumulative damage evolution of serrated joints in original rocks. Zhou et al. [16] studied the shear slip mechanism and failure characteristics of structural planes with different sawtooth sizes, structural-plane roughness, shear rates, and normal stresses through laboratory tests. Okada et al. [17] conducted monotone and cyclic loading tests on rock samples, finding that dynamic shear strength is greater than static shear strength. A mathematical model is proposed to estimate dynamic shear strength. Jafari et al. [18] investigated the variable shear strength of rock joints under cyclic loading. The shear strength of joints is related to the shear velocity, the number of loading cycles, and the stress amplitude. Mirzaghorbanali et al. [19] carried out cyclic shear tests on rock joints in the laboratory under CNS conditions. Shear strength decreases as the loading cycle and shear rate increase. The effect of the shear velocity becomes less pronounced when the normal load is increased. Nguyen et al. [20] conducted cyclic shear tests on slate. The peak shear stress of nodular rocks increases with time under dynamic loading. Aguilar et al. [21] evaluated the macroscopic energy dissipation of concrete under cyclic shear using a punch-through shear test FE model.The aforementioned studies have yielded valuable insights into the comprehensive comprehension of dynamic shear behavior exhibited by rock joints. However, the experimental study of the direct shear still lacks comprehensive tests and analysis under dynamic normal loading conditions due to the limitation of the direct shear equipment. Based on this, sandstone from the same origin was used to fabricate the same plane joints. The RDS-200xl fully automatic servo-controlled rock direct shear equipment was applied to investigate the shear slip characteristics of plane joints under different normal static loads, normal dynamic load amplitudes, and normal dynamic load frequencies.The test program design was based on the controlled single-variable method. The normal static load, normal dynamic load frequency, and normal dynamic load amplitude were used as the research factors for the test. Table 1 shows the test program.The dynamic load frequency test program controlled the normal dynamic load amplitude at a constant value of ±0.50 MPa. Normal static loads were 1, 3, and 5 MPa. The instability characteristics of the plane joint were studied at five normal dynamic loading frequencies (0.10, 0.25, 0.50, 0.75, and 1.00 Hz). The dynamic load amplitude test program controlled the normal dynamic load frequency to a constant value of 0.25 Hz and the normal static load to 1, 3, and 5 MPa. The instability characteristics of plane joints were studied under normal dynamic load amplitudes of ±0.25, ±0.50, ±0.75, ±1.00, and ±1.50 MPa.The specimens are sandstone (100 × 100 × 50 mm). They are labeled in the format of normal static load, normal dynamic load frequency, and normal dynamic load amplitude according to the above test protocol (Figure 1).The mechanical parameters of the sandstone specimens are identified by uniaxial compression tests on rectangular sandstone specimens of 50 × 50 × 100 mm. The uniaxial compressive strength of the specimen is 54.2 MPa, and the elasticity modulus is 3.1 GPa.The test is accomplished with the GCTS RDS-200xl equipment (Figure 2). The original cylindrical shear box is replaced by a square shear box designed in the work. The inner diameter of the shear box is 102 × 102 mm, and the specimen is fixed by squeezing bolts on both sides. The normal driver is installed on the normal rigid support, and the normal rigid support is installed on the plain bearing. The normal device moves with the upper shear box in the shearing process, which can minimize the horizontal friction. The normal load acts vertically on the entire shear displacement surface. Once the sample has been fixed, the upper shear box is positioned onto the platform in alignment with the lower shear box, and it is subsequently secured to the tangential drive using a cotter pin. A gradient is used to check whether the shear direction is horizontal before running the shear driver. Therefore, the shear direction is consistent with the plane joint.Statically loaded direct shear tests are performed on sandstone P (100 × 100 × 50 mm). Figure 3(a) shows the relationship curves of shear stress-shear displacements under three normal static stresses (1, 3, and 5 MPa). The test results are fitted and presented in Figure 3(b), illustrating a linear relationship between the normal stress and peak shear strength according to the Mohr–Coulomb strength criterion. Based on the fitting results of the curve in Figure 3(b), the basic friction angle of the plane joint of the specimen is calculated to be 33.82°.The dynamic load direct shear test process is divided into two stages. The normal stress is applied to the test system using the stress-control normal drive with a normal load rate of 0.05 MPa/s at the first stage. The pressure is kept for 20 seconds after reaching the target value to maintain the accuracy and stability of the normal stress. The second stage commences simultaneously with the application of the normal dynamic load and activation of the shear driver, following the stabilization of the normal static load. The sine wave is superimposed based on the static stress (equation. (1)). The shearing process is subjected to displacement control with a speed of 0.02 mm/s and a displacement of 10 mm.where σsd is the normal superimposed stress; σd is the normal dynamic load amplitude; f is the frequency; t is time.Figure 4 shows the shear-stress and shear-displacement curves at a normal dynamic load amplitude of ±0.5 MPa, a frequency of 0.25 Hz, and normal static loads of 1, 3, and 5 MPa, respectively. The normal stress interval is consistent with the set value, and the shear stress increases nonlinearly with the increased shear displacement. The shear stress remains in a relatively stable interval with the fluctuated normal stress after reaching the peak shear stress, and the median shear stress is unchanged. The fluctuation intervals of shear stresses are 0.33–0.54, 1.05–1.26, and 1.91–2.10 MPa for normal static loads of 1, 3, and 5 MPa, respectively.Figure 5 presents the test results (peak and valley values of normal and shear stresses) of the dynamic load frequency test program. The peak shear stress of the specimens decreases significantly after applying the dynamic loading by comparing the peak shear stress of specimens under static loading. The normal dynamic load is applied with the same frequency of 0.1 Hz. The shear stress of the specimens decreases from 0.66, 2.05, and 3.345 MPa to 0.612, 1.526, and 2.218 MPa (with decreases of 7.3%, 25.6%, and 33.69%) for the normal static loads of 1, 3, and 5 MPa, respectively. The higher the normal static load, the more the decrease in the peak shear stress of specimens after applying the dynamic load.The median shear stresses are 0.47, 1.22, and 2.05 MPa for normal static loads of 1, 3, and 5 MPa, respectively. The variation coefficients (standard deviations divided by the mean value) of the median shear stresses are 0.047 (a static load of 1 MPa), 0.056 (a static load of 3 MPa), and 0.043 (a static load of 5 MPa) at different frequencies for the three static loads, respectively. The variability of the median shear stress is “very low” based on the criteria of the parameter variation in the specification [22]. The increased dynamic loading frequency mainly affects the amplitude of the shear stress and has a low effect on the median shear stress.The peak-to-valley response rate of the shear stress is defined bywhere τmax and τmin are the peak and valley of the shear stress, respectively; σmax and σmin are the peak and valley of the normal stress, respectively.Figure 6 shows the relationship between the peak-to-valley response rate and normal dynamic loading frequency. The peak-to-valley responsivity of the specimen decreases from 32.40%, 44.31%, and 38.66% to 13.13%, 7.56%, and 6.17% for the normal static loads of 1, 3, and 5 MPa, respectively. Fitting the peak-to-valley response rate of shear stress to the normal dynamic loading frequency yields an exponential function. Fitted function curves suggest a threshold value of 0.5 Hz for the attenuating effect of the normal dynamic-loading frequency on the peak-to-valley response rate of the shear stress.The peak-to-valley response rate of the shear stress decreased by 13.8%, 30.3%, and 26.9% for specimens with static loadings of 1, 3, and 5 MPa, respectively, before reaching the dynamic loading frequency threshold (0.5 Hz). The peak–valley response rate of the shear stress decreases significantly with the increased normal frequency. The peak-to-valley response rate of the shear stress decreases by 4.9%, 6.2%, and 6.2% for the specimens with static loadings of 1, 3, and 5 MPa, respectively, after exceeding the dynamic loading frequency threshold (0.5 Hz). The reduction of the shear-stress peak–valley response rate slows down and stabilizes.In conclusion, the normal dynamic load frequency contributes to the reduced shear stress amplitude. The attenuation effect becomes more noticeable as the normal static load increases. The plane joint in the process of dynamic load shear has the characteristics of high normal dynamic load frequency and low-shear-stress amplitude.The change in the normal dynamic load frequency affects the shear slip of the structural plane [23, 24]. Therefore, the evolution of the shear velocity should be studied during the shear process. Figure 7(a) shows the shear displacement–time curves at different normal dynamic load frequencies for the test with a normal static load of 1 MPa and a normal dynamic load amplitude of ±0.5 MPa. The shear displacement curve fluctuates. The shear velocity is obtained by calculating the slopes of each point on the shear displacement–time curve for the whole shearing process. Figure 7(b) presents the shear velocity–time curves for specimens with a normal static load of 1 MPa, a normal dynamic load amplitude of ±0.5 MPa, and a normal dynamic load frequency of 0.1 Hz. The shear velocity is affected by the normal dynamic load, which produces significant fluctuations—floating up and down around the test preset value of 0.02 mm/s. The changing spectrum of the shear velocity in the shear process is analyzed using the short-time Fourier transform (STFT). It verifies the correlation between the shear velocity and the normal dynamic load frequency. Therefore, the spectrum can reflect the time-localized characteristics [25].Figure 8 shows amplitude–vibration frequency data obtained after the STFT. The frequency corresponding to the peak of the amplitude spectrum curve is taken as the dominan frequency (DF) of the shear velocity vibration . Specimens have no significant dominant frequency of the shear-velocity vibration at a normal static load of 1 MPa and a normal dynamic load frequency of less than 0.5 Hz. When the normal dynamic load frequency exceeds 0.5 Hz, the dominant frequency of the shear-velocity vibration becomes apparent and aligns with the set normal dynamic load frequency.A similar situation occurs in the test with a normal static load of 3 MPa. The dominant frequency of the sample’s shear-velocity vibration is obvious with a normal dynamic load frequency greater than 0.5 Hz. It is the same as the set normal dynamic load frequency. The dominant frequency of the shear velocity vibration is not significant at the normal dynamic load frequencies of 0.1, 0.25, 0.5, and 0.75 Hz for the specimens with a normal static load of 5 MPa. Only when the dynamic load frequency is 1 Hz, the dominant frequency is the same as the set normal dynamic load frequency.In summary, the normal dynamic load frequency controls the shear velocity but is limited by the normal static load. The higher the normal static load, the greater the threshold of the dynamic load frequency for controlling the shear velocity. The low-frequency normal dynamic load dominates the shear velocity for specimens under low normal static loads. As the normal static load increases, the low-frequency normal dynamic load no longer controls the shear velocity. The shear velocity of the plane joints, which is controlled by the normal dynamic load frequency of the junction, exhibits regular undulations and is characterized by periodicity. The dominant frequency of the shear velocity vibration is the same as the normal dynamic load frequency.Figure 9 shows the relationship between the shear stress and normal dynamic load amplitude based on test data from the test program of the normal dynamic load amplitude. The shear stress range increases with the increased normal dynamic load amplitude. The peak and valley values of shear stresses decrease significantly, and the peak–valley response rate curve of shear stresses decreases nonlinearly.The peak-to-valley response rate decreases with the increased normal static load for dynamic load amplitudes larger than 0.25 MPa. This is because the plane joints at low normal stresses are more prone to slip when the normal stress amplitude becomes large. The specimen needs more shear displacements at high normal stresses to reach the peak shear stress. However, the peak shear stress is difficult to rise due to the periodic character of the normal dynamic load.It is worth noting that the maximum shear stress does not increase in proportion to the increased maximum normal stress. The conventional Mohr–Coulomb criterion may have overestimated the shear stress in sandstone plane joints under combined dynamic and static loads in direct shear tests.Figure 10 shows the curve fitting of the median shear stresses under three normal static loads. The fitted slopes are −0.3560, −0.3929, and −0.5727 when the normal static loads are 1, 3, and 5 MPa, respectively. The higher the normal static load of the plane joints, the greater the decrease in the median shear stress with the increased normal dynamic load amplitude.The shear stress can be expressed as follows under the Mohr–Coulomb criterion for a combination of dynamic and static loads.According to equation. (3),The median shear stress is defined byequations (5) and (6) are substituted into equation. (7) to obtainThe median shear stress is consistent with the results of the Mohr–Coulomb criterion equation. (9).However, Figures 3 and 9 are compared to obtainTherefore, when a normal dynamic load is added to the direct shear test, the shear stress reduces based on the Mohr–Coulomb criterion. That is, the shear-stress reduction coefficient Kd is obtained byCombined with the fitting results in Figure 10,The equation for the shear-stress reduction coefficient under dynamic loading is defined byFigure 11 shows the relationship between the shear stress and normal stress. A phase shift exists between the normal stress and shear stress. This is because the plane joints require excessive shear displacements to reach the peak shear stress in the dynamic load direct-shear test. The time differences between the peaks of the shear and normal stresses are 0.2, 0.4, and 0.6 seconds under static loads of 1, 3, and 5 MPa, respectively. The shear displacement required to reach the maximum shear stress increases with the increased normal static load.Calculated friction coefficient K=τ/σ and the friction coefficients fluctuate cyclically. Small fluctuations occur at higher normal static loads, indicating the stick–slip motion of the plane joint. The peak friction coefficient lags behind the peak normal stress by approximately half a cycle in all tests. equation. (14) shows the relationship with the normal dynamic load period. The relative time offset is almost constant, which is consistent with Reference 12.where Δt is the time difference between the peak friction coefficient and the peak normal dynamic load; TPeriod is the normal dynamic load period.Figure 12 shows the shear displacement, shear velocity, and normal stress versus time at a normal static load of 3 MPa, a normal dynamic load amplitude of ±0.5 MPa, and a frequency of 0.25 Hz. The base normal stress of 3 MPa and the shear velocity of 0.02 mm/s are marked as the baseline for the fluctuated normal dynamic load and shear velocity. The red datum line represents the normal stress, while the blue datum line corresponds to the shear velocity in the graph.The shear velocity regularly changes with the fluctuated normal stress. The shear velocity is significantly suppressed when the normal stress enters the stage above the datum line. The shear velocity appears to be below the datum line with a minimum value of about 0.005 mm/s. The shear velocity is significantly activated when the normal stress enters the stage below the datum line. The shear velocity is above the datum line, with a maximum value of about 0.04 mm/s.The logarithmic relationship between the joint friction coefficient and the shear velocity is satisfied at the shear rate of 10−4 to 10−1 mm/s (equation. (17)) [26]. The shear stress of the structural plane decreases with the increased shear velocity. Therefore, the macroscopic shear velocity is affected by the normal dynamic load, and high- and low-shear-velocity regions may be the cause of the hysteresis phenomenon of the shear stress in the test.The shear process of the plane joint is simulated with discrete element software PFC2D to further verify the conjecture. The particle contact inside the upper and lower disks is defined as a linear parallel bond contact (the linear parallel bond model), and the contact between the particles is defined as a smooth-joint contact (smooth-joint contact model). Table 2 shows the parameters of particles and contact. The joint plane exhibits a certain level of roughness and undulating heights at the μm scale, which is attributed to the occasional generation and arrangement of particles. It is basically in line with plane-joint roughness [27].The contact of the structural plane will be moved by force during the shear process. If the shear displacement is too large, the contact between the particles will be broken, and a new contact will be created. New contacts are usually assigned as the linear parallel bond model by default in the PFC2D. This mechanism increases the shear resistance and dilatancy of the plane joint. Contact generation is updated in real-time by programming in Fish language. The contact generated by the plane joint is always defined smooth joint.The program code of equation. (1) is written in Fish language to implement the dynamic loading after the basic model is constructed. Note that the numerical model does not precisely replicate the experimental results. It deepens our understanding of the shear stress hysteresis phenomenon by demonstrating the microphysical properties of the plane joint. The complete numerical model is shown in Figure 13(a).Figure 13(b) shows the shear displacement–shear stress curve obtained by the numerical simulation at the normal stress of 3 MPa and the dynamic load amplitude of 0.5 MPa. The numerical simulation results replicate the shear stress curve observed in the test.The shear stresses show similar waveform undulations with the increased normal dynamic load. The shear stress interval is from 1.3 to 1.6 MPa, which is similar to the indoor test results. There is also significant shear-stress hysteresis in the locally amplified normal- and shear-stress curves. The characteristics of the evolution of the physical parameters of the particles on the plane joints are monitored through Fish language (Figure 13(c)). The shear velocity of the particles varies with the normal stress. High- and low-shear-velocity regions also exist.The shear velocity of particles is significantly suppressed, and the peak shear velocity decreases rapidly from 0.06 m/s when the normal stress enters the high normal stress region (greater than the average value of 3 MPa). The shear velocity of the particles is significantly activated—rising rapidly from a valley shear velocity of nearly 0 m/s when the normal stress is less than the average value of 3 MPa.equations. (15) and (16) are used to calculate the smooth joint contact in the numerical model [28]. The increased shear velocity decreases the equivalent cross-sectional area of the contact. Shear stress between the contacts consequently decreases. Eventually, the macroscopic shear strength of the planar joints decreases.where Fs is contact force; ks is the tangential stiffness of the smooth-joint contact; A is the equivalent cross-sectional area of the contact; δs is the shear displacement increment; R(1) and R(2) are the radii of the two particles in contact.In summary, the comparison of the numerical simulation and the data from indoor tests shows that the shear-stress hysteresis of normal dynamic-shear tests may be caused by a combination of two factors. First, the normal dynamic load exerts inhibitory and activating effects on the shear velocity. High- and low-velocity shear regions appear, which affects the shear stress of plane joints. Second, from microscopic contact, the shear velocity of individual contacts increases during the reduction phase of the normal dynamic load. The equivalent cross-sectional area of contact increases with the increased shear stress, which leads to a phase shift between normal and shear stresses.The effects of static load, dynamic load amplitude, and dynamic load frequency on the shear slip characteristics of sandstone plane joints were studied by the RDS-200xl automatic servo-controlled rock direct shear test. The main conclusions are as follows:The normal frequency attenuated the amplitude of shear stress and exhibited an exponential relationship with the peak-to-valley response rate. The threshold value was 0.5 Hz. The peak–valley response rate decreased with the increased frequency within the threshold range. Amplitude attenuation decreased with the increased frequency beyond the threshold range. Besides, the plane joint in the process of dynamic load shear has the characteristics of high normal dynamic load frequency and low-shear-stress amplitude.The shear displacement of the plane joint showed an alternating state of slip-stagnation under the normal dynamic load disturbance. The frequency domain diagram of the shear velocity after the STFT showed that the normal dynamic load frequency controlled the shear velocity. The main frequency of the shear velocity vibration was the same as the normal dynamic load frequency.The shear-stress amplitude increased and the peak-to-valley shear stress decreased with the increased normal dynamic load amplitude. The higher the static load, the lower the peak-to-valley response rate at high dynamic load amplitudes. The median shear stress decreased linearly with the normal dynamic load amplitude. There was a reduction coefficient of the shear stress of the sandstone plane joint under a normal dynamic load direct shear test.There was a phase shift between normal and shear stresses during normal dynamic load direct shear tests. The peak shear stress occurred after the normal stress reached its peak. Plane joints require more shear displacements to reach the peak shear stress. Shear displacements required to reach the maximum shear stress increased with the increased normal static load. The ratio of the time between the peak friction coefficient and the peak normal dynamic load to the normal dynamic load period was approximately 0.5.Based on DEM-simulated stress hysteresis under normal dynamic load shear, the proposed joint dynamic load shear model exhibited smooth behavior. Combined with the indoor tests, the change in the normal stress dominated the evolution of the shear stress. High- and low-shear-velocity regions were the main reason for shear-stress hysteresis.The normal frequency attenuated the amplitude of shear stress and exhibited an exponential relationship with the peak-to-valley response rate. The threshold value was 0.5 Hz. The peak–valley response rate decreased with the increased frequency within the threshold range. Amplitude attenuation decreased with the increased frequency beyond the threshold range. Besides, the plane joint in the process of dynamic load shear has the characteristics of high normal dynamic load frequency and low-shear-stress amplitude.The shear displacement of the plane joint showed an alternating state of slip-stagnation under the normal dynamic load disturbance. The frequency domain diagram of the shear velocity after the STFT showed that the normal dynamic load frequency controlled the shear velocity. The main frequency of the shear velocity vibration was the same as the normal dynamic load frequency.The shear-stress amplitude increased and the peak-to-valley shear stress decreased with the increased normal dynamic load amplitude. The higher the static load, the lower the peak-to-valley response rate at high dynamic load amplitudes. The median shear stress decreased linearly with the normal dynamic load amplitude. There was a reduction coefficient of the shear stress of the sandstone plane joint under a normal dynamic load direct shear test.There was a phase shift between normal and shear stresses during normal dynamic load direct shear tests. The peak shear stress occurred after the normal stress reached its peak. Plane joints require more shear displacements to reach the peak shear stress. Shear displacements required to reach the maximum shear stress increased with the increased normal static load. The ratio of the time between the peak friction coefficient and the peak normal dynamic load to the normal dynamic load period was approximately 0.5.Based on DEM-simulated stress hysteresis under normal dynamic load shear, the proposed joint dynamic load shear model exhibited smooth behavior. Combined with the indoor tests, the change in the normal stress dominated the evolution of the shear stress. High- and low-shear-velocity regions were the main reason for shear-stress hysteresis.The relevant data used in this paper are available from the authors upon request.The authors declare that they have no known personal relationships or competing economic interests that may affect the work reported in this work.The authors gratefully acknowledge the financial support for this work provided by Shandong Energy Group (NO. SNKJ2022BJ03-R28) and the National Natural Science Foundation of China (52104102) and a Project Supported by the Scientific Research Fund of Zhejiang Provincial Education Department (Y202146135).