A plastic explosive source substructure method

Abstract

This invention discloses a method for plastic blast source substructures, belonging to the field of numerical simulation technology for blast resistance in geotechnical engineering. The method includes: small-scale free wave field calculation, establishing a blast source-free field model and extracting displacement time histories; equivalent load calculation, establishing a plastic blast source substructure model with elastic transition elements, fixing internal nodes and applying displacement time histories to obtain boundary node reactions as equivalent blast loads; and large-scale site-structure blast response calculation, applying the equivalent loads to the large-scale model for nonlinear dynamic response analysis. This invention effectively decouples the influence of internal node plastic deformation history on equivalent load calculation by introducing elastic transition elements at the substructure boundaries, solving the problem that traditional substructure methods cannot handle near-field plastic wave propagation. It significantly improves computational efficiency while ensuring computational accuracy and is suitable for large-scale blast resistance numerical simulation in geotechnical engineering.

Description

Technical Field

[0001] This invention relates to the field of numerical simulation technology for blast resistance in geotechnical engineering, and in particular to a method for the substructure of a plastic blast source. Background Technology

[0002] With the rapid construction and large-scale development of large-scale underground infrastructure, geotechnical engineering is increasingly facing severe challenges from extreme dynamic loads such as blasting, rockbursts, and penetrating weapon impacts. These extreme loads are usually accompanied by high-frequency, high-strain-rate stress wave propagation processes, which can easily induce catastrophic engineering accidents such as tunnel collapses and slope instability, posing a serious threat to the stability and service safety of engineering structures. Due to the high cost, high risk, and difficulty in reproducing real explosion tests, numerical simulation has become an indispensable and important tool for studying explosion response problems.

[0003] Currently, mainstream methods for numerical simulation of explosions include the Lagrangian method, the Eulerian method, the arbitrary Lagrangian-Eulerian method, smoothed particle hydrodynamics, and the material point method. These methods have been widely integrated into general-purpose software such as LS-DYNA, AUTODYN, and ABAQUS. However, due to the high frequency and high strain rate characteristics of explosion loads, the mesh size of the numerical model has a significant impact on the accuracy of the calculation results. High-density meshes are usually required to accurately capture the high-frequency components in the explosion response. Furthermore, the stable time step of explicit time-domain integration algorithms is closely related to the smallest element size in the model. For large-scale engineering problems involving the propagation of shock waves in the mid-to-far field, the computational cost of using a comprehensive, refined model covering the explosive, medium, and target is often prohibitive.

[0004] An effective solution is to decouple the load source from the target structure. By constructing an equivalent load, load transfer and mesh transition from the load source to the target structure can be achieved, enabling multi-scale response analysis. Drawing on the concept of internal substructures in earthquake engineering, a multi-scale analysis method based on the blast source substructure is proposed. By transforming the near-field free wave field into an equivalent blast load in a large-scale model, efficient simulation of large-scale site blast response analysis under elastic assumptions is successfully achieved.

[0005] However, existing substructure analysis methods have inherent theoretical limitations when dealing with near-field plastic damage problems. The derivation of traditional substructure methods is based on the linear elastic assumption, requiring the system matrix to satisfy sparsity characteristics to calculate equivalent loads solely based on the boundary displacement field. When considering the plastic behavior of the near-field medium, due to the historical correlation of plastic strain, the stiffness matrix is ​​no longer sparse, making it impossible to decouple the influence of internal nodes in the substructure. Directly applying traditional methods will result in significant errors. This forces existing methods to extrapolate the calculation boundary to the elastic region. However, for high-yield explosions, the radius of the elastic region is typically 80-120 times the proportional radius, severely weakening the regional reduction advantage of substructure methods and limiting their application in strongly nonlinear explosion problems. Summary of the Invention

[0006] The purpose of this invention is to provide a method for plastic burst source substructure, which aims to solve the problem that traditional substructure methods cannot handle the propagation of near-field plastic waves.

[0007] The technical solution of the present invention is as follows:

[0008] A method for constructing a plastic explosive source substructure includes:

[0009] S1. Small-scale free wave field calculation: Establish a small-scale model of the explosion source-free field in the near-explosion source region, set artificial boundary conditions around the model, perform dynamic calculations, and extract displacement time history data of the site motion.

[0010] S2. Equivalent load calculation: Establish a plastic explosive source substructure model containing elastic transition units, fix the internal nodes, apply the motion time history data of the corresponding nodes obtained in step S1 to the boundary nodes, perform dynamic response calculation on the substructure, and obtain the reaction time history data of the boundary nodes as the equivalent explosive load.

[0011] The specific steps are as follows:

[0012] The calculation matrix for the nodal forces in the S2.1 inner zone is expressed as follows:

[0013] (13)

