Integrable and superintegrable quantum mechanical systems with position
dependent masses invariant with respect to one parametric Lie groups. 2.
Systems with dilatation and shift symmetries.

A. G. Nikitin 111E-mail<EMAIL_ADDRESS>Institute of Mathematics,
National Academy of Sciences of Ukraine,
3 Tereshchenkivs’ka Street, Kyiv-4, Ukraine, 01024, and

Universitá del Piemonte Orientale,

Dipartimento di Scienze e Innovazione Tecnologica,

viale T. Michel 11, 15121 Alessandria, Italy

3d quantum mechanical systems with position dependent masses (PDM) admitting
at least one second order integral of motion and symmetries with respect to
dilatation or shift transformations are classified. Twenty seven such systems
are specified and the completeness of the classification results is proved.

##  1 Introduction

Symmetry is one of the most fundamental concept of quantum mechanics, and it
is natural that just the symmetry lies in the focus of this research field. A
systematic search for Lie symmetries of Schrödinger equation started in papers
[1, 2, 3] and [4] where the maximal invariance groups of this equation with
arbitrary scalar potential have been discovered. For symmetries of this
equation with scalar and vector potentials and corrected results of classical
paper [4] see papers [5, 6]. Lie symmetries of Schrödinger equation with
matrix potentials and of quasirelativistic Schrödinger equation are classified
in [7, 8].

More general symmetries, namely, the second order symmetry operators for 2$d$
and 3$d$ Schrödinger equation have been classified in [9], [10] and [11],
[12]. Just such symmetries characterize integrable and superintegrable systems
[13].

One more interesting research field is formed by superintegrable systems with
spin. The first example of such systems presented in paper [14], the
systematic study of them started with papers [15, 16] where the systems with
spin-orbit interaction were classified. Superintegrable systems with Pauli
type interactions were studied in [17, 18] and [20], the examples of
relativistic systems can be found in [20] and [21].

One more research field closely related to the discussed above are symmetries
of Schrödinger equations with position dependent mass. This equation is
requested in many branches of modern theoretical physics, whose list can be
found, e.g., in [22, 24]. These symmetries are studied much less. Nevertheless
symmetries with respect to the continuous groups have been classified in
papers [24] -[26].

The situation with the higher symmetries of the PDM systems is more
complicated. The 2d classical systems with position dependent mass which admit
second order integrals of motion are perfectly classified [27, 28, 29, 30]. In
particular, the two dimensional second-order (maximally) superintegrable
systems for Euclidean 2-space had been classified also algebraic geometrically
[31], see also [32] for the modern proof of the fundamental Koenig theorem
[27]. On the other hand the more interesting for physics 3d systems and their
higher symmetries are not completely classified jet. More exactly, the
maximally superintegrable systems (i.e., admitting the maximal possible number
of integrals of motion) and (or) for the system whose integrals of motion are
supposed to satisfy some special conditions like the functionally linearly
dependence [33] are well studied [34]. The same is true for the so called
nondegenerate systems [35, 36]. A certain progress has been achieved in the
classification of the semidegenerate systems [37].

Surely, the maximally superintegrable and nondegenerate systems are both
important and interesting. In particular, they admit solutions in multi
separated coordinates [38, 39, 40]. On the other hand, there are no reasons to
ignore the PDM systems which admit second order integrals of motion but are
not necessary either superintegrable or nondegenerate. And just such systems
are the subject of our study.

The main stream in studying of superintegrable systems with PDM is the
investigation of classical Hamiltonian systems. Surely there exist the
analogous quantum mechanical systems which in principle can be obtained
starting with the classical ones by the second quantization procedure.
However, the mentioned procedure is not unique, and in general it is possible
to generate few inequivalent quantum systems which have the same classical
limit. In addition, a part of symmetries and integrals of motion can disappear
in the classical limit $h\to 0$ [41]. And it is the reason why we classify
superintegrable quantum systems directly.

In view of the complexity of the total classification problem of integrals of
motion for 3d quantum PDM systems it is reasonable to separate it to well
defined subproblems which can have their own values. The set of such
subproblems can be treated as optimal one if solving them step by step we can
obtain the complete classification.

We choose the optimal set of subproblems in the following way.

It was shown in [24] the PDM Schr¨odinger equation can admit one, two, three,
four, or six parametric Lie symmetry groups. In addition, there are also such
equations which have no Lie symmetry. In other words, there are six well
defined classes of such equations which admit n-parametric Lie groups with
$n=1,2,3,4,6$ or do not have any Lie symmetry. And it is a natural idea to
search for second order integrals of motion consequently for all these
classes.

The systems admitting six- or four-parametric Lie groups are not too
interesting since the related Hamiltonians cannot include non-trivial
potentials. That is why we started our research with the case of three-
parametric groups. The classification of the corresponding PDM systems
admitting second order integrals of motion was obtained in [42]. There were
specified 38 inequivalent PDM systems together with their integrals of motion.
The majority of them are new systems which are not not maximally
superintagrable.

Notice that the superintegrable 3d PDM systems invariant with respect to the
3d rotations have been classified a bit earlier in paper [43] where their
supersymmetric aspects were discussed also. For relativistic aspects of
superintegrability see [21] , [44].

The systems admitting two-parametric Lie groups and second order integrals of
motion have been classified in [45]. We again find a number of new systems in
addition to the known maximally superintegrable ones.

The natural next step is to classify the systems which admit one parametric
Lie groups and second order integrals of motion. As it is shown in [24], up to
equivalence these groups are reduced to dilatations, shifts along the fixed
coordinate axis, rotations around this axis and some specific combinations of
the mentioned transformations.

The integrable and superintegrable PDM admitting the symmetry w.r.t. the
mentioned rotations, i.e., possessing the cylindric symmetry had been
classified in paper [46]. The number of such systems appears to be rather
extended. Namely, the 68 inequivalent systems admitting second order integrals
of motion has been classified in [46]. The majority of these systems are
neither maximally integrable nor degenerate and so are new and cannot be
discovered using the approaches exploited in [28, 29, 30], [35, 36] and [37].

In the present paper we make the next step to the complete classification of
the integrable and superintegrable PDM systems admitting one parametric Lie
symmetry groups. Namely, we present the complete classification of two
subclasses of such systems. The first class include the systems possessing the
dilatation symmetry while the system belonging to the second class are
invariant w.r.t. the shifts along the third coordinate axis.

Notice that the second order integrals of motion which can be admitted by the
PDM systems invariant with respect to the two parametric Lie group including
both the mentioned shifts and dilatation transformations have been classified
in [45] and it is interesting to extend this result to the case of the systems
invariant w.r.t. its one dimensional subgroups. In addition, the second order
integrals of motion of such systems have the following unique property: they
are generated by the specific conformal Killing tensors whose traceless parts
are polynomials of order less than three, and this one more reason to present
both the mention subclasses in one paper.

To classify the mentioned systems we use the optimised algorithm of solution
of the related determining equations proposed in [46]. The number of the
inequivalent systems with the dilatation symmetry appears to be equal to nine,
including four integrable, three superintegrable and two maximally
superintegrable ones. Notice that a special restricted class of such systems
was studied earlier in papers [42] and [47].

We present also eighteen inequivalent systems which admit second order
integrals of motion and are invariant with respect to the shift
transformations. Among them there are seven integrable, seven superintegrable
and four maximally superintegrable ones. This means that the majority of the
found systems is new.

##  2 Main definitions

We will search for superintegrable stationary Schrödinger equations with
position dependent mass of the following generic form:

$\displaystyle\hat{H}\psi=E\psi,$ (1)

where

$\displaystyle
H=\frac{1}{4}(m^{\alpha}p_{a}m^{\beta}p_{a}m^{\gamma}+m^{\gamma}p_{a}m^{\beta}p_{a}m^{\alpha})+\hat{V}(x)$
(2)

where $p_{a}=-i\partial_{a}$, $m=m({\bf x})$ is a function of spatial
variables ${\bf x}=(x_{1},x_{2},x_{3})$, associated with the position
dependent mass, $\alpha,\beta$ and $\gamma$ are the so called ambiguity
parameters satisfying the condition $\alpha+\beta+\gamma=-1$, and and
summation from 1 to 3 is imposed over the repeating index $a$.

There are various physical speculations how to fix the ambiguity parameters in
the particular models based on equation (2) . However, the systems with
different values of these parameters are mathematically equivalent up to
redefinition of potentials $V({\bf x})$. To simplify the following
calculations we will fix them in the following manner:
$\beta=0,\gamma=\alpha=-\frac{1}{2}$ and denote $f=f({\bf x})=\frac{1}{2m({\bf
x})}$. As a result hamiltonian (2) is reduced to the following form

$\displaystyle H=f^{\frac{1}{2}}p_{a}p_{a}f^{\frac{1}{2}}+V({\bf x}).$ (3)

Alternatively, there are some historical reasons to consider the version
$\alpha=\gamma=0,\beta=-1$ which corresponds to the following form of
hamiltonian (2):

$\displaystyle H=p_{a}fp_{a}+\hat{V}.$ (4)

Operators (4) and (3) are equal one to another provided potentials $V=V({\bf
x})$ and $\hat{V}=\hat{V}({\bf x})$ satisfy the following condition:

$\displaystyle V=\hat{V}+V^{k}$ (5)

where

$\displaystyle
V^{(k)}=\frac{1}{4f}((\partial_{1}f)^{2}+(\partial_{2}f)^{2}+(\partial_{3}f)^{2})-\frac{1}{2}\Delta
f$ (6)

and $\Delta$ is the Laplace operator.

Formula (6) represents an example of kinematical potentials which can be
reduced to zero by the rearranging the ambiguity parameters.

In paper [24] all equations (1) with Hamiltonians (4) admitting at least one
parametric Lie symmetry group were classified. It was shown that there are six
inequivalent symmetries which can be accepted by the PDM Schrödinger
equations. They include rotation around the third coordinate axis, shift along
this axis and dilatation groups. In addition, we can fix three combined
symmetries which are superpositions of the mentioned ones.

In paper [46] we find all inequivalent PDM systems invariant with respect to
the mentioned rotations, which admit second order integrals of motion, i.e.,
differential operators of the following generic form

$Q=\partial_{a}\mu^{ab}\partial_{b}+\eta$ (7)

commuting with Hamiltonian (3). The number of such systems appears to be
rather extended. Namely, we specified 68 versions of them defined up to
arbitrary functions and arbitrary parameters.

In the present paper we search for the integrable PDM systems which are
invariant with respect to the dilatation and shift transformations. The
generic form of the related inverse masses $f$ and potentials $V$ can be
represented by the following formulae [24]:

$\displaystyle f=\frac{r^{2}}{F(\varphi,\theta)},\quad V=V(\varphi,\theta)$
(8)

for the case of dilatation symmetry and

$\displaystyle f=F(x_{1},x_{2}),\quad V=V(x_{1},x_{2})$ (9)

for symmetries with respect to the shifts along the third coordinate axis,
where $F(.)$ and $V(.)$ are arbitrary functions whose arguments are fixed in
the brackets,

$\displaystyle
r^{2}=x^{2}_{1}+x_{2}^{2}+x_{3}^{2},\quad\varphi=\arctan\left(\frac{x_{2}}{x_{1}}\right),\
\ \theta=\arctan\left(\frac{\tilde{r}}{x_{3}}\right).$

By definition, operators $Q$ should commute with $H$:

$[H,Q]\equiv HQ-QH=0.$ (10)

Our task is to find all inequivalent PDM Hamiltonians with specific arbitrary
elements $f$ and $V$ whose generic form is fixed in (8) and (9) which admit at
least one integral of motion (7) satisfying (10). As it was noted in [43] it
is reasonable to use representation (2) for the hamiltonian since just this
representation leads to the most simple form of the determining equations for
symmetries. And this is why we will use it in the following calculations.

##  3 Determining equations

Evaluating the commutator in (10) and equating to zero the coefficients for
the linearly independent differential operators
$\partial_{a}\partial_{b}\partial_{c}$ and $\partial_{a}$ we come to the
following determining equations for arbitrary elements $f,V$ of the
Hamiltonian and functions $\mu^{ab}$ , $\eta$ defining integrals of motion (7)
[45]:

$\displaystyle
5\left(\mu^{ab}_{c}+\mu^{ac}_{b}+\mu^{bc}_{a}\right)=\delta^{ab}\left(\mu^{nn}_{c}+2\mu^{cn}_{n}\right)+\delta^{bc}\left(\mu^{nn}_{a}+2\mu^{an}_{n}\right)+\delta^{ac}\left(\mu^{nn}_{b}+2\mu^{bn}_{n}\right),$
(11) $\displaystyle\left(\mu^{nn}_{a}+2\mu^{na}_{n}\right)f-5\mu^{an}f_{n}=0,$
(12) $\displaystyle\mu^{ab}V_{b}-f\eta_{a}+F^{a}=0$ (13)

where $\delta^{bc}$ is the Kronecker delta, $f_{n}=\frac{\partial f}{\partial
x_{n}},\ \mu^{an}_{n}=\frac{\partial\mu^{an}}{\partial x_{n}}$, etc., and
summation is imposed over the repeating indices $n$ over the values $n=1,2,3$.
Moreover,

$\displaystyle F^{a}=(f\mu^{ak}_{km}-\mu^{mn}f_{na})_{m}$ (14)

for the Hamiltonian of form (4) and

$\displaystyle F^{a}=(f\mu^{ak}_{km}-\mu^{mn}f_{na})_{m}+\mu^{ab}V^{(k)}_{b}$
(15)

for the Hamiltonian of form (5).

The term (14) in general is not trivial. However, it generates a kinematical
part of potential. On the other hand the term (14) is reduced to zero provided
equations (11) and (12) are satisfied. And this is why it is reasonable to use
just representation (5) for the Hamiltinian, since in this case equation (13)
is reduced to the maximally simple form, i.e.,

$\displaystyle\mu^{ab}V_{b}-f\eta_{a}=0.$ (16)

Just representation (5) for the PDM Hamiltonian was used in [43] where
rotationally invariant systems were described. The results presented in [42,
45, 46] are valid for the representation (5) also, but not for representation
(4). We present the requested comments in the preprint versions of the
mentioned papers.

