  The study of models with extended Higgs sectors requires to minimize the
corresponding Higgs potentials, which is in general very difficult. Here, we
apply a recently developed method, called numerical polynomial homotopy
continuation (NPHC), which guarantees to find all the stationary points of the
Higgs potentials with polynomial-like nonlinearity. The detection of all
stationary points reveals the structure of the potential with maxima,
metastable minima, saddle points besides the global minimum. We apply the NPHC
method to the most general Higgs potential having two complex Higgs-boson
doublets and up to five real Higgs-boson singlets. Moreover the method is
applicable to even more involved potentials. Hence the NPHC method allows to go
far beyond the limits of the Gr\"obner basis approach.