[0014] In the formula, the subscripts A, B, C, T, and I correspond to the outer node A, the middle node B, the inner node C, the elastic transition element node T, and the internal node I of the substructure. M, C, and K are the mass matrix, damping matrix, and stiffness matrix of the finite element model, respectively. , and These are vectors representing nodal acceleration, velocity, and displacement, respectively. F is the nodal force vector.

[0015] According to equation (13), the nodal forces of the inner nodes are:

[0016] (14)

[0017] S2.2 Assume that node T and all nodes inside it remain stationary, that is, the displacements of nodes T and I are 0:

[0018] (15)

[0019] Substituting equation (15) into equation (14), we get:

[0020] (16)

[0021] In the formula , and This refers to the equivalent value that needs to be calculated.

[0022] The equation of motion for the S2.3 plastic explosion source substructure is expressed as:

[0023] (17)

[0024] In the formula, the superscript PS denotes the plastic substructure model. According to finite element theory, if the mesh size and material properties of the substructure model are completely identical to those at the corresponding locations in the free field model, then the following relationship holds:

[0025] (18)

[0026] Subscript It refers to AB, AC, BA, BB, BC, CA, CB, CC, CT, TC, and TT.

[0027] Formula (19) is obtained:

[0028] (19)

[0029] S3. Calculation of large-scale site-structure explosion response: Establish a large-scale model of site-structure interaction, apply the equivalent load obtained in step S2 to the corresponding nodes, and complete the nonlinear dynamic response analysis.

[0030] Further, the nodes mentioned in step S2 include: substructure outer node A, substructure middle node B, substructure inner node C, elastic transition unit node T, and internal node I, wherein the elastic transition unit is disposed between substructure boundary node C and internal node I.

[0031] Furthermore, the elastic transition unit in step S2 adopts an elastic constitutive structure, and its elastic stiffness is consistent with the elastic stage stiffness of the adjacent plastic material.

[0032] Furthermore, fixing the internal nodes in step S2 specifically means assuming that node T and all nodes inside it remain stationary, i.e., node T is displaced. This decouples the influence of the internal node plastic deformation history on the calculation of equivalent load.

[0033] Furthermore, the equivalent explosive load described in step S2 includes nodal forces at nodes B, C, and T. , and Among them, nodal forces Used to keep internal nodes stationary, while nodal forces , This ensures that the displacement and plastic deformation of the external wave field are consistent with those of the free wave field model.

[0034] Furthermore, the large-scale model described in step S3 maintains consistency between the mesh and substructure in the loading region and includes a layer of elastic transition units.

[0035] Furthermore, the method also includes a multi-scale grid transition step: when extracting motion information at the substructure boundary in the small-scale free field model, according to the grid size ratio of the far-field large-scale model, interval extraction is performed in the dense node group, and the extracted displacement data is applied to the corresponding nodes of the large-scale substructure model.

[0036] Furthermore, the grid size ratio is 1:2 or 1:4.

[0037] This invention also provides a method for constructing a standardized equivalent explosive load library. The method is characterized by pre-calculating and constructing a standardized equivalent explosive load library covering different equivalents and different burial depths based on the above-mentioned plastic explosive source substructure method. When facing multi-condition blast resistance analysis, the corresponding equivalent load time history is called from the database and applied to the structural model as a force boundary condition.

[0038] Beneficial effects:

[0039] 1. Improved computational accuracy: This invention effectively decouples the influence of the internal node plastic deformation history on the equivalent load calculation by introducing elastic transition elements at the substructure boundary, enabling the substructure method to accurately handle near-field plastic wave propagation problems. Two-dimensional and three-dimensional numerical examples demonstrate that the method of this invention can stably and accurately reproduce crack propagation, plastic wave propagation, and rock-structure interaction, with the peak error at each measuring point consistently remaining within 1%, significantly superior to traditional substructure methods.

[0040] 2. Significantly improved computational efficiency: This invention reduces the substructure truncation boundary from approximately 120 times the traditional scale radius to 40 times the scale radius. Under the same accuracy requirements, this reduces the number of mesh elements in the 3D model by approximately 70% to 90%, improving computational efficiency by 20% to 30% compared to traditional methods. Furthermore, by incorporating a multi-scale mesh transition strategy, computational consumption and time costs are further reduced.

[0041] 3. Enhanced Engineering Applicability: The method of this invention can achieve physical and numerical decoupling between the "near-field of the blast source" and the "far-field of the structure." Combined with an experimentally calibrated blast source model, a standardized equivalent blast load library can be pre-established. When facing the blast resistance analysis of large-scale infrastructure with multiple blast sources and multiple working conditions, the complex near-field fluid-structure interaction problem is transformed into a simple pure Lagrangian force boundary loading problem, avoiding repetitive detailed modeling and complex mesh transition processing, significantly improving the operability and robustness of the analysis process. Attached Figure Description

