# Quasi-exactly solvable extensions of the Kepler-Coulomb potential on the
sphere

C. Quesne
Physique Nucléaire Théorique et Physique Mathématique, Université Libre de
Bruxelles,
Campus de la Plaine CP229, Boulevard du Triomphe, B-1050 Brussels, Belgium
E-mail address<EMAIL_ADDRESS>
###### Abstract

We consider a family of extensions of the Kepler-Coulomb potential on a
$d$-dimensional sphere and analyze it in a deformed supersymmetric framework,
wherein the starting potential is known to exhibit a deformed shape invariance
property. We show that the members of the extended family are also endowed
with such a property, provided some constraint conditions relating the
potential parameters are satisfied, in other words they are conditionally
deformed shape invariant. Since, in the second step of the construction of a
partner potential hierarchy, the constraint conditions change, we impose
compatibility conditions between the two sets to build quasi-exactly solvable
potentials with known ground and first-excited states. Some explicit results
are obtained for the first three members of the family. We then use a
generating function method, wherein the first two superpotentials, the first
two partner potentials, and the first two eigenstates of the starting
potential are built from some generating function $W_{+}(r)$ [and its
accompanying function $W_{-}(r)$]. From the results obtained for the latter
for the first three family members, we propose some formulas for $W_{\pm}(r)$
valid for the $m$th family member, depending on $m+1$ constants $a_{0}$,
$a_{1}$, …, $a_{m}$. Such constants satisfy a system of $m+1$ linear
equations. Solving the latter allows us to extend the results up to the
seventh family member and then to formulate a conjecture giving the general
structure of the $a_{i}$ constants in terms of the parameters of the problem.

Keywords: Schrödinger equation, Kepler-Coulomb potential, quasi-exactly
solvable potentials, supersymmetry

PACS Nos.: 03.65.Fd, 03.65.Ge

##  1 Introduction

The quantum oscillator problem in a $d$-dimensional space of constant
curvature, which has arisen as a generalization of the classical nonlinear
oscillator introduced by Mathews and Lakshmanan [1], has led to many studies
[2, 3, 4, 5, 6, 7, 8]. Its exact solvability has been shown [9] to be related
to its deformed shape invariance (DSI) property in a deformed supersymmetric
(DSUSY) framework [10]. Some rational extensions of the quantum oscillator on
the sphere have also been proved to be exactly solvable [9] with bound-state
wavefunctions expressible in terms of exceptional orthogonal polynomials (see,
e.g., Ref. [11] and references quoted therein), instead of classical
orthogonal polynomials for the oscillator alone.

Other types of extensions of the oscillator in a $d$-dimensional constant-
curvature space have been constructed and shown to be quasi-exactly solvable
(QES), i.e., with only a finite number of eigenstates that can be found by
algebraic means, while the remaining ones remain unknown. Whereas the simplest
QES problems are characterized by a hidden sl(2, $\mathbb{R}$) algebraic
structure (see, e.g., Ref. [12] and references quoted therein) and are
connected with polynomial solutions of the Heun equation [13], those
considered for the oscillator in curved space are more complicated QES
problems related to generalizations of the Heun equation. In a first approach,
use has been made [14] of the functional Bethe ansatz method [15]. In a second
approach [16], infinite families of QES extensions of the oscillator in curved
space with known ground and first-excited states have been constructed by
combining a conditionally deformed shape invariance (CDSI) method,
generalizing a procedure known for problems in an Euclidean space [17, 18],
with a generating function method developed both in SUSY [19] and DSUSY [20]
backgrounds.

In addition to the oscillator, another potential considered since the
beginning of quantum mechanics has been the Kepler-Coulomb (KC) one. Its study
on the sphere has indeed a long history, since it dates back to Schrödinger
[21], Infeld [22], and Stevenson [23]. Kalnins et al. discussed the separation
of variables of the Schrödinger equation with a KC potential on the sphere
[24]. Higgs [25] and Leemon [26] analyzed the characteristics of the potential
on the $d$-dimensional sphere. Since then, a lot of works have been devoted to
the subject from several viewpoints, such as superintegrability or addition of
some nonlinear terms (see, e.g., [27, 28, 29, 30, 31, 32] and references
quoted therein).

In particular, the quantum KC problem in a constant-curvature space was shown
to exhibit a DSI property in DSUSY and some exactly-solvable rational
extensions of the KC potential in a hyperbolic space were constructed with
bound-state wavefunctions expressed in terms of exceptional orthogonal
polynomials [9].

The purpose of the present paper is to apply the approaches developed in Ref.
[16] to build extensions of the KC potential on the $d$-dimensional sphere
with known ground and first-excited states. As the problem turns out to be
more complicated than the previous one for the oscillator, we plan to exploit
successively the use of the CDSI method and that of the generating function
one.

The plan of the paper is as follows. In Section 2, the Schrödinger equation
for the KC potential on a $d$-dimensional sphere is reviewed, as well as its
DSUSY description and the corresponding DSI property. In Section 3, the CDSI
method is applied to the first three potentials of a family of extensions. In
Section 4, the results obtained in Section 3 are re-interpreted in the
generating function method language and a generalization to arbitrary members
of the extension family is sketched. Finally, Section 5 contains the
conclusion.

##  2 The Kepler-Coulomb problem on the sphere and deformed supersymmetry

In Ref. [33], it was shown that the Lagrangian of the classical nonlinear
oscillator introduced by Mathews and Lakshmanan [1] can be generalized to a
$d$-dimensional Euclidean space, where it is defined by

$L=\frac{1}{2}\frac{m}{1+\lambda
r^{2}}\left(\sum_{i}v_{i}^{2}+\lambda\sum_{i<j}(x_{i}v_{j}-x_{j}v_{i})^{2}\right)-{\cal
V}(r),\qquad{\cal V}(r)=\frac{1}{2}m\alpha^{2}\frac{r^{2}}{1+\lambda r^{2}}.$
(2.1)

Here all summations run over $i,j=1,2,\ldots,d$, and
$r^{2}\equiv\sum_{i}x_{i}^{2}$, with $r$ running on $(0,+\infty)$ or
$(0,1/\sqrt{|\lambda|})$ according to whether the nonlinearity parameter
$\lambda$ is positive or negative. It was also shown that to $L$ corresponds a
Hamiltonian

$H=\frac{1}{2m}\left(\sum_{i}p_{i}^{2}+\lambda\left(\sum_{i}x_{i}p_{i}\right)^{2}\right)+{\cal
V}(r)=\frac{1}{2m}\left((1+\lambda
r^{2})\sum_{i}p_{i}^{2}-\lambda\sum_{i<j}J_{ij}^{2}\right)+{\cal V}(r),$ (2.2)

where $J_{ij}=x_{i}p_{j}-x_{j}p_{i}$. Such a Hamiltonian in Euclidean space
was then proved to be interpretable as describing the Hamiltonian of an
oscillator in a constant-curvature space with the curvature $\kappa$ related
to the nonlinearity parameter $\lambda$ through the relation
$\kappa=-\lambda$.

A similarr approach was later on applied [27] to the case where the nonlinear
oscillator considered in (2.1) and (2.2) is replaced by a nonlinear KC
potential

${\cal V}(r)=-\frac{Q}{r}\sqrt{1+\lambda r^{2}},\qquad Q>0.$ (2.3)

The latter is also amenable to an interpretation as a KC potential in a space
of constant curvature $\kappa=-\lambda$, generalizing to $d$ dimensions the
one considered in three dimensions by Schrödinger [21], Infeld [22], and
Stevension [23].

To quantize (2.2) with ${\cal V}(r)$ given in either (2.1) or (2.3), one
replaces $\sqrt{1+\lambda r^{2}}\,p_{i}$ and $J_{ij}$ by the operators $-{\rm
i}\sqrt{1+\lambda r^{2}}\partial/\partial x_{i}$ and $\hat{J}_{ij}=-{\rm
i}(x_{i}\partial/\partial x_{j}-x_{j}\partial/\partial x_{i})$, respectively,
thus leading to

$\begin{split}{\hat{H}}&=-\left((1+\lambda r^{2})\hat{\Delta}+\lambda
r\frac{\partial}{\partial r}+\lambda{\hat{J}}^{2}\right)+{\cal V}(r)\\\
&=-\left((1+\lambda r^{2})\frac{\partial^{2}}{\partial r^{2}}+(d-1+d\lambda
r^{2})\frac{\partial}{\partial r}-\frac{{\hat{J}}^{2}}{r^{2}}\right)+{\cal
V}(r),\end{split}.$ (2.4)

where ${\hat{J}}^{2}\equiv\sum_{i<j}\hat{J}_{ij}^{2}$ and $\hat{\Delta}$
denotes the Laplacian in a $d$-dimensional Euclidean space (note that here we
take units wherein $\hbar=2m=1$). In the present paper, we are going to
consider more specifically the KC potential (2.3) on a sphere, which means
that $\kappa=-\lambda>0$, leaving the other case $\kappa=-\lambda<0$ for a
further study.

The Schrödinger equation corresponding to (2.4) is separable in hyperspherical
coordinates and gives rise to the radial equation

$\left(-(1-\kappa r^{2})\frac{d^{2}}{dr^{2}}-(d-1-d\kappa
r^{2})\frac{1}{r}\frac{d}{dr}+\frac{l(l+d-2)}{r^{2}}-\frac{Q}{r}\sqrt{1-\kappa
r^{2}}-{\cal E}\right)R(r)=0,$ (2.5)

where $l=0$, 1, 2, … is the angular momentum quantum number, $Q>0$, and the
radial variable $r$ runs over $(0,1/\sqrt{\kappa})$. The differential operator
in (2.5) being formally self-adjoint with respect to $(1-\kappa
r^{2})^{-1/2}r^{d-1}dr$, bound-state radial wavefunctions $R(r)$ should
therefore be normalizable with respect to $(1-\kappa r^{2})^{-1/2}r^{d-1}dr$.

Equation (2.5) can be rewritten in an alternative form as

$\left(-\sqrt{f(r)}\frac{d}{dr}f(r)\frac{d}{dr}\sqrt{f(r)}+V(r)-E\right)\psi(r)=0,$
(2.6)

where

$\begin{split}&f(r)=\sqrt{1-\kappa r^{2}},\\\
&V(r)=V(r;l,Q)=\frac{L(L+1)}{r^{2}}-\frac{Q}{r}f(r),\quad
L=l+\frac{d-3}{2},\\\ &E={\cal E}+\frac{1}{4}\kappa(d-1)^{2},\\\
&\psi(r)=r^{(d-1)/2}f^{-1/2}(r)R(r),\end{split}.$ (2.7)

and can then be interpreted as a deformed Schrödinger equation, written in
terms of a deformed radial momentum $\hat{\pi}_{r}=\sqrt{f(r)}(-{\rm
i}d/dr)\sqrt{f(r)}$, or as a position-dependent mass (PDM) Schrödinger
equation

$\left(-m^{-1/4}(r)\frac{d}{dr}m^{-1/2}(r)\frac{d}{dr}m^{-1/4}(r)+V(r)-E\right)\psi(r)=0$
(2.8)

with $m(r)=1/f^{2}(r)$ [34], the ordering of the PDM $m(r)$ and the
differential operator $d/dr$ being that of Mustafa and Mazharimousavi [35].
Bound-state wavefunctions correspond to functions $\psi(r)$ normalizable with
respect to $dr$ on the interval $(0,1/\sqrt{\kappa})$.

It turns out [9] that with $f(r)$ and $V(r)$ given in (2.7), Eq. (2.6) has an
infinite number of bound-state eigenvalues, given by

$E_{n_{r}}(l,Q)=-\frac{Q^{2}}{4n^{2}}+\kappa n^{2},\qquad
n=n_{r}+l+\frac{d-1}{2},\qquad n_{r}=0,1,2,\ldots,$ (2.9)

with corresponding wavefunctions111Slightly different, but equivalent, forms
were used in Ref. [9]. expressed in terms of Jacobi polynomials with complex
indices and argument,

$\psi_{n_{r}}(r;l,Q)\propto
r^{n}f^{-1/2}\exp\left(-\frac{Q}{2n\sqrt{\kappa}}\arcsin(\sqrt{\kappa}r)\right)P_{n_{r}}^{\left(-n+\frac{{\rm
i}Q}{2n\sqrt{\kappa}},-n-\frac{{\rm
i}Q}{2n\sqrt{\kappa}}\right)}\left(\frac{{\rm i}f}{\sqrt{\kappa}r}\right),$
(2.10)

or in terms of Romanovski polynomials with real indices and argument [36,
37],222Romanovski polynomials are also termed pseudo-Jacobi polynomials [38].

$\psi_{n_{r}}(r;l,Q)\propto
r^{n}f^{-1/2}\exp\left(-\frac{Q}{2n\sqrt{\kappa}}\arcsin(\sqrt{\kappa}r)\right)R_{n_{r}}^{\left(\frac{Q}{n\sqrt{\kappa}},-n+1\right)}\left(\frac{f}{\sqrt{\kappa}r}\right).$
(2.11)

A deformed Schrödinger equation, such as (2.6), can be discussed in terms of
DSUSY [10]. In the case of unbroken DSUSY, after introducing a rescaled
potential

$V_{1}(r)=V(r)-E_{0},$ (2.12)

defined in terms of the ground-state energy $E_{0}$ of Eq. (2.6), one deals
with a pair of partner Hamiltonians, defined by

$\hat{H}_{1,2}=\hat{\pi}_{r}^{2}+V_{1,2}(r)+E_{0},\qquad
V_{1,2}(r)=W^{2}(r)\mp f(r)\frac{dW}{dr},$ (2.13)

in terms of a superpotential

$W(r)=-f\frac{d}{dr}\log\psi_{0}(r)-\frac{1}{2}\frac{df}{dr}.$ (2.14)

Here $\psi_{0}(r)$ denotes the ground-state wavefunction of (2.6), expressed
in terms of the superpotential as

$\psi_{0}(r)\propto
f^{-1/2}\exp\left(-\int^{r}\frac{W(r^{\prime})}{f(r^{\prime})}dr^{\prime}\right).$
(2.15)

The two partner Hamiltonians can be written as

$\hat{H}_{1}=\hat{A}^{+}\hat{A}^{-}+E_{0},\qquad\hat{H}_{2}=\hat{A}^{-}\hat{A}^{+}+E_{0},$
(2.16)

in terms of two first-order differential opearators

$\hat{A}^{\pm}=\mp\sqrt{f(r)}\frac{d}{dr}\sqrt{f(r)}+W(r).$ (2.17)

They intertwine with $\hat{A}^{+}$ and $\hat{A}^{-}$ as

