The auxiliary field method in quantum mechanics

The auxiliary field method is a new technique to obtain closed formulae for the solutions of eigenequations in quantum mechanics. The idea is to replace a Hamiltonian $H$ for which analytical solutions are not known by another one $\tilde H$, including one or more auxiliary fields. For instance, a potential $V(r)$ not solvable is replaced by another one $P(r)$ more familiar, or a semirelativistic kinetic part is replaced by an equivalent nonrelativistic one. The approximation comes from the replacement of the auxiliary fields by pure real constants. The approximant solutions for $H$, eigenvalues and eigenfunctions, are then obtained by the solutions of $\tilde H$ in which the auxiliary parameters are eliminated by an extremization procedure for the eigenenergies. If $H=T(\bm p)+V(r)$ and if $P(r)$ is a power law, the approximate eigenvalues can be written $T(p_0)+V(r_0)$, where the mean impulsion $p_0$ is a function of the mean distance $r_0$ and where $r_0$ is determined by an equation which is linked to the generalized virial theorem. The general properties of the method are studied and the connections with the envelope theory presented. This method is first applied to nonrelativistic and semirelativistic two-body systems, with a great variety of potentials. Closed formulae are produced for energies, eigenstates, various observables and critical constants, with sometimes a very good accuracy. The method is then used to solve nonrelativistic and semirelativistic many-body systems with one-body and two-body interactions. For such cases, analytical solutions can only be obtained for systems of identical particles, but several systems of interest for atomic and hadronic physics are studied. General results concerning the many-body critical constants are presented, as well as duality relations existing between approximate and exact eigenvalues.

The auxiliary field method is a new technique to obtain closed formulae for the solutions of eigenequations in quantum mechanics. The idea is to replace a Hamiltonian H for which analytical solutions are not known by another oneH, including one or more auxiliary fields, for which they are known. For instance, a potential V (r) not solvable is replaced by another one P(r) more familiar, or a semirelativistic kinetic part is replaced by an equivalent nonrelativistic one. If the auxiliary fields are eliminated by an extremization procedure, the HamiltonianH reduces to Hamiltonian H. The approximation comes from the replacement of these fields by pure real constants. The approximant solutions for H, eigenvalues and eigenfunctions, are then obtained by the solutions ofH in which the auxiliary parameters are eliminated by an extremization procedure for the eigenenergies, which takes the form of a transcendental equation to solve. If H = T (p p p) +V (r) and if P(r) is a power law, the approximate eigenvalues can be written T (p 0 ) +V (r 0 ), where the mean impulsion p 0 is a function of the mean distance r 0 and where r 0 is determined by an equation which is linked to the generalized virial theorem. The general properties of the method are studied and the connections with the envelope theory presented. Its mean field and (anti)variational characters are also discussed. This method is first applied to nonrelativistic and semirelativistic two-body systems, with a great variety of potentials (sum of power laws, logarithm, exponential, square root). Closed formulae are produced for energies, eigenstates, various observables and critical constants (when it is relevant), with sometimes a very good accuracy. The method is then used to solve nonrelativistic and semirelativistic many-body systems with one-body and two-body interactions. For such cases, analytical solutions can only be obtained for systems of identical particles, but several systems of interest for atomic and hadronic physics are studied. General results concerning the many-body critical constants are presented, as well as duality relations existing between approximate and exact eigenvalues. The aim of this work is to present the auxiliary field method (AFM), which is a new and remarkably simple method to find analytical approximate solutions of eigenequations corresponding to Hamiltonians admitting bound states. The quest for exactly solvable problems in quantum physics is actually as old as quantum mechanics itself, and it is worth making general comments about this topic before focusing on the auxiliary field method.
The most famous systems for which the Schrödinger equation admits analytical solutions are probably the harmonic oscillator and the hydrogen atom, i.e. Hamiltonians of the form p p p 2 + r 2 and p p p 2 − 1/r respectively, where p p p and r r r are conjugate variables. The spectrum of the one-dimensional harmonic oscillator is actually known since 1925 thanks to Heisenberg's pioneering work, 1) while the spectrum of the nonrelativistic hydrogen atom was found in Schrödinger's paper 2) as a first and crucial test of his celebrated equation. Notice that the same problem was solved the same year by Pauli. 3) Since those early results, there has been a considerable amount of works devoted to the computation of analytical solutions of the Schrödinger equation, especially in bound state problems. The Schrödinger equation has been found to be exactly solvable in many one-dimensional cases: Dirac comb, exponential and linear potentials, . . . However, only a few analytical three-dimensional solutions are known for any value of the angular momentum (the S-wave channel is very similar to a simpler one-dimensional equation). Besides the aforementioned harmonic oscillator and Coulomb cases, one can mention the Kratzer potential, the particle in a box and the symmetric top, which are of quite limited use in modeling realistic quantum systems. An exhaustive discussion of quantum problems admitting analytical solutions would be out of the scope of the present introduction; therefore we refer the interested reader to the textbooks, 4), 5), 6), 7), 8), 9) in which a wide range of problems is covered.
When exact analytical solutions cannot be found, physicists can either resort to numerical computations or to methods leading to approximate analytical solutions. Even if the computational power that we have at our disposal nowadays allows to solve accurately nearly every eigenequation for few body systems, finding approximate analytical formulae is always useful in physics, not only because of the intrinsic mathematical interest of such a task, but also in view of obtaining informations about the dependence of the observables on the parameters of a model and on the quantum numbers of the states. Calculating an analytical expression being always much less time consuming than solving numerically the corresponding eigenequation, the use of closed formulae can be a great advantage when one tries to fit the parameters of a given model to some experimental data, or even to interpret new data. We find worth mentioning here three approximation schemes that are all "textbook material" and able to lead to analytic results: the WKB or semiclassical approximation, the variational method, and the perturbation theory.
The principle of the WKB approximation has been found by Wentzel, Kramers and Brillouin around 1926, 10) but it is also known that Jeffreys brought important contributions to the field in 1924. 11) It basically consists in a Taylor expansion of the Hamiltonian and the wavefunction in powers ofh. At dominant order for a one-dimensional problem, such an expansion leads for example to the Bohr-Sommerfeld quantization rule a(E) b(E) p(E, x) dx = πh(n + 1/2), where p is the classical momentum expressed as a function of the energy and of the position x, and where a(E) and b(E) are the turning points of the classical trajectory. The energy spectrum can then be extracted from this integral and is generally accurate, but the explicit expression for E(n) is not always possible to find. One has instead often to numerically compute the energy from the Bohr-Sommerfeld quantization rule. Further informations about the semiclassical approximation can be found in e.g. Refs. 4), 7).
While the WKB method is particularly efficient for highly excited states, whereh clearly becomes a small parameter, the variational method due to Ritz 12) allows to gain relevant informations about the ground state (for a modern presentation see Ref. 13)). Indeed, let us consider a Hamiltonian H and a normalized arbitrary function ψ belonging to L 2 space. Then, the Ritz theorem states that E 0 ≤ ψ|H|ψ , where E 0 is the exact ground-state energy. Provided that ψ is chosen such that ψ|H|ψ is analytical, one can then obtain an analytical upper bound on the ground-state energy. Typically, one chooses a gaussian-or exponential-type wavefunction; many examples of analytical upper bounds can be found in Ref. 9). Notice also the MacDonald's theorem stating that, given a set of n orthonormalized functions ψ i belonging to L 2 space, the eigenvalues of the matrix ψ i |H|ψ j provide upper bounds on the n first eigenenergies of H. 14) However, analytical results can hardly be obtained from such a generalization due to the complexity of the underlying calculations.
A last technique to be mentioned is the perturbation theory. 15) Let us consider that the Schrödinger equation with a potential V is exactly solvable. Then the Schrödinger equation with potential V + ξ w can be solved thanks to a Taylor expansion in ξ provided that this parameter is small. At first order in particular, the computation of ξ w with the lowestorder states provides the dominant correction to the whole energy spectrum. Again, a wide range of quantum problems can be studied by using perturbation theory, and we refer to e.g. 5), 9) for explicit examples in which w is analytical.
Those three methods being presented, it appears that a major challenge remains: To find approximate analytical solutions of eigenequations for the complete spectrum (not only the ground state as with the variational method) without using any Taylor expansion as done with the WKB method or the perturbation theory. The AFM, that will be presented in detail in the rest of this work, is an attempt to address that issue. §2. Generalities on the auxiliary field method

Historical aspects
Auxiliary fields, also known as einbein fields, are used in various domains of theoretical physics. Historically, they have been introduced in order to get rid of the square roots typically appearing in relativistic Lagrangians. 16), 17) The most obvious example is the free relativistic particle, described by the Lagrangian L = −m √ẋ 2 , whereẋ is the worldvelocity. An auxiliary field µ can be introduced in this expression so that one is led to the more convenient form L (µ) =ẋ 2 /(2µ) + m 2 µ/2. This last expression is formally simpler than Lagrangian L , but it is equivalent if the equations of motion of the auxiliary field are considered. The auxiliary field is indeed given on shell by µ 0 = − √ẋ 2 /m, and it is readily checked that L (µ 0 ) = L . That feature can be seen as a definition of an auxiliary field: It is a field whose equation of motion is not dynamical but leads to algebraic relations that allow to express it in terms of the other degrees of freedom of the problem.
Besides pointlike particles, auxiliary fields have soon become a cornerstone of string theory because of the systematic replacement of the Nambu-Goto Lagrangian L NG = The potential V in the Hamiltonian H depends on the variable r = |r r r| which is the distance from the center of forces in a one-body problem and the relative distance between the two particles in a two-body problem. The kinetic part T depends on the momentum operator p p p which is the conjugate variable of r r r. Practically, we will consider a nonrelativistic form p p p 2 /(2m) (m is the particle mass for a one-body problem or the reduced mass of a two-body system) as in the Schrödinger equation and a semirelativistic one σ p p p 2 + m 2 (σ = 1 for one-body problems and σ = 2 for two identical particles) as in a spinless Salpeter equation.
An exact analytical expression for all the eigenvalues is known explicitly only for very specific potentials P(r) with a nonrelativistic kinetic parts: the quadratic interaction P(r) = r 2 (harmonic oscillator) and the Coulomb potential P(r) = −1/r (hydrogen-like system) are the most familiar for a practical use. An exact expression is known for a number of other potentials, but only for S-waves (see Ref. 9) for a more detailed discussion) or for one-dimensional problems.
Our goal is the search for approximate analytical solutions for the Hamiltonian of type (2 . 1), relying on Hamiltonians for which solutions are well-known. In other words, we assume that we are able to obtain an analytical solution for the equation h(a)|a = [T (p p p) + aP(r)] |a = e(a) |a , (2 . 2) in which, at this stage, a is a real parameter. This method was formulated and applied for the first time in Ref. 38). We summarize here the principle of the method introducing an auxiliary field, which is a priori an operator. It consists in four steps: 1. We calculate the function K(r) (the prime denotes the derivative with respect to the argument of the function) We denote by J = K −1 the inverse function of K. Thus, one has K(r) = z and J(z) = r. (2 . 4) defined asH (y) = T (p p p) +Ṽ (r, y) = h(y) + B(y).
(2 . 7) The construction ofH seems rather artificial and one can wonder what could be the justification of such a procedure. In fact, there is no mystery. Let us suppose that the y variable is an operator denotedν. It plays the role of an auxiliary field, not present in the original Hamiltonian H, but which is an essential ingredient of the AFM HamiltonianH. Among the infinite number of possibilities for the auxiliary fieldν, let us choose a very specific oneν 0 defined bŷ ν 0 = K(r).
(2 . 8) A more explicit notation would beν 0 (r) indicating clearly that it is an operator depending on the position only. The very important property is thatν 0 , coming from (2 . 8), is an extremum the AFM Hamiltonian, i.e. δH(ν) δν ν=ν 0 = 0. Thus, considering the auxiliary field as an operator and affecting it the value given by (2 . 8) is just an alternative method to solve the original problem. 4. The philosophy of the AFM is to consider the y variable, no longer as an operatorν, but as a pure real number ν. In this case the AFM potentialṼ (r, ν) is of the form νP(r) + B(ν) where ν and B(ν) are no longer operators but must be considered as arbitrary constants. Taking into account (2 . 2), the eigenvalues ofH(ν) are E(ν) = e(ν) + B(ν), (2 . 11) where e(ν) are the eigenvalues of h(ν) which are supposed to be known whatever the radial and orbital quantum numbers n, l. Provided that the J function is calculable, the same property is true for E(ν). Then, we determine the value ν 0 that extremizes E(ν): We propose to consider E(ν 0 ) as the approximate form of the exact energy E of the Hamiltonian H E ≈ E AFM = E(ν 0 ). (2 . 13) In order to obtain an analytical expression for the eigenvalues E AFM , we must then fulfill a second necessary condition: to be able to determine ν 0 and, then, E(ν 0 ) in an analytical way. Denoting |ν the eigenstate ofH(ν) corresponding to the eigenvalue E(ν), we consider |ν 0 as an approximation of the corresponding genuine eigenstate of H.
Presented as such, this procedure appears to be an empirical recipe. Arguments given in the third item show that this prescription makes sense and that the AFM must be considered in essence as a mean field approximation. This point will be developed more deeply below.
Another interesting property concerning the AFM potential (2 . 6) is the following. Let us denote by r 0 the value of the radius defined by ν 0 = K(r 0 ) and r 0 = J(ν 0 ).
(2 . 14) It is just a matter of simple calculation to verify that the AFM potential taken for the special value y = ν 0 can be recast asṼ (r, ν 0 ) = ν 0 [P(r) − P(r 0 )] +V (r 0 ). (2 . 15) Owing to the definition of the K function (2 . 3), one can easily show that the AFM potential coincides both to the exact potential and its derivative at the point r = r 0 , namelỹ V (r 0 , ν 0 ) = V (r 0 ) andṼ ′ (r 0 , ν 0 ) = V ′ (r 0 ). (2 . 16) In other words, at this particular point r 0 , the AFM potential and the exact potential are tangent curves. This property will be exploited in Sect. 2.4 in connection with the envelope theory (see Appendix A). Let us note by |ν an eigenstate ofH(ν). The Hellmann-Feynman theorem 39) states Using this relation, one can show that 38) ν 0 | P(r) |ν 0 = P(r 0 ). (2 . 18) This means that r 0 is a kind of "average point" for the potential P(r). That is why it will be often called "the mean radius" in the following. Using this last relation with the definitions above, we get ν 0 |Z(ν 0 )|ν 0 = Z(ν 0 ) with Z(x) = P(J(x)). (2 . 19) So, our method can actually be considered as a "mean field approximation" with respect to a particular auxiliary field which is introduced to simplify the calculations: ν 0 is the mean value of the operatorν 0 = K(r) through a function Z which can be quite simple. For example, It is in the step of passing from an operator to a constant interpretation for the auxiliary field that lies the approximation. It is clear that the quality of the approximate results strongly depends on the choice of the function P(r). The cleverness of the physicist relies in his ability to guess a form for P(r) as close as possible of V (r) while leading to a manageable form of the corresponding eigenvalues e(ν). A big bulk of this report is devoted to the discussion of the quality of this approximation.
This method is completely general and a priori valid for any potential V (r). Of course, all the four steps mentioned above can be done numerically so that E AFM can be computed numerically. But it seems that the corresponding numerical treatment is as heavy as solving directly the eigenvalue equation to get lastly a poorer result. This is true if we start with an arbitrary function P(r). However we have the complete freedom for this choice.
Let us assume that, for a kinetic part T , we do start with a potential P allowing to obtain analytical results for its eigenvalues e(ν) (e(ν) is a shorthand notation for the most correct expression e(ν; n, l)); the same is true of course for its derivative with respect to ν : e ′ (ν). Using the definition of the K and J functions defined previously, it is easy to show that E ′ (ν) = e ′ (ν) − P(J(ν)).
(2 . 20) The determination of ν 0 results from the condition E ′ (ν 0 ) = 0, that is to say e ′ (ν 0 ) = P(J(ν 0 )). (2 . 21) It appears that obtaining the ν 0 value requires only to solve a transcendental equation. Once this value is obtained, the AFM energy is very easy to calculate E AFM = E(ν 0 ) = e(ν 0 ) + B(ν 0 ). (2 . 22) Thus, provided that we have an analytical expression for e(ν), the AFM results only need solving a transcendental equation, a much easier procedure than solving exactly the genuine eigenvalue equation.
The practical usefulness of the AFM method is to get a final solution which is completely analytical. In order to do that, we must fulfill, as we saw, three conditions: 1. to choose a P function leading to analytical expression for e(ν); 2. to be able to invert relation (2 . 8) in order to have access to the function J defined by (2 . 4) (it appears into the B function); 3. to be able to determine ν 0 (through (2 . 21)) and to calculate the corresponding value E(ν 0 ) in an analytical way. In fact, the second condition is formally not necessary. Indeed, let us introduce the mean radius r 0 = J(ν 0 ). Then ν 0 = K(r 0 ) so that the bridge function is expressed as where, for simplicity, we still maintain the notation B: B = V − KP. The transcendental equation (2 . 21) is transformed into e ′ (K(r 0 )) = P(r 0 ), (2 . 24) while the energy is given by Using the expression (2 . 23) for B(r 0 ) and the definition of e(ν), this last formula can be recast under the form E AFM = T 0 +V (r 0 ), (2 . 26) where, according to (2 . 18), is an average kinetic energy. In this new formulation, any reference to the J function has disappeared so that it is not necessary to try to get it. However, solving analytically the transcendental equation (2 . 24) is more or less of the same difficulty than getting an analytical expression for the J function. Anyhow, this new formulation leads very often to simpler practical calculations and is preferred most of the time.
It is important to stress the following point: if P(r) is chosen to be V (r), it is trivial to check that K(r) = 1, B(r) = 0, r 0 is meaningless since obviously P(r) and V (r) are tangent curves everywhere, and E AFM = e(1) which is precisely the exact value E. In this particular case, AFM recovers the exact result. This property is sometimes used to check special results.

Scaling laws
Scaling laws represent an important property for quantum mechanical systems. They allow to give the expression for the eigenenergies (and wavefunctions) of the most general equation in terms of the corresponding eigenenergies (and wavefunctions) of a reduced equation with less parameters. In fact, the scaling laws are nothing else than a direct consequence of dimensional analysis applied to the various dimensioned parameters entering the problem. Starting from a general Hamiltonian H(α 1 , . . . , α n ) depending on n dimensioned parameters α i , it is generally possible to write where γ has the dimension of an energy and where h is a Hamiltonian expressed in terms of dimensionless conjugate variables and depending on m < n dimensionless parameters β i . The scaling properties of the nonrelativistic Schrödinger equation have been studied and used for a long time, but spinless Salpeter equations also benefit of scaling properties. They will be explicitly described below. In the following, in order to lighten the notations, the AFM will be generally applied to dimensionless Hamiltonians as in Refs. 38), 40), 41). So, the exact and the AFM solutions automatically share the same scaling properties.

Eigenstates and upper/lower bounds
If E(ν 0 ), which is an eigenvalue ofH(ν 0 ), is an approximation of the exact energy E, the corresponding eigenstate |ν 0 is an approximation of a genuine eigenstate of H. The shape of the corresponding wavefunction r r r|ν 0 depends on the quantum numbers via the parameter ν 0 or the mean radius r 0 . Practically, only nonrelativistic harmonic oscillator or nonrelativistic Coulomb wavefunctions can be used. Explicit examples will be presented below.
Since the potentialṼ (r, ν 0 ) is tangent to the potential V (r) at r = r 0 , and since H and H(ν 0 ) have the same kinetic part, the comparison theorem 42), 43), 44) (for both nonrelativistic and relativistic equations) implies that the approximation E(ν 0 ) is an upper (lower) bound on the exact energy ifṼ (r, ν 0 ) ≥ V (r) (Ṽ (r, ν 0 ) ≤ V (r)) for all values of r. Equivalently, a function g(x) can be defined by V (x) = g(P(x)). (2 . 29) It can then be shown that, if g(x) is a concave (convex) function, that is if g ′′ (x) ≤ 0 (g ′′ (x) ≥ 0) ∀x, the approximation E(ν 0 ) is an upper (lower) bound on the exact energy. This property has been demonstrated in the framework of the envelope theory, 45) but can be applied as well to the AFM. 46) Several examples will be presented below. The knowledge of lower and upper bounds on an eigenstate is a first technique to estimate the accuracy of the AFM. It has also been shown that 38) The right-hand side of this equation is the difference between the value of potential V computed at the average point r 0 and the average of this potential for the AFM state |ν 0 considered here as trial state. In some favorable cases (the trial state is a ground state for instance), E ≤ ν 0 |H|ν 0 and a bound on the error can be computed by If the mean value of V can be computed analytically, this constitutes a second procedure to estimate the accuracy of the AFM. Several examples of this calculation are presented in Ref. 38). At last, the eigenstates of a Hamiltonian of type (2 . 1) can be solved numerically with an arbitrary precision. So, as a third possibility, a direct comparison with the AFM results is always possible. So, one can wonder why to use the AFM? Let us recall that the interest of this method is mainly to obtain analytical information about the whole spectra (dependence of eigenenergies on the parameters of the Hamiltonian and on the quantum numbers), without necessarily searching a very high accuracy. Moreover, the AFM approximation can be extended to N-body problems for which exact eigenenergies are not easily reachable, even numerically. In some cases, it is possible to compute analytically the mean value E * (ν 0 ) = ν 0 |H|ν 0 , considering the AFM solution |ν 0 as a trial state (see Sect. 3.3 for two examples). We have then In this case, ∆V ≤ 0 and E * (ν 0 ) ≤ E(ν 0 ). By the Ritz theorem, we know that E * (ν 0 ) ≥ E for the ground state. For this state, it is interesting to compute the mean value E * (ν * ) = ν * |H|ν * , where |ν * is an eigenstate of H(ν * ) with ν * determined in order to minimize the energy E * (ν * ). We have then E * (ν * ) ≤ E * (ν 0 ) since the parameter ν 0 is fixed by the AFM computation.
We can gather these results to obtain: The results mentioned above are directly applicable if H andH(ν 0 ) have the same kinetic part. This is always the case for nonrelativistic Hamiltonians but not for the spinless Salpeter Hamiltonians, whose treatment requires generally the replacement of the square root operator by a nonrelativistic operator. The bounds on these kinds of Hamiltonians are specifically studied in Sect. 5.

