††thanks: ORCiD: 0000-0002-4509-4309††thanks: ORCiD:
0000-0001-8114-0061††thanks: ORCiD: 0000-0002-6286-0722††thanks: ORCiD:
0000-0001-6841-9850

# Improvement of the Envelope Theory for Systems with Different Particles

Cyrille Chevalier<EMAIL_ADDRESS>Cintia T. Willemyns
<EMAIL_ADDRESS>Lorenzo Cimino<EMAIL_ADDRESS>Claude
Semay<EMAIL_ADDRESS>Service de Physique Nucléaire et Subnucléaire,
Université de Mons, UMONS Research Institute for Complex Systems, Place du
Parc 20, 7000 Mons, Belgium

(August 27, 2024)

###### Abstract

Abstract The envelope theory is a method to compute approximate eigensolutions
of quantum $N$-body Hamiltonians with a quite general structure in $D$
dimensions. The advantages of the method are that it is easy to implement and
that $N$ is treated as any other parameters of the Hamiltonian, allowing the
computation for systems of all sizes. If solutions are reliable, they are
generally not very accurate. In the case of systems with identical particles
for $D\geq 2$, it is possible to improve the precision of the eigenvalues by
combining the envelope theory with a generalisation to $N$-body of the
dominantly orbital state method. It is shown that a similar improvement can be
achieved in the case of systems composed of identical particles plus a
different one. The quality of the new procedure is tested with different
systems.

Bound states and Many-body systems and Envelope theory and Dominantly orbital
state method

##  I introduction

Among all the methods to solve quantum many-body problems, the envelope theory
(ET) hall80 ; hall83 ; hall04 , also known as the auxiliary field method
silv10 ; silv12 , occupies a special place, because the computation cost is
independent of the number of particles. It is quite easy to implement if all
the particles are identical sema13  and becomes a little trickier if the
number of different types of particles increases sema20 . Quite general type
of Hamiltonians in $D$ dimensions can be considered sema17 ; sema18b ; sema19
for which the ET can yield analytical bounds for the spectra in favourable
situations. If this is not the case, numerical approximations with or without
variational character can always be computed.

The principle of the ET relies on the fact that the complete solution of a
$N$-body harmonic oscillator Hamiltonian with two-body forces, say
$H_{\textrm{ho}}$, can be obtained by the diagonalisation of a matrix of order
$(N-1)$ hall79 ; silv10 ; cint01 . The idea is to build an auxiliary
Hamiltonian, $\tilde{H}=H_{\textrm{ho}}+B$, where $B$ is a function of the
masses and coupling constants of $H_{\textrm{ho}}$, whose structure is
unequivocally determined from the general $N$-body Hamiltonian under study,
say $H$. The approximate eigenvalues are obtained by an extremisation
procedure performed on the eigenvalues of $\tilde{H}$ sema20 . These
approximate values depend on characteristic global quantum numbers $Q$, whose
the common structure given below implies a strong degeneracy of the levels,
inherited from $H_{\textrm{ho}}$.

In the case of systems composed of identical particles for $D\geq 2$, it is
possible to partly correct this drawback and to improve the accuracy of the ET
approximations by using two results stemming from the study of 2-body systems,
but the price to pay is that the possible variational character of the
approximate eigenvalues cannot longer be guaranteed. This has been done
previously in two steps:

  * •

In loba09 , it is shown that a universal effective quantum number exists for
centrally symmetric 2-body systems. According to this result, a unique
constant $\phi$ is introduced in the quantum number $Q$ of the ET (only one
for identical particles) which yields improvements of the eigenvalues and
observables, very significant for some Hamiltonians sema15a . The optimal
value of $\phi$ for the whole spectrum can be simply deduced from known
accurate numerical results horn14 .

  * •

In the dominantly orbital state method (DOSM), an approximate solution is
found by quantising the radial motion around a semiclassical solution for a
purely orbital motion. Developed at the origin for 2-body systems olss97 , it
has been extended to 3-body systems sema13b . In sema15b , the DOSM is
generalised to $N$-body systems and its approximations are compared with the
solutions from the ET to predict the value of $\phi$. Results are obtained in
perfect agreement with those obtained in sema15a .

The purpose of this paper is to generalise these procedures for systems
containing $N$ identical particles plus a different one, in order to extend
the use of the ET for new systems in different physical contexts. For
instance, in the Large-$N$ approach of QCD, the SU(3) gauge group is replaced
by SU($N$), leading up to baryons made of $N$ quarks witt79 . To treat
hadronic systems within this approach in combination with the constituent
model will16  seems interesting in this energy domain where the usual
expansion methods of the quantum field theory are not relevant. The use of the
ET allows for a systematic study of baryons as a function of $N$ without the
need of an enormous numerical cost, and the accuracy expected is sufficient
for this kind of systems, as it appears in buis11 ; buis12 . So, an extension
of the improved ET to $N+1$ systems will allow for computations for baryons
containing a quark different from the others or for hybrid baryons containing
identical quarks plus a constituent gluon. But, this technique could also be
used in other domains where it is necessary to treat systems with a great
number of particles. For instance, it could be used to estimate the binding
energies of nuclei or clusters of cold atoms for which ab-initio calculations
are already available, as for instance in gatt13 ; gatt11 ; kiev20 . The
method has also allowed for the study of a possible quasi Kepler’s third law
for quantum many-body systems sema21 .

The paper is organised as follows. The improvement of the ET with the DOSM for
$N$ identical particles is recalled in Sec. II. This section is developed to
describe the principle of the calculations in a simpler context, but also to
present new formulas. Indeed, the results given in sema15b  are correct but
uselessly complicated, due to an unnoticed simplification. In Sec. III, the
improvement for the ET is extended to systems with all particles identical
except one. The practical formulas to use are described and two coherence
tests are performed. Different systems are studied in Sec. IV to show the
relevance of this new method but also its limitations. Some concluding remarks
are finally given in Sec. V where some prospects are also presented. In the
following, we work in natural units ($\hbar=c=1$) for $D\geq 2$.

##  II $\bm{N}$ identical particles

For $N$ identical particles, the Hamiltonian considered is written

$H=\sum_{i=1}^{N}T(|\bm{p}_{i}|)+\sum_{i<j=2}^{N}V(|\bm{r}_{i}-\bm{r}_{j}|),$
(1)

where $T$ is a kinetic energy (with some constraints sema18a ) and $V$ a two-
body central potential, and where the centre of mass motion is removed
($\sum_{i=1}^{N}\bm{p}_{i}=\bm{0}$). An ET approximate eigenvalue $E$ is given
by the following compact set of equations for a completely (anti)symmetrised
state of the auxiliary Hamiltonian $\tilde{H}$ associated with (1):

$\displaystyle E=N\,T(p_{0})+C^{2}_{N}\,V(\rho_{0}),$ (2a) $\displaystyle
N\,T^{\prime}(p_{0})\,p_{0}=C^{2}_{N}\,V^{\prime}(\rho_{0})\,\rho_{0},$ (2b)
$\displaystyle Q(N)=\sqrt{C^{2}_{N}}\,\rho_{0}\,p_{0},$ (2c)

where $U^{\prime}(x)=dU/dx$, $C^{2}_{N}=N(N-1)/2$ is the number of pairs,
$\rho_{0}$ and $p_{0}$ are the variables to be determined, and

$Q(N)=\sum_{i=1}^{N-1}\left(2\,n_{i}+l_{i}+\frac{D}{2}\right)$ (3)

is a global quantum number. This set is called compact because all the
relevant variables ($\rho_{0}$ and $p_{0}$, see below) appear in three
equations giving the expression for the energy (2a), the equation of motion
(2b) and the rule for the quantisation (2c). The $\\{n_{i},l_{i}\\}$ are the
quantum numbers associated with the $N-1$ internal Jacobi coordinates silv10 .
According to the nature of the particles, bosons or fermions, only some values
are allowed for $Q(N)$. Their determination is a difficult problem, even for
$N=3$ stan97 . But the ground states are determined for all values of $N$ in
sema20 . Depending on the forms of $T$ and $V$, the approximate value $E$ can
have a variational character sema20 . The method is easy to implement since it
reduces to finding the pair ($\rho_{0}$,$p_{0}$) from a transcendental
equation, and then computing $E$. Let us mention that sema15b

$p_{0}^{2}=\frac{1}{N}\left\langle\sum_{i=1}^{N}\bm{p}_{i}^{2}\right\rangle\quad\textrm{and}\quad\rho_{0}^{2}=\frac{1}{C^{2}_{N}}\left\langle\sum_{i<j=2}^{N}(\bm{r}_{i}-\bm{r}_{j})^{2}\right\rangle,$
(4)

