Algebraic infection of charge correlations at the ionic critical point

We consider a classical two‐component plasma analog of the restricted primitive model of electrolyte, where the hard‐core interaction is replaced by a soft differentiable potential. Within the Born–Green–Yvon hierarchy for the equilibrium distribution functions, we shed light on an infection mechanism where the charge correlations are polluted by the density correlations at the critical point of the liquid–gas transition. This should imply an algebraic decay of critical charge correlations, consistently with the violation of the Stillinger–Lovett sum rule observed in Monte Carlo simulations. In connection with this effect of density correlations on the charge correlations, we also provide a discussion on an extended version of the Debye–Hückel theory that calculates charge correlations in the restricted primitive model for an externally imposed spatial fluctuation of charges.


INTRODUCTION
The liquid-gas transition of a classical electrolyte is a long-standing problem, which has been nicely reviewed by Fisher. [1]his has motivated numerous experimental and theoretical works.On the theoretical side, there still remain difficult questions that have not been satisfactorily answered.A first widely debated problem is about the universality class of the behaviour of thermodynamic properties near the critical point. [2]It was believed, over many years, that the critical exponents are given by mean-field approaches.However, since the end of the last century, careful numerical and experimental inspections at the immediate vicinity of the critical point strongly suggest that these exponents are of the Ising type.There exits now a reasonably general consensus that the transition does belong to the Ising universality class. [3] second controversial question is by concerning the conducting versus dielectric nature of the critical point, in particular for the fully symmetric restricted primitive model (RPM). [1,4]Note that a fully symmetric RPM consists, in a medium having dielectric constant D, of positive and negative spherical ions of equal diameter a, each carrying a charge of equal magnitude q.These charges interact via the long-range Coulomb potential, while a short-range hard-core interaction there prevents overlap between two ions.Away from the critical point, numerical simulations with this model, see, for example, References [4] or [5], have convincingly shown that both the liquid and the gas phases display perfect screening properties that are typical of conducting phases.However, a first suspicion about the breakdown of perfect screening properties, close to the critical point, was pointed out by Caillol. [6]More recently, such breakdown was confirmed through sophisticated Monte Carlo (MC) simulations [4,7] : they show via state-of-the-art finite-size scaling analysis [8] that the second moment Stillinger-Lovett (SL) sum rule [9] for equilibrium charge correlations is violated.This implies that infinitesimal external charges are no longer perfectly screened (see, e.g., Reference [10]).Here note that the internal perfect screening sum rule dictates that the charge-charge structure factor S(k) always vanishes at k = 0, in relation with the charge neutrality in the system.In a conducting phase, for small wave numbers k, it behaves as k 2 times a coefficient, which square root equals the Debye screening length  D , [1] as predicted by the second moment SL sum rule.At the critical point, the MC simulations show that S(k) still behaves as k 2 , but with a multiplicative coefficient different from  2 D .Hence, according to linear response theory, the dielectric constant at k = 0 is finite.
The large-distance decay of charge-charge correlations, S(r), of which S(k) is the Fourier transform, r being the spatial separation between two space points, at the critical point is a central question, directly related to the system's conducting or dielectric nature.If the density-density correlations decay as a power law, [1] namely as 1∕r 1+ with  > 0, various theoretical approaches predict an exponential decay of equilibrium charge correlations for the RPM.In particular, this is observed in a solvable mean-spherical model, [11] which is expected to share common properties with the RPM.In this paper, we consider a fully symmetric two-component plasma (TCP), which is the analog of the RPM where the hard-core interaction is replaced by a smooth differentiable potential.Our main purpose is to show that the critical charge correlations for this TCP are infected by the slow algebraic decay of the critical density correlations.In other words, the charge correlations should also decay in a power law fashion.This is achieved through an inspection of the large-distance behaviour of the equilibrium Born-Green-Yvon (BGY) equations.Our analysis sheds light on the infection mechanism.Note that within the solvable mean-spherical model, [11] an infection mechanism occurs for its charge-asymmetric version but not for the symmetric version analog to the present model and the RPM.
In Section 2, we describe the considered TCP model.Like the RPM, this model is expected to undergo a liquid-vapour transition at low density, general features of which are briefly stated.The BGY hierarchy is introduced in Section 3. We argue that the corresponding BGY equations should remain valid at the critical point.The infection mechanism, which results from the coupling between charge and density correlations, is highlighted in Section 4. We show how this mechanism prevents the exponential decay of S(r).In Section 5, we describe the outline of generalized Debye-Hückel theory, which is a promising formalism for incorporating in a quantitative and tractable fashion the above coupling mechanism.Comments and perspectives are briefly presented in Section 6.

