2. RF argon plasma discharge

A one-dimensional (1-D) RF discharge between two symmetric parallel-plate electrodes separated by a distance of $D = 2.54\ \textrm {cm}$ is simulated. The RF discharge is filled with the argon gas at a pressure of 1 Torr and sustained by a sinusoidal external excitation $V = V_{\textrm {rf}} \sin (2{\rm \pi} f t)$ at $x = 0$ applied to the lower electrode where $f$ is the driving frequency, whereas the electrode at $x= D$ (upper electrode) is grounded. In order to describe the behaviour of an isolated dust particle in the discharge, we adopt for the plasma species a self-consistent fluid approach in drift-diffusion approximation. In addition, in order to examine the dynamics of dust we inject a spherical dust particle into the discharge as shown in figure 1. Assuming that the dust particle has no effects on the discharge characteristics, hence, the system of equations describing the discharge will be solved independently.

Figure 1. Schematic of argon RF capacitive discharge with an injected dust particle.

2.1. RF fluid model

The RF discharge model used in this study is based on the fluid description of the electrons and ions including the continuity equation, momentum balance equation and electron energy balance equation coupled with Poisson's equation. The continuity equations for electrons and positive ions are derived from Boltzmann equation and the particles fluxes are described by the drift-diffusion approximation (Zhao, Liu & Samir Reference Zhao, Liu and Samir2018)

(2.1)\begin{gather} \dfrac{\partial n_{e,i}}{\partial t} + \dfrac{\partial {\varGamma}_{e,i}}{\partial x}- R_{i} =0, \end{gather}

(2.2)\begin{gather} {\varGamma}_{e,i} ={\mp} n_{e,i}\mu_{e,i}{E} - D_{e,i}\dfrac{\partial n_{e,i}}{\partial x}, \end{gather}

where $n_{e,i}$ is the number density of electrons and ions, $R_{i}$ is the source term that contains only the rate of ionization reaction and ${\varGamma }_{e,i}$ represent the fluxes of electrons and ions. In (2.2), $\mu _{e,i}$ and $D_{e,i}$ are the mobilities and the diffusion coefficients of the electrons and the ions. The energy balance of electrons has been also derived from the Boltzmann equation and it is introduced to investigate the electron temperature $T_{e}$ along the inter-electrode discharge. It can read as follows

(2.3)\begin{equation} \dfrac{\partial}{\partial t}\left(\dfrac{3}{2}n_{e}k_{B}T_{e}\right)+ \dfrac{\partial {q}_{e}}{\partial x} +e{\varGamma} _{e} {E} + H_{i}R_{i} = 0, \end{equation}

where $T_{e}$ is the electron temperature, and ${E}$, $H_{i}$ and $k_{B}$ are the electric field, the energy loss coefficient with ionization and the Boltzmann constant, respectively. The heat flux of electron energy ${q}_{e}$ and electron heating rate $P_{\textrm {heat}}$ coming from the electron current are expressed as

(2.4)\begin{gather} {q}_{e} ={-}\kappa_{e}\dfrac{\partial{T}_{e}}{\partial x} + \dfrac{5}{2}k_{B}T_{e} {\varGamma}_{e}, \end{gather}

(2.5)\begin{gather} P_{\textrm{heat}} = e{\varGamma}_{e}E = eD_{e}\dfrac{\partial{n}_{e}}{\partial x}{E} + e\mu_{e}n_{e}{E}^2, \end{gather}

respectively, where $\kappa _{e} = 3(k_{B}D_{e} n_{e})/2$ is the electron thermal conductivity coefficient. Finally, the Poisson equation is used to calculate self-consistently the electric field created by the charged particles in the discharge domain. This equation connects the gradient of the electric field to the densities of charged particles as follows

(2.6)\begin{gather} \dfrac{\partial^{2}V}{\partial x^{2}} = \dfrac{e}{\varepsilon_{0}}(n_{e} - n_{i}), \end{gather}

(2.7)\begin{gather} {E} ={-}\dfrac{\partial V}{\partial x}. \end{gather}

In (2.6), $\varepsilon _{0}$ and $e$ are the permittivity of free space and the elementary charge, respectively.

2.1.1. Boundary conditions

The mathematical problem requires the resolution of the coupled system of (2.1)–(2.3) and (2.6) with a set of adequate boundary conditions. Here we have adopted the same boundary conditions that those used in previous papers (Lymberopoulos & Economou Reference Lymberopoulos and Economou1993; Samir et al. Reference Samir, Liu, Zhao and Zhou2017).

At the bottom electrode ($x =0$)

(2.8a–d)\begin{equation} \varGamma_{i} = n_{i}\mu_{i}E,\quad \varGamma_{e} ={-}k_{s} n_{e} - \gamma \varGamma_{i},\quad T_{e} = 0.5\ (\textrm{eV}) \quad \textrm{and}\quad V = V_{\textrm{rf}}\sin(2{\rm \pi} f t). \end{equation}

