A method for predicting the detonation process of condensed charge based on ZND theory

Abstract

The present application discloses a method for predicting the detonation process of condensed charge based on ZND theory, which relates to the field of detonation technology. The method determines the combustion reaction factor of the current time step according to the detonation velocity and the unit equivalent size of the charge unit, and then obtains the gas phase density of the current time step in combination with the mixed density of the previous time step, and further determines the detonation gas pressure of the current time step by using the JWL state equation in combination with the gas phase unit mass internal energy of the previous time step, and obtains the mixed pressure of the current time step according to the detonation gas pressure and the combustion reaction factor, and substitutes it into the fluid dynamics control equation group to obtain the prediction result and enter the next time step for iterative solution. The method can calculate the detonation process of condensed charge concisely and efficiently, thereby obtaining the prediction results of each charge unit at each time step, and can effectively capture the generation and propagation process of the detonation wave in the detonation process, and improve the prediction accuracy and stability of the detonation process.

Description

Technical Field

[0001] The present application relates to the field of detonation technology, and in particular to a method for predicting the detonation process of condensed charge based on the ZND theory. Background Art

[0002] Whether the condensed charge explodes in water or in the air, its chemical reaction process includes complex phenomena such as flame, deflagration, and detonation. The detonation process mainly includes two basic stages: initiation and propagation. The explosive in the condensed charge loses stability and explodes after being affected by external factors. This process from a stable state to a lost stable state is called the initiation process of the explosive. The propagation of the detonation wave is caused by the impact compression of the shock wave on the explosive. When the intensity of the shock wave exceeds a certain value, it will cause a chemical reaction of the explosive material behind the shock wave front, and the energy released by the chemical reaction is used to support the continued propagation of the shock wave, thus forming the detonation process of the explosive. In the context of studying near-field explosions, it is necessary to understand and master the characteristics of the detonation process in order to provide an accurate and reliable load source for the study of engineering problems.

[0003] The CJ theory proposed by Chapman and Jouguet provides a one-dimensional theoretical description of the detonation phenomenon, which states that after the explosive passes through the shock wave front, the chemical reaction is completed very quickly after a very short time, and the explosive becomes a detonation product after a complete reaction. This theory regards the detonation wave as a strong discontinuity, thereby isolating the explosive and the detonation product for study, laying the foundation for detonation research and having epoch-making significance.

[0004] However, the CJ theory has the disadvantage of ignoring the details of chemical reactions. Therefore, the ZND theory proposed by Zeldovich, Von Neumann, and Doering improved the assumptions of the CJ theory, introduced finite-rate chemical reactions on the wavefront, and correctly revealed the basic structure of the detonation wave for the first time, clarified the triggering mechanism of the detonation chemical reaction, and considered the kinetic process of the chemical reaction, thus obtaining a more complete detonation model. Compared with the CJ theory, the ZND theory is more in line with the actual situation, so it has become another milestone in the development of detonation theory and is also a common theoretical basis for studying the detonation process.

[0005] In the field of detonation process research based on ZND theory, chemical reaction models such as ignition and growth reaction rate model and trinomial ignition model are often used as detonation models for research. Although the calculation results of these models can reflect the structural characteristics of detonation waves to a certain extent, they require many parameters, and the process of obtaining and confirming these parameters lacks sufficiently rigorous verification. In addition, the introduction of these chemical reaction models increases the complexity of the calculation of the entire detonation process and the demand for computing resources, resulting in the high complexity and unsatisfactory accuracy of the detonation process prediction currently achieved using these models, which affects the promotion of engineering applications. Summary of the invention

[0006] In view of the above problems and technical requirements, this application proposes a method for predicting the detonation process of condensed charge based on ZND theory. The technical solution of this application is as follows:

[0007] A method for predicting the detonation process of a condensed charge based on the ZND theory includes: for each charge unit in a condensed charge package at any nth time step:

[0008] According to the detonation velocity D CJ and the unit equivalent size of the drug package unit δ eff Determine the combustion reaction factor λ of the charge unit at the current time step n , integer parameter n≥1;