where the mean value is computed with the eigenstate of the auxiliary
Hamiltonian $\tilde{H}$ corresponding to $Q(N)$. So, $\rho_{0}$ can be
interpreted as the mean distance between two particles and $p_{0}$ as the mean
momentum of each particle. In sema13 ; sema20 , the variable
$r_{0}=\sqrt{C^{2}_{N}}\,\rho_{0}$ is used because one-body potentials were
also considered. It has been checked with the systems studied in sema15a  that
the energies computed by (2) are the same that the ones computed with the
extremisation procedure mentioned in the introduction. Nevertheless, the
compact set (2) is interesting because it yields directly the mean values
$\rho_{0}$ and $p_{0}$ which are meaningful physical quantities. This set has
a nice semiclassical interpretation sema13  and it is a compulsory step of the
procedure explained below for the improvement by the DOSM:

  1. 1.

In order to break the strong degeneracy inherent to the method, the global
quantum number $Q(N)$ is shared into angular and radial contributions by the
introduction of a parameter $\phi$

$\displaystyle Q_{\phi}(N)$
$\displaystyle=\phi\,\nu+\lambda\quad\textrm{where}$ (5a) $\displaystyle\nu$
$\displaystyle=\sum_{i=1}^{N-1}\left(n_{i}+\frac{1}{2}\right)\quad\textrm{and}\quad\lambda=\sum_{i=1}^{N-1}\left(l_{i}+\frac{D-2}{2}\right).$
(5b)

This is inspired from the existence of an effective quantum number $q$ that
determines with high accuracy the level ordering of centrally symmetric 2-body
systems loba09

$q=\phi(n+1/2)+l+(D-2)/2.$ (6)

The parameter $\phi$ in (6) depends on the system, but not on quantum numbers
$n$ and $l$ (for practical purposes, our definition of $\phi$ differs from the
one used in loba09 ). In (5), the same assumption is extended to $N\geq 2$.
Note that the original quantum number $Q(N)$ is recovered with $\phi=2$.

  2. 2.

The second step is to obtain the DOSM solution of $H$. It has been shown in
sema15b  that it can be computed from the ET compact set of equations. The
pair ($\tilde{\rho}_{0}$,$\tilde{p}_{0}$) and the corresponding energy
$\tilde{E}$ are computed for a purely collective orbital motion by solving (2)
but with $Q(N)$ replaced by $\lambda$:

$\displaystyle\tilde{E}=N\,T(\tilde{p}_{0})+C^{2}_{N}\,V(\tilde{\rho}_{0}),$
(7a) $\displaystyle
N\,T^{\prime}(\tilde{p}_{0})\,\tilde{p}_{0}=C^{2}_{N}\,V^{\prime}(\tilde{\rho}_{0})\,\tilde{\rho}_{0},$
(7b) $\displaystyle\lambda=\sqrt{C^{2}_{N}}\,\tilde{\rho}_{0}\,\tilde{p}_{0}.$
(7c)

Under the hypothesis that $\lambda\gg\nu$, a collective radial motion of the
$N$ particles is introduced as a small perturbation by computing the energy
from (7a) with the following replacements

$\tilde{\rho}_{0}\to\tilde{\rho}_{0}+\Delta\rho\quad\textrm{and}\quad\tilde{p}_{0}\to\sqrt{p_{r}^{2}+\frac{\lambda^{2}}{C^{2}_{N}(\tilde{\rho}_{0}+\Delta\rho)^{2}}},$
(8)

inspired by semi-classical considerations developed in olss97 . With the
assumptions given above, $\Delta\rho\ll\tilde{\rho}_{0}$ and
$p_{r}\ll\tilde{p}_{0}=\lambda/(\sqrt{C_{N}^{2}}\tilde{\rho}_{0})$, where
$\Delta\rho$ is a radial displacement and $p_{r}$ is a radial momentum. So, a
power expansion of the radial contribution $\Delta E$ to $\tilde{E}$, limited
to the lowest non-vanishing orders, gives

$\displaystyle\Delta E$
$\displaystyle\approx\frac{1}{2\mu}p_{r}^{2}+\frac{k}{2}\Delta\rho^{2}\quad\textrm{with}$
(9) $\displaystyle\mu=\frac{\tilde{p}_{0}}{N\,T^{\prime}(\tilde{p}_{0})},$
$\displaystyle
k=\frac{2\,N\,\tilde{p}_{0}}{\tilde{\rho}_{0}^{2}}\,T^{\prime}(\tilde{p}_{0})+\frac{N\,\tilde{p}_{0}^{2}}{\tilde{\rho}_{0}^{2}}\,T^{\prime\prime}(\tilde{p}_{0})+C^{2}_{N}\,V^{\prime\prime}(\tilde{\rho}_{0}).$

The term in $\Delta\rho$ is canceled by the equation of motion (7b). This
Hamiltonian must be quantised to compute the radial excitations. As $p_{r}$
and $\Delta\rho$ are effective collective variables, it is not clear if
$\Delta\rho$, or a multiple of it, is the conjugate of $p_{r}$. Nevertheless,
if we assume that $\Delta\rho$ is the conjugate variable of $p_{r}$, the
solution is

$\Delta E\approx\sqrt{\frac{k}{\mu}}\,\left(n+\frac{1}{2}\right).$ (10)

The effective quantum number $n$ must be reinterpreted in terms of the radial
contribution $\nu$. The application of the formulas above to the case of the
$N$-body harmonic oscillator, for which the DOSM gives the exact solution,
shows that

$\sqrt{C^{2}_{N}}\,\left(n+\frac{1}{2}\right)=\nu.$ (11)

This constraint to recover the exact solution for the $N$-body harmonic
oscillator allows to fix the problem about the identification of the correct
pair of conjugate variables. Indeed, the same result for $\Delta E$ is
obtained with another choice of the conjugate variables, but (11) will be
different as in sema15b . The complete solution for the DOSM is then
$\tilde{E}+\Delta E$.

  3. 3.

The system (2) is then solved with $Q(N)\to Q_{\phi}(N)=\lambda(1+\epsilon)$
where $\epsilon=\phi\,\nu/\lambda\ll 1$, in order to be in the same situation
as the one of the previous item. If $\epsilon=0$, the purely orbital solution
$\tilde{E}$ is recovered. Starting from this solution, the energy is expanded
in power of $\epsilon$, giving to the first order,

$\Delta E\approx N\,\tilde{p}_{0}\,T^{\prime}(\tilde{p}_{0})\,\epsilon.$ (12)

This expression is identical to the one found in sema15b  but its form is
simpler due to previously unnoticed vanishing contributions. The comparison of
(10) and (11) with (12) gives directly the value of $\phi$

$\phi=\frac{\lambda}{N\,\tilde{p_{0}}\,T^{\prime}(\tilde{p_{0}})}\sqrt{\frac{k}{C_{N}^{2}\,\mu}},$
(13)

with $\mu$ and $k$ given by (9). This value has been determined in the
specific conditions where the radial motion has a small contribution to the
energy with respect to the orbital motion, because it is computable with the
DOSM. But since the structure of $\phi$ is assumed to be valid for any quantum
numbers, (13) can be used in principle for any state considered.

  4. 4.

With the values $\lambda$ and $\nu$ chosen, and the value of $\phi$ computed
with (13), the system (2) is now solved with $Q(N)$ replaced by $Q_{\phi}(N)$
to yield the improved ET (IET) approximation. The original ET, with its
possible variational solutions, is then recovered with $\phi=2$. For other
values of $\phi$, the possible variational character of the solution cannot be
guaranteed.

As an example of Hamiltonian, let us consider the generic case

$T(x)=F\,x^{\alpha}\quad\textrm{and}\quad
V(x)=\mathrm{sgn}(\beta)\,G\,x^{\beta},$ (14)

where $F$ and $\alpha$ are strictly positive numbers to satisfy the
constraints in sema18a , and where $\beta\neq 0$ and $G$ is a strictly
positive number because of the sign function in $V$. The nonrelativistic case
is covered with $F=1/(2m)$ and $\alpha=2$, and the ultrarelativistic one with
$F=1$ and $\alpha=1$. The application of the above procedure gives

$E=\mathrm{sgn}(\beta)\,(\alpha+\beta)\left(\left(\frac{C_{N}^{2}\,G}{\alpha}\right)^{\alpha}\left(\frac{N\,F}{|\beta|}\right)^{\beta}\left(\frac{Q_{\phi}(N)}{\sqrt{C_{N}^{2}}}\right)^{\alpha\,\beta}\right)^{1/(\alpha+\beta)}$
(15)

with the remarkably simple value

$\phi=\sqrt{\alpha+\beta},$ (16)

independent of $F$, $G$, $N$ and the global quantum number. According to the
meaning of $T$, $E$ can be the mass of the system and then must be a positive
number, or $E$ can be the binding energy and then can be negative or positive
(for a confining potential). In this last case, bound states exist if $E$ and
$\beta$ have the same sign. So, (15) implies that $\alpha+\beta>0$, which is
coherent with (16). If $E$ is the mass (as in the ultrarelativistic case), a
strong attractive potential between the particles does not sound physically
relevant and it seems natural to impose $\beta>0$, and then $\alpha+\beta>0$.
When $\alpha=\beta=2$, the exact result for the $N$-body harmonic oscillator
is recovered, as expected. The quality of (15)-(16) is tested for $N=3$ in
Table 1. The “exact” results from basd90  are obtained with an elaborate
hyperspherical expansion up to a grand orbital momentum $L=8$, which insures a
good convergence of the expansion and a high accuracy of the eigenvalues
rich81 .

