Modeling and calculating method for crystalline rock with bedding structure based on Neper-DEM (Digital Elevation Model)

Abstract

The invention discloses a Neper-DEM-based apparent crystal rock modeling and calculation method with a bedding structure, relates to the field of numerical calculation of a tunnel engineering technology, and solves the problem that a simulation and calculation method for apparent crystal rocks with transverse isotropy such as gneiss and granite gneiss is lacked in the prior art. The method comprises the following steps: constructing a polycrystal structure with crystal characteristics in Neper, and generating a geometric parameter text of crystal grains from the polycrystal structure; constructing a discrete element rock sample model in the DEM; performing a numerical simulation test on the discrete element rock sample model to obtain a stress-strain curve and a sample failure mode, wherein the numerical simulation is realized by adopting a Brazilian splitting test and a triaxial test; according to a numerical simulation result, discrete element sample calculation parameters are calibrated, and the effectiveness of the modeling and calculation method of the discrete element rock sample model is verified; according to the modeling method, construction and assignment of the crystal structure and the bedding structure are achieved at the same time, so that the rock structure is closer to a real rock structure.

Description

Technical Field

[0001] The present invention relates to the field of numerical calculation of tunnel engineering technology, and in particular to a modeling and calculation method of crystalline rock with a bedding structure based on Neper-DEM. Background Art

[0002] As a natural material, rock will produce different degrees of anisotropy and heterogeneity, such as microcracks and defects, under long-term geological conditions. Rock anisotropy refers to the difference in the physical and mechanical properties of rock in different directions, which is mainly affected by the internal structure of the rock, such as bedding, foliation, cracks and mineral arrangement. For example, sedimentary rocks (such as shale and slate) are deposited in different periods and environments to form a bedding structure. Metamorphic rocks (such as diabase and slate) will also form anisotropic characteristics during the metamorphic process.

[0003] Directed pressure during metamorphism causes the mineral grains or layers in the rock to align in parallel, forming structural planes with directional properties. Due to this layered or folded structure, these rocks are considered transversely isotropic rocks, that is, these rocks have the same physical properties in two perpendicular directions, but show different properties in the third direction (usually the direction of bedding planes or major fractures). In addition to the anisotropy caused by the bedding structure, heterogeneous mineral composition also affects the mechanical properties of rocks. As a typical mineral structure rock, granite has a complex microstructure, such as mineral content, shape, size and spatial distribution, as well as defects in the form of voids, microcracks, weak grain boundaries and fractures, which control the generation, propagation and aggregation of cracks and fractures.

[0004] There are many rocks in nature that have both stratification and obvious polycrystalline structures. These rocks are usually metamorphic rocks, and some igneous rocks may also have similar characteristics. Gneiss in metamorphic rocks has obvious gneissic structure, that is, light-colored and dark-colored mineral components present banded and oriented stratification, and the mineral crystals in gneiss are usually large, especially quartz, feldspar and other minerals can form obvious crystals, and biotite and amphibole in the dark bands can also have large crystal particles. Schist also has obvious lamellae stratification, and the oriented arrangement of minerals makes it present a parallel layered structure. Its mineral crystals are usually large, especially mica, chlorite, garnet, etc., which can form obvious lamellae or scaly crystals. Granite, the most typical igneous rock, may also have transverse isotropy. For example, when there are more cracks distributed in a certain direction, the tensile and compressive strength of granite in the direction perpendicular to the cracks may be different from that along the cracks. Metamorphism, sedimentation or intrusive stratification can also form oriented mineral crystals, especially in granite gneiss. These phenomena will make granite show transverse isotropy in micro or macro structure. The crystal structure and bedding structure of rock have important influence on rock stability, stress distribution, deformation characteristics, crack propagation and failure mode. These rocks are also anisotropic and heterogeneous. In terms of deformation, they may usually show larger deformation in the direction perpendicular to the bedding or cracks, while the deformation in the parallel direction is smaller. Due to the presence of non-uniformly distributed minerals, large deformation may occur locally. In terms of failure mode, cracks are often more likely to expand along the direction of bedding, mineral directional arrangement, and weaker mineral structure interfaces, which may lead to non-uniform distribution of cracks and faults, affecting the stability of the rock mass.

[0005] Discrete element software can observe the deformation and failure characteristics of such rocks more clearly and directly, and study their crack propagation mechanism. The discrete element method can simulate the generation, expansion and merging of cracks at the particle scale, which is more advantageous than traditional continuous medium methods (such as finite element method) and is particularly suitable for studying the failure mechanism of rock materials. The discrete element method can also carefully consider the interaction between particles, including adhesion, friction, etc., so that the simulation results are closer to the real mechanical behavior of the material, and can accurately simulate the local failure process of the material, analyze the initiation and propagation path of micro cracks in the material, and reveal the failure mechanism of the material.

