Microscopic numerical simulation method based on mineral crystal model

Abstract

The invention discloses a microcosmic numerical simulation method based on a mineral particle model, which comprises the following steps of: 1, determining the geometric size of a test sample by taking the geometric parameters of an indoor test rock sample as a reference; 2, classifying Voronoi diagram polygons in a geometric model region according to grid structure characteristics of mineral crystals, and coarsely grouping particles; 3, adopting a parallel bonding contact model for the grouped Voronoi polygonal crystals; 4, applying an acting force to the model by moving the upper and lower wall bodies in a loading mode; and 5, continuously debugging mesoscopic parameters based on a trial-and-error method, finally obtaining simulation results of a stress-strain curve, a fracture mode, a cross-crystal tensile crack, a cross-crystal shear crack, an intergranular tensile crack, an intergranular shear crack and the like based on the mineral crystal model numerical sample, and obtaining a simulation test conclusion. According to the method, an effective theoretical basis is provided for exploring a rock fracture failure mechanism, and the expansion process of the cracks in the rock is accurately predicted.

Description

Technical Field

[0001] The invention relates to the technical field of particle flow microscopic numerical simulation experiments, and in particular to a microscopic numerical simulation method based on a mineral crystal model. Background Art

[0002] Natural rock has undergone complex geological structures and mineral formation processes and is a typical heterogeneous material. A large number of indoor microscopic studies have shown that the microscopic structures such as the types of minerals, arrangement and bonding inside the rock jointly determine the mechanical properties of the rock. In addition, through microscopic studies, it is found that intergranular cracks occur at the mineral boundaries and intragranular cracks that penetrate the minerals. However, the results of indoor tests are often obtained through inversion characterization after the test, and it is impossible to obtain the stress field and microscopic crack evolution laws of the rock during loading in real time. With the rapid development of computer technology, numerical calculations have received more and more attention in many engineering fields. Numerical calculation methods have the advantages of strong repeatability, small discreteness of results and low test cost. In addition, they have obvious advantages over in-situ and indoor tests in the creation, calculation and analysis of complex physical models. Therefore, numerical simulation methods have received more and more attention in engineering fields such as mining and geotechnical engineering.

[0003] The particle flow PFC method is used to describe the mechanical behavior of assembled disk or spherical particles. Based on the discrete element method, the microscopic mechanical behavior of the material is characterized according to the interaction and movement between particles, and then the complex macroscopic phenomena and mechanical mechanisms are explained from the microscopic perspective. Most researchers who study the PFC particle flow model use particle matrix units with uniform properties, and the microscopic mechanical parameters of all particle units are also the same. Therefore, the deformation behavior and fracture mechanism of the entire model in terms of microscopic mechanical behavior are relatively consistent. However, the composition of actual rocks is not uniform, but is composed of different types of minerals, which leads to different cementation between different minerals. According to the literature review, most researchers have established mineral particle heterogeneity numerical models mainly for complete specimens. However, there are still few studies on the construction of discrete element microscopic numerical models based on mineral crystals to quantitatively reveal transgranular tensile cracks, transgranular shear cracks, intergranular tensile cracks, and intergranular shear cracks in the process of rock deformation and fracture under different loads.

[0004] In addition, due to the complexity of the internal structure of rocks, the current strength criteria have great difficulties in accurately predicting the propagation process of cracks inside rocks. In order to study the complex crack characteristics of the internal structure of rocks, new technical solutions are needed. Therefore, how to construct a discrete element mesoscopic numerical model that takes into account the heterogeneity of minerals has become a technical problem that needs to be solved urgently. Summary of the invention

[0005] The present invention provides a mesoscopic numerical simulation method based on a mineral crystal model, which adopts a parallel bonding model and a smooth joint model to reveal the microscopic crack characteristics and fracture failure mechanism of a granite sample during the fracture process at a mesoscopic scale, thereby achieving accurate prediction of the expansion process of cracks inside the rock, and solving the difficult problems in the background technology.

[0006] The technical solution of the present invention is as follows: A microscopic numerical simulation method based on a mineral crystal model comprises the following steps:

[0007] Step 1: Based on the geometric parameters of the indoor test rock sample, a two-dimensional discrete element numerical model consistent with the geometric dimensions of the indoor test rock sample is constructed. The two-dimensional discrete element numerical model includes the following parameters:

[0008] The two-dimensional discrete element numerical model includes a parallel bond unit, and the parallel bond unit is calculated using the following formula:

[0009] F c =F l +F d +F b

[0010] M b =M c

[0011] The F l is the linear force, F d is the damping force, F b is the parallel bond force, M b is the parallel bond torque;

[0012] The linear force F l Decomposed into normal and tangential forces:

[0013]

[0014] The F n l is the normal force, F s l is the tangential force, n i is the unit vector of the normal contact surface, n t is the unit vector of the tangential contact surface;

[0015] The parallel key force F b and the parallel key unit moment M b The calculation formula is as follows:

[0016]

[0017]

[0018] The F n b and F s b are normal and tangential force keys, M n b and M s b are torque and bending moment respectively;

[0019] Step 2: Generate a PFC sample, generate boundary walls, generate small radius particles, increase particle radius, and the sample automatically balances;

[0020] Step 3: Classify the Voronoi diagram polygons in the internal area of ​​the boundary wall of the PFC specimen, and the boundary of the Voronoi polygon adopts a smooth joint contact model; the smooth joint direction is not collinear with the contact normal, and the relative displacement increment is decomposed into the joint normal direction and the tangential direction. The parallel bond force and parallel bond unit moment of the parallel bond unit of the Voronoi polygon are calculated based on the following formula:

[0021] Δδ=Δδ n n j +Δδ s

[0022] Δδ n and Δδ s are normal displacement and tangential displacement respectively, Δδ is the relative displacement increment;

[0023] The normal displacement Δδ n and tangential displacement Δδ s Calculate according to the following formula:

[0024] δ n =δ n -Δδ n

[0025] δ s =δ s +Δδ s

[0026] δ n and δ s are the accumulated normal displacement and accumulated tangential displacement respectively;

[0027] The force-displacement relationship of the smooth joint contact model is:

[0028] F c =F,M c ≡0

[0029] The F is the smooth joint contact force;

[0030] The smooth joint contact force F is divided into smooth joint contact force F n and tangential smooth joint contact force F s :

[0031] F=-F n n j +F s

[0032] Step 4: Based on the "trial and error method", the mesoscopic parameters are continuously adjusted, and the obtained simulation results are compared with the experimental results. Finally, a set of mesoscopic mechanical parameters that are close to the indoor test results are formulated. Subsequently, a series of simulation tests are carried out based on this set of correction parameters. Finally, the stress-strain curve, fracture mode, transgranular tensile cracks, transgranular shear cracks, intergranular tensile cracks, intergranular shear cracks and other simulation results of the numerical specimen based on the mineral crystal model are obtained, and the conclusions of the simulation test are drawn.

[0033] Further: it also includes a calculation step of measuring a circle, and the calculation formula of the measuring circle is as follows:

[0034]

[0035] In the formula, σ ij is the average force between crystals; N p is the crystal center of mass; N c is the contact between crystals; V p is the grain area; n is the porosity.

[0036]

[0037] In the formula, σ max is the maximum principal stress; σ x , σ y are the normal stresses in the x and y directions respectively; τ xy is the shear stress in the x direction.

[0038] The beneficial effects of this program:

[0039] 1. This program completed the computer analysis of the digital model of rock samples based on the particle flow method, and compared the analysis results with the experimental data of indoor tests, quantitatively revealing the evolution of transgranular tensile cracks, transgranular shear cracks, intergranular tensile cracks, and intergranular shear cracks in the process of rock fracture from a microscopic scale, and explained and supplemented the fracture mechanism at the macroscopic scale of indoor experiments, thereby providing an effective theoretical basis for exploring the rock fracture failure mechanism and realizing the accurate prediction of the expansion process of internal cracks in rocks.