Pairwise interactions
We consider a two-component classical plasma (TCP) made of two species ( = +, −) of mobile particles, carrying charges ±q, in space dimension d = 3.The particles interact via a sum of pairwise interactions which include the familiar Coulomb potential v C (r) = 1∕r and a short-range repulsive interaction v SR (r).The short-range interaction diverges positively faster than 1∕r when r → 0 in order to avoid the collapse between oppositely charged particles.A possible choice for this is with V 0 > 0. In the following, the analysis will be performed for general forms of v SR (r) and is not specific to the choice (II.2).We restate, the present fully symmetric TCP is quite similar to the celebrated RPM, which is fully symmetric with respect to the charges ±q and the hard-core diameters .Here the hard-core potential is replaced by a soft form that is differentiable everywhere, except at r = 0.This allows us to introduce the BGY hierarchy as described in the next Section 3.

Liquid-vapour phase transition
In this paper, we assume that the thermodynamic limit (TL) of the present model exists for any choice of v SR (r).
Note that this has been proven only for short-range regularizations of the Coulomb interaction [12] (see also the review in Ref. [13]).Overall charge neutrality in the system is imposed.Not at too low temperatures, the system in a fluid phase with a common uniform particle density  for both species is invariant under translations once the TL has been taken.As suggested by numerical simulations for similar systems, [14] like the RPM for instance, [2] the present model is expected to undergo a liquid-vapour phase transition.Similar to what occurs for classical fluids made with neutral particles and Lenard-Jones interactions, [15] at a fixed temperature (T) lower than some critical value (T c ), depending upon the overall density, in equilibrium there will be coexistence of vapour and liquid phases with densities  vap (T) and  liq (T).At the critical point (T c ,  c ),  c being the critical value of overall density,  vap and  liq become identical, that is, In the vapour phase, at low densities, particle correlations should decay exponentially fast at large distances.This is strongly suggested by Debye theory, and also by systematic corrections to this mean-field approach derived within the Abe-Meeron diagrammatic expansions [16,17] (see also the rigorous proof by Brydges and Federbush [18] for similar systems).Perfect screening of external charges is then observed, as encoded in the SL second moment sum rule [9] concerning the charge correlations.When the density increases, the system is expected to remain in a conducting state for both vapour and liquid phases, with perfect screening properties.Note that although particle correlations then might decay slower than an exponential, [19] sufficiently fast power-law decays ensure the validity of the SL rule [20] (see also the review [10] ).Moreover, the numerical simulations show the persistence of free charges, which in turn ensure that screening properties still hold.However, as shown via sophisticated Monte Carlo simulations for the RPM, [7] this picture becomes quite doubtful at the critical point (T c ,  c ).As a consequence of the coupling between fluctuations in particle and charge correlations, the latter is expected to be infected by the slow power-law decay of critical particle correlations, and ultimately the second SL sum rule is violated.The main purpose of the present paper is to analyse the infection mechanism within the BGY hierarchy.

THE BGY HIERARCHY
The BGY hierarchy should be a priori valid for the distribution functions of any infinitely extended equilibrium state, provided that the involved spatial integrals do converge in the infinite space.We first write the second BGY equation for the pair distribution functions   1  2 (r 1 , r 2 ) of the infinite system in an homogeneous fluid phase.Then we show that all terms are well-behaved if we assume weak clustering properties for the particle correlations, which are consistent with their expected large-distance decays in the fluid phase, including the critical point.This confirms that the BGY hierarchy can be safely used, even at the critical point.