Extension to the form aP(r) +V (r)
If the AFM is applicable for some potentials V 1 (r) and V 2 (r) independently, there is no certainty that it still applies for a potential which is their sum: V (r) = V 1 (r) + V 2 (r). However, the method can be used in the case where one of the potentials V i can be identified with the basic P potential. Thus, in this section we consider a Hamiltonian of type H a = T (p p p) + aP(r) +V (r), (2 . 33) whose exact eigenvalues are denoted E a . The extension of AFM to potentials of this type was first presented in Ref. 40). One introduces an auxiliary field ν as before, forgetting about the aP(r) contribution. The first 3 steps of the algorithm remain unchanged. Thus theν 0 field is the same, as is the same the function B(ν). The only difference arises in the expression (2 . 7) ofH and h where νP(r) has to be replaced by (a + ν)P(r). As a consequence, the corresponding energy (2 . 11) has to be replaced by where h is defined by (2 . 2) and where e(a + ν) is an eigenvalue of Hamiltonian h(a + ν). An eigenstate of HamiltoniansH a and h(a + ν) is denoted |a + ν , and we have e(a + ν) = a + ν|h(a + ν)|a + ν . If ν 0 is the value of ν which extremizes (2 . 34), then we could expect that is a good approximation of E a , an exact eigenvalue of Hamiltonian (2 . 33). It seems that (2 . 34) is very similar to (2 . 11), the only difference being the replacement of ν by a + ν in the argument of the function e(ν). Nevertheless, this small difference is important, because, even if the determination of ν 0 from (2 . 11) is technically easy and analytical, it may happen that its determination from (2 . 34) could be much more involved and very often not analytical. The alternative formulation in terms of r 0 is also slightly modified in this case. The transcendental equation now writes e ′ (a + K(r 0 )) = P(r 0 ) (2 . 37) and the AFM energy is given by Using again the Hellmann-Feynman theorem, 39) it can be shown that a + ν 0 |P(r)|a + ν 0 = P (J(ν 0 )) = P(r 0 ). (2 . 39) So, J(ν 0 ) = r 0 is also in this case a kind of "average point" for the potential P(r). Recalling that Z(x) = P(J(x)), we have We confirm also in this case the "mean field approximation" with respect to a particular auxiliary field which is introduced to simplify the calculations.

Recursion procedure
An idea for obtaining analytical expressions of the eigenenergies for an arbitrary potential is the following (see also Ref. 40)). We start with a potential P(r) = P [0] (r) for which the energies of the corresponding Hamiltonian H [0] are exactly known. We then proceed as above to find approximate solutions for the eigenenergies of a Hamiltonian H [1] in which the potential is at present V (r) = P [1] (r). In general, a large class of potentials can be treated in that way. Moreover, as we will see below, by comparison with accurate numerical results, we can even refine the expressions so that analytical forms for the energies are very close to the exact solutions.
Considering now these approximate expressions as the exact ones, we apply once more the AFM with P(r) = P [1] (r) to obtain approximate solutions for the eigenenergies of a Hamiltonian H [2] in which the potential is at present V (r) = P [2] (r). Even if analytical solutions for Hamiltonian H [2] were not attainable directly with P(r) = P [0] (r), it may occur that they indeed are with P(r) = P [1] (r). Pursuing recursively such a procedure, one can imagine to get analytical solutions for increasingly complicated potentials. Presumably, the quality of the analytical expressions deteriorates with the order of the recursion. §3. Schrödinger equation with a power-law potential In this section, we discuss in detail the case of the Schrödinger equation with a powerlaw potential for two reasons. First, it is a simple case to illustrate the AFM. Second, as we will see below, this kind of potentials is in the thick of our method. The Hamiltonian can be written where a is a strictly positive constant and where the sign function is defined by sgn(λ ) = λ /|λ | with λ = 0. The physical values of λ must be such that λ ≥ −2, otherwise the wave equation leads to a collapse. The two unavoidable starting potentials P(r) = r 2 and P(r) = −1/r are indeed two particular cases of power-law potentials. In fact, the only interesting values studied in this paper are those comprised between these extreme values: −1 ≤ λ ≤ 2. They represent most of the physical applications. For example, the linear potential λ = 1 is the traditional form of the confining potential in hadronic physics, 32), 47), 48) the value λ = 2/3 was shown to give the good slope for Regge trajectories in nonrelativistic treatment 49) and λ = 0.1 was considered by Martin. 50) The energy E of (3 . 1) depends on the physical parameters m, a, λ but also on the quantum numbers n (radial) and l (orbital). It is easy to show that the scaling laws allow to write where ε(λ ; n, l) is the eigenvalue of equation This form is chosen in order to match the conventions of Ref. 38). The scaling laws are unable to say anything about the dependence of ε on the power λ and on the quantum numbers (n, l). It is the virtue of the AFM to shed some light on these very important questions. An another interesting potential is the logarithmic potential, discussed for instance in Ref. 51) It is in strong connection with the power-law interaction since it can be rewritten into a similar form Using the the scaling law, we can write where ε(β ) is the eigenenergy of the reduced Schrödinger Hamiltonian We will apply the AFM for these two interactions with the two starting potentials P(r) = r 2 and P(r) = −1/r. Let us note that a number of results presented in the following sections have been already obtained in the framework of the envelope theory. 45) Nevertheless we remind them in order to keep some consistency in our presentation.
Inserting this value in E(r 0 ) (2 . 25) gives the AFM approximation of the eigenenergy It is important to emphasize that this value is an upper bound on the exact energy for the domain of λ we are interested in. Moreover, putting λ = 2 in the previous expression, one recovers the exact result ε AFM (2; n, l) = ε(2; n, l) = Q HO , as expected. Let us consider now the logarithmic potential V (r) = β ln r. We have then K(r) = β /(2r 2 ), e(K(r)) = Q HO β /2/r, e ′ (K(r)) = rQ HO The immediate consequence is that the corresponding spectrum is independent of the mass of the particle. It is remarkable that the basic property (3 . 13) still holds for our approximate expression (3 . 12).
It is worth mentioning that the energy formula (3 . 12) can be understood as a particular limit case of formula (3 . 2), as it was suggested in Ref. 47). Guided by (3 . 5), let us consider the following Hamiltonian (3 . 14) It reduces to the Hamiltonian (3 . 7) in the limit λ → 0. A simple rewriting of formula (3 . 2) for a → β /|λ | and m = 2 gives the eigenenergies of Hamiltonian (3 . 14). We have thus and as expected that is precisely formula (3 . 12), directly obtained from the logarithmic potential. This confirms the idea that the logarithmic potential can be seen as the limit of a power-law potential r λ when λ goes to zero. The same conclusion was obtained with the envelope theory (see Ref. 52)).
The resolution of the Schrödinger equation with the logarithmic potential in (3 . 7) is even simpler with the present choice of P(r). Indeed, one finds in this case K(r) = β r, e(K(r)) = −(rβ /Q C ) 2 , e ′ (K(r)) = −2β r 2 /Q 2 C and B(r) = V (r) + β . The transcendental equation (2 . 24) leads to the value of the mean radius Inserting this value in E(r 0 ) (2 . 25) gives the AFM approximation of the eigenenergy Again, formulae (3 . 12) and (3 . 21) exhibit the same formal expression with just the exchange Q C ↔ Q HO . Consequently, the properties mentioned above concerning the spectrum and the limit of a power-law potential when λ → 0 still hold in this case.

Wavefunctions and observables
In this section, the quality of the wavefunction r r r|ν 0 given by the AFM is tested for the Schrödinger Hamiltonian with the linear potential (results for logarithmic and exponential potentials are presented in Ref. 53)). This Hamiltonian is chosen because it can be solved analytically for a vanishing angular momentum l (see Appendix G). When l = 0, numerical results have been obtained with two different methods. 54), 55) In order to simplify the notation, we will denote |n, l an eigenstate of H and |n = |n, 0 , with r r r|n given by (G . 3). The corresponding energies will be denoted E n,l = n, l|H|n, l and E n = E n,0 = n|H|n . As we need a HamiltonianH with a central potential which is completely solvable to apply the AFM, we can only use in practice a hydrogen-like system (P(r) = −1/r) or a harmonic oscillator (P(r) = r 2 ). As a linear potential seems closer to r 2 than −1/r, we can expect that the use of a harmonic oscillator to start the AFM will give better results. Using the scaling properties, we can consider the following simple Hamiltonian H = p p p 2 + r, (3 . 22) in order to match the conventions of Ref. 53).
As it is shown in the previous section, the approximate AFM energies, denoted here by ε n,l , are given by with Q = 2n + l + 3/2 for P(r) = r 2 and Q = n + l + 1 for P(r) = −1/r. Exact energies, given by (G . 2), reduces to E n = −α n , where α n is the (n + 1) th zero of the Airy functions Ai. A simple approximation of E n can be obtained using the expansion (G . 4) at the first order It is worth noting that the linear potential is not only a toy model to test the AFM method. Effective theories of QCD have proved that it is a good interaction to take into account the confinement of quarks or gluons in potential models of hadronic physics. 32), 47), 48) 3.3.1. AFM with P(r) = −1/r One can ask whether it is possible to obtain good approximations for the solutions of a Schrödinger equation with a linear potential by means of hydrogen-like eigenfunctions. This will be examined in this section. Using the previous results for P(r) = −1/r, we find ν 0 = r 2 0 = 2 2/3 (n + l + 1) 4/3 . Exact eigenstates are then approximated by AFM eigenstates which are, in this case, hydrogen-like states (F . 3) with Such states are denoted |Hy; n, l and |Hy; n = |Hy; n, 0 . Using (F . 6) and results above, it can be shown that (2 . 18) is satisfied, with P(r 0 ) = −1/r 0 = −2 −1/3 (n + l + 1) −2/3 . It is also worth noting that (2 . 19) gives 1/ √ν = 1/ √ ν 0 . We denote ε Hy n,l and ε Hy n = ε Hy n,0 the approximated energies which are given by (3 . 23) with Q = n+ l + 1. SinceṼ (r, ν 0 )−V (r) = −(r − r 0 ) 2 /r ≤ 0, ε Hy n,l are lower bounds on the exact energies. This can also be determined with the function g defined by (2 . 29). In this case, g(y) = −1/y with y < 0. The function g ′′ (y) = −2/y 3 being positive, g is convex as expected.
(3 . 27) The ratio ε Hy n /E n is respectively equal to 0.808, 0.734, 0.712 for n = 0, 1, 2 and tends rapidly toward the asymptotic value 3 1/3 /π 2/3 ≈ 0.672. As expected, these ratios are smaller than 1 since ε Hy n are lower bounds. Two wavefunctions are given in Fig. 1. We can see that the differences between exact r r r|n and AFM r r r|Hy; n wavefunctions can be large. The overlap | n|Hy; n | 2 between these wavefunctions can be computed numerically with a high accuracy. We find respectively the values 0.934, 0.664, 0.298 for n = 0, 1, 2, showing a rapid decrease of the overlap. It is worth noting that, asymptotically, r r r|n ∝ exp − 2 3 r 3/2 while r r r|Hy; n is characterized by an exponential tail. Nevertheless, if an observable is not too sensitive to the large r behavior, this discrepancy will not spoil its mean value. The various observables |ψ n,0 (0)| 2 , r k and p k computed with the AFM states can be obtained using formulae (F . 4), (F . 6) and (F . 7) for l = 0 with the parameter η given by (3 . 25). Results are summed up in Table II. A direct comparison between the structure of exact and AFM observables can be obtained if we remind that the exact ones depend on |α n | which can be well approximated by β n (see (G . 4)). Let us look in detail only at the mean value r . For the exact and AFM solutions, we have respectively , (3 . 28) Hy; n|r|Hy; n = 3 2 2/3 (n + 1) 2/3 ≈ 1.890 (n + 1) 2/3 .
(3 . 29) Some observables like |ψ n,0 (0)| 2 and p 4 are very badly reproduced. Others can be obtained with a quite reasonable accuracy. Despite the fact that exact and AFM wavefunctions differ strongly when n increases, their sizes stay similar. Since it is possible to compute analytically the mean value Hy; n|p p p 2 + r|Hy; n for an arbitrary value of the scale factor η, the relations between the bounds given in Sect. 2.4 can be checked. We can verify that E = −α n ≥ E(ν 0 ) = 3 2 2/3 (n + 1) 2/3 and E * (ν 0 ) = 4 2 2/3 (n + 1) 2/3 ≥ E(ν 0 ) = 3 2 2/3 (n + 1) 2/3 . Moreover, for the ground state, we have E * (ν 0 ) = 4 . The behavior of observables computed with the AFM is similar for values of l = 0, 1, 2. We do not expect strong deviations for larger values of l. This is illustrated with some typical results gathered in Table III. Observables are generally not very well reproduced, but this is expected since eigenstates for a linear potential are very different from eigenstates for a Coulomb potential. Actually, the agreement is not catastrophic, except for |ψ n,0 (0)| 2 and p 4 as mentioned above. It is even surprising that the AFM with P(r) = −1/r could give energies and some observables for a linear potential with a quite reasonable accuracy. Moreover, lower bounds on the energies are obtained. Table III. Ratios between the AFM results (energies and r ) with P(r) = −1/r and the exact results, for several quantum number sets (n, l). l n = 0 n = 1 n = 2 n = 3 n = 4 n = 5 ε Hy n,l /E n,l 0 Since the quadratic potential is closer to the linear potential than a Coulomb one, one can expect better results with harmonic oscillator wavefunctions. This will be examined in this section. Using the previous results for P(r) = r 2 , we find ν 0 = 1/(2r 0 ) = 2 −4/3 (2n + l + 3/2) −2/3 . Exact eigenstates are then approximated by AFM eigenstates which are, in this case, harmonic oscillator states (E . 3) with Such states are denoted |HO; n, l and |HO; n = |HO; n, 0 . Using (E . 6) and results above, it can be shown that (2 . 18) is satisfied, with P(r 0 ) = r 2 0 = 2 2/3 (2n + l + 3/2) 4/3 . It is also worth noting that (2 . 19) gives 1/ν 2 = 1/ν 2 0 . This is in agreement with (20) in Ref. 41) (an auxiliary field φ = 1/ν is used in this last reference). We denote ε HO n,l and ε HO n = ε HO n,0 the approximated energies which are given by (3 . 23) n,l are upper bounds on the exact energies. In this case, g(y) = √ y with y > 0. The function g ′′ (y) = −1/(4y 3/2 ) being negative, g is concave as expected.
Since the scaling parameter λ depends on the quantum numbers, HO; n, l|HO; n, l ′ = δ ll ′ because of the orthogonality of the spherical harmonics, but HO; n, l|HO; n ′ , l = δ nn ′ .
Using the definition (H . 1), we find HO; n, l|HO; n ′ , l = F HO n,n ′ ,l with F HO given by (H . 4). Table IV gives some values of | HO; n, l|HO; n ′ , l | 2 . We can see that the overlap is always small and decreases rapidly with |n − n ′ |. The situation is even better when l increases: for l = 0 → 5, | HO; 0, l|HO; 1, l | 2 = 0.029, 0.023, 0.019, 0.017, 0.014, 0.013. Table IV. Results for P n,n ′ ,l = | HO; n, l|HO; n ′ , l | 2 . Values for l = 0 (l = 1) are given in the lower-left (upperright) triangle of the Table. P n,n ′ ,l = P n ′ ,n,l and P n,n,l = 1. By comparing with (3 . 24), we can see immediately that the situation is more favorable than in the previous case. The ratio ε HO n /E n is respectively equal to 1.059, 1.066, 1.067, for n = 0, 1, 2. The asymptotic value 3 1/3 2 2/3 /π 2/3 ≈ 1.067 is rapidly approached. As expected, these ratios are greater than 1 since ε HO n are upper bounds. Two wavefunctions are given in Fig. 1. We can see that the differences between exact r r r|n and AFM r r r|HO; n wavefunctions are quite small. The overlap | n|HO; n | 2 between these wavefunctions can be computed numerically with a high accuracy. We find respectively the values 0.997, 0.979, 0.951 for n = 0, 1, 2. A value below 0.75 is reached for n = 6 and below 0.25 for n = 14. A wavefunction r r r|HO; n is characterized by an Gaussian tail, while we have asymptotically r r r|n ∝ exp − 2 3 r 3/2 . Again, if an observable is not too sensitive to the large r behavior, this discrepancy will not spoil its mean value.
The various observables |ψ n,0 (0)| 2 , r k and p k computed with the AFM states can be obtained using formulae (E . 4), (E . 6) and (E . 7) for l = 0 with the parameter λ given by (3 . 30). Results are summed up in Table V. Again, a direct comparison between the structure of exact and AFM observables can be obtained by using β n (see (G . 4)) instead of |α n |. Let us look in detail only at the mean value r 2 . For the exact and AFM solutions, we have respectively .

(3 . 34)
In contrast with the previous case, all observables are very well reproduced. This is also the case for n ≫ 1, while the overlap | n|HO; n | 2 tends towards very small values in this limit. This is due to the fact that some observables are not very sensitive to the details of the wavefunctions and that the sizes of exact and AFM states stay similar. The mean value HO; n|p p p 2 + r|HO; n being analytically computable for an arbitrary value of the scale factor λ , the relations between the bounds can also be checked. We can verify that E(ν 0 ) = 3 2 2/3 (2n + 3/2) 2/3 ≥ E = −α n and E(ν 0 ) = 3 2 2/3 (2n + 3/2) 2/3 ≥ E * (ν 0 ) = (n + 3/4) 2/3 + 4 πn! (8n + 6) 1/6 Γ (n + 3/2). Moreover, for the ground state, we have As in the previous case, the behavior of observables computed with the AFM is similar for values of l = 0, 1, 2. We do not expect strong deviations for larger values of l. Some typical results are presented in Table VI. Agreement between AFM and exact results are very good, much better than for the previous case. This is expected, since eigenstates for a quadratic potential are closer to eigenstates for a linear potential than eigenstates for a Coulomb potential. The quantum number dependence of the scaling parameter λ corrects, much better than η, the difference between the shapes of AFM and exact eigenstates. Table VI. Ratios between the AFM results (energies and r ) with P(r) = r 2 and the exact results, for several quantum number sets (n, l). l n = 0 n = 1 n = 2 n = 3 n = 4 n = 5 ε HO n,l /E n,l 0

