An underground engineering blasting simulation method based on ABAQUS secondary development

Abstract

The application relates to the technical field of underground engineering and geotechnical mechanics, and provides an underground engineering blasting simulation method based on secondary development of ABAQUS, which comprises the following steps: step 1, a geometric model is established, and a calculation region is divided into an Euler domain and a Lagrange domain; step 2, a JH-2 constitutive model is called and realized through a UMAT user subroutine interface of ABAQUS; step 3, a JWL state equation is used for modeling in the Euler domain; step 4, a coupled Euler-Lagrange method is applied; step 5, boundary conditions and blasting loading parameters of the geometric model are set; and step 6, transient dynamics calculation is carried out, and results are output for analyzing blasting effects. The scheme not only can provide a high-precision numerical test platform for underground engineering blasting design and safety control, reduces engineering risks, optimizes blasting parameter configuration, but also has remarkable engineering application value and popularization significance.

Description

Technical Field

[0001] This invention relates to the fields of underground engineering and geotechnical mechanics, specifically to an underground engineering blasting simulation method based on ABAQUS secondary development. Background Technology

[0002] In underground engineering construction, blasting excavation is one of the key processes in tunnel excavation, mining, and underground utility tunnel construction. The blasting process involves transient high-pressure impact loads, nonlinear dynamic failure of the rock mass, flow of blast products, and their complex coupling with the surrounding rock. Its mechanical mechanism exhibits significant multi-field coupling and strong nonlinear characteristics. To effectively assess construction safety risks, optimize blasting design parameters, and accurately predict the dynamic response of the surrounding rock, numerical simulation has become an important technical means in research and engineering applications.

[0003] Currently, numerical simulations primarily rely on finite element or finite difference methods. While the traditional Lagrangian method can effectively describe the stress and strain of solid structures, it is prone to mesh distortion under large deformation and material flow conditions, leading to a significant decrease in computational accuracy. The Eulerian method, while capable of handling fluid flow, is insufficient in characterizing solid boundaries and struggles to accommodate complex boundaries and material nonlinearities. In recent years, the Coupled Eulerian-Lagrangian (CEL) method has been increasingly applied to blasting numerical simulations, addressing these issues to some extent, but it still falls short in terms of rock damage evolution and detailed characterization of blast debris. Regarding material constitutive modeling, traditional elastoplastic models or simplified damage models fail to reflect the nonlinear mechanical response of rocks under high strain rate and high-pressure impact conditions. The JH-2 model, designed for brittle materials, can describe strength degradation, equations of state, and damage evolution processes, but it is not directly built into general-purpose finite element software and requires secondary development. Furthermore, in describing blast loads, commonly used ideal gas equations of state or simplified pressure models struggle to accurately simulate the expansion and energy release processes of detonation products, resulting in significant discrepancies between calculated and experimental results. In contrast, the Jones-Wilkins-Lee (JWL) equation of state can more accurately characterize the pressure-volume relationship of explosive detonation products and has been widely verified in explosion mechanics research, but its combined application with complex rock mass material models is still immature.

[0004] In summary, existing numerical simulation techniques for blasting have limitations in handling coupling effects, constitutive description of rock mass damage, and simulation of detonation products, making it difficult to meet the accurate prediction requirements under complex underground engineering conditions. Therefore, it is necessary to propose an underground engineering blasting simulation method based on ABAQUS secondary development and combining the JH-2 constitutive model and the JWL blasting model to achieve high-precision numerical simulation of the dynamic failure mechanism and blasting effects of rock masses. Summary of the Invention

[0005] To address the problems existing in the background technology, this invention proposes an underground engineering blasting simulation method based on ABAQUS secondary development. This method can provide a high-precision numerical test platform for underground engineering blasting design and safety control, reduce engineering risks, optimize blasting parameter configuration, and has significant engineering application value and promotion significance.

[0006] To achieve the above objectives, the present invention adopts the following solution:

[0007] A method for simulating underground engineering blasting based on ABAQUS secondary development includes the following steps:

[0008] Step 1: Establish a geometric model that includes the surrounding rock, support structure and blasting zone, and perform mesh generation. Divide the computational domain into an Eulerian domain for simulating blast products and air and a Lagrange domain for simulating the surrounding rock and support structure.

