Invariant Classification and Limits of Maximally Superintegrable Systems in 3D

The invariant classification of superintegrable systems is reviewed and utilized to construct singular limits between the systems. It is shown, by construction, that all superintegrable systems on conformally flat, 3D complex Riemannian manifolds can be obtained from singular limits of a generic system on the sphere. By using the invariant classification, the limits are geometrically motivated in terms of transformations of roots of the classifying polynomials.


Introduction
The discovery and classification of superintegrable systems is currently undergoing significant activity, see [23] and citations within. Superintegrable systems are dynamical systems, for example classical or quantum Hamiltonian systems, that have more conserved quantities that degrees of freedom. One of the many active areas of research in superintegrable systems is in their classification. On this front, there has been significant progress on maximally superintegrable systems in 2 and 3-D complex Reimannian manifolds with constants of the motion that are at most second-order in the momentum, a necessary condition for multi-separability. In two dimensions, all such systems have been classified [13]. Somewhat surprisingly, all such systems can be obtained through singular limits of a "generic" system defined on the sphere. This fact was first recognized by Bôcher [1] in his search for metrics that admit the maximal number of second-order Killing tensors. Again in the 2D case, the limits between such systems were recently given explicitly in [20] where they were shown to induce contractions on the associated symmetry algebras as well as limits in the representations of such algebras. These representations are given in terms of classical hypergeometric orthogonal polynomials and the limits of the superintegrable systems correspond to the limits within the Askey scheme.
The aim of this paper is to give the corresponding limits of superintegrable systems in 3D. We will show explicitly that all such systems are limiting cases of a "generic" system defined on S 3 and give the coordinate transformation that generates each limit. In order to describe the limits we discuss the classification theory of 3D second-order superintegrable systems [2,3]. The classification makes use of classical invariant theory and associates to each superintegrable system a 7-dimensional representation of the conformal group. By studying the action of the conformal group on this representation, we can identify 10 equivalence classes each of which corresponds to a class of superintegrable systems. From the invariant classification, the necessary limits are abstractly motivated as opposed to the constructions of [20], which are given ad hoc.
The structure of the paper is as follows. We begin in Section 2 with an overview of the results in [2,3] giving the invariant classification of nondegenerate second-order superintegrable systems in 3D. Section 3 contains the results of the paper; namely we prove that all such superintegrable systems are limits of a "generic" system by demonstrating the appropriate limits. Section 4 contains a brief discussion of the results and future applications of this research.