$\hat{A}^{-}\hat{H}_{1}=\hat{H}_{2}\hat{A}^{-},\qquad\hat{A}^{+}\hat{H}_{2}=\hat{H}_{1}\hat{A}^{+}.$
(2.18)

The operator $\hat{A}^{-}$ annihilates the ground-state wavefunction (2.15) of
$\hat{H}_{1}$, whereas $\hat{A}^{+}$ transforms the ground-state wavefunction
$\psi^{\prime}_{0}(r)$ of $\hat{H}_{2}$ into the first-excited state
wavefunction $\psi_{1}(x)$ of $\hat{H}_{1}$.

On considering $\hat{H}_{2}$ as a new starting Hamiltonian, one may in
principle iterate the procedure and obtain a DSUSY pair of partner
Hamiltonians

$\hat{H}^{\prime}_{1,2}=\hat{\pi}_{r}^{2}+V^{\prime}_{1,2}(r)+E^{\prime}_{0},\qquad
V^{\prime}_{1,2}(r)=W^{\prime 2}(r)\mp f(r)\frac{dW^{\prime}}{dr},$ (2.19)

where

$V^{\prime}_{1}(r)+E^{\prime}_{0}=V_{2}(r)+E_{0}.$ (2.20)

The first-excited state wavefunction $\psi_{1}(r)$ of $\hat{H}_{1}$, with
energy $E_{1}=E^{\prime}_{0}$, can then be obtained from the ground-state
wavefunction of $\hat{H}^{\prime}_{1}=\hat{H}_{2}$, given by

$\psi^{\prime}_{0}(r)\propto
f^{-1/2}\exp\left(-\int^{r}\frac{W^{\prime}(r^{\prime})}{f(r^{\prime})}dr^{\prime}\right),$
(2.21)

through the equation

$\psi_{1}(r)\propto\hat{A}^{+}\psi^{\prime}_{0}(r).$ (2.22)

In the case of the KC potential on the sphere of Eq. (2.7), the superpotential
is given by [9]

$W(r)=W(r;l,Q)=-\frac{L+1}{r}f(r)+\frac{Q}{2L+2}$ (2.23)

and the two partner potentials are

$V_{1}(r;l,Q)=\frac{L(L+1)}{r^{2}}-\frac{Q}{r}f(r)-E_{0}$ (2.24)

and

$V_{2}(r;l,Q)=\frac{(L+1)(L+2)}{r^{2}}-\frac{Q}{r}f(r)-E_{0}=V_{1}(r;l+1,Q).$
(2.25)

Hence, the partner is similar in shape and its parameters are obtained by
translation, i.e., $l\to l+1$, $Q\to Q$. This means that the KC potential on
the sphere is DSI, so that a whole hierarchy of Hamiltonians can be
straightforwardly constructed, and this explains the exact solvability of the
starting Schrödinger equation.

In the following Sections, we plan to consider some extensions of the
potential $V(r)$ of Eq. (2.7), namely

$V(r)=\frac{L(L+1)}{r^{2}}-\frac{Q}{r}f+\kappa\sum_{k=1}^{m}\left(B_{2k-1}\frac{r}{f^{2k-1}}+B_{2k}\frac{1}{f^{2k}}\right),$
(2.26)

with $L=l+\frac{d-3}{2}$, $Q>0$, $f(r)=\sqrt{1-\kappa r^{2}}$, and $B_{1}$,
$B_{2}$, …, $B_{2m}$ are $2m$ real parameters, such that $B_{2m}>0$. The
additional terms then break the DSI property of the KC potential, so that the
resulting Schrödinger equation may at most be QES.

##  3 Conditionally deformed shape invariance symmetry method

In the present Section, we plan to consider successivaly the first few members
of the potential family (2.26) and apply the CDSI method of Ref. [16].

###  3.1 First extension

For $m=1$, the potential (2.26) reduces to

$V(r)=\frac{L(L+1)}{r^{2}}-\frac{Q}{r}f+\kappa B_{1}\frac{r}{f}+\kappa
B_{2}\frac{1}{f^{2}},\qquad B_{2}>0,$ (3.1)

and depends on the four parameters $L$, $Q$, $B_{1}$, and $B_{2}$. Let us
assume a superpotential of the type

$W(r)=\frac{\xi}{r}f+\eta+\zeta\frac{r}{f},\qquad\xi\leq 0,\qquad\zeta>0,$
(3.2)

depending on three parameters $\xi$, $\eta$, and $\zeta$. As shown in (2.13),
in DSUSY the rescaled potential (2.12) is represented by
$V_{1}(r)=W^{2}-fdW/dr$. With the extended potential (3.1), this Riccati
equation gives rise to the system of equations

$\begin{split}&\xi(\xi+1)=L(L+1),\\\
&-\kappa\xi^{2}+2\xi\zeta+\eta^{2}-\frac{\zeta^{2}}{\kappa}=-E_{0},\\\
&2\xi\eta=-Q,\\\ &2\eta\zeta=\kappa B_{1},\\\
&\frac{\zeta}{\kappa}\left(\frac{\zeta}{\kappa}-1\right)=B_{2}\end{split}.$
(3.3)

for the three unknowns $\xi$, $\eta$, $\zeta$. The first, third, and fourth of
these equations lead to the values of the latter

$\xi=-L-1,\qquad\eta=\frac{Q}{2(L+1)},\qquad\zeta=\frac{\kappa(L+1)B_{1}}{Q}.$
(3.4)

Since $\zeta>0$, we must have $B_{1}>0$. From the second equation of (3.3),
the ground-state energy $E_{0}$ is obtained in the form

$E_{0}=\kappa(L+1)^{2}-\frac{Q^{2}}{4(L+1)^{2}}+\frac{\kappa(L+1)B_{1}}{Q}\left(\frac{(L+1)B_{1}}{Q}+2L+2\right).$
(3.5)

The last equation of (3.3) then provides a constraint

$B_{2}=\frac{(L+1)B_{1}}{Q}\left(\frac{(L+1)B_{1}}{Q}-1\right),$ (3.6)

connecting the potential parameters. The condition $B_{2}>0$ imposes that
$B_{1}>Q/(L+1)$. To $E_{0}$ corresponds the wavefunction (2.15), which can be
written as

$\psi_{0}(r)\propto
r^{-\xi}\exp\left(-\frac{\eta}{\sqrt{\kappa}}\arcsin(\sqrt{\kappa}r)\right)f^{\frac{\zeta}{\kappa}-\frac{1}{2}},$
(3.7)

with $\xi$, $\eta$, and $\zeta$ given in (3.4).

The partner $V_{2}(r)=W^{2}+fdW/dr$ of $V_{1}(r)$ in DSUSY can be written as

$V_{2}(r)=\frac{L^{\prime}(L^{\prime}+1)}{r^{2}}-\frac{Q^{\prime}}{r}f+\kappa
B^{\prime}_{1}\frac{r}{f}+\kappa B^{\prime}_{2}\frac{1}{f^{2}}+R,$ (3.8)

in terms of some new parameters $L^{\prime}$, $Q^{\prime}$, $B^{\prime}_{1}$,
$B^{\prime}_{2}$, and a constant $R$, which must satisfy the system of
equations

$\begin{split}&L^{\prime}(L^{\prime}+1)=\xi(\xi-1),\\\
&R=-\kappa\xi^{2}+2\xi\zeta+\eta^{2}-\frac{\zeta^{2}}{\kappa},\\\
&-Q^{\prime}=2\xi\eta,\\\ &\kappa B^{\prime}_{1}=2\eta\zeta,\\\
&B^{\prime}_{2}=\frac{\zeta}{\kappa}\left(\frac{\zeta}{\kappa}+1\right).\end{split}$
(3.9)

On using some results obtained above, we get

$L^{\prime}=L+1,\qquad Q^{\prime}=Q,\qquad B^{\prime}_{1}=B_{1},\qquad
B^{\prime}_{2}=B_{2}+2\frac{(L+1)B_{1}}{Q},$ (3.10)

as well as

$R=-E_{0}.$ (3.11)

This shows that the starting potential $V_{1}(r)$ is DSI, but such a deformed
shape invariance is not unconditionally valid since the constraint (3.6) must
be satisfied. The potential is therefore only CDSI.

We may now try to repeat the procedure by taking the partner $V_{2}(r)$ as a
starting potential $V^{\prime}(r)=V_{2}(r)$ with ground-state energy
$E^{\prime}_{0}$. We therefore consider a rescaled potential
$V^{\prime}_{1}(r)=V^{\prime}(r)+E_{0}-E^{\prime}_{0}$ and a new
superpotential

$W^{\prime}(r)=\frac{\xi^{\prime}}{r}f+\eta^{\prime}+\zeta^{\prime}\frac{r}{f},\qquad\xi^{\prime}\leq
0,\qquad\zeta^{\prime}>0.$ (3.12)

By proceeding as in the first step, we obtain for the new parameters

$\xi^{\prime}=-L-2,\qquad\eta^{\prime}=\frac{Q}{2(L+2)},\qquad\zeta^{\prime}=\frac{\kappa(L+2)B_{1}}{Q},$
(3.13)

and for the new ground-state energy

$E^{\prime}_{0}=\kappa(L+2)^{2}-\frac{Q^{2}}{4(L+2)^{2}}+\frac{\kappa(L+2)B_{1}}{Q}\left(\frac{(L+2)B_{1}}{Q}+2L+4\right).$
(3.14)

We also get a new constraint

$B_{2}=\frac{(L+2)B_{1}}{Q}\left(\frac{(L+2)B_{1}}{Q}-1\right)-2\frac{(L+1)B_{1}}{Q},$
(3.15)

relating the potential parameters. The ground-state wavefunction of the
partner reads

$\psi^{\prime}_{0}(r)\propto
r^{-\xi^{\prime}}\exp\left(-\frac{\eta^{\prime}}{\sqrt{\kappa}}\arcsin(\sqrt{\kappa}r)\right)f^{\frac{\zeta^{\prime}}{\kappa}-\frac{1}{2}}.$
(3.16)

The two constraints (3.6) and (3.15) are compatible provided

$B_{1}=Q.$ (3.17)

Then

$B_{2}=L(L+1),\qquad\zeta=\kappa(L+1),$ (3.18)

and the potential

$V(r)=\frac{L(L+1)}{r^{2}}-\frac{Q}{r}f+\kappa
Q\frac{r}{f}+\kappa\frac{L(L+1)}{f^{2}},$ (3.19)

with corresponding superpotentials

$\begin{split}&W(r)=-\frac{L+1}{r}f+\frac{Q}{2(L+1)}+\kappa(L+1)\frac{r}{f},\\\
&W^{\prime}(r)=-\frac{L+2}{r}f+\frac{Q}{2(L+2)}+\kappa(L+2)\frac{r}{f},\end{split}.$
(3.20)

has a ground state and a first-excited state whose energies are given by

$E_{0}=4\kappa(L+1)^{2}-\frac{Q^{2}}{4(L+1)^{2}},\qquad
E_{1}=E^{\prime}_{0}=4\kappa(L+2)^{2}-\frac{Q^{2}}{4(L+2)^{2}}.$ (3.21)

The ground-state wavefunction of potential (3.19) is given by

$\psi_{0}(r)\propto
r^{L+1}\exp\left(-\frac{Q}{2\sqrt{\kappa}(L+1)}\arcsin(\sqrt{\kappa}r)\right)f^{L+\frac{1}{2}},$
(3.22)

which is a normalizable function on $(0,1/\sqrt{\kappa})$. Moreover, its first
excited-state wavefunction $\psi_{1}(r)$ can be obtained from the partner
ground-state wavefunction

$\psi^{\prime}_{0}(r)\propto
r^{L+2}\exp\left(-\frac{Q}{2\sqrt{\kappa}(L+2)}\arcsin(\sqrt{\kappa}r)\right)f^{L+\frac{3}{2}}$
(3.23)

by acting with the operator $\hat{A}^{+}$, defined in (2.17), which can be
rewritten as

$\hat{A}^{+}=-f\frac{d}{dr}-\frac{L+1}{r}f+\frac{Q}{2(L+1)}+\kappa\left(L+\frac{3}{2}\right)\frac{r}{f}.$
(3.24)

The result reads

$\displaystyle\psi_{1}(r)$ $\displaystyle\propto
r^{L+1}\exp\left(-\frac{Q}{2\sqrt{\kappa}(L+2)}\arcsin(\sqrt{\kappa}r)\right)f^{L+\frac{1}{2}}$
$\displaystyle\quad\times\left(-1+2\kappa
r^{2}+\frac{Q}{2(L+1)(L+2)}rf\right),$ (3.25)

which is also a normalizable function on $(0,1/\sqrt{\kappa})$, with a single
zero

$r_{0}=\frac{1}{\sqrt{2\kappa}}\left(1-\frac{Q}{\sqrt{Q^{2}+16\kappa(L+1)^{2}(L+2)^{2}}}\right)^{1/2}$
(3.26)

on this interval, as it should be.

Note that, in the following, it will prove convenient to write the results in
terms of the parameters $L$ and ${\cal Q}=Q/[2(L+1)(L+2)]$, instead of $L$ and
$Q$. With such conventions, we get

$\begin{split}&B_{1}=2(L+1)(L+2){\cal Q},\qquad B_{2}=L(L+1),\\\
&W(r)=-\frac{L+1}{r}f+(L+2){\cal Q}+\kappa(L+1)\frac{r}{f},\\\
&W^{\prime}(r)=-\frac{L+2}{r}f+(L+1){\cal Q}+\kappa(L+2)\frac{r}{f},\\\
&E_{0}=4\kappa(L+1)^{2}-(L+2)^{2}{\cal Q}^{2},\qquad
E_{1}=4\kappa(L+2)^{2}-(L+1)^{2}{\cal Q}^{2}.\end{split}$ (3.27)

In Fig. 1, an example of extended potential (3.19) is plotted and compared
with the starting KC potential. The corresponding (unnormalized) $\psi_{0}(r)$
and $\psi_{1}(r)$ of the former are displayed in Fig. 2.

Figure 1: Plots of potential (3.19) (solid line) and of potential (2.7)
(dashed line), both with $\kappa=L=Q=1$. Figure 2: Plots of ground state
wavefunction $\psi_{0}(r)$ (solid line) and of first-excited state
wavefunction $\psi_{1}(r)$ (dashed line) for potential (3.19) with
$\kappa=L=Q=1$. The corresponding eigenvalues are $E_{0}=15.9375$ and
$E_{1}=35.9722$.

###  3.2 Second extension

For $m=2$, the potential (2.26) reads

