SENAC Vol. 53, No. （1 2020. 1）― 9 ―［共同研究成果］

Simulation of Nanopowder Growth-Transport
in Turbulent Thermal Plasma Flow Field
Masaya Shigeta
Joining and Welding Research InstituteThis paper discusses the theoretical models and numerical methods to simulate a turbulent thermal
plasma flow transporting nanopowder growing collectively. Turbulence models should be capableof expressing not only turbulent but also laminar states because both states co-exist in and aroundthermal plasmas that have large variations of density as well as transport properties. Demonstrative
simulations clearly present that an advanced numerical method with high-order accuracy should
be used for long and robust computation capturing steep gradients of nanopowder concentrationand plasma temperature and 3D dynamic motions of multi-scale vortices which are turbulent
features of thermal plasma flows with low Mach numbers.
1. Introduction
Thermal plasma offers a unique flow with high temperatures of not only electrons but also heavy species,
high chemical reactivity, intensive light emission, and variable properties. Because it is also electricallyconductive fluid, it responds to an electromagnetic field. Therefore, its plasma state can be sustained by Joule
heating and controlled by Lorentz force [1-5]. Thermal plasma is anticipated as a promising device fornanopowder fabrication with remarkably high rates [6].
However, the precise process control is still difficult because the nanopowder’s growth and transport must
be controlled simultaneously. The process involving collective growth by nucleation, condensation, andcoagulation and transport by convection, diffusion, and thermophoresis in/around a plasma is a complicatedmass transfer problem in the microsecond to millisecond time scales. At present, it seems impossible toobserve or measure such a spatially and temporally multi-scale process. Meanwhile, theoretical andnumerical studies are potent approaches to tackle this problem, although they are arduous as well.
Turbulence of thermal plasma flows has been recognized as an important matter. Turbulence offers strong
mixing effects on momentum, energy, and material molecules and particles because of multi-scale vortices.
Pfender et al. [7] proved that a thermal plasma jet forms vortices entraining ambient cold gas at its fringe.
These vortices break into smaller ones. Eventually, the plasma jet transits to a turbulent state with 3D andmulti-scale vortices. For nanopowder fabrication, nanoparticles are generated not in a high-temperatureplasma core but in interfacial regions between plasma and cold gas at low to intermediate temperatures [8].
Lower-temperature regions tend to be more turbulent with multi-scale vortices [9, 10]. The vortices strongly
affect the transport processes of nanopowder. The fundamental mechanisms of turbulence transition ofthermal plasma flows should also be clarified.

Fig. 1. A thermal plasma flow for nanopowder fabrication.

― 10 ―

SENAC Vol. 53, No. （1 2020. 1）

This paper discusses theoretical models and numerical methods to simulate especially a turbulent thermal
plasma jet transporting nanopowder. Demonstrative simulations are performed to highlight the importanceof numerical methods for long and robust computation capturing steep gradients of nanopowder
concentration and plasma temperature and 3D dynamic motions of turbulent multi-scale vortices in/around athermal plasma jet.
2. Model description
Thermal plasma is described by the thermofluid approximation with several typical assumptions: (i) theentire fluid region including plasma and non-ionized gas is in a local thermodynamic equilibrium state, (ii)the plasma is optically thin, and (iii) buoyancy due to density variation is negligible. The governing equations
are given as the conservations of mass, momentum, and energy:    u   0 ,(1) u     uu  P     2 S  2   uI    0 ,t3 (2)

t

and
 h      uh       h   P  u  P  Qrad  Qcon  Φ  0 ,t t CP(3)where  stands for the density of fluid, t represents the time, u is the velocity vector, P denotes the pressure, signifies the viscosity, I is the unit matrix, h is the enthalpy,  is the thermal conductivity, CP is the specificheat at constant pressure, Qrad represents the radiation loss, Qcon denotes the heat generation attributable tocondensation, and  denotes the viscous dissipation. Furthermore, S is the velocity strain tensor defined asS

1u   u tr2(4)where superscript tr denotes the transposition.
In engineering time- and spatial-scales, the aerosol dynamics approach effectively describes the growth
and transport processes of nanopowder produced by thermal plasma with several assumptions: (i)nanopowder consists of spherical nanoparticles, (ii) electric charge effects are neglected, (iii) nanoparticle
temperature is identical to the fluid temperature, and (iv) material vapour is treated as an ideal gas. Extending
the previous model [11] for nanopowder growth through homogeneous nucleation, heterogeneouscondensation and coagulation between nanoparticles, the governing equations that also express thenanopowder’s transport by convection, diffusion, and thermophoresis are written as [5, 8]

np
 n 
n   np / 6 1/ 6
   u   p      D p  p   J  2 2 0 n11f     K th  ln T   0 ,pt      f f f
f
   u        D p    Jg c   0 (nv  ns )n1p/ 3 f 2 / 3     K th  ln T   0 , (6)t     

and(5)

