Stochastic motion of grains with charge gradients in external electric fields

Studies of mechanisms of the energy exchange in systems of interacting particles are of significant interest in various areas of science and engineering (plasma physics, medicine, biology, polymer physics, etc.) [1–8]. A dusty plasma is an ionized gas containing charged grains (dust) of micron sizes. This plasma is widespread in nature and is produced in a number of technological processes [1–3]. The majority of laboratory studies of dusty plasma is carried out in gas discharges of various types [9–13].

In a weakly ionized plasma the massive dust particles efficiently dissipate their kinetic energy through collisions with gas neutrals (atoms/molecules). It is often assumed that as a result of the energy exchange between the dust and neutral components, they are in equilibrium. Nevertheless, experiments show that the stochastic kinetic energy of dust in a laboratory plasma can be considerably higher than the temperature of surrounding gas, $T_{n}$, and can reach the values ${\sim}1\text{--}10\ \text{eV}$. This phenomenon is usually called "anomalous heating". The mechanisms of this "heating" are commonly associated with various temporal and spatial fluctuations of dust charges [14–19]. However, the existing theoretical models do not always allow one to explain the gain of high kinetic energies for typical experimental conditions [20–23].

The charge of a dust particle is not a fixed quantity, it is determined by the local plasma parameters near this particle (for example, by the electron and/or ion number densities and their velocities) [1,2]. Two main factors that can lead to random dust charge fluctuations, that, in turn, can be the reason for the formation of irregular stochastic motion of dust, exist. The first of them is associated with the random nature of the ion/electron currents charging the dust particles [14,15]. The second is determined by the stochastic motion of dust in the volume of a spatially inhomogeneous plasma [24–26].

Since the dust charge is a function of the surrounding plasma parameters, a change in these parameters can lead to changes in the grain charge and the development of various instabilities in the dusty plasma [16–18,26–28]. The formation of regular self-oscillations (large-scale rotations and oscillations of separate particles) in the field of nonelectrostatic forces orthogonal to the grain charge gradient was considered in [16–18]. Besides the development of regular self-oscillations, in the numerical simulations of the problem the "anomalous heating" of the particles (the nature of which has had till now no satisfactory explanation) was also observed. Notice that all types of grain motions mentioned above are observed in laboratory experiments with gas-discharge plasma [17,18,26–30].

Stochastic fluctuations in the dust charge, which are determined by the random motion of grains in a spatially inhomogeneous plasma (due to their thermal or other types of fluctuations) in the presence of a dust charge gradient in the direction of the gravity force, are one of the possible mechanisms for the development of irregular oscillations of grains [26,28]. However, the first of the proposed models [26] describes the mechanism for the excitation of dust oscillations due to the collective fluctuations of interparticle forces in the dust cloud; the second model [28] is based on the assumption of free dust diffusion that is unsuitable for limited (finite) trajectories of dust particles often observed in laboratory plasma [1–3]. Thus, the above-mentioned models do not allow to describe the "heating" of an isolated dust particle due to its charge gradient in the field of a trap.

In this paper we present the mechanism describing the "anomalous heating" of charged grains due to their stochastic motion in the volume of a spatially inhomogeneous plasma. The proposed theoretical model was considered for one and two particles with the gradient of their charge for typical conditions of laboratory experiments in gas-discharge plasma. Analysis of such small-sized systems admits a simple analytical solution of the problem and provides a qualitative picture of the energy exchange in extended systems.

Let us consider the system of two linearized equations of motion that describe the displacements (ξ₁ and ξ₂) of one or two interacting identical particles with mass M from their equilibrium positions in the field of external forces under the action of a random force $F_{b1(2)}$,

$$\begin{eqnarray} &&M\frac{\text{d}^{2}\xi_{1}}{\text{d}t^{2}}=-Mv_{fr} \frac{\text{d}\,\xi_{1}}{\text{d}t}+a_{1} \xi_{1} +b_{1} \xi_{2} +F_{b1}, \end{eqnarray} \tag{ 1a }$$

$$\begin{eqnarray} &&M\frac{\text{d}^{2}\xi_{2}}{\text{d}t^{2}}=-Mv_{fr} \frac{\text{d}\,\xi_{2}}{\text{d}t}+a_{2} \xi_{2} +b_{2} \xi_{1} +F_{b2}, \end{eqnarray} \tag{ 1b }$$

