  A class {\cal R}_p of purely bosonic models is characterized having the
following properties in the Bargmann Hilbert space of analytic functions: (i)
wave function \psi(\epsilon,z)=\sum_{n=0}^\infty \phi_n(\epsilon) z^n is the
{\em generating function} for orthogonal polynomials \phi_n(\epsilon) of a
discrete energy variable \epsilon, (ii) any Hamiltonian \hat{H}_b\in {\cal R}_p
has nondegenerate purely point spectrum that corresponds to infinite discrete
support of measure d\nu(x) in the orthogonality relation of the polynomials
\phi_n, (iii) the support is determined exclusively by the points of
discontinuity of \nu(x), (iv) the spectrum of \hat{H}_b\in {\cal R}_p can be
numerically determined as fixed points of monotonic flows of the zeros of
orthogonal polynomials \phi_n(\upepsilon), (v) one can compute practically an
unlimited number of energy levels (e.g. 2^{53} in double precision). If a model
of {\cal R}_p is exactly solvable, its spectrum can only assume one of four
qualitatively different types. The results are applied to spin-boson quantum
models that are, at least partially, diagonalizable and have at least single
one-dimensional irreducible component in the spin subspace. Examples include
the Rabi model and its various generalizations.
