Solutions of the Bethe-Salpeter equation in Minkowski space and applications to electromagnetic form factors”Relativistic Description of Two- and Three-Body Systems in Nuclear Physics”, ECT*, October 19-23, 2009

# Solutions of the Bethe-Salpeter equation in Minkowski space and applications to electromagnetic form factors††thanks: ”Relativistic Description of Two- and Three-Body Systems in Nuclear Physics”, ECT*, October 19-23, 2009

J. Carbonell and V.A. Karmanov
###### Abstract

We present a new method for solving the two-body Bethe-Salpeter equation in Minkowski space. It is based on the Nakanishi integral representation of the Bethe-Salpeter amplitude and on subsequent projection of the equation on the light-front plane. The method is valid for any kernel given by the irreducible Feynman graphs and for systems of spinless particles or fermions. The Bethe-Salpeter amplitudes in Minkowski space are obtained. The electromagnetic form factors are computed and compared to the Euclidean results.

## 1 Introduction

Bethe-Salpeter (BS) equation for a relativistic bound system was initially formulated in the Minkowski space [1]. It determines the binding energy and the BS amplitude. However, in practice, finding the solution in Minkowski space is made difficult due its singular behaviour. The singularities are integrable, but the standard approaches for solving integral equation fail.

To overcome this difficulty, Wick [2] formulated the BS equation in the Euclidean space, by rotating the relative energy in the complex plane . This ”Wick rotation” led to a well defined integral equation which can be solved by standard methods. Most of practical applications of the BS equation have been achieved using this technique [3, 4] and recent developments make its solution a trivial numerical task [5]. Another method – the variational approach in the configuration Euclidean space – was recently developed in [6].

The binding energy provided by the Euclidean BS equation is the same than the Minkowski one. However, the Euclidean BS amplitude does not allow to calculate some observables, e.g. electromagnetic form factors. The integral providing the form factors contains singularities which are different from those appearing in the BS equation and whose positions depend on the momentum transfer. Their existence invalidates the Wick rotation in the form factor integral. In terms of the Euclidean amplitude, the form factor can be obtained only approximately, in the so called static approximation. The error rapidly increases with the momentum transfer. To avoid this problem, the knowledge of the BS amplitude in Minkowski space is mandatory.

Thus, fifty years after its formulation, finding the BS solutions in the Minkowski space is still a field of active research.

Some attempts have been recently made to obtain the Minkowski BS amplitudes. The approach proposed in [7, 8] is based on the Nakanishi integral representation [9, 10] of the amplitudes and solutions have been found for the ladder scalar case [7, 11] as well as, under some simplifying ansatz, for the fermionic one [12]. Another approach [13] relies on a separable approximation of the kernel which leads to analytic solutions. Recent applications to the system can be found in [14].

In the paper [15] we have proposed a new method to find the BS amplitude in Minkowski space and applied it to the system of two spinless particles. In [16] this approach was generalized to the two-fermion case. Like in the papers [7, 11, 12], our approach is based on Nakanishi representation of the BS amplitude. The main difference between our approach and those followed in [7, 11, 12] is the fact that, in addition to the Nakanishi representation, we use the light-front projection. This eliminates the singularities related to the BS Minkowski amplitudes. Our method is valid for any kernel given by the irreducible Feynman graphs.

In this paper we give brief review of the approach [15] to find the BS amplitude in Minkowski space, of its results and applications.

The plan of the paper is the following. In order to present the method more distinctly, we consider first the case of zero total angular momentum and spinless particles. In sect. 2, corresponding equation for the Nakanishi weight function is given. In sect. 3, it is applied to the ladder kernel and in sect. 4 – to the cross ladder one. In sect. 5 the method is generalized to the two-fermion system. Sect. 6 is devoted to application to the electromagnetic form factor, where advantage of the Minkowski space solution manifests itself in full measure. Sect. 7 contains concluding remarks.

## 2 Equation for the weight function

For a bound state of total momentum and in case of equal mass particles, the BS equation reads

 \mathchar264\relax(k,p)=i2[(p2+k)2−m2+iϵ][(p2−k)2−m2+iϵ]∫d4k′(2π)4iK(k,k′,p)\mathchar264\relax(k′,p), (1)

