Rock nanoindentation numerical simulation method based on PFC discrete element

Abstract

The invention discloses a rock nanoindentation numerical simulation method based on PFC discrete elements, and relates to an unconventional natural gas reservoir fracturing technology. Comprising the following steps: dynamically loading an iterative matching load-displacement curve based on XRD-SEM-EDS mineral spatial distribution of a rock sample and single mineral nanoindentation experimental data microscopic parameters; in combination with mineral space distribution data, constructing a multi-mineral heterogeneous model by stages, and dynamically loading simulation; based on crack number development, force chain distribution and particle displacement vector characteristics, a regulation and control mechanism of a mineral interface coupling effect and pore stress concentration on a crack propagation path is quantitatively analyzed. According to the invention, a heterogeneity evolution mechanism of force chain development and crack propagation in each stage of loading, load holding and unloading can be analyzed, and a dynamic mechanism that mineral components and pore structures cooperatively regulate and control a crack propagation path and nonlinear load response is disclosed; and a solution path is provided for rock microscopic damage mechanism prediction and heterogeneous material discrete element numerical simulation modeling.

Description

Technical Field

[0001] The present invention relates to the technical field of unconventional natural gas reservoir fracturing, and in particular to a rock nanoindentation numerical simulation method based on PFC discrete elements. Background Art

[0002] In the field of unconventional natural gas utilization, its strategic position as a clean energy source is becoming increasingly prominent. Conventional hydraulic fracturing technology relies on the construction of macroscopic fracture networks. However, due to the complex microscopic heterogeneity of unconventional natural gas reservoirs, existing theoretical frameworks struggle to effectively reveal the cross-scale correlations between mineral composition, pore structure, and crack propagation.

[0003] Current research often uses nanoindentation technology to characterize micromechanical parameters, supplemented by continuous medium simulation methods such as finite element methods. Nanoindentation technology uses nanoscale probes to puncture the rock surface and invert micromechanical parameters (such as elastic modulus and hardness) from load-displacement curves. The "pop-in phenomenon" (sudden drop in load) observed in the experiment corresponds to fracture at the mineral particle interface, while the "elbow feature" reflects the creep properties of clay minerals. Traditional methods such as finite element method (FEM) assume that the medium is continuous and have difficulty characterizing contact failure between mineral particles (such as frictional sliding and cementation fracture). Therefore, they face theoretical bottlenecks in analyzing contact failure in multiphase media and the dynamic evolution path of cracks. In particular, for non-continuous mechanical behaviors such as weak cementation effects at mineral interfaces and pore-induced stress concentration, traditional numerical methods have difficulty achieving high-precision simulation of multiphase coupling.

[0004] In recent years, the particle flow discrete element method (PFC) has shown breakthrough advantages. Compared with the continuous medium method, PFC can accurately characterize the contact constitutive relationship between rock mineral phases by constructing a multi-scale particle system. Its explicit dynamics algorithm iteratively calculates the particle motion trajectory through Newton's law, and can dynamically capture the entire process of force chain reconstruction and crack network evolution during loading and unloading. Studies have confirmed that PFC simulation has significant advantages in reproducing key features such as the pop-in phenomenon, elbow and creep displacement response of nanoindentation experiments. However, the traditional GBM modeling method (Grain-Based Model) discretizes rocks into a collection of mineral particles in contact with each other, and simulates the rock failure process by defining parameters such as the bonding strength and friction coefficient between particles. It has technical bottlenecks such as the homogenization of mineral particle size, and in the past, the discrete element simulation results of nanoindentation showed non-physical oscillations, resulting in limited credibility of the simulation results.

[0005] Patent publication number CN116911040A discloses a method for calibrating discrete element parameters of granite using nanoindentation testing. This document constructs a heterogeneous model based on Voronoi polygons and calibrates microscopic parameters using nanoindentation data, improving the scientific accuracy of the parameter calibration. However, this method uses homogenized Voronoi polygon modeling technology and fails to consider the spatial heterogeneity of shale's multi-mineral phases. This results in a deviation between the mineral particle size distribution and the actual geological structure, making it impossible to achieve high-precision simulation of the microscopic damage mechanism of heterogeneous shale. Summary of the Invention