$\nu_{fr}$

$a_{1(2)}$

$b_{1(2)}$

$F_{b1(2)}$

$\langle F_{b1}\rangle = \langle F_{b2}\rangle \equiv 0$

$\langle F_{b1} F_{b2}\rangle = 0$

$\langle F_{b1} V_{2}\rangle = \langle F_{b2} V_{1}\rangle \equiv 0$

$\langle F_{b1} \xi_{2}\rangle = \langle F_{b2} \xi_{1}\rangle \equiv 0$

$\langle F_{b1} \xi_{1}\rangle= \langle F_{b2} \xi_{2}\rangle \equiv 0$

$\langle F_{b1} V_{2}\rangle =\langle F_{b2} V_{1} \rangle \equiv 0$

$V_{1(2)}$

$\langle\ldots\rangle$

$t\to \infty$

The various classes of physical problems (such as the problem of stable positions of one or two particles in the external electric field [31–34], and some problems of stable positions of charged particles in an extended dust cloud [16–18]) can be described by equations similar to eqs. (1a), (1b).

A solution of eqs. (1a), (1b) is always unstable for $\nu_{fr}\le 0$ or $(a_{1}+ a_{2}) > 0$ [24]. Let us investigate the stability of solutions for the system (1a), (1b) with arbitrary coefficients $a_{1(2)}$ and $b_{1(2)}$ for $(a_{1}+ a_{2}) < 0$ and $\nu_{fr} > 0$ in two cases: i) when $s = [(a_{1}-a_{2})^{2}/4 + b_{1}b_{2}]^{1/2}$ is real, and ii) when s is imaginary. In the first case i), the condition for the development of a dissipative instability [7] in the system (1a), (1b) can be represented as

$$\begin{equation} b_{1} b_{2}-a_{1} a_{2} > 0. \end{equation} \tag{ 2a }$$

In the ii) case, the criterion for the development of a dispersive instability [7] in the system can be written as

$$\begin{equation} -(a_{2} + a_{1}) M\nu_{fr}^{2}+ (a_{2}- a_{1})^{2}+ 4b_{1}b_{2} < 0. \end{equation} \tag{ 2b }$$

Let us consider the conditions when the energy of the chaotic motion of grains (their kinetic temperature) ${T=T_{1(2)} \equiv M \langle V_{1(2)}^{2} \rangle}$ can be above their predetermined temperature $(T_{1}^{0}, T_{2}^{0})$, i.e., the possibility of a gain of additional energy $\delta T_{1(2)} = (T_{1(2)}-T_{1(2)}^{0})$ by the dust particles due to the forces acting in the system. Distinction of kinetic temperatures $T_{1}^{0} \ne T_{2}^{0}$ for different degrees of freedom or for different particles can be caused by various mechanisms of the gain of kinetic energy for dust in plasma (for example, fluctuations of dust charges caused by discrete currents of charging [24]). In the absence of sources of additional energy the value of $T_{1}^{0}= T_{2}^{0}\equiv T_{n}$.

In the case of finite trajectories: $\langle \xi_{1} V_{1} \rangle = \langle \xi_{2} V_{2} \rangle \equiv 0$, $\langle V_{1(2)}F_{b 1(2)}\rangle = \nu_{fr} T_{1(2)}^{0}$, and we can write the equations for the correlators of grain velocities and displacements

$$\begin{eqnarray} -M\nu_{fr} \delta T_{1(2)} +b_{1(2)} \left\langle {V_{1(2)} \xi _{2(1)} } \right\rangle &= 0, \end{eqnarray} \tag{ 3a }$$

$$\begin{eqnarray} a_{1(2)} \left\langle {(\xi_{1(2)})^{2}} \right\rangle +b_{1(2)} \left\langle {\xi_{1} \xi_{2}} \right\rangle +T_{1(2)} &= 0, \end{eqnarray} \tag{ 3b }$$

$$\begin{eqnarray} &&-M\nu_{fr} \left\langle {V_{1(2)} \xi_{2(1)} } \right\rangle +a_{1(2)} \left\langle {\xi_{1} \xi_{2} } \right\rangle\nonumber\\ +b_{1(2)} \left\langle {\xi_{1(2)}^{2} } \right\rangle +\left\langle {V_{1} V_{2}} \right\rangle &= 0, \end{eqnarray} \tag{ 3c }$$