At the top electrode ($x =D$)

(2.9a–d)\begin{equation} \varGamma_{i} = n_{i}\mu_{i}E ,\quad \varGamma_{e} = k_{s} n_{e} - \gamma \varGamma_{i},\quad T_{e} = 0.5\ (\textrm{eV})\quad \textrm{and}\quad V =0\ (\textrm{V}). \end{equation}

where $\gamma$ is the secondary electron emission coefficient and $k_{s}$ is the electron recombination coefficient.

2.1.2. Initial conditions

The initial conditions used in this work are similar to those given by Samir et al. (Reference Samir, Liu, Zhao and Zhou2017) and Zhao et al. (Reference Zhao, Liu and Samir2018)

(2.10a–c)\begin{equation} n_{e} = n_{i}= n_{\varepsilon} \left[ { \varepsilon+16 \left( { 1-\dfrac{x}{D} } \right)^2\left( { \dfrac{x}{D} } \right)^{2} }\right], \quad T_{e}= 1\ (\textrm{eV})\quad \textrm{and}\quad V=0\ (\textrm{V}), \end{equation}

where $\varepsilon$ is a small positive number and $n_{\varepsilon }$ is the value of the initial charged particles densities. The remainder of the parameters used in this simulation can be found in table 1.

Table 1. Data used in calculation (Lymberopoulos & Economou Reference Lymberopoulos and Economou1993).

2.2. Discharge characteristics

Figure 2 shows the cycle-averaged profiles of the electron density, electron temperature and electric field at different driving frequencies $6.78$, $13.56$ and 27.12 MHz. As one can see in figure 2(a), the cycle-averaged electron density is higher in the middle of the plasma discharge (bulk region). However, it is lower in the sheath regions close to the electrodes. The results also show that the electron density increases as the driving frequency increases. In the bulk plasma region, the electron density $n_{e}$ increases rapidly with the driving frequency and can reach $3.39 \times 10^{11}\ \textrm {cm}^{-3}$ at 27.12 MHz. However, this increase is less important near the electrodes. This is because with the growth of the driving frequency, more electrons are accelerated toward the opposite electrode during the sheath expansion, leading to high ionization close to the sheath edge.

Figure 2. The cycle-averaged distributions at different frequencies of (a) the electron density, (b) the electron temperature and (c) the electric field. Conditions are 1 Torr pressure and secondary electron emission coefficient $\gamma = 0.01$. The profiles show the discharge characteristics during the 5000th RF cycle.

The electron temperature is given in figure 2(b). It is shown that the profiles are peaked in the sheath regions near the electrodes. These peaks become stronger and narrower when the driving frequency increases from $6.78$ to 27.12 MHz. Hence, the sheath width decreases. This fact is due to the high electric field which also increases in the sheaths with the growth of driving frequency as shown in figure 2(c). Thus, more electrons will gain energy from the electric field in the sheaths to reach the ionization threshold energy which leads to a high electron production rate (Zhao et al. Reference Zhao, Li, Wang, Hu and Wang2016). The maximum electron temperature rises in the sheath regions and changes from $4.91$ eV at 6.78 MHz to 5.38 eV at 27.12 MHz. However, $T_{e}$ remains constant around 1.45 eV in the bulk plasma region. These results agree well with the earlier studies by Samir et al. (Reference Samir, Liu, Zhao and Zhou2017) with small difference in the electron temperature in the pre-sheath region, where a decrease in electron temperature is observed because of the cooling mechanism of electrons in this region. Figure 3 shows the ion mean energy as a function of the normalized distance between electrodes at different frequencies. In the sheath regions, it is shown that the ion mean energy is higher and increases with the driving frequency owing to the strong electric field in this regions (see figure 2c). In the bulk region, the ion mean energy is low and approximately 0.032 eV owing to the weak electric field in this region. Hence, the ion mean energy is mainly determined by the electric field. It should be noted that despite the differences in the operating conditions, the discharge characteristics obtained by the fluid model code show that the calculated data are in satisfactory agreement with the particle-in-cell calculations given by Benyoucef et al. (Reference Benyoucef, Yousfi, Belmadani and Settaouti2010).

Figure 3. The spatial distributions of the ion mean energy $E_{i}$ for three different frequencies $6.78$, $13.56$ and $27.12\ \textrm {MHz}$ in Ar capacitive discharge at 1 Torr.

3. Charging of the dust particle

A dust particle in a RF plasma acts as a small floating object and collects ions and electrons. However, for a spherical dust particle with radius $r_{d}$, the electron and positive ion currents at the dust surface can be calculated using orbital-motion-limited (OML) theory when the condition $r_{d}\ll \lambda _{L}$ applies, where $\lambda _{L}=[({1}/{\lambda _{De}})^{2} + ({1}/{\lambda _{Di}})^{2}]^{-1/2}$ is the linearized Debye length. These orbit-limited currents are given by (Akdim & Goedheer Reference Akdim and Goedheer2003)

