Yuguo CHENG, Gungqing XIA

a PLA 91550 Element 41, Dalian 116023, China

b State Key Laboratory of Structural Analysis for Industrial Equipment, Dalian University of Technology, Dalian 116024, China

c Key Laboratory of Advanced Technology for Aerospace Vehicles of Liaoning Province, Dalian University of Technology,Dalian 116024, China

KEYWORDS

Abstract The pulsed inductive thruster is characterized of no electrode corruption and wide propellant choice.To give insight into the propulsion mechanism of small scale thruster at different propellant mass (m) and energy (E) levels, the transient MagnetoHydroDynamics (MHD) method,completed by high temperature thermodynamic and transport,and plasma electrical models,is developed to study argon plasma response under the excitation of current of high rise rate.By calculating the two-dimensional expansion properties of the thruster with conical pylon,the simulations find that the main energy deposition occurs during the initial pulse rise stage,and the energy density of Joule heat is two magnitudes higher than the deposition in the down side.At propellant mass of 2 mg,average axial velocity of the current sheet increases from about 15 km/s at 750 J to about 21 km/s at 1470 J within the decoupling distance.The velocity variation synchronizes with the pulsed rise in the initial.The monotonically decrease of the temperature along axis results in the growth of low ionization level ions and reducing of high levels.The current sheet maintains the structure formed during the initial pulse rise when moving beyond the decoupling distance.Besides the change in forward velocity,the main difference is the dimension compared with that in the first half period,caused by thermal conduction and particle diffusion. The variations of total impulse It in the range of m from 2 mg to 8 mg and E from 750 J to 1470 J show that It is proportional towhen E is determined.

1. Introduction

The electric thrusters are characterized of high specific impulse and good propellant compatibility,which have substituted the chemical thrusters in orbit keeping, attitude control missions of earth satellites and as the main engines in deep space exploration gradually.1,2Electrostatic thrusters, mainly the ion thruster, can reach an efficiency of 60%,2and based on the evaluation of thrusters already in use,the lifetimes are competent for several years. Although corruption of the electrode exists, the main problem affecting the performance lies in the space charge limit. To avoid this, a promising choice is the pulsed inductive propulsion.3

The planar Pulsed Inductive Thruster (PIT) ionizes and accelerates the neutral particles by transient electromagnetic field,generated by the pulsed current flowing through the specially designed planar coil, and the thrust is produced by the interaction between the azimuthal current in the plasma and the coil,4,5and the experiments found that the energy deposited within the decoupling distance determined the performance. The decoupling distance is a measurement of the interaction between the excitation coil and plasma, which means that the main electromagnetic energy deposits before axial distance of the bulk plasma reaches this value, then the magnetic effect on the plasma decreases quickly when moves far away.6The value increases with greater coil dimension7and once the thruster configuration and materials are established, the factors affecting the performance are the pulsed energy, initial condition and the kind of gas.8

Polzin3generalized from the experiments that a high current rise rate and a uniform gas density over the coil face were needed to improve the performance of PIT, when a uniform current sheet was produced. Dailey9studied the performance of the thruster at static fill and non-zero velocity case in the simulation and found that the static fill case got better performance.

Mikellides et al.10,11numerically studied the impulses and thruster efficiencies of Mark V prototype at energies of 900–1764 J, and propellant mass of about 1.5–10 mg. The simulations employed the MACH2 MHD model,12and the thermodynamics and transport properties employed the SESAME library.13The plasma electrical model was also developed to obtain consistent results,14in which the external circuit was assumed to be of constant resistance and inductance. Monatomic gas,such as Ar and He,and polyatomic gas as NH3were the propellants investigated.10,14The simulations confirmed the constant-efficiency operation of the thruster and the critical mass phenomenon,10and revealed that elevated efficiency of NH3relative to Ar and He owing to the less radiation,suggesting different propellants have similar impulse generation mechanism,11and the main difference is the energy deposited.Due to the similarities, a simplified impulse expression based on one-dimensional model is obtained.10These findings make the in-situ propellant filling possible for high performance can also be realized by using other gas.15