Invariant classification of superintegrable systems
In this section, we describe the invariant classification of nondegenerate second-order superintegrable systems in 3D where the metric is assumed to be conformally flat. The dynamics are given by the following classical Hamiltonian on 6-dimensional phase space (p, x) The conformally flat assumption implies that the metric tensor takes the form g ij = λ(x)δ i j . Any function on the phase space A(p, x) will transform in time as Thus, a function on phase-space will be invariant under the Hamiltonian flow if and only if it Poisson commutes with H, ({A, H} = 0) and is referred to as a constant of the motion. Conversely, the Hamiltonian is also invariant under flow generated by the constant of the motion and so the function A is also often referred to as a symmetry of the Hamiltonian. A Hamiltonian on 2n-dimensional phase space is said to be Liouville integrable if there exist n mutually commuting constants of the motion that are functionally independent. We say that a Hamiltonian is superintegrable if there exist more than n constants of the motion and maximally superinintegrable if there exist 2n − 1 such symmetries. Not all of the symmetries can mutually commute and so we drop the requirement that the integrals commute among themselves. There is a close connection between superintegrability and noncommutative integrability, notably that maximally superintegrable systems are noncommutatively integrebrable. Indeed, the celebrated result of Nekhoroshev [24] concerning foliations of phase space for noncommutative integrable systems includes, as a special case, the result that bounded trajectories of superintegrable systems are closed.
We would like to classify Hamiltonians H on a conformally flat 3D manifold with 5 constants of motion that are second-order in the momentum, i.e. that there exist 5 (including the Hamiltonian) functionally independent constants of the motion such that {A, H} = 0.
If such symmetries exist then a direct computation of the Poisson commutators leads to the following set of differential equations: the Killing equations, which require the second-order part of A to be a Killing tensor on the manifold, and the Bertrand-Darboux (B-D) equations. The B-D equations are the compatibility conditions for the differential equations for W (x) and comprises a set of four linear second-order PDEs for the potential V (x), one for each integral that is not the Hamiltonian. Combining the equations into vector form gives where the subscript f ,k denotes the partial derivative with respect to x k and B and C are matrix functions of the given arguments of dimensions 12 × 5 and 12 × 3, respectively. Suppose in addition that the constants of the motion are functionally, linearly independent. This implies that the B-D equations are of rank 5 and it is possible to solve for second-order derivatives of V as in  Here V ,ee is not a second derivative, but a symmetry adapted variable. The use of V ,ee introduces redundancy in the equations (6 equations instead of 5), which is preserved for the sake of symmetry in the variables. If the compatibility conditions of (4) are assumed to hold identically then the potential depends on 5 parameters, the values of V ,ee , V ,1 , V ,2 , V ,3 and V at a generic point, the last point coming from a trivial additive parameter. Such a system is said to be non-degenrate. For a non-degenerate system, the potential is uniquely determined by the value of the matrix in (4) at a given generic point and V ,11 (x 0 ) are the parameters in the potential. As is clear from the form of the matrix, it depends on 10 functions which are assumed to satisfy the compatibility conditions for (4) identically. The compatibility conditions of (4) are first-order equations in the variables whose compatibility conditions are themselves identically satisfied. Therefore, the values of functions {Q, S, R} at a generic point x 0 are enough to uniquely determine the functions themselves along with a superintegrable system.

Example 1
The harmonic oscillator potential is a non-degenerate second-order superintegrable system. The B-D equations for the potential are  Taking the origin as a generic point, the coefficients are related to V as a = V ,ee (0), b = V ,1 (0), c = V ,2 (0), d = V ,3 (0), and e = V (0).

Equivalence classes
The next step in the classification is to identify the appropriate equivalence classes. Clearly, we would like to consider two potentials as equivalent if they are related by a change of parameters. This is accomplished in the previous section by making a canonical choice of parameters as the value of v(x 0 ). We would also like systems to be equivalent if they are related by translations in the position variables; such a translation corresponds to moving the regular point. We would also like systems to be equivalent if they are related by, possibly complex, rotations. This condition will be studied extensively in this section. Finally, we would like the equivalence classes to include systems that are related via the Stäckel transform [11,12,18,25]. That is, suppose a Hamiltonian can be expressed as H = H 0 + αU then the Stäckel transformed HamiltoniansĤ = 1 U H will also be superintegrable, perhaps on a different conformally flat manifold. Under this transformation, classical trajectories and quantum wave functions as well as the corresponding symmetry algebras are essentially preserved, up to a change of parameters. Thus, we would like to consider two such systems as being equivalent. An analogous classification of equivalence classes for superintegrable systems in 2D has already been performed [21].
In order to consider Stäckel equivalent systems, we focus our attention on conformally superintegrable systems [14]. For conformal superintegrable systems, the Hamiltonian is identically 0 on trajectories. This can be accomplished mathematically by including the trivial added parameter to the system. Thus, along a trajectory the Hamilton-Jacobi equation is expressed as H − E = 0 but with H − E being the Hamiltonian. A function on phase space will then be a constant of the motion whenever its Poisson bracket is given by It is a direct computation to verify (see e.g. Lemma 4.1.4 of [2]) that a conformal integral of the motion A of H will be a conformal integral of the scaled Hamiltonian U (x)H. Thus it is possible to transform any conformal superintegral Hamiltonian on conformally flat space to a conformally superintegrable Hamiltonian on Euclidean space. Furthermore, since the action of the Stäckel transform or couplingconstant metamorphosis corresponds to multiplying the Hamiltonian by a function, it is clear that two Hamiltonians equivalent under this action will be equivalent conformal superintegrable systems, from the perspective of their symmetry algebras.