[0006] To overcome the limitations of existing technologies, this paper proposes a numerical simulation method for rock nanoindentation based on PFC discrete elements. By combining a force servo-controlled dynamic loading algorithm with multi-mineral phase modeling technology, a heterogeneous model system containing pore defects is constructed. By integrating XRD-SEM-EDS multimodal mineral distribution data, the spatial heterogeneity of the rock's main mineral mechanical properties is accurately characterized, revealing the crack network evolution and nonlinear load response characteristics synergistically regulated by the mineral interface coupling effect and pore structure.

[0007] The technical solution adopted in the present invention is as follows:

[0008] A numerical simulation method of rock nanoindentation based on PFC discrete element method includes the following steps:

[0009] Step 1: Based on the XRD-SEM-EDS mineral spatial distribution and single-mineral nanoindentation experimental results of the rock, the microscopic parameters of the single-mineral homogeneous model of the rock are calibrated. The dynamic loading algorithm controlled by force servo is used for loading. The load-displacement curve is iteratively matched through trial-and-error iteration to achieve a high degree of consistency between the simulated and experimental load-displacement curves.

[0010] Step 2: Using the calibrated microscopic parameters and combined with the XRD-SEM-EDS mineral spatial distribution data, a multi-mineral phase separation modeling method is adopted to construct a multi-mineral heterogeneous model, and a force servo-controlled dynamic loading algorithm is used to perform dynamic loading simulation.

[0011] Step 3: Monitor the loading process data and quantitatively analyze the regulation mechanism of the mineral interface coupling effect and pore stress concentration on the crack propagation path based on the crack number development, force chain distribution and particle displacement vector characteristics.

[0012] Furthermore, the step 1 specifically includes:

[0013] In step 1.1, based on the joint characterization results of the mineral spatial distribution using XRD-SEM-EDS, the indentation point data of large-particle minerals in the nanoindentation experiment are selected as the basis for calibrating the parameters of the single mineral homogeneous model. The elastic modulus and hardness of the mineral sample can be calculated based on the load-displacement curve.

[0014] Step 1.2: Establish a PFC 2D single mineral homogeneous model, optimize the size and total number of particles, balance the stress through servo control, and establish a parallel cementation system, including four steps: sampling, pre-pressing, cementation, and pressing head.

[0015] In step 1.3, a force servo-controlled dynamic loading algorithm is used for loading. Loading and unloading are controlled separately. The process follows a one-way control logic to eliminate reverse loading disturbances. The Y-axis data source is switched from the real-time load to the target load, and the target load is recorded synchronously to reproduce the quasi-static response of the test.

[0016] In step 1.4, the microscopic parameters (including particle size, modulus, cementation strength, etc.) are reversely adjusted through the trial-and-error iterative method until the error between the simulated and experimental load-displacement curves is less than the set percentage (e.g., ≤10%), and the simulated and experimental indentation failure characteristics are similar, that is, the crack propagation pattern is consistent with the experiment, and the microscopic parameters of the single mineral homogeneous model are obtained.

[0017] Furthermore, the specific process of step 2 is as follows:

[0018] Step 2.1: Based on the XRD-SEM-EDS joint characterization data, use AutoCAD to perform geometric topological simplification on the spatial distribution structure of the mineral components.

[0019] Step 2.2: Construct five characteristic models, including clay macropores, crack defects, block distribution, layered distribution, and mosaic distribution models. Here, certain modeling guidelines must be followed, considering factors such as model size, pore defects, and mineral distribution, and selecting appropriate areas for modeling.

[0020] Step 2.3: Use the multiphase separation-recombination modeling method to establish a multi-mineral heterogeneous discrete model containing pores.