$$\begin{eqnarray} -2M\nu_{fr} \left\langle {V_{1} V_{2}} \right\rangle +a_{1} \left\langle {\xi_{1} V_{2}} \right\rangle +a_{2} \left\langle {\xi_{2} V_{1}} \right\rangle &= 0. \end{eqnarray} \tag{ 3d }$$

$T_{1}^{0} = T_{2} \equiv T_{n}$

$$\begin{equation} \delta T\equiv \delta T_{1} +\delta T_{2} =\frac{\;2(b_{2} -b_{1} )^{2}T_{n} }{(a_{1} -a_{2} )^{2}+4b_{1} b_{2} -2M\nu_{fr}^{2} (a_{1} +a_{2})}, \end{equation} \tag{ 4 }$$

which is similar to the relation proposed for the case of nonreciprocal interaction of two particles [24,25] and/or for the case of a single grain with charge gradient [24] in a linear electric field. The energy redistribution (between degrees of freedom or between two different grains) in the system being studied can be easily found by taking into account that: $b_{2}\delta T_{1} =-b_{1}\delta T_{2}$. Notice that in the case of $b_{1}=b_{2}$: $(\delta T_{1} + \delta T_{2}) \equiv 0$. In systems with the significant dissipation (for $\nu_{fr} \to \infty$), the increment in kinetic energy $\delta T \quad \to$ 0, while the ratio $T_{2}/T_{1}\to 1$.

Thus, the gain of additional energy $\delta T$ is determined by the stochastic motion of dust in the volume of a spatially inhomogeneous plasma. Since the work of the forces along a closed loop in such systems is $A \propto T_{n}(b_{2}- b_{1})^{2}$, the additional increment in kinetic energy is $\delta T \propto T_{n}$.

Below we will consider the solution of a problem for one and two identical particles with charge $Q^{\ast} = Q^{\ast}(r, z)$ in the field of gravity Mg compensated by the electric field $E(r, z)$ of cylindrical trap with radial component $E_{r} =\beta_{r}r$ and vertical component $E_{z}=E_{z}^{0}+\beta_{z}z$, see fig. 1. Here $r = (x^{2}+y^{2})^{1/2}$ is the radial coordinate, z is the coordinate along an axis parallel to gravity, $\beta_{r}$ and $\beta_{z}$ are the values of gradients of electric field, and $E_{z}^0=Mg/Q$ is defined by balance of the forces acting in the system; here Q is the charge in a balance point.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The system of motion equations for one particle with charge $Q^{\ast}(r,z)$ and its gradient $\alpha(z,r)=\nabla Q^{\ast}$ in the linear electric field $E(r,z)$ can be written in the form (1a), (1b) where the values of $\xi_{1} = z$, $\xi_{2}= r$, and z, r are the displacements of the particle from its equilibrium position, $a_{1} = -Q\beta_{z} + \alpha_{z}E_{z}^{0}$, $b_{1} =\alpha_{r}E_{z}^{0},\ a_{2}= -Q\beta_{r} - \alpha_{r}E_{r}^{0}$, $b_{2}=\alpha_{z}E_{r}^{0},\ \alpha_{r} = \partial Q^{\ast}/\partial r$, $\alpha_{z} = \partial Q^{\ast}/\partial z$, $E_{r}^{0}=0$ and $E_{z}^{0}= Mg/Q$. In this case the value of $\delta T_{r}= 0$, and the solution of (1a), (1b) that describes the gain of additional energy can be presented as

$$\begin{equation} \delta T_{z}= \frac{\;2T_{z}^{0} (\alpha_{r} E_{z}^{0})^{2}} {(Q\beta_{z}\!-\!\alpha_{z} E_{z}^{0}\!-\!Q\beta_{r})^{2}\!+\!2M\nu_{fr}^{2} (Q\beta_{z}\!-\!\alpha_{z} E_{z}^{0}\!+\!Q\beta_{r})}. \end{equation} \tag{ 5 }$$