General considerations
Within the auxiliary field method, if the problem studied is analytically manageable, an eigenstate can then be determined with the same computational effort for any set of quantum numbers (n, l). So this method is very useful if one is interested in obtaining analytical information about the whole spectra, wavefunctions, and observables of a Hamiltonian without necessarily searching a very high accuracy. The selection of the potential P(r) seems crucial to obtain good results for the eigenstates. For a linear potential, P(r) = r 2 is clearly the best choice. In Ref. 53), it is shown that P(r) = −1/r provides much better results for an exponential potential, while the choice of P(r) is not so crucial for the logarithmic interaction.
An eigenvalue equation can also be solved within the variational method by expanding trial states in terms of special basis states. The correct asymptotic tail can be well reproduced if the basis states are well chosen. With this method, a matrix representation of the Hamiltonian is obtained and the solutions are computed by diagonalizing this matrix: M upper bounds on the energies are determined with the corresponding M states, where M is the order of the matrix. A very good accuracy is possible if M is large enough. Even for M = 1, the accuracy can be better than the one provided by the AFM (see Sect. 6.4.3). However, if one is interested in closed-form results, the matrix elements must have an analytical expression and the number of computed states M must be limited to 4. But even for M = 2, the eigenvalues can have very complicated expression, not usable in practice. So, the variational method can only provide, at best, a very limited number of eigenvalues and eigenstates with an analytical form.
The WKB method is also a popular method to solve eigenvalue equations. 7), 10), 11) In principle, it is only valid for high values of the radial quantum number n, but it can sometimes yield very good results for low-lying states. 56) Wavefunctions are not necessary to compute some observables 56) but they can be determined for arbitrary value of n with this method. An advantage is that their asymptotic behavior can be correct but, unfortunately, the WKB method is mainly manageable for S-states. Indeed, for l = 0, the interaction V (r) must be supplemented by the centrifugal potential, which complicates greatly the integrals to compute. Moreover, these wavefunctions are piecewise-defined whose different parts must be connected properly at the turning points. So they are not very practical to use. We have a priori no idea of the dependence of the exact result ε(λ ; n, l) on the parameter λ and on the quantum numbers (n, l). The AFM gives approximate answers: (3 . 10) or (3 . 19). Moreover, a simple glance at both formulae convinces us that they have the same dependence in λ , but they differ in the dependence on (n, l) only through a different expression for the principal quantum number Q which is a remnant of the function P(r) used in the AFM.
Owing to the fact that AFM gives the exact result for λ = −1 and λ = 2, a continuity argument suggests that the exact result should not be too different from the previous formulae, but with a modified expression for the principal quantum number. Thus we propose for the eigenvalues of the Schrödinger equation with power-law potential (3 . 3) the following prescription The choice of the function Q(λ ; n, l) is just a matter of guess but we require that it is a continuous function of λ and that it coincides with Q C (n, l) for λ = −1 and with Q HO (n, l) for λ = 2. Besides, there is some freedom. Asymptotically for large l, both Q C and Q HO are proportional to l with a slope equal to unity. Owing to the fact that Q C and Q HO lead to lower and upper bounds, a unity slope is maintained asymptotically whatever the value of λ . Asymptotically for large n, both Q C and Q HO are proportional to n and we maintain also this characteristic behavior for other values of λ . In this work we propose two forms of Q(λ ; n, l) which improve largely the results as compared to the Q C and Q HO values. Although in a context slightly different, R.L. Hall has considered an improvement of the results based on a function P n,l (λ ) which was shown as monotone increasing. His prescription to improve the results (see Ref. 57)) was different from ours, but the spirit is the same. The first prescription depends on two free functions b, g and is chosen as The preceding conditions impose b(2) = 2, g(2) = 3/2 and b(−1) = 1 = g(−1). The simple form Q s is linear (it is called sometimes the "linear approximation") both in n and l and, thus, looks as simple as Q C and Q HO but allows much better results, as will be shown below.
The reader may use it if he needs a simple formulation giving an average relative precision of order of 10 −2 .
The second prescription depends on five free functions a, c, d, e, f and reads The boundary conditions are in this case (a, c, d, e, f )(λ = 2) = (2, 3, 7/2, 5/2, 3/2) and (a, c, d, e, f )(λ = −1) = (1, 2, 2, 2, 1). The form Q q is a rational fraction with a quadratic numerator (it is called sometimes the "quadratic approximation"). This form is more sophisticated that Q s but it gives results with a very high accuracy of order of 10 −4 . Some examples are given in the next section. Other forms have been tested but formulae (3 . 36) and (3 . 37) appeared very convenient. Contrary to the case Q = Q HO , for which ε AFM is a upper bound and the case Q = Q C for which it is a lower bound, we have no certainty concerning the position of the AFM energy as compared to the exact energy in the cases Q = Q s or Q = Q q ; but in any case, we will show that the approximate value is always very close to the exact value. Except the very special case of the linear potential, which is studied in detail above, no analytical expression for the exact energy ε(λ ; n, l) exists (even in S-state) for potentials with λ = −1 or 2. Therefore to compare our AFM results to the exact ones, one needs very accurate eigenvalues ε num (λ ; n, l) of Hamiltonian (3 . 3). For that purpose, we rely on a very efficient method, called the Lagrange mesh method. 54) Without difficulty, one can get an accuracy of 10 −8 much higher than what could be expected from AFM. Thus, in the following, we identify fully the exact value ε and the numerical value ε num .
The ultimate aim is to determine the approximate AFM energies ε AFM that stick as close as possible to the exact values. Two kinds of approximations, denoted generically (a), are considered; in both of them ε (a) AFM is given by (3 . 35) but in the linear approximation ((a) = s) the principal quantum number is given by (3 . 36) while in the quadratic approximation ((a) = q) this number results from (3 . 37). The free parameters of each approximation are denoted σ (λ ). There are two of them σ = (b, g) for the linear approximation and five of them σ = (a, c, d, e, f ) for the quadratic approximation.

Comparisons to numerical results
We now search for an analytical formula giving σ (λ ). In order to do that, we proceed in the following way. In a first step, we choose a sample {λ 1 = −1, λ 2 , . . . , λ p = 2} of p values of λ in the considered domain −1 ≤ λ ≤ 2. We found that p = 10 is a good compromise between the researched quality and the numerical effort. For each value λ i belonging to the sample, we calculate, using the Lagrange mesh method, 54) an array of "exact" eigenvalues ε(λ i ; n, l) for n max + 1 values of the radial quantum number n (n = 0, . . . , n max ) and for l max + 1 values of the orbital quantum number l (l = 0, . . . , l max ). For the precise study that we want, it is enough to choose n max = 8 and l max = 5. The AFM results depend on the free parameters σ (λ i ) = σ i . The game is to find the values of σ i which make the AFM energies the closest to the exact ones.
A tool to do that consists in the minimization of the chi-square quantity (3 . 38) Other choices are possible but we find this one very convenient. In some very specific situations, this function may lead to bad results and we find convenient to introduce an alternative minimization based on the difference quantity (3 . 39) The minimization of one of these functions provides the values σ i which will be denoted σ num (λ i ) because they rely on a numerical prescription. Generally the chi-square and the difference procedures give very similar results, but sometimes a glance at the results may lead us to favor one against the other. Of course for λ = −1 and λ = 2, we know that the values σ num (λ i ) given by the values discussed earlier lead to the exact results so that the minimization procedure is useless and both χ(σ i ) and ∆ (σ i ) vanish. Over the full range of λ sample (10 values), we found that the maximal value of the chi-square is 4.4 10 −4 and 2.2 10 −7 for the linear and quadratic approximations respectively, while the accuracy given by ∆ is 1.33% and 0.07% respectively. Thus, one sees that the very simple prescription (a) = s gives already a relative accuracy of the order of 10 −2 over the whole spectrum (54 values), while the more sophisticated prescription (a) = q gives practically the exact results. Several formulae exist in literature for approximate values of ε(λ ; n, l) resulting from various techniques. Among them, one of the most precise is reported in Ref. 51) with the use of a WKB approach. It must be stressed that our method, even in its crudest version, looks much simpler and gives much better results. Thus, we are confident that our formulae are really an improvement compared to the previous existing ones. The second step of our study is now the guess of continuous functions σ (λ ) which stick as much as possible to the values σ num (λ i ) for the particular values λ = λ i . The choice is huge because we have freedom on the form of the function, and then on the parameters entering it. Since we research simplicity above all, we find that an hyperbola form for all the parameters is very well suited. The fit of the 3 parameters entering the definition of an hyperbola is done using again a chi-square function χ ′ based on the numerical value σ i and the proposed value σ (λ i ).
We summarize below the definite choice for the various parameters. For the linear approximation (3 . 36), one gets For the quadratic approximation (3 . 37), the corresponding functions look like One verifies that for λ = −1, 2 the exact values are recovered. We choose to use integer numbers below 1000 in order to match at best the real numbers yielded by the minimization procedure.
As a matter of test, we present in detail the special case of the linear potential (λ = 1) studied in Sect. 3.3. ε(1; n, 0) is then just proportional to the opposite of the (n + 1)th zero of the Airy function. An asymptotic solution of Hamiltonian (3 . 3) valid for large n is given by (3 . 24). Equating (3 . 24) and (3 . 36), we find b(1) = π/ √ 3 ≈ 1.814 and g(1) = √ 3π/4 ≈ 1.360. These values are the exact ones to describe the spectrum for l = 0 and large n, but not necessarily for l = 0 or low n values. However the values coming from our fit (3 . 40) are b(1) = 127/71 ≈ 1.789 and g(1) = 22/16 = 1.375, very close to these new values. Assuming that this procedure can be extended for any values of l, the hyperbola can be determined completely by using the exact values b(−1) = 1, g(−1) = 1, b(2) = 2, and g(2) = 3/2. The coefficients are then given by These relations could look quite complicated but it is worth noting the symmetries existing between the absolute values of the numbers present in the formula. When a number appears twice in a coefficient of the numerator (denominator), it is the arithmetic mean of two other corresponding numbers in the denominator (numerator). A number which appears twice in g(λ ) is the double than the corresponding number in b(λ ). In Fig. 2, we give an idea of how good is the fit of σ (λ ) with hyperbolae. We show only the plots for the functions b(λ ) and g(λ ) for the linear approximation. The results look similar for the functions σ (λ ) of the quadratic approximation. The solid lines represent the numerical values σ num (λ i ) obtained by a chi-square procedure; in a sense they are the best that we can find. The short-dashed line represent an hyperbola determined analytically; the corresponding functions are given by (3 . 42). The dashed-dotted line represents the hyperbola which sticks as close as possible to the numerical curve; the corresponding functions are given by (3 . 40). This last prescription appears to be very good.
In Table VII, we present the lower part of the spectrum (n max = l max = 3) obtained for different types of approximations and compare the various results to the exact ones. It Table VII. Values of ε(λ = 1; n, l) for various approximations. For each set (n, l), the first line (bold case) is the exact result obtained by numerical integration, the other lines are obtained with the AFM approximation (3 . 35). For the principal quantum number Q(1; n, l), the second line uses the prescription (3 . 42), the third line the prescription (3 . 36) with (3 . 40) (linear approximation), and the fourth line the prescription (3 . 37) with (3 . 41) (quadratic approximation). is clear that all the considered approximations give very good results. As expected, the approximation with (3 . 42) is very good for l = 0 and large n values, but it becomes worse and worse (although not so bad) as l increases. The linear approximation (3 . 36) is very good everywhere and should be considered as a good compromise between simplicity and accuracy while the quadratic approximation is really excellent everywhere and must be used if we want analytical accurate expression of the exact energies. Other approximate forms for the eigenvalues of Hamiltonian (3 . 3) exist and give good results, 49), 58) but they are not considered here because they do not give the exact result for the cases λ = −1 and λ = 2.
Obviously the functions σ (λ ) depend on the points used for the fit but also on the particular choice of the functions χ and ∆ . Other definitions-relative error instead of absolute error or different summations on quantum numbers-would have given other numbers. Nevertheless, the quality of the results tends to prove that they are very close to the best possible ones. §4. Schrödinger equation with an arbitrary potential

Scaling laws
Let E(m, G, a) be the eigenvalues of a Schrödinger equation corresponding to a system of reduced mass m subject to a potential of intensity G and characteristic inverse length a.
The scaling law gives the relationship between E(m, G, a) and E(m ′ , G ′ , a ′ ). Let us start from the corresponding Schrödinger equations The important point is that it is the same function V (x) which appears in both equations. In (4 . 2), let us make the change of variables r r r = αx x x and multiply it by χ. Now, we choose the arbitrary parameters α and χ in order to fulfill the conditions χ/(m ′ α 2 ) = 1/m and αa ′ = a.
In other words, we impose the following values With these values, (4 . 2) can be recast into the form 2 and a similar relation for the energies. The scaling law is thus expressed in its most general form as In fact, it is always possible to define the function V (x) so that a ′ = 1 (for example in defining a new position operator as a ′ r). In what follows, and without loss of generality, we will apply the scaling law for energies under the form This equality is very powerful since it is valid for the eigenvalues of a nonrelativistic Schrödinger equation with an arbitrary central potential. It allows to express the energy in terms of a dimensionless quantity and some dimensioned factors, as we will see below.
Moreover it is possible to give to the mass m ′ any particular convenient value (for example m ′ = 1) so that E(m, G, a) is expressed through the energies of a very simple reduced equation depending on a single parameter G ′ . Instead, one can choose to impose a given value of the intensity G ′ (for example G ′ = 1) so that E(m, G, a) is expressed through the energies of a very simple reduced equation depending on a single parameter m ′ . To summarize the previous discussion, the scaling laws are very powerful since they allow to express the most general form of the energy E(m, G, a) in terms of the eigenvalues of a reduced equation in which two parameters among the three available can be set at a given fixed value. The chosen parameters as well as their values are determined in the most convenient way.

General formulae for energies and mean radius
In Sect. 3, we remarked that, using either P(r) = r 2 or P(r) = −1/r, we obtain, for power-law potentials, AFM energies which have exactly the same form but differ only by the expression of the principal quantum number Q. In other words, this means that ε AFM (λ ; n, l) appears always in the form F(λ ; Q(n, l)) with a universal form of the function F. The only difference is that, in F, we must use Q HO (n, l) for P(r) = r 2 and Q C (n, l) for P(r) = −1/r. This is obvious from the comparison of (3 . 10) and (3 . 19) (or even from (3 . 12) and (3 . 21)).
This property allowed us to modify the expression of this number Q in order to improve drastically the quality of the results. This universal property seems also valid for the expression of the mean radius, as you can see from (3 . 9) and (3 . 18) (or even for (3 . 11) and (3 . 20)).
We will show in this part that these invariance properties are even more general and are valid not only for a power-law potential but whatever the potential V (r) to be used. 40) This property is even stronger since it persists for any power-law type for the starting function P(r) and not only for a quadratic or Coulomb form. Let us demonstrate this fundamental property.
For beginning, let us recall that the AFM energies for the Hamiltonian (3 . 1) are given by (see (3 . 2) and (3 . 35)) (4 . 7) If we ask that E AFM (m, a, λ ; n, l) = E(m, a, λ ; n, l), we can consider that this formula gives the definition of the principal quantum number Q(λ ; n, l) for all power-law potentials. The exact form is known only in a limited number of cases. Obviously, Q(−1; n, l) = n + l + 1 and Q(2; n, l) = 2n + l + 3/2. Using the notations of Appendix G, we can write In the other cases, we showed that (4 . 7) was able to give the results with a relative accuracy better than 10 −4 if a good choice is made for the value of Q. Following the philosophy of Sect. 2.5, we consider that expression (4 . 7) represents the exact expression of our problem and pursue the general AFM method with the starting function P(r) = P (λ ) (r) = sgn(λ )r λ . In this case, K(r) = V ′ (r)r 1−λ /|λ | and the extremization condition (2 . 24) leads to the transcendental equation It seems that r 0 is really λ -dependent, but replacing K(r) by its expression in terms of the potential, this equation turns out to be For a given potential, the inverse function D(x), is obviously universal and one has r 0 = D(Q 2 /m). The value of the AFM energy follows from (2 . 25) and reads Thus the function is also universal and one has Since it is always possible to choose the same mass for the genuine Schrödinger equation with V (r) and the Schrödinger equation with the AFM potential P (λ ) (r), the previous demonstration shows that both the mean radius and the eigenenergy depends on λ through the principal quantum number Q(λ ; n, l) only. We end up with the very important conclusion that can be stated as a theorem (expressed below for the energy but also valid for the mean radius) If, in the expression E(Q(λ ; n, l)) of the approximate energies resulting from the AFM with P (λ ) (r), one makes the substitution Q(λ ; n, l) → Q(η; n, l) (so that E(Q(λ ; n, l)) → E(Q(η; n, l)) with the same functional form for E), one obtains the approximate eigenenergies resulting from the AFM with P (η) (r). In a sense, as long as we use a power-law potential P (λ ) (r) as starting potential, there is a universality of the approximate AFM expression of the eigenvalue, depending only on the potential V (r). The only remainder of the particular chosen potential P (λ ) (r) is the expression of Q(λ ; n, l), as given by (3 . 36) for instance. This result holds whatever the form chosen for the potential V (r), even if we are unable to obtain analytical expressions for the AFM approximation.
This property was first remarked for the particular case of a power-law potential V (r) = r λ switching from the harmonic oscillator (λ = 2) to the Coulomb potential (λ = −1). We proved here that it is in fact totally general. It is probably related to a well-known property in classical mechanics: one can pass from the motion of a harmonic oscillator to the Kepler motion by a canonical transformation. This universality property will be used extensively in the following applications; it allows us to choose the form P (λ ) (r) which is the most convenient for our particular situation. Nevertheless, the existence of possible lower or upper bounds can only be guaranteed for λ = −1 and 2, since formula (4 . 7) is exact only in these cases.
By rewriting (4 . 10) and recalling that the kinetic energy is denoted by 2m , the AFM formula for the energies can be written into the form where Q stands for Q(λ ; n, l) to lighten the notations and to stress that we have some freedom in the choice of this quantity. Equation (4 . 16) defines a "mean impulsion" p 0 from the mean radius r 0 resulting from (4 . 10). Looking at equation (4 . 17), one can recognize the general form of the virial theorem. 59) Finally, (4 . 15) gives the energy as a sum of the kinetic energy evaluated at the mean impulsion p 0 and the potential energy evaluated at the mean radius r 0 . This makes clearly appear the physical content of the AFM approximation. Let us point out that p 0 is not a linear function of Q, as suggested by (4 . 16), since r 0 has a complicated dependence on Q through (4 . 17). The approximate eigenstate |ν 0 is a solution of h(ν 0 ) = p p p 2 /(2m) + sgn(λ )ν 0 r λ and has a characteristic size depending on ). This quantity can be easily computed once (4 . 17) is solved. The virial theorem applied to HamiltonianH(ν 0 ), given by with ν 0 = K(r 0 ) given by (4 . 9). The use of definition (4 . 16) with property (2 . 18) gives The last equality is another way to write the fundamental property (2 . 18). These equations, as well as the boundary character of the solution, are only applicable when the exact form Q(λ ; n, l) is used. In practice, this occurs when λ = −1 or 2, or λ = 1 for S-states only. Nevertheless, a better accuracy can be obtained by an appropriate choice of Q as shown in Sect. 3.4.
In the rest of this section, several interactions will be studied with the AFM. Approximations for all states in the spectrum will be given. When the potential allows only a finite number of bound states, formulae for corresponding critical constants will be also presented.

Square root potentials
The square root potential that we study in this section is fundamental for the understanding of hybrid mesons, which are exotic mesons actively researched nowadays. Indeed, there are two possible descriptions of hybrid mesons within constituent models: • First, a genuine three-body object made of a quark, an antiquark and a constituent gluon. • Second, a two-body object made of a quark and an antiquark in the potential due to the gluon field in an excited state. It has been shown that these two pictures of the same object are, to a large extent, equivalent. 60), 61), 62) In this part, we focus only on the second aspect of the description.
In general, the string energy, and, therefore, the potential energy between the static quark and antiquark in the excited gluon field is given by 63 where a is the usual string tension while b 2 = 2πaK + C is a term exhibiting the string excitation number K and a constant C. These values depend on the model adopted: a pure string theory 63) or a more phenomenological approach. 62), 64) For the study of heavy hybrid mesons, it is therefore very interesting to calculate the eigenenergies of the Schrödinger equation governed by the Hamiltonian where m is the reduced mass and where the parameters a and b depend on the model adopted, which must be confronted to experimental results ultimately. In this section we give an analytical approximate expression for the eigenvalues of Hamiltonian (4 . 21). In Ref. 65), those results have been exploited to study the mass of a heavy hybrid meson in the picture of an excited color field. With this mass formula, it was possible to predict the general behavior of the mass spectrum as a function of the quantum numbers of the system and search for possible towers of states. These results can then be used as a guide for experimentalists.

Expression of the energies
Using scaling laws (4 . 6), dimensionless variables ε and β can be defined so that the general energy E(m, a, b) is easily calculated Let us apply the AFM in order to find approximate expressions for the eigenvalues of the Hamiltonian (4 . 23). Equation (4 . 10) gives immediately A simple change of variable x = 2r 2 0 /Q 4/3 plus the definition of the parameter allows to transform the transcendental equation in the fourth order reduced equation 4x 4 − 8x − 3Y = 0 whose solution is given in terms of the G − function discussed in Sect. B.2. This allows the determination of ε AFM by (4 . 12).
with Y = Y (β ; n, l) given by (4 . 25). The problem is entirely solved. Equation (4 . 26) is complicated but quite accurate. In order to get a better insight into this formula, it is interesting to calculate the limits when Y ≫ 1 and Y ≪ 1 because they allow an easier comparison with experimental data. This technical point was treated in Ref. 65).

