2 Plasma parameters

The focus of the following is a dusty plasma which, for the sake of simplicity, is assumed to consist of the free electrons with number density $n_{e}$ and electric charge $-e$ , protons with number density $n_{p}$ and electric charge $e$ and dust particles with number density $n_{d}$ and electric charge $-Z_{d}e$ . The presence of neutrals is completely ignored everywhere below since of principal interest is the pure electrostatic interaction between the dust particles, merged into a plasma medium.

To characterize the state of the dusty plasma it is very helpful to introduce a few dimensionless parameters such as the density-ratio parameter

(2.1) $$\begin{eqnarray}\unicode[STIX]{x1D6FD}=\frac{n_{d}}{n_{p}}\end{eqnarray}$$

and the non-isothermicity coefficient

(2.2) $$\begin{eqnarray}\unicode[STIX]{x1D70F}=\frac{T_{e}}{T_{p}}.\end{eqnarray}$$

Here $T_{e}$ and $T_{p}$ stand for the temperatures of the electrons and protons, respectively. The temperature of the dust grains $T_{d}$ is known to significantly vary in magnitude (Avinash, Merlino & Shukla Reference Avinash, Merlino and Shukla2011) but for definiteness it is assumed to be equal to the temperature of the electrons, $T_{d}=T_{e}$ , although the generalization to other practical situations is rather straightforward. As a matter of fact, the temperature of the dust grain matter remains almost equal to that of the neutral gas of the buffer plasma and what is referred to as $T_{d}$ is a quantity that simply determines the average chaotic kinetic energy of the dust particles caused either by the thermally fluctuating electric field in the plasma or by the fluctuating electric charge of dust grains themselves (Quinn & Goree Reference Quinn and Goree2000). Note that at rather low neutral gas pressures ( ${<}200$ mTorr), $T_{d}$ normally lies in the range of 10–300 eV as experimentally evidenced by stereoscopic particle image velocimetry (Fisher & Thomas Reference Fisher and Thomas2010).

To describe the state of the buffer plasma the coupling parameter for the protons is defined as

(2.3) $$\begin{eqnarray}\unicode[STIX]{x1D6E4}_{p}=\frac{e^{2}}{a_{p}k_{B}T_{p}},\end{eqnarray}$$

where $a_{p}=(3/4\unicode[STIX]{x03C0}n_{p})^{1/3}$ denotes the proton sphere radius and $k_{B}$ signifies the Boltzmann constant.

Herein, the finite-size effects of dust particles are under scrutiny, which necessitates a standard definition of the packing fraction as follows

(2.4) $$\begin{eqnarray}\unicode[STIX]{x1D702}={\textstyle \frac{4}{3}}\unicode[STIX]{x03C0}n_{d}R^{3},\end{eqnarray}$$

where $R$ designates the radius of the dust particles, which are viewed as hard spherical balls. It is well established in the literature that the packing fraction cannot exceed the value of $\sqrt{2}\unicode[STIX]{x03C0}/6\approx 0.74$ for the most compact packing of hard balls of the same size.

It is essential for further consideration that the local quasi-neutrality condition is imposed in the form of

(2.5) $$\begin{eqnarray}n_{p}=n_{e}+Z_{d}n_{d}.\end{eqnarray}$$

It has to be strictly emphasized that the electric charge of the dust grains $Z_{d}$ is no longer an independent parameter of the present model. Moreover, all further quantities under study turn out to be functions of the parameters (2.1)–(2.4) defined above.

In closing this section, it is obligatory to mention another two parameters, which are widely used to solely describe the state of the dust component. They are the dust coupling parameter

(2.6) $$\begin{eqnarray}\unicode[STIX]{x1D6E4}_{d}=\frac{Z_{d}^{2}e^{2}}{a_{d}k_{B}T_{d}},\end{eqnarray}$$

and the screening parameter

(2.7) $$\begin{eqnarray}\unicode[STIX]{x1D705}=\frac{a_{d}}{r_{D}}.\end{eqnarray}$$

Here $a_{d}=(3/4\unicode[STIX]{x03C0}n_{d})^{1/3}$ represents the Wigner–Seitz radius of dust particles and $k_{D}=r_{D}^{-1}=(4\unicode[STIX]{x03C0}n_{e}e^{2}/k_{B}T_{e}+4\unicode[STIX]{x03C0}n_{p}e^{2}/k_{B}T_{p})^{1/2}$ denotes the Debye wavenumber inverse to the Debye screening radius $r_{D}$ .

Note that the dust coupling characterizes the non-ideality effects in the intergrain interactions, while the screening parameter estimates the contribution of shielding due to the buffer plasma. However, they both play only a subsidiary role in the following and are needed only for the sake of comparison with the results of other approaches.