[0009] According to the mixture density ρ of the previous time step n-1 Combined combustion reaction factor λ n Determine the gas phase density at the current time step

[0010] According to the gas phase density of the current time step Combined with the gas phase unit mass internal energy e of the previous time step n-1 Determine the detonation gas pressure at the current time step using the JWL state equation

[0011] According to the detonation gas pressure Combined combustion reaction factor λ n Get the mixed pressure p of the current time step n , and substitute it into the fluid dynamics control equations constructed based on the unsteady compressible Euler equations to solve the mixed density ρ at the current time step n and the gas phase unit mass internal energy e at the current time step n And enter the solution of the next time step; where ρ 0 and e 0 All are initial values.

[0012] The beneficial technical effects of this application are:

[0013] The present application discloses a method for predicting the detonation process of condensed charge based on ZND theory. The method uses iterative updating of combustion reaction factors to characterize the evolution of the detonation process. Combined with the JWL state equation and a group of fluid dynamics control equations constructed based on the unsteady compressible Euler equations, the detonation process of condensed charge can be calculated concisely and efficiently, thereby obtaining the prediction results of each charge unit at each time step, and can effectively capture the generation and propagation process of detonation waves during the detonation process, thereby improving the prediction accuracy and stability of the detonation process of underwater explosive charges.

[0014] This application also gives the specific values ​​of the state equation parameters of the JWL state equation under two common types of explosives. There are only 5 state equation parameters and the calculation process has a physical basis. Experiments have confirmed that the accuracy of these state equation parameters is also very high, which can help to obtain calculation results quickly and accurately.

[0015] This method also achieves high-precision solution prediction of the fluid dynamics control equations by adopting high-order numerical discretization methods, improves the calculation accuracy, and provides a high-precision detonation solver for the study of ship structure explosion and shock resistance. BRIEF DESCRIPTION OF THE DRAWINGS

[0016] Figure 1 It is a method flow chart of a method for predicting a condensed charge detonation process according to an embodiment of the present application.

[0017] Figure 2 This is a comparison chart of the mixed density results obtained by the prediction method of this application and the Liu.GR detonation model in the first comparative example.

[0018] Figure 3 This is a comparison chart of the mixed pressure results obtained by the prediction method of the present application and the Liu.GR detonation model in the first comparative example.

[0019] Figure 4 This is a comparison chart of the mixing velocity results obtained by the prediction method of the present application and the Liu.GR detonation model in the first comparative example.

[0020] Figure 5 This is a comparison chart of the results of the gas phase unit mass internal energy obtained by the prediction method of the present application and the Liu.GR detonation model in the first comparative example.

[0021] Figure 6 This is a comparison result of the distribution cloud map of the mixed density obtained by the prediction method of the present application and the AUTODYN commercial software in the second comparative example.

[0022] Figure 7This is a comparison result of the distribution cloud diagram of the mixed pressure obtained by the prediction method of the present application and the AUTODYN commercial software in the second comparative example.

[0023] Figure 8 It is a comparison diagram of the density time history curves at two measuring points obtained by the prediction method of the present application and the AUTODYN commercial software in the second comparative example.

[0024] Fig. 9 It is a comparison chart of the pressure time history curves at two measuring points obtained by the prediction method of the present application and the AUTODYN commercial software in the second comparative example. DETAILED DESCRIPTION

[0025] The specific implementation of the present application is further described below in conjunction with the accompanying drawings.

[0026] The present application discloses a method for predicting the detonation process of condensed charge based on the ZND theory. The method is based on the ZND theory. The method for predicting the detonation process of condensed charge first divides the condensed charge package into a plurality of package units. The division can be customized according to actual needs. After the division, the unit equivalent size δ of each package unit can be determined. eff .

[0027] Then the parameters are initialized, including the charge model type of the condensed charge package, the type of explosive loaded inside, the detonation velocity D CJ , initial density ρ 0 、The initial internal energy of the explosive e 0 , the initial value of the combustion reaction factor λ 0 First, determine the type of explosive. Once the type of explosive is determined, its corresponding detonation velocity and initial density ρ can be found in public literature or manuals. 0 、The initial internal energy of the explosive e 0 The initial value of the combustion reaction factor λ 0 The explosive model type of the condensed charge package is any one of a two-dimensional plane model, a two-dimensional axisymmetric model, a one-dimensional plane model, a one-dimensional cylindrical model, and a one-dimensional spherical symmetric model. In the present application, the type of explosives loaded into the condensed charge package is commonly TNT or PBX9404.

