PD-FEM numerical calculation method and system for full-process simulation of rock mass fracture at engineering scale

By dividing the rock mass calculation model into near-field dynamics zone, finite element zone and coupling zone, combined with the PD-FEM method, the problems of low computational efficiency and insufficient accuracy of rock mass fracture simulation in the existing technology are solved, and efficient engineering is achieved Scaling rock mass fracture simulation.

CN113761760BActive Publication Date: 2024-12-27SHANDONG UNIV
View PDF 0 Cites 0 Cited by

Patent Information

Application Number
CN202110824416.5
Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
Filing Date
2021-07-21
Publication Date
2024-12-27
Estimated Expiration
2041-07-21

AI Technical Summary

Technical Problem

The existing rock mass fracture simulation methods have problems such as low computational efficiency, insufficient accuracy and difficulty in grid division when dealing with discontinuity and crack propagation, especially in the engineering scale, which is difficult to achieve efficient numerical calculations.

Method used

Using PD-FEM numerical calculation method, the rock mass calculation model is divided into near-field dynamics zone, finite element zone and coupling zone. Through the transmission and iterative calculation of force and displacement, engineering scale simulation of rock mass rupture is realized, combining data storage and Post-processing module to improve computing efficiency and accuracy.

Benefits of technology

It effectively reduces the calculation workload, improves the calculation efficiency, ensures the accuracy of the calculation results, and can efficiently simulate the rock mass fracture process on the engineering scale.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN113761760B_ABST
    Figure CN113761760B_ABST
Patent Text Reader

Abstract

The present invention provides a PD-FEM numerical calculation method and system for simulating the whole process of rock mass fracture at an engineering scale. According to whether damage occurs, the rock mass calculation model is divided into a near-field dynamics zone, a finite element zone and a coupling zone; force and displacement are transmitted in the coupling zone, the gravity and tectonic stress of the overlying rock strata are converted into finite element regional force boundary conditions, and the constraints are converted into displacement boundary conditions, and the node displacements of the finite element zone are calculated, which are used as boundary conditions of the near-field dynamics zone; under the above boundary conditions, material point displacements and damage are iteratively calculated to solve the deformation characteristics of the rock mass and realize engineering scale simulation of rock mass fracture. The present invention can effectively reduce the workload, ensure the calculation accuracy and greatly improve the calculation efficiency, which is of great significance for solving the engineering scale problem of rock mass fracture.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] The present invention belongs to the technical field of numerical simulation of the whole process of rock mass fracture, and specifically relates to a PD-FEM numerical calculation method and system for simulating the whole process of rock mass fracture at an engineering scale. Background Art

[0002] The statements in this section merely provide background information related to the present invention and do not necessarily constitute prior art.

[0003] The problem of rock mass fracture and failure has always been one of the important basic scientific problems in the engineering field. In order to simulate the destruction and crack propagation of rock mass, a variety of numerical calculation methods have been proposed. Conventional simulation analysis methods include the finite element method (FEM) based on traditional continuous medium mechanics and the discrete element method (DEM) based on discontinuous medium mechanics. The finite element method assumes that there is a continuous and differentiable displacement field near any point in the structure. However, there are many discontinuities in damage and failure analysis. It is difficult to obtain the partial differential equation of the displacement field in the discontinuous region, and there is a problem of tip singularity. At the same time, it is necessary to preset the crack path and crack propagation criteria, and the mesh needs to be re-divided after the crack propagates, which makes it difficult to solve. Although numerical methods such as the cohesive zone model (CZM) and the extended finite element method (XFEM) can deal with such discontinuous problems, their use also has many limitations. When the cohesive zone model deals with crack propagation problems, cracks can only be distributed on the unit surface; the extended finite element method has difficulties in dealing with crack bifurcation, multi-crack propagation, and three-dimensional crack propagation. The discrete element method cuts the rock mass into discrete blocks as the research object, and uses Newton's second law to describe its movement process. It can simulate the discontinuous deformation and large deformation of the rock mass. Although it overcomes the disadvantage of the finite element method that it requires continuous grid division, when simulating the expansion of rock cracks, this method can only expand along the contact surface between the unit blocks, which is equivalent to artificially determining the expansion path of the crack in advance, affecting the accuracy of the simulation.

