  We use a two-level simulation method to analyse the critical point associated
with demixing of binary hard sphere mixtures. The method exploits an accurate
coarse-grained model with two-body and three-body effective interactions. Using
this model within the two-level methodology allows computation of properties of
the full (fine-grained) mixture. The critical point is located by computing the
probability distribution for the number of large particles in the grand
canonical ensemble, and matching to the universal form for the $3d$ Ising
universality class. The results have a strong and unexpected dependence on the
size ratio between large and small particles, which is related to three-body
effective interactions, and the geometry of the underlying hard sphere
packings.
