# Lagrange-mesh calculations in momentum space

###### Abstract

The Lagrange-mesh method is a powerful method to solve eigenequations written in configuration space. It is very easy to implement and very accurate. Using a Gauss quadrature rule, the method requires only the evaluation of the potential at some mesh points. The eigenfunctions are expanded in terms of regularized Lagrange functions which vanish at all mesh points except one. It is shown that this method can be adapted to solve eigenequations written in momentum space, keeping the convenience and the accuracy of the original technique. In particular, the kinetic operator is a diagonal matrix. Observables in both configuration space and momentum space can also be easily computed with a good accuracy using only eigenfunctions computed in the momentum space. The method is tested with Gaussian and Yukawa potentials, requiring respectively a small or a great mesh to reach convergence.

###### pacs:

02.70.-c,03.65.Ge,03.65.Pm,02.30.Mv^{†}

^{†}thanks: F.R.S.-FNRS Research Fellow

^{†}

^{†}thanks: F.R.S.-FNRS Senior Research Associate

## I Introduction

There are few three-dimensional problems in quantum mechanics which allow a complete analytical solution for any value of the orbital angular momentum (the -wave channel is very similar to a simpler one-dimensional equation) flug99 (). So, numerous methods have been developed to solve numerically with a high accuracy the eigenvalue equations associated with various systems. Among these techniques, the Lagrange-mesh method (LMM), which is especially easy to implement, can produce very accurate results. First created to compute eigenvalues and eigenfunctions of a two-body Schrödinger equation BH-86 (); VMB93 (); Ba-95 (); hess02 (); Baye06 (); baye08 (), it has been extended to treat semirelativistic Hamiltonian sema01 (); brau02 (); buis05 (); lag1 (). The trial eigenstates are developed in a basis of particular functions, the Lagrange functions, which vanish at all mesh points, except one. Once the potential is known in the configuration space, its matrix elements are simply the potential values at the mesh points, if they are computed with an associated Gauss quadrature. At first sight, this method could look like a discrete variational method, but this is absolutely not the case since the eigenfunctions can be computed at any position. Because of the use of the Gauss quadrature scheme, the method is not variational but the results are too a large extent independent of the sole nonlinear parameter of the method fixing the physical scale of the system. Generally, a great accuracy can be reached with a small mesh baye02 ().

At the beginning, the LMM was developed in the configuration space. Recently, it has been shown that the Fourier transform of the eigenfunctions computed in the configuration space can easily be obtained with a good accuracy in the physical domain of the momentum space lacr11 (). Moreover, observables in this space can easily be computed with a good accuracy using only matrix elements and eigenfunctions in the configuration space. But, for some particular problems, it can be preferable to work in the momentum space. This is the case when the potential presents discontinuities in the configuration space karr10 () or when the potential is given in the momentum space. In this last case, if it is possible to use the LMM by computing first the Fourier transform of the potential, we will show here that the LMM can be adapted to solve the eigenequations directly in momentum space. Observables in both configuration and momentum spaces can also be computed with a good accuracy using only eigenfunctions computed in momentum space. Moreover, we will show that the new LMM can provide these types of data very efficiently and very easily, using again the fundamental properties of the Lagrange functions. Let us note that the method presented here relies on a mesh of points built with the zeros of a Laguerre polynomial, but a general procedure for deriving other Lagrange meshes related to orthogonal or non-orthogonal bases has also been developed baye99 ().

Several methods exist to solve eigenequations in momentum space. For instance, iterative procedures have been developed regi84 (); rodr88 (). They are quite accurate but resort finally to numerical integrations on a mesh. Direct computations on a mesh are easier to carry out, but they require very large mesh if a good quadrature rule is not used karr10 (). As we will see, the LMM in momentum space is very easy to implement and can also give accurate results. In order to fix the notations, the eigenequations in momentum space are presented in Sec. II. The LMM adapted in momentum space is described in Sec. III with some details in order that the paper is self-contained. Test calculations are presented in Sec. IV for two potentials with very different properties of convergence. Some concluding remarks are given in Sec. V.