$V(r)=\frac{L(L+1)}{r^{2}}-\frac{Q}{r}f+\kappa B_{1}\frac{r}{f}+\kappa
B_{2}\frac{1}{f^{2}}+\kappa B_{3}\frac{r}{f^{3}}+\kappa
B_{4}\frac{1}{f^{4}},\qquad B_{4}>0,$ (3.28)

and we assume a superpotential of the type

$W(r)=\frac{\xi}{r}f+\eta+\zeta\frac{r}{f}+\frac{\sigma}{f^{2}},\qquad\xi\leq
0,\qquad\sigma>0,$ (3.29)

depending on four parameters $\xi$, $\eta$, $\zeta$, and $\sigma$. From the
Riccati equation for $W(r)$, we get the system of equations

$\begin{split}&\xi(\xi+1)=L(L+1),\\\
&-\kappa\xi^{2}+2\xi\zeta+\eta^{2}-\frac{\zeta^{2}}{\kappa}=-E_{0},\\\
&2\xi(\eta+\sigma)=-Q,\\\ &2(\eta\zeta+\kappa\xi\sigma)=\kappa B_{1},\\\
&2\eta\sigma+\frac{\zeta(\zeta-\kappa)}{\kappa}=\kappa B_{2},\\\
&2\sigma(\zeta-\kappa)=\kappa B_{3},\\\ &\sigma^{2}=\kappa B_{4}.\end{split}$
(3.30)

The latter leads to the values of the superpotential parameters

$\xi=-L-1,\quad\eta=\frac{Q}{2(L+1)}-\sqrt{\kappa
B_{4}},\quad\zeta=\kappa\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+1\right),\quad\sigma=\sqrt{\kappa B_{4}},$ (3.31)

the ground-state energy

$E_{0}=\kappa\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+L+2\right)^{2}-\left(\frac{Q}{2(L+1)}-\sqrt{\kappa B_{4}}\right)^{2},$
(3.32)

and two constraints

$\begin{split}&B_{1}=2\left(\frac{Q}{2(L+1)}-\sqrt{\kappa
B_{4}}\right)\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+1\right)-2(L+1)\sqrt{\kappa B_{4}},\\\
&B_{2}=\frac{1}{\kappa}\left[2\left(\frac{Q}{2(L+1)}-\sqrt{\kappa
B_{4}}\right)\sqrt{\kappa B_{4}}+\kappa\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+1\right)\frac{B_{3}}{2\sqrt{\kappa B_{4}}}\right],\end{split}.$ (3.33)

connecting the potential parameters. To $E_{0}$, given in (3.32), corresponds
the (normalizable) ground-state wavefunction

$\psi_{0}(r)\propto
r^{-\xi}\exp\left(-\frac{\eta}{\sqrt{\kappa}}\arcsin(\sqrt{\kappa
r})-\sigma\frac{r}{f}\right)f^{\frac{\zeta}{\kappa}-\frac{1}{2}}.$ (3.34)

The partner potential can be written as

$V_{2}(r)=\frac{L^{\prime}(L^{\prime}+1)}{r^{2}}-\frac{Q^{\prime}}{r}f+\kappa
B^{\prime}_{1}\frac{r}{f}+\kappa B^{\prime}_{2}\frac{1}{f^{2}}+\kappa
B^{\prime}_{3}\frac{r}{f^{3}}+\kappa B^{\prime}_{4}\frac{1}{f^{4}}+R,$ (3.35)

where

$\begin{split}&L^{\prime}=L+1,\qquad Q^{\prime}=Q,\qquad
B^{\prime}_{1}=B_{1},\qquad B^{\prime}_{2}=B_{2}+\frac{B_{3}}{\sqrt{\kappa
B_{4}}}+2,\\\ &B^{\prime}_{3}=B_{3}+4\sqrt{\kappa B_{4}},\qquad
B^{\prime}_{4}=B_{4},\qquad R=-E_{0}.\end{split}$ (3.36)

Hence, $V_{1}(r)$ is CDSI with the two constraints (3.33).

On taking now a new superpotential in the form

$W^{\prime}(r)=\frac{\xi^{\prime}}{r}f+\eta^{\prime}+\zeta^{\prime}\frac{r}{f}+\frac{\sigma^{\prime}}{f^{2}},\qquad\xi^{\prime}\leq
0,\qquad\sigma^{\prime}>0,$ (3.37)

and proceeding as in the previous subsection, we get the relations

$\xi^{\prime}=-L-2,\quad\eta^{\prime}=\frac{Q}{2(L+2)}-\sqrt{\kappa
B_{4}},\quad\zeta^{\prime}=\kappa\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+3\right),\quad\sigma^{\prime}=\sqrt{\kappa B_{4}},$ (3.38)

and

$E^{\prime}_{0}=\kappa\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+L+5\right)^{2}-\left(\frac{Q}{2(L+2)}-\sqrt{\kappa B_{4}}\right)^{2},$
(3.39)

together with two new constraints

$\begin{split}&B_{1}=2\left(\frac{Q}{2(L+2)}-\sqrt{\kappa
B_{4}}\right)\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+3\right)-2(L+2)\sqrt{\kappa B_{4}},\\\
&B_{2}=\frac{1}{\kappa}\left\\{2\left(\frac{Q}{2(L+2)}-\sqrt{\kappa
B_{4}}\right)\sqrt{\kappa B_{4}}+\kappa\left[\left(\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}\right)^{2}+3\frac{B_{3}}{2\sqrt{\kappa
B_{4}}}+4\right]\right\\},\end{split}$ (3.40)

The two sets of constraints (3.33) and (3.40) are compatible for some specific
values of $B_{3}$ and $B_{4}$, which, in terms of the parameters $L$ and
${\cal Q}$, are given by

$B_{3}=\frac{2\kappa(2L+3){\cal Q}[(2L+1){\cal Q}^{2}-6\kappa]}{({\cal
Q}^{2}+3\kappa)^{2}},\qquad B_{4}=\frac{\kappa(2L+3)^{2}{\cal Q}^{2}}{({\cal
Q}^{2}+3\kappa)^{2}}.$ (3.41)

The resulting values of $B_{1}$ and $B_{2}$ are then obtained as

$\begin{split}&B_{1}=2{\cal Q}\frac{2(L+1)(L+2){\cal Q}^{4}-3\kappa{\cal
Q}^{2}-6\kappa^{2}(L^{2}+3L+3)}{({\cal Q}^{2}+3\kappa)^{2}},\\\
&B_{2}=\frac{2(4L^{2}+10L+7){\cal Q}^{4}+\kappa(4L^{2}+3){\cal
Q}^{2}+18\kappa^{2}}{({\cal Q}^{2}+3\kappa)^{2}}.\end{split}$ (3.42)

The corresponding two superpotentials can be written as

$\begin{split}W(r)&=-\frac{L+1}{r}f+\frac{{\cal Q}[(L+2){\cal
Q}^{2}+\kappa(L+3)]}{{\cal Q}^{2}+3\kappa}+\frac{\kappa[2(L+1){\cal
Q}^{2}-3\kappa]}{{\cal Q}^{2}+3\kappa}\frac{r}{f}\\\
&\quad{}+\frac{\kappa(2L+3){\cal Q}}{{\cal Q}^{2}+3\kappa}\frac{1}{f^{2}},\\\
W^{\prime}(r)&=-\frac{L+2}{r}f+\frac{{\cal Q}[(L+1){\cal Q}^{2}+\kappa
L]}{{\cal Q}^{2}+3\kappa}+\frac{\kappa[2(L+2){\cal Q}^{2}+3\kappa]}{{\cal
Q}^{2}+3\kappa}\frac{r}{f}\\\ &\quad{}+\frac{\kappa(2L+3){\cal Q}}{{\cal
Q}^{2}+3\kappa}\frac{1}{f^{2}}.\end{split}.$ (3.43)

The potential $V(r)$ has a ground state and a first-excited state, whose
energies are given by

$\begin{split}&E_{0}=9\kappa\left(\frac{(L+1){\cal Q}^{2}+\kappa L}{{\cal
Q}^{2}+3\kappa}\right)^{2}-{\cal Q}^{2}\left(\frac{(L+2){\cal
Q}^{2}+\kappa(L+3)}{{\cal Q}^{2}+3\kappa}\right)^{2},\\\
&E_{1}=9\kappa\left(\frac{(L+2){\cal Q}^{2}+\kappa(L+3)}{{\cal
Q}^{2}+3\kappa}\right)^{2}-{\cal Q}^{2}\left(\frac{(L+1){\cal Q}^{2}+\kappa
L}{{\cal Q}^{2}+3\kappa}\right)^{2},\end{split}$ (3.44)

respectively, while their corresponding normalizable wavefunctions are given
by

$\displaystyle\psi_{0}(r)$ $\displaystyle\propto r^{L+1}f^{\frac{(4L+3){\cal
Q}^{2}-9\kappa}{2({\cal Q}^{2}+3\kappa)}}\exp\biggl{(}-\frac{{\cal
Q}}{\sqrt{\kappa}}\frac{(L+2){\cal Q}^{2}+\kappa(L+3)}{{\cal
Q}^{2}+3\kappa}\arcsin(\sqrt{\kappa}r)$
$\displaystyle\quad{}-\frac{\kappa(2L+3){\cal Q}}{{\cal
Q}^{2}+3\kappa}\frac{r}{f}\biggr{)}$ (3.45)

and

$\displaystyle\psi_{1}(r)$ $\displaystyle\propto r^{L+1}f^{\frac{(4L+3){\cal
Q}^{2}-9\kappa}{2({\cal Q}^{2}+3\kappa)}}\exp\biggl{(}-\frac{{\cal
Q}}{\sqrt{\kappa}}\frac{(L+1){\cal Q}^{2}+\kappa L}{{\cal
Q}^{2}+3\kappa}\arcsin(\sqrt{\kappa}r)$
$\displaystyle\quad{}-\frac{\kappa(2L+3){\cal Q}}{{\cal
Q}^{2}+3\kappa}\frac{r}{f}\biggr{)}\left[\left(-1+3\frac{{\cal
Q}^{2}+\kappa}{{\cal Q}^{2}+3\kappa}\kappa r^{2}\right)f+{\cal
Q}r\left(1-\frac{{\cal Q}^{2}+\kappa}{{\cal Q}^{2}+3\kappa}\kappa
r^{2}\right)\right].$ (3.46)

###  3.3 Higher extensions

For higher $m$ values, the calculations are similar. On assuming a
superpotential depending on $m+2$ parameters, the Riccati equation for $V_{1}$
provides a system of $2m+3$ equations, from which one gets the values of the
superpotential parameters and of the ground-state energy, plus a set of $m$
constraints connecting the potential parameters. On repeating the procedure
for the partner potential, one ends up with two sets of $m$ constraints, whose
compatibility leads to the final parameter values. As an example, the results
obtained for $m=3$ are listed in Appendix A.

Since the complexity of the CDSI method considerably increases with the $m$
value, it is useful to explore the possibilities offered by the alternative
method available for building potentials with two known eigenstates. This is
the purpose of Section 4.

##  4 Generating function method

Let us start from two pairs of DSUSY partner Hamiltonians
$(\hat{H}_{1},\hat{H}_{2})$ and $(\hat{H}^{\prime}_{1},\hat{H}^{\prime}_{2})$,
defined as in (2.13) and (2.19) in terms of two superpotentials $W(r)$ and
$W^{\prime}(r)$. The latter are related by the equation

$W^{2}(r)+f(r)\frac{dW(r)}{dr}=W^{\prime
2}(r)-f(r)\frac{dW^{\prime}(r)}{dr}+E_{1}-E_{0},$ (4.1)

which directly follows from (2.20), where we set $E_{1}=E^{\prime}_{0}$. From
$W(r)$ and $W^{\prime}(r)$, we may define the two functions

$W_{+}(r)=W^{\prime}(r)+W(r),\qquad W_{-}(r)=W^{\prime}(r)-W(r),$ (4.2)

which allow us to rewrite Eq. (4.1) as

$f(r)\frac{dW_{+}(r)}{dr}=W_{+}(r)W_{-}(r)+E_{1}-E_{0}.$ (4.3)

Hence, $W_{-}(r)$ can be expressed in terms of $W_{+}(r)$ and the energy
difference $E_{1}-E_{0}$ as

$W_{-}(r)=\frac{f(r)dW_{+}(r)/dr+E_{0}-E_{1}}{W_{+}(r)}.$ (4.4)

The generating function method starts from two functions $W_{+}(r)$ and
$W_{-}(r)$ that are compatible, i.e., such that Eq. (4.3) is satisfied for
some positive constant $E_{1}-E_{0}$. The two superpotentials are then
determined from Eq. (4.2). The starting potential $V_{1}(r)$ [as well as its
partner $V_{2}(r)$] follows from Eq. (2.13) and its ground-state wavefunction
from Eq. (2.15). Its first-excited state wavefunction, given in Eq. (2.22)
with $\hat{A}^{+}$ and $\psi^{\prime}_{0}(r)$ defined in (2.17) and (2.21),
respectively, can be rewritten as

$\psi_{1}(r)\propto
W_{+}(r)f^{-1/2}\exp\left(-\int^{r}\frac{W^{\prime}(r^{\prime})}{f(r^{\prime})}dr^{\prime}\right),$
(4.5)

where use is made of Eq. (4.2). It then remains to check that the functions
$\psi_{0}(r)$ and $\psi_{1}(r)$ are normalizable functions on the definition
interval of $r$.

The difficulty of the generating function method is to guess an appropriate
function $W_{+}(r)$ [and its accompanying function $W_{-}(r)$], such that Eq.
(4.3) is satisfied for some positive constant $E_{1}-E_{0}$ and that the
potential and wavefunctions are well behaved.

###  4.1 Application to extensions of the Kepler-Coulomb potential on the
sphere

From the superpotentials $W(r)$ and $W^{\prime}(r)$ obtained in Eqs. (3.20),
(3.43), and (A.3) for the first three members of the extended KC potential
family (2.26), we obtain for the functions $W_{+}(r)$ and $W_{-}(r)$ the
results

$\begin{split}W_{+}(r)&=(2L+3)\left(-\frac{f}{r}+{\cal
Q}+\kappa\frac{r}{f}\right),\\\ W_{-}(r)&=-\frac{f}{r}-{\cal
Q}+\kappa\frac{r}{f},\end{split}$ (4.6)

if $m=1$,

$\begin{split}W_{+}(r)&=(2L+3)\left(-\frac{f}{r}+{\cal Q}\frac{{\cal
Q}^{2}+\kappa}{{\cal Q}^{2}+3\kappa}+\frac{2\kappa{\cal Q}^{2}}{{\cal
Q}^{2}+3\kappa}\frac{r}{f}+\frac{2\kappa{\cal Q}}{{\cal
Q}^{2}+3\kappa}\frac{1}{f^{2}}\right),\\\ W_{-}(r)&=-\frac{f}{r}-{\cal
Q}+2\kappa\frac{r}{f},\end{split}$ (4.7)

