Numerical simulation method and system for breaking process of ore rock under explosion load based on FDEM
By combining the FDEM method with finite element mesh and unstructured Delaunay triangulation algorithm, setting material and boundary conditions, and using Duvall-type pressure functions and MultiFracS software to calculate contact forces, the problems of large computational load, complex parameter calibration, and insufficient simulation accuracy in the simulation of ore and rock crushing in the existing technology are solved, and accurate simulation and quantitative data provision of the ore and rock crushing process are realized.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- KUNMING UNIV OF SCI & TECH
- Filing Date
- 2026-03-25
- Publication Date
- 2026-05-29
AI Technical Summary
Existing numerical simulation methods suffer from high computational complexity, complex parameter calibration, crack propagation path dependence on mesh generation, and lack of effective explosion pressure functions and boundary conditions when simulating the rock fracturing process under explosive loading, resulting in insufficient simulation accuracy and reliability.
A method based on FDEM was adopted, and the geometric model of the mineral rock specimen was constructed using the finite element meshing software Gmsh. The unstructured Delaunay triangulation algorithm was introduced for mesh generation, and material parameters and explosive load boundary conditions were set. Contact force was calculated by combining Duvall-type pressure function and MultiFracS software. The nodal coordinates and velocities were updated by normal and tangential contact force functions, and contour plots were generated to analyze crack propagation and energy changes.
It achieves accurate simulation of the rock and ore crushing process, provides quantitative data reference, improves the accuracy and efficiency of simulation, and can truly reflect the dynamic response and crack propagation of rock and ore under explosive load, supporting field engineering practice.
Smart Images

