Hui-jun NING*，Nn-peng FENG，Tn-huiWU，Cheng ZHANGHo WANGShool of Energy nd Power Engineering，Nnjing University of Siene nd Tehnology，Nnjing 210094，Jingsu，Chin

bThe Graduate School of Hubei Aerospace Technology Academe，Wuhan 430000，Hubei，China

cDept.of Engineering Mechanics，Shanghai Jiaotong University，Shanghai 200240，China

Received 27 April 2014；revised 24 October 2014；accepted 24 October 2014 Available online 5 January 2015

Semi-analyticalmethod for calculating aeroelastic effect of profiled rod flying at high velocity

Hui-jun NINGa，*，Nan-peng FENGb，Tan-huiWUc，Cheng ZHANGa，Hao WANGa
aSchool of Energy and Power Engineering，Nanjing University of Science and Technology，Nanjing 210094，Jiangsu，China

bThe Graduate School of Hubei Aerospace Technology Academe，Wuhan 430000，Hubei，China

cDept.of Engineering Mechanics，Shanghai Jiaotong University，Shanghai 200240，China

Received 27 April 2014；revised 24 October 2014；accepted 24 October 2014 Available online 5 January 2015

The key technique of a kinetic energy rod（KER）w arhead is to control the flight attitude of rods.The rodsare usually designed to different shapes.A new conceptual KER named profiled rod which has large L/D ratio is described in this paper.The elastic dynam ic equations of this profiled rod flying athigh velocity after detonation are setup on the basisof Euler-Bernoullibeam，and the aeroelastic deformation of profiled rod is calculated by semi-analyticalmethod for calculating the vibration characteristics of variable cross-section beam.In addition，the aeroelastic deformation of theundeformed profiled rod and theaeroelastic deformation of deformed profiled rodwhich is caused by the detonation of explosiveare simulated by computational fluid dynamic and finite elementmethod（CFD/FEM），respectively.A satisfactory agreementof these tw o methods is obtained by the comparison of two methods.The results show that the sem i-analyticalmethod for calculating the vibration characteristics of variable cross-section beam is applied to analyze the aeroelastic deformation of profiled rod flying at high velocity.

Profiled rod；Detonation process；Variable cross-section beam；Aeroelastic deformation；Numerical analysis

1.Introduction

The advanced technologies of kinetic energy rod（KER）warheads have been developed w idely.According to the different threat targets，KER warheadsare usually divided into continuous rod warhead and discrete rod warhead.The advantages of these warhead devices are the high speed and continuous cutting capacity for threat aircraft and cruise missiles［1-5］.

However，the thicker components，which are internal submunition components of missiles，cannot be destroyed by traditional KER warheads in themissile-defense environment. Therefore，it is necessary and urgent for a new warhead mechanism w ith a strong killing effect against the key components of thick targets.

This paper describesa new conceptual profiled rod warhead which is different from the conventional KER warhead.The rods in the warhead are changed into variable cross-section straight rods with large L/D ratio，which we call profiled rods.The conventional rods are constant cross-section cylindrical rods and square rods.The profiled rod is designed to truncated coneshaped rod，as shown in Fig.1.Since the shape of profiled rod is changed，the mass is unevenly distributed along the axial direction of the rod，and the velocity gradient or difference between the head and aft end of the rod is generated，which causes the differentobliquity angles to strike a target.In this process，the rods initially rely on its high velocity to penetrate the target.A lthough the rod velocity decreases under the influence of the air resistance after the detonation，it can provide a better penetration angle that still maintains a high penetrability against targetswhen itmovesfor a long range.On the other hand，the profiled rod warheads also have a continuous cutting capacity for threataircraftand cruisemissilesby setting thedifferentanglesbetween rodsand explosive compared w ith the traditional KER warheads described in Refs.［6-8］.

The formation of profiled rod warhead，the initial dep loym ent velocity，the damage effect and the flight attitude of rod after detonation are analyzed in the developmenteffort.In this paper，we focus on the latter.After the detonation，the rods obtain high deployment velocity and fly against the air resistance until they strike the target.It should be pointed out that the aeroelastic deformation of rod which is caused by overcoming the aerodynamic force at high speed has influence on flight attitude and damage effect.Consequently，an investigation on elastic deformation of profiled rod under the aerodynamic force is very important.

Fig.1.Structural representation of profiled rod.

Fig.2.Initialmotion statemap of a profiled rod.

2.Design princip le of profiled rod