Example 2
The Hamiltonian for the simple Harmonic oscillator p 2 − ω 2 x 2 corresponds to the following conformal Hamiltonian Dividing by x gives Clearly, trajectories that satisfy H = 0 will also satisfy H = 0, however the "energy" of the new Hamiltonian is given by the parameter ω 2 .
Conversely, any conformally superintegrable system that depends linearly on at least one arbitrary constant can be transformed, see e.g. Theorem 4.1.8 [2], into a superintegrable system by coupling constant metamorphosis. Therefore, classifying superintegrable systems on conformally flat manifolds up to Stäckel transform is equivalent to classifying conformally superintegrable systems on Euclidean space. The determining equations for such an A are similar to (4) except that the additive constant is no longer trivial but instead essential for specifying the Hamiltonian. The corresponding structure equations are then  The functions A jk 0 are quadratic functions in the {Q, S, R} given in [2].

Action of the Conformal Group
So far, we have determined that superintegrable systems are uniquely determined by the value of the 10 functions {Q, R, S} at a regular point. We are interested in determining equivalence classes of systems that are stable under the action of the conformal group so we look at the induced action of this group on these 10 functions. Again, the exact formulas and derivations can be found in [2,3] and here we state only the results necessary to understand the limits obtained in Section 3.
A conformal change of variables can be generated by translations, which act trivially on the functions {Q, R, S} and inversions in the spheres of varying radii, which decompose the functions into a 7-dimensional representation {Q, R} and a 3 dimensional representation carried by {S}. Significantly, we can use a group transformation to set the values of {S} to any desired value, described in the following theorem.
Theorem 1 (Re-statement of Thm. 4.2.18 [2]) Given a conformally superintegrable system with function values {Q 0 , R 0 , S 0 } at a regular point x 0 , there exists a conformal group motion that maps the function values to {Q 0 , R 0 , 0} at the transformed regular pointx 0 .
This theorem allows us restrict our attention to the the action of the conformal group on the 7-dimensional space {Q, R}. To understand the action, we consider the continuous generators of the conformal group, namely translations, scaling, rotations and möbius transformations. As discussed before, translations act trivially on our representations. Scaling the coordinates corresponds to a scaling of the functions. More interesting is the action of rotations. To represent this action, we form weight vectors from the functions so that the action of rotations are given by the standard raising and lowering operators for so (3). Using the isomorphism between so 3 (C) and sl 2 (C), we obtain a covariant representation of the action via the polynomial Here the action of SO 3 (C) is represented, via the isomorphism, as the standard action of SL 2 (C) as Möbius transformationsq Thus, the action of the conformal group on the functions {Q, R} can be represented by that action of GL 2 (C) = C × SL 2 (C) on the polynomials q(z) via scaling and (8). Furthermore, the action of Möbius transformations on the coordinates can be represented locally as rotation combined with a scaling and so, if we find invariants that are closed under translation, scaling and rotations they will be automatically closed under the entire action of the conformal group. Finally, we have the following theorem. Therefore, a key to understanding the classification of conformal superintegrable systems is to understand the classification of invariants of the roots of degree 6 polynomials. Furthermore, since we would like the equivalence classes to be stable under translation of the regular point, we require that any invariants obtained be closed under derivations. This significantly reduces the number of equivalence classes to 10, exactly the expected number of equivalent superintegrable systems. For more details on the classification, we refer the reader to the original papers [2,3].

