Continuum theory of electrostatic-elastic coupling interactions in colloidal crystals

Hao Wu; Zhong-Can Ou-Yang; Rudolf Podgornik

doi:10.1088/1572-9494/ad99af

PDF(792 KB)

Communications in Theoretical Physics ›› 2025, Vol. 77 ›› Issue (5) : 55602. DOI: 10.1088/1572-9494/ad99af

Statistical Physics, Soft Matter and Biophysics

Continuum theory of electrostatic-elastic coupling interactions in colloidal crystals

Author information +

History +

Key words

colloids / electrostatics / Poisson and Laplace equations / boundary-value problems / statistical mechanics of classical fluids

Cite this article

EndNote

Ris (Procite)

Bibtex

Download Citations

Hao Wu, Zhong-Can Ou-Yang, Rudolf Podgornik. Continuum theory of electrostatic-elastic coupling interactions in colloidal crystals[J]. Communications in Theoretical Physics, 2025, 77(5): 55602 https://doi.org/10.1088/1572-9494/ad99af

1. Introduction

Colloidal crystals present an interesting model of condensed matter systems where atoms and electrons are substituted by charged colloidal macroions and small, mobile electrolyte ions [1]. They can form in systems with purely repulsive long-range interactions [2–5] or in systems with combined short-range attraction and long-range repulsion [6]. The Deryaguin–Landau–Verwey–Overbeek (DLVO)-type interactions between these constituents can be modified and engineered [7], yielding valuable insights not only into the properties of colloidal crystals themselves but also illuminating the theoretical foundations of atomic and molecular materials such as crystalline elasticity, crystal melting and lattice defects and dynamics, to name just a few [1, 8, 9]. Colloids of various components and different size/charge asymmetries have been assembled into diverse crystalline structures and have served as experimental models to study phase behaviors and self-assembly processes [1, 10]. Of particular interest are the binary colloidal crystals where the smaller, mobile component occupies interstitial space in the lattice, formed by the larger colloid component, and can exhibit a localized-to-delocalized transition [11, 12], as is the case in size-asymmetric DNA-functionalized nanoparticles [1]. A similar type of behavior is also known to occur in a two-dimensional (2D) crystalline array of cylindrical macroions with interstitial linear polyelectrolyte chains exhibiting a behavior analogous to that of electrons in a 2D crystal [13]. In the context of colloidal crystals, this localized-to-delocalized transition has been dubbed the ionic to metallic transition [12].

Attractive interactions between like-charged colloidal particles that cannot be rationalized within the DLVO paradigm [14], first seen in detailed electrostatic double-layer simulations [15], have been directly observed in various experiments [16, 17] and have been studied extensively in the framework of different theoretical approaches [18] without, however, any clearcut consensus emerging as to their underlying universal mechanism(s). This embarras de richesses [19] in theoretical understanding of the non-DLVO colloid interactions stems partially from different microscopic models and partially from different approximations applied to solve these models. Among the variations on the model theme we can refer to the most recent proposal of charge regulation at the interacting colloid surfaces yielding an accurate description of the observed force profiles, including their attractive part, down to a few nm [16, 20, 21]. On the other hand, the departure from the mean-field approximation is still best rationalized in terms of the weak–strong coupling dichotomy of confined Coulomb fluids that yields an attractive electrostatic interaction between identical colloids [15, 18].

As a consequence, discovering new phenomena that display some novel features of the remaining unexamined life of the anomalous attractive forces in charged colloids are thus certainly worth focusing on in more detail [17]. Elastic deformation effects are ubiquitous in crystalline solids [22] where various defects, such as dislocations, interstitials and vacancies, create elastic displacement fields that allow them to interact with long-range interactions [8, 9, 23–25]. These lattice Hookian elasticity-mediated forces between point-defects have been shown to be attractive, leading in the case of interstitials and vacancies to the formation of defect strings despite the purely repulsive forces between the colloidal particles [8, 9, 26, 27]. While the phenomenology of these lattice-mediated attractions bears some similarity to the anomalous attractions between identically charged colloids (see above), the underlying mechanism seems to be quite different and involves an additional (elastic) degree of freedom. Similar lattice-mediated attractions can be surmised to exist between diffusing test charges in the interstitial mobile Coulomb fluid but modified by the ubiquitous electrostatic interactions. This situation is superficially reminiscent of the phonon-mediated attractions between electrons in superconductor Cooper pairing [28] but, of course, lacks its quantum basis, and therefore could be called ‘thermal Cooper pairing' [29]. Motivated by these findings, we will attempt to formulate the simplest possible continuum model of the electrostatic-elastic coupling in a model colloidal crystal and ascertain to what extent it mimicks the behavior of other solid state systems.

In what follows we will present a mean-field continuum model of a background 3D crystalline lattice with an interstitial Coulomb fluid and calculate explicitly the effective interaction between two test interstitial charges. The model is based on continuum elasticity of the crystalline lattice and Poisson–Boltzmann (PB) mean-field electrostatics for the interstitial small, mobile ion component. In addition, we will assume a minimal electrostatic-elastic coupling due to the colloidal lattice–mobile ion interactions, which implies a local lattice dilation/constriction at the location of an interstitial mobile charge. In thermodynamic equilibrium this coupling will lead to two equilibrium equations: a modified Navier equation for the elastic displacement and a modified PB equation for the electrostatic field. Linearizing and solving these two equations explicitly for a test point interstitial charge, we derive an effective interaction between two interstitial test sources that displays a variable attractive component depending on the strength of the bare electrostatic screening as well as the electrostatic-elastic coupling. Notably and contrary to the anomalous attraction in other colloidal systems (see above), the effective attraction emerges already on the mean-field level, the reason being that the electrostatic-elastic coupling leads to a two-field theory and thus straightforwardly implies a decrease of the total interaction free energy.

2. Minimal coupling model: mean-field theory

We consider a minimal model describing the electrostatic-elastic coupling as follows: we have an ideal colloid crystalline lattice of charged macroions, see figure 1, embedded in a bathing Coulomb fluid of mobile univalent positive and negative ions, in thermal equilibrium with a bulk phase. In the spirit of a continuum theory, both of the components, the mobile ions as well as the lattice colloids, are assumed to be point-like. In addition, we assume that the lattice colloids can exhibit only small deviations from their equilibrium lattice positions, while the bathing Coulomb fluid ions are fully mobile.

Figure 1. A colloidal lattice contains an interstitial charged ‘intrusion'. A negative (red circle) or positive (blue circle) interstitial test charge imposes a local elastic dilation (orange arrow) or contraction (brown arrow) on the ideal negatively charged colloidal crystal lattice (magenta circle). The interstitial charge therefore induces a local change in volume/density of the lattice.

Full size|PPT slide

