Time-dependent continuous-discontinuous number method and non-global contact identification method and system

Abstract

The present invention discloses a time-dependent continuous and discontinuous number method for non-global contact identification and a system, the method comprising the following steps: creating a geometric body numerical model, applying boundary conditions, calculating unit forces, identifying units that may generate contact in the entity units at the boundary and the discrete fracture network unit list, storing and updating the potential contact entity unit list, and then performing contact detection in the potential contact entity unit list, updating the state of each unit and repeating the steps S300-S600 in the attached figure until the simulation is completed. The time-dependent continuous and discontinuous number method for non-global contact identification provided by the present invention can reflect the gradual increase of potential contact pairs with crack extension in real time, can effectively inherit the unit forces of failed joint units, maintain the numerical stability of the numerical model, and has sufficient parallel capability, which is of great significance to improving the operational efficiency of continuous and discontinuous numerical models and improving simulation accuracy.

Description

Technical Field

[0001] The present invention belongs to the field of rock mechanics and rock engineering, and more specifically, relates to a time-related continuous-discontinuous number method non-global contact identification method and system. Background Art

[0002] In actual engineering applications, it is usually necessary to design a reasonable construction plan based on specific geological conditions and structural environment, and adjust and optimize the construction plan through field tests. However, due to the limitations of test conditions and the impact of experimental costs, it is difficult to directly carry out large-scale or equal-scale physical model experiments in many cases, such as mining, foundation pit excavation, and tunnel excavation. Numerical simulation can simulate the field test process through numerical calculation methods to reduce costs, and can also simulate boundary conditions that cannot be achieved in physical experiments, which has important reference value. Combining the advantages of continuous and discontinuous methods, the continuous-discontinuous numerical method can truly simulate the evolution of rock structure from continuous to discontinuous, and describe the actual fracture and crushing process of the rock mass. In recent years, it has been a research hotspot and focus of geotechnical engineering calculation methods.

[0003] At present, the contact detection part of most continuous discontinuous number method (FDEM) software is completed by the spatial hashing contact detection method, such as Cundall's FCP method, Munjiza's NBS method and Liu Hongyuan's potential contact pair method. The spatial hashing contact detection method will uniformly grid a certain area (circular or square) around the discrete body, and each discrete body included in the area is regarded as a potential contact pair. The contact detection is realized based on the grid calculation and the distance between the discrete bodies. This method detects the area within a certain range of the discrete unit at each time step, which not only takes up more computing memory, but also cannot respond to the discrete body generated in the current time step in real time. Moreover, since the current discrete body is generated by the rupture of the joint unit in the previous step, the newly detected contact cannot effectively inherit the unit force of the joint unit. Each time the joint unit fails, the contact force needs to be recalculated, which further leads to the instability of the numerical model.

[0004] At the same time, in the FDEM method, the current solutions to improve contact calculations mainly include using CPU or GPU parallel technology to improve contact computing power and using reasonable methods to avoid redundant contact pairs. If the first solution is adopted, although the contact computing power can be significantly improved, it will be subject to the constraints of the parallel mechanism. The contact method does not allow a certain degree of memory interaction, otherwise data competition will occur, resulting in parallel failure. The most suitable method for parallelism is mainly for the traditional direct search method. The complexity of the direct search method is O(N log N), which occupies a large amount of memory and has an insignificant parallel efficiency. Other robust contact search methods, such as the binary tree method and the non-binary tree method, are not suitable for parallelism. The second solution currently only has the potential contact pair method proposed by Liu Hongyuan et al., but because it does not consider the inheritance of the adhesion of closed cracks before and after fracture, it is prone to instability for most calculations and the improvement in computing efficiency is not significant. Summary of the invention

[0005] In view of the above defects of the prior art, the present invention provides a time-dependent continuous-discontinuous number method non-global contact identification method and system, which aims to address the problems of inability to detect contact in real time, large amount of calculation, and difficulty in inheriting contact force in continuous-discontinuous numerical calculations. Under the framework of the finite element-discrete element coupled numerical method (continuous-discontinuous number method), a time-dependent continuous-discontinuous number method non-global contact identification method is provided, and based on the model, a method for calculating the new contact force in the rupture process is invented.