where is the BS amplitude, the interaction kernel, the mass of the constituents and their relative momentum. We will denote by the total mass of the bound state, and by its binding energy.

Our approach consists of two steps. In the first one, the BS amplitude is expressed via the Nakanishi integral representation [9, 10]:

 \mathchar264\relax(k,p)=∫1−1dz′∫∞0%dγ′g(γ′,z′)[k2+p⋅kz′+14M2−m2−γ′+iϵ]3. (2)

Notice that in this representation, the dependence on the two scalar arguments and of the BS amplitude is made explicit by the integrand denominator and that the weight Nakanishi function is non-singular. By inserting the amplitude (2) into the BS equation one finds an integral equation, still singular, for .

In the second step, we apply to both sides of BS equation an integral transform – light-front projection [15] – which eliminates singularities of the BS amplitude. It consists in replacing where is a light-cone four-vector , and integrating over in infinite limits. We obtain in this way, a non-singular equation for the non-singular . After solving it and substituting the solution in eq. (2), the BS amplitude in Minkowski space can be easily computed.

This leads (see ref. [15] for the detail of calculations) to the following equation for the weight function :

 ∫∞0g(γ′,z)dγ′[γ′+γ+z2m2+(1−z2)κ2]2=∫∞0dγ′∫1−1dz′V(γ,z;γ′,z′)g(γ′,z′), (3)

This is just the eigenvalue equation of our method. It is equivalent to the initial BS equation (2). The total mass of the system appears on both sides of equation (2) and is contained in the parameter As calculations [7, 8] show, may be zero in an interval . The exact value where it differs from zero is determined by the equation (2) itself.

The kernel , appearing in the right-hand side of eq. (2), is related to the kernel from the BS equation by

 V(γ,z;γ′,z′)=ω⋅pπ∫∞−∞−iI(k+βω,p)dβ[(p2+k+βω)2−m2+iϵ][(p2−k−βω)2−m2+iϵ], (4)

with

 I(k,p)=∫d4k′(2π)4iK(k,k′,p)[k′2+p⋅k′z′−γ′−κ2+iϵ]3. (5)

For simple kernels given by a Feynman graph the integral (2) over is calculated analytically. The integral (2) over is also calculated analytically and it is expressed via residues. The singularities in the BS equation are removed by the analytical integration over . Equation (2) is valid for an arbitrary kernel , given by a Feynman graph. The particular cases of the ladder kernel and of the Wick-Cutkosky model [2, 17] are detailed in the next section and for the cross ladder kernel – in the paper [18]. Once is known, the BS amplitude can be restored by eq. (2).

The variables () are related to the standard light front (LF) variables as , . The LF wave function (see its definition in [19], for instance) can be easily obtained by [15]:

 ψ(k⊥,x)=1√4π∫∞0x(1−x)g(γ′,1−2x)dγ′[γ′+k2⊥+m2−x(1−x)M2]2. (6)

Eq. (2) can be transformed, in principle, to the equation for the LF wave function , though this requires inverting the kernel in the left-hand side of (2). The initial BS equation (2), projected on the LF plane, can be also approximately transformed to the LF equation:

 (k2⊥+m2x(1−x)−M2)ψ(k⊥,x)=−m22π3∫ψ(k′⊥,x′)VLF(k′⊥,x′;k⊥,x,M2)d2k′⊥dx′2x′(1−x′) (7)

with the LF kernel given, for ladder exchange, in ref. [15].

It is worth noticing that the LF wave function (2) is different from the one obtained by solving the ladder LF equation (2), as it was done e.g. in ref. [20]. The physical reason lies in the fact that the iterations of the ladder BS kernel (Feynman graphs) and the ladder LF kernel (time-ordered graphs) generate different intermediate states. The LF kernel and its iterations contain in the intermediate state only one exchanged particle, whereas the iterations of the ladder Feynman kernel contain also, many-body states with increasing number of exchanged particles (stretched boxes). This leads to a difference in the binding energies, which is however small [20]. Formally, this difference arises because of the approximations which are made in deriving eq. (2) from (2). However, for a kernel given by a finite set of irreducible graphs, both BS (2) and LF (2) equations are already approximate and it is not evident which of them is more ”physical”. The physically transparent interpretation of the LF wave function makes it often more attractive.

