  This is the first paper of a series of our works on the self-similar
orbit-averaged Fokker-Planck (OAFP) equation for distribution function of stars
in dense isotropic star clusters. At the late stage of relaxation evolution of
the clusters, standard stellar dynamics predicts that the clusters evolve in a
self-similar fashion forming collapsing cores. However, the corresponding
mathematical model, the self-similar OAFP equation, has never been solved on
the whole energy domain $(-1< E < 0)$. The existing works based on kinds of
finite difference methods provide solutions only on the truncated domain $-1<
E<-0.2$. To broaden the range of the truncated domain, the present work resorts
to a (highly accurate and efficient) Gauss-Chebyshev pseudo-spectral method. We
provide a spectral solution, whose number of significant figures is four, on
the whole domain. Also, the solution can reduce to a semi-analytical form whose
degree of polynomials is only eighteen holding three significant figures. We
also provide the new eigenvalues; $c_{1}=9.0925\times10^{-4}$,
$c_{2}=1.1118\times10^{-4}$, $c_{3}=7.1975\times10^{-2}$ and
$c_{4}=3.303\times10^{-2}$, corresponding to the core collapse rate
$\xi=3.64\times10^{-3}$, scaled escape energy $\chi_\text{esc}=13.881$ and
power-law exponent $\alpha=2.2305$. Since the solution on the whole domain is
unstable against degree of Chebyshev polynomials, we also provide spectral
solutions on truncated domains ( $-1< E<E_\text{max}$, where
$-0.35<E_\text{max}<-0.03$) to explain how to handle the instability. By
reformulating the OAFP equation in several ways, we improve the accuracy of the
spectral solution and reproduce an existing self-similar solution, which infers
that existing solutions have only one significant figure at best.