3 Dust grain charge

It was recently demonstrated that the plasma electrodynamics in the framework of the density-response formalism provides the following interaction model of dusty plasma particles (Davletov et al. Reference Davletov, Arkhipov and Tkachenko2016):

(3.1) $$\begin{eqnarray}\unicode[STIX]{x1D6F7}_{ab}(r)=\frac{Q_{ab}}{r+R_{ab}}-\frac{Q_{ab}}{r}\left(1-\exp (-rk_{D})-\frac{R_{ab}k_{D}}{2}\text{B}_{ab}(r)\right),\end{eqnarray}$$

where the subscripts $a,b$ take on values of $e$ for electrons, $p$ for protons and $d$ for dust grains, $Q_{ed}=-Q_{pd}=Z_{d}e^{2}$ , $R_{ed}=R_{pd}=R$ and $Q_{dd}=Z_{d}^{2}e^{2}$ , $R_{dd}=2R$ , and

(3.2) $$\begin{eqnarray}\displaystyle \text{B}_{ab}(r) & = & \displaystyle \exp ((R_{ab}+r)k_{D})\text{Ei}((R_{ab}+r)k_{D})-\exp (k_{D}(R_{ab}-r))\text{Ei}(k_{D}R_{ab})\nonumber\\ \displaystyle & & \displaystyle +\exp (-(R_{ab}+r)k_{D})\left[\text{Ei}(-R_{ab}k_{D})-\text{Ei}(-(R_{ab}+r)k_{D})\right],\end{eqnarray}$$

with the exponential integral function $\text{Ei}(x)=\int _{x}^{\infty }\exp (-t)/t\,\text{d}t$ . In formulas (3.1) and (3.2) the distance $r$ is measured between the surfaces of the particles.

The interaction model (3.1) is designed to account for the finite size of dust particles by engaging the parameter $k_{D}R$ which was first proposed by Whipple (Reference Whipple1981). The parameter $k_{D}R$ is known to vary in quite a broad range from very small values up to dozens and the Whipple approximation and its generalizations are widely used in studying various aspects of dusty plasma physics (Tang & Delzanno Reference Tang and Delzanno2014; Delzanno & Tang Reference Delzanno and Tang2015; Momot et al. Reference Momot, Zagorodny and Orel2017).

It is worth mentioning that the Yukawa (Debye–Hückel) interaction potential, $\unicode[STIX]{x1D6F7}_{ab}(r)=Q_{ab}\exp (-k_{D}r)/(r+R_{ab})$ , is simply recovered as a limiting case of expressions (3.1) and (3.2) when $k_{D}R_{ab}\ll 1$ , i.e. the dust particles are quite small. Moreover, it has to be kept in mind that, for the same values of the dust coupling and the screening parameters, the Yukawa potential always underestimates the interparticle interaction potential as compared to formulas (3.1) and (3.2).

Another important limiting case is that of big dust particles $k_{D}R_{ab}\gg 1$ when expressions (3.1) and (3.2) simplify to

(3.3) $$\begin{eqnarray}\unicode[STIX]{x1D6F7}_{ab}(r)=\frac{2Q_{ab}R_{ab}}{k_{D}^{2}r(r+R_{ab})^{3}}-\frac{2Q_{ab}}{k_{D}^{2}R_{ab}^{2}r}\exp (-rk_{D}).\end{eqnarray}$$

It is very interesting to stress that, at large interparticle distances, the exponential decay of the interaction potential (3.1) as well as (3.2) for small dust grains is replaced by the inverse power-law decay for large dust grains, as predicted by the first term on the right-hand side of formula (3.3).

In the remainder of this section the charging process is properly considered on the basis of the orbital motion limited approximation, in which trajectories of plasma particles, say electrons and protons, are considered ballistic such that the interparticle collisions are completely omitted. It is implicitly presumed in the subsequent analysis of the dust particle charging that the material of the dust is a perfect absorber. Under those simplified assumptions, simultaneous application of conservation laws of energy and angular momentum is sufficient to deduce the corresponding absorption cross-sections, whose averaging over the Maxwellian velocity distribution function allows one to derive the corresponding fluxes of plasma particles on the dust grain surface. That the dust grain retains its electric charge over time requires the electron and proton fluxes to be equal, which yields the following general relation, valid for any kind of the interaction potential:

(3.4) $$\begin{eqnarray}\frac{n_{p}}{n_{e}}\sqrt{\frac{m_{e}T_{p}}{m_{p}T_{e}}}\left(1-\frac{\unicode[STIX]{x1D6F7}_{pd}(0)}{k_{B}T_{p}}\right)=\exp \left(-\frac{\unicode[STIX]{x1D6F7}_{ed}(0)}{k_{B}T_{e}}\right).\end{eqnarray}$$

In the case of the Coulomb interaction, $\unicode[STIX]{x1D6F7}_{ab}(r)=Q_{ab}/(r+R_{ab})$ , expression (3.4) is reduced to the well-known equation of the form (Fortov et al. Reference Fortov, Khrapak, Khrapak, Molotkov and Petrov2004)

(3.5) $$\begin{eqnarray}\frac{n_{p}}{n_{e}}\sqrt{\frac{m_{e}T_{p}}{m_{p}T_{e}}}\left(1+\frac{Z_{d}e^{2}}{k_{B}T_{p}R}\right)=\exp \left(-\frac{Z_{d}e^{2}}{k_{B}T_{e}R}\right).\end{eqnarray}$$

For the Yukawa potential with the Whipple approximation (Whipple Reference Whipple1981), $\unicode[STIX]{x1D6F7}_{ab}(r)=Q_{ab}\exp (-k_{D}r)/(1+k_{D}R)(r+R_{ab})$ , expression (3.4) gives rise to

(3.6) $$\begin{eqnarray}\frac{n_{p}}{n_{e}}\sqrt{\frac{m_{e}T_{p}}{m_{p}T_{e}}}\left(1+\frac{Z_{d}e^{2}}{k_{B}T_{p}(1+k_{D}R)R}\right)=\exp \left(-\frac{Z_{d}e^{2}}{k_{B}T_{e}(1+k_{D}R)R}\right)\end{eqnarray}$$

whereas for the interaction potential (3.1), the following equation is obtained

(3.7) $$\begin{eqnarray}\displaystyle & & \displaystyle \frac{n_{p}}{n_{e}}\sqrt{\frac{m_{e}T_{p}}{m_{p}T_{e}}}\left(1+\frac{Z_{d}e^{2}}{k_{B}T_{p}R}\left[1-k_{D}R\left(1+k_{D}R\exp (k_{D}R)\text{Ei}(-k_{D}R)\right)\right]\right)\nonumber\\ \displaystyle & & \displaystyle \quad =\exp \left(-\frac{Z_{d}e^{2}}{k_{B}T_{e}R}\left[1-k_{D}R\left(1+k_{D}R\exp (k_{D}R)\text{Ei}(-k_{D}R)\right)\right]\right).\end{eqnarray}$$

In figure 1 the dust grain charge number $Z_{d}$ is drawn against the non-isothermicity parameter $\unicode[STIX]{x1D70F}$ at fixed values of the density-ratio parameter $\unicode[STIX]{x1D6FD}=1.00\times 10^{-4}$ , the proton coupling $\unicode[STIX]{x1D6E4}_{p}=0.01$ and the packing fraction $\unicode[STIX]{x1D702}=0.1$ . A comparison is provided between the solutions to (3.5)–(3.7) to demonstrate that the grain charge within the Coulomb potential is always underestimated while it is permanently overestimated within the Yukawa potential in comparison with the proposed expression (3.7). The former is prescribed to account for the screening phenomenon whereas the latter is a straightforward consequence of the application of the plasma electrodynamics. It is seen from the same figure 1 that, for all three cases, the dust grain charge number monotonically grows with increasing $\unicode[STIX]{x1D70F}$ and tends to its critical value of $Z_{d}=1/\unicode[STIX]{x1D6FD}$ , as assured by the quasi-neutrality condition (2.5). This is not surprising because the number of protons is kept constant by fixing the density-ratio parameter $\unicode[STIX]{x1D6FD}$ and the negative charge is somehow allocated between the electrons and the dust grains. Then, an increase in $\unicode[STIX]{x1D70F}$ results in significant growth of the electron flux on the grain surface, whose electric charge thus inevitably grows.

Figure 1.

Dust grain charge number $Z_{d}$ as a function of the non-isothermicity parameter $\unicode[STIX]{x1D70F}$ at the fixed values of $\unicode[STIX]{x1D6FD}=1.00\times 10^{-4}$ , $\unicode[STIX]{x1D6E4}_{p}=0.01$ and $\unicode[STIX]{x1D702}=0.1$ . Red line: solution to (3.5); green line: solution to (3.6); blue line: solution to (3.7).

