Strong effect of weak charging in suspensions of anisotropic colloids

Abstract

Suspensions of hard colloidal particles frequently serve as model systems in studies on fundamental aspects of phase transitions. But often colloidal particles that are considered as “hard” are in fact weakly charged. If the colloids are spherical, weak charging has only a weak effect on the structural properties of the suspension, which can be easily corrected for. However, this does not hold for anisotropic particles. We introduce a model for the interaction potential between charged ellipsoids of revolution (spheroids) based on the Derjaguin approximation of Debye–Hückel theory and present a computer simulation study on aspects of the system's structural properties and phase behaviour. In line with previous experimental observations, we find that even a weak surface charge has a strong impact on the correlation functions. A likewise strong impact is seen on the phase behaviour, in particular, we find stable cubatic order in suspensions of oblate ellipsoids.

---

1 Introduction

Colloids are widely used as models to study basic questions of statistical mechanics. In particular, “hard” particles that only interact by volume exclusion have been studied intensively since the 1950s.¹ Hard particle systems are appealing because their phase behaviour is of purely entropic origin and they can easily be treated by computer simulation. For example, in the context of liquid crystals, studies of hard ellipsoids, spherocylinders and platelets have provided valuable insight into the basic phase transition mechanisms.^2–5

In one of the first computer simulation studies of a phase diagram of hard, oblate particles, Veerman and Frenkel⁵ observed that the particles arranged parallely in stacks which in turn formed a suprastructure of perpendicular orientations. They named this phase “cubatic”. The existence of the cubatic phase has since been under heated debate, and recently several simulation studies^6–8 showed that the cubatic is always metastable with respect to either the isotropic or the columnar phase for various round hard platelet models (i.e. platelets with circular cross section). In contrast, the cubatic phase is stable for square plates.⁹ Experimentally cubatic order has recently been detected in dispersions of hexagonal, charged plate-like particles.¹⁰ We will show below that for oblate ellipsoids cubatic order becomes stabilized as soon as there is a small surface charge.

Over the past years, experimental methods to synthesize and characterize suspensions of ellipsoidal colloids have been advanced and theoretical predictions have been tested experimentally.^11–16 In 2011, Cohen et al.¹⁷ measured the structural properties of a PMMA ellipsoid system and showed significant differences when comparing their data to theoretical predictions for hard ellipsoids given by Percus–Yevick theory¹⁸ and simulations.¹⁹ In the following, we will test our theoretical treatment of weakly charged ellipsoids against the results of this experimental study. We will show that weak charging, as it is often present in PMMA-colloid suspensions, changes the pair correlations such that they match those observed experimentally and alters the phase diagram considerably.

The paper is organized as follows. We first derive the interaction potential. Then we present simulation results on the positional correlations in the system and compare them to the experimental results of ref. 17. Finally we present a scan through the phase diagram for changing surface charge density and show that the cubatic phase is stabilized.

2 Derivation of the interaction potential

The interaction of weakly charged colloids in an electrolyte suspension can be treated in an adiabatic fashion: one assumes that co- and counterions instantaneously readjust upon a change in the colloidal positions, giving rise to an effective interaction between the colloids, possibly of multi-body nature. If the Debye–Hückel screening length is smaller than the extensions of the colloid, then this effective potential can be well approximated by a sum of two-body terms. For the interaction between two colloids, we use the Debye–Hückel approximation (linearized Poisson–Boltzmann (PB) theory). For two spheroids (ellipsoids with one rotational symmetry axis), the effective potential depends on four variables (the center-to-center distance and three angles) such that an explicit tabulation of the PB solutions, let alone the determination of PB solutions “on the fly” in a simulation code appears forbidding.^20,21 Here, we resort to the venerable Derjaguin approximation which has been often used to calculate effective colloid-wall or colloid–colloid interactions in the literature but we are not aware of its practical use in further simulation or theoretical studies of concentrated solutions involving anisotropic particles. The Derjaguin approximation rests upon the following argument: suppose the free energy of the interaction between two planar walls is known, and its density will be denoted by f(h) where h is the distance between the walls. The interaction potential between two convex bodies (of the same type as the walls) can be approximated by just integrating over geometrically opposing area elements dA (at distance h) where the free energy of interaction between the area elements is given by the wall free energy f(h)dA. The mathematical elaboration of this approximation is given in Appendix A and results in the following expression for the free energy of interaction F(H₀) between two convex bodies with the minimal distance H₀ between their surfacesHere, the product εε′ is given by

