2nd-order conformal superintegrable systems in $n$ dimensions are Laplace equations on a manifold with an added scalar potential and $2n - 1$ independent 2nd order conformal symmetry operators. They encode all the information about Helmholtz (eigenvalue) superintegrable systems in an efficient manner: there is a 1-1 correspondence between Laplace superintegrable systems and Stackel equivalence classes of Helmholtz superintegrable systems. In this paper we focus on superintegrable systems in two dimensions, $n = 2$, where there are 44 Helmholtz systems, corresponding to 12 Laplace systems. For each Laplace equation we determine the possible 2-variate polynomial subspaces that are invariant under the action of the Laplace operator, thus leading to families of polynomial eigenfunctions. We also study the behavior of the polynomial invariant subspaces under a Stackel transform. The principal new results are the details of the polynomial variables and the conditions on parameters of the potential corresponding to polynomial solutions. The hidden gl_3-algebraic structure is exhibited for the exact and quasi-exact systems. For physically meaningful solutions, the orthogonality properties and normalizability of the polynomials are presented as well. Finally, for all Helmholtz superintegrable solvable systems we give a unified construction of 1D and 2D quasi-exactly solvable potentials possessing polynomial solutions, and a construction of new 2D PT-symmetric potentials is established.
Download