The second BGY equations
In a fluid phase, the particle distribution functions are invariant under translations, and the second BGY equations for the pair distribution functions , can be written by fixing one particle at the origin, r 1 = 0, and by taking the gradient with respect to the position r 2 = r of the second particle.This provides ∇ ++ (0, r) =  ++ (0, r)F ++ (0, r) and In these equations, is the force exerted on a particle with species  2 and position r, by a particle with species  1 and position r ′ .This force can be decomposed as with the short-range part and the Coulomb part (III.5) For further purposes, it is useful to express the two-and three-body distribution functions in terms of the corresponding particle correlations, whose dimensionless counterparts are the Ursell functions, namely Thanks to both the translational and rotational invariance of the fluid phase, the two-body Ursell functions h  1  2 only depend on the relative distance between the fixed particles.Similarly, the three-body Ursell function h (3) , that is, the three sides of the triangle formed by the particles.Using the decompositions (III.3) and (III.6) of the force and of the distribution functions, we recast the BGY Equations (III.1) and (III.2) as and In these equations, N ( r ′ ) is the correlation between particle densities at points 0 and r ′ , is the correlation between charge densities at points 0 and r ′ , The BGY equations for the density-density and charge-charge correlations are readily obtained by combining Equations (III.7) and (III.8), and they become for r ≠ 0. and (III.12) The various three-body correlations H (3) are defined as linear combinations of the three-body Ursell functions, and they are related to the equilibrium averages of products of three microscopic particle-density or charge-density operators.As shown by the structure of Equations (III.11) and (III.12), the density and charge correlations are coupled together, as it can be a priori expected.

Validity of the BGY equations at the critical point
The various integrals over r ′ involved in the BGY Equations (III.11) and (III.12)do converge under rather weak clustering assumptions on the decay of two-and three-body particle correlations.Indeed, since the short-range force decays as an exponential at large distances |r ′ − r|, the integrals upon r ′ of F SR ( r ′ − r ) times particle correlations are always well behaved.This is not the case of the integrals with the Coulomb force since F C ( r ′ − r ) decays as 1∕|r ′ − r| 2 when r ′ is separated from a fixed r by infinite distance.In order to ensure the (absolute) convergence of the related integrals, the correlations S ( r ′ ) , H (3) dc ( 0, r, r ′ ) and H (3) cc ( 0, r, r ′ ) have to decay faster than 1∕|r ′ | 1+ϵ with ϵ > 0 when |r ′ | → ∞.For any equilibrium state where such weak algebraic decays hold, all terms in the BGY equations are finite: this strongly suggests that these equations are indeed satisfied by the corresponding equilibrium particle correlations.
At the critical point, one expects a slow algebraic decay of all n-body Ursell functions, with n = 2, 3, … , typically as 1∕r 1+ with a strictly positive exponent  > 0. Hence, correlations S ( r ′ ) , H (3) dc ( 0, r, r ′ ) and H (3) cc ( 0, r, r ′ ) decay at least as 1∕|r ′ | 1+ when |r ′ | → ∞ since they are linear combinations of two-and three-body Ursell functions.According to the previous analysis, this implies that the BGY equations remain valid at the critical point.

THE INFECTION MECHANISM
In order to extract constraints from the BGY equations for N(r) and S(r), we first introduce weak assumptions for the respective decays of two-and three-body Ursell.Such assumptions are shown to be consistent with the internal charge sum rules, which are expected to hold in any phase and at the critical point (Section 4.1).Then, we highlight an infection mechanism that prevents the exponential decay of S(r) (Section 4.2).