## Ii Eigenequations in position and momentum spaces

Let us consider the following eigenequation ()

(1) |

in which the kinetic part and the potential depend respectively on the relative square momentum and on the radial distance between the particles. The wavefunctions in configuration space and in momentum space can be written using the spherical representation

(2) | ||||

(3) |

where is a modified spherical harmonic var (), and is the number of nodes at finite distance. These functions are linked by the following Fourier transforms saku93 ()

(4) | ||||

(5) |

These equations lead to lacr11 ()

(6) | ||||

(7) |

where is a spherical Bessel function abra65 ().

Written in the momentum space, (1) takes the following form

(8) |

with , the Fourier transform of , given by

(9) |

This potential is a continuous function of the momentum, even if parts of the interaction in configuration space present discontinuities. One can think of square well or Dirac delta function (repulsive only in a 3D-space). As the potential depends only on , we have and (9) becomes grad80 ()

(10) |

Using the standard decomposition of a radial function var (), the eigenvalue equation (8) takes the form of an integral equation for the wavefunction

(11) |

with the partial potentials

(12) |

The Legendre polynomial depends on the variable .

For a Schrödinger equation, the kinetic operator is given by

(13) |

and the eigenvalue is the binding energy of the system. In a spinless Salpeter Hamiltonian, the kinetic operator takes the following form

(14) |

and the eigenvalue is the mass of the system. This kind of Hamiltonian is sometimes denoted semirelativistic since it is not a covariant formulation. This equation can be considered as a Schrödinger equation with its nonrelativistic kinetic part replaced by a relativistic counterpart. More rigorously, it is obtained from the covariant Bethe-Salpeter equation salp51 () with the following approximations: Elimination of any dependences on timelike variables and neglect of particle spin degrees of freedom as well as negative energy solutions grei94 (). The spinless Salpeter Hamiltonian is often used in hadronic physics to study bound states of quarks or gluons godf85 (); fulc94 (); buis07 (); math08 ().

Within this formulation, the action of the kinetic operator is just an ordinary multiplication. So, nonrelativistic and semirelativistic systems are computed with the same manner. Moreover, more complicated kinetic parts, with momentum-dependent masses szcz96 (); llan00 (); agui11 (), can be equally treated. Though the formulations in the configuration and momentum spaces are completely equivalent, this does not mean that the technical difficulties to solve the eigenequations are the same in both spaces.

If the potential is known in the configuration space, (10) and (12) allows the computations of the partial potentials . can also be directly obtained from a physical theory, like a field theory naturally written in the momentum space. One can think of effective potentials obtained from Feynman diagrams, for instance. It is not always possible to obtain an analytical form for (12), but such a numerical integration can be rapidly and accurately performed.

## Iii Method in momentum space

### iii.1 Lagrange functions

The LMM relies on the existence of a -point mesh , which is associated with an orthonormal set of indefinitely derivable functions , called the Lagrange function BH-86 (); VMB93 (); Ba-95 (). Each function satisfies the Lagrange conditions,

(15) |

that is to say it vanishes at all mesh points except one. The and are respectively the abscissae and the weights of a Gauss quadrature formula

(16) |

Several quadratures are possible, but, as we work with wavefunctions depending only on the module of a variable, we consider the case of the Gauss-Laguerre quadrature whose domain of interest is . The Gauss formula (16) is exact when is a polynomial of degree at most, multiplied by . The mesh points are the zeros of a Laguerre polynomial of degree : BH-86 (). These zeros can be determined with a high precision with usual methods to find the roots of a polynomial numrec () (the Mathematica expression Root does efficiently the job) or as the eigenvalues of a particular tridiagonal matrix golu69 (). The weights can be computed by the following formula baye02 ()

(17) |

### iii.2 Matrix equation