[0009] Step 2: Through the UMAT user subroutine interface of ABAQUS, call and implement the JH-2 constitutive model to describe the strength change, volume compression and damage evolution process of rock under blasting load.

[0010] Step 3: In the Eulerian domain, the JWL equation of state is used to model the pressure-volume relationship and energy release law of the explosive detonation products;

[0011] Step 4: Apply the coupled Eulerian-Lagrange method to simulate the coupled interaction between the detonation products and the surrounding rock;

[0012] Step 5: Set the boundary conditions and blasting loading parameters of the geometric model, wherein the boundary conditions include applying non-reflective boundary conditions to the far-field boundary of the geometric model to eliminate stress wave reflection;

[0013] Step 6: Perform transient dynamic calculations and output the results for analyzing the blasting effect. The output results include the surrounding rock stress field, damage variable distribution, velocity field, and energy distribution.

[0014] Optionally, in step 2, the JH-2 constitutive model includes a material strength model, an equation of state (EOS), and a damage model, and is implemented by embedding the expressions of the material strength model, the equation of state (EOS), and the damage model into a user subroutine of ABAQUS.

[0015] Optionally, the material strength model expression is:

[0016] ,

[0017] In the formula, To normalize the complete strength, To normalize the fracture strength, Let be the damage variable, where 0 ≤ D ≤1;

[0018] The expression for the state equation EOS is:

[0019] ,

[0020] In the formula, m = r / r 0-1, r For the current density, r 0 represents the initial density. K 1, K 2, K 3 represents the bulk modulus parameter of the material;

[0021] The evolution law of the damage variable D in the damage model is expressed as follows:

[0022] ,in, ,

[0023] In the formula, D 1 and D 2 is a material constant, Δ e p For each calculation step, the equivalent plastic strain increment, e f To disrupt the equivalent plastic strain, P ∗ =P / P HEL For normalization pressure, T ∗ =T / P HEL Normalized tensile pressure.

[0024] Optionally, in step 3, the mathematical expression of the JWL state equation is:

[0025] ,

[0026] In the formula, The pressure of the detonation products; V / V 0 represents the ratio of the current volume to the initial volume; A , B , R 1, R 2, oh This is an empirical constant for explosives; E The specific internal energy of the detonation products.

[0027] Optionally, in step 4, the Eulerian domain and the Lagrange domain achieve bidirectional coupling through an interface contact algorithm; the solution process is based on the mass, momentum and energy conservation equations, and the Lagrange-Euler algorithm is used to ensure numerical stability.

[0028] Optionally, the mass conservation equation is expressed as:

[0029] ,

[0030] In the formula, r For fluid density, v Let velocity be the vector; the momentum conservation equation is:

[0031] ,

[0032] In the formula, σ is the stress tensor; the energy conservation equation is:

[0033] ,

[0034] In the formula, E For total energy density, P For pressure.

[0035] Optionally, in step 5, the mathematical expression for the non-reflective boundary condition is:

[0036] ,

[0037] In the formula, s n For boundary normal stress, r For the density of the medium, c For wave speed, v n The boundary normal velocity.

[0038] Optionally, in step 5, setting the detonation loading parameters includes: defining the initial pressure field of the detonation products in the Eulerian domain based on the JWL equation of state, and applying Chapman-Jouguet state parameters, including the initial density, at the detonation center. r 0. Detonation speed D CJ pressure P CJ and specific energy E 0.

[0039] Optionally, the initial CJ ​​condition applied at the explosion center is expressed by the formula:

[0040] ,

[0041] In the formula, c The thermal index is the thermal index of the explosive product.

[0042] Optionally, the transient dynamics calculation is performed using the ABAQUS / Explicit solver based on an explicit integration algorithm through small time step iterations; the output results incorporate an energy balance relationship:

[0043] ,

[0044] In the formula, The total energy released by the detonation. The internal energy of the detonation products. The kinetic energy of the surrounding rock and detonation gases. The strain energy absorbed by the surrounding rock.

[0045] The beneficial effects of this invention are as follows: First, this scheme achieves high-precision and high-fidelity simulation of the entire blasting process. Specifically, through the secondary development of the JH-2 constitutive model, it can realistically reproduce the strength degradation, volume compression, and damage evolution of rocks under high strain rate and high pressure conditions, thereby significantly improving the accuracy and reliability of underground engineering blasting simulation. Through the application of the JWL equation of state, accurate modeling of the expansion of detonation products, shock wave propagation, and their interaction with the surrounding rock is achieved, making the numerical simulation of blasting loads closer to the real working conditions, thereby effectively improving the reliability and accuracy of blasting process prediction.