(3.1)\begin{gather} I_{e} = 4{\rm \pi} r_{d}^{2}e n_{e} \sqrt{\dfrac{k_{B}T_{e}}{2 {\rm \pi} m_{e}}} \exp \left( \dfrac{e V_{{fl}}}{k_{B}T_{e}} \right), \end{gather}

(3.2)\begin{gather} I_{i}=4 {\rm \pi} r_{d}^{2} e n_{i} \sqrt{\frac{k_{B} T_{i}}{2 {\rm \pi} m_{i}}}\left(1-\frac{e V_{f l}}{k_{B} T_{i}}\right), \end{gather}

where $m_{i}$ and $m_{e}$ are the ion and electron masses, $T_{e}$ and $T_{i}$ are the ion and electron temperatures and $V_{{fl}}$ is the floating potential at the surface of the dust.

Recent works by Matthews et al. (Reference Matthews, Sanford, Kostadinova, Ashrafi, Guay and Hyde2020), Khrapak et al. (Reference Khrapak, Ratynskaia, Zobnin, Usachev, Yaroshenko, Thoma, Kretschmer, Höfner, Morfill and Petrov2005) and Land et al. (Reference Land, Matthews, Hyde and Bolser2010) have shown the importance of ion–neutral collisions for the charging of dust in plasma. In order to take into account the collision of ions with the neutral particles, a modified formula for ion current is used to replace (3.2) and is given by Land et al. (Reference Land, Matthews, Hyde and Bolser2010) based on the previous study by Khrapak et al. (Reference Khrapak, Ratynskaia, Zobnin, Usachev, Yaroshenko, Thoma, Kretschmer, Höfner, Morfill and Petrov2005):

(3.3)\begin{equation} I_{i}= 4 {\rm \pi} r_{d}^{2} e n_{i} \sqrt{\frac{E_{i}}{2 m_{i}}}\left[1-\frac{e V_{{fl}}}{E_{i}}+0.1\left(\frac{e V_{{fl}}}{E_{i}}\right)^{2} \frac{\lambda_{L}}{l_{m f p}}\right]. \end{equation}

The mean energy $E_{i}$ for ions is adopted to take into account the drift component of the ion velocity (Akdim & Goedheer Reference Akdim and Goedheer2003; Bleecker et al. Reference Bleecker, Bogaerts and Goedheer2006; Goedheer & Land Reference Goedheer and Land2008)

(3.4)\begin{equation} E_{i} = \dfrac{1}{2}m_{i} v_{s}^{2} = \dfrac{4k_{B}T_{g}}{{\rm \pi}} + \dfrac{1}{2} m_{i} v_{i}^{2}, \end{equation}

where the first part is the ion thermal energy and the second part is a kinetic energy from the drift velocity.

The ions current in (3.3) depends on the mean free path and the linearized Debye length around the dust particle because ions near the dust can lose their angular momentum when they collide with neutral particles. If the ions are close enough to the dust, they can be collected, which increases the ion current and reducing the negative floating potential. Making (3.1) and (3.3) equal, the floating potential $V_{{fl}}$ of the dust particle is obtained. Then, the dust charge is determined by using

(3.5)\begin{equation} Q_{d} = 4{\rm \pi}\epsilon_0 r_{d}V_{{fl}}. \end{equation}

Figure 4 represents the distribution of the dust charge as a function of the distance between electrodes at different driving frequencies and with particle radius $1.5\ \mathrm {\mu }\textrm {m}$. As we can see, the charge accumulated on the particle is larger and negative near the electrodes because of the high electron temperature in these two sheath regions. However, in the bulk plasma region, the charge is lower because of the weak electric field and low electron temperature. Moreover, the maximums of the dust charge change and move toward the electrodes when the driving frequency increases from 6.78 to 27.12 MHz in the sheath regions, because of the growth in the electric field with the driving frequency. However, the charge slightly increases in the bulk region.

Figure 4. The spatial distributions of the dust charge for three different frequencies $6.78$, $13.56$ and 27.12 MHz corresponding to $r_{d}= 1.5\ \mathrm {\mu }\textrm {m}$ in Ar capacitive discharge at 1 Torr.

4. Forces acting on the dust particle

In the present model, we consider the gravitational, electric, neutral drag and ion drag forces, which are expected to be important, whereas the thermophoretic force is neglected. The gravitational force is given by

(4.1)\begin{equation} F_{g} = \tfrac{4}{3}{\rm \pi} r_{d}^{3}\rho_{d}g, \end{equation}

where $g$ is the gravitational acceleration, $r_{d}$ is the dust particle radius and $\rho _{d}$ is the mass density of the dust particle. The electric force is given by

(4.2)\begin{equation} F_{e} = Q_{d}E, \end{equation}

where $E$ is the electric field and $Q_{d}$ is the dust particle electric charge.

The neutral drag force due to the collision of dust particle with neutral gas is taken to be of the following form (Epstein Reference Epstein1924; Nitter Reference Nitter1996)