εε′ = ε₁ε^′₁ + ε₂ε^′₂ + (ε₁ε₂ + ε^′₁ε^′₂)sin² ω + (ε₁ε^′₂ + ε^′₁ε₂)cos² ω. (2)

It involves the principal curvatures ε_i, ε^′_i of the surface of body i = 1, 2 in the planes tangential to the distance vector between the bodies, and also the angle ω between the coordinate systems in the tangential planes with coordinate axes given by the directions of the principal curvatures. For these geometric definitions, see Fig. 1.

Fig. 1 Two surfaces, separated by the minimal distance H₀ = dist(OO′). Around the points O and O′, the surfaces can be approximated by the quadratic forms z_i = (ε_i/2)x_i² + (ε^′_i/2)y_i² where ε_i, ε^′_i are the principal curvatures of surface i at point O resp. O′ and x_i, y_i are coordinate axes in the direction of the principal curvatures. Most generally, the coordinate axes x₁ and x₂ include an angle ω.

2.1 Charged ellipsoids

We will apply these ideas to the interaction between hard, charged ellipsoids (spheroids) with main axes a and b where b is the main axis in the plane perpendicular to the rotational symmetry axis. The aspect ratio is given by t = a/b. In Debye–Hückel approximation, the electrostatic potential ψ fulfills

Δψ − κ²ψ = 0, (3)

where κ⁻¹ is the Debye–Hückel screening length. For a charged wall with charge density σ, the solution is

ϕ_w(z) = ϕ₀ exp(−κz) (4)

with the wall contact potential ϕ₀ = σ/(ε_sκ) (ε_s is the dielectric constant of the solvent). We approximate the solution for two charged walls at distance h by

ϕ_2w(z) ≈ ϕ_w(z) + ϕ_w(h − z) (5)

and the pressure (force density per unit area between the plates) f̃_2w is obtained most easily by evaluating the stress tensor at the midplane z = h/2 which has there only a contribution from the ion osmotic pressure:

The free energy density f is then found through integration

Using this, the Derjaguin free energy (1) becomes

The Derjaguin free energy decays exponentially with H₀, as expected for the Debye–Hückel approximation. This decay is fast enough that the approximation is accurate enough for practical purposes, see Appendix B for a discussion. Note that the anisotropy in the free energy has two sources: H₀ depends on the different orientations as well as the curvature term . The latter one has a strong influence on the interaction of oblate ellipsoids (see below). We write the Derjaguin free energy as

Note that V(H₀) is dimensionless in the last equation. The prefactor (with dimension of energy) in eqn (9), V₀ = 2σ²b/(ε_sκ²), contains the charge density σ as a parameter. It is advantageous to introduce the dimensionless charge density [small sigma, Greek, tilde] = σeβ/(κε_s) (where e is the elementary charge and β = 1/(k_BT) is the inverse temperature) as well as the Bjerrum length of the solvent, λ_B = βe²/(4πε_s). With these definitions the prefactor becomes

Below, in the comparison with the experiments of ref. 17, we will treat [small sigma, Greek, tilde] as a fitting parameter. Regarding the interpretation of this value one should keep in mind that it is an effective or renormalized charge density. For [small sigma, Greek, tilde] ≲ 1, the bare and the renormalized charge density are approximately the same, whereas in the limit of very large bare charge densities the renormalized charge density approaches a constant, [small sigma, Greek, tilde] → 4.²² Thus the approach is consistent only for [small sigma, Greek, tilde] < 4.

2.2 Numerical implementation

We computed numerically the dimensionless, exponentiated free energy exp(−V(H₀)) (eqn (9)) on a 4-dimensional grid with axes characterizing the relative configurational state of two ellipsoids (center-to-center distance and three angles). For each configuration, the minimal distance H₀ was determined by a conjugate gradient routine and the radii and directions of principal curvature were determined through the first and second fundamental form of the ellipsoid surfaces at the two points O and O′ (whose distance is H₀, see Fig. 1). In the Monte Carlo simulations (see below), the such tabulated free energy was used as the acting potential between pairs of ellipsoids, together with linear interpolation to determine the potential at off-grid values of the variables characterizing the relative configurational state.

