APP下载

Numerical simulation and experimental validation of multiphysics field coupling mechanisms for a high power ICP wind tunnel∗

2021-06-26MingHaoYu喻明浩ZheWang王哲ZeYangQiu邱泽洋BoLv吕博andBoRuiZheng郑博睿

Chinese Physics B 2021年6期

Ming-Hao Yu(喻明浩) Zhe Wang(王哲) Ze-Yang Qiu(邱泽洋)Bo Lv(吕博) and Bo-Rui Zheng(郑博睿)

1Faculty of Mechanical and Precision Instrument Engineering,Xi’an University of Technology,Xi’an 710048,China 2School of Automation and Information Engineering,Xi’an University of Technology,Xi’an 710048,China

Keywords: inductively coupled plasma,multiphysics field,coupling mechanism,simulation and experiment

1. Introduction

Reconnaissance satellites,manned spacecrafts,and space shuttles should return to the ground in time after they have completed their orbital missions. However, because of their high speed flight, they encounter severe aerodynamic heating problems when they re-enter the dense atmosphere. Accordingly,a perfect thermal protection system is required to ensure the safety of the spacecraft and astronauts. To improve the thermal protection system of re-entry vehicles,an inductively coupled plasma (ICP) wind tunnel is built on the ground for investigating the thermal protection materials for these vehicles.

The ICP wind tunnel can provide flawless long-term stable operation of a high enthalpy plasma jet and accurately simulate the flight environment of re-entry aircraft.It is,therefore,a piece of ideal equipment for investigating the plasma sheath of hypersonic aircraft. The ICP wind tunnel is mainly composed of an intake part, a high frequency discharge part, and a test part. The structural diagram of an ICP wind tunnel is shown in Fig.1. The high frequency discharge part is mainly comprised of the high frequency power supply,induction coil,and quartz tube.The electrified induction coil is wound around the high temperature resistant quartz tube. After the working gas enters the quartz tube from the air inlet,it is rapidly heated in the induction coil area to form a high temperature high enthalpy induction coupling plasma. The plasma produced in the quartz tube (also called the plasma torch) usually continues to flow forward into the large test chamber. In this chamber,the ICP can be used to test the aerodynamic performance of the aircraft,develop thermal protection materials[1]and investigate communication blackouts and electromagnetic wave transmissions.[2]

The ICP has been widely utilized in aerospace, materials, chemical industry, and other fields since the 1960s and 1970s. Presently,with rapid industrial development,many researchers have performed experiments and numerical simulations on ICPs.[1–5]In an experimental research,Moonet al.[4]measured the plasma density and coil voltage in ICP by using floating and high voltage probes,respectively,and studied the influence of coil turns on the ICP.They found that the plasma density increased with the number of coil turns because of the enhancement of inductive coupling. Moreover, the coil voltage was affected by both the plasma density and the coil.Their research can considerably aid the design of future ICP sources.Kimet al.[5]performed experiments on a CO2ICP that had a pressure of 100 mTorr (1 Torr=1.33322×102Pa) and power of 13.56 MHz. Electron energy distribution functions were measured using a Langmuir probe, and other plasma parameters, such as electron number density and electron temperature,were also obtained. According to the electrostatic probe technology, Panet al.[6]developed an electrostatic probe diagnosis system suitable for the test environment of a high frequency plasma wind tunnel. This system was employed to diagnose the flow field parameters of the ICP wind tunnel. Fujitaet al.[7]measured the temperature parameters of plasma in an ICP wind tunnel by means of a radiation spectrum.

Fig.1. Structural diagram of inductively coupled plasma wind tunnel.