(4.3)\begin{equation} F_{n} ={-}\frac{16}{3} \sqrt{{\rm \pi}} r_{d}^{2} n_{n} k T_{g} s\left(1+\frac{{\rm \pi}}{8}\right), \end{equation}

where $n_{n}$ is the density of the neutral gas atoms. Here $s = v_{d}/v_{t h, n}$, with $v_{d}$ is the dust velocity and $v_{t h, n}$ is the neutral thermal velocity. In our discharge plasma conditions with argon, $P = 1$ Torr, $m_{n}= 6.64\times 10^{-26}\ \textrm {kg}$ and $T_{g} = 300\ \textrm {K}$, the dust particle velocity is always much lower than the neutral thermal velocity ($v_{d}\ll v_{t h, n}$). Thus, in this case the Epstein's expression for the neutral drag force is appropriate (Nitter Reference Nitter1996).

The total ion drag force resulting from the momentum transfer from the ions following the dust particle is the sum of the collection force $F_{i}^{\textrm {coll}}$ due to the direct hitting by ions and the Coulomb force $F_{i}^{\textrm {coul}}$ due to Coulomb scattering of ions by the charged dust. Different approaches have been used to calculate the ion drag force (Khrapak et al. Reference Khrapak, Ivlev, Morfill and Thomas2002; Ivlev et al. Reference Ivlev, Khrapak, Zhdanov, Morfill and Joyce2004; Hutchinson Reference Hutchinson2006; Goedheer, Land & Venema Reference Goedheer, Land and Venema2009; Khrapak Reference Khrapak2013). In our case, we adopt the same expression as Goedheer et al. (Reference Goedheer, Land and Venema2009) who considered the collisions between ions and neutrals by incorporating the results of previous studies (Khrapak et al. Reference Khrapak, Ivlev, Morfill and Thomas2002; Ivlev et al. Reference Ivlev, Khrapak, Zhdanov, Morfill and Joyce2004; Hutchinson Reference Hutchinson2006). The total ion drag force is given by

(4.4)\begin{align} F_{i} &=F_{i}^{\textrm{coll}} + F_{i}^{\textrm{coul}} \nonumber\\ &=n_{i} m_{i} v_{s} {v}_{i} \times\left(\sigma_{c}(v_{s})+{\rm \pi} \rho_{0}^{2}(v_{s})\left[\varGamma(v_{s})+\mathcal{K}\left(\frac{\lambda_{L}(v_{s})}{l_{m f p}}\right)\right]\right), \end{align}

where

(4.5)\begin{equation} \mathcal{K}(x)=x \arctan (x)+\left(\sqrt{\frac{{\rm \pi}}{2}}-1\right) \frac{x^{2}}{1+x^{2}}-\sqrt{\frac{{\rm \pi}}{2}} \ln \left(1+x^{2}\right) \end{equation}

is a collisional function to take into account loss of angular momentum in collisions of ions and background neutrals (Ivlev et al. Reference Ivlev, Zhdanov, Khrapak and Morfill2005), $\sigma _{c}$ is the ion capture cross section given by

(4.6)\begin{equation} \sigma_{c}\left(v_{s}\right)={\rm \pi} r_{d}^{2}\left(1+\frac{\rho_{0}\left(v_{s}\right)}{r_{d}}\right) \end{equation}

and $\rho _{0}=Z e^{2} / 2 {\rm \pi} \epsilon _{0} m_{i} v_{s}^{2}$ is the Coulomb radius where $Z$ is the number of electron charges on the dust particle. To include ions scattered beyond the screening length, the expression of the Coulomb logarithm is given by (Khrapak et al. Reference Khrapak, Ivlev, Morfill and Thomas2002)

(4.7)\begin{equation} \varGamma = \ln \left[\frac{\rho_{0}(v_{s})+\lambda_{L}(v_{s})}{\rho_{0}(v_{s})+r_{d}}\right]. \end{equation}

However, this equation is not correct for larger ion drift velocity, hence, an effective ion velocity suggested by Hutchinson (Reference Hutchinson2005) can give a good agreement in the calculation of the Coulomb logarithm. In our case, we adopt the ion mean velocity using the same approach as Goedheer et al. (Reference Goedheer, Land and Venema2009) and Land et al. (Reference Land, Matthews, Hyde and Bolser2010), based on the fitting results from previous studies by Hutchinson (Reference Hutchinson2005)

(4.8)\begin{equation} v_{s}^{2}=\frac{8 k_{B} T_{g}}{{\rm \pi} m_{i}}+v_{i}^{2} \times\left[1+\left(\frac{v_{i}}{u_{b}} /\left(0.5+0.05 \ln \left(\tilde{m}_{i} / Z\right)+\left(T_{i} / T_{e}\right)^{1 / 2}\right)\right)^{3}\right], \end{equation}