The use of the LMM in configuration space is described in BH-86 (); VMB93 (); Ba-95 (); Baye06 (); baye08 (); hess02 (); sema01 (); lacr11 (). We will present here the formulation for the integral equation (11) within the LMM. The idea is to expand the wavefunction with the regularized Lagrange functions in such a way that a trial state is written

(19) |

This formula is identical to formula (6) in lacr11 () with just the replacement of the variable by . The coefficients are linear variational parameters and the scale factor is a nonlinear parameter aimed at adjusting the mesh to the domain of physical interest. The parameter has here the dimension of a momentum ( plays the same role in lacr11 () but has the dimension of a distance). Contrary to some other mesh methods, the wavefunction is also defined between mesh points by (18) and (19). For a good value of and a sufficiently high value of (see sec. IV), the function

(20) |

can be a good approximation of the exact function . Note that an advantage of the LMM is that the value for has not to be known with a great accuracy. It is sufficient that this value be located within a given interval, as we will check in the following.

Basis states built with the regularized Lagrange functions are not exactly orthogonal. But, at the Gauss approximation, we have . So, in the following, all mean values will be performed using the Gauss quadrature formula (16). Inserting expansion (20) in (11) gives

(21) |

where is a dimensionless variable. We can now multiply this equation by and integrate on with again the Gauss quadrature formula (16). Finally, we obtain

(22) |

The Hamiltonian matrix is symmetric since . A similar expression is obtained for calculations with the LMM in the configuration space for a nonlocal potential hess02 (). With the LMM in momentum space, the solution of a quantum equation reduces (as it is often the case) to the determination of eigensolutions of a given matrix. So, once the partial potentials are known, (22) shows that this method is very easy to implement. The computations of the Fourier transform of the wavefunction and of observables is no more complicated, as presented in the following.

In Sec. IV, we will study the convergence of the method as a function of the scale parameter and the number of mesh points for nonrelativistic an semirelativistic kinematics. Let us note that an automatic determination of has been developed for the LLM in the configuration space brau98 (); lacr11 (). But, the generalization of such a technique to the LLM in the momentum space is very difficult due to the nonlocal nature of the interaction in (11). We will cross-check our results by comparing eigenvalues and mean values of observables computed with the LMM in both configuration and momentum spaces. Moreover, a wavefunction computed in the momentum space will be compared with the Fourier transform of the corresponding wavefunction computed in the configuration space, and vice versa.

### iii.3 Fourier transform

It has been shown in lacr11 () that a good approximation of a function can be obtained from the Fourier transform of the corresponding solution computed in the configuration space by the LLM. This can be performed by taking benefit of the very special properties of the regularized Lagrange functions. Similarly, a good approximation of the function can be obtained from the Fourier transform of a solution computed in the momentum space by the LLM. Using the Gauss quadrature rule (16) with the Lagrange condition (15), the integral (7) for a given trial state (20) simply reduces to

(23) |

This formula is identical to formula (22) in lacr11 () with just the replacement of the variable by . This results from the similarity of (4) and (5), and the choice of expansion (19). For a sufficiently high value of (see sec. IV), can be a very good approximation of the genuine function in the configuration space. Above a critical value of , generally in the asymptotic tail of , this function can present large unphysical rapid oscillations. These oscillations do not develop in , because they are damped by the rapid decrease of the regularized Lagrange functions.

### iii.4 Mean values of momentum-dependent observables

The mean value of the operator for a trial state is given by

(24) |

Using the Lagrange condition (15) and the Gauss quadrature (16), this integral reduces to

(25) |

If is the identity, we recover the normalization condition as expected. As we will see in sec. IV, a very good accuracy can be reached for the mean values .

### iii.5 Mean values of radial observables

The mean value of the operator for a trial state is given by

(26) |

The method to compute matrix elements relies on the fact that in the momentum space luch90 (). Let us first look at the matrix whose elements are . With (16), these matrix elements are given by

(27) |

where

(28) |