if $m=2$, and

$\begin{split}W_{+}(r)&=(2L+3)\biggl{(}-\frac{f}{r}+{\cal Q}\frac{{\cal
Q}^{2}+4\kappa}{{\cal Q}^{2}+10\kappa}+3\kappa\frac{{\cal
Q}^{2}+2\kappa}{{\cal Q}^{2}+10\kappa}\frac{r}{f}+\frac{6\kappa{\cal Q}}{{\cal
Q}^{2}+10\kappa}\frac{1}{f^{2}}\\\ &\quad{}+\frac{6\kappa^{2}}{{\cal
Q}^{2}+10\kappa}\frac{r}{f^{3}}\biggr{)},\\\ W_{-}(r)&=-\frac{f}{r}-{\cal
Q}+3\kappa\frac{r}{f},\end{split}$ (4.8)

if $m=3$. These functions satisfy Eq. (4.3) with the constant $E_{1}-E_{0}$
derived from Eqs. (3.21), (3.44), and (A.4), and given by

$\begin{split}E_{1}-E_{0}&=(2L+3)({\cal Q}^{2}+4\kappa),\\\
E_{1}-E_{0}&=(2L+3)\frac{({\cal Q}^{2}+\kappa)({\cal Q}^{2}+9\kappa)}{{\cal
Q}^{2}+3\kappa},\\\ E_{1}-E_{0}&=(2L+3)\frac{({\cal Q}^{2}+4\kappa)({\cal
Q}^{2}+16\kappa)}{{\cal Q}^{2}+10\kappa},\end{split}$ (4.9)

respectively.

As a generalization for any $m\in\mathbb{N}^{+}$, we propose to consider

$\begin{split}&W_{+}(r)=(2L+3)\left(-\frac{f}{r}+{\cal
Q}a_{0}+\kappa\sum_{k=0}^{[(m-1)/2]}a_{2k+1}\frac{r}{f^{2k+1}}+\sum_{k=1}^{[m/2]}a_{2k}\frac{1}{f_{2k}}\right),\\\
&W_{-}(r)=-\frac{f}{r}-{\cal Q}+m\kappa\frac{r}{f},\\\
&E_{1}-E_{0}=(2L+3)a_{0}[{\cal Q}^{2}+(m+1)^{2}\kappa],\end{split}.$ (4.10)

where $a_{0}$, $a_{1}$, …, $a_{m}$ are $m+1$ constants to be determined in
such a way that Eq. (4.3) is satified. These generalized formulas agree with
the results previously obtained for $m=1$, 2, and 3 if the parameters $a_{0}$,
$a_{1}$, …, $a_{m}$ are chosen as follows:

$\begin{split}a_{0}&=1,\quad a_{1}=1,\quad\text{if $m=1$},\\\
a_{0}&=\frac{{\cal Q}^{2}+\kappa}{{\cal Q}^{2}+3\kappa},\quad
a_{1}=\frac{2{\cal Q}^{2}}{{\cal Q}^{2}+3\kappa},\quad
a_{2}=\frac{2\kappa{\cal Q}}{{\cal Q}^{2}+3\kappa},\quad\text{if $m=2$},\\\
a_{0}&=\frac{{\cal Q}^{2}+4\kappa}{{\cal Q}^{2}+10\kappa},\quad
a_{1}=\frac{3({\cal Q}^{2}+2\kappa)}{{\cal Q}^{2}+10\kappa},\quad
a_{2}=\frac{6\kappa{\cal Q}}{{\cal Q}^{2}+10\kappa},\quad
a_{3}=\frac{6\kappa}{{\cal Q}^{2}+10\kappa},\\\ &\qquad{}\text{if
$m=3$}.\end{split}.$ (4.11)

###  4.2 Results for higher $m$ values

To find the values of the $a_{i}$ coefficients in (4.10), let us insert the
explicit expressions for $W_{+}(x)$, $W_{-}(x)$, and $E_{1}-E_{0}$ in Eq.
(4.3) and use the following identities:

$\begin{split}&\frac{df}{dr}=-\kappa\frac{r}{f},\\\ &\kappa r^{2}=1-f^{2},\\\
&\frac{1}{rf^{2k}}=\frac{1}{r}+\kappa\left(\frac{r}{f^{2}}+\frac{r}{f^{4}}+\cdots+\frac{r}{f^{2k}}\right),\quad
k=1,2,\ldots,\\\
&\frac{1}{rf^{2k+1}}=\frac{f}{r}+\kappa\left(\frac{r}{f}+\frac{r}{f^{3}}+\cdots+\frac{r}{f^{2k+1}}\right),\quad
k=0,1,2,\ldots.\end{split}$ (4.12)

On both sides of (4.3), we then obtain a linear combination of 1, $r/f$,
$r/f^{2k+1}$ ($1\leq k\leq[(m-1)/2]$), $1/f^{2k}$ ($1\leq k\leq[m/2]$), and
$f/r$. On equating their coefficients, we get the relations

$\displaystyle(m+1)a_{0}-a_{1}=1,$ (4.13) $\displaystyle m{\cal Q}a_{0}-{\cal
Q}a_{1}-\sum_{k=1}^{[m/2]}a_{2k}=0,$ (4.14) $\displaystyle(m-2k)a_{2k}-{\cal
Q}a_{2k+1}-\sum_{l=k+1}^{[m/2]}a_{2l}=0,\quad
k=1,2,\ldots,\left[\frac{m-1}{2}\right],$ (4.15)
$\displaystyle\kappa(m+1-2k)a_{2k-1}-{\cal
Q}a_{2k}-\kappa(m+1-2k)a_{2k+1}=0,\quad
k=1,2,\ldots,\left[\frac{m}{2}\right]-1,$ (4.16)
$\displaystyle\frac{1}{2}\bigl{(}1-(-1)^{m}\bigr{)}(2\kappa a_{m-2}-{\cal
Q}a_{m-1}-2\kappa a_{m})+\frac{1}{2}\bigl{(}1+(-1)^{m}\bigr{)}(\kappa
a_{m-1}-{\cal Q}a_{m})=0,$ (4.17) $\displaystyle{\cal Q}-{\cal
Q}a_{0}-\sum_{k=1}^{[m/2]}a_{2k}=0,$ (4.18)

which form a system of $m+1$ linear equations in the $m+1$ unknowns $a_{0}$,
$a_{1}$, …, $a_{m}$ because the last equation (4.18) turns out to be a
combination of Eqs. (4.13) and (4.14). The solutions of this system of
equations for some values of $m$ higher than three are given in Appendix B.

From the results obtained for $m\leq 7$ in (4.11) and in Appendix B, we feel
justified in formulating the following conjecture.

Conjecture For any $m>1$, the solution of the system of equations
(4.13)–(4.17) can be written as

$\displaystyle a_{0}$ $\displaystyle=\frac{({\cal Q}^{2}+4\kappa)({\cal
Q}^{2}+16\kappa)\ldots({\cal Q}^{2}+4\mu^{2}\kappa)}{{\cal
Q}^{2\mu}+b_{1}\kappa{\cal Q}^{2\mu-2}+b_{2}\kappa^{2}{\cal
Q}^{2\mu-4}+\cdots+b_{\mu}\kappa^{\mu}},$ (4.19) $\displaystyle a_{2k}$
$\displaystyle=\frac{(2\mu+1)!}{(2\mu-2k+1)!}\frac{\kappa^{k}{\cal Q}({\cal
Q}^{2}+4\kappa)({\cal Q}^{2}+16\kappa)\ldots({\cal
Q}^{2}+4(\mu-k)^{2}\kappa)}{{\cal Q}^{2\mu}+b_{1}\kappa{\cal
Q}^{2\mu-2}+b_{2}\kappa^{2}{\cal Q}^{2\mu-4}+\cdots+b_{\mu}\kappa^{\mu}},$
$\displaystyle\quad{}k=1,2,\ldots,\mu,$ (4.20) $\displaystyle a_{2k+1}$
$\displaystyle=\frac{(2\mu+1)!}{(2\mu-2k)!}\frac{\kappa^{k}({\cal
Q}^{2\mu-2k}+c^{(k)}_{1}\kappa{\cal Q}^{2\mu-2k-2}+c^{(k)}_{2}\kappa^{2}{\cal
Q}^{2\mu-2k-4}+\cdots+c^{(k)}_{\mu-k}\kappa^{\mu-k})}{{\cal
Q}^{2\mu}+b_{1}\kappa{\cal Q}^{2\mu-2}+b_{2}\kappa^{2}{\cal
Q}^{2\mu-4}+\cdots+b_{\mu}\kappa^{\mu}},$
$\displaystyle\quad{}k=0,1,\ldots,\mu,$ (4.21)

if $m=2\mu+1$, and

