2.1. Particle dynamics

Following Newton's second law of motion, the dynamic equation of a single particle can be written as

(2.1)\begin{equation} m_p \frac{{\rm d} u_{p,i}}{{\rm d}t} = F_i, \end{equation}

where $m_p$ is the mass of the particle, $u_{p,i}$ is the velocity in the $i$th direction and $F_i$ is the $i$th component of the total force exerted on the particle.

In the frame of kinetic theory of rapid granular flows used for gas–solid flow modelling, the inter-particle forces are treated as instantaneous collisions. Equation (2.1) represents the particle dynamic equation between collisions and the force $F_i$, in the framework of the point-particle approximation, is written as (Gatignol Reference Gatignol1983; Maxey & Riley Reference Maxey and Riley1983)

(2.2)\begin{equation} F_i ={-}V_p \frac{\partial P_{g@p}}{\partial x_i} - \frac{m_p}{\tau_p}(u_{p,i} -u_{g@p,i}) + m_pg_i + q_pE_i, \end{equation}

where $V_p$ is the particle volume, $P_{g@p}$ is the pressure of the undisturbed flow at the particle position, $\tau _p$ is the particle relaxation time, $u_{g@p,i}$ is the $i$th component of the undisturbed gas flow velocity, $g_i$ is the gravitational acceleration, $q_p$ is the particle electric charge and $E_i$ is the $i$th component of electric field obtained by solving Gauss's flux law

(2.3)\begin{equation} E_i ={-}\boldsymbol{\nabla} \phi, \end{equation}

and the electric potential $\phi$ is given by

(2.4)\begin{equation} \boldsymbol{\nabla} \left(\varepsilon \boldsymbol{\nabla} \phi \right) ={-}\varrho, \end{equation}

where $\varepsilon$ is the medium permittivity and $\varrho$ is the electric charge density. Equations (2.2), (2.3) and (2.4) are exact if the electric potential is solved using the local instantaneous charge density. However, in this work, we will follow the strategy used in all previous studies in this subject, and we will solve the electric potential using the mean electric charge distribution $\varrho = n_p Q_p$ (where $n_p$ is the particle number density and $Q_p$ the mean particle electric charge, which will be clearly defined later). Therefore, the calculated electric field $E_i$ is an average, or macroscopic, electric field that does not take into account the instantaneous local distribution of charged particles. The effect of short-range electric interaction forces on the particle dynamics has been studied in very dilute turbulent flows (see e.g. Boutsikakis, Fede & Simonin Reference Boutsikakis, Fede and Simonin2022), but is generally completely neglected in dense flows and should be the subject of future studies using computational fluid dynamics/discrete element method (CFD/DEM) (Kolehmainen et al. Reference Kolehmainen, Ozel, Boyce and Sundaresan2016). On the other hand, the short-range electrical interaction between particles has a very important effect on the charge exchange between colliding particles and it is taken into account using the electrical collision model presented in the next section.

2.2. Particle dynamic and electrical collision model

In this work, we assume binary instantaneous collisions of frictionless inelastic spheres in translation, but the extension to inelastic and frictional collisions with rotation of the particles is feasible within the proposed approach (Yang, Padding & Kuipers Reference Yang, Padding and Kuipers2016). During these collisions, the impacting particles will only exchange linear momentum and electric charge.

Let us consider two colliding particles, $p_1$ and $p_2$, whose centres are located at $\boldsymbol {x_{p_1}}$ and $\boldsymbol {x_{p_2}}$. Before the collision, the particles have given velocities $\boldsymbol {c_{p_1}}$ and $\boldsymbol {c_{p_2}}$, and given electric charges $\xi _{p_1}$ and $\xi _{p_2}$. We define the vector $k_i$ as the unit vector going from the centre of $p_1$ to the centre of $p_2$. We also define the vector $\boldsymbol {g_{r}}$ as the relative velocity between $p_1$ and $p_2$: $\boldsymbol {g_{r}}=\boldsymbol {c_{p_1}} - \boldsymbol {c_{p_2}}$.

Following previous works on the Eulerian modelling of electrostatic forces, we are going to neglect the effect of any electric field (mean or fluctuating) during the particle–particle collisions. Hence, the velocities after the collision are given by $\boldsymbol {c_{p_1}^{+}}$ and $\boldsymbol {c_{p_2}^{+}}$

(2.5)$$\begin{gather} c_{p_1,i}^+= c_{p_1,i} - \tfrac{1}{2} \left( 1+e_c \right)g_{r,j} k_j k_i, \end{gather}$$

(2.6)$$\begin{gather}c_{p_2,i}^+= c_{p_2,i} + \tfrac{1}{2} \left( 1+e_c \right)g_{r,j} k_j k_i, \end{gather}$$

where $e_c$ is the collision restitution coefficient.