We also tested a further approximation to the Derjaguin free energy in which the curvature term is replaced by a constant, the average radius of the ellipsoid. Then V(H₀) only depends on the minimal distance H₀, which can be well approximated by an extension of the Perram–Wertheim routine²³ frequently used for checking overlap of hard ellipsoids (see Appendix C). In this way, the potential can be determined “on the fly”, and it works reasonably well for aspect ratios 0.8 ≲ t ≲ 2. Note, however, that the short-range anisotropy of the potential increases rapidly with the aspect ratio becoming small. For oblate ellipsoids (t < 1, disk-like particles), the ratio between the potential at contact in side–side configuration (flat sides of the disk touching) and in edge–edge configuration (rims of the disk are touching) is (t² + 1)/(2t³) and thus scales for small t as 1/t³. Therefore, this further approximation to the Derjaguin free energy is not applicable to flat oblates.

3 Simulation results

We have carried out Monte Carlo simulations at constant temperature T, constant number of particles N and volume V with periodic boundary conditions, and computed equilibrium structural properties of the system as a function of the effective surface charge density and the packing fraction. The particle number N ranged from 3000 to 3840, T was set to 300 K.

First, we discuss the structure of the isotropic phase. This part of our work has been inspired by recent experimental measurements of the radial, orientation-averaged pair correlation function g(r) in suspensions of prolate ellipsoids (aspect ratio t = 1.6, a = 3.2 μm, b = 2.0 μm) and oblate ellipsoids (t ≈ 0.25, a ≈ 0.96 μm b ≈ 3.8 μm, with a considerable experimental uncertainty on the polydispersity and thus on the value of t).¹⁷

The structural correlations that are presented in ref. 17 are much stronger than one would expect for a system of hard ellipsoids, and this was taken as an indication that on the theory side, the correlations in suspensions of ellipsoids are not sufficiently well understood. However, the experimental suspension was additionally stabilized by a surfactant which introduced a small amount of charge on the particles. In a later study,¹² the authors investigated the influence of charge on g(r) for the prolate particles by simulation and found it non-negligible: with a small charge density of σ ≈ 9e μm⁻² (where e = 1.6 × 10⁻¹⁹ C is the elementary charge), distributed on particles modelled by an assembly of three cut spheres to approximate the shape of the ellipsoids, the experimental g(r) could be reproduced. No corresponding results for the oblate particles have been reported, though.

In order to model the experiment of ref. 17, we set λ_B ≈ 22 nm appropriate for a solvent with an average dielectric constant of ε_s = 2.5 and κ⁻¹ = 0.3 μm for the Debye–Hückel screening length. We only vary the effective surface charge density to reproduce the experimental data. Fig. 2(top) shows the radial positional distribution function g(r) of ellipsoids (circles) for an aspect ratio t = 1.6 and a packing fraction ϕ = 0.31. The simulation data perfectly match the experimental data (triangles, from Fig. 2 in ref. 17). The corresponding dimensionless charge density is given by [small sigma, Greek, tilde] = 0.83, i.e. the effective charge density is σ = 10e μm⁻². This value is reasonable for the experimental system used in ref. 17, it is in good agreement with the deduced effective charge on the colloids of ref. 12, and it justifies in retrospect the assumption of the Debye–Hückel approximation used to derive the interaction potential. (Note that the modelling of the electrostatic particle interactions in ref. 12 approximately corresponds to the Derjaguin free energy with curvature neglected (see Section 2.2 above)) which works well for the moderate aspect ratio of 1.6 but not for particles with larger curvatures. This rationalizes the good agreement between our results and those of ref. 12 for σ and g(r).

Fig. 2 Radial positional distribution function of ellipsoids with aspect ratio t = 1.6 and packing fractions ϕ = 0.26 (top) and ϕ = 0.31 (bottom), simulation data (circles), experimental data (triangles) and data for hard prolate ellipsoids (HPE, dashed line) taken from Fig. 2(a) and (b) in ref. 17.