$\displaystyle a_{0}$ $\displaystyle=\frac{({\cal Q}^{2}+\kappa)({\cal
Q}^{2}+9\kappa)\ldots({\cal Q}^{2}+(2\mu-1^{2}\kappa)}{{\cal
Q}^{2\mu}+b_{1}\kappa{\cal Q}^{2\mu-2}+b_{2}\kappa^{2}{\cal
Q}^{2\mu-4}+\cdots+b_{\mu}\kappa^{\mu}},$ (4.22) $\displaystyle a_{2k}$
$\displaystyle=\frac{(2\mu)!}{(2\mu-2k)!}\frac{\kappa^{k}{\cal Q}({\cal
Q}^{2}+\kappa)({\cal Q}^{2}+9\kappa)\ldots({\cal
Q}^{2}+(2\mu-2k-1)^{2}\kappa)}{{\cal Q}^{2\mu}+b_{1}\kappa{\cal
Q}^{2\mu-2}+b_{2}\kappa^{2}{\cal Q}^{2\mu-4}+\cdots+b_{\mu}\kappa^{\mu}},$
$\displaystyle\quad{}k=1,2,\ldots,\mu,$ (4.23) $\displaystyle a_{2k+1}$
$\displaystyle=\frac{(2\mu)!}{(2\mu-2k-1)!}$
$\displaystyle\quad{}\times\frac{\kappa^{k}{\cal Q}^{2}({\cal
Q}^{2\mu-2k-2}+c^{(k)}_{1}\kappa{\cal
Q}^{2\mu-2k-4}+c^{(k)}_{2}\kappa^{2}{\cal
Q}^{2\mu-2k-6}+\cdots+c^{(k)}_{\mu-k-1}\kappa^{\mu-k-1}}{{\cal
Q}^{2\mu}+b_{1}\kappa{\cal Q}^{2\mu-2}+b_{2}\kappa^{2}{\cal
Q}^{2\mu-4}+\cdots+b_{\mu}\kappa^{\mu}},$
$\displaystyle\quad{}k=0,1,\ldots,\mu-1,$ (4.24)

if $m=2\mu$, where $b_{1}$, $b_{2}$, …, $b_{\mu}$, and $c^{(k)}_{1}$,
$c^{(k)}_{2}$, …, $c^{(k)}_{\mu-k}$ or $c^{(k)}_{\mu-k-1}$ are some positive
constants.

From (4.2) and (4.10), the superpotential $W(r)$ assumes the form

$\displaystyle W(r)$
$\displaystyle=\frac{1}{2}\biggl{\\{}-(2L+2)\frac{f}{r}+{\cal
Q}[(2L+3)a_{0}+1]+\kappa[(2L+3)a_{1}-m]\frac{r}{f}$
$\displaystyle\quad{}+\kappa(2L+3)\sum_{k=1}^{[(m-1)/2]}a_{2k+1}\frac{r}{f^{2k+1}}+(2L+3)\sum_{k=1}^{[m/2]}a_{2k}\frac{1}{f^{2k}}\biggr{\\}}$
(4.25)

and from the equation $W^{2}-fdW/dr=V-E_{0}$, resulting from (2.12) and
(2.13), one can easily check that $V(r)$ is an extended KC potential, in other
words it contains a term $L(L+1)/r^{2}$, as well as a term

$\displaystyle\frac{1}{4}\biggl{\\{}-2(2L+2){\cal
Q}[(2L+3)a_{0}+1]-2(2L+2)(2L+3)\sum_{k=1}^{[m/2]}a_{2k}\biggr{\\}}\frac{f}{r}$
$\displaystyle\quad=-2(L+1)(L+2){\cal Q}\frac{f}{r}=-Q\frac{f}{r},$ (4.26)

where use has been made of Eq. (4.18). The other parameters $B_{1}$, $B_{2}$,
…, $B_{2m}$ of Eq. (2.26) can in principle be calculated from the same
equation, but their explicit expressions in terms of the parameters $a_{i}$
being very complicated, we will not give them here. In contrast, the constant
term in $W^{2}-fdW/dr$ easily provides the value of $-E_{0}$ in terms of
$a_{0}$ and $a_{1}$,

$\displaystyle-E_{0}$
$\displaystyle=\frac{1}{4}\biggl{\\{}-(2L+2)^{2}\kappa+{\cal
Q}^{2}[(2L+3)a_{0}+1]^{2}-[(2L+3)a_{1}-m^{2}]\kappa$
$\displaystyle\quad{}-2(2L+2)[(2L+3)a_{1}-m]\kappa\biggr{\\}}.$ (4.27)

On expressing $a_{1}$ in terms of $a_{0}$ through Eq. (4.13), we get

$E_{0}=\frac{1}{4}\left\\{-{\cal
Q}^{2}[(2L+3)a_{0}+1]^{2}+(m+1)^{2}\kappa[(2L+3)a_{0}-1]^{2}\right\\}.$ (4.28)

Furthermore, on combining this result with Eq. (4.10), we obtain for the
first-excited state energy

$E_{1}=\frac{1}{4}\left\\{-{\cal
Q}^{2}[(2L+3)a_{0}-1]^{2}+(m+1)^{2}\kappa[(2L+3)a_{0}+1]^{2}\right\\}.$ (4.29)

We conclude that the first two eigenvalues of $V(r)$ are expressed in terms of
the single parameter $a_{0}$.

The corresponding wavefunctions, determined from Eqs. (2.15) and (4.5), are
given by

$\displaystyle\psi_{0}(r)$ $\displaystyle\propto
r^{L+1}f^{\frac{1}{2}[(2L+3)a_{1}-m-1]}$
$\displaystyle\quad{}\times\exp\biggl{\\{}-\frac{{\cal
Q}}{2\sqrt{\kappa}}[(2L+3)a_{0}+1]\arcsin(\sqrt{\kappa}r)-\frac{1}{4}(2L+3)\sum_{k=1}^{[(m-1)/2]}\frac{a_{2k+1}}{kf^{2k}}$
$\displaystyle\quad{}-\frac{1}{2}(2L+3)\sum_{k=1}^{[m/2]}\frac{(2k-3)!!}{(2k-2)!!}\frac{r}{f^{2k-1}}\sum_{l=k}^{[m/2]}\frac{(2l-2)!!}{(2l-1)!!}a_{2l}\biggr{\\}}$
(4.30)

and

$\displaystyle\psi_{1}(r)$ $\displaystyle\propto
W_{+}(r)r^{L+2}f^{\frac{1}{2}[(2L+3)a_{1}+m-1]}$
$\displaystyle\quad{}\times\exp\biggl{\\{}-\frac{{\cal
Q}}{2\sqrt{\kappa}}[(2L+3)a_{0}-1]\arcsin(\sqrt{\kappa}r)-\frac{1}{4}(2L+3)\sum_{k=1}^{[(m-1)/2]}\frac{a_{2k+1}}{kf^{2k}}$
$\displaystyle\quad{}-\frac{1}{2}(2L+3)\sum_{k=1}^{[m/2]}\frac{(2k-3)!!}{(2k-2)!!}\frac{r}{f^{2k-1}}\sum_{l=k}^{[m/2]}\frac{(2l-2)!!}{(2l-1)!!}a_{2l}\biggr{\\}},$
(4.31)

respectively. They are normalizable on the interval $(0,1/\sqrt{\kappa})$
provided $a_{m}>0$, which is the case for all $m$ values for which explicit
calculations have been carried out.

##  5 Conclusion

In the present paper, we have shown that the combination of methods previously
used to build families of extended oscillator potentials in a curved space
with known ground and first-excited states can be applied to extended KC
potentials on the sphere containing $2m$ additional terms, where $m$ may take
any value in $\mathbb{N}^{+}$.

We have first extended the DSI symmetry, known to be valid for the KC
potential alone, by completing it with $m$ constraints relating the parameters
of the extended potential, which therefore turns out to be CDSI. The second
step in the construction of an extended potential hierarchy has then provided
us with a second set of $m$ constraints among the extended potential
parameters. Compatibility conditions between the two sets of $m$ constraints
have given rise to QES extensions of the KC potential with two known
eigenstates. In this way, some explicit results have been obtained for the
first three members of the extended KC potential family.

To generalize those outcomes to other members of the extended KC potential
family, we have then turned ourselves to the generating function method,
wherein the first two DSUSY superpotentials (and hence the potentials and
their first two eigenstates) are expressed in terms of a function $W_{+}(r)$
[and its accompanying function $W_{-}(r)$]. From the results obtained for
$W_{\pm}(r)$ for the first three members of the family, we have proposed some
general formulas for $W_{\pm}(r)$, depending on $m+1$ constants $a_{0}$,
$a_{1}$, …, $a_{m}$, and we have shown that such constants must satisfy a
system of $m+1$ linear equations. On solving such a system, we have obtained
some explicit results for $m=4$, 5, 6, and 7. We have then formulated a
conjecture giving the general structure of the $a_{i}$ constants in terms of
the parameters of the problem. All the aspects of this conjecture and the
consequences drawn from it for the first two eigenstates of $V(r)$ have been
confirmed for $m\leq 7$.

The problem that remains open is to prove our conjecture by finding the
explicit values of the constants appearing there. Solving it would be an
interesting topic for future investigation. Applying the methods presented
here to extensions of the KC potential in a hyperbolic space would also be an
open question for future work.

Taking into account the large number of physical problems wherein the KC
potential (either in curved space or in a PDM background) provides a first
crude approximation and some additional terms are needed to reproduce
experimental data, it is obvious that the exact results presented here for
potentials including some additional terms may find some interesting
applications to practical problems. The search for such applications would be
another interesting topic for future investigation.

## Acknowledgments

The author was supported by the Fonds de la Recherche Scientifique - FNRS
under Grant Number 4.45.10.08.

## Author declarations

### Conflicts of interest

The author declares no conflict of interest.

## Data availability

All data supporting the findings of this study are included in the article.

## Appendix A: Results for the third extension of the Kepler-Coulomb potential

The third extension of the KC potential reads

$V(r)=\frac{L(L+1)}{r^{2}}-\frac{Q}{r}f+\kappa B_{1}\frac{r}{f}+\kappa
B_{2}\frac{1}{f^{2}}+\kappa B_{3}\frac{r}{f^{3}}+\kappa
B_{4}\frac{1}{f^{4}}+\kappa B_{5}\frac{r}{f^{5}}+\kappa B_{6}\frac{1}{f^{6}}.$
(A.1)

The CDSI method leads to the following values of the coefficients:

$\begin{split}B_{1}&=6{\cal Q}\frac{(L+1)(L+2){\cal
Q}^{4}+4\kappa(L^{2}+3L+1){\cal Q}^{2}-4\kappa^{2}(3L^{2}+9L+13)}{({\cal
Q}^{2}+10\kappa)^{2}},\\\ B_{2}&=3\frac{(7L^{2}+19L+14){\cal
Q}^{4}+4\kappa(3L^{2}+5L+7){\cal Q}^{2}-4\kappa^{2}(13L^{2}+29L-17)}{({\cal
Q}^{2}+10\kappa)^{2}},\\\ B_{3}&=\frac{6\kappa(2L+3){\cal Q}[4(L+1){\cal
Q}^{2}+5\kappa(2L-1)]}{({\cal Q}^{2}+10\kappa)^{2}},\\\
B_{4}&=\frac{9\kappa(2L+3)[4(L+1){\cal Q}^{2}+\kappa(2L-17)]}{({\cal
Q}^{2}+10\kappa)^{2}},\\\ B_{5}&=\frac{18\kappa^{2}(2L+3)^{2}{\cal Q}}{({\cal
Q}^{2}+10\kappa)^{2}},\\\ B_{6}&=\frac{9\kappa^{2}(2L+3)^{2}}{({\cal
Q}^{2}+10\kappa)^{2}}.\end{split}.$ (A.2)

The corresponding two superpotentials are

$\begin{split}W(r)&=-\frac{L+1}{r}f+\frac{{\cal Q}[(L+2){\cal
Q}^{2}+\kappa(4L+11)]}{{\cal Q}^{2}+10\kappa}+\frac{3\kappa[(L+1){\cal
Q}^{2}+2\kappa(L-1)]}{{\cal Q}^{2}+10\kappa}\frac{r}{f}\\\
&\quad+\frac{3\kappa(2L+3){\cal Q}}{{\cal
Q}^{2}+10\kappa}\frac{1}{f^{2}}+\frac{3\kappa^{2}(2L+3)}{{\cal
Q}^{2}+10\kappa}\frac{r}{f^{3}},\\\ W^{\prime}(r)&=-\frac{L+2}{r}f+\frac{{\cal
Q}[(L+1){\cal Q}^{2}+\kappa(4L+1)]}{{\cal
Q}^{2}+10\kappa}+\frac{3\kappa[(L+2){\cal Q}^{2}+2\kappa(L+4)]}{{\cal
Q}^{2}+10\kappa}\frac{r}{f}\\\ &\quad+\frac{3\kappa(2L+3){\cal Q}}{{\cal
Q}^{2}+10\kappa}\frac{1}{f^{2}}+\frac{3\kappa^{2}(2L+3)}{{\cal
Q}^{2}+10\kappa}\frac{r}{f^{3}}.\end{split}.$ (A.3)

The first two bound states of the potential defined in (A.1) and (A.2) are
characterized by the energies

$\begin{split}E_{0}&=16\kappa\left(\frac{(L+1){\cal Q}^{2}+\kappa(4L+1)}{{\cal
Q}^{2}+10\kappa}\right)^{2}-{\cal Q}^{2}\left(\frac{(L+2){\cal
Q}^{2}+\kappa(4L+11)}{{\cal Q}^{2}+10\kappa}\right)^{2},\\\
E_{1}&=16\kappa\left(\frac{(L+2){\cal Q}^{2}+\kappa(4L+11)}{{\cal
Q}^{2}+10\kappa}\right)^{2}-{\cal Q}^{2}\left(\frac{(L+1){\cal
Q}^{2}+\kappa(4L+1)}{{\cal Q}^{2}+10\kappa}\right)^{2},\end{split}$ (A.4)

and their wavefunctions can be written as

$\begin{split}\psi_{0}(r)&\propto r^{L+1}f^{\frac{(6L+5){\cal
Q}^{2}+2\kappa(6L-11)}{2({\cal Q}^{2}+10\kappa)}}\exp{\biggl{(}-\frac{{\cal
Q}}{\sqrt{\kappa}}\frac{(L+2){\cal Q}^{2}+\kappa(4L+11)}{{\cal
Q}^{2}+10\kappa}\biggr{)}}\\\
&\quad{}\times\exp{\biggl{(}-\frac{3\kappa(2L+3){\cal Q}}{{\cal
Q}^{2}+10\kappa}\frac{r}{f}-\frac{3\kappa(2L+3)}{2({\cal
Q}^{2}+10\kappa)}\frac{1}{f^{2}}\biggr{)}},\\\ \psi_{1}(r)&\propto
r^{L+1}f^{\frac{(6L+5){\cal Q}^{2}+2\kappa(6L-11)}{2({\cal
Q}^{2}+10\kappa)}}\exp{\biggl{(}-\frac{{\cal
Q}}{\sqrt{\kappa}}\frac{(L+1){\cal Q}^{2}+\kappa(4L+1)}{{\cal
Q}^{2}+10\kappa}\biggr{)}}\\\
&\quad{}\times\exp{\biggl{(}-\frac{3\kappa(2L+3){\cal Q}}{{\cal
Q}^{2}+10\kappa}\frac{r}{f}-\frac{3\kappa(2L+3)}{2({\cal
Q}^{2}+10\kappa)}\frac{1}{f^{2}}\biggr{)}}\\\ &\quad{}\times\biggl{[}{\cal
Q}\biggl{(}1-\frac{{\cal Q}^{2}+4\kappa}{{\cal Q}^{2}+10\kappa}\kappa
r^{2}\biggr{)}rf-1+\frac{5{\cal Q}^{2}+32\kappa}{{\cal Q}^{2}+10\kappa}\kappa
r^{2}-4\frac{{\cal Q}^{2}+4\kappa}{{\cal
Q}^{2}+10\kappa}\kappa^{2}r^{4}\biggr{]},\end{split}$ (A.5)

respectively.

## Appendix B: Solution of the system of equations (4.13)–(4.17) and
eigenvalues $E_{0}$ and $E_{1}$ for $m=4$, 5, 6, 7

The solution of the system of $m+1$ equations (4.13)–(4.17) is given by

$\displaystyle a_{0}$ $\displaystyle=\frac{({\cal Q}^{2}+\kappa)({\cal
Q}^{2}+9\kappa)}{{\cal Q}^{4}+22\kappa{\cal Q}^{2}+45\kappa^{2}},\quad
a_{1}=\frac{4{\cal Q}^{2}({\cal Q}^{2}+7\kappa)}{{\cal Q}^{4}+22\kappa{\cal
Q}^{2}+45\kappa^{2}},\quad a_{2}=\frac{12\kappa{\cal Q}({\cal
Q}^{2}+\kappa)}{{\cal Q}^{4}+22\kappa{\cal Q}^{2}+45\kappa^{2}},$
$\displaystyle a_{3}$ $\displaystyle=\frac{24\kappa{\cal Q}^{2}}{{\cal
Q}^{4}+22\kappa{\cal Q}^{2}+45\kappa^{2}},\quad a_{4}=\frac{24\kappa^{2}{\cal
Q}}{{\cal Q}^{4}+22\kappa{\cal Q}^{2}+45\kappa^{2}},$ (B.1)

if $m=4$,

$\displaystyle a_{0}$ $\displaystyle=\frac{({\cal Q}^{2}+4\kappa)({\cal
Q}^{2}+16\kappa)}{{\cal Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}},\quad
a_{1}=\frac{5({\cal Q}^{4}+16\kappa{\cal Q}^{2}+24\kappa^{2})}{{\cal
Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}},$ $\displaystyle a_{2}$
$\displaystyle=\frac{20\kappa{\cal Q}({\cal Q}^{2}+4\kappa)}{{\cal
Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}},\quad a_{3}=\frac{60\kappa({\cal
Q}^{2}+2\kappa)}{{\cal Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}},$
$\displaystyle a_{4}$ $\displaystyle=.\frac{120\kappa^{2}{\cal Q}}{{\cal
Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}},\quad
a_{5}=\frac{120\kappa^{2}}{{\cal Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}},$
(B.2)

if $m=5$,

$\displaystyle a_{0}$ $\displaystyle=\frac{({\cal Q}^{2}+\kappa)({\cal
Q}^{2}+9\kappa)({\cal Q}^{2}+25\kappa)}{{\cal Q}^{6}+65\kappa{\cal
Q}^{4}+919\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}},$ $\displaystyle a_{1}$
$\displaystyle=\frac{6{\cal Q}^{2}({\cal Q}^{4}+30\kappa{\cal
Q}^{2}+149\kappa^{2})}{{\cal Q}^{6}+65\kappa{\cal Q}^{4}+919\kappa^{2}{\cal
Q}^{2}+1575\kappa^{3}},$ $\displaystyle a_{2}$
$\displaystyle=\frac{30\kappa{\cal Q}({\cal Q}^{2}+\kappa)({\cal
Q}^{2}+9\kappa)}{{\cal Q}^{6}+65\kappa{\cal Q}^{4}+919\kappa^{2}{\cal
Q}^{2}+1575\kappa^{3}},$ $\displaystyle a_{3}$
$\displaystyle=\frac{120\kappa{\cal Q}^{2}({\cal Q}^{2}+7\kappa)}{{\cal
Q}^{6}+65\kappa{\cal Q}^{4}+919\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}},$
$\displaystyle a_{4}$ $\displaystyle=\frac{360\kappa^{2}{\cal Q}({\cal
Q}^{2}+\kappa)}{{\cal Q}^{6}+65\kappa{\cal Q}^{4}+919\kappa^{2}{\cal
Q}^{2}+1575\kappa^{3}},$ $\displaystyle a_{5}$
$\displaystyle=\frac{720\kappa^{2}{\cal Q}^{2}}{{\cal Q}^{6}+65\kappa{\cal
Q}^{4}+919\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}},$ $\displaystyle a_{6}$
$\displaystyle=\frac{720\kappa^{3}{\cal Q}}{{\cal Q}^{6}+65\kappa{\cal
Q}^{4}+919\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}},$ (B.3)