[0028] After completing the parameter initialization, each drug package unit can be used as the minimum calculation unit for iterative calculation in time steps. The calculation method for each drug package unit at each time step is the same. Taking the calculation process of any drug package unit in the condensed matter drug package at any nth time step as an example, it includes the following steps, please refer to Figure 1 The flowchart shown:

[0029] Step 1: According to the detonation velocity D CJand the unit equivalent size of the drug package unit δ eff Determine the combustion reaction factor λ of the charge unit at the current time step n , integer parameter n≥1.

[0030] As the detonation process progresses, the detonation wave will undergo an evolution process such as generation and development. The explosives at each charge unit will gradually react and produce detonation gas. In the initial stage of the reaction, the charge unit is completely in a solid explosive state, and in the complete reaction stage, the charge unit is completely in a detonation gas state. It can be considered that each charge unit is in a mixed state of solid explosive and detonation gas, and this mixed state is constantly changing with the evolution of the detonation wave. The present application characterizes the evolution of the mixed state by the change of the combustion reaction factor, so the combustion reaction factor is iteratively updated at each time step.

[0031] In one embodiment, the combustion reaction factor λ at the current nth time step is n for:

[0032]

[0033] Among them, λ 0 is the initial value of the combustion reaction factor, and Δt is the duration between two adjacent time steps.

[0034] Step 2: According to the mixture density ρ of the previous time step n-1 Combined combustion reaction factor λ n Determine the gas phase density at the current time step

[0035] The combustion reaction factor characterizes the solid-gas mixing state. Therefore, the mixed density of the charge unit is the density of the solid explosive and the detonation gas mixed according to the weighted calculation of the combustion reaction factor. Therefore, combined with the combustion reaction factor λ n To determine the gas phase density In one embodiment, the gas phase density at the current nth time step is The calculation formula is:

[0036]

[0037] Among them, ρ s is the density of the explosive in the condensed charge, which can be determined after the type of explosive is determined. When the type of explosive in the condensed charge is TNT, the density of the explosive is s =1630kg / m 3 When the type of explosive in the condensed charge is PBX9404, the density of the explosive is s =1842kg / m 3 .

[0038] The known ρs , n-1 And the λ obtained in the previous step n Substituting into formula (2), we can get the gas phase density at the current nth time step:

[0039] Step 3: Based on the gas phase density at the current time step Combined with the gas phase unit mass internal energy e of the previous time step n-1 Determine the detonation gas pressure at the current time step using the JWL state equation

[0040] Determine the detonation gas pressure at the current time step using the JWL state equation for:

[0041]

[0042] in, is a relative specific volume and

[0043] The form of the JWL state equation was proposed by Lee based on the work of Jones and Wilkins, when the state equation includes unknown state equation parameters A, B, R1, R2, and ω.

[0044] In this application, the above-mentioned state equation parameters are related to the type of explosive in the condensed charge. This application determines the specific values ​​of the state equation parameters for two commonly used types of explosives by combining cylinder explosion tests and two-dimensional fluid elastic-plastic numerical calculation programs:

[0045] When the type of explosive in the condensed charge is TNT, the state equation parameters A = 371.2 GPa, B = 3.21 GPa, R1 = 4.15, R2 = 0.95, ω = 0.3, and the density of TNT ρ s =1630kg / m 3 .

[0046] When the type of explosive in the condensed charge is PBX9404, the state equation parameters A = 852.4 GPa, B = 18.02 GPa, R1 = 4.6, R2 = 1.3, ω = 0.38, and the density ρ of PBX9404 s =1842kg / m 3 .

[0047] The values ​​of the state equation parameters for two different types of explosives are shown in the following table:

[0048] State equation parameters TNT PBX9404 A(GPa) 371.2 852.4 B(GPa) 3.21 18.02 <![CDATA[R1]]> 4.15 4.6 <![CDATA[R2]]> 0.95 1.3 ω 0.3 0.38

[0049] When the type of explosive is TNT or PBX9404, the above determined state equation parameters are directly substituted into the JWL state equation for calculation.

[0050] Step 4: According to the detonation gas pressure Combined combustion reaction factor λ n Get the mixed pressure p of the current time step n Similarly, since each charge unit can be regarded as a mixture of solid explosive and detonation gas, the detonation gas pressure is determined by Finally, the pressure of the solid explosive must be combined with the combustion reaction factor λ n Perform weighted calculation to obtain the mixed pressure p of the mixed state n .

[0051] The mixed pressure p of the drug package unit at the current nth time step n for:

[0052]

[0053] Among them, p s is the pressure of the explosive in the condensed charge, which can generally be directly taken as p s =0.

[0054] Step 5: The obtained mixed pressure p n Substitute it into the fluid dynamics control equations constructed based on the unsteady compressible Euler equations and solve it to obtain the mixing density ρ at the current time step n and the gas phase unit mass internal energy e at the current time step n And enter the solution of the next time step.

[0055] Since the detonation process is extremely short, the transport effects such as viscosity, heat conduction and diffusion in the detonation process can be ignored within the error range. The fluid dynamics control equations are constructed based on the unsteady compressible Euler equations. The calculation formula of the fluid dynamics control equations at any nth time step is:

[0056]

[0057] in,

[0058] In the above formula, is the mixing velocity u of the drug package unit at the current time step n The component in the r direction, is the mixing velocity u of the drug package unit at the current time step n The component in the z direction. The r direction is the radial direction of the condensed charge, the z direction is the axial direction of the condensed charge, and t is time.

[0059] N is a parameter related to the charge model type of the condensed charge package. When the charge model type is a two-dimensional plane model, N = 0 and When the charge model type is a two-dimensional axisymmetric model, N = 1 and When the charge model type is a one-dimensional plane model, N = 0 and When the charge model type is a one-dimensional cylindrical model, N = 1 and When the charge model type is a one-dimensional spherically symmetric model, N = 2 and

[0060] So the resulting mixed pressure p n Substituting into formula (5), the mixing density ρ at the nth time step can be obtained by solving the fluid dynamics control equations: n , total energy per unit mass E n and the mixing speed u n It is generally believed that the speed of the gas phase and the solid phase are both the mixing speed u n , and the initial value of the mixing speed is generally zero. Then further follow Calculate the gas phase unit mass internal energy e at the nth time step n .

[0061] When solving the fluid dynamics control equations, high-order WENO reconstruction and HLLC approximate Riemann solver are used to spatially discretize the fluid dynamics control equations, and the third-order Runge-Kutta method is used to time discretize the fluid dynamics control equations.

[0062] Therefore, for each drug package unit, the mixing pressure p of the drug package unit can be determined and output at each time step. n , mixture density ρ n , total energy per unit mass E n , gas phase unit mass internal energy e n and the mixing speed u n , so that the state parameters of the detonation mixing zone at the charge unit can be determined. These state parameters are used to enter the iterative solution of the next time step. In this way, the iterative solution can output the state parameters of the detonation mixing zone at each charge unit at each time step, realize the prediction of the detonation process, and determine the propagation process of the detonation wave.

[0063] In addition, the detonation time t0 and the detonation speed D CJ Determine the detonation time t of each charge unit ig :

[0064]

[0065] Where x is the unit position vector of the charge unit, x d is the detonation point position vector, t0 is the detonation time of the charge, these are all known parameters, so combined with the detonation velocity D CJ The detonation time t of each charge unit can be determined ig .

[0066] Then we can calculate the detonation time t of each charge unit. ig Determine the detonation state of the charge unit at the current time step: When the nth time step reaches the detonation time t of the charge unit ig When the detonation state of the charge unit at the nth time step indicates that the charge unit has been detonated. ig When the detonation state of the charge unit at the nth time step is determined to indicate that the charge unit has not detonated, the prediction result obtained at each time step includes not only the state parameters of the detonation mixing zone at each charge unit, but also the detonation state of each charge unit.