The second data set discussed in ref. 17 has been measured in a more dilute suspension of prolate ellipsoids at a packing fraction ϕ = 0.26. Using the same system parameters for λ_B, κ⁻¹ and even σ as for ϕ = 0.31, we again obtain very good agreement of the radial distribution functions, see Fig. 2(bottom). The data set for hard ellipsoids from ref. 17 is presented as well for comparison.

The last radial distribution function that is presented in ref. 17, was measured in a suspension of oblate ellipsoids of an aspect ratio of “t ≈ 0.25” (with a larger polydispersity than in the prolate case). As explained in Section 2.2, the curvature around the rim of the particles in this case is important for the electrostatic interactions, hence the approach of ref. 12 could not be applied here. Fig. 3 shows our simulation results for oblate ellipsoids. We set again the same value for λ_B, κ⁻¹ and σ. The experimental and theoretical data agree reasonably well given the uncertainty of the aspect ratio of the experimental system.

Fig. 3 Radial positional distribution function of oblate ellipsoids with a packing fraction of ϕ = 0.35, simulation data (circles) for t = 0.35, experimental data for “t ≈ 0.25” (triangles) and data for hard oblate ellipsoids (HOE) taken from Fig. 3 in ref. 17 (dashed line, data on the abscissa is multiplied by a factor 1.4 with respect to ref. 17 to undo the rescaling and recover units of b.).

To conclude this section, we validated the Derjaguin approximation for the electrostatic interaction of charged ellipsoids. A small amount of surface charge has a strong influence on the pair correlations, due to the small dielectric constant (large Bjerrum length) of the solvent. The short-range anisotropy of the electrostatic interaction is especially important for oblate ellipsoids.

4 Impact of the surface charges on the nematic phase

We now consider oblate ellipsoids of aspect ratio t = 0.25 at a packing fraction of ϕ = 0.48. In suspensions of hard ellipsoids the nematic phase is located at packing fractions ϕ > 0.4². We study the effect of increasing surface charge density, ranging from σ = 0.00002e μm⁻² to σ = 2.0e μm⁻², on the structural properties of the liquid. Note that these surface charge densities are even lower than the value that was needed to reproduce the experimental findings of ref. 17. The Bjerrum length λ_B and the Debye–Hückel screening length κ⁻¹ are not modified with respect to the previous section.

As the initial configuration of a first set of simulation runs we used an fcc crystal in which the ellipsoids were oriented in parallel. We let the system relax until its energy had reached a stable value. In the case of perfect charge screening (i.e. almost hard ellipsoids) the system relaxed into the expected nematic phase, see Fig. 4. A similar degree of nematic ordering formed for a surface charge density of σ = 0.02e μm⁻². In contrast we observe qualitatively different behaviour for surface charge densities σ = 0.2e μm⁻² and σ = 2.0e μm⁻², see snapshots in Fig. 5 and Fig. 6. Note that these configurations evolved from an initial fcc configuration with parallel orientation of the ellipsoids.

Fig. 4 Snapshot of the suspension of almost hard oblates (σ = 0.00002e μm⁻²) (aspect ratio t = 0.25) at a packing fraction ϕ = 0.48 in the nematic phase. Colour code according to orientation.

Fig. 5 Snapshot of the suspension of charged oblates at surface charge density σ = 0.2e μm⁻² (aspect ratio t = 0.25) at a packing fraction ϕ = 0.48. Colour code according to orientation.

Fig. 6 Snapshot of the suspension of charged oblates at surface charge density σ = 2.0e μm⁻² (aspect ratio t = 0.25) at a packing fraction ϕ = 0.48. Colour code according to orientation.

Fig. 7 shows g(r) for different surface charge densities σ at ϕ = 0.48. In the nematic phase, g(r) has a cusp at a distance r that corresponds to in-plane rim–rim configurations. In contrast, at σ = 0.2e μm⁻² there is pronounced positional order with peak positions at multiples of the length of the small axis a plus a small distance to account for the electrostatic repulsion.

Fig. 7 Radial positional distribution, aspect ratio t = 0.25, for different surface charge densities σ at a packing fraction ϕ = 0.48.

Fig. 8 shows the orientational distribution function

