Paper received: 2011-12-20, paper accepted: 2012-03-21 © 2012 Journal of Mechanical Engineering. All rights reserved. A Model for Shaped Charge Warhead Design Slobodan Jaramaz1 - Dejan Mickovic1.* - Predrag Elek1 - Dragana Jaramaz2 - Dušan Mickovic1 1 University of Belgrade, Faculty of Mechanical Engineering, Serbia ^University UNION-Nikola Tesla, Faculty of Civil Management, Serbia A model for shaped charge warhead design was developed. The model is incorporated in the computer code - CUMUL. The code includes detonation wave profile estimation, liner collapse, arrival of collapsed liner to the centerline of shaped charge, jet creation and jet breakup. The penetration phenomena are discussed and governing equations are presented. Two cases dealing with the target type are included: homogenous and non-homogeneous targets. For the purpose of verifying CUMUL, a set of 20 specimens of shaped charges was tested. The tests were directed to investigate the effect of cone apex angle and stand-off distance on the performance of shaped charge. From the comparison between experiments and CUMUL results, it was concluded that CUMUL program shows a good agreement with the experiments. That enables it to be a powerful tool for shaped charge warhead design. Keywords: shaped charge warhead, shaped charge jet, warhead design, computer code, penetration 0 INTRODUCTION Shaped charges are extremely useful when an intense, localized force is required for the purpose of piercing a barrier. The main application is in the military arena, for high explosive antitank (HEAT) rounds including hand-held (bazooka type) rounds, gun-launched rounds (e.g., rifle grenades), cannon-launched rounds, and various bombs. The targets are armors, bunkers, concrete or geological fortifications, and vehicles [1] and [2]. The shaped charge was analysed using an analytical approach for preliminary analysis and parametric studies to determine an approximate design that could satisfy technical requirements. For this purpose we developed models for the following phases: a) estimation of explosive properties, b) detonation wave properties and profile, c) calculation of liner driven velocity, d) calculation of liner collapse velocity and angle, e) jet length determination, f) estimation of target penetration. These models are included in CUMUL computer code. CUMUL program calculations are compared with experimental results that include the study of liner apex angle and stand-off distance influence on the penetration of 64 mm anti-tank rocket with shaped charge. 1 STRUCTURE OF CUMUL PROGRAM The main outline and different approaches adopted in CUMUL are shown in block diagram in Fig. 1. Performing the shaped charge design phases by utilizing the mostly known approaches, make CUMUL a very powerful design tool. In addition, it also provides the designer with a great chance to have a wide range of calculated results which indicate the expected performance for the designed shaped charge. 1.1 Input Data Input data needed for running CUMUL can be classified as follows: • explosive input data, • shaped charge liner shape and dimensions, • target data, • options for approaches to be used for calculations. 1.2 Explosive Properties Determination The explosive properties such as Chapman-Jouguet pressure, detonation velocity and Gurney constant are calculated by the use of empirical equations. The explosive name or names if it is a mixture, densities and percentage are given in the input file. Another data base file contains basic characteristics of 28 explosives. 1.3 Estimation of Detonation Wave Profile Determination of detonation wave attack angle and jet mass calculations are fully depended on the liner shape used. The model of logarithmic spiral is used as a detonation profile technique. Due to the variety of liner shapes, separated subroutines were created, so each liner shape has its own subroutine. There are four subroutines for four different liner shapes: (1) conical, (2) parabolic, (3) biconic, (4) Gaussian. 1.4 Evaluation of Initial Driven Liner Velocity When the detonation wave arrives at the liner, the liner element will accelerate to an initial velocity V0. In seeking all the approved ideas for this phase, CUMUL provides four approaches for V0 calculations: 1. Asymmetric sandwich. Fig. 1. Block diagram for CUMUL program V 42Ë M N (Xn - XL ) C + C X2 1 X 3 X. 1 + (( - Xl )3 x3 -1/2 (1) where: ■JlE Gurney constant, M liner mass per unit area, C explosive charge mass per unit area, N confinement mass per unit area, Xm location of liner surface, Xn location of confinement surface, Xx location of stationary surface: Xs = Xn (1/2 + N/C)(M/C + N/C +1). 2. Gurney formula for imploding cylinder Vo IM + N (Rn - Rs) + 1 c c (R - Rm) + 1 (( - Rm )(3Rm + Rs ) + 6 R - R 2 +1 (Rn - Rs )3 (3R + Rs ) 6 ( -Rl )( -Rm ) -12 (2) where Rm radial location of liner surface (Fig. 2) Rn radial location of confinement surface Rs Lagrangian radial position of an assumed stationary surface: R; + 3R„ (( + Rm)| CR + NR -3RmRn (Rn + Rm ) 3. Chanteret formula [1]: V 4ÏE '( 02 2 A R - Rm v rs - Rm j M N_ C C M 1 — + — C 6 -RR = 0. -V2 (3) 4. Hirsch formula [2]: Vn 42Ë where: A ( N + 3ß +1 ï_ A + M + ß + 3 C 6ß + 6 J 3 C 6ß + < M ß+ß2 A = CP 3ß+3 ß= R. A N + 2ß +1 ' P Rm • -1/2 (4) C 3ß + 3 Fig. 2. Control mass for concentric cylinder liner 1.5 Evaluation of Liner Collapse A diagram of the general charge geometry and liner collapse is shown in Fig. 3. The explosive is detonated at B. As the detonation wave passes P, the liner element originally at P begins to collapse in a direction that makes an angle 5 with the normal to the original liner at P. The detonation wave speed D is considered normal to the wave front. However, the velocity with which the wave sweeps the liner is not constant , however, and is given by the formula: u (x ) = d cosy (5) where y = y(x) is the angle between the normal to the detonation wave at P and the tangent to the liner at P (see Fig. 3). Fig. 3. Charge geometry and collapse The projection angle for the unsteady case is: V 1,1, 5= ■V0- --tV + -tV , 2D 2 4 (6) where the prime indicates differentiation with respect to the Lagrangian liner coordinate x. The symbol t is time parameter related to the acceleration of the liner. When the liner element collapses, the driven velocity history is assumed to follow one of the following profiles (Fig. 4): T t C. Exponential Acceleration Fig. 4. Liner collapse velocity profile; a) instantaneous acceleration, b) constant acceleration, c) exponential acceleration The liner elements in the classical theory were assumed to reach collapse velocity instantaneously (Vc = V0). The first level of refinement assumes that velocity increases linearly over a short period until it reaches final velocity V0 or collapses on the axis. The velocity history that assumes an exponential form was proposed in [5]. The time base used is t = 0 when the detonation wave is at x = 0 on the liner. T = T(x) is the time when the detonation wave reaches the element x on the liner. The three velocity profiles were involved in CUMUL. For example, for each initial driven liner velocity (i.e. asymmetric sandwich, Gurney, Chanteret or Hirsch), the three categories are applied to cover all the assumed possibilities that describe the way of movementand traveling path, which the liner element had followed [6]. In order to develop the contour of the collapsing liner we write the coordinates of general point P' in Fig. 3, as it collapses at time t to the point M. The coordinates of point M are: z = x +1 (x, t )sin (a +5), (7) r = R (x)-/(x,i)cos(a + 5), (8) where z is the axial coordinate, r is the radial coordinate, R is the original liner radius, and l(x,t) is the distance the element has travelled from P' to M. The angle of impact of the liner with the axis is given by: tan ß = , dl (x, tc ) ( . R---—— cos (a + 5) _dx 1 J 1 sin (a + 5) + R (a+ 5) R (a + 5)tan (a + 5) 1 +Ôl(^ sin (a+5) + R (a+ 5) + , (9) where tc = tc(x) is the time of impact and the prime indicates differentiation with respect to x. The derivative dl/ dx is evaluated at a proper impact time. Once fi has been calculated the velocities of each element of jet and slug are calculated from: v. = _v_ J • ß sin — V = V ß cos— 2 cos sin a + ö-ß |, ß (10) (11) where V = Vc(x,tc). A ring element of the liner of mass dm splits into an element of jet of mass dmj and an element of slug dms. These masses are defined by: dm dm dms dm J = sin2—, ß 2 = cos2 ß. (12) (13) 1.6 Evaluation of Jet Forming and Stretching During this phase, the first step is to consider the jet forming case, that i, at the moment when the collapsed liner elements arrive to center line of shaped charge. This consideration does not include stretching and particulation of the jet. At the moment of jet formation the initial jet length, mass, diameter, arrival time and arrival velocity are calculated. On the other hand, the same procedure is applied to calculate the slug properties as well, but, because the slug part does not contribute in penetration phenomena, it is coming out of the picture of interest. In order to trace the location of the jet elements and study their stretching after formation the coordinate system conventions shown in Fig. 5 are adopted. Fig. 5. Description of coordinates Lagrangian coordinate x defines the original position of the liner elements along the axis. The position coordinate in the jet is measured from the original position of the liner apex and is designated as £(x,t).The position of an element at the moment it just reaches the axis is denoted by z (x). From the basic collapse geometry, (Fig. 3), the position is given by: z (x) = x + R tan (a +5). (14) At any time t > tc a portion of the element originally at x will be in the jet. Assuming that each element of the jet travels at a constant velocity Vj(x) immediately after it is formed, the position of the x element at time t is: ^(x,t) = z + (t - tc )Vj, t > tc. (15) The one-dimensional extension of the jet elements may now be defined in the following manner. Consider two points xx and x2 on the liner separated by a distance Dx as shown in Fig. 6. Fig. 6. Extension of a jet element During the collapse process the point x: reaches the axis and then proceeds to jet along the axis. + Meanwhile, x2 has begun to collapse and reaches the axis at time t0 = tc(x2). At this moment t0, jetting material originally at xl is located at £(xbt0) and x2 has just reached the axis and is located at £(x2,t0). At some arbitrary later time t, their locations are £(xbt) and £(x2,t), respectively. One-dimensional jet extension is defined as the ratio of the increase in length of a jet element to its length when first formed: AI-A4 E = lim.. A^o , )/g(xl5.0 ) (16) -1, where to = tc(xj). The second step is to consider jet stretching and particulation phenomena. The main outlines of these phenomena can be summarized as following: A. Jet Tip Calculation. It is a start point for jet forming calculation. Jet tip mass, length, velocity and number of liner elements involved in jet tip creation are founded, keeping in mind, that, the tip is not exposed to stretching [7]. In many cases the liner has a region where the elements do not reach the final collapse velocity. In this region, close to the charge axis, an element with jetting velocity Vjl may be followed by an element with greater jetting velocity (Vj2>Vji). This inverse velocity profile usually continues throughout this region and the mass piles up forming the jet tip. Each element is considered to impact until the first jetting element whose velocity is less than the velocity of the combined tip projectile. Then, conventional jetting ensures. Assuming a pefectly plastic impact of elements the conservation of linear momentum leads to the following expression: V (p )= ■"tip i V (x)■ dx -dx ' dm j dx (17) dx where Vj is velocity of the combined tip particle, and dmj/dx mass of jet element per unit lenght of original cone. Eq. (17) is integrated step by step until a point xtip is found such that: Vj Up )< V (x ). (18) The value of xtip is considered as the point on the liner that distinguishes where the formation of the tip stops and where normal jetting begins. This value of xtip is presented in Fig. 7. Fig. 7. Typical jet velocity distribution curve B. Breaking up Time. Breaking up time is the maximum time spent in jet stretching operation before jet particulating take place. CUMUL provides two formulas for estimating the possible breaking up time for shaped charge [8] to [10]. C. Jet Stretching. After finding break up time, the maximum jet stretching is calculated, and therefore, the final or the total jet length is founded [11]. 1.7 Penetration Referring to Fig. 1, it can be seen that the penetration phase is created according to the following classification: 1.7.1 Homogenous Target When the target consists of one material it is called a homogenous target. In this case, there are two possible ways of calculations: a) No virtual origin is applied. In this case the stand off distance is taken as it is provided by the user. b) Applying virtual origin approach [12]. A virtual origin is an important issue for determining the stand off distance. By finding the location of virtual origin, the stand off distance can be easily achieved by simple addition operation. The penetration value for both cases can be calculated by the following techniques (this can be decided by the user) [13]: 1. density law formula (DL), 2. minimum jet velocity (Vmin), 3. minimum penetration velocity (Umin). 1. For a jet of constant velocity, assuming that the penetration stops when the jet length is consumed, the penetration is given by the density law formula: (19) where Lj is jet length, pj liner density, and pt target density. 2. For a jet of non-uniform velocity distribution, the jet length is not constant but increases with time. Three cases are considered: a) Penetration before jet break-up P = S v.. Jtp V ■ V mm / -1 (20) where S is stand-off distance, V]tip jet tip velocity, and Vmin minimum jet velocity capable to penetrate the target material y = (pT / Pj)1/2. In this case S is bounded by: 0 < ^ < VmJb K. V.. V j"p / where tb is breaking up time of the jet. b) Jet breaks during penetration 1 Y P _(1 + b S1+Y- Vmntb - (21) Y ' where S is bounded by: V ■ tb min b V V V Jt'p y < ^ < Vjtiptb. c) Jet breaks before reaching the target ((jtip — ^min P = (22) for stand-off in the range Vjtip tb < S < . 3. Formulae for penetration based on minimum penetration velocity (Umin) are similar to those given for minimum jet velocity (Vmin). 1.7.2 Non-Homogenous Target A non-homogenous target is defined as the target which consists of many layers of different material. Due to the target non-homogeneity a target resistance factor was defined. By using this factor the expected penetration value can be calculated using the formula [14]: P = 1 - Pt 2R Pj i - PT P j j (23) PT 2R p/ j j'ip 1 - Pt_ P j J Pr_ Pj where pt is target density given by pH /pj , ph hydrodynamic density, and R target resistance factor. 2 EXPERIMENTAL WORK AND MODEL VERIFICATION The main aim from the present experimental work is to measure the penetration caused by specified shaped charge. This includes the study of changing cone apex angle and stand off distance on the resulting penetration. Fig. 8. Conical shaped charge (apex angle = 50°) Two models of shaped charge with conical liner were implemented in experimental work. The first model with apex angle 2a = 50°, second model has an apex angle 2a = 60°. The defined shaped charge having a 64 mm diameter and HMX explosive material was chosen. The two models are shown in Figs. 8 and 9, respectively. A hard cylindrical paper was used to control the standoff distance. Three steel plates with 300, 10 and 25 mm thickness were used as a target. Three shaped charges with different standoff distance are shown in Fig. 10. A set of 20 conical shaped charges were divided into 4 groups. r r Y Y A complete analysis and comparison of the theoretical and experimental results were performed. Three different acceleration histories were implemented. Penetration results were calculated using density law, minimum jet velocity and minimum penetration velocity approaches. Four approaches for calculating liner collapse velocity were involved. 30 10^ >1.4 R 8.34 V1.2 Q o