Equation (11) defines the conformal Killing tensor. Its particular solution is
$\mu^{ab}=\mu_{0}^{ab}$ where

$\displaystyle\mu_{0}^{ab}=\delta^{ab}g({\bf x})$ (17)

with arbitrary function $g({\bf x})$. In addition the conformal Killing
tensors include fourth order polynomials in $\bf x$.

The determining equations (12) represent the coupled system of three nonlinear
partial differential equation equations for two unknowns $g({\bf x})$ and
$f({\bf x})$. However this system can be linearizing by introduction the new
dependent variables

$\displaystyle M=\frac{1}{f},\ N=\frac{g}{f}$ (18)

which reduces (12) to the following form:

$\displaystyle\left(\mu^{nn}_{a}+2\mu^{na}_{n}\right)M+5(\mu^{an}M_{n}+N_{a})=0.$
(19)

Equation (16) in its turn can be effectively linearised by introducing the
following new dependent variables $\tilde{M}$ and $R$:

$\displaystyle\tilde{M}=MV,\ R=\frac{N\tilde{M}}{M}-\eta$ (20)

which reduce (16) to the following equation:

$\displaystyle\left(\mu^{nn}_{a}+2\mu^{na}_{n}\right)\tilde{M}+5(\mu^{an}\tilde{M}_{n}+R_{a})=0$
(21)

which is analogous to (19).

Notice that the related inverse mass potential can be expressed via $M$ and
$\tilde{M}$ in the following form:

$\displaystyle f=\frac{1}{M},\ V=\frac{\tilde{M}}{M}$ (22)

while the corresponding functions $g$ and $\eta$ are expressed via $M$ and
$\tilde{M}$ in the following manner:

$\displaystyle g=\frac{N}{M},\ \eta=-\frac{N\tilde{M}}{M}-R.$ (23)

Just the linearised determining equation (18) will be used in the following to
solve our classification problem.

##  4 PDM systems invariant with respect to dilatation transformations

In this section we solve the first part of our classification problem and
consider the systems which are invariant with respect to the one parametric
Lie group of dilatation transformations. The corresponding determining
equations can be effectively decoupled and solved.

###  4.1 Equivalence relations

Any classification problem is solved up to equivalence relations, and we need
a clear definition of them in our particular case.

Non degenerated changes of dependent and independent variables of equations
(1), (2) are equivalence transformations provided they keep their generic form
up to the explicit expressions of arbitrary elements $f$ and $V$. They have a
structure of a continuous group which is extended by some discrete elements.

In accordance with [24] the maximal continuous equivalence group of equation
(1) is the group of conformal transformations of the 3d Euclidean space which
we denote as C(3). The corresponding Lie algebra c(3) is a linear span of the
following differential operators:

$\displaystyle\begin{split}&P^{a}=p^{a}=-i\frac{\partial}{\partial
x_{a}},\quad L^{a}=\varepsilon^{abc}x^{b}p^{c},\\\
&D=x_{n}p^{n}-\frac{3{\mathrm{i}}}{2},\quad
K^{a}=r^{2}p^{a}-2x^{a}D\end{split}$ (24)

where $r^{2}=x_{1}^{2}+x_{2}^{2}+x_{3}^{2}$. Operators $P^{a},$ $L^{a},$ $D$
and $K^{a}$ generate shifts, rotations, dilatations and pure conformal
transformations respectively.

In addition, equation (1) is form invariant with respect to the reflections of
spatial variables and exactly invariant w.r.t. the following discrete
transformations:

$\displaystyle x_{a}\to\tilde{x}_{a}=\frac{x_{a}}{r^{2}},\quad\psi({\bf
x})\to\tilde{x}^{3}\psi(\tilde{\bf x})$ (25)

where
$\tilde{x}=\sqrt{\tilde{x}_{1}^{2}+\tilde{x}_{2}^{2}+\tilde{x}_{3}^{2}}.$

The presented speculations are valid for an abstract system (1). However in
our case the system should be invariant w.r.t. the conformal transformations,
and to keep this property we have to reduce the equivalence group to its
subgroup which does not affect the generator of these transformation, i.e.,
$D$. In other words, the equivalence algebra whose base is presented in (24)
is reduced to its subalgebra whose elements commute with $D$. And this
subalgebra is nothing but so(3) $\oplus e_{1}$ spanned on $L_{1},L_{2},L_{3}$
and $D$. In addition, the mentioned discrete transformations are valid too.

###  4.2 Decoupling of the determining equations

The autonomous subsystem (11) defines the conformal Killing tensor which is a
linear combination of the term presented in (17) and the following terms (see,
e.g., [48]) :

$\displaystyle\mu^{ab}_{1}=\lambda_{1}^{ab},$ (26)
$\displaystyle\mu^{ab}_{2}=\lambda_{2}^{a}x_{b}+\lambda_{2}^{b}x_{a}-2\delta^{ab}\lambda_{2}^{c}x_{c},$
(27)
$\displaystyle\mu^{ab}_{3}=(\varepsilon^{acd}\lambda_{3}^{cb}+\varepsilon^{bcd}\lambda_{3}^{ca})x_{d},$
(28)
$\displaystyle\mu^{ab}_{4}=(x_{a}\varepsilon^{bcd}+x_{b}\varepsilon^{acd})x_{c}\lambda^{d}_{4},$
(29)
$\displaystyle\mu^{ab}_{5}=\delta^{ab}r^{2}+k(x_{a}x_{b}-\delta^{ab}r^{2}),$
(30)
$\displaystyle\mu^{ab}_{6}=\lambda_{6}^{ab}r^{2}-(x_{a}\lambda_{6}^{bc}+x_{b}\lambda_{6}^{ac})x_{c}-\delta^{ab}\lambda_{6}^{cd}x_{c}x_{d},$
(31)
$\displaystyle\mu^{ab}_{7}=(x_{a}\lambda_{7}^{b}+x_{b}\lambda_{7}^{a})r^{2}-4x_{a}x_{b}\lambda_{7}^{c}x_{c}+\delta^{ab}\lambda_{7}^{c}x_{c}r^{2},$
(32)
$\displaystyle\mu^{ab}_{8}=2(x_{a}\varepsilon^{bcd}+x_{b}\varepsilon^{acd})\lambda_{8}^{dn}x_{c}x_{n}-(\varepsilon^{ack}\lambda_{8}^{bk}+\varepsilon^{bck}\lambda_{8}^{ak})x_{c}r^{2},$
(33)
$\displaystyle\mu^{ab}_{9}=\lambda_{9}^{ab}r^{4}-2(x_{a}\lambda_{9}^{bc}+x_{b}\lambda_{9}^{ac})x_{c}r^{2}+(4x_{a}x_{b}+\delta^{ab}r^{2})\lambda_{9}^{cd}x_{c}x_{d}+\delta^{ab}\lambda_{9}^{cd}x_{c}x_{d}r^{2}$
(34)

where $r=\sqrt{x_{1}^{2}+x_{2}^{2}+x_{3}^{2}}$, $\lambda_{n}^{a}$ and
$\lambda_{n}^{ab}$ are arbitrary parameters, satisfying the conditions
$\lambda_{n}^{ab}=\lambda_{n}^{ba},\lambda_{n}^{bb}=0$. Thus we can represent
the generic integral of motion (7) in the following form:

$\displaystyle Q=Q_{(0)}+\sum Q_{(n)}+\eta$ (35)

where

$\displaystyle Q_{(0)}=P_{a}g({\bf x})P_{a},\ Q_{(n)}=P_{a}\mu_{m}^{ab}P_{b},\
m=1,2,..9$ (36)

Formulae (17) and (26)-(34) include an arbitrary function and 35 arbitrary
parameters, so the system of determining equations (19) is rather complicated.
Happily the considered PDM systems by definition should be invariant w.r.t.
the scaling transformations. It means that the related equations (19) cannot
include linear combinations of all polynomials listed in (26)-(34) but only
homogeneous ones. Indeed, integral of motion (7) including inhomogeneous terms
by construction is not invariant w.r.t. the mentioned transformations and so
can be transformed to linear combination of integrals whose coefficients are
homogeneous polynomials with arbitrary multipliers. In other words, the
determining equations (19) are reduced to the five decoupled subsystems
corresponding to the Killing tensors which are $n$-order homogeneous
polynomials with $n=0,1,2,3,4$. Moreover, since our Hamiltonians (3) with
arbitrary elements (8) are invariant with respect to the inverse
transformation (25) we can restrict ourselves to the polynomials of order
$n<3$, since symmetries with $n$=3 and $n$=4 appears to be equivalent to ones
with $n=1$ and $n=0$ correspondingly. The corresponding Killing tensors are
$\mu_{0}^{ab}+\mu_{1}^{ab},\mu_{0}^{ab}+\mu_{2}^{ab}+\mu_{3}^{ab}$ and
$\mu_{0}^{ab}+\mu_{5}^{ab}+\mu_{6}^{ab}$ which we represent in (17) and the
following formulae:

$\displaystyle\mu^{ab}=\lambda^{ab}+\delta^{ab}g({\bf r}),$ (37)
$\displaystyle\mu^{ab}=\lambda^{a}x_{b}+\lambda^{b}x_{a}-(2-\kappa)\delta^{ab}\lambda^{c}x_{c}+(\varepsilon^{acd}\lambda^{cb}+\varepsilon^{bcd}\lambda^{ca})x_{d}+\delta^{ab}g({\bf
r}),$ (38)
$\displaystyle\mu^{ab}=\mu_{(1)}^{ab}+\mu_{(2)}^{ab}+\delta^{ab}(g({\bf r})$
(39)

where

$\displaystyle\begin{split}&\mu_{(1)}^{ab}=(x_{a}\varepsilon^{bcd}+x_{b}\varepsilon^{acd})\lambda^{d}x_{c},\\\
&\mu_{(2)}^{ab}=\lambda^{ab}x^{2}-(x_{a}\lambda^{bc}+x_{b}\lambda^{ac})x_{c}+\delta^{ab}(g({\bf
r})-2\lambda^{cd}x_{c}x_{d}).\end{split}$ (40)

We omit the subindices present in the related formulae (26)-(32) and ignore
the tensor $\mu_{4}^{ab}$ which corresponds to the squared dilatation operator
being the transparent symmetry.

Notice that functions $g({\bf x})$ should be homogeneous functions of
independent variables. Moreover their homogeneity degree should be zero in
case (37), one in case (38) and two in case (39).

To classify the scale invariant PDM systems admitting second order integrals
of motion it is sufficient to solve the determining equations (18) where
$\mu^{ab}$ are special Killing tensors fixed in (38), (39). In other words,
our classification problem is decoupled to three subproblems corresponding to
symmetries whose differential terms are independent on $x_{a}$ or are linear
and quadratic in these variables.

###  4.3 Symmetries with differential terms $Q_{(n)}$ (36) independent on
spatial variables

Consider step by step all equations (18) including Killing tensors (37), (38),
and (39).

In the case presented in (37) operator (7) and equations (18) take the
following form

$\displaystyle Q=P_{a}(\lambda^{ab}+\delta^{ab}g({\bf x}))P_{b}+\eta$ (41)

and

$\displaystyle\lambda^{ab}M_{b}+N_{a}=0.$ (42)

Let us note that in fact we have to more equations additional to system (42).
By definition function $M=\frac{1}{f}$ satisfies the condition (refer to (8))

$\displaystyle x_{a}M_{a}=-2M.$ (43)

On the other hand, since function $g=g({\bf x})$ in (36) in fact is a function
of the Euler angles otherwise the dilatation transformation will change such
integral of motion. It follows from the above that function $N=gM$ satisfies
the same additional condition as $M$, i.e., $x_{a}N_{a}=-2N$.

Multiplying (42) by $x_{a}$, summing up with respect to the repeating index
$a$ and using (43) we can express the unknown function $N$ via the bilinear
combination variables $x_{a}$ and derivatives $M_{b}$:

$\displaystyle N=\lambda^{ab}x_{a}M_{b}$ (44)

Up to rotation transformations which belong to the equivalence group symmetric
tensor $\lambda^{ab}$ can be reduced to the diagonal form. Moreover, we can
restrict ourselves to the following nonzero entries

$\displaystyle\lambda^{11}=c_{1},\ \lambda^{33}=c_{2}$ (45)

since it is possible to include nonzero $\lambda^{22}$ to $g({\bf x})$ by
adding the term $P_{a}\lambda^{22}P_{a}$ to $Q$.

The corresponding equation (42) is qualitatively different in two cases:
$c_{1}c_{2}=0$ and $c_{1}c_{2}\neq 0$. In the first case we have only one
nonzero coefficient, say, $c_{3}$. Setting $c_{3}=1$ we come to the following
version of equation (42):

$\displaystyle M_{3}+N_{3}=0,\ N_{1}=N_{2}=0.$ (46)

The system of equations (43), (46) is solved by the following functions

$\displaystyle M=\frac{F(\varphi)}{\tilde{r}^{2}}+\frac{c_{1}}{x_{3}^{2}},\
N=-\frac{c_{1}}{x_{3}^{2}}$ (47)

where $F(\varphi)$ is an arbitrary function and $c_{1}$ is the integration
constant.

Solutions of equations (21) are analogous to (47) but include another
arbitrary function and integration constant. Substituting these solutions into
formulae (22) we come to the inverse mass and potential presented in Item 1 of
Table 2 where two integrals of motion are indicated, namely,

$\displaystyle Q=P_{3}^{2}-\left(\frac{c_{1}}{x_{3}^{2}}\cdot
H\right)+\frac{c_{2}}{x_{3}^{2}}$ (48)

and

$\displaystyle\tilde{Q}=K_{3}^{2}-\left(\frac{c_{1}r^{4}}{x_{3}^{2}}\cdot
H\right)+\frac{c_{2}r^{4}}{x_{3}^{2}}$ (49)

where we use the notation $(\Phi\cdot H)=P_{a}\Phi fP_{a}+\Phi V$. These
integrals can be transformed one into another via the symmetry transformation
(25).

In the case when both coefficients $c_{1}$ and $c_{2}$ are nontrivial equation
(42) is reduced to the following system

$\displaystyle c_{1}M_{1}+N_{1}=0,\ c_{2}M_{3}+N_{3}=0,\ N_{2}=0$ (50)

which is solved by the following functions

$\displaystyle
M=\frac{c_{3}}{x_{1}^{2}}+\frac{c_{4}}{x_{2}^{2}}+\frac{c_{5}}{x_{3}^{2}},\
N=\frac{c_{3}c_{1}}{x_{1}^{2}}+\frac{c_{2}c_{4}}{x_{3}^{2}}.$ (51)

Comparing (51) with (47) it is not difficult to understand that in fact we
have two symmetries corresponding to $c_{1}=0,c_{3}\neq 0$ and $c_{1}\neq
0,c_{3}=0$ and deal with their linear combination. Moreover, since all
variables $x_{1},x_{2}$ and $x_{3}$ are involved into $M$ presented in (51) in
the completely symmetric form we can predict one more symmetry whose
differential part is equal to $P_{2}^{2}$. Just these results are presented in
Item 5 of Table 2.

###  4.4 Symmetries with differential terms $Q_{(n)}$ (35) linear in spatial
variables

Thus we find the generic form of PDM Hamiltonians admitting integrals of
motion whose differential terms do not depend on spatial variables. The next
step is to classify the systems whose integrals of motion are generated by the
Killing tensors (38) which are linear in these variables.

Substituting (38) into (19) we come to the following equation

$\displaystyle\mu^{ab}M_{b}+N_{a}=0$ (52)

since the first term in (19) appears to be identical zero.

By definition function $M$ and $N$ satisfy condition (43) and the following
one:

$\displaystyle x_{a}N_{a}=-N$ (53)

which can be used to solve equation (52) for $N$. Indeed, multiplying the
latter equation by $x_{a}$, summing up with respect to the repeating index $a$
and using (43), (53) we obtain the following expression for $N$:

$\displaystyle N=r^{2}\lambda^{b}M_{b}+2\lambda^{b}x_{b}M$ (54)

or, using notations (8),

$\displaystyle N=\lambda^{b}F_{b}.$ (55)

Equations (52), (55) look rather gentle but in fact it is a rather complicated
system including eight arbitrary parameters. Fortunately the number of these
parameters can be reduced to five using the equivalence transformations which
include rotations and dilatations. To do it let us mention that symmetry
operator generated by the Killing tensor (38) can be represented as a bilinear
combination of first order differential operators (24) in the following
manner:

$\displaystyle
Q=\lambda^{ab}\\{P_{a},L_{b}\\}+\lambda^{a}\\{P_{a},D\\}+\delta^{ab}g({\bf
x}).$ (56)

Thus up to rotation transformations and the freedom in definition of their
trace we can restrict ourselves to the following nonzero entries of tensor
coefficients $\lambda^{ab}$:

$\displaystyle\lambda^{12}=\lambda^{21}=\mu,\ \lambda^{33}=\nu$ (57)

while the remaining entries can be reduced to zero.

The next reduction of the number of arbitrary parameters can be obtained using
the following speculations. System (52) includes three linear equations for
function $F$ dependent on two variables, and so to be consistent it should be
degenerated. In other words, the determinant of the matrix whose entries are
$\tilde{\mu}^{ab}$ should be equal to zero. This condition is satisfied in the
following cases:

$\displaystyle\mu=0,\ \nu=0,\ \kappa=1,$ (58) $\displaystyle\lambda_{1}=0,\
\lambda_{2}=0,\ \mu^{2}+\lambda_{3}^{2}-\nu^{2}=0,\ \kappa=0$ (59)
$\displaystyle\lambda_{2}=\pm\lambda_{1},\ \lambda_{3}=0,\
\nu^{2}\mp\nu\mu+\lambda_{1}^{2}=0,\ \kappa=0$ (60)
$\displaystyle\lambda_{1}=0,\ \lambda_{2}=0,\ \nu=0,\ \kappa=1,$ (61)
$\displaystyle\lambda_{1}=\pm\lambda_{2},\ \mu=\mp\nu,\
\lambda_{3}=0,\kappa=1$ (62)

while the parameters which are not included into (58)-(62) are arbitrary. In
other words effectively we have four versions of system (54) any of which
include as maximum two arbitrary parameters since up to rotation
transformations in the case (58) we can restrict ourselves to the only nonzero
parameter $\lambda^{3}$.

Let us consider versions (58)-(62) step by step. In the case (58) the only
nonzero parameter is $\lambda^{3}$ and equation (52) is reduced to the
following system:

$\displaystyle x_{a}M_{3}+N_{a}=0,a=1,2,3.$ (63)

In addition, unknowns $M$ and $N$ should satisfy the following conditions:

$\displaystyle x_{a}M_{a}+2M=0,$ (64) $\displaystyle x_{a}N_{a}+2N=0.$ (65)

System (63), (64) is easy solvable, its solutions are presented in Item 2 of
Table 2, but in terms of functions (18).

Considering version (59) we find that the related equations (52) and (64) have
nontrivial solutions only for the cases $\lambda_{3}=0,\mu=\pm\nu$ and
$\lambda_{3}=\pm\nu,\mu=0$.

In the first case these equations are solved by the following functions:

$\displaystyle M=\frac{c_{1}}{(x_{1}\pm x_{2})^{2}+2x_{3}^{2}},N=0.$ (66)

The corresponding PDM Schrödinger equations admit two parametric Lie symmetry
group including dilatations and simultaneous shifts of variables $x_{1}$ and
$x_{2}$. Such systems are completely classified in [45].

In the second case we obtain the solution

$\displaystyle M=\frac{c_{1}}{x_{1}^{2}+x_{2}^{2}},\ N=0$ (67)

while for the case (69) we again have solutions (66) which also admits two Lie
symmetries, i.e., dilatations and shifts along the third coordinate axis.

Version (60) generates a solvable system (52) and (64) in two cases:

$\displaystyle\lambda_{1}=0,\ \lambda_{2}=0,\ \lambda_{3}=0,\ \nu=0,\
\kappa=0,$ (68) $\displaystyle\lambda_{1}=0,\ \lambda_{2}=0,\ \lambda_{3}=0,\
\mu=\pm\nu,\ \kappa=0$ (69)

In the case (68) the related system (52) and (64) is solved by functions (67)
while for the case (69) we again have solutions (66). In other words we again
have the systems with two parametric Lie symmetry groups.

The version (60) generates solutions (67). However, for the special cases
$\lambda_{3}=\mu$ and $\lambda_{3}=0$ we obtain nice solutions presented in
Items 3 of Table 1 and 2 of Table 2.

Table 1. Inverse masses, potentials and integrals of motion for integrable
scale invariant systems.

No | $f$ | $V$ | Integrals of motion
---|---|---|---
1 | $\frac{r^{2}\tilde{r}^{2}}{r^{2}F(\theta)+c_{1}\tilde{r}^{2}\varphi}$ | $\frac{r^{2}G(\theta)+c_{2}\tilde{r}^{2}\varphi}{r^{2}F(\theta)+c_{1}\tilde{r}^{2}\varphi}$ | $\begin{array}[]{c}\\{L_{3},D\\}-(c_{1}\ln(r)\cdot H)\\\ +c_{2}\ln(r)\end{array}$
2 | $\frac{\tilde{r}^{2}}{G(\theta)+F(\varphi)}$ | $\frac{\tilde{G}(\theta)+\tilde{F}(\varphi)}{G(\theta)+F(\varphi)}$ | $\begin{array}[]{c}L_{3}^{2}-(F(\varphi)\cdot H)+\tilde{F}(\varphi)\end{array}$
3 | $\frac{\tilde{r}^{4}}{c_{1}(x_{3}^{2}+r^{2})e^{\pm 2\varphi}+c_{2}x_{3}\tilde{r}e^{\pm\varphi}+c_{3}\tilde{r}^{2}}$ | $\frac{c_{4}(x_{3}^{2}+r^{2})e^{\pm 2\varphi}+c_{5}x_{3}\tilde{r}e^{\pm\varphi}+c_{6}\tilde{r}^{2}}{c_{1}(x_{3}^{2}+r^{2})e^{\pm 2\varphi}+c_{2}x_{3}\tilde{r}e^{\pm\varphi}+c_{3}\tilde{r}^{2}}$ | $\begin{array}[]{c}\\{P_{3},(D\mp L_{3})\\}\\\ +2\left((\frac{c_{2}e^{\pm\varphi}}{\tilde{r}}+\frac{c_{2}x_{3}e^{\pm 2\varphi}}{\tilde{r}^{2}})\cdot H\right)\\\ -2(\frac{c_{5}e^{\pm\varphi}}{\tilde{r}}+\frac{c_{6}x_{3}e^{\pm 2\varphi}}{\tilde{r}^{2}}),\\\ \\{K_{3},(D\mp L_{3})\\}-15x_{3}\\\ +2\left(r^{2}(\frac{c_{2}e^{\pm\varphi}}{\tilde{r}}+\frac{c_{2}x_{3}e^{\pm 2\varphi}}{\tilde{r}^{2}})\cdot H\right)\\\ -2r^{2}(\frac{c_{5}e^{\pm\varphi}}{\tilde{r}}+\frac{c_{6}x_{3}e^{\pm 2\varphi}}{\tilde{r}^{2}})\end{array}$
4 | $\frac{x_{3}^{2}\omega}{F({\omega+z})+G(\omega-z)}$ | $\frac{\tilde{F}({\omega+z})+\tilde{G}(\omega-z)}{F({\omega+z})+G(\omega-z)}$ | $\begin{array}[]{c}\mu\\{P_{1},K_{1}\\}+\nu\\{P_{2},K_{2}\\}\\\ +\left(g(F,G)\cdot H\right)+\eta(\tilde{F},\tilde{G})\end{array}$
5 | $\begin{array}[]{c}\sum_{m}c_{m}\Phi_{m},\end{array}$ | $\begin{array}[]{c}\frac{\sum_{m}\tilde{c}_{m}\Phi_{m}}{\sum_{m}c_{m}\Phi_{m}}\end{array}$ | $\begin{array}[]{c}\mu\\{P_{1},K_{1}\\}+\kappa\\{P_{2},K_{2}\\}\\\ -\left(\rho\sum_{m}c_{m}\Phi_{m-1}\cdot H\right)\\\ +\rho\sum_{m}\tilde{c}_{m}\Phi_{m-1}-\hat{\eta}\end{array}$

In Table 1 and the following Tables 2–5 $F(.),G(.),\tilde{F}(.),\tilde{G}(.)$
are arbitrary functions whose arguments are fixed in the brackets,
$\omega,z,\rho$ and $\Phi_{m}$ are given by relations (81) and (83),
$c_{1},c_{2},...$ are arbitrary real parameters,
$\hat{\eta}=3r^{2}+15x_{3}^{2}$.

###  4.5 Symmetries with differential terms $Q_{(n)}$ (36) quadratic in
spatial variables

The last class of symmetries which we are supposed to find are integrals of
motion (7) generated by Killing tensors (39). The specificity of this class is
that it includes integrals of motions with different transformation properties
with respect to the discrete symmetry (25). Namely, the vector integrals of
motion generated by $\mu_{(1)}^{ab}$ are invariant with respect to
transformations (25) while the tensor integrals of motion generated by
$\mu_{(2)}^{ab}$ change their sign under this transformation. It means that in
fact we have two independent subclasses of integrals of motion since the
linear combinations of integrals of motion with different transformations
properties with respect to symmetry transformations are forbidden.

Let us start with vector integrals of motion. Up to rotation transformations
the only nonzero parameter in $\mu_{(1)}^{ab}$ (39) is $\lambda^{3}$, and the
related equations (52) are reduced to the following form:

$\displaystyle\begin{split}&(x_{1}^{2}-x_{2}^{2})M_{2}-2x_{1}x_{2}M_{1}-x_{3}x_{2}M_{3}+N_{1}=0,\\\
&(x_{1}^{2}-x_{2}^{2})M_{1}+2x_{1}x_{2}M_{2}+x_{3}x_{1}M_{3}+N_{2}=0,\\\
&x_{3}x_{1}M_{2}-x_{3}x_{2}M_{1}+N_{3}=0.\end{split}$ (70)

Using relation (64) and the following condition

$\displaystyle x_{a}N_{a}=0$ (71)

(which is a consequence of (70)) equations (70) can be reduced to the
following system including two equations only:

$\displaystyle\begin{split}&x_{1}^{2}M_{2}-x_{1}x_{2}M_{1}+N_{1}=0,\\\
&x_{1}x_{2}M_{2}-x_{2}^{2}M_{1}+N_{2}=0.\end{split}$ (72)

Solutions of the system including equations (72) and (64) are presented in
Item 1 of Table 1.

Consider now the tensor integrals of motion generated by $\mu_{(2)}^{ab}$
(39). Up to rotation transformations and arbitrary element $g({\bf x})$ we can
restrict ourselves to the following version of nonzero parameters
$\lambda^{ab}$:

$\displaystyle\lambda^{33}=1$ (73)

or

$\displaystyle\lambda^{11}=\mu,\ \lambda^{22}=\nu$ (74)

where $\mu$ and $\nu$ are arbitrary real constants.

In the case (73) integral of motion (7) and determining equations (52) take
the following forms:

$\displaystyle Q=L_{3}^{2}+Q_{0}+\eta$ (75)

and

$\displaystyle\begin{split}&x_{3}x_{\alpha}M_{3}=N_{\alpha},\ \ \
\alpha=1,2,\\\ &\tilde{r}^{2}M_{3}-x_{3}x_{a}M_{a}=N_{3}.\end{split}$ (76)

The system (64), (76) is solved by the following functions

$\displaystyle M=\frac{G(\theta)+F(\varphi)}{\tilde{r}^{2}},\ N=F(\varphi)$

where $G(\theta)$ and $F(\varphi)$ are arbitrary functions of Euler angles,
refer to Item 2 of Table 1 where the corresponding potential is represented
also.

For some special form of the mentioned arbitrary functions represented in Item
2 of Table 1 the related PDM system admits more integrals of motion and is
superintrgrable or even maximally superintegrable, see Items 1–3, of Table 2.
One more superintegrable system which admit two integrals of motion of type
(75) includes the arbitrary elements represented in the following formulae:

$\displaystyle\begin{split}&f=\frac{x_{1}^{2}x_{2}^{2}x_{3}^{2}r^{2}}{x_{3}^{2}(x_{3}^{2}(c_{1}x_{2}^{2}+c_{2}x_{1}^{2})+c_{1}x_{2}^{4}+c_{2}x_{1}^{4})+c_{4}\tilde{r}^{2}x_{1}^{2}x_{2}^{2}},\\\
&V=\frac{x_{3}^{2}(x_{3}^{2}(c_{5}x_{2}^{2}+c_{6}x_{1}^{2})+c_{5}x_{2}^{4}+c_{6}x_{1}^{4})+c_{8}\tilde{r}^{2}x_{1}^{2}x_{2}^{2}}{x_{3}^{2}(x_{3}^{2}(c_{1}x_{2}^{2}+c_{2}x_{1}^{2})+c_{1}x_{2}^{4}+c_{2}x_{1}^{4})+c_{4}\tilde{r}^{2}x_{1}^{2}x_{2}^{2}},\end{split}$
(77)

while the related integrals of motion look as:

$\displaystyle\begin{split}&Q_{1}=L_{1}^{2}-\left(\frac{c_{2}x_{3}^{2}+c_{4}x_{2}^{2}}{x_{2}^{2}x_{3}^{2}}\cdot
H\right)+\frac{c_{6}x_{3}^{2}+c_{8}x_{2}^{2}}{x_{2}^{2}x_{3}^{2}},\\\
&Q_{2}=L_{2}^{2}-\left(\frac{c_{1}x_{3}^{4}+c_{2}x_{1}^{4}}{x_{1}^{2}x_{3}^{2}}\cdot
H\right)+\frac{c_{5}x_{3}^{4}+c_{6}x_{1}^{4}}{x_{1}^{2}x_{3}^{2}}.\end{split}$
(78)

The case (74) is the most complicated and interesting. The related determining
equations (52) is reduced to the following system:

$\displaystyle\begin{split}&\mu(x_{2}^{2}+x_{3}^{2})M_{1}-\nu
x_{1}x_{2}M_{2}+N_{1}=0,\\\ &\nu(x_{1}^{2}+x_{3}^{2})M_{2}-\mu
x_{1}x_{2}M_{1}+N_{2}\end{split}$ (79)

where we use relations (64) and (71) to simplify the system.

Excluding unknown variable $N$ we come to the following second order equation
for $M$:

$\displaystyle\left(\mu(x_{2}^{2}+x_{3}^{2})-\nu(x_{1}^{2}+x_{3}^{2})\right)M_{12}+x_{1}x_{2}(\mu
M_{11}-\nu M_{22}+3(\mu x_{2}M_{1}-\nu x_{1}M_{2})=0$ (80)

which looks as drastically complicated one. However, changing dependent and
independent variables in the following way:

$\displaystyle M=\frac{F(\omega,z)}{\omega x_{3}^{2}}$

where

$\displaystyle\begin{split}&\omega=\sqrt{z^{2}+4\rho},\ \rho=\frac{\mu
x_{1}^{2}-\nu x_{2}^{2}}{x_{3}^{2}}+\mu\nu(\mu-\nu),\\\ &z=\frac{\nu
x_{1}^{2}+\mu
x_{2}^{2}+(\mu-\nu)(\mu^{2}-\nu^{2})x_{3}^{2}}{x_{3}^{2}\sqrt{\mu\nu(\mu-\nu))}}\end{split}$
(81)

we can reduce it to the standard D’Alembert equation for function
$F(\omega,z)$

$\displaystyle(\partial_{\omega\omega}-\partial_{zz})F(\omega,z)=0$

and so

$\displaystyle M=\frac{F(\omega+z)+G(\omega-z)}{\omega x_{3}^{2}}.$ (82)

For any fixed $F(\omega+z)$ and $G(\omega-z)$ we can solve the system (79),
(71) and find functions $N=N(F,G)$ corresponding to $M$ defined in (82).
Unfortunately, it is seemed be impossible to represent functions $N=N(F,G)$ in
closed form for $F$ and $G$ arbitrary. However, it is possible do it for
rather extended classes of functions $F$ and $G$. In particular it can be done
if the sum $F(.)+G(.)$ is a homogeneous function of $z$ and $\omega$, say, if
$F(.)+G(.)=\Phi_{m}$ where

$\displaystyle\Phi_{m}=\frac{(\omega+z)^{m}+(-1)^{m+1}(\omega-z)^{m}}{\omega}$
(83)

Substituting (83) into (82) and integrating the corresponding equations (79),
(71) we find the related functions $M$ and $N$ in closed form for $m$
arbitrary:

$\displaystyle M=\frac{\Phi_{m}}{\omega x_{3}^{2}},\ N=-\rho\Phi_{m-1}.$ (84)

Notice that linear combinations of solutions (84) also solve equations (79),
(71). Just these results are presented in Item 5 of Table 1.

Table 2. Inverse masses, potentials and integrals of motion for
superintegrable scale invariant systems.

No | $f$ | $V$ | Integrals of motion
---|---|---|---
1 | $\frac{x_{3}^{2}\tilde{r}^{2}}{c_{1}\tilde{r}^{2}+F(\varphi)x_{3}^{2}}$ | $\frac{G(\varphi)x_{3}^{2}+c_{2}\tilde{r}^{2}}{c_{1}\tilde{r}^{2}+F(\varphi)x_{3}^{2}}$ | $\begin{array}[]{c}P_{3}^{2}+\frac{c_{2}}{x_{3}^{2}}-\left(\frac{c_{1}}{x_{3}^{2}}\cdot H\right),\\\ K_{3}^{2}+\frac{c_{2}r^{4}}{x_{3}^{2}}-\left(\frac{c_{1}r^{4}}{x_{3}^{2}}\cdot H\right)+\hat{\eta},\\\ L_{3}^{2}+G(\varphi)-(F(\varphi)\cdot H)\end{array}$
2 | $\frac{\tilde{r}^{2}r}{F(\varphi)r+c_{1}x_{3}}$ | $\frac{G(\varphi)+c_{2}x_{3}}{F(\varphi)r+c_{1}x_{3}}$ | $\begin{array}[]{c}\\{P_{3},D\\}-\frac{c_{2}}{r}+(\frac{c_{1}}{r}\cdot H),\\\ \\{K_{3},D\\}+(r\cdot H)-c_{2}r-15x_{3},\\\ L_{3}^{2}+G(\varphi)-\left(F(\varphi)\cdot H\right)\end{array}$
3 | $\frac{x_{3}^{2}x_{1}^{2}\tilde{r}}{c_{1}\tilde{r}x_{3}^{2}+c_{2}x_{2}x_{3}^{2}+c_{3}x_{1}^{2}\tilde{r}}$ | $\frac{c_{4}\tilde{r}x_{3}^{2}+c_{5}x_{2}x_{3}^{2}+c_{6}x_{1}^{2}\tilde{r}}{c_{1}\tilde{r}x_{3}^{2}+c_{2}x_{2}x_{3}^{2}+c_{3}x_{1}^{2}\tilde{r}}$ | $\begin{array}[]{c}\\{L_{3},P_{1}\\}+\frac{2c_{4}x_{2}}{x_{1}^{2}}+\frac{c_{5}(\tilde{r}^{2}+x_{2}^{2})}{\tilde{r}x_{1}^{2}}\\\ -\left(H\cdot(\frac{2c_{1}x_{2}}{x_{1}^{2}}+\frac{c_{2}(\tilde{r}^{2}+x_{2}^{2})}{\tilde{r}x_{1}^{2}})\right),\\\ \\{L_{3},K_{1}\\}+\frac{2c_{4}x_{2}r^{2}}{x_{1}^{2}}+\frac{c_{5}r^{2}(\tilde{r}^{2}+x_{2}^{2})}{\tilde{r}x_{1}^{2}}\\\ -\left(H\cdot(\frac{2r^{2}c_{1}x_{2}}{x_{1}^{2}}+\frac{c_{2}r^{2}(\tilde{r}^{2}+x_{2}^{2})}{\tilde{r}x_{1}^{2}})\right),\\\ L_{3}^{2}-\left(\frac{c_{1}\tilde{r}^{2}+c_{2}x_{2}\tilde{r}}{x_{1}^{2}}\cdot H\right)+\frac{c_{4}\tilde{r}^{2}+c_{5}x_{2}\tilde{r}}{x_{1}^{2}},\\\ P_{3}^{2}-\left(\frac{c_{3}}{x_{3}^{2}}\cdot H\right)+\frac{c_{6}}{x_{3}^{2}},\\\ K_{3}^{2}-\left(\frac{c_{3}r^{4}}{x_{3}^{2}}\cdot H\right)+\frac{c_{6}r^{4}}{x_{3}^{2}}+\hat{\eta}\par\end{array}$
4 | $\frac{x_{1}^{2}x_{2}^{2}x_{3}^{2}}{c_{1}x_{2}^{2}x_{3}^{2}+c_{2}x_{1}^{2}x_{3}^{2}+c_{3}x_{1}^{2}x_{2}^{2}}$ | $\frac{c_{4}x_{2}^{2}x_{3}^{2}+c_{5}x_{1}^{2}x_{3}^{2}+c_{6}x_{1}^{2}x_{2}^{2}}{c_{1}x_{2}^{2}x_{3}^{2}+c_{2}x_{1}^{2}x_{3}^{2}+c_{3}x_{1}^{2}x_{2}^{2}}$ | $\begin{array}[]{c}P_{2}^{2}-\left(\frac{c_{2}}{x_{2}^{2}}\cdot H\right)+\frac{c_{5}}{x_{2}^{2}},\\\ P_{1}^{2}-\left(\frac{c_{1}}{x_{1}^{2}}\cdot H\right)+\frac{c_{4}}{x_{1}^{2}},\\\ K_{2}^{2}-\left(\frac{c_{2}r^{4}}{x_{2}^{2}}\cdot H\right)+\frac{c_{5}r^{4}}{x_{2}^{2}}\\\ +3r^{2}+15x_{2}^{2},\\\ K_{1}^{2}-\left(\frac{c_{1}r^{4}}{x_{1}^{2}}\cdot H\right)+\frac{c_{4}r^{4}}{x_{1}^{2}},\\\ +3r^{2}+15x_{1}^{2},\\\ L_{1}^{2}-\left(\frac{c_{3}x_{2}^{4}+c_{2}x_{3}^{4}}{x_{2}^{2}x_{3}^{2}}\cdot H\right)+\frac{c_{6}x_{2}^{4}+c_{5}x_{3}^{4}}{x_{2}^{2}x_{3}^{2}},\\\ L_{2}^{2}-\left(\frac{c_{3}x_{1}^{4}+c_{1}x_{3}^{4}}{x_{1}^{2}x_{3}^{2}}\cdot H\right)+\frac{c_{6}x_{1}^{4}+c_{3}x_{3}^{4}}{x_{1}^{2}x_{3}^{2}}\end{array}$

Thus we find all inequivalent PDM systems which are scale invariant and admit
at least one second order integral of motion. They are represented in Tables 1
and 2 where $F(.),G(.),\tilde{F}(.),\tilde{G}(.)$ are arbitrary functions
whose arguments are fixed in the brackets, $\omega,z,\rho$ and $\Phi_{m}$ are
given by relations (81) and (83), $c_{1},c_{2},...$ are arbitrary real
parameters. One more scale invariant and superintegrable PDM system is
presented in equations (77) and (LABEL:equa) since the related arbitrary
elements are too cumbersome to be placed in the table.

Notice that systems whose arbitrary elements are present in Table 1 are
integrable while the superintegrable systems are represented in Table 2 and
equation ones are superintegrable. Moreover, Items 3 and 4 ot this table
represent maximally superintegrable systems.

Notice that systems whose arbitrary elements are present in Table 1 are
integrable while the superintegrable systems are represented in Table 2 and
equation ones are superintegrable. Moreover, Items 3 and 4 ot this table
represent maximally superintegrable systems.

##  5 PDM systems invariant with respect to shifts

The last class of PDM system we discus in the present paper are those ones
which by definition admit Lie symmetry whose generator is $P_{3}$. Since the
related arbitrary elements $M$ and $V$ do not depend on $x_{3}$ it is possible
a priori reduce the number of admissible second order integrals of motion.

Let Q be an integral of motion (7) admitted by equation (1) with arbitrary
elements (9). By definition $P_{3}$ is the integral of motion too, the same is
true for the commutators $[P_{3},Q],[P_{3},[P_{3},Q]$ and
$[P_{3},[P_{3},[P_{3},Q]].$ Thus any second order symmetry induces the
symmetry generated by $\mu^{ab}$ independent on $x_{3}$ which are itemized
below together with one additional special tensor:

  * •

$\mu_{1}^{ab}$ with arbitrary $\lambda_{1}^{ab}$,

  * •

$\mu_{2}^{ab}$ with $\lambda_{2}^{3}=0$,

  * •

$\mu_{3}^{ab}$ and $\mu_{6}^{ab}$ where the only nonzero parameters are
$\lambda_{3}^{33}$ and $\lambda_{6}^{33},$

  * •

$\mu_{2}^{ab}$ with $\lambda_{2}^{\alpha}=0,\ \lambda_{2}^{3}\neq 0$

refer to (26)-(34).

The integrals of motion (7) generated by $\mu^{ab}$ specified in the first
three items commute with $P_{3}$ and so do not induce additional symmetries.
However, the analogous property belongs to integral of motion generated by
$\mu_{2}^{ab}$ with $\lambda_{2}^{3}\neq 0$. Its commutator with $P_{3}$ is
not trivial but proportional to $P_{3}^{2}$ which is by now means an
additional symmetry for the considered class of PDM systems.

Classifying PDM systems admitting these induced symmetries we can obtain the
complete list of the systems invariant with respect to shifts and admitting
second order integrals of motion. Then it will be necessary to specify the
additional integrals of motion which can be admitted by the found systems. As
a result we obtain the complete list of the systems presented in the following
Tables 3-5. The detailed calculations are presented in Sections 5.1- 5.4.

###  5.1 Symmetries whose differential part $Q_{(n)}$ (36) does not depend on
$x_{1},x_{2}$ and $x_{3}$

The first step of our analysis is to find the systems with the special
property specified in the subsection title. The corresponding Killing tensors
are specified in the items presented in the above.

Let us start with the constant Killing tensors
$\mu^{ab}=\mu^{ab}_{1}=\lambda_{1}^{ab}$. The related integrals of motion look
as follows:

$\displaystyle Q_{1}=\lambda^{ab}P_{a}P_{b}+\left(N\cdot H\right)+\eta$ (85)

where we omit the subindex 1 for $\lambda_{1}^{ab}$.

Up to rotations around the third coordinate axis and the arbitrariness of
function $g({\bf x})$ we can restrict ourselves to the following nonzero
entries: $\lambda^{3\alpha},\alpha=1,2$, $\lambda^{11}=-\lambda^{22}=\mu$ and
$\lambda^{12}=\lambda^{21}=\nu$ where $\mu$ and $\nu$ are arbitrary real
numbers satisfying the condition $\mu\nu=0$. The corresponding operator (85)
and determining equations (19) take the following forms:

$\displaystyle Q_{2}=\mu(P_{1}^{2}-P_{2}^{2})+\nu
P_{1}P_{2}+\lambda^{31}P_{3}P_{1}+\lambda^{32}P_{3}P_{2}+\left(N\cdot
H\right)+\eta$ (86)

and

$\displaystyle\begin{split}&\mu M_{1}+\nu M_{2}+N_{1}=0,\\\ &\nu M_{1}-\mu
M_{2}+N_{2}=0,\\\ &\lambda^{31}M_{1}+\lambda^{32}M_{2}+N_{3}=0.\end{split}$
(87)

Let $\lambda^{31}=\lambda^{32}=\nu=0$ than solutions of equations (87) are

$\displaystyle M=F(x_{1})+G(x_{2}),\ N=\mu(-F(x_{1})+G(x_{2}))$ (88)

where $F(.)$ and $G(.)$ are arbitrary functions. The generic form of the
related potential $V$ is the same as it has $M$ given in the above, but with
another arbitrary functions, i.e.,

$\displaystyle V=\tilde{F}(x_{1})+\tilde{G}(x_{2})$ (89)

while the corresponding symmetry is reduced to the following form:

$\displaystyle Q_{2}=P_{1}^{2}-P_{2}^{2}+\left((-F(x_{1})+G(x_{2}))\cdot
H\right)+\tilde{F}(x_{1})-\tilde{G}(x_{2}).$ (90)

Just this system is represented in Item 1 of Table 3.

In the case of trivial $\mu$ and $\nu$ and nontrivial $\lambda^{3\alpha}$ we
can apply the rotations again and reduce $\lambda^{32}$ to zero. As a result
we come to the following solution

$\displaystyle M=c_{1}x_{1}+G(x_{2}),\ V=c_{2}x_{1}+\tilde{G}(x_{2}),\
N=-\lambda^{31}c_{1}x_{3}.$ (91)

The corresponding symmetry (87) looks as:

$\displaystyle Q_{3}=2P_{3}P_{1}-(c_{1}x_{3}\cdot H)+c_{2}x_{3}.$ (92)

Notice that function $M$ in (91) is nothing but a particular case of the
function presented in (88). It means that the PDM system with the mass and
potential given in (91) admits both symmetries (89) and (92). This speculation
is supported by the following observation: functions $M$ and $V$ fixed by
equation (91) solve equations (87) where $\nu=\lambda^{33}=0$ but both $\mu$
and $\lambda^{31}$ are not trivial. The obtained solutions are presented in
Item 2 of Table 4 where the additional symmetry of PDM hamiltonian with
arbitrary elements (91) is represented also.

Finally, if all parameters $\mu,\lambda^{31}$ and $\lambda^{32}$ are
nontrivial, the related solutions of (87) for the mass and potential are:

$\displaystyle M=c_{1}x_{1}+c_{2}x_{2},\ V=c_{3}x_{1}+c_{4}x_{2}+c_{5}$ (93)

and we have again a particular case of $M$ represented in (88) (and in (91)).
In other words, the PDM Hamiltonian with the mass and potential given in (93)
admits three independent integrals of motion which we represent in Item 2 of
Table 4 together with the additional linearly independent but functionally
dependent ones.

Table 3. Inverse masses, potentials and integrals of motion for integrable
shift invariant systems

No | $f$ | $V$ | Integrals of motion
---|---|---|---
1 | $\frac{1}{F(x_{1})+G(x_{2})}$ | $\frac{\tilde{F}(x_{1})+\tilde{G}(x_{2})}{F(x_{1})+G(x_{2})}$ | $\begin{array}[]{c}P_{1}^{2}-P_{2}^{2}+\tilde{F}(x_{1})-\tilde{G}(x_{2})\\\ +\left((G(x_{2})-F(x_{1}))\cdot H\right)\end{array}$
2 | $\frac{\tilde{r}^{2}}{F(\tilde{r})-G(\varphi)}$ | $\frac{N(\varphi)+M(\tilde{r})}{G(\varphi)-F(\tilde{r})}$ | $L_{3}^{2}+(G(\varphi)\cdot H)+N(\varphi)$
3 | $\begin{array}[]{c}\frac{\tilde{r}}{F(z_{+})+G(z_{-})},\\\ z_{\pm}=\tilde{r}\pm x_{2}\end{array}$ | $\frac{W(z_{+})+K(z_{-})}{F(z_{+})+G(z_{-})}$ | $\begin{array}[]{c}\\{L_{3},P_{1}\\}+P_{a}x_{2}P_{a}\\\ -((F(z_{+})-G(z_{-}))\cdot H)\\\ +W(z_{+})-K(z_{-})\end{array}$
4 | $\frac{1}{c_{1}\varphi+F(\tilde{r})}$ | $\frac{G(\tilde{r})+c_{2}\varphi}{c_{1}\varphi+F(\tilde{r})}$ | $\\{P_{3},L_{3}\\}-c_{1}(x_{3}\cdot H)+c_{2}x_{3}$
5 | $\frac{\tilde{r}^{2}}{\lambda\tilde{r}^{2}+F(\varphi-\nu\ln(\tilde{r}))}$ | $\frac{\kappa\tilde{r}^{2}+G(\varphi-\nu\ln(\tilde{r}))}{\lambda\tilde{r}^{2}+F(\varphi-\nu\ln(\tilde{r}))}$ | $\begin{array}[]{c}\\{P_{3},D\\}+\nu\\{P_{3},L_{3}\\}\\\ -\left(2\lambda x_{3}\cdot H\right)+2\kappa x_{3}\end{array}$
6 | $\frac{x_{2}^{2}}{x_{2}^{2}F(x_{1})+c_{1}}$ | $\frac{x_{2}^{2}\tilde{F}(x_{1})+c_{2}}{x_{2}^{2}F(x_{1})+c_{1}}$ | $\begin{array}[]{l}L_{1}^{2}-\left(\frac{c_{1}x_{3}^{2}}{x_{2}^{2}}\cdot H\right)+\frac{c_{2}x_{3}^{2}}{x_{2}^{2}},\\\ \\{P_{2},L_{1}\\}+2\left(\frac{c_{1}x_{3}}{x_{2}^{2}}\cdot H\right)-2\frac{c_{2}x_{3}}{x_{2}^{2}},\\\ P_{1}^{2}-P_{2}^{2}+\left((\frac{c_{1}}{x_{2}^{2}}-F(x_{1}))\cdot H\right)\\\ -\frac{c_{2}}{x_{2}^{2}}+\tilde{F}(x_{1})\end{array}$
7 | $\frac{\tilde{\omega}}{F(\tilde{\omega}+y)+G(\tilde{\omega}-y)}$ | $\frac{\tilde{F}(\tilde{\omega}+y)+\tilde{G}(\tilde{\omega}-y)}{F(\tilde{\omega}+y)+G(\tilde{\omega}-y)}$ | $\begin{array}[]{c}L_{3}^{2}\pm(P_{1}^{2}-P_{2}^{2})\\\ +\left(g(F,G)\cdot H\right)-g(\tilde{F},\tilde{G})\end{array}$
8 | $\begin{array}[]{c}\sum_{m}c_{m}\Omega_{m},\end{array}$ | $\begin{array}[]{c}\frac{\sum_{m}\tilde{c}_{m}\Omega_{m}}{\sum_{m}c_{m}\Omega_{m}}\end{array}$ | $\begin{array}[]{c}L_{3}^{2}\pm(P_{1}^{2}-P_{2}^{2})\\\ -\left(y\sum_{m}c_{m}\Phi_{m-1}\cdot H\right)\\\ +y\sum_{m}\tilde{c}_{m}\Omega_{m-1}\end{array}$

In Table 3 we use the following notations:
$\Omega_{m}=\frac{(\tilde{\omega}+y)^{m}+(-1)^{m+1}(\tilde{\omega}-y)^{m}}{\tilde{\omega}},$
$y=x_{1}^{2}-x_{2}^{2}\mp 1,\ \tilde{\omega}=\sqrt{\tilde{r}^{4}-4y},c_{m}$
with $m=1,2,..$ are arbitrary parameters. This table represents all
inequivalent PDM systems which admit two integrals of motion, one of which is
the generator of shift transformations along the third coordinate axis and the
other is a second order integral of motion. The two additional integrals of
motion presented in Item 6 can be expressed via commutators of $P_{3}$ and
$L_{1}^{2}$ and so are not independent. In other words the systems presented
in Table 3 are integrable while superintegrable systems are represented in
Table 4.

###  5.2 Integrals of motion whose differential part $Q_{(n)}$ (36) is
independent on $x_{3}$ but depends on $x_{1}$ and $x_{2}$

The second step of our programm is to classify the systems whose integrals of
motion are generated by the Killing tensors specified at the beginning of
Section 5, i.e., by the following linear combinations:

$\displaystyle\begin{split}&\mu^{\alpha\beta}=\lambda^{\alpha\beta}_{1}+\lambda_{2}^{\alpha}x_{\beta}+\lambda_{2}^{\beta}x_{\alpha}+\lambda_{6}^{33}(\delta^{\alpha\beta}\tilde{r}^{2}-x_{\alpha}x_{\beta}),\\\
&\mu^{3\alpha}=\lambda^{33}_{3}\varepsilon^{3\alpha
c}x_{c}+\lambda_{1}^{3\alpha}+\lambda_{2}^{3}x_{\alpha}\end{split}$ (94)

where $\alpha$ and $\beta$ can take the values 1, 2.

Using the rotations around the third coordinate axis we reduce
$\lambda_{2}^{2}$ to zero, and denote
$\lambda_{6}^{33}=\sigma,\lambda_{1}^{11}=-\lambda_{1}^{22}=\omega,\lambda_{1}^{12}=\lambda_{1}^{21}=\kappa,\lambda^{3}_{2}=\alpha,\lambda_{3}^{33}=\mu$.
The corresponding integrals of motion (7) take the following form:

$\displaystyle\begin{array}[]{l}Q=\omega(P_{1}^{2}-P_{2}^{2})+\kappa
P_{1}P_{2}+\lambda^{3\alpha}P_{3}P_{\alpha}+\mu
L_{3}P_{3}+\nu\\{P_{1},L_{3}\\}+\sigma L_{3}^{2}+\alpha\\{P_{3},D\\}\\\
+(N\cdot H)+\eta\end{array}$ (97)

while the related determining equations (19) are reduced to the following
system:

$\displaystyle\begin{split}&(x_{2}(\sigma
x_{2}-2\nu)+\omega)M_{1}+((x_{1}(\nu-\sigma x_{2})+\kappa)M_{2}+N_{1}=0,\\\
&(x_{1}(\nu-\sigma x_{2})+\kappa)M_{1}+(\sigma
x_{1}^{2}+\omega)M_{2}+N_{2}=0,\\\ &(\lambda^{31}-\mu x_{2}+\alpha
x_{1})M_{1}+(\mu x_{1}+\lambda^{32}+\alpha x_{2})M_{2}+2\alpha
M+N_{3}=0\end{split}$ (98)

where we used the shifts of variables $x_{1}$ and $x_{2}$ to nullify the
maximally possible number of parameters $\lambda_{1}^{ab}$.

System (98) has non-constant and nonlinear solutions if parameters
$\lambda^{3\alpha}$ are zero or can be reduced to zero by shifts. Linear
solutions admit shift transformations and so have to be ignored. Thus without
loss of generality these parameters can be nullified.

Equations (98) with $\lambda^{3\alpha}$ have consistent (i.e., non-invariant
with respect to additional Lie groups) solutions provided we impose one of the
following conditions to the remaining parameters:

$\displaystyle\mu\neq 0,\ \omega=0,\ \kappa=0,\ \nu=0$ (99)
$\displaystyle\mu=0,\ \sigma\neq 0,\ \kappa=0,\ \omega\nu=0,$ (100)
$\displaystyle\mu=\sigma=0,\ \nu\neq 0,\ \omega=\kappa=0,$ (101)

where we used shift transformations to nullify as many arbitrary parameters as
possible.

Consider versions (99)-(101) consequently. Let arbitrary parameters satisfy
conditions (99) then equations (98) are reduced to the following form:

$\displaystyle\begin{split}&\sigma x_{2}^{2}M_{1}-\sigma
x_{1}x_{2}M_{2}+N_{1}=0,\\\ &\sigma x_{1}^{2}M_{2}-\sigma
x_{1}x_{2}M_{1}+N_{2}=0,\\\ &(\alpha x_{1}-\mu x_{2})M_{1}+(\mu x_{1}+\alpha
x_{2})M_{2}+2\alpha M+N_{3}=0\end{split}$ (102)

Solutions of system (102) for the cases $\alpha=0$ and $\alpha\neq 0$ are
represented in Items 4 and 5 of Table 3 respectively.

Table 4. Inverse masses, potentials and integrals of motion for shift
invariant superintegrable systems

No | $f$ | $V$ | Integrals of motion
---|---|---|---
1 | $\frac{1}{c_{1}x_{1}+c_{2}x_{2}+2c_{3}\tilde{r}^{2}+c_{4}}$ | $\frac{c_{5}x_{1}+c_{6}x_{2}+2c_{7}\tilde{r}^{2}+c_{8}}{c_{1}x_{1}+c_{2}x_{2}+2c_{3}\tilde{r}^{2}+c_{4}}$ | $\begin{array}[]{c}P_{1}^{2}-P_{2}^{2}\\\ +\left((c_{2}x_{2}-c_{1}x_{1}+c_{3}(x_{2}^{2}-x_{1}^{2}))\cdot H\right)\\\ +c_{5}x_{1}-c_{6}x_{2}+c_{7}(x_{1}^{2}-x_{2}^{2}),\\\ P_{1}P_{2}-\frac{1}{2}\left(((2c_{3}+c_{2})x_{1}+c_{1}x_{2})\cdot H\right)\\\ +\frac{1}{2}((2c_{7}+c_{6})x_{1}+c_{5}x_{2})\end{array}$
2 | $\frac{1}{c_{1}x_{1}+F(x_{2})}$ | $\frac{c_{2}x_{1}+\tilde{F}(x_{2})}{c_{1}x_{1}+F(x_{2})}$ | $\begin{array}[]{c}P_{1}^{2}-P_{2}^{2}+c_{1}x_{1}-\tilde{F}(x_{2})\\\ +\left((F(x_{2})-c_{1}x_{1})\cdot H\right),\\\ 2P_{3}P_{1}-(c_{1}x_{3}\cdot H)+c_{2}x_{3},\\\ \\{P_{3},L_{2}\\}-\frac{c_{1}}{2}(x_{3}^{2}\cdot H)+\frac{1}{2}c_{2}x_{3}^{2}\end{array}$
3 | $\frac{x_{1}^{2}x_{2}^{2}}{c_{1}x_{1}^{2}+c_{2}x_{2}^{2}+c_{3}x_{1}^{2}x_{2}^{2}\tilde{r}^{2}+c_{4}x_{1}^{2}x_{2}^{2}}$ | $\frac{c_{5}x_{1}^{2}+c_{6}x_{2}^{2}+c_{7}x_{1}^{2}x_{2}^{2}\tilde{r}^{2}+c_{8}x_{1}^{2}x_{2}^{2}}{c_{1}x_{1}^{2}+c_{2}x_{2}^{2}+c_{3}x_{1}^{2}x_{2}^{2}\tilde{r}^{2}+c_{4}x_{1}^{2}x_{2}^{2}}$ | $\begin{array}[]{c}L_{3}^{2}-\left(\frac{c_{1}x_{1}^{4}+c_{2}x_{2}^{4}}{x_{1}^{2}x_{2}^{2}}\cdot H\right)+\frac{c_{5}x_{1}^{4}+c_{6}x_{2}^{4}}{x_{1}^{2}x_{2}^{2}},\\\ P_{1}^{2}-P_{2}^{2}\\\ +\left(\frac{c_{1}^{2}x_{2}^{2}-c_{2}^{2}x_{2}^{2}-c_{3}(x_{1}^{2}-x_{2}^{2})x_{1}^{2}x_{2}^{2}}{x_{1}^{2}x_{2}^{2}}\cdot H\right)\\\ -\frac{c_{5}^{2}x_{2}^{2}-c_{6}^{2}x_{2}^{2}-c_{7}(x_{1}^{2}-x_{2}^{2})x_{1}^{2}x_{2}^{2}}{x_{1}^{2}x_{2}^{2}}\end{array}$
4 | $\frac{x_{1}^{2}}{c_{1}(x_{1}^{2}+4x_{2}^{2})x_{1}^{2}+c_{2}x_{1}^{2}x_{2}+c_{3}+c_{4}x_{1}^{2}}$ | $\frac{c_{5}(x_{1}^{2}+4x_{2}^{2})x_{1}^{2}+c_{6}x_{1}^{2}x_{2}+c_{7}+c_{8}x_{1}^{2}}{c_{1}(x_{1}^{2}+4x_{2}^{2})x_{1}^{2}+c_{2}x_{1}^{2}x_{2}+c_{3}+c_{4}x_{1}^{2}}$ | $\begin{array}[]{c}P_{1}^{2}-P_{2}^{2}-c_{6}x_{2}+\frac{c_{7}}{x_{1}^{2}}\\\ +\left((c_{1}(4x_{2}^{2}-x_{1}^{2})+c_{2}x_{2})\cdot H\right)\\\ -\left(\frac{c_{3}}{x_{1}^{2}}\cdot H\right)-c_{5}(4x_{2}^{2}-x_{1}^{2}),\\\ \\{P_{1},L_{3}\\}+\frac{4x_{2}(c_{5}x_{1}^{4}-c_{7})+c_{5}x_{1}^{4}}{2x_{1}^{2}}\\\ -\left(\frac{4x_{2}(c_{1}x_{1}^{4}-c_{3})+c_{1}x_{1}^{4}}{2x_{1}^{2}}\cdot H\right)\end{array}$
5 | $\frac{\tilde{r}^{2}}{c_{1}\tilde{r}^{2}+F(\varphi)}$ | $\frac{c_{2}\tilde{r}^{2}+G(\varphi)}{c_{1}\tilde{r}^{2}+F(\varphi)}$ | $\begin{array}[]{c}\\{P_{3},D\\}-2c_{1}(x_{3}\cdot H)+2c_{2}x_{3},\\\ L_{3}^{2}-(F(\varphi)\cdot H)-G(\varphi),\\\ D^{2}-c_{1}(r^{2}\cdot H)+c_{2}r^{2}\end{array}$
6 | $\frac{\tilde{r}^{2}}{c_{1}(\mu\ln(\tilde{r})-\nu\varphi)+c_{2}\tilde{r}^{2}+c_{3}}$ | $\frac{c_{4}\mu(\ln(\tilde{r})-\nu\varphi)+c_{5}\tilde{r}^{2}+c_{6}}{c_{1}\mu(\ln(\tilde{r})-\nu\varphi)+c_{2}\tilde{r}^{2}+c_{3}}$ | $\begin{array}[]{c}\mu P_{3}L_{3}+\nu\\{P_{3},D\\}\\\ -\left(2c_{2}\nu x_{3}\cdot H\right)+2c_{5}\nu x_{3},\\\ L_{3}^{2}+2(c_{1}\varphi\cdot H)-2c_{4}\varphi\end{array}$
7 | $\frac{\tilde{r}x_{1}^{2}}{c_{1}x_{1}^{2}+c_{2}x_{2}+c_{3}\tilde{r}x_{1}^{2}+c_{4}\tilde{r}}$ | $\frac{c_{5}x_{1}^{2}+c_{6}x_{2}+c_{7}\tilde{r}x_{1}^{2}+c_{8}\tilde{r}}{c_{1}x_{1}^{2}+c_{2}x_{2}+c_{3}\tilde{r}x_{1}^{2}+c_{4}\tilde{r}}$ | $\begin{array}[]{c}L_{3}^{2}-\left(\frac{c_{2}\tilde{r}x_{2}+c_{4}x_{2}^{2}}{x_{1}^{2}}\cdot H\right)-\frac{c_{6}\tilde{r}x_{2}+c_{8}x_{2}^{2}}{x_{1}^{2}},\\\ \\{P_{1},L_{3}\\}-\frac{2c_{8}x_{2}\tilde{r}+c_{6}(x_{1}^{2}+2x_{2}^{2})+c_{5}x_{2}x_{1}^{2}}{x_{1}^{2}\tilde{r}}\\\ +\left(\frac{2c_{4}x_{2}\tilde{r}+c_{2}(x_{1}^{2}+2x_{2}^{2})+c_{1}x_{2}x_{1}^{2}}{x_{1}^{2}\tilde{r}}\cdot H\right)\end{array}$

The next versions to be considered correspond to the trivial parameter $\mu$.
In the case (100) equations (98) look as follows:

$\displaystyle\begin{split}&x_{2}(\sigma x_{2}-2\nu)M_{1}+x_{1}(\nu-\sigma
x_{2})M_{2}+N_{1}=0,\\\ &x_{1}(\nu-\sigma x_{2})M_{1}+\sigma
x_{1}^{2}M_{2}+N_{2}=0,\\\ &\alpha(x_{1}M_{1}+x_{2}M_{2})+2\alpha
M+N_{3}=0\end{split}$ (103)

for $\omega=0$, and

$\displaystyle\begin{split}&x_{2}(\sigma x_{2}+\omega)M_{1}-\sigma
x_{1}x_{2}M_{2}+N_{1}=0,\\\ &\sigma(x_{1}^{2}M_{2}-x_{1}x_{2}M_{1})-\omega
M_{2}+N_{2}=0,\\\ &\alpha(x_{1}M_{1}+x_{2}M_{2})+2\alpha M+N_{3}=0\end{split}$
(104)

for $\nu=0$.

Solution of equations (103) with $\nu=0$ and $\alpha=0$ are presented in Item
2 of Table 2. If $\nu\neq 0$ and $\alpha=0$ we obtain the following solutions:

$\displaystyle M=c_{1}+\frac{F(\varphi)}{\tilde{r}^{2}},\ N=-2\nu
c_{1}x_{3}-\sigma F(\varphi).$ (105)

Function $M$ (108) solves the related equations (103) with arbitrary
parameters $\nu$ and $\sigma$. In particular we can set $\nu=1,\sigma=0$ or
$\nu=0,\sigma=1$, and so the PDM Schrödinger equation with the mass presented
in (108) admits as minimum two integrals of motion, i.e,

$\displaystyle\\{P_{3},D\\}-2c_{1}(x_{3}\cdot H)+2c_{2}x_{3}$ (106)

and

$\displaystyle L_{3}^{2}-(F(\varphi)\cdot H)-G(\varphi),$ (107)

see Item 6 of Table 3.

Analogously, for $\nu\neq 0$ and $\alpha\neq 0$ we come to solutions presented
in Item 7 of Table 3. The corresponding mass function is a particular case of
function $M$ presented in (108).

Consider now equations (104). If parameters $\sigma,\omega$ and $\alpha$ are
nontrivial we easy obtain solutions presented in Item 4 of Table 3, if
$\alpha=0$ we have solutions presented in Item 3 of the same table.

The most complicated case is $\alpha=0$ and nontrivial $\sigma$ and $\omega$.
Up to scaling of variables $x_{a}$ we can set $\sigma=\pm\omega=1$ and reduce
(104) to the following form:

$\displaystyle\begin{split}&(x_{2}^{2}\pm 1)M_{1}-x_{1}x_{2}M_{2}+N_{1}=0,\\\
&(x_{1}^{2}\mp 1)M_{2}-x_{1}x_{2}M_{1}+N_{2}=0,\\\ &N_{3}=0\end{split}$ (108)

Excluding variable $N$ we come to the following consistency condition of this
system:

$\displaystyle(x_{2}^{2}-x_{1}^{2}\pm
2)M_{12}+x_{1}x_{2}(M_{11}-M_{22})+3(x_{2}F_{1}-x_{1}F_{2})=0$ (109)

which looks rather cumbersome. However, choosing new dependent and independent
variables, namely,

$\displaystyle\tilde{M}=\tilde{\omega}M,\
\tilde{x}_{1}=\tilde{\omega}+x_{1}^{2}-x_{2}^{2}\pm 1,\
\tilde{x}_{2}=\tilde{\omega}+x_{2}^{2}-x_{1}^{2}\mp 1$ (110)

where $\tilde{\omega}=\sqrt{\tilde{r}^{4}+4(x_{2}^{2}-x_{1}^{2}\mp 1)}$ we can
reduce (117) to the standard D’Aalembert equation for
$\tilde{M}=\tilde{M}(\tilde{x}_{1},\tilde{x}_{2})$ whose solutions are

$\displaystyle\tilde{M}=F_{1}(\tilde{x}_{1})+F_{2}(\tilde{x}_{2})$ (111)

where $F_{1}(.)$ and $F_{2}(.)$ are arbitrary functions. Just these solutions
are presented in Item 8 of Table 3.

For any fixed functions $F_{1}(.)$ and $F_{2}(.)$ we can solve equation (108)
and find the corresponding function $N$. Unfortunately it is impossible to
find $N$ in closed form for $F_{1}(.)$ and $F_{2}(.)$ arbitrary. However, it
can be done for the special class of these functions presented in Item 8 of
Table 2.

Let us represent the system of determining equations (98) corresponding to the
remaining version (101):

$\displaystyle\begin{array}[]{l}\nu(x_{1}M_{2}-2x_{2}M_{1})+N_{1}=0,\\\ \nu
x_{1}M_{1}+N_{2}=0,\\\ \alpha(x_{1}M_{1}+x_{2}M_{2}+2M)+N_{3}=0\end{array}$
(115)

Their solutions for $\alpha=0,\nu\neq 0$, $\alpha\neq 0,\nu=0$ and
$\alpha\nu\neq 0$ are represented in Item 3 of Table 3, Item 6 of Table 4 and
Item 3 of Table 5 respectively.

###  5.3 Systematic search for superintegrable systems

In the previous section we searched for integrable systems which admit one
second order integral of motion. Just such systems are collected in Table 3.

However some of the found systems accidentally admit two independent second
order integrals of motion and so are superintegrable. Such systems are
collected in Table 4, see Items 2-7 there.

But it is necessary to provide a systematic search for superintegrable and
maximally superintegrable systems. To do it we will start with the solutions
found in Sections 5.1 and 5.2, and test their compatibility with equations
(98) for all admissible sets of arbitrary parameters
$\alpha,\mu,\nu,\sigma,\omega,\kappa,\lambda^{31},\lambda^{32}$ discussed in
the above.

Let us start with the system whose arbitrary parameters are represented in
Item 1 of Table 3, refer to equations (88) and (89). They satisfy equations
(87) with $\lambda^{31}=\lambda^{32}=\nu=0$ and $\mu=1$. Our goal is to verify
if they can satisfy these equations with another values of these parameters.

We already know that it is the case if the generic potential (89) is reduced
to the forms presented in (91) and (93). The related equations (87) include
nontrivial parameters $\lambda^{3\alpha}$ while parameter $\nu$ was proposed
to be trivial. Notice that the triviality of this parameter is supported by
the possibility to reduce it to zero using the rotation transformation.

But what has happen wether the system admits a pair of integrals of motion of
generic form (85)? In general we can nullify the coefficient $\nu$ only in one
of them, while the second integral can include the term $\nu P_{1}P_{2}$. In
other words, the system admitting integral of motion (90) and so including the
mass function $M$ presented in (88) can admit one more integral of motion(85)
provided this function satisfies equations (87) with $\nu\neq 0$ but $\mu=0$.
Solving these equations for all possible values of arbitrary parameters we
come to the following versions.

$\displaystyle\lambda^{31}=\lambda^{32}=0,\
M=c_{1}x_{1}+c_{2}x_{2}+c_{3}\tilde{r}^{2}+c_{4},\
V=c_{5}x_{1}+c_{6}x_{2}+c_{7}\tilde{r}^{2}+c_{8},$ (116)
$\displaystyle(\lambda^{31})^{2}+(\lambda^{32})^{2}\neq 0,\
M=c_{1}x_{1}+c_{2}x_{2}+c_{3},\ V=c_{4}x_{1}+c_{5}x_{2}+c_{6}.$ (117)

Version (117) had been already discussed in Section 5.1. Integrals of motion
for the PDM systems with arbitrary elements (116) are represented in Item 2 of
Table 4.

The next step is to verify wether the PDM system with arbitrary elements (88)
and (89) can admit additional integral of motion of generic form (97). To do
it we are supposed to solve the determining equations (98) for function $M$
represented in (88). Substituting (88) into (98) we come to the following
system:

$\displaystyle\begin{split}&(\sigma x_{2}^{2}-\nu x_{2})F_{1}+(\nu
x_{1}+\kappa-\sigma x_{1}x_{2})G_{2}+N_{1}=0,\\\ &\sigma
x_{1}^{2}G_{2}+((\nu-\sigma x_{2})x_{1}+\kappa)F_{1}+N_{2}=0,\\\ &(\alpha
x_{1}-\mu x_{2})F_{1}+(\mu x_{1}+\alpha x_{2})G_{2}+N_{3}=0.\end{split}$ (118)

Notice that whenever parameter $\sigma$ is not trivial it is possible to
reduce parameter $\nu$ to zero making a shift of variable $x_{2}$, so we can
impose the condition $\sigma\nu=0$.

Equations (118) have nonconstant solutions for $M$ provided their parameters
satisfy one of the following conditions:

$\displaystyle\nu=0,\alpha=0,\mu\neq 0,\sigma\neq 0$ (119)
$\displaystyle\nu=0,\alpha=0,\mu=0,\kappa=0,\sigma\neq 0,$ (120)
$\displaystyle\nu=0,\alpha\neq 0,\mu=0,\kappa=0,\sigma\neq 0,$ (121)
$\displaystyle\sigma=0,\mu=0,\kappa=0,\alpha=0,\nu\neq 0,$ (122)
$\displaystyle\sigma=0,\mu=0,\kappa=0,\alpha\neq 0,\nu\neq 0$ (123)

Considering equations (118) step by step for the values of arbitrary
parameters presented in (119)-(122) we find the following solutions:

$\displaystyle F=c_{1}x_{1}^{2}+c_{2},\ G=c_{1}x_{2}^{2}+c_{3}$ (124)

for the case (119),

$\displaystyle F=c_{1}x_{1}^{2}+c_{2}+\frac{c_{3}}{x_{1}^{2}},\
G=c_{4}x_{2}^{2}+c_{5}+\frac{c_{6}}{x_{2}^{2}}$ (125)

for the case (120),

$\displaystyle F=c_{1}+\frac{c_{2}}{x_{1}^{2}},\
F=c_{3}+\frac{c_{4}}{x_{1}^{2}}$ (126)

for the case (121),

$\displaystyle F=c_{1}x_{1}^{2}+c_{2}+\frac{c_{3}}{x_{1}^{2}},\
G=4c_{1}x_{2}^{2}+c_{5}+{c_{6}}{x_{2}}$ (127)

for the case (122), and

$\displaystyle F=\frac{c_{1}}{x_{1}^{2}}+c_{2},\ G=c_{3}+{c_{4}}{x_{2}}$ (128)

for the case (123).

Version (130) correspond to the system invariant w.r.t. two parametric Lie
group including shifts and rotations around the third coordinate axis and so
can be omitted as considered before in in paper [45]. Versions (132), (127),
and (133) (128) are presented in Items 3, 4 of Table 4 and Items 2, 4 of Table
5 correspondingly.

We have analysed all possible superintegrable systems admitting integral of
motion (90) and find seven versions of such systems. In analogous way it is
possible to classify the superintegrable systems admitting integrals of motion
(106) and find four more superintegrable systems presented in Items 5, 6, 7 of
Table 4 and Item 3 of Table 5. Notice that there are two superintegrable
systems admitting integrals of motion (90) and (106) both, see Item 4 of Table
1 and Item 2 of Table 5, and there are no superintegrable systems which do not
admit either integral (90) or integral (106).

Table 5. Inverse masses, potentials and integrals of motion for shift
invariant maximally superintegrable systems

No | $f$ | $V$ | Integrals of motion
---|---|---|---
1 | $\frac{1}{c_{1}x_{1}+c_{2}x_{2}}$ | $\frac{c_{3}x_{1}+c_{4}x_{2}+c_{5}}{c_{1}x_{1}+c_{2}x_{2}}$ | $\begin{array}[]{l}P_{1}^{2}-P_{2}^{2}+c_{3}x_{1}-c_{4}x_{2}\\\ +\left((c_{2}x_{2}-c_{1}x_{1})\cdot H\right),\\\ 2P_{1}P_{2}+c_{3}x_{2}+c_{4}x_{1}\\\ -((c_{1}x_{2}+c_{2}x_{1})\cdot H),\\\ 2P_{3}P_{1}-(c_{1}x_{3}\cdot H)+c_{3}x_{3},\\\ 2P_{3}P_{2}-(c_{2}x_{3}\cdot H)+c_{4}x_{3},\\\ \\{P_{3},L_{2}\\}-\frac{c_{1}}{2}(x_{3}^{2}\cdot H)+\frac{1}{2}c_{3}x_{3}^{2},\\\ \\{P_{3},L_{1}\\}+\frac{c_{2}}{2}(x_{3}^{2}\cdot H)-\frac{1}{2}c_{4}x_{3}^{2}\end{array}$
2 | $\frac{x_{1}^{2}x_{2}^{2}}{c_{1}x_{1}^{2}+c_{2}x_{2}^{2}+c_{3}x_{1}^{2}x_{2}^{2}}$ | $\frac{c_{4}x_{1}^{2}+c_{5}x_{2}^{2}+c_{6}x_{1}^{2}x_{2}^{2}}{c_{1}x_{1}^{2}+c_{2}x_{2}^{2}+c_{3}x_{1}^{2}x_{2}^{2}}$ | $\begin{array}[]{l}L_{3}^{2}-\left(\frac{c_{1}x_{1}^{4}+c_{2}x_{2}^{4}}{x_{1}^{2}x_{2}^{2}}\cdot H\right)+\frac{c_{5}x_{1}^{4}+c_{6}x_{2}^{4}}{x_{1}^{2}x_{2}^{2}},\\\ L_{1}^{2}-\left(\frac{c_{1}x_{3}^{2}}{x_{2}^{2}}\cdot H\right)+\frac{c_{4}x_{3}^{2}}{x_{2}^{2}},\\\ \\{P_{2},L_{1}\\}-\left(\frac{2c_{1}x_{3}}{x_{2}^{2}}\cdot H\right)+\frac{2c_{4}x_{3}}{x_{2}^{2}},\\\ P_{1}^{2}-P_{2}^{2}+\left(\frac{c_{1}x_{1}^{2}-c_{2}x_{2}^{2}}{x_{1}^{2}x_{2}^{2}}\cdot H\right)-\frac{c_{4}x_{1}^{2}-c_{5}x_{2}^{2}}{x_{1}^{2}x_{2}^{2}},\\\ \\{P_{3},D\\}-2(c_{3}x_{3}\cdot H)+2c_{6}x_{3},\\\ D^{2}-\left(c_{3}x_{3}r^{2}\cdot H\right)+c_{6}x_{3}r^{2}\end{array}$
3 | $\frac{\tilde{r}x_{1}^{2}}{c_{1}x_{2}+c_{2}\tilde{r}x_{1}^{2}+c_{3}\tilde{r}}$ | $\frac{c_{4}x_{2}+c_{5}\tilde{r}x_{1}^{2}+c_{6}\tilde{r}}{c_{1}x_{2}+c_{2}\tilde{r}x_{1}^{2}+c_{3}\tilde{r}}$ | $\begin{array}[]{l}L_{3}^{2}-\left(\frac{c_{1}\tilde{r}x_{2}+c_{3}x_{2}^{2}}{x_{1}^{2}}\cdot H\right)-\frac{c_{4}\tilde{r}x_{2}+c_{6}x_{2}^{2}}{x_{1}^{2}},\\\ \\{P_{1},L_{3}\\}+\left(\frac{2c_{3}x_{2}\tilde{r}+c_{1}(x_{1}^{2}+2x_{2}^{2})}{x_{1}^{2}\tilde{r}}\cdot H\right)\\\ -\frac{2c_{6}x_{2}\tilde{r}+c_{4}(x_{1}^{2}+2x_{2}^{2})}{x_{1}^{2}\tilde{r}},\\\ \\{P_{3},D\\}-\left(2c_{2}x_{3}\cdot H\right)+2c_{5}x_{3},\\\ D^{2}-\left(c_{2}r^{2}\cdot H\right)+c_{5}r^{2},\\\ \\{L_{2},L_{3}\\}-\frac{2c_{6}x_{3}\tilde{r}+c_{4}x_{3}(x_{1}^{2}+2x_{2}^{2})}{x_{1}^{2}\tilde{r}}\\\ +\left(\frac{2c_{3}x_{3}\tilde{r}+c_{1}x_{3}(x_{1}^{2}+2x_{2}^{2})}{x_{1}^{2}\tilde{r}}\cdot H\right)\end{array}$
4 | $\frac{x_{1}^{2}}{x_{2}x_{1}^{2}+c_{1}}$ | $\frac{c_{2}x_{2}x_{1}^{2}+c_{3}+c_{4}x_{1}^{2}}{x_{2}x_{1}^{2}+c_{1}}$ | $\begin{array}[]{l}L_{2}^{2}-\left(\frac{c_{1}x_{3}^{2}}{x_{1}^{2}}\cdot H\right)+\frac{c_{3}x_{3}^{2}}{x_{1}^{2}},\\\ \\{P_{1},L_{2}\\}-\left(\frac{2c_{1}x_{3}}{x_{1}^{2}}\cdot H\right)+\frac{2c_{3}x_{3}}{x_{1}^{2}},\\\ \\{P_{2},D\\}-\frac{1}{2}\left((x_{1}^{2}+4x_{2}^{2}+x_{3}^{2})\cdot H\right)\\\ +\frac{c_{2}}{2}(x_{1}^{2}+4x_{2}^{2}+x_{3}^{2})+2c_{4}x_{2},\\\ P_{1}^{2}-P_{2}^{2}+\left(\frac{x_{2}x_{1}^{2}-c_{1}}{x_{1}^{2}}\cdot H\right)-\frac{c_{2}x_{2}x_{1}^{2}-c_{3}}{x_{1}^{2}},\\\ P_{2}P_{3}-(c_{1}x_{3}\cdot H)+c_{3}x_{3},\\\ \\{P_{3},L_{1}\\}+\left(\frac{c_{1}x_{3}^{2}}{2}\cdot H\right)-\frac{1}{2}c_{3}x_{3}^{2}\end{array}$

###  5.4 Integrals of motion dependent on $x_{3}$

In the previous sections we have found all inequivalent PDM systems which are
shift invariant and admit second order integrals of motion of the generic
forms (85) and (97). The specificity of such integrals of motion is that their
commutators with symmetry $P_{3}$ are trivial or proportional to $P_{3}$.
Moreover, relations (85) and (97) represent the most general forms of second
order integrals of motion with possess this property. On the other hand,
commuting the other possible integrals of motion with $P_{3}$ $n$ times, for
some $n$ we have to come to a symmetry of form (85) or (97). In other words,
we have found all inequivalent shift invariant systems which admit at least
one second order integral of motion. However, we did not present additional
integrals of motion of more general form. Just these integrals are classified
in the previous section.

The first step in the promised classification is to specify the bilinear
combinations of generators (24) whose commutators with $P_{3}$ have the form
presented in (85) or (97). To achieve this goal we will use the following
commutation relations

$\displaystyle\begin{split}&[P_{3},P_{\alpha}]=0,\
[P_{3},D]=-{\mathrm{i}}P_{3},\\\
&[P_{3},L_{1}]={\mathrm{i}}P_{2},[P_{3},L_{2}]=-{\mathrm{i}}P_{1},\
[P_{3},L_{3}]=0,\\\ &[P_{3},K_{3}]=2{\mathrm{i}}D,\
[P_{3},K_{1}]=-2{\mathrm{i}}L_{2},\
[P_{3},K_{2}]=2{\mathrm{i}}L_{1}.\end{split}$ (129)

Let us present the generic form of operators $\tilde{Q}_{2}$ and
$\tilde{Q}_{3}$ whose commutator with $P_{3}$ is proportional to $Q_{2}$ and
$Q_{3}$ (97) respectively:

$\displaystyle\begin{split}&\tilde{Q}_{2}=\frac{\omega}{2}(\\{P_{1},L_{2}\\}+\\{P_{2},L_{1}\\})+\kappa(P_{1}L_{1}-P_{2}L_{2})\\\
&+\frac{1}{2}(\lambda^{32}\\{P_{3},L_{1}\\}-\lambda^{31}\\{P_{3},L_{2}\\})+\left(N_{3}\cdot
H\right)+\eta_{3}\end{split},$ (130)
$\displaystyle\tilde{Q}_{3}=\tilde{Q}_{2}-\mu
L_{3}D-\nu\\{L_{2},L_{3}\\}-\frac{\alpha}{2}\\{P_{3},K_{3}\\}$ (131)

In accordance with (129) operators $Q_{2},Q_{3}$ and $\tilde{Q}_{2},Q_{3}$
satisfy the following relation

$\displaystyle[P_{3},\tilde{Q}_{\alpha}]={\mathrm{i}}Q_{c},c=2,3$ (132)

where $Q_{2}$ is operator (86) and $Q_{3}$ is operator (97) with $\sigma=0$.
Notice that equation (132) has no solutions for $c=3,\sigma\neq 0$.

Let operator $\tilde{Q}_{c}$ with some combinations of arbitrary parameters
$\omega,\kappa,\lambda^{31},\lambda^{32},\mu,\nu,\alpha$ be an integral of
motion of equation (1) with arbitrary elements (9). Then operator $Q_{c}$ with
the same values of the mentioned parameters is the integrals of motion too.
Thus the necessary condition of the commutativity $\tilde{Q}_{c}$ with the
Hamiltonian $H$ is the commutativity of $Q_{c}$ with $H$. Since we already
found such operators $Q_{c}$, we know all admissible combinations of arbitrary
parameters in integrals of motion (130). However, the mentioned necessary
condition is not sufficient one, and so we are supposed to verify wether
operators $\tilde{Q}_{c}$ commute with the found Hamiltonians. Making this
routine job for the admissible values of arbitrary parameters fixed in
(99)-(101) we recover the $x_{3}$ dependent integrals of motion presented in
Items 2, 5 of Table 3 and Items 1-4 of Table 4.

A bit more sophisticated speculations are requested to find the integral of
motion $\hat{Q}=\\{P_{2},L_{1}\\}+2\left(\frac{c_{1}x_{3}}{x_{2}^{2}}\cdot
H\right)-2\frac{c_{2}x_{3}}{x_{2}^{2}}$ presented in Item 6 of Table 2. It
satisfied the following relation:

$\displaystyle[P_{3},\tilde{Q}]={\mathrm{i}}(Q+P_{a}P_{a}-P_{3}^{2})$ (133)

where $Q$ is a particular case of the integral of motion presented in Item 1
of Table 3:

$\displaystyle
Q=P_{1}^{2}-P_{2}^{2}+\left(\left(\frac{c_{1}}{x_{2}^{2}}-F(x_{1})\right)\cdot
H\right)-\frac{c_{2}}{x_{2}^{2}}+\tilde{F}(x_{1}).$ (134)

Relation (133) is more general than (132). Its right hand size includes two
additional terms one of which, i.e., $P_{3}^{2}$, by definition commutes with
the Hamiltonian, while $P_{a}P_{a}$ can be included into $Q_{0}$ by adding the
unity to $g({\bf x})$.

Thus we already know the integrals of motion $\tilde{Q}$ linear in variable
$x_{3}$. The next step is to find such integrals $\hat{Q}$ which are quadratic
polynomials in this variable. By definition they satisfy the following
relation

$\displaystyle[P_{3},\hat{Q}]=\tilde{Q}.$ (135)

Repeating the speculations presented in the above we find integrals of motion

$\displaystyle\begin{split}&\hat{Q}_{1}=L_{1}^{2}-\left(\frac{c_{1}x_{3}^{2}}{x_{2}^{2}}\cdot
H\right)+\frac{c_{2}x_{3}^{2}}{x_{2}^{2}},\\\
&\hat{Q}_{2}=\\{L_{2},L_{3}\\}-\frac{2c_{6}x_{3}\tilde{r}+c_{4}x_{3}(x_{1}^{2}+2x_{2}^{2})}{x_{1}^{2}\tilde{r}}\\\
&+\left(\frac{2c_{3}x_{3}\tilde{r}+c_{1}x_{3}(x_{1}^{2}+2x_{2}^{2})}{x_{1}^{2}\tilde{r}}\cdot
H\right)\end{split}$ (136)

first of which is presented in Items 6 of Table 3 and Item 2 of Table 5 while
the second one can be found in Item 3 of Table 5.

Notice that integrals of motion whose differential part $Q_{(n)}$ (36) being
third and second order polynomials in $x_{3}$ are not admitted by the
considered shift invariant PDM systems.

##  6 Discussion

We classify inequivalent quantum mechanical systems with position dependent
masses which admit second order integrals of motion and one out of two one
parametric Lie groups, i.e., dilatation and shift ones. The total number of
such systems is equal to twenty seven, including nine ones invariant w.r.t.
the dilatation transformations and eighteen shift invariant systems. The
classification results are presented in five tables and include eleven
integrable, ten superintegrable and six maximally superintegrable systems. Let
us note that this statement is a bit conventional since the obtained systems
are defined up to arbitrary functions and (or) arbitrary parameters, and in
fact we enumerate the subclasses of the systems under study including the
infinite sets of them. One more note is that the systems represented in Item 5
of Table 1 and Item 8 of Table 3 form the subclasses of the systems
represented in Item 4 of Table 1 and Item 7 of Table 3 correspondingly.

The presented results can be treated as the next step in our programm of the
complete classification of the 3d PDM systems admitting second order integrals
of motion. Such systems possessing the symmetry w.r.t. three and two
parametric Lie groups were classified in paper [42] and [45] respectively. The
next step presupposes the classification of the systems admitting at least one
parametric Lie symmetry group. In accordance with [24], up to equivalence it
is sufficient to restrict ourselves to the two subclassess of the groups. The
first of them includes rotations around the fixed axis, shifts along this axis
and dilatations, we call these symmetries natural. The second subclass
includes the so called combined transformations and includes superpositions of
shifts and rotations, dilatations and rotations, and shift, rotations and
conformal transformations. The main result of the present paper consists in
completing the classification of the systems admitting the natural symmetries
started in paper [47].

Thus to complete the classification of the 3d PDM systems admitting second
order integrals of motion and a one parametric Lie symmetry group it is
sufficient to classify the the systems admitting the combined symmetries
mentioned in the above. This work is in progress.

Acknowledgement. I am indebted with Universitá del Piemonte Orientale and and
Dipartimento di Scienze e Innovazione Tecnologica for the extended stay as
Visiting Professor.

## References

  * [1] C. R. Hagen, Scale and conformal transformations in Galilean-invariant conformal field theory, Phys. Rev. D 5, 377–388 (1972).
  * [2] U. Niederer, The maximal kinematical invariance group of the free Schrödinger equation, Helv. Phys. Acta, 45, 802–810 (1972).
  * [3] R. L. Anderson, S. Kumei, C. E. Wulfman, Invariants of the equations of wave mechanics. I., Rev. Mex. Fis., 21, 1–33 (1972).
  * [4] C. P. Boyer, The maximal kinematical invariance group for an arbitrary potential, Helv. Phys. Acta, 47, 450–605 (1974).
  * [5] A. G. Nikitin, The maximal ”kinematical” invariance group for an arbitrary potential revised, Journal of Mathematical Physics, Analysis, Geometry 14, 519-531 (2018).
  * [6] A. G. Nikitin, Symmetries of Schrödinger equation with scalar and vector potentials, J. Phys. A: 53, 455202 (2020).
  * [7] A. G. Nikitin, Symmetries of the Schrödinger-Pauli equation for neutral particles, J. Math. Phys. 62, 083509 (2021).
  * [8] A. G, Nikitin Symmetries of the Schrödinger-Pauli equations for charged particles and quasirelativistic Schrödinger equations, J. Phys. A 55, 115202 (2022).
  * [9] P. Winternitz, J. Smorodinsky, M. Uhlir̆, and I. Fris̆. Symmetry groups in classical and quantum mechanics. Yad. Fiz. 4, 625-635 (1966) (English translation: Sov. J. Nucl. Phys. 4, 444-450 (1967)).
  * [10] A. Makarov, J. Smorodinsky, Kh. Valiev and P. Winternitz. A systematic search for non-relativistic systems with dynamical symmetries. Nuovo Cim. A 52, 1061-1084 (1967).
  * [11] N. Evans. Group theory of the Smorodinsky-Winternitz system. J. Math. Phys. 32, 3369-3375 (1991).
  * [12] N. W. Evans. Super-integrability of the Winternitz system. Phys. Lett. 147, 483-486 (1990).
  * [13] W. Miller, Jr, S. Post and P. Winternitz Classical and Quantum Superintegrability with Applications J. Phys. A: Math. Theor. 46 423001 (2013).
  * [14] G.P. Pron’ko and Y.G. Stroganov, ”New example of quantum mechanical problem with hidden symmetry”, Sov. Phys. JETP 45, 1075–1077 (1977).
  * [15] P. Winternitz and I. Yurdusen. Integrable and superintegrable systems with spin. J. Math. Phys. 47 103509 (2006).
  * [16] J. -F. Désilets, P. Winternitz and I. Yurdusen. Superintegrable systems with spin and second-order integrals of motion. Phys. A: Math. Theor. 45 475201 (2012).
  * [17] Nikitin A G 2012 Matrix superpotentials and superintegrable systems for arbitrary spin J. Phys. A: Math. Theor. 45 225205
  * [18] Nikitin A G 2012 New exactly solvable systems with Fock symmetry J. Phys. A: Math. Theor. 45 485204
  * [19] Nikitin A G 2013 Laplace-Runge-Lenz vector for arbitrary spin J. Math. Phys. 54 123506
  * [20] Nikitin A G 2014 Laplace-Runge-Lenz vector with spin in any dimension J. Phys. A: Math. Theor. 47 375201
  * [21] E. Ferraro, N. Messina and A.G. Nikitin, ”Exactly solvable relativistic model with the anomalous interaction”, Phys. Rev. A 81, 042108 (2010).
  * [22] O. Rosas-Ortiz, Position-dependent mass systems: Classical and quantum pictures. In Geometric Methods in Physics XXXVIII, Birkhäuser, Cham, pp. 351-361 (2020).
  * [23] O. von Roos Position-dependent effective masses in semiconductor theory Phys. Rev. B 27, 7547–7552 ( 1983).
  * [24] A. G. Nikitin and T. M. Zasadko, Superintegrable systems with position dependent mass, J. Math. Phys. 56, 042101 (2015).
  * [25] A.G Nikitin and T. M. Zasadko, Group classification of Schrödinger equations with position dependent mass, J. Phys. A: Math. Theor. 49, 365204 (2016).
  * [26] A. G. Nikitin, Kinematical invariance groups of the 3d Schrödinger equations with position dependent masses, J. Math. Phys. 58, 083508 (2017).
  * [27] G. Koenigs, Sur les géodésiques $\grave{a}$ intégrales quadratiques, in Leçons sur la théorie générale des surfaces, Vol. 4, Editor J.G. Darboux, Chelsea Publishing, New York, 1972, 368–404
  * [28] E.G. Kalnins, J.M. Kress, G.S. Pogosyan, W. Miller, Jr., Completeness of superintegrability in two-dimensional constant-curvature spaces, J. Phys. A: Math. Gen. 34, 4705–4720 (2001).
  * [29] E. G. Kalnins, J. M. Kress and W. Miller, Jr., Second order superintegrable systems in conformally flat spaces. II. The classical two-dimensional Stäckel transform, J. Math. Phys. 46 (2005), 053510.
  * [30] J. M. Kress, Equivalence of superintegrable systems in two dimensions, Phys. Atomic Nuclei 70, 560–566 (2007).
  * [31] J. Kress and K. Schöbel, An algebraic geometric classification of superintegrable systems in the Euclidean plane. Journal of Pure and Applied Algebra 223, 1728-1752 (2019).
  * [32] G. Valent, Koenigs Theorem and Superintegrable Liouville Metrics, SIGMA 19, 048 (2023).
  * [33] B. K. Berndson, E. G. Kalnins and W. Miller, Jr. Toward Classification of 2nd Order Superintegrable Systems in 3-Dimensional Conformally Flat Spaces with Functionally Linearly Dependent Symmetry Operators, SIGMA 16 (2020).
  * [34] E. G. Kalnins, J. M. Kress and W. Miller, Jr., Second order superintegrable systems in conformally flat spaces. V. Two- and three-dimensional quantum systems, J. Math. Phys. 47 , 093501 (2006).
  * [35] J. J. Capel and J. M. Kress, Invariant Classification of Second-order Conformally Flat Superintegrable Systems, J. Phys. A: Math. Theor. 47, 495202 (2014).
  * [36] A. Vollmer, Stäckel equivalence of non-degenerate superintegrable systems and invariant quadrics, SIGMA 17, 015 (2021).
  * [37] M. A. Escobar-Ruiz and W. Miller, Toward a classification of semidegenerate 3D superintegrable systems. J. Phys. A: Math. Theor., 50(9), 095203, (2017).
  * [38] J. F. Cariñena, M. F. Rañada and M. Santander, Superintegrability of three-dimensional Hamiltonian systems with conformally Euclidean metrics. Oscillator-related and Kepler-related systems, J. Phys. A: Math. Theor. 54, 105201 (2021).
  * [39] A. Ballesteros, A. Enciso, F. J. Herranz, O. Ragnisco and D. Riglioni, Superintegrable Oscillator and Kepler Systems on Spaces of Nonconstant Curvature via the Stäckel Transform, SIGMA 7, 048 (2011).
  * [40] O. Ragnisco and D. Riglioni, A Family of Exactly Solvable Radial Quantum Systems on Space of Non-Constant Curvature with Accidental Degeneracy in the Spectrum, SIGMA 6 097 (2010).
  * [41] J. Hietarinta, Pure quantum integrability Phys. Lett. A 246 97-104 (1998).
  * [42] A. G. Nikitin, Superintegrable quantum mechanical systems with position dependent masses invariant with respect to three parametric Lie groups, J. Math. Phys. 64, 113505 (2023).
  * [43] A. G. Nikitin. Superintegrable and shape invariant systems with position dependent mass. J. Phys. A: Math. Theor. 48 (2015) 335201
  * [44] J. Niederle, A. G. Nikitin. Relativistic wave equations for interacting, massive particles with arbitrary half-integer spins. Phys. Rev. D 64, 125013 (2001).
  * [45] A. G. Nikitin, Superintegrable quantum mechanical systems with position dependent masses invariant with respect to two parametric Lie groups. J. Phys. A: Math. Theor. 56, 395203 (2023).
  * [46] A. G. Nikitin. Integrable and superintegrable quantum mechanical systems with position dependent masses invariant with respect to one parametric Lie groups. 1. Systems with cylindric symmetry. J. Phys. A: Math. Theor. 57 (2024) 265201 (35pp).
  * [47] A. G. Nikitin,. Superintegrable and Scale-Invariant Quantum Systems with Position-Dependent Mass. Ukrainian Mathematical Journal, 74(3) (2022).
  * [48] A. G. Nikitin, Generalized Killing tensors of arbitrary valence and order, Ukrainian Mathematical Journal 43, 734–743 (1991).