[0040] 2. Based on the constructed mineral crystal model, this scheme creatively provides a method based on the measuring circle method to calculate the maximum principal stress and shear stress of all particles in the sample, and obtains the distribution law of the maximum principal stress field and shear stress field of the rock under different loading stages, so as to further reveal the deformation and fracture mechanism of the rock. BRIEF DESCRIPTION OF THE DRAWINGS

[0041] Figure 1 The present invention is a flowchart of a specific implementation of the microscopic numerical simulation method. DETAILED DESCRIPTION

[0042] The present invention will be further described below in conjunction with the accompanying drawings and embodiments:

[0043] like Figure 1 As shown, the present invention discloses a microscopic numerical simulation method based on a mineral crystal model, comprising the following steps:

[0044] 1. Establishing numerical model with particle flow method

[0045] The contact models of the PFC particle flow method mainly include the contact bonding model and the parallel bonding model. The characteristic of the contact bonding model is that the particles can only transmit forces but not torques. For the parallel bonding model, the contact of the particle bond runs through the entire unit body. No matter where it breaks, the entire parallel bond will be completely removed. Therefore, the parallel bonding model is used for the crystal inside the Voronoi polygon in this step.

[0046] The parallel bonding model is mainly composed of linear units and parallel bonds. The linear units only allow the transmission of elastic effects between particles, and cannot provide stretching and rotational effects. The parallel bond unit can transmit forces and torques between particles until the relative motion force between them exceeds the bond strength of the bond. Therefore, this scheme is based on the parallel bond unit for research.

[0047] F c =F l +F d +F b

[0048] M b =M c

[0049] In the formula, F l is the linear force, F d is the damping force, F b is the parallel bond force, M b is the parallel bond torque.

[0050] Linear force F l It can be decomposed into normal and tangential forces:

[0051]

[0052] In the formula, F n l is the normal force, F s l is the tangential force, n i is the unit vector of the normal contact surface, n t is the unit vector of the tangential contact surface.

[0053] The calculation formulas for the parallel key force and parallel key unit moment are as follows:

[0054]

[0055]

[0056] In the formula, F n b and F s b are normal and tangential force keys, M n b and M s b are torque and bending moment respectively.

[0057] 2. The Voronoi polygon boundary adopts a smooth joint contact model; the smooth joint direction is not collinear with the contact normal, and the relative displacement increment is decomposed into the joint normal direction and the tangential direction. The parallel bond force and parallel bond unit moment of the parallel bond unit are calculated based on the following formula:

[0058] Δδ=Δδ n n j +Δδ s

[0059] In the formula, Δδ n and Δδ s are the normal displacement and tangential displacement respectively, and Δδ is the relative displacement increment.

[0060] The normal displacement and shear displacement are calculated as follows:

[0061] δ n =δ n -Δδ n

[0062] δ s =δ s +Δδ s

[0063] In the formula, δ n and δ s are the accumulated normal displacement and the accumulated tangential displacement, respectively.

[0064] The force-displacement relationship of the smooth joint contact model is:

[0065] F c =F,M c ≡0

[0066] Where F is the contact force of smooth joints;

[0067] The contact force of smooth joints is divided into normal and tangential forces:

[0068] F=-F n n j +F s

[0069] In the formula, F n and F s are the normal smooth joint contact force and the tangential smooth joint contact force;

[0070] 3. Based on the geometric parameters of the rock samples tested indoors, a two-dimensional discrete element numerical model consistent with the geometric dimensions of the rock samples tested indoors is constructed through PFC. In addition, the spatial position and geometric dimensions of the prefabricated cracks in the model are consistent with those of the sandstone samples tested indoors. The number and size of particles have a great influence on the macroscopic mechanical parameters of the model. Considering the accuracy and efficiency of the calculation results, the minimum particle size in the two-dimensional discrete element numerical model is set to 0.2mm and the maximum particle size is 0.4mm, which is close to the average particle size of sandstone measured indoors of 0.15mm. In the specific implementation of this plan, three types of granite samples are constructed, namely complete samples, non-filled samples and gypsum-filled samples.