[0006] In order to achieve the above object, according to a first aspect of the present invention, a time-related continuous-discontinuous number method non-global contact identification method is provided, comprising the following steps:

[0007] S100, creating a geometric body, setting a mesh size, dividing the model into entity units through mesh division, and inserting joint units that play a bonding role between the entity units;

[0008] S200, setting calibrated material parameters, applying boundary conditions and performing balance;

[0009] S300, creating an initial entity unit list, a discrete fracture network unit list, and a joint unit state list, performing numerical simulation calculations on the balanced model, and in the first time step, creating a potential contact entity unit list, identifying entity units at the boundary and units in the discrete fracture network unit list that may come into contact, and storing them in the potential contact entity unit list;

[0010] S400, counting the calculation results of each contact, solid unit, and joint unit in each time step, updating the joint unit state list in real time, and then storing the solid unit sharing the joint unit in the potential contact solid unit list according to the destroyed joint unit, and inheriting the joint unit force;

[0011] S500, creating a contact unit pair list, performing contact detection on the entity units in the potential contact entity unit list by using a contact detection method, calculating the contact force and storing it in the contact unit pair list;

[0012] S600, returning to step S100 to update the states of each node, entity unit, and joint unit, update the list of potential contact entity units in S200 and S300, update the list of contact unit pairs and contact forces in S400, and complete the calculation of one time step;

[0013] S700, repeat steps S300-S600 until a certain number of steps are cycled or the model reaches equilibrium.

[0014] Furthermore, the calculation results in step S300 include contact force, nodal force, nodal displacement, nodal velocity, solid unit stress, solid unit strain, joint unit force, joint unit opening, joint unit sliding displacement and joint unit state.

[0015] Furthermore, step S400 provides a method for calculating the new contact force inherited by the solid element during the fracture of the joint element:

[0016] S401, calculating the stress conditions of each entity unit and joint unit in the current time step, and calculating the connection state of the joint unit according to the equilibrium equation and the yield criterion;

[0017] S402, updating the entity units linked to the joint units that failed in the current time step into the potential contact entity unit list, and inheriting the unit internal force of the joint unit by the linked entity unit according to the failure form of the joint unit, that is, converting the internal force of the joint unit into an external force on the entity unit;

[0018] S403. Recalculate the internal force increment of the solid unit that inherits the bond force of the joint unit, and update the unit stress state according to the constitutive relationship and the internal force increment.

[0019] Furthermore, in step S402, the calculation formula of the unit internal force of the joint unit that inherits the failure of the solid unit is:

[0020] If the failure mode is type I tensile failure, the external force on the solid element can be expressed as formula (1):

[0021]

[0022] If the failure mode is type II shear failure, the external force on the solid element can be expressed as formula (2):

[0023]

[0024] Among them, F n and F s F represents the normal force increment vector and tangential force increment vector on the side of the solid element connected to the failed joint element; jojntn and F joints Represents the normal force vector and tangential force vector of the failed joint unit before failure.

[0025] Furthermore, in step S403, the calculation formula for the internal force increment of the solid unit that inherits the bonding force of the joint unit is:

[0026] Δp'=F n +F s +Δp

[0027] Among them, Δp and Δp' represent the internal force increment vectors before and after the solid element inherits the bond force of the joint element.

[0028] Furthermore, in step S403, according to the constitutive relationship and the internal force increment, the calculation formula for calculating the unit stress is:

[0029]

[0030] Among them, δσ and δτ represent the normal stress increment and tangential stress increment of the entity element, n and s represent the normal vector and tangential vector, respectively. n and δ s denote the normal displacement and tangential displacement of the entity element respectively. represents the stress function.

[0031] According to a second aspect of the present invention, there is provided a time-correlated continuous-discontinuous number method non-global contact recognition system, comprising:

[0032] Model creation module, used to create geometric numerical models, set grid size, divide the model into solid units through grid division and insert bonding units between solid units, set calibrated material parameters, apply boundary conditions and balance;

[0033] The simulation calculation module creates the initial solid unit list, discrete fracture network unit list and joint unit state list, and then performs numerical simulation calculations on the model;