[0006] The existing sample modeling and calculation technology has the following shortcomings: (1) It can only simulate and calculate single bedding structures or rocks with crystal structures, but there is little research on the interaction mechanism and mechanical properties of the two; (2) The research on rock samples is concentrated on the two-dimensional plane, ignoring the three-dimensional ductility and joint crack tendency; (3) The discrete element calibration method is mainly based on trial and error, with huge calculation volume and low accuracy; (4) The grain size is relatively uniform and the distribution rules are quite different from the distribution of mineral crystals in rocks in nature. (5) There is little analysis on the crack generation and expansion process and the proportion of various types of cracks, which cannot reflect the development of rock failure process. Summary of the invention

[0007] In order to solve the problems existing in the above-mentioned prior art, the present invention provides a modeling and calculation method of phaneritic rocks with bedding structure based on Neper-DEM, which solves the problem that the prior art lacks simulation and calculation methods for phaneritic rocks with transverse isotropy such as gneiss and granite gneiss.

[0008] Modeling and calculation methods of crystalline rocks with bedding structure based on Neper-DEM, including:

[0009] Construct a polycrystalline structure with crystal characteristics in Neper, and then generate the polycrystalline structure into a geometric parameter text of the grains;

[0010] A discrete element rock sample model is constructed in DEM, including constructing discrete element samples and importing geometric parameter text as geometric information, grouping discrete element samples according to the geometric information, constructing multiple mutually parallel bedding planes and inter-group joint planes whose directions are consistent with the grain intersection plane in the discrete element samples, and then assigning contact models and parameters of three types of contact modes to the discrete element samples, the three types of contact modes including intra-group contact of discrete element samples and indirect contact of discrete element samples with joint planes;

[0011] A numerical simulation test is performed on a discrete element rock sample model to obtain a stress-strain curve and a sample failure mode, wherein the numerical simulation is implemented using a Brazilian splitting test and a triaxial test;

[0012] The calculation parameters of the discrete element sample are calibrated according to the numerical simulation results to verify the effectiveness of the modeling and calculation method of the discrete element rock sample model.

[0013] Furthermore, the polycrystalline structure creation process includes:

[0014] The lens is constructed using the 3D Thiessen polygon tessellation method based on the crystal characteristics, namely the number of grains, grain shape, and grain distribution. The number of grains is specified in the Neper command, the grain shape is reflected by the normal distribution of equivalent diameter and sphericity, and the grain distribution is controlled by random positions and seed generation weights. Under the influence of different seed weights, the boundaries between grains become non-equidistant between seeds. The greater the seed weight, the larger the associated grain. After the Thiessen polygon tessellation is completed, the various parameters of the lens are output as an STL file.

[0015] Furthermore, grouping discrete element samples according to geometric information includes:

[0016] Importing geometric information, the geometric information includes normal vectors of each face of the convex grain and node coordinates of each face, and multiple faces together constitute a grain entity; using DEM to read the geometric information, generate the same node domain surface in the DEM, and construct a single Geometry based on the convex grain composed of multiple faces and name it according to the import order;

[0017] The existing discrete element samples are grouped according to the named Geometry. When multiple particles fall within the Geometry polyhedron, they are identified as a separate group (Group), and its name is consistent with the Geometry.

[0018] Repeat the above steps until all the grains constituting the polycrystalline body are introduced and grouped.

[0019] Furthermore, the construction of a plurality of mutually parallel bedding planes and inter-group joint planes whose directions are consistent with the grain intersection planes comprises:

[0020] A certain number of parallel joint planes are added to the polycrystalline discrete element sample, and the bedding direction is controlled by the dip and inclination angle to be consistent with the original sample and the spacing between each bedding is kept consistent;

[0021] Construct inter-group joint surfaces according to the directions of the imported grain surfaces, and determine whether the two contacting grains are in the same group. If they are in the same group, do not process them. If they are not in the same group, read the group names respectively, and the group names correspond to the Geometry with the same name.

[0022] Find the two corresponding Geometries according to the name, determine the common face of the Geometries, obtain the direction of the face and assign the same direction to the two contacting particle joint surfaces.

[0023] Furthermore, the contact models and parameters for assigning the three types of contact modes include:

[0024] According to the established bedding planes and inter-group joint planes, the contact model and parameter assignment are performed on the discrete element model;