Figure 2 portrays the dependence of the dust grain charge number $Z_{d}$ on the logarithm of the density-ratio parameter $\unicode[STIX]{x1D6FD}$ at fixed values of the non-isothermicity parameter $\unicode[STIX]{x1D70F}=1$ , the proton coupling $\unicode[STIX]{x1D6E4}_{p}=0.01$ and the packing fraction $\unicode[STIX]{x1D702}=0.1$ . It is obvious from the physical point of view that the increase in the number density of dust grains gives rise to a decrease of the dust grain charge number since the corresponding number of absorbing centres grows, which is clearly proved by figure 2.

Figure 2.

Dust grain charge number $Z_{d}$ as a function of the logarithm of the density-ratio parameter $\unicode[STIX]{x1D6FD}$ at fixed values of $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6E4}_{p}=0.01$ and $\unicode[STIX]{x1D702}=0.1$ . Red line: solution to (3.5); green line: solution to (3.6); blue line: solution to (3.7).

Figure 3 is specifically shown to demonstrate the influence that the finite-size effects exert on the electric charge of the dust grains. Once again, it clearly shows that the grain charge number is underestimated within the Coulomb interaction and overestimated within the Yukawa model and in all cases the charge number grows with $\unicode[STIX]{x1D702}$ because of the resulting increase in the dust particle radius.

Figure 3.

Dust grain charge number $Z_{d}$ as a function of the packing fraction $\unicode[STIX]{x1D702}$ at the fixed values of $\unicode[STIX]{x1D6FD}=1.00\times 10^{-4}$ , $\unicode[STIX]{x1D6E4}_{p}=0.01$ and $\unicode[STIX]{x1D70F}=1$ . Red line: solution to (3.5); green line: solution to (3.6); blue line: solution to (3.7).

4 Thermodynamics and dust-acoustic waves

It is well known from the statistical physics of equilibrium systems that the complete many-body distribution function is exactly expressed by the Gibbs distribution. It still, however, carries a great deal of unnecessary information, which is not needed in practice. From the point of view of calculating the thermodynamic characteristics of the medium, it is sufficient to know only the pair correlation function $h(\boldsymbol{r})$ or the corresponding radial distribution function $g(\boldsymbol{r})$ , which simply represents the probability density of finding two particles at a certain distance $\boldsymbol{r}$ from each other.

For a system of single particle species, interacting via the potential $\unicode[STIX]{x1D6F7}(\boldsymbol{r})$ , the direct correlation function $c(\boldsymbol{r})$ is expressed in terms of the pair correlation function $h(\boldsymbol{r})$ with the help of the Ornstein–Zernike relation

(4.1) $$\begin{eqnarray}h(\boldsymbol{r})=c(\boldsymbol{r})+n\int c(\boldsymbol{r}-\boldsymbol{r}^{\prime })h(\boldsymbol{r}^{\prime })\,\text{d}\boldsymbol{r}^{\prime },\end{eqnarray}$$

where

(4.2) $$\begin{eqnarray}h(\boldsymbol{r})=g(\boldsymbol{r})-1\end{eqnarray}$$

with $n$ being the particle number density.

Equations (4.1) and (4.2) stay unclosed from the viewpoint of mathematics since they still contain three unknown functions. For neutral systems like ordinary liquids, an important closure is provided by the Percus–Yevick relation

(4.3) $$\begin{eqnarray}c(\boldsymbol{r})=g(\boldsymbol{r})\left[1-\exp \left(\frac{\unicode[STIX]{x1D6F7}(\boldsymbol{r})}{k_{B}T}\right)\right],\end{eqnarray}$$

where $T$ stands for the particle temperature.

The set of equations (4.1)–(4.3) was especially successful in describing various properties of hard sphere (HS) systems since it can be explicitly solved for a wide range of packing fraction and a very good agreement was found with the results of numerical simulations.

Another general expression, appropriate for charged systems at any coupling, appears in the diagrammatic expansion in the following form

(4.4) $$\begin{eqnarray}g(\boldsymbol{r})=\exp \left(-\frac{\unicode[STIX]{x1D6F7}(\boldsymbol{r})}{k_{B}T}+h(\boldsymbol{r})-c(\boldsymbol{r})+B(\boldsymbol{r})\right),\end{eqnarray}$$

where $B(\boldsymbol{r})$ denotes the bridge function.

Equation (4.4) is of no help unless it is directly stated how to calculate the bridge function $B(\boldsymbol{r})$ . The simplest possible method is the so-called hypernetted-chain (HNC) approximation, which completely drops the bridge function by assuming

(4.5) $$\begin{eqnarray}B_{\text{HNC}}(\boldsymbol{r})=0.\end{eqnarray}$$