A profiled rod w ith diametersof 3mm and 5mm，length of 98 mm and length-to-diameter ratio of 24.5 is proposed.Its model is shown in Fig.1.

An initialmotion stateof theprofiled rod after detonation is assumed，as shown in Fig.2.In an analysis on the elementof one end of the profiled rod，v0is the radial velocity of profiled rod，vris the linear velocity around the y-axis，and vc0is the center ofmass velocity.

Assume thatαis the incidence angle，Lris the length of profiled rod，ω0is the initial transverse angular velocity，ω1is the initial turning angular velocity，and r and R represent the radii of the top and bottom of profiled rod，respectively.The relation betweenαand vris deduced by

whereθisw ritten as

here，

Thus the initial transverse angular velocityω0can be w ritten as

andω1can bew ritten as

Substituting Eq.（2）into Eq.（4），ω1can be w ritten as

3.Elastic dynam ic equations of p rofiled rod

3.1.Assumptions

According to Refs.［9，10］，we take the follow ing assumptions.

1）Rod bending due to detonation is not taken into account.

The profiled rod is considered as elastic rod，which accords w ith the assumptions of Euler-Bernoulli beam belonging to variable cross-section Euler-Bernoullibeam. 2）Gravity is not taken into account，and the effect on acceleration of gravity is ignored.

3）The fixed-axis rotation of the profiled rod hasoccurs in the initial state.

4）The effectson tip resistance and viscous friction of the rod are neglected，and the aerodynam ic drag is perpendicular to the axis of the profiled rod.

5）The elastic deformation ispresumed not to bring about the changes in aerodynam ic drag in order to resolve the corresponding vibration used rod aerodynamic force as external-exciting force.The elastic deformation is loaded to the elastic rod，and then the aerodynamic drag can be achieved.Ultimately，the flight trajectory of the rod is obtained by repeating above-mentioned process.3.2.Establishment of elastic dynamic equations

The flightmotion of the profiled rod is made up of the fixed-axis rotation about the center ofmassand the translation about the center of mass，as shown in Figs.3 and 4.The differential equation for transverse vibration on the rotational plane and translation plane and the differential equation for longitudinalvibration on the translation plane can be setup by New ton's second law ofmotion combined w ith the theory of Euler-Bernoullibeam.The elastic dynam ic equations can be w ritten as

where u（x，t），v（x，t）are the longitudinal and transverse displacements deviating from the origin of the cross-section x at themoment t，respectively；ρis the density of rod；A（x）is the area of the cross-section；E is Young's modulus；ω0is the angular velocity of the rod on rotational plane；q（x，t）is the transverse force（aerodynam ic drag）of the rod in a unit length on rotation plane；I（x）is the second moment of area w ith respect to the neutral axis of the cross-section；w（x，t）is the transverse displacement from the origin of the cross-section x at the moment t；and q（x，t）is the transverse force（aerodynam ic drag）of the rod in a unit length on translation plane.

3.3.Semi-analyticalmethod for calculating transverse

vibration characteristics of profiled rod

Eqs.（6）-（8）are used to calculate the force vibration of variable cross-section beam，and the vibration theory isused to acquire u（x，t），v（x，t），w（x，t）［11，12］.

Take an example of equation for transverse vibration of profiled rod on the rotational plane.Obviously，Eq.（7）is the force vibration about lateraldeformation v（x，t）under the force q（x，t）.Based on the theory of Euler-Bernoulli beam，the equation for free vibration of the rod w ith length of L，direction of transverse vibration is

where v（x，t）is the function of deformation；EI（x）andρA（x）are the transverse bending stiffness and linear density，respectively.

The profiled rod is divided into a finite number of homogeneous stepped beam s w ith constant cross-section by using the finiteelementmethod，both theendsof thebeam are called node.Therefore，the profiled rod is changed into an aggregation with finite beam elements interconnected at the node，and the deflection of arbitrary element can be described as linear combinationscomposed by nodaldisplacementso thatwe can obtain the approximate solution of the profiled rod while the segments are more enough.The profiled rod in Fig.5 is divided into a finite number of homogeneous stepped beams w ith constant cross-section by using the finite elementmethod，as shown in Fig.6.

Introduce（EI）iand（ρA）ias follows

where liis the length of the i th section beam；（EI）iis the equivalent longitudinal bending stiffness；and（ρA）iis the equivalent linear density.