Table 1: Bosonic ground state energies for the case (14) with $N=3$, $\alpha=2$, $F=G=0.5$ and several values of $\beta$ (arbitrary units). The results from the ET (upper bound for $\beta<2$ and lower bound for $\beta>2$) and the IET with (16) are compared with the “exact” ones basd90 . The relative errors in % are indicated between square brackets. In the case of $\beta=2$, ET and IET give the exact result as expected. $\beta$ | “exact” | ET | IET
---|---|---|---
$-1$ | $-0.26675$ | $-0.12500\,[53]$ | $-0.28125\,[5.4]$
$-0.5$ | $-0.59173$ | $-0.49139\,[17]$ | $-0.59977\,[1.4]$
0.1 | 1.88019 | 1.91406 [1.8] | 1.87743 [0.15]
0.5 | 2.91654 | 3.08203 [5.7] | 2.90211 [0.49]
1 | 3.86309 | 4.08852 [5.8] | 3.84130 [0.56]
2 | 5.19615 | 5.19615 [0] | 5.19615 [0]
3 | 6.15591 | 5.68394 [7.7] | 6.22479 [1.1]

The IET is tested on various systems in sema15b  with a dramatic improvement
in some cases. That is the reason why it is interesting to generalise the
combination of the ET with the DOSM to systems with different particles. This
is done in the next section.

##  III $\bm{N_{a}+1}$ systems

The ET for Hamiltonians with different particles is detailed in sema20  and
tested with systems containing $N_{a}$ particles of type $a$ and only one
different particle of type $b$. The general Hamiltonian for such systems is
written

$H=\sum_{i=1}^{N_{a}}T_{a}(|\bm{p}_{i}|)+T_{b}(|\bm{p}_{N}|)+\sum_{i<j=2}^{N_{a}}V_{aa}(|\bm{r}_{i}-\bm{r}_{j}|)+\sum_{i=1}^{N_{a}}V_{ab}(|\bm{r}_{i}-\bm{r}_{N}|),$
(17)

with $N=N_{a}+1$. The corresponding compact set of five ET equations, for an
eigenstate completely (anti)symmetrised for the $N_{a}$ particles, is given in
cimi21

$\displaystyle
E=N_{a}\,T_{a}(p_{a}^{\prime})+T_{b}(P_{0})+C_{N_{a}}^{2}\,V_{aa}(r_{aa})+N_{a}\,V_{ab}(r_{0}^{\prime}),$
(18a) $\displaystyle
N_{a}\,T_{a}^{\prime}(p_{a}^{\prime})\,\frac{p_{a}^{2}}{p_{a}^{\prime}}=C_{N_{a}}^{2}\,V_{aa}^{\prime}(r_{aa})\,r_{aa}+\frac{N_{a}-1}{2}V^{\prime}_{ab}(r_{0}^{\prime})\frac{r_{aa}^{2}}{r_{0}^{\prime}},$
(18b)
$\displaystyle\frac{1}{N_{a}}\,T^{\prime}_{a}(p_{a}^{\prime})\,\frac{P_{0}^{2}}{p_{a}^{\prime}}+T_{b}^{\prime}(P_{0})\,P_{0}=N_{a}\,V^{\prime}_{ab}(r_{0}^{\prime})\,\frac{R_{0}^{2}}{r_{0}^{\prime}},$
(18c) $\displaystyle Q(N_{a})=\sqrt{C_{N_{a}}^{2}}\,p_{a}\,r_{aa},$ (18d)
$\displaystyle Q(2)=P_{0}\,R_{0},$ (18e) with
$p_{a}^{\prime\,2}=p_{a}^{2}+\frac{P^{2}_{0}}{N_{a}^{2}}\quad\textrm{and}\quad
r_{0}^{\prime\,2}=R_{0}^{2}+\frac{N_{a}-1}{2\,N_{a}}\,r_{aa}^{2},$ (18f)

and with $Q(N)$ given by (3). $Q(N_{a})$ is associated with the internal
motions of the $N_{a}$ particles and its values are constrained by their
bosonic or fermionic nature. But any value is allowed for $Q(2)$ which is
associated with $R_{0}$, $R_{0}^{2}$ being the mean quadratic distance between
the particle $b$ and the centre of mass of the $N_{a}$ particles sema20 ;
cimi21 . Checks of (18) are provided in cimi21 .

To partially break the strong degeneracy implied by the global quantum numbers
$Q(N_{a})$ and $Q(2)$ in (18), and thus to improve the results of the ET, the
four-step procedure described above for $N$ identical particles is generalised
to the $N_{a}+1$ system in the four following subsections.

###  III.1 Separation between angular and radial dynamics

According to the situation for systems with all identical particles, the two
global quantum numbers are splitted into an angular part and a radial part:

$\displaystyle Q_{\phi}(N_{a})$
$\displaystyle=\phi_{a}\,\nu_{a}+\lambda_{a}\quad\textrm{where}$ (19a)
$\displaystyle\nu_{a}=\sum_{i=1}^{N_{a}-1}\left(n_{i}+\frac{1}{2}\right)\quad\textrm{and}\quad\lambda_{a}=\sum_{i=1}^{N_{a}-1}\left(l_{i}+\frac{D-2}{2}\right),$
(19b) $\displaystyle Q_{\phi}(2)$
$\displaystyle=\phi_{b}\,\nu_{b}+\lambda_{b}\quad\textrm{where}$ (19c)
$\displaystyle\nu_{b}=n_{b}+\frac{1}{2}\quad\textrm{and}\quad\lambda_{b}=l_{b}+\frac{D-2}{2}.$
(19d)

Two new parameters $\phi_{a}$ and $\phi_{b}$ are now introduced. With the
nature of $Q(2)$ mentioned above, the quantum numbers
$\\{\nu_{b},\lambda_{b}\\}$ are linked to the relative motion between the
particle $b$ and the centre of mass of the $N_{a}$ particles. The quantum
numbers $\\{\nu_{a},\lambda_{a}\\}$ are associated with a global excitation of
the $N_{a}$ identical particles.

###  III.2 DOSM for $\bm{N_{a}+1}$ systems

The DOSM solution for a $N_{a}+1$ system is obtained in solving (18) with
$Q(N_{a})$ replaced by $\lambda_{a}$ and $Q(2)$ replaced by $\lambda_{b}$,
that is to say in considering only orbital motions without radial excitations.
An energy $\tilde{E}$ is then computed with the associated variables
$(\tilde{p}_{a},\tilde{r}_{aa},\tilde{P}_{0},\tilde{R}_{0})$. Under the
hypothesis that $\lambda_{a}\gg\nu_{a}$ and that $\lambda_{b}\gg\nu_{b}$,
radial motions are introduced as small perturbations by computing the energy
from (18a) with the following replacements

$\tilde{r}_{aa}\to\tilde{r}_{aa}+\Delta
r,\quad\tilde{p}_{a}\to\sqrt{p_{r}^{2}+\frac{\lambda_{a}^{2}}{C^{2}_{N_{a}}(\tilde{r}_{aa}+\Delta
r)^{2}}}$ (20)

and

$\tilde{R}_{0}\to\tilde{R}_{0}+\Delta
R,\quad\tilde{P}_{0}\to\sqrt{P_{r}^{2}+\frac{\lambda_{b}^{2}}{(\tilde{R}_{0}+\Delta
R)^{2}}},$ (21)

in analogy with the replacement (8) made for all identical particles. With the
assumptions given above, $\Delta r\ll\tilde{r}_{aa}$, $p_{r}\ll\tilde{p}_{a}$,
$\Delta R\ll\tilde{R}_{0}$ and $P_{r}\ll\tilde{P}_{0}$. So, a power expansion
of the sum $\Delta E$ of the radial contributions to $\tilde{E}$, limited to
the lowest non-vanishing orders, gives after a lengthy calculation