if $m=6$, and

$\displaystyle a_{0}$ $\displaystyle=\frac{({\cal Q}^{2}+4\kappa)({\cal
Q}^{2}+16\kappa)({\cal Q}^{2}+36\kappa)}{{\cal Q}^{6}+98\kappa{\cal
Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}},$ $\displaystyle a_{1}$
$\displaystyle=\frac{7({\cal Q}^{6}+50\kappa{\cal Q}^{4}+544\kappa^{2}{\cal
Q}^{2}+720\kappa^{3})}{{\cal Q}^{6}+98\kappa{\cal Q}^{4}+2464\kappa^{2}{\cal
Q}^{2}+13392\kappa^{3}},$ $\displaystyle a_{2}$
$\displaystyle=\frac{42\kappa{\cal Q}({\cal Q}^{2}+4\kappa)({\cal
Q}^{2}+16\kappa)}{{\cal Q}^{6}+98\kappa{\cal Q}^{4}+2464\kappa^{2}{\cal
Q}^{2}+13392\kappa^{3}},$ $\displaystyle a_{3}$
$\displaystyle=\frac{210\kappa({\cal Q}^{4}+16\kappa{\cal
Q}^{2}+24\kappa^{2})}{{\cal Q}^{6}+98\kappa{\cal Q}^{4}+2464\kappa^{2}{\cal
Q}^{2}+13392\kappa^{3}},$ $\displaystyle a_{4}$
$\displaystyle=\frac{840\kappa^{2}{\cal Q}({\cal Q}^{2}+4\kappa)}{{\cal
Q}^{6}+98\kappa{\cal Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}},$
$\displaystyle a_{5}$ $\displaystyle=\frac{2520\kappa^{2}({\cal
Q}^{2}+2\kappa)}{{\cal Q}^{6}+98\kappa{\cal Q}^{4}+2464\kappa^{2}{\cal
Q}^{2}+13392\kappa^{3}},$ $\displaystyle a_{6}$
$\displaystyle=\frac{5040\kappa^{3}{\cal Q}}{{\cal Q}^{6}+98\kappa{\cal
Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}},$ $\displaystyle a_{7}$
$\displaystyle=\frac{5040\kappa^{3}}{{\cal Q}^{6}+98\kappa{\cal
Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}},$ (B.4)

if $m=7$.

The corresponding eigenvalues of $V(r)$ are obtained in the form

$\displaystyle E_{0}$ $\displaystyle=25\kappa\biggl{(}\frac{(L+1){\cal
Q}^{4}+2\kappa(5L+2){\cal Q}^{2}+9\kappa^{2}(L-1)}{{\cal Q}^{4}+22\kappa{\cal
Q}^{2}+45\kappa^{2}}\biggr{)}^{2}$ $\displaystyle\quad{}-{\cal
Q}^{2}\biggl{(}\frac{(L+2){\cal Q}^{4}+2\kappa(5L+13){\cal
Q}^{2}+9\kappa^{2}(L+4)}{{\cal Q}^{4}+22\kappa{\cal
Q}^{2}+45\kappa^{2}}\biggr{)}^{2},$ $\displaystyle E_{1}$
$\displaystyle=25\kappa\biggl{(}\frac{(L+2){\cal Q}^{4}+2\kappa(5L+13){\cal
Q}^{2}+9\kappa^{2}(L+4)}{{\cal Q}^{4}+22\kappa{\cal
Q}^{2}+45\kappa^{2}}\biggr{)}^{2}$ $\displaystyle\quad{}-{\cal
Q}^{2}\biggl{(}\frac{(L+1){\cal Q}^{4}+2\kappa(5L+2){\cal
Q}^{2}+9\kappa^{2}(L-1)}{{\cal Q}^{4}+22\kappa{\cal
Q}^{2}+45\kappa^{2}}\biggr{)}^{2},$ (B.5)

if $m=4$,

$\displaystyle E_{0}$ $\displaystyle=36\kappa\biggl{(}\frac{(L+1){\cal
Q}^{4}+10\kappa(2L+1){\cal Q}^{2}+32\kappa^{2}(2L-1)}{{\cal
Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}}\biggr{)}^{2}$
$\displaystyle\quad-{\cal Q}^{2}\biggl{(}\frac{(L+2){\cal
Q}^{4}+10\kappa(2L+5){\cal Q}^{2}+32\kappa^{2}(2L+7)}{{\cal
Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}}\biggr{)}^{2},$ $\displaystyle
E_{1}$ $\displaystyle=36\kappa\biggl{(}\frac{(L+2){\cal
Q}^{4}+10\kappa(2L+5){\cal Q}^{2}+32\kappa^{2}(2L+7)}{{\cal
Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}}\biggr{)}^{2}$
$\displaystyle\quad{}-{\cal Q}^{2}\biggl{(}\frac{(L+1){\cal
Q}^{4}+10\kappa(2L+1){\cal Q}^{2}+32\kappa^{2}(2L-1)}{{\cal
Q}^{4}+40\kappa{\cal Q}^{2}+264\kappa^{2}}\biggr{)}^{2},$ (B.6)

if $m=5$,

$\displaystyle E_{0}$ $\displaystyle=49\kappa\biggl{(}\frac{(L+1){\cal
Q}^{6}+5\kappa(7L+4){\cal Q}^{4}+\kappa^{2}(259L-71){\cal
Q}^{2}+225\kappa^{3}(L-2)}{{\cal Q}^{6}+65\kappa{\cal
Q}^{4}+912\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}}\biggr{)}^{2}$
$\displaystyle\quad-{\cal Q}^{2}\biggl{(}\frac{(L+2){\cal
Q}^{6}+5\kappa(7L+17){\cal Q}^{4}+\kappa^{2}(259L+848){\cal
Q}^{2}+225\kappa^{3}(L+5)}{{\cal Q}^{6}+65\kappa{\cal
Q}^{4}+912\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}}\biggr{)}^{2},$ $\displaystyle
E_{1}$ $\displaystyle=49\kappa\biggl{(}\frac{(L+2){\cal
Q}^{6}+5\kappa(7L+17){\cal Q}^{4}+\kappa^{2}(259L+848){\cal
Q}^{2}+225\kappa^{3}(L+5)}{{\cal Q}^{6}+65\kappa{\cal
Q}^{4}+912\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}}\biggr{)}^{2}$
$\displaystyle\quad{}-{\cal Q}^{2}\biggl{(}\frac{(L+1){\cal
Q}^{6}+5\kappa(7L+4){\cal Q}^{4}+\kappa^{2}(259L-71){\cal
Q}^{2}+225\kappa^{3}(L-2)}{{\cal Q}^{6}+65\kappa{\cal
Q}^{4}+912\kappa^{2}{\cal Q}^{2}+1575\kappa^{3}}\biggr{)}^{2},$ (B.7)

if $m=6$, and

$\displaystyle E_{0}$ $\displaystyle=64\kappa\biggl{(}\frac{(L+1){\cal
Q}^{6}+7\kappa(8L+5){\cal Q}^{4}+56\kappa^{2}(14L-1){\cal
Q}^{2}+72\kappa^{3}(32L-45)}{{\cal Q}^{6}+98\kappa{\cal
Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}}\biggr{)}^{2}$
$\displaystyle\quad-{\cal Q}^{2}\biggl{(}\frac{(L+2){\cal
Q}^{6}+7\kappa(8L+19){\cal Q}^{4}+56\kappa^{2}(14L+43){\cal
Q}^{2}+72\kappa^{3}(32L+141)}{{\cal Q}^{6}+98\kappa{\cal
Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}}\biggr{)}^{2},$
$\displaystyle E_{1}$ $\displaystyle=64\kappa\biggl{(}\frac{(L+2){\cal
Q}^{6}+7\kappa(8L+19){\cal Q}^{4}+56\kappa^{2}(14L+43){\cal
Q}^{2}+72\kappa^{3}(32L+141)}{{\cal Q}^{6}+98\kappa{\cal
Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}}\biggr{)}^{2}$
$\displaystyle\quad-{\cal Q}^{2}\biggl{(}\frac{(L+1){\cal
Q}^{6}+7\kappa(8L+5){\cal Q}^{4}+56\kappa^{2}(14L-1){\cal
Q}^{2}+72\kappa^{3}(32L-45)}{{\cal Q}^{6}+98\kappa{\cal
Q}^{4}+2464\kappa^{2}{\cal Q}^{2}+13392\kappa^{3}}\biggr{)}^{2},$ (B.8)

if $m=7$.

## References

  * [1] P.M. Mathews, M. Lakshmanan, Q. Appl. Math. 32 (1974) 215.
  * [2] J.F. Cariñena, M.F. Rañada, M. Santander, Rep. Math. Phys. 54 (2004) 285.
  * [3] J.F. Cariñena, M.F. Rañada, M. Santander, Ann. Phys., NY 322 (2007) 434.
  * [4] J.F. Cariñena, M.F. Rañada, M. Santander, Ann. Phys., NY 322 (2007) 2249.
  * [5] J.F. Cariñena, M.F. Rañada, M. Santander, J. Math. Phys. 48 (2007) 102106.
  * [6] J.F. Cariñena, M.F. Rañada, M. Santander, J. Phys. A 45 (2012) 265303.
  * [7] A. Schulze-Halberg, J.R. Morris, J. Phys. A 45 (2012) 305301.
  * [8] C. Quesne, Phys. Lett. A 379 (2015) 1589.
  * [9] C. Quesne, J. Math. Phys. 57 (2016) 102101.
  * [10] B. Bagchi, A. Banerjee, C. Quesne, V.M. Tkachuk, J. Phys. A 38 (2005) 2929.
  * [11] D. Gómez-Ullate, Y. Grandati, R. Milson, J. Math. Phys. 55 (2014) 043510.
  * [12] A.V. Turbiner, Phys. Rep. 642 (2016) 1.
  * [13] A. Ronveaux, Heun Differential Equations, Oxford University Press, Oxford, 1995.
  * [14] C. Quesne, J. Math. Phys. 58 (2017) 052104.
  * [15] Y.-Z. Zhang, J. Phys. A 45 (2012) 065206.
  * [16] C. Quesne, J. Math. Phys. 59 (2018) 042104.
  * [17] B. Chakrabarti, J. Phys. A 41 (2008) 405301.
  * [18] S. Bera, B. Chakrabarti, T.K. Das, Phys. Lett. A 381 (2017) 1356.
  * [19] V.M. Tkachuk, Phys. Lett. A 245 (1998) 177.
  * [20] O. Voznyak, V.M. Tkachuk, J. Phys. Stud. 16 (2012) 1003 (in Ukrainian).
  * [21] E. Schrödinger, Proc. R. Ir. Acad., Sec. A 46 (1940) 9.
  * [22] L. Infeld, Phys. Rev. 59 (1941) 737.
  * [23] A.F. Stevenson, Phys. Rev. 59 (1941) 842.
  * [24] E.G. Kalnins, W. Miller, Jr., P. Winternitz, SIAM J. Appl. Math. 30 (1976) 630.
  * [25] P.W. Higgs, J. Phys. A 12 (1979) 309.
  * [26] H.I. Leemon, J. Phys. A 12 (1979) 489.
  * [27] J.F. Cariñena, M.F. Rañada, M. Santander, J. Math. Phys. 46 (2005) 052702.
  * [28] A. Ballesteros, F.J. Herranz, J. Phys. A 42 (2009) 245203.
  * [29] M.F. Rañada, Phys. Lett. A 380 (2016) 2204.
  * [30] T. Hakobyan, A. Nersessian, H. Shmavonyan, Phys. Rev. D 95 (2017) 025014.
  * [31] J.F. Cariñena, M.F. Rañada, M. Santander, J. Phys. A 54 (2021) 105201.
  * [32] J.F. Cariñena, M.F. Rañada, M. Santander, J. Phys. A 54 (2021) 365201.
  * [33] J.F. Cariñena, M.F. Rañada, M. Santander, M. Senthilvelan, Nonlinearity 17 (2004) 1941.
  * [34] C. Quesne, V.M. Tkachuk, J. Phys. A 37 (2004) 4267.
  * [35] O. Mustafa, S.H. Mazharimousavi, Int. J. Theor. Phys. 46 (2007) 1786.
  * [36] V. Romanovski, C. R. Acad. Sci. Paris 188 (1929) 1023.
  * [37] A.P. Raposo, H.J. Weber, D. Alvarez-Castillo, M. Kirchbach, Cent. Eur. J. Phys. 5 (2007) 253.
  * [38] R. Koekoek, P.A. Lesky, R.F. Swarttouw, Hypergeometric Orthogonal Polynomials and Their $q$-Analogues, Springer, Berlin, 2010.


11institutetext: Itec, imec research group at KU Leuven, Kortrijk, Belgium
22institutetext: KU Leuven, Campus KULAK, Department of Public Health and
Primary Care, Kortrijk, Belgium 33institutetext: KU Leuven, Department of
Management, Strategy and Innovation, Leuven, Belgium 33email:
<EMAIL_ADDRESS>
# HypeRS: Building a Hypergraph-driven ensemble Recommender System

Alireza Gharahighehi 1122 Celine Vens 1122 Konstantinos Pliakos 33

###### Abstract

Recommender systems are designed to predict user preferences over collections
of items. These systems process users’ previous interactions to decide which
items should be ranked higher to satisfy their desires. An ensemble
recommender system can achieve great recommendation performance by effectively
combining the decisions generated by individual models. In this paper, we
propose a novel ensemble recommender system that combines predictions made by
different models into a unified hypergraph ranking framework. This is the
first time that hypergraph ranking has been employed to model an ensemble of
recommender systems. Hypergraphs are generalizations of graphs where multiple
vertices can be connected via hyperedges, efficiently modeling high-order
relations. We differentiate real and predicted connections between users and
items by assigning different hyperedge weights to individual recommender
systems. We perform experiments using four datasets from the fields of movie,
music and news media recommendation. The obtained results show that the
ensemble hypergraph ranking method generates more accurate recommendations
compared to the individual models and a weighted hybrid approach. The
assignment of different hyperedge weights to the ensemble hypergraph further
improves the performance compared to a setting with identical hyperedge
weights.

###### Keywords:

Recommender systems Hypergraph learning Ensemble methods.

##  1 Introduction