[0004] Peridynamics (PD) is a numerical calculation method that builds models based on the idea of ​​non-local action and describes the mechanical behavior of materials by solving spatial integral equations. It combines the advantages of molecular dynamics and network-free methods, breaks through the limitations of traditional methods based on continuity assumptions in solving discontinuous problems, avoids the singularity of solving differential equations at the crack tip, and has unique advantages in describing material crack generation, initiation, expansion, and destruction. It has been widely used in the field of solid mechanics. Compared with the finite element method based on traditional continuous medium theory, peridynamics has lower computational efficiency than the finite element method due to its non-local effects and the use of meshless particle method for solution, and has higher requirements for computer memory. Summary of the invention

[0005] In order to solve the above problems, the present invention proposes a PD-FEM numerical calculation method and system for simulating the whole process of rock fracture at an engineering scale. The present invention can effectively reduce the workload, ensure the calculation accuracy and greatly improve the calculation efficiency, which is of great significance for solving the engineering scale problem of rock fracture.

[0006] According to some embodiments, the present invention adopts the following technical solutions:

[0007] A PD-FEM numerical calculation method for simulating the whole process of rock mass fracture at an engineering scale includes the following steps:

[0008] According to whether damage occurs, the rock mass calculation model is divided into the peridynamic area, the finite element area and the coupling area;

[0009] The force and displacement are transferred in the coupling zone, the gravity and tectonic stress of the overburden are converted into the force boundary conditions of the finite element zone, the constraints are converted into displacement boundary conditions, the node displacements of the finite element zone are calculated, and this is used as the boundary conditions of the peridynamic zone;

[0010] Under the above boundary conditions, the material point displacement and damage are iteratively calculated to solve the deformation characteristics of the rock mass and realize the engineering scale simulation of rock mass fracture.

[0011] As an optional implementation method, the specific process of dividing the rock mass calculation model into a near-field dynamics zone, a finite element zone and a coupling zone includes: establishing a corresponding size model according to the actual project, dividing the model into three areas, the area where the rock mass may break is the near-field dynamics zone, the area far away from the crack is the finite element zone, and the part where the two areas are connected is the coupling zone.

[0012] As a further limitation, the coupling zone is composed of connection units, which are finite element units embedded with peridynamic material points. The embedded material points can be used to transmit force and displacement, and can also serve as virtual boundary layers in the peridynamic zone, as objects to which boundary conditions are applied, thereby weakening the influence of boundary effects on calculation results.

[0013] As an optional implementation, the node displacements of the finite element zone are calculated and used as the boundary conditions of the peridynamic zone. The peridynamic zone and the finite element zone are connected by transferring forces and displacements through the coupling zone. In the connection unit, the displacements of the embedded peridynamic material points are obtained by interpolating the node displacements of the finite element unit through shape functions. At the same time, the force generated by the interaction between the embedded peridynamic material points and the material points outside the connection unit is called the coupling force, which is distributed to the nodes of the connection unit located at the interface through the shape function as a part of the internal force received by the nodes.

[0014] As an optional implementation method, the specific process of calculating the node displacement of the finite element zone includes: iteratively solving through the central difference method, first solving the finite element equation to obtain the displacement of the finite element node, using the obtained finite element node displacement as the boundary condition of the near-field dynamics zone, applying it to the material points in the connection unit, and solving the velocity and displacement of the material point.

[0015] As an optional implementation method, the specific process of iteratively calculating the displacement and damage of material points includes: judging whether the bonds of material points in the near-field dynamics zone meet the destruction conditions, assigning corresponding constants to the bonds of corresponding material points according to different judgment results, recording the local damage conditions respectively, and then describing the fracture of the rock mass, and finally realizing the engineering scale simulation of rock fracture.

[0016] A PD-FEM numerical calculation system for simulating the whole process of rock mass fracture at an engineering scale, including:

[0017] A model division module is configured to divide the rock mass computational model into a peridynamics area, a finite element area, and a coupling area according to whether damage occurs;

[0018] The boundary condition determination module is configured to transfer forces and displacements in the coupling zone, convert the overburden gravity and tectonic stress into the force boundary conditions of the finite element zone, convert the constraints into displacement boundary conditions, calculate the node displacements of the finite element zone, and use them as the boundary conditions of the peridynamic zone;

[0019] The solution judgment module is configured to iteratively calculate the displacement and damage of material points under the above boundary conditions, solve the deformation characteristics of the rock mass, and realize the engineering scale simulation of rock mass fracture.

[0020] It also includes a data storage and post-processing module, which is configured to store the calculation results and intuitively display the calculation results in the form of a cloud map.

[0021] A computer-readable storage medium is used to store computer instructions, and when the computer instructions are executed by a processor, the steps in the above method are completed.

[0022] An electronic device comprises a memory and a processor, and computer instructions stored in the memory and executed on the processor. When the computer instructions are executed by the processor, the steps in the above method are completed.