[0034] An update module identifies the entities at the boundary and the entities in the discrete fracture network entity list that may come into contact, stores and updates the potential contact entity list, updates the joint entity status list in real time using the calculation results, stores the entity unit pairs that share the failed joint unit in the potential contact entity list according to the failed joint unit, and updates them;

[0035] A detection and calculation module performs contact detection on the entity elements in the potential contact entity element list by using a contact detection method to form contact pairs, calculates the contact force and stores it in the contact element pair list;

[0036] The loop module is used to update the state of each unit after completing the calculation of a time step, and repeat the update module, detection calculation module and loop module in a loop until all calculations are completed.

[0037] According to a third aspect of the present invention, there is provided a terminal device comprising:

[0038] Memory for storing computer programs;

[0039] A processor is used to implement the steps of the above-mentioned time-related continuous-discontinuous number method non-global contact identification method when executing the computer program.

[0040] According to a fourth aspect of the present invention, there is provided a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, is used to execute the steps of a time-related continuous-discontinuous number method non-global contact identification method.

[0041] In general, the above technical solutions conceived by the present invention can achieve the following beneficial effects compared with the prior art:

[0042] 1. The present invention proposes a non-global contact identification method based on the continuous-discontinuous number method of time evolution. Compared with the traditional spatial hashing contact detection method, it occupies less computer memory, can respond to the discrete body generated by the current time step in real time, and can effectively inherit the unit force of the joint unit to maintain the numerical stability of the numerical model.

[0043] 2. The present invention proposes a non-global contact identification method based on the continuous-discontinuous number method of time evolution, which has a certain parallel capability and is simpler and more efficient in calculation compared with other algorithms with parallel capability.

[0044] 3. The present invention proposes a time-related continuous-discontinuous number method non-global contact identification method. Compared with other methods, fewer potential contact units and contact pairs need to be calculated, which greatly simplifies the model and improves the operating efficiency and calculation speed.

[0045] Additional aspects and advantages of the present application will be partially given in the following description, which will become apparent from the following description, or will be understood through the practice of the present application. BRIEF DESCRIPTION OF THE DRAWINGS

[0046] Figure 1 It is a flow chart of a time-related continuous-discontinuous number method for non-global contact identification in an embodiment of the present invention;

[0047] Figure 2 Schematic diagram of a direct shear test model with a central prefabricated crack in an embodiment of the present invention;

[0048] Figure 3 It is a list of potential contact units and a schematic diagram of the evolution of contact pairs before the rupture of a direct shear test model containing a central prefabricated crack in an embodiment of the present invention;

[0049] Figure 4 It is a list of potential contact units and a schematic diagram of the evolution of contact pairs after the direct shear test model with a central prefabricated crack in an embodiment of the present invention is ruptured;

[0050] Figure 5 It is a schematic diagram of the stability of the calculation results of the direct shear test model containing a central prefabricated crack in the time-dependent contact mode in an embodiment of the present invention;

[0051] Figure 6 A comparison chart of the contact search complexity between the method in the embodiment of the present invention and other traditional methods;

[0052] Figure 7 It is a comparison chart of the contact numbers of the method in the embodiment of the present invention and other traditional methods.

[0053] Figure 8 4 is a diagram of a non-global contact recognition system using a time-related continuous-discontinuous number method in an embodiment of the present invention. DETAILED DESCRIPTION

[0054] In order to make the purpose, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below in conjunction with the accompanying drawings and embodiments.

[0055] It should be understood that the specific embodiments described herein are only used to explain the present invention, and are not used to limit the present invention.

[0056] In addition, the technical features involved in the various embodiments of the present invention described below can be combined with each other as long as they do not conflict with each other.

[0057] Those skilled in the art will understand that, unless otherwise defined, all terms (including technical and scientific terms) used herein have the same meaning as generally understood by those skilled in the art to which this application belongs.

[0058] It should also be understood that terms such as those defined in common dictionaries (including technical and scientific terms) have the same meaning as generally understood by one of ordinary skill in the art to which this application belongs.

[0059] It should also be understood that terms such as those defined in general dictionaries should be understood to have meanings consistent with the meanings in the context of the prior art, and will not be interpreted with idealized or overly formal meanings unless specifically defined as in the embodiments of the present application.