where $\tilde {m}_{i}$ is the ion mass in amu ($\tilde {m}_{i} / Z= 40$ for argon ions) and $u_{b}=(k T_{e} / m_{i})^{1 / 2}$ is the Bohm velocity (Hutchinson Reference Hutchinson2005).

Having specified the forces, the equation of motion for the dust particle in the whole discharge can be written in the form

(4.9)\begin{equation} m_{d}\ddot{x}= F_{t}(x, \dot{x}) = {F}_{2}(x)+ {F}_{i}(x, \dot{x})+{F}_{n}(\dot{x}), \end{equation}

where $F_{t}$ represents the total force acting on the dust particle and $\dot {x}$ is the dust velocity. Here, $F_{2}$ combines the force of gravity and the electric force.

By linearizing around the equilibrium position $x=x_{0}$ and $\dot {x}=0$, equation (4.9) can be expressed as

(4.10)\begin{equation} m_{d} \ddot{x}_{1}=\left(\dfrac{\partial F_{2}}{\partial x}\right)_{x_{0}} x_{1} +\left(\dfrac{ \partial F_{i}}{\partial x}\right)_{x_{0}, \dot{x}=0} x_{1} +\left(\dfrac{\partial F_{i}}{\partial \dot{x}}\right)_{x_{0}, \dot{x}=0} \dot{x}_{1}+\left(\dfrac{\partial F_{n}}{\partial \dot{x}}\right)_{x_{0}, \dot{x}=0} \dot{x}_{1}, \end{equation}

where $x_{1}=x-x_{0}$. As the force of gravity remains constant, equation (4.10) can be written as an equation of the damped harmonic oscillator

(4.11)\begin{equation} m_{d} \ddot{x}_{1}=k_{2} x_{1}+k_{3} \dot{x}_{1}, \end{equation}

where $k_{2}=k_{2 e}+k_{2 i}$ and $k_{3}=k_{3 i}+k_{3 n}$ describe the coefficient of oscillation and the coefficient of damping caused by the neutral drag and the ion drag forces, respectively. Here, $k_{2 e}$ is the coefficient describing the contribution of the electric force, $k_{2 i}$ and $k_{3 i}$ combine the contribution of the collection and the Coulomb parts of the ion drag force as

(4.12a,b)\begin{equation} k_{2 i} = k_{2i, \textrm{coll}} + k_{2 i, \textrm{coul}} ,\quad k_{3 i} = k_{3 i, \textrm{coll}} + k_{3 i, \textrm{coul}}. \end{equation}

where the coefficients $k_{2 e}$, $k_{2 i, \textrm {coll}}$ and $k_{2 i, \textrm {coul}}$ are calculated from (4.2) and (4.4) and given by

(4.13)\begin{align} k_{2 e}&={-}4 {\rm \pi} \varepsilon_{0} r_{d} \left[\dfrac{\partial V_{{fl}}}{\partial x} \dfrac{\partial V}{\partial x}+V_{{fl}}\left(\dfrac{\partial^{2} V}{\partial x^{2}}\right)\right], \end{align}

(4.14)\begin{align} k_{2 i, \textrm{coll}} &= m_{i} n_{i} v_{i} \left[\left(\dfrac{1}{v_{i}} \dfrac{\partial v_{i}}{\partial x}+\dfrac{1}{n_{i}} \dfrac{\partial n_{i}}{\partial x}\right) v_{s} \sigma_{c}+\sigma_{c} \dfrac{\partial v_{s}}{\partial x}+\dfrac{2 {\rm \pi} e r_{d}}{m_{i} v_{s}^{2}}\right.\nonumber\\ &\quad \times\left.\left(v_{s} \dfrac{\partial V_{f l}}{\partial x}-2 V_{f l} \dfrac{\partial v_{s}}{\partial x}\right)\right] \end{align}

(4.15)\begin{align} k_{2 i, \textrm{coul}} &= {\rm \pi} m_{i} n_{i} v_{i} \rho_{0}^{2} \left[\left(\varGamma+\mathcal{K}(\lambda_{n})\right)\left(\left(\dfrac{1}{v_{i}} \dfrac{\partial v_{i}}{\partial x}+\dfrac{1}{n_{i}} \dfrac{\partial n_{i}}{\partial x}\right) v_{s} + \dfrac{\partial v_{s}}{\partial x}+\dfrac{2}{V_{f l}} \right)\right. \nonumber\\ &\quad +\left.\dfrac{v_{s}}{\rho_{0}+r_{d}} \left(\dfrac{e r_{d}}{m_{i} v_{s}^{3}}\left(\dfrac{r_{d}-\lambda_{L}}{\rho_{0}+r_{d}} \right) \left(v_{s} \dfrac{\partial V_{f l}}{\partial x}-2 V_{f l} \dfrac{\partial v_{s}}{\partial x}\right)+\dfrac{\partial \lambda_{{L}}}{\partial x}\right)\right. \nonumber\\ &\quad \left.+\,v_{s}\left(\dfrac{\partial \lambda_{n}}{\partial x} \arctan \left(\lambda_{n}\right) - \dfrac{\lambda_{n}}{1+\lambda_{n}^{2}}-2\left(\sqrt{\dfrac{{\rm \pi}}{2}}-1\right) \dfrac{\lambda_{n}^{3}}{\left(1+\lambda_{n}^{2}\right)^{2}}\right) \right], \end{align}