Although MACH2 simulations reproduce the experiment results for a range of energy levels and propellant mass values,several improvements can be done to give more insight into the plasma behavior of PIT: the first is the revise of spatial deviation scheme,so that high resolution flow field can be obtained;the second is the development of thermodynamics and transport models, for the solution of density of individual species is not provided by SESAME library; the last is the twodimensional expansion characteristics of the plasma, in the first half period of the pulsed current and thereafter, which are less reported.

In this paper, the transient MHD method, including the high temperature thermodynamic and transport models, is developed. Numerical simulations, in which the twodimensional characteristics of flow field are analyzed in detail,are performed on a smaller scale than Mark V,with the aim to penetrate into the operation characteristics of smaller thruster and to make use of PITs in attitude control,orbit keeping missions,and so on,not only as the main engine. The discussions in this work are expected to be applicable to general cases.

2. Model description

2.1. Governing equation

The current sheet forms during the initial pulse rise,which is a high density region moving forward,the collision frequency of the particles in the sheet are high enough so that temperature difference between the electrons and heavy particles is within 1%, as calculated in Ref.11. In this work, a thermodynamic equilibrium MHD model are employed to depict the flow of transient inductive plasma, and the governing equation are revised by the hyperbolic cleaning algebra equation to maintain the conservation of magnetic divergence in the calculation,16the continuum, momentum, energy, magnetic inductive and hyperbolic cleaning equations are as follows:

whereL(T)is the volume radiation loss rate,curve fitted using the model in Ref.18.

2.2. Thermodynamic and transport properties

2.2.1. Thermodynamic property

Under the excitation of pulsed current, temperature of the plasma can be greater than 2.0 eV (1 eV=11605 K) at low propellant mass or high pulsed energy case, and ionization up to the forth level is considered to include the main species that may be generated. By solving Eqs. (1)–(5), the internal energy and density are obtained directly, which is related to the number density and temperature by Eq. (7):

In Fig.1, the following equations are solved together to obtain the analytical expression off(p,ρ) at different pressure and propellant mass levels.

The first is the perfect gas law:

The second is the charge conservation law:

and the last are the Saha equations:

In this work,argon is used as the working gas.The number density (n) and specific heat ratio (γ) at different temperature and pressure are shown in Fig.2.Comparison of Fig.2(a)with the result in Ref.19shows good consistency,demonstrating the validity of thermodynamic calculation. And in Fig.2(b) γ decreases sharply with increasing temperature,indicating more energy will deposit with less kinetic energy increment.

2.2.2. Transport property

The electric conductivity is one of the most important transport parameters affecting the plasma discharge,especially during the initial pulse rise, when the number density varies rapidly. In order to obtain the electric conductivity at varied temperature and pressure, the first order approximation solution of the Boltzmann equation is adopted,and the electrical conductivity is expressed as20:

Fig.1 Numerical procedure calculating p, T and number densities of species.

Fig.2 Variation of the number densities at 104 Pa and specific heat ratios at different pressure with temperature of argon plasma.

whereNis the number of charged species,mArj+is particle mass of Arj+,Dijis the multi-component diffusion coefficient calculated using the method in Ref.20, and the Sonine polynomial ofDijis solved in terms of the collision integral Ω （l，s）ij:

the deflection angle χ in Eq. (13) is:

To resolve the singular and orbiting problems occurring in calculating the collision integrals, a convergence algorithm is developed, the comparison of the deflection angle and cross section with the solutions in Refs.21,22, where orbiting phenomenon of Morse potential for different energies and piecewise characteristics of the cross section of Lennard-Jones (L-J) potential are accurately investigated, are given in Fig.3, showing well agreement.

The following potential are taken into account for different collisions:

(1) The HFDTCS2 potential is employed for Ar-Ar collision23.

(2) The elastic and charge exchange collisions are considered for Ar+-Ar collision.

(3) The values in Ref.20are used to evaluate the e--Ar cross section.

(4) The screened coulomb potential is adopted for collision between charged particles.

Variations of diffusion coefficient of Ar+-e-and electrical conductivity in thermal equilibrium state with temperature are shown in Fig.4. The curves have the same form as those in Ref.19, which verify correctness of the algorithm.

2.3. Calculation configuration and boundary conditions

2.3.1. Calculation configuration

Fig.3 Verification of algorithm of collision integral.

Fig.4 Diffusion coefficient DAr +- e- and σe as T is varied at 1 atm.