We also need to account for the transfer of electric charge between the particles during contact. For this, we use the model proposed by Kolehmainen et al. (Reference Kolehmainen, Ozel, Boyce and Sundaresan2017) following the work of Laurentie, Traoré & Dascalescu (Reference Laurentie, Traoré and Dascalescu2013). According to this model, the transfer of electric charge is written as a function of the interparticle contact area due to the elastic deformation of the particles during collisions. Some charge might be transferred when the particles are moving away from each other, however, this is considered negligible in their modelling approach. Therefore, when two particles $p_1$ and $p_2$ collide, the charge evolution of the particle $p_1$ during the collision, can be written as

(2.7)\begin{equation} \frac{{\rm d} q_{p_1}}{{\rm d}t} = \frac{{\rm d}\mathcal{A}}{{\rm d}t} (\kappa_1 - \kappa_2 q_{p_1}) \quad \text{for } \ \frac{{\rm d}\mathcal{A}}{{\rm d}t} > 0, \end{equation}

where $\kappa _1$ is a coefficient that depends on the collision type and the pre-collision particle charge, $\kappa _2$ is a geometrical coefficient that depends only on the particle diameter and $\mathcal {A}$ is the contact area during the collision. This contact area can be written as a function of the overlap distance $\delta$

(2.8)\begin{equation} \mathcal{A} = \frac{8 {\rm \pi}\delta}{d_p}, \end{equation}

where

(2.9)\begin{equation} \delta = \frac{d_{p_1} + d_{p_2}}{2} - \left| \boldsymbol{x_{p_1}} - \boldsymbol{x_{p_2}} \right|. \end{equation}

Finally, the total charge transferred during the collision is given by the solution of (2.7) when the contact area is maximum. Using these hypotheses, Kolehmainen et al. (Reference Kolehmainen, Ozel, Boyce and Sundaresan2017) showed the particle charge after the collisions are given by

(2.10)$$\begin{gather} \xi_{p_1}^+= \xi_{p_1} - \varepsilon_0 \mathcal{A}_{max} E_i^{+}k_i, \end{gather}$$

(2.11)$$\begin{gather}\xi_{p_2}^+= \xi_{p_2} + \varepsilon_0 \mathcal{A}_{max} E_i^{+}k_i. \end{gather}$$

Also, $\mathcal {A}_{max}$ is given by

(2.12)\begin{equation} \mathcal{A}_{max} = {\rm \pi}\frac{d_p}{2}\left( \frac{30m_p ( 1 - \nu^2 )}{32Y\sqrt{d_p}} \right)^{{2}/{5}} (g_{r,i}k_i)^{{4}/{5}}, \end{equation}

where $d_p$ is the particle diameter, $Y$ is the particle Young's modulus, $\nu$ is the particle Poisson's ratio and $\varepsilon _0$ is the vacuum permittivity.

In this case, we do decompose the local electric field $E_i^{+}$ into its two components, the mean macroscopic electric field $E_i$ and the local instantaneous electric field at the contact point generated by the two colliding particles

(2.13)\begin{equation} E_{i}^{+} = E_i - \frac{\xi_{p_2} - \xi_{p_1}}{{\rm \pi} \varepsilon_0 d_p^2}k_i. \end{equation}

According to this tribocharging model, we can write the electric charge transferred during a particle–particle collision as a contribution of two separate terms. One contribution is proportional to the projection of the mean electric field on the vector $\boldsymbol {k}$ and another contribution proportional to the difference on the electric charge of the colliding particles

(2.14)$$\begin{gather} \xi_{p_1}^{+} = \xi_{p_1} + [-\beta E_i k_i + \gamma (\xi_{p_2} - \xi_{p_1})] {(g_{r,i} k_i)}^{{4}/{5}}, \end{gather}$$

(2.15)$$\begin{gather}\xi_{p_2}^{+} = \xi_{p_2} - [-\beta E_i k_i + \gamma (\xi_{p_2} - \xi_{p_1})] (g_{r,i} k_i)^{{4}/{5}}, \end{gather}$$

where $\beta$ and $\gamma$ are given by

(2.16)$$\begin{gather} \beta = \varepsilon_0 {\rm \pi}\frac{d_p}{2}\left( \frac{30m_p (1 - \nu^2)}{32Y\sqrt{d_p}} \right)^{{2}/{5}}, \end{gather}$$

(2.17)$$\begin{gather}\gamma = \frac{1}{2\varepsilon_0 d_p}\left( \frac{30m_p(1-\nu^2)}{32Y\sqrt{d_p}} \right)^{{2}/{5}}. \end{gather}$$

2.3. Transport equation for the solid phase mean properties