[0067] The present application provides a concise and fast method for predicting the detonation process. Combined with the state equation parameters of two commonly used explosive types provided by the present application, the prediction results can be obtained quickly and accurately. In order to verify the effectiveness of the method of the present application, the detonation process prediction method of the present application is compared with the detonation calculation results of Liu.GR (Liu.GR, Liu MB, Smoothed particle hydrodynamics-a meshfree particle method. World Scientific Publishing, Singapore, 2010, 215-218.). When the charge model type is a one-dimensional plane model and the type of explosive is TNT, Liu.GR uses the SPH method to simulate the flow field distribution at a typical time within 10μs after detonation. The mixed density predicted by the method of the present application and the mixed density at different distances x from the detonation point obtained by the Liu.GR detonation model at 2μs, 4μs, 5μs, 6μs, and 10μs after detonation are compared as shown in the figure below. Figure 2 The mixed pressure predicted by the method of the present application and the mixed pressure at different distances x from the detonation point obtained by the Liu.GR detonation model at 2μs, 4μs, 5μs, 6μs, and 10μs after detonation are compared as shown in (a), (b), (c), (d), and (e). Figure 3 The comparison of the mixing velocity predicted by the method of the present application and the mixing velocity at different distances x from the detonation point obtained by the Liu.GR detonation model at 2μs, 4μs, 5μs, 6μs, and 10μs after detonation is shown in (a), (b), (c), (d), and (e). Figure 4 The results of the comparison of the gas phase unit mass internal energy predicted by the method of the present application and the gas phase unit mass internal energy at different distances x from the detonation point obtained by the Liu.GR detonation model at 2μs, 4μs, 5μs, 6μs, and 10μs after detonation are shown in (a), (b), (c), (d), and (e). Figure 5 As shown in (a), (b), (c), (d), and (e) in the method of the present application, the size of the drug package unit is 0.1 mm, and the size of the drug package unit in the Liu.GR calculation method is 0.025 mm. Figure 2 to Figure 5 It can be seen that the prediction results calculated by the present application are very close to the SPH calculation results of Liu.GR, which shows the accuracy of the prediction results obtained by the present application. In addition, the method of the present application reduces the calculation time required for 10μs after detonation by about 30% on a single computing node compared with the calculation time of the Liu.GR detonation model, and the accuracy is higher. The same is true for calculations at other times. Therefore, the method of the present application has good performance in both calculation speed and calculation accuracy.

[0068] In another example, when the charge model type is a two-dimensional plane model and the explosive type is TNT, the calculation domain is [-0.1, 0.1] × [-0.1, 0.1], and the surrounding solid wall boundary conditions are used. The prediction results of this application are compared with the commercial software AUTODYN. The distribution cloud diagrams of the mixed density obtained by the method of this application and AUTODYN at 5 typical moments are compared as follows: Figure 6 As shown, the comparison results of the distribution cloud diagram of the mixed pressure at 5 typical moments obtained by the present application method and AUTODYN are as follows Figure 7 As shown by Figure 6 and Figure 7 It can be seen that the prediction results of this application are in good agreement with the results of AUTODYN. In addition, the density time history curves of the two measuring points (0.04, 0) and (0.06, 0) in the calculation domain are as follows Figure 8 As shown in Figure 2, the pressure-time curves of these two measuring points are as follows: Fig. 9 As shown by Figure 8 and Fig. 9 From the comparison of the results, it can be seen that the AUTODYN software has obvious non-physical oscillations in the reflected wave part, and the method of the present application is superior to the AUTODYN software in terms of calculation accuracy and stability. In addition, the method of the present application reduces the calculation time by about 20% compared with the AUTODYN commercial software, and also has better performance in calculation speed.

[0069] The above is only a preferred embodiment of the present application, and the present application is not limited to the above embodiments. It is understood that other improvements and changes directly derived or associated by those skilled in the art without departing from the spirit and concept of the present application should be considered to be included in the protection scope of the present application.

Claims