[0021] Step 2.4: Import the single mineral microscopic parameters obtained by calibration in step 1, use the force servo to control the indenter to perform dynamic loading, and obtain the simulation experimental results by loading.

[0022] Furthermore, in step 3, during the loading process, the load-displacement curve of the model, the number of different types of cracks, occurrence locations, occurrence angles, occurrence frequencies, force chain size and evolution laws, particle displacement vectors and other data are monitored and analyzed.

[0023] The beneficial technical effects brought about by the present invention are as follows.

[0024] This paper proposes a numerical simulation method of rock nanoindentation based on PFC discrete elements, constructs a three-dimensional quantitative analysis system for the synergistic effect of mineral spatial distribution, crack dynamic evolution and pore structure, breaking through the limitations of traditional experimental two-dimensional observation.

[0025] This study utilizes a force servo-controlled dynamic loading algorithm, combined with XRD-SEM-EDS multimodal mineral characterization data for multi-mineral phase modeling, to construct a nanoindentation discrete element model of heterogeneous rocks containing pore defects. This model breaks through the traditional GBM modeling method's assumption of homogeneous mineral particle size. It accurately reveals the regulatory mechanisms of mineral particle size heterogeneity and pore distribution on crack initiation and propagation, achieving high-precision simulation of the microscopic damage mechanism of heterogeneous rocks. The force servo-controlled dynamic loading algorithm effectively suppresses the abnormal oscillation of the load-displacement curve seen in previous nanoindentation PFC discrete element numerical simulations, enabling high-precision replication of nanoindentation experiments.

[0026] The present invention can analyze the heterogeneous evolution mechanism of force chain development and crack propagation in the loading, load holding and unloading stages, and reveal the dynamic mechanism of the coordinated regulation of crack propagation path and nonlinear load response by mineral components and pore structure, providing a solution for the prediction of rock micro-damage mechanism and discrete element numerical simulation modeling of heterogeneous materials. BRIEF DESCRIPTION OF THE DRAWINGS

[0027] Figure 1 is a flow chart of the simulation method of the present invention;

[0028] Figure 2 It is a flow chart for constructing a single mineral homogeneous model;

[0029] Figure 3 This is a flow chart for establishing a PFC 2D single mineral homogeneous model;

[0030] Figure 4 It is a flow chart for constructing a multi-mineral heterogeneous model;

[0031] Figure 5 It is a flow chart for establishing the PFC 2D multi-mineral heterogeneous model;

[0032] Figure 6 It is the load-displacement curve and failure mode diagram of the numerical simulation of the single mineral quartz homogeneous model;

[0033] Figure 7 This is a numerical simulation failure analysis diagram of a multi-mineral heterogeneous model. DETAILED DESCRIPTION

[0034] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments:

[0035] This example proposes a numerical simulation method for rock nanoindentation based on PFC discrete elements. Through cross-scale coupling analysis of nanoindentation experiments and discrete element simulations, combined with XRD-SEM-EDS multimodal mineral characterization data, a multi-mineral heterogeneity model considering pore defects is constructed. Using a force servo-controlled dynamic loading algorithm, the regulation mechanism of rock microheterogeneity on mechanical behavior is clarified, revealing the intrinsic relationship between the mechanical properties of mineral components and the dynamic evolution of cracks.

[0036] See also Figure 1 This embodiment is based on the rock nanoindentation numerical simulation method using PFC discrete elements, taking shale as an example, and mainly includes the following three steps:

[0037] Step 1: Based on the XRD-SEM-EDS and single-mineral nanoindentation experimental results of the shale sample, the microscopic parameters of quartz, calcite, and clay were calibrated. Dynamic loading was performed using a force servo-controlled dynamic loading algorithm. Through iterative matching using a trial-and-error iterative method, a high degree of consistency between the simulated and experimental load-displacement curves was achieved, and a single-mineral homogeneous model of shale nanoindentation was constructed.