Clustering assumptions and charge sum rules
At the critical point, according to their respective definitions of particle density (III.9) and charge (III.10)correlations, N(r) and S(r) decay at least as 1∕r 1+ at large distances r.Such a decay should hold for the density correlations N(r), in agreement with the divergence of its integral over r, implied by the compressibility sum rule (note that here the total density is 2) Indeed, the isothermal compressibility  T = [P∕] −1 , where P is the pressure of the system, diverges at the critical point.For the charge correlations S(r), one expects a decay faster than 1∕r 1+ in order to satisfy the internal perfect screening rule which requires the integrability of S(r) over the whole space.Discarding oscillatory behaviours, this implies that S(r) decays at least as 1∕r 3+ϵ with ϵ > 0 when r → ∞.This leads us to infer where the common amplitude A does not depend on the charges carried by the particles.Note that the charge sum rule (IV.2) is crucial for the consistency of the present picture: it guarantees a minimal screening of Coulomb interactions, which in turn do not affect the leading critical tails.The perfect screening rule (IV.2) means that the total charge carried by the polarization cloud surrounding a given fixed particle exactly cancels its charge.If now two particles with charges e  1 and e  2 are fixed at positions 0 and r, the total charge carried by the corresponding polarization cloud should exactly reduce to − ( e  1 + e  2 ) .Hence, the three-body Ursell functions are expected to satisfy the sum rules. [10]dr ′  [ h (3) +++ ( 0, r, r ′ ) − h (3) ++− ( 0, r, r ′ ) ] = −2h ++ (0, r), (IV.4) and In order to pursue our analysis, we introduce a quite plausible ansatz about the decay of correlations between a given pair of particles at fixed positions r 1 and r 2 on the one hand, and a third particle with position r 3 sent to infinity on the other hand.Similarly to the ansatz (IV.3) for the decay of two-particle correlations, we assume where R 12 = (r 1 + r 2 ) ∕2 is the barycentre of the pair of fixed particles, and n 12,3 = (r 3 − R 12 ) ∕|r 3 − R 12 | is the unit vector fixing the direction in which r 3 is sent to infinity.This decay does not depend on the charge e  3 of the third particle, so the difference ] decays as least at 1∕|r 3 − R 12 | 3+ϵ as S(r).Moreover, the amplitudes where n is a fixed unit vector are expected to decay as the two-particle correlations, namely when (r 2 − r 1 ) is sent to infinity in a direction, making an angle  with respect to n.As for particle correlations, the difference is expected to decay as least at 1∕|r 2 − r 1 | 3+ϵ .Note that, since the present TCP is fully symmetric, the amplitude functions satisfy the symmetry relations A ++ = A −− and A +− = A −+ .Now, let us consider the three-body sum rule (IV.4) in the limit r → ∞.In the integral in the l.h.s., the leading contributions arise from the two regions (i) r ′ close to the origin (ii) r ′ close to r.This estimation can be justified by assuming uniform bounds on three-body correlations for large triangular configurations consistent with the decays IV.6.According to the decays (IV.6), both regions (i) and (ii) give identical contributions, which lead to Comparing the behaviour (IV.8) with the large-r decay of the r.h.s of (IV.4) inferred from (IV.3), we find Analogous manipulations can be repeated for the sum rule (IV.5).The leading contributions of regions r ′ close to the origin, and r ′ close to r, then exactly cancel out by virtue of the symmetry relations A ++ = A −− and A +− = A −+ .Hence, no additional constraints on A ++ and A +− are imposed by sum rule (IV.5).