[0023] Compared with the prior art, the present invention has the following beneficial effects:

[0024] (1) One or more embodiments of the present invention introduce the PD-FEM coupling calculation theory into the field of underground engineering, and divide the engineering-scale rock mass calculation model into a peridynamic zone, a finite element zone, and a coupling zone according to whether damage occurs. Compared with the pure peridynamic method, the number of material points is greatly reduced, and the calculation workload is reduced;

[0025] (2) One or more embodiments of the present invention only need to model the core calculation area, and convert the gravity load of the overlying rock strata on the calculation area and the horizontal tectonic stress into stress boundary conditions on the calculation area, and convert the constraints into displacement boundary conditions, thereby improving the calculation efficiency and ensuring the accuracy of the calculation results;

[0026] (3) One or more embodiments of the present invention realize coupled calculation of the near-field dynamics zone and the finite element zone through force and displacement transmission between the two, while giving full play to the advantages of both, ensuring the calculation accuracy while greatly improving the calculation efficiency, thereby realizing engineering-scale simulation of rock fracture.

[0027] In order to make the above-mentioned objects, features and advantages of the present invention more obvious and easy to understand, preferred embodiments are given below and described in detail with reference to the accompanying drawings. BRIEF DESCRIPTION OF THE DRAWINGS

[0028] The accompanying drawings in the specification, which constitute a part of the present invention, are used to provide a further understanding of the present invention. The exemplary embodiments of the present invention and their descriptions are used to explain the present invention and do not constitute improper limitations on the present invention.

[0029] Figure 1 is a flow chart of Embodiment 1 of the present invention;

[0030] Figure 2 is a schematic diagram of a coupling method according to a first embodiment of the present invention;

[0031] Figure 3 It is a schematic diagram of a large underground cavern group excavation model according to the second embodiment of the present invention. DETAILED DESCRIPTION

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

[0033] It should be noted that the following detailed descriptions are all illustrative and intended to provide further explanation of the present invention. Unless otherwise specified, all technical and scientific terms used herein have the same meanings as those commonly understood by those skilled in the art to which the present invention belongs.

[0034] It should be noted that the terms used herein are only for describing specific embodiments and are not intended to limit exemplary embodiments according to the present invention. As used herein, unless the context clearly indicates otherwise, the singular form is also intended to include the plural form. In addition, it should be understood that when the terms "comprising" and / or "including" are used in this specification, it indicates the presence of features, steps, operations, devices, components and / or combinations thereof.

[0035] The following is combined with Figure 1 The present invention is described in detail, specifically, as follows:

[0036] This embodiment provides a PD-FEM numerical calculation method for simulating the entire process of rock mass fracture at an engineering scale, including the following steps:

[0037] (1) According to the actual project, a corresponding rock mass calculation model is established and divided into three regions. The area where fracture may occur is the peridynamic area, the area far away from the crack is the finite element area, and the part where the two areas are connected is the coupling area, which is composed of connecting units.

[0038] (2) The finite element area and the coupling area are meshed. The two-dimensional model uses plane four-node isoparametric elements, and the three-dimensional model uses hexahedron node isoparametric elements. The near-field dynamics area is discretized into a series of material points in space with material physical and mechanical information, and a certain number of material points are embedded in the finite element unit of the coupling area to form a connection unit. The discretization method of the embedded material points is the same as that of the near-field dynamics area, and the number is related to the range of the near-field dynamics neighborhood. See the schematic diagram of the coupling model for details. Figure 2 .

[0039] The material points embedded in the connection unit can be used to transmit force and displacement, and also serve as a virtual boundary layer in the near-field dynamics zone. As the object of boundary conditions, it weakens the influence of boundary effects on the calculation results.

[0040] The neighborhood refers to the near-field range with which a material point interacts:

[0041] H x ={x'∈Ω:||x'-x||≤δ};

[0042] Wherein, Ω is the calculation area, x is any material point in the calculation area, and x' is any other material point within a certain spatial range of material point x. If the distance between two points is not greater than a given constant δ, there is a certain interaction relationship between the two points, and the range size δ is the neighborhood size.

[0043] The interaction relationship between material points can be expressed by the concept of bonds.

[0044] (3) Calculate the unit stiffness matrix and unit mass matrix of each unit, and integrate the overall stiffness matrix and overall mass matrix through the unit node degree of freedom transformation matrix.

[0045] Element stiffness matrix:

[0046] Among them, Ω (e) is the unit area, C (e) is the unit elastic constant, B (e) is the unit strain matrix;

[0047] Element mass matrix:

[0048] Among them, ρ (e) is the unit density, H (e) is the displacement interpolation matrix in the finite element;

[0049] Overall stiffness matrix: K = G T K e G

[0050] Overall mass matrix: M = G T M e G

[0051] Among them, the conversion matrix G is specifically:

[0052] Number the nodes.

[0053] (4) According to the rock mass in-situ stress conditions, equivalent in-situ stress loads are applied to the boundaries of the finite element region to form a node equivalent external load array F.

[0054] Specifically,

[0055]

[0056] in, is the equivalent surface force at the unit node, is the unit equivalent strength.

[0057] (5) Calculate the internal force matrix F of the node int , the node displacements are solved according to the finite element equations of motion.

[0058] The internal force of the unit nodes on the interface between the finite zone and the coupling zone and the finite element zone is:

[0059] F I int =K I U I

[0060] Among them, K Iis the stiffness of node I, U I is the displacement of node I;

[0061] The coupling force refers to the near-field force generated by the interaction between the embedded material points in the connecting unit and the material points in the peridynamic zone. The interaction between the material points in the connecting unit is not considered as the coupling force. The coupling force borne by the embedded material points in the unit is distributed to the nodes on the interface between the coupling zone and the peridynamic zone through shape function interpolation:

[0062]

[0063] Among them, φ I is the shape function of node I in the bridge unit, (ξ,η) is the natural coordinate of the embedded material point in the bridge unit projected on the interface between the coupling region and the peridynamic region, and f i cp is the coupling force on embedded material point i.

[0064] The internal force of the node on the interface between the coupling zone and the peridynamic zone is the sum of the finite element internal force and the coupling force, and its motion equation is:

[0065]

[0066] Among them, M I is the total mass of node I, is the acceleration of node I, External force, internal force

[0067] (6) The obtained node displacements of the connecting unit are used as the boundary conditions of the peridynamic region, and the displacements of the material points within the connecting unit are calculated by shape function interpolation.

[0068] Specifically:

[0069] Among them, φ I is the shape function of node I in the bridge unit, is the natural coordinate of the embedded material point in the bridge unit projected on the interface between the coupling region and the peridynamic region, U I is the displacement of node I.

[0070] (7) Calculate the bond forces between material points. The displacement of material points in the peridynamic zone is obtained by solving the peridynamic control equations. The displacement of material points in the coupling zone is not included in the calculation.

[0071] The equation of motion for material point i is:

[0072]

[0073] Where ρ is the density of material point i, ü is the acceleration of material point i, is the neighborhood range of material point i, the interaction force function f between material points is called the “constitutive force function”, and b is the body force exerted on the material point.

[0074] Since the interaction force function between material points is consistent with the distribution direction of the material points at the bond ends, the magnitude of the peridynamic constitutive force function can be expressed as a scalar function:

[0075]

[0076] (8) During the iterative solution process, determine whether the bonds of the material points in the near-field dynamics zone meet the destruction conditions and record the local damage conditions.

[0077] The failure condition is the judgment of the integrity of the material point bond expressed by the critical elongation:

[0078]

[0079] Among them, s 0 is the elongation of a given material point bond, s is the elongation of a material point bond, expressed as Among them, η is the relative displacement between any two material points, and ξ is the relative elongation between any two material points. That is, when the tensile deformation s of the material point bond exceeds the given limit value s 0 When , the bond breaks, and at this time, there is no longer an interaction relationship between the two interacting material points.

[0080] Local damage is defined as the ratio of the number of remaining intact bonds to the initial number of bonds after the material point bonds are broken, expressed as:

[0081]

[0082] Among them, V ξ is the volume of material point i. Note that Among them, 0 represents an intact state, and 1 represents a completely damaged state. The value between 0 and 1 is a quantitative representation of the degree of local damage.

[0083] Iterative solution refers to the use of central difference to solve the velocity and displacement of unit nodes and material points in each time step. If the specified time step is not completed, the velocity and displacement of the next time step are iteratively solved.

[0084] Embodiment 2:

[0085] This embodiment provides a PD-FEM numerical calculation method for simulating the entire process of rock mass fracture at an engineering scale, including the following steps:

[0086] (1) Initialize model parameters

[0087] This example simulates the excavation process of a large underground cavern group. The model length and width are 120m and 90m, and three underground caverns of different sizes are excavated in the middle of the model. The rock mass elastic modulus is 6GPa, Poisson's ratio is 0.33, and density is 2650kg / m 3 The burial depth is 120m and the lateral pressure coefficient is 1.5.

[0088] (2) Model region division and discretization

