Using the concept of conformable fractional derivative, we study the properties of fractional [Formula: see text]-dimensional Schrödinger equation for the potential [Formula: see text]. The extended Nikiforov–Uvarov method is generalized to the fractional domain and then employed to obtain the analytic exact energy eigenvalues and eigenfunctions and their dependence on the fractional order [Formula: see text] and the dimension [Formula: see text]. To test its applicability, we apply the method on heavy quarkonia systems, and reproduce their mass spectra and fractional radial probabilities at different values of [Formula: see text] and [Formula: see text]. Comparing the mass spectra with the experimental data, we discuss to what extent fractional models can account for some features in the description of heavy quarkonia at certain dimensional space.