[0071] 4. Loading method: Generate PFC specimens, generate walls in sequence, generate small radius particles, increase particle radius, and the specimens are automatically balanced;

[0072] In order to prevent the stress-strain curve from fluctuating greatly and to ensure that the entire force is in a quasi-static loading state, the loading rate of the numerical experiment is set to 0.1 mm / s, and the calculation termination condition is 0.4σ after the peak. c . It should be noted that the prefabricated cracks in the specimen are achieved by deleting particles;

[0073] The loading method is to apply force to the model by moving the upper and lower walls, and use the built-in Fish language to monitor the mechanical parameters during the entire loading process.

[0074] To ensure that the specimen does not move in the axial direction during the lateral pressure application, an initial force is first applied in the axial direction (upper and lower walls). Then, the lateral force is applied until the lateral pressure increases to the target value. Next, the lateral pressure is kept constant and a constant speed of 0.05 mm / s is applied to the upper and lower walls in the axial direction until the stress drops to 40% after the peak, and the loading is terminated.

[0075] 5. Calibration of micro-model parameters

[0076] The initial microscopic parameters were input, and uniaxial compression tests were first carried out to calibrate the mesoscopic Young's modulus and Poisson's ratio, and then uniaxial tension tests were carried out to calibrate the mesoscopic Young's modulus and Poisson's ratio, respectively. Next, the mesoscopic parameters were continuously debugged based on the "trial and error method", and the obtained simulation results were compared with the experimental results. Finally, a set of mesoscopic mechanical parameters that were close to the indoor test results were proposed. Subsequently, a series of simulation tests were carried out using this set of calibrated mesoscopic parameters, and finally the stress-strain curves of the intact sandstone under different loading conditions and the final failure and fracture mode results of the specimens were obtained.

[0077] Comparing the numerical simulation and test results, we can see that although the simulation results entered the elastic stage from the beginning of loading, that is, there was no pore crack compression deformation stage of the indoor sample results, it was mainly because the pores and particle skeletons of the sample were evenly distributed under the initial hydrostatic pressure, and the particles were already in a stress equilibrium state before the axial force was applied. However, from the mechanical strength and failure mode of the two, we can see that except for the initial loading stage, other mechanical parameters such as peak stress and strain corresponding to peak stress are basically the same, which shows the rationality of the selection of micromechanical parameters, and further shows that the numerical simulation results have excellent reference value.

[0078] Special note: The evolution characteristics of the full-field stress of the sample cannot be directly obtained from the existing testing methods. It is well known that the stress field is crucial to revealing the failure mechanism of rock fracture, and the existing PFC particle flow module cannot directly obtain the evolution characteristics of the maximum principal stress field and the shear stress field. In order to further explore the evolution law of the stress field during rock fracture, the scheme uses the method of arranging measurement circles for the entire sample to reveal the evolution law of the stress field at different loading stages. The calculation formulas for the average stress, maximum principal stress and maximum shear stress of the particles in the measurement circle are shown below.

[0079]

[0080] In the formula, σ ij is the average force between particles; N p is the particle mass center; N c is the particle contact; V p is the particle area; n is the porosity.

[0081]

[0082] In the formula, σ max is the maximum principal stress; σ x , σ y are the normal stresses in the x and y directions respectively; τ xy is the shear stress in the x direction.

[0083] In addition, the existing indoor monitoring technology has carried out detailed research on the deformation behavior and crack rupture process of rock during loading. However, for the fractured specimens, the stress distribution around the fracture and the influence of the filling on the stress evolution around the fracture cannot be directly obtained from the existing testing methods. Based on the measuring circle method, the maximum shear stress between particles in the typical non-filled fractured sandstone specimen is calculated and processed to obtain the evolution law of the maximum shear stress field under different stress stages. Consistent with the variation law of the maximum principal stress field, high stress accumulation is generated in the rock bridge and fracture tip areas, while low stress accumulation areas are nurtured around the fractures. In general, under different fracture geometry conditions, the maximum shear stress fields of the three typical specimens are both similar and different. In common, high shear stress accumulation areas and low accumulation areas appear alternately, and the high shear stress accumulation areas are mainly distributed in the rock bridge area and the fracture tip. On the contrary, the low shear stress areas are mainly distributed around the prefabricated fractures. The difference is that the range of the lower shear stress accumulation zone decreases with the increase of the crack inclination angle, and it can be seen from different loading stages that the stress peak in the low shear stress accumulation zone gradually increases. As the crack inclination angle increases, the peak shear stress around the prefabricated crack, the rock bridge area and the crack tip at different stress stages gradually increase.