where u_i is the unit vector along the axis of particle i and the average is over all pairs of particles and the canonical ensemble. At small σ g₂(r) decays smoothly to a non-vanishing value at large distances, which is characteristic for nematic ordering. At σ = 0.2e μm⁻² and above there is parallel order at short distances and random orientation at large distances, i.e. g₂(r) decays to zero. The first “perpendicular peak” with g₂ < 0 appears at the distance that corresponds to a configuration in which the rim of one ellipsoid points to the pole of the other. This peak is superposed with the second layer of parallel stacking in g(r). We conclude that stacks of ellipsoids are arranged perpendicular to each other to form cubatic order.

Fig. 8 Radial orientational distribution, aspect ratio t = 0.25, for different surface charge densities σ at a packing fraction ϕ = 0.48.

To test whether this phase is metastable, we then initialized simulations in the nematic phase, the columnar phase and a perfect long-range cubatic phase. All runs equilibrated into the phase discussed above, thus we conclude that it is the most stable.

5 Conclusion

We have added surface charges to hard ellipsoids and treated them numerically using the Derjaguin approximation. With the Derjaguin approximation, the short-range anisotropy of the electrostatic interactions is captured quantitatively correctly as long as the screening length is smaller than the extensions of the particles. Even very small surface charges, as they often are present in experiments on “hard” particles have a strong effect on the structure of the suspension.

We showed for charged oblate ellipsoids that a phase of perpendicularly oriented short stacks (a cubatic) is thermodynamically more stable than the nematic phase, which is stable for uncharged ellipsoids at the same packing fraction. The cubatic phase is also more stable than the crystal and the columnar phase. It should be accessible to experiments on suspensions of PMMA ellipsoids.

Appendix

A Derjaguin approximation: mathematical derivation

The Derjaguin approximation²⁴ involves the following steps: suppose the free energy of the interaction between two planar walls is known, and its density will be denoted by f(h) where h is the distance between the walls. For calculating the interaction potential between two convex bodies (of the same type as the walls) one determines the minimal distance H₀ between the two surfaces and the tangential planes (see Fig. 1). The free energy of interaction between the two bodies is approximated by

F = ∫∫dxdyf(H₀ + z₁ + z₂), (12)

where the integral runs over (one of) the tangential planes and z₁, z₂ are the distances of the point on surface i described by z_i(x, y) to their respective tangential plane (see Fig. 1). If f(h) quickly decays to zero, it is safe to integrate over the whole plane.

This integral is greatly simplified if we approximate the surfaces by quadratic forms around the points O and O′ respectively:

where ε₁, ε^′₁ are the principal curvatures of surface 1 at point O and x₁, y₁ are coordinates in the tangential plane in the direction of the principal curvatures. Likewise for surface 2. The directions of the principal curvatures of surface 1 and 2 do not agree but include an angle ω, thus:

The distance z₁ + z₂ becomes

A = ε₁ + ε₂ cos² ω + ε^′₂ sin² ω (17)

B = ε^′₁ + ε^′₂ cos² ω + ε₂ sin² ω (18)

C = (ε₂ − ε^′₂)cos ωsin ω . (19)

This distance is a quadratic form. We may perform a rotation to another coordinate system x, y where the off-diagonal matrix elements become zero, i.e.

Using that result, the free energy becomes

F = ∫∫dxdyf(H₀ + εx²/2 + ε′y²/2). (21)

Introducing new coordinates r, ϕ via and :

where the second line (which is eqn (1)) follows from the substitution h = H₀ + r²/2. The force (in direction of OO′) between the two bodies is just given by

The product εε′ is the determinant of the matrix in eqn (20) and must be equal to the determinant of the matrix in eqn (16). Thus we find:

εε′ = ε₁ε^′₁ + ε₂ε^′₂ + (ε₁ε₂ + ε^′₁ε^′₂)sin² ω + (ε₁ε^′₂ + ε^′₁ε₂)cos² ω (25)

which is eqn (2).

B Validity of the Derjaguin approximation

We can estimate the validity of the Derjaguin approximation by considering the example of two interacting, charged spheres with radius r₀ and charge density σ at center distance d for which we can compare the “exact” Debye–Hückel result with the corresponding Derjaguin approximated result. In Debye–Hückel approximation, the potential of a single sphere is given byand the charge Q_eff is determined through the boundary condition ∂Φ/∂r|_r=r0 = −σ/ε_s, giving