In this case, the set of equations (4.1), (4.2), (4.4) and (4.5) can be solved numerically by using an iterative procedure and straightforward comparison with the results of numerical simulations shows that the HNC approximation works very well for weakly and moderately coupled systems like one-component plasmas. For a strongly coupled regime, however, a discrepancy with Monte Carlo data was discovered for both the radial distribution function and the static structure factor, which provoked further studies on the topic.

On the one hand, it was immediately realized that for very strong couplings the shape of the radial distribution function was practically reminiscent of that observed for a HS system. On the other hand, the tremendous success of analytical predictions for a HS model stimulated further advance by invoking a number of extensions to systems with an arbitrary interaction potential. The original idea was to incorporate the hard sphere system as a reference by presenting a real interaction potential as a sum of two parts such that the interaction potential varies from that of the reference system to that of the system under investigation by adjusting a free parameter. All these techniques somehow rely on the supposition that the bridge function is weakly dependent on the interaction potential and can, therefore, be taken from the solution of the Percus–Yevick relation for hard spheres (Wertheim Reference Wertheim1963). Among most fruitful approaches that deserve to be mentioned here are the perturbative HNC (PHNC) (Kang & Ree Reference Kang and Ree1995), the reference HNC (RHNC) (Lado Reference Lado1973, Reference Lado1976), the modified HNC (MHNC) (Lado Reference Lado1982; Lado, Foiles & Ashcroft Reference Lado, Foiles and Ashcroft1983) and, finally, the variationally modified HNC (VMHNC) (Rosenfeld & Ashcroft Reference Rosenfeld and Ashcroft1979; Rosenfeld Reference Rosenfeld1986; Faussurier Reference Faussurier2004). All of which are primarily designed for systems with point-like particles and an adjustable packing fraction of a HS reference model is chosen to satisfy some additional conditions. For instance, in MHNC, the hard sphere diameter is chosen to ensure thermodynamic consistency, whereas, in VMHNC, the packing fraction serves as a parameter for minimizing the system free energy.

Figure 4.

Radial distribution function as a function of the normalized intergrain spacing $r/a_{d}$ at fixed values of $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=1.00\times 10^{-4}$ , $\unicode[STIX]{x1D6E4}_{p}=0.001$ and $\unicode[STIX]{x1D702}=1.00\times 10^{-6}$ , which correspond to $\unicode[STIX]{x1D6E4}_{d}=10.6$ and $\unicode[STIX]{x1D705}=1.5$ . Red line: HNC scheme with the bridge correction for point-like particles; blue line: RHNC approximation for particles with finite dimensions.

Figure 5.

Radial distribution function as a function of the normalized intergrain spacing $r/a_{d}$ at fixed values of $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=0.061$ , $\unicode[STIX]{x1D6E4}_{p}=0.109$ and $\unicode[STIX]{x1D702}=0.3$ , which correspond to $\unicode[STIX]{x1D6E4}_{d}=10.05$ and $\unicode[STIX]{x1D705}=1.5$ . Red line: HNC scheme with the bridge correction for point-like particles; blue line: RHNC approximation for particles with finite dimensions.

Figure 6.

Radial distribution function as a function of the normalized intergrain spacing $r/a_{d}$ at the fixed values of $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=1.00\times 10^{-5}$ , $\unicode[STIX]{x1D6E4}_{p}=2.20\times 10^{-4}$ and $\unicode[STIX]{x1D702}=1.00\times 10^{-6}$ , which correspond to $\unicode[STIX]{x1D6E4}_{d}=132$ and $\unicode[STIX]{x1D705}=1.66$ . Red line: HNC scheme with the bridge correction for point-like particles; blue line: RHNC approximation for particles with finite dimensions.

Contrarily, it is crucial for the present consideration that the dust grains have their own size, hence, the packing fraction cannot be treated as an adjustable parameter. The aforementioned logically necessitates further engagement of the RHNC calculation scheme with the fixed packing fraction. Thus, figures 4–7 show a comparison of the radial distribution function obtained within the RHNC method for the finite-size dust particles interacting via potential (3.1) with that derived from the HNC scheme with the bridge correction (Daughton, Murillo & Thode Reference Daughton, Murillo and Thode2000) for point-like particles interacting via the Yukawa potential. In-depth analysis of the obtained graphic data allows us to draw the following conclusions. For rather small values of the dust coupling parameter and low packing fractions, both methods give practically the same results, as figure 4 clearly displays. Increase in the packing fraction results in strong oscillations of the radial distribution function obtained within the RHNC method, whereas a single weakly expressed maximum is observed for the HNC scheme, as evidenced by figure 5. In the strong dust coupling regime and for rather low packing fractions, the above two methods differ significantly from each other, see figure 6, namely, the HNC scheme gives a higher first peak located to the right of that in the RHNC method. It is interesting to note that, at high packing fractions and strong dust couplings, the opposite picture is observed, as exemplified by figure 7. It is also seen from figures 5 and 7 that there is a discontinuity in the radial distribution function at high values of the packing fraction, which is a rudiment of the HS model involved into the RHNC scheme.