[0038] Step 2: Using the microscopic parameters calibrated in step 1 and combined with the XRD-SEM-EDS mineral spatial distribution data, multi-mineral phase modeling is performed to construct a shale nanoindentation multi-mineral heterogeneous model and perform dynamic loading simulation.

[0039] Step 3: Based on the crack number and development, force chain distribution, and particle displacement vector characteristics, quantitatively analyze the mechanism by which mineral interface coupling and pore stress concentration regulate crack propagation paths. During the loading process, monitor and analyze the model's load-displacement curve, the number of different crack types, their locations, inclination angles, and frequencies, the size and evolution of force chains, and particle displacement vectors.

[0040] like Figure 2 As shown, the single mineral homogeneous model construction process in step 1 specifically includes the following steps:

[0041] Step 1.1: Based on the XRD-SEM-EDS combined characterization results, the spatial distribution of the three main minerals, quartz, calcite, and clay, was divided. The indentation point data of large-particle minerals from the nanoindentation experiment was selected as the basis for calibrating the parameters of the single-mineral homogeneous model. Based on the load-displacement curve, the elastic modulus and hardness of the mineral sample can be calculated. The specific equation is as follows:

[0042]

[0043] Where H is the hardness, Pmax is the maximum load, Ac is the contact area, Er is the reduction modulus, S = dP / dh is the contact stiffness, ɑ = 1.034 is the shape factor of the Bohli indenter, E and v are the elastic modulus and Poisson's ratio of the shale sample, respectively, Ei and Vi are the elastic modulus and Poisson's ratio of the indenter, respectively, and Ei = 1140 GPa and vi = 0.07 for the Bohli indenter.

[0044] Step 1.2: Construct a PFC 2D single mineral homogeneous model, optimize the size and total number of particles, balance the stress through servo control and establish a parallel cementation system, including sampling, pre-pressing, cementation, press head, and loading. The specific process is as follows: Figure 3 As shown:

[0045] A two-dimensional discrete element model of the shale specimen was established. Based on the principle of experimental comparability, an array of 10μm equally spaced indentations was constructed to eliminate boundary effects. While maintaining computational accuracy, the total number of particles was optimized by setting the specimen aspect ratio to 0.5 (10μm × 5μm), effectively improving computational efficiency.

[0046] A servo control mechanism was used to adjust the confining pressure baffle displacement, achieving initial stress equilibrium in the specimen based on a particle elastic modulus threshold. This equilibrium stress threshold, determined based on the particle elastic modulus, allows for a small amount of elastic deformation in the specimen, which has been shown to have no impact on subsequent experimental results.

[0047] A linear parallel bond system between particles is constructed. When constructing, the contact between particles adopts the parallel bond (ParallelBond) model, and the particle-baffle contact maintains a linear relationship. The calculation formulas for the total force and total torque are as follows:

[0048]

[0049] Among them, F c is the total force, F l is the linear force, F d is the damping force, is the bonding force on the interface. c is the total torque, is the torque.

[0050] The upper baffle was removed, and the remaining three baffles were fixed. After the specimen was balanced, a two-dimensional model of the Berkovich indenter was constructed. Based on previous research, the three-dimensional structure of the Berkovich indenter was converted to a two-dimensional conical indenter with a 70.3° half-apex angle. This geometric transformation maintains a constant projected area-to-depth ratio, ensuring equivalent load-displacement responses.

[0051] During the loading phase, since the concept of time step is built into PFC, the loading speed in the simulation should be much higher than the actual experimental loading speed to shorten the time required for a single simulation without affecting the authenticity of the results.

[0052] Step 1.3: Use the force servo control dynamic loading algorithm to perform dynamic loading. The specific process is as follows:

[0053] Establish time-varying target load function F target1 =K×t. K is the proportional coefficient, representing the rate of load change. When the real-time load F actual <F target When the load is accelerated, the loading rate v1=K×(F target -F actual ); when F actual ≥F target When the loading rate is immediately reset to zero, the loading time parameter is set to obtain the final load F final .

