Rarefied flow regime of an underexpanded supersonic jet

. The results of a numerical study of a supersonic underexpanded jet flowing from a conical nozzle into a rarefied environment are presented. The modeling was performed by the direct simulation Monte Carlo method. The range of parameters corresponding to the expansion ratio 40 < n < 240, the hydrodynamic regime of the flow in the nozzle and the rarefied regime of interaction of the jet with the environment, characterized by Knudsen numbers in the range 0.03 < Kn L < 0.2, is considered. It is shown that in the specified range of Kn L a drastic rearrangement of the flow structure occurs. For the lower limit of the Kn L range, a shock-wave structure typical of a highly underexpanded jet is observed. For Kn = 0.2, it completely degrades. Data have been obtained on the process of dimer formation in an expanding jet. An increase in the mole fraction of dimers with increasing distance from the nozzle throat is shown both in the conical nozzle and in the initial section of the jet expansion region.


Introduction
Interest in the study of rarefied regimes of interaction of supersonic jets with the environment is associated with the development of aerospace applications and vacuum technologies.Evaluation of the parameters of rocket engine jets flowing into a rarefied space is necessary when designing spacecraft in order to exclude the harmful effects of the jet on the spacecraft surfaces [1].In vacuum gas jet technologies for the synthesis of nanomaterials, the parameters of the outflow of a highly underexpanded jet predetermine the conditions for the formation of clusters in the volume of the gas phase, and in the case of the growth of a nanostructured film directly on the sprayed surface -the rate of deposition of the material on the substrate [2].
Underexpanded jets are the subject of many experimental and theoretical studies [3][4][5][6][7][8][9][10][11][12][13][14].The parameters that determine the flow pattern of a supersonic jet in a submerged space include: the expansion ratio n = pa /p∞ (pa is the gas pressure at the nozzle exit, p∞ is the ambient pressure), the Mach number at the nozzle exit Ma, the angle of inclination of the nozzle contour in the exit section a, the Reynolds number determined from the parameters of the gas at the nozzle exit Rea, the specific heat ratio , the Prandtl number Pr and the Schmidt number Sc for the gases of the jet and the surrounding space [4].Instead of Rea, the Reynolds number Re* is often used, determined from the parameters of the gas in the throat.The mode of interaction between the outflowing jet and the environment can be characterized by an additional parameter, the Reynolds number Re, determined from the jet exit velocity, the density and viscosity of the environment, and the characteristic size of the jet.The distance to the Mach disk is usually chosen as the characteristic dimension.
The modes n >> 1 and Re  1000 are typical for vacuum technologies.These regimes are relatively well studied experimentally [5,8].However, numerical studies of such jets are few.The latter circumstance is associated with the nonequilibrium nature of the flow in the far field of the jet and the need to use kinetic approaches to simulate the flow.In the papers [13,14], kinetic by nature the direct simulation Monte Carlo method (DSMC) was used to model the underexpanded jet.The papers consider the case of sonic outflow and the ratio N = p0 /p∞ <5000 (p0 is the pressure in the stagnation chamber).
The regularities that manifest themselves when the Re number decreases below 1000 include, firstly, the transition from the Mach to the X-shaped configuration of shock waves in the initial section of the jet flowed from supersonic nozzle [5,8].Secondly, with a decrease in the number Re below approximately 100, thickening and blurring of the characteristic shocks in the flow region are observed [9].
It should be noted that the outflow of a highly underexpanded jet is accompanied by a significant drop in temperature in the initial section and may be accompanied by gas condensation.In turn, gas condensation can lead to a change in the characteristic dimensions of the jet [7,8].
This paper is devoted to the numerical study of a highly underexpanded low-density supersonic jet outflowed from conical nozzle in the range of parameters n >> 1 (N > 6500), 800 < Re* < 3000, and Re < 100.The work has two goals.First, to study the influence of Re (or the related Knudsen number) on the structure and parameters of an underexpanded jet.
Second, to study the processes of cluster formation in the jet and the possible influence of the condensation process on the gas dynamics of the flow.