Fig.5 Calculation configuration of pulsed inductive plasma discharge.

The two-dimensional axial symmetric calculation configuration inrOzcoordinate system is shown in Fig.5, in which the thruster is comprised of the coil,conical pylon and confining cuff to emulate the main physical structure of experiment prototype.

The calculation geometry, as is noted in Section 2.1 by domainD,is enclosed by polygonABCDEFGin Fig.5,whereABrepresents the planar excitation coil transmitting the capacitor energy and is separated from the plasma by nonconducting transparent medium;BCis the confining cuff and is considered as non-conducting wall boundary;CD,DEare the virtual outlet boundaries;EFis the symmetric axis;FGAis the conical pylon and is considered as non-conducting wall boundary.

The initial static filling assumption is implemented in the simulation, and the gas is uniformly distributed within the axial length of segmentBC, meanwhile, the background is assumed to be of a density 1/10 of the surface. Low density assumption in the near field of the background may be proximate in practical situation due to the diffusion of denser gas to the vacuum.

Spatial deviation employs the M-AUSMPW+ scheme and the third TVD-RK is adopted in the time advancement.

2.3.2. Boundary conditions

For the aerodynamic conditions, no slip assumption is implemented in wall boundaries and axial symmetric condition is applied in the symmetric axis; for the magnetic conditions,the axial and radial magnetic field intensities at nonconducting walls are in accordance with those introduced in Ref.12, radial magnetic intensity at the coil-plasma interface is calculated according to the Ampere’s law:

In Eq. (15), the coil currentJcis determined by the plasma electrical model.The circuit model is shown in Fig.6,in which the plasma is considered as a load of resistanceRpand inductanceLpand the circuit a load of resistanceReand inductanceLe.

The voltage conservation equation of the circuit in Fig.6 is10:

Fig.6 Circuit model of discharge.

2.4. Model demonstration

To verify the ability of the code in resolving strong discontinuity of the MHD flow and the model’s accuracy in simulating the transient inductive discharge, the following two cases are taken to test the calculation results.

Case 1. Orszag-Tang vortex problem24

The flow is developed within a square region inxoycoordinate, as shown in Fig.7. The flow parameters are dimensionless in this case and the non-dimensional initial condition is:

whereu,vis the velocity inx,ydirection, respectively,Bx,Byare the magnetic field intensity inx,ydirections, respectively,the specific heat ratio γ=5/3. Periodic boundary condition is implemented in the flow domainx,y∈[0,2π].The initial sinusoidal waves are discontinuous as time advances, the density contour lines are at different time are shown in Fig.7.

From Fig.7, the flow becomes complex due to the nonlinear interference among waves, and the flow turns into MHD turbulence in the end. The ‘thick solid lines’ magnetic shock waves generated from the non-contact discontinuities are distinguished by the code and the structure is in accordance with the Galerkin algorithm25.

Case 2. Impulse of Ar of PIT atE=1764 J10

In this case,the impulse of Ar at pulsed energyE=1764 J is calculated and current of the form in Ref.10is used, and the inner and outer radius are 0.2 m,0.5 m,respectively,geometric parameters are in accordance with the experiment prototype5.

Fig.8 gives the comparison of calculated impulses with the experiment at varied propellant mass, where the value gap reduces with increased mass is distinct whenm＜6 mg and the relative errors are about 4.1%, 3.1%, 1.3%, 2.5% atm≈3, 4, 5.6 and 9.2 mg, respectively, which are acceptable in predicting the performance of thruster. Since impulse is an integral value, the result also indicates correctness of predicting critical flow parameters.

Two factors may lead to the disparity between the results:the first is the initialization condition; the second is the background settings. In this work, the initial temperature is set to be 0.1 eV, different from the 0.025 eV in MACH2, which avoids the initial breakdown process and may result in more energy deposition, and the background is assumed to be vacuum in MACH2 simulation,slightly different from the settings here.

The above two cases demonstrate reliability of the model in simulating the operation characteristics of the mass range in Fig.8 and medium pulsed energy.

3. Results and discussion

In this section, geometric parameters of simulated PIT are:ri=0.1 m(LOA),ro=0.3 m(LOB),LBC=0.15 m,LOF=0.2 m,LGF=0.05 m,LOE=0.6 m. The propellant mass and pulsed energy of the following ranges are calculated:

(1) Propellant mass from 2 to 8 mg;

(2) The capacitor parameters are:C=15 μF,V0=10–14 kV, with the corresponding pulsed energy ranges from 750 J to 1470 J.

The above ranges enable the excitation of high ionization levels in the low mass and high energy case, for increasing the energy absorbed per particle, and the thermal and transport properties will deviate from quasilinear variations at higher temperature, as shown in Figs. 2 and 4.

The form of excitation current and temporal evolution of the transient electromagnetic force and impulse at propellant mass (m=4 mg) and pulsed energy (V0=13 kV,E=1267.5 J) are shown in Fig.9. The axial transient forceFis calculated by volume integration ofjθBr(jθ= ?×B)and thermal expansion velocity, and the total impulse is the time integration of the transient force.

Fig.7 Flow field of Orszag-Tang problem at different non-dimensional time t-.

Fig.8 Comparison of calculated total impulse with experiment at different propellant mass of Mark V prototype.

The evolution in Fig.9(b) gives insight into the operation process, it can be seen that peak transient force appears in the first half period,during which the values rise sharply,indicating large azimuthal current is generated and main ionization occurs. The plasma is expelled by the repulsion force and the interaction within the decoupling distance determines peak value of the transient force. The main excitation current pulse lasts about 14 μs,shown in Fig.9(a),with the maximum radial magnetic intensity of 0.49 T and appears at 1.6 μs; nevertheless, the time needed for getting stable total impulse is about 100 μs, showing that self-induced current in the sheet plays an important role in the expansion after the first half period.

3.1. Two-dimensional expansion characteristics of flow

Similar evolution characteristics of the fields in the mass and energy range mentioned before are obtained, and in this section, flow field ofm=4 mg andE=1267.5 J (V0=13 kV)are presented to analyze the flow properties at different stage.

As can be inferred from Fig.9(b) and also generated from the simulations within the mass and energy range stated, two different processes occur sequentially based on the interaction of the coil current and induced azimuthal current in the plasma: the first is the ionization and acceleration within the decoupling distance mainly in the pulse rise of first half cycle;the second is the acceleration in the second half cycle and thereafter. The former is dominated by the intense transient electromagnetic action,during which the species of the mixture and intensity of the current sheet and coil current vary rapidly;the latter is controlled by weaker electromagnetic action, and thermal conduction and particle diffusion take effect as time advances, meanwhile, the high density sheet maintains similar structure to that in the first half period when moving forward.

3.1.1. Flow evolution in the first half period

The two-dimensional density evolution in the first half period is shown in Fig.10,where a sharp high density front is apparent,and the neutral gas is ionized by the transient electric field generated by the high rise pulsed current, forming the current sheet. Peak density within the sheet, larger than 2×10-4kg/m3in the front, is one order magnitude higher than residual plasma close to the coil surface, indicating efficient deposition of the pulsed energy in the first half period.The peak density rises from about 2.2×10-4kg/m3at 2 μs to about 3.5×10-4kg/m3at 5 μs, compared with 6×10-5kg/m3in the initial to 2.4×10-4kg/m3at the rise side, the result demonstrates slower increment of denstiy in the downside of the curve in Fig.9 (a). The synchronous change of the pulsed current and density shows the time dimension of the field change is determined by the current during this period, which reflects the importance of high current rise rate in the initial pulse rise.

The inductive magnetic field intensities at 2 μs and 5 μs are shown in Fig.11. Maximum intensity is about 0.46 T at 2 μs,with uniform radial and monotonically decreasing axial distribution. When the peak of current is reached, the field close to the coil surface is disturbed by reversing of change rate of the excitation field. Besides changing of the values, the distinction is the increased nonlinearity at 5 μs, both axially and radially.The intensities along axial direction at different radial locations of the two moments are shown in Fig.12,which suggests exponentially decreasing at 2 μs and the curve can be fitted by:

Fig.9 Excitation magnetic field intensity and performance of thruster.

Fig.10 Distributions of density in the first half period.

Fig.11 Inductive magnetic field intensities in the first half period.

Fig.12 Radial magnetic intensities along axis.