Regarding the numerical studies of the ICPs,Ivanov and Zverev[8]numerically studied the influence of hydrogen on the parameters of plasma in the ICP torch with considering chemical reactions of Ar with H2. They found that as the hydrogen percentage increases over 10% at a constant total Ar+H2mass flow rate, the electrical conductivity of plasma increases slightly, but the enthalpy of plasma increases significantly due to high value of the dissociation energy of hydrogen molecules. Degrezet al.[9]also numerically studied the air ICP under a local thermodynamic equilibrium and a chemical non-equilibrium. They provided the calculation results of the streamline, temperature, and chemical concentration field in the air ICP torch. Alavi and Mostaghimi[10]numerically and experimentally studied a novel ICP torch which possesses a conical nozzle. The inductive coil were wound around the conical nozzle rather than around the conventional quartz tube. They found that the new torch has higher power density, better plasma stability, and better resistance against extinguishing factors than the traditional torch. Namet al.[11]conducted a numerical analysis of the power distribution and heat flow characteristics of Ar–N2ICP with nitrogen content of 0 mol%–50 mol%at a plasma power level of 50 kW.They obtained relevant data, such as heat loss and load resistance.They also employed numerical simulation and equivalent circuit analysis to optimize the design and analysis of a high power ICP torch system with a vacuum tube oscillator. Lu and Feng[12]investigated the effects of the operating parameters on the ICP properties by a non-equilibrium magnetohydrodynamic model. Their numerical model can aid in determining the optimal conditions for producing the ideal plasma.Sunet al.[13,14]studied ICP plasma characteristics by using two-dimensional self-consistent fluid model. It was found that the spatial distribution of plasma characteristics was controlled by changing the bias power of single frequency bias or high frequency bias power ratio of double frequency offset. The plasma density and ion species ratio were adjusted by changing bias voltage and pulse frequency. Gaoet al.[15]found that different RF bias voltage has different effect on plasma characteristics: in the low RF voltage range, the electron number density decreases and the effective electron temperature increases; in the high RF voltage range the electron number density increases and the effective electronic temperature decreases.

In addition, with the rapid development of commercial softwares such as ANSYS, COMSOL,etc., many research groups prefer to simulate the ICP using these as simulation tools. For example, Zhuet al.[16]used the ANSYS FLUENT software to simulate a two-dimensional ICP torch and derived the electromagnetic field distribution,temperature distribution,and velocity distribution in the torch. Chenet al.[17]employed the ANSYS FLUENT software to calculate the temperature of pure argon thermal plasma and the flow field distribution in the plasma generator device. They found that the temperature distribution and flow field distribution can be controlled by changing the corresponding device structure and working frequency. Zhao[18]established a fluid model of argon ICP by using the COMSOL software.It was observed that the temperature of the argon ICP changed non-monotonically with power. To explain this phenomenon, a new equation for the average electron energy was proposed. Leiet al.[19]employed the COMSOL software to simulate the ICP discharge and studied the influences of different coil designs and gas compositions on the discharge characteristics of a twodimensional ICP. The best discharge characteristic is considered to be afforded by a 20-mm coil gap, with the first coil positioned 130 mm away from the top of the device, and the best discharge mode is considered to be argon discharge mode.This numerical model is a typical ICP discharge,and its results are important for further experimental studies.

Although the above-mentioned research groups have conducted in-depth research on the numerical simulation of ICP,the majority of their simulation results have not been validated experimentally. Moreover,most of researches focused on low power(1 kW–50 kW)ICP wind tunnels. There is a lack of detailed studies to reveal the flow characteristics and the multiphysics coupling mechanism of the high power(such as a 1.2-MW class) ICP wind tunnel. For example, on the one hand,for a low power ICP wind tunnel,[20–23]its coil current,working pressure, plasma density are usually small;[20–22]hence,its electromagnetic field is small,[22,23]and the coupling effect between the electromagnetic field and the flow field is also relatively weak. One the other hand, for a high power ICP wind tunnel,[24,25]its working parameters such as the input power,operating pressure and the discharge frequency are quite different from those for a low power ICP wind tunnel.Its plasma density, electric-field intensity, plasma temperature, and heat flux are all higher than those of the low power one.[26,27]Thus,for the high power ICP wind tunnel,the coupling effects among the electromagnetic, flow, chemical, and thermodynamic fields are much stronger and more complex.Therefore,further researches on the high power ICP wind tunnel should be conducted to make clear such complex coupling effects. While so far there has been few such work available.