[0042] Figure 1 This is a schematic diagram of the equivalent load input for the explosion source substructure method of the present invention;

[0043] Figure 2 This is a schematic diagram of a traditional explosion source substructure;

[0044] Figure 3 This is a schematic diagram of the node division of the free wave field model in the plastic substructure method of this invention;

[0045] Figure 4 This is a schematic diagram of the plastic explosive source substructure of the present invention;

[0046] Figure 5 This is a comparative schematic diagram of the plastic substructure of the present invention and the traditional substructure;

[0047] Figure 6 This is a schematic diagram of the calculation process for the plastic explosive source substructure method of the present invention;

[0048] Figure 7 This is a schematic diagram of the jump node input of the present invention;

[0049] Figure 8 This is a schematic diagram of a two-dimensional analysis model in an embodiment of the present invention;

[0050] Figure 9 This is a schematic diagram of the three-dimensional analysis model in an embodiment of the present invention;

[0051] Figure 10 These are comparison diagrams of crack propagation using different methods in embodiments of the present invention;

[0052] Figure 11 This is a comparison diagram of numerical simulation and experimental damage in an embodiment of the present invention;

[0053] Figure 12 This is a comparison diagram of overpressure propagation results of different methods in the embodiments of the present invention;

[0054] Figure 13 This is a comparison chart of overpressure results using different methods in the embodiments of the present invention;

[0055] Figure 14 This is a comparison diagram of acceleration response using different methods in embodiments of the present invention;

[0056] Figure 15 This is a comparison chart of the reduction ratio of the number of units in different methods in the embodiments of the present invention;

[0057] Figure 16 This is a schematic diagram illustrating the application scenario of the explosion source substructure model of this invention.

[0058] Figure 17 This is a schematic diagram of the method flow of the present invention. Detailed Implementation

[0059] To make the technical solution of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0060] Example 1: Verification of Two-Dimensional Crack Propagation Morphology

[0061] This embodiment uses a two-dimensional analysis model to compare and verify the differences in simulation accuracy of the plastic region between the traditional explosion source substructure method and the plastic explosion source substructure method proposed in this invention.

[0062] Step 1: Establish a two-dimensional refined finite element model

[0063] Experimental study on crack propagation induced by rock blasting, establishing such as Figure 8 The model shown is a two-dimensional finite element model. It consists of granite, air, and TNT explosive, with a geometric dimension of 5000 mm × 5000 mm. The air domain is 1500 mm × 1500 mm, and the explosive radius is 12 mm. To control crack propagation direction, a 100 mm hole is left in the center of the model. The mesh size is 10 mm, and a single-layer solid element model is used. Normal symmetric boundary conditions are applied to both the front and back planes of the model to simulate two-dimensional plane strain. Non-reflective boundary conditions are set on the surrounding rock sides to eliminate the reflection effect of the blast shock wave at the computational boundaries.

[0064] Numerical calculations employed the arbitrary Lagrange-Euler method, with the air domain and explosive domain using the Euler algorithm, and the surrounding rock material using the Lagrange algorithm. Fluid-structure interaction was achieved using the keyword CONSTRAINED_LAGRANGE_IN_SOLID. The surrounding rock was modeled using the MAT_RHT constitutive model, the explosive domain using the MAT_HIGH_EXPLOSIVE_BURN model, and the air domain using the MAT_NULL model. Material parameters were calibrated based on explosion test results, and specific parameters are shown in Tables 1 and 2.

[0065] Table 1. Parameters of the constitutive model for granite

[0066]

[0067] Table 2 Constitutive Model Parameters for Air and Explosives

[0068]

[0069] Step 2: Substructure Setup and Calculation

[0070] Rectangular substructures with side lengths of 1400 mm, 1600 mm, and 1800 mm were respectively set around the explosion source. Numerical simulations were conducted using the overall model method, the traditional explosion source substructure method, and the plastic explosion source substructure method proposed in this invention, respectively.

[0071] For the traditional burst source substructure method, according to Figure 2 The diagram shows a substructure consisting of three layers of nodes A, B, and C, and two layers of internal units. Node C is fixed, and the corresponding free-field motion is applied to nodes A and B. The equivalent force is calculated using equation (5). and This is used as an equivalent explosive load input to the site-structure interaction model.