Breakdown of the exponential decay of charge correlations
Within the previous clustering assumptions, we have seen that some combinations of correlations decay faster than 1∕r 1+ because of cancellations.In particular such mechanism arises for the charge correlations S(r), which should decay at least as 1∕r 3+ϵ .Let us assume a priori that S(r) decays exponentially fast at the critical point.Consistently, we then also assume that the combinations of three-body correlations, similar to r), where cancellations of the critical 1∕r 1+ -tails occur, also decay exponentially fast.Using the behaviours (IV.6) for the three-body Ursell functions, the corresponding exponential-decay scenario (EDS) reads.
• EDS3: Then, the strategy consists in showing that this exponential-decay scenario is not consistent with the large-r behaviour of the BGY Equation (III.12).
Let us analyse, within EDS, the large-distance behaviour of the various terms in Equation (III.12).In the l.h.s.∇S(r) decays exponentially fast by virtue of EDS1.In the r.h.s., we first consider the two terms involving the short-range force F SR .The direct short-range term, obviously decays exponentially fast.Because of the exponential decay of F SR , the sole contributions in the three-body short-range term which might decay slower than an exponential arise from the region where r ′ is close to r.However, because of EDS2, H cd ( 0, r, r ′ ) decays exponentially fast for such configurations.Hence, R (3) SR (r) also decays exponentially fast.In the second step, we study the three terms, which involve the Coulomb force F C .The mean-field term decays exponentially fast, by virtue of the charge sum rule (IV.2) and of the rotational invariance of S(r) = S(r).The direct Coulomb term decays algebraically, namely discarding exponentially fast decaying corrections.In the three-body Coulomb term there are exponentially decaying contributions from the region r ′ close to r as a consequence of EDS3.However, there are algebraic contributions from the region r ′ close to the origin 0, which arise from the large-distance behaviour discarding exponentially decaying terms.In (IV.16), the barycentre position is R = r ′ ∕2 and the unit vector is n R = (r − R)∕|r − R|.Hence, we find discarding exponentially decaying terms.The previous analysis shows that all terms in the BGY Equation (III.12)decay exponentially fast, except the sum [ R C (r) + R (3) C ( ] , which, according to Equations (IV.14) and (IV.17) provides the algebraic contribution At large distances r, its asymptotic representation in power series of 1∕r is generated by the expansion of the integrand in powers of r ′ ∕r.The a priori leading term is obtained by setting n R = n = r∕r, while |r − r ′ ∕2| 1+ and F C ( r − r ′ ) can be replaced by 1∕r 1+ and F C (r) respectively.Hence it reads so it decays as 1∕r 3+ .However, the corresponding amplitude vanishes by virtue of the sum rule (IV.9).Thus, the leading behaviour of Equation (IV.18) is provided by the next terms in its multipole-like expansion, and it decays at least as 1∕r 4+𝜂 .
It is quite unlikely that, in the multipole-like expansion of Equation (IV.18), all the amplitudes of the successive terms of order n (n ≥ 1), which decay as 1∕r 3+n+ , identically vanish.Such statement is well supported by the following argument.If one discards the anisotropic effects associated with the occurrence of n R in both A ++ and A +− , the corresponding integral (IV.18) does not depend entirely on the total like-charge inside the sphere with radius r because the effective force F C ( r − r ′ ) ∕|r − r ′ ∕2| 1+ no longer satisfies Gauss theorem.Hence, and contrarily to the case of the total Coulomb force created by a spherically localized charge distribution with zero net charge, the total effective force here decays algebraically, even if the like-charge distribution is exponentially localized.The presence of anisotropic effects enforces this argument, i.e., the existence of non-vanishing amplitudes of order n ≥ 1 in the above multipole-like expansion.
Since all the other terms than ] in Equation (III.12)decay exponentially fast, we conclude that the assumed EDS is not consistent with the BGY hierarchy.Thus, at the critical point, the charge correlations cannot decay exponentially fast, and they are polluted by the algebraic tails in the density-density correlations.The infection mechanism arises from the contributions of three-body correlations, where the critical 1∕r 1+ -tails are coupled to the Coulomb 1/r 2 -force.