[0060] This embodiment provides a time-dependent continuous-discontinuous number method for non-global contact identification, which is applied to establish a physical model of field conditions in geotechnical engineering such as mining, foundation pit and tunnel excavation.

[0061] The flow chart of the time-related continuous-discontinuous number method non-global contact identification method in this embodiment is as follows: Figure 1 As shown, the steps include:

[0062] S100, create a geometric model, set the mesh size, divide the model into solid units through mesh division, and insert joint units that act as bonding between the solid units

[0063] In one embodiment, the geometric model is as follows: Figure 2 shown.

[0064] Among them, a polygonal specimen with a length of 75 mm and a width of 60 mm needs to be set in the middle of the geometric body, a prefabricated crack is inserted in the middle of the specimen, and a shear box with a thickness of 5 mm is set outside the specimen, and boundary conditions are directly applied to the shear box.

[0065] Furthermore, after meshing, joint units that play a bonding role need to be inserted between the geometric solid units of the specimen to pre-create cracks.

[0066] It will be appreciated that no joint elements are inserted into the shear box blocks.

[0067] S200, assign parameters to the geometric model, apply boundary conditions and perform balancing.

[0068] Furthermore, the calibrated material parameters are set according to different groups.

[0069] In one embodiment, the boundary blocks used to apply boundary conditions are set to two groups of upper shear boxes and lower shear boxes, and the specimens are divided into two groups of upper plates and lower plates according to the two sides of the prefabricated cracks as models.

[0070] Furthermore, if Figure 2As shown, boundary conditions are added, that is, a downward normal load is applied to the upper shear box, and the velocity of the bottom of the model in the vertical direction is limited to 0, and a horizontal velocity v is applied to the right side of the lower shear box, waiting for the next simulation and preliminary balance.

[0071] S300, creating an initial list of entity units, a list of discrete fracture network units, and a list of joint unit states, performing numerical simulation calculations on the balanced model, creating a list of potential contact entity units in the first calculation, identifying entity units at the boundary and units in the list of discrete fracture network units that may come into contact, and storing them in the list of potential contact entity units.

[0072] Further, such as Figure 3 As shown in the figure, the triangular units marked in red are potential contact entity units. The contact detection algorithm will only detect in these units, which are represented by the number 1 in the program; the remaining blue-marked units are non-potential contact entity units, which are represented by the number 0 in the program.

[0073] S400, counting the calculation results of each contact, solid unit, and joint unit in each time step, updating the joint unit status list in real time, and storing the solid units sharing the joint unit in the potential contact solid unit list according to the failed joint unit.

[0074] Furthermore, the simulation calculation results include contact force, nodal force, nodal displacement, nodal velocity, solid unit stress, solid unit strain, joint unit force, joint unit opening, joint unit sliding displacement and joint unit state.

[0075] Furthermore, if Figure 4 As shown in the figure, the crack is generated after the joint unit fails. The corresponding solid units on both sides of the crack inherit the force of the failed joint unit, and are then stored in the list of potential contact solid units (red units) as the objects of the next contact detection algorithm calculation.

[0076] Furthermore, the specific steps of the entity unit inheriting the force of the joint unit are as follows:

[0077] S401, calculating the stress conditions of each entity unit and joint unit in the current time step, and calculating the connection state of the joint unit according to the equilibrium equation and the yield criterion;

[0078] S402, updating the potential contact entity unit list according to the entity units linked to the joint units that failed in the current time step, and inheriting the unit internal force of the joint unit by the linked entity units according to the failure form of the joint units, that is, converting the internal force of the joint unit into the external force on the entity unit;

[0079] Furthermore, if the failure mode of the joint unit is type I tensile failure, the external force on the solid unit can be expressed as formula (1):

[0080]

[0081] If the failure mode of the joint unit is type II shear failure, the external force on the solid unit can be expressed as formula (2):

[0082]

[0083] Among them, F n and F s F represents the normal force increment vector and tangential force increment vector on the side of the solid element connected to the failed joint element; jointn and F joints It represents the normal force vector and tangential force vector of the failed joint unit before failure;

[0084] S403. Recalculate the internal force increment of the solid unit that inherits the bonding force of the joint unit. The calculation formula is:

[0085] Δp'=F n +F s +Δp,