Simulation of Nanopowder Growth-Transport
in Turbulent Thermal Plasma Flow Field n 
n   nv    u   v      Dv  v   Jg c   0 (nv  ns )n1p/ 3 f 2 / 3 .
t      ― 11 ―(7)where n is the number density, D is the diffusion coefficient, and T is the temperature. The subscripts p, v,and s denote particle, vapour, and saturated state, respectively. The variable f is defined as f = npg, where g is
the average monomer number in a nanoparticle. 0 is the parameter related to collision frequency [11]. J is
the homogeneous nucleation rate and gc is the number of monomers composing a nanoparticle in a criticalstate [12]. Kth is the thermophoresis coefficient [13].
3. To express turbulent effects
A turbulent flow is composed of multi-scale vortices. A high-temperature argon thermal plasma flow
cannot have sub-millimetre vortices, while a low-temperature argon flow can have vortices larger than submillimetres [9, 10]. Simulation of a turbulent flow of thermal plasma interacting with non-ionized cold gas
requires a theoretical model and a numerical method that can express such a field.
This study adopts a Large Eddy Simulation (LES) approach with Coherent Structure (CS) model [14, 15]
to deal with laminar (non-turbulent) regions and turbulent regions in and around plasma at the same time.
The detailed discussion for the mathematical expression with the transformed equations is found in Refs. [8,10].
The widely ranging temperature causes large variations of the transport properties and density with several
orders of magnitude. Meanwhile, the Mach number is on the order of 10-3~10-1. Therefore, to obtain thesolution in a practically sufficient time-scale, a thermal plasma flow is calculated as an incompressible flow
with temperature- dependent density. Therefore, the following requirements should be satisfied on numericalmethods.
(i) Small vortices are resolved with a large grid element size to reduce computational costs.
(ii) Steep gradients of variables are also captured by the computational grid element width.
(iii) The discretization errors are small even with the large grid element width.
(iv) The time-integration errors are small even with a large time increment.
Table 1. Conventional method (Method-I in Ref. [10]).
Algorithm

SIMPLE [16]Convection terms

Upwind differencing with 1st-order accuracy

Transient terms

Backward Euler scheme with 1st-order accuracy

Diffusion terms

Central differencing with 2nd-order accuracy

Source terms

Central differencing with 2nd-order accuracy

Table 2. Advanced method (Method-III in Ref. [10]).
Algorithm

IPISO [17]Convection terms

Hybrid-upwind K-K differencing
with mainly 3rd-order accuracy [18]Transient terms

Adams-Bashforth-Moulton schemewith 3rd-order accuracy

Diffusion terms

Central differencing with 2nd-order accuracy

Source terms

Central differencing with 2nd-order accuracy

― 12 ―

SENAC Vol. 53, No. （1 2020. 1）

Fig. 2. Instantaneous thermal flow fields in/around a thermal plasma jet obtained by the conventional method.

Fig. 3. Instantaneous thermal flow fields in/around a thermal plasma jet obtained by the advanced method.

(v) A solution algorithm for incompressible flows is preferable even with large density variation.
(vi) Computation is robust even with large variations of density and transport properties of plasma andextremely large values of the source terms by nanopowder formation.
Tables 1 and 2 introduce two numerical methods for the demonstrative simulations. The conventional method(Method-I [10]) has often been used for thermal plasma flow simulations because it obtains a numericalsolution relatively easily. Meanwhile, the advanced method (Method-III [10]) is expected not only to express
turbulent features but also to perform numerically stable computation even with a large increment of timesteps. More discussion about the differencing schemes and solution algorithms is found in Refs. [8, 10].
4. Results and discussion
Figure 1 shows a schematic illustration of the computation domain. A non-transferred argon thermal plasmajet is ejected from a circular nozzle with the diameter of 8.0 mm. At the nozzle exit, the plasma jet has steadyprofiles of temperature with the maximum 12000 K and velocity with the maximum 400 m/s. These
conditions are obtained from the approximated curves [10] with the experiment data [7]. Silicon vapour is
supplied at 0.1 g/min with the plasma jet. The geometries of the domain are (x, y, z) = (0.0~255.8 mm, 51.2~51.1 mm, -51.2~51.1 mm). Computations are performed using a Cartesian staggered grid system withthe uniform spatial interval ofx =y =z = 0.1 mm and the time increment of t = 0.1 ms on asupercomputer SX-ACE (NEC) at Tohoku University, Japan. The thermodynamic and transport propertiesof argon thermal plasma and the material properties of silicon are obtained from Refs. [19] and [20],respectively.
Figures 2 and 3 show the snapshots of the instantaneous temperature and velocity fields and vortex
structures obtained by the conventional and advanced methods, respectively. It is noted that the region near
the jet nozzle (x, y, z) = (0.0~125.0 mm, -31.5~31.5 mm, -31.5~31.5 mm) is presented here. It is apparent
that the advanced method simulates a complex flow field with widely spreading multi-scale vortices. Tovisualize the vortices, the coherent structures of vortices were extracted from the velocity data by"Qcriterion" [21]. Here, the vortices are presented by the isosurfaces of the second invariant of the velocity
gradient tensor of the value 0.25 which was normalized by the mean velocity of 160 m/s and theSimulation of Nanopowder Growth-Transportin Turbulent Thermal Plasma Flow Field― 13 ―Fig. 4. Instantaneous distributions of nanopowder obtained by the conventional method.