Figure 7.

Radial distribution function as a function of the normalized intergrain spacing $r/a_{d}$ at the fixed values of $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=8.00\times 10^{-3}$ , $\unicode[STIX]{x1D6E4}_{p}=0.05$ and $\unicode[STIX]{x1D702}=0.3$ , which correspond to $\unicode[STIX]{x1D6E4}_{d}=124$ and $\unicode[STIX]{x1D705}=2.03$ . Red line: HNC scheme with the bridge correction for point-like particles; blue line: RHNC approximation for particles with finite dimensions.

It has to be remarked at this point that, in an ordinary dusty plasma encountered in real experimental installations, the packing fraction is usually very small such that the HNC approximation retains its credibility. However, this is not the case for colloidal suspensions, to which the constructed theoretical approach is undoubtedly applicable.

The main thermodynamic quantity studied below is the pressure $P$ , which is written below in reduced units of the excess pressure $p_{ex}=P/n_{d}k_{B}T-1$ . It is well known from the statistical physics of many-body systems that the reduced excess pressure is expressed in terms of the radial distribution function $g(r)$ as

(4.6) $$\begin{eqnarray}p_{ex}=-\frac{2\unicode[STIX]{x03C0}n_{d}}{3k_{B}T}\int _{0}^{\infty }r^{3}\unicode[STIX]{x1D6F7}_{dd}^{\prime }(r)g(r)\,\text{d}r.\end{eqnarray}$$

In table 1 the isothermal compressibility factor $Z=1+p_{ex}$ is shown as $Z_{\text{present}}$ calculated within the RHNC scheme for $\unicode[STIX]{x1D702}=1.00\times 10^{-15}$ to formally disregard finite-size effects. In the same table the comparison is made with the data of the Monte Carlo simulations (Meijer & Frenkel Reference Meijer and Frenkel1991) and the results of the approach proposed by Tejero et al. (Reference Tejero, Lutsko, Colot and Baus1992), which also benefits from merits of the HNC and Percus–Yevick integral equation techniques. A very good agreement is observed, even for large values of the dust coupling, which completely verifies the validity of the present approach for the point-like dust grains corresponding to vanishing packing fractions. In figure 8 comparison of the excess pressure at $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=1.00\times 10^{-3}$ and $\unicode[STIX]{x1D6E4}_{p}=0.001$ is provided with the practical expressions (Khrapak et al. Reference Khrapak, Khrapak, Ivlev and Morfill2014; Khrapak & Thomas Reference Khrapak and Thomas2015b ) obtained for point-like dust grains. It is visibly inferred from figure 8 that, when the packing fraction grows, the hard sphere contribution begins to dominate over the like-charge repulsion. Namely, for a packing fraction of $\unicode[STIX]{x1D702}\sim 10^{-3}$ the correction to the excess pressure due to the finite-size effects reaches the value of almost 20 %.

Table 1.

Compressibility factor $Z=1+p_{ex}$ for a wide range of the dust coupling and screening parameters. $Z_{\text{MC}}$ : Monte Carlo simulations (Meijer & Frenkel Reference Meijer and Frenkel1991);  $Z_{\text{DRY}}$ : DRY approach by (Tejero et al. Reference Tejero, Lutsko, Colot and Baus1992);  $Z_{\text{present}}$ : present results for $\unicode[STIX]{x1D702}=1.00\times 10^{-15}$ .

Figure 8.

Excess pressure (4.6) as a function of the logarithm of the packing fraction $\log _{10}\unicode[STIX]{x1D702}$ at fixed values of $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=1.00\times 10^{-3}$ and $\unicode[STIX]{x1D6E4}_{p}=0.001$ . Blue line: results of the present approach; red line: practical expressions (Khrapak et al. Reference Khrapak, Khrapak, Ivlev and Morfill2014; Khrapak & Thomas Reference Khrapak and Thomas2015b ).

Finally, the dust-acoustic wave dispersion is studied below in the strongly coupled regime. Various approaches have been developed here (Rosenberg & Kalman Reference Rosenberg and Kalman1997; Kaw & Sen Reference Kaw and Sen1998; Murillo Reference Murillo1998; Filippov et al. Reference Filippov, Starostin, Tkachenko, Fortov, Ballester and Conde2010, Reference Filippov, Starostin, Tkachenko and Fortov2011) but the simplest hydrodynamic model of the dust component with a neutralizing background of electrons and ions yields (Khrapak et al. Reference Khrapak, Khrapak, Ivlev and Morfill2014)