wherea1=0.57,a2=-34.71 anda3=-2.6×10-3at 2 μs in this case.The trend is in accordance with the one-dimensional model concerning radial magnetic field profile proposed in Ref.10.The result shows that the deposited energy mainly converts into the axial directed kinetic energy and radial diffusion is insignificant in the current rise side,as seen in Fig.13(a).After peak value is reached, changing rate of the excitation filed is reversed,meanwhile,self-induced magnetic intensity in the high density region will exceed the region close to boundaryAB,leading to divergence flow on the coil surface, as shown in Fig.13(b)and magnetic profiles deviating from the exponential trend at 5 μs, with local maximums exist within the sheet, displayed in Fig.12(b).

The number densities of the mixtures at 5 μs are shown in Fig.14,where a narrow high number density region is obvious.Maximum density ofnAr3+is about 4.3×1021m-3, which exceed those ofnAr+,nAr2+,nAr4+, showing high temperature of the region, as illustrated in Fig.15. As can be seen from Fig.14, spatial separation of heavy particles is apparent. The higher ionization levels species locate on the back of the sheet and the lower levels on the front, indicating axial decreasing characteristics of the temperature. Under the assumption of electric neutrality in the model, it can be concluded that ambipolar diffusion exists in the movement of inductive plasma.

Fig.13 Two-dimension distributions of axial velocities Vz.

Fig.14 Number densities of species of different ionization levels at 5 μs.

Fig.15 Temperature distributions at 5 μs.

The Joule heating as time is varied at points (0.025 m,0.1 m),(0.025 m,0.2 m)and(0.025 m,0.3 m)near the coil surface is plotted in Fig.16,calculated by ｜J｜2/σe.The main energy deposition is completed before 2.5 μs,and peak energy density is about 1010W/m3, and the calculated thermal conductivity?qis about 106–108W/m3. The magnitude comparison shows that Joule heating dominates the ionization and acceleration in the rise side of the excitation current,thermal conduction takes effect when the electromagnetic heating is weakened as time advances, as discussed in Section 3.1.2.

3.1.2. Flow evolution after the first half period

Fig.16 Energy distribution as time is varied.

In the second half period and afterwards, the excitation intensity is nearly one order magnitude lower than before;furthermore, abscissa value of the sheet is greater than the decoupling distance, implying the weakening of electromagnetic effect. Two factors affect the thermal movement thereafter: the first is expansion of the sheet to larger space when crossing the corner of conical pylon; the second is the interaction between the coil and self-induced current in the plasma sheet in a low level.

The density and velocity evolutions at moments 10 μs,20 μs and 35 μs are shown in Fig.17. Density in the sheet varies slowly, with the maximum values around 1.8×10-4kg/m3in central radial region, lower than that at the end of first half period. Comparing the density image at 35 μs with previous time, the high density region shifts to axis when moving forward, and the value, larger than 2.4×10-4kg/m3, is greater than those of the central region. The convergence can be explained as follows: the expanding configuration of conical shape guides the flow to the center axis when crossing the turning point. The reason results in the high density effect toward the axis. The convergence characteristic is beneficial to the improvement of the performance,for divergence plume usually means greater radial kinetic loss, which leads to lower energy utilization efficiency. When crossing the turning point, the plasma diffuses immediately downward to the low density background, the abrupt turn of the domain leads to large radial velocity, which is visible in the streamlines of Fig.18,forming a low density region that is nearly circular shaped at outer contour.

Judging from the values of parameters in Fig.17,impact of the turning on the structure of the sheet is small, and velocity in the sheet decrease monotonically from about 1.6×104m/s at 10 μs to about 9×103m/s at 35 μs at the core. Different from the density, the velocity reaches local maximum around the corner, which can be two or more times the values of the sheet, due to the diffusion of low density plasma.

Fig.17 Two-dimensional distributions of density and axial velocity of m=4 mg and E=1267.5 J.

Fig.18 Streamlines of flow at 15 μs when crossing the turning of conical pylon.

