Direct Measurement of Thermodynamic Properties of Colloidal Hard Spheres Measurement of Thermodynamic Properties of Colloidal Hard Spheres

— Direct


INTRODUCTION
The intimate relation between colloids and statistical mechanics was pioneered by Perrin in his determination of Boltzmann's constant k B from microscopy images of sedimenting colloids at low densities [1].At these densities, the equation of state is given by van't Hoff's law for the ideal (osmotic) pressure P = ρk B T , with ρ the colloid number density and T the absolute temperature.The measurement of the equation of state for interacting systems is somewhat more elaborate.In most interacting systems, the structure is dominated by the repulsive part of the potential [2].The simplest model liquid with a purely repulsive potential is a fluid consisting of hard spherical particles.Since the properties of hard spheres are determined by excluded volume effects, geometry plays an important role [3,4].Moreover, the available volume to insert an additional sphere and the surface area of that volume are geometric quantities that are directly related to important thermodynamic quantities [3,5,6].In particular, geometry provides a microscopic and straightforward route to the equation of state, the chemical potential and thus the free energy of the hard sphere system.
The equation of state of colloidal hard spheres, widely accepted as an experimental realization of hard spheres [7], has been measured using several techniques.The static structure factor can be extrapolated to zero wave vector in scattering methods [8][9][10].In this limit, the structure factor is proportional to the compressibility.Furthermore, the pressure was also measured by integrating the sedimentationdiffusion profile of colloidal hard spheres [11][12][13].To directly determine thermodynamic quantities of hard spheres via the microscopic route, the accurate acquirement of experimental hard sphere configurations in three dimensions is essential.Several in colloid science commonly used experimental techniques, like light scattering, are unable to provide such microscopic information.However, the colloidal length scale is also very suitable for direct microscopy as pioneered by Perrin [1].In particular, the combination of model colloids and confocal microscopy has proven to be very powerful in analyzing colloidal systems in three dimensions and on a single-particle level.For example, issues like crystallization and glass formation in concentrated colloidal suspensions have been addressed in great detail using confocal microscopy [14][15][16][17][18][19].Nevertheless, it is non-trivial to measure accurate hard sphere configurations in more diluted suspensions.
Here, we combine fast confocal microscopy and model colloids to obtain three dimensional hard sphere configurations throughout almost the whole fluid region of the hard sphere phase diagram.From these snapshots the available volume to insert an additional sphere and the surface area of that volume are determined.Applying modern statistical mechanics enables us to directly measure the pressure, the A 0 V 0 Figure 1 2D illustration of V 0 (grey area) and A 0 (black lines).The white spheres represent the excluded volume spheres of the dark grey hard spheres.chemical potential and the free energy density of a homogeneous colloidal hard sphere suspension as was first shown by us in [20].Here, we provide further theoretical and experimental details.