[0025] The bedding plane is given a smooth joint model, which has the same direction as the bedding plane and a small contact parameter to simulate a weak structural plane;

[0026] The inter-group joint surface is also given a smooth joint model to simulate the interface contact between crystals, and the contact parameter value is larger than that of the layer;

[0027] The particles that are not covered by inter-group joint planes and bedding planes are assigned a flat joint model.

[0028] Furthermore, the experimental process of the numerical simulation experiment is as follows:

[0029] Create wall units at the top and bottom of the discrete element rock specimen model, use the wall units to simulate the pressure head of the Brazilian splitting test, set the modulus and width of the wall units, and fix the wall loading speed;

[0030] When the wall starts to be loaded, the displacement of the wall compression, the force during the compression process, the crack distribution and number of the polycrystalline model are recorded. The data and pictures are automatically exported through the DEM command to obtain the stress-strain curve, elastic modulus, failure mode and failure strength.

[0031] Further, obtaining the stress-strain curve includes:

[0032] The axial force, axial strain and displacement of the wall are monitored in the discrete element simulation of the Brazilian splitting test. The maximum peak force is continuously updated during the monitoring process. If the newly generated force is less than the threshold of the peak force, the specimen is judged to be damaged.

[0033] At this point, the calculation is completed, and the axial force and axial strain of the wall are output. The axial force is converted into stress using the formula, and the stress-strain curve is plotted. The stress σ t The calculation formula is as follows:

[0034]

[0035] where σ t is the Brazilian splitting stress, P is the wall axial force, D is the specimen diameter, and t is the specimen thickness.

[0036] Furthermore, the obtaining of the stress-strain curve further includes: monitoring the axial force, axial strain and displacement of the wall in the discrete element simulation in the triaxial test, continuously updating the maximum peak force during the monitoring process, and determining that the sample is damaged when the current force is less than a threshold value of the peak force;

[0037] At this point, the calculation is completed, and the axial force and axial strain of the wall are output. The axial force is converted into stress using the formula, and the stress-strain curve is plotted. The stress σ c The calculation formula is as follows:

[0038]

[0039] where σ c is the axial stress, P 1 is the bottom wall pressure, P 2 is the top wall pressure, and D is the diameter of the cylindrical specimen.

[0040] Furthermore, the calibrating of discrete element sample calculation parameters based on numerical simulation results includes:

[0041] The deformation parameters of the discrete element model are calibrated by the triaxial test results. When the bedding direction is perpendicular to the loading direction, the bedding plane has the least effect on the elastic modulus of the specimen. At this time, the stiffness parameters of the smooth joint model and the flat joint model between grains are calibrated.

[0042] When the bedding direction coincides with the loading direction, the bedding plane stiffness has the greatest influence, and the stiffness parameters of the smooth joint model of the bedding plane are calibrated;

[0043] The strength parameters of the discrete element model are calibrated according to the Brazilian splitting test results. When the bedding plane direction is perpendicular to the loading direction, the strength parameters of the smooth joint model and the flat joint model are calibrated.

[0044] When the bedding direction coincides with the loading direction, the stiffness of the bedding plane has the greatest influence, and the strength parameters of the smooth joint model of the bedding plane are calibrated.

[0045] Furthermore, the process of obtaining the crack distribution of the polycrystalline model includes:

[0046] According to the different characteristics and failure modes of the three types of contact modes, the crack types are divided into 6 categories, namely, flat joint model tensile crack, flat joint shear crack, intergranular smooth joint model tensile crack, intergranular smooth joint model shear crack, bedding plane smooth joint model tensile crack, bedding plane smooth joint model shear crack;

[0047] During the test loading process, the value of the increase in the number of cracks as the strain increases is recorded, and a crack change picture is drawn to determine which type of damage accounts for the largest proportion in the test.

[0048] The beneficial effects of the present invention include:

[0049] The present invention improves the existing numerical simulation modeling methods for rocks such as gneiss, granite gneiss and schist, while taking into account the influence of the heterogeneity of the crystal structure and the transverse isotropy of the bedding structure. In the modeling method, the construction and assignment of the crystal structure and the bedding structure are simultaneously realized, making it closer to the real rock structure.

[0050] The present invention improves the polycrystal generation method in the existing discrete element software, introduces Neper software to construct a polycrystal model inlaid with three-dimensional Thiessen polygons, generates a more complex numerical model that conforms to the actual crystal distribution, and realizes the import and grouping of the crystal model generated in Neper in the DEM software, which effectively solves the problem that the crystal structure generation in the discrete element model is difficult and the crystals are single and uniform.