[0089] The area within one times the diameter of the excavated cavern is divided into the near-field dynamics zone, which is 60m long and 30m wide respectively. Outwards are the coupling zone and the finite element zone.

[0090] The size of each lattice in the near-field dynamics zone is 0.15m×0.15m×0.15m, the material point is located at the center of the lattice, and the near-field range of the material point is 0.4725m (about 3.15 times the lattice spacing). The finite element zone unit uses a plane four-node isoparametric element, and the size of each unit is 1.5m×1.5m. The unit division of the coupling zone is the same as that of the finite element zone, and three layers of lattices with a size of 0.15m×0.15m×0.15m are embedded in the units of the coupling zone. A total of 4000 units and 83600 material points.

[0091] (3) Coupling calculation

[0092] (a) Initialize the bonds of all material points:

[0093] Retrieve the numbers of other material points within a given neighborhood range (||x'-x||≤0.4725m) of each material point and store them in the matrix to initialize the scalar coefficient matrix Each element of μ is 1, Each element is 0, that is, in the initial case, the bonds of all material points are complete and there is no local damage.

[0094] (b) Calculation of overall stiffness matrix and overall mass matrix:

[0095] Element stiffness matrix:

[0096] Among them, Ω (e) is the unit area, C (e) is the unit elastic constant, B (e) is the unit strain matrix;

[0097] Element mass matrix:

[0098] Among them, ρ (e) is the unit density, H (e) is the displacement interpolation matrix in the finite element;

[0099] Overall stiffness matrix: K = G T K e G

[0100] Overall mass matrix: M = G T M e G

[0101] Among them, the conversion matrix G is specifically:

[0102] Number the nodes.

[0103] (c) Applying global boundary conditions

[0104] The vertical geostress generated by the overlying rock mass under the action of gravity is converted into the equivalent node force load of the upper boundary of the finite element, and the horizontal geostress generated by the tectonic action is converted into the equivalent node force load of the right boundary of the finite element. Normal fixed constraints are applied to the left boundary and the lower boundary, that is, the model cannot undergo rigid body displacement in the spatial coordinate system, forming a node equivalent external load array F.

[0105] Specifically,

[0106]

[0107] in, is the equivalent surface force at the unit node, is the unit equivalent strength.

[0108] (d) Calculation of internal forces at finite element nodes and solution of motion equations:

[0109] The internal force of the unit nodes on the interface between the finite zone and the coupling zone and the finite element zone is:

[0110] F I int =K I U I

[0111] Among them, K I is the stiffness of node I, U I is the displacement of node I;

[0112] The coupling force refers to the near-field force generated by the interaction between the embedded material points in the connecting unit and the material points in the peridynamic zone. The interaction between the material points in the connecting unit is not considered as the coupling force. The coupling force borne by the embedded material points in the unit is distributed to the nodes on the interface between the coupling zone and the peridynamic zone through shape function interpolation:

[0113]

[0114] Among them, φ Iis the shape function of node I in the bridge unit, (ξ,η) is the natural coordinate of the embedded material point in the bridge unit projected on the interface between the coupling region and the peridynamic region, is the coupling force on embedded material point i.

[0115] The internal force of the node on the interface between the coupling zone and the peridynamic zone is the sum of the finite element internal force and the coupling force, and its motion equation is:

[0116]

[0117] Among them, M I is the total mass of node I, is the acceleration of node I, External force, internal force

[0118] (e) Boundary conditions applied in the peridynamic region:

[0119] The obtained node displacements of the connecting unit are used as the boundary conditions of the peridynamic region, and the displacements of the material points in the connecting unit are calculated by shape function interpolation.

[0120] Specifically:

[0121] Among them, φ I is the shape function of node I in the bridge unit, is the natural coordinate of the embedded material point in the bridge unit projected on the interface between the coupling region and the peridynamic region, U I is the displacement of node I.

[0122] (f) Calculation of material point displacement in the peridynamic region:

[0123] The equation of motion for material point i is:

[0124]

[0125] Where ρ is the density of material point i, ü is the acceleration of material point i, is the neighborhood range of material point i, the interaction force function f between material points is called the “constitutive force function”, and b is the body force exerted on the material point.

[0126] Since the interaction force function between material points is consistent with the distribution direction of the material points at the bond ends, the magnitude of the peridynamic constitutive force function can be expressed as a scalar function:

[0127]

[0128] (4) Determine the initial equilibrium state

[0129] The displacement changes of specific points in the model are monitored, and the ground stress is considered to be balanced when the displacement is stable.

[0130] (5) Underground cavern excavation

[0131] An underground cavern is excavated using the material point dormancy method. A monitoring point is arranged above, below, left and right of the excavated cavern. Once the displacement of the monitoring point is stable, the surrounding rock of the excavated cavern is considered stable, and the corresponding iteration step is determined as the excavation step (i.e., the excavation is performed once after a certain number of iteration steps).

[0132] The underground caverns are excavated step by step at fixed intervals of iterative steps (excavation steps). In this embodiment, there are 12 excavation steps in total.

[0133] (6) Damage judgment

[0134] The excavation of underground caverns breaks the original ground stress balance. At this time, the surrounding rock produces a large deformation under the action of unloading, that is, the material point produces a large displacement. The critical elongation is used to judge the fracture condition of each material point bond. If the material point elongation s exceeds the critical elongation s 0 , then the bond constant μ of the corresponding material point is 0; if the elongation s of the material point does not exceed the critical elongation s 0 , then the bond constant μ of the corresponding material point is 1; and by integration, we can get

[0135] (7) Post-processing

[0136] After the calculation is completed, the displacement and damage data of the finite element nodes and material points in the model can be obtained. The deformation and damage cloud map of the model during excavation is drawn based on the data, and the deformation law and damage state of the surrounding rock are obtained.

[0137] Embodiment three:

[0138] This embodiment provides a PD-FEM numerical calculation system for simulating the entire process of rock mass fracture at an engineering scale, including:

[0139] Model building and area division module: establish a computational model of corresponding size according to the actual project, divide it into peridynamics area, finite element area and coupling area, discretize the peridynamics area into material points, divide the finite element area and coupling area into unit grids, and embed a certain number of material points in the units of the coupling area. The number of embedded material points depends on the size of the peridynamic neighborhood.

[0140] The ground stress application module converts the gravity load of the overlying rock strata on the calculation area and the horizontal tectonic stress into stress boundary conditions on the calculation area, and converts the constraints into displacement boundary conditions.

[0141] The overlying rock strata refer to the rock mass that is buried at a great depth. Due to the constraints of computational efficiency and computational cost, it is impossible to discretize all the strata for modeling. Therefore, certain assumptions are needed. Only a certain calculation area is discretized using a limited number of finite element units and peridynamic material points, and the external loads such as the gravity load of the overlying rock strata in the calculation model and the horizontal tectonic stress are converted into stress boundary conditions on the boundary of the calculation area.

[0142] At the same time, the space where the rock mass is located is mostly semi-infinite space, and it is necessary to apply stress boundary conditions and displacement boundary conditions in the calculation area to constrain the displacement of the calculation area and eliminate the influence of stiffness displacement. Stress boundary conditions and displacement boundary conditions need to be applied to the finite element boundary.

[0143] Force and displacement transmission module: the displacement of the embedded peridynamic material point in the connection unit is obtained by interpolating the displacement of the finite element node through the shape function. At the same time, the coupling force generated by the interaction between the embedded peridynamic material point and the material point outside the connection unit is distributed to the nodes of the connection unit at the interface through the shape function.

[0144] The solution judgment model is used to solve the velocity and displacement of the finite element nodes and peridynamic material points, and to judge whether the bonds between the material points in the peridynamic zone meet the failure conditions and record the local damage.

[0145] The data storage and post-processing module is used to store the calculation results (node ​​displacement, material point displacement and damage) of each iteration step, and to intuitively display the calculation results in the form of cloud maps.

[0146] Embodiment 4:

[0147] An embodiment of the present invention also provides an electronic device, including a memory, a processor, and a computer program stored in the memory and executable on the processor. When the processor executes the program, the PD-FEM numerical calculation method for simulating the entire process of rock mass fracture at an engineering scale is implemented.

[0148] Embodiment five:

[0149] This embodiment provides a computer-readable storage medium on which a computer program is stored. When the program is executed by a processor, the PD-FEM numerical calculation method for simulating the entire process of rock mass fracture at an engineering scale is implemented.

[0150] The steps involved in the above embodiments 3, 4, and 5 correspond to the method embodiment 1. For the specific implementation, please refer to the relevant description part of embodiment 1. The term "computer-readable storage medium" should be understood as a single medium or multiple media including one or more instruction sets; it should also be understood to include any medium that can store, encode or carry an instruction set for execution by a processor and enable the processor to execute any method in the present invention.

[0151] Those skilled in the art will appreciate that embodiments of the present invention may be provided as methods, systems, or computer program products. Therefore, the present invention may take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware. Moreover, the present invention may take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0152] The present invention is described with reference to flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. It should be understood that each process and / or block in the flowchart and / or block diagram, as well as the combination of processes and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to produce a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0153] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1 A function specified in one or more boxes.