Comparisons to numerical results
As we discussed in Sect. 4.2, it is possible to improve drastically the results given by (4 . 26) if, instead of using the natural forms Q = Q HO or Q = Q C , we switch to a more convenient expression. The only relevant dynamical parameter in this case being β , it is natural to propose an expression depending on this parameter: Q(β ; n, l). For power-law potentials we proposed a "linear" (in terms of n, l) approximation Q s and a "quadratic" approximation Q q . For our special potential, one can adopt the same prescriptions, or choose a new better one, if one wishes. Since our aim is not to obtain the most accurate possible analytic expressions at the price of complicated formulae, but rather to show that the AFM is able to give already very good results even with very simple expressions, we adopt here a linear expression, so that we choose Q(β ; n, l) under the form We allow the β value to vary between β = 0 for which the potential is purely linear and β = ∞ for which it is harmonic. We are again in a situation for which Q = Q HO gives an upper bound while Q = Q C provides a lower bound. The procedure used is similar to the one presented in Sect. 3.4 and is described in detail in Ref. 65). In order to obtain functions which are as simple as possible, continuous in β , and which reproduce at best the exact values, we choose hyperbolic forms for A(β ) and C(β ). Explicitly, we find The integers appearing in A and C are rounded numbers whose magnitude is chosen in order to not exceed too much 100. The A and C functions have been constrained to exhibit the right behavior A → 2 and C → 3/2 for very large values of β . formulae (4 . Our results are exact for β → ∞, and the error is maximal for small values of β . But, over the whole range of β values, the results given by our analytical expression can be considered as excellent. Just to exhibit a quantitative comparison, we report in Table VIII the exact ε(β ; n, l) and approximate ε AFM (β ; n, l) values obtained for β = 1, a value for which the corresponding potential is neither well approximated by a linear one nor a harmonic one. As can be seen, our approximate expressions are better than 1% for any value of n and l quantum numbers. Such a good description is general and valid whatever the parameter β chosen.
The upper bounds obtained with P(r) = r 2 are far better than the lower bounds computed with P(r) = −1/r. This is expected since the potential √ a 2 r 2 + b 2 is closer to a harmonic interaction than to a Coulomb one. Better lower bounds could be obtained with P(r) = r. But, the exact form of Q is not known for this potential, except for l = 0 for which Q is given by (4 . 8). With the approximate form Q = (π/ √ 3)n + l + √ 3π/4, 38), 40) we have checked that results obtained are good but the variational character cannot be guaranteed.  In this section, we apply the AFM to find approximate closed analytical formulae for central potentials of exponential form, that is −α e −(β r) η where α, β , η are positive real number. We thus start with the following Hamiltonian Two particular cases are of current use in physics: • The first one is the Gaussian potential (η = 2) which was intensively used in molecular and nuclear physics because it allows very often an exact analytical expression of various matrix elements. • The second one is the pure exponential (η = 1), used also in nuclear physics and for which exact analytical expression for the energies in S-state are known.

Expression of the energies
Using the scaling laws, the eigenenergies E(m, α, β , η; n, l) can be written where ε(η, g; n, l) is en eigenvalue of the dimensionless Hamiltonian An interesting feature of exponential-type potentials is that they all admit a finite number of bound states that depends on the dimensionless parameter g, ruling the potential depth, and defined by (4 . 32). There exists thus "critical constants": Potential depths beyond which new bound states appear. We refer the reader to Ref. 66) for detailed explanations about how to compute critical constants in a given potential. The definition and the properties of the Lambert function, that will frequently appear in our calculations, are given in Appendix C. Some of these potentials have been studied with the AFM in Ref. 67).
The application of (4 . 10) gives Introducing the parameter the solution of (4 . 33) is obtainable using the properties of the Lambert function Lastly the AFM energy is given by (4 . 12). A straightforward calculation gives the final result with Y η given by (4 . 34). The energy ε AFM must be negative, which implies that W (−Y η ) ≥ −2/(η + 2). This is only possible for the branch W 0 with a negative argument, which is the case since −Y η < 0. The quality of this formula for the Gaussian potential (η = 2) is discussed in Ref. 68) and for the pure exponential (η = 1) in Refs. 53), 67).
Let us just mention that eigenenergies of Hamiltonian (4 . 31) for the pure exponential (η = 1) can be analytically computed for l = 0 only. But even in this case, the expression of ε(g) is not very tractable since it is formally defined by the relation 9) where J ρ is a Bessel function of the first kind. Consequently, the AFM result is interesting since it yields an analytical formula for the energy levels of the pure exponential potential that is of simpler use than (4 . 37) for l = 0 and, above all, that remains valid for arbitrary n and l quantum numbers.

Critical constants
By definition, the critical constants g η;nl of the potential are such that ε AFM (η, g = g η;nl ; n, l) = 0, with ε given by (4 . 36). Of course this is the AFM approximation of these critical constants, but according to this last equation, one finds that the energy vanishes for Y η = 2/(η + 2) e −2/(η+2) , a value which is lower than 1/e as required by (4 . 35). Equivalently (equating this particular value of Y η with the general value (4 . 34)), one can say that the critical constants are given by They are such that, if g > g n 0 l 0 , the potential admits a bound state with the given quantum numbers n 0 , l 0 . It is remarkable that our approximation scheme, based on the potential P(x) = x λ for which an infinite number of bound states is present, is able to predict that only a finite number of bound states will be present in exponential potentials. This is another test of the ability of the AFM to analytically reproduce the qualitative features of a given eigenvalue problem. Let us note that the AFM does not give any information about the optimal power-law potential to determine the form of the number Q(n, l). For that, we can rely on direct comparison with numerical solution. It is also very remarkable that whatever the power η, the critical constants are always proportional to the square of the principal quantum number Q. The critical constants for the Gaussian potential are discussed in Ref. 68) and for the pure exponential in Ref. 67). Let us just discuss a little bit this last case. Equation (4 . 37) leads to a determination of the exact critical constants g n0 (η = 1 is here always assumed and is no longer indicated); let us denote them as g * n0 . Indeed these critical constants are such that they correspond to a energy level ε(g * n0 ) = 0. We are thus led to the equation where j n is the (n + 1) th zero of the Bessel function J 0 . At large n, these are given by j n ≈ π(n + 3/4), 69) leading to the asymptotic expression This result is qualitatively similar to (4 . 38), stating that g n0 ∝ Q(n, 0) 2 .
We can try to use formula (4 . 40) to improve the result (4 . 38). Assuming that Q nl = b n + l + c with no constraint on parameters b and c, it is easy to see that formulae (4 . 38) and (4 . 40) coincide for g n0 with b = π/e and c = 3π/(4e). We can then try a new relation to compute the critical constants This formula is in good agreement with exact results: For n ∈ [0, 5] and l ∈ [0, 5], the minimal, maximal and mean relative errors are respectively 0.03%, 8.7% and 3.5%. 67)

Yukawa potential
Among all the central interactions, the Yukawa potential is widely used in atomic physics (effective interaction), nuclear physics (long range behavior of the nucleon-nucleon interaction, due to one pion exchange), hadronic physics (screened Coulomb force). It is thus interesting to apply the AFM to this very important case. We start with the following Hamiltonian 4.5.1. Expression of the energies Using the scaling laws, the eigenenergies E(m, α, β ; n, l) can be written where ε(g; n, l) is en eigenvalue of the dimensionless Hamiltonian This potentials admits also a finite number of bound states that depends on the dimensionless parameter g, ruling the potential depth, and defined by (4 . 45). The application of the AFM to this case leads to the transcendental equation (see (4 . 10)) As far as we know, there does not exist an analytical expression giving r 0 in terms of T . Let us introduce formally the inverse function Ω = t −1 of the function t(u) = e u u(u − 1). By definition, we have the properties Thus one has r 0 = −Ω (T ). A calculation giving the AFM energy through (4 . 12) leads to the expression ε AFM (g; n, l) = g 2 . (4 . 48) Ω (u) is a multi-valued function composed of 3 monotonic branches (see Fig. 3). However, we have the constraints that r 0 > 0 and ε AFM < 0. This imposes the condition −1 ≤ Ω (T ) < 0, with T > 0. This condition automatically selects the branches, passing through the origin, defined in where φ is the golden ratio. So, there is no ambiguity. It is possible to obtain an approximate AFM energy formula in terms of usual simple functions, but it is then necessary to use an approximate form for the solution of the transcendental equation (4 . 46). 67) Since With the choice Q = Q C , this corresponds to the Coulomb energy. It is expected since g → ∞ is equivalent to β → 0. Upper bounds can be obtained by choosing Q = Q C = n + l + 1 in (4 . 48). This seems quite natural since a Yukawa potential is similar in form to a Coulomb one. It is possible to improve the accuracy of the AFM formula by using an appropriate form Q = A(g)n + l + C(g) 67) but the variational character of the approximation cannot then guaranteed. An accurate energy formula has been found in Ref. 70) from a fit of the numerically computed energy levels of the Yukawa potential. In our notations, it reads where A ′ = 1.9875, B ′ = 1.2464 and σ = 0.003951, and where This formula is rather different from ours. However, they coincide at the limit g → ∞ when Q = Q c (see (4 . 49)). Equation (4 . 51) gives the critical constants of the Yukawa potential with a relative accuracy around 0.4%. Table IX. Eigenvalues ε(g; n, l) for a Yukawa potential with g = 30 as a function of (n, l) sets. First line: exact value; Second line: formula (4 . 48) with Q = Q C = n + l + 1; Third line: formula (4 . 50). A * indicates a non real or a non negative value.
The quality of the two approximations presented above can be appraised in Tables IX for a particular value of g. Not all bound states can be found with (4 . 48). Let us note that formula (4 . 50) gives generally the correct number of bound states, but not in all cases. 67) Formula (4 . 48) with Q = Q C gives quite good results for the lowest eigenvalues. The quality of the fit for formula (4 . 50) is better but comparable to the quality of AFM formula (4 . 48) with a more sophisticated parameterization of Q. 67) Nevertheless, formula (4 . 50) is empirical, while AFM result (4 . 48) is obtained from an explicit analytical resolution of the Schrödinger equation.

Critical constants
Let us discuss below more deeply the critical constants relative to the Yukawa potential. It is easy to find from (4 . 48) that the T value which cancels the AFM energy is equal to T = 2/e. Using the definition (4 . 46), one obtains the AFM critical constants for the Yukawa potential g yuk;nl = eQ 2 .
(4 . 53) It is again remarkable that the critical constants for the Yukawa potential are also proportional to the square of the principal quantum number, a property that seems universal for exponential-like potentials. Due to the smallness of γ and δ in (4 . 51), we can replace S l by S 0 in (4 . 51). It clearly ap- asymptotically. This is precisely the quadratic behavior predicted by our analytical results. It is important to stress that formula (4 . 51) was specially designed for fitting the Yukawa critical constants and that it needs several free parameters, while our AFM result (4 . 53) follows a general algorithm valid in any circumstance. The fact that both give the same asymptotic behavior is again a very strong argument in favor of the relevance of the AFM. Introducing a constant which allows to obtain the exact result for n = l = 0 and considering the asymptotic expansion for large values of l, (4 . 51) is approximately given by g G;nl ≈ (1.248 n + 1.652 l + 1.296) 2 .
(4 . 54) This is not too different from our result (4 . 53) with Q = Q C = n + l + 1 649. This formula is unfortunately not very accurate: For n ∈ [0, 4] and l ∈ [0, 4], the minimal, maximal and mean relative errors are respectively 6%, 42% and 26%. The minimal, maximal and mean relative errors are respectively 0.01%, 0.9% and 0.3% with formula (4 . 51). The very good quality of these last results is due to the use of a complicated formula fitted on exact results. Our formula is simpler and the general behavior is predicted by the AFM. Let us note that it is possible to improve the values of AFM critical constants for the Yukawa potential using an appropriate form Q = An + Bl +C, with A, B and C fixed, and keeping the simplicity of the original formula (4 . 53). 67)

Sum of power-law potentials
In Sect. 2.4.1 we suggested that potentials of type aP(r) + V (r) could be treated with the AFM in some occasions. Since a power-law potential is very often a good starting point, it is interesting to see whether an analytical AFM solution exists for a second contribution which is itself a power-law potential. The resulting total potential is thus the sum of two power-law potentials. Those types of potential are frequently used in various domains of quantum mechanics. This is the subject of this section, whose calculations are based on (2 . 37) and (2 . 38).

Expression of the energies
It is very instructive to start with the potential P(r) = sgn(η)r η = P (η) (r). In this case, the extremization condition (2 . 37) is written (see also (4 . 9)) where Q η is a priori the optimal principal quantum number for P (η) (r). The AFM energy (2 . 38) is given by With this in mind, let us introduce the potential V (r) in the form of a power-law potential sgn(λ )br λ so that the total potential for our system is given by The above considerations show that the extremization condition and the AFM energy look like Should have we taken instead P(r) = P (λ ) (r) and V (r) = P (η) (r), we would have recovered the same physical problem with the potential W (r) given by (4 . 58). Due to the form of basic equations (2 . 37) and (2 . 38) the only change would have been Q η → Q λ . Because of the symmetry of the problem there is no reason to prefer the choice Q η , as in the previous equations, instead of the choice Q λ ; moreover in any application, we already saw that it is convenient to change the value of the principal quantum number. In consequence, in the following, we will use the more neutral notation Q. This result can also be directly obtained by using (4 . 15)-(4 . 17), but it is interesting to illustrate that the AFM depends only on the potential P(r) through the principal quantum number Q.
In general, one does not have an analytical expression for the root of the extremization equation (4 . 59). This is possible only for very specific values of the powers η and λ .
In this section, we will study such a problem for the most favorable cases, where η is chosen to give an exact expression for the eigenvalues (in practice η = 2 and η = −1) and where λ is chosen in order to have an analytical result. In the following, it is assumed that a > 0, b > 0 and −2 ≤ λ ≤ 2.

Solvable potentials
Let us examine first the case P(r) = r 2 . The extremization equation is a transcendental equation for which an analytical solution does not exist automatically. The only cases for which we are sure that an analytical solution exists is when it can be transformed into a polynomial of degree less than or equal to 4. In order to investigate this condition, let us put (λ + 2)/4 = p/q (p and q are relatively prime integers), and define the new variable X = r 4/q 0 . The corresponding extremization condition is transformed into All the solvable potentials should verify the conditions 0 ≤ p ≤ 4, 1 ≤ q ≤ 4. An exhaustive research of all the solvable potentials in our domain for λ leads to the following, non trivial, values of the power λ : We will study in detail the cases: • λ = −2 because it corresponds to a centrifugal term with a real parameter instead of the usual l(l + 1) term; • λ = −1 because it corresponds to a simplified potential for hadronic systems with a short range Coulomb potential and a quadratic confinement; • λ = 1 because this anharmonic potential is sometimes used in molecular physics. Now we investigate the case P(r) = −1/r. The same type of argument with λ +2 = p/q and X = r 1/q 0 leads to the following polynomial equation An exhaustive list of the non trivial solvable potentials is given below Among them we will study: • λ = 2 because it corresponds to a potential −a/r + br 2 which is already studied with P(r) = r 2 . This is the only potential that can be described with either P(r) = r 2 or P(r) = −1/r and this property allows very fruitful comparisons. • λ = 1 because it corresponds to the funnel potential (Coulomb + linear) which is widely used in hadron spectroscopy. 32) Finding approximate analytical values for the energies corresponding to this potential is thus a very interesting question. To our knowledge, such formulae have not been proposed in the literature. It could be also interesting to introduce potentials with λ ≤ −2 but with the restriction that it is repulsive at the origin (for instance, van der Waals forces or Lennard-Jones type of potentials). In this case, the transcendental equation is modified a little bit since λ + 2 ≤ 0 implying a negative fraction p/q. An exhaustive research of the different values of λ gives the following results: In this list, we will just consider the case: • λ = −2 because the corresponding potential, known as the Kratzer potential, exhibits its spectrum under an analytical form for all values of radial quantum number n and orbital quantum number l.

Kratzer potential
The Kratzer potential 9) is defined in its simpler form as It can mimic the interaction of two atoms in a diatomic molecule. Moreover, it presents some interest as a benchmark since it is one of the rare potentials for which ones knows an exact analytical expression of the energies valid for any n and l quantum numbers. In fact, one can define the Kratzer potential in a more general formulation (see Ref. 71)) as One recovers the simplest form putting G = a 2 and f = 2/a. It can be shown (see for example the hint given in Ref. 9)) that the exact eigenenergies read Applying the previous AFM expressions to this potential leads to the following mean radius (η = −1, λ = −2, a = G f and b = −G in (4 . 58)) which, inserted in (4 . 60), gives the desired result .
(4 . 70) It is natural to set Q = Q C = n + l + 1, considering that P(r) = −1/r. We thus remark that the approximate value E (K) presents the correct asymptotic behavior for large n and for large l. Just to have an idea of the quality of this approximation, let us calculate the difference δ between the terms in brackets appearing in the denominators of E and E (K) . It is just a matter of simple algebra to find (4 . 71) One always has δ > 0 so that E(n, l) > E (K) (n, l) meaning that the AFM approximation is a lower bound. Moreover, for small intensity and/or mass, mG ≪ 1, or for large angular momentum, l ≫ 1, the approximate value tends to the exact one and we have more explicitly For these limit conditions, the AFM approximation tends towards the exact result. This behavior is easily understandable because, under those conditions, the contribution due to 1/r is predominant as compared to the contribution of 1/r 2 , and both expressions tend towards the same exact Coulomb result.

Quadratic + centrifugal potential
We consider now the potential (for an attractive centrifugal potential, not all values of b are relevant 72) ) The incorporation of the term ±b/r 2 into the l(l + 1)/r 2 term already present in p p p 2 allows to get the exact eigenvalue using the same kind of arguments than those developed in the harmonic oscillator case. 9) Explicitly, we obtain E(n, l) = a 2m 2(2n + 1) + (2l + 1) 2 ± 8mb . (4 . 74) Using the AFM with P(r) = r 2 and setting Y = 2mb/Q 2 , the value r 0 that extremizes this energy comes from a first degree equation and reads Substituting this value into the energy (4 . 60), one obtains a very simple expression for the approximate energy E = 2 ab(1 ±Y )/Y or, alternatively This quantity and the corresponding exact one depend on three parameters m, a, b but we know that the general scaling law properties allow us to write them in a more pleasant form E(m, a, b; n, l) = a 2m ε(β = 2mb; n, l), (4 . 77) where ε(β ; n, l) is an eigenvalue of the reduced Schrödinger equation for which the Hamiltonian depends now on a single dimensionless parameter β h = p p p 2 + r 2 ± β r 2 .
(4 . 78) The exact eigenvalues of this Hamiltonian are given by ε(β ; n, l) = 2(2n + 1) The approximate values immediately come from (4 . 76) It seems natural to choose Q = Q HO , considering that P(r) = r 2 . In this case, one can check that the relative error between ε and ε (qc) decreases as l −3 for a fixed value of n and large l; it decreases as n −1 for a fixed value of l and large n. This behavior is easily understandable because, for large values of the quantum numbers, the contribution due to r 2 is predominant as compared to the contribution of 1/r 2 and both expressions tends towards the same exact harmonic oscillator result. Let us assume that β ≪ 1 and let us choose Q = Q HO ; a Taylor expansion truncated to first order leads to It is easy to check that this expression can also be obtained by perturbation theory. In particular for β = 0, one recovers the exact value 2Q HO , as expected. ). Substituting this value into the energy (4 . 60) and making a little algebra leads to the desired approximate energy with Y given by (4 . 84). It is natural to use the value Q = Q HO , considering that P(r) = r 2 .
As in the previous case, this quantity and the corresponding exact one depend on three parameters m, a, b but the general scaling law allows us to write them in terms of a reduced quantity depending on a single parameter β E (m, a, b; The parameter β could also be associated with the quadratic potential, but this less interesting case is not considered here. Let us assume that β ≪ 1 and let us choose Q = Q HO ; a Taylor expansion truncated to first order leads to In particular for β = 0, one recovers the exact value √ 3Q HO , as it should be. This result comes also from the perturbation theory.
The limit β → ∞ is not physically relevant, but it is interesting to consider it in order to check the formula. In this limit, we find ε (an) (β ; n, l) = 3(4β The dominant term is the result expected for a pure linear potential, as given by (4 . 7).

Quadratic + Coulomb potential
We study now the quadratic + Coulomb potential defined as The equation giving the mean radius looks like  7)). Substituting this value into the energy (4 . 60), one is led, after some manipulations, to the desired approximate energy with Y given by (4 . 93). As usual with the AFM, the quality of the approximation depends on the choice of Q. For instance, the use of Q = Q C and Q = Q HO in the expression for Y gives respectively a lower bound and an upper bound on the true energy. The quantity (4 . 94) and the corresponding exact one depend again on three parameters m, a, b but the general scaling law allows us to write them in terms of a reduced quantity depending on a single dimensionless parameter β . One can imagine two formulations depending on whether β is part of the quadratic contribution or of the Coulomb contribution: ; n, l , ; n, l .
(4 . 96) The ε and η energies are the eigenvalues of the reduced Schrödinger equations for the respective Hamiltonians h ε et h η : The approximate values corresponding to these reduced Hamiltonians follow from (4 . 94): Let us assume that β ≪ 1 or β ′ ≪ 1 (it is equivalent to the formulation V (r) ≪ P(r)): • A Taylor expansion truncated to first order for the formulation based on the ε form leads to In particular for β = 0, one recovers the exact value √ 3Q/4 if Q = Q HO . • The Taylor expansion for the formulation based on the η form gives In particular for β ′ = 0, one recovers the exact value In all these cases, the approximate formulae resulting from AFM agree with the result of perturbation theory. The limits β → ∞ and β ′ → ∞ are interesting to consider in order to check the formulae: • The ε-based formulation gives If Q = Q C , the dominant term is the exact result for a Coulomb potential. • The η-based formulation gives It is worth mentioning that analytical solutions of the Schrödinger equation with quadratic plus Coulomb potentials are presented in Ref. 73), within the framework of a two-dimensional system of two-interacting electrons (−1/r) in a confining magnetic field (r 2 ). Closed-form solutions are found for particular values of magnetic field and spatial confinement length. It is argued that a generalization to a three-dimensional space is possible. But, as the corresponding formulae are not given explicitly, a comparison with our results is not available. More generally, the Schrödinger equation with a −a/r + b r + c r 2 potential (particular cases are studied in Sect. 4.6.6 and 4.6.7) is directly linked with the biconfluent Heun's equation. 74) 4.6.7. Funnel potential In this section, we study the funnel potential defined as This potential is particularly important and its interest is discussed in Sect. 4.6.2. The extremization condition writes in this case with Y given by (4 . 108). Let us mention that another form of this equation can be found thanks to the following relation with the change of variables Y = sinh(3θ ).
The quantity (4 . 109) and the corresponding exact one depend again on three parameters m, a, b but the general scaling law allows us to write them in terms of a reduced quantity depending on a single dimensionless parameter β . In hadronic physics, the dominant interaction between a quark and an antiquark is a confining linear potential. 32) So, although the linear potential has no analytical exact solution for all values of the quantum numbers, it is also interesting to consider two formulations depending on whether β is part of the linear contribution or of the Coulomb contribution: ; n, l .
The approximate values corresponding to these reduced Hamiltonian follow from (4 . 109): Let us assume that β ≪ 1 or β ′ ≪ 1: • The Taylor expansion for the formulation based on the ε form gives In particular for β = 0, one recovers the value expected for a pure linear potential, as given by (4 . 7).
• The Taylor expansion for the formulation based on the η form gives In particular for β ′ = 0, one recovers the exact value −3 5/3 /(4Q 2 ) for the Coulomb potential if Q = Q C . In all these cases, the approximate formulae resulting from AFM agree with the result of perturbation theory.
The limits β → ∞ and β ′ → ∞ are interesting to consider in order to check the formulae: • The ε-based formulation gives The dominant term is the exact result for a Coulomb potential if Q = Q C . • The η-based formulation gives The dominant term is the result expected for a pure linear potential, as given by (4 . 7).