In superposition approximation (as before), the interaction free energy of two such spheres is given by

where H₀ = d − 2r₀ is the minimal surface-to-surface distance. The Derjaguin approximated free energy follows from eqn (8) using εε′ = 4/r₀²:

It is the limit κr₀ ≫ 1, H₀ ≪ r₀ of eqn (28). The ratio is given by

From this equation one sees that the Derjaguin approximation produces a large relative error for H₀ ≳ 2r₀. Appropriate for the examples studied in Section 3, we set r₀ = 1 μm and the screening length κ⁻¹ = 300 nm. Thus, κr₀ ≈ 1/3, and the Derjaguin approximation overestimates the free energy by a factor 3.4 at H₀ = 2r₀ = 2 μm. However, at that distance the potential has dropped by a factor exp(−2κr₀) ≈ 0.0013 compared to its value at contact, so the free energy itself at that distance is small and likewise the absolute error the Derjaguin approximation produces is small. The effect is presumably negligible in our simulations which are sensitive to the short-range behaviour of the effective free energy between the particles.

C Perram–Wertheim approximation for the minimal distance

Perram and Wertheim²³ developed a criterion to check for overlap of two ellipsoids. The algorithm can also be used to compute an approximate minimal distance. We summarize this approach as follows: suppose we have two ellipsoids with half axes a₁, a₂, a₃ and corresponding axis orientation vectors (of unit length) in a lab-fixed coordinate system u₁, u₂, u₃ (ellipsoid 1) and v₁,v₂, v₃ (ellipsoid 2). One defines two matrices:

Let r_a and r_b denote the center positions of ellipsoid 1 and 2, respectively. We define quadratic forms

F_A = (r − r_a)^TA(r − r_a), (33)

F_B = (r − r_b)^TB(r − r_b), (34)

and the ellipsoid surfaces are given by the solutions to the equations F_A = 1 and F_B = 1. For points inside the ellipsoid, F_A[B] < 1, for points outside F_A[B] > 1. Further we define a quadratic form

F(r, λ) = λF_A + (1 − λ)F_B. (35)

and its minimum, depending on λ ∈ [0, 1]:

F(r(λ), λ) = min_rF(r, λ). (36)

For λ = 0, r(0) = r_b and for λ = 1, r(1) = r_a. Thus as λ varies from 0 to 1, then r(λ) moves from the center of ellipsoid 2 to the center of ellipsoid 1. If the two ellipsoids do not overlap, then there exists a particular λ for sure for which r(λ) is outside both ellipsoids and F(r(λ), λ) > 1 there. If the ellipsoids overlap, then F(r, λ) < 1 in the overlap region and thus for each λ the minimal point r(λ) can not lie outside both ellipsoids since there F(r, λ) > 1. Therefore a useful overlap criterion is formulated with introducing

s = max_λ∈(0,1)F(r(λ), λ) (37)

which fulfills

This criterion is convenient to use due to the explicit form for F(r(λ), λ) which Perram and Wertheim provide:²³

F(r(λ), λ) = λ(1 − λ)(r_b − r_a)^TC(r_b − r_a), (39)

C = (λB⁻¹ + (1 − λ)A⁻¹)⁻¹. (40)

The maximization needed in eqn (37) has to be done numerically, though.

The value of s can also be used to calculate an approximative minimal distance through

The interpretation, according to Paramonov and Yaliraki,²⁵ is as follows. First, the geometrical meaning of the λ maximization in eqn (37) becomes clear by looking at

The term ∇F is zero by virtue of the definition of F(r(λ), λ) in eqn (36) and thus the derivative above is zero when s = F_A(r(λ_max), λ_max) = F_B(r(λ_max), λ_max). This, however, describes the tangential contact between ellipsoids with scaled half-axes . (r(λ_max) is the tangential contact point since ∇F_A || ∇F_B there.) Then the points s_a, s_b defined by