In the present study,a 1.2-MW high power ICP wind tunnel is studied as a research object. Furthermore,the coupling mechanisms of the electromagnetic field,flow field,chemical field,and thermal field of the wind tunnel are explored and discussed by combining numerical simulation with experimental validation. The electron number density, temperature, pressure,electric field intensity,Lorentz force,and other important parameters are discussed and investigated through numerical simulation. The plasma parameters under different working conditions are compared with each other to reveal the distribution characteristics of the flow field, electromagnetic field,chemical field, and thermal field of the high power wind tunnel. The accuracy of the numerical simulation is verified by comparing the simulation results with experimental data.[24]

In the following sections, the computational grid, operating conditions, control equations, numerical solution methods,and boundary conditions of the high power ICP wind tunnel are presented. Thereafter,by analyzing the distribution of flow field parameters(e.g.,temperature,electron density,and pressure) and electromagnetic field parameters (e.g., electric field intensity and Lorentz force) of the ICP wind tunnel under two different power values, the interaction and influence mechanism between the flow field and electromagnetic field in the high power ICP wind tunnel are determined. Finally,the accuracy of the research results is verified by comparing the experimental results with the numerical results.

2. Research methods

The research object of this study is an ICP wind tunnel with a power of up to 1.2 MW. The plasma torch is made of a quartz tube with an inner diameter and a wall thickness of

80 mm and 3 mm,respectively;six electrified coils are wound around this tube. As shown in Fig.2, the calculation domain of the electromagnetic field is−100 mm≤x ≤570 mm and 0 mm≤y ≤450 mm,with each interval consisting of 161×97 grid nodes;the ICP torch consists of 111×41 grid nodes. The whole flow field grids including the vacuum chamber are composed of 329×41 nodes as shown in Fig.3.

To calculate the electromagnetic field near the coil and air velocity at the air inlet accurately,the calculation grid is encrypted in both the area of the induction coil and the area of inlet.In this study,it is assumed that the plasma is steady,electrically quasineutral and two-dimensionally axisymmetric. The fluid micro-element is in a thermochemical non-equilibrium state. In the simulation calculation, the Lorentz force and induced magnetic field generated by the plasma are considered.

Fig.2. Computational mesh of electromagnetic-field and ICP torch.

Fig.3. Computational mesh of the whole flow field.

Usually, the plasma discharge is ignited by an appropriate argon injection and thereafter,gradually replaced by an air mixture (16 g/s standard dust-free atmospheric air) until stable working conditions are reached. To measure the temperature of plasma particles,spectrum analyzers are located at 195,320,and 445 mm away from the plasma torch outlet. For the simulation study,two typical working conditions: cases 1 and 2 for the 1.2-MW ICP wind tunnel are selected. The specific working conditions are listed in Table 1.The influences of different working conditions on the plasma temperature,velocity,pressure,electron density,electric field intensity,and other parameters are compared to demonstrate the flow field and electromagnetic field distribution characteristics of the high power ICP wind tunnel. These comparisons are based on the multifield coupling solutions of the flow field,electromagnetic field,thermodynamic field,and chemical field.

Table 1. Working conditions of computational cases.

2.1. Governing equations

The vector form of the system of non-equilibrium magnetohydrodynamic equations is as follows:

The gas temperature consists of the translational temperature(Ttr),rotational temperature(Trot),vibration temperature(Tvib),and electronic temperature(Te). The conservation vectorQ, flux vectorsFandFv, and source term vectorWin Eq.(1)are expressed as follows:

whereδi,j,FLi,Sjoule,Sint, ∆h0s, andΘvib,sdenote the Kronecker function,Lorentz force,Joule heating efficiency,internal energy exchange rate,component formation enthalpy,and vibration characteristic temperature, respectively; subscriptsrepresents the chemical species;Mis the molecular number;˙ωsis the mass generation rate of substancesand is the source term of the mass conservation equation of the particle at line 4 inW,Xs,andDsare the mole fraction and diffusion coefficient ofs,respectively;ui(oruj)andxjrepresent the velocity and tensor in the rectangular coordinate system, respectively;τijandqjare the stress tensor and the heat flux vector,respectively,and expressed as

The total internal energy(E)of the gas is defined as follows:

The translational energy(Etr),rotational energy(Erot),vibration energy(Evib),and electronic internal energy(Ee)of gas particle are calculated by the corresponding equations given in the work of Yuet al.[22]

The gas viscosityµ, and translationalλtr, vibrationalλvib, rotational thermal conductivityλrotin the above equations are all computed from Yos’ formulas.[28,29]The diffusion coefficientDsis calculated from the formula of Curtiss and Hirschfelder.[30]The conductivity of plasmaσand the electron thermal conductivityλeare evaluated from the thirdorder Sonine polynomials of air.[31–33]

The Joule heating rateSjoule, and the axial and radial Lorenz force,i.e.,FLxandFLycan be calculated from the following equations,respectively:

Here,Eθ=ER+ iEI.ERandEIare the real and imaginary part of the electric field intensity, respectively.They are obtained by solving the magnetic vector-potential equations[22]at each time step.

In addition,in Eq.(2),the rotational energy exchange rate(Sint,rot),vibration energy exchange rate(Sint,vib),and electron energy exchange rate(Sint,e)are calculated from the following equations:

Here,νb,rsandνf,rsare the stoichiometric coefficient of componentsbefore and after therchemical reaction, respectively;Msrepresents the molar mass of substances;sandjare the subscripts representing the component numbers of chemical species. The mass production rate of speciess,i.e.,˙ωsacts as a mass source term in the mass conservation equation of chemical species in the fourth row of the source termWin Eq. (2). It participates in iterative computations of the non-equilibrium magnetohydrodynamic equations mentioned above,and updates solutions at each time step.

2.2. Iterative computational methods

To iteratively solve the governing equations mentioned above, the flow field equations are discretized by the finitevolume method. The viscous term is calculated by the secondorder center difference method. To compute the inviscid fluxes of the flow field equations at the cell interfaces, the simple low dissipation advection upwind[36](SLAU)splitting method, which is a parameter-free and simple compressible numerical flux function of the AUSM family method,is used.The monotonic upstream-centered scheme for conservation law(MUSCL).[37,38]is adopted to avoid computational oscillations and spurious solutions in the region of high gradients.In the presently-used MUSCL, the physical properties at the cell interface are determined by linearly interpolating properties between neighboring cells. Thus, the accuracy of the inviscid numerical flux is kept to the second order. To overcome the numerical stiffness for a chemical nonequilibrium plasma,the point-implicit method[39]is utilized in this study. In addition, the point-implicit lower and upper symmetric Gauss–Seidel(LU-SGS)method[40]is used to update the solutions of the flow field results at each time step.

The finite-difference method is adopted to discretize the magnetic vector-potential equations to figure out the distributions of the electric field intensitiesERandEI. The secondorder, central-difference method is used to evaluate their numerical flux. The sub-relaxation iterative method is used to solve the discretized electromagnetic-field equations stably.The relaxation factor is set to be 0.5. When the relative error of the magnetic-vector potential betweennandn+1 steps is less than 10−5,the calculation is considered to be converged.