[0072] Regarding the plastic explosion source substructure method proposed in this invention, according to Figure 3 and Figure 4 As shown, an elastic transition element is added between the traditional substructure boundary node C and the internal node I. The interface node between the transition element and the internal element is denoted as T. A plastic blast source substructure containing nodes A to T is established, with the elastic stiffness of the elastic transition element consistent with the elastic stage stiffness of the adjacent plastic material. Node T is fixed, and the displacement time history data of the corresponding nodes obtained in step 1 are applied to nodes A, B, and C. Dynamic response calculations are performed on the substructure to obtain the reaction time history data of nodes B, C, and T as equivalent blast loads. This equivalent load is applied to a large-scale site-structure interaction model to complete the dynamic response analysis.

[0073] Crack propagation morphologies obtained by different methods are as follows: Figure 10As shown in the figure. The results indicate that the traditional substructure method has poor applicability in the plastic zone, and its crack propagation characteristics deviate significantly from the overall model results. Moreover, the error increases sharply as the substructure approaches the explosion source. For the extremely small substructure with a side length of 1400 mm, the crack area difference of the traditional substructure method can reach a maximum of 34.42%, and obvious spurious derived cracks also appear.

[0074] In contrast, the plastic burst source substructure method proposed in this invention effectively shields the interference of the plastic strain history of inner nodes on the nonlinear response of external elements by leveraging the physical isolation effect of elastic transition units. Its crack propagation morphology is highly consistent with the overall model results, and the error is insensitive to the size of the substructure, with the maximum error consistently remaining within a small range (less than 4%). Given that crack evolution is extremely sensitive to initial conditions, accurately reproducing the crack propagation process is a stringent verification criterion; this result strongly demonstrates the high fidelity of this method in characterizing complex plastic failure processes.

[0075] Example 2: Three-dimensional plastic wave field propagation and structural response verification

[0076] This embodiment uses a three-dimensional analysis model to evaluate the computational accuracy of the proposed plastic burst source substructure method in three-dimensional plastic wave field propagation and surrounding rock-structure interaction problems.

[0077] Step 1: Establish a three-dimensional refined finite element model

[0078] Established based on relevant experimental research, such as Figure 9 The three-dimensional finite element model shown is composed of air, explosives, surrounding rock, filling sand, soil, and a reinforced concrete underground structure. C20 concrete with similar density and P-wave velocity is used to simulate Class IV surrounding rock, while C40 concrete is used for the reinforced concrete tunnel structure. The explosive equivalent is 13 kg, with a side length of 200 mm. The geometric dimensions of the air zone and surrounding rock are 4500 mm × 4000 mm × 3000 mm. The tunnel has a longitudinal length of 3000 mm, a wall thickness of 20 mm, and a span of 400 mm. Soil surrounds the surrounding rock, extending the overall model dimensions to 6000 mm × 5000 mm × 5000 mm. The charging holes are filled with sand to a depth of 1430 mm. The overall mesh size of the model is controlled at 10 mm. Consistent with experimental results, pressure measuring points P-1 to P-3 are placed at 900 mm, 600 mm, and 300 mm above the tunnel, and acceleration measuring points A-1 and A-2 are placed at the crown and 45° angle to the arch.

[0079] The numerical analysis employed the arbitrary Lagrange-Eulerian method, using the Lagrange algorithm for the surrounding rock and soil layers, and simulating the interaction between the soil and surrounding rock using a surface-to-surface contact algorithm. The tunnel structure consisted of a steel mesh and concrete; the concrete was modeled using Lagrange solid elements, and the steel reinforcement using Beam elements. The tunnel structure was connected to the surrounding rock via shared nodes. Non-reflective boundary conditions were applied to all four sides and the bottom of the surrounding rock.

[0080] Material constitutive models: Concrete uses MAT_RHT, steel reinforcement uses MAT_PLASTIC_KINEMATIC, sand uses MAT_DRUCKER_PRAGER, and air and explosive materials are consistent with the two-dimensional model. Specific parameters are shown in Tables 3 to 5.

[0081] Table 3 Concrete Constitutive Model Parameters

[0082]

[0083] Table 4 Parameters of the Reinforcement Constitutive Model

[0084]

[0085] Table 5 Constitutive model parameters for soil and filling sand

[0086]

[0087] Step 2: Model Validation and Substructure Calculation

[0088] First, to ensure the reliability of the baseline solution, the simulation results of the three-dimensional overall model are verified against the field test data. For example... Figure 11 As shown, the macroscopic damage phenomena of the surface and underground structures are in good agreement with the experimental observations, and the overall model has good accuracy.

[0089] Based on the overall model, a cubic substructure with a side length of 1000 mm was set around the explosion source for wave field propagation analysis. Calculations were performed using the overall model, the traditional substructure method, and the plastic substructure method proposed in this invention, following the same method as in Example 1.

[0090] Step 3: Result Comparison and Analysis