[0154] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process in the computer or other programmable device. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.

[0155] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. For those skilled in the art, the present invention may have various modifications and variations. Any modification, equivalent replacement, improvement, etc. made by those skilled in the art within the spirit and principle of the present invention without creative labor shall be included in the protection scope of the present invention.

Claims

1. A PD-FEM numerical calculation method for simulating the whole process of rock mass fracture at an engineering scale, characterized by: The following steps are involved: According to whether damage occurs, the rock mass calculation model is divided into the peridynamic area, the finite element area and the coupling area; The specific process includes: establishing a model of corresponding size according to the actual project, dividing the model into three areas, the area where the rock mass may break is the peridynamic area, the area far away from the crack is the finite element area, and the part where the two areas are connected is the coupling area; specifically, excavating underground caverns of different sizes in the middle of the model, setting the initialization model parameters, and dividing the excavated caverns within the set cave diameter range outward into the peridynamic area, and outward into the coupling area and the finite element area in turn, discretizing the peridynamic area into material points, and dividing the finite element area and the coupling area into unit grids, and at the same time embedding a certain number of material points in the units of the coupling area, and the number of embedded material points depends on the size of the peridynamic neighborhood range; The underground caverns of different sizes excavated in the middle of the model include: the excavated caverns are divided into the peridynamics area within the range of one times the diameter of the cavern, with a length and width of 60m and 30m respectively, and the coupling area and the finite element area are arranged outwards; the size of each lattice in the peridynamics area is 0.15m×0.15m×0.15m, the material point is located at the center of the lattice, the near-field range of the material point is 0.4725m, the unit of the finite element area adopts a plane four-node isoparametric element, and the size of each unit is 1.5m×1.5m, the unit division of the coupling area is the same as that of the finite element area, and three layers of lattices with a size of 0.15m×0.15m×0.15m are embedded in the units of the coupling area; The force and displacement are transmitted in the coupling zone, the gravity and tectonic stress of the overlying rock strata are converted into the force boundary conditions of the finite element region, the constraints are converted into displacement boundary conditions, the node displacements of the finite element region are calculated, and this is used as the boundary conditions of the near-field dynamics zone; including: applying overall boundary conditions: converting the vertical ground stress generated by the overlying rock mass under the action of gravity into the equivalent node force load of the upper boundary of the finite element, converting the horizontal ground stress generated by the tectonic action into the equivalent node force load of the right boundary of the finite element, and applying normal fixed constraints to the left and lower boundaries, that is, the model cannot undergo rigid body displacement in the spatial coordinate system, forming a node equivalent external load array F; Specifically, in, is the equivalent surface force at the unit node, is the unit equivalent physical strength, Ω (e) is the unit area, H (e) is the displacement interpolation matrix in the finite element; Under the above boundary conditions, material point displacement and damage are iteratively calculated to solve the deformation characteristics of the rock mass and realize engineering scale simulation of rock mass fracture; The coupling force refers to the near-field force generated by the interaction between the embedded material points in the connecting unit and the material points in the near-field dynamics zone. The interaction between the material points in the connecting unit is not considered as the coupling force. The coupling force borne by the embedded material points in the unit is distributed to the nodes on the interface between the coupling zone and the near-field dynamics zone through shape function interpolation.