In this study,the local time stepping scheme is adopted to update the flow field results. When the residual of the mass,the momentum and the energy conservative equation are less than 10−5, 10−5, and 10−3, respectively, then the simulation is thought to be converged. In the present study, because the chemical nonequilibrium process of the chemical reaction of air is considered and modeled in the simulation,the residual of the energy conservative equation can hardly reach 10−5.Thus,when the residual of the energy conservative equation reaches 10−3and the flow parameters such as the temperature,velocity and pressure are no more changed after 1.2 million iterative computational steps,then we think the simulation is converged. In addition, all the equations and models mentioned above are programmed by the computer language Fortran 90 by ourselves. The computational mesh of the ICP wind tunnel is divided by Gridgen V15. All simulations are carried out by using the Tianhe-2 Supercomputer.

2.3. Boundary conditions

The working gas is injected into the quartz tube; the mass flow rate and initial gas temperature arem=16 g/s andT0=300 K, respectively. At the exit of the ICP wind tunnel laboratory,the working pressure is set to be a fixed value,and other flow properties are evaluated through the linear interpolation of adjacent internal grid points. For the central axis, the axisymmetric boundary condition (Qi,j=Qi,j+1) is adopted, whereQdenotes the conservation variable. For the magnetic vector-potential equations,the outer boundary of the electromagnetic field calculation area is set to bex=−10,x=570, andy=450 mm. The electric field intensitiesERandEIare set to be zero at these outer boundaries:Ei,0=Ei,1.It is assumed that there is no slippage nor catalytic effect nor pressure gradient on the wall of the plasma torch. The wall temperature is determined using the heat conduction equation,∂Tw/∂n=qjmax, whereqjmaxis the heat flux density at the inner wall. For the wall in the test room,a constant wall temperature boundary condition is adopted,i.e.,Tw=300 K.

3. Results and discussion

3.1. Numerical results

In this subsection, we will show the simulated results of flow,electromagnetic,and chemical fields,respectively. As is well known,when the input power is different,the flow properties such as the pressure, temperature, velocity also change and may demonstrate different physical results. So in order to make clear this issue for the high power ICPWT, we take case 1 and case 2 listed in Table 1 as baseline cases to show the flow field differences under the different input power conditions. The simulated results of these two cases will be discussed comparatively below.

Fig. 4. Distribution of axial Lorentz force (upper half) and radial Lorentz force(lower half)for case 1.

Figures 4 and 5 show the distributions of the axial and radial Lorentz forces in the plasma torches for case 1 and case 2,respectively. First,it can be observed that the maximum axial and radial Lorentz force position are under the third coil and fourth coil,respectively.The axial Lorentz forces in both cases are positive and negative, whereas the radial Lorentz forces are both always negative. The radial Lorentz force direction is from the wall to the central axis. Second, the maximum values of the radial Lorentz force are 2.23 and 3.05 times greater than the maximum values of the axial Lorentz force for case 1 and case 2, respectively. This proves that the radial momentum transfer of the plasma is stronger than the axial momentum transfer. Finally,because the greater the input power,the greater the electric field intensity and Lorentz force are.Hence, the Lorentz force for case 2 is greater than that for case 1.

Fig. 5. Distribution of axial Lorentz force (upper half) and radial Lorentz force(lower half)for case 2.

Figures 6 and 7 show the distribution diagrams of the electronic number density and imaginary part of the electricfield intensity in the plasma torch for case 1 and case 2,respectively. The electronic number density and electric field intensity considerably vary in the plasma torch,and their maximum value positions are under the third and fourth coil,respectively.The upper half of the figure shows that the electron number density is maximum at the central axis and minimum at the wall. This is because the radial Lorentz force is negative,and its direction points from the wall to the central axis. The electrons are gathered at the central axis under the action of the radial Lorentz force. The maximum electron number densities for cases 1 and 2 are 1.01×1022m−3and 9.05×1021m−3,respectively. The lower half for each of the figures shows that the maximum negative electric-field intensity is on the plasma torch wall and gradually decreases away from the wall. The maximum electric-field intensity in case 2 is greater than that in case 1. This is because a negative electric field is generated by the current of the external induction coil; it is maximum near the wall and minimum near the central axis. The positive electric field generated by the plasma is the largest near the central axis;it can offset part of the negative electric field.The electric field intensity in case 2 is greater than that in case 1,because the input power of the coil in case 2 is greater than that in case 1.