Nowadays, people use digital services more and more to fulfill their needs.
The owners of these services monitor users’ behavior and utilize users’
interactions with provided items, such as movies, songs, commercial products,
to predict users’ preferences. This enables the personalization of digital
services and the rise of effective recommender systems (RSs) which learn from
users’ preferences and provide them with accurate recommendations. Generally,
there are two main categories in RSs: content-based filtering and
collaborative filtering approaches. Content-based RSs use the features that
describe the items for computing similarities between the items and the user
interaction profile. Next, they recommend items that are more similar to this
user profile. Upon a recommendation query for a target user, these RSs do not
consider the interactions of the other users in generating the recommendation
list. In contrast to that, collaborative filtering approaches infer the users’
preferences by processing the collaborative information between users or
items. In many applications, collaborative filtering RSs generate more
accurate [1] and less obvious [17] recommendations compared to content-based
approaches.

Each type of RS processes the information based on different assumptions to
decide which items should be ranked higher among many available ones. For
instance, memory-based collaborative filtering approaches (user-based and
item-based) assume that users (items) with similar interactions have similar
interests. Therefore, these approaches form neighborhoods to generate
recommendations. Model-based collaborative filtering approaches assume that
users and items can be represented in a common feature space and they use
different learning methods to learn these latent features. While these
approaches might vary in prediction power, they convey relevant information
from different perspectives, following practically different learning
strategies for the same recommendation task. Ensemble methods include multiple
learning methods and integrate their predictive power into a single system,
achieving superior predictive performance to individual models. Examples of
ensembles in machine learning are bagging and boosting. In recommendation
tasks a hybrid RS can be applied to exploit several data sources or the
prediction power of different RSs to generate more relevant recommendations.
An ensemble RS is a hybrid model that employs the ranking lists of multiple
RSs to decide which items should be recommended to each user [2].

In this paper we propose an ensemble hypergraph learning framework for
recommendation. This way we integrate the predictive power of several models
into a unified RS powered by hypergraph ranking. Unlike regular graphs, where
edges connect pairs of nodes, in hypergraphs multiple nodes can be connected
via hyperedges. These higher order relations in hyperedges empower hypergraphs
to cast more reliable information in the model [23]. Furthermore, hypergraph
learning can inherently model the complex relations between different types of
entities in a unified framework. It is therefore a deliberate choice for the
construction of an ensemble of individual RSs driven by different types of
information. Moreover, as was shown in [13], hypergraph ranking-based methods
can mitigate popularity bias, enhance fairness and coverage as well as act as
innate multi-stakeholder RSs. The main contribution of this paper is to
construct a hypergraph as an ensemble framework for recommendation tasks.
Despite its capability to stack multiple connections in a unified model, to
the best of our knowledge hypergraphs have not been employed to form ensembles
of RSs.

This article is an extension of our previous conference paper [12] in the
following main directions:

  * •

We have extended the hypergraph-based ensemble RS by further incorporating two
effective individual RSs, namely Weighted Approximate-Rank Pairwise (WARP)
[25] and Multi-variational autoencoder (MVAE) [16]. This way, we managed to
substantially enhance the performance of our RS. We also demonstrated the
potential of our model, as it has the capacity to integrate more effective
individual RSs that might appear in the future, boosting its recommendation
performance even further.

  * •

In the previous paper we considered equal weights for all the hyperedges in
the ensemble hypergraph. Here, we have differentiated the actual links between
users and items from the ones predicted by the individual RSs. We also
assigned different weights to hyperedges that relate to different individual
RSs based on the performance of these RSs.

  * •

We have extended the experimental study by including more evaluation metrics
as well as by comparing our method to two additional baseline models, namely
Weighted Approximate-Rank Pairwise (WARP) [25] and Multi-variational
autoencoder (MVAE) [16]. Especially the MVAE is a state-of-the-art method that
emerged among the winning models in a recent comparison study of RSs [7].

The structure of this paper is as follows: Studies about applications of
hypergraph learning in RSs are presented in Section 2. Next, in Section 3, we
show how a unified hypergraph can be formed as an RS (Section 3.1) and how it
can formulate an ensemble of RSs (Section 3.2). In Section 4, four
recommendation datasets are described and the experimental setup in designing
and testing the proposed model is described. Next, the obtained results of
comparing the proposed ensemble model against other methods on these four
datasets are presented and discussed in Section 5. Finally, we draw
conclusions and outline some directions for future research in Section 6.

##  2 Related work

Hypergraph learning has been applied to generate recommendation lists in
several applications. For instance in the music domain, Bu et al. [3] used
hypergraph learning to recommend music tracks where the relations between
users, tracks, albums and artists were modeled using a unified hypergraph.
Hypergraph ranking has been also used in news recommendation tasks [15, 13].
News usually contains very rich features such as text, tags and named
entities. Therefore, hypergraph learning can effectively model the relations
between these entities. Moreover, Pliakos et al. [21] used hypergraph ranking
for a tag recommendation task. They built a hypergraph ranking model to
capture the complex relations between different entities in the system, such
as users, images, tags, and geo-tags. Hypergraph-based RSs have been also used
in e-commerce applications [18, 24]. For instance in [18], a multipartite
hypergraph is used to model the relations between users, restaurants and
attributes in a multi-objective setting. In such applications, item attributes
and sequences of user-item interactions are effectively modeled in
hypergraphs.

Hypergraph learning has been employed to address various issues in RSs. A
hypergraph can model the relations between different types of stakeholders and
objects and therefore, it can be intrinsically used as a multi-stakeholder RS
[13]. Additionally, it can be used to burst the filter bubble around the user
by querying a more diverse recommendation list based on the user history [15,
13]. Moreover, hypergraph learning has been used to address fairness [13], the
cold-start problem [27] as well as context-awareness [26] in recommendation
tasks.

An ensemble RS is a type of hybrid RSs that integrates the recommendations of
multiple individual RSs. Aggarwal [2] categorized hybrid RSs to monolithic,
ensembles, and mixed RSs. Burke et al. [4] provided another categorization
where hybrid models are categorized into weighted, switching, cascade, feature
augmentation, feature combination, meta-level and mixed RSs. A weighted hybrid
RS uses the weighted average of the scores from individual RSs to generate the
recommendation list. For instance, Do et al. [6] applied a weighted hybrid RS
based on collaborative and content-based filtering approaches on Movielens
dataset and showed that it is more effective compared to the individual
collaborative and content-based RSs. Here, we employ a unified hypergraph as
an ensemble RS. Although hypergraph learning is very promising and effective
in addressing many problems in RSs, to the best of our knowledge, it has never
been studied as an ensemble RS.

##  3 Methodology

###  3.1 Hypergraphs as recommender systems

Hereafter, uppercase bold letters are used for matrices, lowercase bold
letters represent vectors, uppercase non-bold letters are used for sets and
lowercase non-bold letters represent constants. The element in $i^{th}$ row
and $j^{th}$ column of matrix X is denoted as X(i,j).

A hypergraph consists of a set of nodes (vertices)
$N:\\{n_{1},n_{2}\cdots,n_{|N|}\\}$ and a set of hyperedges
$E:\\{e_{1},e_{2}\cdots,e_{|E|}\\}$ that connect the nodes. Each hyperedge can
connect multiple nodes in the hypergraph. Based on the application, different
types of hyperedges can be defined that capture different forms/sources of
information. We define these hyperedge types in Section 3.2. In a typical
collaborative filtering setting there are two types of entities in a
hypergraph: users $U:\\{u_{1},u_{2}\cdots,u_{|U|}\\}$ and items
$I:\\{i_{1},i_{2}\cdots,i_{|I|}\\}$. Therefore, the set of nodes $N$ in a
hypergraph is formed based on users and items ($N:\\{U\cup I\\}$).

Let $\mathbf{H}$ of size $|N|\times|E|$ be the incidence matrix of the
hypergraph, where $H(n,e)=1$, if node $n$ is in hyperedge $e$ and zero
otherwise. Based on $\mathbf{H}$, the symmetric matrix $\mathbf{A}$ can be
formed using Eq.1:

$\mathbf{A}=\mathbf{D_{n}}^{-1/2}\mathbf{H}\mathbf{W}\mathbf{D_{e}}^{-1}\mathbf{H}^{T}\mathbf{D_{n}}^{-1/2}$
(1)

where $\mathbf{D}_{n}$ and $\mathbf{D}_{e}$ are the diagonal matrices that
contain the node and hyperedge degrees and $\mathbf{W}$ is the diagonal
hyperedge weight matrix. Each element $\mathbf{A}(i,j)$ reflects the
relatedness between nodes $i$ and $j$. Higher values indicate stronger
relations between the corresponding nodes. Then, the recommendation problem is
formulated as finding a ranking (score) vector $\mathbf{f}\in{\rm
I\\!R}^{|N|}$ that minimizes the following loss function [3]:

$Q(\mathbf{f})=\frac{1}{2}\mathbf{f}^{T}\mathbf{L}\mathbf{f}+\vartheta||\mathbf{f}-\mathbf{y}||_{2}^{2}$
(2)

where $\mathbf{L}$ is the hypergraph Laplacian matrix (i.e.
$\mathbf{L}=\mathbf{I}-\mathbf{A}$), $\vartheta$ is a regularizing parameter
and $\mathbf{y}\in{\rm I\\!R}^{|N|}$ is the query vector. Every item of the
ranking vector $\mathbf{f}$ or query vector $\mathbf{y}$ corresponds to a
node. Typically, to generate the recommendation list for user $u$ in a regular
recommendation task, one can query the hypergraph for user $u$ by setting the
corresponding value in the query vector to one ($\mathbf{y}(u)=1$) and all the
other values that correspond to other nodes to zero. By solving the
optimization problem in Eq.2, the optimal score (ranking) vector can be
calculated using Eq.3:

${\mathbf{f}}^{\ast}=\frac{\vartheta}{1+\vartheta}\bigl{(}\mathbf{I}-\frac{1}{1+\vartheta}\mathbf{A}\bigr{)}^{-1}\mathbf{y}.$
(3)

Finally, the top k items that have the highest scores in $\mathbf{{f}^{\ast}}$
are recommended to the user $u$.

###  3.2 A hypergraph-based ensemble recommender system, HypeRS

An ensemble RS111The source code is available at
https://github.com/alirezagharahi/ensemble_hypergraph. utilizes the decisions
of multiple individual RSs to decide which items should be ranked higher in
the final recommendation lists. Let $M:{\\{m_{1},m_{2}\cdots,m_{|M|}\\}}$ be
the set of individual methods that we want to incorporate in our ensemble RS.
Each of these individual methods $m_{i}$ can generate its own top $k$ rankings
$\mathbf{R}_{i}\in{\rm I\\!R}^{|U|\times k}$ where each row in
$\mathbf{R}_{i}$ is the top $k$ ranked items for the corresponding user. Then,
based on the recommendation lists of each RS, hyperedges are formed to connect
users to their top $k$ recommendations.

As is mentioned previously, the hypergraph consists of multiple types of
hyperedges. We consider three types of hyperedges, which are defined in Table
1. The $E_{UI}$ hyperedges connect the users with the items that they have
interacted with. To make the relations between users with similar tastes more
explicit, the $E_{UU}$ hyperedges connect users to their k nearest neighbors.
To find these neighbors we use the user-item interaction matrix $\mathbf{Z}$,
where $\mathbf{Z}(i,j)\in\\{0,1\\}$. The k nearest neighbors of user u are
users that have the highest cosine similarity with $u^{th}$ row of matrix
$\mathbf{Z}$. The $E_{M}$ hyperedges are considered to integrate the
recommendations of multiple RSs in the hypergraph. These RSs can be from
different families such as collaborative filtering or content-based
approaches. The fact that recommendations from any type of RS can be directly
modeled as hyperedges in our system is a vital advantage of the proposed
method.

We constructed the $E_{M}$ hyperedge set using four well-established and
powerful matrix completion-based recommendation methods, namely BPR [22], WARP
[25], WRMF [14, 20] and MVAE [16]. BPR is a learning-to-rank matrix completion
approach which uses user-specific relative pair-wise preferences between
observed and unobserved items to learn items’ and users’ low rank matrices.
Similar to BPR, WARP is also a pair-wise learning-to-rank approach but with a
different objective function. BPR is optimized for approximation of AUC
whereas WARP is optimized for precision. WRMF is a matrix factorization
approach for implicit feedback datasets that uses the alternating-least-
squares optimization process to learn items’ and users’ parameters. MVAE is a
CF model for implicit feedback that assumes that user logs are from a
multinomial distribution and uses variational autoencoders to learn users’ and
items’ parameters.

Table 1: Hyperedge definitions hyperedge |  definition | # of hyperedges
---|---|---
$E_{UI}$ |  Each user is connected to the items that the user has interacted with | $|U|$
$E_{UU}$ |  Each user is connected to the k most similar users | $|U|$
$E_{M}$ |  Each user is connected to top k recommended items by a RS | $|M|\times|U|$

The hypergraph and its incidence matrix $\mathbf{H}$ are constructed using the
hyperedge sets of Table 1. Following that, the affinity matrix $\mathbf{A}$ is
computed and the recommendation task is addressed as was described in Section
3.1.

The weight of a hyperedge reflects the relative importance of the hyperedge
compared to the other hyperedges in the hypergraph. We propose to assign lower
weights for $E_{M}$ hyperedges (0.5) compared to $E_{UI}$ ones (1.0), as the
$E_{M}$ are based on predicted links, whereas the $E_{UI}$ are based on real
links between users and items. One should also consider that the individual
RSs in the ensemble have different predictive power. Therefore we assign a
decay weight based on their performance in the validation set. Hyperedges that
are associated with the top ranked RS receive no decay weight while hyperedges
that are related to lower ranked RSs receive lower weights (according to their
ranks). In this paper we consider a linear decay (10%) per rank.

##  4 Experimental Setup

To evaluate the performance of the proposed approach we use four datasets from
news, music and movie application domains. These datasets are described in
Table 2. AOTM is a publicly available dataset collected from the Art-of-the-
Mix platform that is based on user playlists [19].
Movielens222http://www.grouplens.org is a publicly available movie rating
dataset [5]. As we only encode interactions in the hypergraph for this dataset
we transform ratings to binary feedback. Globo333http://www.globo.com and
Roularta444http://www.roularta.be are news datasets that contain readers’
interactions with news articles.

Table 2: Datasets descriptions | AOTM | Movielens | Globo | Roularta
---|---|---|---|---
item type | music track | movie | news article | news article
# users | 1,605 | 1,573 | 3,903 | 5.082
# items | 2,199 | 2,053 | 1,246 | 2,739
sparsity | 3.8% | 19.9% | 5.7% | 8.5%

In our experiments we consider the following eight approaches555For BPR and
WRMF we used implicit library.:

  * •

BPR: Bayesian Personalized Ranking (BPR) [22] is a pair-wise learning-to-rank
matrix completion approach as presented in the previous section.

  * •

WARP: Weighted Approximate-Rank Pairwise (WARP) [25] is a another pair-wise
learning-to-rank method as explained in the previous section.

  * •