This expression is exact for some Lagrange meshes, but this is not the case for the regularized Laguerre mesh. An exact expression can easily be obtained (see the Appendix in VMB93 ()). However, as shown in Ba-95 (), it is preferable to use the approximation (27)-(28). The quantities are then easy to obtain and read Ba-95 ()

(29) |

If is the diagonal matrix formed by the eigenvalues of , we have

(30) |

where is the transformation matrix composed of the normalized eigenvectors. Let us call the diagonal matrix obtained by taking the function of all diagonal elements of (remember that is linked to the matrix elements of , not ). The numbers are well approximated by the elements of the matrix obtained by using the transformation (30): . As we will see in sec. IV, a very good accuracy can be reached for the mean values .

## Iv Numerical results

### iv.1 Gaussian potential

We first consider a Gaussian potential whose Fourier transform is well known grad80 ()

(31) |

Let us note that the asymptotic behavior of the potential is similar in the configuration and momentum spaces. Using (12), we can compute

(32) |

where is an exponential integral abra65 () and means the integer part of . This type of potentials is the prototype for a short-range interaction and can be used in many field of physics.

With a nonrelativistic kinematics, we can write

(33) |

where is the solution of the dimensionless Hamiltonian

(34) |

We choose to study the performance of the LMM in momentum space with the value . In this case, two bound states exist: and .

Let us look at the results obtained for (34) by solving (22). The variation of the eigenvalue as a function of the nonlinear parameter is presented for the two states in Fig. 1. For a value of as small as 20, a plateau is present with abrupt variations of the eigenvalue at the borders. The length of the plateau increases with the number of mesh points. For the excited states, the convergence is usually slower than for the ground states: The plateau is less extended and the variations around the plateau are larger. Even with , very good results can be obtained provided is taken in the small corresponding quasi-plateau (the fluctuation in the rapid variation). Note that the flatness of the plateau is a criterion obviously depending on the accuracy one wants to reach. It is shown on Fig. 2 that eigenvalues, corresponding to different values of taken in plateaus, converge towards the same value as increases. This convergence can be achieved from above or from below. All these results clearly show the non variational nature of the LMM.

It is remarkable that, even for small values of , the value of must not be determined with a high precision. This is a great advantage for the LMM. Nevertheless, if is too small, a significant part of the wavefunction is not covered by the points of the mesh. When is too large, all points of the mesh are located in the asymptotic tail of the wavefunctions and it is then impossible to obtain good eigenvalues. Let us note that the same behaviors are observed for the LMM in configuration space BH-86 (); Ba-95 (); sema01 ().

Eigenvalues and some observables have been calculated with the LMM in both configuration and momentum spaces. Results are presented and compared in Table 1 for the ground state only, but results are similar for the excited state. A very good accuracy can be obtained for the computation in momentum space, even for small values of . For most of the physical problems, the required precision can probably be reached with . It seems that position-dependent observables converge more slowly than momentum-dependent observables. We have checked that the situation is the opposite for LMM calculations in configuration space.

In principle, we must have . This is not the case since momentum-dependent and radial observables are not computed with the same method (see Sec. III.4 and III.5). We can see in Table 1 that the difference between the two quantities is around the accuracy of the eigenvalues .

Conf. | Mom. | |||
---|---|---|---|---|

The wavefunctions produced by the LMM in both spaces have also been compared. A typical example is shown on Fig. 3. On the left part, the wavefunction obtained directly in the momentum space and the wavefunction computed by the Fourier transform of the wavefunction obtained directly in the configuration space are superposed. The situation is the opposite on the right part. In both cases, the agreement is very good for low values of the arguments, but the Fourier transform wavefunctions present large unphysical oscillation for larger values of the argument. The starting points of these oscillations can be rejected to high values of the argument in the asymptotic tail, but it is then necessary to work with greater values of , typically lacr11 (). A good representation of the wavefunction in momentum space can be obtained with a small mesh only by working directly in this space. This an advantage of the LMM in momentum space.