[0091] Pressure propagation wave field snapshots from different methods, such as Figure 12 As shown in the figure. The results indicate that the traditional substructure method has significant errors when applied to plastic wave field analysis, and the wave field exhibits obvious disordered characteristics. In contrast, the plastic substructure method proposed in this invention can more accurately capture the propagation characteristics of the plastic wave field, and its pressure contour map shows good consistency with the overall model results.

[0092] To quantify the error analysis, a comparative analysis was conducted on the pressure time history (P1 to P3) at measuring points in the rock below the blast source and the acceleration response at key structural locations (A1 at the arch crown and A2 at the 45° angle to the circular arch). Figure 13 and Figure 14 As shown in Table 6, the error results are as follows.

[0093] Table 6. Error Calculation Results for Different Methods

[0094]

[0095] The overall model numerical simulation results agree well with the experimental data, with errors all controlled within 10%. Based on the overall model results, the traditional substructure method shows significant errors in pressure and acceleration response, with a maximum error reaching 200%, indicating severe distortion. In contrast, the plastic explosive source substructure method proposed in this invention exhibits time history curves that almost coincide with the benchmark solution, with peak errors at each measurement point consistently remaining within 1%.

[0096] Comprehensive analysis shows that the plastic burst source substructure method proposed in this invention successfully breaks through the bottleneck of traditional methods being unable to handle near-field plasticity, and can accurately characterize the plastic crack evolution process, the plastic wave field propagation law, and the dynamic interaction effect between the surrounding rock and the structure.

[0097] Example 3: Computational Efficiency Analysis

[0098] This embodiment quantitatively evaluates the computational efficiency of the method of the present invention. The percentage reduction in the overall number of cells compared to a fine mesh is used as the metric. The free field radius is set to 100 times the charge radius, and the substructure mesh size is defined as x. The outer mesh size is taken as 2 times and 4 times that of the inner mesh, respectively, denoted as 2x and 4x. Figure 15 The number of elements used by different calculation methods under different mesh sizes and substructure region sizes was compared.

[0099] Depend on Figure 15 As can be seen, compared to the overall model, introducing the substructure method into the two-dimensional analysis model reduces the total number of elements by 50% when the outer mesh size is twice that of the inner mesh, and by more than 60% when the mesh size is increased by four times. In the three-dimensional analysis model, the substructure method demonstrates even greater advantages in computational efficiency, reducing the total number of elements by approximately 70% to 90%. Furthermore, the settling time step in explicit dynamic analysis is controlled by the minimum mesh size; increasing the size of the outer mesh not only reduces the number of elements but also allows for larger time steps locally, thus improving computational efficiency in both spatial and temporal dimensions.

[0100] Further comparison of the traditional substructure method and the plastic substructure method of this invention reveals that the plastic substructure method only adds one layer of elastic transition unit to the traditional substructure, resulting in a very limited additional computational cost and minimal impact on the overall computational scale. Conversely, this method achieves further efficiency improvements through the physical isolation of the inner and outer regions: the traditional substructure method, limited by theoretical assumptions, must extrapolate the truncation boundary to the far-field elastic region (approximately 120 times the proportional radius) to avoid plasticity errors; while the method of this invention successfully shrinks the applicable boundary to the interface between the near-field crushing region and the plastic region (approximately 40 times the proportional radius). The reduction in the computational domain radius means a sharp decrease in the geometric progression of the high-density mesh region volume. Data shows that, under the same accuracy requirements, the computational efficiency of the plastic substructure method of this invention is 20% to 30% higher than that of the traditional method.

[0101] Example 4: Multi-scale grid transition application

[0102] This embodiment illustrates the application of the method of the present invention in multi-scale mesh transition. To achieve an efficient transition from a fine near-field mesh to a sparse far-field mesh, the present invention employs a wavefield transfer strategy of "interval node mapping" (i.e., "skip" input), such as... Figure 7 As shown.

[0103] In practice, a high-precision wavefield is first calculated in a small-scale near-field free-field model. When extracting motion information at the substructure boundary, it is not necessary to traverse all fine nodes. Instead, extraction is performed at intervals within the fine node group based on the mesh size ratio of the large-scale far-field model (e.g., 1:2 or 1:4). Subsequently, the extracted displacement data is applied to the corresponding nodes of the large-scale substructure model and converted into equivalent nodal forces through dynamic calculations.

[0104] This strategy has two advantages: 1) Avoiding mesh size effects: By transferring displacement and back-calculating nodal forces from substructures, the mesh mismatch error caused by directly mapping nodal forces is cleverly avoided, ensuring the accuracy of multi-scale models; 2) Reducing data processing scale: By "jumping" acquisition, the amount of data interaction at the interface is significantly reduced, avoiding the problems of mesh quality degradation and modeling complexity caused by the use of transition units in traditional methods.

[0105] Example 5: Construction and Application of Standardized Equivalent Explosive Load Library

