  Certain models with purely repulsive pair interactions can form cluster
crystals with multiply-occupied lattice sites. Simulating these models'
equilibrium properties is, however, quite challenging. Here, we develop an
expanded isothermal-isobaric $[N]pT$ ensemble that surmounts this problem by
allowing both particle number and lattice spacing to fluctuate. We apply the
method with a Monte Carlo simulation scheme to solve the phase diagram of a
prototypical cluster-crystal former, the penetrable sphere model (PSM), and
compare the results with earlier theoretical predictions. At high temperatures
and densities, the equilibrium occupancy $n_{\mathrm{c}}^{\mathrm{eq}}$ of
face-centered cubic (FCC) crystal increases linearly. At low temperatures,
although $n_{\mathrm{c}}^{\mathrm{eq}}$ plateaus at integer values, the crystal
behavior changes continuously with density. The previously ambiguous crossover
around $T\sim0.1$ is resolved.