[0086] Among them, Δp and Δp' represent the internal force increment vectors before and after the solid element inherits the bond force of the joint element.

[0087] The element stress state is updated according to the constitutive relationship and internal force increment calculation. The calculation formula is:

[0088]

[0089] Among them, δσ and δτ represent the normal stress increment and tangential stress increment of the solid element, n and s represent the normal vector and tangential vector, respectively. n and δ s denote the normal displacement and tangential displacement of the entity element respectively. represents the stress function.

[0090] S500, create a contact unit pair list, perform contact detection on the entity units in the potential contact entity unit list through a contact detection method, and form contact pairs, calculate the contact force and store it in the contact unit pair list

[0091] Furthermore, the entity element in the potential contact entity element is Figure 3 and Figure 4 The red cell in .

[0092] S600, return to steps S2 and S3 to update the list of potential contact entity units, update the contact unit pair list and contact force in step S4, update the status of each node, entity unit, and joint unit in step S100, and complete the calculation of one time step.

[0093] Furthermore, after the calculation of the first time step is completed, steps S300-S600 are repeated until a certain number of steps are cycled or the simulation reaches equilibrium.

[0094] Furthermore, the stability of the model simulated by the algorithm was evaluated based on the stress cloud diagram and crack morphology during the simulation process.

[0095] It is understandable that if Figure 4 , Figure 5 As shown in the figure, the stress is evenly distributed in other positions of the model and changes in a gradient, mainly concentrated at the crack tip and accompanied by crack breakage, which verifies the stability of the algorithm; the crack extends from the prefabricated crack to the specimen, and the crack extension direction is not unique and is random, which proves the feasibility of the algorithm.

[0096] Furthermore, the complexity of the algorithm is compared with the computational complexity of other algorithms such as direct search method, binary method, NBS method, etc.

[0097] It is understandable that if Figure 6 As shown, the complexity of this method increases slowly with the increase of the number of partition units, which is much lower than the direct search method and has lower complexity than other methods.

[0098] Furthermore, the number of detection units and the number of contact pairs of this method are compared with those of the direct search method.

[0099] It is understandable that if Figure 7 As shown, the potential contact units and contact pairs required to be calculated by this method are much smaller than those of the direct search method, which greatly simplifies the model and improves the operating efficiency and calculation speed.

[0100] This embodiment also provides a time-related continuous-discontinuous number method non-global contact recognition system, Figure 8 This is a schematic diagram of the time-dependent continuous-discontinuous number method non-global contact identification system.

[0101] like Figure 8 As shown, the time-related continuous-discontinuous number method non-global contact identification system includes a model creation module, a simulation calculation module, an update module, a detection calculation module and a loop module;

[0102] Model creation module, used to create geometric numerical models, set grid size, divide the model into solid units through grid division and insert bonding units between solid units, set calibrated material parameters, apply boundary conditions and balance;

[0103] The simulation calculation module creates the initial solid unit list, discrete fracture network unit list and joint unit state list, and then performs numerical simulation calculations on the model;

[0104] The update module is used to identify the units that may come into contact in the list of entity units at the boundary and discrete fracture network units, store and update the list of potential contact entity units, update the list of joint unit status in real time using the calculation results, store the entity unit pairs that share the joint unit in the list of potential contact entity units according to the failed joint unit, and update the list of potential contact entity units again.

[0105] The detection and calculation module is used to create a contact unit pair list, perform contact detection on the entity units in the potential contact entity unit list through a contact detection method, create contact unit pairs, calculate contact forces and store them in the contact unit pair list.

[0106] The loop module is used to update the state of each unit after completing the calculation of a time step, and repeat the update module, detection calculation module and loop module in a loop until all calculations are completed.

[0107] It should be noted that the time-correlated continuous-discontinuous number method non-global contact identification system provided in this embodiment can also be a computer program (including program code) running in a computer device. For example, the time-correlated continuous-discontinuous number method non-global contact identification system as an application can be used to execute the corresponding steps in the time-correlated continuous-discontinuous number method non-global contact identification method provided in this embodiment.