Statement of the problem and simulation method
An axisymmetric outflow of gas from a chamber through a conical nozzle into a submerged space is considered.The outflowing gas and the ambient gas is argon.The pressure and temperature of the gas in the chamber are p0 and T0, respectively, in the surrounding space are p∞ and T∞, T0/ T∞  1.The nozzle is a truncated diverging cone with a height of 6 mm, a minimum section (throat) diameter adjoining the chamber d* = 1 mm, and an exit section diameter da = 3 mm.The dimensions of the considered part of the chamber before entering the nozzle are 2d* 1.35d*.The dimensions of the rectangular region of jet expansion outside the nozzle are 300d*40d*.The boundary of the expansion region adjacent to the nozzle exit is a solid wall.The origin of coordinates corresponds to the center of the throat.
Simulations were carried out for three cases, which differ in the expansion ratio n, the Re* and Re numbers.The Re* number is defined as where *, v* and T* are the gas density, velocity and temperature in the throat, respectively, (T*) is gas viscosity corresponded to the temperature T*.In the case T0/T∞  1, instead of the Re parameter (see section Introduction), the parameter ReL is usually used where N = p0 / p∞.Case 1 corresponds to the parameters n = 239.6,N = 27920, Re* = 830 and ReL = 5; case 2 corresponds to n = 60.2,N = 6980, Re* = 830 and ReL = 9.9; case 3 corresponds to n = 40, N = 6980, Re* = 2503 and ReL = 30.The Knudsen number, which characterizes the mode of interaction between the jet and the ambient gas, can be introduced into consideration.The Knudsen number is defined as KnL = λ∞/L  1/ReL (λ∞ is the mean free path of molecules in the environment) [3].For calculation case 1 KnL = 0.2, for case 2 KnL = 0.1, for calculation case 3 KnL = 0.033.
The study was carried out by the direct simulation Monte Carlo method (DSMC) [15,16].The NTC (no time counter) collision scheme was used in the calculations.The simulations were carried out using a parallel algorithm based on the decomposition of the computational domain using MPI technology on the Polytechnic-Tornado cluster.The number of computational particles in the region ranged from 4•10 8 to 10 9 , and the number of cells from 10 6 to 3.5•10 6 , depending on the calculation case.
The boundary conditions at the remote boundaries inside the chamber corresponded to a gas at rest with parameters p0 and T0.At the outer boundaries of the expansion region, particles reaching the boundary were removed from the calculation, and a particle flux was generated inside the region according to the semi-Maxwellian velocity distribution function with parameters p∞ and T∞.Particles are diffusely reflected from the walls of the chamber, the nozzle, and the boundary of the expansion region adjacent to the nozzle exit with total energy accommodation and temperature T0.At the initial time, there were no particles in the expansion region and nozzle.The calculation continued until a steady state was established, at which the number of particles in the computational domain did not change.
Elastic collisions of argon atoms were described in terms of the VHS model (variable hard sphere model) with following parameters: viscosity index  = 0.81, reference diameter dref = 4.17•10 -10 m and reference temperature Tref = 273K [15,16].
For case 2, inelastic processes were also considered.The forward process of formation of an argon dimer by a three-particle collision of atoms These reactions were described in terms of the TCE model (total collision energy model) [15][16][17].The kinetic model assumes that the reaction occurs upon collision of particles with the probability p1(3) for the dimerization process (3) or p2c for the collisional decay of the dimer (4).The probabilities are a function of the total energy of the colliding particles (kinetic energy of relative motion and internal energy).This functional dependence includes parameters of macroscopic reaction rate constants written in Arrhenius form.For reaction (3), the forward rate constant had the form [15] , where A = 3.863•10 -43 m 6 K 0.504 s −1 , b = -0.504.The macroscopic rate of reverse reaction (4) was determined based on the principle of detailed balance, taking into account the rate of the forward reaction and the equilibrium constant where C = 3.254•10 26 m −3 K -0.872 , α = 0.872, ε2 = 1.98•10 -21 J. Parameters of reaction constants are adapted taking into account the data [18].
The argon dimer has two rotational degrees of freedom and one vibrational, characterized by double energy capacity.Quantum effects were not taken into account, and the internal energy distribution function was assumed to be continuous.The Larsen-Borgnakke model was used to simulate the energy exchange between internal and translational degrees of freedom in dimer-monomer collisions.The probability of VRT exchange during collisions of the dimer with argon atoms was assumed to be 0.1.
The used DSMC code has been verified in a number of papers.In the current paper, for the purpose of additional verification, a comparison was made with the results of the numerical solution of the Navier-Stokes equations for the flow through the nozzle in the described statement (cases 1 and 2).The flow regime through the nozzle is close to hydrodynamic, and the solutions based on the DSMC method and the Navier-Stokes equations should coincide.The flow in the nozzle is quite complex and includes an oblique shock, which is formed due to the appearance of a region with a positive pressure gradient along the surface of the nozzle behind the throat [19].Figure 1 shows the density and velocity distributions on the nozzle axis (0 < X/d* < 6).It can be seen that the results agree well.