To analyze the effect of reversed flow, as can been seen in the velocity distributions at 10 μs (Fig.17(a)) and 15 μs(Fig.18), which may emerge due to the diffusion of high density current sheet to the low density region (ρ ＜0.2ρmax) on the coil surface, volume mass integration∫ρdVof the two regions and the ratio of the high density region to total mass are calculated,the results are listed in Table 1,and two magnitude order larger of the former (high) than the latter (low),when the flow is fully developed (t≥20 μs), is obtained. The comparison of mass shows well confinement of the plasma within the sheet formed in the first half period when moving forward. It is noteworthy that, mass in the sheet may be over evaluated here by underestimating random diffusion of the plasma into the extreme low density background in practice,nevertheless, the transient thrust and total impulse are little effected for the electromagnetic force is dominant.

The inductive field intensities at 20 μs and 35 μs are shown in Fig.19, different from the images in Fig.11, where the values are positive in the main stream,the nonlinearity at the two moments is strengthened, especially around the turning point of conical pylon, where opposite signs of radial magnetic field exist simultaneously.In Eq.(4),the inductive intensity is determined by the electrical conductivity and convection velocity of the flow, calculation results in Fig.17 have shown that large velocity gradient exists around the turning, thus the magnetic field shown.

At 20 μs,the radial intensity is on the magnitude of 10-2T and even lower in the far field.Sincejθ∝Br,the Lorentz force is two or three orders of magnitude lower that in the first half period forBrof the range 4×10-3–3.5×10-2T, as can be seen in Fig.9, where the peak transient force is 56 kN, and is 0.8 kN at 20 μs. The contrast demonstrates the way by changing amplitude of the excitation current to adjust the thrust, and subsequently the impulse, as required, while improvement of propulsion efficiency by increasing the current will be uncertain compared with the case here, for energy utilization in the second half period is low caused by the movement beyond decoupling distance of the current sheet. The intensity within the sheet maintains the positive values and with relatively regular distribution during the period investi-gated,meanwhile,high value region shifts to low ordinate side,as the main stream moves toward the axis.

Table 1 Mass comparison of two regions at different time.

The number densities of the mixture at 15 μs are shown in Fig.20,where the higher level ions dominate the flow,and the lower in the narrow front.Maximum density is 2.5×1021m-3of Ar2+, indicating decreasing of the temperature with time.Compared with the densities in Fig.14, the main differences are the location,magnitude and affected area,relative position of different levels are not changed,and the slowly varied magnitudes and enlarged affected area in the post half period shows that thermal conduction is prominent in the flow,whose characteristic time is longer than the pulse period. The profile of the high and residual density interface, which is distinct in Fig.20, is related to the angle of conical pylon.

In general, the flow characteristics within the current sheet is maintained in the post first half period,and the nonlinearity mainly exists around the corner of conical pylon, where the abrupt turn leads to large velocity gradient and unregularly distributed magnetic inductive field. Compared with the magnitudes of thrust in the first half period,the effect of the conical pylon on the performance can be neglected.

3.2. Flow characteristics and performances at different propellant mass and energy levels

3.2.1.Comparison of density and inductive magnetic intensity at different energy levels

In this section, the two-dimensional expansion characteristics,mainly the mass density and magnetic field at different propellant mass and energies, are compared. Images att=50 μs under propellant mass of 8 mg, energies of 750, 1080, 1470 J are given in Fig.21 and Fig.22,The excitation magnetic intensity approaches zero at 50 μs, and the flow is fully developed.For the first ionization level of argon, the collision mainly comprised of ionization, excitation and elastic processes, and accurate reaction rate constants of the processes are available.The collision energy loss εcper electron-ion created is calculated by26

As calculated in Ref.26, the collision energy loss is about 100 eV at electron temperature of 2 eV for Ar plasma. AtE=750 J and taking into consideration the gas close to the surface, e.g. a third of the propellant mass,m1/3≈0.7 mg form=2 mg andm1/3≈2.7 mg form=8 mg, the average energies deposited per neutral particle are about 465 eV,116 eV,which can assure first level ionization of the surface gas at the lowest energy and mass ratio (E/m) considered. Besides affecting the ionization level, the propellant mass influences the skin depth,which is closely related to the electrical conductivity, and consequently the kinetic energy of the flow.

Fig.19 Inductive magnetic field intensities after the first half period.

Fig.20 Number densities of species of different ionization levels at 15 μs.

Fig.21 Density distributions at different pulsed energies of m=8 mg and t=50 μs.

Fig.22 Magnetic field distributions at different pulsed energies of m=8 mg and t=50 μs.