Fig.6. Distribution of electron number density(upper half)and imaginary part of electric-field intensity(lower half)for case 1.

Fig.7. Distributionof electron numberdensity(upper half)and imaginary partof electric-field(lower half)for case2.

Figures 8 and 9 show the distribution of the Joule heating rate and radial Lorentz force in the plasma torch for case 1 and case 2,respectively. The electromagnetic field and flow field of the wind tunnel are coupled and connected through the Joule heating rate and Lorentz force. First,it is observed that the distribution shape of the Joule heating rate is similar to that of radial Lorentz force,and the maximum value positions of Joule heating rate are also below the third and fourth coil.This proves that the position of the Joule heating phenomenon is mainly governed by the radial Lorentz force. Second, according to the formula, the Joule heating rate is directly proportional to the electric field intensity. Figures 6 and 7 show that the maximum electric field intensity is close to the wall,

whereas the maximum Joule heating rate is not on the wall but approximately 40 mm away from the wall. This is because a cooling device is used for the wall of the plasma torch,thereby reducing the Joule heating rate. The power in case 2 is greater than that in case 1. According to the proportional relationship between power and electric field intensity,the greater the electric power,the greater the electric field intensity and the Joule heating rate are. The maximum Joule heating rate in cases 1 and 2 are 7.0×107W/m3and 8.2×107W/m3,respectively.

Fig. 8. Distribution of Joule heating rate (upper half) and radial Lorentz force(lower half)for case 1.

Fig. 9. Distribution of Joule heating rate (upper half) and radial Lorentz force(lower half)for case 2.

Figures 10 and 11 show the streamline and velocity distribution in the entire flow field for case 1 and case 2,respectively;it can be observed that the maximum velocity in case 1 and case 2 are 110 m/s and 129 m/s, respectively. The maximum plasma flow velocity is found to be at the torch outlet.The maximum of flow velocity is on the central axis, and its minimum is on the test chamber wall. In the plasma torch, a strong vortex can be observed from the gas inlet to the second coils in Fig. 11. Such one or two recirculation flows are frequently found in the ICP simulations.[22,26,41]It has been reported that this vortex results mainly from the effect of the Lorentz force and the negative pressure gradient[22]generated in the coil region (see Figs. 12 and 14 below). On the one hand,the radial and axial component of the Lorentz force can affect corresponding momentum transfer of charged species significantly for a plasma flow.The strong radial Lorentz force strengthens the pinch effect,which can narrow the radial electrical channel of electrons and squeeze the plasma flow in the radial direction in the coil region.[22]On the other hand,both the negative axial Lorentz force and the negative pressure gradient beneath the inductive coils act as resistive forces to prevent the plasma flow from moving forward smoothly. As a result, some of the plasma elements move backwards and form a vortex near the center axis in the coil region.

Fig. 10. Distribution of streamlines (upper half) and velocity (lower half)for case 1.

Fig. 11. Distribution of streamlines (upper half) and velocity (lower half)for case 2.

In addition, near the wall of the test chamber, there are several notable vortexes. Their formation and intensity are affected mainly by the pressure characteristics of the spacious vacuum chamber. First,due to the high speed plasma flowing along the center axis from the torch outlet to the test chamber outlet, it will results in pressure differences between the chamber wall and the outer boundary of the plasma column in the spacious test chamber. The local pressure difference in a spacious chamber easily results in recirculation flow in the chamber. While such recirculation flows near the chamber wall hardly affect the flow properties of the plasma column near the center axis due to the far distance.