(likewise for a → b) lie on the surface of ellipsoid 1 and 2, respectively. We see that s_b, s_a is parallel to the center distance vector r_b − r_a and that d_PW = |s_b, s_a|. Thus d_PW is the minimal directional distance between the ellipsoids (i.e. minimal distance between two points on the surfaces of ellipsoid 1 and 2 in the direction of the center distance vector.

This approximation is not exact but as presented in Fig. 9 the error is small and hence the approach using the simple determination of the Perram–Wertheim distance justified for the calculation.

Fig. 9 For given configurations of two prolates of aspect ratio t = 1.6, we plot the minimal distance relative to the Perram–Wertheim distance.

The relation is not one to one because of the different relative orientations but there is agreement for moderate aspect ratios.

Acknowledgements

This project has been financially supported by the DFG (SFB Tr6 and SPP1296) and by the National Research Fund, Luxembourg co-funded under the Marie Curie Actions of the European Commission (FP7-COFUND) and under the project FRPTECD. Computer simulations presented in this paper were carried out using the HPC facility of the University of Luxembourg.

References

1. G. Gompper and M. Schick, Complex Colloidal Suspensions, Wiley-VCH Verlag GmbH & Co. KGaA, 2007, vol. 2, pp. 285–287.
2. G. Odriozola, J. Chem. Phys., 2012, 136, 134505.
3. G. Bautista-Carbajal, A. Moncho-Jordá and G. Odriozola, J. Chem. Phys., 2013, 138, 064501.
4. L. Onsager, Ann. N. Y. Acad. Sci., 1949, 51, 627–659.
5. J. A. C. Veerman and D. Frenkel, Phys. Rev. A: At., Mol., Opt. Phys., 1992, 45, 5632–5648.
6. P. D. Duncan, A. J. Masters and M. R. Wilson, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2011, 84, 011702.
7. M. Marechal, A. Patti, M. Dennison and M. Dijkstra, Phys. Rev. Lett., 2012, 108, 206101.
8. P. D. Duncan, M. Dennison, A. J. Masters and M. R. Wilson, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2009, 79, 031702.
9. M. R. Wilson, P. D. Duncan, M. Dennison and A. J. Masters, Soft Matter, 2012, 8, 3348–3356.
10. S. J. S. Qazi, G. Karlsson and A. R. Rennie, J. Colloid Interface Sci., 2010, 348, 80–84.
11. Z. Zhang, P. Pfleiderer, A. B. Schofield, C. Clasen and J. Vermant, J. Am. Chem. Soc., 2011, 133, 392–395.
12. A. P. Cohen, E. Janai, D. C. Rapaport, A. B. Schofield and E. Sloutskin, J. Chem. Phys., 2012, 137, 184505.
13. J.-W. Kim, R. J. Larsen and D. A. Weitz, Adv. Mater., 2007, 19, 2005–2009.
14. J. J. Crassous, H. Dietsch, P. Pfleiderer, V. Malik, A. Diaz, L. A. Hirshi, M. Drechsler and P. Schurtenberger, Soft Matter, 2012, 8, 3538–3548.
15. I. Martchenko, H. Dietsch, C. Moitzi and P. Schurtenberger, J. Phys. Chem. B, 2011, 115, 14838–14845.
16. A. Mohraz and M. J. Solomon, Langmuir, 2005, 21, 5298–5306.
17. A. P. Cohen, E. Janai, E. Mogilko, A. B. Schofield and E. Sloutskin, Phys. Rev. Lett., 2011, 107, 238301.
18. M. Letz and A. Latz, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1999, 60, 5865–5871.
19. J. Talbot, D. Kivelson, M. P. Allen, G. T. Evans and D. Frenkel, J. Chem. Phys., 1990, 92, 3048–3057.
20. D. Chapot, L. Bocquet and E. Trizac, J. Chem. Phys., 2004, 120, 3969.
21. C. Álvarez and G. Téllez, J. Chem. Phys., 2010, 133, 144908.
22. L. Bocquet, E. Trizac and M. Aubouy, J. Chem. Phys., 2002, 117, 8138–8152.
23. J. W. Perram and M. Wertheim, J. Comput. Phys., 1985, 58, 409–416.
24. B. Derjaguin, Kolloid-Z., 1934, 69, 155–164.
25. L. Paramonov and S. N. Yaliraki, J. Chem. Phys., 2005, 123, 194111.

---