Themodal function of the i th section beam can be obtained on the basis of vibration mode for constant cross-section as where Xi（x）=βi（x-xi-1），xi-1≤x≤xi，（i=1，2，…，N），x0=0；Ai，Bi，Ci，Diare the undefined coefficients of the i th section beam；besides，

Fig.3.Rotationmotion of profiled rod.

whereωis the natural lateral vibration frequency of variable cross-section beam.

Under the continuous conditions of deflection，rotation，moment and shear force at the point of xibetween the i th section beam and the（i+1）th section beam.We can acquire the relationships as follows

Substituting Eqs.（11）and（13）into Eq.（14），we have

and

Z（i）can bew ritten as

where

These recursive relations are used to establish a new algorithm for A（1）and A（N）

The relationship of undeterm ined coefficient between the first equivalent beam and the N th equivalent beam is built by matrix Z where all the elements are the function of natural frequency.Consequently，the solution of modal function can be obtained from Eqs.（11）and（18）by solving the natural frequency w ith four boundaries at both ends of the beam.

After the detonation，the profiled rods fly rapidly in parallel to the axes of the warhead，as the result of the geometric asymmetry，themotion performances of the rods are complicated.The different flight attitudes of the rods were demonstrated in forward flight.Two boundary conditions were presented to illustrate the aeroelastic response by terms of different flight attitudes.

Fig.5.Variable circular cross-section beam.

Fig.4.Translationmotion of profiled rod.

3.4.Boundary conditions

3.4.1.Clamped-clamped beam

For instance，the casewhere thebeam isclamped atboth its ends is considered，as shown in Fig.7.The governing boundary conditions are as follows

Eqs.（19）and（20）are rew ritten as the form ofmatrix

where

The determinants of coefficient matrix must be equal to zero if Eq.（16）has a solution.

W hen characteristic Eq.（22）is a nonlinear functionω，it can be solved by New ton-Raphson iteration extensions to obtain the natural lateral vibration frequency in responding boundary conditions.Thenωis substituted into Eqs.（11）and（15）to obtain the principalmode of the beam.Furthermore，the elastic deformation can be also solved by amode superpositionmethod［13，14］.

3.4.2.Free-clamped beam

Here，the differences between this condition，as shown in Fig.8，and the above-mentioned com putation are the boundary conditions.

where

Fig.6.Stepped beam w ith N segments.

Fig.7.Clamped-clamper beam.

Fig.8.Free-clamper beam.

4.Num erical sim ulation method

A two-way coupling approach isproposed in order to solve the static aeroelastic problem of a profiled rod to conduct this static aeroelastic analysis，ANSYSWorkbench multi-physics coupling platform，which includes fluid，structural solvers and couplingmodule，namely FLUENT，ANSYSMechanical and System coup ling，is used to solve the disp lacem ents associated w ith the aerodynam ic pressure loading and find out the static aeroelastic behavior of this variable cross-section rod.CFD grid is deformed to obtain the aerodynamic solutionsof the deformed geometry using dynamicmeshwhich is called diffusion-based smoothing method.For the structural analysis，the ANSYSMechanical inertia relief option［15，16］，which isbased on d'Alembertprinciple，isused to simulate the unconstrained rockets in flightandmakes sure the rocket has no rigid body displacement，is used w ith the linear elastic solver.The aerodynam ic coefficients distribution of the profiled rod and the elastic deformation are calculated and compared w ith program results.The steps of static aeroelastic calculation are described and the coupling procedure is explained in the follow ing section.

4.1.CFD governing equations

For CFD problems，the directnumerical simulation（DNS）is used to solve the Navier-Stokes（N-S）equations w ithout the turbulentmodelwhich needshigh speed and largememory computer due to their high nonlinearity and complexity.It is impossible to adopt this method in practical engineering. Decom posing the N-S equations into the RANS equations m akes it possible to simulate the engineering fluid dynam ic problems.

Reynolds averaged N-S equations（RANS）can be expressed as

where

whereρis the fluid density；u，v，w are the three com ponentsof Cartesian coordinate system；p is pressure；and E is total energy of unitmass.

For the closure of above equations，the shear stress transport（SST）k-ωturbulence model developed by Menter［17，18］is used.SST model combines the advantages of standard k-εand standard k-ωmodels.k andωtransport equations can bew ritten as

whereτijis shear stress

Blending function F1

where

The eddy viscosity is defined as

whereΩis the absolute value of the vorticity.F2is given by

where

The constantφof SSTmodel are calculated from the constantsφ1，φ2as follow s