Comparison with numerical results
The potentials studied in this section are rather sophisticated and, presumably, the AFM should lead to less good results than power-law potentials. Thus, it is interesting to test this method in these less favorable cases. A detailed discussion can be found in Ref. 40).
Here we focus our analysis on the funnel potential only, because of its physical importance, and rely on the reduced Hamiltonian h ε (4 . 113). The AFM energies ε (f) (β ; n, l) are given by (4 . 115). For most of the hadronic systems, the physical values of β vary from 0 to about 1.5.
As usual, the simplest prescription to improve harmonic oscillator or Coulomb-like results is to adopt a principal quantum number of the form (4 . 122) The protocol to determine the best values for the b and c coefficients follows essentially the same steps than those explained in Sect. 3.4.2 and is not repeated here. The only difference is the form of the functions b(β ) and c(β ) which are better fitted with Gaussian functions than hyperbolae. Explicitly, the parameterization retained is  Allowing a β -dependence for the coefficients b and c improves greatly the approximate eigenvalues since one gains a factor 10 −2 to 10 −3 on the chi-square values. 40) The quality of the fits can be appraised for β = 0.5 by examining the values of approximate results compared with exact ones presented in Table X. The values Q = Q HO and Q = Q C give respectively an upper and a lower bound on the exact energy. However, as can be seen from Table X for the case Q = Q C , these bounds are not very close to the exact value and, thus not very interesting (the error is of the order of 10% or more). The value Q = Q(β ) gives more accurate results (the error is less than 5% and is often of order of 2-3%), but we have no certainty concerning the (anti)variational character of the approximation. The global accuracy is less good than for power-law potentials but remains quite acceptable. §5. Two-body spinless Salpeter equation In this section, we will apply the AFM to the spinless Salpeter Hamiltonian for which the kinetic energy has a semirelativistic form H = σ p p p 2 + m 2 +V (r).
(5 . 1) The term "relativistic" is generally not used because such a Hamiltonian is not manifestly covariant. The parameter σ is equal to 1 for one-body systems (r is then the distance of the particle from the origin of the force field) and to 2 for two-body systems with equal masses (r is then the distance between the two particles). However, it is very convenient to choose it as an arbitrary positive real parameter. We will prefer the notation M instead of E for the eigenvalue of this kind of Hamiltonian because it does include the rest mass of the particles. Obviously, there is much less material concerning the properties of eigenvalues for semirelativistic equations as compared to the huge bulk of results for nonrelativistic treatments. Nevertheless the quadratic potential has been studied taking the Fourier transform of of the potential. 75) Moreover, working in momentum space, Boukraa and Basdevant were able to obtain interesting properties concerning a class of potentials. 76) For two-body systems, the presence of two different masses cannot be dealt with a reduced mass, so that we need the introduction of two different auxiliary fields to treat the kinetic energy term. This situation makes the problem much more complex. It is briefly studied in Sect. 5.6. Thus, in the following, we focus on Hamiltonian (5 . 1).

Scaling laws
Let us consider the spinless Salpeter equation where G and a are two real numbers defining the potential, a being related to its extension and G to its intensity. The eigenvalues M depend a priori on the four parameters m, G, a, σ . Scaling laws allow to obtain the general form of M(m, G, a, σ ) in terms of the energies of another system depending on less free parameters. Let us examine this point. Defining ρ ρ ρ = ar r r, π π π = p p p/a, one has π π π 2 + (m/a) 2 This is just the spinless Salpeter equation We are thus led to the following scaling law for the mass spectrum It means that both one parameter of the potential (here the extension of the potential a) and the number of constituents in the system σ can always be set equal to 1 in order to simplify the computations without loss of generality. The general energy formula will then be recovered thanks to (5 . 6).

Generalities
To deal with the worrying square root, it is interesting to introduce one auxiliary field µ in order to build the AFM Hamiltonian corresponding to (5 . 1) The semirelativistic operator is then replaced by a nonrelativistic counterpart with a parameter µ to be determined. This technique to get rid of the square root operator is not new and has been used previously to study hadronic physics. 28), 29), 30), 31), 32) Let us remark that, asH(µ) and H have not the same kinetic part, most of the results about lower and upper bounds developed in Sect. 2.4 are not applicable. At this stage, one can imagine two cases which are discussed below.

Potential V (r) is solvable
Let us assume that the analytical solution of the Schrödinger equation with the potential V (r) is known. In particular, this approach is interesting if V (r) = r 2 or V (r) = −1/r. One can then write the eigenvalues ofH(µ) as This can be also demonstrated using a simpler version of the AFM which can be used only when it is interesting to replace an operator by its square. 41) Using the Hellmann-Feynman theorem, 39) it has been shown that the value of µ 0 is given by 32), 41) µ 0 = 2 σ p p p 2 + m 2 , where the mean value is computed with an eigenstate ofH(µ 0 ). An estimation of the difference δ between the exact and approximate eigenenergies is given by 41) whereμ 0 is the auxiliary field obtained with the AFM applied to the Hamiltonian H ur .

Potential V (r) is not solvable
If the eigenvalues of h = p p p 2 /µ +V (r) are not known, it is tempting to treat the kinetic energy term and the potential term on equal footing, in the spirit of what is done in Sect. 2.2, introducing the auxiliary field µ to deal with T and the auxiliary field ν to deal with V . The case of a general function P(r) cannot be solved analytically, but restricting to the important power-law potential P(r) = sgn(λ )r λ , very interesting expressions can be obtained.
The eigenenergies e(µ, ν) of the nonrelativistic Hamiltonian h(µ, ν) = p p p 2 /µ + νP(r) are given by (3 . 2) and (3 . 35) and the AFM mass by (the K and J functions are defined as usual) The extremization of M needs to solve two coupled equations with respect to µ and ν which provide the values µ 0 and ν 0 . Introducing the new variable it can be shown that both µ 0 and ν 0 can be expressed in terms of x 0 only. Namely Reporting these values into the definition (5 . 16) allows the determination of x 0 through the transcendental equation One deduces easily the expression of the AFM mass It is possible to go further in the simplification using the expression of the K function and defining the new mean radius The transcendental equation (5 . 19) and the value of the AFM mass (5 . 20) are given now by simpler expressions: Despite the fact that two auxiliary fields were introduced, only one transcendental equation is necessary to solve the problem. This great simplification is due to the use of a power-law expression for the basic function P(r). It is also clear from (5 . 22) and (5 . 23) that, very much in the same way than the nonrelativistic case, the expression of the mean radius and of the mass depend on the λ variable only through the value of the principal quantum number Q(λ ; n, l), which can thus be modified at will in order to improve the results. So, in the following, this quantum number will be simply denoted Q.
As in the nonrelativistic case, the AFM also yields approximations for the eigenstates. The approximant states, which are eigenstates of p p p 2 /µ 0 + ν 0 sgn(λ )r λ , are characterized by a size parameter depending on the product µ 0 ν 0 (see (E . 3

), (F . 3) and (G . 3)). This quantity is given by
(5 . 24) A lot of observables can then be computed in terms of r 0 only. Another procedure to solve the problem is to start from formula (5 . 8) and replace the exact value e(µ) by a corresponding approximation given by an AFM solution obtained from the Hamiltonian h(µ, ν 0 ) = p p p 2 /µ + ν 0 P(r) with P(r) = sgn(λ )r λ . After this first step, the optimal value ν 0 depends not only on V (r) but also on µ. One can then show that the subsequent extremization on µ gives exactly the results (5 . 22) and (5 . 23) presented above.
This procedure in two steps is nevertheless interesting to obtain informations about the bounds on exact energies. The replacement of the square root operator by a nonrelativistic counterpart without modifying the potential, H →H(µ) (see (5 . 1) and (5 . 7)), yields upper bounds. If the exact potential V (r) is then replaced by an AFM approximantṼ (r) ≥ V (r), the upper bounds will be further increased. On the contrary, if the new potential is such that V (r) ≤ V (r), the upper bounds will be decreased. As a consequence, the energies given by (5 . 22) are upper bounds on the exact energies ifṼ (r) ≥ V (r), or equivalently if the function is convex, the antivariational character of the approximation cannot be guaranteed.
The previous equations take a particularly simple form for massless particles. Indeed, for m = 0, they reduce to Let us define the function D(x) as the inverse function of and the related function F(x) defined by Once the potential V is given, these functions are universal in the sense that they do not depend on the mass m (except if the potential is mass-dependent) and can be computed once and for all. The mass for a system in the ultrarelativistic limit is thus written in a very simple form, namely M AFM = F(σ Q).
which is a kind of "total energy of a free particle". After some algebra, one can check that the virial theorem applied to Hamiltoniañ H(µ 0 , ν 0 ), whose solutions are M(µ 0 , ν 0 ) given by (5 . 15), implies that where |µ 0 , ν 0 is an eigenstate ofH(µ 0 , ν 0 ). These equations, as well as the boundary character of the solution, are only applicable when the exact form Q(λ ; n, l) is used. In practice, this occurs when λ = −1 or 2, or λ = 1 for S-states only. Nevertheless, a better accuracy can be obtained by an appropriate choice of Q as we will see below.

General case
The power-law potential was considered as a prototype for testing the AFM in the framework of the Schrödinger equation. It was shown that an analytical expression for the eigenenergy is obtained whatever the power λ under consideration. We apply in this section the AFM to find approximate analytical energy formulae for the spinless Salpeter Hamiltonian (5 . 1) with the potential V (r) = sgn(λ )a r λ where a is positive. We consider that λ > −2, as needed in the nonrelativistic case to have bound states. In this cases, (5 . 22) and (5 . 23) become With the change of variable a form which will be very convenient to use. With this convention, the mass (5 . 36) is also written in a simple form In order to obtain analytical energy formulae, one should be able to analytically solve (5 . 38). Let us now examine the cases for which this is possible. First, we set λ + 2 = p/q and X 0 = x 1/q 0 , where p and q are co-prime integers. Then, (5 . 38) becomes σ 2 4 A polynomial possesses analytical roots if its order is less or equal to 4; therefore all the solvable potentials should verify the conditions 0 ≤ p ≤ 4, 1 ≤ q ≤ 4. An exhaustive research of all the solvable potentials leads to the following values for the power λ : Among these allowed values, three are of particular interest: λ = −1 and 2 because the Coulomb problem and the harmonic oscillator play a central role in theoretical physics, and λ = 1 since a linearly rising potential is generally considered to be a relevant approximation of the confining potential in QCD. 47) These three cases are explicitly solved in the following sections. Thus, in contrast to the nonrelativistic case where an analytical expression exists for an arbitrary value of λ , the spinless Salpeter equation leads to a restricted list of favorable values (5 . 41). Fortunately, the most interesting potentials belong to this list. Let us briefly discuss the existence of physical solutions, as function of the parameters. Starting from (5 . 38) and (5 . 39), it is possible to prove the following properties: • For λ > 0, there exists always a physical solution with M(λ ) > 0.
• For −1 ≤ λ < 0 there exists a physical solution with M(λ ) > 0 as long as the parameter a is less than a critical value a c (λ ) given by Let us remark that a c (λ ) depends on m, except a c (−1) = σ Q. • For −2 ≤ λ < −1 either there is no solution for (5 . 38), or there exists two solutions but, in the latter case, these solutions are not compatible with the nonrelativistic expressions when m tends toward infinity. In both cases, the corresponding solutions are not physical so that these types of potential must be discarded from our study. In the following, we will denote Q λ = Q(λ ; n, l) the best possible form of the principal quantum number for the power-law potential V (r) = sgn(λ )a r λ in a Schrödinger equation. We have then Q 2 = Q HO , Q −1 = Q C and Q 1 given by (4 . 8). For other values of λ , a formula of type (3 . 36) can be chosen for instance.
With Q = Q 2 these formulae yield an upper bound on the energy, since this setting gives the exact solution for the corresponding nonrelativistic Hamiltonian. In the nonrelativistic limit (Y → ∞), the last equations reduce to as expected for a nonrelativistic harmonic oscillator. For large value of m, the choice Q = Q 2 is clearly optimal. For small values of m, another choice could give better results.

Linear potential
The resolution of the case λ = 1 is rather similar to the one of the harmonic oscillator. Defining Following (B . 2), the solution is given by X (Y 1 ) = F − (Y 1 ). After a rearrangement of (5 . 39), the mass spectrum reads (under one of the equivalent forms) The energy spectrum (5 . 39) is finally given by It is obvious from this last equation that bound states exist only when a < σ Q, in agreement with (5 . 42). The most stringent upper bound for a is actually found for the ground state, that is a < σ Q| n=l=0 . With Q = Q −1 , (5 . 56) yields an upper bound on the energy since this setting gives the exact solution for the corresponding nonrelativistic Hamiltonian. This is confirmed by the results of Refs. 77), 78) in which formula (5 . 56) was already given.
As the coupling constant a is dimensionless in this case, the only mass scale of the Hamiltonian is m. So the nonrelativistic limit cannot be obtained by setting m → ∞, as usual. It is well known that the mean speed of the particle is independent of m and proportional to a in the Coulomb case. So, the nonrelativistic limit is achieved for a → 0. When a ≪ 1, (5 . 56) reduces to

Ultrarelativistic limit
In the ultrarelativistic limit, that is to say m = 0, the various elimination equations for the auxiliary field and the energy formulae become simpler. Equation (5 . 38) is then written as and the final energy spectrum is given by (using the non trivial solution of (5 . 62)) Let us stress that, in this special case, an analytical expression is obtained whatever the power λ of the potential, in contrast to the particular values (5 . 41) resulting from the general case. One can check that (5 . 45) and (5 . 51) reduce to (5 . 63) when m → 0. The quantity M ur λ is physically a mass. It appears that M ur λ < 0 when −1 < λ < 0 and M ur −1 = 0. No bound state of massless constituent particles can be found in these cases. Although (5 . 63) is positive for −2 < λ < −1, this range of values for λ has been proved to be unphysical.
For the particular case λ = 1, this last formula becomes M ur 1 = 2 √ σ a Q, and one has consequently (M ur 1 ) 2 ∝ Q. Such a linear behavior between the squared masses and the quantum numbers of the different states (with Q = Q 1 for instance) is a well-known property of the spinless Salpeter Hamiltonian with a linear potential and massless constituents. Such a Hamiltonian is one of the simplest ways to describe a light meson in a potential approach of QCD, and it is experimentally checked that the square masses of the light mesons mainly grow linearly with their angular momentum (Regge trajectories). See for example Refs. 47), 82) for a discussion of that point.
It is worth looking at the ultrarelativistic harmonic oscillator. The energy spectrum of the Hamiltonian h ho = 2 p p p 2 + a r 2 (5 . 64) can be analytically computed for l = 0 and reads 83) With (5 . 64) as a starting point, a natural choice for Q is Q 2 , but it is clear that both M ho and M ur 2 are identical provided that Q = Q(1; n, 0), the principal quantum number giving the exact energy formula in the case of a nonrelativistic kinetic energy with linear potential. The explanation of this fact is that the Fourier transform of Hamiltonian (5 . 64) is a nonrelativistic Hamiltonian with linear potential. This is a supplementary indication that Q should be different in the semirelativistic and nonrelativistic cases for the same potential.

Improved formulae for ultrarelativistic power-law potentials
Exactly as it was the case in the nonrelativistic framework, it is possible to improve drastically the energy formulae of the spinless Salpeter simply by changing the value of the principal quantum number Q. To study this problem, we employ the scaling law properties and consider the dimensionless Hamiltonian valid for massless particles H = 2 p p p 2 + r λ (5 . 68) with λ > 0. Very accurate numerical values for the eigensolutions of this Hamiltonian can be obtained with the Lagrange mesh method. 54) The approximate energy spectrum is given by (5 . 63) with a = 1 and σ = 2. With the choice Q = Q 2 = 2 n + l + 3/2, upper bounds are obtained. As mentioned before, another choice for the n-and l-dependences of Q can greatly improve the results. By using the form  Table XI.  The approximate energy spectrum is given by (5 . 56) with m = 1 and σ = 2. With the choice Q = Q −1 = n + l + 1, upper bounds are obtained and the results are exact in the limit a → 0. As mentioned before, another choice for the n-and l-dependences of Q can greatly improve the results. By using the form  Table XII.
One can see that Q = Q −1 = n + l + 1 is a better choice for large values of n or l. This can be understood as a kind of nonrelativistic behavior since the limits a → 0 and Q → ∞ have similar effects.

Square root potential
In Sect. 4.3, we introduced the square root potential V (r) = √ a 2 r 2 + b 2 and commented its important role in the theoretical description of hybrid mesons in a nonrelativistic framework. Nevertheless, such a type of potential allied with a relativistic kinetic energy operator is important to describe hybrid mesons with low mass quarks as well. The general formula with a non vanishing quark mass looks rather difficult and we prefer here to present the two extreme limits and show the close connection with the previous section through a Fourier transform.

Nonrelativistic kinematics
The nonrelativistic limit follows from the Hamiltonian H = p p p 2 2m + √ a 2 r 2 + b 2 and as been studied in Sect. 4.3. It is easy to show that the Fourier transform of this Hamiltonian, denoted as H FT , reads This last Hamiltonian is nothing else than a spinless Salpeter one (5 . 1) with a harmonic potential. Provided that the proper substitution rules are taken into account, the energy spectrum in this case is the same than the one computed in Sect. 5.3.2. We have thus This potential was studied in more detail in Sect. 4.3, and we recovered here the results presented there (see (4 . 26)). This remarkable result proves that AFM is an approximation consistent with the Fourier transform.

Ultrarelativistic limit
In the ultrarelativistic limit (m = 0), the corresponding Salpeter Hamiltonian is H = σ p p p 2 + √ a 2 r 2 + b 2 . To get the eigenenergies, it is not necessary to start a new calculation. In fact, the Fourier transform of H is simply H FT = σ p p p 2 + M 2 + a r, (5 . 80) with M = b/σ . But this last Hamiltonian appears to be a spinless Salpeter one (5 . 1) with a linear potential. Provided that the proper substitution rules are taken into account, the energy spectrum in this case is the same than the one computed in Sect. 5.3.3. We have thus When b → 0, these last equations logically reduce to (5 . 63) with λ = 1.

Funnel potential
As we stressed in Sect. 4.6.7, the funnel potential V (r) = ar − b/r, with a and b both positive, is of crucial importance in hadronic physics. It was shown that an analytical expression is available in a nonrelativistic framework. However the glueball resonances or the spectroscopy of very light mesons need a semirelativistic description and a spinless Salpeter equation governed with a funnel potential is the most economic way to consider this situation.