[0046] Moreover, by setting the boundary conditions and detonation loading parameters of the geometric model, non-physical reflection of waves in finite domain simulation is avoided, and the authenticity of the detonation product loading process is ensured. By monitoring and verifying the calculation process through the energy balance relationship, it is proved that the simulation results conform to the law of conservation of energy, and the rationality of the calculation process and the credibility of the results are ensured from the physical principle.

[0047] In addition, this scheme effectively solves the coupling problem between the large deformation flow of detonation products and the nonlinear failure of the surrounding rock structure: the coupled Eulerian-Lagrange (CEL) method is adopted, which not only utilizes the advantages of the Eulerian method in handling the large deformation flow of detonation products, but also leverages the strengths of the Lagrange method in characterizing solid stress damage. Furthermore, the interaction is realized through the interface contact algorithm, which ensures that the dynamic loading effect of the detonation gas on the surrounding rock can be truly transmitted, effectively avoiding the problems of mesh distortion or unclear boundary characterization under a single method. Attached Figure Description

[0048] Figure 1 This is a flowchart of the method of the present invention;

[0049] Figure 2 This is a flowchart illustrating the implementation of the JH-2 constitutive model in ABAQUS in this embodiment of the invention.

[0050] Figure 3The diagrams shown are model diagrams in embodiments of the present invention, where a is a three-dimensional numerical calculation model diagram and b is a typical cross-sectional calculation model diagram. Detailed Implementation

[0051] To make the present invention clearer and more understandable, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the given embodiments are merely one implementation method and do not represent all embodiments.

[0052] Example 1

[0053] Combination Figure 1-Figure 2 This embodiment provides a method for simulating underground engineering blasting based on ABAQUS secondary development, including the following steps:

[0054] Step 1: Establish a geometric model including the surrounding rock, support structure, and blasting zone, and mesh the computational domain into Eulerian and Lagrangian domains. The Eulerian domain is mainly used to simulate explosion products and flowable substances such as air, and can handle high-speed expansion and large-scale flow problems; the Lagrangian domain is mainly used to simulate the surrounding rock and support structure, and is suitable for characterizing the stress-strain relationship and damage process of solid materials.

[0055] Through reasonable domain partitioning and grid settings, specifically, based on the JWL pressure decay and wave mechanics criteria, a "three-domain, one-zone" partition is established: the range where P ≤ 0.1·PCJ or Ma ≤ 1 is defined as the Eulerian kernel domain, and a value of 0.5–1.0 r is set. core A thick coupling transition zone; the outer edge of the Lagrange domain is at least c units away from the explosion center. p ·T win It also includes a non-reflective absorption band. The grid scale follows Δx. E ≤min{0.1·d hole , r core / 40}、Δx L ≤λ s Constraints such as / 10 and CFL≈0.3 were applied, and ALE adaptive smoothing and local refinement were performed under the |▽P| / P threshold triggering. This ensured the accuracy of the detonation product flow resolution while precisely capturing the JH-2 damage evolution of the surrounding rock. The established geometric model can accurately simulate the dynamic behavior of the detonation products and realistically reflect the mechanical response of the surrounding rock structure under blasting loads, laying the foundation for subsequent coupled calculations.

[0056] Step 2: Through the UMAT user subroutine interface of ABAQUS, call and implement the JH-2 constitutive model to describe the strength change, volume compression and damage evolution process of rock under blasting load.

[0057] Since the ABAQUS software does not have a built-in JH-2 model suitable for brittle rocks, this embodiment uses the UMAT user subroutine interface of ABAQUS to independently write code and call the JH-2 constitutive model. The JH-2 constitutive model is divided into a material strength model, an equation of state (EOS), and a damage model. The expressions of the material strength model, the equation of state (EOS), and the damage model are embedded into the ABAQUS user subroutine to comprehensively describe the dynamic response of brittle rocks under blasting.

[0058] Specifically, the material strength model is used to characterize the strength change process of rocks under high-pressure impact conditions. The normalized equivalent stress in the JH-2 constitutive model can be expressed as:

[0059] ,