[0108] In some feasible implementations, the time-related continuous and non-continuous number method non-global contact recognition system provided in this embodiment can be implemented in a combination of software and hardware. As an example, the time-related continuous and non-continuous number method non-global contact recognition system of the embodiment of the present application can be a processor in the form of a hardware decoding processor, which is programmed to execute the time-related continuous and non-continuous number method non-global contact recognition system provided in the embodiment of the present application. For example, the processor in the form of a hardware decoding processor can adopt one or more application-specific integrated circuits (ASICs), digital signal processors (DSPs), programmable logic devices (PLDs), complex programmable logic devices (CPLDs), field programmable gate arrays (FPGAs) or other electronic components.

[0109] In some feasible implementations, the time-related continuous and non-continuous number method non-global contact identification system provided in this embodiment can be implemented in software, which can be software in the form of programs and plug-ins, and includes a series of modules to implement the time-related continuous and non-continuous number method non-global contact identification method of the embodiment of the present invention.

[0110] This embodiment also provides a terminal device including:

[0111] Memory for storing computer programs;

[0112] The processor is used to implement the various steps of the above-mentioned time-related continuous-discontinuous number method non-global contact identification method when executing the computer program. For details, please refer to the implementation methods provided in the above-mentioned various steps, which will not be repeated here.

[0113] This embodiment also provides a computer-readable storage medium on which a computer program is stored. When the program is executed by a processor, it is used to execute various steps of a time-related continuous-discontinuous number method non-global contact identification method.

[0114] The computer-readable storage medium provided in the present embodiment creates a geometric body, sets a grid size on the geometric body, divides the model into solid units through grid division, inserts joint units that play a bonding role between the solid units, assigns parameters, applies boundary conditions and performs balance to complete the model construction; creates an initial list of solid units, a list of discrete fracture network units and a list of joint unit states, performs numerical simulation calculations on the balanced model, and at the first time step, it is also necessary to create a list of potential contact solid units, identify the solid units at the boundary and the units that may come into contact in the list of discrete fracture network units, and store them in the list of potential contact solid units; counts the calculation results of each contact, solid unit, and joint unit at each time step, updates the list of joint unit states in real time, and then stores the solid units that share the joint unit in the list of potential contact solid units according to the failed joint units; creates a list of contact unit pairs, performs contact detection on the solid units in the list of potential contact solid units through a contact detection method, calculates the contact force and stores it in the list of contact unit pairs;

[0115] It should be understood that, although the steps in the flowchart of the accompanying drawings are displayed sequentially as indicated by the arrows, these steps are not necessarily executed sequentially in the order indicated by the arrows.

[0116] Unless otherwise specified herein, there is no strict order limitation for the execution of these steps and they may be executed in other orders.

[0117] Moreover, at least part of the steps in the flowchart of the accompanying drawings may include multiple sub-steps or multiple stages. These sub-steps or stages are not necessarily executed at the same time, but can be executed at different times. The order of execution is not necessarily sequential, but can be executed in turn or alternately with other steps or at least part of the sub-steps or stages of other steps.

[0118] The above description is only a partial implementation method of the present application. It should be pointed out that for ordinary technicians in this technical field, several improvements and modifications can be made without departing from the principles of the present application. These improvements and modifications should also be regarded as the scope of protection of the present application.

[0119] It will be easily understood by those skilled in the art that the above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the protection scope of the present invention.

Claims

1. A time-dependent continuous-discontinuous number method for non-global contact identification, characterized in that: The following steps are involved: S100, creating a geometric model, setting a mesh size, dividing the model into entity units through mesh division, and inserting joint units that play a bonding role between the entity units; S200, setting calibrated material parameters, applying boundary conditions and performing balance; S300, creating an initial entity unit list, a discrete fracture network unit list, and a joint unit state list, performing numerical simulation calculations on the balanced model, and in the first time step, creating a potential contact entity unit list, identifying entity units at the boundary and units in the discrete fracture network unit list that may come into contact, and storing them in the potential contact entity unit list; S400, counting the calculation results of each contact pair, entity unit, and joint unit in each time step, updating the joint unit state list in real time, and then storing the entity unit pairs sharing the joint unit in the potential contact entity unit list according to the failed joint unit; S500, creating a contact unit pair list, performing contact detection on entity unit pairs in the potential contact entity unit list by using a contact detection method, calculating contact forces and storing them in the contact unit pair list; S600, returning to step S100 to update the states of each node, entity unit, and joint unit, update the list of potential contact entity units in S200 and S300, update the list of contact unit pairs and contact forces in S400, and complete the calculation of one time step; S700, repeat steps S300-S600 until a certain number of steps are cycled or the model reaches equilibrium.