[0106] This embodiment illustrates the application of constructing a standardized equivalent explosive load library based on the method of this invention. For example... Figure 16 As shown, based on a high-precision blast source model that has been experimentally calibrated, a "standardized equivalent blast load library" covering different yields and burial depths can be pre-calculated and constructed.

[0107] When dealing with multi-condition blast-resistant design for large and complex engineering projects, there is no need to repeatedly perform time-consuming near-field blast simulations. Instead, the corresponding equivalent load time histories can be retrieved from the database and applied as force boundary conditions to the structural model. In this case, the overall calculation is transformed into a structural dynamics problem within a pure Lagrangian framework, completely eliminating the costly fluid-structure interaction (FSI) calculations. This "one-time calculation, multiple reuse" strategy significantly improves the efficiency of multi-blast source combinations and parameter sensitivity analysis.

[0108] This method can be applied to the evaluation of the protective effectiveness of earth-penetrating weapons, the safety analysis of major engineering facilities such as bridges, dams and nuclear power plants under extreme explosive loads, and the study of the disturbance effect of tunnel blasting on the surrounding rock mass and structure. It has good engineering applicability and broad application prospects.

[0109] Traditional explosive source substructure methods and their limitations

[0110] To calculate the equivalent load, the traditional blast source substructure method first divides the blast source-free field model into nodes. For example... Figure 1 As shown, the model nodes are divided into internal nodes I, three layers of substructure nodes (denoted as A, B, and C from the outside in), and external nodes E. Based on the above node division, the motion equations of the explosion source-free field model can be expressed as:

[0111] (1)

[0112] In the formula, the subscripts A, B, C, E, and I correspond to the five types of nodes mentioned above. M, C, and K are the mass matrix, damping matrix, and stiffness matrix of the finite element model, respectively. , and These are the vectors of nodal acceleration, velocity, and displacement, respectively; F is the nodal force vector.

[0113] Node A divides the computational domain into an inner and outer region. To ensure that the wavefield in the outer region is consistent with the free wavefield of the overall model, the nodal forces in the inner region should satisfy the following equation:

[0114] (2)

[0115] The superscript "0" indicates the free wave field response. According to equation (2), the nodal forces of the inner nodes are:

[0116] (3)

[0117] From equation (3), it can be seen that the nodal force of node I at this time is This means that to achieve accurate wavefield input, equivalent nodal forces must be applied to all inner nodes, without substantially reducing computational complexity. Therefore, it is necessary to further simplify the calculation process and decouple the influence of the inner nodes. According to wave propagation theory, the outer traveling wave field can be uniquely determined by the displacement field of the boundary node B. Based on this, it can be assumed that all nodes inside node B remain stationary, i.e., the displacements of nodes C and I are 0.

[0118] (4)

[0119] Substituting equation (4) into equation (3), we get:

[0120] (5)

[0121] at this time This eliminates the influence of nodal forces on internal nodes. The nodal forces of nodes B and C in the formula are... and This is the equivalent force that needs to be calculated. As shown in equation (5), thanks to the sparsity of the system matrix, the calculation of the equivalent force at nodes A and B is only related to the free-field motion of the corresponding nodes and the sub-matrix formed by adjacent nodes A, B, and C. Therefore, the original model can be extracted to contain only the three-layer node structure of nodes A, B, and C and its two internal layers of units (e.g., ...). Figure 2 (As shown). Fix node C, and then apply the corresponding free field motion to nodes A and B of the substructure model. Calculate the equivalent force using equation (5). and This can be used as an equivalent explosive load input into the site-structure interaction (SSI) model.

[0122] The physical basis for the efficient computation of traditional substructure methods lies in the banded sparsity of the finite element mass matrix, damping matrix, and stiffness matrix. For the mass matrix, the sparsity stems from the lumped mass assumption; for the damping matrix, the Rayleigh damping assumption is commonly used in engineering, and its expression is as follows:

[0123] (6)

[0124] In the formula, α and β are Rayleigh damping coefficients. and These are the initial mass matrix and the initial stiffness matrix, respectively. Since the initial elastic stiffness matrix is ​​based on the locality assumption that "nodal forces are only related to the strain of adjacent elements," it is inherently sparsity. The damping matrix, as a linear combination of the mass and stiffness matrices, also maintains a banded sparsity characteristic.

[0125] However, once the material enters the plastic stage, due to the historical correlation of plastic deformation, the overall stiffness matrix no longer satisfies the sparsity assumption. For plastic materials, incremental calculations are typically used, and the load increment caused by the stiffness matrix can be expressed as:

[0126] (7)