$\displaystyle\Delta E$
$\displaystyle\approx\frac{1}{2}\left(\frac{1}{\mu_{a}}p_{r}^{2}+\frac{1}{\mu_{b}}P_{r}^{2}+k_{a}\,\Delta
r^{2}+k_{b}\,\Delta R^{2}+k_{c}\,\Delta r\Delta R\right)\quad\textrm{with }$
(22)
$\displaystyle\mu_{a}=\frac{\tilde{p}_{a}^{\prime}}{N_{a}T_{a}^{\prime}(\tilde{p}_{a}^{\prime})},$
$\displaystyle\mu_{b}=\left(\frac{T_{a}^{\prime}(\tilde{p}_{a}^{\prime})}{N_{a}\tilde{p}_{a}^{\prime}}+\frac{T_{b}^{\prime}(\tilde{P}_{0})}{\tilde{P}_{0}}\right)^{-1},$
$\displaystyle
k_{a}=\frac{N_{a}T^{\prime\prime}_{a}(\tilde{p}_{a}^{\prime})\tilde{p}_{a}^{4}}{\tilde{r}_{aa}^{2}\tilde{p}_{a}^{\prime\,2}}+\frac{N_{a}T_{a}^{\prime}(\tilde{p}_{a}^{\prime})\tilde{p}_{a}^{2}}{\tilde{r}_{aa}^{2}}\left(\frac{3}{\tilde{p}_{a}^{\prime}}-\frac{\tilde{p}_{a}^{2}}{\tilde{p}_{a}^{\prime\,3}}\right)+C_{N_{a}}^{2}V^{\prime\prime}_{aa}(\tilde{r}_{aa})$
$\displaystyle\text{\phantom{imm}}+\frac{(N_{a}-1)^{2}\tilde{r}_{aa}^{2}}{4N_{a}\tilde{r}_{0}^{\prime\,2}}V^{\prime\prime}_{ab}(\tilde{r}_{0}^{\prime})+\frac{(N_{a}-1)}{2}\left(\frac{1}{\tilde{r}_{0}^{\prime}}-\frac{(N_{a}-1)\tilde{r}_{aa}^{2}}{2N_{a}\tilde{r}_{0}^{\prime\,3}}\right)V_{ab}^{\prime}(\tilde{r}_{0}^{\prime}),$
$\displaystyle
k_{b}=\frac{T^{\prime\prime}_{a}(\tilde{p}_{a}^{\prime})\tilde{P}_{0}^{4}}{N_{a}^{3}\tilde{R}_{0}^{2}\tilde{p}_{a}^{\prime\,2}}+\frac{T_{b}^{\prime\prime}(\tilde{P}_{0})\tilde{P}_{0}^{2}}{\tilde{R}_{0}^{2}}+\frac{T_{a}^{\prime}(\tilde{p}_{a}^{\prime})\tilde{P}_{0}^{2}}{N_{a}\tilde{R}_{0}^{2}}\left(\frac{3}{\tilde{p}_{a}^{\prime}}-\frac{\tilde{P}_{0}^{2}}{N_{a}^{2}\tilde{p}_{a}^{\prime\,3}}\right)$
$\displaystyle\text{\phantom{imm}}+\frac{2T^{\prime}_{b}(\tilde{P}_{0})\tilde{P}_{0}}{\tilde{R}_{0}^{2}}+\frac{N_{a}\tilde{R}_{0}^{2}}{\tilde{r}_{0}^{\prime\,2}}V^{\prime\prime}_{ab}(\tilde{r}_{0}^{\prime})+N_{a}\left(\frac{1}{\tilde{r}_{0}^{\prime}}-\frac{\tilde{R}_{0}^{2}}{\tilde{r}_{0}^{\prime\,3}}\right)V^{\prime}_{ab}(\tilde{r}_{0}^{\prime}),$
$\displaystyle
k_{c}=\frac{2\tilde{p}_{a}^{2}\tilde{P}_{0}^{2}}{N_{a}\tilde{p}_{a}^{\prime\,2}\tilde{r}_{aa}\tilde{R}_{0}}\left(T_{a}^{\prime\prime}(\tilde{p}_{a}^{\prime})-\frac{T_{a}^{\prime}(\tilde{p}_{a}^{\prime})}{\tilde{p}_{a}^{\prime}}\right)+\frac{(N_{a}-1)\tilde{r}_{aa}\tilde{R}_{0}}{\tilde{r}_{0}^{\prime\,2}}\left(V_{ab}^{\prime\prime}(\tilde{r}_{0}^{\prime})-\frac{V_{ab}^{\prime}(\tilde{r}_{0}^{\prime})}{\tilde{r}_{0}^{\prime}}\right).$

There is no term in $\Delta r$ nor in $\Delta R$, but a term $\Delta r\Delta
R$ couples the two radial motions. This energy $\Delta E$ must be quantised.
From the nature of the variable $R_{0}$, it is clear that $\Delta R$ is the
conjugate variable of $P_{r}$. For $\Delta r$ and $p_{r}$, the situation is
similar to the case of all identical particles: $\Delta r$ is taken as the
conjugate variable of $p_{r}$. The solution of the two one-dimensional coupled
oscillators is given by (35) in Appendix A. In order to interpret correctly
the two quantum numbers $n$ and $n^{\prime}$ associated with (35), the
formulas above are applied to the equivalent harmonic oscillator whose exact
solution is given by formula (37) in sema20  without one-body potentials. The
comparison between the DOSM solution and the exact one shows that

$\sqrt{C^{2}_{N_{a}}}\,\left(n+\frac{1}{2}\right)=\nu_{a}\quad\textrm{and}\quad\left(n^{\prime}+\frac{1}{2}\right)=\nu_{b}.$
(23)

The complete solution for the DOSM is then

$E_{\textrm{DOSM}}=\tilde{E}+\sqrt{\frac{A}{C^{2}_{N_{a}}\mu}}\sum_{i=1}^{N_{a}-1}\left(n_{i}+\frac{1}{2}\right)+\sqrt{\frac{B}{\mu}}\left(n_{b}+\frac{1}{2}\right),$
(24)

with $\mu$, $A$ and $B$ given by (35).

In order to check these formulas, we can consider the limit in which all
particles are identical, that is to say $T_{a}(x)=T_{b}(x)=T(x)$ and
$V_{aa}(x)=V_{ab}(x)=V(x)$. In this case, an eigenstate must be completely
(anti)symmetrised. So, we must impose $r_{aa}=r_{0}^{\prime}$ and
$p_{a}^{\prime}=P_{0}$, and as a consequence,
$R_{0}^{2}=(N_{a}+1)r_{aa}^{2}/(2N_{a})$ and
$p_{a}^{2}=(N_{a}^{2}-1)P_{0}^{2}/N_{a}^{2}$ cimi21 . These constraints must
then be imposed on both orbital and radial motions. Indeed, it can be shown
that (22) reduces to (9) when it is imposed that
$\tilde{R}_{0}^{2}=(N_{a}+1)\tilde{r}_{aa}^{2}/(2N_{a})$, $\Delta
R^{2}=(N_{a}+1)\Delta r^{2}/(2N_{a})$,
$\tilde{p}_{a}^{2}=(N_{a}^{2}-1)\tilde{P}_{0}^{2}/N_{a}^{2}$ and
$p_{r}^{2}=(N_{a}^{2}-1)P_{r}^{2}/N_{a}^{2}$. Let us note that this limit
cannot be checked in (24), since no constraint has been imposed on the the
radial motions controlled by $\Delta R$ and $\Delta r$ to obtain this formula.

###  III.3 ET for $\bm{N_{a}+1}$ systems with small radial motions

To compute the solutions of (18) for small radial excitations, $Q(N_{a})$ is
replaced by $Q_{\phi}(N_{a})=\lambda_{a}(1+\epsilon_{a})$ where
$\epsilon_{a}=\phi_{a}\,\nu_{a}/\lambda_{a}\ll 1$, and $Q(2)$ by
$Q_{\phi}(2)=\lambda_{b}(1+\epsilon_{b})$ where
$\epsilon_{b}=\phi_{b}\,\nu_{b}/\lambda_{b}\ll 1$. The case
$\epsilon_{a}=\epsilon_{b}=0$ corresponds to the purely orbital solution of
Sec. III.2. So, we can write $r_{aa}=\tilde{r}_{aa}+\Delta r$ with $\Delta
r\ll\tilde{r}_{aa}$, $p_{a}=\tilde{p}_{a}+\Delta p$ with $\Delta
p\ll\tilde{p}_{a}$, $R_{0}=\tilde{R}_{0}+\Delta R$ with $\Delta
R\ll\tilde{R}_{0}$ and $P_{0}=\tilde{P}_{0}+\Delta P$ with $\Delta
P\ll\tilde{P}_{0}$. Equations (18d) and (18e) give a link between $\Delta r$,
$\Delta p$ and $\epsilon_{a}$, and between $\Delta R$, $\Delta P$ and
$\epsilon_{b}$. Then, (18a) can be developed at first order in $\epsilon_{a}$
and $\epsilon_{b}$. All these lengthy calculations give

$E\approx\tilde{E}+T^{\prime}_{a}(\tilde{p}_{a}^{\prime})\frac{N_{a}\tilde{p}_{a}^{2}}{\tilde{p}_{a}^{\prime}}\epsilon_{a}+\left(T^{\prime}_{a}(\tilde{p}_{a}^{\prime})\frac{\tilde{P}_{0}^{2}}{N_{a}\tilde{p}_{a}^{\prime}}+T_{b}^{\prime}(\tilde{P}_{0})\tilde{P}_{0}\right)\epsilon_{b}.$
(25)