## 3 Spinless particles. Ladder kernel

As illustration, we give here the kernel of equation (2) for the ladder BS kernel, which reads:

 iK(L)(k,k′,p)=i(−ig)2(k−k′)2−μ2+iϵ. (8)

We substitute it in eq. (2), then substitute (2) in (2) and calculate the integrals. The details of these calculations are given in ref. [15] and one obtains:

 V(γ,z;γ′,z′)={W(γ,z;γ′,z′),if −1≤z′≤z≤1W(γ,−z;γ′,−z′),if −1≤z≤z′≤1 (9)

where has the form:

 W(γ,z,γ′,z′)=αm22π(1−z)2D0∫10v2D2dv (10)

with and

 D0 = γ+m2z2+(1−z2)κ2 D = v(1−v)(1−z′)γ+(1−z)[(1−v)μ2+vγ′] + vm2[(1−v)(1−z′)z2+vz′2(1−z)]+vκ2(1−z)(1−z′)[1+z−v(z−z′)]

In the particular case we are considering here, the angular integral over can be performed analytically and obtains the simple expression:

 W(γ,z;γ′,z′) = αm22π(1−z)2γ+z2m2+(1−z2)κ2 × 1b22(b+−b−)3[(b+−b−)(2b+b−−b+−b−)(1−b+)(1−b−)+2b+b−logb+(1−b−)b−(1−b+)]
 b± = −12b2(b1±√b21−4b0b2),b0=(1−z)μ2, b1 = γ+γ′−(1−z)μ2−γ′z−γz′+(1−z′)[z2m2+(1−z2)κ2], b2 = −γ(1−z′)−(z−z′)[(1−z)(1−z′)κ2+(z+z′−zz′)m2].

One can show [15] that in the case , which constitutes the original Wick-Cutkosky model [2, 17], the solution has the form . The -dependence disappears from the equation which, after that, exactly coincides with the Wick-Cutkosky equation [2, 3, 17].

Equation (2) has been solved by using the method explained in Appendix A, i.e. by expanding the solution on a bicubic cubic spline basis. By keeping fixed and varying the grid parameters to ensure four digits accuracy, we obtain for and and unit constituent mass () the values displayed in table 1. They correspond to , , . With all shown digits, they are in full agreement with the results we have obtained, similarly to [20], by using the Wick rotation and the method of [5]. Increasing to changes at most one unit in the last digit. This demonstrates the validity of our approach.

The weight function for a system with and is plotted in fig. 1. It has been obtained with and the same grid parameters than in table 1. Its -dependence is not monotonous and has a nodal structure; the -variation is also non trivial. We have remarked a strong dependence of relative to values of the parameter smaller than , in contrast to high stability of corresponding eigenvalues. However the corresponding BS amplitude and LF wave function , obtained from by the integrals (2) and (2), show the same strong stability as the eigenvalues.

The BS amplitude in Minkowski space in the rest frame is shown in fig. 2. The -dependence is rather smooth but the -dependence, due to poles of the propagators in (2), exhibits a singular behaviour at , i.e. moving with and .

Note that our solution gives also the BS amplitude in Euclidean space, by substituting in (2) . It is indistinguishable from the one obtained by a direct solution of the Wick-rotated BS equation [15].

## 4 Spinless particles. Cross ladder kernel

Non-ladder effects, within the same model, using Feynman-Schwinger representation, were considered in ref. [22]. In this work the full set of all irreducible cross-ladder graphs in a bound state calculation was included. Refs. [23, 24] contain results on the binding energy found by solving the BS equation for (L+CL) in Euclidean space. In the LFD framework, binding energy for (L+SB) kernel is calculated in [25], whereas the (L+CL+SB) contribution is incorporated in [24]. In [26] the effect of the cross-ladder graphs in the BS framework was estimated with the kernel represented through a dispersion relation. Non-ladder self-energy effects have been incorporated in [27, 28, 11].

As mentioned, the equation (2) is valid for any kernel. Derivation of in this equation for the cross-ladder BS kernel, shown in fig. 3, is quite similar but more lengthy, since the kernel itself is more complicated. We have first to calculate to the kernel corresponding to the diagram in fig. 3, substitute the result in (2), then in (2) and find in this way the cross-ladder contribution to the kernel in equation (2). The full kernel – including ladder and cross-ladder graphs – will be written in the form:

 V(γ,z;γ′,z′)=V(L)(γ,z;γ′,z′)+V(CL)(γ,z;γ′,z′).

The ladder kernel is given in refs. [15, 19]. The cross-ladder contribution was calculated in the paper [18].

We compare the results obtained in the BS approach with the equivalent ones found in Light-Front Dynamics (LFD). The latter incorporates the graphs with two bosons in flight forming six cross boxes and two stretched boxes. They were also calculated in [18],

The binding energy as a function of the coupling constant is shown in figures 4 for exchange masses and respectively and .

We see that for the same kernel – ladder or (ladder + cross-ladder) – and exchange mass – or – the binding energies obtained by BS and LFD approaches are very close to each other. The BS equation is slightly more attractive than LFD. At the same time, the results for ladder and (ladder +cross-ladder) kernels considerably differ from each other. The effect of the cross-ladder is strongly attractive. Though the stretched box graphs are included, its contribution to the binding energy is smaller than 2% and attractive as well.

The zero binding limit of fig. 4 deserves some comments. It was found in [20] that for massive exchange, the relativistic (BS and LF) ladder results do not coincide with those provided by the Schrödinger equation and the corresponding non relativistic kernel (Yukawa potential) even at very small binding energies. Their differences increase with the exchanged mass and do not vanish in the limit . We have displayed in fig. 5 a zoom of fig. 4 (right) for small values of . The cross ladder and stretched box diagrams reduce the differences but are not enough to cancel it.

We see that the cross-ladder contribution, relative to the ladder one, results in a strong attractive effect. The BS and LFD approaches give very close results for any kernel, with BS equation being always more attractive. These approaches differ from each other by the stretched-box diagrams with higher numbers of intermediate mesons. Our results indicate that the higher order stretched box contributions are small. This agrees with direct calculations in LFD of stretched box kernel with two-meson states [29] and with calculations of the higher Fock sector contributions [30] in the Wick-Cutkosky model. Calculation in LFD of binding energy with the stretched box contribution (L+CL+SB) and its comparison with (L+CL) also shows that the stretched box contribution is attractive but small.

The comparison of our results with those obtained in [22], evaluating the binding energy for the complete set of all irreducible diagrams, shows that the effect of the considered cross ladder graphs, though being very important, represent only a small part of the total correction. Thus for and the corresponding binding energies obtained with BS equation are , and .

Qualitatively, our large attractive effect of CL (for and ) is in agreement with [23] (where ). For , our numerical results for binding energy, both for the BS equation in Minkowski space and in the LFD approach, are close to those found in [24] (within accuracy of extracting the data from plots). Our results are smaller than the ones found in [26] by a factor 3.

## 5 Two fermions

We have considered the following fermion () - meson () interaction Lagrangians:

(i) Scalar coupling for which

(ii) Pseudoscalar coupling for which .

(ii) Massless vector exchange with and as vector propagator.

Each interaction vertex has been regularized with a vertex form factor by the replacement and we have chosen in the form:

 F(q)=μ2−\mathchar259\relax2q2−\mathchar259\relax2+iϵ. (12)

Let us first consider the case of the scalar coupling and the corresponding ladder kernel. The BS equation for the amplitude reads:

 \mathchar264\relax(k,p)=i(m+12^p+^k)(12p+k)2−m2+iϵ[∫d4k′(2π)4\mathchar264\relax(k′,p)(−ig2)F2(k−k′)(k−k′)2−μ2+iϵ]i(m−12^p+^k)(12p−k)2−m2+iϵ. (13)

In the case of state, the BS amplitude has the following general form:

 \mathchar264\relax(k,p)=S1ϕ1+S2ϕ2+S3ϕ3+S4ϕ4 (14)

where are independent spin structures ( matrices) and are scalar functions of and . The choice of is to some extent arbitrary. To benefit from useful orthogonality properties we have taken

 S1=γ5,S2=1M^pγ5,S3=k⋅pM3^pγ5−1M^kγ5,S4=iM2σμνpμkνγ5,

where . The antisymmetry of the amplitude (14) with respect to the permutation implies for the scalar functions: , . A decomposition similar to (14) was used in [12] to solve the BS equation for a quark-antiquark system but the solution was approximated keeping only the first term .

We substitute (14) in eq. (13), multiply it by and take traces. As we will see, the kernel in the resulted equation, in contrast to the spinless case, is still singular. These singularities are integrable numerically. They do not prevent from finding numerical solution, but they reduce its precision. This can be avoided by a proper regularization of equation, multiplying both sides of it by the factor

 (15)

This factor has the form of a product of two scalar propagators with mass . It plays the role of form factor suppressing the high off-mass shell values of the constituent four-momenta and tends to 1 when . In this way, we get the following system of equations for the invariant functions :

 η(k,p)ϕi(k,p) = η(k,p)[(p2+k)2−m2+iϵ][(p2−k)2−m2+iϵ] (16) × ∫d4k′(2π)4ig2F2(k−k′)(k−k′)2−μ2+iϵ4∑j=1cij(k,k′,p)ϕj(k′,p),

Since , the equation thus obtained is strictly equivalent to one where is cancelled. We will see however that, due to the presence of the factor, the LF projection modifies the resulting kernels which become less singular functions.

The coefficients are determined by traces. We do not give here their explicit form, which can be found in [16].

Then we represent each of the BS components by means of the Nakanishi integral (2) and, similarly to the scalar case, apply the light-front projection to the set of coupled equations for the corresponding weight functions . As explained in sect. 2, this projection, which is an essential ingredient of our works [15, 18, 31], consists in replacing in eq. (16) and integrating over in all the real domain.

The technical details of the light-front projection are similar to those given in ref. [15]. We obtain in this way a set of coupled two-dimensional integral equations:

 ∫∞0dγ′∫1−1dz′Vg(γ,z;γ′,z′)gi(γ′,z′)=∑j∫∞0dγ′∫1−1dz′Vdij(γ,z;γ′,z′)gj(γ′,z′) (17)

The kernel and also for all types of couplings and states are given in [16]. These kernels depend on the parameter . Closer is to , smoother is the kernel and more stable are the numerical solutions. However the weight functions as well as binding energies provided by (17) are independent of . To avoid spurious singularities in (17) due to the factor (15), must be larger than , what is fulfilled for . In practical calculations we have taken and .

The kernel is finite and it vanishes for . For a fixed values of and , is a continuous function of with a discontinuous derivative at .

As already mentioned, without -factor, most of the kernel matrix elements are singular. Namely, they are discontinuous at . In some cases – like e.g. displayed in fig. 6 (left) – the value of the discontinuity, although being finite at fixed value of , diverges when . This creates some numerical difficulties when computing the solutions in the vicinity of . They can be in principle solved by properly taking into account the particular type of divergence. However the latter may depend on the particular matrix element, on the type of coupling, the quantum number of the state and other details of the calculation.

For , i.e., for a finite value of in (15), the -dependence of the regularized kernels is much more smooth and therefore better adapted for obtaining accurate numerical solutions. In fig. 6 (right) we plotted the regularized kernel as a function of for the same arguments and parameters than in fig. 6 (left), where it was calculated without the factor. As one can see, the kernel is now a continuous function of . A discontinuity of the first derivative however remains at .

We would like to emphasize again that despite the fact that the non-regularized and regularized kernels are very different from each other (compare e.g. figs. 6 left and right) and that the regularized ones strongly depends on the value of , they provide – up to numerical errors – the same binding energies and weight functions . We construct in this way a family of equivalent kernels.

The solutions of eq. (17) have been obtained using the techniques detailed in Appendix A. We have computed the binding energies and BS amplitudes, for the two fermion system interacting with massive scalar (S) and pseudoscalar (PS) exchange kernels and for the fermion-antifermion system interacting with massless vector exchange in Feynman gauge. In the limit of an infinite vertex form factor parameter , the later case would correspond to positronium with an arbitrary value of the coupling constant. All the results presented in this section are given in the constituent mass units () and with .

The binding energies obtained with the form factor parameter are given in left table 2. For the scalar and pseudoscalar cases, we present the results for and boson masses. They have been compared to those obtained in a previous calculation in Euclidean space [32] using a slightly different form factor, which differs from our one by a factor. Once taken into account this correction, our scalar results are in full agreement (four digits) with [32]. The pseudoscalar ones show small discrepancies (). We have also computed the binding energies by directly solving the fermion BS equation the Euclidean space using a method independent of the one used in [32]. Our Euclidean results are in full agreement with those given in the table 2.

The dependence for the scalar and pseudoscalar couplings is plotted in figs. 7. Notice the different scales of both dependences. The pseudoscalar binding energies are fast increasing functions of and thus more sensitive to the accuracy of numerical methods. This sharp behaviour was also exhibit when solving the corresponding LF equation [33].

It is worth noticing that the stability properties of the BS solutions for the scalar coupling are very similar to the LF ones. In the latter case, we have shown [34, 35] the existence of a critical coupling constant below which the system is stable without vertex form factor while for , the system ”collapses”, i.e. the spectrum is unbounded from below. The numerical value was found to be [34, 35], which corresponds to . Performing the same analysis than in our previous work – eq. (71) from [34, 35] – we found that for BS equation the critical coupling constant is , in agreement with [32]. The difference between the numerical values of is apparently due to the different contents of the intermediate states in the two approaches. The ladder BS equation incorporates effectively the so-called stretch-boxes diagrams which are not taken into account in the ladder LF results.

The positronium case deserves some comments. First we would like to notice that in our formalism, the singularity of the Coulomb-like kernels in terms of the momentum transfer is absent. This is a combined consequence of the Nakanishi transform (2) – which allows to integrate over analytically in the right hand side of the BS equation (16) – and of the consecutive LF projection integral. After this integration, the Coulomb singularity does not anymore manifest itself in the kernel. This can be explicitly seen in the kernel of the Wick-Cutkosky model obtained in eq. (22) of our previous work [15].

A second remark concerns the dependence of the positronium results. Using the methods developed in [34, 35] we found that in the BS approach with ladder kernel there also exists a critical value of the coupling constant . Note that, as in the scalar coupling, the very existence and the value of this critical coupling constant is independent on the constituent () and exchange masses but depends on the quantum number of the state.

The ground state positronium binding energies without vertex form factor are given in the left table 2 for values of the coupling below , Nonrelativistic results are included for comparison. One can see that the relativistic effects are repulsive.

These results are displayed in fig. 8 (left), black solid line, and compared to the binding energies obtained with two values of the form factor parameter (dashed) and (dot-dashed). The stability region is limited by a vertical dotted line at . Beyond this value the binding energy without form factor becomes infinite and we have found .

The inclusion of the form form factor has a repulsive effect, i.e. for a fixed value of the coupling constant it provides a binding energy of the system which is smaller than in the limit (no cut-off). This is also illustrated in fig. 8 (right) where we have plotted the dependence of for two different energies. One can see that the value of the coupling constant to produce a bound state is a decreasing function of . The size of the effect depends strongly on the binding energy but for both energies the asymptotics is reached at . This behaviour is understandable in terms of regularizing the short range singularity of interactions.

We present in figs. 9 some examples of the Nakanishi weigh functions . They correspond to a state with the scalar coupling and the same parameters , than in table 2. In the left figure is shown the -dependence for a fixed value of and in the right figure the -dependence for a fixed . Notice the regular behaviour of these functions as well as their well defined parity with respect to are even and is odd. As in the scalar case, the -dependence of is more important than for the binding energy.

Corresponding BS amplitudes are displayed in figs. 10. The figure 10 (left) represents the dependence of for a fixed value of . They exhibit a singular behaviour which corresponds to the pole of free propagators , like in the scalar case. The figure 10 (right) represents the dependence of the amplitudes for a fixed value . For this choice of arguments, the amplitudes are smooth functions of , though they will be also singular for . Notice that all the infinities in the BS amplitudes come from the free propagator. However the BS amplitude has other non-analytical points. For instance it obtains an imaginary part for .

## 6 Electromagnetic form factor

We demonstrate in this section the advantages of using the BS amplitude in Minkowski space for computing the electromagnetic (em) form factor shown in fig. 11.

For simplicity, we will restrict ourselves to the spinless case for which one has the following expression in terms of the Minkowski BS amplitude

 (p+p′)μFM(Q2)=−i∫d4k(2π)4(p+p′−2k)μ(m2−k2)\mathchar264\relaxM(p2−k,p)\mathchar264\relaxM(p′2−k,p′) (18)

As mentioned in the Introduction, the Euclidean BS amplitude cannot provide the right em form factor. This fact has been discussed in some detail in our previous work [37, 38] where two kinds of reasons were advocated.

The first reason was the impossibility for performing the Wick rotation in the integral (18). Indeed this transformation, represented schematically in figure 12, consists in replacing the integral over the variable along the real axis by the integral along the imaginary axis . The result is unchanged provided no any singularity of the integrand is crossed when rotating anticlockwise the real into the imaginary integration interval, that is when there are no singularities in the first or third quadrants. This is the case, for instance, of those poles associated with the free propagators (denoted by in fig. 12) which are located in the second and fourth quadrants. It turns out however that the integrand of (18) has singularities in the variable (denoted by X in fig. 12) which, for some values of , are located in the first or third quadrants. That prevents from performing a naive Wick rotation, i.e. by simply replacing in equation (18) and , the Euclidean BS amplitudes. Taking properly into account the contribution of these additional singularities in the contour integral would require a careful analysis on the complex plane of the form factor integrand, a task which is hardly doable even by assuming that and are known numerically.

This problem is already manifested at zero momentum transfer in the simplest model case of a BS amplitude given by the product of free propagators corresponding to the external lines:

 \mathchar264\relaxM=1[k2−m2+iϵ][(p−k)2−m2+iϵ] (19)

In this case, the integral (18) corresponds to the Feynman amplitude represented in fig. 11 with constant vertices (). The form factor is obtained by multiplying its result (18) by . For it reads:

 4M2FM(0)=i∫d4k(2π)4(4M2−4k\makebox[2.27622pt]$⋅$p)[k2−m2+iϵ][(p−k)2−m2+iϵ]2 (20)

Since the BS amplitude (19) – and consequently – are not normalized, eq. (20) determines the normalization factor. Using the Feynman parametrization

 1ab2=∫102xdx[(1−x)a+xb]3

the integral (20) can be calculated analytically and reads:

 FM(0)=116π2M3⎛⎜ ⎜⎝4m2arctanM√4m2−M2√4m2−M2−M⎞⎟ ⎟⎠ (21)

which for and takes the numerical value .

In the rest frame the form factor (20) has the form:

 FM(0)=i(2π)4M∫∞04πk2dk∫∞−∞(M−k0)dk0[k20−k2−m2+iϵ][(k0−M)2−k2−m2+iϵ]2 (22)

In order to check the validity of the ”naive Wick rotation” in the form factor integral (22), we replace there and integrate over in real infinite limits. We obtain in this way:

 FNW(0)=−i(2π)4M∫∞04πk2dk∫∞−∞(M−ik4)idk4[k24+k2+m2][(k4+iM)2+k2+m2]2 (23)

We have omitted in denominator, since the integrand in (23) is now non-singular at real . For particular values of the parameters this integral can be computed numerically. For , we get: which differs from the value obtained above, using Minkowski BS amplitudes. This shows, by a simple example, that performing the Wick rotation to the variable – the argument of the BS amplitude – in the expression for the em form factor is not allowed.

The reason for this difference is the existence of a pole in the integrand of (22) created by the zeros of the the denominator rightest factor which are given by . If , the pole corresponding to

 k0=M−√k2+m2+iϵ,

lies in the first quadrant and it is crossed by the Wick rotation. The residue at this pole is

 Res=−i(M−√k2+m2)25π2M3√k2+m2(2√k2+m2−M)2

and its contribution to the form factor is given by:

 FRes(0)=∫√M2−m20(2πiRes)4πk2dk (24)

Numerical calculation for and gives and the sum of both contributions is in perfect agreement with the value .

 FNW(0)+FRes(0)=3.15404⋅10−4+1.83837⋅10−4=4.99241⋅10−4

This example illustrates, hopefully clearly, how the singularities in the first quadrant indeed appear and that they prevent from the Wick rotation. Taking properly into account their contribution restores the form factor value. However, this cannot be done if one knows only the numerical values of the Euclidean BS amplitude.

The second reason results from the fact that even neglecting the above explained impossibility of the Wick rotation in the variable , i.e., neglecting the contribution of the singularities denoted by X in fig. 12, we still cannot compute the form factor in terms of the Euclidean BS amplitude at rest. Indeed, the amplitude depends on the scalars and . After Wick rotation , the first scalar becomes . Suppose that we are working at the rest frame . Then the second scalar turns into . This purely imaginary value is just the argument of the Euclidean BS amplitude and we do not cause any problem.

However, for non-zero momentum transfer we have . The argument of the amplitude – right vertex in fig 11 – transforms in the Wick rotation as , i.e. a complex value. The amplitude with complex – not purely imaginary – value of its argument is not the Euclidean one, which enters the left vertex. Therefore, even after performing the invalid Wick rotation, the form factor integrand, cannot be expressed via BS amplitude obtained by real boost from the Euclidean BS amplitude at rest.

To avoid the complex boost one can solve the Euclidean BS equation in a moving frame with and obtain in this way the BS amplitude for any value of . This approach was followed in [39, 40] in the framework of the quark model. The form factor was computed in terms of the Euclidean amplitudes at but without taking into account the contribution of the singularities discussed above.

Another approximate way to compute the form factors, widely used in the literature, is the so called static approximation [36]. It consists in replacing the BS amplitude in a moving frame by the rest frame amplitude boosting only the spatial part of its argument , i.e. letting unchanged the time-component . This means in practice that the left vertex of fig. 11 involves the Euclidean BS amplitude with , , while the right vertex involves . In this way, the form factor is approximately expressed via Euclidean BS amplitude. The accuracy of static approximation was estimated in [36] perturbatively, at relatively small .

When using the Minkowski BS amplitude, represented via Nakanishi integral (2), we can calculate the integral of analytically and express the form factor in terms of the computed Nakanishi weight function . The procedure is now theoretically safe and the precision of the result depends only on the accuracy of the numerical solution of the functions.

The exact results [37, 38], together with the static approximation, are shown in fig. 13 (left). The difference between the exact calculation and static approximation is small at small momentum transfer but it strongly increases with increase of .

In fig. 13 (right) we compare the exact (Minkowski) form factor with the one found through the LF wave function (2). At all they are almost indistinguishable from each other.

## 7 Conclusions

We presented a new method for obtaining the solutions of the Bethe-Salpeter equation in Minkowski space, both for two spinless bosons and for two fermions. It is based on a Nakanishi integral representation of the BS amplitude and on subsequent Light-Front projection.

The binding energies for these systems are calculated. They coincide with the ones found via Euclidean space solution, thus providing a validity test of our method. For two fermions, the solutions for the scalar exchange and positronium states without vertex form factor () are stable below some critical value of the coupling constant, respectively (scalar exchange) and (positronium).

The BS amplitudes are obtained in terms of the computed Nakanishi weight functions. They exhibit a singular behaviour due, on one hand, to the poles of the free propagators and, on the other hand, to non analytic branch cuts which are responsible for the appearance of an imaginary part.

Some applications to em form factors are also presented. We demonstrate that the Wick rotation in the variable – argument of the BS amplitude – is not allowed. It is prevented by the singularities which appear in the first quadrant of the complex -plane and whose contributions are hardly evaluable in practice. Neglecting these contributions, one obtains an approximate expression which supposes the knowledge of Euclidean BS amplitude in a non-zero momentum frame. The latter requires to solve numerically a much more involved - though solvable - equation than at rest. On the contrary, in terms of the BS amplitude in Minkowski space, represented in the Nakanishi form (2), the form factors can be calculated without any problem.

## Appendix A: Numerical methods

The numerical methods wil be illustrated by considering the general case of coupled two-dimensional integral equations for the Nakanishi weight functions in the form:

 ∫∞0dγ′∫+1−1dz′Vga(γ,z,γ′,z′)ga(γ′,z′)=na∑a,a′=1∫∞0dγ′∫+1−1dz′Vdaa′(γ,z,γ