The kinetic theory of rapid granular flow is based on the analogy between the motion of solid particles in a gas–solid flow and the thermal motion of molecules in a gas. In this approach, we define $f( \boldsymbol {x}, \boldsymbol {c_p}, \xi _p ; t )\delta \boldsymbol {x}\delta \boldsymbol {c_p}\delta \xi _p$ as the mean probable number of particles at time $t$ with their centre $\boldsymbol {x_p} ( t )$ in the volume element $[\boldsymbol {x}, \boldsymbol {x} + \delta \boldsymbol {x}[$, with a velocity $\boldsymbol {u_p} ( t )$ in the range $[ \boldsymbol {c_{p}}, \boldsymbol {c_p} + \delta \boldsymbol {c_p}[$ and an electric charge $q_p(t)$ in the range $[\xi _{p}, \xi _p + \delta \xi _p[$, respectively. Where $\boldsymbol {x}$, $\boldsymbol {c_p}$ and $\xi _p$ are the phase space coordinates for the particle position, velocity and electric charge. Also, $\boldsymbol {x_{p}} ( t )$, $\boldsymbol {u_p} ( t )$ and $q_p ( t )$ are the short form notation for any $n$-particle position, velocity and electric charge for a given realization $r$ of the ensemble averaging $\boldsymbol {x_p} ( t ) = \boldsymbol {x_p}^{( n, r )} ( t )$, $\boldsymbol {u_p} ( t ) = \boldsymbol {u_p}^{( n, r )} ( t )$ and $q_p ( t ) = q_p^{( n, r )} ( t )$, where $n = 1, 2, 3,\ldots, N_p$ with $N_p$ being the total number of particles.

This probability density function allows us to define the particle number density

(2.18)\begin{equation} n_p \left( \boldsymbol{x}, t \right) = \int_{\mathbb{R}} \int_{\mathbb{R}^{3}} f ( \boldsymbol{x}, \boldsymbol{c_p}, \xi_p, t) {\rm d}\boldsymbol{c_p}\,{\rm d}\xi_p. \end{equation}

We can also define the mean value for any particle property $\phi (t, \boldsymbol {x}, \boldsymbol {c_p}, \xi _p)$

(2.19)\begin{equation} \langle \phi \rangle \left( \boldsymbol{x}, t \right) = \frac{1}{n_p} \int_{\mathbb{R}} \int_{\mathbb{R}^3} \phi f (\boldsymbol{x}, \boldsymbol{c_p}, \xi_p, t) {\rm d}\boldsymbol{c_p}\,{\rm d}\xi_p. \end{equation}

This mean property $\langle \phi \rangle$ allows us to define the fluctuant value $\phi ^{\prime }$

(2.20)\begin{equation} \phi^{\prime} = \phi - \langle \phi \rangle. \end{equation}

By definition, we have that $\langle \phi ^{\prime } \rangle = 0$.

Using the Grad–Boltzmann limit from the Liouville equation, we can deduce the dynamic equation for $f$

(2.21)\begin{align} &\frac{\partial f}{\partial t}+\frac{\partial}{\partial x_{i}}\left[c_{p,i}\,f\right]+\frac{\partial}{\partial c_{p,i}}\left[\left\langle \frac{{\rm d}u_{p,i}}{{\rm d}t} \mid \boldsymbol{x},\boldsymbol{c_{p}}, \xi_{p}\right\rangle f\right]\nonumber\\ &\quad + \frac{\partial}{\partial\xi_{p}}\left[\left\langle\frac{{\rm d}q_{p}}{{\rm d}t}\mid\boldsymbol{x},\boldsymbol{c_{p}}, \xi_{p}\right\rangle f\right]=\left(\frac{\partial f}{\partial t}\right)_{coll} , \end{align}

where the notation $\langle \varPhi \mid \boldsymbol {x}, \boldsymbol {c_p}, \xi _p \rangle$ is a short form for the conditional mean $\langle \varPhi \mid \boldsymbol {x_p}(t) = \boldsymbol {x}, \boldsymbol {u_p}(t)=\boldsymbol {c_p}, q_p(t)=\xi _p ; t \rangle$.

Experimental results have proven that the gas–particle contact does not charge the particles (Mehrani, Bi & Grace Reference Mehrani, Bi and Grace2005). The only charging mechanisms are the particle–particle or particle–wall collisions. Therefore, the last term on the left-hand side of (2.21) can be discarded

(2.22)\begin{equation} \frac{{\rm d}q_p}{{\rm d}t} = 0. \end{equation}

If we multiply (2.21) by $\phi \,{\rm d}\boldsymbol {c_p^{\prime }}\,{\rm d}\xi _p^{\prime }$, and then integrate over the whole phase space, we can derive the transport equation for the mean value $\left \langle \phi \right \rangle$

(2.23) \begin{align} &\frac{{\rm D}n_{p}\left\langle \phi\right\rangle }{{\rm D}t}+n_{p}\left\langle \phi\right\rangle \frac{\partial U_{p,i}}{\partial x_{i}}+\frac{\partial n_{p}\langle \phi c'_{p,i}\rangle }{\partial x_{i}}-n_{p}\left\langle \frac{{\rm D}\phi}{{\rm D}t}\right\rangle -n_{p}\left\langle c'_{p,i}\frac{\partial\phi}{\partial x_{i}}\right\rangle\nonumber\\ &\quad -n_{p}\left\langle \left\langle \frac{{\rm d} u_{p,i}}{{\rm d}t} | \boldsymbol{x}, \boldsymbol{c_p}, \xi_p \right\rangle \frac{\partial\phi}{\partial c'_{p,i}}\right\rangle +n_{p}\frac{{\rm D}U_{p,i}}{{\rm D}t}\left\langle \frac{\partial\phi}{\partial c'_{p,i}}\right\rangle +n_{p}\left\langle c'_{p,j} \frac{\partial\phi}{\partial c'_{p,i}}\right\rangle \frac{\partial U_{p,i}}{\partial x_{j}}\nonumber\\ &\quad +n_{p}\frac{{\rm D}Q_{p}}{{\rm D}t}\left\langle \frac{\partial\phi}{\partial\xi_{p}^{\prime}}\right\rangle -n_{p}\left\langle \left\langle \frac{{\rm d}q_{p}}{{\rm d}t} | \boldsymbol{x}, \boldsymbol{c_p}, \xi_p \right\rangle \frac{\partial\phi}{\partial\xi_{p}^{\prime}}\right\rangle +n_{p}\left\langle \frac{\partial\phi}{\partial\xi_{p}^{\prime}}c_{p,i}^{\prime}\right\rangle \frac{\partial Q_{p}}{\partial x_{i}} = \mathcal{C}\left( \phi \right), \end{align}

where $U_{p,i} = \langle c_{p,i} \rangle$ and $Q_p = \langle \xi _p \rangle$.

The term on the right-hand side of (2.23) accounts for the mean transfer rate of property $\phi$, during particle–particle collisions (Jenkins & Savage Reference Jenkins and Savage1983; Jenkins & Richman Reference Jenkins and Richman1985). The details of how this term is calculated are given in Montilla et al. (Reference Montilla, Ansart and Simonin2020). From this right-hand side term, we will obtain all the collisional transport terms, such as the collisional triboconductivity and the collisional dispersion terms, in contrast with what we will call the kinetic terms that are associated with the transport of the property $\phi$ due to the particle-velocity fluctuations (such as the third term in (2.23)).

2.4. Mean electric charge transport equation

Kolehmainen et al. (Reference Kolehmainen, Ozel and Sundaresan2018b) were the first to derive a transport equation for the mean particle electric charge $Q_p$. Assuming that the particle velocity and electric charge are not correlated, they were able to propose a first approximation to the mean electric charge transport equation. Later, Montilla et al. (Reference Montilla, Ansart and Simonin2020) extended this work, proposing a linear model for the conditional mean $\langle \xi _p \mid \boldsymbol {c_p} \rangle$. This led them to the more general equation

(2.24) \begin{align} n_{p}\frac{\partial Q_{p}}{\partial t}+n_{p}U_{p,i} \frac{\partial Q_{p}}{\partial x_{i}}+\frac{\partial}{\partial x_i} ( n_{p} ( 1 + \eta_{coll} ) \langle \xi'_{p}c'_{p,i}\rangle) ={-} \frac{\partial}{\partial x_i}( \sigma_p E_i) + \frac{\partial}{\partial x_i}\left( n_p D_p \frac{\partial Q_p}{\partial x_i} \right) , \end{align}

where the triboconductivity coefficient $\sigma _p$, the collisional dispersion coefficient $D_p$ and $\eta _{coll}$ are given by

(2.25)$$\begin{gather} \sigma_p = \lambda_{1.1} d_{p}^{3}\beta g_{0}n_{p}^{2}(\varTheta_{p})^{{9}/{10}} , \end{gather}$$

(2.26)$$\begin{gather}D_p = \lambda_{1.2} d_{p}^{4}\gamma g_{0}n_{p}(\varTheta_{p})^{{9}/{10}} , \end{gather}$$

(2.27)$$\begin{gather}\eta_{coll} = \lambda_{1.3} d_p^3\gamma g_0 n_p \varTheta_p^{{2}/{5}}, \end{gather}$$

(2.28)$$\begin{gather}\varTheta_p = \frac{R_{p,ii}}{3}, \end{gather}$$

(2.29)$$\begin{gather}R_{p,ij} = \langle c_{p,i}^{\prime} c_{p,j}^{\prime} \rangle, \end{gather}$$

and $\lambda _{1.1} = \lambda _{1.2} \approx 1.825$, and $\lambda _{1.3} \approx 5.936$.

In order to close (2.24), the second-order moment $\langle \xi _p^{\prime } c_{p,i}^{\prime } \rangle$ needs to be modelled in terms of computed variables. The existing literature has focused on algebraic closure laws for this particle charge–velocity covariance term, either by analogy with the kinetic dispersion of the particle temperature (Kolehmainen et al. Reference Kolehmainen, Ozel and Sundaresan2018b) or by simplifying the covariance transport equation (Ray et al. Reference Ray, Chowdhury, Sowinski, Mehrani and Passalacqua2019; Montilla et al. Reference Montilla, Ansart and Simonin2020). Their results have shown that this modelling approach leads to extra dispersion and triboconductivity effects due to the random motion of particles. In this study, we will analyse a more complex approach by keeping the full transport equation for the particle charge–velocity covariance vector and we will derive a transport equation for the charge variance $\langle \xi _p^{\prime } \xi _p^{\prime } \rangle$. Additionally, we will propose algebraic closure laws for the third-order moments appearing in those transport equations.

5.1. Dimensionless analysis

Given the simplicity of this system, we can rewrite the governing equations in a simpler dimensionless form. We can choose $L_{{ref}} = L$ as reference length, $Q_{p,{ref}} = Q_{p,0}$ as reference electric charge, $U_{p,{ref}} = \sqrt {\varTheta _p}$ as reference velocity and $E_{{ref}} = {n_pQ_{p,{ref}}L}/{\varepsilon _0}$ as reference electric field. These choices lead to a reference time equal to $t_{{ref}} = {L}/{\sqrt {\varTheta _p}}$. Using these characteristic scales, we can express the governing equations as:

1. (i) The dimensionless mean particle charge transport equation

   (5.2) \begin{equation} \frac{\partial Q_p^*}{\partial t^*} + \underbrace{\left( 1+\eta_{coll} \right)\frac{\partial \langle c_p^{*\prime} \xi_p^{* \prime}\rangle}{\partial x^*}}_{\textit{Kinetic flux}} ={-} \underbrace{ \frac{1}{\tau_{\sigma}^*} \frac{\partial E^*}{\partial x^*} }_{\textit{Triboconductivity}} + \underbrace {\frac{1}{Pe}\frac{\partial ^ 2 Q_p^*}{\partial x^{*2}}}_{\textit{Collisional dispersion}}. \end{equation}

2. (ii) The dimensionless particle charge–velocity covariance transport equation

   (5.3) \begin{align} & \frac{\partial \langle c_p^{*\prime} \xi_p^{*\prime}\rangle}{\partial t^*} + \underbrace{ \left( 1 + \lambda_{2.2}e_c\frac{1}{Pe}\frac{L}{d_p} \right)\frac{\partial Q_p^*}{\partial x^*}}_{\textit{Production}} - \underbrace{\lambda_{2.1}e_c\frac{1}{\tau_{\sigma}^*}\frac{L}{d_p}E^*}_{\textit{Production}} - \underbrace{ \frac{u_e}{u_k}\langle \xi_p^{*\prime} \xi_p^{*\prime} \rangle E^* }_{\textit{Elec. force}} \nonumber\\ &\quad ={-} \underbrace{ \left( \frac{1+e_c}{3}\frac{1}{\tau_c^*} + \frac{1}{\tau_p^*} + \frac{3-e_c}{5}\frac{1}{\tau_{\xi}^*} \right)\langle c_p^{*\prime} \xi_p^{*\prime} \rangle }_{\textit{Destruction}} + \left[ 2 + \frac{1}{\tau_{\xi}^*}\frac{L}{d_p} \left( 2\lambda_{2.9} +\lambda_{2.10} \right)\left( 1+e_c \right) \right. \nonumber\\ & \quad \left. \vphantom{\frac{1}{\tau_{\xi}^*}} + 3\lambda_{2.11}\frac{1}{\tau_{\xi}^*}{\frac{L}{d_p}}\left( 1+e_c \right)^2 \right] \frac{1}{K_1 + K_2} \underbrace{ \frac{\partial^2 \langle \xi_p^{*\prime} c_{p}^{*\prime} \rangle}{\partial x^{*2}}}_{\textit{Collisional + Kinetic dispersion}}. \end{align}

3. (iii) The dimensionless particle charge variance transport equation

   (5.4) \begin{align} & \frac{\partial\langle \xi_{p}^{*\prime}\xi_{p}^{*\prime}\rangle}{\partial t^*}+ \underbrace{ 2\langle \xi_{p}^{*\prime}c_{p}^{*\prime}\rangle\frac{\partial Q^{*}_{p}}{\partial x^{*}} }_{\textit{Production}} ={-}\underbrace{\left(\frac{1}{\tau_{\xi}^{*}}-\lambda_{7}\frac{1}{{Pe}^2}\tau_c^* \left( \frac{L}{d_p} \right)^{4} \right)\langle \xi_{p}^{*\prime}\xi_{p}^{*\prime}\rangle}_{\textit{Destruction}} \nonumber\\ &\quad + \underbrace{\lambda_{8} \frac{\tau_c^*}{\tau_{\sigma}^2} \left( \frac{L}{d_p} \right)^2 \left| E^* \right|^2}_{\textit{Production}} + \underbrace{ \lambda_{9} \frac{1}{Pe^2}\tau_c^*\left( \frac{L}{d_p} \right)^2 \frac{\partial^2\langle \xi_{p}^{*\prime}\xi_{p}^{*\prime}\rangle }{\partial x^{*2}}}_{\textit{Collisional dispersion}} \nonumber\\ &\quad + \underbrace{\left[ \frac{1+e_c}{3}\frac{1}{\tau_c^*} + \frac{1}{\tau_p^*} + \frac{1}{\tau_{\xi}^*}\left( 2.8 - 20.39\frac{1}{Pe}\tau_c^* \left( \frac{L}{d_p} \right)^2 \right) \right]^{{-}1} \frac{\partial^2\langle \xi_{p}^{*\prime}\xi_{p}^{*\prime}\rangle }{\partial x^{*2}}}_{\textit{Kinetic dispersion}}. \end{align}

And the dimensionless Maxwell equations can be written as

(5.5)$$\begin{gather} E^* ={-}\boldsymbol{\nabla} \phi^*, \end{gather}$$

(5.6)$$\begin{gather}\nabla^2\phi^* ={-} Q_p^*, \end{gather}$$

where the dimensionless variables $Q_p^{*}$, $\langle c_{p,i}^{\prime *} \xi _{p}^{\prime *} \rangle$, $\langle \xi _p^{\prime *} \xi _p^{\prime *} \rangle \ E^{*}$, $x^*$, $t^*$ are given by

(5.7)$$\begin{gather} Q_p^{*} = \frac{Q_p}{Q_{p,0}} \end{gather}$$

(5.8)$$\begin{gather}E^{*} = \frac{E}{E_{{ref}}} \end{gather}$$

(5.9)$$\begin{gather}\langle c_{p,i}^{\prime *} \xi_{p}^{\prime *} \rangle = \frac{\langle c_{p,i}^{\prime } \xi_{p}^{\prime } \rangle}{Q_{p,0}\varTheta_{p}} \end{gather}$$

(5.10)$$\begin{gather}\langle \xi_{p}^{\prime *} \xi_{p}^{\prime *} \rangle = \frac{\langle \xi_{p}^{\prime } \xi_{p}^{\prime } \rangle}{Q_{p,0}Q_{p,0}} \end{gather}$$

(5.11)$$\begin{gather}t^{*} = t \frac{\sqrt{\varTheta_p}}{L} \end{gather}$$

(5.12)$$\begin{gather}x^{*} = \frac{x}{L}. \end{gather}$$

The seven dimensionless parameters are expressed as

(5.13)$$\begin{gather} Pe = \frac{L \sqrt{\varTheta_p}}{D_p} \end{gather}$$

(5.14)$$\begin{gather}\tau_{\sigma}^* = \frac{\varepsilon_0}{\sigma_p} \frac{\sqrt{\varTheta_p}}{L} \end{gather}$$

(5.15)$$\begin{gather}\tau_c^{*} = \tau_c \frac{\sqrt{\varTheta_p}}{L} \end{gather}$$

(5.16)$$\begin{gather}\tau_p^{*} = \tau_p \frac{\sqrt{\varTheta_p}}{L} \end{gather}$$

(5.17)$$\begin{gather}\frac{u_e}{u_k} = \frac{\dfrac{1}{2} \varepsilon_0 E_{ref}^2}{\dfrac{1}{2} n_p m_p \varTheta_p} \end{gather}$$

(5.18)$$\begin{gather}\frac{L}{d_p} \end{gather}$$

(5.19)$$\begin{gather}e_c . \end{gather}$$

This non-dimensionalization process has revealed that the original system of equations can be expressed as a function of seven independent dimensionless parameters: $Pe$, $\tau _c^{*}$, $\tau _p^{*}$, $\tau _{\sigma }^*$, ${L}/{d_p}$, ${u_e}/{u_k}$ and $e_c$. The dimensionless parameters $\tau _{\xi }^{*}$ and $\eta _{coll}$ can be written as a function of the previous dimensionless parameters

(5.20)$$\begin{gather} \tau_{\xi}^{*} = \lambda^{\prime} Pe \left( \frac{d_p}{L} \right)^{2} \end{gather}$$

(5.21)$$\begin{gather}\eta_{coll} = \lambda^{\prime\prime} \frac{1}{Pe}\frac{L}{d_p}, \end{gather}$$

with $\lambda ^{\prime } = 12$ and $\lambda ^{\prime \prime } \approx 3.256$.

The dimensionless analysis of this simple system has allowed us to highlight some of the characteristics of the equations derived in the previous sections. First of all, the dimensionless mean electric charge transport equations are mainly controlled by 2 parameters: $Pe$ and $\tau _{\sigma }^{*}$. The first one corresponds to the inverse dimensionless dispersion coefficient while, $\tau _{\sigma }^*$ characterizes the strength of the triboconductivity effect. It is also worth noting that, in our case, the dimensionless collision time ($\tau _c^*$) is also equivalent to a Knudsen number, as it represents the ratio between a collisional mean free path and the system characteristic length scale. The dimensionless covariance transport equations show that the relaxation characteristic time of this variable is a function of the collision time ($\tau _c^*$), the particle relaxation time ($\tau _p^*$) and the variance relaxation time ($\tau _{\xi }^*$), and this relaxation characteristic time is guaranteed to be at least as small as the smallest of the three previously mentioned times. Finally, analysing the dimensionless variance transport equation, we note that its relaxation characteristic time will always be greater than $\tau _{\xi }^*$. This point is important because it shows that the charge variance will always have a larger characteristic relaxation time than the charge–velocity covariance.

These equations allow us to derive the conditions in which some mechanisms are more dominant than others. For example, one of the most important distinctions is the difference between the dense or collisional regime and the dilute or kinetic regime. The former is characterized for low kinetic transport and high collisional charge transfer and, in contrast, the latter is associated with high kinetic transport and low collisional charge transfer. Looking at the dimensionless form of the transport equation, a dense regime can be achieved in configurations where the product between $Pe$ and the covariance relaxation time is small. Conversely, kinetic regimes correspond to configuration where this product is large.

7.1. Full coupled algebraic model

The simplest approach could be to try to extend the algebraic model proposed by Montilla et al. (Reference Montilla, Ansart and Simonin2020) to take into account the electric charge variance. This will reduce the model to a single partial differential equation with two algebraic expressions for the charge–velocity covariance vector and the charge variance. In order to obtain this more general algebraic model, we can follow the methodology used to derive the previous algebraic model. First, we simplify the charge–velocity covariance and the charge variance equations neglecting the contributions of the Lagrangian derivative, the velocity gradient, the third-order moments and the dispersion terms. This reduces the covariance transport equation (3.1) and the variance transport equation (4.1) to the following system of algebraic equations:

(7.1) \begin{equation} \left\{\begin{array}{@{}l} D_1 \langle c_{p,i}^{\prime} \xi_p^{\prime} \rangle + D_{2,i} \langle \xi_p^{\prime} \xi_p^{\prime} \rangle = D_{3,i} \\ F_{1,i} \langle c_{p,i}^{\prime} \xi_p^{\prime} \rangle + F_2 \langle \xi_p^{\prime} \xi_p^{\prime} \rangle = F_3 \end{array} \right. \end{equation}

where

(7.2)$$\begin{gather} D_1 = \frac{1+e_c}{3}\frac{1}{\tau_c} + \frac{1}{\tau_p} + \frac{(3-e_c)}{5} \frac{1}{\tau_{\xi}} \end{gather}$$

(7.3)$$\begin{gather}D_{2,i} ={-}\frac{n_p}{m_p} E_i \end{gather}$$

(7.4)$$\begin{gather}D_{3,i} ={-}\left( 1 + \lambda_{2.2} e_c\frac{\sqrt{\varTheta_p}}{d_p} D_p \right)n_p\frac{\partial Q_p}{\partial x_i} + \lambda_{2.1} e_c \frac{\sqrt{\varTheta_p}}{d_p}\sigma_p E_i \end{gather}$$

(7.5)$$\begin{gather}F_{1,i} = 2n_p\frac{\partial Q_p}{\partial x_i} \end{gather}$$

(7.6)$$\begin{gather}F_{2} = n_p \left( \frac{1}{\tau_{\xi}} -\lambda_{3.1}(D_p)^2 \frac{\tau_c}{d_p^4} \right) \end{gather}$$

(7.7)$$\begin{gather}F_{3} = \lambda_{3.2}(\sigma_p)^2\frac{\tau_c}{n_p d_p^2} E_iE_i. \end{gather}$$

This system of equations can be solved to find an algebraic model coupling the particle charge–velocity covariance and the charge variance

(7.8)$$\begin{gather} \langle c_{p,i}^{\prime} \xi_p^{\prime} \rangle = \frac{D_{3,i}}{D_{1}} - \frac{D_{2,i}}{D_1}\langle \xi_p^{\prime} \xi_p^{\prime} \rangle, \end{gather}$$

(7.9)$$\begin{gather}\langle \xi_p^{\prime} \xi_p^{\prime} \rangle = \frac{F_3D_1 - F_{1,i}D_{3,i}}{D_1F_2 - D_{2,j}F_{1,j}}. \end{gather}$$

Equations (7.8) and (7.9) give us a set of algebraic expressions that can be used to close the mean particle electric charge transport equation. With these expressions, we no longer have to solve a system of multiple partial differential equations (PDEs). We would also like to remark that this model reduces to the aforementioned simpler algebraic model of Montilla et al. (Reference Montilla, Ansart and Simonin2020) if $\langle \xi _p^{\prime } \xi _p^{\prime } \rangle = 0$. So this model is indeed a more general algebraic model including the effects of the charge variance.

However, we would like to note that this algebraic closure law could produce non-physical results. In particular, the positive sign in (7.9) cannot be guaranteed. In addition to this, the denominator of this equation could also approach zero, leading to indeterminate variance values in some regions. A clear example of this can be shown if we look at the variance prediction for the initial condition of the two cases studied earlier (figure 6a,b). These figures show that this coupled algebraic model could lead to an unphysical prediction even for valid initial conditions. Given this behaviour, we cannot use this approach to simulate these two test cases. In our tests, we remarked that this problem is more prone to occur when the variance term is not negligible.

Figure 6.

Initial variance profile predicted by the coupled algebraic model for the two cases studied; (a) $\tau _c^*=10^{-4}$, (b) $\tau _c^*=10^{-1}$.

7.2. Semi-algebraic model

Another possible reduced-order model can be derived by only simplifying one of the second-order moment transport equations. Here, we choose to simplify the charge–velocity covariance transport equation to an algebraic expression. Therefore, the system is governed by the mean charge transport equation (2.24), the variance transport equation (4.1) and the algebraic simplification for the covariance (7.8). The reason behind this choice is that, as shown in the dimensionless analysis, the variance characteristic relaxation time will always be larger than the covariance characteristic relaxation time. Therefore, the quasi-steady state hypothesis can be reached first for the covariance but it is harder to satisfy for the variance.

We tested this approach for the same two configurations we have been studying in this paper. The first thing we would like to point out is that, unlike the previous model, this semi-algebraic model produced meaningful physical results at all times. Figures 7 and 8 show a comparison between this semi-algebraic model and the solution of the full second-order transport equation model. We observe that, in the smaller collision time configuration, the proposed semi-algebraic model matches perfectly the result obtained by solving the full model. However, for the larger collision time configuration, the semi-algebraic model fails to correctly reproduce the dynamics.

Figure 7.

Comparison between the predictions of the full transport equation model and the semi-algebraic model in the test case with $\tau _c^*=10^{-4}$. (a) Dimensionless mean particle electric charge. (b) Dimensionless electric field. (c) Dimensionless particle charge–velocity covariance. (d) Dimensionless particle charge variance.

Figure 8.

Comparison between the predictions of the full transport equations model and the semi-algebraic model in the test case $\tau _c^*=10^{-1}$. (a) Dimensionless mean particle electric charge. (b) Dimensionless electric field. (c) Dimensionless particle charge–velocity covariance. (d) Dimensionless particle charge variance.

In order to explain the success and failure of this model, we need to look at the hypotheses taken to reduce the covariance transport equation to an algebraic expression. The algebraic expression is obtained by assuming that the temporal derivative and any dispersion term can be neglected. Therefore, the semi-algebraic model should produce good results when these two assumptions hold simultaneously. This is exactly the case for the small collision time configuration. The term-by-term analysis of the covariance terms (figure 4b) showed that the transient and the kinetic dispersion terms are negligible compared with the production and destruction terms. In this case, the small collision time imposes also a very small relaxation characteristic time for the covariance, which ensures the quasi-steady state, corresponding to equilibrium between production and destruction terms. Also, the covariance dispersion coefficient is limited by the smallest of the particle collision time, the particle relaxation time and the variance destruction time. In our case, the very small particle collision time ensures a very weak covariance dispersion coefficient that could be neglected.

However, in the larger collision time configuration, the semi-algebraic model is not a suitable alternative. In this case, the assumptions mentioned above are no longer valid and therefore none of the terms in the equation can be neglected. The large collision time also means that the covariance characteristic destruction time is very large which means that no quasi-steady state can be reached. Also, the kinetic dispersion term is stronger as the collision times becomes larger. This is confirmed by looking at the term-by-term analysis of this case (figure 5b). We can clearly see that both the dispersion and temporal derivative terms are relevant to correctly characterize the dynamics of the covariance. In our example, we see that the semi-algebraic model has overestimated the mean charge transport. This is due to the overestimation of the kinetic transport imposed by the quasi-steady-state hypothesis.

Given these results, we wanted to determine the validity region in which the semi-algebraic model is a valid alternative to model the dynamics of a gas flow system with electrostatic charge. In order to do this, we performed several simulations varying the different dimensionless numbers to see the impact on the quality of the semi-algebraic model. According to our results, two of the most important parameters to determine the suitability of the semi-algebraic model are $\tau _c^*$ and $Pe$. Using our simulations, we were able to distinguish the conditions in which the semi-algebraic model is a valid approach to model these kinds of systems, as shown by figure 9. We have also marked the two test cases we have studied in this article. As we can see, the boundary is a nonlinear curve. For high $Pe$ values we only need to ensure that the collision time is small enough to guarantee rapid relaxation towards the steady state and a weak kinetic dispersion for the charge covariance. However, as $Pe$ decreases, the electric charge transport due to collision increases. So the collision time has to be small enough to continue to ensure the quasi-steady-state hypothesis in the charge–velocity covariance transport equation.

Figure 9.

Suitability regions of the semi-algebraic model as a function $Pe$ and $\tau _c^*$. The two marked points correspond to the test cases studied in this work.