In order to check this formula, we again consider the limit in which all
particles are identical as in the previous subsection. Considerations about
the solution $(\tilde{p}_{a},\tilde{r}_{aa},\tilde{P}_{0},\tilde{R}_{0})$ are
unchanged, and we have to deal with the quantities $\epsilon_{a}$ and
$\epsilon_{b}$. When all particles are identical, only the sums
$Q(N_{a})+Q(2)=Q(N_{a}+1)$ and $\lambda_{a}+\lambda_{b}$ are relevant, which
implies that $\epsilon_{a}=\epsilon_{b}$. This is not surprising since these
two numbers are ratios between radial and orbital excitations. Under these
constraints, (25) reduces to (12), as expected.

Making explicit $\epsilon_{a}$ and $\epsilon_{b}$ in terms of $\phi_{a}$ and
$\phi_{b}$, (25) is written

$\displaystyle
E\approx\tilde{E}+D_{a}\frac{\phi_{a}}{\lambda_{a}}\sum_{i=1}^{N_{a}-1}\left(n_{i}+\frac{1}{2}\right)+D_{b}\frac{\phi_{b}}{\lambda_{b}}\left(n_{b}+\frac{1}{2}\right)\quad\textrm{with}$
(26) $\displaystyle
D_{a}=T^{\prime}_{a}(\tilde{p}_{a}^{\prime})\frac{N_{a}\tilde{p}_{a}^{2}}{\tilde{p}_{a}^{\prime}}\quad\textrm{and}\quad
D_{b}=T^{\prime}_{a}(\tilde{p}_{a}^{\prime})\frac{\tilde{P}_{0}^{2}}{N_{a}\tilde{p}_{a}^{\prime}}+T_{b}^{\prime}(\tilde{P}_{0})\tilde{P}_{0}.$

The comparison of (24) with (26) gives

$\phi_{a}=\frac{\lambda_{a}}{D_{a}}\sqrt{\frac{A}{C_{N_{a}}^{2}\mu}}\quad\textrm{and}\quad\phi_{b}=\frac{\lambda_{b}}{D_{b}}\sqrt{\frac{B}{\mu}}.$
(27)

The quantities $\mu$, $A$, $B$, $D_{a}$ and $D_{b}$ can be computed from the
values of the solution
$(\tilde{p}_{a},\tilde{r}_{aa},\tilde{P}_{0},\tilde{R}_{0})$ fixed by the
choice of the collective orbital quantum numbers $\lambda_{a}$ and
$\lambda_{b}$. If the complete expressions for $\phi_{a}$ and $\phi_{b}$ are
quite complicated, their computation is straightforward.

###  III.4 Summary of the procedure

Once the values $\lambda_{a}$, $\lambda_{b}$, $\nu_{a}$ and $\nu_{b}$ are
chosen, the value of $\phi_{a}$ and $\phi_{b}$ can be computed with (27). The
system (18) is then solved with $Q(N_{a})$ replaced by $Q_{\phi}(N_{a})$ and
$Q(2)$ by $Q_{\phi}(2)$ given by (19) to yield the IET. If a solution of the
original ET has a variational character, this cannot be guaranteed for the
corresponding one obtained by the IET.

By comparing the procedures described in Sec. II and Sec. III, one can
appreciate the increase in complexity by simply adding one different particle
to a set of identical ones. The system (18) is so complex that no analytical
solution was found, except in the case of harmonic oscillators. For this last
system, it has been checked that the exact solution (formula (37) in sema20
without one-body potentials) is recovered by this procedure and that
$\phi_{a}=\phi_{b}=2$ in this case, as expected.

##  IV Numerical tests

For systems with identical particles, it was possible to test the ET ($D=1$
and 3) and the IET ($D=3$) sema15a ; sema15b ; sema19  thanks to the existence
of very accurate results on simple many-body systems in the literature horn14
; timo12 ; timo17 . Accurate results for simple systems with one particle
different from all others are more scarce. In the following, the many-body
systems considered in sema20  with $D=3$ are studied here with the ET and the
IET. As no analytical solutions were found for the set (18) for these systems,
all calculations are performed numerically.

###  IV.1 Three-body ultrarelativistic oscillators

The first system considered is composed of 3 ultrarelativistic bosons with a
vanishing mass, interacting via harmonic oscillator potentials. The
Hamiltonian in arbitrary units is written as

$H=\sum_{i=1}^{3}|\bm{p}_{i}|+(\bm{r}_{1}-\bm{r}_{2})^{2}+\kappa\sum_{i=1}^{2}(\bm{r}_{i}-\bm{r}_{3})^{2}.$
(28)

When $\kappa\neq 1$, the 3rd particle is different since it has a different
interaction with the two other ones. Accurate numerical solutions for this
Hamiltonian do not seem to be available in the literature. Fortunately, it is
shown in salo20  that such solutions can be obtained from a rescaling of
solutions for a nonrelativistic 3-body systems with a linear potential. This
latter property can also be seen on the ET and IET solutions (15) and (16)
where it is clear that the energy is invariant under the exchanges
$\alpha\leftrightarrow\beta$ and $F\leftrightarrow G$ when $\beta>0$ and
$N=3$. The accurate solutions of (28) presented in Table 2 are properly
rescaled nonrelativistic accurate results provided by J.M. Richard. He
computed them with the previously mentioned method developed in rich81 ;
basd90 . They are compared in Table 2 with the upper bounds obtained with the
ET and the results computed with the IET. For two values of $\kappa$, the
three states presented have a positive parity with a vanishing total orbital
angular momentum, the lowest one being the ground state. The quantum numbers
$\nu_{a}$, $\lambda_{a}$, $\nu_{b}$ and $\lambda_{b}$ are assigned thanks to
the knowledge of the global quantum numbers and the symmetry between the two
identical particles of the exact solutions.

Table 2: Positive parity spherical states of Hamiltonian (28) for two values of $\kappa$. The results from the ET (upper bounds) and the IET with ($\phi_{a},\phi_{b}$) and the global numbers $\langle\nu_{a},\lambda_{a},\nu_{b},\lambda_{b}\rangle$ (see text) are compared with the “exact” ones (see text). The relative errors in % are indicated between square brackets. $\kappa=0.1$ |  | $\kappa=10$
---|---|---
“exact” | ET | IET | $(\phi_{a},\phi_{b})$ | $\langle\nu_{a},\lambda_{a},\nu_{b},\lambda_{b}\rangle$ |  | “exact” | ET | IET | $(\phi_{a},\phi_{b})$ | $\langle\nu_{a},\lambda_{a},\nu_{b},\lambda_{b}\rangle$
5.288 | 5.597 [5.8] | 5.307 [0.4] | (1.76,1.79) | $\langle 1/2,1/2,1/2,1/2\rangle$ |  | 14.506 | 15.353 [5.8] | 14.699 [1.3] | (1.82,1.80) | $\langle 1/2,1/2,1/2,1/2\rangle$
6.570 | 6.970 [5.9] | 6.571 [0.0] | (1.76,1.79) | $\langle 1/2,1/2,3/2,1/2\rangle$ |  | 19.134 | 20.272 [5.9] | 19.291 [0.8] | (1.82,1.80) | $\langle 3/2,1/2,1/2,1/2\rangle$
7.515 | 7.868 [4.7] | 7.625 [1.5] | (1.76,1.79) | $\langle 1/2,3/2,1/2,3/2\rangle$ |  | 20.340 | 21.580 [6.1] | 21.032 [3.4] | (1.82,1.80) | $\langle 1/2,3/2,1/2,3/2\rangle$

The accuracy obtained by the ET upper bounds are around 4-6%, while the
relative error can be divided by more than 10 for the IET results. Even if its
magnitude seems quite unpredictable, there is always an improvement with the
IET. For both values of $\kappa$, the parameters $\phi_{a}$ and $\phi_{b}$ are
close to $\sqrt{3}\approx 1.732$, the value (16) for the case of identical
particles ($\kappa=1$). The values of $(\phi_{a},\phi_{b})$ depend on $\kappa$
but are independent of the quantum numbers of the states considered. This
situation is similar to the one for all identical particles. Supplementary
calculations would be desirable to see if this property is preserved for all
the spectrum and for an arbitrary number of particles. Further comments about
the accuracy of these results are given in the following section.

###  IV.2 Nonrelativistic three-body systems with a power law potential

The second system is composed of three bosons interacting via different power
law potentials, the 3rd one having a different mass. The Hamiltonian in
arbitrary units is written as

$H=\sum_{i=1}^{2}\frac{\bm{p}_{i}^{2}}{2}+\frac{\bm{p}_{3}^{2}}{2m}+\frac{1}{2}\,\textrm{sgn}(\beta)\sum_{i<j=2}^{3}\bm{r}_{ij}^{\beta}.$
(29)

Variational bounds from the ET and the results computed with the IET for the
bosonic ground states are given in Table 3. The very accurate (“exact”)
results are taken from basd90 , and have been computed with the previously
mentioned method developed in rich81 .