Dependences of $T_{z}/T_{z}^{0} = (T_{z}^{0} + \delta T_{z})/T_{z}^{0}$ on the ratio of $\omega_{r}/\nu_{fr}$ for one charged grain in the systems with $\alpha_{z} = 0$, $\alpha_{r}/Q = 0.05\ \text{cm}^{-1}$ and $\alpha_{r}/Q = 0.1\ \text{cm}^{-1}$ for various parameters of electric trap $\omega_{z}$ are shown in fig. 2; here $\omega_{z}^{2} = Q \beta_{z}/M$, $\omega_{r}^{2}= Q \beta_{r}/M$. Detailed study of an influence of plasma densities on the gradients of grain charges was presented in ref. [18] where it has been shown that in the laboratory plasma of capacitive RF discharges the value of the relative gradients of dust charges can reach $|\alpha_{r}/Q| \sim0.1\ \text{cm}^{-1}$, $|\alpha_{z}/Q|\sim0.5\ \text{cm}^{-1}$.

It is easy to see (fig. 2), that the gain of additional energy by the grains due to the gradients of their charges is comparable in magnitude to the dust heating caused by discrete currents of charging [14,15,24,35]. Because the latter mechanism (associated with the random nature of the ion/electron currents charging of dust is inherent for all types of plasma, it can be assumed that the efficiency of the anomalous heating of grains due to the dust charge gradients may turn out to be considerably higher because the gain of additional energy (5) will be defined by the value of $T_{z}^{0}> T_{n}$.

Here we note that for the case $\alpha_{z} =0$ with $\omega_{z} = \omega_{r}$ the energy $\delta T \propto \nu_{fr}^{-2}$ (5), and it can be assumed that $\delta T$ can increase indefinitely with the friction forces decrease (for $\nu_{fr}\to 0$). However, in real physical systems (such as, for example, a dusty plasma), a dynamic equilibrium is often established between the incoming and scattering energies. The increase in the kinetic energy of particles in such systems can be restricted owing to the existence of certain boundary conditions or peculiarities of the spatial distribution of parameters of the medium, which ensure a progressive growth of the dissipative losses, and due to various dispersion effects, which phenomenologically play a dissipative role [7,16]. In addition, various nonlinear effects can manifest themselves during the anomalous heating of particles [7,16–18,34].

We consider the system of two identical particles located at a distance d that are interacting with the force $F_{d}$, see fig. 1(b). For this case the motion equations may be presented as

$$\begin{eqnarray} &&M\frac{\text{d}^{2}z_{+}}{\text{d}t^{2}}=-Mv_{fr} \frac{\text{d}\,z_{+}}{\text{d}t}\!+\!a_{1+} z_{+}\!+\!b_{1+} r_{+}\!+\!F_{b1}^{z}\!+\!F_{b2}^{z}, \end{eqnarray} \tag{ 6a }$$

$$\begin{eqnarray} &&M\frac{\text{d}^{2}r_{+}}{\text{d}t^{2}}=-Mv_{fr} \frac{\text{d}\,r_{+}}{\text{d}t}\!+\!a_{2+} r_{+}\!+\!b_{2+} z_{+}\!+\!F_{b1}^{r}\!+\!F_{b2}^{r}, \end{eqnarray} \tag{ 6b }$$

$$\begin{eqnarray} &&M\frac{\text{d}^{2}z_{-}}{\text{d}t^{2}}=-Mv_{fr} \frac{\text{d}\,z_{-}}{\text{d}t}\!+\!a_{1-} z_{-}\!+\!b_{1-} r_{-}\!+\!F_{b1}^{z}\!-\!F_{b2}^{z}, \end{eqnarray} \tag{ 6c }$$

$$\begin{eqnarray} &&M\frac{\text{d}^{2}r_{-}}{\text{d}t^{2}}=-Mv_{fr} \frac{\text{d}\,r_{-}}{\text{d}t}\!+\!a_{2-} r_{-}\!+\!b_{2-} z_{-}\!+\!F_{b1}^{r}\!-\!F_{b2}^{r}. \end{eqnarray} \tag{ 6d }$$

$z_{1}+z_{2}= z_{+}$

$z_{1}-z_{2}=z_{-}$

$r_{1}+ r_{2} = r_{+}$

$r_{1}- r_{2} = r_{-}$

$z_{1}$

$z_{2}$

$r_{1}$

$r_{2}$

$a_{1+} = -Q\beta_{z} + \alpha_{z}E_{z}^{0}$

$b_{1+} \equiv b_{1-}=\alpha_{r}E_{z}^{0}$

$a_{2+}= - Q\beta_{r} - \alpha_{r}E_{r}^{0} +2\partial F_{d}/\partial r$

$b_{2+}= - \alpha_{z}E_{r}^{0} +2\partial F_{d}/\partial z \equiv - \alpha_{z}E_{r}^{0}+ 4 F_{d}\alpha_{z}/Q$

$a_{1-}=a_{1+}+2F_{d}/d$

$a_{2-} = - Q\beta_{r} - \alpha_{r}E_{r}^{0}$

$b_{2-}=\alpha_{z}E_{r}^{0}$

$Q\beta_{r}d = 2F_{d}$

$E_{r}^{0}=\beta_{r}d/2$

$E_{z}^{0}= Mg/Q$

Zoom In Zoom Out Reset image size

The conditions of instabilities for a solution of (6a)–(6d) with $\nu_{fr} > 0$ and $\{a_{1+(-)}+ a_{2+(-)}\} < 0$ can be represented as

$$\begin{eqnarray} b_{1+} b_{2+} - a_{1+} a_{2+} &> 0, \end{eqnarray} \tag{ 7a }$$

$$\begin{eqnarray} b_{1-} b_{2-} - a_{1-} a_{2-} &> 0, \end{eqnarray} \tag{ 7b }$$

$$\begin{eqnarray} - (a_{2+} + a_{1+}) M\nu_{fr}^{2} + (a_{2+} -a_{1+})^{2}+4b_{1+}b_{2+} &< 0, \quad \end{eqnarray} \tag{ 7c }$$

$$\begin{eqnarray} - (a_{2-} + a_{1-}) M\nu_{fr}^{2} + (a_{2-} -a_{1-})^{2}+4b_{1-}b_{2-} &< 0.\qquad \end{eqnarray} \tag{ 7d }$$

Regions of the steady solutions of a problem for two grains with Coulomb interaction $F_{d} = Q^{2}/d^{2}$, $\partial F_{d}/\partial r = - 2Q^{2}/d^{3} +2\alpha_{r}F_{d}/Q$ depending on the gradients of their charge, including the lines of development of dissipative (7b) and dispersive (7c) instabilities, are shown in a fig. 3 for the systems with parameters $d =750\,\mu\text{ m}$, $\omega_{r}=10\,\text{s}^{-1}$ and $\omega_{z} \cong 14\,\text{s}^{-1}$. In the systems without the gradients of charges (at $\alpha_{r} =\alpha_{z} \equiv 0$), the conditions (7a), (7b) completely correspond to the criteria of infringement of stability of a horizontal configuration of two particles interacting with various types of potentials [31–34].

For two identical grains: $T_{z 1}^{0}= T_{z 2}^{0} \equiv T_{z}^{0}$, ${T_{z 1}= T_{z 2} \equiv T_{z}}$, $T_{r 1}^{0}= T_{r 2}^{0} \equiv T_{r}^{0}$ and $T_{r 1}= T_{r 2} \equiv T_{r}$. Then, in case of $T_{z}^{0} = T_{r}^{0}\equiv T_{n}$, the solution of (6a)–(6d) can be written as

$$\begin{equation} \delta T \equiv \delta T_{r} +\delta T_{z} = D_{+}+ D_{-}, \end{equation} \tag{ 8 }$$

where

$$\begin{eqnarray} & \hspace{-10pt}D_{+(-)}=\nonumber\\ & \hspace{-10pt}\frac{T_{n} (b_{2+(-)}-b_{1+\!(-)})^{2}} {(a_{1+(-)}\!\!-\!a_{2+(-)})^{2}\!+\!4b_{1+(-)} b_{2+(-)} \!\!-\!\!2M \nu_{\!fr}^{2} (\!a_{1+(-)} \!\!+\! a_{2+(-)})}.\nonumber\\ \end{eqnarray} \tag{ 9 }$$

For the correlation functions of the velocities $\langle V_{z1}V_{z2} \rangle$ and $\langle V_{r1}V_{r 2} \rangle$, where $V_{z 1(2)}= \text{d}z_{1(2)}/\text{d}t$, $V_{r 1(2)} = \text{d}r_{1(2)}/\text{d}t$, we obtain

$$\begin{equation} M (\langle V_{z1}V_{z2} \rangle + \langle V_{r1}V_{r2} \rangle) = D_{+} - D_{-}, \end{equation} \tag{ 10 }$$

$$\begin{eqnarray} b_{2+} (\delta T_{z} +M\langle V_{z1} V_{z2} \rangle)&= -b_{1+} (\delta T_{r} +M \langle V_{r1}V_{r2} \rangle),\,\qquad \end{eqnarray} \tag{ 11a }$$

$$\begin{eqnarray} b_{2-} (\delta T_{z} -M\langle V_{z1} V_{z2} \rangle)&= -b_{1-} (\delta T_{r} -M \langle V_{r1}V_{r2} \rangle).\, \end{eqnarray} \tag{ 11b }$$

$$\begin{eqnarray} &&\delta T_{r} =\frac{\;b_{2+} D_{+} }{b_{2+} -b_{1+} }+\frac{\;b_{2-} D_{-}}{b_{2-} -b_{1-}}, \end{eqnarray} \tag{ 12a }$$

$$\begin{eqnarray} &&\delta T_{z} =\frac{\;b_{1+} D_{+} }{b_{1+} -b_{2+} }+\frac{\;b_{1-} D_{-}}{b_{1-} -b_{2-}}. \end{eqnarray} \tag{ 12b }$$

$\alpha_{z} =0$

$\delta T_{r} = 0$

$\delta T_{z} = D_{+}+D_-$

$\alpha_{r} =0$

$\delta T_{z} = 0$

$\delta T_{r} = D_{+}+ D_{-}$

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Notice that in most cases $E_{r}^{0} = \beta_{r}d/2 \ll E_{z}^{0} = Mg/Q$, and, accordingly, for comparable values of $\alpha_{z}$ and $\alpha_{r}$ the magnitude $\delta T_{r} \ll \delta T_{z}$. Dependences of $T_{z}/T_{n}$ on $\alpha_{r}/Q$ for various parameters of the system are shown in fig. 4, and the dependence of $T_{z}/T_{n}$ on $\alpha_{z}/Q$ is presented in fig. 5.

The dynamics of one $(N_{p}=1)$ and two $(N_{p}=2)$ grains with charge $Q^{\ast} = Q^{\ast} (r,z)$ in the field of gravity $|\textbf{F}_{g}|= Mg$ balanced by the electric field E of the cylindrical trap with components $E_{r} =\beta_{r}r$ and $E_{z}=E_{z}^{0} + \beta_{z}z$ (see figs. 1(a), (b)) was studied by the Langevin molecular dynamic method. This method is based on the solution of the system of $N_{p}$ differential equations of motions with random force $F_{b}$ which is a source of stochastic ("thermal") motions of particles with the kinetic temperature $T_{n}$ [1,2]. In our case the equations of motions can be presented as

$$\begin{equation} M\frac{\text{d}^{2}{\rm {\bf l}}_{k}}{\text{d}t^{2}}=\sum\limits_j {F_{d} (l)\frac{{\rm {\bf l}}_{k} -{\rm {\bf l}}_{j} }{\left| {{\rm {\bf l}}_{k} -{\rm {\bf l}}_{j} } \right|}-M\nu_{fr} \frac{\text{d}{\rm {\bf l}}_{k} }{\text{d}t}+Q^{\ast }{\rm {\bf E}}+{\rm {\bf F}}_{g} +{\rm {\bf F}}_{b}}, \end{equation} \tag{ 13 }$$

where $l=|\textbf{l}_{k} - \textbf{l}_{j}|$ is the inter-grain distance, $\textbf{l}_{k}$ is the vector of position of the mass centre for the k-grain, and $F_{d}$ is the force of interaction between two particles ${(F_{d} =0~\text{for}~N_{p}=1)}$.

The analytical model described above is valid for pair interactions of various types. Here for an examination of this model the simulations of the systems consisting of two particles were carried out for Coulomb pair interaction, $F_{d} = (Q^{\ast}/l)^{2}$. Notice that the Coulomb asymptotes $\propto l^{-2}$ for inter-grain interactions were observed in a laboratory RF discharge plasma [36].

The time integration step $\Delta t$ for solution of eq. (13) was from ${\sim}1/40\ \Omega$ to ${\sim}1/100\ \Omega$ where $\Omega=\max\{\nu_{fr}; \omega_{z}; \omega_{r}\}$ [1,2]. The temperature of grains was set as $T_{z}^{0}= T_{r}^{0}\equiv T_{n}$, where the magnitude of $T_{n}$ was changed from ${\sim}0.025\ \text{eV}$ to ${\sim}0.1\ \text{eV}$. Relative values of charge gradients $|\alpha_{r}/Q|$ and $|\alpha_{z}/Q|$ were varied from 0 to ${\sim}1.5\ \text{cm}^{-1}$. The magnitude of the vertical electric field of the trap $\beta_{z}$ was from ${\sim}1.5\beta_{r}$ to $3\beta_{r}$, that provided the steady existence of a horizontal configuration for the two particles ($\beta_{z} > \beta_{r}$ [31–34]). The characteristic frequency of the trap $\omega_{r} = (eZ \beta_{r}/M)^{1/2}$ was changed from $10\,\text{s}^{-1}$ to $100\,\text{s}^{-1}$. The value of the scaling parameter $\omega_{r}/\nu_{fr}$ was from ${\sim}1$ to ${\sim}25$, that is typical for the conditions of laboratory experiments in a gas discharge plasma [1–3].

Under conditions corresponding to the steady solution of a problem an additional kinetic energy $\delta T$ was registered in the simulated systems. In all cases the observed distributions of grain velocities were close to Maxwellian functions, with $\delta T_{z} \gg\delta T_{r}$, and $T_{z} \ge T_{r} \approx T_{n}$. The value of additional energy was proportional to the originally predetermined temperature of particles, $\delta T \propto T_{n}$. Thus, the "anomalous heating" systems which was accompanied by a nonuniform distribution of energy on the degrees of freedom was observed. Under a steady condition of the simulated system an additional kinetic energy $\delta T$ was registered. In all cases the observed distributions of grain velocities were close to Maxwellian functions, with $\delta T_{z} \gg \delta T_{r}$, and $T_{z} \ge T_{r} \approx T_{n}$. The value of the additional energy was proportional to the predetermined temperature of particles, $\delta T \propto T_{n}$. Thus, the "anomalous heating" systems which was accompanied by a nonuniform distribution of energy on the degrees of freedom was observed.

The results of numerical simulations for the case of one particle are presented in fig. 2 together with analytical solutions of this problem. Comparison of the results of numerical simulations with theoretical curves for two-particle systems is shown in figs. 4, 5. It is easy to notice the good accordance of numerical data with the proposed analytical relations. Deviations between the theoretical and numerical data were less ${\sim}5{\%}$ (up to a point of formation of instabilities in the two particle-systems). In the case of development of dispersive instability the motion of grains represented their synchronous regular oscillations which are considerably different from stochastic fluctuations. And in the case of development of dissipative instability we observed the change of horizontal configuration in the arrangement of particles on their vertical orientation.

In conclusion, we underline that in this paper we present the results of numerical and theoretical investigations of the dynamics of grains with charge gradients in external electric fields. The mechanism describing the "anomalous heating" of charged grains due to their stochastic motion in the volume of a spatially inhomogeneous plasma was considered, and new analytical relations for conditions of the grain heating were proposed. It has been shown that the presented mechanism is an effective source of additional stochastic energy for charged dust particles in laboratory plasmas along with alternative mechanisms such as the discrete currents of charging and/or the nonreciprocal interactions (caused by ion flows), that were proposed earlier for an explanation of the nature of "anomalous heating" of grains [14,15,24,25]. The estimation of the influence of the presented mechanism on the behavior of grains in plasma can be obtained for concrete conditions (with additional estimations of gradients of dust charges [18]) on the basis of proposed relations.

The results of the presented study can be useful for analyzing the energy exchange in inhomogeneous systems, which are of interest in plasma physics, medicine, biology, and the physics of polymers and colloidal systems.

We acknowledge financial support by the Russian Foundation for Basic Research (Grants Nos. 16-08-00594, 15-32-21159).