Distribution of electrons, ions, and fine (dust) particles in cylindrical fine particle (dusty) plasmas: drift-diffusion analysis

Our system is composed of electrons, ions and fine particles (simply 'particles' in what follows). We denote the physical quantities of these components by suffixes $e,~i$, and α, respectively, α being the species of the particles. For densities ${{n}_{e,~i,~\alpha}}$ and flux densities ${{\mathbf{\Gamma}}_{e,~i,~\alpha}}$, we have the equations of continuity in the general form,

$$\begin{eqnarray}&&\frac{\partial {{n}_{e}}}{\partial t}+\nabla \centerdot {{\mathbf{\Gamma}}_{e}}=\frac{\delta {{n}_{e}}}{\delta t},\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&\frac{\partial {{n}_{i}}}{\partial t}+\nabla \centerdot {{\mathbf{\Gamma}}_{i}}=\frac{\delta {{n}_{i}}}{\delta t},\end{eqnarray} \tag{ 2 }$$

and

$$\begin{eqnarray}&&\frac{\partial {{n}_{\alpha}}}{\partial t}+\nabla \centerdot {{\mathbf{\Gamma}}_{\alpha}}=0.\end{eqnarray} \tag{ 3 }$$

On the right-hand sides, $\delta {{n}_{e}}/\delta t$ and $\delta {{n}_{i}}/\delta t$ express the generation of electrons and ions by the ionisation and their loss from the system, respectively. In this article, we consider only the stationary case where $\partial {{n}_{e,~i,~\alpha}}/\partial t=0$.

We assume that flux densities are expressed by diffusion coefficients ${{D}_{e,~i,~\alpha}}$ and mobilities ${{\mu}_{e,~i,~\alpha}}$ in the forms,

$$\begin{eqnarray}&&{{\mathbf{\Gamma}}_{e}}=-{{D}_{e}}\nabla {{n}_{e}}-{{n}_{e}}{{\mu}_{e}}\left(\frac{{{\mathbf{F}}_{e}}}{-e}\right),\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&{{\mathbf{\Gamma}}_{i}}=-{{D}_{i}}\nabla {{n}_{i}}+{{n}_{i}}{{\mu}_{i}}\left(\frac{{{\mathbf{F}}_{i}}}{e}\right),\end{eqnarray} \tag{ 5 }$$

and

$$\begin{eqnarray}&&{{\mathbf{\Gamma}}_{\alpha}}=-{{D}_{\alpha}}\nabla {{n}_{\alpha}}-{{n}_{\alpha}}{{\mu}_{\alpha}}\left(\frac{{{\mathbf{F}}_{\alpha}}}{-{{Q}_{\alpha}}e}\right).\end{eqnarray} \tag{ 6 }$$

Here  −e, e, $-{{Q}_{\alpha}}e~\left({{Q}_{\alpha}}>0\right)$ and ${{\mathbf{F}}_{e,~i,~\alpha}}$ are the charges and the forces for an electron, an ion and a particle of species α, respectively. We also assume that the diffusion coefficients and mobilities are related by Einstein relations as

$$\begin{eqnarray}&&\frac{{{D}_{e}}}{{{\mu}_{e}}}=\frac{{{k}_{\text{B}}}{{T}_{e}}}{e},~~~\frac{{{D}_{i}}}{{{\mu}_{i}}}=\frac{{{k}_{\text{B}}}{{T}_{i}}}{e},~~~\frac{{{D}_{\alpha}}}{{{\mu}_{\alpha}}}=\frac{{{k}_{\text{B}}}{{T}_{\alpha}}}{{{Q}_{\alpha}}e}\end{eqnarray} \tag{ 7 }$$

with the temperature of the corresponding components, ${{T}_{e,~i,~\alpha}}$. Note that, since we have no generation/loss of particles, the stationary state of our system depends, within the drift-diffusion approximation, only on ${{D}_{\alpha}}/{{\mu}_{\alpha}}$ and we need not specify ${{D}_{\alpha}}$ nor ${{\mu}_{\alpha}}$ when ${{T}_{\alpha}}$ is given.

The electric field $\mathbf{E}$ is related to the net charge density by the Poisson's equation

$$\begin{eqnarray}&&{{\varepsilon}_{0}}\nabla \centerdot \mathbf{E}=e\left({{n}_{i}}-{{n}_{e}}-\underset{\alpha}{\mathop \sum}\,{{Q}_{\alpha}}{{n}_{\alpha}}\right).\end{eqnarray} \tag{ 8 }$$

The force acting on a particle of species α is given by

$$\begin{eqnarray}&&{{\mathbf{F}}_{\alpha}}=\left(-{{Q}_{\alpha}}e\right)\mathbf{E}+\mathbf{F}_{\alpha}^{id},\end{eqnarray} \tag{ 9 }$$

where the first term is due to the electric field $\mathbf{E}$ and the second term is the ion drag force. The force on an ion is given by

$$\begin{eqnarray}&&{{\mathbf{F}}_{i}}=e\mathbf{E}-\underset{\alpha}{\mathop \sum}\,\frac{{{n}_{\alpha}}}{{{n}_{i}}}\mathbf{F}_{\alpha}^{id},\end{eqnarray} \tag{ 10 }$$

where the second term is the reaction of $\mathbf{F}_{\alpha}^{id}$ in (9). The force on an electron is written as

$$\begin{eqnarray}&&{{\mathbf{F}}_{e}}=-e\mathbf{E}.\end{eqnarray} \tag{ 11 }$$

The ion drag force has been investigated extensively for usual positive ions, as in our system [2, 10–12] and also for negative ions [13]. It is along the direction of the average ion flow ${{\mathbf{u}}_{i}}$ and here we estimate the magnitude by the expression [10] (see appendix A)