Massless particle
The ultrarelativistic limit (m = 0), is particularly important for the description of a glueball composed of two massless gluons or a meson composed of two massless quarks. Then, the extremization condition (5 . 22) becomes ar 2 0 + b = σ Q whose solution is Inserting this value in the expression of the mass (5 . 23) leads to the mass formula It is astonishing that such a complicated potential used in conjunction with a semirelativistic kinetic energy admits analytical approximate eigenenergies of such simple form. For b = 0, (5 . 86) reduces to 2 √ aσ Q, that is the expected expression from (5 . 63) in the case λ = 1. For a = 0, M ur → 0 as shown in Sect. 5.3.5. If the prescription Q = Q 2 is chosen, then (5 . 86) is an upper bound on the exact result. As already mentioned, other choice for Q could improve the accuracy of the formula but without guarantee about the variational character of the result. This point is studied in the next section.

Improved formulae for the ultrarelativistic limit
Written with dimensionless variables, the ultrarelativistic Hamiltonian for two massless particles interacting via the funnel potential is given by with β ≥ 0. The approximate energy spectrum is given by (5 . 86) with σ = 2, a = 1 and b replaced by β (to avoid confusion with the coefficient of n). This kind of Hamiltonian is often used in hadronic physics with typical values for β ≈ 0.4. Very accurate numerical values for the eigensolutions of this Hamiltonian can be obtained with the Lagrange mesh method. 54) With the choice Q = Q 2 = 2 n + l + 3/2, upper bounds are obtained. As shown in previous cases, another choice for the n-and l-dependences of Q can greatly improve the results. By using the form Defining a new variable x and a new parameter α by the previous equation can be recast as a third order polynomial equation: In consequence, it is analytically solvable and can be put under the reduced form (B . 1).
Reporting the obtained x value in the expression (5 . 23) for the mass leads to the final expression. Explicitly, the procedure to get the analytical AFM expression for the spinless Salpeter equation based on the Hamiltonian (5 . 91) needs the following steps. Calculate the α value from (5 . 93), then introduce the parameter The x value is obtained from the expression (5 . 96) while the AFM mass is given by Once again, one sees that the AFM is very powerful to get analytical expressions, even for quite sophisticated problems. We have checked that (5 . 86) [with σ = 1, a = 1 and b = β ] is recovered in the limit χ → 0. It can also be checked that (5 . 97) and (4 . 109) [with σ = 1, m = χ, a = 1 and b = β ] tend towards the same limit when χ → ∞.

Low mass expression
By a "small mass", we mean the condition m ≪ √ a. For m = 0, we saw already that r 0 = (σ Q − b)/a. The principle of the method is based on a limited expansion of the equations at first order in m 2 . Thus, we set (1 + ε). Inserting these values in the expression (5 . 23) for the mass, we obtain the final result Here again, it is amazing that the approximate eigenvalues of a so sophisticated Hamiltonian take such a simple form.

Unequal masses
Some particular problems require to deal with a system of two particles with unequal masses. In this case, a general spinless Salpeter Hamiltonian is given in the rest frame by H = p p p 2 + m 2 1 + p p p 2 + m 2 2 +V (r).

(5 . 101)
A general result can be obtained for Hamiltonians with two different masses. Let us consider the following two-body Hamiltonians H = T 1 + T 2 + V , H 1 = 2T 1 + V , and H 2 = 2T 2 + V whose ground-state energies are respectively M = φ |H|φ , M 1 = φ 1 |H 1 |φ 1 , and The square roots appearing in the kinetic terms can still be treated by resorting to the AFM. But this time, two auxiliary fields, µ 1 and µ 2 , have to be introduced. One is led to the HamiltonianH playing the role of a reduced mass. It has been shown in Ref. 85) that the AFM leads to the system (5 . 30)-(5 . 32), but with T (x) = x 2 + m 2 1 + x 2 + m 2 2 . It is then easy to obtain numerical upper bounds, but one can wonder if it is possible to obtain analytical solutions? Due to the presence of the square roots, it is only conceivable with very special conditions: one mass is vanishing (m 1 = 0 and m 2 = m) and the potential is well chosen. Analytical solutions can be obtained for the funnel potential, often used in hadronic physics. 85) Let us note that an analytical upper bound for a logarithmic potential (also relevant in hadronic physics) can also be computed.
The solution for the general funnel potential is very complicated. As an illustration, the case of the Coulomb interaction, −a/r, is only presented here. The solution for the associated system (5 . 30)-(5 . 32) gives This formula is an upper bound if Q = n + l + 1, but, as we pointed out several times, it is justified to take a more sophisticated expression in order to get a better accuracy. It can be seen from (5 . 107) that lim m→0 M = 0, as expected. From equations above, a solution M exists only when the following condition is verified: The right part of this inequality corresponds to the cancellation of the binding (M → m when Q → a), and the left part corresponds to a collapse (M → 0 when Q → a/2). Indeed, the mean radius r 0 → ∞ when Q → a, and r 0 → 0 when Q → a/2. A sufficiently strong interaction must exist to bind the system, but unphysical values of the mass could appear if it is too strong. With the Coulomb potential, the nonrelativistic limit is not defined by m → ∞ but by a → 0. So, (5 . 108) implies that this limit is irrelevant for the solution (5 . 107).

Semiclassical interpretation
A semiclassical interpretation of the system (5 . 30)-(5 . 32) is also possible. Let us assume a classical circular motion for the two particles, as illustrated in Fig. 4. In this case, the force F i acting on particle i is given by Taking into account that both particles are characterized by the same momentum p 0 , This equation becomes The rigid rotation constraint, v 1 /r 1 = v 2 /r 2 , implies that If the force acting on i comes from the potential V (r) generated by j, then F 1 = F 2 = V ′ (r 0 ). (5 . 110) and (5 . 111) can be recast onto the form (5 . 32), and it is obvious than (5 . 30) gives the mass of the system. The total orbital angular momentum is r 1 p 0 + r 2 p 0 = r 0 p 0 . A semiclassical quantification gives thus r 0 p 0 = L + 1/2, and we obtain a system very similar to (5 . 30)-(5 . 32). Nevertheless, the AFM produces more general results. Up to now, the AFM was applied essentially to one-and two-body systems with either a nonrelativistic or a relativistic kinetic energy operator. It was shown that, in many situations, the method was able to provide an accurate analytical expression for the corresponding eigenenergies. In this section we address the important point of whether AFM is also able to give informations concerning the eigenenergies of a system composed of N particles. Except the very special case of particles whose dynamics is governed by quadratic interactions, which is treated in detail in Appendix J, no analytical solution is known for N-body systems. Nevertheless a number of very interesting physical problems need to solve a N-body equation. As in the 2-body case, the AFM cannot pretend to give a very high accuracy in the results, but it exhibits the interesting virtue to clearly show the dependence of the results in terms of the physical parameters and the various quantum numbers.
Indeed, the envelope theory was used in the N-body problem to find bounds for the nonrelativistic binding energies if the two-body potential fulfills certain restrictive conditions. 86) In this section, we apply AFM to the N-body problem in a more general context • we consider the possibility of dealing with a semirelativistic kinetic energy potential; • we allow the presence of one-body terms in addition to two-body interactions. Moreover, we are not interested only by the existence of bounds which can sometimes give rather poor values, but the purpose of the present section is to show that the AFM can be successfully applied to find approximate analytical mass formulae for general relativistic

N-body Hamiltonians of the form
V i j (|r r r i − r r r j |), (6 . 1) with ∑ N i=1 p p p i = P P P = 0 0 0, where r r r i and p p p i are respectively the position and momentum of particle i with a mass m i , and where R R R is a global variable defined in Appendix I. As stated, we consider a relativistic kinematics and allow for both pairwise V and one-body U interactions. Particular systems exhibiting one-body potentials of the previous form will be studied below. Most of the material of this section can be found in Ref. 68).
A quantity which appears very often in this section is the number of interacting particles. In order to simplify the notations, we use a special symbol to designate this number and set 6.1. Main mass formulae 6.1.1. General case As stated, the AFM can be straightforwardly generalized to the case of N-body Hamiltonians of the form (6 . 1). The basic idea is to introduce auxiliary fields so that this Hamiltonian is formally replaced by a Hamiltonian for which an analytical solution can eventually be found. The auxiliary fields are denoted as µ i , ν i , and ρ i j , and are introduced as follows where (notice the definition of r i = |r r r i |) r i = |r r r i − R R R|, r i j = |r r r i − r r r j |, (6 . 4) and where where the prime denotes the derivative with respect to the argument. It is sufficient to set µ k = m k ∀k to treat a nonrelativistic kinematics. We stress that, in order for the AFM to apply, it is assumed that the mass spectrum of (6 . 3) is analytically computable. In practice, this is only possible for the choice P(x) = S(x) = x 2 under some conditions (see Appendix J).
The reason of such a definition is that Hamiltonians (6 . 1) and (6 . 3) are equivalent if the auxiliary fields are properly eliminated. It can indeed be checked that and finally that H(μ k ,ν k ,ρ kl ) = H. The approximation underlying the AFM is now that, as we have done in previous sections, the auxiliary fields will be seen as real variational parameters. Let us note where µ k,0 , ν k,0 , ρ kl,0 are the real values of the auxiliary fields such that M 0 is extremal. |ϕ 0 is an eigenstate with fixed quantum numbers and a given symmetry. Using the same mathematical techniques than in Sect. 2, one can show that µ k,0 = ϕ 0 |p p p 2 i + m 2 i |ϕ 0 with P P P = 0 0 0, (6 . 10) ϕ 0 |P(r i )|ϕ 0 = P(r i,0 ) with r i,0 = I i (ν i,0 ), (6 . 11) ϕ 0 |S(r i j )|ϕ 0 = S(r i j,0 ) with r i j,0 = J i j (ρ i j,0 ).
(6 . 12) Using the Hellmann-Feynman theorem, 39) it can be shown that ϕ 0 |μ 2 i |ϕ 0 = µ 2 i,0 , ϕ 0 |P(I i (ν i ))|ϕ 0 = P(I i (ν i,0 )) and ϕ 0 |S(J i j (ρ i j ))|ϕ 0 = S(J i j (ρ i j,0 )). Let us rewrite the Hamiltonian (6 . 3) under the form where T (µ k ) stands for the kinetic part. One can also check that This means that the approximate potentialŨ i (Ṽ i j ) and the corresponding genuine potential . From now on, this choice is adopted for the rest of this section. Second, the extremization of the eigenvalues with respect to the real auxiliary fields is a classical numerical problem which can be solved with a high accuracy.

Identical particles
To obtain analytical closed solutions of the eigenvalue and extremization problems associated to Hamiltonian (6 . 3), it is necessary to simplify the system. First of all, we will only consider systems with identical particles, that is with m i = m. In this case, it is reasonable to consider identical interactions between them, namely U k (x) = U (x) and V kl (x) = V (x). This means that K k (x) = K(x), I k (x) = I(x), L kl (x) = L(x), and J kl (x) = J(x). The global variable R R R (see (I . 1)) is then the center of mass of the system, even for relativistic kinematics.
Let us denoteP i j the permutation operator exchanging particles i and j. We can writê P i j |ϕ 0 = ±|ϕ 0 if this state |ϕ 0 is completely (anti)symmetrical. Then P(r j,0 ) = ϕ 0 |P(r j )|ϕ 0 = ϕ 0 |P i j P(r i )P i j |ϕ 0 = ϕ 0 |P(r i )|ϕ 0 = P(r i,0 ). (6 . 16) If P(x) is monotonic, which is always the case in practice, then r j,0 = r i,0 . Finally, if I(x) is invertible, which must be the case to solve the problem, we have ν i,0 = K(r i,0 ) = K(r j,0 ) = ν j,0 . So all optimal values ν i,0 are the same. Using the same reasoning, we can draw the same conclusion for other auxiliary fields. Under these conditions, we can set µ i = µ, ν i = ν, and ρ i j = ρ in the expression of the eigenenergies. We are thus led to replace the original Hamiltonian (6 . 3) by the following simpler Hamiltonian, which now depends on only 3 auxiliary fields (remember that P( which, by virtue of (J . 16), has the following mass spectrum (P P P = 0 0 0) which is characterized by a high degeneracy due to the form of the principal quantum number Q. The last step needed to get the final mass formula is to find the optimal values µ 0 , ν 0 , and ρ 0 from the extremization conditions Writing explicitly (6 . 20) and after some algebra, we arrive at the mass formula where the auxiliary fields are solutions of At this stage, functions I and J are not known. They depend on the specific forms of U and V . Optimal values of µ 0 , ν 0 , and ρ 0 are not known anymore. Moreover, they also depend on the state considered through the variable Q. Nevertheless, formula (6 . 21) makes clearly appear the mean-field nature of the AFM, because µ 0 can be interpreted as the average kinetic energy of one particle, while U (I(ν 0 )) and V (J(ρ 0 )) can be respectively seen as the average potential energy of one particle in the potential U (r) and of a pair in the potential V (r).
A general remark about the AFM can now be done: The analytical results obtained for a given potential U (r) approximated by a potential P(r) can be used as a starting point for finding analytical results for an other potential W (r), approximated this time by U (r). The same procedure can be applied for the potential V (r). This considerably enlarges the domain of applicability of this method.
Notice that the idea of rewriting a N-body Hamiltonian with pairwise interactions of the form g(r 2 i j ) as a N-body harmonic oscillator has already been investigated in Ref. 87) within the framework of the envelope theory. That method shares many similarities with the auxiliary field method, as shown in Appendix A, but the results presented hereafter have, to our knowledge, be presented for the first time in Ref. 68).

Simplified formulations
The four equations (6 . 21)-(6 . 24) give the AFM approximation for the mass spectrum of a N-body system where the constituent particles are identical. They can be even simplified by defining X 0 = 2µ 0 (ν 0 + Nρ 0 ). (6 . 25) We have indeed where the remaining auxiliary field is a solution of and where K(x) = U ′ (x)/(2x) and L(x) = V ′ (x)/(2x) (see (6 . 5)). This condition ensures that M(X 0 ) is extremal. In so doing, we have simplified a lot the original formulation. The set (6 . 22)-(6 . 24) is a complicated system of non linear coupled equations, while the equivalent formulation (6 . 27) is just a transcendental equation which is a priori easier to be solved. As we shall see in the next sections, this task can be achieved in several cases of interest. It is worth mentioning that the AFM results (6 . 26) and (6 . 27) are valid not only for the ground state but also for excited states.
For particles with null mass (ultrarelativistic limit), one obtains even simpler formulae by simply setting m = 0 in formulae (6 . 26) and (6 . 27). In the nonrelativistic limit, the auxiliary field µ tends towards m. In this case also the various formulae look simpler but cautions must be taken in the limit. Explicitly, (6 . 26) reduces to with m t = Nm and where X 0 is the solution of Let us introduce the distance r 0 = NQ/X 0 and the momentum p 0 = Q/r 0 , and let us define T by T (x x x) = √ x x x 2 + m 2 or by its nonrelativistic counterpart m + x x x 2 2m . It is a simple algebra exercise to show that formulae (6 . 26)-(6 . 27) can be written as 88) With this formulation, an AFM eigenvalue given by (6 . 30) is simply the kinetic operator evaluated at the mean momentum p 0 plus the potential energy computed at some mean radius depending on r 0 . This generalizes the relations found previously for two-body systems. As one could expect, the kinetic energy and the one-body potential energy are proportional to the number of particles and the two-body potential energy is proportional to the number of pairs. Formula (6 . 30) looks like a semiclassical approximation but this is absolutely not the case. The AFM yields an approximate N-body wavefunction, 68), 89), 90) and the relation (6 . 31) between p 0 and r 0 is a full quantum link, function of the quantum numbers of the system. Lastly, the value of r 0 (and thus of p 0 ) is the solution of a transcendental equation (6 . 32) which is the translation into the AFM variables of the generalized virial theorem 59) which comes from very general properties of quantum mechanics. These considerations prove that the AFM really relies on very sound physical basis. Once the system (6 . 30)-(6 . 32) is written, it can appear finally quite natural to obtain such a result. The problem is to find a relevant link between the mean values r 0 and p 0 . This is solved by the AFM. The eigenstates of H(µ 0 , ν 0 , ρ 0 ) are built with harmonic oscillator states. A N-body state with a mass M(X 0 ) is written where ψ n j l j (λ j , x x x j ) is a three-dimensional harmonic oscillator wavefunction (see Appendix E), depending on the Jacobi coordinate x x x j (see Appendix I) and decreasing asymptotically like e −λ 2 j x x x 2 j /2 (the magnetic quantum number is omitted). The quantum numbers are such that ∑ N−1 j=1 (2n j + l j ) = B = Q − 3 2 (N − 1) and the scale parameters λ j are given by The quantum number B is generally called the band number. The state (6 . 33) has neither a defined total angular momentum nor a good symmetry, but its is characterized by a parity (−1) B . By combining states (6 . 33) with the same value of B (or Q), it is generally possible to build a physical state with good quantum numbers and good symmetry properties, but the task can be technically very complicated. 91) Some observables can be easily computed: Since these results only depend on the quantum numbers via Q, and since a physical state must be a combination of eigenstates with the same value of Q, formulae (6 . 35) and (6 . 36) are also valid for a physical state. This shows that r 0 can be considered as a mean radius for the system and p 0 as a mean momentum per particle. Indeed, (6 . 35) and (6 . 36) imply that for arbitrary i = j since the mean values are taken with completely symmetrized states. These results can also be obtained using the more general relations (6 . 10)-(6 . 12) relevant for P(x) and S(x) different from x 2 . Identifying (6 . 21) and (6 . 30), it appears that µ 0 = m 2 + p 2 0 . (6 . 38) In the above expression, Q is the principal quantum number defined in (6 . 19). However, it has been shown many times for the two-body problem that a much better approximation of the exact energies can be obtained with a slight modification of the principal quantum number. A particularly simple form which seems to work quite well for N-body systems is given by For instance, such a formula is tested numerically for baryons in Sect. 6.4.3; it improves substantially the results as compared to the expression (6 . 19). It is worth noting that there is no systematic procedure to determine the values of parameters α, β and γ, which depend on both the interaction and the kinematics.

Critical constants
Some interactions, as the exponential or the Yukawa potentials, admit only a finite number of bound states (see Sects. 4.4.2 and 4.5.2). Let us assume that such an interaction can be written as W (x) = −κ w(x), where κ is a positive quantity which has the dimension of an energy and w(x) a "globally positive" dimensionless function such that lim x→∞ w(x) = 0. The critical constant κ({θ }), where {θ } stands for a set of quantum numbers, is such that, if κ > κ({θ }), the potential admits a bound state with the quantum numbers {θ }. The interaction energy for the state with quantum numbers {θ } is then just vanishing for κ = κ({θ }). 66) Let us consider a nonrelativistic N-body system (no manageable calculation can be performed for a semirelativistic kinematics) with one-body potentials U (x) = −k u(x) and two-body potentials V (x) = −g v(x), both independent of the particle mass and both admitting only a finite number of bound states. The system (6 . 30)-(6 . 32) for a vanishing energy gives 41) where k N and g N are the critical constants for the system with N particles. The elimination of the ratio N Q 2 /(m r 2 0 ) from both equations yields the equality When potentials u and v are both taken into account, nothing interesting can be said. So let us consider one type of potential at once. Assuming that only two-body forces are present, (6 . 42) reduces to 88) 2 where the parameter g N has disappeared. Introducing the new variable y 0 = r 0 / √ C N , we can rewrite (6 . 40) and (6 . 41) as The variable y 0 , determined by (6 . 45), is independent of N, Q and m, and depends only on the form of the function v(x). So, (6 . 44) gives precise information about the dependence of the many-body critical constant g N as a function of all the characteristics of the system. With the system (6 . 44)-(6 . 45), it is easy to recover previous AFM results obtained in Sects. 4.4.2 and 4.5.2, and in Refs. 67), 68). Within the AFM approximation, the ground state (GS) of a boson-like system is characterized by Q = 3 2 (N − 1). We obtain in this case the following very general relation valid, at the AFM approximation, for all pairwise potentials with a finite number of bound states This ratio has previously been obtained and numerically checked for several exponentialtype potentials. 92), 93) Similarly, in the same general situation, g N (GS) = 2 N g 2 (GS), (6 . 47) indicating that in order to bind a N-body system, a coupling N/2 times smaller than the coupling for a two-body problem is sufficient. 92), 93) Assuming that only one-body forces are present, a similar calculation gives 88) 2 u(y 0 ) + y 0 u ′ (y 0 ) = 0, (6 . 49) where the change of variable y 0 = r 0 /N has been used. Again, (6 . 48) gives precise information about the dependence of the one-body critical constant k N as a function of all the characteristics of the system. These results are strongly different from those for pairwise forces.
If the AFM gives upper (lower) bounds for the exact eigenvalues, the critical constants predicted by formulae above are upper (lower) bounds for the exact critical constants. In the favorable situation where the AFM gives both upper and lower bounds for the eigenvalues, it is possible to approximate from above and below the exact critical constants.