As in the nonrelativistic case, we will not try to study a “realistic” relativistic system. We will consider quite arbitrary, but convenient, values for the parameters,

(35) |

for which only one bound state exists (). As we can see by examining Fig. 4, the behaviors of the eigenvalues as a function of and is similar for nonrelativistic and semirelativistic Hamiltonians. But, the convergence is less good for semirelativistic systems, with shorter plateaus and larger variations around the plateaus. This situation is similar in configuration space sema01 (); lacr11 (). Various observables have been computed and are presented in Table 2. A reasonable accuracy can be obtained even with a small mesh.

Conf. | Mom. | |||
---|---|---|---|---|

### iv.2 Yukawa potential

Some numerical tests are also performed using the Yukawa potential whose asymptotic behavior is very different in the configuration and momentum spaces grad80 ()

(36) |

Using (12), we can compute

(37) |

where is a Legendre function of the second kind magn66 (). This type of potentials is used in many field of physics, even in hadronic physics. If the interaction in a quark-antiquark system or between two gluons in the vacuum can be simulated by a funnel potential (Coulomb+linear), the confinement vanishes inside a quark-gluon plasma above a critical temperature and the interaction could turn into a Yukawa potential brau07 ().

With a nonrelativistic kinematics, we can write

(38) |

where is the solution of the dimensionless Hamiltonian

(39) |

We choose to study the performance of the LMM in momentum space with the value . In this case, three bound states exist. A phenomenological approximate formula gree82 () gives , and .

On Fig. 5, we can see the behaviors of the ground state eigenvalue of (39), obtained by solving (22), as a function of and . The situation is at first sight similar to the case of the Gaussian potential, but the convergence is much slower for the Yukawa potential. For a mesh of 50 points, a plateau does not appear clearly, just a slowing down of the variation of . For computations in the configuration space, we have checked that a small plateau is already present for . Consequently, eigenvalues computed in momentum space for two different values of will only coincide for large values of . About ten times more points are necessary for computations in momentum space to reach an accuracy similar to the one of computations in configuration space. It is illustrated in Table 3 where eigenvalues and some observables are calculated with the LMM in both configuration and momentum spaces. Again, the mean value , which must be equal to , is computed. The disagreement for the results in momentum space shows that the convergence is not so good than for the results in configuration space. We think that this peculiarity is linked to the asymptotic behavior of the Yukawa potential in the momentum space: Its decrease like a power-law (36) is quite slow. The wavefunction in this space is then very extended and requires its computation at high values of to be well described on its physical domain.

Conf. | Mom. | Conf. | Mom. | Conf. | Mom. | |

() | () | () | () | () | () | |

16.340426 | 16.340415 | 0.6053933 | 0.6053975 | 0.205082327 | 0.205082331 | |

23.788977 | 23.788942 | 2.95238 | 2.95241 | 2.70792857 | 2.70792862 | |

40.1294 | 40.1200 | 3.55778 | 3.55743 | 2.913010896 | 2.913010877 | |

16.340426 | 16.331047 | 0.6053933 | 0.6050217 | 0.205082327 | 0.205082257 |

The wavefunctions produced by the LMM in both spaces have also been compared for the Yukawa potential on Fig. 6. With a small mesh of 20 points, the relative error on the ground state eigenvalues is around for the computation in configuration space, and around for the computation in momentum space. Curiously, the wavefunction computed in momentum space is already very well reproduced and is even better than the Fourier transform computed from the computation in configuration space with the same mesh. By increasing , the beginning of the unphysical oscillations in the Fourier transforms can be rejected far in the asymptotic tail.

We will briefly comment the Yukawa interaction with a semirelativistic kinematics. As in the nonrelativistic case, we did not try to study a “realistic” system. Quite arbitrary but convenient values are considered for the parameters,

(40) |