Within this model we assume the free energy of the system as being decomposable into an elastic part, describing the lattice colloids, and an electrostatic part, describing the mobile Coulomb fluid, as well as a minimal coupling part that connects the electrostatic field with elastic deformation. Our aim is then to calculate an effective electrostatic-elastic interaction between two point-like interstitial test particles (‘intrusions') immersed in the system that interact with the Coulomb fluid components electrostatically as well as with the background colloidal lattice elastically.

The electrostatic free energy of the Coulomb fluid is furthermore assumed to coincide with the standard PB free energy [14], while the elastic energy of the colloidal lattice is of a standard Hookeian form, allowing a further simplification if one assumes that the elastic medium is isotropic [22]. The elastic constants can be expressed as lattice sums of the interaction potential between colloid macroions in an undistorted lattice [9], which is furthermore assumed to be three-dimensional [30]. Finally, the electrostatic-elastic coupling term embodies the effect that the interstitial charges have on the dilation (or contraction) of the background colloidal lattice. It features in two different contexts: the first is the coupling between the Coulomb fluid charge density and the elastic displacement, and the second one is the coupling between the test charge particles and the neighboring background colloidal lattice. The two characteristics have both been included in our recent letter [29].

The electrostatic free energy is based on the mean-field PB description and is composed of an electrostatic energy and an ideal gas entropy of the mobile ions. It can be written in different equivalent forms, either in the form of a density functional theory [31] or a field theory [32]. We choose here a hybrid description of the form [14]

\begin{array}{rcl} \begin{array}{l} F_{E S} [n_{\pm}, ϕ] = \int_{V} d^{3} x [- \frac{1}{2} ε {(\nabla ϕ)}^{2} + e_{0} (n_{+} - n_{-}) ϕ \\ + k_{B} T \sum_{\pm} (n_{i} l n n_{i} a^{3} - n_{i}) - \sum_{\pm} μ_{i} n_{i} + e_{0} ρ_{e} ϕ], \end{array} \end{array}

(1)

where ϵ is proportional to the dielectric permittivity, e₀ is the elementary charge, n_± are the densities of the diffusing mobile ions, μ_± = μ₀ are their chemical potentials, which are identical since the mobile components are electroneutral in the bulk, φ is the mean electrostatic potential, while ρ_e here and below corresponds to the source term pertaining to the point-like test charge source at the origin, with the density ρ_e = δ(x). Just as in the case of the (linearized) PB theory, the test charge will allow us to calculate the electrostatic potential and the charge–charge interaction. The minimization of this free energy yields the PB equation and the chemical potentials are expressed by the concentrations of the ions in the bulk.

For the elastic part, we first of all rewrite the standard Hookeian elastic free energy for an isotropic elastic medium that is usually written in the form [22]

\begin{array}{rcl} F_{E L} [u_{i k}] = \int_{V} d^{3} x [\frac{1}{2} λ {(T r u_{i k})}^{2} + μ T r u_{i k}^{2} - f_{i} u_{i}], \end{array}

(2)

where u_i is the elastic displacement vector, f_i is the external force density, u_ik is the symmetrized deformation tensor, while λ and μ are Lamé coefficients, which can be expressed as lattice sums of the undistorted lattice [9]. The last term is the source term, pertaining to the point-like interstitial elastic test source at the origin. The source term in equation (1) and the source term in equation (2) both correspond to a test particle at the origin, but the former describes its electrostatic interaction and the latter describes its elastic interaction with the system.

The elastic constants can be obtained straightforwardly in the absence of thermal effects, starting from the total energy of the colloidal crystal [9]

\begin{array}{rcl} E = \frac{1}{2} \sum_{i \neq j} υ (| r_{i} - r_{j} |), \end{array}

(3)

where

υ

(∣r_i − r_j∣) is the pair potential between colloid particles i and j in the undistorted lattice. Its form is often assumed to coincide with the Debye–Hückel (DH) interaction with an effective charge of the colloid [2–5] that includes the effects of distorted double layers at high macroion volume fractions, but it should also include the details of the dissociation process, which are more difficult to incorporate. From here, in the case of an isotropic elastic material, one can derive the following identities for the Lamé elastic constants λ and μ [33, 34]

\begin{array}{rcl} \begin{array}{rcl} λ + \frac{2}{3} μ & = & \frac{ρ}{18} \sum_{i}^{'} (υ^{''} (r_{i}) r_{i}^{2} - 2 υ^{'} (r_{i}) r_{i}) \\ μ & = & \frac{ρ}{30} \sum_{i}^{'} (υ^{''} (r_{i}) r_{i}^{2} + 4 υ^{'} (r_{i}) r_{i}), \end{array} \end{array}

(4)

where