[0051] The present invention improves the problem of parameter calibration in existing discrete element software, and calibrates different parameters by considering the influence of bedding plane on rock strength respectively, thus solving the problem of time-consuming, low efficiency and inaccurate calibration of "trial and error method".

[0052] The present invention improves the problems of crack generation and statistics, realizes accurate grouping of cracks in discrete element software, effectively distinguishes tensile cracks and shear cracks generated in crystals, crystal interfaces and bedding planes, and monitors the development and proportion of the number of each crack, thereby better realizing the analysis and research of rock failure modes. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] Figure 1 This is a flow chart of the modeling and calculation method of crystalline rocks with bedding structure based on Neper-DEM involved in the embodiments of the present application.

[0054] Figure 2 Create grains for the Neper software Thiessen polygon tessellation method involved in this application.

[0055] Figure 3 The initial sample model constructed by the DEM software involved in this application.

[0056] Figure 4 The polycrystalline model of the polycrystalline structure in the Neper software involved in this application is imported into the DEM software.

[0057] Figure 5 A model with completed contact assignment for the discrete element model of the layered polycrystal involved in this application.

[0058] Figure 6 This is a flow chart of discrete element layered polycrystal parameter calibration involved in this application.

[0059] Figure 7 (a) is a comparison chart of the discrete element calibration strength and the test results when the inclination is 30° involved in this application.

[0060] Figure 7 (b) is a comparison chart of the discrete element calibration strength and the test results when the inclination is 60° involved in this application.

[0061] Figure 8 This is a comparison chart between the discrete element calibration failure mode and the test results involved in this application.

[0062] Fig. 9 (a) is a graph showing the variation of the failure strength with the inclination and the dip angle in the sensitivity analysis involved in this application.

[0063] Fig. 9 (b) is a diagram showing the change of failure mode with inclination and dip angle in the sensitivity analysis involved in this application.

[0064] Fig. 9 (c)-(i) are graphs showing the variation of crack proportion with inclination and dip angle in the sensitivity analysis involved in this application. DETAILED DESCRIPTION

[0065] In order to make the purpose, technical scheme and advantages of the embodiments of the present application clearer, the technical scheme in the embodiments of the present application will be clearly and completely described below in conjunction with the drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all of the embodiments. Therefore, the following detailed description of the embodiments of the present application provided in the drawings is not intended to limit the scope of the application for protection, but merely represents the selected embodiments of the present application. Based on the embodiments of the present application, all other embodiments obtained by those skilled in the art without making creative work belong to the scope of protection of the present application.

[0066] The following is combined with Figure 1-9 The specific embodiments of the present invention are described in detail;

[0067] Modeling and calculation methods of crystalline rocks with bedding structure based on Neper-DEM, such as Figure 1 As shown, including:

[0068] Construct a polycrystalline structure with crystal characteristics in Neper, and then generate the polycrystalline structure into a geometric parameter text of the grains;

[0069] A discrete element rock sample model is constructed in DEM, including constructing discrete element samples and importing geometric parameter text as geometric information, grouping discrete element samples according to the geometric information, constructing multiple mutually parallel bedding planes and inter-group joint planes whose directions are consistent with the grain intersection plane in the discrete element samples, and then assigning contact models and parameters of three types of contact modes to the discrete element samples, the three types of contact modes including intra-group contact of discrete element samples and indirect contact of discrete element samples with joint planes;

[0070] A discrete element rock sample model is numerically simulated to obtain a stress-strain curve and a sample failure mode, wherein the numerical simulation is implemented using a Brazilian splitting test and a triaxial test;

[0071] The calculation parameters of the discrete element sample are calibrated according to the numerical simulation results, the effectiveness of the modeling and calculation method of the discrete element rock sample model is verified, and this method is used to carry out sensitivity analysis.

[0072] In another embodiment, the polycrystalline structure creation process includes:

[0073] The lens is constructed using the 3D Thiessen polygon tessellation method based on the crystal characteristics, namely the number of grains, grain shape, and grain distribution. The number of grains is specified in the Neper command, the grain shape is reflected by the normal distribution of equivalent diameter and sphericity, and the grain distribution is controlled by random positions and seed generation weights. Under the influence of different seed weights, the boundaries between grains become non-equidistant between seeds. The greater the seed weight, the larger the associated grain. After the Thiessen polygon tessellation is completed, the various parameters of the lens are output as an STL file.