respectively, where $\lambda _{n} = \lambda _{L}(v_{s})/l_{m f p}$. Moreover, the coefficients $k_{3 i, \textrm {coll}}$, $k_{3 i, \textrm {coul}}$ and $k_{3 n}$ are calculated from (4.3) and (4.4) and given by

(4.16)\begin{align} k_{3 i, \textrm{coll}} &={-} m_{i} n_{i}{\left[v_{s} \sigma_{c}+v_{i}^{2}\left(5\left(\frac{v_{i}}{u_{b}} /\left(0.5+0.05 \ln \left(\tilde{m}_{i} / Z_{i}\right)+\left(T_{i} / T_{e}\right)^{1 / 2}\right)\right)^{3}+2\right)\right.} \nonumber\\ &\quad \left.\times\left(\dfrac{\sigma_{c}}{2 v_{s}}-\dfrac{2 {\rm \pi} e r_{d} V_{f l}}{m_{i} v_{s}^{3}}\right)\right] \end{align}

(4.17)\begin{align} k_{3 i, \textrm{coul}} &= {\rm \pi} m_{i} n_{i} \rho_{0}^{2} v_{i}^{2}\left[\vphantom{\frac{v_{s}}{v_{i}^{2}}}\left(\varGamma+\mathcal{K}(\lambda_{n})\right)\right.\nonumber\\ &\quad \times\left(\frac{v_{s}}{v_{i}^{2}}-\frac{3}{2v_{s}}\left(5\left(\frac{v_{i}}{u_{b}} /\left(0.5+0.05 \ln \left(\tilde{m}_{i} / Z_{i}\right)+\left(T_{i} / T_{e}\right)^{1 / 2}\right)\right)^{3}+2\right)\right)\nonumber\\ &\quad \left.-\frac{v_{s}}{v_{i}}\left(\varGamma^{\prime \prime}+\mathcal{K}^{\prime \prime}\left(\lambda_{n}\right)\right)\right], \end{align}

respectively, where $\varGamma ^{\prime \prime }$ and $\mathcal {K}^{\prime \prime }(\lambda _{n})$ denote the derivatives with respect to dust velocity of the Coulomb logarithm and of the collisional function, respectively,

(4.18)\begin{align} &\varGamma^{\prime \prime} = \frac{1}{\left(\rho_{0}+\lambda_{L}\right)\left(\rho_{0}+r_{d}\right)}\left[\vphantom{\frac{\lambda_{L}^{3} v_{i}}{\lambda_{D_{i}}^{2} v_{t h, i}^{2}\left(1+\dfrac{v_{i}^{2}}{v_{t h,i}^{2}}\right)}}\left(r_{d}-\lambda_{L}\right)\right.\nonumber\\ &\qquad \times \left(\frac{e r_{d} V_{f l} }{m_{i} v_{s}^{4}} v_{i}\left(5\left(\frac{v_{i}}{u_{b}} /\left(0.5+0.05 \ln \left(\tilde{m}_{i} / Z_{i}\right)+\left(T_{i} / T_{e}\right)^{1 / 2}\right)\right)^{3}+2\right)\right)\nonumber\\ &\qquad \left. - \frac{\lambda_{L}^{3} v_{i}}{\lambda_{D_{i}}^{2} v_{t h, i}^{2}\left(1+\dfrac{v_{i}^{2}}{v_{t h,i}^{2}}\right)}\right], \end{align}

(4.19)\begin{align} \mathcal{K}^{\prime \prime} (\lambda_{n}) &={-}\dfrac{\lambda_{L}^{3} v_{i}}{l_{mfp} \lambda_{D i}^{2} v_{t h,i}^{2}\left(1+\dfrac{v_{i}^{2}}{v_{t h,i}^{2}}\right)}\nonumber\\ &\quad \times\left[\arctan \left(\lambda_{n}\right) - \dfrac{\lambda_{n}}{1+\lambda_{n}^{2}}-2\left(\sqrt{\dfrac{{\rm \pi}}{2}}-1\right) \dfrac{\lambda_{n}^{3}}{\left(1+\lambda_{n}^{2}\right)^{2}}\right] \end{align}

and

(4.20)\begin{equation} k_{3 n} = \left. -\frac{16}{3} \sqrt{{\rm \pi}} r_{d}^{2} n_{n} k_{B} T_{g} \left(1+\frac{{\rm \pi}}{8}\right)\right/v_{\textrm{th},n}. \end{equation}

From (4.11), it can be shown that the dust particle acts as a damped harmonic oscillator. Hence, the damping time $\tau _{d}$, the oscillation frequency $\omega$ and the oscillation time $\tau$ are defined as