1. A method for predicting the detonation process of condensed charge based on ZND theory, characterized in that: The method for predicting the detonation process of condensed charge comprises: time steps for each charge unit in the condensed charge: According to the detonation speed and the unit equivalent size of the drug package unit Determine the combustion reaction factor of the package unit at the current time step , integer parameter , is the initial value of the combustion reaction factor, is the duration of two adjacent time steps; According to the mixture density of the previous time step Combined combustion reaction factor Determine the gas phase density at the current time step ; According to the gas phase density of the current time step , combined with the gas phase unit mass internal energy of the previous time step Determine the detonation gas pressure at the current time step using the JWL state equation ; According to the detonation gas pressure Combined combustion reaction factor Get the mixed pressure at the current time step , and substitute it into the fluid dynamics control equations based on the unsteady compressible Euler equations to solve the mixed density of the current time step , the gas phase unit mass internal energy at the current time step and the mixed speed at the current time step And enter the solution of the next time step; among them, and All are initial values.

2. The method for predicting the detonation process of condensed charge according to claim 1, characterized in that: Determine the detonation gas pressure at the current time step using the JWL state equation for: in, is a relative specific volume and , is the density of the explosive in the condensed charge; , , , , They are all state equation parameters and are related to the type of explosive in the condensed charge.

3. The method for predicting the detonation process of condensed charge according to claim 2, characterized in that: When the type of explosive in the condensed charge is TNT, the state equation parameter , , , , , and the density of TNT .

4. The method for predicting the detonation process of condensed charge according to claim 2, characterized in that: When the type of explosive in the condensed charge is PBX9404, the state equation parameter , , , , , and the density of PBX9404 .

5. The method for predicting the detonation process of condensed charge according to claim 1, characterized in that: According to the mixture density of the previous time step Combined combustion reaction factor Determine the gas phase density at the current time step Including according to the following calculation formula: in, is the density of the explosive in the condensed charge.

6. The method for predicting the detonation process of condensed charge according to claim 1, characterized in that: According to the detonation gas pressure Combined combustion reaction factor Get the mixed pressure at the current time step Including according to the following calculation formula: in, It is the pressure of the explosive in the condensed charge.

7. The method for predicting the detonation process of condensed charge according to claim 1, characterized in that: The fluid dynamics control equations based on the unsteady compressible Euler equations are The calculation formula for a time step is: in, , , , ; is the mixing speed of the drug package unit at the current time step exist The weight in direction, is the mixing speed of the drug package unit at the current time step exist Directional weight; The direction is the radial direction of the condensed charge package. The direction is the axial direction of the condensed charge package. It’s time; It is a parameter related to the charge model type of the condensed charge package. When the charge model type is a two-dimensional plane model and , when the charge model type is a 2D axisymmetric model and , when the charge model type is a one-dimensional plane model and , when the charge model type is a one-dimensional cylindrical model and , when the charge model type is a one-dimensional spherically symmetric model and ; The first Mixed pressure at each time step Substituting into the fluid dynamics control equations, we get The mixture density at time steps , total energy per unit mass and mixing speed ,according to Calculate the first The internal energy per unit mass of the gas phase in time steps .

8. The method for predicting the detonation process of condensed charge according to claim 1, characterized in that: Solving the governing equations for fluid dynamics includes: The fluid dynamics control equations are spatially discretized using high-order WENO reconstruction and HLLC approximate Riemann solver, and the fluid dynamics control equations are temporally discretized using the third-order Runge-Kutta method to solve them.

9. The method for predicting the detonation process of condensed charge according to claim 1, characterized in that: The method for predicting the detonation process of condensed charge further comprises: According to the detonation time of the explosive pack and explosive speed Determine the detonation time of each charge unit ;in, is the unit position vector of the charge unit, is the detonation point position vector, It's explosive speed, It is the moment when the explosive charge is detonated; According to the detonation time of each charge unit Determine the detonation state of the charge unit at the current time step; time steps to reach the detonation moment of the charge unit When the charge unit is at The detonation state of the time step indicates that the charge unit has been detonated; The detonation time of the charge unit has not been reached by time steps When the charge unit is at The detonation state of the time step indicates that the charge unit has not detonated.

Images