[0127] In the formula, , and These represent the nodal displacement increment, nodal load increment, and nonlinear strain increment at the Gaussian integration point, respectively. The initial elastic stiffness matrix of the overall structure remains unchanged during the calculation. and This is the stiffness matrix associated with incremental nonlinear deformation. Related to the location of the incremental nonlinear deformation region, It is related not only to the location of the region where nonlinear deformation occurs, but also to the tangential modulus or tangential constitutive matrix of the material within these regions. The determination of this matrix during the calculation depends on the nonlinear constitutive model of the material. From equation (7), we can obtain:

[0128] (8)

[0129] After simplification, the calculation of nodal forces satisfies the following formula:

[0130] (9)

[0131] In the formula This is the integrated plastic stiffness matrix. Existing research has demonstrated... The matrix is ​​a full matrix. The matrix is ​​a diagonal matrix, which leads to The matrix is ​​also full. At this point, the matrix for calculating the forces at the inner nodes is:

[0132] (10)

[0133] According to equation (10), we can obtain:

[0134] (11)

[0135] At this point, the nodal forces in the inner region are globally coupled, even if relevant constraints are applied to cause displacement at node C. , It is still not equal to 0, and the influence of internal nodes cannot be effectively decoupled:

[0136] (12)

[0137] This indicates that the calculation of nodal forces in the plastic state depends not only on the nodal displacement at the current moment, but also on the plastic deformation history of the entire element. The degrees of freedom of internal nodes cannot be simply "constrained" and thus decoupled as in the elastic state. This is the fundamental reason why the traditional substructure method fails when dealing with the near-field plastic wave propagation problem.

[0138] This invention relates to a method for constructing a plastic explosive source substructure.

[0139] To address the aforementioned problems, this invention proposes an improved strategy: adding an elastic transition element between the traditional substructure boundary node C and the internal node I to decouple the influence of the internal node. Let T be the interface node between this transition element and the internal element. The node division of the improved plastic substructure free-field model is as follows: Figure 3 As shown.

[0140] In this model, since the transition elements adopt an elastic constitutive model, their corresponding stiffness matrix recovers its sparse characteristics. Therefore, the calculation matrix for the nodal forces in the inner region can be expressed as:

[0141] (13)

[0142] According to equation (13), the nodal forces of the inner nodes are:

[0143] (14)

[0144] To achieve decoupling, a constraint assumption similar to that in traditional methods is introduced to decouple the influence of internal nodes. It is assumed that node T and all nodes inside it remain stationary, i.e., the displacements of nodes T and I are zero.

[0145] (15)

[0146] Substituting equation (15) into equation (14), we get:

[0147] (16)

[0148] At this time, the nodal force of node I The nodal forces at the internal nodes are eliminated. The nodal forces at nodes B, C, and T are shown in the equation. , and This is the equivalent force that needs to be calculated. From equation (16), it can be seen that the calculation of the equivalent force at nodes A, B, C, and T is only related to the free-field motion of nodes A, B, C, and T, and the corresponding element sub-matrices. Therefore, the plastic burst source substructure containing nodes A to T (such as...) can be extracted from the original model. Figure 4 As shown), calculations are performed independently.

[0149] The equation of motion for the substructure of the plastic explosion source can be expressed as:

[0150] (17)

[0151] In the formula, the superscript PS denotes the plastic substructure model. According to finite element theory, if the mesh size and material properties of the substructure model are completely identical to those at the corresponding locations in the free field model, then the following relationship holds:

[0152] (18)

[0153] Subscript Refers to AB, AC, BA, BB, BC, CA, CB, CC, CT, TC, and TT. Substitute equation (21) into equation (20) and fix node T, i.e., the displacement of node T. =0, and then apply the corresponding free field motion to nodes A, B, and C, that is = , = , = We obtain the following formula:

[0154] (19)

[0155] Comparing equation (19) and equation (16), it can be seen that the nodal forces at nodes B, C, and T calculated using the plastic explosive source substructure model are consistent with the results calculated using the free wave field model. Therefore, only the equivalent forces calculated using the substructure composed of four layers of substructure units are needed. , and This allows for the accurate input of the equivalent explosive load. The nodal force at node T is... Used to keep the internal nodes stationary, while the nodal forces of nodes B and C... , This ensures that the displacement and plastic deformation of the external wave field are consistent with those of the free wave field model. The elastic stiffness of the elastic transition element should be consistent with the elastic stage stiffness of the adjacent plastic material to ensure that the initial displacement and damping matrix are consistent with the overall model during plasticity calculations, thus avoiding errors in determining the plasticity state.

[0156] Compared to traditional explosive source substructure methods, the core improvement of the plastic explosive source substructure method of this invention lies in the introduction of an elastic transition unit (such as...). Figure 5(As shown). The elastic transition element physically isolates the internal plastic region from the substructure region, and "cuts off" the transmission path of the internal plastic deformation history to the substructure region by utilizing the sparsity of the elastic stiffness matrix. For LS-DYNA software, adding elastic transition elements can be achieved by creating a new part, which is convenient. When extracting nodal forces, only one more layer of nodal forces needs to be extracted, and the increase in computational complexity is negligible.