[0060] In the formula, To normalize the complete strength, To normalize the fracture strength, For damage variables (0≤D≤1); when D=0 The materials were complete at the time. D=1 The material completely fails at that time. Normalized strength is generally expressed in the form:

[0061] ,

[0062] In the formula, P ∗ =P / P HEL For normalization pressure, T ∗ =T / P HEL To normalize the tensile stress, Here, A, C, and N represent the normalized strain rate and material constants.

[0063] The two formulas above represent two levels of the same material strength model (JH-2), jointly characterizing the changes in strength with pressure / tension, strain rate, and damage evolution. The normalized equivalent stress formula provides the damage weighting for the current (normalized) equivalent strength, while the normalized strength formula provides the parameterized form of the strength surface. The two are used in a coupled manner. The normalized strength formula defines the strength surface, and the normalized equivalent stress formula interpolates damage between the two strength surfaces, jointly describing the evolution of strength with load state and damage.

[0064] The equation of state (EOS) is used to reflect the hydrostatic pressure change of the rock under volumetric strain, and its mathematical expression is:

[0065] ,

[0066] In the formula, m= r / r 0-1, r For the current density, r 0 represents the initial density. K 1, K 2, K 3 represents the bulk modulus parameter of the material. When the material is under tension, take... P=K 1 m .

[0067] The damage model is used to describe the process of rock gradually evolving from an intact state to failure during plastic deformation, where the evolution law of damage variable D can be expressed as:

[0068] ,in, ,

[0069] In the formula, D 1 and D 2 is a material constant, Δ e p For each calculation step, the equivalent plastic strain increment, e f To disrupt the equivalent plastic strain. As the equivalent plastic strain accumulates, the damage gradually increases until... D=1 The material has completely failed.

[0070] By embedding the above theoretical formulas into the user subroutine of ABAQUS, the secondary development of the JH-2 constitutive model was realized. It can realistically reproduce the strength degradation, volume compression and damage evolution of rocks under high strain rate and high pressure conditions, thereby greatly improving the accuracy and reliability of underground engineering blasting simulation.

[0071] Step 3: In the Eulerian domain, the JWL equation of state is used to model the pressure-volume relationship and energy release law of the explosive detonation products. The JWL equation of state is an empirical equation that can characterize the expansion characteristics of detonation products under high pressure and high temperature. Its mathematical expression is:

[0072] ,

[0073] In the formula, The pressure of the detonation products; V / V 0 represents the ratio of the current volume to the initial volume; A , B , R 1, R 2, oh This is an empirical constant for explosives; E Let be the specific internal energy of the detonation products. Among them, the parameters are... A , B, R 1, R 2, oh This can be determined through cylindrical expansion experiments or hydrodynamic calculations. The equation simultaneously reflects the mechanical behavior of detonation products under initial high pressure and after expansion into an approximately ideal gas state. Specifically, in the initial stage of detonation (i.e., near the Chapman–Jouguet point), the first exponential term… Second index term The dominant pressure of the detonation products is characterized by a rapid decay process under high pressure. However, after the detonation products expand significantly, the third term ωE / V plays a major role, at which point the detonation products exhibit expansion characteristics similar to those of an ideal gas.

[0074] In ABAQUS / Explicit simulations, by inputting the parameters and initial conditions (such as detonation velocity, initial density, and specific energy) of the above equations into the Eulerian domain, accurate modeling of the expansion of detonation products, shock wave propagation, and their interaction with the surrounding rock can be achieved. This embodiment utilizes this advantage of the JWL equation of state to make the numerical simulation of blasting loads closer to real working conditions, thereby effectively improving the reliability and accuracy of blasting process prediction.

[0075] Step 4: Apply the Coupled Eulerian-Lagrange (CEL) method to simulate the coupled interaction between detonation products and surrounding rock.

[0076] Within the CEL framework, the Lagrangian domain is used to simulate solid materials such as surrounding rock and support structures. Its nodes move with the material's motion, characterizing the stress-strain evolution and damage process of solids. The Eulerian domain is used to simulate detonation products and flowing materials such as air. Its mesh is fixed, while the material can flow freely within the mesh, thus avoiding mesh instability under high-speed flow and large deformation. The two domains interact through an interface contact algorithm, ensuring that the dynamic loading effect of the detonation gas on the surrounding rock is realistically transmitted.