WRMF: Weighted Regularized Matrix Factorization (WRMF) [14, 20] is a MF
approach using the alternating-least-squares optimization process to learn
items and users’ parameters as presented in the previous section.

  * •

MVAE: Multi-variational autoencoder (MVAE) [16] is a CF approach based on
variational autoencoders as explained in the previous section.

  * •

Hybrid: A weighted hybrid model that uses scores of BPR, WARP, MVAE and WRMF,
and then considers the weighted average of these scores to generate the final
ranking lists.

  * •

H: A hypergraph-based RS explained in Section 3.1 that only contains the
hyperedge types of $E_{UI}$ and $E_{UU}$ from Table 1.

  * •

$\mathbf{HypeRS}$: The proposed hypergraph-based ensemble RS explained in
Section 3.2.

  * •

$\mathbf{HypeRS_{W}}$: The proposed hypergraph-based ensemble RS with proposed
hyperedge weights.

To validate the performance of the proposed method against the compared
methods we randomly hide ten interactions of each user from training and then
measure the ability of the methods in predicting these hidden interactions
666Users with few interactions are omitted from experiments.. We use
precision@10, recall@10 and F1-score@10 to measure the accuracy of
predictions. Precision and Recall are standard information retrieval accuracy
measures that reflect the proportion of the recommendation list which is
relevant and the proportion of relevant items that is recommended
respectively. To have a balanced measure, one can calculate F1-score that
combines the precision and recall into a single measure by taking their
harmonic mean. These accuracy measures can be calculated as follows:

$Precision=\frac{1}{|U|}\sum_{u\in U}\frac{|T_{u}\cap R_{u}|}{|R_{u}|},$ (4)
$Recall=\frac{1}{|U|}\sum_{u\in U}\frac{|T_{u}\cap R_{u}|}{|T_{u}|}.$ (5)

$F1-score=2\times\frac{precision\times recall}{precision+recall}.$ (6)

where $U$ is the set of users, $T_{u}$ is the test items for user $u$, and
$R_{u}$ is the recommendation list for user $u$. As the number of relevant
items and the length of the recommendation lists are the same (10 items), the
reported results from precision, recall and F1-score are eventually all the
same.

The compared methods have some hyperparameters to be tuned. BPR and WARP have
number of latent features, number of iterations, regularizing parameter and
learning rate, WRMF has number of latent features, number of iterations,
regularizing parameter, MVAE has batch size, number of iterations and number
of anneal steps, Hybrid model has hybridization weights and H as well as
$HypeRS$ have a regularizer as a hyperparameter. To tune these hyperparameters
we form a validation set for each dataset by randomly drawing five
interactions of each user from the training set as the validation set. The
final tuned hyperparameter values are based on precision@10 and are reported
in Table 3.

Table 3: Hyperparameters |  | range | AOTM | Movielens | Globo | Roularta
---|---|---|---|---|---|---
BPR | # iterations | [1000,2000] | 1645 | 1984 | 1598 | 1984
# latent features | [100,250] | 179 | 129 | 168 | 129
regularizing parameter | [0.01,0.05] | 0.0194 | 0.0412 | 0.0374 | 0.0412
learning rate | [0.001,0.07] | 0.0284 | 0.0092 | 0.0174 | 0.0092
WARP | # iterations | [200, 850] | 809 | 810 | 810 | 650
# latent features | [15,40] | 21 | 21 | 21 | 26
regularizing parameter | [10e-3,10e-6] | 4.7e-06 | 3.4e-6 | 8.1e-6 | 7.4e-6
learning rate | [0.001,0.1] | 0.0472 | 0.0165 | 0.0191 | 0.0191
WRMF | # iterations | [1000,2000] | 1276 | 1393 | 29 | 1288
# latent features | [100,250] | 201 | 107 | 493 | 109
regularizing parameter | [0.01,0.05] | 0.0374 | 0.0225 | 0.0432 | 0.0315
MVAE | # iterations | [10,250] | 29 | 18 | 29 | 28
batch size | [25,500] | 469 | 34 | 152 | 127
# anneal steps | [100000,300000] | 127692 | 100212 | 244065 | 223544
Hybrid | BPR_weight | [0.01,0.99] | 0.181 | 0.152 | 0.443 | 0.285
WARP_weight | [0.01,0.99] | 0.422 | 0.193 | 0.290 | 0.184
WRMF_weight | [0.01,0.99] | 0.0588 | 0.490 | 0.458 | 0.168
MVAE_weight | [0.01,0.99] | 0.329 | 0.255 | 0.469 | 0.427
H | regularizing parameter | [0.01,0.99] | 0.2414 | 0.2414 | 0.0656 | 0.0616
$HypeRS$ | regularizing parameter | [0.01,0.99] | 0.4554 | 0.4554 | 0.8301 | 0.6325

##  5 Results and Discussion

The results of the proposed hypergraph-based ensemble RS and the selected
approaches on the four datasets are reported in Table 4. The reported values
are in terms of average precision@10 of the recommendation lists generated by
the compared approaches. As is shown in Table 4, the proposed hypergraph-based
ensemble RS ($HypeRS$) has superior predictive performance compared to all the
competitor approaches including the weighted hybrid model in all datasets. The
competitor methods have different performance rankings in the four datasets.
Each of these methods processes the information based on different assumptions
and learning approaches. The effectiveness of these assumptions and learning
approaches differs across different applications. An ensemble RS exploits the
combined predictive power of the individual methods. It considers all
assumptions and decisions of various independent RSs and achieves overall
superior performance regardless of the application domain of the
recommendation task. The weighted ensemble RS ($HypeRS_{W}$) performs better
compared to the ensemble RS with identical hyperedge weights ($HypeRS$) in all
datasets.

Table 4: Results w.r.t. precision@10/recall@10/F1-score@10 | AOTM | Movielens | Globo | Roularta
---|---|---|---|---
BPR | 0.0373 | 0.1716 | 0.0937 | 0.0704
WRMF | 0.0402 | 0.1718 | 0.0921 | 0.0918
MVAE | 0.0426 | 0.1690 | 0.1487 | 0.0837
WARP | 0.0360 | 0.1870 | 0.1375 | 0.0907
H | 0.0338 | 0.1503 | 0.1125 | 0.0657
Hybrid | 0.0412 | 0.1861 | 0.1012 | 0.0944
$HypeRS$ | 0.0479 | 0.1957 | 0.1499 | 0.0950
$HypeRS_{W}$ | 0.0483 | 0.1995 | 0.1504 | 0.0955

In this study we keep the experiments simple by only using the collaborative
information, i.e. user-item interactions, to make them applicable on available
datasets and various application fields (i.e. movies, music, news).
Nevertheless, in cases where side information is available for users or items,
content-based approaches can be included in the ensemble RS. Hypergraph
learning has the natural capability of modeling the complex relations between
different types of entities in a unified hypergraph and therefore is a
deliberate choice to construct an ensemble of RSs with different types of
information.

##  6 Conclusion

We proposed a new ensemble hypergraph learning-based RS. A unified hypergraph
can integrate multiple connections between entities (here users and items) and
therefore can combine the predictive power of various individual RSs boosting
the precision of final recommendation lists. We empirically tested this method
on four datasets from different application domains, such as news, music, and
movies. The obtained results showed that the hypergraph-based ensemble RS
achieves superior performance compared to all the individual models, as well
as compared to a weighted hybrid approach that averages individual scores to
produce final rankings, in all datasets.

For future work we outline the following directions:

  * •

Beyond accuracy evaluation: In this paper we only used user-item interactions.
Future approaches could include additional information and relevant
stakeholders so that fairness [13] and diversity [10] are also taken into
account.

  * •

Consumption level: We only captured the binary feedback between users and
items. In real applications usually the user feedback is graded [8] which
shows to what extend the user is interested in the item. This graded feedback
could be reflected in the hypergraph to model user preferences more precisely.

  * •

Long-term vs short-term preferences: In some applications such as news [9] and
music [11] recommendation tasks, users’ short-term preferences play important
roles. Session-based RSs have been used to model such user short-term
preferences. An ensemble RS could include models for both long-term and short-
term preferences.

#### Acknowledgments

This work was executed within the imec.icon project NewsButler, a research
project bringing together academic researchers (KU Leuven, VUB) and industry
partners (Roularta Media Group, Bothrs, ML6). The NewsButler project is co-
financed by imec and receives project support from Flanders Innovation &
Entrepreneurship (project nr. HBC.2017.0628). The authors also acknowledge
support from the Flemish Government (AI Research Program).

## References

  * [1] Adomavicius, G., Tuzhilin, A.: Toward the next generation of recommender systems: A survey of the state-of-the-art and possible extensions. IEEE transactions on knowledge and data engineering 17(6), 734–749 (2005)
  * [2] Aggarwal, C.C.: Ensemble-based and hybrid recommender systems. In: Recommender Systems: The Textbook. pp. 199–224. Springer International Publishing, Cham (2016)
  * [3] Bu, J., Tan, S., Chen, C., Wang, C., Wu, H., Zhang, L., He, X.: Music recommendation by unified hypergraph: combining social media information and music content. In: Proceedings of the 18th ACM international conference on Multimedia. pp. 391–400 (2010)
  * [4] Burke, R.: Hybrid recommender systems: Survey and experiments. User modeling and user-adapted interaction 12(4), 331–370 (2002)
  * [5] Cantador, I., Brusilovsky, P., Kuflik, T.: 2nd workshop on information heterogeneity and fusion in recommender systems (hetrec 2011). In: Proceedings of the 5th ACM conference on Recommender systems. RecSys 2011, ACM, New York, NY, USA (2011)
  * [6] Do, H., Le, T., Yoon, B.: Dynamic weighted hybrid recommender systems. In: 2020 22nd International Conference on Advanced Communication Technology (ICACT). pp. 644–650 (2020). https://doi.org/10.23919/ICACT48636.2020.9061465
  * [7] Ferrari Dacrema, M., Cremonesi, P., Jannach, D.: Are we really making much progress? a worrying analysis of recent neural recommendation approaches. In: Proceedings of the 13th ACM Conference on Recommender Systems. p. 101–109. RecSys ’19, ACM, New York, NY, USA (2019). https://doi.org/10.1145/3298689.3347058
  * [8] Gharahighehi, A., Vens, C.: Extended bayesian personalized ranking based on consumption behavior. In: Artificial Intelligence and Machine Learning, pp. 152–164. Springer (2020)
  * [9] Gharahighehi, A., Vens, C.: Making session-based news recommenders diversity-aware. In: Proceedings of the Workshop on Online Misinformation- and Harm-Aware Recommender Systems (2020)
  * [10] Gharahighehi, A., Vens, C.: Diversification in session-based news recommender systems. Personal and Ubiquitous Computing (2021). https://doi.org/10.1007/s00779-021-01606-4
  * [11] Gharahighehi, A., Vens, C.: Personalizing diversity versus accuracy in session-based recommender systems. SN Computer Science 2 (2021). https://doi.org/10.1007/s42979-020-00399-2
  * [12] Gharahighehi, A., Vens, C., Pliakos, K.: An ensemble hypergraph learning framework for recommendation. In: International Conference on Discovery Science. pp. 295–304. Springer (2021)
  * [13] Gharahighehi, A., Vens, C., Pliakos, K.: Fair multi-stakeholder news recommender system with hypergraph ranking. Information Processing and Management 58(5), 102663 (2021). https://doi.org/10.1016/j.ipm.2021.102663
  * [14] Hu, Y., Koren, Y., Volinsky, C.: Collaborative filtering for implicit feedback datasets. In: 2008 Eighth IEEE International Conference on Data Mining. pp. 263–272. Ieee (2008)
  * [15] Li, L., Li, T.: News recommendation via hypergraph learning: encapsulation of user behavior and news content. In: Proceedings of the sixth ACM international conference on Web search and data mining. pp. 305–314 (2013)
  * [16] Liang, D., Krishnan, R.G., Hoffman, M.D., Jebara, T.: Variational autoencoders for collaborative filtering. In: Proceedings of the 2018 world wide web conference. pp. 689–698 (2018)
  * [17] Lops, P., de Gemmis, M., Semeraro, G.: Content-based Recommender Systems: State of the Art and Trends, pp. 73–105. Springer US, Boston, MA (2011).
  * [18] Mao, M., Lu, J., Han, J., Zhang, G.: Multiobjective e-commerce recommendations based on hypergraph ranking. Information Sciences 471, 269 – 287 (2019). https://doi.org/10.1016/j.ins.2018.07.029
  * [19] McFee, B., Lanckriet, G.R.: Hypergraph models of playlist dialects. In: ISMIR. vol. 12, pp. 343–348. Citeseer (2012)
  * [20] Pan, R., Zhou, Y., Cao, B., Liu, N.N., Lukose, R., Scholz, M., Yang, Q.: One-class collaborative filtering. In: 2008 Eighth IEEE International Conference on Data Mining. pp. 502–511. IEEE (2008)
  * [21] Pliakos, K., Kotropoulos, C.: Simultaneous image tagging and geo-location prediction within hypergraph ranking framework. In: 2014 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). pp. 6894–6898. IEEE (2014)
  * [22] Rendle, S., Freudenthaler, C., Gantner, Z., Schmidt-Thieme, L.: Bpr: Bayesian personalized ranking from implicit feedback. arXiv preprint arXiv:1205.2618 (2012)
  * [23] Vinayak, R.K., Hassibi, B.: Crowdsourced clustering: Querying edges vs triangles. In: Advances in Neural Information Processing Systems. pp. 1316–1324 (2016)
  * [24] Wang, J., Ding, K., Hong, L., Liu, H., Caverlee, J.: Next-item recommendation with sequential hypergraphs. In: Proceedings of the 43rd International ACM SIGIR Conference on Research and Development in Information Retrieval. p. 1101–1110. SIGIR ’20, Association for Computing Machinery, New York, NY, USA (2020). https://doi.org/10.1145/3397271.3401133
  * [25] Weston, J., Bengio, S., Usunier, N.: Wsabie: Scaling up to large vocabulary image annotation. In: Twenty-Second International Joint Conference on Artificial Intelligence (2011)
  * [26] Yu, C.A., Tai, C.L., Chan, T.S., Yang, Y.H.: Modeling multi-way relations with hypergraph embedding. In: Proceedings of the 27th ACM International Conference on Information and Knowledge Management. pp. 1707–1710 (2018)
  * [27] Zheng, X., Luo, Y., Sun, L., Ding, X., Zhang, J.: A novel social network hybrid recommender system based on hypergraph topologic structure. World Wide Web 21(4), 985–1013 (2018)