How dense can one pack spheres of arbitrary size distribution?

We present the first systematic algorithm to estimate the maximum packing density of spheres when the grain sizes are drawn from an arbitrary size distribution. With an Apollonian filling rule, we implement our technique for disks in 2d and spheres in 3d. As expected, the densest packing is achieved with power-law size distributions. We also test the method on homogeneous and on empirical real distributions, and we propose a scheme to obtain experimentally accessible distributions of grain sizes with low porosity. Our method should be helpful in the development of ultra-strong ceramics and high performance concrete.

💡 Research Summary

The paper addresses the long‑standing geometric problem of determining the maximal packing density of spheres (or disks) when the particle sizes follow an arbitrary distribution, a question of direct relevance to the production of ultra‑strong ceramics and ultra‑high‑performance concrete (UHPC). Existing analytical results are limited to monodisperse packings or to a few idealized bimodal or multimodal distributions; no general method exists for a continuous, possibly complex, size distribution. The authors therefore develop a systematic algorithm that yields an upper bound on the achievable packing density for any given size distribution.

The core of the method is a hierarchical, Apollonian‑type filling scheme. First, the size distribution is discretized into bins, each bin containing particles of similar radius. Within each bin the particles are placed in the densest possible monodisperse arrangement: in two dimensions the hexagonal close‑packing (hcp) with density ρ_hcp≈0.9069, and in three dimensions a tetrahedral arrangement of four mutually tangent spheres with density ≈0.3633. Starting from the bin of the largest particles, the algorithm proceeds downward: the voids left by the hcp (or tetrahedral) arrangement are identified, and the required volume of smaller particles needed to fill those voids is calculated using the positive solution of the Soddy‑Gossett (Descartes) equation. If the immediate smaller bin does not contain enough material, the algorithm draws the necessary volume from still smaller bins, continuing until the voids are completely filled or the material is exhausted. After each step the effective volume V_eff (including both particles and voids) and the actual particle volume V are updated for every bin. The final net density is then ρ_net = Σ V / Σ V_eff.

A key theoretical insight is that when the particle‑volume distribution follows a power law V(r) ∝ r^{−α}, the optimal exponent α_opt that maximizes ρ_net is directly linked to the fractal dimension d_f of the Apollonian packing: α_opt = d – d_f, where d is the Euclidean dimension (d = 2 or 3). In 2D the Apollonian fractal dimension is d_f ≈ 1.306, giving α_opt ≈ 0.694; in 3D d_f ≈ 2.474, giving α_opt ≈ 0.526. Numerical experiments confirm that as the lower cutoff ε (minimum particle radius) is reduced, the optimal α approaches these theoretical values and the net density approaches unity. For example, with α = 0.71 and ε = 10^{−5} in 2D the algorithm yields ρ_net ≈ 0.9997, while in 3D with α = 0.55 and ε = 3×10^{−3} it yields ρ_net ≈ 0.9325, a porosity roughly one‑third of that of a monodisperse close‑packed sphere assembly.

The authors also test the algorithm on a realistic UHPC particle size distribution comprising four components (crushed quartz, cement, sand, silica fume). The naïve application of the method gives an upper bound ρ_net ≈ 0.8203, consistent with typical experimental values (~0.8). To improve upon this, they propose an optimization scheme in which each component is modeled as a Gaussian distribution characterized by mean size μ_i, standard deviation σ_i, and peak height H_i. By fixing the volume fractions (e.g., equal for two components) and scanning the (μ,σ) parameter space, the algorithm identifies combinations that minimize porosity. Adding further Gaussian components iteratively refines the mixture. Applying this procedure to the UHPC example raises the achievable density to ρ_opt ≈ 0.9501, a three‑fold reduction in porosity relative to the original mixture. The authors suggest that such optimized mixtures could be realized experimentally through controlled sieving, filtration, or other size‑selection techniques.

The paper concludes by emphasizing the generality of the approach: it provides (1) a practical computational tool for estimating the densest possible packing of arbitrarily distributed spheres, (2) a theoretical justification for why power‑law size distributions are optimal, and (3) a concrete, step‑by‑step recipe for engineering low‑porosity granular composites. Limitations include the assumption of spherical particles and the neglect of particle deformation, friction, or dynamic compaction effects. Future work is suggested on extending the method to non‑spherical shapes, mixed‑shape ensembles, and incorporating realistic processing dynamics. Overall, the study bridges a gap between abstract packing theory and the practical design of high‑performance granular materials.