$$\begin{eqnarray}&&F_{\alpha}^{id}\sim {{c}_{1}}\frac{{{\left({{Q}_{\alpha}}e\right)}^{2}}}{4\pi {{\varepsilon}_{0}}}\frac{{{n}_{i}}{{e}^{2}}}{{{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{i}}}\frac{{{u}_{i}}}{{{v}_{i}}}~~~~~\left(\text{when}~\frac{{{u}_{i}}}{{{v}_{i}}}<1\right),\end{eqnarray} \tag{ 12 }$$

where

$$\begin{eqnarray}&&{{c}_{1}}=\frac{1}{3}{{\left(\frac{2}{\pi}\right)}^{1/2}}\Lambda\end{eqnarray} \tag{ 13 }$$

and $\Lambda$ is the generalised Coulomb logarithm

$$\begin{eqnarray}&&\Lambda\left({{\beta}_{T}},{{r}_{\alpha}}/\lambda \right)=-{{e}^{{{\beta}_{T}}/2}}\text{Ei}\left(-{{\beta}_{T}}/2\right)+{{e}^{\left[\left({{\beta}_{T}}/2\right)\left(\lambda /{{r}_{\alpha}}\right)\right]}}\text{Ei}\left[-\left({{\beta}_{T}}/2\right)\left(\lambda /{{r}_{\alpha}}\right)\right],\end{eqnarray} \tag{ 14 }$$

$\text{Ei}(x)$ being the exponential integral (${{\beta}_{T}}\sim \left[|{{Q}_{\alpha}}|{{e}^{2}}/4\pi {{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{i}}\right]/$${{\left({{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{i}}/{{e}^{2}}{{n}_{i}}\right)}^{1/2}}$ as defined in appendix A). The thermal velocity of the ions is denoted by v_i with the definition ${{v}_{i}}={{\left({{k}_{\text{B}}}{{T}_{i}}/{{m}_{i}}\right)}^{1/2}}$ (m_i is the ion mass) and we have ${{u}_{i}}<{{v}_{i}}$ in all cases shown below: the particles are distributed near the symmetry axis and they do not approach the domain near the wall, where we have the opposite condition ${{u}_{i}}>{{v}_{i}}$. We also note that the average ion flux ${{n}_{i}}{{\mathbf{u}}_{i}}$ is given by

$$\begin{eqnarray}&&{{n}_{i}}{{\mathbf{u}}_{i}}=-{{D}_{i}}\left[\nabla {{n}_{i}}-\frac{e{{n}_{i}}}{{{k}_{\text{B}}}{{T}_{i}}}\left(\mathbf{E}-\underset{\alpha}{\mathop \sum}\,\frac{{{n}_{\alpha}}}{{{n}_{i}}}\frac{\mathbf{F}_{\alpha}^{id}}{e}\right)\right].\end{eqnarray} \tag{ 15 }$$

We assume that the axial electric field maintains the discharge. This field induces a discharge current composed mainly of electrons, and they eventually have the temperature of the order of a few electron volts. These electrons are considered to generate the plasma (electrons and ions) by the impact ionisation, and the generation rate is proportional to the neutral atom density and the electron density; the proportionality constant is dependent on the electron temperature [14]. We regard factors other than the electron density as uniform, and express the rate per volume by ${{c}_{g}}{{n}_{e}}$ with the coefficient c_g, which is constant in the system (we assume uniform neutral gas density).

We assume that the plasma is lost by recombination after diffusing onto the outer boundary of the system (wall) or the surface of fine particles, neglecting the recombination in the bulk of the plasma. The charge of the particle α, $-{{Q}_{\alpha}}e$, is determined by the balance between the ion current ${{I}_{i,~\alpha}}$ and the electron current ${{I}_{e,~\alpha}}$ onto the surface of particle α,

$$\begin{eqnarray}&&{{I}_{i,~\alpha}}={{I}_{e,~\alpha}},\end{eqnarray} \tag{ 16 }$$

or explicitly,

$$\begin{eqnarray}&&{{Q}_{\alpha}}e=z\frac{{{k}_{\text{B}}}{{T}_{e}}}{e/4\pi {{\varepsilon}_{0}}{{r}_{\alpha}}},\end{eqnarray} \tag{ 17 }$$

$$\begin{eqnarray}&&{{I}_{i,~\alpha}}={{(8\pi )}^{1/2}}r_{\alpha}^{2}{{n}_{i}}{{v}_{i}}\left[1+\frac{{{T}_{e}}}{{{T}_{i}}}z+0.1{{\left(\frac{{{T}_{e}}}{{{T}_{i}}}\right)}^{2}}{{z}^{2}}\frac{\lambda}{{{\ell}_{i}}}\right],\end{eqnarray} \tag{ 18 }$$

and

$$\begin{eqnarray}&&{{I}_{e,~\alpha}}={{(8\pi )}^{1/2}}r_{\alpha}^{2}{{n}_{e}}{{v}_{e}}\exp (-z),\end{eqnarray} \tag{ 19 }$$

where ${{r}_{\alpha}}$ is the particle radius and v_e is the the electron thermal velocity defined by ${{v}_{e}}={{\left({{k}_{\text{B}}}{{T}_{e}}/{{m}_{e}}\right)}^{1/2}}$ (m_e is the electron mass). In (18), the effect of ion-neutral collisions is expressed by the third term in square brackets (| |) [15] where ${{\ell}_{i}}$ is the ion mean free path and λ is the screening length defined by

$$\begin{eqnarray}&&\frac{1}{{{\lambda}^{2}}}=k_{i}^{2}+k_{e}^{2},\end{eqnarray} \tag{ 20 }$$

k_i and k_e being the ion and electron Debye wave numbers, respectively,

$$\begin{eqnarray}&&k_{i}^{2}=\frac{{{e}^{2}}{{n}_{i}}}{{{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{i}}},~~~k_{e}^{2}=\frac{{{e}^{2}}{{n}_{e}}}{{{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{e}}}.\end{eqnarray} \tag{ 21 }$$

Since we assume ${{T}_{e}}\gg {{T}_{i}}$, $\lambda \sim {{\left({{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{i}}/{{e}^{2}}{{n}_{i}}\right)}^{1/2}}$. Denoting the balanced current ${{I}_{i,~\alpha}}={{I}_{e,~\alpha}}$ by ${{c}_{\alpha}}$, we have $\delta {{n}_{e}}/\delta t$ and $\delta {{n}_{i}}/\delta t$ in (1) and (2) as

$$\begin{eqnarray}&&\frac{\delta {{n}_{i}}}{\delta t}=\frac{\delta {{n}_{e}}}{\delta t}={{c}_{g}}{{n}_{e}}-\underset{\alpha}{\mathop \sum}\,{{c}_{\alpha}}{{n}_{\alpha}}.\end{eqnarray} \tag{ 22 }$$

We assume that our system is stationary and cylindrically symmetric. In this case, $\mathbf{E}$, $\mathbf{F}_{\alpha}^{id}$, and ${{\mathbf{u}}_{i}}$ has only the radial component E, $F_{\alpha}^{id}$, and u_i, respectively, and we define the ratio of the magnitudes of the ion drag force to the electric field force for a particle α, ${{f}_{\alpha}}$, by

$$\begin{eqnarray}&&{{f}_{\alpha}}=\frac{F_{\alpha}^{id}}{{{Q}_{\alpha}}eE}.\end{eqnarray} \tag{ 23 }$$

In terms of s  =  R² where R is the radius, equations (1)–(6) and (8) reduce to

$$\begin{eqnarray}&&2{{D}_{e}}\frac{\text{d}}{\text{d}s}s\left[-~2\frac{\text{d}}{\text{d}s}{{n}_{e}}-\frac{eE}{{{s}^{1/2}}{{k}_{\text{B}}}{{T}_{e}}}{{n}_{e}}\right]={{c}_{g}}\left({{n}_{e}}-\underset{\alpha}{\mathop \sum}\,\frac{{{c}_{\alpha}}}{{{c}_{g}}}{{n}_{\alpha}}\right),\end{eqnarray} \tag{ 24 }$$

$$\begin{eqnarray}&&2{{D}_{i}}\frac{\text{d}}{\text{d}s}s\left[-~2\frac{\text{d}}{\text{d}s}{{n}_{i}}+\left(1-\underset{\alpha}{\mathop \sum}\,{{f}_{\alpha}}\frac{{{Q}_{\alpha}}{{n}_{\alpha}}}{{{n}_{i}}}\right)\frac{eE}{{{s}^{1/2}}{{k}_{\text{B}}}{{T}_{i}}}{{n}_{i}}\right]={{c}_{g}}\left({{n}_{e}}-\underset{\alpha}{\mathop \sum}\,\frac{{{c}_{\alpha}}}{{{c}_{g}}}{{n}_{\alpha}}\right),\end{eqnarray} \tag{ 25 }$$

$$\begin{eqnarray}&&2{{D}_{p}}\frac{\text{d}}{\text{d}s}s\left[-~2\frac{\text{d}}{\text{d}s}{{n}_{\alpha}}-\left(1-{{f}_{\alpha}}\right){{Q}_{\alpha}}\frac{eE}{{{s}^{1/2}}{{k}_{\text{B}}}{{T}_{p}}}{{n}_{\alpha}}\right]=0,\end{eqnarray} \tag{ 26 }$$

and

$$\begin{eqnarray}&&2{{\varepsilon}_{0}}\frac{\text{d}}{\text{d}s}{{s}^{1/2}}E=e\left({{n}_{i}}-{{n}_{e}}-\underset{\alpha}{\mathop \sum}\,{{Q}_{\alpha}}{{n}_{\alpha}}\right).\end{eqnarray} \tag{ 27 }$$

These equations, together with (16) and (17), determine n_e, n_i, ${{n}_{\alpha}}$, ${{Q}_{\alpha}}$, and E. Equation (15) reduces to

$$\begin{eqnarray}&&2\frac{\text{d}}{\text{d}s}{{s}^{1/2}}{{n}_{i}}{{u}_{i}}={{c}_{g}}\left({{n}_{e}}-\underset{\alpha}{\mathop \sum}\,\frac{{{c}_{\alpha}}}{{{c}_{g}}}{{n}_{\alpha}}\right).\end{eqnarray} \tag{ 28 }$$

Noting that the densities and E/s^1/2 are finite at s  =  0, we integrate both sides of these equations with respect to s and introduce ${{N}_{e,~i,~\alpha}}$ defined, respectively, by

$$\begin{eqnarray}&&{{N}_{e,~i,~\alpha}}(s)\equiv \int_{0}^{s}\text{d}s{{n}_{e,~i,~\alpha}}(s)\end{eqnarray} \tag{ 29 }$$

to have the equations as shown in appendix B, which are used in the numerical analysis.

Let us first consider the case where the effect of the ion drag force (expressed by ${{f}_{\alpha}}$) can be neglected and K₁ reduces to $\left({{T}_{e}}/{{n}_{e}}\right)\left(Q_{\alpha}^{2}{{n}_{\alpha}}/{{T}_{\alpha}}+{{n}_{i}}/{{T}_{i}}+{{n}_{e}}/{{T}_{e}}\right)$. We then have

$$\begin{eqnarray}&&\frac{{{a}_{1}}}{A}\sim {{\left[\frac{\left({{T}_{e}}/{{T}_{i}}\right)\left(1-{{c}_{\alpha}}{{n}_{\alpha}}/{{c}_{g}}{{n}_{e}}\right)}{R_{a}^{2}\left({{e}^{2}}/{{\varepsilon}_{0}}{{k}_{\text{B}}}\right)\left(Q_{\alpha}^{2}{{n}_{\alpha}}/{{T}_{\alpha}}+{{n}_{i}}/{{T}_{i}}+{{n}_{e}}/{{T}_{e}}\right)}\right]}_{0}}\\ &&\sim {{\left[\frac{\left({{T}_{e}}/{{T}_{i}}\right)\left(1-{{c}_{\alpha}}{{n}_{\alpha}}/{{c}_{g}}{{n}_{e}}\right)}{R_{a}^{2}\left({{e}^{2}}/{{\varepsilon}_{0}}{{k}_{\text{B}}}\right)\left(Q_{\alpha}^{2}{{n}_{\alpha}}/{{T}_{\alpha}}+{{n}_{i}}/{{T}_{i}}\right)}\right]}_{0}}.\end{eqnarray} \tag{ 53 }$$

Here, in the denominator, fine particles also contribute to the factor

$$\begin{eqnarray}&&\left(Q_{\alpha}^{2}{{n}_{\alpha}}/{{T}_{\alpha}}+{{n}_{i}}/{{T}_{i}}+{{n}_{e}}/{{T}_{e}}\right)\sim \left(Q_{\alpha}^{2}{{n}_{\alpha}}/{{T}_{\alpha}}+{{n}_{i}}/{{T}_{i}}\right)\end{eqnarray} \tag{ 54 }$$

in contrast to the case without particles (36), where the factor reduces to ${{n}_{i}}/{{T}_{i}}$ (equation (53) includes the latter case with ${{n}_{\alpha}}=0$).

We assume that the charge density of the fine particles is less than a few tens % of the total negative charges. However, since ${{Q}_{\alpha}}$ is usually as large as 10³–10⁴, we expect the inequality holds (at R  =  0)

$$\begin{eqnarray}&&Q_{\alpha}^{2}{{n}_{\alpha}}\gg \left\{{{n}_{i}},~{{n}_{e}}\right\}>{{Q}_{\alpha}}{{n}_{\alpha}},\end{eqnarray} \tag{ 55 }$$

unless ${{n}_{\alpha}}(0)$ is very small. Then the first terms in () of (54) and in [ ] of K₁ are dominant (${{T}_{\alpha}}\sim {{T}_{i}}$) and a₁ become much smaller than unity (the value without fine particles);

$$\begin{eqnarray}&&{{a}_{1}}\sim {{\left[\frac{\left({{T}_{e}}/{{T}_{i}}\right)\left(1-{{c}_{\alpha}}{{n}_{\alpha}}/{{c}_{g}}{{n}_{e}}\right)}{\left({{T}_{e}}/{{n}_{e}}\right)\left(Q_{\alpha}^{2}~{{n}_{\alpha}}/{{T}_{\alpha}}\right)}\right]}_{0}}<\frac{{{T}_{\alpha}}}{{{T}_{i}}}{{\left(\frac{{{n}_{e}}}{Q_{\alpha}^{2}{{n}_{\alpha}}}\right)}_{0}}\ll 1.\end{eqnarray} \tag{ 56 }$$

The densities around the axis are thus given by

$$\begin{eqnarray}&&\frac{{{n}_{e}}}{{{n}_{e}}(0)}=1-\frac{{{a}_{1}}}{4}{{\left(\frac{R}{{{R}_{a}}}\right)}^{2}}+\ldots,\end{eqnarray} \tag{ 57 }$$

$$\begin{eqnarray}&&\frac{{{n}_{i}}}{{{n}_{e}}(0)}=\frac{{{n}_{i}}(0)}{{{n}_{e}}(0)}-\frac{{{D}_{a}}}{4{{D}_{i}}}{{\left(1-\frac{{{c}_{\alpha}}{{n}_{\alpha}}}{{{c}_{g}}{{n}_{e}}}\right)}_{0}}{{\left(\frac{R}{{{R}_{a}}}\right)}^{2}}+\ldots,\end{eqnarray} \tag{ 58 }$$

and

$$\begin{eqnarray}&&\frac{{{Q}_{\alpha}}{{n}_{\alpha}}}{{{n}_{e}}(0)}=\frac{{{Q}_{\alpha}}(0){{n}_{\alpha}}(0)}{{{n}_{e}}(0)}-\frac{{{D}_{a}}}{4{{D}_{i}}}{{\left(1-\frac{{{c}_{\alpha}}{{n}_{\alpha}}}{{{c}_{g}}{{n}_{e}}}\right)}_{0}}{{\left(\frac{R}{{{R}_{a}}}\right)}^{2}}+\ldots.\end{eqnarray} \tag{ 59 }$$

For n_e, the coefficient of ${{R}^{2}}/R_{a}^{2}$ is much smaller than unity and also than those for n_i and ${{Q}_{\alpha}}{{n}_{\alpha}}$, the latter two being the same and giving charge densities compensating each other. From (44), we see that the charge neutrality near the axis is satisfied to much higher accuracy than the case without particles, where ${{a}_{1}}\approx 1$.

When (55) is satisfied, the behaviour of the solutions near the axis is thus characterised by (a) much enhanced charge neutrality and almost flat potential, (b) almost flat electron density and (c) almost compensating (parallel) the charge densities of the ions and particles. As shown in (53) and (54), they are the result of the participation of the particles in the screening of the space charge. These will also be confirmed by the numerical solutions in section V.C. Though the above analysis is for the case of cylindrical symmetry, we may qualitatively apply the conclusions to the general case where we have appreciable amount of fine particle charges in the sense of (55). It should be noted that the dependency on the square of the charge number also appeared in (42) in the case without particles.

We emphasise again that (53) describes the behaviour of the space charge around the axis based on the expansion of solutions of the drift-diffusion equations (24)–(27). Although (54) takes the form of the Debye wave number squared when multiplied by ${{e}^{2}}/{{\varepsilon}_{0}}{{k}_{\text{B}}}$, the usual Debye theory for the screening of point charge is not applied, as in the case without particles (36) and (42).

In the experiments, one of controllable parameters may be the linear density of the fine particles along the axis N_z. We now relate N_z to the central density and the radius of distribution. From (59), we have the approximate radius of particle distribution R₀

$$\begin{eqnarray}&&{{\left(\frac{{{R}_{0}}}{{{R}_{a}}}\right)}^{2}}\sim 4\frac{{{Q}_{\alpha}}(0){{n}_{\alpha}}(0)}{{{n}_{e}}(0)}\frac{{{D}_{i}}}{{{D}_{a}}}\left(1-\frac{{{c}_{\alpha}}{{n}_{\alpha}}}{{{c}_{g}}{{n}_{e}}}\right)_{0}^{-1}\\ &&\sim 4\frac{{{Q}_{\alpha}}(0){{n}_{\alpha}}(0)}{{{n}_{e}}(0)}\frac{{{T}_{i}}}{{{T}_{e}}}\left(1-\frac{{{c}_{\alpha}}{{n}_{\alpha}}}{{{c}_{g}}{{n}_{e}}}\right)_{0}^{-1},\end{eqnarray} \tag{ 60 }$$

which is small compared with R_a, ${{R}_{0}}/{{R}_{a}}\ll 1$. This estimation will be compared with the numerical solutions in section V, and shown to be in reasonable agreement. (The factor ${{\left[1-{{c}_{\alpha}}{{n}_{\alpha}}(0)/{{c}_{g}}{{n}_{e}}(0)\right]}^{-1}}$ could be much larger than unity, but in the examples given there, ${{c}_{\alpha}}{{n}_{\alpha}}(0)/{{c}_{g}}{{n}_{e}}(0)<0.4$.) Then we have

$$\begin{eqnarray}&&{{N}_{z}}\sim 2\pi {{n}_{\alpha}}(0)\int_{0}^{{{R}_{0}}}RdR\left[1-{{\left(\frac{R}{{{R}_{0}}}\right)}^{2}}\right]=\frac{\pi}{2}{{n}_{\alpha}}(0)R_{0}^{2}.\end{eqnarray} \tag{ 61 }$$

In terms of N_z, we have

$$\begin{eqnarray}&&{{\left(\frac{{{Q}_{\alpha}}{{n}_{\alpha}}}{{{n}_{e}}}\right)}_{0}}\sim {{\left[\frac{1}{2}\frac{{{T}_{e}}}{{{T}_{i}}}{{\left(1-\frac{{{c}_{\alpha}}{{n}_{\alpha}}}{{{c}_{g}}{{n}_{e}}}\right)}_{0}}\frac{{{Q}_{\alpha}}(0){{N}_{z}}}{\pi R_{a}^{2}{{n}_{e}}(0)}\right]}^{1/2}},\end{eqnarray} \tag{ 62 }$$

$$\begin{eqnarray}&&\frac{{{R}_{0}}}{{{R}_{a}}}\sim {{2}^{3/4}}{{\left(\frac{{{T}_{i}}}{{{T}_{e}}}\right)}^{1/4}}\left(1-\frac{{{c}_{\alpha}}{{n}_{\alpha}}}{{{c}_{g}}{{n}_{e}}}\right)_{0}^{-1/4}{{\left(\frac{{{Q}_{\alpha}}(0){{N}_{z}}}{\pi R_{a}^{2}{{n}_{e}}(0)}\right)}^{1/4}}.\end{eqnarray} \tag{ 63 }$$

These expressions determine the shape of the particle distribution when the linear density N_z is given. The central density and the distribution radius increase in proportion to $N_{z}^{1/2}$ and $N_{z}^{1/4}$, respectively.

Let us now consider the effects of the ion drag force on the particles and discuss the condition for the formation of voids in cylindrical systems. We first assume that we have fine particles of one species α.

According to (48), the distributions have the same characteristics, (57)–(59) and therefore (a), (b) and (c), as far as we have the inequality at R  =  0

$$\begin{eqnarray}&&\left(1-{{f}_{\alpha}}\right)Q_{\alpha}^{2}{{n}_{\alpha}}\gg \left\{{{n}_{i}},~{{n}_{e}}\right\}>{{Q}_{\alpha}}{{n}_{\alpha}},\end{eqnarray} \tag{ 64 }$$

even if ${{f}_{\alpha}}$ increases and approaches unity. The characteristics really change when (at R  =  0)

$$\begin{eqnarray}&&\left(1-{{f}_{\alpha}}\frac{{{Q}_{\alpha}}{{n}_{\alpha}}}{{{n}_{i}}}\right)\frac{{{n}_{i}}}{{{T}_{i}}}\sim \frac{{{n}_{i}}}{{{T}_{i}}}>\left(1-{{f}_{\alpha}}\right)\frac{Q_{\alpha}^{2}{{n}_{\alpha}}}{{{T}_{\alpha}}}.\end{eqnarray} \tag{ 65 }$$

Then a₁ reduces to the value without particles

$$\begin{eqnarray}&&{{a}_{1}}\sim {{\left\{\frac{\left({{D}_{a}}/{{D}_{i}}\right)\left[1-{{c}_{\alpha}}{{n}_{\alpha}}/{{c}_{g}}{{n}_{e}}\right]}{\left[1-{{f}_{\alpha}}{{Q}_{\alpha}}{{n}_{\alpha}}/{{n}_{i}}\right]\left({{T}_{e}}/{{T}_{i}}\right)\left({{n}_{i}}/{{n}_{e}}\right)}\right\}}_{0}}\\ &&\sim {{\left\{\frac{\left({{D}_{a}}/{{D}_{i}}\right)\left[1-{{c}_{\alpha}}{{n}_{\alpha}}/{{c}_{g}}{{n}_{e}}\right]}{\left({{T}_{e}}/{{T}_{i}}\right)\left({{n}_{i}}/{{n}_{e}}\right)}\right\}}_{0}}\sim 1\end{eqnarray}$$

and we have

$$\begin{eqnarray}&&\frac{{{n}_{e}}}{{{n}_{e}}(0)}\sim 1-\frac{1}{4}{{\left(\frac{R}{{{R}_{a}}}\right)}^{2}},\end{eqnarray}$$

$$\begin{eqnarray}&&\frac{{{n}_{i}}}{{{n}_{e}}(0)}=\frac{{{n}_{i}}(0)}{{{n}_{e}}(0)}-\frac{1}{4}\frac{{{D}_{a}}}{{{D}_{e}}}{{\left(1-\frac{{{c}_{\alpha}}{{n}_{\alpha}}}{{{c}_{g}}{{n}_{e}}}\right)}_{0}}{{\left(\frac{R}{{{R}_{a}}}\right)}^{2}}+\ldots,\end{eqnarray}$$

$$\begin{eqnarray}&&\frac{{{n}_{\alpha}}}{{{n}_{e}}(0)}=\frac{{{n}_{\alpha}}(0)}{{{n}_{e}}(0)}-\frac{1}{4}\frac{{{T}_{e}}}{{{T}_{\alpha}}}{{\left[\left(1-{{f}_{\alpha}}\right)\frac{{{Q}_{\alpha}}{{n}_{\alpha}}}{{{n}_{e}}}\right]}_{0}}{{\left(\frac{R}{{{R}_{a}}}\right)}^{2}}+\ldots,\end{eqnarray}$$

where (by (65))

$$\begin{eqnarray}&&\frac{{{T}_{e}}}{{{T}_{\alpha}}}{{\left[\left(1-{{f}_{\alpha}}\right)\frac{{{Q}_{\alpha}}{{n}_{\alpha}}}{{{n}_{e}}}\right]}_{0}}<\frac{{{T}_{e}}}{{{T}_{i}}}{{\left[\frac{{{n}_{i}}}{{{Q}_{\alpha}}{{n}_{e}}}\right]}_{0}}\ll 1.\end{eqnarray}$$

Thus, the electron distribution reduces to the case of no particles and the charge density distributions of the particles and ions are still parallel (compensating), both being nearly flat (${{D}_{a}}/{{D}_{e}}\approx \left({{T}_{e}}/{{T}_{i}}\right)\left({{D}_{i}}/{{D}_{e}}\right)\sim \left({{T}_{e}}/{{T}_{i}}\right)\left({{v}_{i}}/{{v}_{e}}\right)\left({{\sigma}_{e}}/{{\sigma}_{i}}\right)\ll 1$, v_e being the electron thermal velocity and ${{\sigma}_{en}}$ and ${{\sigma}_{in}}$, the electron and the ion-neutral collision cross sections).

When $1-{{f}_{\alpha}}=0$, the particle distribution becomes flat. We expect this occurs at R  =  0 first (with respect to parameters): in order to have $1-{{f}_{\alpha}}=0$ first at some finite values of ${{R}^{\prime}}>0$, we have to satisfy two apparently independent conditions at ${{R}^{\prime}}$, $1-{{f}_{\alpha}}\left({{R}^{\prime}}\right)=0$ and ${{\left[\text{d}{{f}_{\alpha}}(R)/\text{d}R\right]}_{{{R}^{\prime}}}}=0$ and this is not realistic.

With a further increase in ${{f}_{\alpha}}$, $1-{{f}_{\alpha}}$ becomes negative. However, K₁ cannot be negative: a negative K₁ gives a₁  <  0 and therefore E  <  0, the inward electric field, which makes $\text{d}{{n}_{\alpha}}/\text{d}R$ and $\text{d}{{n}_{e}}/\text{d}R$ positive and $\text{d}{{n}_{i}}/\text{d}R$ negative, causing the increase in $|E|$ with the radius and its divergence. In order for K₁ to be still positive, ${{n}_{\alpha}}$ should satisfy (at R  =  0)

$$\begin{eqnarray}&&1\gg \frac{\left[1-{{f}_{\alpha}}{{Q}_{\alpha}}{{n}_{\alpha}}/{{n}_{i}}\right]}{{{Q}_{\alpha}}}\frac{{{n}_{i}}}{{{n}_{e}}}\sim \frac{1}{{{Q}_{\alpha}}}>\left(\,{{f}_{\alpha}}-1\right)\frac{{{Q}_{\alpha}}{{n}_{\alpha}}}{{{n}_{e}}}>0,\end{eqnarray}$$

because ${{Q}_{\alpha}}$ is 10³–10⁴. This indicates that ${{Q}_{\alpha}}{{n}_{\alpha}}/{{n}_{e}}$ is almost zero (formation of a void), except for the case where ${{f}_{\alpha}}$ is very close to unity. Thus, the void formation condition is approximately given by

$$\begin{eqnarray}&&{{f}_{\alpha}}(0)={{\left[\frac{F_{\alpha}^{id}}{{{Q}_{\alpha}}eE}\right]}_{R\to 0}}>1~~~\text{when}~~~{{n}_{\alpha}}(0)=0.\end{eqnarray} \tag{ 66 }$$

We note that this condition is naturally expected from the competition between the confining electric force and the ion drag force.

By (12) for $F_{\alpha}^{id}$ and noting ${{u}_{i}}\sim \left(e/{{m}_{i}}\right)\left({{\ell}_{i}}/{{v}_{i}}\right)E$, we can rewrite the condition (66) into the form

$$\begin{eqnarray}&&{{n}_{e}}(0)>n_{e}^{(c)}(0),\end{eqnarray} \tag{ 67 }$$

where the critical value $n_{e}^{(c)}(0)$ is determined by

$$\begin{eqnarray}&&n_{e}^{(c)}(0)\sim \frac{1}{4\pi {{c}_{1}}}{{\left[{{Q}_{\alpha}}(0){{\left(\frac{{{e}^{2}}}{4\pi {{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{i}}}\right)}^{2}}{{\ell}_{i}}\right]}^{-1}}.\end{eqnarray} \tag{ 68 }$$

On the right-hand side, c₁, ${{Q}_{\alpha}}$ and ${{\ell}_{i}}$ depend on parameters including n_e as in (14), (17), and

$$\begin{eqnarray}&&{{\ell}_{i}}\sim \frac{1}{{{n}_{n}}{{\sigma}_{in}}}=\frac{{{k}_{\text{B}}}{{T}_{n}}}{{{\sigma}_{in}}{{p}_{n}}}\end{eqnarray}$$

(${{\sigma}_{in}}$ is the ion-neutral gas collision frequency and n_n, p_n, and T_n are the neutral gas number density, pressure and temperature, respectively). Using the interpolation for $\Lambda$, (A.3) in appendix A, and the approximate T_e-dependence of z, $z\propto T_{e}^{-0.5}$, we have the parameter dependency of $n_{e}^{(c)}(0)$ from (68) in the form

$$\begin{eqnarray}&&n_{e}^{c}(0)\propto \frac{T_{i}^{1.7}}{T_{e}^{0.36}}{{\left(\frac{{{p}_{n}}}{{{T}_{n}}}\right)}^{1.43}}\frac{1}{r_{\alpha}^{0.57}}.\end{eqnarray} \tag{ 69 }$$

This will be compared with the results of the numerical analysis in section V.D.

It should be noted that when $1-{{f}_{\alpha}}>0$ is satisfied by all species, (48) leads to features (a), (b) and (c), even in the case of multiple species of particles. We now consider the case of two species of particles, α and β, with radii ${{r}_{\alpha}},~{{r}_{\beta}},~{{r}_{\alpha}}>{{r}_{\beta}}$, and therefore ${{f}_{\alpha}}(0)/{{f}_{\beta}}(0)={{Q}_{\alpha}}(0)/{{Q}_{\beta}}(0)>1$. When ${{f}_{\alpha}}(0)>1>{{f}_{\beta}}(0)$, we have a distribution of particles β in the void formed by particles α. The reason is as follows. Without particles β, particles α have a distribution with a void at the centre, as determined by the one-species case. Since the solutions in the void reduce to the ones without particles, we may apply, at least around the centre, the case of one species to the distribution of particles β. We thus have a structure where smaller particles distribute inside the void of larger particles. Details of the structure are difficult to have analytically and are left to the numerical treatment of the drift-diffusion equations in section V.E. As a distribution of the species α, there is a void at the centre. However, being different from the simple void, this void is filled with smaller fine particles of species β.

In numerical analyses, it is necessary to specify some parameters of the system. We assume the gas is Ar. The adopted values of the basic parameters and the corresponding values of the typical characteristic parameters are listed in tables 1 and 2, respectively. The radii 0.5 and 1 μm are within the applicability of (12)–(14) and are included in the range of particle size in experiments by PK-4.

Table 1. Basic parameters.

${{k}_{\text{B}}}{{T}_{e}}$	3 eV, 1 eV
Ti	300 K
${{T}_{\alpha}}$	300 K
ne(R  =  0)	$1\centerdot {{10}^{8}}$ $\text{c}{{\text{m}}^{-3}}$
${{r}_{\alpha}}$	0.5 μm, 1 μm
Rw	1.5 cm
pn	$10\sim 100$ Pa
${{\sigma}_{in}}$	$1\centerdot {{10}^{-14}}~\text{c}{{\text{m}}^{2}}$
${{\sigma}_{en}}$	$2\centerdot {{10}^{-16}}~\text{c}{{\text{m}}^{2}}$

Table 2. Typical values of the characteristic parameters (for T_e-dependent ones, the left and right values being for T_e  =  3 eV and T_e  =  1 eV, respectively).

${{R}_{a}}={{R}_{w}}/2.40$	0.624 cm
$1/{{k}_{e}}(0)={{\left[{{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{e}}/{{n}_{e}}(0){{e}^{2}}\right]}^{1/2}}$	$1.29\centerdot {{10}^{-1}},7.43\centerdot {{10}^{-2}}\,\text{cm}$
$1/{{k}_{i}}(0)\approx {{\left[{{\varepsilon}_{0}}{{k}_{\text{B}}}{{T}_{i}}/{{n}_{e}}(0){{e}^{2}}\right]}^{1/2}}$	$1.20\centerdot {{10}^{-2}}\,\text{cm}$
$\lambda (0)={{\left[{{k}_{e}}{{(0)}^{2}}+{{k}_{i}}{{(0)}^{2}}\right]}^{-1/2}}$	$1.19\centerdot {{10}^{-2}},1.18\centerdot {{10}^{-2}}\,\text{cm}$
$A=\left({{T}_{i}}/{{T}_{e}}\right){{\left[{{k}_{i}}(0){{R}_{a}}\right]}^{2}}$	23.5, 70.4

We numerically solve (24)–(27) introducing ${{N}_{e,~i,~\alpha}}$, as described in appendix B. The densities normalised by n_e(0) are integrated outwards from R  =  0 by the sixth-order Runge–Kutta method with some assumed value of $\Delta$. We adjust $\Delta$ so as to have asymptotically vanishing densities near the wall: starting integration with inappropriate values of $\Delta$ leads to a divergence of densities before reaching the wall. The solutions are quite sensitive to $\Delta$ and numerical computations with very high accuracy are necessary.

Examples of the set of solutions without fine particles are shown in figure 1 and compared with the analytical values obtained by the ambipolar diffusion equation, (33), (34) and (35). Starting from the quasi-neutral domain around the centre, the positive space charge increases with the radius, connecting to the domain of the sheath near the outer boundary. We observe that the known solutions given in section III are reproduced for R  <  R_a. It should be noted that, in our approach described in section II, the densities near the wall are not intended to be accurate, similar to the description by the ambiploar diffusion equation (30). In this sense, the comparison of our numerical result with the analytical one in section III has a physical meaning only within the domain R  <  R_a, both being inaccurate otherwise. We observe, however, that, our numerical solution gives a gradual increase in the space charge smoothly connecting to the domain of the sheath, in contrast to the diverging of the space charge given by the extrapolation of the analytic solution (35).

Zoom In Zoom Out Reset image size

In principle, a separate treatment is necessary for detailed analyses of the domain of the sheath. However, since our interest is mainly in the behaviour of the fine particles and the plasma around the axis and not in the domain of the sheath, we may adopt our approach for our purpose. We also emphasise that R_w is an outer boundary extrapolated from the domain described by the ambipolar diffusion equation and, in all the cases analysed in what follows, the domain of the particles' existence is limited well within R_w. Therefore, this separate treatment is not expected to be influenced by the existence of particles in our examples.

Based on the equations for the velocity-moments of the distribution functions, we can approximately describe the average ion and electron velocities as functions of radius including the sheath domain near the wall [14]. In figure 2, we show the ratio of the ion drift velocity u to the ion sound velocity c_s, obtained by moment equations in the cases of T_e  =  1 eV and 3 eV with T_i  =  300 K. This ratio is larger for higher electron temperatures and increases with the radius, starting from 0 and reaching unity around the edge of the sheath; for T_e  =  3 eV, $v/{{c}_{s}}\sim 0.1R/{{R}_{a}}$ in R/R_a  <  0.5. In our results shown below, the radius of the particle distribution is mostly well within 0.5 R_a, except for the case of and near the void formation, where it is of the order of 0.5 R_a. Although the ion thermal velocity is smaller than the ion sound velocity by the factor ${{\left({{T}_{i}}/{{T}_{e}}\right)}^{1/2}}$, the ratio of the ion drift velocity to the ion thermal velocity is less than 0.5 in R/R_a  <  0.5. This may partly justify our approach, where the ion flux is described by the diffusion coefficients and mobilities, as in (4) and (5).

Zoom In Zoom Out Reset image size

In section III we assumed that the outer boundary of the plasma is related to the radius of the apparatus R_w by ${{R}_{w}}=2.40{{R}_{a}}$ based on the behaviour of the solution of (30). Since the electron density nearly vanishes at $R\sim 3{{R}_{a}}$ in figure 2, it might seem more appropriate to assume ${{R}_{w}}>2.40{{R}_{a}}$. On the other hand, the decrease in the potential of the solution of (30) in figure 1 is closer to the result of the moment equation, especially for T_e  =  3 eV. Therefore, we still adopt ${{R}_{w}}=2.40{{R}_{a}}$ in this article keeping in mind that, if we should take ${{R}_{w}}=\xi {{R}_{a}}$ with $\xi >2.40$, the parameter R_w given in table 1 should be regarded as ${{R}_{w}}=(\xi /2.40)1.5$ cm.

Typical examples of densities in the case of fine particles of one species at T_e  =  1 eV and T_e  =  3 eV are shown in figure 3, where the potential and the net charge density are also shown. The radius of the particle distribution increases with the increase in the central particle density ${{n}_{\alpha}}(0)$. In these results, we observe the features (a), (b) and (c) as described in section IV.A.

Zoom In Zoom Out Reset image size

In the bottom panels in figure 3, we see that the charge neutrality is largely enhanced only in the domain where particles exist and it returns to values without fine particles outside that domain. In the domain of the particles' existence, the electrostatic potential and electron distribution become almost flat, and the charge density of the ions and the (absolute value of) the particle charge density are almost parallel, in accordance with the enhanced charge neutrality.

At the radius where the particle density vanishes, the drop in the ion density is slower than particles, and we have overshooting in the net positive charge density. While the potential is almost flat in the domain of the particles' existence, the potential changes at the edge of particle distribution. There, the decay of the particle density is sharper than that of the electrons or ions due to large charges on the particles, causing the overshooting of the positive charge.

The radius of particle distribution is compared with the estimation (60) in table 3. We observe that they are in good agreement, (60) overestimating the radius only slightly in most cases. We also show the ratio of the plasma loss onto the surface of particles ${{c}_{\alpha}}{{n}_{\alpha}}$ to the plasma generation ${{c}_{g}}{{n}_{e}}$ at the centre. Although the loss increases with the increase in the particle density, it is still less than half even when the central charge density of particles is about 20% of the electrons. The rest of the generated plasma is transported to the wall by the ambipolar diffusion.

Table 3. Distribution radius and the ratio ${{\left({{c}_{\alpha}}{{n}_{\alpha}}/{{c}_{g}}{{n}_{e}}\right)}_{R=0}}$ (${{r}_{\alpha}}=1$ μm).

${{k}_{\text{B}}}{{T}_{e}}$ (eV)	${{\left({{Q}_{\alpha}}{{n}_{\alpha}}/{{n}_{e}}\right)}_{R=0}}$	${{R}_{0}}/{{R}_{a}}$	${{R}_{0}}/{{R}_{a}}$ by (60)	${{\left({{c}_{\alpha}}{{n}_{\alpha}}/{{c}_{g}}{{n}_{e}}\right)}_{R=0}}$
1	0.19	0.17	0.18	0.36
1	0.097	0.10	0.11	0.18
1	0.049	0.07	0.075	0.090
3	0.17	0.075	0.086	0.19
3	0.088	0.045	0.058	0.095
3	0.035	0.015	0.036	0.038

With the increase in n_e(0), the particle distribution increases its radius and becomes flatter, as shown by the examples in figure 4. The corresponding behaviour of the potential and the net charge density is also shown.

Zoom In Zoom Out Reset image size

When n_e(0) exceeds some critical value $n_{e}^{(c)}(0)$, the particle density increases with R, even starting from a very small value at R  =  0: This signals the appearance of the void at the centre. The critical value $n_{e}^{(c)}(0)$ for the void formation decreases either by the increase in the particle radius ${{r}_{\alpha}}$, the increase in the electron temperature T_e or the decrease in the neutral gas pressure p_n. On the right-hand side of (67) in section IV.B, the increase in ${{r}_{\alpha}}$ or T_e gives the increase for ${{Q}_{\alpha}}(0)$, and the decrease in p_n gives the increase in ${{\ell}_{i}}$, all giving the decrease in $n_{e}^{(c)}(0)$. We have analysed the domain $0.5<{{r}_{\alpha}}\left[\mu \text{m}\right]<1$, $10<{{p}_{n}}\left[\text{Pa}\right]<100$, and $1<{{T}_{e}}\left[\text{eV}\right]<3$, and the results are plotted in figure 5. The dependency on ${{r}_{\alpha}}$, T_e, and p_n is approximately expressed by

$$\begin{eqnarray}&&n_{e}^{c}(0)\left[{{10}^{8}}~\text{c}{{\text{m}}^{-3}}\right]\sim 0.040\frac{p_{n}^{1.56}\left[\text{Pa}\right]}{r_{\alpha}^{0.42}\left[\mu \text{m}\right]T_{e}^{0.24}\left[\text{eV}\right]}\end{eqnarray} \tag{ 70 }$$

to within a few %. Although the above dependence on parameters is not exactly the same but consistent with our expectation (69), the latter being obtained using approximate interpolations for the parameter-dependence of $\Lambda$ and z. The main point of our results may be that, when ${{p}_{n}}\geqslant 40$ Pa, the critical density at the centre is sufficiently larger than 10⁹ cm⁻³ and rapidly increases with the increase in p_n.

Zoom In Zoom Out Reset image size

Examples of the distributions, the potential and the net charge density with the void are shown in figure 6 in the cases where n_e(0) slightly exceeds the critical value $n_{e}^{(c)}(0)$ for T_e  =  3 eV and ${{r}_{\alpha}}=1$ μm. Starting from very small values, the particle charge density reaches the maximum and decays with a further increase in R. In the domain of the void, we observe the increase in the net charge density or the weakened charge neutrality: This confirms that the enhanced charge neutrality (feature (a) in section IV.A) is due to the existence of the particles, as described there. (Since (B.4) indicates that ${{n}_{\alpha}}(0)=0$ means ${{n}_{\alpha}}(R)=0$ everywhere, we have to start with some finite but very small values for ${{n}_{\alpha}}(0)$. The solutions, however, are not so sensitive as the starting value ${{n}_{\alpha}}(0)$. The results shown here are obtained by assuming ${{Q}_{\alpha}}(0){{n}_{\alpha}}(0)/{{n}_{e}}(0)\sim {{10}^{-4}}$ and the behaviour of the solutions is almost the same, at least when the starting value is changed within factor 2.)

Zoom In Zoom Out Reset image size

Here, we consider the case where we have particles of two species or two radii. The void formation condition depends on the radius of the particles ${{r}_{\alpha}}$ through the charge on a particle ${{Q}_{\alpha}}$, which is approximately proportional to ${{r}_{\alpha}}$ and through the coefficient c₁.

When larger particles satisfy the no-void condition, the smaller particles are automatically under the same condition. Examples of distributions in this case are shown in figure 7. Even if ${{f}_{\alpha}}$ is smaller for smaller particles, it does not necessarily mean that the arithmetic sum of the inward electric force and the outward ion drag is larger for smaller particles. The results of the numerical analysis shown here indicate that the smaller particles are located inside the larger particles and, with R, the density of the smaller particles decays fast and is taken over by the larger particles. Since we also have features (a), (b) and (c) in multicomponent cases, the total charge density of the particles is kept almost constant. As a result, the addition of smaller particles (which only distribute very near to the axis) induces a higher density of larger particles in the domain where no smaller particles exist: this can be clearly seen in the lower figures in figure 7, where the density of the larger particles at the centre is kept constant.

Zoom In Zoom Out Reset image size

When the larger particles are under the void condition, the smaller particles are with either the void condition or the no-void condition. Here, we consider the latter case. An example of the distributions, the potential and the net charge density is shown in figure 8. Without the smaller particles, the larger particles form a void at the centre, as shown in figure 6. When we have both larger and smaller particles, the smaller ones distribute around the centre as if to fill the void of the larger particles. As shown in figure 8, due to the existence of smaller particles, the charge neutrality becomes satisfied to a higher accuracy in the domain very near to the centre, which is originally in the void of the larger particles. We also notice that the takeover of the particle charge density by the larger particles is not so smooth as the case shown in figure 7, as indicated by the overshooting of the positive charge density.

Zoom In Zoom Out Reset image size