[0157] The embodiments described above are merely illustrative of several implementations of the present invention, and while the descriptions are specific and detailed, they should not be construed as limiting the scope of the present invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and these modifications and improvements all fall within the scope of protection of the present invention. Therefore, the scope of protection of this patent should be determined by the appended claims.

Claims

1. A method for constructing a plastic explosive source substructure, characterized in that, include: S1. Small-scale free wave field calculation: Establish a small-scale model of the explosion source-free field in the near-explosion source region, set artificial boundary conditions around the model, perform dynamic calculations, and extract displacement time history data of the site motion. S2. Equivalent load calculation: Establish a plastic explosive source substructure model containing elastic transition units, fix the internal nodes, apply the motion time history data of the corresponding nodes at the substructure location obtained in step S1 to the boundary nodes, perform dynamic response calculation on the substructure, and obtain the reaction time history data of the boundary nodes as the equivalent explosive load. The calculation matrix for the nodal forces in the S2.1 inner zone is expressed as follows: （13） In the formula, the subscripts A, B, C, T and I correspond to the outer node A of the substructure, the middle node B of the substructure, the inner node C of the substructure, the elastic transition element node T and the internal node I; M, C and K are the mass matrix, damping matrix and stiffness matrix of the finite element model, respectively. , and These are vectors representing nodal acceleration, velocity, and displacement, respectively. F is the nodal force vector; According to equation (13), the nodal forces of the inner nodes are: （14） S2.2 Assume that node T and all nodes inside it remain stationary, that is, the displacements of nodes T and I are 0: （15） Substituting equation (15) into equation (14), we get: （16） In the formula , and That is, the equivalent value that needs to be calculated; The equation of motion for the S2.3 plastic explosion source substructure is expressed as: （17） In the formula, the superscript PS denotes the plastic substructure model; according to the finite element theory, if the mesh size and material properties of the substructure model are completely consistent with those of the corresponding location in the free field model, then the following relationship holds: （18） Subscript Refers to AB, AC, BA, BB, BC, CA, CB, CC, CT, TC, and TT; Formula (19) is obtained: （19） S3. Calculation of large-scale site-structure explosion response: Establish a large-scale model of site-structure interaction, apply the equivalent load obtained in step S2 to the corresponding nodes, and complete the nonlinear dynamic response analysis.

2. The method for constructing a plastic explosive source substructure according to claim 1, characterized in that, The nodes mentioned in step S2 include: substructure outer node A, substructure middle node B, substructure inner node C, elastic transition unit node T, and internal node I; wherein the elastic transition unit is disposed between substructure boundary node C and internal node I.

3. The method for constructing a plastic explosive source substructure according to claim 1, characterized in that, The elastic transition unit in step S2 adopts an elastic constitutive model, and its elastic stiffness is consistent with the elastic stage stiffness of the adjacent plastic material.

4. The method for constructing a plastic explosive source substructure according to claim 1, characterized in that, The step S2, which involves fixing the internal nodes, specifically assumes that node T and all nodes inside it remain stationary, i.e., node T is displaced. .

5. The method for constructing a plastic explosive source substructure according to claim 1, characterized in that, The equivalent explosive load mentioned in step S2 includes the nodal forces at nodes B, C, and T. , and Among them, nodal forces Used to keep internal nodes stationary, while nodal forces , This ensures that the displacement and plastic deformation of the external wave field are consistent with those of the free wave field model.

6. The method for constructing a plastic explosive source substructure according to claim 1, characterized in that, The large-scale model described in step S3 maintains consistency between the mesh and substructure in the loading region and includes a layer of elastic transition units.

7. The method for constructing a plastic explosive source substructure according to claim 1, characterized in that, It also includes a multi-scale grid transition step: when extracting motion information at the substructure boundary in the small-scale free field model, according to the grid size ratio of the far-field large-scale model, interval extraction is performed in the dense node group, and the extracted displacement data is applied to the corresponding nodes of the large-scale substructure model.

8. The method for constructing a plastic explosive source substructure according to claim 7, characterized in that, The grid size ratio is 1:2 or 1:

4.

9. A method for constructing a standardized equivalent explosive load library, characterized in that, Based on the plastic explosive source substructure method according to any one of claims 1 to 8, a standardized equivalent explosive load library covering different equivalents and different burial depths is pre-calculated and constructed. When facing multi-condition blast resistance analysis, the corresponding equivalent load time history is called from the database and applied to the structural model as force boundary conditions.

Images