Table 3: Bosonic ground state energies of Hamiltonian (29) for several values of $m$ and $\beta$. The results from the ET (upper bound for $\beta<2$ and lower bound for $\beta>2$) and the IET with ($\phi_{a},\phi_{b}$) are compared with the “exact” ones basd90 . The relative errors in % are indicated between square brackets. In the case of $\beta=2$, ET and IET give the exact result as expected. Results for $m=1$ are given in Table 1. |  | $m=0.2$ |  | $m=5$
---|---|---|---|---
$\beta$ |  | “exact” | ET | IET | $(\phi_{a},\phi_{b})$ |  | “exact” | ET | IET | $(\phi_{a},\phi_{b})$
$-1$ |  | $-0.1398$ | $-0.0645\,[54]$ | $-0.1316\,[5.9]$ | (1.07,1.14) |  | $-0.3848$ | $-0.1797\,[53]$ | $-0.3029\,[21]$ | (1.05,1.64)
$0.1$ |  | 1.9452 | 1.9804 [1.8] | 1.9489 [0.2] | (1.55,1.53) |  | 1.8486 | 1.8820 [1.8] | 1.8568 [0.4] | (1.48,1.77)
$1$ |  | 4.9392 | 5.2278 [5.8] | 4.9687 [0.6] | (1.79,1.77) |  | 3.4379 | 3.6386 [5.8] | 3.4753 [1.1] | (1.74,1.88)
$2$ |  | 7.5730 | 7.5730 [0] | 7.5730 [0] | (2,2) |  | 4.3729 | 4.3729 [0] | 4.3729 [0] | (2,2)
$3$ |  | 9.7389 | 8.9925 [7.7] | 9.6703 [0.7] | (2.16,2.20) |  | 5.0166 | 4.6320 [7.7] | 4.9693 [0.9] | (2.20,2.15)

The accuracy obtained by the ET bounds is better than 8%, except for the
Coulomb case. Except for one case, the Coulomb interaction with $m=5$, the
improvement is again quite large with the IET, even if its magnitude seems
quite unpredictable. The parameters $\phi_{a}$ and $\phi_{b}$ have a smooth
variation with $\beta$ and have values close to $\sqrt{2+\beta}$, the value
(16) for all identical particles. It is a little bit frustrating to note that
the results obtained here for the IET are slightly less good that for the
naive choice $\phi_{a}=\phi_{b}=\sqrt{2+\beta}$ made in sema20 . But this
trial had no a firm justification, contrary to the calculations developed
here.

Hamiltonians (28) and (29) have both interacting parts of power law form.
Results from Tables 1, 2 and 3 shows that very good results can be obtained
when both exponent are positive. When the potential is singular, the accuracy
is much lower. In Table 1, the comparison between $\beta=-1$ and $\beta=-0.5$
indicates that the accuracy for systems with identical particles is
deteriorated when the singularity increases. Nevertheless, the results are
remarkably good in some situations. Further comments are given in Sec. V.

###  IV.3 Atoms

In atomic units, the Hamiltonian for $N_{e}$ electrons in an atom of charge
$Z$ is written as

$H=\frac{1}{2}\sum_{i=1}^{N_{e}}\bm{p}_{i}^{2}+\frac{1}{2m}\bm{p}_{N}^{2}-Z\sum_{i=1}^{N_{e}}\frac{1}{|\bm{r}_{iN}|}+\sum_{i<j=2}^{N_{e}}\frac{1}{|\bm{r}_{ij}|},$
(30)

where $N=N_{e}+1$ is the number of the nucleus with a mass $m$. Energies in eV
are obtained by multiplying the eigenvalues of $H$ by the usual factor
$\alpha^{2}\,m_{e}=27.21$ eV. Hamiltonian (30) contains the main contributions
to the binding energy in an atom. So, the approximate results are compared
with the experimental data about ionisation energies nist  which are very
close to the eigenvalues of (30). Due to the mixing of attractive and
repulsive parts, the approximate eigenvalues have no variational character.
Taking into account the fermionic nature of the electrons (see Appendix B),
ground state binding energies for some atoms computed with the ET and the IET
are presented in Table 4 with the experimental values. Note that the two
electrons in 4He, 6Li+, 12C4+ and 16O6+ can be treated as bosons because of
the double-degeneracy due to the spin.

Table 4: Ground state binding energies (in eV) of Hamiltonian (30) for some atoms with two or more electrons. Results from the ET and the IET with ($\phi_{a},\phi_{b}$) are compared with the experimental values (see text). | Exp. | ET | IET | $(\phi_{a},\phi_{b})$
---|---|---|---|---
4He | 79 | 33 | 47 | $(1.21,1.78)$
6Li+ | 198 | 85 | 123 | $(1.18,1.77)$
6Li | 203 | 66 | 95 | $(1.03,1.99)$
12C4+ | 882 | 386 | 568 | $(1.16,1.77)$
12C | 1030 | 321 | 496 | $(0.96,2.10)$
16O6+ | 1611 | 707 | 1047 | $(1.15,1.77)$
16O | 2044 | 672 | 1062 | $(0.94,2.14)$

It is clear that the poor quality of the ET energies are not significantly
improved by the IET. Curiously, the naive choice $\phi_{a}=\phi_{b}=1$ made in
sema20  improves very much the binding energies but only when there is only
two electrons in the atom. Even if it produces better results in some cases,
this trial had no a firm justification, contrary to the calculations developed
here. We identify three possible sources for this problem: the action of the
Pauli principle; the mixing of attractive and repulsive parts in the
potentials, forbidding the existence of a definite variational character for
the solutions; the singularity of the Coulomb interaction.

The problem could arise from the fermionic nature of the electrons, but
results are not better when the two only electrons are treated as bosons.
Moreover, good results are obtained for the fermionic ground state of the one-
dimensional Calogero model with linear forces supplemented by inverse-cube
forces sema19 . So, the origin of the problem should be found elsewhere.

Some preliminary calculations with different three-body systems have been
undertaken to clarify the situation. Nonrelativistic systems with different
power law potentials have been considered. The results from the ET and the IET
have been compared not with very accurate numerical calculations but with the
HOB0 approximation described in sema20 , that is to say a variational method
based on a trial three-body state built with two Gaussian functions with
adjustable sizes. We know from our previous results that this method can give
very satisfactory upper bounds. It appears that the mixing of attractive and
negative interactions is not the cause of the poor results but instead the
presence of singular potentials. In this last situation, the binding energy
computed by the ET is then always of bad quality and the IET seems to make up
the situation only for identical particles sema15a ; sema15b . Supplementary
studies are necessary to check the relevance of these first results.

##  V concluding remarks and prospects

The envelope theory (ET) is a remarkably easy to use method to obtain
approximate but reliable eigensolutions of some quantum many-body systems
hall80 ; hall83 ; hall04 ; silv10 ; silv12 ; sema13 ; sema20 ; sema17 ;
sema18b ; sema19 . It is particularly useful when analytical or numerical
results with only a reasonable accuracy are asked for sema21 , as simple tests
for numerical codes sema15a ; char15 , or when a systematic study as a
function of the number of particles is necessary as in the Large-$N$ approach
of QCD buis11 ; buis12 . When all particles are identical for $D\geq 2$, the
non-physical strong degeneracy inherent to this method can be partly removed
by combining it with the dominantly orbital state method leading to what we
called the improved envelope theory (IET). This IET can then provide improved
eigenvalues sema15b , with large effects in some cases. That is the reason why
this IET is generalised in this paper to systems with identical particles plus
a different one to open new domains of applicability of this method.

The new equations to solve are explicitly written for these systems, for
general kinetic parts and two-body central potentials. They have more complex
structures than the ones for all identical particles, which does not allow to
find analytical solutions, but they are not more complicated to solve
numerically. The procedure is quite straightforward and can be easily
implemented. Tests of accuracy are performed on three different systems. The
main finding is that the accuracy of the eigenvalues is always improved, but
the variational character can be lost. When the ET solution is very poor, the
IET cannot provide an enough large improvement leading to a reasonable
accuracy in some cases. Nevertheless, it seems that results obtained here
already show that the IET can be safely used for a large variety of systems.

The IET relies on the generalisation to $N$-body systems of the effective
quantum number $q$ (6) which gives good spectra for 2-body systems loba09 .
The results obtained here show that this idea works well, at least for some
Hamiltonians. But this is not a rigorous demonstration. As the determination
of $q$ is performed within the WKB approximation, the same kind of study could
be attempted for $N$-body systems. WKB computations are possible for such
systems azbe90 ; niss10 . A promising idea could be to perform an hyperradial
approximation of the considered Hamiltonian as in timo12 ; timo17  and solve
it with the WKB approach. In this way, it could be possible to obtain
analytical information about the quality of an effective quantum number for
$N$-body systems, and possibly improve the ET in another way.

It is mentioned in Sec. I, that the ET and the IET can be used for various
systems. A special type of three-body forces is used in clusters of cold atoms
to produce a strong short-range repulsion gatt11 ; gatt13 ; kiev20  or in
hadrons to take into account peculiarities of the colour interaction between
quarks ferr95 ; dmit01 ; pepi02 . As this kind of three-body force can be
handled by the ET sema18b , it could be interesting to try an extension within
the IET.