for which only one bound state exists (). Results about the convergence are presented in Fig. 7. Graphics are similar to the ones produced in the nonrelativistic case: Existence of plateaus for energy as a function of ; increase of the length of the plateau with . But, the convergence is slower than in the nonrelativistic computations: must be greater to achieve a quasi flat plateau (Fig. 7, left) and to reach a convergence between eigenvalues computed with different values of (Fig. 7, right). We have checked that it is also the case to obtain an agreement between observables computed in both spaces and to obtain better Fourier transforms. With a semirelativistic Hamiltonian, the kinetic energy increases as not as . So, higher values of the momentum can be reached by the wavefunction, and it can be more extended than with a nonrelativistic kinematics. This amplifies the problem mentioned in the previous section about the asymptotic behavior of the potential.

## V Concluding remarks

The Lagrange-mesh method is a procedure to compute accurate eigenvalues and eigenfunctions of quantum equations. Implemented at the origin in the configuration space, this technique requires only the computation of the potential at some mesh points BH-86 (); sema01 (). This is due to the use of a Gauss quadrature rule with the fact that the basis functions satisfy the Lagrange conditions, that is to say they vanish at all mesh points except one.

Using this very special property, we have shown that the Lagrange-mesh method can be adapted to be used directly in momentum space with a nonlocal potential (but local in the configuration space). In this case, nonrelativistic and semirelativistic kinetic parts are treated with the same manner. Moreover, it is possible to treat interactions which present discontinuities in configuration space. The wavefunction in momentum space is directly obtained under the form of expansion coefficients for the Lagrange functions. Using only these coefficients, it is easy to compute the wavefunction in configuration space by Fourier transform, as well as mean values of -dependent or -dependent operators. Though the technique is not variational, the convergence can be checked: Eigenvalues present a plateau as a function of the sole nonlinear parameter of the method; eigenvalues computed with different nonlinear parameters in a plateau tends towards the same value as the number of mesh points increases.

The method has been tested with two potentials. With a Gaussian interaction and a nonrelativistic kinematics, a good accuracy can be obtained for eigenvalues and observables with a very small mesh, as for the method in configuration space: 20 points are probably sufficient for most of the physical applications. For a Yukawa interaction, it is necessary to use more points than for the Lagrange-mesh method in the configuration space, typically 10 times more. For both potentials: i) The convergence is slowed down for a semirelativistic kinematics, as for the method in configuration space. But a good accuracy can nevertheless be reached sema01 (); ii) It is necessary to use a large mesh to obtain a correct Fourier transform of the wavefunctions, as for the computations in configuration space lacr11 (). iii) A good wavefunction in momentum space can obtained with a small mesh by a direct computation in momentum space. It does not present the unphysical oscillations of the Fourier transform of the wavefunctions in configuration space, which can only be eliminated by using a large mesh.

The purpose of this paper is to test the feasibility and the reliability of the Lagrange-mesh method in momentum space. We think that it can be safely applied to physical problems, since convergence tests exists and a good accuracy can always be obtained. Even if several hundreds of points are necessary to treat some potentials, eigenvalues can be rapidly computed. The method is very easy to implement once the interaction is known in momentum space, and it can be particularly useful if the kinetic part is an unusual function of the momentum, since its corresponding matrix is diagonal and easy to compute. This kind of situation appears in hadronic physics where quarks or gluons can be considered with a momentum dependent mass which can be very complicated to define in the configuration space szcz96 (); llan00 (); agui11 (), e.g. . For these systems, the dominant interaction can be simulated by a linear potential in the configuration space. This type of potential is a priori highly singular in the momentum space, but it exists in a distributional sense and can be used in this context hers93 ().

## Acknowledgments

G.L. and C.S. would thank the F.R.S.-FNRS for the financial support.

## References