where set1（φ1）is the constantof k-ωmodel，and set2（φ2）is the constant of k-εmodel.The constants of set 1（φ1）are

All other parameters are given in Ref.［18］.

Density，velocity，pressure，etc.atevery grid node could be computed by solving theaboveequations.Then，normal，axial，lift，drag，pitching moment and center of the pressure coefficients can be obtained by the follow ing equations，respectively

where q∞，α，S，xcpand xcgare the dynamic pressure，angle of attack，referencearea，center of pressure location and center of gravity location，respectively；l is the length of rod.

4.2.Static structure analysis equations

A static structural analysis determines the displacements，stresses，strains，and forces in the structures or components caused by loads that do not induce significant inertia and damping effects.Steady loading and response conditions are assumed，that is，the loads and the structure's response vary slow ly w ith respect to time.This is an isotropic linear elastic problem.The static structure analysis equation can bew ritten as

where［K］is stiffnessmatrix；｛δ｝is displacement vector；and｛F｝is force vector.

For these linear elastic structuralm echanics problems，K is a constant，and the aerodynam ic force F is calculated by CFD code.Then the deformationδcan be obtained.

4.3.Coupling boundary conditions

In the fluid-structure interaction work，there is a boundary where the interactionsof fluid and structural domainsoccur at this boundary or interface.Deformation compatibility and force equilibrium conditions should be satisfied on the fluidstructure coupling interface［19］：

where d，q，T，andτare displacement，heat flux，temperature and stress field on the fluid-structure coupling interface，respectively；n is the normal direction of interface，and subscripts f and s are the fluid and the solid，respectively.

Fig.9.Themaximum deformation of undeformed profiled rod w ith transverse vibration.

Fig.10.Profiled rod after detonation.

Fig.11.Deformation of deformed rod w ith bending vibration.

5.Discussions on the comparison of two methods

The profiled rod is shown in Fig.5.Here the initial deformation caused by detonation is not considered，the constraint conditions can be used as clamped-clamped based on the initial conditions，and the aerodynamic load and the aerodynam ic coefficient obtained by using quasi-steady aerodynam ics are 890 N and 1.16，respectively.The aeroelastic deformation can be solved by the equations for elastic vibration.Maximum deformation was obtained at the left side（where D=3 mm）of the rod at 39 mm from the origin.

In order to test the accuracy of the numerical results，we have developed a special low-cost numerical method，computational fluid dynamic finite element method（CFD/ FEM），useful for studies，based on the actualmodel.

The numerical resultsare summarized in Fig.9.As shown in Fig.9，themaximum elastic deformation of the numerical results is 0.9 mm，the maximum elastic deformation of the CFD results is0.7mm，and themaximum elastic deformation appearsat the leftside（where D=3mm）of the rod at39mm from the origin.This gave a good agreementw ith results on the undeformed profiled rod.The calculated aeroelastic deformation is substituted into Eq.（12）to solve the aerodynamics and the drag aerodynamic coefficient.The drag aerodynam ic coefficient is 1.16 that indicates the elastic deformation almost has no influence on the flight performance.

It is well known that the rods，which are explosively deployed from its warhead mechanism，generate obvious plastic deformation（Fig.10）.Therefore，the discussion on the flight attitude of profiled rod should be taken account for initial deform ation generated from previous explosion wave impacts.

Fig.11 shows the elastic deformation of deformed rod. Here，it isnoticed that themaximum elastic deformation of the undeformed profiled rod issignificantly higher than thatof the deformed profiled rod.It is considered that this difference is due to the initial deformation of the profiled rod.Namely，when the initial deformation of unreformed rod is considered，the resultsw ill be consistent.On the other hand，the relative elastic deformation is so small thatwe can neglect its impact on the flight performance.Therefore，the sem i-analytical m ethod can be used to calculate the vibration characteristics of variable cross-section beam，which solves the elastic deformation for elastic dynam ic equations of profiled rod based on the theory of Euler-Bernoulli beam model.Moreover，the rigid bodymotion can beused to describe themotion of the profiled rod flying againstair resistance athigh velocity.

6.Conclusions

A new conceptual KER warhead named profiled rod warhead was proposed firstly.The design principle of profiled rod was described and the deployment velocity of rod was obtained.

Regarding the profiled rod as Euler-Bernoulli beam model，the elastic dynamic equations of the profiled rod were established to analyze the aeroelastic effect of rod flying at high velocity，and sem i-analyticalmethod was used to calculate the aeroelastic deformation.