The ET can give estimations of critical coupling constants in nonrelativistic
$N$-body systems sema13 . Let us consider a set of identical particles with a
mass $m$ interacting via a two-body potential $V(r)=-g\ v(r)$ admitting only a
finite number of bound states, where $g$ is a positive quantity which has the
dimension of an energy and $v(x)$ a globally positive dimensionless function
vanishing at infinity. The approximate value of $g$ allowing for a bound state
with quantum numbers $\\{n_{i},l_{i}\\}$ is given by sema13

$g=\frac{1}{u^{2}\,v(u)}\frac{2}{N(N-1)^{2}}\frac{Q(N)^{2}}{m},$ (31)

where $Q(N)$ is given by (3) and where the number $u$ is the solution of
$2\,v(u)+u\,v^{\prime}(u)=0$. It is then possible to link the values of the
critical coupling constants between two sets with $N$ and $N+1$ particles. For
a bosonic ground state (39), the universal relation is

$g(N+1)=\frac{N}{N+1}g(N).$ (32)

For a fermionic ground state with a large number of particles (42), the
corresponding relation is

$g(N+1)=\left(\frac{N}{N+1}\right)^{(D-2)/D}g(N).$ (33)

These results can be connected with the Efimov physics naid17 . At the origin,
the Efimov effect is the existence of three-body bound states whereas the
interaction between the particles cannot bind two of them. But, it can be
extended to a larger number of particles, and could then be studied by the ET
and the IET for sets of identical or different particles. It is worth
mentioning that both methods could have a very limited interest for a small
number of particles for the Efimov physics. Indeed, the ET cannot find the
bound state of a few identical bosons interacting via a weak Gaussian
potential both at $D=1$ sema19  and $D=3$ sema15a , whereas the bound states
can be computed with a good accuracy for a large number of particles.

As a next step, we plan to perform a systematic study of the accuracy of the
ET and the IET for 3-body systems containing bosons and fermions with an
expansion of the eigenstates in large harmonic oscillator basis nunb77 ;
silv96 ; silv00 . This will allow for a careful analysis of spectra as the one
performed in sema15a  to better precise the validity of the ET and the IET. As
already mentioned in sema20 , it would be desirable to have available accurate
computations of eigenvalues for large many-body systems containing different
particles, bosons or fermions, as the ones performed in horn14  for identical
bosons, to test the ET and the IET on a larger scale.

###### Acknowledgements.

C.C. and L.C. would like to thank the Fonds de la Recherche Scientifique -
FNRS for the financial support. This work was also supported under Grant
Number 4.45.10.08. All authors would like to thank Jean-Marc Richard for
providing the accurate results in Table 2 and the anonymous referee for the
useful suggestions.

##  Appendix A Two coupled oscillators

Let us consider the Hamiltonian for two one-dimensional coupled oscillators

$H=\frac{1}{2}\left(\frac{1}{\mu_{a}}p_{1}^{2}+\frac{1}{\mu_{b}}p_{2}^{2}+k_{a}\,x_{1}^{2}+k_{b}\,x_{2}^{2}+k_{c}\,x_{1}\,x_{2}\right),$
(34)

where $p_{j}=-i\partial/\partial x_{j}$ with $j=1,2$. Using the procedure
described in maka18 , one can find

$E_{nn^{\prime}}=\sqrt{\frac{A}{\mu}}\left(n+\frac{1}{2}\right)+\sqrt{\frac{B}{\mu}}\left(n^{\prime}+\frac{1}{2}\right),\\\
$ (35a) where $\displaystyle\mu=\sqrt{\mu_{a}\,\mu_{b}},$
$\displaystyle\begin{cases}A=\sqrt{\frac{\mu_{b}}{\mu_{a}}}k_{a}\\\
B=\sqrt{\frac{\mu_{a}}{\mu_{b}}}k_{b}\end{cases}\quad\text{if}\quad k_{c}=0,$
$\displaystyle\begin{cases}A=\sqrt{\frac{\mu_{b}}{\mu_{a}}}k_{a}-\frac{k_{c}}{2}\\\
B=\sqrt{\frac{\mu_{b}}{\mu_{a}}}k_{a}+\frac{k_{c}}{2}\end{cases}\quad\text{if}\quad\epsilon=\frac{1}{k_{c}}\left(\sqrt{\frac{\mu_{a}}{\mu_{b}}}k_{b}-\sqrt{\frac{\mu_{b}}{\mu_{a}}}k_{a}\right)=0,$
(35b)
$\displaystyle\begin{cases}A=\sqrt{\frac{\mu_{b}}{\mu_{a}}}k_{a}-\frac{k_{c}}{2}\left(\textrm{sgn}(\epsilon)\sqrt{1+\epsilon^{2}}-\epsilon\right)\\\
B=\sqrt{\frac{\mu_{a}}{\mu_{b}}}k_{b}+\frac{k_{c}}{2}\left(\textrm{sgn}(\epsilon)\sqrt{1+\epsilon^{2}}-\epsilon\right)\end{cases}\quad\text{if}\quad\epsilon\neq
0,$

with

$\lim_{\epsilon\,\to\,0^{\pm}}\left(\textrm{sgn}(\epsilon)\sqrt{1+\epsilon^{2}}-\epsilon\right)=\pm
1.$ (36)

Let us note that the equivalent formulas given in maka18  are only valid if
$\mu_{a}=\mu_{b}$.

##  Appendix B Bosonic and fermionic Ground states

An approximate eigenstate obtained with the ET is an eigenstate of a many-body
harmonic oscillator Hamiltonian. So the symmetry and other characteristics of
the approximate solutions are given by the properties of the solutions of this
Hamiltonian. If we consider only the spatial symmetry of the many-body states,
it is shown that the internal energy of a system of $N$ identical oscillators,
interacting via two-body harmonic forces, is equivalent to the energy of a
system of decoupled oscillators in a common central field, provided the
contribution of the centre of mass is removed levy68 . Using this result, the
bosonic ground state (BGS) for the ET solution for $N$ identical particles is
trivially given by

$Q_{\textrm{BGS}}(N)=(N-1)\frac{D}{2},$ (37)

since all quantum numbers are vanishing. For fermions with degeneracy $d$, the
situation is more complicated. In this case, the fermionic ground state (FGS)
can be computed with the two following relations

$\displaystyle Q_{\textrm{FGS}}(N)$
$\displaystyle=d\,D\binom{q+D-1}{D+1}+q\,r+(N-1)\frac{D}{2}\quad\textrm{with}$
(38a) $\displaystyle N$ $\displaystyle=d\binom{q+D-1}{D}+r,$ (38b)

where $q$ is the greatest natural number such that $r\geq 0$ in (38b) sema20 .

How to modify these formulas when the number $\phi$ is introduced in the
global quantum number $Q$? We consider that the energy is still obtained using
the rules defined in levy68 . But in piling the particles on levels with
higher and higher energies, we must take into account that the presence of
$\phi$ modifies the order of the levels with respect to the usual harmonic
oscillator levels. For bosons, the problem is still trivial since all quantum
numbers can be put to zero

$Q_{\phi\,\textrm{BGS}}(N)=(N-1)\frac{D+\phi-2}{2}.$ (39)

For fermions, when $\phi=1$, it is possible to perform a calculation similar
to the one giving (38). The first step is to determine the degeneracy
$W_{l}^{D}$ of a level of orbital angular momentum $l$ in $D$ dimensions
($D\geq 2$) yane94

$W_{l}^{D}=\sum_{m_{1}=0}^{l}\sum_{m_{2}=0}^{m_{1}}...\sum_{m_{D-2}=-m_{D-3}}^{m_{D-3}}d=d\
\frac{2l+D-2}{D-2}\binom{l+D-3}{D-3},$ (40)

where the $\\{m_{i}\\}$ are the magnetic quantum numbers. This formula can be
demonstrated by induction. One can check that $W_{l}^{2}=d(2-\delta_{l0})$ and
that $W_{l}^{3}=d(2l+1)$, as expected. Then, the fermionic ground state
$Q_{\phi\,\textrm{FGS}}(N)$ for $\phi=1$ can be computed with the two
following relations

$\displaystyle Q_{1\,\textrm{FGS}}(N)$ $\displaystyle=d\
\frac{(2\,q\,D-2\,D+D^{2}+1)}{D+1}\binom{q+D-2}{D}+q\,r+(N-1)\frac{D-1}{2}\quad\textrm{with}$
(41a) $\displaystyle N$ $\displaystyle=d\
\frac{2\,q+D-2}{D}\binom{q+D-2}{D-1}+r,$ (41b)

where $q$ is the greatest natural number such that $r\geq 0$ in (41b).