(4.21)\begin{gather} \tau_{d}=\dfrac{2 m_{d}}{\left|k_{3}\right|}, \end{gather}

(4.22)\begin{gather} \omega=\left[-\left(\dfrac{k_{2}}{m_{d} }\right)_{x_{0}}\right]^{1 / 2}, \end{gather}

(4.23)\begin{gather} \tau= \dfrac{2{\rm \pi}}{\omega}, \end{gather}

respectively. Calculation of the oscillation and damping parameters described previously yields data on the dust motion in the discharge.

In figure 5, we show the oscillation time in seconds and the damping time compared with the oscillation time (damping time/oscillation time) for the ion collection drag, ion Coulomb drag and neutral drag. The results are obtained for a micrometre-sized dust particle as a function of the normalized inter-electrode distance at three different frequencies $6.78$, $13.56$ and 27.12 MHz. As we can see, there are some critical points in the discharge domain where the pulsation $\omega$ tends to be zero when the dust particle approaches the critical point corresponding to zero oscillation coefficient ($k_{2} =0$), which means that the oscillation time depends on the position in the discharge. In addition, the damping effect due to the neutral drag is more important than that due to ion drag, which indicates that the damping is mainly caused by the neutral drag force. When the driving frequency increases from $6.78$ to 27.12 MHz, the oscillation time decreases in the sheath regions from $0.021$ to $8.1\times$ $10^{-3}\ \textrm {s}$, which means that the oscillation frequency becomes important and the particle oscillates rapidly in these regions.

Figure 5. The oscillation time in seconds and the damping time compared with the oscillation time due to the ion collection drag (labelled ion collection), ion Coulomb (labelled ion Coulomb) and neutral labelled (neutral drag) for a micrometre-sized $r_{d}=1.5\ \mathrm {\mu }\textrm {m}$ dust particle versus inter-electrode distance at different frequencies (a) 6.78 MHz, (b) 13.56 MHz and (c) 27.12 MHz. Here $\rho _{d} = 2000\ \textrm {kg}\,\textrm {m}^{-3}$. The vertical dashed lines indicate the position of the critical points in the discharge.

5. Dust dynamics

Note that we have determined the plasma discharge characteristics in § 2 and the charge and the forces acting on the dust in §§ 3 and 4, we now turn to the dust particle to examine its motion in plasma discharge under the influence of driving frequency by using (4.11).

Figure 6 shows the velocities of a dust particle with $r_{d} = 1.5\ \mathrm {\mu }\textrm {m}$ released from the lower electrode given as a function of the position and time for different frequencies. It is seen that this dust particle follows a damped oscillatory motion around its equilibrium position with higher oscillation frequency when the driving frequency increases, because of the competition between the parameters describing the amplification and the damping of the dust oscillations. In our conditions and for all frequencies, the parameter $k_{3i,\textrm {coul}}$ is positive around the dust equilibrium position, which means that the contribution of the Coulomb part of the ion drag force is to increase the amplitude of oscillation. However, the parameter $k_{3i,\textrm {coll}}$ is negative, i.e. the contribution of the collection part of the ion drag force is in the damping of dust oscillation in this position. Meanwhile, $k_{3n}$ is negative, which means that the neutral drag force also acts as a damping force on the dust motion. On the other hand, the location of the dust particle changes and becomes closer to the powered electrode with the growth of the driving frequency because of the reduction in the sheath width and the height from the electrode is reduced from $3.25$ to 1.72 mm when the frequency increases from $6.78$ to 27.12 MHz. In addition, the oscillation frequency increases because the oscillation time is small near the electrodes compared with that away from electrodes which can be clearly seen in figure 5.

Figure 6. The spatiotemporal distributions of a dust particle velocity with $r_{d}=1.5\ \mathrm {\mu }\textrm {m}$. The dust particle is released from the lower electrode with zero initial velocity $v_{d} = 0\ \textrm {m}\,\textrm {s}^{-1}$ at (a) $6.78\ \textrm {MHz}$, (b) $13.56\ \textrm {MHz}$ and (c) $27.12\ \textrm {MHz}$, with $\rho _{d} = 2000\ \textrm {kg}\,\textrm {m}^{-3}$.

Note that in figure 6, we considered a dust particle with $r_{d} = 1.5\ \mathrm {\mu }\textrm {m}$. To visualize the effect of the dust particle radius on the movement and the equilibrium position in the discharge as a function of the driving frequency, we inject into the discharge an other dust particle with small radius $r_{d} = 1.1\ \mathrm {\mu }\textrm {m}$. Figure 7 gives the velocity as a function of position and time for a dust particle released from the lower electrode without initial velocity, these results are given for three different frequencies $6.78$, $13.56$ and 27.12 MHz. In all cases, the figures 7(a), 7(b) and 7(c) show that this dust particle also follows a damped oscillatory motion around its equilibrium position with different oscillation frequencies. Moreover, the height from the electrode is reduced from $3.40$ to 1.75 mm when the driving frequency increases from $6.78$ to 27.12 MHz, which means that this small particle with $r_{d} = 1.1\ \mathrm {\mu }\textrm {m}$ levitates relatively away from the electrode with lower oscillation frequency compared with the large dust particle with $r_{d} = 1.5\ \mathrm {\mu }\textrm {m}$ (see figure 6). These results agree with the experimental data given by Tomme et al. (Reference Tomme, Annaratone and Allen2000) which show that the particle with a large size levitates closer to the electrode than the particle with a small radius.