υ^{'} (r)

and

υ

″(r) are the first and second derivative of the pair potential, respectively, ρ is the colloid number density and r_i is the distance of particle i from the origin. In the summation the floating prime indicates that we pick one colloid at the center and sum over all the neighbors. The colloid–colloid interaction potential is—at least in part—connected with electrostatic interactions [35] and thus follows from the electrostatic free energy with colloidal charges as external sources [36].

We now cast the Hookeian elastic energy into a different form by inserting the symmetrized deformation tensor

u_{i k} = \frac{1}{2} (\partial_{i} u_{k} + \partial_{k} u_{i})

into equation (2), while using the Gauss theorem for the terms that can be cast into the form of a divergence. We thus obtain the elastic free energy as a sum of a volume integral and a surface integral, where the former contains only the scalar invariants composed of ∇ · u and ∇ × u, plus some irrelevant surface terms⁶

6 $μ \oint_{S} n \cdot [\frac{1}{2} \nabla u^{2} - u \times (\nabla \times u) - u (\nabla \cdot u)] d S$

. This is of course possible only in the case of an isotropic elastic medium. Since the colloid crystal is an anisotropic elastic medium, this is an approximation that ignores the directional aspects of the elastically mediated interaction but, on the other hand, allows for a simplified derivation of the salient features of the positional dependence.

In addition, the source term of equation (2), corresponding to a point-like interstitial test particle in the elastic medium, can be modified by taking into account that the point-like body force can be written as f = v₀(λ + 2μ)∇ρ_e [9, 37], where v₀ is the local volume change induced by the point-like test particle, acting as a local dilation/contraction center, see figure 1. In fact, v₀ of a point-like interstitial test particle plays the same role in elastic interactions as e₀ of a point-like charge plays in the electrostatic interactions, and the relaxation of the deformation away from the test particle is described by the Navier equation and the electrostatic potential by the PB equation.

Furthermore, it is clear that for an infinite sample volume one can derive the equality

\begin{array}{rcl} \int_{V} d^{3} x \nabla ρ_{e} \cdot u ⟶ - \int_{V} d^{3} x ρ_{e} \nabla \cdot u, \end{array}

(5)

which allows us to write the elastic source term in an analogous form as in the electrostatic free energy equation (1). Effectively the same form of the elastic coupling, based on a 2D point interstitial defect modelled by two orthogonal pairs of forces, has been used in [9]. The couplings of the elastic displacement field and the electrostatic potential to the density of the interstitial test particles ρ_e are thus completely analogous and given by

\begin{array}{rcl} v_{0} (λ + 2 μ) \int_{V} d^{3} x ρ_{e} \nabla \cdot u a n d e_{0} \int_{V} d^{3} x ρ_{e} ϕ, \end{array}

(6)

for the elastic and the electrostatic case, respectively. Clearly, v₀(λ + 2μ) plays the role of an effective elastic ‘charge'.

With these provisos the elastic energy of an isotropic elastic medium then assumes a much simplified final form

\begin{array}{rcl} \begin{array}{rcl} F_{E L} [u] & = & \int_{V} d^{3} x [\frac{1}{2} (λ + 2 μ) {(\nabla \cdot u)}^{2} + \frac{1}{2} μ {(\nabla \times u)}^{2} \\ + v (λ + 2 μ) ρ_{e} \nabla \cdot u] . \end{array} \end{array}

(7)

Note the analogy but also the fundamental difference between the electrostatic and the elastic free energies: while the former one is formulated in terms of the scalar field and, after an expansion to the second order, exhibits (Debye) screening, the latter is formulated in terms of a vector field and does not exhibit any screening, remaining ‘critical' (i.e. having an inifinite screening length) for all values of the parameters. This is a general property of the elastic Hamiltonians in soft condensed matter that further implies Casimir-like phenomena in e.g. nematic, smectic and cholesteric liquid crystals [38, 39].

To the electrostatic and elastic free energies we now add an explicit minimal coupling term, corresponding to an energy cost associated with dilatation/contraction of the colloidal lattice due to a local excess of positive or negative mobile charge. Since in an ideal crystal there is a one-to-one mapping of particles to lattice positions [40], the relation between the density change, δc, of the colloidal lattice with mean density c, and the divergence of the displacement field can be written as [41]

\begin{array}{rcl} δ c = - c \nabla \cdot u . \end{array}

(8)

In our model, we assume that the local density change of the colloidal lattice is a consequence of the local mobile charge density mismatch proportional to (n₊ − n₋), and therefore the electrostatic-elastic minimal coupling term can be assumed to have the form

\begin{array}{rcl} F_{C} [n_{\pm}, u] = v_{0} (λ + 2 μ) \int_{V} d^{3} x (n_{+} - n_{-}) \nabla \cdot u, \end{array}

(9)

consistent with equation (6) and based upon the assumption that the density change for positive and negative ions is the same but of opposite sign. We additionally assumed that the elastic coupling is the same—except for the sign—for both types of mobile charges, so there is only a single coupling constant. The above coupling term is just the first term in a series of higher-order symmetry-permitted scalar invariants that we will not consider explicitly in this work.

Gathering all the contributions to the free energy in this minimal model of electrostatic-elastic coupling, the total free energy is then given by

\begin{array}{rcl} \begin{array}{rcl} F & = & F_{E S} [n_{\pm}, ϕ] + F_{E L} [u] + F_{C} [n_{\pm}, u] \\ = & \int_{V} d^{3} x g (n_{\pm}, ϕ, u), \end{array} \end{array}

(10)

and the thermodynamic equilibrium is described by the corresponding Euler–Lagrange (EL) variational equations that are of the form

\begin{array}{rcl} \frac{\partial g}{\partial n_{\pm}} = 0 a n d \frac{\partial g}{\partial ϕ} - \nabla \cdot (\frac{\partial g}{\partial \nabla ϕ}) = 0, \end{array}

(11)

as well as

\begin{array}{rcl} \frac{\partial g}{\partial u_{i}} - \nabla_{k} (\frac{\partial g}{\partial \nabla_{k} u_{i}}) = 0. \end{array}

(12)

The explicit form of these equations is as follows. For the mobile charges

\begin{array}{rcl} \begin{array}{l} n_{+} : e_{0} ϕ + v_{0} (λ + 2 μ) \nabla \cdot u + k_{B} T l n n_{+} a^{3} - μ_{0} = 0, \\ n_{-} : - e_{0} ϕ - v_{0} (λ + 2 μ) \nabla \cdot u + k_{B} T l n n_{-} a^{3} - μ_{0} = 0. \end{array} \end{array}

(13)

For the electrostatic potential we derive a modified Poisson equation, and for the elastic displacement we derive a modified Navier equation

\begin{array}{rcl} \begin{array}{l} ϕ : e_{0} ρ_{e} + e_{0} (n_{+} - n_{-}) - \nabla \cdot (- ε \nabla ϕ) = 0, \\ u : v_{0} (λ + 2 μ) (\nabla ρ_{e} + \nabla (n_{+} - n_{-})) \\ + μ \nabla^{2} u + (λ + μ) \nabla (\nabla \cdot u) = 0. \end{array} \end{array}

(14)

These equations can be cast into a much simpler form by introducing the charge density of the mobile ions

\begin{array}{rcl} \begin{array}{l} e_{0} (n_{+} - n_{-}) \equiv ρ (ϕ, u) \\ = - \frac{2 e_{0} e^{β μ_{0}}}{a^{3}} \sinh β [e_{0} ϕ + v_{0} (λ + 2 μ) \nabla \cdot u] . \end{array} \end{array}

(15)

Note the difference: ρ is the mean-field charge density, while ρ_e is the external test particle source density.

The first equation in equation (14) presents a generalization of the PB equation, while the second one generalizes the Navier equation. Without the source term, the mean-field solution corresponds to the bulk state of vanishing electrostatic potential and vanishing elastic displacement.

We now invoke the standard Helmholtz decomposition [9] for the elastic displacement field in terms of the scalar (Φ) and vector (A) elastic potentials

\begin{array}{rcl} u = \nabla Φ + \nabla \times A, \end{array}

(16)

representing the curl-free longitudinal and the divergence-free transversal part of the elastic deformation. Since all the electrostatic-elastic couplings in the free energy equation (10) are in terms of ∇ · u, we can assume the trivial solution

A = c o n s t .

[9] for an infinite sample. In addition, in the spirit of the DH approximation [14], we linearize the charge density to obtain

\begin{array}{rcl} ρ (ϕ, u) ≃ - 2 e_{0} λ_{0} (e_{0} ϕ + v_{0} (λ + 2 μ) \nabla \cdot u), \end{array}

(17)

with

λ_{0} = \frac{β e^{β μ_{0}}}{a^{3}} = β n_{0}

, where n₀ is the density in the bulk reservoir. Just as in the case of the standard DH theory of Coulomb fluids, the linearization ansatz is consistent if

\begin{array}{rcl} β e_{0} ϕ ≪ 1 a n d β v_{0} (λ + 2 μ) \nabla \cdot u ≪ 1. \end{array}

(18)

Since for point-like test particles both the electrostatic and the elastic potentials formally diverge at the source, see next section for details, one would need to introduce a finite-size core inside which the theory breaks down, which actually means one should use effective values for both e₀ and v₀.

This linearization furthermore implies that the electrostatic free energy equation (1) is expanded to the second order in the charge density and/or potential, effectively leading to Debye screening of electrostatics. Since, as we already noted, the elastic energy equation (7) is ‘critical' and implies no separate elastic screening, one can expect that the coupled description will also yield only one type of electrostatic screening.

Taking into account the linearization and the Helmholtz decomposition of the elastic displacement, we end up with the following coupled set of linear equations

\begin{array}{rcl} \begin{array}{rcl} e_{0} ρ_{e} & = & 2 λ_{0} {e_{0}}^{2} ϕ - ε \nabla^{2} ϕ + 2 λ_{0} e_{0} v_{0} (λ + 2 μ) \nabla^{2} Φ, \\ v_{0} ρ_{e} & = & 2 λ_{0} v_{0} e_{0} ϕ + [2 λ_{0} v_{0}^{2} (λ + 2 μ) - 1] \nabla^{2} Φ, \end{array} \end{array}

(19)

where the electrostatic-elastic coupling terms are characterized by the product e₀v₀. The solution of these two equations gives the electrostatic potential and the elastic potential in the system up to the lowest order electrostatic-elastic coupling terms.

3. Electrostatic and elastic deformation fields

The EL equations, equation (19), can now be solved explicitly for a point-like test particle with an electrostatic charge e₀, and an elastic ‘charge' v₀(λ + 2μ), described by a delta function-like density ρ_e = δ(r) located at the origin. The solution yields the spatial relaxation of the mobile ion charge density, ∇ · E = − ϵ∇²φ, and the colloid lattice relative density change , δc/c = ∇ · u, away from the origin.

The solution of this set of linear equations can be obtained straightforwardly and we omit the unnecessary details. The electrostatic potential, proportional to the electrostatic charge of the test particle e₀, has the form

\begin{array}{rcl} ϕ (r) = \frac{e_{0} α}{4 π ε} \frac{e^{- κ r}}{r}, \end{array}

(20)

with

\begin{array}{rcl} α = \frac{1}{1 - ξ} a n d κ^{2} = α κ_{0}^{2}, \end{array}

(21)

with

ξ = 2 β n_{0} v_{0}^{2} (λ + 2 μ)

the electrostatic-elastic coupling constant and

κ_{0}^{2} = 2 λ_{0} e_{0}^{2} / ε

the inverse square of the standard Debye screening length. Since ξ ≥ 0 it follows that the screening length with electrostatic-elastic coupling decreases compared to the standard Debye screening length. The electrostatic field is then obtained from

\begin{array}{rcl} E = - \nabla ϕ = \frac{e_{0} α}{4 π ε} (1 + κ r) \frac{e^{- κ r}}{r^{2}} \frac{r}{r}, \end{array}

(22)

and describes the spatial relaxation of the mobile ion charge density proportional to ∇ · E. The elastic deformation potential, proportional to the elastic ‘charge' v₀(λ + 2μ), has the form

\begin{array}{rcl} Φ (r) = \frac{v_{0} α}{4 π} \frac{e^{- κ r}}{r}, \end{array}

(23)

and taking into account u = ∇Φ, we have the elastic deformation vector in the medium as

\begin{array}{rcl} u \cdot \frac{r}{r} = - \frac{v_{0} α (1 + κ r)}{4 π} \frac{e^{- κ r}}{r^{2}}, \end{array}

(24)

describing the spatial relaxation of the colloid lattice relative density change, ∇ · u, as

\begin{array}{rcl} \nabla \cdot u = \nabla^{2} Φ = v_{0} α \frac{κ^{2} e^{- κ r}}{4 π r} . \end{array}

(25)

Recalling the discussion on the limits of linearization in the previous section, we can now write the solutions for the point-like test particles as

\begin{array}{rcl} \begin{array}{l} β e_{0} ϕ (r) = 2 (κ ℓ_{B}) α \frac{e^{- κ r}}{κ r}, \\ β v_{0} (λ + 2 μ) \nabla \cdot u (r) = 2 {(κ ℓ_{E})}^{3} α \frac{e^{- κ r}}{κ r}, \end{array} \end{array}

(26)

where we introduced two characteristic length scales, the electrostatic Bjerrum length

\begin{array}{rcl} ℓ_{B} = β e_{0}^{2} / (8 π ε), \end{array}

(27)

and an elastic characteristic length given by

\begin{array}{rcl} ℓ_{E}^{3} = β v_{0}^{2} (λ + 2 μ) / (8 π) . \end{array}

(28)

Note the electrostatic-elastic analogy in the definition of both lengthscales. We get an estimate for both lengths by assuming an aqueous solvent ϵ = ϵ_rϵ₀ ∼ 78ϵ₀ with ℓ_B ∼ .7 nm, and a colloidal crystal bulk modulus (λ + 2μ) ∼ 100 dyne/cm² implying that

ℓ_{E} \sim {(10^{- 5} v_{0}^{2})}^{1 / 3}

if v₀ is in nm. Consequently ℓ_E and ℓ_B would be comparable for v₀ being a fraction of a colloid lattice unit cell of ∼μm lattice spacing. In this case, the validity of the linearization ansatz can then be expressed as (κℓ_B)α(e^−κa)/(κa) < 1, which implies also

{(κ ℓ_{E})}^{3} α (e^{- κ a}) / (κ a) < 1

, where a would correspond to an effective core radius of the interstitial test particle.

We next calculate the total interaction energy between two otherwise identical point-like test sources with the same electrostatic charge and elastic ‘charge', one located at x₁ and the other one at x₂, with ∣x₁ − x₂∣ = D. From the total free energy and taking into account the linearized mean-field equations, equation (19), we obtain

\begin{array}{rcl} F = \frac{1}{2} \int_{V} d^{3} x (e_{0} ρ_{e} ϕ + v_{0} (λ + 2 μ) ρ_{e} \nabla \cdot u), \end{array}

(29)

which can be evaluated explicitly. Note that ρ_e, φ and u refer to the sum of all the source(s) and all the field(s) and thus necessarily contain also the self energies; however, we will be interested only in the part that depends on the separation D between the sources. Using the second equation in equation (19) we derive equation (29) in the final form

\begin{array}{rcl} F (D) = \frac{{e_{0}}^{2} α^{2}}{8 π ε} \frac{e^{- κ D}}{D} = \frac{{e_{0}}^{2} α^{2}}{8 π ε} \frac{e^{- \sqrt{α} κ_{0} D}}{D} . \end{array}

(30)

For non-zero ξ, it follows from its definition that α ≥ 1, while the decay length is obviously smaller (shorter-ranged) than the bare Debye length, κ ≥ κ₀. This follows from the fact (see section 4) that the inverse square of the screening length is proportional to the charge density response function. Consequently, any terms additive in free energy, like the elastic terms in our case, add to the charge density response function, thus decreasing the effective screening length.

The interaction energy

F (D)

can be either larger or smaller then the standard (no background colloidal lattice) DH screened electrostatic interaction

\begin{array}{rcl} F_{D H} (D) = \frac{{e_{0}}^{2}}{8 π ε} \frac{e^{- κ_{0} D}}{D}, \end{array}

(31)

depending on the interplay of the increased magnitude and increased screening. The electrostatic-elastic coupling thus creates a non-trivial modification of the interactions between interstitial inclusions.

4. Fluctuations around mean-field solution

Notably, the solutions of the linearized mean-field equations, as derived above, can in principle be extended to a negative value of the coupling constant α, or equivalently to imaginary screening length leading to oscillatory field solutions. Would this be real? An important issue that needs to be addressed at this point is therefore the range of validity of the mean-field solution that we derived above. The standard approach is to investigate the fluctuations around this solution, governed by the field Hessian, and to ascertain that all the eigenvalues of the Hessian are positive.

In order to assess the role of fluctuations around the mean-field (MF) solution, δu and δφ, we need to analyze the deviations from the mean-field solution quantified by u = u_MF + δu and φ = φ_MF + δφ, where u_MF, φ_MF are the solutions of the EL equations derived above. In order to proceed, we first note that the mean-field free energy equation (10) can be cast in the form

\begin{array}{rcl} \begin{array}{rcl} F [ϕ, u] & = & \int_{V} d^{3} x (- \frac{1}{2} ε {(\nabla ϕ)}^{2} + \frac{1}{2} (λ + 2 μ) {(\nabla \cdot u)}^{2} \\ - 2 k_{B} {T n}_{0} \cosh β (e_{0} ϕ + v_{0} (λ + 2 μ) \nabla \cdot u)), \end{array} \end{array}

(32)

which now contains only the field variables. Without the elastic displacement terms this free energy would correspond to the well-known alternative form of the PB free energy that contains only the mean-field electrostatic potential [42].

The exact partition function now has to be written as a functional integral over the imaginary field iφ [43, 44] and the field u. Since we are interested only in the limit of stability of the mean-field solution, we expand this field partition function to the second order in δφ and δu around the saddle point of the field action, corresponding to the evaluation of the field Hessian, which will fortunately turn out to be rather simple to analyze, and will consequently answer our original question on the range of validity of the mean-field theory.

Expanding the total field free energy to the second order we are then led to a deviation from the mean-field value of the free energy for the electrostatic potential and the longitudinal component of the deformation field

\begin{array}{rcl} \begin{array}{l} F [ϕ_{M F} + δ ϕ, u_{M F} + δ u] - F [ϕ_{M F}, u_{M F}] \\ = \int_{V} d^{3} x (- \frac{1}{2} ε {(\nabla δ ϕ)}^{2} {+ \frac{\partial ρ}{\partial ϕ} |}_{M F} δ ϕ^{2} \\ + \frac{1}{2} {(λ + 2 μ) (\nabla \cdot δ u)}^{2} + \frac{1}{2} v_{0} (λ + 2 μ) {\frac{\partial ρ}{\partial (\nabla \cdot δ u)} |}_{M F} \\ \times {(\nabla \cdot δ u)}^{2} \\ + i v_{0} (λ + 2 μ) {\frac{\partial ρ}{\partial ϕ} |}_{M F} δ ϕ (\nabla \cdot δ u)), \end{array} \end{array}

(33)

since the transverse component of δu is not coupled to electrostatics. Above we introduced the derivatives of the mean-field charge density, equation (15). As already noted, the i in front of the term linear in δφ stems from the functional integral representation of the Coulomb fluid partition function [43, 44]. In addition, because of equation (15), we have

\begin{array}{rcl} e_{0} \frac{\partial ρ}{\partial (\nabla \cdot δ u)} = v_{0} (λ + 2 μ) \frac{\partial ρ}{\partial ϕ} . \end{array}

(34)

The derivative of the charge density,

\frac{\partial ρ (ϕ)}{\partial ϕ}

, by definition the capacitance density response function, is in fact proportional to the square of the inverse Debye screening length [45]

\begin{array}{rcl} κ_{0}^{2} = {- \frac{1}{ε} \frac{\partial ρ (ϕ)}{\partial ϕ} |}_{M F} = \frac{2 λ_{0} e_{0}^{2}}{ε} . \end{array}

(35)

The longitudinal component of the elastic displacement vector derivatives yields a Gaussian integral and can be integrated out from the partition function, yielding a purely electrostatic expression with a renormalized screening length

\begin{array}{rcl} κ^{2} ⟶ κ_{0}^{2} α = κ_{0}^{2} \frac{1}{1 - ξ} . \end{array}

(36)

This of course tallies with the screening length derived from the mean-field equations.

Since the condition that the eigenvalues of the Hessian be positive is thus reduced to the condition that the screening length is real, the stability of the mean-field solution is defined by 1 − ξ ≥ 0, with strict identity corresponding to the limit of stability.

Let us investigate then the consequences of this result. From the definition of the elastic constants, equation (4), it follows that

ξ \sim (n_{0} a^{3}) {(v_{0} / a^{3})}^{2} (β υ^{''} (a) a^{2})

, where a is the colloidal lattice spacing and the colloid–colloid interaction potential was assumed to be short-ranged. The strength of this interaction potential depends on the colloidal charge, which cannot be measured directly [35]. Taking instead the measured colloidal crystal bulk modulus (λ + 2μ) ∼ 100 dyne/cm², the regime of ξ ≤ 1 holds if β

υ

″(a) ≪ 1/a² and n₀ ≪ 1/a³. If, on the other hand, these criteria are not satisfied and the colloidal crystal is strongly charged, the mean-field would cease to be a valid description and, following the general analogy with Coulomb fluids [46], one would need to formulate a strong coupling theory for this problem.

Based on the above analysis, the applicability of the mean field theory is therefore limited by the requirement that ξ ≤ 1, which can be reduced to

\begin{array}{rcl} ξ = (n_{0} a^{3}) {(\frac{v_{0}}{a^{3}})}^{2} {(κ a)}^{2} β v (a) ≪ 1, \end{array}

(37)

while the validity of the linearization ansatz, equation (18), together with the solution of the linearized equations, equations (25) and (20), can be expressed as

\begin{array}{rcl} ξ α (\frac{a}{b}) e^{- κ b} ≪ 1, \end{array}

(38)

where b would correspond to an effective core radius of the interstitial test particle. Clearly, assuming ξ ≪ 1, so that mean field is valid, implies that the linearization of the mean field equations is also allowed if the screening length is much larger then the core radius of the colloid. Assuming the values relevant to experiment [47], b ∼ .1 μm, a ∼ 1 − 4 μm, and the effective salt concentration estimated from the dissociation of the colloidal particles to be below ∼μM, both above conditions are satisfied.

As we noted above, outside the stable regime of parameter values the MF solutions would formally acquire a periodic component, perhaps indicating a delocalized-to-localized transition in the mobile ion density. It remains to be explored either within a strong-coupling framework or in detailed simulations that fully include the electrostatic interactions.

5. Results and discussion

We now analyze the general results derived above for the electrostatic field and elastic displacement field of a single test particle, and also the interaction between two test particles located at a finite separation, in terms of the parameters of the model. The point-like test particles are characterized by an electrostatic charge e₀ and an elastic ‘charge' v₀(λ + 2μ) that are both taken as parameters of the model.

The radial component of the electrostatic field of a single point source can be cast in the form

\begin{array}{rcl} \begin{array}{rcl} E_{r} (R, ξ) & = & - \frac{d (β e_{0} ϕ (r))}{d (κ_{0} r)} \\ = & \frac{{e_{0}}^{2} e^{- κ_{0} r}}{4 π ε r} \frac{1}{1 - ξ} (\frac{1 + \frac{κ_{0} r}{\sqrt{1 - ξ}}}{κ_{0} r}) e^{- (\frac{1}{\sqrt{1 - ξ}} - 1) κ_{0} r} \\ = & F_{D H} (r) {\bar{E}}_{r} (R, ξ), \end{array} \end{array}

(39)

where the radial electrostatic field E_r(R, ξ) can be written in a dimensionless form

{\bar{E}}_{r} (R, ξ)

as a function of variables R = κ₀r and

ξ = 2 β n_{0} v_{0}^{2} (λ + 2 μ)

. Its dependence E_r(R, ξ) is presented in figure 2.

Figure 2. The behaviors of the dimensionless electrostatic field $\bar{E} (R, ξ)$ , equation (39), as a function of both relevant variables R = κ₀r and $ξ = 2 β n_{0} v_{0}^{2} (λ + 2 μ)$ (inset). Clearly, the dependence on the coupling parameter ξ is non-monotonic in the inset, but the dependence on the dimensionless radial distance from the point source R is monotonic, and the decay is the fastest for the pure electrostatic system ξ = 0.

Full size|PPT slide

Next, we characterize the radial component of the elastic displacement of the colloidal lattice as follows

\begin{array}{rcl} \begin{array}{rcl} u_{r} (R, ξ) & = & u (r) \cdot \frac{r}{r} \\ = & - \frac{v_{0}}{4 π r^{2}} (\frac{1}{1 - ξ}) (1 + \frac{κ_{0} r}{\sqrt{1 - ξ}}) e^{- \frac{κ_{0} r}{\sqrt{1 - ξ}}} \\ = & v_{0} κ_{0}^{2} {\bar{u}}_{r} (R, ξ), \end{array} \end{array}

(40)

where

{\bar{u}}_{r} (R, ξ)

is a dimensionless form as a function of variables R = κ₀r and

ξ = 2 β n_{0} v_{0}^{2} (λ + 2 μ)

and quantifies the elastic displacement of the colloidal lattice induced by the electrostatic and elastic interactions. Its dependence

{\bar{u}}_{r} (R, ξ)

is presented in figure 3.

Figure 3. The behaviors of the dimensionless radial elastic displacement field ${\bar{u}}_{r} (R, ξ)$ , equation (40), as a function of both variables R = κ₀r and $ξ = 2 β n_{0} v_{0}^{2} (λ + 2 μ)$ (inset). Clearly, the dependence on the coupling parameter ξ is non-monotonic, while the dependence on the dimensionless radial distance from the point source is monotonic, and the increase is fastest for the pure electrostatic system.

Full size|PPT slide

We finally analyze the consequences of the electrostatic-elastic coupling on the interaction free energy, equation (29), between two test particles characterized by the same e₀, v₀ but located D apart. Comparing the lattice-mediated electrostatic interaction with the standard DH interaction, we write the result in the following form

\begin{array}{rcl} \begin{array}{rcl} Δ F (D) & \equiv & F (D) - F_{D H} (D) \\ = & \frac{{e_{0}}^{2} e^{- κ_{0} D}}{8 π ε D} [{(\frac{1}{1 - ξ})}^{2} e^{- (κ - κ_{0}) D} - 1] \\ = & F_{D H} (D) \bar{W} (\tilde{R}, ξ), \end{array} \end{array}

(41)

where

\tilde{R} = κ_{0} D

, and

\bar{W} (\tilde{R}, ξ)

, the relative change with respect to the DH interaction, can be cast into a dimensionless form as a function of dimensionless variables

\tilde{R}

and ξ, quantifying the relative change in the DH interaction energy due to electrostatic-elastic coupling. It can be either positive,

F (D) ⩾ F_{D H} (D)

, corresponding to repulsive elastic coupling-induced interactions, or negative,

F (D) ⩽ F_{D H} (D)

, corresponding to attractive elastic coupling-induced interactions, dependent on both the coupling strength ξ and the distance D.

In order to understand the phase diagram better, we take the logarithm of both sides of

\bar{W} (\tilde{R}, ξ) = 0

to further obtain

\begin{array}{rcl} - 2 l n (1 - ξ) = (\frac{1}{\sqrt{1 - ξ}} - 1) \tilde{R} . \end{array}

(42)

Let

1 / \sqrt{1 - ξ} = 1 + δ

, then the above equation can be rewritten as

\begin{array}{rcl} 4 l n (1 + δ) = δ \tilde{R} . \end{array}

(43)

With ξ approaching 0, δ also approaches 0 and

\tilde{R}

thus goes to a universal ‘critical' constant

{\tilde{R}}^{*} = 4

, which does not depend on any parameters or variables (see the green dot indicated by an arrow in figure 4).

Figure 4. Behaviors of the dimensionless total electrostatic-elastic coupling free energy $\bar{W} (\tilde{R}, ξ)$ , equation (41), as a function of both variables $\tilde{R} = κ_{0} D$ and $ξ = 2 β n_{0} v_{0}^{2} (λ + 2 μ)$ . Under a moderate dimensionless distance $\tilde{R} < 4.0$ , the dimensionless total energy shows a non-monotonic dependence on the electrostatic-elastic coupling parameter ξ, but at large enough distance $\tilde{R} > 4.0$ , it changes to a monotonic decaying dependence on the coupling parameter ξ and the interaction is purely attractive (Inset (a)). The smallest attainable value of $W (\tilde{R}, ξ)$ is obviously −1, which corresponds to the vanishing of the electrostatic-elastic coupling interaction, with the attractive contribution compensating the repulsive contribution. But the dependence on the dimensionless radial distance from the point source is decreasing, and decaying is slower as the electrostatic-elastic coupling constant is weaker. Finally, it vanishes for the pure electrostatic system ξ = 0 (Inset (b)), and the ‘phase diagram' in the parameter space delineates a repulsive region in an otherwise attractive sea. The magenta dotted line shows our theory goes back to the standard DH theory for the coupling strength ξ = 0. The blue dashed line indicates a ‘critical line' ${\tilde{R}}^{}$ as a function of ξ separating the repulsive and attractive regions. A universal ‘critical' constant ${\tilde{R}}^{} = 4$ is indicated by a arrow for ξ = 0.

Full size|PPT slide

Inset (a) in figure 4 displays the dependence of

\bar{W} (\tilde{R}, ξ)

on ξ for different values of

\tilde{R}

. Clearly

\bar{W} (\tilde{R}, ξ)

corresponds to additional electrostatic-elastic coupling-generated interactions, apart from and on top of the usual DH interactions.

\bar{W} (\tilde{R}, ξ) = - 1

indicates that there are no interactions, so the DH repulsion is completely compensated by the additional electrostatic-elastic coupling interaction.

\bar{W} (\tilde{R}, ξ) ⩾ 0

implies additional interactions that are repulsive and

- 1 ⩽ \bar{W} (\tilde{R}, ξ) ⩽ 0

corresponds to additional attractive electrostatic-elastic coupling interactions that, however, do not overpower the underlying DH repulsion.

\bar{W} (\tilde{R}, ξ)

starts off at zero at ξ = 0, and then depending on

\tilde{R}

either develops a local negative minimum

\bar{W} = - 1

, for

\tilde{R} < 4.0

, or spikes into a local positive maximum for

\tilde{R} > 4.0

until again settling down at

\bar{W} = - 1

for ξ = 1. Furthermore, as can be seen from figure 4,

\bar{W} (\tilde{R}, ξ) ⩾ 0

is possible only for

\tilde{R} ⩽ {\tilde{R}}_{c} = 4.0

, which represents in some sense a ‘critical line'

{\tilde{R}}^{*}

as a function of ξ of the electrostatic-elastic coupling interaction, separating a regime of

F (D) ⩽ F_{D H} (D)

from

F (D) ⩾ F_{D H} (D)

. For

\tilde{R} > {\tilde{R}}_{c}

- 1 ⩽ \bar{W} (\tilde{R}, ξ) ⩽ 0

corresponds to a regime where the electrostatic-elastic coupling contributes a purely attractive additional interaction for all values of ξ.

Inset (b) in figure 4 displays a monotonic dependence of

\bar{W} (\tilde{R}, ξ)

\tilde{R}

for different values of ξ. Clearly,

\bar{W} (\tilde{R}, ξ)

decays exponentially, but depending on ξ it can show a regime where the electrostatic-elastic coupling contributes a purely repulsive additional interaction for small

\tilde{R}

, reverting invariably to an attractive regime of additional interaction for large enough

\tilde{R}

. In addition, the larger ξ is, the faster is the decay of the additional electrostatic-elastic coupling-generated interaction. The screening length can be in fact much smaller then the Debye length, so the electrostatic-elastic coupling effect can be short-range.

A simple continuum model of electrostatic-elastic coupling between small, mobile ions and a fixed colloidal lattice, based on the mean-field electrostatic energy and Hookeian elasticity, coupled by minimal coupling terms and solved in the linearized mean-field framework, can thus lead to an attractive additional component in the overall interaction between two test particles in the colloidal lattice. While this lattice-mediated interaction evaluated on the mean-field level cannot overpower the underlying electrostatic repulsion between charges, it can effectively eliminate it, i.e. when

\bar{W} (\tilde{R}, ξ) = - 1

. On the level of this theory, the electrostatic-elastic coupling can thus substantially diminish the repulsion between two mobile charges but cannot completely reverse its sign.

Our approach is necessarily approximate and model-dependent, as such being open to generalizations in different directions. While the continuum elastic lattice approximation allows for straightforward analytical developments, the connection between the colloid–colloid interaction potential

υ

(∣r_i − r_j∣), see equation (3), and the interstitial Coulomb fluid is not direct, requiring an additional, even if very reasonable, assumption that it is of a screened Coulomb type. Note also that the colloid–colloid interaction potential depends on the assumed value of the dissociated charge of the colloids, which can include the effects of distorted double layers at high macroion volume fractions, but it is also difficult to disentangle from experiments and even more from a detailed calculation based on explicit charge dissociation models [45]. Using the phenomenological elastic constants as opposed to phenomenological effective colloidal charges seems in our view more appropriate and going beyond would necessitate a different formulation that would loose the intuitive clarity of the Hookeian elasticity.

There are two obvious generalizations that one could pursue in a more detailed description of this system. The first one would be to alter the model and bypass the continuum Hookeian lattice approximation. Instead one would consider explicit finite-size colloids on a lattice and solve the PB equation in the interstitial space with the proper boundary conditions [11, 36], possibly including also charge regulation [45] at the surface of the macroions. While this could allow for a more detailed inclusion of the finite-size effects [48], the theory would have to be formulated differently possibly lacking a direct connection with Hookeian elasticity.

The next step could be to go beyond the mean-field approximation. As the colloids are strongly charged, the colloidal crystal can certainly be conceived as a strongly coupled electrostatic system in the strict sense of the weak–strong coupling dichotomy [15]. The strong coupling in this context refers to the coupling between mobile ions and the lattice colloids, and the putative extension of our theory would then follow the very same pathway as in the standard strong-coupling formalism [18]. It seems reasonable to expect that on the level of the strong-coupling theory one would end up with a net attraction between test particles, but with possibly even more complicated condensed structures of mobile charges within a fixed, charged colloidal lattice. Only after the electrostatics is properly included and modeled, the simulation as the comparison with our theory makes sense. In addition, it would be interesting to compare not only the relative probability of tracer pairs in the background colloid lattice, but an actual interaction potential that could then be compared directly with equation (41).

Finally, our theory proposes a new mechanism for generating effective attractive interactions between equally charged ions in a colloid crystalline lattice that stems from the coupling between electrostatics and lattice elasticity, proposing a new mechanism that could prove to be important for colloid physics.

References

List( Publishing order | Descend order by publishing year | Descend order by cited within ) Chart analysis

1	Girard M, Wang S, Du J S, Das A, Huang Z, Dravid V P, Lee B, Mirkin C A, de la Cruz M O 2019 Particle analogs of electrons in colloidal crystals Science 364 1174 1178 https://doi.org/10.1126/science.aaw8237 Cited in this article [4]

2	Williams R, Crandall R S 1974 The structure of crystallized suspensions of polystyrene spheres Phys. Lett. A 48 225 226 https://doi.org/10.1016/0375-9601(74)90555-6 Cited in this article [2]

3	Leunissen M E, van Blaaderen A 2008 Concentrating colloids with electric field gradients. II. Phase transitions and crystal buckling of long-ranged repulsive charged spheres in an electric bottle J. Chem. Phys. 128 164509 https://doi.org/10.1063/1.2909200

4	Sirota E B, Ou-Yang H D, Sinha S K, Chaikin P M, Axe J D, Fujii Y 1989 Complete phase diagram of a charged colloidal system: a synchrotron x-ray scattering study Phys. Rev. Lett. 62 1524 1527 https://doi.org/10.1103/PhysRevLett.62.1524

5	Russell E R, Spaepen F, Weitz D A 2015 Anisotropic elasticity of experimental colloidal wigner crystals Phys. Rev. E 91 032310 https://doi.org/10.1103/PhysRevE.91.032310 Cited in this article [2]

6	Dotera T, Oshiro T, Ziherl P 2014 Mosaic two-lengthscale quasicrystals Nature 506 208 211 https://doi.org/10.1038/nature12938 Cited in this article [1]

7	Groenewold J, Kegel W K 2001 Anomalously large equilibrium clusters of colloids J. Phys. Chem. B 105 11702 11709 https://doi.org/10.1021/jp011646w Cited in this article [1]

8	Lechner W, Dellago C 2009 Defect interactions in two-dimensional colloidal crystals: vacancy and interstitial strings Soft Matter 5 2752 2758 https://doi.org/10.1039/b903091c Cited in this article [3]

9	Lechner W, Dellago C 2009 Point defects in two-dimensional colloidal crystals: simulation versus elasticity theory Soft Matter 5 646 659 https://doi.org/10.1039/B810554E Cited in this article [10]

10	Dinsmore A D, Crocker J C, Yodh A G 1998 Self-assembly of colloidal crystals Curr. Opin. Colloid Interface Sci. 3 5 11 https://doi.org/10.1016/S1359-0294(98)80035-6 Cited in this article [1]

11	Schmitz K S 1999 Many-bodied effects and the structure of colloidal crystals Phys. Chem. Chem. Phys. 1 2109 2117 https://doi.org/10.1039/a808229d Cited in this article [2]

12	Lin Y, de la Cruz M O 2022 Superionic colloidal crystals: Ionic to metallic bonding transitions J. Phys. Chem. B 126 6740 6749 https://doi.org/10.1021/acs.jpcb.2c04041 Cited in this article [2]

13	Podgornik R, Saslow W M 2005 Long-range many-body polyelectrolyte bridging interactions J. Chem. Phys. 122 204902 https://doi.org/10.1063/1.1908870 Cited in this article [1]

14	Markovich T, Andelman D, Podgornik R 2021 Charged membranes: Poisson–Boltzmann theory, the DLVO paradigm, and beyond Handbook of Lipid Membranes Safynia C R, Raedler J O London Taylor & Francis Cited in this article [4]

15	Boroudjerdi H, Kim Y-W, Naji A, Netz R R, Schlagberger X, Serr A 2005 Statics and dynamics of strongly charged soft matter Phys. Rep. 416 129 199 Cited in this article [3]

16	Smith A M, Maroni P, Borkovec M 2018 Attractive non-DLVO forces induced by adsorption of monovalent organic ions Phys. Chem. Chem. Phys. 20 158 164 https://doi.org/10.1039/C7CP06383K Cited in this article [2]

17	Ludwig M, von Klitzing R 2020 Recent progress in measurements of oscillatory forces and liquid properties under confinement Curr. Opin. Colloid Interface Sci. 47 137 152 https://doi.org/10.1016/j.cocis.2020.02.002 Cited in this article [2]

18	Naji A, Kanduč M, Forsman J, Podgornik R 2013 Perspective: Coulomb fluids—weak coupling, strong coupling, in between and beyond J. Chem. Phys. 139 150901 https://doi.org/10.1063/1.4824681 Cited in this article [3]

19	Podgornik R, Andelman D 2020 Embarras de richesses in non-DLVO colloid interactions Journal Club for Condensed Matter Physics https://doi.org/10.36471/JCCMseptember202002 Cited in this article [1]

20	Trefalt G, Palberg T, Borkovec M 2017 Forces between colloidal particles in aqueous solutions containing monovalent and multivalent ions Curr. Opin. Colloid Interface Sci. 27 9 17 https://doi.org/10.1016/j.cocis.2016.09.008 Cited in this article [1]

21	Moazzami-Gudarzi M, Adam P, Smith A M, Trefalt G, Szilágyi I, Maroni P, Borkovec M 2018 Interactions between similar and dissimilar charged interfaces in the presence of multivalent anions Phys. Chem. Chem. Phys. 20 9436 9448 https://doi.org/10.1039/C8CP00679B Cited in this article [1]

22	Sutton A P 2020 Physics of Elasticity and Crystal Defects Oxford Oxford University Press Cited in this article [3]

23	Lechner W, Polster D, Maret G, Keim P, Dellago C 2013 Self-organized defect strings in two-dimensional crystals Phys. Rev. E 88 060402 https://doi.org/10.1103/PhysRevE.88.060402 Cited in this article [1]

24	Eisenmann C, Gasser U, Keim P, Maret G, von Grünberg H-H 2005 Pair interaction of dislocations in two-dimensional crystals Phys. Rev. Lett. 95 185502 https://doi.org/10.1103/PhysRevLett.95.185502

25	He B, Chen Y 2013 Interaction energy and point-defect configurations in two-dimensional colloidal crystals Solid State Commun. 159 60 64 https://doi.org/10.1016/j.ssc.2013.01.024 Cited in this article [1]

26	Pertsinidis A, Ling X S 2001 Equilibrium configurations and energetics of point defects in two-dimensional colloidal crystals Phys. Rev. Lett. 87 098303 https://doi.org/10.1103/PhysRevLett.87.098303 Cited in this article [1]

27	Lechner W, Büchler H-P, Zoller P 2014 Role of quantum fluctuations in the hexatic phase of cold polar molecules Phys. Rev. Lett. 112 255301 https://doi.org/10.1103/PhysRevLett.112.255301 Cited in this article [1]

28	Bardeen J, Cooper L N, Schrieffer J R 1957 Theory of superconductivity Phys. Rev. 108 1175 https://doi.org/10.1103/PhysRev.108.1175 Cited in this article [1]

29	Wu H, Ou-Yang Z-C, Podgornik R 2024 Electrostatic-elastic coupling in colloidal crystals Europhys. Lett. 148 47001 https://doi.org/10.1209/0295-5075/ad8f5d Cited in this article [2]

30	Xie Q, Huang M 1994 Application of lattice inversion method to embedded-atom method Phys. Status Solidi (b) 186 393 402 https://doi.org/10.1002/pssb.2221860207 Cited in this article [1]

31	Löwen H, Hansen J-P, Madden P A 1993 Nonlinear counterion screening in colloidal suspensions J. Chem. Phys. 98 3275 3289 https://doi.org/10.1063/1.464099 Cited in this article [1]

32	Maggs A C, Podgornik R 2016 General theory of asymmetric steric interactions in electrostatic double layers Soft Matter 12 1219 1229 https://doi.org/10.1039/C5SM01757B Cited in this article [1]

33	Schwerdtfeger P, Burrows A, Smits O R 2021 The Lennard-Jones potential revisited: analytical expressions for vibrational effects in cubic and hexagonal close-packed lattices J. Phys. Chem. A 125 3037 3057 https://doi.org/10.1021/acs.jpca.1c00012 Cited in this article [1]

34	Ashcroft N W, Mermin N D 1976 Solid State Phys. Cengage Learning Cited in this article [1]

35	Kung W 2009 Geometry and Phase Transitions in Colloids and Polymers World Scientific Cited in this article [2]

36	Dobnikar J, Chen Y, Rzehak R, von Grünberg H H 2002 Many-body interactions in colloidal suspensions J. Phys. Condens. Matter 15 S263 https://doi.org/10.1088/0953-8984/15/1/335 Cited in this article [2]

37	Teodosiu C 1982 The Elastic Field of Point Defects Springer Cited in this article [1]

38	Kardar M, Golestanian R 1999 The ‘friction' of vacuum, and other fluctuation-induced forces Rev. Mod. Phys. 71 1233 1245 https://doi.org/10.1103/RevModPhys.71.1233 Cited in this article [1]

39	Karimi F, Haddadan P, Naji A, Podgornik R 2019 Casimir-like interactions and surface anchoring duality in bookshelf geometry of smectic-a liquid crystals Soft Matter 15 2216 2222 https://doi.org/10.1039/C8SM02328J Cited in this article [1]

40	Walz C, Fuchs M 2010 Displacement field and elastic constants in nonideal crystals Phys. Rev. B 81 134110 https://doi.org/10.1103/PhysRevB.81.134110 Cited in this article [1]

41	Landau L D, Lifshitz E M 1970 Theory of Elasticity Course of Theoretical Physics Vol. 7 Pergamon Cited in this article [1]

42	Verwey E J W, Overbeek J T G 1948 Theory of the Stability of Lyophobic Colloids Elsevier Cited in this article [1]

43	Wiegel F W 1975 Path integral methods in statistical mechanics Phys. Rep. 16 57 114 https://doi.org/10.1016/0370-1573(75)90030-7 Cited in this article [2]

44	Blossey R, Maggs A C 2018 A fluctuation-corrected functional of convex Poisson–Boltzmann theory J. Phys. A: Math. Theor. 51 385001 https://doi.org/10.1088/1751-8121/aad352 Cited in this article [2]

45	Avni Y, Andelman D, Podgornik R 2019 Charge regulation with fixed and mobile charged macromolecules Curr. Opin. Electrochem. 13 70 77 https://doi.org/10.1016/j.coelec.2018.10.014 Cited in this article [3]

46	Naji A, Kanduč M, Netz R R, Podgornik R 2011 Exotic electrostatics: unusual features of electrostatic interactions between macroions Understanding Soft Condensed Matter via Modeling and Computation Singapore World Scientific Cited in this article [1]

47	Liu X, Stenhammar J, Wennerström H, Sparr E 2022 Vesicles balance osmotic stress with bending energy that can be released to form daughter vesicles J. Phys. Chem. Lett. 13 498 507 https://doi.org/10.1021/acs.jpclett.1c03369 Cited in this article [1]

48	Kornyshev A A 2007 Double-layer in ionic liquids:paradigm change J. Phys. Chem. B 111 5545 5557 https://doi.org/10.1021/jp067857o Cited in this article [1]

Acknowledgments

HW is partially supported by the open research fund of Songshan Lake Materials Laboratory No. 2023SLABFN20, the General Program of National Natural Science Foundation of China (NSFC) under Grant No. 12374210 and the startup fund under Grant No. WIUCASQD2022005 from Wenzhou Institute University of Chinese Academy of Sciences (WIUCAS). Z-CO-Y was supported by the Major Program of the NSFC under Grant No. 22193032. RP acknowledges the support of UCAS and funding from the Key Program of NSFC under Grant No. 12034019.