Advanced model of electromagnetic launcher

POWER ENGINEERING AND ELECTRICAL ENGINEERING VOLUME: 13 |NUMBER: 3 |2015 |SEPTEMBER Advanced Model of Electromagnetic Launcher Karel LEUBNER, Radim LAGA, Ivo DOLEZEL Department of Electrical Power Engineering, Faculty of Electrical Engineering, Czech Technical University, Technicka 2, 166 27 Prague, Czech Republic [email protected], [email protected], [email protected] DOI: 10.15598/aeee.v13i3.1419 Abstract. An advanced 2D model of electromagnetic launcher is presented respecting the influence of eddy currents induced in the accelerated ferromagnetic body. The time evolution of electromagnetic field in the system, corresponding forces acting on the projectile and time evolutions of its velocity and current in the field circuit are solved numerically using own application Agros2d. The results are then processed and evaluated in Wolfram Mathematica. The methodology is illustrated with an example whose results are discussed. Keywords Coil gun, electromagnetic acceleration, electromagnetic field, electromagnetic force, numerical modeling. 1. Introduction Devices working on the principle of electromagnetic launching have been known more than 150 years. We can mention, for example, electromagnetic actuators of various constructions, parts of switching gears, electromagnetic guns, and many others. But despite their relatively simple construction, their design and modeling is not an easy business because of the presence of nonlinear structural elements; some of them, moreover, move. Even when the physics of the system is known, the problem has not been solved completely so far because of strongly nonlinear interactions of involved magnetic field and two circuits (electric circuit and mechanical circuit). For example, the first attempts to design and build a magnetic coil gun were only based on more or less sophisticated experiments, which is the way able to provide a functioning device, but mostly with poor operation parameters. This part of a rather old history is well summarized and documented in [1]. A more extensive development in modeling such devices began in eighties of the last century, which was possible due to achievements in hardware and software capabilities. But still, there appeared serious problems with the numerical solution due to strong nonlinearities involved in the task. Therefore, the first models were solved in the weakly coupled formulations, i.e., the particular equations describing the system were calculated separately. This approach could be considered only orientative and the results often exhibited unacceptable errors compared with experiments. More reliable methodologies and also the first attempts to optimize selected structural parts started appearing in nineties of the last century and at the beginning of this century. Worth mentioning are several works solving the problem in the quasi-coupled formulation, which means that the mathematical model was solved iteratively as a whole, but specific nonlinearities were only respected in selected time steps (in the remaining time steps the system was considered linear). Information about these more advanced approaches can be found in [2], [3], [4], [5], [6], [7], [8] and others. The presented paper represents an attempt to solve the problem respecting all nonlinearities in the system together with the influence of eddy currents due to both transformation and motion. 2. Formulation of Technical Problem Consider the arrangement in Fig. 1. The system consists of a field coil carrying current of an appropriate time evolution, a ceramic or plastic leading barrel and a ferromagnetic body (projectile) that is to be accelerated. Geometry and physical properties of the above c 2015 ADVANCES IN ELECTRICAL AND ELECTRONIC ENGINEERING 223 POWER ENGINEERING AND ELECTRICAL ENGINEERING VOLUME: 13 |NUMBER: 3 |2015 |SEPTEMBER Fig. 1: Basic arrangement of electromagnetic launcher. parts are supposed to be known, as well as the parameters of the feeding electric circuit. The time variable magnetic field generated by the field coil produces magnetic force acting on the projectile that starts moving in the direction indicated by the arrow. The field coil must be switched off at the moment when the projectile is approximately in the middle of the coil. Then it continues moving just by its own inertia. The aim of this work is to model the movement of the projectile (particularly its acceleration and velocity) with the help of computational software Agros2D [9], [10] and Mathematica 7 [11]. 3. Mathematical Model The basic mathematical model of electromagnetic launcher consists of a nonlinear partial differential equation describing the magnetic field in the system (taking into account the eddy currents due to the velocity of motion) and two highly nonlinear ordinary differential equations for the time evolution of the field current and velocity of motion of the projectile. The corresponding equations are presented and in the following text. The electric circuit (see Fig. 2) obeys the equation Ri +d dt(Li) + 1 CZt 0 idτ=UC0,(1) Fig. 2: Feeding circuit. with the initial conditions i(0) = 0,di(0) dt=UC0 L0 . Here, Rdenotes the overall resistance of the circuit, Lis the inductance of the field coil that is a nonlinear function of the position sof the projectile, its instantaneous velocity ~v, and instantaneous value iof the field current, Cstands for the capacitance of the battery, and UC0represents its initial voltage. Finally, L0=L(~sstart,~v = 0, i = 0) is the initial inductance of the system, ~sstart denoting the starting position of the projectile. The time evolutions of velocity and trajectory of the projectile are governed by the equations md~v dt=~ Fem −~ Fdr, ~v(0) = 0, ~v =d~s dt, ~s(0) = ~sstart, (2) where mis the mass of the projectile, ~v is its velocity in, ~ Fem is the electromagnetic force acting on it and ~ Fdr denotes the sum of the drag forces (that are given by friction in the leading barrel and aerodynamic force). The quantity ~ Fem, analogously to the inductance L, is also a nonlinear function of the instantaneous position ~s of the projectile, its velocity ~v and instantaneous value iof the field current. Both the above quantities, ~ Fem and Lhave to be calculated from the current distribution of magnetic field in the system. The distribution of magnetic field may be described using several formulations based on various potentials. In our case, the magnetic vector potential ~ Ais used as the principal quantity. The full equation for this quantity may be written in the form [12] curl µ−1curl ~ A+γ ∂~ A ∂t −~v ×curl ~ A!=~ Jext,(3) where µdenotes the magnetic permeability, γis the electrical conductivity, and ~ Jext is the field current density in the field coil. The boundary condition along a sufficiently distant boundary is of the Dirichlet type. The term −γ ∂~ A ∂t −~v ×curl ~ A!represents the total eddy current densities induced in the system. The first term −γ∂~ A ∂t is eddy current densities due to transformation, the second term γ~v ×curl ~ Ais eddy current densities due to the movement. Their effects act against the magnetic field and lead to a slight deceleration of the projectile. The magnetic force ~ Fem acting on the projectile is determined using the magnetic Maxwell stress tensor c 2015 ADVANCES IN ELECTRICAL AND ELECTRONIC ENGINEERING 224 POWER ENGINEERING AND ELECTRICAL ENGINEERING VOLUME: 13 |NUMBER: 3 |2015 |SEPTEMBER Tem. The vector of the force is in our case given by the integral [13] ~ Fem =ZZN D TmdS, (4) where the integration is performed over the whole boundary of the projectile. The magnetic Maxwell stress tensor itself is given by the relation Tm=~ H⊗~ B−1 2~ H·~ BI,(5) where ~ Hand ~ Bare the vectors of magnetic field strength and magnetic flux density, respectively ( ~ B= curl ~ A,~ B=µ~ H), Istands for the unit matrix and symbol ⊗represents the dyadic product. The inductance Loccurring in Eq. (1) can be calculated from the expression L=2Wm i2,(6) where Wmis the energy of magnetic field of the system that may be calculated using the formula Wm=1 2Z V ~ J·~ AdV , (7) where Vis the volume of the definition area. As the arrangement may practically be considered axi-symmetric, both quantities ~ Aand ~ Jext have only one nonzero component in the tangential direction, while the trajectory ~s, velocity ~v and both forces ~ Fem and ~ Fdr have only one nonzero component in the xdirection (see Fig. 1). 4. Numerical Solution The task was solved numerically by means of our own application Agros2D and commercial software Mathematica 7. The code Agros2D (that is completely free) is based on a fully adaptive higher-order finite element method that is intended for numerical solution of 2D nonlinear and nonstationary multiphysics problems described by a set of partial differential equations. It is characterized by a number of quite unique features such as finite elements up to the 10th order, efficient multimesh technology, hanging nodes of any level, combination of various elements including curvilinear elements (for approximation of curvilinear boundaries and interfaces) etc. The evaluation of the results and some auxiliary computations were carried out in Mathematica 7 using a lot of own procedures and scripts. The whole algorithm consists of the following steps. 4.1. Inductances and Magnetic Forces Computation of the dependencies of the inductance L and x-component of the force Fem,x on the field current iand position xfor particular velocities ~v. In this way we will obtain two sets of nomograms (for particular velocity) from which we can easily find the values of both above quantities for values i,xand vx(and, where necessary, also their derivatives). The final values of L and Fem,x are for any current, position of the projectile and its velocity calculated using linear interpolation in the above nomograms. 4.2. Circuit Equations Equation (1) and Eq. (2) are solved simultaneously but successively, in time steps. •Solution of Eq. (1). First it is necessary to modify the term d(Li)/dt. d dt(Li) = Ldi dt+idL dt=Ldi dt+idL dx dx dt=(8) =Ldi dt+i~v dL dx. For further processing, we determine the second derivative of term Li with respect to time: d2 dt2(Li) = Ld2i dt2+di dt dL dx dx dt+di dt dL dx dx dt+ +idL dx d2x dt2+i~v d2L dxdt=Ld2i dt2+ 2~v di dt dL dx+ +id~v dx dL dx+i~v d2L dxdt. (9) After another differentiating Eq. (1) with respect to time tand substituting for d2(Li)/dt2we obtain Rdi dt+Ld2i dt2+ 2~v di dt dL dx+id~v dx dL dx+i1 C= 0.(10) Introducing the differences instead of the derivatives and substitution g=di/dtleads to a system of difference equations in the form: gk+1 =gk1−R∆tk Lk −2~vk LkdL dxk−hk,(11) hk=1 Lkikd~vx dtkdL dxk −ik 1 C∆tk,(12) and ik+1 =gk∆tk+ik,(13) c 2015 ADVANCES IN ELECTRICAL AND ELECTRONIC ENGINEERING 225 POWER ENGINEERING AND ELECTRICAL ENGINEERING VOLUME: 13 |NUMBER: 3 |2015 |SEPTEMBER with the initial conditions i0= 0, g0=UC0 L0 ,(14) Here, adenotes the acceleration, kis the number of step and ∆tkstands for the k-th time step. •Solution to Eq. (2). No modification of this equation is necessary. •The values L,dL dxand Fem,x can easily be determined from the above nomograms using an appropriate interpolation or extrapolation technique. On the other hand, the nomograms must be appropriately dense. 4.3. Numerical Problems Connected with Solution The principal complication in the numerical solution is computation of the nomograms for higher values of the field current. High field currents bring about either partial or complete oversaturation of the ferromagnetic projectile. But this oversaturation may be evaluated only approximately, because the magnetization curve of material in the relevant domain can only be estimated. This leads to non-estimable errors of results, and, moreover, the related iterative processes are often accompanied by various undesirable phenomena such as oscillations, and generally, a very poor convergence rate, which means a long time of computations. 5. Illustrative Example The mathematical model was tested on a portable three-stage electromagnetic accelerator depicted in Fig. 3 that was completely developed and built by the second author. Both the model and measurements were realized just on the first stage of the device. The principal aim of the task was to compare the measured and calculated velocity and trajectory of the projectile. 5.1. Input Data The principal dimensions of the coil of the first stage, leading barrel and projectile are shown in Fig. 4. Their values are: D1= 36 mm, d1= 8 mm, L1= 50 mm, D2= 6.75 mm, L2= 52 mm and initial position of the end of the projectile xstart =−38 mm. The voltage UC0of the capacitor battery is 350 V. The parameters of the electric circuit were measured for frequency 156 Hz, because the shape of the corresponding current corresponds to the shape of real current in the series RLC electric circuit. For the above frequency, the capacitance C= 7.11 mF, the total resistance of the circuit R= 0.145 Ω and inductance of the coil L= 0.220 mH. The coil contains 203 turns wound in 7 layers, each having 29 turns. The conductor is of circular cross section and its diameter is 1.6 mm. The leading barrel is made of brass with a longitudinal gap, so that it has no effect on electromagnetic quantities. The projectile is made of material Vacoflux 48 (produced by German company Vacuumschmelze GmbH & Co). After manufacturing, the surface of the projectile was hardened, but the conditions of hardening are not known. That is why we used the curve µr=µr(B)of its primary magnetization that is depicted in Fig. 5. Electric conductivity of this material is approximately 2.5 MS·m−1and its mass density is 8120 kg·m−3. 5.2. Results For the accuracy of results the decisive role is played by the correctness of the time evolution of current i in the field coil. This current was measured and also modelled using Eq. (1). The results are presented in Fig. 6. The blue line (measured) shows the measured current that was switched off at the moment when the projectile reached the middle of the coil. The black line (calculated 1) shows the evolution for circuit parameters given in subsection 5.1. (here, however, the current was not switched off). And just for comparison, the red line (calculated 2) shows another calculated evolution for circuit parameters determined by their measurement at the industrial frequency, i.e., not quite correctly. Obviously, the accordance of the blue and black lines is extremely good. Important are the nomograms of the inductance L, its derivative with respect to coordinate x, and magFig. 3: A three-stage electromagnetic accelerator. c 2015 ADVANCES IN ELECTRICAL AND ELECTRONIC ENGINEERING 226 POWER ENGINEERING AND ELECTRICAL ENGINEERING VOLUME: 13 |NUMBER: 3 |2015 |SEPTEMBER Fig. 4: Principal dimensions of the first-stage coil and projectile. netic force Fem,x, that are functions of the field current i, and position xand velocity ~v of the projectile. Figure 7 contains the corresponding nomogram of Lfor velocity ~v = 0 m·s−1. In fact, the figure contains two nomograms. The transparent one shows the case when the relative permeability of the projectile corresponds to Fig. 5, but ends at the value µr= 4.6for B= 1.5T (for lower values of µrthere appeared severe problems with the convergence). The real nomogram, however, is the lower one. For its construction, a special restricting function had to be proposed, being able to simulate even lower values of the relative permeability. In case of higher values of the field current i, the projectile becomes strongly oversaturated and its relative permeability drops to 1. The analogous nomograms for the derivative of inductance Land magnetic force Fem,x are depicted in Fig. 8 and Fig. 9. The nomogram containing the derivative dL/dxis not as smooth as it should be. The visible oscillations are caused by numerical errors due to problems with the convergence of the iterative processes in the domain of a high saturation of the projectile. A better level of smoothness could be achieved by using interpoFig. 5: Magnetization curve of material Vacoflux 48. Fig. 6: Time evolution of current in the field coil. Fig. 7: Inductance of the coil. Fig. 8: Derivative of the real inductance with respect to x. lation polynomials of higher orders, but at a cost of a substantially longer time of computation. On the other hand, many tests proved that these oscillations do not affect the resultant evolution of the field current too significantly. The time evolution of motion of the projectile was calculated using Eq. (2), where the drag force is represented by the aerodynamic resistance. Its value is given by the formula Fdr,x =1 2CρSv2,(15) where ρis the density of air, Sis the area of the cross section of the projectile and Cis a constant that depends on the shape of the projectile. In our case c 2015 ADVANCES IN ELECTRICAL AND ELECTRONIC ENGINEERING 227 POWER ENGINEERING AND ELECTRICAL ENGINEERING VOLUME: 13 |NUMBER: 3 |2015 |SEPTEMBER Fig. 9: Force acting on the projectile in the direction of x: transparent nomogram - constant permeability of the projectile, low (blue) nomogram - real permeability. C= 0.4, but some variations of this coefficient proved to have only a small influence on the time evolution of velocity. Figure 10 depicts the time evolution of the position of the rear face of the projectile and Fig. 11 shows the time evolution of its velocity. The time step ∆t= 10−6s. The red lines are the dependences without respecting eddy currents in the projectile, the blue ones represent them. The differences are, however, very Fig. 10: Time dependence of trajectory of the projectile. Fig. 11: Time dependence of velocity of the projectile. small. The eddy currents decelerate the velocity of the motion, but only very slightly. The average muzzle velocity measured on the device was 28.9 m·s−1. The difference from the calculated values is about 1 m·s−1(over 3 %), which can be considered quite excellent. 6. Conclusion The paper is aimed at numerical modeling of operation characteristics and parameters of an electromagnetic accelerator. The proposed mathematical model consisting of one partial and two ordinary strongly nonlinear differential equations was solved numerically in the quasi-coupled formulation, with pre-calculated nomograms of the inductance of the system, its derivative and magnetic force acting on the projectile. The computations had to be carried out extremely carefully because of poor convergence due to oversaturation of projectile at high values of the field current. The model proved to provide realistic results whose agreement with the measured data is excellent. Next work will be focused on modeling all three stages of the electromagnetic gun. Acknowledgment The financial support of the Grant project GACR P102/11/0498 is gratefully acknowledged. References [1] MCNAB, I. R. Early electric gun research. IEEE Transactions on Magnetics. 1999, vol. 35, iss. 1, pp. 250–261. ISSN 0018-9464. DOI: 10.1109/20.738413. [2] BENGUI, Z., C. YANJIE, W. JIE, W. HUIJIN and C. XUEHUI. Magnetic-Structural Coupling Analysis of Armature in Induction Coilgun. IEEE Transactions on Plasma Science. 2011, vol. 39, iss. 1, pp. 65–70. ISSN 0093-3813. DOI: 10.1109/TPS.2010.2077742. [3] LIU, W., C. YANJIE, Y. ZHANG, J. WANG and D. YANG. Parameters Optimization of Synchronous Induction Coilgun Based on Ant Colony Algorithm. IEEE Transactions on Plasma Science. 2010, vol. 39, iss. 1, pp. 100–104. ISSN 00933813. DOI: 10.1109/TPS.2010.2076315. c 2015 ADVANCES IN ELECTRICAL AND ELECTRONIC ENGINEERING 228 POWER ENGINEERING AND ELECTRICAL ENGINEERING VOLUME: 13 |NUMBER: 3 |2015 |SEPTEMBER [4] AKSOY, S., A. BALIKCI, Z. ZABAR and L. BIRENBAUM. Numerical Investigation of the Effect of a Longitudinally Layered Armature on Coilgun Performance. IEEE Transactions on Plasma Science. 2010, vol. 39, iss. 1, pp. 5–8. ISSN 0093-3813. DOI: 10.1109/TPS.2010.2055585. [5] HAGHMARAM, R. and A. SHOULAIE. Study of Traveling Wave Tubular Linear Induction Motors. In: International Conference on Power System Technology. Singapore: IEEE, 2004, pp. 288–293. ISBN 0-7803-86108. DOI: 10.1109/ICPST.2004.1460008. [6] TZENG, J. T. and E. M. SCHMIDT. Comparison of Electromagnetic and Conventional Guns from a Mechanics and Material Aspect. In: 23rd International Conference on Ballistics. Tarragona: Academia, 2007, pp. 597–604. ISBN 978-84-7493379-6. [7] VOKOUN, D., M. BELEGGIA and L. HELLER. Magnetic Guns with Cylindrical Permanent Magnets. Journal of Magnetism and Magnetic Materials. 2012, vol. 324, iss. 9, pp. 1715–1719. ISSN 0304-8853. DOI: 10.1016/j.jmmm.2011.12.026. [8] GUO, L., N. GUO, S. WANG, J. QIU, J. G. ZHU, Y. GUO and Y. WANG. Optimization for Capacitor-Driven Coilgun Based on Equivalent Circuit Model and Genetic Algorithm. In: IEEE Energy Conversion Congress and Exposition. San Jose: IEEE, 2009, pp. 234–239. ISBN 978-1-42442893-9. DOI: 10.1109/ECCE.2009.5316199. [9] KARBAN, P., F. MACH, P. KUS, D. PANEK and I. DOLEZEL. Numerical Solution of Coupled Problems Using Code Agros2D. Computing. 2013, vol. 95, iss. 1, pp. 381–408. ISSN 1436-5507. DOI: 10.1007/s00607-013-0294-4. [10] Agros2D. Agros2D [online]. 2015. Available at: www.agros2d.org. [11] Wolfram Mathematica: The world’s definitive system for modern technical computing. Wolfram [online]. 2015. Available at: http://www. wolfram.com/mathematica. [12] KUCZMANN, M. and A. IVANYI. The finite element method in magnetics. Budapest: Akademiai Kiado, 2008. ISBN 978-963-0586-498. [13] STRATTON, J. A. and A. IVANYI. Electromagnetic theory. Hoboken: John Wiley, 2007. ISBN 978-047-0131-534. About Authors Karel LEUBNER was born in 1982. He received his M.Sc. from the Czech Technical University in Prague in 2006. Now he is finishing his doctoral studies. His research interests include computational electromagnetics and modeling of multiphysics problems. Author or co-author of about 10 papers and conference presentations. Radim LAGA was born in 1990. He is a student of the second course of magister studies at the Czech Technical University in Prague. His research interests include electrical engineering, electronics and experimental techniques. Ivo DOLEZEL was born in 1949. He received his M.Sc. from the Czech Technical University in Prague in 1973. His research interests are in the mathematical and computer modeling of electromagnetic fields and coupled problems in heavy current and power applications. Author or co-author of two monographs (USA), more than 400 research papers, and several large codes. c 2015 ADVANCES IN ELECTRICAL AND ELECTRONIC ENGINEERING 229