(4.7) $$\begin{eqnarray}\frac{\unicode[STIX]{x1D714}^{2}}{\unicode[STIX]{x1D714}_{p}^{2}}=\frac{q^{2}}{q^{2}+\unicode[STIX]{x1D705}^{2}}+\frac{q^{2}}{3\unicode[STIX]{x1D6E4}_{d}}\unicode[STIX]{x1D6FE}\unicode[STIX]{x03BC}_{p},\end{eqnarray}$$

where $\unicode[STIX]{x1D714}_{p}^{2}=4\unicode[STIX]{x03C0}Z_{d}^{2}e^{2}n_{d}/m_{d}$ designates the dust plasma frequency with $m_{d}$ being the mass of the dust particles, $\unicode[STIX]{x1D707}_{p}=(1/k_{B}T)(\unicode[STIX]{x2202}P/\unicode[STIX]{x2202}n_{d})_{T}$ is the inverse reduced isothermal compressibility, $q=ka_{d}$ stands for the dimensionless wavenumber, $\unicode[STIX]{x1D6FE}=C_{p}/C_{V}$ refers to the adiabatic index.

In its turn, the inverse reduced isothermal compressibility is expressed via the reduced excess pressure with a neutralizing background

(4.8) $$\begin{eqnarray}p_{ex}^{0}=-\frac{2\unicode[STIX]{x03C0}n_{d}}{3k_{B}T}\int _{0}^{\infty }r^{3}\unicode[STIX]{x1D6F7}_{dd}^{\prime }(r)\left(g(r)-1\right)\,\text{d}r,\end{eqnarray}$$

in the following form

(4.9) $$\begin{eqnarray}\unicode[STIX]{x1D707}_{p}=1+p_{ex}^{0}+\unicode[STIX]{x1D6FD}\frac{\unicode[STIX]{x2202}p_{ex}^{0}}{\unicode[STIX]{x2202}\unicode[STIX]{x1D6FD}}+\unicode[STIX]{x1D702}\frac{\unicode[STIX]{x2202}p_{ex}^{0}}{\unicode[STIX]{x2202}\unicode[STIX]{x1D702}},\end{eqnarray}$$

in which it is employed that

(4.10a,b ) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6FD}}{\unicode[STIX]{x2202}n_{d}}=\frac{\unicode[STIX]{x1D6FD}}{n_{d}},\qquad \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D702}}{\unicode[STIX]{x2202}n_{d}}=\frac{\unicode[STIX]{x1D702}}{n_{d}}.\end{eqnarray}$$

Expressions (4.9) and (4.10) are quite similar to those obtained in Khrapak et al. (Reference Khrapak, Khrapak, Ivlev and Morfill2014) but bear new original meaning. Indeed, one has to be very careful in determining the isothermal compressibility of the dust component. In a series of papers (Mamun et al. Reference Mamun, Shukla and Farid2000; Xie & Yu Reference Xie and Yu2000a ,Reference Xie and Yu b ) the variation of dust particle charge in DAW propagation was treated for the strongly coupled regime but the isothermal compressibility was inaccurately assessed. This inference is simply justified because, for example, the following relation

(4.11) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6E4}_{d}}{\unicode[STIX]{x2202}n_{d}}=\frac{\unicode[STIX]{x1D6E4}_{d}}{3n_{d}}\end{eqnarray}$$

has been implicitly used, which is only valid for $Z_{d}=\text{const.}$ as evidenced by formula (2.6).

Indeed, the isothermal compressibility of the dust component is the derivative of the pressure over the number density, which, in the strongly coupled regime, severely depends on the dust grain charge. Thus, the derivation of the dust grain charge on the dust grain number density should appear in the isothermal compressibility, which was completely overlooked in the literature. Moreover, it is seen from formula (2.6) that when the dust component strongly affects the local plasma quasi-neutrality, expression (4.11) should be rewritten as

(4.12) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6E4}_{d}}{\unicode[STIX]{x2202}n_{d}}=\frac{\unicode[STIX]{x1D6E4}_{d}}{3n_{d}}+\frac{2\unicode[STIX]{x1D6E4}_{d}}{Z_{d}}\frac{\unicode[STIX]{x2202}Z_{d}}{\unicode[STIX]{x2202}n_{d}}.\end{eqnarray}$$