Fig.21 gives the densities at different energies ofm=8 mg.The images show different magnitudes and moving distances caused by varied energies, and the greater the energy, the larger the abscissa value of the sheet. In these cases, the fields have similar distribution: a high density front where the neutral gas is compressed and ionized, and a low density background, mainly the residual plasma diffusing randomly. The results indicate no distinct difference of the propulsion mechanism although the initial condition is varied,which is in accordance with studies on He and NH310,11.The property indicates the use of wider monotonic or polynomial gas is possible,with no significant change on the thrust generation mechanism within the range studied.

The magnetic field at different energies ofm=8 mg is shown in Fig.22. At the moments studied, the induced field is mainly generated by the thermal expansion of plasma, and similar to the density distributions, the magnetic fields at different regions,especially the conical turning and sheet regions,have similar distributions.

3.2.2. Comparison of velocity and total impulse at different propellant mass and energy levels

The average axial velocity of the high density sheet is shown in Fig.23, the velocity is calculated by:

The density range is chosen to incorporate the majority positive moving particles and to exclude the particles at the edge of the sheet, which may diffuse to the low density region and with negative velocity. The velocity increases with higher energy and lower propellant mass, as maximum velocity increases from about 15 km/s at 750 J to about 21 km/s at 1470 J atm=2 mg, within the decoupling distance. The trends show that once the decoupling distance is reached, the plasma moves nearly without power source and mainly thermal expansion exists, and the velocity decreases subsequently.

In the down sides of the velocity curves, local peak values appear in these cases, and the corresponding time decreases with lower propellant mass. The moments ofm=2 mg are about 14.7, 10.3, 8.5 μs at 750, 1080, 1470 J, respectively.The local peaks show increment of axial velocity of the flow when passing the turning of conical pylon, as also can be seen in Fig.17, where the main stream approaches the turning at 10 μs ofm=4 mg, and the more time needed as mass is increased is a result of lower kinetic obtained within the decoupling distance.

The total impulses,calculated by the time integration of the transient force,at 2–8 mg and 750–1470 J are given in Fig.24,the mass interval calculated is 0.5 mg at each energy level.The impulse rises from 1.3×10-2, 1.8×10-2, 2.2×10-2N·s at 2 mg to 2.6×10-2, 3.6×10-2, 4.4×10-2N·s at 8 mg ofE=750, 1080, 1470 J, respectively. The monotonically increase of the impulse as propellant mass grows is in accordance with the calculation of Ar,He and NH310,11at other coil radius, showing that the operation mechanism is similar at different dimension. The trend suggests the evaluation of performance in a predictable way at certain energy level.In Ref.4,the impulse versus mass curve was fitted by analytical expression at about 0.5–3.0 mg for NH3.The discrete points in Fig.8 can be approximate as:

Fig.23 Average axial velocity of current sheet.

Fig.24 Total impulses variation at different mass and energy levels.

The similar variations of impulse with mass on different propellant make the performance prediction of other propellant that can be used possible, at least in the range studied here, and an effective way which can reduce the experiment efforts is practical: by curve fitting the impulses at discrete mass, the values at certain energy level are obtained preliminarily.

4. Conclusions

In this paper, the transient MHD method, including the high temperature thermodynamic and transport, and plasma electrical models, is developed, which is verified by MHD vortex problem and the experiment results of PIT. The twodimensional ionization and acceleration characteristics of the inductive plasma are then studied in the range ofmfrom 2 mg to 8 mg andEfrom 750 J to 1470 J , and the variations of density, velocity, magnetic field, energy deposition and impulse with time are analyzed, the calculations of the flow field find that:

(1) The Joule heat dominates the energy deposition in the pulse rise side, and the value maintains a high level in this period, which is 1010W/m3from 1 μs to 2 μs atm=4 mg,E=1267.5 J, and is two magnitudes higher than the deposition in the down side.

(2) The velocity rapidly rises to maximum value, synchronizing with the initial pulsed current rise, and radial velocity is negligible compared with axial velocity.Radial component grows mainly after the first half period.

(3) The current sheet maintains the structure during the whole period studied. The core density of sheet is one or two magnitudes larger than that near the coil. The decrease of temperature in axial direction leads to the increase of low ionization level ions and reducing of high levels.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Nos. 11675040 and 11702319).