When $\phi$ is arbitrary, we must resort to a numerical routine to determine
$Q_{\phi\,\textrm{FGS}}(N)$. We have checked that this routine is in agreement
with formulas (38) and (41). We have noticed that the following approximate
formula

$Q_{\phi\,\textrm{FGS}}(N)\approx\frac{D}{D+1}\left(\frac{\phi\,D!}{2\,d}\right)^{1/D}N^{(D+1)/D}$
(42)

is valid for small values of $D$, $d$ and $\phi$ when $N\gg 1$. Note that the
bosonic and fermionic ground states are the same if $d\geq N$.

## References

  * (1) R.L. Hall, Energy trajectories for the $N$-boson problem by the method of potential envelopes. Phys. Rev. D 22, 2062 (1980)
  * (2) R.L. Hall, A geometrical theory of energy trajectories in quantum mechanics. J. Math. Phys. 24, 324 (1983)
  * (3) R.L. Hall, W. Lucha, F.F. Schöberl, Relativistic $N$-boson systems bound by pair potentials $V(r_{ij})=g(r^{2}_{ij})$. J. Math. Phys. 45, 3086 (2004)
  * (4) B. Silvestre-Brac, C. Semay, F. Buisseret, F. Brau, The quantum $\mathcal{N}$-body problem and the auxiliary field method. J. Math. Phys. 51, 032104 (2010)
  * (5) B. Silvestre-Brac, C. Semay, F. Buisseret, The Auxiliary Field Method in Quantum Mechanics. J. Phys. Math. 4, P120601 (2012)
  * (6) C. Semay, C. Roland, Approximate solutions for $N$-body Hamiltonians with identical particles in $D$ dimensions. Results Phys. 3, 231 (2013)
  * (7) C. Semay, L. Cimino, C. Willemyns, Envelope Theory for Systems with Different Particles. Few-Body Syst. 61, 19 (2020)
  * (8) C. Semay, F. Buisseret, Bound Cyclic Systems with the Envelope Theory. Few-Body Syst. 58, 151 (2017)
  * (9) C. Semay, G. Sicorello, Many-Body Forces with the Envelope Theory. Few-Body Syst. 59, 119 (2018)
  * (10) C. Semay, L. Cimino, Tests of the Envelope Theory in One Dimension. Few-Body Syst. 60, 64 (2019)
  * (11) R.L. Hall, B. Schwesinger, The complete exact solution to the translation-invariant $N$-body harmonic oscillator problem. J. Math. Phys. 20, 2481 (1979)
  * (12) C.T. Willemyns, C. Semay, Some specific solutions to the translation-invariant $N$-body harmonic oscillator Hamiltonian. J. Phys. Commun. 5, 115002 (2021)
  * (13) A.A. Lobashev, N.N. Trunov, A universal effective quantum number for centrally symmetric problems. J. Phys. A 42, 345202 (2009)
  * (14) C. Semay, Numerical Tests of the Envelope Theory for Few-Boson Systems. Few-Body Syst. 56, 149 (2015)
  * (15) J. Horne, J.A. Salas, K. Varga, Energy and Structure of Few-Boson Systems. Few-Body Syst. 55, 1245 (2014)
  * (16) M.G. Olsson, Universal behavior in excited heavy-light and light-light mesons. Phys. Rev. D 55, 5479 (1997)
  * (17) C. Semay, F. Buisseret, Two- and three-body calculations within the dominantly orbital state method. Phy. Lett. A 377, 1826 (2013)
  * (18) C. Semay, Improvement of the envelope theory with the dominantly orbital state method. Eur. Phys. J. Plus 130, 156 (2015)
  * (19) E. Witten, Baryons in the $1/N$ Expansion. Nucl. Phys. B 160, 57 (1979)
  * (20) C. Willemyns, C. Schat, Operator analysis of effective spin-flavor interactions for $L=1$ excited baryons. Phys. Rev. D 93, 034007 (2016)
  * (21) F. Buisseret, C. Semay, Light baryon masses in different large-$N_{c}$ limits. Phys. Rev. D 82, 056008 (2010)
  * (22) F. Buisseret, N. Matagne, C. Semay, Spin contribution to light baryons in different large-$N$ limits. Phys. Rev. D 85, 036010 (2012)
  * (23) M. Gattobigio, A. Kievsky, M. Viviani, Spectra of helium clusters with up to six atoms using soft-core potentials. Phys. Rev. A 84, 052503 (2011)
  * (24) M. Gattobigio, A. Kievsky, M. Viviani, Six-Bodies Calculations Using the Hyperspherical Harmonics Method. Few-Body Syst. 54, 657 (2013)
  * (25) A. Kievsky, A. Polls, B. Juliá-Díaz, N.K. Timofeyuk, M. Gattobigio, Few bosons to many bosons inside the unitary window: A transition between universal and nonuniversal behavior. Phys. Rev. A 102, 063320 (2020)
  * (26) C. Semay, C.T. Willemyns, Quasi Kepler’s third law for quantum many-body systems. Eur. Phys. J. Plus 136, 342 (2021)
  * (27) C. Semay, Three theorems of quantum mechanics and their classical counterparts. Eur. J. Phys. 39, 055401 (2018)
  * (28) Fl. Stancu, _Group Theory in Subnuclear Physics_. Oxford: Oxford University Press 1996
  * (29) J.L. Basdevant, A. Martin, J.M. Richard, Improved bounds on many-body Hamiltonians (II). Baryons from mesons in the quark model. Nucl. Phys. B 343, 69 (1990)
  * (30) J.M. Richard, From mesons to baryons. Phys. Lett. B 100, 515 (1981)
  * (31) L. Cimino, C. Semay, Compact Equations for the Envelope Theory. Braz. J. Phys. 52, 45 (2022)
  * (32) N.K. Timofeyuk, Convergence of the hyperspherical-harmonics expansion with increasing number of particles for bosonic systems. Phys. Rev. A 86, 032507 (2012)
  * (33) N.K. Timofeyuk, D. Baye, Hyperspherical Harmonics Expansion on Lagrange Meshes for Bosonic Systems in One Dimension. Few-Body Syst. 58, 157 (2017)
  * (34) I. Salom, V. Dmitrašinović, Relativistic Three-Body Harmonic Oscillator. In: V. Dobrev (ed) _Lie Theory and Its Applications in Physics_ 473, Springer Proceedings in Mathematics & Statistics 335, Springer Singapore (2020)
  * (35) A. Kramida, Y. Ralchenko, J. Reader, NIST ASD Team, _NIST Atomic Spectra Database (ver. 5.7.1)_. (National Institute of Standards and Technology, Gaithersburg, MD, 2019), https://physics.nist.gov/asd. Accessed 1 Nov 2021
  * (36) Y. Chargui, A. Dhahbi, A. Trabelsi, Exact analytical treatment of the asymmetrical spinless Salpeter equation with a Coulomb-type potential. Phys. Scr. 90, 015201 (2015)
  * (37) M. Ya Azbel, P.M. Platzman, Many-Body Aspects of the Tunneling of Electrons from a Helium Surface. Phys. Rev. Lett. 65, 1376 (1990)
  * (38) F. Nissen, J. Keeling, Wentzel-Kramers-Brillouin approach and quantum corrections to classical dynamics in the Josephson problem. Phys. Rev. A 81, 063628 (2010)
  * (39) M. Ferraris, M.M. Giannini, M. Pizzo, E. Santopinto, L. Tiator, A three-body force model for the baryon spectrum. Phys. Lett. B 364, 231 (1995)
  * (40) V. Dmitrašinović, Cubic Casimir operator of SU${}_{\textrm{C}}$(3) and confinement in the nonrelativistic quark model. Phys. Lett. B 499, 135 (2001)
  * (41) S. Pepin, Fl. Stancu, Three-body confinement force in hadron spectroscopy. Phys. Rev. D 65, 054032 (2002)
  * (42) P. Naidon, S. Endo, Efimov physics: a review. Rep. Prog. Phys. 80, 056001 (2017)
  * (43) P. Nunberg, D. Prosperi, E. Pace, An application of a new harmonic-oscillator basis to the calculation of trinucleon ground-state observables. Nucl. Phys. A 285, 58 (1977)
  * (44) B. Silvestre-Brac, Spectrum and Static Properties of Heavy Baryons. Few-Body Syst. 20, 1 (1996)
  * (45) B. Silvestre-Brac, R. Bonnaz, C. Semay, F. Brau, Quantum three-body problems using harmonic oscillator bases with different sizes. Internal Report ISN-00-66 (2000) [arXiv:2003.11028]
  * (46) D.N. Makarov, Coupled harmonic oscillators and their quantum entanglement. Phys. Rev. E 97, 042203 (2018)
  * (47) J.M. Lévy-Leblond, Generalized uncertainty relations for many-fermion system. Phys. Lett. A 26, 540 (1968)
  * (48) R.J. Yáñez, W. Van Assche, J.S. Dehesa, Position and momentum information entropies of the $D$-dimensional harmonic oscillator and hydrogen atom. Phys. Rev. A 50, 3065 (1994)