To check the correctness of the velocity simulated in our study,we consult several relevant references.[12,25,42]Approximate values of the simulated velocity are found for the high power ICPs.[25,42]For example, Vasil’evskii and Kolesnikov[25]simulated an air ICP flow under similar working conditions: the nominal input powerP= 90 kW,f=0.44 MHz,p=10 kPa,m=2.8 g/s. Their simulated maximum velocity at the torch exit is about 150 m/s.[25]Although their input power is lower thanP=160 kW in our case 1,their working pressure and mass flow rate are also lower thanp=15 kPa,m=16 g/s in case 1,respectively.Lu and Feng[12]has reported that the plasma velocity decreases as the pressure and the mass flow rate increase. Thus, for case 1 the simulated maximum velocity (110 m/s) has the same order of magnitude as the one (150 m/s) reported by Vasil’evskii and Kolesnikov.[25]Thus,the maximum plasma velocities obtained in our simulations are reasonable.

Figures 12–14 show the distributions of the static pressure for case 1 and case 2. The maximum pressure positions are both below the third and fourth coil. The maximum pressures for cases 1 and 2 are 15027.9 Pa and 15076.7 Pa,respectively; thus, the maximum pressure for case 2 is greater than that for case 1. From Fig.14 it follows that the pressure variation trends of these two cases are similar along the center axis,Along the positive direction of thexaxis, the pressure value under the induction coil first increases and then decreases until the pressure value tend to be an identical constant of approximately 1.5×104Pa in the test chamber for these two computational cases. However, for these two cases, their maximum pressure drops in the whole flow field are 27.9 Pa and 76.7 Pa,respectively. They are both small.

Fig.12. Distribution of static pressure in the whole flow field for case 1.

Fig.13. Distribution of static pressure in the whole flow field for case 2.

Fig.14. Comparison of simulated static pressure along center axis between case 1 and case 2.

Fig.15.Radial distribution of mole fractions for air species at the coil center x=225 mm for case 1.

Figure 15 shows the radial distributions of mole fractions for air chemical species at the coil centerx=225 mm for case 1. As can be seen from the figure,the atomic N and O are the dominant chemical species beforey ≤40 mm. It indicates that in the vicinity of the central axis,after a series of dissociation and ionization reactions, nitrogen molecules decompose into atomic nitrogen N and ionic nitrogen N+, oxygen molecules also decompose into oxygen atomic O and ionic oxygen O+.The mole fractions of electrons e−and ionic N+are about 0.1,and almost the same. That is to say, the negative charge of electrons just offsets the positive charge of ionic atomic nitrogen,so the simulated plasma flow is electrically quasi-neutral.With the increase of the distance from the central axis, the mole fractions of species N+2,O+,N+,e−,NO+decrease fast aftery ≥15 mm.Molecular nitrogen becomes the main chemical species aftery ≥40 mm.

3.2. Experimental validation of simulated results

To validate the correctness of the numerical methods and the simulated results presented in this study, the quantitative comparisons between the simulated and the measured temperatures are performed. Four simulations listed in Table 1 are carried out according to the corresponding experimental study in the 1.2-MW ICP wind tunnel for comparison.

First, following the spectroscopic experimental study of Cipulloet al.,[24]the air inductive plasma jet produced by the 1.2-MW ICPWT is investigated. Intrusive and nonintrusive measurement technique are combined to characterize the plasma flow at a broad range of pressure and input power. The high speed camera imaging technique is used to study the flow features of the plasma,and a frequency-spatial-domain elaboration procedure is adopted to analyze the measured data. In the experiment, the high speed camera imaging is performed by using a complementary metal–oxide–semiconductor camera to ensure the complete visualization of the plasma flow field. Three spectrometers are used to observe the plasma at three different distances from the torch inlet (x1=665 mm,x2=790 mm,x3=915 mm)along the central axis. The temperatures at these three positions are estimated by fitting the measured spectrum with the theoretical spectra.

Second, the simulated and measured temperatures are compared at these three positions for cases 1 and 3 (P=160 kW), and cases 2 and 4 (P=200 kW), respectively, in Figs.16 and 17. In these two figures,the curves and the symbols denote the numerical and experimental results, respectively. Additionally,because the error range for cases 1 and 3 is not obtained in the experimental study,[24]none of the error bars of the experimental points are depicted in Fig.16.

