1 Introduction
For the past few decades a great deal of the attention of plasma physics researchers has been paid to complex (dusty) plasmas that are frequently encountered in various settings both in nature and the laboratory (Fortov & Morfill Reference Fortov and Morfill2010). In particular, a dusty plasma appears in various astrophysical contexts (Forsberg et al. Reference Forsberg, Brodin, Marklund, Shukla and Moortgat2006; Malmrose et al. Reference Malmrose, Marscher, Jorstad, Nikutta and Elitzur2011; Seok, Koo & Hirashita Reference Seok, Koo and Hirashita2015), space and Earth experiments (Heidemann et al. Reference Heidemann, Couëdel, Zhdanov, Sütterlin, Schwabe, Thomas, Ivlev, Hagl, Morfill and Fortov2011; Fedoseev et al. Reference Fedoseev, Sukhinin, Abdirakhmanov, Dosbolayev and Ramazanov2016; Izvekova & Popel Reference Izvekova and Popel2016), nanotechnology (Szetsen, Hsiu-Feng & Chien-Ju Reference Szetsen, Hsiu-Feng and Chien-Ju2007; Kundrapu & Keidar Reference Kundrapu and Keidar2012), cancer therapy in medicine (Keidar et al. Reference Keidar, Shashurin, Volotskova, Stepp, Srinivasan, Sandler and Trink2013; Walk et al. Reference Walk, Snyder, Scrinivasan, Kirch, Diaz, Blanco, Shashurin, Keidar and Sandler2013), etc. Furthermore, it is nowadays a working substance in controlled nuclear fusion investigations (Castaldo et al. Reference Castaldo, Ratynskaia, Pericoli, de Angelis, Rypdal, Pieroni, Giovannozzi, Maddaluno, Marmolino and Rufoloni2007; Tolias et al. Reference Tolias, Ratynskaia, De Angeli, De Temmerman, Ripamonti, Riva, Bykov, Shalpegin, Vignitchouk and Brochard2016) and plasma etching in electronics (Kokura et al. Reference Kokura, Yoneda, Nakamura, Mitsuhira, Nakamura and Sugai1999; Kersten et al. Reference Kersten, Deutsch, Stoffels, Stoffels, Kroesen and Hippler2001) because solid micron-sized particles readily penetrate into the plasma medium as a consequence of its contact with electrodes and chamber walls, which substantially modifies surface properties of the confining material and the local plasma characteristics as well.
What makes complex plasma a unique object for intensive investigations is the presence of micron-sized dust particles, called grains, which are capable of acquiring a high, mostly negative, electric charge (Khrapak & Morfill Reference Khrapak and Morfill2009; Shukla & Eliasson Reference Shukla and Eliasson2009), thereby invoking diverse manifestations of strong coupling effects (Bonitz, Henning & Block Reference Bonitz, Henning and Block2010). For instance, experiments with dusty plasmas clearly demonstrate that, under certain external conditions, strong electrostatic interactions between grains prevail over their thermal kinetic energy resulting in the formation of so-called plasma crystals (Kählert & Bonitz Reference Kählert and Bonitz2010; Dietz & Thoma Reference Dietz and Thoma2016; Piel Reference Piel2017). The latter are quite similar, in physical properties, to ordered structures in liquids and solids, such that phase transitions of the first and second orders were theoretically predicted and experimentally observed (Vaulina et al. Reference Vaulina, Vladimirov, Petrov and Fortov2002; Kundu et al. Reference Kundu, Avinash, Sen and Ganesh2014).
In order to reliably describe the properties of strongly coupled systems it is crucial to somehow establish the form of interaction energy between the constituent elements (Momot, Zagorodny & Orel Reference Momot, Zagorodny and Orel2017). As for dusty plasmas, the interaction potential between grains is conventionally taken in the form of the screened Coulomb (Yukawa) potential (Kalman et al. Reference Kalman, Hartmann, Donkó, Golden and Kyrkos2013; Khrapak & Thomas Reference Khrapak and Thomas2015a ,Reference Khrapak and Thomas b ), which is exceptionally valid for point-like dust particles immersed into the buffer plasma, whose role is virtually reduced to shielding of the electric field of grains. Such an idealized one-component model remains intrinsically unphysical in essence because a dust grain must have an appropriate size in order to be able to absorb the buffer plasma particles and thus to be able to acquire the electric charge. To overcome this deficiency a convenient interaction model was proposed in the framework of the density-response formalism (Davletov, Arkhipov & Tkachenko Reference Davletov, Arkhipov and Tkachenko2016) and its logical and proper application is the core motivation of the present consideration.
In most dusty plasma simulations treating the strongly coupled regime, the charge of dust particles is assumed to be invariable, which is only reasonable if the charging time is much larger than the characteristic time of the processes under study. Otherwise, the dust grain charge can no longer be considered a constant, which is true, for example, for dust-acoustic wave (DAW) propagation when local variations of the dust number density in the wave can seriously affect the local quasi-neutrality balance in the medium.
The influence of dust particle charge variation on wave phenomena in cold and ideal dusty plasmas was studied in many papers, see for instance Ivlev et al. (Reference Ivlev, Samsonov, Goree, Morfill and Fortov1999), Ivlev & Morfill (Reference Ivlev and Morfill2000), Ostrikov et al. (Reference Ostrikov, Vladimirov, Yu and Morfill2000), Khrapak & Morfill (Reference Khrapak and Morfill2001), by incorporating the particle charging dynamics. It is advocated herein that such a detailed consideration is irrelevant for DAWs because the recharging time in ordinary dusty plasmas is much lower in magnitude than the inversed frequency of DAWs, which means that, in the course of DAW propagation, dust grains are almost instantly recharged. Of interest in the following is the strongly coupled regime in which the dust charge variation for DAWs must play a much more crucial role because the thermal pressure is predominated by the mutual interaction of likely charged dust particles. DAWs in strongly coupled dusty plasmas with varying grain charge were intensively studied in the literature (Mamun, Shukla & Farid Reference Mamun, Shukla and Farid2000; Xie & Yu Reference Xie and Yu2000a ,Reference Xie and Yu b ; Kaw Reference Kaw2001) but unfortunately the isothermal compressibility of the dust component was calculated under the assumption of the constant grain charge. The present consideration is intended to amend this situation.
With the purpose of correcting the determination of the isothermal compressibility of the dust component, one has to use unconventional parametrization in which the dust particle charge is no longer an independent quantity and is governed by the ambient plasma characteristics as well as by the dust grain size and number density. That is why we start from the examination of the charging process and proceed to the evaluation of the thermodynamic properties of the strongly coupled dust component. The final goal pursued is to show that, in contrast to some previous investigations, the isothermal compressibility found in such a way takes into account the variation of dust particle charge, which has a decisive impact on the DAW spectrum in strongly coupled dusty plasmas.
The rest of the paper is organized as follows. The following § 2 introduces the dimensionless parameters relevant to the description of the state of dusty plasmas. Section 3 is devoted to the determination of the dust grain charge within the orbital motion limited approach. Thermodynamic quantities and the dust-acoustic wave dispersion are extensively dealt with in § 4 starting with the derivation of the correlation functions in the reference hypernetted-chain (RHNC) approximation initially proposed by Lado (Reference Lado1973, Reference Lado1976). Main inferences and general conclusions are summarized in § 5.
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