Results
The results of calculating the flow parameters in the considered range of Reynolds numbers ReL: 5  ReL  30 demonstrate a drastic change in the flow structure.Figure 2 shows the density distributions on the flow axis.For case 3 (ReL = 30) at a distance of X/d*  75, a shock is observed.The density of the gas increases by about an order of magnitude.Further, after the shock downstream, the density decreases nonmonotonically, approaching the density of the surrounding space at a distance X/d* = 300.Such a picture is typical for relatively small Reynolds numbers corresponding to the presence of an X-shaped configuration of shock waves in the initial section of the supersonic jet [5].The region of free expansion is limited by a blurred oblique shock that comes to the axis in the region X/d*  75 (Fig. 3).The presence of a reflected shock leads to a turn of the streamlines.The presence of a characteristic shockwave structure is also clearly seen in the field of the Mach number (Fig. 4).Beyond the region of interaction of the incident shock with the axis, the flow is supersonic.The Mach number on the axis exceeds 20, then drops after the shock to approximately the value M = 3, then increases to 5, and then decreases downstream.For case 2 (ReL = 9.9), the shock-wave structure in the expansion region degrades.There is no jump in the axial density distribution (Fig. 2).However, in the region of 40 < X/d* < 60, a plateau is observed, which is the "remnant" of the diffuse shock-wave structure.The Mach number reaches a maximum value of about 13 on the axis in the region X/d*  50.In the rest of the region, the Mach number decreases monotonically both in the axial and in the radial direction.In this case, the field of Mach numbers is somewhat deformed (Fig. 5).
For case 1 (ReL = 5), the shocks are completely blurred.The density of the gas decreases monotonically with distance from the nozzle exit, reaching at infinity its value in the surrounding gas (Fig. 2).In contrast to the density, the Mach number initially increases, then the outflowing gas is decelerated by the environment (Fig. 6).
From the Reynolds number Re* one can pass to the Knudsen number Kn*.The Knudsen number Kn* ~ M*/Re* for cases 1 and 2 is 1.5•10 -3 , for case 3 -5•10 -4 .Thus, the flow regime in the nozzle is collisional and close to hydrodynamic.The temperature in the nozzle and in the initial section of the jet expansion drops by fifty times (Fig. 7).Thus, the appearance of gas condensation effects can be expected.For case 2, an additional simulation was carried out taking into account the process of formation and decay of dimers according to equations ( 3) and ( 4).The dimer density field is shown in Figure 8, and the axial distributions of the density and mole fraction of dimers are shown in Figure 9.The maximum numerical density of dimers exceeding 2•10 -4 n0 (n0 is the stagnation numerical density) is observed in the chamber region.Further, the expansion process leads to a decrease in the numerical density of dimers.The sharpest drop in density is observed in the region of the throat before the shock inside the conical nozzle.After the shock, the density and temperature of the monomers increase (see Fig. 1).For case 2, the change in these parameters before and after the shock does not exceed 1.5 times.Further downstream, the gradients of these parameters turn out to be significantly smaller with respect to the region of the throat.The number of triple collisions is proportional to the monomer numerical density to the third power, and the dependence of the collision cross section on temperature is weaker.After the shock (X/d* > 1), a nearly linear dependence of the dimer density on the coordinate is observed.The mole fraction of dimers in the nozzle and the beginning of the jet expansion region increases with increasing distance from the throat.At a distance of the order of X/d*  20, the mole fraction freezes, reaching a value of the order of 10 -3 .The increase in the mole fraction is associated with two factors.Basically, with the predominance of the process of formation of dimers over the process of decay in the jet.Also, an additional contribution is associated with the effect of focusing the heavier species (dimers) to the flow axis during the expansion of the gas mixture [20].
Thus, the case under consideration is characterized by a small mole fraction of dimers in the flow field.The dimerization process practically does not affect the distribution of gas-dynamic parameters in the jet.

Conclusion
The paper presents the results of simulations of the outflow of a supersonic jet from a conical nozzle into a rarefied ambient gas.The range of parameters corresponding to the hydrodynamic flow regime in the nozzle and the transitional flow regime in terms of the Knudsen number in the region of interaction of the outflowing jet with the surrounding gas is considered.It is shown that the drastic restructuring of the flow occurs in a narrow range of KnL (ReL) numbers.The structure characteristic of a highly underexpanded jet flowing into a rarefied environment is observed at ReL = 30.Significant degradation of the structure occurs at ReL = 10, and at ReL = 5, there are no shocks in the region of jet expansion and a monotonous variation in density is observed.
For one of the cases considered in the work, the flow field was simulated taking into account the dimerization process.Data on variations in the density and mole fraction of dimers in the flow region are presented.The mole fraction of dimers in the flow field does not exceed 0.1 % and does not significantly affect the gas-dynamic parameters of the jet.

Fig. 1 .
Fig. 1.Argon density and velocity distribution in the nozzle for case 2. ρ0 is the stagnation density in the chamber, a0 is the sound velocity in the chamber.

Fig. 3 .
Fig. 3. Density field and isolines for case 3 in the vicinity of shock interaction with axis.

Fig. 9 .
Fig. 9. Dimer density and mole fraction variation along flow axis for case 2.