[0077] The fundamental governing equations of the CEL method are the mass conservation equation, momentum conservation equation, and energy conservation equation. In the Eulerian domain, the mass conservation equation can be expressed as:

[0078] ,

[0079] In the formula, r For fluid density, v Let velocity be the vector. The momentum conservation equation is:

[0080] ,

[0081] in, s Let be the stress tensor. The energy conservation equation is:

[0082] ,

[0083] In the formula, E For total energy density, P For pressure.

[0084] During the calculation, ABAQUS / Explicit employs the Arbitrary Lagrangian–Eulerian (ALE) algorithm to achieve relative slip between the mesh and the material, ensuring numerical stability. This method allows the high-speed expansion flow of detonation products to evolve freely in the Eulerian domain, while the surrounding rock in contact with it responds to stress and damage in the Lagrangian domain.

[0085] This embodiment utilizes the CEL method to achieve realistic coupling between detonation products and surrounding rock: on the one hand, the impact pressure of the detonation gas can accurately act on the surface of the surrounding rock and induce its damage evolution; on the other hand, the deformation and failure of the surrounding rock, in turn, changes the flow path and pressure distribution of the detonation products. Through bidirectional coupling simulation, this embodiment can reproduce the complete physical process of "detonation product flow - rock crack propagation - stress wave propagation" during underground engineering blasting, thereby significantly improving the accuracy and reliability of numerical simulation.

[0086] Step 5: Set the boundary conditions and blasting loading parameters of the geometric model, wherein the boundary conditions include applying non-reflective boundary conditions to the far-field boundary of the geometric model to eliminate stress wave reflection.

[0087] This step is crucial for ensuring the accuracy and stability of the numerical simulation. Because high-speed shock waves and stress waves propagate during underground engineering blasting, improper boundary conditions can easily lead to non-physical reflections at the computational domain boundaries, thus interfering with the numerical results of the blasting process. Therefore, this embodiment eliminates wave reflection effects by applying appropriate boundary constraints and introducing absorbing boundary conditions when establishing the finite element model.

[0088] Specifically, for the model boundary far from the explosion source, this embodiment uses a non-reflecting boundary (NRB). This means that by introducing a damping term on the boundary that matches the incident wave, the outgoing wave can leave the computational domain in the same way as the infinite domain. Its mathematical expression can be written as:

[0089] ,

[0090] In the formula, s n For boundary normal stress, r For the density of the medium, c For wave speed, v nThe velocity is the normal velocity at the boundary. This non-reflective boundary condition ensures the conservation of momentum for both the incident and outgoing waves at the boundary, thus effectively absorbing stress waves.

[0091] For the input of the blast load, this embodiment introduces the initial pressure field of the explosive detonation products in the Eulerian domain through the JWL equation of state, and simultaneously applies Chapman-Jouguet (CJ) state parameters, including the initial density, at the detonation center. r 0. Detonation speed D CJ pressure P CJ and specific energy E 0. The initial CJ ​​condition applied at the explosion center is expressed as:

[0092] ,

[0093] In the formula, c The thermal index is the thermal index of the explosive product.

[0094] By setting the boundary conditions and loading parameters as described above, non-physical reflection of waves in finite domain simulation is avoided, and the authenticity of the detonation product loading process is ensured, thereby significantly improving the stability and accuracy of underground engineering blasting simulation.

[0095] Step 6 involves performing transient dynamic calculations and outputting the results for analyzing the blasting effect. The output includes the surrounding rock stress field, damage variable distribution, velocity field, and energy distribution. This step, by combining stress wave propagation theory and the principle of energy conservation, achieves high-precision numerical reproduction of the entire underground engineering blasting process, significantly improving the scientific rigor and operability of blasting design and safety assessment.

[0096] Specifically, after setting the boundary conditions and loading parameters, the ABAQUS / Explicit solver is used to perform transient dynamic calculations on the model, obtaining the stress field distribution, strain evolution, crack propagation, and detonation product flow patterns during the underground engineering blasting process. During the calculation, the solver uses an explicit integration algorithm, employing small time-step iterations to achieve stable solutions for the high strain rate dynamic process, thus ensuring the accuracy and convergence of the blasting simulation.

[0097] Regarding the output results, this is based on the propagation laws of stress waves. The propagation speed of stress waves in a medium is related to the mechanical parameters of the medium, and the basic relationship is as follows:

[0098] ,