Figure CN121920157B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the fields of rock mechanics and explosive impact technology, specifically to a numerical simulation method and system for the rock and ore fracturing process under FDEM explosive loading. Background Technology
[0002] In the study of rock fracture under explosive loading, numerical simulation is one of the key research methods. Existing numerical simulation methods are mainly divided into two categories: continuous methods and discontinuous methods. The former is represented by the finite element method (FEM), while the latter is typically represented by the discrete element method (DEM). The finite element method (FEM) is mainly based on the assumption of a continuous medium. However, when it comes to simulating discontinuous phenomena such as crack propagation and block separation, this method has inherent limitations. The discrete element method (DEM) focuses on simulating the motion and interaction of blocks, but it cannot accurately describe the stress propagation and deformation behavior in a continuous medium.
[0003] To overcome the limitations of single methods, a method combining finite element and discrete element methods is proposed. This method can simulate the stress-strain behavior of continuous media and effectively handle the contact and separation problems of discontinuous bodies. However, it still faces many challenges when using numerical simulation to solve the rock and ore fracturing process under explosive loading: First, traditional numerical simulation uses continuous media structures and solves by sharing node elements, resulting in a huge computational load; second, the numerical simulation results are extremely sensitive to the values of joint element penalty parameters, making the parameter calibration process complex and time-consuming; third, the crack propagation path is greatly affected by the mesh generation method, reducing the accuracy and reliability of the simulation; finally, there is a lack of effective methods for applying explosive pressure functions and boundary conditions. These problems further limit the application of numerical simulation in the study of rock and ore fracturing processes.
[0004] Therefore, how to accurately simulate the crushing process of ore and rock under explosive load, explore the intrinsic relationship between ore and rock crushing characteristics and energy consumption, and provide quantitative data reference for subsequent field practice is an important problem that urgently needs to be solved. Summary of the Invention
[0005] To address the shortcomings of existing methods and their limitations in practical applications, and in order to achieve accurate simulation of the rock-breaking process under explosive loading, this invention provides a numerical simulation method for the rock-breaking process under FDEM explosive loading. This method aims to deeply analyze the intrinsic relationship between rock-breaking characteristics and energy consumption, providing quantitatively based data references for field engineering practice. The method includes the following steps: constructing a geometric model of the rock-breaking specimen using the finite element meshing software Gmsh; meshing the geometric model of the rock-breaking specimen to obtain a triangular mesh; and setting materials based on the geometric model of the rock-breaking specimen and the triangular mesh. The parameters and boundary conditions of the explosive load are used to simulate the explosive load application process to obtain explosive load simulation information. Based on the explosive load simulation information, the normal contact force calculation function and the tangential contact force calculation function are derived and established. The node coordinates and velocities of the triangular element are updated through the normal contact force calculation function and the tangential contact force calculation function to obtain the updated coordinate and velocity information of the triangular element nodes. Based on the coordinate and velocity information, a cloud map is generated, the crack situation is statistically analyzed, and the energy change is analyzed to obtain the crack propagation, stress field distribution, and energy evolution results of the ore and rock under the explosive load.
[0006] The results obtained by this invention regarding crack propagation, stress field distribution, and energy evolution of ore and rock under explosive loading can provide quantitative data references for field practice.
[0007] Optionally, the step of constructing a geometric model of the ore and rock specimen using the finite element meshing software Gmsh, and then meshing the geometric model of the ore and rock specimen to obtain a triangular mesh, includes: introducing an unstructured Delaunay triangulation algorithm; using the unstructured Delaunay triangulation algorithm to mesh the geometric model of the ore and rock specimen and inserting joint elements at the boundaries of adjacent elements to obtain a triangular mesh of the geometric model of the ore and rock specimen.
[0008] The geometry of mineral and rock specimens may contain irregular boundaries, holes, grooves, and other features. The unstructured Delaunay triangulation algorithm of this invention has strong adaptability and can fit complex geometric boundaries well to generate high-quality triangular meshes, ensuring the consistency between the mesh and the geometry of the mineral and rock specimens.
[0009] Optionally, setting material parameters and explosive load boundary conditions based on the geometric model of the rock specimen and the triangular mesh includes: obtaining the mechanical parameters of the triangular elements and joint elements based on the geometric model of the rock specimen and the triangular mesh; determining material parameters and explosive load boundary conditions based on the mechanical parameters, wherein the material parameters include elastic modulus, Poisson's ratio, density, tensile strength, cohesion, internal friction angle, fracture energy release rate of joint elements, normal penalty parameter, tangential penalty parameter, and penalty parameter of joint elements.
[0010] The material parameters and boundary conditions of this invention help determine key parameters, provide a theoretical basis for parameter selection and optimization in practical engineering, make the analysis results more accurate and reliable, and provide an information foundation for practical engineering.
[0011] Optionally, the step of simulating the application process of the explosive load based on the material parameters and the explosive load boundary conditions to obtain explosive load simulation information includes: calibrating the fracture energy according to the material parameters and the explosive load boundary conditions, wherein the fracture energy includes the Type I fracture energy release rate. and Type II fracture energy release rate A Duvall-type pressure function is introduced; the Duvall-type pressure function simulates the application process of the explosive load based on the fracture energy, the material parameters and the explosive load boundary conditions, and obtains the explosive load simulation information.
[0012] The simulation process of applying explosive loads in this invention can effectively analyze the initiation, propagation, and penetration of cracks inside the ore and rock, which helps to reveal the fragmentation mechanism of ore and rock under explosive loads, thereby understanding the process of ore and rock gradually breaking into blocks from a continuous medium.
[0013] Optionally, the explosion load simulation information satisfies the following relationship:
[0014] ,
[0015] in, For time Pressure at the location, The peak explosion pressure acting on the borehole wall. As a constant one, The constant is two. For time variables, for The normalized time parameter.
[0016] ,
[0017] in, for The normalized time parameter, The constant is two. It is a constant.
[0018] The formula of this invention not only considers the stage of pressure rising from zero to the peak, but also covers the process of pressure gradually decaying from the peak. By adjusting the constant value, the shape of the pressure rise and decay curves can be flexibly changed to adapt to the changing characteristics of explosion pressure under different types of explosives, different charge structures, and different environmental conditions, making the simulation results more comprehensive and accurate.
[0019] Optionally, deriving and establishing the normal contact force calculation function and the tangential contact force calculation function based on the explosion load simulation information includes: using MultiFracS software for FDEM solving and using the NBS contact algorithm for contact determination; dividing the solution region into triangular meshes based on the MultiFracS software, the NBS contact algorithm, and the explosion load simulation information, and determining the contact situation of triangular elements in adjacent meshes; and establishing the normal contact force calculation function and the tangential contact force calculation function based on the triangular meshes and the contact situation.
[0020] The present invention establishes a calculation function based on explosion load simulation information, which can fully consider dynamic characteristics and calculate the corresponding contact force according to the explosion pressure value at different times, so that the simulation results can more realistically reflect the dynamic response of the ore and rock under the explosion impact, including instantaneous changes in contact force, vibration and fluctuation, etc.
[0021] Optionally, establishing the normal contact force calculation function and the tangential contact force calculation function based on the triangular mesh and the contact condition includes:
[0022] The function for calculating the normal contact force satisfies the following relationship:
[0023] ,
[0024] in, For normal contact force, For normal penalty parameters, One of two adjacent triangular units It is the other one of two adjacent triangular units. For gradient, for Dot at The momentum, For overlapping regions and At a point within the intersection, for Dot at The momentum, For overlapping regions and At a point within the intersection, Let be the area of the overlapping portion of the triangular units. Let be the vector representing the direction of the outward normal to the overlapping portion. For contact triangle The outer boundary of the overlapping portion, For unit Potential function on, For unit Potential function on.
[0025] In practical engineering, contact boundaries are not regular and may have complex shapes such as curves and polygons. The function of this invention is based on triangular mesh for calculation, which can adapt well to various irregular boundaries, accurately describe the geometry of the contact boundary, and handle various complex contact scenarios.
[0026] Optionally, establishing the normal contact force calculation function and the tangential contact force calculation function based on the triangular mesh and the contact condition includes: setting the relative misalignment condition between two adjacent triangular elements based on Coulomb's friction law; and establishing the tangential contact force calculation function based on the relative misalignment condition, the triangular mesh, and the contact condition.
[0027] The tangential contact force calculation function satisfies the following relationship:
[0028] ,
[0029] in, for Tangential contact force at moment, For time variables, For time step, For the damping of the element node, It is the normal contact force.
[0030] The function established by this invention based on the relative misalignment condition can take into account the relative motion between two adjacent triangular units, and can more accurately describe the influence of complex relative motion on tangential contact force.
[0031] Optionally, updating the node coordinates and velocities of the triangular element using the normal contact force calculation function and the tangential contact force calculation function to obtain the updated coordinate and velocity information of the triangular element nodes includes: constructing a coordinate and velocity update model for the triangular element nodes based on Newton's second law, the normal contact force calculation function, and the tangential contact force calculation function; updating the node coordinates and velocities of the triangular element using the coordinate and velocity update model, and obtaining the updated coordinate and velocity information of the triangular element nodes.
[0032] The coordinate and velocity update model satisfies the following relationship:
[0033] ,
[0034] in, for The speed of time for The speed of time For total nodal force, For time step, For the quality of the node, for The distance of time, for The distance of time.
[0035] The model of this invention explicitly considers the influence of normal and tangential contact forces on node motion, and can comprehensively and accurately describe the force situation of nodes during the contact process, thereby more realistically simulating the dynamic behavior of objects.
[0036] Secondly, this invention also provides a numerical simulation system for the rock and ore fracturing process under FDEM explosive loading, which can efficiently execute the numerical simulation method for the rock and ore fracturing process under FDEM explosive loading provided by this invention. The system includes an input device, a processor, an output device, and a memory, wherein the input device, processor, output device, and memory are interconnected. The memory includes a computer-readable storage medium as described in the first aspect of this invention, and is used to store a computer program. The computer program includes program instructions, and the processor is configured to call the program instructions. The numerical simulation system for the rock and ore fracturing process under FDEM explosive loading provided by this invention has a compact structure, strong applicability, and greatly improves operating efficiency. Attached Figure Description
[0037] Figure 1 This is a flowchart of the numerical simulation method for the rock and ore fracturing process under FDEM explosive loading based on the present invention.
[0038] Figure 2 This is a schematic diagram of the FDEM principle of the present invention;
[0039] Figure 3 This is a schematic diagram of the geometric model of the mineral and rock specimen of the present invention;
[0040] Figure 4 This is a schematic diagram of the joint element fracture constitutive model of the present invention;
[0041] Figure 5 This is a schematic diagram illustrating the contact force conditions of the present invention;
[0042] Figure 6 This is a cloud map showing the displacement changes during the rock and ore fracturing process under explosive loading, as described in this invention.
[0043] Figure 7 This is a schematic diagram of the crack rose under explosive loading according to the present invention;
[0044] Figure 8 This is a schematic diagram of the numerical simulation system for the rock and ore fracturing process under FDEM explosive loading based on the present invention. Detailed Implementation
[0045] Specific embodiments of the present invention will now be described in detail. It should be noted that the embodiments described herein are for illustrative purposes only and are not intended to limit the invention. In the following description, numerous specific details are set forth in order to provide a thorough understanding of the invention. However, it will be apparent to those skilled in the art that these specific details are not necessary to practice the invention. In other instances, well-known circuits, software, or methods have not been specifically described to avoid obscuring the invention.
[0046] Throughout this specification, references to "an embodiment," "an embodiment," "an example," or "an example" mean that a particular feature, structure, or characteristic described in connection with that embodiment or example is included in at least one embodiment of the invention. Therefore, the phrases "in an embodiment," "in an embodiment," "an example," or "an example" appearing in various places throughout the specification do not necessarily refer to the same embodiment or example. Furthermore, specific features, structures, or characteristics can be combined in one or more embodiments or examples in any suitable combination and / or sub-combination. Moreover, those skilled in the art will understand that the illustrations provided herein are for illustrative purposes and are not necessarily drawn to scale.
[0047] Please see Figure 1 To accurately simulate the rock-breaking process of ore under explosive loading, deeply explore the intrinsic relationship between ore fracturing characteristics and energy consumption, and provide quantitative technology and data reference for field engineering practice, this invention provides a numerical simulation method for ore fracturing under explosive loading based on FDEM. The method includes the following steps:
[0048] S1. Construct a geometric model of the ore and rock specimen using the finite element meshing software Gmsh. Mesh the geometric model of the ore and rock specimen to obtain a triangular mesh. The implementation details are as follows:
[0049] The core idea of the FDEM method in this embodiment is to discretize the continuous medium structure into a finite element network composed of triangular elements. Please refer to the schematic diagram of the basic principle of FDEM. Figure 2 ,based on Figure 2 It is evident that this includes the topological structure of triangular elements and joint elements. For subsequent numerical simulation studies based on FDEM, the open-source finite element meshing software Gmsh was used in this embodiment to construct the geometric model of the mineral rock specimen.
[0050] In an optional embodiment, a standard cylinder is used as the geometric model of the ore specimen. The cylinder has a diameter of 50 mm and a height of 50 mm. To effectively simulate an actual blasting scenario, a borehole with a diameter of 6 mm and a depth of 30 mm is pre-drilled in the middle of the cylinder. Since subsequent simulations need to focus on the ore fragmentation process in a two-dimensional plane, a two-dimensional finite discrete element coupled model is ultimately established. Please refer to the schematic diagram of the model. Figure 3 .
[0051] To achieve mesh generation for the geometric model of the rock specimen while meeting the requirements for simulating crack initiation and propagation, the unstructured Delaunay triangulation algorithm was introduced in this embodiment. This algorithm can generate relatively regular triangular elements to ensure mesh quality and can also adapt to complex geometries. It has good adaptability to rock specimen models with special structures such as blast holes.
[0052] When meshing the geometric model of the ore and rock specimen using the unstructured Delaunay triangulation algorithm, joint elements are simultaneously inserted at the boundaries of adjacent elements. The introduction of these joint elements is crucial for the FDEM method to handle discontinuities in ore and rock, allowing the model to more accurately simulate the dynamic evolution of cracks, including but not limited to crack initiation, propagation, and interactions. An example of a joint element fracture constitutive model is provided in the documentation. Figure 4 .
[0053] Taking into account both simulation accuracy and computational efficiency, in an optional embodiment, the mesh size is set to 0.0006 mm. After being divided by the unstructured Delaunay triangulation algorithm, the entire geometric model of the ore specimen is discretized into 17,422 triangular elements, generating 8,931 nodes. The above mesh division result not only ensures that the micromechanical behavior of the ore under explosive loading can be effectively analyzed, but also controls the computational scale to a certain extent, thus improving the feasibility of the simulation.
[0054] S2. Based on the geometric model of the rock specimen and the triangular mesh, set the material parameters and explosive load boundary conditions. The implementation details are as follows:
[0055] The embodiments are mainly based on the constructed geometric model of the mineral rock specimen and the divided triangular mesh to obtain the mechanical parameters of the triangular elements and joint elements. These parameters are the basis for subsequent determination of material parameters and explosive load boundary conditions, and are crucial for accurately simulating the mechanical behavior of mineral rock under explosive load.
[0056] Material parameters and explosive load boundary conditions are determined based on mechanical parameters. In this embodiment, the material parameters mainly include elastic modulus, Poisson's ratio, density, tensile strength, cohesion, internal friction angle, fracture energy release rate of joint elements, normal penalty parameter, tangential penalty parameter, and penalty parameter of joint elements.
[0057] Furthermore, the aforementioned material parameters cover multiple aspects, which can be further categorized into the following three types in the FDEM method:
[0058] The mechanical parameters of a triangular element include its elastic modulus. Poisson's ratio and density These parameters can reflect the elastic properties and material density of the triangular unit itself, and play a key role in simulating the overall stiffness and deformation characteristics of the ore and rock.
[0059] The microscopic parameters of joint elements include tensile strength. Cohesion internal friction angle and the rate of energy release during joint unit fracture and The aforementioned microscopic parameters are mainly used to describe the mechanical properties of joint elements, and can reflect the influence of joint surfaces on the strength and deformation of ore and rock. They play an important role in simulating the crack initiation and propagation process.
[0060] The penalty parameters include the normal penalty parameters of the triangular elements. Tangential penalty parameters Penalty parameters of joint elements In the FDEM method, the penalty parameter is mainly used to handle contact and constraint problems, which can ensure the stability and accuracy of the model during the calculation process.
[0061] In this embodiment, the mechanical parameters of the triangular and joint elements can be obtained through relevant laboratory experiments. For example, parameters such as elastic modulus, Poisson's ratio, density, tensile strength, cohesion, and internal friction angle can be accurately measured using standard rock mechanics testing methods. Obtaining these parameters through laboratory experiments ensures their authenticity and reliability, providing a solid foundation for subsequent numerical simulations.
[0062] In one embodiment, the microscopic parameters in the FDEM model are determined based on the physical and mechanical parameters of the rock mass. Taking a certain ore rock as an example, the average saturated uniaxial compressive strength of the rock is known. It has a strength of 50.87 MPa and a density of Elastic modulus Poisson's ratio ,tensile strength Cohesion internal friction angle =22.30.
[0063] To ensure a high degree of agreement between the numerical simulation results and the mechanical test results, extensive trial and error and repeated verification were conducted to finally determine when... and The simulation yields the best results. Simultaneously, the normal penalty parameter should be set. With tangential penalty parameters All were 39.26 GPa, joint unit penalty parameters The Pa value is 3926 GPa. These parameter settings allow for a more accurate simulation of the rock and ore fracture process under explosive loading, providing a reliable reference for subsequent research and engineering practice.
[0064] S3. Simulate the application process of the explosive load based on material parameters and boundary conditions to obtain explosive load simulation information. Derive and establish the calculation functions for normal contact force and tangential contact force based on the explosive load simulation information. The specific implementation details are as follows:
[0065] First, the process of applying the explosive load is simulated based on material parameters and explosive load boundary conditions to obtain explosive load simulation information.
[0066] The fracture energy is determined based on material parameters and explosive load boundary conditions. The fracture energy includes the Type I fracture energy release rate. and Type II fracture energy release rate .
[0067] In this embodiment, the fracture energy is calibrated based on the determined material parameters and explosive load boundary conditions. The fracture energy primarily includes the Type I fracture energy release rate. and Type II fracture energy release rate The embodiment employs an adaptive FDEM method, which, without prior calibration to obtain optimal joint penalty parameter values, yields simulation results equivalent to those with optimal parameter values, significantly saving calibration time. However, the fracture energy still needs calibration. The main control is the tensile failure of the joint elements. The main control is the shear failure of the joint elements, and the two work together to determine the failure mode of the entire continuum element. Based on the above, it can be seen that after continuous trial and error and repeated verification, when... and At that time, the numerical simulation results and the mechanical test results showed a high degree of agreement.
[0068] A Duvall-type pressure function is introduced. This Duvall-type pressure function simulates the application process of explosive load based on fracture energy, material parameters, and explosive load boundary conditions, and obtains the explosive load simulation information.
[0069] Because applying boundary conditions to the explosion model requires secondary development of the software, the specific steps are as follows: Use Microsoft Visual Studio to rewrite the explosion load code in the native SBoundary.dll file of the MultiFracS software, regenerate the dll file, and then replace the old file with the newly generated dll file. This allows the modified explosion load to be recognized and successfully applied to the boundary conditions. Simultaneously, using a Duvall-type pressure function to simulate the explosion load application process accurately reflects the peak value and pressure decay process of the explosion load.
[0070] The above-mentioned explosion load simulation information, namely the peak explosion load and the pressure decay process, satisfies the following relationship:
[0071] ,
[0072] in, For time Pressure at the location, The peak explosion pressure acting on the borehole wall. It is a constant one. The constant is two. For time variables, for The normalized time parameter.
[0073] Furthermore, the normalized time parameter satisfies the following relationship:
[0074] ,
[0075] in, for The normalized time parameter, The constant is two. The constant is 1. This formula primarily changes the decay time ratio. Adjusting the explosion pressure rise and fall time, as well as the pressure peak, can simplify the complexity of the explosion pressure curve.
[0076] Then, based on the explosion load simulation information, the calculation functions for normal contact force and tangential contact force are derived and established.
[0077] By introducing MultiFracS software for FDEM solving and using the NBS contact algorithm for contact determination, the entire solution domain can be divided into triangular meshes, and the contact situation of triangular elements within adjacent meshes can be determined.
[0078] By dividing the solution domain into triangular meshes using MultiFracS software, NBS contact algorithm, and blast load simulation information, the contact situation of triangular elements within adjacent meshes can be determined.
[0079] For a more detailed diagram of the contact forces, please refer to [link / reference]. Figure 5 .
[0080] The calculation functions for normal contact force and tangential contact force are mainly established based on the triangular mesh and contact conditions.
[0081] The above-mentioned normal contact force is obtained by the surface integral of the gradient of the potential function in the overlapping region, that is, the above-mentioned normal contact force calculation function satisfies the following relationship:
[0082] ,
[0083] in, For normal contact force, For normal penalty parameters, One of two adjacent triangular units It is the other one of two adjacent triangular units. For gradient, for Dot at The momentum, For overlapping regions and At a point within the intersection, for Dot at The momentum, For overlapping regions and At a point within the intersection, Let be the area of the overlapping portion of the triangular units. Let be the vector representing the direction of the outward normal to the overlapping portion. For contact triangle The outer boundary of the overlapping portion, For unit Potential function on, For unit Potential function on.
[0084] The tangential contact force is calculated using the penalty function principle, and its formula is as follows:
[0085] ,
[0086] in, for Tangential contact force at moment, for Tangential contact force at moment, For tangential penalty parameters, for The displacement increment within a given time period This represents the current time step.
[0087] Based on Coulomb's law of friction, a relative misalignment condition is set between two adjacent triangular elements. That is, when the following conditions are met, two adjacent contacting triangular elements will be relatively misaligned and generate relative displacement.
[0088] ,
[0089] in, for Tangential contact force at moment, For normal contact force, The relevant friction angle.
[0090] Based on the aforementioned relative misalignment conditions, triangular mesh, and contact situation, a tangential contact force calculation function is established. When the relative misalignment conditions are met, the tangential contact force can be calculated using the following function:
[0091] ,
[0092] in, for Tangential contact force at moment, For time variables, For time step, For the damping of the element node, It is the normal contact force.
[0093] S4. Update the nodal coordinates and velocities of the triangular element using the normal contact force calculation function and the tangential contact force calculation function to obtain the updated coordinate and velocity information of the triangular element nodes. The specific implementation details are as follows:
[0094] In this embodiment, Newton's second law is used to solve the system control equations. Newton's second law can accurately describe the changes in the motion state of an object under the action of force, providing a theoretical basis for constructing a coordinate and velocity update model. The expression of the above control equations is as follows:
[0095] ,
[0096] Where M is the mass of the element node. Let x be the second partial derivative with respect to t. Let x be the first-order partial derivative with respect to t. For the damping of the element node, For total node force.
[0097] The aforementioned element node mass reflects the inertial influence of the node on motion; the second-order partial derivative represents the node's acceleration over time, reflecting the rate of change of the node's velocity; the first-order partial derivative is the node's velocity, describing the speed of the node's motion in space; element node damping reflects the resistance experienced by the node during motion, causing the node's motion to gradually decrease; the total node force is the comprehensive manifestation of various forces acting on the node, mainly including the contact forces on the node. Nodal forces in the deformation of triangular and jointed elements Nodal forces caused by external loads The adhesive force of joint units .
[0098] Furthermore, the aforementioned total nodal force consists of multiple parts, the specific analysis of which is as follows:
[0099] Contact force at the node Contact forces arise from the contact between nodes and other structures or elements. Under complex mechanical environments such as explosive loads, the presence of contact forces has a significant impact on the movement of nodes. Under the impact of an explosion, the nodes of adjacent triangular elements squeeze and rub against each other, generating contact forces.
[0100] When triangular and jointed elements deform, stress is generated within them. This stress is transmitted through the nodes, forming nodal forces. The greater the degree of deformation, the greater the nodal force is usually generated.
[0101] In practice, structures are often subjected to various external loads, such as gravity and explosive loads. These external loads act directly on the nodes, generating corresponding nodal forces. The shock wave from an explosion exerts enormous pressure on the structural nodes, thus creating nodal forces induced by the external load. .
[0102] Jointed elements are a special type of unit that may exist in a structure. They contain adhesive forces that generate bonding forces at the nodes, thus constraining the movement of the nodes. This is the adhesive force of the jointed element. .
[0103] Based on the above system control equations and the composition of total nodal forces, an explicit method was used to further construct a coordinate and velocity update model for triangular unit nodes. The explicit method has the advantages of high computational efficiency and ease of implementation, and is suitable for real-time updating of node coordinates and velocities.
[0104] The coordinate and velocity update model in this embodiment satisfies the following relationship:
[0105] ,
[0106] in, for The speed of time for The speed of time For total nodal force, For time step, For the quality of the node, for The distance of time, for The distance of time.
[0107] The above The speed at a given moment reflects the node's position. After a time step The speed of movement afterwards; The speed at a given moment is the speed at the node. The state of motion at any given moment; the total nodal force acts on the nodes. The resultant force of all forces on the node; the time step determines the update interval, and the smaller the time step, the higher the update accuracy; node quality is one of the basic attributes of a node; The distance at any given moment reflects the node's position in time. After a time step The subsequent positional change; The distance at time is the distance between nodes. The position and state at any given moment.
[0108] Finally, the coordinate and velocity update model is used to update the node coordinates and velocities of the triangular element, and the updated coordinate and velocity information of the triangular element nodes is obtained.
[0109] The coordinate and velocity update model described above is used to update the node coordinates and velocities of the triangular element. During the update process, the current time is considered. The node velocity, total node force, node mass, and time step are used to calculate the... The node velocity at time point. Then, based on... Calculate the node coordinates, node velocity, and time step at time t. The node coordinates at each moment are determined by progressively updating the coordinates and velocity information of the triangular element nodes. This information accurately reflects the changes in the motion state of the triangular element nodes under mechanical actions such as explosive loads, providing an important basis for subsequent structural analysis and research.
[0110] S5. Generate cloud maps based on coordinate and velocity information, statistically analyze crack conditions and energy changes to obtain the results of crack propagation, stress field distribution and energy evolution in the ore under explosive loading. The specific implementation content is as follows:
[0111] The visualization analysis content is as follows;
[0112] Import the VTK result file containing the calculation results into the visualization software ParaView. The VTK format has good versatility and data compatibility, and can completely save various information in the calculation process, providing accurate data support for subsequent visualization analysis.
[0113] I. Generating various cloud maps and vector fields:
[0114] Displacement field cloud map: The imported data is processed by ParaView software to generate a displacement field cloud map. The displacement field cloud map can intuitively show the magnitude and direction of displacement of the ore and rock at various locations under the action of explosive load. Different colors and gray areas represent different displacement values. By observing the cloud map, the deformation of the ore and rock can be clearly seen, and further information such as displacement and displacement direction can be obtained.
[0115] The displacement change cloud map of the ore and rock fracturing process under explosive load in this embodiment can be found in the following example. Figure 6
[0116] Stress cloud map: Stress cloud maps are generated using ParaView software. Stress cloud maps can reflect the stress distribution state inside the ore and rock, including but not limited to normal stress and shear stress. Different colored areas correspond to different stress levels. By analyzing the stress cloud map, high stress areas and low stress areas in the ore and rock can be identified, thereby judging the dangerous areas of ore and rock failure, and providing an important basis for the analysis of crack initiation and propagation.
[0117] Velocity vector field: Generating a velocity vector field can intuitively show the magnitude and direction of the velocity at each point in the ore under the action of an explosive load. The velocity vector can be represented in the form of an arrow. The length of the arrow can represent the magnitude of the velocity, and the direction of the arrow can represent the direction of the velocity. By observing the velocity vector field, we can understand the propagation of the explosive shock wave in the ore and the movement trend of the ore.
[0118] II. Extraction and Statistical Analysis of Crack Evolution Patterns
[0119] Extracting crack evolution patterns: Based on the visualization results and raw data generated above, the evolution patterns of cracks can be extracted. By observing the changes in the shape and location of cracks at different times, the initiation, propagation and penetration of cracks can be analyzed. It is possible to observe how cracks gradually develop from initial micro-defects into macro-cracks, and how cracks connect with other cracks during the propagation process.
[0120] Crack count: The extracted cracks are counted to accurately record the total number of cracks in the ore at different times. The change in the number of cracks can reflect the degree of damage to the ore by the explosion load. As the explosion load acts, the number of cracks usually increases gradually.
[0121] Statistical analysis of crack types and orientations: Cracks are classified and statistically analyzed, primarily into tensile cracks and shear cracks. Tensile cracks are caused by tensile stress on the ore, while shear cracks are caused by shear stress. By statistically analyzing the number and proportion of different types of cracks, we can further understand the main failure modes of ore under explosive loading. Simultaneously, the orientation distribution of cracks is statistically analyzed to determine whether the crack propagation direction exhibits a certain regularity or whether it concentrates along a specific direction.
[0122] The crack rose diagram in this embodiment illustrates the crack direction distribution characteristics. Please refer to [link / reference]. Figure 7 .
[0123] III. Analysis of Energy Evolution Laws
[0124] Analyzing total energy means analyzing how the total energy of the ore and rock changes under explosive loads. Total energy mainly includes the sum of various energy forms such as strain energy and kinetic energy. By plotting the curve of total energy change over time, we can understand the accumulation and release of energy during the explosion. In the early stages of the explosion, the total energy increases rapidly, and as the explosive load decays and the ore and rock are destroyed, the total energy gradually decreases.
[0125] Strain energy is the energy stored in ore and rock due to deformation. Analyzing the changes in strain energy can help us understand the degree of deformation and energy storage of ore and rock under explosive loads. When ore and rock undergo large deformation, strain energy will increase significantly, while some strain energy will be released when cracks initiate and propagate.
[0126] Kinetic energy is the energy possessed by rocks due to their movement. Under the action of explosive loads, rocks will move, and kinetic energy will be generated and change accordingly. Analyzing the changes in kinetic energy can help us understand the motion state of rocks and the influence of the propagation of explosive shock waves on the motion of rocks.
[0127] Strain energy specificity is the strain energy per unit volume that reflects the energy density inside the ore or rock. By analyzing the changes in strain energy specificity, we can gain a deeper understanding of the energy distribution inside the ore or rock and identify areas of concentrated energy. These areas may be high-risk areas for crack initiation and propagation.
[0128] Furthermore, the numerical simulation method of rock and ore fracturing process based on FDEM explosive loading is compared with existing technologies.
[0129] The method of this invention adopts an adaptive joint element insertion mechanism, which differs from the traditional FDEM method of pre-inserting joint elements globally. The traditional method pre-inserts a large number of joint elements, which consumes a lot of computing resources during the calculation process, resulting in excessively long calculation time. However, the method of this invention can dynamically insert joint elements according to the actual stress and failure characteristics of the ore and rock, avoiding unnecessary computational waste, thereby significantly reducing the calculation time and improving the calculation efficiency.
[0130] The adaptive FDEM method of this invention eliminates the need for pre-obtaining optimal joint penalty parameter values through complex calibration. In traditional methods, the selection of joint penalty parameter values significantly impacts the accuracy of simulation results, requiring extensive experimentation and calibration to determine the optimal values—a cumbersome and time-consuming process. The method of this invention can obtain simulation results equivalent to those with optimal parameter values without prior acquisition, simplifying the parameter acquisition process while ensuring the accuracy of the simulation results.
[0131] This invention combines an improved joint element fracture constitutive model with the NBS contact algorithm, enabling more accurate simulation of crack initiation, propagation, and penetration. The improved joint element fracture constitutive model can more accurately describe the mechanical behavior of joint elements under stress, while the NBS contact algorithm can better handle the contact and interaction between joint elements, thus making the crack simulation more consistent with reality and providing more reliable data support for in-depth research on the failure mechanism of minerals and rocks.
[0132] The method of this invention employs a Duvall-type pressure function and a secondary development interface, which can more realistically simulate the rise and decay process of explosion pressure. The aforementioned Duvall-type pressure function can more accurately describe the change law of explosion pressure over time, while the secondary development interface can further adjust and optimize the pressure function according to actual needs, thereby improving the realism of explosion pressure simulation and providing a guarantee for more accurate analysis of the effect of explosion load on ore and rock.
[0133] The method of this invention provides a comprehensive analysis function from displacement field, stress field to crack statistics and energy evolution. Based on the analysis function, it is possible to fully understand the various responses of ore and rock under explosive load, including but not limited to deformation, stress distribution, crack propagation and energy change. The above analysis results can provide a quantitative basis for the selection of blasting parameters on site, help optimize blasting schemes, improve blasting efficiency and accuracy, and reduce the impact of blasting on the surrounding environment.
[0134] Please see Figure 8In an optional embodiment, the present invention also provides a numerical simulation system for the rock and ore fracturing process under FDEM explosive loading. This system includes a processor, an input device, an output device, and a memory, all interconnected. The memory stores a computer program, which includes program instructions. The processor is configured to call the program instructions and execute the specific steps of the numerical simulation method and related embodiments of the rock and ore fracturing process under FDEM explosive loading provided by the present invention. The numerical simulation system for the rock and ore fracturing process under FDEM explosive loading of the present invention is structurally complete and objectively stable.
[0135] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some or all of the technical features therein. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention, and they should all be covered within the scope of the claims and specification of the present invention.
Claims
1. A numerical simulation method for the rock and ore fracturing process under FDEM explosive loading, characterized in that, Includes the following steps: The geometric model of the mineral and rock specimen was constructed using the finite element meshing software Gmsh, and the geometric model of the mineral and rock specimen was meshed to obtain a triangular mesh of the geometric model of the mineral and rock specimen. Material parameters and explosive load boundary conditions are set based on the geometric model of the mineral and rock specimen and the triangular mesh; Based on the material parameters and the explosion load boundary conditions, the explosion load application process is simulated to obtain explosion load simulation information. Based on the explosion load simulation information, the normal contact force calculation function and the tangential contact force calculation function are derived and established. The node coordinates and velocities of the triangular element are updated using the normal contact force calculation function and the tangential contact force calculation function to obtain the updated coordinate and velocity information of the triangular element nodes. Based on the coordinate and velocity information, a cloud map is generated, the crack situation is statistically analyzed, and the energy changes are analyzed to obtain the crack propagation, stress field distribution, and energy evolution results of the ore and rock under the action of explosive load. The process of constructing a geometric model of the mineral and rock specimen using the finite element meshing software Gmsh, and then meshing the geometric model of the mineral and rock specimen to obtain a triangular mesh, includes: Introduce the unstructured Delaunay triangulation algorithm; The unstructured Delaunay triangulation algorithm is used to mesh the geometric model of the mineral rock specimen and joint elements are inserted at the boundaries of adjacent elements to obtain a triangular mesh of the geometric model of the mineral rock specimen. The process of simulating the application of explosive load based on the material parameters and the explosive load boundary conditions to obtain explosive load simulation information includes: The fracture energy is calibrated based on the material parameters and the explosive load boundary conditions. The fracture energy includes the Type I fracture energy release rate. and Type II fracture energy release rate ; Introduce a Duvall-type pressure function; The Duvall-type pressure function simulates the application process of the explosive load based on the fracture energy, the material parameters, and the explosive load boundary conditions, and obtains the explosive load simulation information. Based on the triangular mesh and contact conditions, the calculation functions for normal contact force and tangential contact force are established, including: The relative misalignment condition between two adjacent triangular elements is set based on Coulomb's friction law; A tangential contact force calculation function is established based on the relative misalignment condition, the triangular mesh, and the contact situation; The tangential contact force calculation function satisfies the following relationship: , in, for Tangential contact force at moment, For time variables, For time step, For the damping of the element node, It is the normal contact force.
2. The numerical simulation method for ore and rock fracturing process under FDEM explosive loading according to claim 1, characterized in that, The setting of material parameters and explosive load boundary conditions based on the geometric model of the rock specimen and the triangular mesh includes: The mechanical parameters of the triangular elements and joint elements are obtained based on the geometric model of the mineral rock specimen and the triangular mesh. The material parameters and explosive load boundary conditions are determined based on the mechanical parameters. The material parameters include elastic modulus, Poisson's ratio, density, tensile strength, cohesion, internal friction angle, fracture energy release rate of joint elements, normal penalty parameter, tangential penalty parameter, and penalty parameter of joint elements.
3. The numerical simulation method for rock and ore fracturing process under FDEM explosive loading according to claim 1, characterized in that, The explosion load simulation information satisfies the following relationship: , in, For time Pressure at the location, The peak explosion pressure acting on the borehole wall. As a constant one, The constant is two. For time variables, for The normalized time parameter, , in, for The normalized time parameter, The constant is two. It is a constant.
4. The numerical simulation method for rock and ore fracturing process under FDEM explosive loading according to claim 1, characterized in that, The process of deriving and establishing the normal contact force calculation function and the tangential contact force calculation function based on the explosion load simulation information includes: MultiFracS software was introduced to solve FDEM, and the NBS contact algorithm was used for contact determination. Based on the MultiFracS software, the NBS contact algorithm, and the explosion load simulation information, the solution domain is divided into triangular meshes, and the contact situation of triangular elements in adjacent meshes is determined. Based on the triangular mesh and the contact conditions, establish the normal contact force calculation function and the tangential contact force calculation function.
5. The numerical simulation method for rock and ore fracturing process under FDEM explosive loading according to claim 4, characterized in that, The step of establishing the normal contact force calculation function and the tangential contact force calculation function based on the triangular mesh and the contact condition includes: The function for calculating the normal contact force satisfies the following relationship: , in, For normal contact force, For normal penalty parameters, One of two adjacent triangular units It is the other one of two adjacent triangular units. For gradient, for Dot at The momentum, For overlapping regions and At a point within the intersection, for Dot at The momentum, For overlapping regions and At a point within the intersection, Let be the area of the overlapping portion of the triangular units. Let be the vector representing the direction of the outward normal to the overlapping portion. For contact triangle The outer boundary of the overlapping portion, For unit Potential function on, For unit Potential function on.
6. The numerical simulation method for rock and ore fracturing process under FDEM explosive loading according to claim 1, characterized in that, The process of updating the node coordinates and velocities of the triangular element using the normal contact force calculation function and the tangential contact force calculation function to obtain the updated coordinate and velocity information of the triangular element nodes includes: Based on Newton's second law, the normal contact force calculation function, and the tangential contact force calculation function, a coordinate and velocity update model for the triangular element nodes is constructed. The coordinates and velocities of the triangular element nodes are updated using the coordinate and velocity update model, and the updated coordinate and velocity information of the triangular element nodes is obtained. The coordinate and velocity update model satisfies the following relationship: , in, for The speed of time, for The speed of time, For total nodal force, For time step, For the quality of the node, for The distance of time, for The distance of time.
7. A numerical simulation system for the rock and ore fracturing process under FDEM explosive loading, characterized in that, The system includes a processor, an input device, an output device, and a memory, which are interconnected. The memory is used to store a computer program, which includes program instructions. The processor is configured to call the program instructions to execute the numerical simulation method for the rock and ore fracturing process under FDEM explosive loading as described in any one of claims 1-6.