- (1) S. Flügge, Practical Quantum Mechanics (Springer, New York, 1994).
- (2) D. Baye and P.-H. Heenen, J. Phys. A 19, 2041 (1986).
- (3) M. Vincke, L. Malegat, and D. Baye, J. Phys. B 26, 811 (1993).
- (4) D. Baye, J. Phys. B 28, 4399 (1995).
- (5) M. Hesse, J. Roland, and D. Baye, Nucl. Phys. A 709, 184 (2002).
- (6) D. Baye, Phys. Stat. Sol. (b) 243, 1095 (2006).
- (7) D. Baye and K. D. Sen, Phys. Rev. E 78, 026701 (2008).
- (8) C. Semay, D. Baye, M. Hesse, and B. Silvestre-Brac, Phys. Rev. E 64, 016703 (2001).
- (9) F. Brau and C. Semay, J. Phys. G: Nucl. Part. Phys. 28, 2771 (2002).
- (10) F. Buisseret and C. Semay, Phys. Rev. E 71, 026705 (2005).
- (11) F. Buisseret and C. Semay, Phys. Rev. E 75, 026705 (2007).
- (12) D. Baye, M. Hesse, and M. Vincke, Phys. Rev. E 65, 026701 (2002).
- (13) G. Lacroix and C. Semay, Phys. Rev. E 84, 036705 (2011).
- (14) W. A. Karr, C. R. Jamell, and Y. N. Joglekar, Am. J. Phys. 78, 407 (2010).
- (15) D. Baye and M. Vincke, Phys. Rev. E 59, 7195 (1999).
- (16) P. E. Regier and A. J. Thakkar, Phys. Rev. A 30, 30 (1984).
- (17) W. Rodriguez and Y. Ishikawa, Chem. Phys. Lett. 146, 515 (1988).
- (18) D. A. Varshalovich, A. N. Moskalev, and V. K. Khersonskii, Quantum Theory of Angular Momentum (World Scientific, Singapore, 1988).
- (19) J. J. Sakurai, Modern Quantum Mechanics (Addison Wesley, 1993).
- (20) M. Abramowitz and I. A. Stegun, Handbook of Mathematical Functions (Dover, New York, 1965).
- (21) I. S. Gradshteyn and I. M. Ryzhik, Tables of Integrals, Series, and Products (Academic Press, New York, 2007).
- (22) E. E. Salpeter and H. A. Bethe, Phys. Rev. 84, 1232 (1951).
- (23) W. Greiner and J. Reinhardt, Quantum Electrodynamics (Springer, Berlin, 1994).
- (24) S. Godfrey and N. Isgur, Phys. Rev. D 32, 189 (1985).
- (25) L. P. Fulcher, Phys. Rev. D 50, 447 (1994).
- (26) F. Buisseret, C. Semay, V. Mathieu, and B. Silvestre-Brac, Eur. Phys. J. A 32, 123 (2007).
- (27) V. Mathieu, F. Buisseret, and C. Semay, Phys. Rev. D 77, 114022 (2008).
- (28) A. Szczepaniak, E. S. Swanson, C.-R. Ji, and S. R. Cotanch, Phys. Rev. Lett. 76, 2011 (1996).
- (29) F. J. Llanes-Estrada and S. R. Cotanch, Phys. Rev. Lett. 84, 1102 (2000).
- (30) A. C. Aguilar, D. Binosi, and J. Papavassiliou, Phys. Rev. D 84, 085026 (2011).
- (31) W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, Numerical Recipes (Cambridge University Press, 2007).
- (32) G. H. Golub and J. H. Welsch, Math. Comput. 23, 221 (1969).
- (33) F. Brau and C. Semay, J. Comput. Phys. 139, 127 (1998).
- (34) W. Lucha, Mod. Phys. Lett. A 5, 2473 (1990).
- (35) W. Magnus, F. Oberhettinger, and R. P. Soni, Formulas and Theorems for the Special Functions of mathematical Physics (Springer, New York, 1966).
- (36) F. Brau and F. Buisseret, Phys. Rev. C 76, 065212 (2007).
- (37) A. E. S. Green, Phys. Rev. A 26, 1759 (1982).
- (38) H. Hersbach, Phys. Rev. D 47, 3027 (1993).