DISCUSSION ON CALCULATION OF CHARGE CORRELATIONS IN RPM VIA A GENERALIZED DEBYE-HÜCKEL THEORY
The discussion above provides important information about the influence of density fluctuations on the behaviour of charge correlations that can lead to the violation of the SL second moment sum rule at the ionic critical point.Here we briefly describe, following Reference [21], a generalized version of the Debye-Hückel (DH) theory [22] that calculates charge-charge correlations in the RPM.How the implementation of the above fact in this theory concerning a coupling between density and charge fluctuations can help obtain a quantitative picture of the violation of the above as well as other sum rules?This generalized version of the theory was provided as a step forward from homogeneous to the case of inhomogeneous spatial distributions of charges.
For inhomogeneous particle densities   (r) ( = ±), Debye charging process allows the calculation of Helmholtz free energy, F, as where F 0 is typically the hard-core contribution that applies to the RPM, q  (0 ≤  ≤ 1) is the charge on a particle of species  and   is the average electrostatic potential experienced by an -type ion, sitting at r ′ , that is, The potential   ( r; r ′ ) at r, originates from a charge fixed at r ′ and the induced charge distribution, and it satisfies the Poisson equation Here g  (,  = ±) is the radial distribution function of -type ions about a central ion of type .One splits   ( r; r ′ ) into φ ( r; r ′ ) , the local induced potential, plus the imposed overall electrostatic potential, V C (r) created by the fixed charge distribution ∑  q    (r).Incorporating the DH approximation, which connects g  with φ via a Boltzmann factor for |r − r ′ | > a, one has the Poisson-Boltzmann equation Lee and Fisher linearized the right-hand side of Equation (V.4) to write where κD (r) = ( 4 ∑  q 2    (r)∕D ) 1∕2 is the inhomogeneous Debye wavenumber in the imposed inhomogeneous situation.
Also one has, for |r − r ′ | < a, The last two are the generalized DH (GDH) equations.Note that V C (r), which is in the left-hand side of Equation (V.3), obeys having no contribution in the Boltzmann factor, which determines g  .
Let us consider the imposed particle-density fluctuations for each species  = ±, where Δ is a small given dimensionless parameter.Note that the corresponding total charge density ∑  q    (r) = 2qΔ cos(k ⋅ r) oscillates with the wavenumber k, while the total particle density ∑    (r) = 2 remains uniform.Using (V.8) in the GDH equations, assisted also by the functional relation between free energy density and charge-charge structure factor, Lee and Fisher then obtained an expression for S(k).While this GDH expression works very well in certain regions, it does not pick up neither the violation of the SL second moment sum rule, nor the divergence of the fourth moment. [7]A possible way to improve is to go beyond the linearization approximation used above.On the other hand, it should be noted that the charge fluctuation corresponding to Equation (V.8) keeps the total density space-independent.Recall that fluctuations in the latter are important near the critical point.Thus, incorporation of a different charge fluctuation, such that the fluctuations in the total density field are not washed out, may help obtain a quantitative understanding of the behaviour of the various moments of S(r) at the critical point via this theory.

CONCLUDING COMMENTS
The infection mechanism described in the previous section leads to an algebraic decay of charge correlations.The key ingredients of our derivations are plausible a priori assumptions on the decay of three-body Ursell functions.In a forthcoming paper, [23] we will proceed to a more refined analysis of the BGY equations, inspired by the numerical results obtained within sophisticated Monte Carlo simulations. [7]Such simulations show the violation of the Stillinger-Lovett second moment sum rule on the one hand, and the divergence of the fourth moment of S(r) on the other hand: these facts will be taken as starting inputs in our improved BGY analysis.This provides the exponent s of the 1∕r s -decay of S(r) with 5 +  ≤ s ≤ 6 + , in agreement with heuristic estimations [24] inferred within a completely different approach from ours.
Besides the exact analysis, we will also discuss the predictions derived from approximate theories.A first approximation consists in representing three-body correlations by an ad-hoc truncation (see, e.g., Reference [25]) of systematic diagrammatic series where the bonds are the two-body Ursell functions (see e.g., Reference [26]).Another path, as briefly sketched above in Section 5, relies on an extension of the generalized Debye-Hückel theory of Lee and Fisher, [21] via the introduction of charge fluctuation that complements the density fluctuations: an analogous coupling between the charge-charge and density-density correlations can possibly be realized.