  The construction of analytic solutions for quasi-exactly solvable systems is
an interesting problem. We revisit a class of models for which the odd
solutions were largely missed previously in the literature: the anharmonic
oscillator, the singular anharmonic oscillator, the generalized quantum
isotonic oscillator, non-polynomially deformed oscillator, and the
Schr\"odinger system from the kink stability analysis of $\phi^6$-type field
theory. We present a systematic and unified treatment for the odd and even
sectors of these models. We find generic closed-form expressions for
constraints to the allowed model parameters for quasi-exact solvability, the
corresponding energies and wavefunctions. We also make progress in the analysis
of solutions to the Bethe ansatz equations in the spaces of model parameters
and provide insight into the curves/surfaces of the allowed parameters in the
parameter spaces. Most previous analyses in this aspect were on a case-by-case
basis and restricted to the first excited states. We present analysis of the
solutions (i.e. roots) of the Bethe ansatz equations for higher excited states
(up to levels $n$=30 or 50). The shapes of the root distributions change
drastically across different regions of model parameters, illustrating
phenomena analogous to phase transition in context of integrable models.
Furthermore, we also obtain the $sl(2)$ algebraization for the class of models
in their respective even and odd sectors in a unified way.