[0084] Therefore, based on the measuring circle method, the maximum principal stress of all particles in the sample is calculated and processed to obtain the distribution law of the maximum principal stress field of the rock sample under different loading stages.

[0085] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the protection scope of the present invention.

Claims

1. A mesoscopic numerical simulation method based on a mineral crystal model, characterized in that: The steps include: Step 1: Based on the geometric parameters of the indoor test rock sample, a two-dimensional discrete element numerical model consistent with the geometric dimensions of the indoor test rock sample is constructed. The two-dimensional discrete element numerical model includes the following parameters: The two-dimensional discrete element numerical model includes a parallel bond unit, and the parallel bond unit is calculated using the following formula: F c =F l +F d +F b M b =M c The F l is the linear force, F d is the damping force, F b is the parallel bond force, M b is the parallel bond torque; The linear force F l Decomposed into normal and tangential forces: The F n l is the normal force, F s l is the tangential force, n i is the unit vector of the normal contact surface, n t is the unit vector of the tangential contact surface; The parallel key force F b and the parallel key unit moment M b The calculation formula is as follows: The F n b and F s b are normal and tangential force keys, M n b and M s b They are torque and bending moment respectively; Step 2: Generate PFC specimen, generate boundary wall in sequence, generate small radius particles, increase particle radius, and the specimen is automatically balanced; Step 3: Classify the Voronoi diagram polygons in the internal area of ​​the boundary wall of the PFC specimen, and the boundary of the Voronoi polygon adopts a smooth joint contact model; the smooth joint direction is not collinear with the contact normal, and the relative displacement increment is decomposed into the joint normal direction and the tangential direction. The parallel bond force and parallel bond unit moment of the parallel bond unit of the Voronoi polygon are calculated based on the following formula: Δδ=Δδ n n j +Dd s Δδ n and Δδ s are normal displacement and tangential displacement respectively, Δδ is the relative displacement increment; The normal displacement Δδ n and tangential displacement Δδ s Calculate according to the following formula: d n =d n -Dd n d s =d s +Dd s δ n and δ s are the accumulated normal displacement and accumulated tangential displacement respectively; The force-displacement relationship of the smooth joint contact model is: F c =F,M c ≡0 The F is the smooth joint contact force; The smooth joint contact force F is divided into smooth joint contact force F n and tangential smooth joint contact force F s : F=-F n n j +F s Step 4: Based on the "trial and error method", the mesoscopic parameters are continuously adjusted, and the obtained simulation results are compared with the experimental results. Finally, a set of mesoscopic mechanical parameters that are close to the indoor test results are formulated. Subsequently, a series of simulation tests are carried out based on this set of correction parameters. Finally, the stress-strain curve, fracture mode, transgranular tensile cracks, transgranular shear cracks, intergranular tensile cracks, intergranular shear cracks and other simulation results of the numerical specimen based on the mineral crystal model are obtained, and the conclusions of the simulation test are drawn.

2. The microscopic numerical simulation method based on a mineral crystal model according to claim 1, characterized in that: The method also includes a step of calculating a measuring circle, wherein the calculation formula of the measuring circle is as follows: In the formula, σ ij is the average force between crystals; N p is the crystal center of mass; N c is the contact between crystals; V (p) is the grain area; n is the porosity; In the formula, σ max is the maximum principal stress; σ x , σ y are the normal stresses in the x and y directions respectively; τ xy is the shear stress in the y direction.

Images