[0074] Specifically in this embodiment, the size and crystal characteristics of the model are determined, and a Thiessen polygon mesh is constructed in Neper to create a polycrystalline structure with the desired crystal characteristics. Neper generates specific graphics such as Figure 2 The constructed polycrystalline structure is generated as a text of the geometric parameters of the grains; a discrete element rock sample model is constructed in the DEM software, and the initial parameters of the discrete element model are adjusted to achieve the initial balance. The initial discrete element model is as follows Figure 3 As shown;

[0075] The crystal characteristics include three parameters: generation quantity, equivalent diameter and sphericity; the discrete element initial parameters include particle effective modulus, stiffness ratio, friction coefficient, particle size ratio, minimum particle size and particle density.

[0076] Some parameters of the Neper modeling and DEM modeling are as follows: (the parameters in the table are actual values ​​in the selected embodiment)

[0077] Table 1 Parameter definition table

[0078]

[0079]

[0080] In another embodiment, grouping discrete element samples according to geometric information comprises:

[0081] Importing geometric information, the geometric information includes normal vectors of each face of the convex grain and node coordinates of each face, and multiple faces together constitute a grain entity; using DEM to read the geometric information, generate the same node domain surface in the DEM, and construct a single Geometry based on the convex grain composed of multiple faces and name it according to the import order;

[0082] Existing discrete element samples are grouped according to the named Geometry. When multiple particles fall within the Geometry polyhedron, they are identified as a separate group (Group), and its naming is consistent with the Geometry.

[0083] Repeat the above steps until all the grains constituting the polycrystalline body are introduced and grouped.

[0084] Specifically in this embodiment, the STL file exported from Neper contains the normal vectors of each face of the convex grain and the node coordinates that make up each face. Multiple faces together constitute a grain entity, with a total of 2442 grains; the DEM software can read the command of the STL file, regenerate the same nodes and faces in the DEM, and construct a single Geometry based on the convex grains composed of multiple faces and name them according to the import order; the existing discrete element samples are grouped according to the named Geometry, and when multiple particles fall within the range of the Geometry polyhedron, they are identified as a separate group (Group), and their naming is consistent with the Geometry. Repeat the above steps until all 2442 grains that constitute the polycrystal are imported and divided into 2442 groups. The imported polycrystal discrete elements are as follows: Figure 4 shown.

[0085] In another embodiment, constructing a plurality of mutually parallel bedding planes and inter-group joint planes whose directions are consistent with the grain intersection planes comprises:

[0086] A certain number of parallel joint planes are added to the polycrystalline discrete element sample, and the bedding direction is controlled by the dip and inclination angle to be consistent with the original sample and the spacing between each bedding is kept consistent;

[0087] Construct inter-group joint surfaces according to the directions of the imported grain surfaces, and determine whether the two contacting grains are in the same group. If they are in the same group, do not process them. If they are not in the same group, read the group names respectively, and the group names correspond to the Geometry with the same name.

[0088] Find the two corresponding Geometries according to the name, determine the common face of the Geometries, obtain the direction of the face and assign the same direction to the two contacting particle joint surfaces.

[0089] Specifically, in this embodiment, 5 parallel bedding planes are added to the polycrystalline discrete element sample, and the bedding direction is controlled by the inclination and dip angle to be consistent with the original sample and the spacing between each bedding is consistent, and the spacing is controlled at 0.75mm; the inter-group joint plane is constructed according to the direction of the imported grain plane, and it is determined whether the two contacting particles are the same group. If they are the same group, they are not processed. If they are not the same group, the group names are read separately. It can be seen from the above that the group names correspond to the Geometries with the same names; the corresponding two Geometries are found according to the names, the common faces of the Geometries are determined, the direction of the face is obtained, and the joint faces of the two contacting particles are given the same direction.

[0090] In another embodiment, the contact models and parameters for assigning three types of contact modes include:

[0091] According to the established bedding planes and inter-group joint planes, the contact model and parameter assignment are performed on the discrete element model;

[0092] The bedding plane, i.e. the joint plane contact, is assigned a smooth joint model, which has the same direction as the bedding plane and a small contact parameter to simulate a weak structural plane;

[0093] The inter-group joint surface, i.e. the contact between discrete element specimens, is also given a smooth joint model to simulate the interface contact between crystals, and the contact parameter value is larger than that of the layer.

[0094] A flat joint model is assigned to particles that are not within the coverage of inter-group joint planes and bedding planes, i.e., particles that are in contact within discrete element sample groups.

[0095] Specifically, in this embodiment, the discrete element model after the three types of contact assignment are completed is as follows: Figure 5 shown.

[0096] In another embodiment, the test process of the numerical simulation test is:

[0097] Create wall units at the top and bottom of the discrete element rock specimen model, use the wall units to simulate the pressure head of the Brazilian splitting test, set the modulus and width of the wall units, and fix the wall loading speed;

[0098] When the wall starts to be loaded, the displacement of the wall compression, the force during the compression process, the crack distribution and number of the polycrystalline model are recorded. The data and pictures are automatically exported through the DEM command to obtain the stress-strain curve, elastic modulus, failure mode and failure strength.

[0099] In another embodiment, obtaining the stress-strain curve comprises:

[0100] The axial force, axial strain and displacement of the wall are monitored in the discrete element simulation of the Brazilian splitting test. The maximum peak force is continuously updated during the monitoring process. If the newly generated force is less than the threshold of the peak force, the specimen is judged to be damaged.

[0101] At this point, the calculation is completed, and the axial force and axial strain of the wall are output. The axial force is converted into stress using the formula, and the stress-strain curve is plotted. The stress σ t The calculation formula is as follows:

[0102]

[0103] where σ t is the Brazilian splitting stress, P is the wall axial force, D is the specimen diameter, and t is the specimen thickness.

[0104] The axial force, axial strain and displacement of the wall are monitored in the discrete element simulation of the triaxial test. The maximum peak force is continuously updated during the monitoring process. When the current force is less than the threshold of the peak force, the specimen is judged to be damaged.

[0105] At this point, the calculation is completed, and the axial force and axial strain of the wall are output. The axial force is converted into stress using the formula, and the stress-strain curve is plotted. The stress σ c The calculation formula is as follows:

[0106]

[0107] where σ c is the axial stress, P 1 is the bottom wall pressure, P 2 is the top wall pressure, and D is the diameter of the cylindrical specimen.

[0108] In this embodiment, the threshold is 80%.

[0109] In another embodiment, the calibrating the discrete element sample calculation parameters based on the numerical simulation results includes:

[0110] The deformation parameters of the discrete element model are calibrated by the triaxial test results. When the bedding direction is perpendicular to the loading direction, the bedding plane has the least effect on the elastic modulus of the specimen. At this time, the stiffness parameters of the smooth joint model and the flat joint model between grains are calibrated.

[0111] When the bedding direction coincides with the loading direction, the bedding plane stiffness has the greatest influence, and the stiffness parameters of the smooth joint model of the bedding plane are calibrated;

[0112] The strength parameters of the discrete element model are calibrated according to the Brazilian splitting test results. When the bedding plane direction is perpendicular to the loading direction, the strength parameters of the smooth joint model and the flat joint model are calibrated.

[0113] When the bedding direction coincides with the loading direction, the stiffness of the bedding plane has the greatest influence, and the strength parameters of the smooth joint model of the bedding plane are calibrated.

[0114] Specifically in this embodiment, the main calibrated parameters are the flat joint model parameters of the grains inside the rock and the smooth joint model parameters between groups (between grains) and weak bedding planes; the deformation parameters of the discrete element model are calibrated by the triaxial test results. When the bedding direction is perpendicular to the loading direction, the bedding plane has the least effect on the elastic modulus of the specimen. At this time, the stiffness parameters of the smooth joint model and the flat joint model between the grains are calibrated; when the bedding direction coincides with the loading direction, the bedding plane stiffness has the greatest effect, and the stiffness parameters of the smooth joint model of the bedding plane are calibrated; the strength parameters of the discrete element model are calibrated according to the Brazilian splitting test results. When the bedding plane direction is perpendicular to the loading direction, the strength parameters of the smooth joint model and the flat joint model are calibrated; when the bedding direction coincides with the loading direction, the bedding plane stiffness has the greatest effect, and the strength parameters of the smooth joint model of the bedding plane are calibrated. The calibration process is as follows: Figure 6 As shown, the deformation parameters include the effective modulus of the flat joint model Flat joint model stiffness ratio Crystal

[0115] Normal stiffness of smooth joint model Tangential stiffness of intergranular smooth joint model Normal stiffness of smooth joint model of bedding plane Shear stiffness of smooth joint model of bedding plane The strength parameters include the flat joint model tensile strength Flat joint model cohesion Friction angle of flat joint model Tensile strength of smooth intergranular joint model Intergranular smooth joint model cohesion Friction angle of intergranular smooth joint model Tensile strength of smooth joint model on bedding plane Bedding plane smooth joint model cohesion Friction angle of smooth joint model of bedding plane The calibration results are shown in Table 2, and compared with the test strength Figure 7 As shown, the failure mode comparison is Figure 8 shown.