and the non-isothermicity coefficient

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

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

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

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

and the screening parameter

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):

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

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

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:

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)

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

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

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

where

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

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

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

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

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)

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

in the following form

in which it is employed that

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

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

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.
5 Conclusions
This paper has been generally focused on estimating the impacts that finite-size effects and the screening phenomenon have on the electric charge of dust particles, its thermodynamics and dust-acoustic wave dispersion.
The starting point of the whole consideration is an original intergrain interaction model that stems from the plasma electrodynamics formulated within the linear density-response formalism.
In the framework of the orbital motion limited approximation, the charge number of the dust grain in a plasma has been evaluated as a function of the non-isothermicity, the density-ratio and the packing fraction parameters to show that the screening effects are responsible for the increase of the charge number as compared to the pure Coulomb interaction. It has also been verified that the Yukawa potential always predicts higher values of the dust charge than the intergrain interaction model that takes advantage of the plasma electrodynamics.
The cornerstone of the present consideration has been the self-consistent evaluation of the thermodynamic quantities of the dust component starting from the determination of the grain charge. The radial distribution function has been calculated within the RHNC scheme to take into account the finite dimension of dust particles and it has been proved that the HNC method (Davletov et al. Reference Davletov, Yerimbetova, Mukhametkarimov and Ospanova2014; Tolias, Ratynskaia & De Angelis Reference Tolias, Ratynskaia and De Angelis2014; Yazdi et al. Reference Yazdi, Ivlev, Khrapak, Thomas, Morfill, Löwen, Wysocki and Sperl2014) works very well for small values of the packing fraction. Such a consistent determination of the thermodynamic characteristics has allowed us to study the dispersion of the dust-acoustic waves that accounts for variation of the dust grain charge in the dust isothermal compressibility.
In this regard, the entire analysis of the DAW propagation with variable dust charge (Xie & Yu Reference Xie and Yu2000a ,Reference Xie and Yu b ; Kaw Reference Kaw2001) must be revisited because in the strongly coupled regime almost all transport coefficients begin to depend on the dust grain charge, which is to be implemented elsewhere.
It is now rather timely to discuss a few limitations of the present consideration imposed by the approximations applied. First of all, the neutral gas pressure is assumed to be rather low such that, on the one hand, atoms exert absolutely no influence on dust particles motion in the course of DAW propagation, and, on the other hand, the fluctuating electric field heats up the dust component to a rather high temperature as mentioned in § 2. Another restriction is shared with the orbital motion limited approximation stating that the free flight paths of electrons and ions should be much larger than the dust grain radius, which is again assured by the low neutral gas pressure. In a more accurate model of DAW one has also take into account the dynamics of ions and electrons in the buffer plasma, thereby necessitating strict treatment of the plasma background contribution to the thermodynamic quantities. That the material of the dust grain is a perfect absorber is much harder to handle but is still possible in the framework of the chemical model of dusty plasmas proposed in Davletov et al. (Reference Davletov, Arkhipov and Tkachenko2016).
One of the provisions for future work is the analytic determination of the static correlation functions of § 4 in the framework of the generalized Poisson–Boltzmann equation (Arkhipov et al. Reference Arkhipov, Baimbetov, Davletov and Ramazanov1999), which was shown to be exceptionally valid for weakly and moderately coupled regimes (Arkhipov, Baimbetov & Davletov Reference Arkhipov, Baimbetov and Davletov2003, Reference Arkhipov, Baimbetov and Davletov2005). This will allow us to study the dynamic properties of the dust component (Arkhipov et al. Reference Arkhipov, Askaruly, Ballester, Davletov, Tkachenko and Zwicknagel2010; Dzhumagulova et al. Reference Dzhumagulova, Masheeva, Ramazanov and Donkó2014; Ott et al. Reference Ott, Bonitz, Stanton and Murillo2014) in order to accurately derive the spectrum and decrement of damping of DAWs within the method of moments (Arkhipov et al. Reference Arkhipov, Askaruly, Davletov, Dubovtsev, Donkó, Hartmann, Korolov, Conde and Tkachenko2017) and to make instructive comparison with the results of the quasi-localized charge approximation (Kalman et al. Reference Kalman, Hartmann, Donkó and Rosenberg2004; Donkó, Kalman & Hartmann Reference Donkó, Kalman and Hartmann2008).
Acknowledgements
The authors are indebted to the Ministry of Education and Science of the Republic of Kazakhstan for financial support through the state grants AP05132677, AP05132333 and the program BR05236730.