Connection with the perturbation theory
It has been shown in Ref. 40) that, for one-and two-body nonrelativistic systems, the AFM and the perturbation theory give similar results when the potential is an exactly solvable one plus a small perturbation. This result is extended here for the general Hamiltonian (6 . 1). 88) Let us first assume that each pairwise potential V (r i j ) is supplemented by a term ε v(r i j ), with ε ≪ 1 in order that ε v(x) ≪ V (x) in the physical domain of interest. In the system (6 . 30)-(6 . 32), the potential V (x) is replaced by V (x) + ε v(x). In this case, new values r 1 and p 1 for the mean radius and momentum will be the solution of the new system 50) Writing r 1 = (1 + δ )r 0 , we can expect δ ≪ 1 since ε ≪ 1. In this case, power expansions at first order can be computed. We have p 1 ≈ (1 − δ )p 0 from (6 . 51), and we can write This result could seem quite obvious, but it demonstrates that the knowledge of r 0 is sufficient to obtain the contribution of the perturbation at the first order.
Let us now assume too that each [one-body potential U (|s s s i |) / kinetic operator T (|p p p i |)] is supplemented by a term [η u(|s in the physical domain of interest. With similar calculations, we finally find The parameter δ is determined at the same order by the following relation Perturbed observables and wavefunctions can then be computed at first order, since r 1 = (1 + δ )r 0 and p 1 = (1 − δ )p 0 at this order. The contribution of a perturbation at the first order can thus be very easily computed within the AFM once the unperturbed problem is solved. In order to check the quality of this approximation, let us consider a case in which the unperturbed Hamiltonian H can be solved exactly by the AFM, that is M 0 is the exact solution. If the small perturbation potential is written ε ∑ N i< j=1 v(r i j ), the quantum perturbation theory says that the solution M * is given by for any pair (i j). The mean value is taken with a completely symmetrized eigenstate of the unperturbed Hamiltonian H. The comparison of (6 . 56) with (6 . 53) shows that v(r i j ) is replaced by v r 0 / √ C N within the AFM. This is to be compared with the exact relation S(r i j ) = S r 0 / √ C N for the auxiliary potential (see (6 . 12)). So, the AFM does not give the same result as the perturbation theory. But the agreement can be very good, as shown with several examples calculated above. Similar discussions can be made for small one-body perturbation potentials or small perturbations of the kinematics.

General results
The first explicit example that will be considered below is the case of power-law potentials, i.e. the Hamiltonian where λ , η ≥ −1 (in the nonrelativistic case, one can further consider λ , η > −2). When only a one-body or a two-body interaction is present, parameters a or b must be positive. If both types of potentials are present, there are less constraints on the sign provided that a bound state can exist. Following the definitions (6 . 5), it is readily computed that K(x) = a |λ | x λ −2 /2 and L(x) = b |η| x η−2 /2. Then, by defining The sufficient condition to get a closed analytical formula for M(X 0 ) is to solve (6 . 60) analytically. This is possible if this last equation can be rewritten as a polynomial equation of the fourth degree at most. If λ = η, it can be computed that the following couples will lead to such an analytically solvable equation Note that the last three cases are only allowed with a nonrelativistic kinematics. The problem actually becomes particularly simple when λ = η. In this case, the following values lead to an analytical solution The last six values are only allowed in the case of a nonrelativistic kinematics. The solution of (6 . 60) reads 63) and the mass formula (6 . 59) becomes Among the values (6 . 62), three cases are of obvious physical interest: λ = −1, 1, and 2. The case λ = 2 corresponds to the harmonic oscillator and is solved in Appendix J. The Coulomb problem, i.e. λ = −1, will be specifically considered in Sect. 6.4.4, while we will discuss the case λ = 1 in Sect. 6.4.2. Notice that closed mass formulae for any λ can be obtained in the nonrelativistic and ultrarelativistic limits, as we also show in the following. In the case m = 0, (6 . 63) allows us to extract the X 0 value for arbitrary power (6 . 65) Inserting this value in (6 . 64) provides us with the expression of the energy in the ultrarelativistic limit We already know from Sect. 5.3.4, that a bound state of massless particles exists only if λ > 0. Setting N = 2 and a = 0, this equation reduces to (5 . 63) with σ = 2 (and a renamed b).
On the contrary, assuming a large mass, m t ≈ M(X 0 ), it is possible to extract X 0 from (6 . 63) Except the Coulomb case (λ = −1) which will be treated subsequently, the term m 2 is dominant with respect to X 0 . Expanding the expression of the energy (6 . 64) at lower order in 1/m leads to the nonrelativistic value of the energy (6 . 68) Setting N = 2 and a = 0, this equation reduces to (4 . 7) with the reduced mass m replaced by m/2 (and a renamed b).

The linear potential
The case λ = 1 corresponds to a linearly rising confining potential. It is of great interest in hadronic physics since a linearly rising potential appears to be the best way of modeling the QCD confining interaction within potential approaches, see for example Ref. 94) for more details. For the three-body problem, the one-body potential (term in a) corresponds to the so-called Y-junction, while the two-body potential (term in b) corresponds to the so-called ∆ approximation. For physical problems, one has the choice of retaining one approximation or the other or both. Here for sake of generality we present the general case including both. In the real QCD world, one must take N = 3. But a general N value is interesting to study alternative approaches of QCD. 89), 90) For λ = 1, one can write The solution of (6 . 63) when λ = 1 is given by where F − (Y ) is given by (B . 2) and where Introducing this value in the expression of the energy (6 . 64), a simple calculation allows us to obtain the mass under the form . (6 . 72) In the ultrarelativistic limit, where m = 0, one obtains, from (6 . 66), the very simple relation M(X 0 ) 2 = 4Nc Q. (6 . 73) Such a linear behavior of the square mass versus the principal quantum number is a wellknown fact in hadronic physics, where light mesons and baryons are known from experimental data to exhibit Regge trajectories.
In the nonrelativistic limit, the mass formula (6 . 68) can be recast under the form The second term is the "binding energy", which is a positive quantity in case of a positive linear potential.
In the special case of a linear potential, there is a peculiar relationship between the mass of the two-body problem and the mass of the N-body problem. More precisely, let us call M (2) (σ , m, b, Q 2 ) the AFM eigenvalues of the Hamiltonian H (2) = σ p p p 2 + m 2 + br (6 . 75) calculated with the natural quantum number Q 2 = 2n + l + 3/2. Then, using (6 . 73) and its counterpart for N = 2 with the value of the quantum number Q defined by (6 . 19), the mass M (N) (N, m, a, b, Q) of the original Hamiltonian is given by

Baryonic case
In this section, we want to discuss a special situation which can have an immediate application for hadronic systems, especially baryons or glueballs. 89), 90) Let us consider a one-body potential with a linear shape and a two-body potential of Coulomb type (λ = 1 and η = −1). We noticed that an analytical solution does exist, but one needs to solve a general polynomial of degree 3 and the corresponding solution is quite involved. We prefer to give here the solution in the ultrarelativistic limit (m = 0). In this special case the solution of (6 . 60) is quite simple and (6 . 59) can be simplified as follows Since the argument of the square root must be positive, it is readily checked that, either there exists a maximal allowed number of particles for a fixed value of b, either there exists a maximal allowed value of b for a fixed number of particles. This kind of mass formula has been applied to the computation of light baryon masses for various theories of QCD with a large number of colors. 89), 90) In order to test the relevance of our method, we apply it for a particular relativistic three-body system: the light baryon composed of three massless quarks. In the framework of constituent models, the Hamiltonian for such a system is given by 95), 96) The dominant interaction is a confinement by three strings with density σ s meeting at the center of mass. The short-range part is given by pairwise interactions of Coulomb nature with a strong coupling constant α S . In order to allow a more realistic calculation of the correct baryon masses, this Hamiltonian must be completed with terms that can be computed in perturbation. 95) Since these contributions represent no interest for the AFM, they are not considered in this study. Let us denoteH B the AFM counterpart of H B . Approximations of eigenmasses of the genuine Hamiltonian are given by the application of (6 . 78) with N = 3, a = σ s and b = 2α S /3, where B = 2(n 1 + n 2 ) + (l 1 + l 2 ) is the band number (see (6 . 19)). One can notice the strong degeneracy due to the particular form of B. This formula generalizes a result obtained for a small value of α S in Ref. 95). Following above considerations, these masses are upper bounds on the exact masses.
In order to test the accuracy of this formula, it is necessary to compute accurately eigenvalues of H B . It is possible, for instance, to use a variational method relying on the expansion of trial states with a harmonic oscillator basis. 97) We can write |φ (B, q(B)) , (6 . 81) where B characterizes the band number of the basis state and where q summarizes all the quantum numbers of the state (which can depend on B). This procedure is specially interesting since the eigenstates of H B are expanded in terms of eigenstates ofH B (up to a length scale factor). In practice, a relative accuracy better than 10 −4 is reached with B max = 20. Such results are denoted "exact" in the following. Before comparing exact masses with formula (6 . 80), it is significant to describe the baryon wavefunctions. Quarks are fermions of spin 1/2 with isospin and color degrees of freedom. The global color function is unique and completely antisymmetrical. For total spin S = 3/2 (isospin T = 3/2) the spin (isospin) wavefunction is completely symmetrical. The corresponding wave function are of mixed symmetry for S = 1/2 or T = 1/2. Using degenerate eigenstates ofH B , it is always possible to build three-quark states completely antisymmetrical, 91), 97) which are characterized by the same mass. For instance, when S = T = 3/2 or S = T = 1/2, the baryon states possess a spatial wavefunction completely symmetrical. In the real world, the degeneracies are removed by spin-or isospin-dependent operators which can be computed as perturbations. 95) When an eigenstate of H B is computed for given values of S, T , and the total angular momentum L, we determine the band number of its dominant components. It can then be compared with the eigenstate ofH B with the same spin-isospin quantum numbers and the same band number.
In Table XIV, exact masses are compared with the predictions of formula (6 . 80) for usual values of the parameters σ s and α S . The relative error can be large but the quality of the approximation is quite reasonable for a so simple formula. Due to the strong degeneracy of the harmonic oscillator, all AFM eigenstates with the same band number have the same energy. This approximation is better for high values of L. Let us note that the relative error can be reduced by a factor around 2 when α S → 0. When α S = 0, the relative error becomes independent of the energy scale factor √ σ s .   We have shown many times that it is possible to improve two-body mass formulae by changing the structure of the principal quantum numbers Q. By fitting another form on exact eigenvalues, a very high accuracy can sometimes be reached. Here, we will proceed differently and will try to use analytical results to find the best shape for the mass formula. Using a good upper bound on the ground state given by the trial state of (6 . 81) at B max = 0, a modified AFM formula was proposed One can see in Table XIV that the masses are greatly improved, but the variational character of the formula (6 . 82) cannot longer be guaranteed: some masses are now below the exact ones. At last, using information coming from the WKB method, 56) the problem of the strong degeneracy has been improved by a new formula Details about these procedures are given in Ref. 68). One can see in Table XIV, that the relative error is now around 1% and generally below. Despite its simplicity and its nonvariational character, (6 . 83) is then a very good mass formula for the eigenstates of the Hamiltonian (6 . 79). It is not sure that the procedure used here to improve the mass formula for baryons could work so well for other Hamiltonians. But, this shows that an improvement is possible, at least in some particular cases.

Atomic systems
We turn now our attention to atom-like systems, i.e. systems described by the following Hamiltonian For this section, we remain completely general, the only hypothesis being that the parameters α andᾱ should be real positive numbers. The obvious application could correspond in a first approximation to N identical electrons, with charge e, feeling the attraction of a central static source (the nucleus) and their own repulsion. In this particular situation, one must make the identification α = Ne 2 andᾱ = e 2 . At this stage, we point out that our formalism is spin-independent and thus that this model cannot be directly applied to a real physical atom. Moreover, one should be very careful in applying the formulae already given. Two cautions are in order: • In (6 . 84), r i must be the distance between the particle i and the center of mass of the system of N electrons. In the atomic interpretation, this center of mass must coincide with the nucleus. This is not always the case but this prescription should be valid in case of a spherical atom. • Electrons are fermions and the total wavefunction (space and spin) must be completely antisymmetrical. Since the spatial wavefunction cannot be completely symmetrical, all possible values are not allowed for Q given by (J . 17). For closed shell atoms, Q AGS (J . 19) should be used to estimate the ground-state energy only. Nevertheless, solving Hamiltonian (6 . 84) has an intrinsic interest since, to our knowledge, no corresponding analytical mass formula is known so far. A general solution can moreover be found with the relativistic kinematics, starting from (6 . 63) and (6 . 64). Applying them to the case λ = −1, a = α, b = −ᾱ and defining the solution of equation (6 . 63) is given by The mass formula (6 . 64) then reads In the ultrarelativistic limit, this expression gives the value M = 0. This property is well known in the two-particle system and prolongates to the N-body problem. It is quite obvious that a system of massless particles cannot be characterized by a (negative) binding energy since the resulting mass would be negative. Moreover, with a Coulomb-type potential, the only energy scale of the problem is the particle mass. So, if this mass vanishes, no bound state can exist. In the nonrelativistic regime, one must assume that the D quantity is small in (6 . 87). In this case It is clear from formula (6 . 87) that the value of the parameters cannot be arbitrary. One must satisfy the condition D < 1. Explicitly this means Let us emphasize that this condition is independent of the mass of the particle. 98) Let us notice that, for α = Ne 2 andᾱ = e 2 , a maximal value of N exists if Q increases less rapidly that N 2 .
6.5. Duality relations 6.5.1. Generalities In this section, we will show that AFM is a very powerful method to give relationships between the energies of states (excited or ground states) for two different systems. We call the relations between the energies of both systems "duality relations". The case for which one of the two systems contains 2 particles is specially interesting since the corresponding eigenenergies are rather easy to obtain. These relations are exact if we consider the energies obtained in the AFM; they are no longer exact for the true eigenstates, but if AFM gives a good precision on the true states, we hope that the duality relations are good approximations of the physical situations and gives, at least, general trends for the exact levels. A rather complete treatment of duality relations can be found in Ref. 99). We present here a simplified version where we restrict our discussion to systems of N identical particles interacting only through a two-body potential V (r). Readers interested in more specific details may have a glance to Ref. 99) for further informations.
Our original Hamiltonian is thus V (|r r r i − r r r j |). (6 . 90) In this expression r r r i is the position operator for particle i, p p p i its conjugate momentum and m the common mass of all particles. The AFM mass of the system is given by (6 . 26) which is written now as while the X 0 quantity is solution of the transcendental equation (6 . 27) which turns out to be now 92) and where L(x) = V ′ (x)/(2x) as usual (see (6 . 5)).
In general, the potential V (r) depends on physical parameters {τ 1 , τ 2 , . . . , τ p } ≡ {τ} and should be noted more precisely V ({τ}; r). The other parameters of the problem are the number of particles N and the mass m of the particles. Lastly, we want to have a description of the whole spectrum, so that the principal quantum number Q also enters the game. In consequence, a complete notation for the eigenmasses would be M (N) ({τ}; m, Q).
Let us stress now a very important point: for realistic potentials, the parameters {τ} could depend on N or/and m. We do not consider such a particular behavior here. Thus we assume that the parameters of the potentials are independent of N and m. This means that the m dependence of H is only through the kinetic energy and its N dependence through the numbers of terms in the summation. We suppose that the potentials are given once for all and consider that their form and parameters do not vary for all studied systems. Therefore we are finally interested only by the N, m, Q dependence of the eigenmasses and use the simplified notation M (N) (m, Q). What we call a "duality relation" is just a relation between M (N) (m, Q) and M (p) (m ′ , Q ′ ).

General case
Instead of the X 0 quantity, let us define the new variable Then, the transcendental equation (6 . 92) and the mass expression (6 . 91) are replaced by the following ones A little algebra on those equations allows to arrive at the following duality relation between the N-body and the p-body systems with particles interacting via the same twobody interaction In this case, the spectrum of the N-body system is the same as the p-body system (with the same two-body potential) provided we consider different particle masses in both situations and different excitation states.

Ultrarelativistic limit
The case of ultrarelativistic systems, characterized by a vanishing mass m = 0, presents some very specific and interesting features. The case of systems composed of gluons or/and light quarks can be well represented in this scheme. The Hamiltonian of the system is then V (|r r r i − r r r j |).
while the corresponding AFM eigenmass is given by The eigenmass, depending now only on N and Q, will be noted M (N) u (Q) (the index u stands for "ultrarelativistic").
With the same definition (6 . 93) for s 0 , the transcendental equation (6 . 98) reduces to Defining the D(x) and F(x) function as in (5 . 27) and (5 . 28) calculated with the potential V (x), one has As we already mentioned, the D and F functions are universal. From expression (6 . 102), one deduces immediately the duality relation which is the special case of (6 . 96) with m = 0.
6.5.4. Nonrelativistic limit Another interesting limit of the theory is the nonrelativistic one, valid when the mass of the particles is large compared to the mean potential. In this case, the considered Hamiltonian is simply V (|r r r i − r r r j |). (6 . 104) Instead of dealing with the total mass M (N) , it is better to consider the binding energy obtained by removing the total rest mass: E (N) = M (N) − Nm. The AFM approximation E (N) (m, Q) of the binding energy is given by the following equations, and We will see that, in this nonrelativistic limit, one obtains additional interesting properties. Introducing s 0 by (6 . 93) as above, the transcendental equation (6 . 105) reduces to Defining the D(x) and F (x) function as in (4 . 11) and (4 . 13) calculated with the potential V (x), one has As the case studied in the previous section, the D and F functions are universal. For a nonrelativistic system, an additional property appears. Remarking that E (N) depends only on the ratio Q 2 /m, one has the general duality relation valid for any value of the real parameter β . From expression (6 . 109), one deduces immediately the duality relation which is identical to the general case (6 . 96). One can use the property (6 . 110) to obtain many other possibilities. For example, choosing the value β = p(N − 1)/(N(p − 1)), one has the alternative simpler duality relation In this last relation, let us choose β = N/p, one arrives at the relation In this expression, we decide to keep the same mass for systems with N and p particles. The duality relation leads to a link between different excited states of both systems. Choosing the value β = C N /C p in equation (6 . 112), one obtains the alternative expression In this expression we decide to maintain a one to one correspondence in the spectrum but for systems with different particle masses.
6.5.5. Passing from nonrelativistic to ultrarelativistic limits We showed that both the ultrarelativistic limit and the nonrelativistic limit for the eigenmasses share the property of being expressed in terms of universal functions. The F function for the ultrarelativistic case and the F function for the nonrelativistic case are independent of the system and of the excitation quantum numbers, but depends only on the form of the potential under consideration.
If the same form of potential is used in both situations, the expressions of the F function and of the F function are not identical so that the corresponding spectra are quite different. However, we showed in Ref. 99) that if the potentials are different but linked by a certain relationship, one can arrive at very interesting conclusions.
We report here the main conclusions, skipping all the rigorous proofs that can be found in the previous reference. Let us consider a system described by a nonrelativistic treatment based on a two-body potential V (r) whose binding energy is E (N) (m, Q). One knows that If one uses rather the potential W (r) = V (α √ r) , one has If we take the special value S = 2Q/((N − 1) √ C N ), the value C N F(S) represents the ultrarelativistic mass M (N) u (Q) of the same system but obtained with the potential W (r). The condition on S and the link between S and R leads to the condition defining the value of α, namely . From its definition, the potential W (r) depends on m so that the notation M (N) u (Q) which we have used up to now, depends indirectly on m, as imposed by the theorem. 6.5.6. Applying duality relations to exact levels All the duality relations presented above are exact for the AFM solutions of quantum systems. One can wonder to what extent these constraints are satisfied for the corresponding exact solutions. We discuss briefly this point and, again, refer the reader to Ref. 99) for additional information. We are interested essentially to the case of nonrelativistic systems.
Let us call ε (N) (m; {n i , l i }) the exact eigenenergy for state of the Hamiltonian (6 . 90) with {n i , l i } = {n 1 , l 1 , n 2 , l 2 , . . . , n N−1 , l N−1 }, while ε(m; n, l) = ε (2) (m; {n, l}) gives the exact spectrum of the corresponding two-body problem. The first thing to do is to connect the exact level ε (N) (m; {n i , l i }) to the corresponding AFM approximation E (N) (m, Q). Of course, one can always impose an equality between both This equality defines in practice the value of the principal quantum number Q that must be used in the AFM solution. The problem is that the value obtained in that way does depend explicitly on the mass m, in addition to the basic quantum numbers {n i , l i }. This is contrary to the philosophy of the AFM approach for which Q depends on {n i , l i } only. To apply AFM formulas, one must give up the dependence on m for Q and choose a dependence on quantum numbers only, Q({n i , l i }). With such a constraint, the equality (6 . 118) does not hold anymore and is replaced by an approximate value The cleverness of the physicist is to guess the form for Q (N) ({n i , l i }) which makes the previous approximation as precise as possible. The choice (6 . 19) is the most natural one since it follows directly from the construction of the AFM levels. However, the corresponding approximation is not always good and, moreover, it is plagued by degeneracies which are absent in the real physical spectrum. A choice like the one proposed above (6 . 39) could be much better in this respect, but very often the determination of the parameters entering it needs a preliminary study of the problem. Let us suppose that we have the technical ability to choose a practical and correct form for the principal quantum number. The idea of the method is to approximate the exact values by the AFM ones and to take the exact duality relations on the AFM results to extend them to the true states. This can be done in two steps 1. the link between the excited state of the N-body problem to the ground state of another N-body problem; 2. the link between the ground state of the N-body problem to the ground state of a twobody problem which can be solved easily. These two steps are the consequence of the duality relations given above and their proofs can be found in Ref. 99). Here we focus essentially on the practical result.
The two-body problem can be solved exactly (at least numerically) for two particles of mass m interacting through the potential V (r). Thus, given the potential V , the ground sate energy is function of the mass only. Let us define the f function by f (m) = ε(m; 0, 0). (6 . 120) This function is universal in the sense that it depends only on the form of the potential V (r) and it can be computed once for all. An example of the function f (m) for the Hamiltonian H = p p p 2 /m + r is given in Fig. 5. The conclusion of the two steps indicated previously can be summarized by the following approximation with the following expression for the mass M entering as the argument in the f function: (6 . 122) The conclusion of the relation (6 . 121) is very strong. It proves that the whole spectrum of all systems can be obtained approximatively by the calculation of a universal function, f (m), corresponding to the ground state of the 2-body system with the same potential, and for arguments M(m; {n i , l i }) given by (6 . 122). Getting f (m) is a very easy task. Solving the 2-body system can be performed with a great accuracy for any potential (let us recall that this potential must not depend on m and N); moreover, obtaining the ground state energy is free from possible numerical complications arising for excited states. Let us note also that the duality relation (6 . 121), which is obviously an approximation, concerns the exact eigenvalues only; the AFM values, which were very convenient intermediate quantities in our demonstration, have completely disappeared. The accuracy of this duality relation for exact states has been tested with success in Ref. 99) for some systems. §7.