[0054] Reach F final After that, maintain the target load F target =F final Perform closed-loop force control and set load holding time parameters.

[0055] Build time-varying uninstall function F target2 =F final -K×t. When F actual >F target2 When the uplink unloading is accelerated, the unloading rate v2=K×(F actual -F target2 ); when F actual ≤F target2 The unloading time and loading time are designed to be symmetrical, and the final load converges to zero.

[0056] The target load and real-time displacement data are recorded synchronously to construct the load-displacement characteristic curve.

[0057] Here, a force servo-controlled dynamic loading algorithm is used to strictly follow the one-way control logic during the loading and unloading process. During the loading process, when the real-time load is lower than the target load, the load is forward loaded, and when it is higher than the target load, the speed becomes zero, and the unloading process is reversed. The reverse loading disturbance is eliminated, and the mechanical path of the nanoindentation experiment is accurately reproduced; the Y-axis data source is switched from the real-time load to the target load, which effectively suppresses the signal oscillation caused by the elastic collision of discrete particles, and is more in line with the data output logic of the nanoindenter. Compared with the traditional loading mode, the advantages are: the uniform speed loading mode does not conform to the force control loading logic of the nanoindenter; in the ordinary force servo control mode, when the sample particles are separated from the cementation constraint, the bidirectional movement of the indenter causes repeated collisions between the detached particles and the indenter, resulting in high-frequency fluctuations in the measured load. The force servo-controlled dynamic loading algorithm in the present invention adopts a one-way force servo control strategy in the loading / unloading process to avoid reverse motion disturbances, and replaces the Y-axis data source with the target load value instead of the real-time load measurement value to avoid interference from discrete particle collisions and successfully eliminate abnormal fluctuations. Numerical verification based on quartz specimens shows that the improved loading mode successfully reproduces the quasi-static load response characteristics of indoor nanoindentation experiments. The choice of loading mode has a decisive influence on the PFC simulated nanoindentation experiment.

[0058] Step 1.4: Adjust the main microscopic parameters such as particle size, elastic modulus, and cementation strength through trial and error iteration. Through the parameter inverse identification method, the relative error between the simulated load-displacement curve and the experimental data is ≤10%. At the same time, the crack propagation pattern is highly consistent with the physical experiment, verifying the validity of the model and obtaining the microscopic parameters of the single mineral homogeneous model. Figure 6 shown.

[0059] See further Figure 4 The specific process of step 2 for constructing a multi-mineral heterogeneous model is as follows:

[0060] Step 2.1: Based on the combined XRD-SEM-EDS characterization data, use AutoCAD to simplify the two-dimensional spatial distribution structure of the mineral components. XRD reveals that the shale sample contains three major minerals: quartz, calcite, and clay. The distribution of Al, Si, and Ca in the EDS analysis results can be used to briefly delineate the distribution areas of clay, quartz, and calcite. Combined with the SEM results, a simplified mineral composition diagram can be generated using AutoCAD.

[0061] Here, XRD-SEM-EDS joint characterization data is a mineralogical data set obtained collaboratively through three analytical techniques, including XRD (X-ray diffraction), SEM (scanning electron microscopy), and EDS (energy dispersive spectrometry). It aims to achieve a comprehensive characterization of mineral composition, morphology, and elemental distribution, breaking through the limitations of a single technology and revealing mineral heterogeneity (such as the distribution patterns of clay, quartz, and pyrite in shale) through data fusion.

[0062] Step 2.2: Construct five types of characteristic models: clay macropores, crack defects, block distribution, layered distribution, and mosaic distribution models.

[0063] In order for the model to reflect the heterogeneity characteristics of shale, the following modeling guidelines must be followed:

[0064] First, since the longitudinal and transverse distribution characteristics of shale materials are similar at the same point, the XY plane cross-section SEM data can be used as the basis for longitudinal cross-section modeling; second, since the distance between indentation points needs to be greater than three times the indentation radius to eliminate the mutual influence of indentation points, the distance between the indentation center and the boundary should be ≥5μm to eliminate the boundary effect; third, the selected area needs to contain ≥2 mineral phases and consider the pore distribution to meet the research needs of multi-mineral heterogeneity of shale.

[0065] According to the above modeling criteria and the simplified diagram of mineral components obtained by XRD-SEM-EDS results, five types of characteristic models can be selected: clay macropores, crack defects, massive distribution, layered distribution, and mosaic distribution models.

[0066] Step 2.3: In order to overcome the particle size homogenization defect of the traditional GBM method in micron-scale modeling, the present invention proposes a multiphase separation-recombination modeling method to construct a multi-mineral heterogeneous model containing voids, which mainly includes phase construction, particle mosaic and weak cementation treatment, such as Figure 5 shown.

[0067] The specific method is:

[0068] According to the mineral composition of the sample, single mineral rectangular models of quartz, calcite and clay were constructed by phase separation. The particles in each single mineral distribution area were retained, and the particles in the redundant areas were deleted. Multi-phase particles were mosaicked and merged to obtain a shale nanoindentation multi-mineral heterogeneous model. All parameters of quartz, calcite and clay are controllable, showing multi-scale particle size distribution characteristics.

[0069] A servo control mechanism was used to adjust the displacement of the confining pressure baffle to achieve initial stress equilibrium in the specimen. The equilibrium stress threshold was determined based on the particle elastic modulus. This process allowed for a small amount of elastic deformation in the specimen, which was verified to have no impact on subsequent experimental results.

[0070] The LPBM model is used uniformly for intracrystalline and intercrystalline cementation, and the intercrystalline cementation strength is appropriately reduced to simulate the weak cementation characteristics between different minerals. The linear contact model is used for mineral boundary contact, and the pore system is selectively deleted by regional particles after cementation is completed.

[0071] The upper baffle was removed, and the remaining three baffles were fixed. After the specimen was balanced, a two-dimensional model of the Berkovich indenter was constructed. The three-dimensional structure of the Berkovich indenter was converted to a two-dimensional conical indenter with a 70.3° half-apex angle. This geometric transformation maintained a constant projected area-to-depth ratio, ensuring equivalent load-displacement responses.

[0072] Step 2.4: Based on the multi-mineral heterogeneous model obtained in the previous step, the microscopic parameters of each single mineral obtained in step 1 are imported. The indentation position can be controlled by moving the indenter left and right or deleting particles in the surface area. The force servo control mode is used in the loading stage. The results of the multi-mineral heterogeneous model with an 8mN load are as follows: Figure 7 As shown, the upper figure represents the distribution of cracks and particle displacements caused by the model being compressed, and the lower figure represents the distribution of the pressure chain.

[0073] Step 3: Monitor and obtain data during the loading process. Based on the crack number development, force chain distribution, and particle displacement vector characteristics, quantitatively analyze the regulation mechanism of the mineral interface coupling effect and pore stress concentration on the crack propagation path.

[0074] The data that needs to be monitored include: the load-displacement curve of the model, the number of different types of cracks, their occurrence locations, inclination angles, and frequencies, the size and evolution of the force chain, the particle displacement vector, etc. The specific process is as follows:

[0075] Monitor the target load and indenter displacement in real time and construct a load-displacement characteristic curve.

[0076] According to different contact types, the servo function is used to classify the cracks and monitor their respective quantities, including mineral intracrystalline cracks and intercrystalline cracks, tensile cracks and shear cracks.

[0077] According to the occurrence of different types of cracks, the servo function is used to statistically monitor the number of various types of cracks generated in this short period of time at intervals from the start to the end of loading to obtain the frequency of occurrence of different types of cracks.

[0078] By utilizing PFC 2D’s built-in monitoring functions such as ball, fracture dfn name, fracture rosette, contact force, and ball displacement, the crack location, crack inclination, force chain size and evolution, and particle displacement vector can be monitored in real time.

