  It is well known that using the conventional non-Hermitian but ${\cal
PT}-$symmetric Bose-Hubbard Hamiltonian with real spectrum one can realize the
Bose-Einstein condensation (BEC) process in an exceptional-point limit of order
$N$. Such an exactly solvable simulation of the BEC-type phase transition is,
unfortunately, incomplete because the standard version of the model only offers
an extreme form of the limit characterized by a minimal geometric multiplicity
$K=1$. In our paper we describe a rescaled and partitioned direct-sum
modification of the linear version of the Bose-Hubbard model which remains
exactly solvable while admitting any value of $K\geq 1$. It offers a complete
menu of benchmark models numbered by a specific combinatorial scheme. In this
manner, an exhaustive classification of the general BEC patterns with any
geometric multiplicity is obtained and realized in terms of an exactly solvable
generalized Bose-Hubbard model.