[0099] In the formula, c For the propagation speed of longitudinal waves, E The elastic modulus of the material. rThe density of the material is given. When a detonation wave acts on the surrounding rock, the wave velocity determines the propagation path of the stress wavefront and the rate of energy transfer.

[0100] The transient response of the surrounding rock under blasting is solved by using ABAQUS / Explicit three-dimensional explicit dynamics and CEL coupling. Specifically, the numerical calculation is performed by ABAQUS / Explicit three-dimensional CEL explicit dynamics. This includes applying the JWL equation of state to the detonation products and the JH-2 constitutive model to the surrounding rock. Through constitutive coupling, non-reflective boundaries, and adaptive time steps, the principal stress / equivalent stress and damage field are output at each time step, thereby obtaining stress contour lines and determining stress concentration areas and crack initiation locations.

[0101] Simultaneously, the energy conservation principle is introduced into the output results to verify the rationality of the calculation. The detonation energy can be decomposed into three parts during propagation: internal energy, kinetic energy, and strain energy. Its energy balance relationship can be expressed as:

[0102] ,

[0103] In the formula, The total energy released by the detonation. The internal energy of the detonation products. The kinetic energy of the surrounding rock and detonation gases. The strain energy absorbed by the surrounding rock. The energy evolution curve over time can intuitively reflect the law of energy transfer from detonation products to the surrounding rock structure during the blasting process.

[0104] Furthermore, the method described in this embodiment can output results such as crack distribution, velocity field, damage variable cloud map, and stress-strain curve, forming a complete dataset for blasting dynamics analysis. By comparing the numerical simulation results with field test data, the reliability of the model can be further verified, and it can be used to optimize the design of blasting parameters, such as charge quantity, borehole layout, and detonation sequence, thereby providing a scientific basis for underground engineering blasting.

[0105] Simulation Example:

[0106] To verify the effectiveness and accuracy of the underground engineering blasting simulation method based on ABAQUS secondary development proposed in this invention, this embodiment uses this method to carry out finite element numerical simulation analysis step by step, as follows:

[0107] First, a model is established. Specifically, a large underground railway station and surrounding rock geometry model are created in the ABAQUS finite element platform. The underground station is 18m wide, 18m high, and 20m deep. The surrounding rock model has dimensions of 60m × 60m × 110m to avoid boundary effects. The computational domain is meshed, dividing the region into Eulerian and Lagrange domains. For example... Figure 3 As shown.

[0108] Second, in ABAQUS, the JH-2 constitutive model is called via the UMAT subroutine. The execution logic of the UMAT subroutine is as follows: Figure 2 As shown.

[0109] Third, the detonation products are modeled using the JWL equation of state:

[0110] ,

[0111] Taking TNT explosives as an example, take the parameters A =524GPa, B =9.5GPa, R 1 = 4.2 R 2 = 1.1, oh =0.3. Detonation velocity D =6930m / s, initial density r 0 = 1630 kg / m 3 .

[0112] Fourth, the Coupled Eulerian-Lagrange (CEL) method is used to solve the interaction between the Eulerian and Lagrange domains. The Eulerian domain has a fixed mesh within which material flows freely; the Lagrange domain moves with material deformation. A contact algorithm is used to dynamically load the detonation products onto the surrounding rock, and the feedback effect of surrounding rock failure on the flow path of the detonation products is considered.

[0113] Fifth, to eliminate boundary stress wave reflection, non-reflective boundary conditions were applied to the far-field boundary of the model, and accurate initial detonation parameters were set at the detonation center based on Chapman-Jouguet (CJ) theory.

[0114] Sixth, use the ABAQUS / Explicit solver to perform transient dynamic calculations and output the surrounding rock stress field, damage variable distribution, velocity field, and energy distribution.

[0115] Simulation results show that the maximum principal stress in the surrounding rock after blasting reached 145 MPa, and the cracks were distributed radially. Comparison shows that the numerical results agree well with the field test data, effectively verifying the accuracy of the invention.

[0116] In summary, the application process and numerical results of this simulation embodiment in tunnel blasting simulation verify its superiority in accurately describing the blasting damage evolution, detonation product flow, and fluid-structure interaction of brittle rock masses.