[0116] Table 2 Calibration parameters

[0117]

[0118] In another embodiment, the process of acquiring the crack distribution of the polycrystalline model includes:

[0119] According to the different characteristics and failure modes of the three types of contact modes, the crack types are divided into 6 categories, namely, flat joint model tensile crack, flat joint shear crack, intergranular smooth joint model tensile crack, intergranular smooth joint model shear crack, bedding plane smooth joint model tensile crack, bedding plane smooth joint model shear crack;

[0120] During the test loading process, the value of the increase in the number of cracks as the strain increases is recorded, and a crack change picture is drawn to determine which type of damage accounts for the largest proportion in the test.

[0121] In another embodiment, a sensitivity analysis is performed on the model to analyze the effect of the change in the bedding direction on the strength, failure mode and crack propagation of the polycrystalline model. The inclination direction changes from 15° to 90°, and the inclination angle changes from 0° to 90°. The strength change results are as follows: Fig. 9 shown. Fig. 9 a represents the destructive strength of the specimens at various inclinations and inclination angles under the Brazilian splitting test. Fig. 9 b Distribution of the six types of cracks after the specimen loading is completed: flat joint model tensile crack, flat joint shear crack, intergranular smooth joint model tensile crack, intergranular smooth joint model shear crack, bedding plane smooth joint model tensile crack, bedding plane smooth joint model shear crack. Fig. 9 c- Fig. 9 i is the change diagram of crack proportion with inclination and inclination in sensitivity analysis. Specifically, Fig. 9 c- Fig. 9 e are the specific crack distribution diagrams at inclinations of 15°30°45°60°75°90°, Fig. 9 i means Fig. 9 c- Fig. 9 The diagram of the six types of cracks involved in e, GBSC refers to intra-group shear cracks, GBTC refers to intra-group tensile cracks, BPSC refers to joint shear cracks, BPTC refers to intra-group tensile cracks, IGSC refers to inter-group shear cracks, and IGTC refers to inter-group tensile cracks.

[0122] The above-mentioned embodiments only express the specific implementation methods of the present application, and the descriptions thereof are relatively specific and detailed, but they cannot be understood as limiting the protection scope of the present application. It should be pointed out that, for ordinary technicians in this field, several variations and improvements can be made without departing from the technical solution concept of the present application, and these all belong to the protection scope of the present application.

Claims

1. Modeling and calculation method of crystalline rock with bedding structure based on Neper-DEM, characterized by: include: Construct a polycrystalline structure with crystal characteristics in Neper, and then generate the polycrystalline structure into a geometric parameter text of the grains; A discrete element rock sample model is constructed in DEM, including constructing discrete element samples and importing geometric parameter text as geometric information, grouping discrete element samples according to the geometric information, constructing multiple mutually parallel bedding planes and inter-group joint planes whose directions are consistent with the grain intersection plane in the discrete element samples, and then assigning contact models and parameters of three types of contact modes to the discrete element samples, the three types of contact modes including intra-group contact of discrete element samples and indirect contact of discrete element samples with joint planes; A numerical simulation test is performed on a discrete element rock sample model to obtain a stress-strain curve and a sample failure mode, wherein the numerical simulation is implemented using a Brazilian splitting test and a triaxial test; The calculation parameters of the discrete element sample are calibrated according to the numerical simulation results to verify the effectiveness of the modeling and calculation method of the discrete element rock sample model.

2. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The polycrystalline structure creation process includes: The lens is constructed using the 3D Thiessen polygon tessellation method based on the crystal characteristics, namely the number of grains, grain shape, and grain distribution. The number of grains is specified in the Neper command, the grain shape is reflected by the normal distribution of equivalent diameter and sphericity, and the grain distribution is controlled by random positions and seed generation weights. Under the influence of different seed weights, the boundaries between grains become non-equidistant between seeds. The greater the seed weight, the larger the associated grain. After the Thiessen polygon tessellation is completed, the various parameters of the lens are output as an STL file.

3. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The grouping of discrete element samples according to geometric information comprises: Importing geometric information, the geometric information includes normal vectors of each face of the convex grain and node coordinates of each face, and multiple faces together constitute a grain entity; using DEM to read the geometric information, generate the same node domain surface in the DEM, and construct a single Geometry based on the convex grain composed of multiple faces and name it according to the import order; The existing discrete element samples are grouped according to the named Geometry. When multiple particles fall within the Geometry polyhedron, they are identified as a separate group, and their naming is consistent with the Geometry. Repeat the above steps until all the grains constituting the polycrystalline body are introduced and grouped.

4. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The method of constructing a plurality of mutually parallel bedding planes and inter-group joint planes whose directions are consistent with the grain intersection planes comprises: A certain number of parallel joint planes are added to the polycrystalline discrete element sample, and the bedding direction is controlled by the dip and inclination angle to be consistent with the original sample and the spacing between each bedding is kept consistent; Construct inter-group joint surfaces according to the directions of the imported grain surfaces, and determine whether the two contacting grains are in the same group. If they are in the same group, do not process them. If they are not in the same group, read the group names respectively, and the group names correspond to the Geometry with the same name. Find the two corresponding Geometries according to the name, determine the common face of the Geometries, obtain the direction of the face and assign the same direction to the two contacting particle joint surfaces.

5. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The contact models and parameters for assigning the three types of contact modes include: According to the established bedding planes and inter-group joint planes, the contact model and parameter assignment are performed on the discrete element model; The bedding plane is given a smooth joint model, which has the same direction as the bedding plane and a small contact parameter to simulate a weak structural plane; The inter-group joint surface is also given a smooth joint model to simulate the interface contact between crystals, and the contact parameter value is larger than that of the layer; The particles that are not covered by inter-group joint planes and bedding planes are assigned a flat joint model.

6. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The experimental process of the numerical simulation test is as follows: Create wall units at the top and bottom of the discrete element rock specimen model, use the wall units to simulate the pressure head of the Brazilian splitting test, set the modulus and width of the wall units, and fix the wall loading speed; When the wall starts to be loaded, the displacement of the wall compression, the force during the compression process, the crack distribution and number of the polycrystalline model are recorded. The data and pictures are automatically exported through the DEM command to obtain the stress-strain curve, elastic modulus, failure mode and failure strength.

7. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The obtaining of the stress-strain curve comprises: The axial force, axial strain and displacement of the wall are monitored in the discrete element simulation of the Brazilian splitting test. The maximum peak force is continuously updated during the monitoring process. If the newly generated force is less than the threshold of the peak force, the specimen is judged to be damaged. At this point, the calculation is completed, and the axial force and axial strain of the wall are output. The axial force is converted into stress using the formula, and the stress-strain curve is plotted. The stress σ t The calculation formula is as follows: where σ t is the Brazilian splitting stress, P is the wall axial force, D is the specimen diameter, and t is the specimen thickness.

8. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The obtaining of the stress-strain curve further comprises: The axial force, axial strain and displacement of the wall are monitored in the discrete element simulation of the triaxial test. The maximum peak force is continuously updated during the monitoring process. When the current force is less than the threshold of the peak force, the specimen is judged to be damaged. At this point, the calculation is completed, and the axial force and axial strain of the wall are output. The axial force is converted into stress using the formula, and the stress-strain curve is plotted. The stress σ c The calculation formula is as follows: where σ c is the axial stress, P1 is the bottom wall pressure, P2 is the top wall pressure, and D is the diameter of the cylindrical specimen.

9. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 1, characterized in that: The calibrating of discrete element sample calculation parameters according to the numerical simulation results includes: The deformation parameters of the discrete element model are calibrated by the triaxial test results. When the bedding direction is perpendicular to the loading direction, the bedding plane has the least effect on the elastic modulus of the specimen. At this time, the stiffness parameters of the smooth joint model and the flat joint model between grains are calibrated. When the bedding direction coincides with the loading direction, the bedding plane stiffness has the greatest influence, and the stiffness parameters of the smooth joint model of the bedding plane are calibrated; The strength parameters of the discrete element model are calibrated according to the Brazilian splitting test results. When the bedding plane direction is perpendicular to the loading direction, the strength parameters of the smooth joint model and the flat joint model are calibrated. When the bedding direction coincides with the loading direction, the stiffness of the bedding plane has the greatest influence, and the strength parameters of the smooth joint model of the bedding plane are calibrated.

10. The method for modeling and calculating crystalline rocks with bedding structure based on Neper-DEM according to claim 6, characterized in that: The process of obtaining the crack distribution of the polycrystalline model includes: According to the different characteristics and failure modes of the three types of contact modes, the crack types are divided into 6 categories, namely, flat joint model tensile crack, flat joint shear crack, intergranular smooth joint model tensile crack, intergranular smooth joint model shear crack, bedding plane smooth joint model tensile crack, bedding plane smooth joint model shear crack; During the test loading process, the value of the increase in the number of cracks as the strain increases is recorded, and a crack change picture is drawn to determine which type of damage accounts for the largest proportion in the test.

Images