Contractions
For each of the contractions, we need only record the following data. The initial regular point x 0 , the rotation angles t 1 , t 2 , t 3 , the scaling parameter c and the final regular point y 0 . We note that the rotations can be represented in SL 2 (C) by the following matrices: The scaling c can be encoded in the same manner by the matrix In the following sections, unless otherwise mentioned, the angles are assumed to be set to 0 and the scale factor c to be set to 1. The action of the matrices on the polynomials are given by as in (8), which transform the roots asâ The action of the SL 2 (C) matrices can also be interpreted in terms of their action on a stereographic projection of the complex plane onto the unit sphere using the formula (X, Y, Z) = 2x 1 + x 2 + y 2 , For real t 1 , the action of R 1 (t 1 ) is to rotate the roots around the X axis by an angle t 1 counter-clockwise (i.e. from Y to Z). Similarly, the action of R 2 (t 2 ) is to rotate the roots around the Y axis clockwise (i.e. from X to Z). Finally, the action of R 3 (t 3 ) is to rotate the roots around the Z axis clockwise (i.e from Y to X).
The action of these rotations on R 3 can be recovered as In total, the change of coordinates for the potential is given by Notice that moving the regular point corresponds to translating the coordinates. For the remainder of this section, we give the limits between the systems. The first few case involve three or fewer roots, or (as in the case of [3111b]) four roots with a fixed cross ratio. The limits between these cases can be achieved solely with the action of GL 2 (C), that is, without appealing to translation the regular point. For most of these case the limits can be achieved by fixing one root and collapsing the rest together at infinity (or another prescribed point in C * ). Several of the limits are elaborated to give a more complete understanding of the process.

[6] V
The covariant polynomial for the [0] equivalence class is simply 0, so to contract down to this equivalence class, we need only scale by a factor of ǫ −1 and move the regular point appropriately. For V A , the regular point is already x 0 = (0, 0, 0) so it is unaffected by the scaling. The contraction is then x 0 = (0, 0, 0), y 0 = (0, 0, 0).
We start with the potential The conformal change of coordinates y = ǫ −1 x introduces a factor of ǫ −2 to the second order terms in the Hamiltonian (due to the conformal scaling of the metric) and so we consider the rescaled potential V A = ǫ 2 V A . Adjusting the parameters to bê gives the potential which tends in the ǫ → 0 limit to the potential V O =â 2 (y 2 1 + y 2 2 + y 2 3 ) +by 1 +ĉy 2 +dy 3 +ê.

[51] V V II to [6] V A
The regular points for V V II and V A are the origins, x 0 = 0 and y 0 = 0. The change of variables is given by which transforms the polynomial q V II = −36iz → q A = −iz 6 .
Consider first the roots of the (projective) sextic q V II = −36iz. There are five at z = ∞ and one at z = 0.
In terms of a stereographic projection of the complex plane z = x + iy (with the plane of projection intersecting the equator) the action of R 2 (t 2 ) is to rotate the sphere around the y-axis (i.e. the imaginary axis which is directed into the plane) by the angle t 2 . So consider Under the action described above the covariant polynomial becomes Which has brought the fifth root down from infinity down to z = −1 and moved the root at zero up to z = 1. Re(z) Setting t 3 = i 2 ln (ǫ) gives the matrix On the polynomial this induces the change So finally, scaling by c = 9 2 ǫ −3 gives This limit can be visualised as all the points on the sphere (except ∞) being drawn down towards the origin.
Or equivalently the roots can remain stationary while the sphere descends into the plane of projection. From this point of view, the image of any point on the sphere (expect ∞) end up eventually being projected inside the circle where the plane and sphere intersect, which is shrinking to the origin in the limit.
Here, we would like to take the two roots of q OO , namely 0 and ∞, and move them both to 0. This is going to be essentially the same contraction as the previous one. The contraction is The regular points are x 0 = (0, 0, 1), y 0 = (0, 0, 0). The covariant polynomial is transformed as