[0117] The specific embodiments of the present invention have been described in detail above with reference to the accompanying drawings, but the present invention is not limited to the described embodiments. For those skilled in the art, various changes, modifications, substitutions, and variations can be made to these embodiments without departing from the principles and spirit of the present invention, and these variations still fall within the protection scope of the present invention.

Claims

1. A method for simulating underground engineering blasting based on secondary development of ABAQUS, characterized in that, The method comprises the following steps: Step 1, a geometric model comprising surrounding rock, support structure and blasting action area is established, and meshing is performed, so that the calculation area is divided into Euler domain and Lagrangian domain; Step 2, a JH-2 constitutive model is called and implemented through a UMAT user subroutine interface of ABAQUS, so as to describe the strength change, volume compression and damage evolution process of rock under blasting load; Step 3, in the Euler domain, a JWL state equation is used to model the pressure-volume relationship and energy release law of explosive detonation products; Step 4, a coupled Euler-Lagrangian method is applied to realize the coupled interaction simulation between detonation products and surrounding rock; Step 5, boundary conditions and blasting loading parameters of the geometric model are set, wherein the boundary conditions comprise a non-reflective boundary condition applied to a far-field boundary of the geometric model to eliminate stress wave reflection; Step 6, transient dynamics calculation is performed, and results are output for analyzing blasting effect, wherein the output results comprise surrounding rock stress field, damage variable distribution, velocity field and energy distribution.

2. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 1, characterized in that: In step 2, the JH-2 constitutive model comprises a material strength model, a state equation EOS and a damage model, and expressions of the material strength model, the state equation EOS and the damage model are embedded into a user subroutine of ABAQUS to realize the JH-2 constitutive model.

3. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 2, characterized in that: The expression of the material strength model is as follows: ， wherein is the normalized full intensity, is the normalized broken intensity, is a damage variable, where 0≤ D ≤1; The expression of the state equation EOS is as follows: ， wherein The evolution law of a damage variable D of the damage model is expressed as follows: = 0.5 The mathematical expression of the JWL state equation in step 3 is as follows: / The evolution law of a damage variable D of the damage model is expressed as follows: 0-1, In step 4, the Euler domain and the Lagrangian domain realize bidirectional coupling through an interface contact algorithm; a solving process is based on mass, momentum and energy conservation equations, and a Lagrangian-Euler algorithm is used to ensure numerical stability. is the current density, The mass conservation equation is expressed as follows: 0 is the initial density, K 1, K 2, K 3 is the material bulk modulus parameter; The energy conservation equation is as follows: wherein , wherein D 1 and D 2 are material constants, Δ In step 5, the mathematical expression of the non-reflective boundary condition is as follows: p is the equivalent plastic strain increment for each calculation step, The CJ initial condition applied to the blast center is expressed by a formula as follows: f is the failure equivalent plastic strain, P ∗ The calculation of the transient dynamics is realized through small time step iteration based on an explicit integration algorithm of an ABAQUS / Explicit solver; and the output results introduce an energy balance relationship as follows: HEL is the normalized pressure, T ∗ ​ HEL is the normalized tensile pressure.

4. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 1, characterized in that: ​ ， wherein P is the pressure of the detonation products; V / V 0 is the ratio of the current volume to the initial volume; A , B , R 1, R 2, ​ K is an empirical constant for the explosive; E is the specific internal energy of the detonation products.

5. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 1, characterized in that: ​ 6. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 5, characterized in that, ​ ， wherein ​ is the fluid density, v is the velocity vector; the momentum conservation equation is: ， ​ ， wherein E is the total energy density, P is the pressure.

7. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 4, characterized in that: ​ ， wherein ​ n is the boundary normal stress, ​ is the medium density, c is the wave velocity, v n is the boundary normal velocity.

8. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 1, characterized in that: In step 5, the setting of the blast loading parameters includes defining the initial pressure field of the detonation products in the Euler domain based on the JWL state equation, imposing the Chapman-Jouguet state parameters at the blast center, including the initial density ​ 0, Detonation velocity D , CJ pressure P CJ and specific energy E 0.

9. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 8, characterized in that: ​ ， wherein ​ is the adiabatic index of the explosive product.

10. The underground engineering blasting simulation method based on ABAQUS secondary development according to claim 1, characterized in that: ​ ， wherein is the total energy released by the explosion, is the internal energy of the detonation products, is the kinetic energy of the surrounding rock and detonation gas, is the strain energy absorbed by the surrounding rock.

Images