2. A time-dependent continuous-discontinuous number method for non-global contact identification according to claim 1, characterized in that: The numerical simulation calculation results in step S300 include contact force, nodal force, nodal displacement, nodal velocity, solid unit stress, solid unit strain, joint unit force, joint unit opening, joint unit sliding displacement and joint unit state.

3. A time-dependent continuous-discontinuous number method for non-global contact identification according to claim 1 or 2, characterized in that: In step S400, according to the failed joint unit, the entity unit pair sharing the joint unit is stored in the potential contact entity unit list, including: S401, calculating the stress conditions of each entity unit and joint unit in the current time step, and calculating the connection state of the joint unit according to the equilibrium equation and the yield criterion; S402, updating the entity units linked to the joint units that failed in the current time step into the potential contact entity unit list, and inheriting the unit internal force of the joint unit by the linked entity unit according to the failure form of the joint unit, that is, converting the internal force of the joint unit into an external force on the entity unit; S403. Recalculate the internal force increment of the solid unit that inherits the bond force of the joint unit, and update the unit stress state according to the constitutive relationship and the internal force increment.

4. The time-related continuous-discontinuous number method non-global contact identification method according to claim 3, characterized in that: In the step S402, the calculation method of the unit internal force of the joint unit that inherits the failure of the solid unit is: If the failure mode of the failed joint unit is type I tensile failure, the external force on the solid unit can be expressed as formula (1): Among them, F n and F s Represents the normal force increment vector and tangential force increment vector on the side of the solid element connected to the failed joint element.

5. The time-related continuous-discontinuous number method non-global contact identification method according to claim 4, characterized in that: In the step S402, If the failure mode of the failed joint unit is type II shear failure, the external force on the solid unit can be expressed as formula (2): Among them, F jointn and F joints Represents the normal force vector and tangential force vector of the failed joint unit before failure.

6. The time-related continuous-discontinuous number method non-global contact identification method according to claim 5, characterized in that: In step S403, the calculation method of the internal force increment of the solid unit that inherits the bonding force of the joint unit is: Δp'=F n +F s +Δp Among them, Δp and Δp' represent the internal force increment vectors before and after the solid element inherits the bond force of the joint element.

7. The time-related continuous-discontinuous number method non-global contact identification method according to claim 6, characterized in that: In step S403, according to the constitutive relationship and the internal force increment, the calculation method of the unit stress is: Among them, δσ and δτ represent the normal stress increment and tangential stress increment of the entity element, n and s represent the normal vector and tangential vector, respectively. n and δ s denote the normal displacement and tangential displacement of the entity element respectively. represents the stress function.

8. A time-dependent continuous-discontinuous number method non-global contact recognition system, comprising: Model creation module, used to create geometric numerical models, set grid size, divide the model into solid units through grid division and insert bonding units between solid units, set calibrated material parameters, apply boundary conditions and balance; The simulation calculation module creates the initial solid unit list, discrete fracture network unit list and joint unit state list, and then performs numerical simulation calculations on the model; An updating module is used to identify the entities at the boundary and the entities in the discrete fracture network entity list that may come into contact, store and update the potential contact entity list, update the joint unit status list in real time using the calculation results, store the entity unit pairs that share the joint unit in the potential contact entity list according to the failed joint unit, and update them; A detection and calculation module performs contact detection on the entity elements in the potential contact entity element list by using a contact detection method to form contact pairs, calculates the contact force and stores it in the contact element pair list; The loop module is used to update the state of each unit after completing the calculation of a time step, and repeat the update module, detection calculation module and loop module in a loop until all calculations are completed.

9. A terminal device, characterized in that: include: Memory for storing computer programs; A processor, configured to implement the steps of the time-related continuous-discontinuous number method non-global contact identification method according to any one of claims 1 to 7 when executing the computer program.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the program is executed by a processor, the steps of the time-related continuous-discontinuous number method non-global contact identification method according to any one of claims 1 to 7 are implemented.

Images