[0079] Based on the above embodiments, it can be seen that the present invention constructs a multi-mineral heterogeneous model considering pore defects through cross-scale coupling analysis of nanoindentation experiments and discrete element simulations, combined with XRD-SEM-EDS multimodal mineral characterization data, and uses force servo control dynamic loading algorithm to reveal the regulatory mechanism of rock micro-heterogeneity on crack evolution; it breaks through the limitations of traditional experimental two-dimensional observation and the bottleneck of nanoindentation numerical simulation, and can systematically analyze the heterogeneous evolution law of the force chain in the full cycle of loading-holding-unloading; establishes a quantitative mapping relationship between the mechanical properties of mineral components and the crack propagation path, and clarifies the synergistic damage mechanism of pore structure-induced stress concentration and mineral interface cementation strength; realizes accurate simulation of the three-dimensional dynamic evolution of the multiphase medium crack network, and reveals the energy conversion mechanism of elastic strain energy release driving the initiation of circumferential secondary cracks; accurately reproduces the differential deformation characteristics of quartz elasticity-dominated and clay plasticity-dominated, so that the simulation results are highly consistent with the actual nanoindentation experimental conditions, providing a theoretical basis and method support for the optimization of unconventional natural gas reservoir fracturing parameters and discrete element modeling of heterogeneous materials.

Claims

1. A numerical simulation method of rock nanoindentation based on PFC discrete element, characterized by: include: Step 1: Based on the XRD-SEM-EDS mineral spatial distribution of the rock sample and the single-mineral nanoindentation experimental results, the microscopic parameters of the single-mineral homogeneous model were calibrated. The loading was performed using a force servo-controlled dynamic loading algorithm. The load-displacement curve was matched through a trial-and-error iterative method to achieve a high degree of agreement between the simulated and experimental load-displacement curves. Step 2: Using the calibrated microscopic parameters and combined with XRD-SEM-EDS mineral spatial distribution data, multi-mineral phase modeling is performed to construct a multi-mineral heterogeneous model, and dynamic loading simulation is performed using a force servo-controlled dynamic loading algorithm; Step 3: Monitor the loading process data and quantitatively analyze the regulation mechanism of the mineral interface coupling effect and pore stress concentration on the crack propagation path based on the crack number development, force chain distribution and particle displacement vector characteristics.

2. The rock nanoindentation numerical simulation method based on PFC discrete element according to claim 1 is characterized in that: The step 1 comprises: Step 1.1: Based on the combined characterization results of mineral spatial distribution using XRD-SEM-EDS, the indentation point data of large-particle minerals in the nanoindentation experiment are selected as the basis for calibrating the microscopic parameters of the single mineral homogeneous model, and the elastic modulus and hardness of the single mineral are calibrated; Step 1.2, establish a PFC 2D single mineral homogeneous model; Step 1.3: Loading is performed using a force servo-controlled dynamic loading algorithm. The loading and unloading processes follow a one-way control logic to eliminate reverse loading disturbances. The Y-axis data source is switched from the real-time load to the target load, and the target load is recorded synchronously to reproduce the quasi-static response of the test. In step 1.4, the microscopic parameters are adjusted by trial-and-error iteration until the error between the simulated and experimental load-displacement curves is less than the set percentage, so that the simulated and experimental indentation failure characteristics are consistent, and the microscopic parameters of the single mineral homogeneous model are obtained.

3. The rock nanoindentation numerical simulation method based on PFC discrete element according to claim 1 is characterized in that: In step 1.1, the elastic modulus and hardness of the rock sample are calculated based on the load-displacement curve. The specific equations are as follows: Where H is the hardness, Pmax is the maximum load, Ac is the contact area, Er is the reduction modulus, S = dP / dh is the contact stiffness, ɑ is the shape factor of the Bohler indenter, E and v are the elastic modulus and Poisson's ratio of the rock specimen, respectively, and Ei and vi are the elastic modulus and Poisson's ratio of the indenter, respectively.