Besides，some efforts are applied to test the numerical resultsby using CFD/FEM.A satisfactory agreement is obtained by comparing the numerical results and CFD/FEM results.

Theproposedmodel hasbeen proven to be able to simulate the actual phenomena that the deformed rod after detonation flying at high velocity.The model can be modified for the parametric studies of intricate design change details，such as diameteror length of profiled rod，initialvelocity ofdetonation or other similar changes in the design of the profiled rod warhead.In addition，themodel built in this papermay be the foundation for the research on flight perform ance of the profiled rod deployed by detonation，and provides a theoretical basis of describing the aerodynam ic and the bomb fall characteristics.Moreover，themodelmay give great convenience to the design of the profiled rod warhead.

［1］L loyd Richard M，Sebeny JL.Physics of kinetic energy rod warheads against TBM submunition payloads.AIAA 1999：1-20.

［2］Lloyd Richard M，Sebeny JL.Novel penetrator study for defeatof ballistic m issile payloads.Int J Impact Eng 2006；33：380-9.

［3］Richard M.L loyd.Kinetic energy rod warhead w ith optimal penetrator［P］.US，6，779，462 B2，2004.

［4］Richard M.Lloyd.Warhead w ith aligned projectiles［P］.WO，03/042624 A 2，2003.

［5］Richard M.Lloyd.Kinetic energy rod warhead w ith lower deployment angles［P］.WO，05/022074 A2，2005.

［6］Qian JP，Cao B，Shen PH.Overview of ammunition.Beijing：National Defence Industry Press；2004.p.267-70［in Chinese］.

［7］Sun CJ，Lu ZH，Lu YG，Wen Y，Sun Y.Motion analysisof controllable rotation discrete rod.Explos Shock Waves 2008；28（4）：379-83［in Chinese］.

［8］Liu ZJ，Xiao C，Wang QH.The design analysis on the controllable rotation discrete rod.JExplos 2000；23（4）：48-50［in Chinese］.

［9］Cai ZB.Dynamic characteristics of high-speed flight discrete rod and penetration effect to the target.Ph.D.Thesis.Nanjing：Nanjing University of Science and Technology；2005［in Chinese］.

［10］Wang LM.An analysison the flexibility in flightof projectiles or rockets having high L/D ratios.Acta Armamentarii 2000；21（2）：109-11［in Chinese］.

［11］Cui C，Jiang H，li YH.Semi-analyticalmethod for calculating vibration characteristics of variable cross-section beam.J Vib Shock 2012；31（14）：85-8［in Chinese］.

［12］Ni ZH.Vibration mechanics.Xi'an：Xi'an Jiaotong University Press；1989［in Chinese］.

［13］Shabana AA.Resonance conditions and deformable body coordinate systems.JSound Vib 1996；192（1）：389-98.

［14］Escalona JL，Valverde J，Mayo J，Dominguez.Reference motion in deformable bodies under rigid body motion and vibration.Part II：evaluation of the coefficient of restitution for impacts.J Sound Vib 2003；264：1057-72.

［15］Lin L.A study of inertia relief analysis.In：52nd AIAA/ASME/ASCE/ AHS/ASC structures，structural dynam ics and materials conference，Denver，Colorado；2011.AIAA 2011-2002.

［16］Chen ZT，Sun Q.Applications of inertia reliefmethod in aircraft static aero-elasticity.Flight Dyn 2008；26（5）：71-4［in Chinese］.

［17］M enter FR.Zonal two equation k-ωturbulencemodels for aerodynam ic flow s.In：24th fluid dynam ics conference，Orlando，F(xiàn)lorida；1993. p.1-21.AIAA 93-2906.

［18］Menter FR.Two-equation eddy-viscosity models for engineering applications.AIAA 1994；32（8）：1598-605.

［19］Han ZZ.FLUENT-examples and analysis of flow filed calculation of engineering.Beijing：Beijing Institute of Technology Press；2009. p.13-9［in Chinese］.

.

E-mail address：ninghui85@163.com（H.J.NING）.

Peer review under responsibility of China Ordnance Society.

http：//dx.doi.org/10.1016/j.dt.2014.10.002

2214-9147/Copyright?2015，China Ordnance Society.Production and hosting by Elsevier B.V.All rights reserved.

Copyright?2015，China Ordnance Society.Production and hosting by Elsevier B.V.All rights reserved.