Figure 7. The spatiotemporal distributions of a dust particle velocity with $r_{d}=1.1\ \mathrm {\mu }\textrm {m}$. The dust particle is released from the lower electrode with zero initial velocity $v_{d} = 0\ \textrm {m}\,\textrm {s}^{-1}$ at (a) $6.78\ \textrm {MHz}$, (b) $13.56\ \textrm {MHz}$ and (c) $27.12\ \textrm {MHz}$, with $\rho _{d} = 2000\ \textrm {kg}\,\textrm {m}^{-3}$.

Many parameters can influence the possibility of suspension of the injected dust particle, for example, plasma parameters, initial particle position and velocity, particle shape and mass density. Figure 8 illustrates the position as a function of time for a dust particle with $r_{d} = 1.1\ \mathrm {\mu }\textrm {m}$ released from two different positions (from the middle and from the upper electrode) without initial velocity. In all cases, when the release is from the upper electrode, the dust particle can be levitated near the sheath edge with higher oscillation frequency when the driving frequency increases. On the other hand, when the release is from the middle of the discharge, the dust particle oscillates near the sheath region of the powered electrode at $6.78\ \textrm {MHz}$, while its movement ends by hitting the powered electrode at $13.56$ and $27.12\ \textrm {MHz}$, because the dust particle accelerates on its trajectory from the bulk region, and the amplitude of the oscillation initially increases. In figure 9, the trajectories of a micrometre-sized dust particle with $r_{d} = 1.1\ \mathrm {\mu }\textrm {m}$ are shown for three frequencies $6.78$, $13.56$ and $27.12\ \textrm {MHz}$. The particle is injected into the plasma from the middle of the discharge (central bulk plasma region) with two initial velocities $v_{d} = 0.25$ and $v_{d} = -0.25\ \textrm {m}\,\textrm {s}^{-1}$. As we can see, the dust particle performs a damped oscillation around its equilibrium position close to the grounded sheath edge when the initial velocity is positive. In this case, the dust rises and can reach 9.2 mm from its initial position at $6.78\ \textrm {MHz}$, whereas at $13.56\ \textrm {MHz}$, the dust particle hits the grounded electrode. The same scenario is observed at $27.12\ \textrm {MHz}$ with rapid hitting even if we give the same positive initial velocity to the dust particle. On the other hand, the dust particle at $6.78\ \textrm {MHz}$ can also be suspended near the powered sheath edge despite the negative initial velocity, whereas the dust particle hits rapidly the powered electrode at $13.56$ and $27.12\ \textrm {MHz}$, because the electric force becomes more important with the growth in the driving frequency and it is anti-parallel to the dust velocity. Note that if we examine the motion of a large dust particle with $r_{d} = 1.5\ \mathrm {\mu }\textrm {m}$ in the same initial conditions taken in figure 9, the particle with positive initial velocity at $6.78\ \textrm {MHz}$ rises but with a low height from the initial position compared with the small particle, whereas with negative initial velocity, the dust particle hits the powered electrode due to the increase in the gravitational force. On the other hand, if the driving frequency increases, this dust particle cannot rise in the sheath regions and rapidly hits the electrodes because the equilibrium position becomes closer to the electrodes, and the amplitude of oscillation initially becomes higher.

Figure 8. The trajectory of the dust particle with $r_{d}=1.1\ \mathrm {\mu }\textrm {m}$ and $v_{d} = 0\ \textrm {m}\,\textrm {s}^{-1}$ released from the middle of the discharge and from the upper electrode at (a) $6.78\ \textrm {MHz}$, (b) $13.56\ \textrm {MHz}$ and (c) $27.12\ \textrm {MHz}$, with $\rho _{d} = 2000\ \textrm {kg}\,\textrm {m}^{-3}$.

Figure 9. The trajectory of the dust particle with $r_{d}=1.1\ \mathrm {\mu }\textrm {m}$ released from the middle of the discharge for different initial velocities $v_{d} = 0.25\ \textrm {m}\,\textrm {s}^{-1}$ and $v_{d} = -0.25\ \textrm {m}\,\textrm {s}^{-1}$ at (a) $6.78\ \textrm {MHz}$, (b) $13.56\ \textrm {MHz}$ and (c) $27.12\ \textrm {MHz}$, with $\rho _{d} = 2000\ \textrm {kg}\,\textrm {m}^{-3}$.