STATISTICAL GEOMETRY OF HARD SPHERES
The available volume to insert an additional identical sphere into a fluid of N hard spheres, V 0 , and the surface area of that volume, A 0 , are the key geometric quantities to consider.The geometrical interpretation of V 0 and A 0 is illustrated in Figure 1.For hard sphere fluids, statistical geometry provides exact relations between V 0 and A 0 on the one hand and the chemical potential μ and the pressure P on the other hand [5,6]: Here, d is the particle diameter, V is the total system volume and μ 0 the ideal gas chemical potential.In this section we briefly reproduce the derivation of Equations 1 and 2 within the framework of scaled particle theory following [3].We start with an equilibrium configuration of N hard spheres with diameter d.The radius r of the excluded volume spheres is given by r = (d+d )/2, where d is the diameter of an additional hard sphere we want to add to the system.To accommodate this additional sphere we have to find a cavity that is large enough.The probability of finding such a cavity is given by This probability can also be expressed in terms of the reversible work to create the cavity W(r) as given by fluctuation theory [21]: Combining Equations 3 and 4 yields an expression for the reversible work to create the cavity, in other words, to insert the additional hard sphere of diameter d : If we set the diameter of the N hard spheres and of the additional sphere to zero, the ideal gas chemical potential μ 0 = k B T ln Λ 3 N/V is recovered, which represents the work required to insert a point particle into a fluid of point particles.Here, Λ is the thermal de Broglie wavelength.
If the diameter of the additional sphere is set equal to the diameter of the other spheres, i.e. d = d, the work to insert this sphere, in addition to μ 0 is given by W(d), which is Equation 5 with r = d.Thus, the total reversible work μ = μ 0 + W(d) involved in adding the extra sphere is: This is exactly Equation To derive Equation 2 we consider a spherical cavity of radius r from which the centers of the hard spheres are excluded.The local number density of spheres on the surface of the cavity is given by ρG(r, ρ).G(r, ρ) is the conditional probability density for a center to be at r when it is known that no spheres are within the cavity.It is straightforward to see that G(r, ρ) is nothing else than the radial distribution function g(r, ρ) at contact (r = d): As the momentum transfer between a hard cavity and a hard sphere is purely impulsive, the equilibrium force per unit area on the surface of the cavity is f (r) = k B T ρG(r, ρ) and on the whole surface it is 4πr 2 f (r).Now, the reversible work to increase the radius of the cavity by an amount dr is given by Combining Equations 5 and 8 results in an expression for G(r, ρ) in terms of V 0 (r): From Figure 1 and purely geometrical considerations it is clear that dV 0 = −A 0 (r)dr.Using this and subsequently substituting Equation 9 in Equation 7 yields the following expression for the radial distribution function at contact (r = d): Substituting Equation 10into the equation of state for a hard sphere fluid [22] given by leads to the pressure in terms of V 0 and A 0 : which is just Equation 2. Once the pressure and the chemical potential are known, the central thermodynamic quantity, the free energy density F/V is directly given by:

EXPERIMENTAL SECTION
We used monodisperse polymethyl(methacrylate) spheres consisting of a fluorescent core and a non-fluorescent shell [23,24].The particles were dispersed in a mixture of tetralin (Merck, for synthesis), cis-decalin (Merck, for synthesis) and carbon tetrachloride (Merck, for spectroscopy) which matches the refractive index and the mass density of the colloids.The volume ratios were 36% tetralin, 31.5% cis-decalin and 32.5% carbon tetrachloride.In this mixture, the particles behave as hard spheres [25].The particle diameter d is 1.3 μm.
A sequence of samples with target volume fractions φ ≡ ρv (with v the particle's volume) ranging from 0.05 to 0.45 were prepared and contained in small homemade vials.For that purpose, a stock dispersion of φ ≈ 0.30 was centrifugated in the sample cell until the particles were completely sedimented.The supernatant was removed and the volume fraction at this point was taken to be 0.64, corresponding to random close packing [26].Subsequently, all the different samples were diluted to their target volume fraction.Before acquiring the microscopy data, the samples were thoroughly homogenized and it was verified that there was no flow in the system while imaging.During the experiment, the samples were carefully weighed to check for possible evaporation of the solvent.Furthermore, the mass fractions of the samples were determined afterwards.These were in good agreement with the volume fractions directly obtained from the microscopy data.
The particles were imaged using a Leica DM IRBE inverted microscope equipped with a Perkin Elmer Ultraview RS confocal imaging system in combination with a ArKr laser (568 nm, Melles Griot) and an oil-immersion lens (100×, NA 1.4, Leica).The frame rate of this Nipkow spinning disk confocal system was 19.6 images per second, which corresponds to imaging a volume of 25 × 25 × 25 μm 3 by successive two-dimensional (2D) slices in about 6 seconds.This time interval is about 4 times the Brownian time of the particles at infinite dilution.The resolution of the 2D images was 384 × 384 pixels and the stepsize in z between the frames was 0.20 μm.For every sample, i.e. at every volume fraction, typically series of 11 z-stacks were taken with a delay time of 12 seconds between the sequential stacks.The delay time being on the order of the Brownian time of the particles, ensured the statistical independence of the stacks.The data were acquired at least 15 μm from the bottom wall to exclude possible effects of the presence of the container wall.The combination of the core-shell character of the colloids and the high-speed confocal microscope enabled us to obtain the particle coordinates in full 3D [27].Hence, the number of particles in known volume directly gave the volume fraction.
The experiment was mimicked by computer simulations.Hard sphere configurations were generated using Monte Carlo simulations with periodic boundary conditions [28].Firstly, the initial configurations consisting of 500 up to 5000 particles were generated in a random fashion.Subsequently, the step length during equilibration was adjusted such that the Monte Carlo acceptance rate was at all volume fractions in the range 0.3-0.5.The obtained hard sphere configurations were analyzed in exactly the same way as the experimental data.
The determination of V 0 and A 0 from the particle configurations is a purely geometrical problem.Although an exact method is reported in [2,29], we use a simple approximate, yet accurate and fast method to measure V 0 and A 0 .First of all, the particle positions are mapped onto a fine lattice.For convenience, here one lattice site is called a pixel.Subsequently, three types of lattice sites are distinguished, from which V 0 and A 0 are directly determined.As illustrated in Figure 2A, V 0 is given by the number of pixels outside the excluded volume spheres (the white pixels).In principle, V 0 can be determined by the random insertion of test particles as well [5,28].However, at higher densities this procedure is very time consuming.The surface area A 0 is given by the number of pixels at the edges of the excluded volume spheres (the dark grey pixels in Fig. 2).For low densities A 0 can also be obtained by taking the numerical derivative of V 0 with respect to differently sized excluded volume spheres.The remaining pixels (light grey in Fig. 2) correspond to the interior of the excluded volume spheres.
For the production of Figure 2A, the dark grey pixels corresponding to A 0 were assigned using a radial distance criterion.However, due to the inefficient scanning of the 3D numerical lattice, this approach turned out to be extremely time consuming when applied to the experimental data sets.(A) 2D illustration of an excluded volume sphere corresponding to a hard sphere with a diameter of 19 pixels, mapped onto the lattice.The white pixels contribute to V 0 , the dark grey pixels are constituting A 0 and the light grey pixels represent the interior of the particles.The pixels have been assigned using a radial distance criterion.(B) Determination of A 0 by sequential scanning the y-pixels for each successive x-pixel.Only the pixels located at the outer y-edges, indicated by the white dot, contribute to A 0 .
The computing time was significantly shortened by sequential scanning the lattice as illustrated in Figure 2B (in 2D).For each successive pixel in the x-direction, the pixels in the y-direction were scanned.Only those pixels located at the outer y-edges of the excluded volume sphere contributed to A 0 .These pixels are labeled by a white dot in Figure 2B.Note that A 0 is underestimated since not all pixels constituting A 0 are taken into account (the dark grey pixels without a white dot in Fig. 2B).In addition, there is a fundamental conflict between the 3D character of a pixel and the 2D character of A 0 : the surface area of a cut through a unit volume does not immediately correspond to a unit surface area.To correct for these issues, we compared the exact surface area of a sphere to the surface area of an identical sphere as measured in a simulation using the lattice method.The ratio between the exact and measured surface areas provided a correction factor for A 0 .This correction factor can be evaluated for all sphere diameters.In this work, the hard sphere diameter equals 19 pixels corresponding to a correction factor for A 0 of 2.075.This correction factor becomes 2 for infinitesimally small lattice sites.
In the case of V 0 , we also compared the exact volume of a sphere to the volume of an identical sphere as determined by the lattice method in a simulation.Evaluating the ratio between the exact and measured volumes for a sphere diameter of 19 pixels, yields a correction factor of 1.047 for V 0 .Obviously, the determination of V 0 and A 0 becomes more accurate upon decreasing the size of a pixel with respect to the particle size.In this work, the number of pixels per particle is more than 2.5 • 10 4 .We checked that this finite number of pixels did not affect the accuracy of the measurements.The lattice method to measure V 0 and A 0 has been verified by comparing the results of the hard sphere Monte Carlo simulations, without further adjustment, to predictions that follow from the Carnahan-Starling equation of state [30].Here, N z is the number of particles between z ± Δz bin /2 with Δz bin the binsize; L x and L y are the x and y dimensions of the 2D cross-section images.The density profiles for the system at φ = 0.05, 0.19, 0.43 are shown in Figure 4.At all volume fractions, the density as a function of the height is fairly constant indicating that the system is isotropic in the plane of gravity.The structure of the hard sphere fluid at different densities can also be quantified by the (three dimensional) radial distribution function: g(r) = ρ −2 i j i δ(r i )δ(r j − r) .The indices i and j run over all particles.Figure 5A shows this correlation function of the system at all volume fractions.The contributions to g(r) at distances smaller than one particle diameter are due to errors in the particle tracking.Note that this effect becomes stronger upon diluting the system -hence -increasing particle movement.At lower volume fractions the noise in the correlation function increases for statistical reasons as well.The effect of the increasing volume fraction is clearly reflected by the enhanced structure in g(r).The gradual increase of the peak heights exactly follows the volume fraction sequence.The fact that even at the highest volume fraction the second peak shows no splitting, points towards the absence of local hexagonal order.Furthermore, the position of the first peak of g(r) being at the same distance for all volume fractions is a clear indication of the hard-sphere like interparticle potential.In addition, the radial distribution functions are in agreement with g(r)'s calculated using Percus-Yevick theory, which is known to be accurate for hard spheres [22] (see Fig. 5B).Thus all together, the confocal images, the density profiles and the radial distribution functions show that the obtained snapshots represent an equilibrium hard sphere fluid at the various volume fractions.

RESULTS AND DISCUSSION
From the hard sphere configuration the available volume V 0 and the accompanying surface area A 0 were measured as explained above.Figure 6 shows V 0 and A 0 as a function of the volume fraction measured from the experiment (dots), from the simulations (crosses) and calculated from the Carnahan-Starling equation of state (solid line) [30].First of all, the agreement between the simulations and Carnahan-Starling for V 0 as well as A 0 demonstrates the accuracy of the lattice approach, and thus eliminates the possibility that discrepancies are due to the lattice method.In other words, it verifies the observed agreement between the experimentally obtained V 0 and A 0 and those deduced from Carnahan-Starling.The available volume and the surface area show completely different trends.While V 0 decreases rapidly upon increasing volume fraction, A 0 increases up to a maximum at φ = 0.10 before it starts to decrease.This is probably correlated to a change in topography of the available volume [2].At low densities the excluded volume spheres are not overlapping and with every additional sphere, V 0 decreases while A 0 increases.At some density the excluded volume spheres start to overlap and eventually the available volume becomes disconnected [2,31,32].At  V 0 (black) and A 0 (grey) in number of pixels as a function of volume fraction.The solid lines follow from the Carnahan-Starling equation of state [30].Dots denote measurements and crosses represent Monte Carlo simulations.The inset shows V 0 and A 0 on a log scale for 0.23 ≤ φ ≤ 0.43.
this point, the surface area starts to decrease upon increasing density.Note that below φ = 0.10, V 0 is approximately one order of magnitude larger than A 0 .Upon further densification, both V 0 and A 0 decrease rapidly.Moreover, at φ ≈ 0.30 both quantities occupy less than 1% of all pixels.
Using Equations 1 and 2, the measured V 0 and A 0 directly yield the chemical potential and the pressure.Figure 7A and B respectively show, as a function of the volume fraction, the experimentally measured chemical potential and pressure (dots) compared to the Carnahan-Starling equation of state (solid lines) [30].The simulation data points are also shown (crosses).The observed agreement demonstrates that both the pressure and chemical potential can be measured from pictures alone.In fact, it is the first time that the chemical potential is directly measured in an experimental system.At higher volume fractions the chemical potential starts to deviate from the Carnahan-Starling equation of state.The pressure shows similar behavior: at low volume fractions hardly any deviation is observed, while upon increasing concentration the measured pressure is slightly below the Carnahan-Starling equation of state.
Ideally, many more snapshots than 11 must be analyzed to assure proper statistical averaging.Since the available volume and corresponding surface area become exceedingly small, they are very difficult to determine accurately at high densities.Slight deviations may have a relatively large effect on the extremely small V 0 and A 0 and thus also on the chemical potential and the pressure.In addition, the small polydispersity of the particles (∼ 4%) could induce deviations in V 0 and A 0 at high densities.In fact, the pressure of a (A) The chemical potential μ and (B) the pressure P as a function of volume fraction φ.Dots denote the measurements, crosses the Monte Carlo simulations and the full curves follow from the Carnahan-Starling equation of state [30].polydisperse fluid is always lower than the corresponding monodisperse system at the same density, since polydispersity induces a more efficient packing and thereby increases the available volume [33].
Combining the pressure and chemical potential yields the free energy density (Eq.13) of the system, which is presented in Figure 8.The experimental free energy density shows a reasonable agreement with the Carnahan-Starling equation of state.Note that the errors in the pressure and chemical potential add up to induce larger errors in the free energy density, especially at intermediate densities.At φ = 0.43 the experimental free energy density even shows a decrease.Both the measured pressure and The free energy density F/V as a function of volume fraction φ.Dots denote the measurements, crosses the Monte Carlo simulations and the full curves follow from the Carnahan-Starling equation of state [30].
chemical potential are too small, but the contribution from the chemical potential in Equation 13dominates leading to an overall decrease in the free energy density.Nevertheless, it is intriguing to see that the free energy can be measured directly from images alone.

CONCLUSIONS
We have taken three-dimensional snapshots of a colloidal hard sphere suspension throughout almost the whole fluid region using confocal microscopy.From these particle configurations the available volume to insert an additional sphere and the surface area of that volume were determined.Statistical geometry of hard spheres provides direct relations between these quantities and thermodynamics.Using this concept, we measured in a direct and non-interfering way -in principle -all thermodynamic properties, here demonstrated for the pressure, the chemical potential and the free energy density of a homogeneous hard sphere suspension by optical microscopy.Because real liquids are only a perturbation away from a hard-sphere fluid, it now becomes feasible to explore the relation between local geometrical quantities and thermodynamics beyond hard spheres as well [34][35][36][37][38].For example, the presented method can be applied to a suspension of bidisperse colloids [39,40].Another possible extension lies in the use of nonspherical hard colloids, for example, fluorescent rods [41,42], or polyhedral colloids [43] where exactly the same concepts can be applied.Constructing the available volume and its surface area then becomes more cumbersome, but in principle possible.

First
of all, a representative 2D confocal image of the system at φ = 0.25 is shown in Figure 3A.The distribution of particles is reminiscent of a homogeneous fluid.The 3D computer reconstruction of the system presented in Figure 3B confirms the isotropic structure in three dimensions.Although the mass density of the particles is almost matched by the mass density of the solvent, we have to rule out that gravity induces a gradient in the concentration within the imaged volume.Therefore, we calculated the number density profiles along the z-direction ρ(z) = N z /(L x L y Δz bin ).

Figure 5 (
Figure 5 (A) Radial distribution functions of the system at all volume fractions (as indicated).(B) The experimental (circles) and Percus-Yevick (solid line) radial distribution functions at φ = 0.43, 0.30, 0.19, 0.11, 0.05.The curves in A and B have been shifted vertically for clarity.