Conclusions
In this paper, we presented a detailed review of the auxiliary field method (AFM) to solve eigenequations in the context of quantum mechanics. We showed that this method is very powerful to obtain closed analytical approximate expressions for the properties of many physically interesting systems.
The idea of the method is to replace a Hamiltonian H for which analytical solutions are not known by another oneH for which they are known. For instance, a potential V (r) not solvable is replaced by another one P(r) more familiar, or a semirelativistic kinetic part is replaced by an equivalent nonrelativistic one. The bridge between both Hamiltonians is ensured by a special function including one or more auxiliary fields. These fields can be determined by an extremization procedure and HamiltonianH reduces to Hamiltonian H. In so doing the original problem is completely equivalent to the new one. The approximation comes from the replacement of the auxiliary fields by pure real constants. The approximant solutions for H are then obtained by the solutions ofH in which the auxiliary parameters are eliminated by an extremization procedure for the eigenenergies.
The great advantage of this method is that it is able to predict the behavior of the observables in terms of the various parameters entering into the problem and also in terms of the quantum numbers. We presented the principles of the method, its characteristics, its major properties and we discussed the quality of the results for many relevant situations appearing in various domains of physics. The only restriction is the fact that, for relativistic systems or many-body problems, the computation is manageable only for systems composed of identical particles.
Although we focused essentially on the search for closed formulae of the eigenenergies, we showed that the method is also able to provide good results for the eigenstates and for different types of observables. In this context the harmonic oscillator wavefunctions and the hydrogen-like wavefunctions are particularly important. As an example, we proved that the AFM is very simple to implement for the Schrödinger equation with a linear potential and gives very satisfactory results.
The main part of the paper deals with the approximate values coming from the AFM concerning the binding energies of the Schrödinger equation (nonrelativistic) and the eigenmasses of the spinless Salpeter equation (semirelativistic). We proved a number of very important and pleasant general features of the method: • The AFM is essentially a kind of mean field approximation. In particular, using a power-law form for P(r), one can define a "mean radius" r 0 and a "mean momentum" p 0 which allows to express the final AFM eigenvalue in the very simple form T (p 0 ) + V (r 0 ), where T (|p p p|) and V (|r r r|) are respectively the kinetic and potential part of H. A very simple link exists between p 0 and r 0 . This last parameter is the solution of a transcendental equation which looks like a semiclassical version of a generalized virial theorem. • Even in complicated situations (for example N-body systems with relativistic kinematics, and presence of both one and two-body interactions), the AFM can always be brought to the resolution of a transcendental equation, a procedure quite easy from the numerical point of view. In many cases, this transcendental equation can be solved analytically and the corresponding results have been presented in this paper. • We discovered, and demonstrated, an amazing property of the AFM solutions that we called "universality of the form". This means that, provided a power-law potential is chosen for the function P(r), the AFM gives an universal expression of the resulting energies for any potential V (r) under consideration, whatever the exponent chosen for power-law potential. The only remnant of this exponent is entirely contained in the form of a principal quantum number Q(n, l). This property of universality is valid not only for the eigenvalues but also for the parameter r 0 and is maintained passing from a nonrelativistic treatment to a semirelativistic one. • The AFM is fully equivalent to the envelope theory (ET) and, as such, it can benefit of all the properties already demonstrated in the framework of this method which was introduced a long time ago. In particular, in many occasions, one can deduce that the AFM energies are upper or lower bounds on the exact energy depending on the convexity of a certain function linked to the potential. As in the ET, the energy function (potential and/or kinetics) appearing inH is, in a sense, the best possible energy function tangent to the exact one. But, in the AFM, we also pushed the theory in tracks that have not been explored by ET. • If we consider only the AFM approximation for the energy, we were able to prove several exact duality relations between different types of problems. For instance, a nonrelativistic problem with a potential V (r) can be reduced to the ultrarelativistic (m = 0) problem with a related potential W (r) = V (α √ r), and the mass of a N-body system with identical particles is related to the mass of a two-body system with a simple renormalization of the parameters. Moreover, these duality relations are shown approximately verified for exact solutions. The 2-body Schrödinger equation with a power-law potential was the prototype for the starting potential P(r). The resulting AFM energies are particularly simple. Using the "universality of the form" we proposed expressions for the principal quantum numbers Q(n, l) which allow to improve drastically the quality of the results. In some cases, a sophisticated expression for Q gives a relative accuracy of the order of 10 −4 over tens of the lowest sates of the spectrum, while even a very simple prescription is already able to give an accuracy of the order of 10 −2 . Since the AFM can be used recursively, the power-law potential is very often chosen as the basic potential P(r).
In the framework of a nonrelativistic approach, we applied the AFM to a great number of 2-body problems with potentials of various forms: sum of powers, square root, exponential potentials. In each case the AFM reproduces in a very simple way the basic properties of the solutions, even when the starting interaction P(r) is very far from the genuine potential V (r). We were able to give a closed form for the energies of two particles interacting via a funnel potential. This situation is quite realistic in hadronic systems and, to our knowledge, it is the first time that an analytical expression for the corresponding eigenenergies is given.
The AFM gives a very general formula for the critical coupling constants of nonrelativistic Hamiltonians with a finite number of bound states. The dependence on the quantum numbers, the mass m of the particles, the number N of particles, and the structure of the potential are predicted. Different N behaviours are obtained depending on the one-body or pairwise character of the interaction. If the AFM gives upper (lower) bounds for the exact eigenvalues, the critical coupling constants predicted are upper (lower) bounds for the exact critical coupling constants.
The 2-body Salpeter equation is also tested for a variety of potentials: power laws, square root, funnel. We showed that, at the limit of a large mass, the nonrelativistic expression is recovered in each case. The ultrarelativistic limit leads to particularly simple expressions. The expression for the funnel potential is really astonishing of simplicity in the case of massless particles.
Once a problem is solved within the AFM (quantities p 0 and r 0 found), it is very easy to compute the contribution of a small perturbation at the first order. It is given by the perturbation Hamiltonian evaluated at the mean momentum p 0 for a kinetic energy or at a function of the mean radius r 0 for a potential. The result does not coincide with the one obtained by the quantum perturbation theory, but the agreement can be very good.
Lastly the AFM is able to provide analytical results for the N-body problems even for quite sophisticated types of potential: Coulomb, linear, funnel interactions are presented. Although the accuracy would probably become worse and worse with increasing values of N, we showed that the spectrum of a realistic 3-body system can be reproduced with a 1% accuracy if a good definition of the principal quantum number Q is adopted. Let us mention also that the AFM N-body results were recently applied to the computation of light baryon masses for various theories of QCD with a large number of colors. 89), 90) In summary, the auxiliary field method is easy to use, can be applied to a great varieties of problems in quantum mechanics regardless the number of particles, exhibits remarkable properties and is able to provide closed analytical expressions of sometimes very impressive accuracy.

A.1. Presentation of the envelope theory
As the AFM, the envelope theory (ET) 45) is a method aiming to get approximate analytical energy formulae from an arbitrary Hamiltonian. Initially, it has been introduced to obtain bounds on the eigenenergies of the N-body problem (see Ref. 86)). Later, it has been explored and refined to simpler systems and to other purposes. We only present here its key features and refer the reader to Refs. 45), 100), 101), 102) for a detailed discussion about the basis and applications of ET.
Let us set V (r) = v f (r) in Hamiltonian (2 . 1). Then the energy spectrum of this Hamiltonian is formally given by E = F(v), where the dependence on the usual quantum numbers n and l will be dropped for simplicity. The function F(v) is concave but not necessarily monotonic. This allows to define a so-called kinetic potential k(S) by using the Legendre transformation (here the prime denotes the derivative with respect to v) This transformation can be understood as follows. |Ψ being the eigenstates of Hamiltonian (2 . 1), one can define S = Ψ | T (p p p 2 ) |Ψ and rewrite formally the energy spectrum as The transformation (A . 1) follows from these relations. One is consequently led to the exact formula What can now be done to go a step further in ET is to assume that V (r) = g(P(r)), where P(r) is a potential for which the solution of the eigenequation is analytically known. Then, can be analytically computed. It can moreover be shown that the kinetic potential corresponding to V (r), namely K(s), is given approximately by where k A (s) is the kinetic potential associated to P(r). One then obtains an approximate form for the eigenenergies that reads 45) E ≈ E = min and another approximate energy formula coming from the rewriting of (A . 6) is This last formula is called the principal envelope formula in Refs. 101), 102). It is possible to understand (A . 8) as follows. If V (r) = g(P(r)), with g(x) a smooth function of x, then we can define the "tangential potential" V t (r) at the point r = t as V t (r) = a(t) P(r) + g(P(t)) − a(t) P(t) with a(t) = V ′ (t) P ′ (t) = g ′ (P(t)). (A . 9) Such a particular form is obtained by demanding that V t (r) and its derivative agree with V (r) and V ′ (r) at the point of contact r = t. If ε ≪ 1, one has indeed V (t + ε) −V t (t + ε) = ε 2 2 P ′ (t) 2 g ′′ (P(t)) + O(ε 3 ).
(A . 10) The eigenenergies of Hamiltonian H t = T (p p p 2 ) +V t (r), denoted by E (t), are given by E (t) = ε A (a(t)) + g(P(t)) − a(t) P(t). The final energy spectrum has to be extremized with respect to v, so we have . 14) and the physical energy reads that is nothing else than the principal envelope formula (A . 8).
We have just shown that ET can lead to analytical approximate energy formulae, namely (A . 6) and (A . 8), which are both equivalent. Moreover, it has been shown in Ref. 45), as suggested by (A . 10), that E is a lower (upper) bound on the exact energy if the function g is convex (concave), that is if g ′′ > 0 (g ′′ < 0). Let us note that H and H t have the same kinetic part. The tangential potential indeed always underestimates (overestimates) the exact potential in this case. A clear interest of ET is thus that it allows to know the variational or antivariational nature of the approximation that is performed. In practice, this information can be obtained only for a nonrelativistic kinematics, since it is necessary to know exactly ε A (v) in (A . 3).

A.2. Equivalence between AFM and ET
The similarity of the starting points of ET and the AFM is obvious: In both cases, a potential for which no analytical solution is known is "approximated" by an other potential for which analytical solutions exist. It suggests that a connection between both approaches should exist; and it will indeed be established in this section. Let us apply the AFM as described above with V (r) = g(P(r)). We find the following expression for the energy (2 . 11) E(ν) = ε A (ν) + g (P(J(ν))) − ν P(J(ν)), (A . 16) the function J(x) = K −1 (x) being computed from the relation (2 . 3). Remarkably, this AFM formula is equal to the ET one (A . 13) since J(x) = P −1 (g ′−1 (x)). Consequently, the AFM and the ET lead to the same final energy formula (A . 6). The link between both approaches is given by ν = a(t). It is then easy to see thatṼ (r 0 ) = V (r 0 ) and thatṼ ′ (r 0 ) = V ′ (r 0 ). So, the potentialṼ (r) is tangent to the potential V (r). An explicit example is presented in Ref. 46). The function J(x) can be defined if the function K(x) can be inverted. In order to fulfill this condition, it is sufficient that K(x) is monotonic, that is to say that K ′ (x) has a constant sign. But, from the definitions above, we have K(x) = g ′ (P(x)), which implies that K ′ (x) = g ′′ (P(x)) P ′ (x). Since K(x) must be monotonic, the convexity of the function g is well defined if P(x) is also monotonic. This is the case if P(x) is a power-law potential, for instance. In these conditions, the convexity of the function g can also be used to determine the variational character of the AFM.
Let us summarize our results. The auxiliary field ν can be introduced as an operator in the Hamiltonian (2 . 1), and leads to an equivalent formulation of this Hamiltonian. If one considers it as a variational parameter rather than an operator, as in the AFM, the results are approximate but can be analytical. We have shown in this section that the auxiliary field, when seen as a variational parameter, is nothing else than the function a(t) generating the tangential potential in ET. This shows that, although obtained in different ways, the AFM and the ET lead to the same results. In this way, some formulae about the power-law potentials obtained in Ref. 103) by the ET were rediscovered with the AFM in Ref. 38), but supplementary results are given in this last reference.
Taking this equivalence into account, we can now better understand the meaning of the variational parameter v in the ET: Its optimal value can be seen as given by a mean field approximation since Z(v 0 ) = Z(g ′ (P(r 0 ))) (see (2 . 19)). Moreover, the properties of the AFM that have been proven in Refs. 38),40) also hold for ET. Finally, we can now have an a priori knowledge of the (anti)variational nature of the AFM energy formulae provided that we express V (r) as g(P(r)) and compute whether g is convex or concave. Equivalently, since the potentialṼ (r, ν 0 ) is tangent to the potential V (r) at r = r 0 , the approximation E(ν 0 ) is an upper (lower) bound on the exact energy ifṼ (r, ν 0 ) ≥ V (r) (Ṽ (r, ν 0 ) ≤ V (r)) for all values of r. 46) All these results holds provided that H andH(ν 0 ) have the same kinetic part. Several examples are presented above.

Appendix B Reduced equations
Finding analytical energy formulae for the potentials that we study in this work requires an analytical knowledge of the roots of particular cubic and quartic polynomial equations. In each case of interest, algebraic manipulations allow to transform the original equations to one of the following reduced equations. We sum up these equations in this appendix and put their roots in a form that is as convenient as possible to deal with.

B.1. Third order equation
We begin by the cubic equation (Y ≥ 0) So F − (Y ) is well defined for all positive values of its argument. It can be checked that the following approximate forms hold The following approximate expressions can also be useful if Y ≫ 1. (B . 10) Let us briefly recall some points concerning the Lambert function (also called Omega function or product-log), that we will denote W (z). Complements of information can be found in Ref. 104). First of all, W (z) is defined as the inverse function of z e z . Consequently, it has the following properties:

Appendix C Lambert function
the Lambert function thus exist * ) , respectively denoted as W 0 (z), defined for z ≥ −1/e, and W −1 (z), defined for −1/e ≤ z ≤ 0 They are plotted in figure 6. Obviously, both branches of the Lambert function share the same properties (C . 1) and (C . 2) as W (z); they meet in W −1 (−1/e) = W 0 (−1/e) = −1. It can moreover be checked that W 0 (|x| ≪ 1) ≈ x, lim x→0 W −1 (x) = −∞, and lim x→∞ W 0 (x) = ∞. For our purpose, it is worth mentioning that the following equation in which a, b, ρ are real numbers such that a z + b is positive is analytically solvable with the Lambert function. One has where either W 0 or W −1 has to be chosen following the range of θ and z that the function z(θ ) has to cover.

Appendix D Extended virial theorem
In this section and in the following ones, we will denote O = n, l|O|n, l and n|O|n = n, 0|O|n, 0 , |n, l being an eigenstate of H with a radial quantum number n and an orbital quantum number l. If A is an arbitrary operator, it follows that This recurrence relation is particularly useful to compute r k mean values for potentials with integer power, choosing s as an integer. and the normalized eigenvectors by ψ nlµ (r r r) = (2γ nl ) 3/2 n! 2(n + l + 1)(n + 2l + 1)! (2γ nl r) l e −γ nl r L 2l+1 n (2γ nl r)Y l µ (r r r), (F . 3) with γ nl = η/(n + l + 1) and η = mν. At the origin, the S-states are such that |ψ n,0 (0)| 2 = η 3 π(n + 1) 3 .

Appendix H Overlap with dilated functions
The scalar product of two radial functions R n,l (r) and R n ′ ,l (r) of a set of orthonormal states is simply given by δ nn ′ . When one of these functions is scaled by a positive factor a, the overlap The first relation stems from the definition (H . 1), the second one from the Schwarz inequality, and the others are due to scaling properties. Using the dilation properties of the Laguerre polynomials and the various existing recurrence relations, 106) it is possible to compute analytically the formula F n,n ′ ,l (a) for hydrogenlike systems and harmonic oscillators.
In the first case (hydrogen-like systems), one obtains Denoting p p p i and π π π i the conjugate momenta associated respectively to r r r i and x x x i , it is easy to prove that π π π i = N ∑ j=1 B ji p p p j and p p p i = N ∑ j=1 U ji π π π j . (I . 4) π π π N = P P P = p p p 1 + p p p 2 + . . . + p p p N is the total momentum of the system.

Appendix J The N-body harmonic oscillator
The case of a quadratic potential (the N-body nonrelativistic harmonic oscillator) is practically the only one for which an exact solution is reachable in three dimensions, at least formally. This section is devoted to this problem, since it is the basic ingredient of the N-body AFM.

J.1. Nonrelativistic general case
Let us start with the most general harmonic-oscillator-like Hamiltonian, corresponding to N particles of arbitrary masses, a one-body and a two-body quadratic potentials with arbitrary spring constants, so that the Hamiltonian looks like The kinetic energy operator, expressed in terms of Jacobi variables, allows the correct separation of the center of mass motion and appears decoupled in the various variables (this is in fact the justification of the form of Jacobi coordinates) where m t = αm = m 1 + m 2 + · · · m N is the total mass of the system. The kinematical quantities λ i are calculated as For further convenience, it is judicious to switch from the standard Jacobi variables to renormalized conjugate Jacobi coordinates defined as y y y i = x x x i λ i ; ρ ρ ρ i = λ i π π π i . (J . 4) Working in the center of mass frame (P P P = 0 0 0) and using these new variables, the kinetic energy operator has a very simple form With these new variables, the one-body operator is written F lm y y y l · y y y m , (J . 6) where the symmetric definite positive matrix F is defined by the matrix B being given by (I . 3). In the very same way, the two-body operator is written G lm y y y l · y y y m (J . 8) where the symmetric definite positive matrix G is defined by Introducing the matrix J = F + G, the total potential V = V 1 +V 2 is expressed as and consequently the energy of the system is given by ω i (2n i + l i + 3/2), (J . 14) where n i and l i are respectively the radial and orbital quantum numbers associated to the coordinate z z z i . The problem is now completely solved. Moreover, this result is valid for any excited state. Even if the expression of ω i is in general not analytical, the result (J . 14) is exact and could be calculated with a high accuracy. Notice that formula (J . 14) extends a previous result, 107), 108) where an equivalent mass formula is obtained in the case a i = 0. An explicit result for the case N = 3 is given in Ref. 68).

J.2. Case of identical particles
It is of interest to rewrite the solution in the case where all the particles are identical, implying that they have the same mass m i = m and the same spring constants a i = a, b i j = b. Hamiltonian (J . 1) indeed reads in this case (r r r i − r r r j ) 2 , (J . 15) and it is easy to see that the J matrix is already of diagonal from the very beginning, so that its eigenvalues are analytically known (they all read a + Nb). Consequently, the eigenenergies of the system are also analytical. Explicitly they are given by where Q is the total principal number Wavefunction with good global quantum numbers and various symmetries can be built by appropriate linear combinations of states characterized by the same value of the principal quantum number Q. 91) All symmetrized states have then the same energy than the nonsymmetrized states. For spatial wave functions completely symmetrical, the ground state is obtained for the values n i = l i = 0 ∀i so that the principal quantum number is simply For mixed symmetry or completely antisymmetrical spatial wavefunctions, the situation is much more involved. An estimation of the ground-state energy for the completely antisymmetrical case can be computed by choosing different values for the quantum numbers and piling the states (d identical values of the same n j and l j per state in order to take care of a possible degeneracy due to internal degrees of freedom) up to the Fermi level. By considering only particle number insuring a saturated Fermi level (closed shell), one obtains