Third, as shown in Fig. 16, it can be observed that the temperature value sharply increases to a maximum and then slightly decreases in the plasma torch. Due to the strong external electromagnetic field generated by the electrified coil,the working gas entering from the inlet is heated rapidly by the electrified coil. So the gas temperature rises rapidly in the coil region. The simulated maximum temperatures for cases 1 and 3 are 10514 K and 10693 K, respectively. In addition,the measured (red circle and blue delta) and the simulated(curves) temperatures in Fig. 16 decrease with the increase of the torch outlet distance due to the interaction between the plasma and its surrounding cold gas. In the vacuum chamber,the simulated temperatures accord well with the temperatures measured at the three positions. For example, at the positionx3=915 mm, the experimental and numerical temperatures for case 1 are 6108 K and 5941 K, respectively. The relative error between them is 2.7%.

Fig. 16. Comparison between simulated and the measured temperatures along center axis for cases 1 and 3(P=160 kW).

Additionally,as can be seen in Fig.16,the simulated temperatures show a slightly increasing tendency near the chamber outlet and in a region of 800 mm≤x ≤880 mm for case 1 (p=15 kPa) and case 3 (p=10 kPa), respectively. This abrupt temperature fluctuation might be caused by uncertainties of chemical reactions of air in the simulation,because similar abrupt temperature increment and fluctuation were also found in the arc-heated air plasma simulation.[43]However,the fluctuation of the temperature simulated in our study is about 200 K. The relative error of the simulated temperature caused by this fluctuation is within 3.5%. So the effect of this temperature fluctuation might be ignorable. To make clear the definite reasons for this temperature fluctuation,this issue will be studied further in the future.

Fig. 17. Comparison between simulated and measured temperatures along center axis for cases 2 and 4(P=200 kW).

Finally,as shown in Fig.17,generally the simulated and measured temperatures also accord well with each other in the test chamber. Although a few deviations can be seen at the positionx3=915 mm for the 10-kPa case(blue curve and blue symbol delta), the maximum relative error is 6.8%, because the measured and the simulated temperatures are 6465 K and

6024 K respectively.In summary,the maximum relative errors between the simulated and measured results for cases 1–4 are all within 6.8%, and thus indicating that the numerical methods and results presented in this study are reliable and suitable.

4. Conclusions

In this study, the basic flow and electromagnetic properties of a 1.2-MW high power ICP wind tunnel,such as its distribution of the temperature,pressure,electronic density,electric field intensity, Lorentz force, Joule heating rate,etc., are obtained by solving the multiphysics coupling fluid dynamic equations. The simulated temperatures are compared with the experimentally measured ones,and they accord well with each other. So the numerical methods used in this study are reliable and suitable.

Through the present study, it is also found that the maximum values of the temperature,pressure,electric field intensity, Lorentz force, Joule heating rate, and speed are all generated near the center of the induction coil. The temperature increases to the maximum value in the plasma torch,and then decreases with the increase in the distance from the torch outlet. The pressure value is kept basically the same in the test chamber,but slightly changes in the plasma torch. The change trend of the pressure is similar to that of the translational temperature. In addition, it is also observed that the maximum value of the radial Lorentz force is practically twice or thrice bigger than that that of the axial the Lorentz force. As a result,the momentum transfer in the radial direction is stronger.Moreover,the movement of electrons is also strongly affected by the radial Lorentz force pointing to the central axis. The electron number density and positive electric field intensity generated by the electron are largest at the central axis. Consequently,the negative electric-field intensity generated by the coil current is partially offset by the positive electric field intensity generated by electrons. At the entrance of the torch and near the wall of the test chamber,it is found that the eddy flows are easily produced.

Acknowledgment

All simulations were performed by using the Tianhe-2 Supercomputer at the National Supercomputer Center in Guangzhou,China.