It is thus clear from formula (4.12) that the dust particle charge variation must be also taken into account in the isothermal compressibility of the dust component. A similar procedure must be applied for the screening parameter (2.7), in which the number density of electrons is a function dependent on the dust grain charge and number density.

To strictly evaluate the second term on the right-hand side of expression (4.12) it is necessary to calculate the grain charge as a function of its number density and, as is shown below, the result strongly depends on the dust particle size. This is why the finiteness of dust particle dimension must be appropriately treated right from the charging process up to determination of the corresponding thermodynamic quantities.

In the present consideration, another strategy is proposed stemming from the parametrization of § 1. Namely, the dust grain charge $Z_{d}$ is evaluated as a function of parameters (2.1)–(2.4), and is then used to find the radial distribution function and the reduced excess pressure of the dust component in terms of the same parameters (2.1)–(2.4). Since the dust grain charge is, thus, not an independent parameter, its variation with the dust number density is automatically handled in the evaluation of the dust isothermal compressibility, as is witnessed by expression (4.9). To validate such an interpretation, the simultaneous fulfilment of two conditions is required. First, for the dust grains to be almost instantly recharged and to avoid the direct consideration of the recharging dynamics, the characteristic time of the charging process must be much less in magnitude than the inverse dust-acoustic wave frequency, which is usually the case under ordinary experimental conditions (Fortov et al. Reference Fortov, Khrapak, Khrapak, Molotkov and Petrov2004). Secondly, the dust component must exert a strong influence on the local quasi-neutrality of the plasma medium, which normally means that the Havnes parameter $P=n_{d}Z_{d}/n_{e}$ should be of the order or larger than 1. This is usually the case in the strongly coupled regime of interest herein such that the local variation of the dust particle number density in the dust-acoustic wave leads to a respective change in the local electron number density, thereby affecting the dust grain charge.

Figure 9.

Frequency of the dust-acoustic waves (4.7) as a function of the reduced wavenumber $q=ka_{d}$ at fixed values of $\unicode[STIX]{x1D705}=0.3$ and $\unicode[STIX]{x1D6E4}_{d}=144$ . Red line: $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=1.37\times 10^{-4}$ , $\unicode[STIX]{x1D6E4}_{p}=7.00\times 10^{-5}$ and $\unicode[STIX]{x1D702}=1.00\times 10^{-5}$ ; blue line: $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=3.07\times 10^{-5}$ , $\unicode[STIX]{x1D6E4}_{p}=1.95\times 10^{-5}$ and $\unicode[STIX]{x1D702}=1.00\times 10^{-7}$ ; green line: $\unicode[STIX]{x1D70F}=1$ , $\unicode[STIX]{x1D6FD}=2.30\times 10^{-6}$ , $\unicode[STIX]{x1D6E4}_{p}=3.00\times 10^{-6}$ and $\unicode[STIX]{x1D702}=1.00\times 10^{-9}$ ; black circles: MD data with error bars (Ohta & Hamaguchi Reference Ohta and Hamaguchi2000; Hamaguchi & Ohta Reference Hamaguchi and Ohta2001).

To demonstrate the influence of finite-size effects, figure 9 shows the dust-acoustic wave dispersion (4.7) at three different sets of values of the dimensionless parameters (2.1)–(2.4), which all correspond to the same values of $\unicode[STIX]{x1D705}=0.3$ and $\unicode[STIX]{x1D6E4}_{d}=144$ . In the same figure the results of the numerical experiments (Ohta & Hamaguchi Reference Ohta and Hamaguchi2000; Hamaguchi & Ohta Reference Hamaguchi and Ohta2001) are drawn as circles with the corresponding error bars. It is clearly seen in figure 9 that, when the charge variation is properly treated in the dust isothermal compressibility, even qualitative behaviour of the DAW dispersion starts to strongly depend on the packing fraction, which still remains very small. An accidental coincidence with the MD data is found for $\unicode[STIX]{x1D702}=1.00\times 10^{-7}$ , whereas for $\unicode[STIX]{x1D702}=1.00\times 10^{-5}$ the DAW dispersion remains positive and for $\unicode[STIX]{x1D702}=1.00\times 10^{-9}$ a severe cutoff in the wavenumber is discovered since the DAW frequency cannot turn negative. Such a picture is prescribed to the exact evaluation of the dust isothermal compressibility that accounts for the dust grain charge variation in the course of DAW propagation. In this respect, it would be interesting to work out the molecular dynamics, in which the grain charge depends on its local number density so that the total energy of the grain subsystem would no longer be conserved. In our opinion, such a computer simulation technique would certainly take a model of Yukawa one-component plasma much closer to real dusty plasma experiments.