[411] V III to [51] V V II
The covariant polynomial for V III is given by q III = −3(1 + 3z 2 ). In order to transform this polynomial into q V II = −36iz we need to take one of the finite roots and send it to zero and the other to infinity. The contraction that accomplishes this is The regular points for each system do not move and they are x 0 = (0, 1, 0), y 0 = (0, 0, 0).
The limit of the polynomial is then

[3111b] V V I to [51] V V II
The necessary contraction is given by x 0 = (0, 0, 2i), y 0 = (0, 0, 0). In this limit, all of the roots coalesce to infinity except for a single root at zero. The limit of the polynomial is then

[3111b] V V I to [33] V OO
For this contraction, we would like to move the roots of the polynomial q V I by sending the three nonzero roots to infinity and leaving the 0 root fixed. This can be accomplished by choosing an imaginary, singular value for t 3 which has the effect of scaling the variable z. The required contraction is then x 0 = (0, 0, 2), y 0 = (0, 0, 1).
The polynomials transform as

[3111a] V II to [411] V V
From this case onward the polynomials involve four or more roots and require additional information about the cross ratios of the roots in order to differentiate between equivalence classes. Translation of the regular point allows the cross ratios of these roots to change. Without translation of the regular point the only way singular limits could be taken is by collapsing all, save one, of the roots together. This means that, starting at a point which has represents root multiplicities [3111] the only possible limits without translation of the regular point are [51], [33] or their degenerations. In all of the examples that follow, translation of the regular point will be necessary to modify the cross ratios (i.e. move the roots in a way which the Möbius transformations cannot) and, in effect, rescue some of the roots in the singular limit.
To contract the covariant polynomial of V II , we observe that the polynomial evaluated at ( Specialising the first two coordinates to (x 1 , x 2 ) = (0, −i) gives and so by moving the regular point to ∞, we obtain the desired contraction. The regular points and scaling parameter is then chosen as This contraction gives the required limits of the polynomial, There are two distinct equivalence classes of systems whose covariant polynomials have the root structure [311]; these equivalence classes are distinguished by the characteristic that the cross ratio of the root be exp(iπ/3) (case [3111b]) or not (case [311a]) [2,3]. For the limit between these two systems, we begin with a polynomial with a triple root and three simple roots with a cross ratio of the roots away from 0, 1, ∞ or exp(±iπ/3). We would like to find a limit that will contract this polynomial to one with the same root structure but with a cross ratio of exp(iπ/3). We begin with the covariant polynomial for V II which has a triple root at infinity and move this root to 0 by rotating around the x 2 axis by π. We then move the x 1 coordinate of the regular point to infinity and perform a singular rotation about x 3 via This limit transforms the polynomials as (from ǫ = 1 to ǫ = 0) For this contraction, we take the equivalence class described by a polynomial with 6 simple roots satisfying the multi-ratio condition, a discriminating condition for the [111111] root structure [2,3], to one with two simple finite roots and the remaining at infinity. The coordinate limit is given by x 0 = (0, 1 The covariant polynomial contracts as

[111111c] V IV to [3111b] V V I
The covariant for V IV is given by An obvious easy way to achieve a triple root is to let the second coordinate x 2 tend to infinity. However this has the unintended side effect of collapsing the other three roots at infinity. The singular rotation R 3 (i ln(ǫ)) can be used to counteract this growth, provided the parameter ǫ correctly balances the growth of x 2 . This action only serves to accelerate the roots tending to zero, and hence gives us the limit we require. This contraction is accomplished by a singular rotation in around the x 3 -axis and a singular translation of the regular point. The contraction is The covariant polynomials contract as

[111111b] V I to [3111a] V II
Here we would like to send 3 of the finite roots to infinity. The covariant polynomial for V I is Rotations around the z axis have the effect of multiplying the variable z by e −2it3 and the entire polynomial by a factor of e 6it3 . Thus if we take a purely imaginary angle, the z 3 variable will be unchanged and it will be possible to send the coefficient of z 4 and z 6 to 0. In order to retain the lower order terms, it is necessary to move the regular point so as to offset the singular limit. The required contraction is The polynomials transform as 3.11.1 "Geometric inspired limit" This contraction can be carried out in a manner similar to the previous one. By moving the regular point correctly we can collapse three roots at z = 0 and three roots are z = ∞. Using the same singular rotation as before, the collapse of one of these triplets of roots can halted, and the other accelerated.
The covariant polynomial for V I is given by By choosing to work on the hypersurface x 2 = −ix 1 we ensure that q I has a double root at z = 0. The choice x 1 = η, x = −iη gives Sending η to infinity gives the limit collapses a single root to join the double root at zero, and sends the other three roots to infinity. Using a singular rotation around the x 3 axis the root tending to zero can be saved, and the roots tending infinity accelerated. The required contraction is therefore The polynomials transform as To obtain this contraction, we move the regular point and scale the potential. The required contraction is c = 1 ǫ , This transforms the covariant polynomials as The case [111111a] covers the most generic forms for the covariant polynomials. A particularly nice example of the polynomials in this class is However, this polynomial doesn't show the full regular point dependence of q SW (z) The first contraction is a scaling by ǫ −1 and translation of the regular point to the origin c = 1 ǫ , x 0 = (ǫ, ǫ, ǫ), y 0 = (1, 1, 1).
Following the regular point (ǫ, ǫ, ǫ) the covariant polynomial for the is given by Here ǫ must follow a path from 1 to 0 in the complex plane, avoiding points where ǫ 4 = 1 9 . The limit of this contraction transforms the covariant polynomial to

Conclusion
By using the invariant classification of conformally superintegrable systems on conformally flat 3D manifolds, we have been able to give explicitly the coordinate limits between the 10 equivalence classes. We have therefore proven that all such systems are limits of the generic system, V SW .
On of the most interesting avenues of this research is the connection between superintegrable systems, orthogonal polynomials (OPs) and special functions. From the beginning of the study of superintegrable systems, going back to the harmonic oscillator, Kepler Coulomb system and the discovery in '65 of multiseparable systems by Smorodinsky, Winternitz and collaborators [4,5,22], there has been a clear connection between superintegrable systems and the special functions of mathematical physics. In particular, the wave functions of superintegrable systems have been expressed in terms of OPs [16] and all superintegrable systems are conjectured to have this exactly-solvable nature [27]. More recently, this connection has been extended by studying the representation theory of the symmetry algebras associated with superintegable systems and interbasis expansion coefficients [6][7][8][9][10]17,19,26]. In [20], the Askey scheme of classical hypergeometric OPs in one variable as well as their limits were related to superintegrable systems in 2D and the contractions of their symmetry algebras. Work is currently ongoing to apply the limits obtained in this paper to find contractions of the associated quadratic algebras and their polynomial representations.
Finally, we would like to mention that as is evident from the classification theory, the superintegrable systems are uniquely determined by the functionally linearly independent set of second-order Killing tensors. Therefore we expect, as in the 2D case [15], that the quadratic algebras and their contractions will be completely determined by contractions of the underlying Lie algebra. However, this is still an open question and one which we anticipate will be aided by understanding the connection between the invariant classification and the quadratic algebras.

A Potentials
Here we list the potentials discussed above. Each of the limits above has been chosen so that the limit potential is in the form given below. Note that the parameters are chosen is for aesthetic reasons, to give the potentials in the most simple form. They are related by a linear transformation to the parameters {V (x 0 ), V 1 (x 0 ), V 2 (x 0 ), V 3 (x 0 ), V ee (x 0 )}, which of course depend on the regular point x 0 .
Note that V V is Stäckel equivalent to V III by choice of the Stäckel multiplier U = (x 1 + ix 2 ) −2 .
[111111a] I [111111b] x x r r r r r r r r r r r r r IV x x r r r r r r r r r r r r r r r r III/V