Smearing techniques are widely used in first-principles calculations of metallic and magnetic materials where they improve the accuracy of Brillouin-zone sampling and lessen the impact of level-crossing instabilities. Smearing introduces a fictitious electronic temperature that smooths the discontinuities of the integrands; consequently, a corresponding fictitious entropic term arises, and needs to be considered in the total freeenergy functional. Advanced smearing techniques-such as Methfessel-Paxton and cold smearing-have been introduced to guarantee that the system's total free energy remains independent of the smearing temperature, at least, up to the second order. In doing so, they give rise to nonmonotonic occupation functions (and, for Methfessel-Paxton, nonpositive definite), which can result in the chemical potential not being uniquely defined. We explore this shortcoming in detail and introduce a numerical protocol utilizing Newton's minimization method that is able to identify the desired Fermi energy. We validate the method by calculating the Fermi energy of similar to 20 000 materials and comparing it with the results of standard bisection approaches. In passing, we also highlight how traditional approaches, based on Fermi-Dirac or Gaussian smearing, are actually equivalent for all practical purposes, provided the smearing width is appropriately renormalized by a factor of similar to 2.565.