Fig. 5. Instantaneous distributions of nanopowder obtained by the advanced method.

diameter of 8.0 mm at the nozzle exit. Many vortices are generated even far from the plasma jet cores. Thesevisualized vortex structures are similar to the Schlieren photograph by Pfender et al. [7].
Large vortices have higher temperatures whereas small vortices have lower temperatures as predicted on
the basis of Kolmogorov theory [9, 10]. The magnified snapshot Fig. 3(d) clearly shows a double-layer
structure of high-temperature thicker vortex rings surrounded by low-temperature thinner vortex rings in theupstream region. As they flow downstream, the high-temperature thicker vortex rings deform largely whereas
the low-temperature thinner vortex rings break up into smaller vortices.
Compared with the advanced method, the conventional method does not simulate the features of a thermal
plasma jet as shown in Fig. 2. Although vortices are generated in the upstream region, they disappear in thedownstream region because the upwind differencing scheme with the first-order accuracy for the convection
terms inherently includes a large discretization error that unphysically flattens the spatial gradients ofvariables. Therefore, the temperature and velocity fields are calm and the vortices attenuate by the nonphysical diffusional effects.
Figures 4 and 5 show the snapshots of the instantaneous number density and mean diameter distributionsof silicon nanopowders obtained by the conventional and advanced methods, respectively. The silicon vapour
in the plasma jet is transported by convection and diffusion. At the plasma fringe with lower temperature, thesilicon vapour becomes supersaturated and changes its phase into nanopowder by homogeneous nucleationand heterogeneous condensation. Because nanopowder collectively grows up and diffuses outside the plasma
region, nanopowder distributes widely in the field. The larger size regions coincide with smaller number
density regions. This result indicates that simultaneous coagulation decreasing particle number plays animportant role for nanopowder growth as well.
It is also obvious that the conventional method and the advanced method obtain very different results
although the governing equations, the boundary conditions, and the given parameters are the same. As
discussed above, the conventional method blurs gradients of variables. As a result, the vapour and thenanopowder experience unphysically large diffusions. Note, those methods obtains the very different
temperature and velocity fields, which also cause the differences of the vapour and nanopowder distributions.
The more detailed discussion is found in Ref. [22].
Using the advanced method, a numerical study has also presented that the turbulence of thermal plasma
and the growth of nanopowder can be controlled by an external magnetic field [5].

― 14 ―

SENAC Vol. 53, No. （1 2020. 1）

5. Acknowledgements
This work was partially supported by JSPS KAKENHI (Grant No. 16K13737). The numerical results wereobtained using supercomputing resources at Cyberscience Center, Tohoku University.
6. References[1] T. Sato et al., Int. J. Thermal Sciences, 40, 273 (2001).
[2] M. Shigeta et al., Int. J. Heat Mass Transfer, 47, 707 (2004).
[3] M. Shigeta and H. Nishiyama, Transactions of ASME, Journal of Heat Transfer, 127, 1222 (2005).
[4] Y. Tanaka et al., Applied Physics Letters, 89, 031501 (2006).
[5] M. Shigeta, Journal of Flow Control, Measurement & Visualization, 6, 107 (2018).
[6] M. Shigeta and A.B. Murphy, Journal of Physics D: Applied Physics, 44, 174025 (2011).
[7] E. Pfender et al., Plasma Chemistry and Plasma Processing, 11, 529 (1991).
[8] M. Shigeta, IEEJ Transactions on Electrical and Electronic Engineering, 14, 16 (2019).
[9] M. Shigeta, Plasma Sources Science and Technology, 21, 055029 (2012).
[10] M. Shigeta, Journal of Physics D: Applied Physics, 49, 493001 (2016).
[11] V.A. Nemchinsky and M. Shigeta, Modelling and Simulation in Materials Science and Engineering, 20,045017 (2012).
[12] S.L. Girshick et al., Aerosol Science and Technology, 13, 465 (1990).
[13] L. Talbot et al., Journal of Fluid Mechanics, 101, 737 (1980).
[14] H. Kobayashi, Physics of Fluids, 17, 045104 (2005).
[15] H. Kobayashi, Physics of Fluids, 18, 045107 (2006).
[16] S. V. Patanker and D.B. Spalding, Int. J. Heat Mass Transfer, 15, 1787 (1972).
[17] P.J. Oliveira and R.I. Issa, Numerical Heat Transfer B, 40, 473 (2001).
[18] S. Komurasaki, Proceedings of 7th Int. Conf. on Computational Fluid Dynamics, ICCFD7–3001, Kohala
Coast, USA (2012).
[19] M.I. Boulos et al., Thermal Plasmas Fundamentals and Applications, 1, Plenum Press, New York (1994).
[20] Japan Institute of Metals, Metal Data Book, Maruzen, Tokyo (1993). (in Japanese)[21] J.C.R. Hunt et al., Center for Turbulence Research Proceedings of the Summer Program, 193 (1988).
[22] M. Shigeta, Plasma Chemistry and Plasma Processing, (2020),(Available online https://doi.org/10.1007/s11090-020-10060-8).