4. The rock nanoindentation numerical simulation method based on PFC discrete element according to claim 2 is characterized by: The construction process of step 1.2 includes: A two-dimensional discrete element model of the rock sample was established, and the sample aspect ratio was set to 10 μm × 5 μm; A servo control mechanism is used to adjust the displacement of the confining pressure baffle to achieve initial stress balance of the specimen; Construct a linear parallel cementation system between particles; Delete the upper baffle, fix the remaining three side baffles, and after the sample is balanced, perform two-dimensional modeling of the glass indenter.

5. The rock nanoindentation numerical simulation method based on PFC discrete element according to claim 3 is characterized by: In step 1.2, a linear parallel bond model is used to construct a linear parallel bond system, and the calculation formulas for its total force and total moment are as follows: Among them, F c is the total force, F l is the linear force, F d is the damping force, is the bonding force on the interface. c is the total torque, is the torque.

6. The rock nanoindentation numerical simulation method based on PFC discrete element method according to any one of claims 1 to 5, characterized in that: The step 2 includes: Step 2.1, based on the XRD-SEM-EDS joint characterization data, perform geometric topological simplification on the spatial distribution structure of the mineral components; Step 2.2: Construct five types of heterogeneous characteristic models, including clay macropores, crack defects, block distribution, layered distribution, and mosaic distribution models; Step 2.3: Use the multiphase separation-recombination modeling method to establish a multi-mineral heterogeneous discrete model containing pores; Step 2.4: Import the microscopic parameters obtained by calibration in step 1, use the force servo to control the indenter to perform dynamic loading, and obtain the simulation experimental results.

7. The rock nanoindentation numerical simulation method based on PFC discrete element according to claim 6 is characterized by: The model of step 2.2 follows three criteria: based on the principle of orthotropic material, XY plane cross-section data is selected; the distance between the indentation center and the boundary is ≥5μm to avoid boundary effects; The selected area must contain ≥2 mineral phases and consider the pore distribution.

8. The rock nanoindentation numerical simulation method based on PFC discrete element according to claim 6 is characterized by: The step 2.3 specifically includes: The multiphase separation-recombination modeling method is used to construct a single mineral model by phase separation, and a multi-mineral heterogeneous void model is generated by particle mosaic merging to achieve multi-scale particle size distribution parameter control; Based on the particle elastic modulus threshold, a servo control mechanism is used to adjust the displacement of the confining pressure baffle to achieve initial stress balance of the specimen. Unify the LPBM models for intracrystalline and intercrystalline cementation, appropriately reduce the intercrystalline cementation strength, use a linear contact model for mineral boundary contacts, and generate pore systems through selective particle deletion. Delete the upper baffle, fix the remaining three side baffles, and after the sample is balanced, perform two-dimensional modeling of the glass indenter.

9. The rock nanoindentation numerical simulation method based on PFC discrete element method according to any one of claims 1 to 5, characterized in that: The monitoring of step 3 includes: monitoring the load-displacement curve of the model during loading, monitoring the number, occurrence location, occurrence angle, and occurrence frequency of different types of cracks, monitoring the size and evolution law of the force chain, and the particle displacement vector.

10. The rock nanoindentation numerical simulation method based on PFC discrete element according to claim 9 is characterized in that: During the monitoring of step 3: Monitor the target load and indenter displacement in real time and construct a load-displacement characteristic curve. According to different contact types, the servo function is used to classify the cracks and monitor their respective quantities, including mineral intracrystalline cracks and intercrystalline cracks, tensile cracks and shear cracks. According to the occurrence of different types of cracks, the servo function is used to statistically monitor the number of various types of cracks generated in this short period of time at intervals from the start to the end of loading to obtain the frequency of occurrence of different types of cracks. The monitoring function of PFC 2D is used to monitor the crack location, crack inclination, force chain size and evolution, particle displacement vector, etc. in real time.

Images