2. The PD-FEM numerical calculation method for simulating the whole process of rock mass fracture at an engineering scale as claimed in claim 1 is characterized by: The coupling zone is composed of connection units, which are finite element units embedded with peridynamic material points. The embedded material points can be used to transmit force and displacement, and can also serve as virtual boundary layers in the peridynamic zone, as objects to which boundary conditions are applied, thereby weakening the influence of boundary effects on calculation results.

3. The PD-FEM numerical calculation method for simulating the whole process of rock mass fracture at an engineering scale as claimed in claim 1 is characterized by: In the process of calculating the node displacements of the finite element zone and using them as the boundary conditions of the peridynamic zone, the peridynamic zone and the finite element zone are connected by transferring forces and displacements in the coupling zone. In the connecting unit, the displacements of the embedded peridynamic material points are obtained by interpolating the node displacements of the finite element unit through shape functions. At the same time, the force generated by the interaction between the embedded peridynamic material points and the material points outside the connecting unit is called the coupling force, which is distributed to the nodes of the connecting unit located at the interface through shape functions as part of the internal force acting on the nodes.

4. The PD-FEM numerical calculation method for simulating the whole process of rock mass fracture at an engineering scale as claimed in claim 1 is characterized by: The specific process of calculating the node displacement of the finite element zone includes: iterative solution through the central difference method, first solving the finite element equation to obtain the displacement of the finite element node, using the obtained finite element node displacement as the boundary condition of the near-field dynamics zone, applying it to the material points in the connection unit, and solving the velocity and displacement of the material point.

5. The PD-FEM numerical calculation method for simulating the whole process of rock mass fracture at an engineering scale as claimed in claim 1, characterized in that: The specific process of iterative calculation of material point displacement and damage includes: judging whether the bonds of material points in the near-field dynamics zone meet the failure conditions, assigning corresponding constants to the bonds of corresponding material points according to different judgment results, recording the local damage conditions respectively, and then describing the fracture of the rock mass, and finally realizing the engineering scale simulation of rock mass fracture.

6. A PD-FEM numerical calculation system for simulating the whole process of rock mass fracture at an engineering scale, characterized by: include: A model division module is configured to divide the rock mass computational model into a peridynamics region, a finite element region, and a coupling region according to whether damage occurs; The specific process includes: establishing a corresponding size model according to the actual project, dividing the model into three areas, the area where the rock mass may break is the peridynamic area, the area far away from the crack is the finite element area, and the part where the two areas connect is the coupling area; specifically, excavating underground caverns of different sizes in the middle of the model, setting the initialization model parameters, and dividing the excavated caverns into a peridynamic area within a set cave diameter range, and then dividing them into a coupling area and a finite element area outwards, discretizing the peridynamic area into material points, and dividing the finite element area and the coupling area into unit grids. At the same time, a certain number of material points are embedded in the units of the coupling area, and the number of embedded material points depends on the range of the peridynamic neighborhood. Size; The underground caverns of different sizes excavated in the middle of the model include: the excavated caverns are divided into the peridynamics area within the range of one times the diameter of the cavern, with a length and width of 60m and 30m respectively, and the coupling area and the finite element area are arranged outwards; the size of each lattice in the peridynamics area is 0.15m×0.15m×0.15m, the material point is located at the center of the lattice, and the near-field range of the material point is 0.4725m. The unit of the finite element area adopts a plane four-node isoparametric element, and the size of each unit is 1.5m×1.5m. The unit division of the coupling area is the same as that of the finite element area, and three layers of lattices with a size of 0.15m×0.15m×0.15m are embedded in the units of the coupling area; The boundary condition determination module is configured to transfer force and displacement in the coupling zone, convert the gravity and tectonic stress of the overlying rock strata into the force boundary conditions of the finite element region, convert the constraints into displacement boundary conditions, calculate the node displacement of the finite element region, and use this as the boundary condition of the near-field dynamics region; including: applying overall boundary conditions: converting the vertical ground stress generated by the overlying rock mass under the action of gravity into the equivalent node force load of the upper boundary of the finite element, converting the horizontal ground stress generated by the tectonic action into the equivalent node force load of the right boundary of the finite element, and applying normal fixed constraints to the left boundary and the lower boundary, that is, the model cannot have rigid body displacement in the spatial coordinate system, forming a node equivalent external load array F; Specifically, in, is the equivalent surface force at the unit node, is the unit equivalent physical strength, Ω (e) is the unit area, H (e) is the displacement interpolation matrix in the finite element; The solution judgment module is configured to iteratively calculate the displacement and damage of material points under the above boundary conditions, solve the deformation characteristics of the rock mass, and realize the engineering scale simulation of rock mass fracture; The coupling force refers to the near-field force generated by the interaction between the embedded material points in the connecting unit and the material points in the near-field dynamics zone. The interaction between the material points in the connecting unit is not considered as the coupling force. The coupling force borne by the embedded material points in the unit is distributed to the nodes on the interface between the coupling zone and the near-field dynamics zone through shape function interpolation.

7. The PD-FEM numerical calculation system for simulating the whole process of rock mass fracture at an engineering scale as claimed in claim 6, characterized in that: It also includes a data storage and post-processing module, which is configured to store the calculation results and intuitively display the calculation results in the form of a cloud map.

8. A computer-readable storage medium for storing computer instructions, which, when executed by a processor, complete the steps of the PD-FEM numerical calculation method for simulating the entire process of rock mass fracture at an engineering scale as described in any one of claims 1 to 5.

9. An electronic device comprising a memory and a processor and computer instructions stored in the memory and running on the processor, wherein when the computer instructions are run by the processor, the steps of the PD-FEM numerical calculation method for simulating the entire process of rock mass fracture at an engineering scale described in any one of claims 1 to 5 are completed.