# Unitarization of electron scattering with an external potential at NLO in
QED

J. A. Oller<EMAIL_ADDRESS>Departamento de Física, Universidad de Murcia, E-30071
Murcia, Spain Marcela Peláez<EMAIL_ADDRESS>Instituto de Física, Facultad
de Ingeniería, Universidad de la República,
J. H. y Reissig 565, 11000 Montevideo, Uruguay

###### Abstract

We have calculated the one-loop scattering amplitude of an electron by an
external Coulomb potential in QED free of infrared divergences. This feature
is achieved by applying the Faddeev-Kulish formalism, which implies a
redefinition of both the asymptotic electronic states and of the $S$ matrix.
Additionally, we have also derived the infrared-finite one-loop partial-wave
amplitudes for this process by applying a recent method in the literature.
Next, these partial-wave amplitudes are unitarized based on analyticity and
unitarity by employing three different methods of unitarization: the algebraic
$N/D$ method, the Inverse Amplitude Method and the first-iterated $N/D$
method. Then, we have studied several partial waves both for physical momentum
and for complex ones to look for bound-state poles. The binding momentum for
the fundamental bound state in $S$ wave is discussed with special detail. This
is a wide-ranging method for calculating nonperturbative partial-wave
amplitudes for infinite-range interactions that could be applied to many other
systems.

###### Contents

  1. 1 Introduction
  2. 2 Standard Diagrams
    1. 2.1 Leading-order amplitude
    2. 2.2 Vacuum polarization contribution
    3. 2.3 Vertex correction
    4. 2.4 The iterated diagram
  3. 3 Infrared divergences
    1. 3.1 Redefinition of the asymptotic electron states
    2. 3.2 Redefinition of the $S$ matrix
  4. 4 Partial-wave amplitudes
    1. 4.1 Finite PWAs after redefinition of the $S$ matrix and unitarity
    2. 4.2 Handling different angular-projection integrals
  5. 5 Unitarization methods
    1. 5.1 The $g$-method
    2. 5.2 The Inverse Amplitude Method
    3. 5.3 First-iterated $N/D$ method
  6. 6 Scattering in partial waves
    1. 6.1 $J=1/2$ in $S$-wave
    2. 6.2 $J=3/2$ and $J=1/2$ in $P$-waves
  7. 7 The bound-state equation for the different unitarization methods
  8. 8 Binding momentum for the ground state in the $J=1/2$ and $\ell=0$ PWA
    1. 8.1 Bound state from the unitarized LO QED PWA
    2. 8.2 Bound state from the unitarized NLO QED PWA
  9. 9 Conclusions and outlook
  10. A QED one-loop results
    1. A.1 Electron self energy
    2. A.2 Vacuum polarization
    3. A.3 Vertex modification of the electron-photon coupling
    4. A.4 Iterated diagram
    5. A.5 Master Integrals
      1. A.5.1 Master integrals in Euclidean three dimensions
  11. B No infrared-divergent phase for $n=m$ in Ref. [1]
  12. C Partial-wave amplitudes in the $\ell SJ$ basis

##  1 Introduction

The scattering of particles by a Coulomb field is a fundamental topic because
it represents the interaction between charged particles, such as electrons and
nuclei. In nonrelativistic scenarios, the scattering can be computed by
solving the Schrödinger equation, yielding closed-form solutions for the
scattering amplitude or differential cross section. As it is well known, the
cross section for nonrelativistic electron scattering off a central Coulomb
potential is described by the classical Rutherford formula. The wave function
for the scattering of a Dirac electron by a Coulomb potential was obtained in
Ref. [2] as a series of spherical harmonics. Based on this calculation, the
lowest-order correction to the Rutherford cross section for a Dirac electron
(positron) scattering off a heavy, point-like target with charge Z was first
derived, though the correct result had to await until Dalitz’s calculation in
Ref. [3] by evaluating the second Born approximation.

To obtain more precise and accurate predictions, it is necessary to consider
higher-order corrections in QED, typically calculated as a power series
expansion in $\alpha=e^{2}/4\pi$, the fine structure constant, with $e$ the
electron charge. Higher-order perturbation theory for this problem has been
explored in several articles, such as [2, 4, 3, 5, 6]. Additionally, two-loop
corrections of the $g$-factor in a nonperturbative Coulomb field are discussed
in [7]. Similar two-loop calculations for processes like Möller scattering [8,
9], Bhabha scattering [10, 11], and electron-muon scattering [12] can be
found. Of particular note, Refs. [13, 14] compute the high-energy elastic
scattering cross section for an electron off a classical point source to next-
to-leading order (NLO) accuracy. They employ the MS renormalization scheme to
manage ultraviolet divergences, while the Kinoshita-Lee-Nauenberg theorem
addresses infrared singularities at the level of the cross section.

In this work, we are interested in the scattering of an electron by an
external Coulomb potential, as an example of an infinite-range potential. Our
central objects are the partial-wave amplitudes (PWAs), first evaluated
perturbatively in QED up to NLO, and then unitarized by implementing non-
perturbative unitarization techniques [15, 16, 17, 18]. These results are also
applied to computing the bound states of an electron in a Coulomb potential.
By considering corrections to one loop in QED, in addition to the once-time
iterated potential, already provided by the second Born approximation of the
Dirac electron [3], we take into account the effects of the particle self-
energy and vertex-modified corrections, like the anomalous magnetic moment.
The latter ones modify the energy of the bound state compared to the
prediction based on the Dirac equation.

Despite the lack of convergence of a partial-wave expansion for the Coulomb
potential, they are often used when considering hadron and nuclear physics.
For instance, Coulomb effects [19] are important since they typically
contribute to the low-energy cross section near the threshold of a reaction
involving charged particles. Furthermore, the study of the modifications of
Coulomb atoms due to strong interactions is also an important technique for
hadron physics. And vice versa, an accurateness calculation of the hadronic
spectrum requires to take into account electromagnetic effects.

This article builds upon the work of Ref. [18], which investigates
nonrelativistic Coulomb scattering as a benchmark for the unitarization method
of Refs. [16, 17] for the case of infinite-range potentials, due to the null
mass of photons and gravitons. It was also shown as a clear counterexample
that the unitarization method used in Ref. [20] is unsuitable for unitarizing
Coulomb scattering, because the impact of infrared divergences on scattering
amplitudes is not properly taken into account. Ref. [18] shows that the exact
Coulomb $S$-wave cannot be approximated using the approach described in Ref.
[20], because it is typically mismatched by an order of magnitude at least.
Furthermore, it lacks of bound states with non-zero binding momentum.

Our present study extends the work of [18] by unitarizing the scattering
amplitude of an electron by an external Coulomb potential that we calculate up
to NLO in perturbative QED. The perturbative QED scattering amplitudes are
affected by two types of infinities, infrared (IR) and ultraviolet (UV)
divergences. The latter can be treated according to the general theory of
renormalization, that we apply here to end with UV-finite scattering
amplitudes. For the IR divergences one has to account for two main effects in
order to finally have scattering amplitudes free of IR divergences. At this
point it is important to differentiate between the standard formalism to treat
IR divergences at the level of transition rates (cross sections) [21, 1], and
the less used Faddeev-Kulish (FK) formalism [22], which allows to obtain
finite $S$-matrix elements.

On the one hand, there are IR divergences at the one-loop level that in the
formalism based on the finiteness of transition rates requires to account for
the radiation of low-energy photons. In the FK formalism the IR divergences of
this type are taken into account by a redefinition of the asymptotic
electronic states based on introducing a proper cloud of real soft photons.
Such a definition of asymptotic electron states (without hard photons) was
firstly proposed in the literature by Chung in Ref. [23]. This step to cure IR
divergences was not needed in Refs. [16, 17, 18] because only tree-level
amplitudes were unitarized. Therefore, they are accounted here for the first
time in order to unitarize infinite range interactions.

On the other hand, there is still the issue of the well-known divergent
Coulomb phase that already appears in the nonrelativistic case [3, 1, 22]. To
account for this type of divergence one needs to redefine the $S$-matrix
operator itself by a phase factor depending on the initial and final states
[22]. The partial-wave projection of this new $S$ matrix is considered, and in
terms of it one can define IR-finite PWAs, as firstly derived in Refs. [16,
17, 18].

By applying this scheme of work we start by calculating perturbatively the IR-
finite PWAs up to NLO in QED. Next, we unitarize them by employing three
different unitarization methods: the algebraic $N/D$ method, $g$-method for
brevity [24, 25], the Inverse Amplitude Method (IAM) [26, 27, 28], and the
first-iterated $N/D$ method [29, 30, 15]. Several PWAs are calculate with
these unitarization methods, and we compare between them and with the
perturbative PWAs from NLO QED. The resulting PWAs from the different methods
are typically in agreement in a broad range of momentum, $m\alpha\lesssim
p\lesssim m/\alpha$. In the lower end nonperturbative effects become very
important, since one is close to the hydrogen atom, whose reproduction can
only be accomplished after unitarization, so that the perturbative QED PWAs
depart from the unitarized ones. For the upper limit, it occurs that terms
suppressed by higher powers of $\alpha$ become more important than those with
a lower power due to their momentum dependence. This implies that the standard
IAM method can fail at such high momenta (as it does, in fact, for some
partial waves), similarly as it is deficient in meson-meson scattering around
an Adler zero, and variants of the method should be used [31, 32, 33].

We also pay special attention to the bound-state equations that results from
the three unitarization methods, and their application to calculate the
binding momenta. We argue that because at NLO the perturbative PWAs have a
left-hand cut (LHC) for $p^{2}<0$ in the complex $p^{2}$ plane, it is only the
first-iterated $N/D$ method that is robust enough to deal with this situation.
The special difficulty stems from the fact that the pole associated to a bound
state in the PWA lies on top of the LHC. Importantly, the application of our
unitarization scheme offers a novel procedure to study bound states in QED
[34, 35], that is improvable order by order as the perturbative PWAs are
calculated at higher orders in QED. Further application to higher orders and
to different systems would be valuable in the future.

The article is organized as follows. After this introduction we elaborate in
Sec. 2 on the perturbative calculation of QED scattering amplitudes up to NLO,
with more technical details and a compendium of formulas given in appendix A.
In order to calculate the IR-finite PWAs of Sec. 4 we first discuss the
removal of these divergences in Sec. 3, with its associated appendix B. The
unitarization methods are introduced in Sec. 5, with Sec. 6 dedicated to their
application for studying scattering amplitudes. The binding equation is the
object of study of Sec. 7, and its detailed application to the ground state is
discussed in Sec. 8. Finally, concluding remarks and a brief outlook are given
in Sec. 9. The appendix C is added to give more technical details for deriving
the final formula used to calculate PWAs.

##  2 Standard Diagrams

The calculation up to NLO in QED of the scattering amplitude of an electron by
an external potential produced by an infinitely heavy particle of charge
$e^{\prime}$ can be obtained from QED Lagrangian:

$\displaystyle\mathcal{L}_{\text{QED}}=\bar{\psi}(i\gamma^{\mu}\partial_{\mu}-m_{0})\psi-\frac{1}{4}F^{\mu\nu}F_{\mu\nu}-e_{0}\bar{\psi}\gamma^{\mu}A_{\mu}\psi$
(2.1)

where $m_{0}$ and $e_{0}$ represent the bare mass and electron charge,
respectively.

As one loop diagrams are plagued with UV divergences it is necessary to
introduce the renormalization factors defined as:

$\displaystyle\psi$ $\displaystyle=\sqrt{Z_{2}}\psi^{R}~{},$ (2.2)
$\displaystyle A_{\mu}$ $\displaystyle=\sqrt{Z_{3}}A_{\mu}^{R}~{},$
$\displaystyle e_{0}$ $\displaystyle=Z_{e}e~{},$ $\displaystyle m_{0}$
$\displaystyle=Z_{m}m~{}.$

The Lagrangian expressed in term of counterterms takes the form (now the
constants and fields are the renormalized ones. The latter are expressed again
as $\psi$ and $A_{\mu}$ for ease of writing),

$\displaystyle\mathcal{L}_{\text{QED}}$
$\displaystyle=\bar{\psi}(i\gamma^{\mu}\partial_{\mu}-m)\psi-\frac{1}{4}F^{\mu\nu}F_{\mu\nu}-e\bar{\psi}\gamma^{\mu}A_{\mu}\psi$
(2.3) $\displaystyle+\delta
Z_{2}\bar{\psi}(i\gamma^{\mu}\partial_{\mu})\psi-(Z_{2}Z_{m}-1)\bar{\psi}m\psi-\delta
Z_{3}\frac{1}{4}F^{\mu\nu}F_{\mu\nu}-\delta
Z_{1}e\bar{\psi}\gamma^{\mu}A_{\mu}\psi~{},$

where $\delta Z_{2}=Z_{2}-1$, $\delta Z_{3}=Z_{3}-1$ and $\delta
Z_{1}=Z_{1}-1$ being $Z_{1}=Z_{2}Z_{e}\sqrt{Z_{3}}$.

With the QED Lagrangian we can compute Feynman diagrams up to one loop, or
equivalently up to ${\cal O}(\alpha^{2})$, by following the standard
techniques [36, 37]. The first part of the Lagrangian gives the standard
Feynman rules in terms of the renormalized mass and charge, while the other
terms give additional rules for the counterterms. In our chosen
renormalization scheme, $Z_{2}$ also absorbs the IR divergences of the
electron self-energy and is renormalized to ensure that the electron
propagator has a residue equal to one. In this scenario, loop corrections to
external legs are nullified with the counterterm of the external leg, allowing
us to disregard diagrams modifying external legs.111This procedure is also
followed for the two-loop calculation of Babbha scattering in Ref. [11]. The
counterterm $Z_{3}$ is similarly fixed by requiring that the residue of the
photon propagator is $-g_{\alpha\beta}$, where
$g_{\alpha\beta}={\text{diag}}(1,-1,-1,-1)$ is the Minkowski metric. The set
of Feynman diagrams at LO, diagram (a), and NLO in QED, diagrams (b)–(f), are
shown in Fig. 1. The vertex-modification diagram [diagram (b)] encompasses
both IR and UV divergences, which are absorbed in $Z_{1}$. As it is well
known, the Ward identity dictates that $Z_{1}=Z_{2}$, and therefore
$\frac{1}{2}\delta Z_{3}+\delta Z_{e}=0$ is fixed by the definition of the
charge
$\bar{u}(\boldsymbol{p})\Lambda^{\mu}(p,p)u(\boldsymbol{p})=\bar{u}(\boldsymbol{p})\gamma^{\mu}u(\boldsymbol{p})$.
A detailed account of our procedure for the calculation of the renormalized
scattering amplitude for the scattering of an electron by an external
potential is given in the appendix A.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 1: All the one-loop diagrams contributing to the electron scattering by
a Coulomb potential.

###  2.1 Leading-order amplitude

The leading order (LO) amplitude of the electron scattering by an external
potential due to a charge $e^{\prime}$, represented by diagram (a) of Fig.1,
is

$\displaystyle A_{0}$
$\displaystyle=-ee^{\prime}\,\bar{u}_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\not{A}({\boldsymbol{q}}^{2})u_{\sigma}({\boldsymbol{p}})~{},$
(2.4)

where $e$ is the charge of the electron ($e<0$),
$u_{\sigma}({\boldsymbol{p}})$ is the Dirac spinor with linear momentum
${\boldsymbol{p}}$ and third component of spin $\sigma$ along the $z$ axis
(the Dirac spinors are always associated to on-shell electrons). For the case
of an infinitely heavy charged point particle one has the Coulomb potential
$\displaystyle{A_{\mu}({\boldsymbol{q}}^{2})=\frac{\delta_{\mu}^{0}}{{\boldsymbol{q}}^{2}+\lambda^{2}}}$,
where a vanishing small photon mass $\lambda$ has been included in order to
handle IR divergences.

###  2.2 Vacuum polarization contribution

We calculate the correction arising from the diagrams (d) and (e) of Fig. 1,
properly renormalized. The final result is

$\displaystyle A_{1}$
$\displaystyle=-(ee^{\prime})\bar{u}_{\sigma^{\prime}}\left({\boldsymbol{p}}^{\prime}\right)\not{A}({\boldsymbol{q}}^{2})u_{\sigma}({\boldsymbol{p}})\left\\{1+\frac{\alpha}{9\pi
t}\left[3\left(t+2m^{2}\right)\sigma\left(t\right)\ln\frac{\sigma\left(t\right)+1}{\sigma\left(t\right)-1}-5t-12m^{2}\right]\right\\}~{},$
(2.5)

which has no IR divergences. In the expression, the function
$\sigma\left(t\right)=\sqrt{1-4m^{2}/t}$ and $t=q^{2}$, where
$q=p^{\prime}-p$.

###  2.3 Vertex correction

The diagrams (b) and (c) of Fig. 1 account for this effect. The former diagram
is usually written using two form factors:

$\displaystyle\bar{u}_{\sigma^{\prime}}\left(\boldsymbol{p}^{\prime}\right)(-ie\Gamma^{\mu})u(\boldsymbol{p})=-ie\bar{u}_{\sigma^{\prime}}\left(\boldsymbol{p}^{\prime}\right)\left(F_{1}(q)\gamma^{\mu}+\frac{i\sigma^{\mu\nu}}{2m}q_{\nu}F_{2}(q)\right)u_{\sigma}(\boldsymbol{p})$
(2.6)

with $\sigma^{\mu\nu}=\frac{i}{2}[\gamma^{\mu},\gamma^{\nu}]$ and $q$ the
transfer momentum $q=p^{\prime}-p$. Standard techniques for computing one-loop
Feynman diagrams are used in order to obtain the final expression [37]. The
renormalized vertex is expressed in terms of the operator
$-ie\Lambda^{\mu}(p^{\prime},p)\equiv-
ie\left(\gamma^{\mu}+\Gamma^{\mu}(p^{\prime},p)+\delta
Z_{1}\gamma^{\mu}\right)$. In order to renormalize, we impose the following
renormalization condition:

$\displaystyle\bar{u}_{\sigma^{\prime}}\left(\boldsymbol{p}^{\prime}\right)\Lambda^{\mu}(p^{\prime},p)u_{\sigma}(\boldsymbol{p})\Big{|}_{\not{p^{\prime}}=\not{p}=m}=\bar{u}_{\sigma^{\prime}}\left(\boldsymbol{p}^{\prime}\right)\gamma^{\mu}u_{\sigma}(\boldsymbol{p})\Big{|}_{\not{p^{\prime}}=\not{p}=m}$
(2.7)

This gives the definition of the $\delta Z_{1}$ renormalization factor, being
$\delta Z_{1}=-ieF_{1}(0)$. Finally, the vertex contribution can be summarized
as :

$\displaystyle A_{2}$
$\displaystyle=-(ee^{\prime})A^{\mu}\left(\boldsymbol{q}^{2}\right)\bar{u}_{\sigma^{\prime}}\left(\boldsymbol{p}^{\prime}\right)\left[\frac{i\alpha
m\left(p^{\prime}-p\right)^{\nu}\sigma_{\mu\nu}}{2\pi\left(t-4m^{2}\right)}\sigma(t)\ln\frac{\sigma(t)-1}{\sigma(t)+1}\right.$
$\displaystyle+\frac{\alpha\gamma_{\mu}}{4\pi
t\sigma(t)}\left\\{4\left(t\sigma(t)+\left(t-2m^{2}\right)\ln\frac{\sigma(t)-1}{\sigma(t)+1}\right)\ln\frac{m}{\lambda}\right.$
$\displaystyle-4t\sigma(t)+\left(8m^{2}-3t+\left(t-2m^{2}\right)\ln\left(4-\frac{t}{m^{2}}\right)\right)\ln\frac{\sigma(t)-1}{\sigma(t)+1}$
$\displaystyle\left.\left.+2\left(2m^{2}-t\right)\left(\operatorname{Li}_{2}\left(\frac{1}{2}-\frac{1}{2\sigma(t)}\right)-\operatorname{Li}_{2}\left(\frac{1}{2}+\frac{1}{2\sigma(t)}\right)\right)\right\\}\right]u_{\sigma}(\boldsymbol{p})\,,$
(2.8)

and it has an infrared divergent term.

###  2.4 The iterated diagram

The diagram (f) of Fig.1 for a Coulomb external potential was first correctly
calculated in Ref. [3]. We have also calculated it independently and the final
expression for $\lambda\to 0$ is,

$\displaystyle A_{3}$
$\displaystyle=-(ee^{\prime})^{2}\bar{u}_{\sigma^{\prime}}\left({\boldsymbol{p}}^{\prime}\right)\left(E\gamma^{0}\left[\frac{1-\csc\frac{\theta}{2}}{32p^{3}\cos^{2}\frac{\theta}{2}}+i\frac{\ln\frac{2p}{\lambda}}{16\pi
p^{3}\cos^{2}\frac{\theta}{2}}-i\frac{1+\cos^{2}\frac{\theta}{2}}{32\pi
p^{3}\sin^{2}\frac{\theta}{2}\,\cos^{2}\frac{\theta}{2}\,\sqrt{1+\frac{\lambda^{2}}{p^{2}\sin^{2}\frac{\theta}{2}}}}\ln\frac{\sqrt{1+\frac{\lambda^{2}}{p^{2}\sin^{2}\theta/2}}+1}{\sqrt{1+\frac{\lambda^{2}}{p^{2}\sin^{2}\theta/2}}-1}\right]\right.$
$\displaystyle\left.+m\left[-\frac{1-\csc\frac{\theta}{2}}{32p^{3}\cos^{2}\frac{\theta}{2}}+i\frac{\ln\sin\frac{\theta}{2}}{16\pi
p^{3}\cos^{2}\frac{\theta}{2}}\right]\right)u_{\sigma}({\boldsymbol{p}}),$
(2.9)

Let us notice that $A_{3}$ is affected by IR divergences. Extra $\lambda$
dependence is kept in the denominators and argument of the log as directly
given by the loop integrals. They will be needed when performing the partial-
wave projection of this amplitude, which requires an angular integration with
$\theta\in[0,\pi]$. Note that for $\theta\to 0$ the factor
$\sqrt{1+\frac{\lambda^{2}}{p^{2}\sin^{2}\theta/2}}$ cannot be expanded in
powers of $\lambda$.

##  3 Infrared divergences

When computing scattering amplitudes in gauge theories, one encounters
infrared (IR) divergences linked to massless gauge fields [38]. Infrared
divergences are a common issue in quantum electrodynamics (QED), arising due
to the massless character of the photon, which results in the long-range
nature of electromagnetic interactions. These divergences occur when
considering the radiative contributions due to soft photons, i.e. photons
whose energy and momenta are much smaller than the typical energies and masses
in the process, when integrating over all possible photon momenta in loop
diagrams. Furthermore, the radiation of real soft photons also lead to the
same type of divergences when integrating over their phase space. In the
calculation of observables it is well known that these divergences cancel
between each other to all orders in perturbation theory leading to transitions
rates free of IR divergences [21, 39, 1].

However, this approach does not resolve the issue of IR divergences appearing
in transition amplitudes, c.f. Eqs.(2.3)-(2.4). These divergences stem from
the methodology of constructing a quantum field theory from free fields,
assuming that incoming and outgoing states do not interact in the asymptotic
limit (far before and after a scattering event). Conceptually one must take
into account that in quantum electrodynamics (QED) interactions are
significant even in the asymptotic limit due to the infinite range of massless
force-carrying bosons [22]. By considering these asymptotic interactions, one
can identify an operator that transforms free Fock states into asymptotically
interacting Faddeev-Kulish (FK) states [22]222See, e.g. [40] for a pedagogical
summary.. This type of states was first introduced by Chung in Ref. [23]. For
asymptotic states not containing hard photons but only charge particles, the
Chung’s and FK states are equivalent. A further point is that Ref. [22] also
considered the case when hard photons are part of asymptotic states, so that
FK states should be employed then. Heuristically, these FK states could be
pictured as consisting of charged particles and hard photons surrounded by a
cloud of real soft photons, making all scattering processes infrared finite at
the $S$-matrix level. The full FK procedure implies not only a redefinition of
the asymptotic states but also a redefinition of the $S$-matrix itself,
getting rid of a generalized Coulomb phase that does not affect the transition
rates. Directly solving the problem of infrared divergences at the level of
the $S$ matrix allows one to perform meaningful partial-wave expansions of QED
scattering amplitudes, unitarize them by employing techniques stemming from
general principles of $S$-matrix theory (for instance, commonly used in hadron
physics [41, 29, 42, 15]), and, in terms of them, to study bound states and
resonances. Previous studies along these lines are Refs. [16, 17, 18].

###  3.1 Redefinition of the asymptotic electron states

In order to introduce the FK states we need some preliminary notation. The
photon annihilation/creation operators
$a({\boldsymbol{k}})/a({\boldsymbol{k}})^{\dagger}$, correspond to

$\displaystyle
a^{(h)}({\boldsymbol{k}})=\varepsilon^{(h)}_{\mu}a^{\mu}({\boldsymbol{k}})~{},$
(3.1)

where $\varepsilon^{(h)}({\boldsymbol{k}})\,$ $h=1,2$ are the polarization
vectors; they are transversal, $\varepsilon^{(h)}({\boldsymbol{k}})\cdot k=0$,
and orthogonal
$\varepsilon^{(h)}({\boldsymbol{k}})^{*}\cdot\varepsilon^{(h^{\prime})}({\boldsymbol{k}})=-\delta_{hh^{\prime}}$.
The $a_{\mu}({\boldsymbol{k}})$ and
$a_{\nu}^{\dagger}({\boldsymbol{k}}^{\prime})$ satisfy the usual creation and
annihilation commutation relations

$\displaystyle[a_{\mu}({\boldsymbol{k}}),a_{\nu}({\boldsymbol{k}}^{\prime})^{\dagger}]$
$\displaystyle=-g_{\mu\nu}2k^{0}(2\pi)^{3}\delta({\boldsymbol{k}}-{\boldsymbol{k}}^{\prime})~{},$
(3.2)
$\displaystyle[a_{\mu}({\boldsymbol{k}}),a_{\nu}({\boldsymbol{k}}^{\prime})]$
$\displaystyle=[a_{\mu}({\boldsymbol{k}})^{\dagger},a_{\nu}({\boldsymbol{k}}^{\prime})^{\dagger}]=0~{}.$

In turn, the annihilation/creation operators for electrons are
$c_{\sigma}({\boldsymbol{p}})/c_{\sigma}({\boldsymbol{p}})^{\dagger}$, with
$\sigma=\pm 1/2$ the third component of the electron spin at rest.

The FK initial state may be expressed as [23, 22]

$\displaystyle\ket{{\boldsymbol{p}},\sigma}_{i}=e^{R_{f}(p_{i})}c_{\sigma}({\boldsymbol{p}})^{\dagger}\ket{0}$
(3.3)

where $R_{f}$ takes the form:

$\displaystyle
R_{f}(p)=\sum_{h=1}^{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}\left(S_{i}^{(h)}({\boldsymbol{k}})a^{(h)}({\boldsymbol{k}})^{\dagger}-S_{i}^{(h)}({\boldsymbol{k}})^{*}a^{(h)}({\boldsymbol{k}})\right)~{},$
(3.4)

The operator $\exp(R_{f}(p))$ of Eq. (3.4) can be interpreted as generating a
surrounding cloud of real soft photons around the external electrons. The
integral over ${\boldsymbol{k}}$ is cut at $|{\boldsymbol{k}}|<\Lambda$ to
preserve the soft character of the photons involved, being $\Lambda$ much
smaller than the characteristic momentum of the scattering process and the
electron mass. The function $S^{(h)}_{i}({\boldsymbol{k}})$ is derived from
the asymptotic potential, and in the case of QED [22] it can be written as
[23, 22]

$\displaystyle S_{i}^{(h)}({\boldsymbol{k}})$
$\displaystyle=-e\frac{p\cdot\varepsilon^{(h)}({\boldsymbol{k}})}{p\cdot
k}~{},$ (3.5)

where $p$ is the momentum of the incoming electron. Similarly, we define the
analogous factors $S_{f}^{(h)}({\boldsymbol{k}})$ associated to the final
electronic state of momentum $p^{\prime}$. Let us notice that that $R_{f}(p)$
respects the Gupta-Bleuler condition because it only involves transversal
polarization vectors [22].

(a)

(b)

(c)

(d)

(e)

Figure 2:  In this figure we show the Chung’s diagrams contribution to the
second-order scattering amplitude. The depicted photons not corresponding to
the external potential are real. Diagram (a): Soft photons from the clouds of
the initial and final electrons propagate without interacting. Diagrams (b)
and (c): Soft photons from the initial electron cloud are absorbed by the
initial and final electronic lines, respectively. Diagrams (d) and (e): Soft
photons are emitted from the initial and final electrons, respectively, to the
cloud of soft photons of the final state.

Next, we move to the right the annihilation operators in Eq. (3.4), by making
use of the commutation relations in Eq. (3.2). Then, one can rewrite the
states $\ket{{\boldsymbol{p}},\sigma}_{i}$ as

$\displaystyle\ket{{\boldsymbol{p}},\sigma}_{i}=\exp\left(-\frac{1}{2}\sum_{h=1}^{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}|S_{i}^{(h)}({\boldsymbol{k}})|^{2}\right)\exp\left(\sum_{h=1}^{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}S_{i}^{(h)}({\boldsymbol{k}})a^{(h)}({\boldsymbol{k}})^{\dagger}\right)c_{\sigma}({\boldsymbol{p}})^{\dagger}\ket{0}~{}.$
(3.6)

In our present calculation up to one loop in QED we only need to keep terms up
to ${\cal O}(e^{2})$ in the definition of $\ket{{\boldsymbol{p}},\sigma}_{i}$.
Therefore,

$\displaystyle\ket{{\boldsymbol{p}},\sigma}_{i}$
$\displaystyle=\left(1-\frac{1}{2}\sum_{h=1}^{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}|S_{i}^{(h)}({\boldsymbol{k}})|^{2}+\sum_{h=1}^{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}S_{i}^{(h)}({\boldsymbol{k}})a^{(h)}({\boldsymbol{k}})^{\dagger}+{\cal
O}(e^{2})\right)c_{\sigma}({\boldsymbol{p}})^{\dagger}\ket{0}$ (3.7)

Analogously, the final FK state is

$\displaystyle_{f}\\!\bra{{\boldsymbol{p}}^{\prime},\sigma^{\prime}}=\bra{0}c_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\left(1-\frac{1}{2}\sum_{h=1}^{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}|S_{f}^{(h)}({\boldsymbol{k}})|^{2}+\sum_{h=1}^{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}S_{f}^{(h)}({\boldsymbol{k}})^{*}a^{(h)}({\boldsymbol{k}})+{\cal
O}(e^{2})\right)~{}.$ (3.8)

The redefinition of the asymptotic states in Eqs. (3.7) and (3.8) generates
new contributions to the scattering amplitude. The corresponding diagrams up
to ${\cal O}(e^{2})$ are shown in Fig. 2. There, diagrams (b) and (c)
correspond to the absorption of real soft photons sourced by the initial-
electron cloud, while in diagrams (e) and (f) real soft photons are emitted to
the final-electron cloud. Non-interacting real soft photons that propagate
from the initial- to the final-electron clouds are shown in Diagram (a).

We refer to the diagrams in Fig. 2 as Chung’s diagrams, as they were
introduced by him for the first time in [23]. In this reference, Chung studied
in detail the scattering of a charged particle by an external field. It is
explicitly shown that the IR divergences from these diagrams cancel the IR
divergences from the vertex-modification diagram. The cancellation of IR
divergences was shown also in Ref. [23] to all orders in perturbation
theory.333No iteration diagrams were considered in Ref. [23]. Calculations are
done including a small photon mass $\lambda$ in order to explicitly handle
infrared divergences. Other approaches without a photon mass can be found in
[43].

To conclude, in order to take into account the FK states, together with the
standard diagrams in Fig. 1, we must consider all the possible interactions
involving real photons in the clouds of the electronic states in Fig. 2. The
explicit calculation of the latter diagrams give

$\displaystyle A_{\text{4}}$
$\displaystyle=-A_{0}\frac{e^{2}}{2}\int^{\Lambda}\frac{d^{3}k}{(2\pi)^{3}2k^{0}}\left(\frac{p}{p\cdot
k}-\frac{p^{\prime}}{p^{\prime}\cdot k}\right)^{2}~{}.$ (3.9)

For the evaluation of the integral we proceed as in Section V of Ref. [1], and
write,

$\displaystyle A_{4}$
$\displaystyle=-A_{0}\frac{e^{2}}{2(p^{0})^{2}}\int_{0}^{\Lambda}\frac{|{\boldsymbol{k}}|^{2}d|{\boldsymbol{k}}|}{(2\pi)^{3}2(k^{0})^{3}}\int
d\Omega_{{\boldsymbol{u}}}\left(\frac{p}{1-{\boldsymbol{v}}\cdot{\boldsymbol{u}}}-\frac{p^{\prime}}{1-{\boldsymbol{v}}^{\prime}\cdot{\boldsymbol{u}}}\right)^{2}~{},$
(3.10)

where ${\boldsymbol{u}}={\boldsymbol{k}}/|{\boldsymbol{k}}|$,
${\boldsymbol{v}}={\boldsymbol{p}}/p^{0}$, and
${\boldsymbol{v}}^{\prime}={\boldsymbol{p}}^{\prime}/p^{0}$ (let us recall
that $p^{0}={p^{\prime}\,}\\!^{0}$). In this way, all the IR divergences are
comprised in the first integral over the modulus of ${\boldsymbol{k}}$, while
the angular integration is just finite. The former gives

$\displaystyle\int_{0}^{\Lambda}\frac{|{\boldsymbol{k}}|^{2}d|{\boldsymbol{k}}|}{({\boldsymbol{k}}^{2}+\lambda^{2})^{3/2}}=\ln\frac{2\Lambda}{\lambda}-1+{\cal
O}(\lambda^{2})~{}.$ (3.11)

After performing the integration of the solid angle, the result from the
Chung’s diagrams of Fig. 2 is

$\displaystyle A_{4}$
$\displaystyle=\frac{A_{0}e^{2}}{4\pi^{2}}\left(-1+\frac{2m^{2}-t}{t\sigma(t)}\ln\frac{\sigma(t)-1}{\sigma(t)+1}\right)\left(\ln\frac{2\Lambda}{\lambda}-1\right)$
(3.12)

As said above, when added to the vertex-modification diagram (b) in Fig. 1 the
IR divergences mutually cancel. The final finite sum of amplitudes of the
vertex-modification and Chung’s diagrams, $A_{2}$ and $A_{4}$, respectively,
gives:

$\displaystyle A_{2}+A_{4}$
$\displaystyle=-(ee^{\prime})A^{\mu}({\boldsymbol{q}}^{2})\bar{u}_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\Bigg{[}\frac{i\alpha
m(p^{\prime}-p)^{\nu}\sigma_{\mu\nu}}{2\pi(t-4m^{2})}\sigma(t)\ln\frac{\sigma(t)-1}{\sigma(t)+1}$
(3.13) $\displaystyle+\frac{\alpha\gamma_{\mu}}{4\pi
t\sigma(t)}\Bigg{\\{}4\bigg{(}t\sigma(t)+(t-2m^{2})\ln\frac{\sigma(t)-1}{\sigma(t)+1}\bigg{)}\bigg{(}\ln\frac{m}{2\Lambda}+1\bigg{)}$
$\displaystyle-4t\sigma(t)+\bigg{(}8m^{2}-3t+(t-2m^{2})\ln(4-\frac{t}{m^{2}})\bigg{)}\ln\frac{\sigma(t)-1}{\sigma(t)+1}$
$\displaystyle+2(2m^{2}-t)\bigg{(}{\rm
Li}_{2}(\frac{1}{2}-\frac{1}{2\sigma(t)})-{\rm
Li}_{2}(\frac{1}{2}+\frac{1}{2\sigma(t)})\bigg{)}\Bigg{\\}}\Bigg{]}u_{\sigma}({\boldsymbol{p}})~{}.$

###  3.2 Redefinition of the $S$ matrix

The $S$ matrix is standardly defined as the transition amplitude between
asymptotic states, describing the scattering of particles in the far past to
the far future. However, when dealing with QED, the presence of soft photons
leads to infrared divergences, which pose a challenge to the standard
$S$-matrix approach. We have discussed in the previous section the
redefinition of the asymptotic states. But, on top of that, it is still
necessary to take care of a divergent phase factor affecting the $S$ matrix,
which was first established by Dalitz in Ref. [3]. There, the Born series for
the scattering of a Dirac electron by an external Coulomb potential of an
electron is studied up to second order, the diagram (f) in Fig. 1. It was
conjectured that the IR divergent terms affecting this diagram can be resumed
in a factor phase

$\displaystyle\exp\left(i\frac{ee^{\prime}p}{2\pi
E(p)}\ln\frac{\lambda}{\Lambda}\right)~{}.$ (3.14)

Let us note that there is no effect on the cross section of this factor, since
it has unit modulus. The Born series was studied up to third order in the same
reference for the nonrelativistic case, and the expansion of this phase factor
also accounts for its IR divergences. The validity to all orders of this
conjecture was demonstrated by Weinberg [1].

Faddeev and Kulish addressed this issue in Ref. [22] and, in addition to
redefining the asymptotic states, as discussed in Sec. 3.1, a related
redefinition of the $S$ matrix was derived which removes this type of phase
factors. The general form of this phase factor is [1, 22]

$\displaystyle S_{W;f,i}$ $\displaystyle=S_{f}\cdot S_{i}~{},$ (3.15)
$\displaystyle S_{i}$
$\displaystyle=\\!\\!\\!\\!\\!\\!\\!\\!\prod_{\,\,\,\,{\small\begin{tabular}[]{l}$n\neq
m$\cr$n,m\in\ket{i}$\end{tabular}}}\displaystyle{\\!\\!\\!\\!\\!\\!\\!\\!\exp\left(i\frac{e_{n}e_{m}}{4\pi\beta_{nm}}\ln\frac{\lambda}{\Lambda}\right)}~{},$
(3.18) $\displaystyle S_{f}$
$\displaystyle=\\!\\!\\!\\!\\!\\!\\!\\!\prod_{\,\,\,\,{\small\begin{tabular}[]{l}$n\neq
m$\cr$n,m\in\ket{f}$\end{tabular}}}\displaystyle{\\!\\!\\!\\!\\!\\!\\!\\!\exp\left(i\frac{e_{n}e_{m}}{4\pi\beta_{nm}}\ln\frac{\lambda}{\Lambda}\right)}~{}.$
(3.21)

Where $\ket{i}$ and $\ket{f}$ refer to the initial and final states,
respectively, and $\beta_{nm}$ is the Lorentz-invariant relative velocity
between the particles $n$ and $m$,

$\displaystyle\beta_{nm}=\sqrt{1-\frac{m_{n}^{2}m_{m}^{2}}{(p_{n}\cdot
p_{m})^{2}}}~{},$ (3.22)

with $m_{n}$ and $m_{m}$ the masses of the particles. We notice that Ref. [22]
explicitly states that the natural numbers $n$ and $m$ labelling the charged
states in the initial (final) state must be different, as we do here, but this
is not indicated in [1]. This fact is a source of confusion from this latter
reference, because the exponent of the exponential becomes infinite for
$\beta_{nn}=0$. Indeed, the real part of the same loop for $n=m$, both
belonging to the initial or final state, must be accounted for so as to end
with IR-finite transition rates (this loop has a zero imaginary part if
calculated explicitly). In the appendix B we show that the loop integral of
Ref. [1] called $J_{nm}$, whose imaginary part gives rise to the phases in Eq.
(3.15), is purely real for $n=m$. These diagrams correspond to the
contributions from the wave-function renormalization of the electron.

The $S$ matrix for the transition $\alpha\to\beta$ free of IR divergences is
then [22]

$\displaystyle S_{\beta,\alpha}$
$\displaystyle=\widetilde{S}_{\beta,\alpha}S_{W;\beta,\alpha}^{-1}~{},$ (3.23)

with $\widetilde{S}_{\beta\alpha}$ the $S$ matrix that results from the QED
calculation with the asymptotic states dressed by the cloud of real soft
photons, c.f. Sec. 3.1. In our present study we employ the IR free $S$ matrix
$S_{{\boldsymbol{p}}^{\prime}s^{\prime},{\boldsymbol{p}}s}$. The correct phase
factor is Eq. (3.14), and it can be obtained from Eq. (3.15) by taking that
the source of the external potential is a infinitely heavy particle of charge
$e^{\prime}$, so that $\beta=p/E(p)$.

The IR divergent contribution of the iterated diagram in Eq. (2.4) for
$\theta\neq 0$ is

$\displaystyle
i(ee^{\prime})^{2}\bar{u}_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\gamma^{0}u_{\sigma}({\boldsymbol{p}})\frac{E(p)\ln
2p/\lambda}{8\pi p^{3}\sin^{2}\frac{\theta}{2}}~{}.$ (3.24)

The IR divergent term proportional to $\ln\lambda$ is equal to $A_{0}$, Eq.
(2.4), times the phase factor in Eq. (3.14) expanded up to ${\cal
O}(ee^{\prime})$ and removing the one [3].

##  4 Partial-wave amplitudes

The expansion of a scattering amplitude in PWAs is familiar from elementary
Quantum Mechanics. This is a technique that typically simplifies the analysis
of a particle interaction with a potential of finite range, as unitarity is
imposed more simply. The key idea is to take advantage of rotational symmetry
and decompose the scattering amplitude into PWAs, each characterized by having
a well defined total-angular momentum quantum number $J$.

In our case, the interaction is of infinite range, which makes that the PWA
expansion for the Coulomb interaction does not converge.444Despite this lack
of convergence, the PWA expansion is typically used in nuclear physics to
disentangle the effects of Coulomb scattering because of the interference with
strong interactions [19]. The point here is to evaluate separately those
low-$\ell$ PWAs affected by the strong interactions, and then to resum the
rest of Coulomb PWAs by knowing the Coulomb scattering amplitude without
strong interactions. In this case, the Lehmann ellipse collapses because the
coalescence of the LHC and RHC to a point at $|{\boldsymbol{p}}|=0$ [44, 45,
46, 47]. Nonetheless, as stressed in Ref. [16], due to rotational invariance
the angular momentum operator $\boldsymbol{J}$ and the Hamiltonian $H$
commute, $[H,\boldsymbol{J}]=0$, in this case as well. Therefore, it is
meaningful to study scattering in the basis of states that simultaneously
diagonalize $\boldsymbol{J}$ and $H$, whose scattering amplitudes are the
PWAs.

Additionally, the unitarity requirement of the $S$ matrix can be expressed in
a particularly simple way via PWAs, which is one of their main appeals.
Related to this, the study of PWAs also allows one to directly address the
spectrum of the theory and find bound states and resonances with the selected
quantum numbers of the PWAs coupled. These points will be treated in detail
below in Secs. 5-7.

We project the scattering amplitudes for an electron scattered by a external
potential in the so-called $\ell sJ$ basis. Every state of this basis is
characterized by having orbital angular momentum $\ell$, spin $s$, and angular
momentum $J$. It is denoted by $|J\mu,\ell s\rangle$, with $\mu$ the third
component of the total angular momentum $\boldsymbol{J}$ along the $z$ axis.
Although it is not indicated in the set of quantum numbers inside a ket, a
$\ell sJ$ state has well defined modulus of momentum $|{\boldsymbol{p}}|$. Of
course, because of rotational invariance and the Wigner-Eckart theorem the
PWAs do not depend on $\mu$.

The partial-wave expansion in the $\ell sJ$ scheme is developed in detail in
Ref. [30], to which we refer for more details, and here we only review briefly
the main points. The one-particle states with well-defined momentum are
defined as

$\displaystyle|{\boldsymbol{p}},\sigma s\rangle$
$\displaystyle=U({\boldsymbol{p}})|\boldsymbol{0},\sigma s\rangle~{},$ (4.1)

with the standard Lorentz transformation directly taking the state at rest to
its final momentum ${\boldsymbol{p}}$. It can be decomposed as

$\displaystyle
U({\boldsymbol{p}})=R({\hat{\boldsymbol{p}}})B_{z}(p)R({\hat{\boldsymbol{p}}})^{\dagger}~{},$
(4.2)

with $B_{z}(p)$ a Lorentz boost along the $z$ axis of velocity
$-{\boldsymbol{p}}/E$, and the rotation $R({\hat{\boldsymbol{p}}})$ takes the
$z$ axis to ${\boldsymbol{p}}$. In terms of two rotations around the $y$ and
$z$ axis, $R_{y}$ and $R_{z}$, it can be expressed as
$R({\hat{\boldsymbol{p}}})=R_{z}(\phi)R_{y}(\theta)$, with $\theta$ and $\phi$
the polar and azimuthal angles of ${\boldsymbol{p}}$, in this order. In
passing, we mention that if
$U({\boldsymbol{p}})=R({\hat{\boldsymbol{p}}})B_{z}(p)$ were used in Eq. (4.1)
the state would have helicity $\sigma$.

Within our normalization the one particle states have the following Lorentz
invariant scalar product,

$\displaystyle\langle{\boldsymbol{p}}^{\prime},\sigma^{\prime}|{\boldsymbol{p}},\sigma\rangle$
$\displaystyle=(2\pi)^{3}2E(p)\delta({\boldsymbol{p}}^{\prime}-{\boldsymbol{p}})\delta_{\sigma^{\prime}\sigma}~{},$
(4.3) $\displaystyle E(p)$
$\displaystyle=\sqrt{{\boldsymbol{p}}^{2}+m^{2}}~{}.$

Let us stress that the composition
$\boldsymbol{{s}}+\boldsymbol{\ell}=\boldsymbol{J}$ is compatible with a
relativistic formalism. This is so because the states

$\displaystyle|\ell m,\sigma s\rangle=\frac{1}{\sqrt{4\pi}}\int
d\hat{{\boldsymbol{p}}}\,Y_{\ell}^{m}(\hat{{\boldsymbol{p}}})|{\boldsymbol{p}},\sigma
s\rangle~{},$ (4.4)

transform under rotations as the direct product of the representations of
angular momentum $\ell$ and $s$ [30]. Therefore, the states with definite
angular momentum are

$\displaystyle|J\mu,\ell s\rangle$
$\displaystyle=\sum_{\sigma,m}C(m\sigma\mu\mid\ell sJ)\left|\ell m,\sigma
s\right\rangle=\frac{1}{\sqrt{4\pi}}\sum_{m=-\ell}^{\ell}\sum_{\sigma=-s}^{s}C(m\sigma\mu\mid\ell
sJ)\int
d\hat{{\boldsymbol{p}}}\,Y_{\ell}^{m}(\hat{{\boldsymbol{p}}})\left|{\boldsymbol{p}},\sigma
s\right\rangle~{}.$ (4.5)

where $C(m\sigma\mu\mid lsJ)$ refers to Clebsh-Gordon coefficients of the
composition of spins $\ell$ and $s$ to give $J$, with third components $m$,
$\sigma$ and $\mu$, respectively. The normalization of the $|J\mu,\ell
s\rangle$ is a consequence of Eqs. (4.3) and (4.5). It corresponds to

$\displaystyle\langle J^{\prime}\mu^{\prime},\ell^{\prime}s|J\mu,\ell
s\rangle$
$\displaystyle=\frac{4\pi^{2}}{p}\delta(E^{\prime}-E)\delta_{\ell^{\prime}\ell}\delta_{\mu^{\prime}\mu}\delta_{J^{\prime}J}~{}.$
(4.6)

When computing the PWAs for scattering processes, we have to project the
initial and final states onto $\ell sJ$ states. This projection involves
integrating over the solid angle to extract the angular momentum components of
the initial and final states, according to Eq. (4.5). Employing rotational
invariance, one can get rid of the solid angle integration over the initial
state and to keep its momentum fixed along the $z$ axis [30]. One can
furthermore take an extra rotation to restrict the solid angle integration
over the final state to a single angle integration over $\theta\in[0,\pi]$, by
restricting the final momentum to belong to the $xz$ plane. This last step is
explained in Appendix C. Restricting ourselves to the case of interest for the
calculation of the PWAs of an electron scattered by an external potential,
$T^{(J,\ell)}$, the final equation is

$\displaystyle T^{(J,\ell)}\equiv\langle J\mu,\ell s|T|J\mu,\ell s\rangle$
$\displaystyle=\frac{1}{2}\int_{-1}^{1}d\cos\theta\,P_{\ell}(\cos\theta)\langle{\boldsymbol{p}}_{xz},\sigma^{\prime}=\frac{1}{2}s|T|p\hat{\boldsymbol{z}},\sigma=\frac{1}{2}s\rangle$
(4.7)
$\displaystyle+\frac{(-1)^{\ell+\frac{1}{2}-J}}{2(J+\frac{1}{2})}\int_{-1}^{1}d\cos\theta\,P_{\ell}^{1}(\cos\theta)\,\langle{\boldsymbol{p}}_{xz},\sigma^{\prime}=-\frac{1}{2}s|T|p\hat{\boldsymbol{z}},\sigma=\frac{1}{2}s\rangle~{},$

and $s=1/2$. Here we are also denoting the PWAs as $\langle J\mu,\ell
s|T|J\mu,\ell s\rangle$. In turn, the scattering amplitudes between states of
definite linear momentum are designed as
$\langle{\boldsymbol{p}}_{xz},\sigma^{\prime}s|T|p\hat{\boldsymbol{z}},\sigma
s\rangle$, with the momentum of the incoming electron fixed in the
$\hat{\boldsymbol{z}}$ direction. In Eq. (4.7) we have taken into account that
since $J$ is conserved because rotation symmetry, it follows that
$\displaystyle{\ell=J\pm\frac{1}{2}}$, as $s=1/2$. Since the parity of an
electron state of the type $|J\mu,\ell s\rangle$ is $(-1)^{\ell}$ it follows
then that $\ell$ is conserved, and PWAs with different $\ell$ do not mix.

###  4.1 Finite PWAs after redefinition of the $S$ matrix and unitarity

The redefinition of the $S$ matrix in PWAs is a straightforward consequence of
Eq. (3.23) because the phase factor $S_{W}$ does not depend on spin or on the
scattering angle (for one-body scattering in the presence of an external
potential and also for two-body scattering of a closed system). Thus, it
factorizes out in the process of calculating the PWAs, Eq. (4.7).

We denote by $\bar{S}^{(J,\ell)}\equiv 1+2i\rho\bar{T}^{(J,\ell)}$ the $S$
matrix in PWAs that follows from the direct calculation within perturbative
QED up to some order. The phase space factor $\rho$ is fixed by the
normalization of the $|J\mu,\ell s\rangle$ states, Eq. (4.6), by definition is

$\displaystyle\rho=\frac{p}{4\pi}~{}.$ (4.8)

We have indicated by $\bar{T}^{(J,\ell)}$ the PWA between asymptotic
electronic states corrected by the soft photon cloud. The PWAs satisfy
unitarity,

$\displaystyle\Im\bar{T}^{(J,\ell)}$
$\displaystyle=\rho\left|\bar{T}^{(J,\ell)}\right|^{2}$ (4.9)

The redefined IR-finite $S$ matrix in PWAs $S^{(J,\ell)}$ up to ${\cal
O}(\alpha^{2})$ is

$\displaystyle S^{(J,\ell)}$
$\displaystyle=\bar{S}^{(J,\ell)}S_{W}^{-1}=1+2i\rho\bar{T}_{1}^{(J,\ell)}+\frac{iee^{\prime}}{2\pi
v}\ln\frac{\lambda}{\Lambda}+2i\rho\bar{T}_{2}^{(J,\ell)}-\rho\frac{ee^{\prime}}{\pi
v}\ln\frac{\Lambda}{\lambda}\bar{T}_{1}^{(J,\ell)}-\frac{1}{8\pi^{2}}\frac{(ee^{\prime})^{2}}{v^{2}}\ln^{2}\frac{\Lambda}{\lambda}+\mathcal{O}(e^{3})$
(4.10)

In this equation the subscripts indicate the order in the calculation within
the QED expansion,
$\bar{T}^{(J,\ell)}=\sum_{n}\bar{T}^{(J,\ell)}_{n}\alpha^{n}$.

The following step was firstly introduced in Ref. [17] in order to define IR
finite PWAs $T^{(J,\ell)}$ out of the IR-finite $S^{(J,\ell)}$. The procedure
is based on identifying the PWAs free of IR divergences $T^{(J,\ell)}$ from
the relation

$\displaystyle S^{(J,\ell)}$ $\displaystyle=1+2i\rho T^{(J,\ell)}~{}.$ (4.11)

Then, from Eq. (4.10) it follows that the IR-finite PWAs at leading and next-
to-leading order, $T^{(J,\ell)}_{1}$ and $T^{(J,\ell)}_{2}$, respectively, are
given by

$\displaystyle T_{1}^{(J,\ell)}$
$\displaystyle=\bar{T}_{1}^{(J,\ell)}+\frac{ee^{\prime}E}{p^{2}}\ln\frac{\Lambda}{\lambda}~{},$
(4.12) $\displaystyle T_{2}^{(J,\ell)}$
$\displaystyle=\bar{T}_{2}^{(J,\ell)}+i\frac{ee^{\prime}E}{2\pi
p}\ln\frac{\Lambda}{\lambda}\,\bar{T}_{1}^{(J,\ell)}+i\frac{(ee^{\prime})^{2}E^{2}}{4\pi
p^{3}}\ln^{2}\frac{\Lambda}{\lambda}~{}.$ (4.13)

This redefinition only affects the imaginary of $T^{(J,\ell)}_{2}$. This fact
is in agreement with our NLO calculation of the iteration diagram (f) in Fig.
1, whose IR divergences affect only its imaginary, cf. Eq. (2.4).

The redefinition of the $S$-matrix elements $\bar{S}^{(J,\ell)}$ to end with
$S^{(J,\ell)}$ only involves the phase factor $S_{W}$ and unitarity of the $S$
matrix is preserved. Therefore,

$\displaystyle S^{(J,\ell)}{S^{(J,\ell)}}^{*}=I~{}.$ (4.14)

Taking into account the relation between $S^{(J,\ell)}$ and $T^{(J,\ell)}$,
Eq. (4.11), it also follows the unitarity relation for the PWAs,

$\displaystyle\Im T^{(J,\ell)}$
$\displaystyle=\rho\left|T^{(J,\ell)}\right|^{2}~{},~{}p^{2}>0~{}.$ (4.15)

This also settles the unitarity or right-hand cut (RHC) for physical values of
the momentum.

###  4.2 Handling different angular-projection integrals

Let us first notice that for the term with $P_{\ell}^{1}(\cos\theta)$ in Eq.
(4.7) the angular integration is straightforwardly convergent because of the
extra factor of $\sin\theta$ that is part of the generalized Legendre
Polynomials $P_{\ell}^{m}(\theta)$ with $m>0$,

$\displaystyle
P_{\ell}^{m}(\cos\theta)=(-1)^{m}(1-\cos^{2}\theta)^{m/2}\frac{d^{m}}{d\cos^{m}\theta}P_{\ell}(\cos\theta)~{}.$
(4.16)

For this reason, when calculating the partial-wave projection of the iteration
diagram, for the term with $P_{\ell}^{1}(\cos\theta)$ in Eq. (4.7) one may use
instead of Eq. (2.4) a simpler expression given by

$\displaystyle
A_{3}=-(ee^{\prime})^{2}\bar{u}_{\sigma^{\prime}}\left(\boldsymbol{p}^{\prime}\right)$
$\displaystyle\left(E\gamma^{0}\left[\frac{1-\csc\frac{\theta}{2}}{32p^{3}\cos^{2}\frac{\theta}{2}}-i\frac{\ln\frac{2p}{\Lambda}}{8\pi
p^{3}\sin^{2}\frac{\theta}{2}}-i\frac{\left(1+\cos^{2}\frac{\theta}{2}\right)\,\ln\left(\sin\frac{\theta}{2}\right)}{16\pi
p^{3}\sin^{2}\frac{\theta}{2}\,\cos^{2}\frac{\theta}{2}}\right]\right.$ (4.17)
$\displaystyle\left.+m\left[-\frac{1-\csc\frac{\theta}{2}}{32p^{3}\cos^{2}\frac{\theta}{2}}+i\frac{\ln\sin\frac{\theta}{2}}{16\pi
p^{3}\cos^{2}\frac{\theta}{2}}\right]\right)u_{\sigma}(p\boldsymbol{z}),$

Actually, the term proportional to $\ln 2p/\Lambda$ in the previous equation
results by the cancellation of an IR divergent term proportional to $\ln
2p/\lambda$ with another term that stems from the redefinition of the $S$
matrix, from the operation,

$\displaystyle S_{W}^{-1}A_{0}\to i\frac{ee^{\prime}E}{2\pi
p}\ln\left(\frac{\Lambda}{\lambda}\right)A_{0}~{},$ (4.18)

and then $\ln 2p/\lambda\to\ln 2p/\Lambda$ in Eq. (4.17).

Let us consider now the term $P_{\ell}(\cos\theta)$ in Eq. (4.7), which
requires a more careful treatment in order to cancel the $\ln\lambda$ terms.
The IR divergent angular integrals of $A_{3}$ occur for $\theta\to 0$. On the
one hand, one has the terms in Eq. (2.4) that do not drive to IR divergent
integrals, and their angular integrals for the partial-wave projection can be
done in a straightforward manner. Namely, these terms are

$\displaystyle-(ee^{\prime})^{2}\bar{u}_{\sigma^{\prime}}\left({\boldsymbol{p}}^{\prime}\right)\left(E\gamma^{0}\frac{1-\csc\frac{\theta}{2}}{32p^{3}\cos^{2}\frac{\theta}{2}}+m\left[-\frac{1-\csc\frac{\theta}{2}}{32p^{3}\cos^{2}\frac{\theta}{2}}+i\frac{\ln\sin\frac{\theta}{2}}{16\pi
p^{3}\cos^{2}\frac{\theta}{2}}\right]\right)u_{\sigma}(p\boldsymbol{z})~{}.$
(4.19)

On the other hand, the terms in $A_{3}$, Eq. (2.4), that give rise to IR
divergent angular integrals when multiplied by $P_{\ell}(\cos\theta)$ when
taken to Eq. (4.7), are the last two terms in between the square brackets of
this equation, namely:

$\displaystyle A_{3}^{\text{IR}}$
$\displaystyle=-\frac{i(ee^{\prime})^{2}u_{\frac{1}{2}}^{\dagger}\left({\boldsymbol{p}}^{\prime}_{xz}\right)u_{\frac{1}{2}}(p\hat{\boldsymbol{z}})E}{16\pi
p^{3}\cos^{2}\frac{\theta}{2}}\left(\ln\frac{2p}{\lambda}-\frac{\cos^{2}\frac{\theta}{4}}{\sin^{2}\frac{\theta}{2}\,\sqrt{1+\frac{\lambda^{2}}{p^{2}\sin^{2}\frac{\theta}{2}}}}\ln\frac{\sqrt{1+\frac{\lambda^{2}}{p^{2}\sin^{2}\theta/2}}+1}{\sqrt{1+\frac{\lambda^{2}}{p^{2}\sin^{2}\theta/2}}-1}\right).$
(4.20)

The angular integration of this term times $P_{\ell}(\cos\theta)$ is not
algebraic, but we need an algebraic expression in order to handle with the
removal of the IR divergence factor $\ln\lambda$ after including terms from
the $S$-matrix redefinition. In order to accomplish this purpose, we split the
angular integral over $\theta\in[0,\pi]$ in two intervals. The first one is
from 0 up to $\delta$, with $1\gg\delta\gg\lambda^{2}/p^{2}$, and the other
from $[\delta,\pi]$, such that power-law dependence on $\delta$ are neglected
since one takes at the end $\delta\to 0^{+}$.

For the former integration, one can simplify the angular integral in Eq. (4.7)
as

$\displaystyle\frac{1}{2}\int_{0}^{\delta}d\theta\sin\theta
A_{3}^{\text{IR}}P_{\ell}(\cos\theta)$
$\displaystyle=i(ee^{\prime})^{2}E^{2}\int_{0}^{\delta}\frac{\theta
d\theta}{2\pi
p^{3}\theta^{2}\sqrt{1+\frac{4\lambda^{2}}{p^{2}\theta^{2}}}}\ln\frac{1+\sqrt{1+\frac{4\lambda^{2}}{p^{2}\theta^{2}}}}{\frac{2\lambda}{p\theta}}+\ldots=\frac{i(ee^{\prime})^{2}E^{2}}{4\pi
p^{3}}\ln^{2}\frac{\lambda}{p\delta}+\ldots$ (4.21)

where the ellipsis involve terms suppressed by powers of $\delta$. In the
previous equation we have taken into the matrix elements in Eq. (C.16) from
the Dirac algebra, and that $P_{\ell}(\cos\theta)\to 1$ for $\theta\ll 1$, so
that this result is valid for all partial waves. Notice that the term
proportional to $\ln(2p/\lambda)$ in Eq. (4.20) is not included in Eq. (4.21)
because its integral vanishes linearly in $\delta$ for $\delta\to 0$.

For the integral with $\theta\in[\delta,\pi]$, since
$\delta\gg\lambda^{2}/p^{2}$, we can expand $A_{3}^{\text{IR}}$ in powers of
$\lambda/p^{2}$. The integral gets simpler, as

$\displaystyle\frac{1}{2}\int_{\delta}^{\pi}d\theta\sin\theta
A_{3}^{\text{IR}}P_{\ell}(\cos\theta)$
$\displaystyle=\frac{i(ee^{\prime})^{2}E}{4\pi
p^{3}}\int_{\delta}^{\pi}d\theta\cot\frac{\theta}{2}\left(E+m+(E-m)\cos\theta\right)\,P_{\ell}(\cos\theta)\left(\ln\frac{2p}{\lambda}+\frac{\left(1+\cos^{2}\frac{\theta}{2}\right)\ln\sin\frac{\theta}{2}}{2\cos^{2}\frac{\theta}{2}}\right).$
(4.22)

Of course, when summing the results of the integrals in Eq. (4.21) and (4.22)
the dependence on $\delta$ vanishes for $\delta\to 0$. For instance, for
$\ell=0$ the result from the sum of these two contributions is

$\displaystyle A_{3}^{\text{IR};\ell=0}$
$\displaystyle=\frac{i(ee^{\prime})^{2}E}{4\pi
p^{3}}\left[E\ln^{2}\lambda-\left(m+E(-1+\ln 4)+2E\ln
p\right)\ln\lambda+\frac{1}{24}\left(-m(6+\pi^{2})+6E+24E\ln^{2}2p-24(E-m)\ln
2p\right)\right],$ (4.23)

and the $\ln\delta$ dependence has disappeared.

Additionally, the $\ln\lambda$ dependence present in Eq. (4.21) and (4.22)
precisely cancels with the one from the $S$-matrix redefinition. Again for
illustration let us take $\ell=0$. The IR divergent piece from the $S$-matrix
redefinition, cf. Eq. (4.13), is

$\displaystyle i\frac{ee^{\prime}E}{2\pi
p}\ln\frac{\Lambda}{\lambda}\bar{T}_{1}^{(\frac{1}{2},0)}+i\frac{(ee^{\prime})^{2}E^{2}}{4\pi
p^{3}}\ln^{2}\frac{\Lambda}{\lambda}$
$\displaystyle=\frac{i(ee^{\prime})^{2}E}{4\pi
p^{3}}\left[-E\ln^{2}\lambda+\left(m+E\left(-1+\ln 4\right)+2E\ln
2p\right)\ln\lambda\right.$
$\displaystyle\left.+\left(E-m+E\left(\ln\Lambda-2\ln
2p\right)\right)\ln\Lambda\right]$ (4.24)

and its dependence on $\ln\lambda$ cancels with that from
$A_{3}^{\text{IR};\ell=0}$ in Eq. (4.23), such that the combination $\ln
2p/\Lambda$ is the one remaining.

The calculation of the angular integrals in Eq. (4.7) for the partial-wave
projection of the NLO QED amplitudes due to the vacuum polarization, $A_{1}$,
and vertex modification, $A_{2}$, can be done numerically in a straightforward
way. Nonetheless, for small momenta $|{\boldsymbol{p}}|\ll m$ one can derive
algebraic expressions by making an expansion in powers of $t$ of $A_{1}$ and
$A_{2}$ in Eqs. (2.5) and (2.3), respectively. In terms of this expansion the
angular integrals of Eq. (4.7) can be done algebraically.

For instance, applying this procedure to the vacuum polarization diagram for
$\ell=0$ we have

$\displaystyle
A_{1}^{\ell=0}=-ee^{\prime}\alpha\frac{\left(E+m\right)\left(8E^{4}+16E^{3}m-59E^{2}m^{2}-43Em^{3}+267m^{4}\right)}{2835\pi
m^{6}}+{\cal O}\left((p/m)^{6}\right)~{},$ (4.25)

and for the vertex-modification diagram the result is

$\displaystyle A_{2}^{\ell=0}$
$\displaystyle=\frac{ee^{\prime}\alpha}{20160\pi
m^{6}}\left(-96E^{5}+80E^{3}m^{2}+4216Em^{4}+m^{5}(832-945\pi)+20E^{2}m^{3}(68+21\pi)-E^{4}m(512+105\pi)\right.$
(4.26)
$\displaystyle\left.+64(E+m)(8E^{4}+16E^{3}m-53E^{2}m^{2}-37Em^{3}+171m^{4})\ln\frac{m}{2\Lambda}\right)+{\cal
O}\left((p/m)^{6}\right)~{}.$

The partial-wave projection of the Born amplitude of an electron by a Coulomb
potential was already discussed in Refs. [16, 17, 18] for the nonrelativistic
case. The situation is analogous for the LO QED scattering amplitude $A_{0}$,
Eq. (2.4). Its PWA projection according to Eq. (4.7) can be evaluated
algebraically. It is IR divergent, but these divergences are cancelled out
when adding the ${\cal O}(ee^{\prime})$ term after the redefinition the $S$
matrix, following Eq. (4.12).

##  5 Unitarization methods

Partial wave unitarity in quantum field theory serves as a cornerstone
principle, guaranteeing that the combined probability of all conceivable
outcomes in a scattering event remains constrained to 1. This foundational
concept is essential for upholding the coherence of quantum field theories.

For a bound states $|p|\sim m\alpha\ll m$, and this makes that the chain of
unitarity diagrams comprising the Lippmann-Schwinger or Schrödinger equations
are very much enhanced because of unitarity,

$\displaystyle\Im
T^{(J,\ell)}=\frac{p}{4\pi}|T^{(J,\ell)}|^{2}~{},~{}p^{2}>0~{},$ (5.1)

as deduced in Eq. (4.15). Since $T^{(J,\ell)}$ at leading order in powers of
$\alpha$ goes like $p^{-2}$ for $|p|\ll m$, it is clear that then unitarity
produces ever increasing powers of $p^{-1}$, which correspondingly require
higher powers of $m$ in the numerator so as to keep right the dimensions of
the amplitude resulting. Of course, this is to be expected since the
generation of bound states requires non-perturbative physics. This enhancement
can also be seen by analyzing the pinch singularity in the unitarity diagrams,
similarly as originally discussed in [48] for nucleon-nucleon scattering
(there one has the Deuteron among other non-perturbative effects). This same
effect was also unveiled for $\bar{K}N$ scattering and the $\Lambda(1405)$
together with other non-perturbative phenomena and near resonances in [49, 50,
51, 42].

Next, we briefly summarize the three unitarization methods that we consider
along the manuscript: the $g$-method, the Inverse Amplitude Method (IAM) and
the first-iterated $N/D$ method. We apply and compare the results obtained
between them. For scattering along physical values of momentum, that is, real
and positive ones, the $g$-method and the 1st-iterated $N/D$ method are
mutually compatible. For the IAM this is usually also the case, except for
some qualifications at relatively large momentum, as we discussed in detail
below. However, for the study of the fundamental state of the hydrogen atom,
due to the overlap between the pole position and the left-hand cut (LHC), only
the first-iterated $N/D$ method is the one with the analytical structure sound
enough to be applied in such demanding circumstance.

Common to all these three unitarization methods is to use that for the inverse
of a PWA the unitarity relation in partial waves, Eq. (5.1), implies

$\displaystyle\Im\frac{1}{T^{(J,\ell)}}$
$\displaystyle=-\frac{p}{4\pi}~{},~{}p^{2}>0~{}.$ (5.2)

For more details in the presentation of the different unitarization methods
here employed, including extra references, the interested reader can consult
the text book [30], and the reviews [29, 15].

###  5.1 The $g$-method

Based on Eq. (5.2) one can isolate the branch-point singularity at $p^{2}=0$
due to unitarity by writing

$\displaystyle\frac{1}{T^{(J,\ell)}(p^{2})}$
$\displaystyle=\frac{1}{V^{(J,\ell)}(p^{2})}+g(p^{2})~{},$ (5.3)

with the unitarity loop function $g(p^{2})$ having a RHC for $p^{2}>0$. The
branch point singularity at $p^{2}=0$ also has another source due to the LHC
as photons have null mass. For instance, at NLO QED one has a LHC contribution
for $p^{2}<0$ from iterating the LO QED result along the RHC. This fact poses
some extra difficulties that are properly treated in Sec. 5.3. From Eq. (5.3)
it follows that the unitarized PWA then reads

$\displaystyle T^{(J,\ell)}(p^{2})$
$\displaystyle=\frac{1}{V^{(J,\ell)}(p^{2})^{-1}+g(p^{2})}=\frac{V^{(J,\ell)}(p^{2})}{1+V^{(J,\ell)}(p^{2})g(p^{2})}~{}.$
(5.4)

The unitarity loop function, accounting for the propagation of the electron in
the scattering process, obeys the following once-subtracted dispersion
relation in the complex $p^{2}$ plane

$\displaystyle g(p^{2})$
$\displaystyle=g(0)-\frac{p^{2}}{4\pi^{2}}\int_{0}^{\infty}\frac{dk^{2}k}{(k^{2}-p^{2})k^{2}}=g(0)-i\frac{\sqrt{p^{2}}}{4\pi}~{}.$
(5.5)

Next, for $p^{2}\to 0$ we match the unitarity loop function with its
nonrelativistic counterpart [18], which requires that $g(0)=0$. Thus, our
final formula for $g(p^{2})$ is

$\displaystyle g(p^{2})$ $\displaystyle=-i\frac{\sqrt{p^{2}}}{4\pi}~{}.$ (5.6)

As a matter of fact, it then coincides with its nonrelativistic expression for
all $p^{2}$ [18].

Once we have calculated $g(p^{2})$ we can apply Eq. (5.3) to calculate
$V^{(J,\ell)}(p^{2})$ in a series expansion in powers of $\alpha$,

$\displaystyle V^{(J,\ell)}(p^{2})=\sum_{n=1}V_{n}^{(J,\ell)}(p^{2})~{},$
(5.7)

where $V^{(J,\ell)}_{n}$ is ${\cal O}(\alpha^{n})$. To fix them we match the
perturbative QED PWA and its unitarized formula order by order in an expansion
in powers of $\alpha$. We denote by $T_{n}^{(J,\ell)}(p^{2})$ the ${\cal
O}(\alpha^{n})$ contribution to the scattering amplitude calculated
perturbatively in QED. Then, by expanding Eq. (5.4) up to NLO one has that

$\displaystyle V_{1}^{(J,\ell)}(p^{2})$
$\displaystyle=T_{1}^{(J,\ell)}(p^{2})~{},$ (5.8) $\displaystyle
V_{2}^{(J,\ell)}(p^{2})$
$\displaystyle=T_{2}^{(J,\ell)}(p^{2})+g(p^{2})V_{1}^{(J,\ell)}(p^{2})^{2}~{}.$
(5.9)

###  5.2 The Inverse Amplitude Method

The IAM [26, 27, 52, 28, 53] can be derived by an expansion of the inverse of
a PWA up to NLO,

$\displaystyle\frac{1}{T_{1}^{(J,\ell)}+T_{2}^{(J,\ell)}}$
$\displaystyle=\frac{1}{T_{1}^{(J,\ell)}}-\frac{T_{2}^{(J,\ell)}}{{T_{1}^{(J,\ell)}}^{2}}+{\cal
O}(\alpha)~{}.$ (5.10)

Then, this result is directly inverted such that the unitarized PWA is

$\displaystyle T^{(J,\ell)}$
$\displaystyle=\frac{T_{1}^{(J,\ell)}(p^{2})^{2}}{T_{1}^{(J,\ell)}(p^{2})-T_{2}^{(J,\ell)}(p^{2})}~{}.$
(5.11)

Unitarity is satisfied because of the perturbative unitarity fulfilled by
$T_{2}^{(J,\ell)}(p^{2})$, such that

$\displaystyle\Im\,T_{2}^{(J,\ell)}(p^{2})=\frac{p}{4\pi}T_{1}^{(J,\ell)}(p^{2})^{2}~{},~{}p^{2}>0~{}.$
(5.12)

An issue that raises when applying the IAM for $p\gg m$ is that
$T_{2}^{(J,\ell)}$ decreases more slowly for large $p$ than
$T_{1}^{(J,\ell)}$. While the latter does it as $1/p$, the former is typically
ultraviolet enhanced by logarithmic factors of the type
$\displaystyle{\frac{\ln^{2}(p^{2}/m^{2})}{p}}$, due to the dressed-vertex
contribution, as well as by other less strong ones as
$\displaystyle{\frac{\ln(p^{2}/m^{2})}{p}}$, that stem from the rest of
contributions to $T_{2}^{(J,\ell)}$. As a result, for large $p\gtrsim{\cal
O}(m/\alpha)$, the expansion of
$1/(T_{1}+T_{2})=1/T_{1}-T_{2}/T_{1}^{2}+\ldots$ loses its meaning. This
manifests in the appearance of a marked peak when $T_{1}^{(J,\ell)}-\Re
T_{2}^{(J,\ell)}=0$, wherever these perturbative amplitudes have the same sign
for large enough $p$. We show below some examples in this respect concerning
the PWAS with $(J,\ell)=(\frac{1}{2},0)$ and $(\frac{1}{2},1)$. However, these
peaks do not correspond to any resonance pole in the second Riemann sheet and
are artefacts of the approach.

This behavior is in fact similar to the one that appears when applying the IAM
in Chiral Perturbation Theory for meson-meson scattering at around an Adler
zero [53, 31, 32, 33]. Around this region, the leading $T_{1}$ crosses zero
and the subleading $T_{2}$ is then larger, such that $T_{1}-T_{2}$ becomes
zero at some point in those energies. This translates in a deficient behavior
of the IAM in a region around the Adler zero until $|T_{1}|$ becomes clearly
larger again than $|T_{2}|$. However, there is a difference between this
situation and our problem at hand on Coulomb scattering. Now,
$|T_{1}|<|T_{2}|$ remains so for large enough $p$ up to infinity, as just
discussed, so that a proper behavior for the unitarized PWA by employing IAM
cannot be restored for large momenta when $T_{1}$ and $\Re T_{2}$ have the
same sign.

###  5.3 First-iterated $N/D$ method

To end having unitarized PWAs with the right analytical properties, we apply
the first iterated $N/D$ method [29, 30, 15]. The issue with the $g$-method
and IAM is that they drive to structures in which terms with RHC and LHC are
directly multiplied.555For instance, these are not the right analytical
properties of a triangular loop diagram [54]. These deficiencies manifest when
looking for the poles in the unitarized amplitudes following these methods. In
particular, complex poles with momentum having positive imaginary part and
nonvanishing real part are found in the first or physical Riemann sheet, which
is against the due properties of a PWA, cf. Eq. (8.18) below. The right
secular equation for looking for poles in the PWAs should only have RHC. In
the $N/D$ method they correspond to zeros of the $D$ function, while in the
general theory of integral equations for the Lippmann-Schwinger equation one
has to look for the zeroes of the Fredholm determinant [55, 56].

In the first-iterated $N/D$ method the function ${\cal V}^{(J,\ell)}(p^{2})$,
having only LHC and to be fixed below, is taken as the $N$ function. Then, the
unitarized PWA and the resulting $D$ function, $D^{(J,\ell)}(p^{2})$, are
given by the expressions,

$\displaystyle T^{(J,\ell)}(p^{2})$ $\displaystyle=\frac{{\cal
V}^{(J,\ell)}(p^{2})}{D^{(J,\ell)}(p^{2})}~{},$ (5.13) $\displaystyle
D^{(J,\ell)}(p^{2})$
$\displaystyle=1-\frac{1}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k{\cal
V}^{(J,\ell)}(k^{2})}{k^{2}-p^{2}}~{}.$ (5.14)

Notice that $D^{(J,\ell)}(p^{2})$ has only RHC. Equation (5.14) resembles to
the denominator in the last term to the right of Eq. (5.4), the difference is
that now the analogous to the unitarization kernel $V^{(J,\ell)}(p^{2})$ is
integrated.

The function ${\cal V}^{(J,\ell)}(p^{2})$ is expressed as a series expansion
in powers of $\alpha$ as

$\displaystyle{\cal V}^{(J,\ell)}(p^{2})$
$\displaystyle=\sum_{n=1}^{\infty}{\cal
V}^{(J,\ell)}_{n}(p^{2})\alpha^{n}~{},$ (5.15)

where the ${\cal V}_{n}^{(J,\ell)}$ are determined by requiring the matching
of Eq. (5.13) order by order with the perturbative calculation in QED of the
PWA. Therefore,

$\displaystyle{\cal V}_{1}^{(J,\ell)}(p^{2})$
$\displaystyle=T_{1}^{(J,\ell)}(p^{2})~{},$ (5.16) $\displaystyle{\cal
V}_{2}^{(J,\ell)}(p^{2})$
$\displaystyle=T_{2}^{(J,\ell)}(p^{2})-\frac{T_{1}^{(J,\ell)}(p^{2})}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{kT_{1}^{(J,\ell)}(k^{2})}{k^{2}-p^{2}}~{}.$
(5.17)

We notice that because of perturbative unitarity, cf. Eq. (5.12), ${\cal
V}_{2}^{(J,\ell)}$ is real for $p^{2}>0$, cancelling mutually the imaginary
parts in $T_{2}^{(J,\ell)}$ and the integral in Eq. (5.17).

At this stage Eq. (5.17) is rather formal since further discussions are
necessary because the RHC integration for calculating $D^{(J,\ell)}(p^{2})$ is
divergent and, at least, one subtraction is needed to make it convergent.
Thus, the unitarization of the LO QED amplitude takes the form:

$\displaystyle T^{(J,\ell)}(p^{2})$ $\displaystyle=\frac{{\cal
V}_{1}^{(J,\ell)}(p^{2})}{D_{{\rm LO}}^{(J,\ell)}(p^{2})}~{},$ (5.18)
$\displaystyle D^{(J,\ell)}_{{\rm LO}}(p^{2})$ $\displaystyle=1+C_{\rm
LO}^{(J,\ell)}-\frac{p^{2}+p_{0}^{2}}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k{\cal
V}_{1}^{(J,\ell)}(k^{2})}{(k^{2}-p^{2})(k^{2}+p_{0}^{2})}~{}.$ (5.19)

Within the iterated $N/D$ method we can keep fixed to 1 the constant term for
$p\to 0$ in $D^{(J,\ell)}(p^{2})$. Then, we remove the extra constant term
arising from the real part of the integral for $p\to 0$ by properly choosing
the constant $C_{\rm LO}^{(J,\ell)}$,

$\displaystyle C_{\rm LO}$ $\displaystyle=\lim_{p^{2}\to
0}\frac{p_{0}^{2}+p^{2}}{4\pi^{2}}\mathchoice{{\vbox{\hbox{$\textstyle-$}}\kern-4.86108pt}}{{\vbox{\hbox{$\scriptstyle-$}}\kern-3.25pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-2.29166pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-1.875pt}}\\!\int_{0}^{\infty}dk^{2}\frac{k{\cal
V}_{1}^{(J,\ell)}(k^{2})}{(k^{2}-p^{2})(k^{2}+p_{0}^{2})}~{},$ (5.20)

with
$\mathchoice{{\vbox{\hbox{$\textstyle-$}}\kern-4.86108pt}}{{\vbox{\hbox{$\scriptstyle-$}}\kern-3.25pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-2.29166pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-1.875pt}}\\!\int$
denoting the Cauchy principal value of the integral. As a result the
$p_{0}^{2}$ dependence in $D_{\rm LO}^{(J,\ell)}(p^{2})$ is fully removed.
This is a legitimate process because this reshuffling in the subtraction
constant is finally reabsorbed in the LHC function ${\cal V}^{(J,\ell)}$
during the matching procedure in order to reproduce the perturbative QED PWA
calculated in higher orders. Our results are more accurate the higher the
order in which $T^{(J,\ell)}$ is calculated in QED.

This treatment, namely, to keep 1 the constant term in $D^{(J,\ell)}$ for
$p\to 0$, is also applied to the contributions of $D^{(J,\ell)}(p^{2})$
stemming from the integration of ${\cal V}_{2}^{(J,\ell)}(p^{2})$. Once this
constant term is removed there is only $p_{0}^{2}$ dependence left from the
RHC integration involving the real part of $T_{2}^{(J,\ell)}$ that stems from
the QED diagram representing the first iteration of LO QED amplitude, diagram
(f) in Fig. 1. These are pure relativistic terms which are zero in the
nonrelativistic limit [18], where the NLO QED amplitude is purely
imaginary.666Notice that if the leading nonrelativistic contribution to the
first iteration amplitude is ${\cal O}(\alpha^{2}p^{-3})$ (as required by
unitarity) then its relativistic correction is expected to be an order
$p^{2}/m^{2}$ higher, so that it would start at ${\cal O}(\alpha^{2}p^{-1})$,
as it is indeed the case. Now, the point we want to comment upon here is that
this term, ${\cal V}_{\rm it}^{(J,1/2)}(p^{2})$, is given by
$\alpha^{2}\pi^{2}C_{J\ell}/p$, with $C_{J\ell}$ a PWA dependent constant
($C_{\frac{1}{2}0}=1$), so that its RHC integration contributing to
$D^{(J,\ell)}$ is

$\displaystyle\delta D_{\rm it}^{(J,\ell)}(p^{2})$
$\displaystyle=-\frac{(p^{2}+p_{0}^{2})}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k{\cal
V}_{\rm
it}^{(J,\ell)}(k^{2})}{(k^{2}-p^{2})(k^{2}+p_{0}^{2})}=\frac{\alpha^{2}C_{J\ell}}{4}\ln\frac{-p^{2}}{p_{0}^{2}}~{}.$
(5.21)

For such a case we take $p_{0}^{2}=m^{2}$, since any change of subtraction
point $p_{0}^{2}$ just adds a constant that is again removed by the
redefinition of ${\cal V}^{(J,\ell)}$ at higher orders.777The factor in front
of $\ln(-p^{2}/p_{0}^{2})$ in Eq. (5.21) is fixed by unitarity. This value for
the subtraction point is taken because the scale that settles the relativistic
expansion is the electron mass $m$. We show the variation in the results due
to this $p_{0}^{2}$ dependence left by considering values of
$p_{0}^{2}\lesssim m^{2}$, namely, for $p_{0}^{2}=m^{2}/4$ and $m^{2}$. It is
remarkable that for scattering this variation will be shown to be very small,
or even negligible for momenta clearly below or above $p=m$.

A remark is in order for the unitarization of the other terms in ${\cal
V}_{2}^{(J,\ell)}(p^{2})$ that arise from the vacuum polarization (${\cal
V}_{\rm vc}^{(J,\ell)}$), the dressed vertex (${\cal V}_{\rm vx}^{(J,\ell)}$),
and the principal value part of the integration present in the right-hand side
of Eq. (5.17) (${\cal V}_{D}^{(J,\ell)}$). An interesting fact is that these
functions are smooth in the threshold region, and admit a Taylor expansion in
powers of $p^{2}$ around threshold, with higher orders suppressed by negative
powers of $m^{2}$. Then, for $|p^{2}|\ll m^{2}$ the integration

$\displaystyle-\frac{p^{2}}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k\left[{\cal
V}_{\rm vc}^{(J,\ell)}(k^{2})+{\cal V}_{\rm vx}^{(J,\ell)}(k^{2})+{\cal
V}_{\rm D}^{(J,\ell)}(k^{2})\right]}{(k^{2}-p^{2})k^{2}}$
$\displaystyle\approx-\left[{\cal V}_{\rm vc}^{(J,\ell)}(0)+{\cal V}_{\rm
vx}^{(J,\ell)}(0)+{\cal
V}_{D}^{(J,\ell)}(0)\right]\frac{p^{2}}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k}{k^{2}(k^{2}-p^{2})}$
(5.22) $\displaystyle\to\left[{\cal V}_{\rm vc}^{(J,\ell)}(0)+{\cal V}_{\rm
vx}^{(J,\ell)}(0)+{\cal V}_{D}^{(J,\ell)}(0)\right]g(p^{2})~{}.$

This makes that for the study of the contributions of these terms to the
binding momenta ($p_{b}$) of the hydrogen atom one gets the same results
either from the first-iterated $N/D$ method or from the $g$-method.

Then, we have the following final expression for these contributions to
$D^{(J,\ell)}(p^{2})$,

$\displaystyle\delta D_{\rm vc}^{(J,\ell)}+\delta D_{\rm vx}^{(J,\ell)}+\delta
D_{D}^{(J,\ell)}$
$\displaystyle=-\frac{p^{2}}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k\left[{\cal
V}_{\rm vc}^{(J,\ell)}(k^{2})+{\cal V}_{\rm vx}^{(J,\ell)}(k^{2})+{\cal
V}_{\rm D}^{(J,\ell)}(k^{2})\right]}{(k^{2}-p^{2})k^{2}}~{}.$ (5.23)

All in all, the unitarized PWA obtained from the first-iterated $N/D$ method
applied to the perturbative QED amplitude calculated up to NLO is

$\displaystyle T^{(J,\ell)}(p^{2})$ $\displaystyle=\frac{{\cal
V}_{1}^{(J,\ell)}(p^{2})+{\cal V}_{2}^{(J,\ell)}(p^{2})}{D_{{\rm
LO}}^{(J,\ell)}(p^{2})+\delta D_{{\rm it}}^{(J,\ell)}(p^{2})+\delta D_{{\rm
vc}}^{(J,\ell)}(p^{2})+\delta D_{{\rm vx}}^{(J,\ell)}(p^{2})+\delta
D_{D}^{(J,\ell)}(p^{2})}~{},$ (5.24)

where

$\displaystyle{\cal V}_{2}{(J,\ell)}(p^{2})$
$\displaystyle=T_{2}^{(J,\ell)}-\frac{p^{2}+p_{0}^{2}}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k{\cal
V}_{1}^{(J,\ell)}(k^{2})}{(k^{2}-p^{2})(k^{2}+p_{0}^{2})}+\displaystyle{\lim_{q^{2}\to
0}}\frac{q^{2}+p_{0}^{2}}{4\pi^{2}}\mathchoice{{\vbox{\hbox{$\textstyle-$}}\kern-4.86108pt}}{{\vbox{\hbox{$\scriptstyle-$}}\kern-3.25pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-2.29166pt}}{{\vbox{\hbox{$\scriptscriptstyle-$}}\kern-1.875pt}}\\!\int_{0}^{\infty}dk^{2}\frac{k{\cal
V}_{1}^{(J,\ell)}(k^{2})}{(k^{2}-q^{2})(k^{2}+p_{0}^{2})}~{}.$ (5.25)

The function $D^{(J,\ell)}(p^{2})$ is $D_{\rm LO}^{(J,\ell)}+D_{{\rm
NLO}}^{(J,\ell)}$, with

$\displaystyle D_{{\rm NLO}}^{(J,\ell)}(p^{2})$ $\displaystyle=\delta D_{{\rm
it}}^{(J,\ell)}(p^{2})+\delta D_{{\rm vc}}^{(J,\ell)}(p^{2})+\delta D_{{\rm
vx}}^{(J,\ell)}(p^{2})+\delta D_{D}^{(J,\ell)}(p^{2})~{}.$ (5.26)

To simplify the writing we introduce the following notation: The unitarized
PWAs obtained by applying Eqs. (5.4), (5.11), and (5.13) are denoted by
$T^{(J,\ell)}_{\rm{on}}(p^{2})$, $T^{(J,\ell)}_{\rm{IAM}}(p^{2})$ and
$T^{(J,\ell)}_{\rm{nd}}(p^{2})$, respectively. When needed, the order of the
perturbative QED PWA unitarized is indicated as a subscript preceded by a
semicolon.

##  6 Scattering in partial waves

###  6.1 $J=1/2$ in $S$-wave

Figure 3: $|T^{(\frac{1}{2},0)}(p^{2})|$ as a function of $p$ in units of
$m^{-1}$ and $m$, respectively, in log scale. We show the results from the LO
and NLO $g$-method and 1st-iterated $N/D$ methods, the IAM, and the LO and NLO
QED perturbative PWAs. The types of lines are explained in the keys of the
figure. The IAM result shows a fake resonant structure for large $p$, which
does not correspond to any pole in the 2nd RS.

We show in Fig. 3 the modulus of $T^{(\frac{1}{2},0)}(p^{2})$ as a function of
$p>0$ up to $100m$ in log scale. Energy units are taken such that $m=1$. There
we show both the unitarized PWA with different methods and the perturbative
QED calculations at LO and NLO. The curve corresponding to each case is
indicated by the keys in the figure. Namely, the 1st-iterated $N/D$ method for
unitarizing the LO and NLO QED amplitudes correspond to the magenta and green
lines, respectively. The $g$-method for unitarizing the LO and NLO QED
amplitudes are the light blue and orange lines, in this order. The IAM is the
yellow line, and the LO and NLO QED perturbative PWAs are the blue and red
lines, respectively.

Regarding the perturbative amplitudes we see that they depart from each other
for small momenta, $p\lesssim 0.1~{}m$, and for large ones, $p\gtrsim 6~{}m$.
It is interesting to discuss the reasons involved. For the small momenta the
point is that for $|p|={\cal O}(m\alpha)$ Coulomb scattering becomes
nonperturbative, with $im\alpha$ the binding momentum of the hydrogen atom,
and the perturbative expansion breaks down as $p\to 0$. The interested reader
can consult Fig. 1 of Ref. [18] where this point is explicitly shown for
nonrelativistic scattering. For the large $p$ part the issue there is, as
already pointed out in Sec. 5.2, that $T_{1}^{(\frac{1}{2},0)}(p^{2})$
vanishes as $1/p$ for $p\to\infty$, while $T_{2}^{(\frac{1}{2},0)}$ is
asymptotically dominated in this limit by the contribution from the dressed
electromagnetic vertex of the electron, because despite it also decreases like
$1/p$ its size is relatively enhanced by factors $\ln^{2}(p/m)^{2}$ [coming
from the function $Li_{2}$ in Eq.(3.13)]. This is not a particularity of
$J=1/2$ in $S$-wave, but it applies to any PWA. In this way, for $p\sim
m/\alpha$ this ultraviolet enhancement is large enough to compensate the extra
factor of $\alpha$ of $T_{2}^{(J,\ell)}$ relative to $T_{1}^{(J,\ell)}$, and
$T_{2}^{(J,\ell)}$ starts becoming larger in absolute value than
$|T_{1}^{(J,\ell)}|$ as $p$ increases.

Concerning the unitarization methods we would like to stress two points: i)
The different unitarization results follow closely the perturbative QED
amplitudes, such that they cannot typically be distinguished in the figure,
except for small $p\lesssim 5m\alpha$. It is worth noting that there the
results from the three unitarization methods overlap each other, both by
unitarizing the LO and NLO perturbative QED amplitudes, and they cannot be
distinguished in the scale of the figure. However, the blue and red lines
corresponding to the LO and NLO QED perturbative amplitudes, respectively, can
be distinguished between them and with respect to the unitarized results for
the small momentum region. ii) The figure clearly shows the problem of the IAM
at large momenta, as already anticipated in Sec. 5.2. The peak occurs at
$p\simeq 38.5~{}m$ because of the cancelation of $T_{1}^{(\frac{1}{2},0)}-\Re
T_{2}^{(\frac{1}{2},0)}$ at this crossover point. The peak in Fig. 3 does not
correspond to any resonance manifested as a pole in the complex $p^{2}$ plane
in the second Riemann sheet, but it is just an artifact of the method.

Figure 4: We plot in log scale $R_{i}$ defined in Eq. (6.1) as a function of
$p$ (in units of $m$) for the unitarization with the 1st-iterated $N/D$ method
of the LO and NLO QED amplitudes. For the latter case we plot the results for
$p_{0}^{2}=m^{2}$ and $m^{2}/4$. The types of lines are explained in the keys
of the figure.

In order to ascertain the impact of the unitarization of the perturbative QED
PWAs $T^{(\frac{1}{2},0)}_{i}(p^{2})$, we define the ratio $R_{i}^{(J,\ell)}$
as

$\displaystyle R^{(J,\ell)}_{i}$ $\displaystyle=\left|\frac{T^{(J,\ell)}_{{\rm
nd};i}}{\displaystyle{\sum_{j=1}^{i}T^{(J,\ell)}_{j}}}-1\right|~{}.$ (6.1)

This is shown in Fig. 4 for the unitarization of the LO ($i=1$) and NLO
($i=2$) perturbative QED $S$-wave by employing the first-iterated $N/D$
method. We also show the results for two values of the subtraction point
$p_{0}^{2}$ at $m^{2}$ and $m^{2}/4$. Unitarization corrections to the LO QED
PWA at typical momenta $p={\cal O}(m)$ are ${\cal O}(\alpha^{2})$ in absolute
terms, and ${\cal O}(\alpha)$ in relative ones. This sets the size of the
magenta line in Fig. 4 for momenta $p={\cal O}(m)$. For the unitarization of
the NLO QED PWA the corrections with respect to the NLO QED PWA are of ${\cal
O}(\alpha^{3})$ and ${\cal O}(\alpha^{2})$ in absolute and relative terms,
respectively. This sets the size of the unitarization corrections for the blue
and orange lines in Fig. 4 for typical momenta of ${\cal O}(m)$. In more
detail, for the NLO case we can write

$\displaystyle R_{2}$ $\displaystyle=\left|\frac{{\cal V}_{1}^{(J,\ell)}+{\cal
V}_{2}^{(J,\ell)}-(T_{1}^{(J,\ell)}+T_{2}^{(J,\ell)})D^{(J,\ell)}}{(T_{1}^{(J,\ell)}+T_{2}^{(J,\ell)})D^{(J,\ell)}}\right|={\cal
O}(\alpha^{3})~{},$ (6.2)

because of the matching procedure in Eq. (5.14).

However, the corrections are much larger for $p={\cal O}(m\alpha)$ because QED
is nonperturbative there and, therefore, the largest corrections are precisely
expected at such momenta.888Notice that in addition to the nominal global
power in $\alpha$ corresponding to $T_{1}^{(J,\ell)}$ and $T_{2}^{(J,\ell)}$,
the momenta themselves carry powers of $\alpha$ when $p={\cal O}(m\alpha)$.
On the other side, for large momenta $p={\cal O}(m/\alpha)$, the unitarity
corrections to the NLO QED PWAs grow and become comparable to the corrections
of the LO QED PWA. This is due to the aforementioned faster decrease for
$p\to\infty$ of the LO amplitude at higher momentum compared with the NLO QED
one. It is also worth stressing that the left dependence of the results on
$p_{0}^{2}$ due to the unitarization of ${\cal V}_{\rm it}^{(\frac{1}{2},0)}$
is already barely visible in Fig. 4, and it can be neglected in front of the
rest of contributions.

###  6.2 $J=3/2$ and $J=1/2$ in $P$-waves

For $P$-wave scattering there are two possible values for the total angular
momentum, $J=3/2$ and $J=1/2$. We proceed along similar lines as for the
$S$-wave. First, we show in Fig. 5 the modulus of $T^{(J,1)}(p^{2})$ for
$J=3/2$ (left panel) and $J=1/2$ (right panel), as resulting from the
application of the unitarization methods and the pure perturbative
calculations. We follow the same convention for the type of lines as in Fig.
3. The basic features that one can observe in the Fig. 5 for the $P$ waves are
similar to those for the $S$ wave. The $g$-method and the 1st-iterated $N/D$
method provide unitarized results that cannot be distinguished from the
perturbative PWAs in the scale of the figures. The LO and NLO QED PWAs can
only be distinguished in the large momentum region because of the faster
vanishing for large of $p$ of the LO result. For $J=1/2$ the difference
between the red and blue lines is more visible in the log scale used in the
figure because the displacement in the position of the zero of the LO QED PWA
at $6.4m$ and the one in the real part of the NLO QED PWA at $5.4m$.999Let us
indicate that its imaginary part does not vanish because it is given by
unitarity, and $T_{1}^{(\frac{1}{2},1)}$ is not zero at $p\approx 5.4m$ but at
$p\approx 6.4m$. It is however very small because $T_{1}^{(\frac{1}{2},1)}$ is
suppressed due to its nearby zero. The application of the IAM in the $J=1/2$
case is handicapped because $T_{1}^{(\frac{1}{2},0)}-\Re
T_{2}^{(\frac{1}{2},0)}$ vanishes at two points with momenta equal to $8.7m$
and $29.6m$. This is not the case for $J=3/2$ since $T_{1}^{(\frac{3}{2},1)}$
and $T_{2}^{(\frac{3}{2},1)}$ have opposite signs, and
$T_{2}^{(\frac{3}{2},1)}$ remains very small in comparison with
$T_{1}^{(\frac{3}{2},1)}$ in the range of momentum shown.

|
---|---
Figure 5: $|T^{(3/2,1)}(p^{2})|$ (left panel) and $|T^{(1/2,1)}(p^{2})|$ (right panel) are plotted in log scale as a function of $p$ in units of $m^{-1}$ and $m$, respectively. We show the LO and NLO $g$\- and 1st-iterated $N/D$ methods, the IAM and the LO and NLO QED perturbative PWAs. The types of lines are explained in the keys of the figure, being the same as in Fig. 3.  |
---|---
Figure 6: We plot in log scale $R_{i}^{(3/2,1)}$ (left panel) and
$R_{i}^{(1/2,1)}$ (right panel) as functions of $p$ (in units of $m$) for the
unitarization with the 1st-iterated $N/D$ method of the LO and NLO QED
amplitudes. For every case we plot the results for two subtractions points
$p_{0}^{2}=m^{2}$ and $m^{2}/4$. The types of lines are explained in the keys
of the figure, being the same as in Fig. 4.

The ratios $R_{i}^{(J,1)}$, cf. Eq. (6.1), are shown in Fig. 6, with the lines
corresponding to the same casuistry as in Fig. 4. From Fig. 6 we see that the
relative size of the unitarization corrections for typical momenta of ${\cal
O}(m)$ are as expected, namely, of ${\cal O}(\alpha)$ and ${\cal
O}(\alpha^{2})$ for the unitarization of the LO and NLO QED PWAs, in this
order. The peak of the relative unitarization corrections in the $J=1/2$ NLO
QED PWA just reflects the almost zero of $T_{2}^{(\frac{1}{2},1)}$ at
$p\approx 5.4m$. However, $R^{(\frac{1}{2},1)}_{1}$ is smooth when
$T_{1}^{(\frac{1}{2},1)}$ crosses its zero at $p\approx 5.4m$ because ${\cal
V}_{1}^{(J,\ell)}=T_{1}^{(J,\ell)}$, and the zero at $T_{1}^{(\frac{1}{2},1)}$
drops in the ratio between the unitarized and the perturbative LO QED PWA. As
in the $S$-wave case, the change in the results from the two values of
$p_{0}^{2}$ is very small and it can only be disentangled in the figure for
$p$ close to $m$.

##  7 The bound-state equation for the different unitarization methods

For nonrelativistic Coulomb scattering one knows the full $S$ matrix. In terms
of it Ref. [18] gave the exact expression for $V_{\rm nr}^{(\ell)}(p)$, the
nonrelativistic counterpart of $V^{(J,\ell)}$ of Eq. (5.3), such that

$\displaystyle V_{\rm nr}^{(\ell)}(p)$
$\displaystyle=-\frac{4i\pi}{p}\frac{\Gamma(1+\ell-i\frac{m\alpha}{p})-\Gamma(1+\ell+i\frac{m\alpha}{p})}{\Gamma(1+\ell-i\frac{m\alpha}{p})+\Gamma(1+\ell+i\frac{m\alpha}{p})}~{},$
(7.1)

where we have suppressed the superscript $J$ because in the nonrelativistic
case $V^{(J,\ell)}$ only depends on $\ell$.

The expansion of $V_{\rm nr}^{(\ell)}(p)$ of Eq. (7.1) in powers of $\alpha$
breaks down for $|p|<p_{\rm critical}$, with $p_{\rm critical}$ the largest
momentum at which the denominator in Eq. (7.1) vanishes. The values of such
critical momenta are given in Table 1 for small $\ell\leq 5$. For large values
of $\ell$ an approximate expression for $p_{\rm critical}$ can be obtained by
expanding the denominator in Eq. (7.1) in powers of $\ell^{-1}$. Its
asymptotic behavior for $\ell\gg 1$ is

$\ell$ | 0 | 1 | 2 | 3 | 4 | 5
---|---|---|---|---|---|---
$\frac{p_{\rm critical}}{m\alpha}$ | 0.18 | 0.43 | 0.62 | 0.81 | 0.96 | 1.09
Table 1: Largest root of the denominator in Eq. (7.1) for $\ell\leq 5$.

$\displaystyle\Gamma(1+\ell-i\frac{m\alpha}{p})+\Gamma(1+\ell+i\frac{m\alpha}{p})\longrightarrow\sqrt{2\pi
l}\exp\left(\ell[-1+\ln\ell]\right)\ell^{-\frac{im\alpha}{p}}\left(1+\ell^{\frac{2im\alpha}{p}}\right)~{},$
(7.2)

from where we deduce that

$\displaystyle p_{\rm critical}\to\frac{2m\alpha\ln\ell}{\pi}~{}.$ (7.3)

Let us notice that for $\ell=5$ this formula gives the approximate result of
$p_{\rm critical}\approx 1.02m\alpha$, already rather close to the numerical
result in Table 1. We should compare $p_{\rm critical}$ with the modulus of
the binding momentum of the nonrelativistic hydrogen atom with principal
quantum number $n$, $im\alpha/n$, such that the angular momentum $\ell$ is
within the range $0\leq\ell\leq n-1$. Because the binding momentum decreases
in modulus with $\ell$ while $p_{\rm critical}$ increases, as it is clear from
Table 1 and Eq. (7.3), it follows that $|p_{b}|>p_{\rm critical}$ only for
$\ell<2$.101010In this statement we take into account that QED corrections to
the nonrelativistic Quantum Mechanics values of the binding momentum are
suppressed by extra powers of $\alpha$ and are small for the physical value of
the fine structure constant.

For $P$ wave one needs at least $n=2$, and $1/2$ is quite close to the radius
of convergence $p_{\rm critical}/m\alpha=0.43$ for $\ell=1$ in Table 1, such
that one really needs high orders in the expansion of $V^{(1)}_{\rm nr}$. It
is straightforward to check from the expansion of Eq. (7.1) that $V^{(1)}_{\rm
nr}$ becomes negative for $p^{2}<0$ and $|p|\ll m$ (and then able to give rise
to a bound state) only when including contributions of ${\cal O}(\alpha^{5})$,
though higher orders are still necessary to reach convergence. This would
require more than 4 loops in QED, certainly beyond the scope of this study
based on unitarized NLO QED. For $S$ wave, in principle, one could study
$n\leq 5$ as $p_{\rm critical}=0.18$, cf. Table 1. However, the unitarized LO
QED amplitude only gives rise to the fundamental state of the hydrogen atom,
and higher-order QED calculations are needed to reach accurateness in the
result. Indeed, $p_{b}$ for the fundamental state is reproduced with a
precision of a few $\permil$ only after including the ${\cal O}(\alpha^{3})$
contribution in the expansion of $V_{{\rm nr}}^{(0)}$ (see Table 1 of Ref.
[18]). The calculation is more and more demanding as $n$ increases.
Nevertheless, it is a compelling prospect for future studies the possibility
of studying within this new formalism the Lamb shift effect, corresponding to
the energy difference between the ${}^{2}S_{1/2}$ and ${}^{2}P_{1/2}$ energy
levels by unitarizing high-loop calculations in QED. In this respect, a two-
loop calculation would be of interest since it is at this order when the self-
energy of the electron inside the iteration diagram first appears. This is a
key effect for the Lamb shift, as it is well known [57, 58, 38].

Next, we explain a procedure to solve the bound-state equation for the
$g$-method and the IAM. The $g$-method can be applied to unitarize the LO and
NLO QED PWAs, while the IAM for being applied requires at least the NLO QED
PWAs. We will not treat separately the bound-state equation from the IAM
because, as shown below, the bound-state equation for the $g$-method with a
perturbative treatment of the NLO contributions encompasses that of the IAM.

The bound-state equation from the knowledge of the LO QED PWAs is

$\displaystyle\frac{1}{V_{1}^{(J,\ell)}(p^{2})}-i\frac{p}{4\pi}=0~{}.$ (7.4)

In order to handle with transcendental equations we typically proceed
perturbatively expanding the solution in powers of $\alpha$, being the
nonrelativistic result the leading contribution. This process is illustrated
below with explicit examples.

The total result at NLO QED contains the vacuum polarization, vertex diagrams
plus the Chung’s diagrams, and the once-time iteration, that also includes the
redefinition of the $S$ matrix up to ${\cal O}(\alpha^{2})$ to make it IR
finite. Since the NLO QED result is ${\cal O}(\alpha^{2})$, its contribution
to $V^{(J,\ell)}(p^{2})$ can be calculated by applying Eq. (5.9) and is also
${\cal O}(\alpha^{2})$. The partial-wave projection of the vacuum polarization
and dressed vertex diagram, including Chung’s diagrams, has no algebraic
expression, though it can be done numerically. However, their expansion at low
momentum can be algebraically worked out order by order in a power series in
$p^{2}$, which is precisely what we need in order to implement a perturbative
treatment for solving the bound-state equation. All these terms, $V_{{\rm
vp}}^{(J,\ell)}$ and $V_{\rm vx}^{(J,\ell)}$ start at ${\cal
O}(\alpha^{2}p^{0})$, cf. Eq. (5.22). In contrast, the LO QED amplitude is
${\cal O}(\alpha/p^{2})$ and for $p={\cal O}(m\alpha)$, as it is the case for
the binding momentum, it counts as ${\cal O}(\alpha^{-1})$. Therefore, the
contributions from the vacuum polarization and dressed vertex diagram to the
binding momentum count three extra powers of $\alpha$, and they are ${\cal
O}(\alpha^{4})$. Higher order terms in the low-momentum expansion would
contribute with extra powers of $\alpha^{2}$, giving rise to ${\cal
O}(\alpha^{6})$ contributions that we do not discuss.

The once-iterated diagram $V_{\rm it}^{(J,\ell)}$ depends in momentum as
${\cal O}(\alpha^{2}p^{-1})$. This is a pure relativistic contribution because
the leading nonrelativistic contribution to the first-iterated scattering
amplitude is purely imaginary for real $p$ and it is removed when calculating
$V_{2}^{(J,\ell)}$, as shown in [18]. For momentum of ${\cal O}(m\alpha)$ it
then counts as ${\cal O}(\alpha)$, two powers higher than the LO QED amplitude
that counts as ${\cal O}(\alpha^{-1})$ in the same momentum region. Therefore,
the contribution to the binding momentum from the first-iterated diagram is
${\cal O}(\alpha^{3})$. Thus, from the unitarization of the NLO QED scattering
amplitude, the leading corrections to the nonrelativistic result for the
binding momentum are ${\cal O}(\alpha^{3})$, and stem from the relativistic
corrections to the nonrelativistic LO amplitude and from the once-iterated
amplitude.

Since, as discussed, $V_{2}^{(J,\ell)}$ contributes perturbatively to the
binding momentum, the bound-state equation Eq. (7.4) with the NLO QED
contribution included can be expanded on $V_{2}^{(J,\ell)}$ and it becomes
then

$\displaystyle\frac{1}{V_{1}^{(J,\ell)}(p^{2})+V_{2}^{(J,\ell)}(p^{2})}-i\frac{p}{4\pi}=\frac{1}{V_{1}^{(J,\ell)}}-\frac{V_{2}^{(J,\ell)}}{(V_{1}^{(J,\ell)})^{2}}+{\cal
O}(\alpha)-i\frac{p}{4\pi}=0~{}.$ (7.5)

This extra contribution proportional to $V_{2}^{(J,\ell)}$ starts at ${\cal
O}(\alpha^{0})$. If the previous equation is multiplied by
$V_{1}^{(J,\ell)}(p^{2})^{2}$ one realizes, after taking into account Eq.
(5.9) to determine $V_{2}^{(J,\ell)}$ and the expression of $g(p)$ in Eq.
(5.6), that it can also be written as $T_{1}^{(J,\ell)}-T_{2}^{(J,\ell)}=0$.
Thus, Eq. (7.5) is the same bound-state equation as for the IAM case (and this
is why we do not treat it separately).

However, Eq. (7.5) cannot be so straightforwardly applied to the bound-state
problem and it requires some amends. The problem is that the once-iterated
amplitude contributing to $V^{(J,\ell)}_{2}(p^{2})$ is proportional to
$1/\sqrt{p^{2}}$, so that it has a LHC for $p^{2}<0$ that overlaps with the
bound-state binding momentum $p_{b}$. As a result, the direct application of
Eq. (7.5) drives to a pole in the first Riemann sheet with $p_{b}$ having a
positive imaginary part but also a non-vanishing real part. This last feature
is unacceptable from the analytical point of view. This issue is solved by the
1st-iterated $N/D$ method that keeps the right analytical structure of a PWA,
as already discussed in Sec. 5.3. The general bound-state equation for
calculating the binding momentum $p_{b}$ within the 1st-iterated $N/D$ method
is

$\displaystyle D^{(J,\ell)}(p_{b}^{2})$ $\displaystyle=0~{}.$ (7.6)

This equation is then applied to $D_{{\rm LO}}^{(J,\ell)}$ and $D_{{\rm
LO}}^{(J,\ell)}+D_{{\rm NLO}}^{(J,\ell)}$, Eqs. (5.19) and (5.26),
respectively. Since $D^{(J,\ell)}(p^{2})$ only has RHC, the problem in the
secular equation due to the overlapping of the binding momentum with the LHC
is avoided by construction.

Finally, it is important to realize that for the bound-state problem not all
the contributions to $V_{n}^{(J,\ell)}$, with $n\geq 3$, can be treated
perturbatively in comparison with $V_{1}^{(J,\ell)}$, as discussed in Ref.
[18]. To show it let us consider the expansion $V_{\rm nr}^{(\ell)}$, Eq.
(7.1) in powers of $\alpha$, namely, $V_{{\rm nr}}^{(\ell)}=\sum_{j=0}V_{{\rm
nr};2j+1}\alpha^{2j+1}$, since it only involves odd powers of $\alpha$.
Therefore, for $p={\cal O}(m\alpha)$ the functions $V_{{\rm
nr;}2j+1}^{(\ell)}(p^{2})$, with $j\geq 1$, are not suppressed compared to
$V_{{\rm nr};1}^{(\ell)}(p^{2})$ (the expansion in powers of $\alpha$ of Eq.
(7.1) is actually an expansion in powers of $\alpha m/p$, the so-called
Sommerfeld parameter).

##  8 Binding momentum for the ground state in the $J=1/2$ and $\ell=0$ PWA

In this section $J=1/2$ and $\ell=0$. Firstly, we analyze the results for the
binding momentum from the unitarized LO QED PWA, and then we move on and
consider the NLO QED PWA.

###  8.1 Bound state from the unitarized LO QED PWA

The LO QED result plus the contribution from the $S$ matrix redefinition, that
removes the infrared divergences, is

$\displaystyle T^{(\frac{1}{2},0)}_{1}(p^{2})$
$\displaystyle=\frac{2\pi\alpha(m-E(1-2\ln a))}{p^{2}}~{}.$ (8.1)

The bound-state equation Eq. (7.4) is now

$\displaystyle\frac{p^{2}}{2\pi\alpha(m-(1-2\ln a)E)}-i\frac{p}{4\pi}=0~{}.$
(8.2)

Its solution for the binding momentum $p_{b}$ is

$\displaystyle p_{b}$ $\displaystyle=im\alpha\frac{1+|1-2\ln
a|\sqrt{1-\alpha^{2}(1-\ln a)\ln a}}{2\left(1-\alpha^{2}(1-\ln a)\ln
a\right)+\alpha^{2}/2}$ (8.3) $\displaystyle=im\alpha\ln
a+\frac{1}{4}im\alpha^{3}(1-2\ln
a)\ln^{2}a+\frac{im\alpha^{5}\ln^{3}a}{16}(2+\ln a(-7+6\ln a))+{\cal
O}(\alpha^{7})~{}.$

In the pure nonrelativistic case the previous result reduces to
$p_{b}=im\alpha\ln a$ [18].

Next, we consider the application of the first-iterated $N/D$ method for
unitarizing $V_{1}^{(\frac{1}{2},0)}$. For that we evaluate algebraically the
function $D^{(\frac{1}{2},0)}(p^{2})$, Eq. (5.19), that explicitly reads

$\displaystyle D^{(\frac{1}{2},0)}_{{\rm LO}}(p^{2})$
$\displaystyle=1+C^{(\frac{1}{2},0)}-\frac{\alpha(p^{2}+p_{0}^{2})}{\pi}\int_{0}^{\infty}dk\frac{m-(1-2\ln
a)\sqrt{m^{2}+k^{2}}}{(k^{2}-p^{2})(k^{2}+p_{0}^{2})}$ (8.4)
$\displaystyle=1+C^{(\frac{1}{2},0)}-i\alpha\frac{m-(1-2\ln
a)E}{2p}-\frac{\alpha(1-2\ln a)}{2\pi p}E\ln\frac{E+p}{E-p}$
$\displaystyle+\frac{\alpha}{2\pi p_{0}}\left(\pi m+2(1-2\ln
a)\sqrt{p_{0}^{2}-m^{2}}{\rm
arctanh}\sqrt{1-\frac{m^{2}}{p_{0}^{2}}}\right)~{}.$

Here $C^{(J,\ell)}$ is a subtraction constant that we expand in powers of
$\alpha$ as

$\displaystyle C^{(J,\ell)}$
$\displaystyle=\sum_{n=1}C_{n}^{(J,\ell)}\alpha^{n}~{},$ (8.5)

where $C_{1}^{(J,\ell)}$ and $C_{2}^{(J,\ell)}$ are going to be fixed in the
unitarization process of the LO and NLO QED amplitudes, respectively, by
following Sec. 5.3. The expansion in powers of $p$ of Eq. (8.4) is

$\displaystyle D^{(\frac{1}{2},0)}_{{\rm LO}}(p^{2})$
$\displaystyle=1-i\frac{m\alpha\ln a}{p}+C^{(\frac{1}{2},0)}+\frac{\alpha}{\pi
p_{0}}\left(\frac{m\pi}{2}+(1-2\ln a)\left[\sqrt{p_{0}^{2}-m^{2}}{\rm
arctanh}\sqrt{1-\frac{m^{2}}{p_{0}^{2}}}-p_{0}\right]\right)$ (8.6)
$\displaystyle+i\frac{\alpha(1-2\ln a)p}{4m}+{\cal O}(p^{2})~{}.$

By requiring that its real part goes to 1 for $p\to 0$ we then have that

$\displaystyle C_{1}^{(\frac{1}{2},0)}$ $\displaystyle=-\frac{1}{\pi
p_{0}}\left(\frac{m\pi}{2}+(1-2\ln a)\left[\sqrt{p_{0}^{2}-m^{2}}{\rm
arctanh}\sqrt{1-\frac{m^{2}}{p_{0}^{2}}}-p_{0}\right]\right)~{}.$ (8.7)

Taking this expression into $D_{{\rm LO}}^{(\frac{1}{2},0)}(p^{2})$ it follows
that

$\displaystyle D_{{\rm LO}}^{(\frac{1}{2},0)}(p^{2})$
$\displaystyle=1-i\alpha\frac{m-(1-2\ln a)E}{2p}+\frac{\alpha(1-2\ln a)}{2\pi
p}\left[2p-E\ln\frac{E+p}{E-p}\right]~{}.$ (8.8)

The $D_{{\rm LO;nr}}^{(\frac{1}{2},0)}$ function calculated by employing the
nonrelativistic version of $V_{1}^{(0,\frac{1}{2})}$, denoted above by
$V_{{\rm nr}}^{(0)}$, is

$\displaystyle D_{{\rm LO;nr}}^{(0)}(p^{2})$
$\displaystyle=1-\frac{2m\alpha\ln
a}{\pi}\int_{0}^{\infty}\frac{dk}{k^{2}-p^{2}}=1-i\frac{m\alpha\ln a}{p}~{},$
(8.9)

in agreement with [18]. Notice that the integration in Eq. (8.9) is convergent
and it does not require of any subtraction. In the limit $p\to 0$ we have from
Eqs. (8.8) and (8.9) that

$\displaystyle D_{{\rm LO}}^{(\frac{1}{2},0)}(p^{2})$ $\displaystyle\to
D_{{\rm LO;nr}}^{(0)}(p^{2})+{\cal O}(p)~{}.$ (8.10)

The bound-state equation $D_{\rm LO}^{(\frac{1}{2},0)}(p_{b}^{2})=0$ is a
transcendental equation and we proceed to solve it in a power series expansion
in $\alpha$. For that, it is convenient to multiply Eq. (8.8) by
$1/V_{1}^{(\frac{1}{2},0)}$, and the resulting secular equation is

$\displaystyle-i\frac{p}{4\pi}+\frac{p^{2}}{2\pi\alpha(m-E[1-2\ln
a])}+\frac{(1-2\ln a)p}{4\pi^{2}(m-E[1-2\ln
a])}\left(2p-E\ln\frac{E+p}{E-p}\right)=0~{}.$ (8.11)

Except for the term $-ip/4\pi$ this equation only contains even powers of $p$
if expanded around $p=0$. The leading order contribution is
$p_{b;1}=im\alpha\ln a$, by constructions the same result as in the
nonrelativistic treatment. Next, to get the ${\cal O}(\alpha^{2})$ correction
we write the binding momentum $p_{b}=p_{1}+\alpha^{2}\kappa_{2}$, and expand
Eq. (8.11) up to ${\cal O}(\alpha^{2})$. One then obtains that $\kappa_{2}=0$.
For the first non-vanishing corrections we write
$p_{b}=p_{1}+\alpha^{3}\kappa_{3}$, expand Eq. (8.11) up to ${\cal
O}(\alpha^{3})$ and the results is

$\displaystyle\kappa_{3}$ $\displaystyle=i\frac{m(1-2\ln a)\ln^{2}a}{4}~{}.$
(8.12)

Therefore,

$\displaystyle p_{b}=im\alpha\ln a+im\alpha^{3}\frac{(1-2\ln
a)\ln^{2}a}{4}+{\cal O}(\alpha^{4})~{}.$ (8.13)

In connection with the Dirac equation for the binding momentum $p_{D}^{(J,n)}$
for the hydrogen atom,

$\displaystyle p^{(J,n)}_{D}$
$\displaystyle=im\alpha\left[\alpha^{2}+\left(n-J-\frac{1}{2}+\sqrt{(J+\frac{1}{2})^{2}-\alpha^{2}}\right)^{2}\right]^{-\frac{1}{2}}~{},$
(8.14)

with $n$ the principal quantum number with $J+\frac{1}{2}\leq n$, we notice
that the iteration of the LO QED amplitude should represent the same physics.
For the case of the lowest energy level of the hydrogen atom $n=1$ and one has
that $p_{D}^{(\frac{1}{2},1)}=im\alpha$, as in the nonrelativistic case. This
result differs from Eq. (8.13) by ${\cal O}(\alpha^{3})$. The difference, of
course, stems from the need to know $V^{(J,\ell)}$ to all orders, similarly as
for the nonrelativistic case, treated in detail in Ref. [18]. There, it was
shown that once higher orders terms in the expansion in powers of $\alpha$ of
$V_{\rm nr}^{(0)}$ are taken into account from Eq. (7.1), the exact value
$p_{b}=im\alpha$ is obtained, despite that at LO one has $p_{b}=im\alpha\ln a$
with $\ln a=\gamma_{E}$. One can see from Table 1 of [18] that the largest
correction to $p_{b}$ comes from the ${\cal O}(\alpha^{3})$ contribution to
$V_{\rm nr}^{(0)}$, while higher order terms in its expansion gives
contributions less than a few per mil. It is certainly worth stressing that
the suppression of these higher orders in the binding momentum is purely
numerical, and does not arise because higher powers of $\alpha$ are involved.

The ${\cal O}(\alpha^{3})$ contribution found in Eq. (8.13) to
$p_{D}^{(\frac{1}{2},1)}=im\alpha$ is connected to the relativistic correction
to $D_{{\rm LO;nr}}^{(0)}$ for $p={\cal O}(m\alpha)$ in Eq. (8.10), since the
latter is also two orders higher in the $\alpha$ counting than the LO
contribution. Our freedom to choose $C_{1}^{(\frac{1}{2},0)}$ as in Sec. 5.3,
cf. Eq. (8.7), defines our renormalization scheme when applying the once-
iterated $N/D$ method. Of course, one could choose other scheme and, as a
result, ${\cal V}_{2}^{(\frac{1}{2},0)}$ would be different in order to
guarantee that the perturbative QED amplitude is reproduced order by order. In
this respect, we can choose $C_{1}^{(\frac{1}{2},0)}$ so as to reproduce the
exact Dirac result for the binding momentum by adding to $\alpha
C_{1}^{(\frac{1}{2},0)}$, Eq. (8.7), the amount

$\displaystyle\delta C^{(\frac{1}{2},0)}$
$\displaystyle=\frac{m\alpha}{4\pi^{2}\left(1-\sqrt{1-\alpha^{2}}(1-2\ln
a)\right)}\left(\pi+(1-2\ln
a)2\left[\alpha+\sqrt{1-\alpha^{2}}\left\\{\frac{\pi}{2}-\arcsin\alpha\right\\}\right]\right)~{}.$
(8.15)

###  8.2 Bound state from the unitarized NLO QED PWA

Consider first the application of the $g$-method, such that
$V_{2}^{(\frac{1}{2},0)}(p^{2})$ according to Eq. (5.9) is given by

$\displaystyle
V_{2}^{(\frac{1}{2},0)}=T_{2}^{(\frac{1}{2},0)}-i\alpha^{2}\pi\frac{(m-(1-2\ln
a)E)^{2}}{p^{3}}~{}.$ (8.16)

The partial contribution to Eq. (8.16) from the once-iterated diagram plus the
$S$-matrix renormalization at ${\cal O}(\alpha^{2})$ is

$\displaystyle V_{\rm it}^{(\frac{1}{2},0)}(p^{2})$
$\displaystyle=\frac{\alpha^{2}\pi^{2}}{\sqrt{p^{2}}}~{}.$ (8.17)

Let us notice that this contribution has a branch-point singularity at
$p^{2}=0$ and a LHC for $p^{2}<0$. As a result, its use in Eq. (5.4) drives to
incorrect analytical properties of a PWA (in comparison with the right
analytical properties of a triangular loop [54]). This issue manifests in the
calculation of the solution to the bound-state equation,

$\displaystyle 1-i\alpha\frac{(m-E_{b}(1-2\ln
a))}{2p_{b}}-i\frac{\pi\alpha^{2}}{4}=0~{},$ (8.18)

since then $\displaystyle{p_{b}=im\alpha\ln a\left(1+\frac{\alpha^{2}(1-2\ln
a)\ln a}{4}\right)-\alpha^{3}m\pi\frac{\ln a}{4}+{\cal O}(\alpha^{5})}$ has
also a real together a positive imaginary part, which is unacceptable.

We now proceed to apply the first-iterated $N/D$ method to unitarize the NLO
QED PWA. We have to fix ${\cal V}_{2}^{(\frac{1}{2},0)}$ by employing Eq.
(5.17) which implies that

$\displaystyle{\cal V}_{2}^{(\frac{1}{2},0)}$
$\displaystyle=T_{2}^{(\frac{1}{2},0)}-i\alpha^{2}\pi\frac{(m-(1-2\ln
a)E)^{2}}{p^{3}}+\alpha^{2}\frac{(1-2\ln a)(m-(1-2\ln
a)E)}{p^{3}}\left[2p-E\ln\frac{E-p}{E+p}\right]~{}.$ (8.19)

$V_{{\rm it}}^{(\frac{1}{2},0)}(p^{2})$ is also part of ${\cal
V}_{2}^{(\frac{1}{2},0)}$ and now it contributes to
$D^{(\frac{1}{2},0)}(p^{2})$ as, cf. Eq. (5.21),

$\displaystyle-\frac{\alpha^{2}(p^{2}+p_{0}^{2})}{4}\int_{0}^{\infty}\frac{dk^{2}}{(k^{2}-p^{2})(k^{2}+p_{0}^{2})}=\frac{\alpha^{2}}{4}\ln\frac{-p^{2}}{p_{0}^{2}}~{},$
(8.20)

so that it has only RHC and then does not pose any problem on analyticity of
the $T$ matrix.

The contribution from Eq. (8.20) to the binding momentum is obtained by
multiplying it by $-p_{b}$ and evaluate the result with $p\to im\alpha\ln a$
[this is clear from the secular equation Eq. (7.6)]. It is an ${\cal
O}(\alpha^{3})$ contribution given by

$\displaystyle\delta p_{b;{\rm it}}$ $\displaystyle=-\frac{i}{2}m\alpha^{3}\ln
a\,\ln\\!\left(\\!\frac{m\alpha\ln a}{p_{0}}\\!\right)+{\cal
O}(\alpha^{5})~{}.$ (8.21)

We set $p_{0}=m$, and vary $p_{0}^{2}$ between $m^{2}/4$ and $m^{2}$, which
implies a relative uncertainty of $\frac{\alpha^{2}}{2}\ln 2$ in ${\cal
V}_{2}^{(\frac{1}{2},0)}$. For the seeking of the binding momentum this source
of uncertainty is suppressed in front of the infrared enhanced log-factor
$\frac{\alpha^{2}}{2}\ln\frac{|p_{b}|}{m}\sim\frac{\alpha^{2}}{2}\ln\alpha$
with $p_{b}\approx im\alpha$. The relative error here is around a 14%.

We next consider the contributions from the one-loop vacuum polarization
(${\cal V}_{\rm vc}^{(\frac{1}{2},0)}$), the vertex corrections plus Chung’s
diagrams (${\cal V}_{\rm vx}^{(\frac{1}{2},0)}$), and the last term to the
right in Eq. (8.19) (${\cal V}_{D}^{(\frac{1}{2},0)}$), equal to
$V_{1}^{(\frac{1}{2},0)}(\Re D_{{\rm LO}}^{(\frac{1}{2},0)}-1)$. ${\cal
V}_{\rm vc}^{(\frac{1}{2},0)}$ and ${\cal V}_{D}^{(\frac{1}{2},0)}$ admit a
low-momentum expansion in powers of $p^{2}$, and tend to a constant for $p\to
0$. For ${\cal V}_{\rm vx}^{(\frac{1}{2},0)}$ an expansion in powers of
$p^{2}$ is also possible modulo $\ln p$. One then has,

$\displaystyle{\cal V}_{\rm vc}^{(\frac{1}{2},0)}(p^{2})$
$\displaystyle=\frac{8\alpha^{2}}{15m}+{\cal O}(p^{2})~{},$ (8.22)
$\displaystyle{\cal V}_{\rm vx}^{(\frac{1}{2},0)}(p^{2})$
$\displaystyle=-\frac{\alpha^{2}}{24m}\left(28-3\pi+16\ln\frac{m^{4}a^{4}}{256p^{4}}\right)+{\cal
O}(p^{2})~{},$ $\displaystyle{\cal V}_{D}^{(\frac{1}{2},0)}$
$\displaystyle=-\frac{4\alpha^{2}}{3m}(1-2\ln a)\ln a+{\cal O}(p^{2})~{}.$

When inserted in the calculation of the function $D^{(\frac{1}{2},0)}(p^{2})$
a subtraction is needed because they vanish as $1/p$(modulo powers of $\ln p$)
for $p\to\infty$. Within the iterated $N/D$ method we can proceed as we did
concerning $D_{LO}^{(\frac{1}{2},0)}(p^{2})$ and Eq. (8.20) from ${\cal
V}_{{\rm it}}^{(\frac{1}{2},0)}$, so that the constant term that adds to the 1
in $D^{(J,\ell)}(p^{2})$, cf. Eq. (8.7), can be removed. This is then
reshuffled as a contribution to the ${\cal V}^{(\frac{1}{2},0)}(p^{2})$
function by the matching procedure explained in Sec. 5.3, such that the
perturbative QED PWA can be reproduced at any order in the expansion in powers
of $\alpha$. The removal of this constant term is equivalent to take the
subtraction at $p_{0}=0$ with zero subtraction constant when calculating their
contribution to $D^{(\frac{1}{2},0)}(p^{2})$, denoted by $\delta
D^{(\frac{1}{2},0)}_{{\rm vc;vx;}D}$. Thus, we have

$\displaystyle\delta D^{(\frac{1}{2},0)}_{{\rm vc;vx;}D}(p^{2})$
$\displaystyle=-\frac{p^{2}}{4\pi^{2}}\int_{0}^{\infty}dk^{2}\frac{k\\{{\cal
V}_{{\rm vc}}^{(\frac{1}{2},0)}(k^{2})+{\cal V}_{{\rm
vx}}^{(\frac{1}{2},0)}(k^{2})+{\cal
V}_{D}^{(\frac{1}{2},0)}(k^{2})\\}}{k^{2}(k^{2}-p^{2})}\,,$ (8.23)

and it contributes to the binding momentum by an amount

$\displaystyle\delta p_{b}$ $\displaystyle=-p_{b;1}\delta
D^{(\frac{1}{2},0)}_{{\rm vc;vx;}D}(p_{b;1}^{2})+{\cal O}(\alpha^{6})~{},$
(8.24)

and it is ${\cal O}(\alpha^{4})$. This can be manifestly shown by introducing
the dimensionless variable $y$ such that $k=m\alpha\ln a\,y$, together with
the dimensionless function $W^{(\frac{1}{2},0)}=m{\cal
V}^{(\frac{1}{2},0)}/\alpha^{2}$. Then, we can rewrite Eq. (8.23) for
$p=im\alpha\ln a$ as

$\displaystyle\delta p_{b}$
$\displaystyle=-i\frac{m\alpha^{4}\ln^{2}a}{2\pi^{2}}\int_{0}^{\infty}dy\frac{\left(W_{{\rm
vc}}^{(\frac{1}{2},0)}+W_{{\rm
vx}}^{(\frac{1}{2},0)}+W_{D}^{(\frac{1}{2},0)}\right)}{y^{2}+1}~{}.$ (8.25)

These integrals are calculated numerically, and we have checked that they
perfectly agree with the algebraic results that one obtains by applying the
$g$-method

$\displaystyle\delta p_{b;{\rm vc}}$
$\displaystyle=-im\alpha^{4}\gamma_{E}^{2}\frac{2}{15\pi}=-im\alpha^{4}\gamma_{E}^{2}\,4.2\,10^{-2}~{},$
(8.26) $\displaystyle\delta p_{b;{\rm vx}}$
$\displaystyle=im\alpha^{4}\gamma_{E}^{2}\frac{1}{96\pi}\left(28-3\pi+64\ln\frac{a}{4\alpha\ln
a}\right)=im\alpha^{4}\gamma_{E}^{2}\,1.0~{},$ $\displaystyle\delta p_{b;D}$
$\displaystyle=-im\alpha^{4}\gamma_{E}^{3}\frac{1-2\gamma_{E}}{3\pi}=-im\alpha^{4}\gamma_{E}^{2}\,9.4\,10^{-2}~{},$

plus ${\cal O}(\alpha^{6})$ for all cases.

##  9 Conclusions and outlook

In this study we have unitarized in partial wave amplitudes the one-loop
scattering amplitude in QED of a relativistic electron by an external Coulomb
potential. This work is an extension of Refs. [16, 18] that unitarized
nonrelativistic Coulomb scattering. The approach of Refs. [16, 17] has been
extended to deal with infrared divergences that go beyond the generalized
divergent Coulomb phase already accounted for by the method of these
references. We have then applied the redefinition of the electron asymptotic
states by accounting for a cloud of real soft photons, first introduced by
Chung in Ref. [23], and afterwards generalized by Faddeev-Kulish [22]. The
resulting QED scattering amplitudes that we have calculated up to NLO are then
infrared finite.

Based on them, we have calculated infrared finite partial-wave amplitudes,
first calculated perturbatively and then unitarized by applying three
different unitarization methods: the algebraic $N/D$ method, renamed
$g$-method for brevity, the Inverse Amplitude Method (IAM), and the first-
iterated $N/D$ method. In this way, we have presented explicit results for the
$S$\- and $P$-wave scattering. The results from the different unitarization
methods typically agree, though for $|p|\gtrsim m/\alpha$ the IAM fails in
those partial-wave amplitudes for which the real part of the $T$ matrix
crosses zeroes. This deficiency at high momentum is similar to the well-known
issue for the straightforward application of this method around an Adler zero
in meson-meson scattering [31, 32, 33].

We have also discussed the secular equation for the calculation of the binding
momenta by applying the different unitarization methods. However, it is only
the first-iterated $N/D$ method that is suitable for looking for bound states
which position lie on top of the left-hand cut, as it is the case here with
the perturbative QED scattering amplitude calculated to the one-loop level. We
have then evaluated the binding momentum for the fundamental state from the
unitarized perturbative QED partial-wave amplitudes calculated at leading and
up to next-to-leading orders.

As an outlook, we propose that this method could be applied to study
properties of the hydrogen bound states by unitarizing higher-order
perturbative QED scattering amplitudes. This is a novel method that does not
require to use quantum mechanical wave functions for the electron, and it is
based entirely on the knowledge of the $S$ matrix. It can be applied for bound
states with $\ell\leq 1$ and principal quantum number up to $n=5$ in $S$ wave
and $n=2$ for $P$ waves. Similarly, this unitarization method could be applied
to study the interplay between strong and electromagnetic interactions in
processes of interest, like hadronic or nuclear low-energy scattering of
charged particles, electromagnetic corrections to hadronic and nuclear bound
states and resonances, etc. For the near future, we are interested in
extending the results of Refs. [16, 17] by unitarizing graviton-graviton
scattering from a one-loop calculation in the low-energy effective-field
theory of gravity [59]. Of special interest will be to look for the graviball
resonance pole found in these references from the unitarized leading-order
graviton-graviton scattering amplitudes.

## Acknowledgements

JAO would like to acknowledge partial financial support to the Grant
PID2022-136510NB-C32 funded by
MCIN/AEI/10.13039/501100011033/ and FEDER, UE, and to the EU Horizon 2020
research and innovation program, STRONG-2020 project, under grant agreement
no. 824093. MP acknowledges financial support from the ANII-FCE-175902
project.

##  Appendix A QED one-loop results

In this appendix we provide a concise summary of QED results up to one loop.
The results are first expressed as Feynman integrals, then simplified to
Passarino-Veltman integrals, and finally presented analytically. All the
integrals are evaluated in dimensional regularization. A photon mass, denoted
by $\lambda$, appears explicitly only in the terms that give rise to infrared
divergences, otherwise is set to zero. Some definitions are repeated here
again for consistency in the appendix. In the following $q=p^{\prime}-p$,
$t=q^{2}$ and $\sigma=\sigma(t)=\sqrt{1-4m^{2}/t}$.

###  A.1 Electron self energy

Figure 7:  One-loop diagram to the electron self energy.

The one-loop diagram for the electron self energy is depicted in Fig. 7, which
expression is denoted as $-i\Sigma(\not{p})$. It takes the form,

$\displaystyle-i\Sigma(\not{p})$
$\displaystyle=-e^{2}\mu^{-2\epsilon}\int\frac{d^{d}k}{(2\pi)^{d}}\frac{(2-d)(\not{p}-\not{k})+dm}{\left(k^{2}-\lambda^{2}+i\varepsilon\right)\left((p-k)^{2}-m^{2}+i\varepsilon\right)}~{},$
(A.1)

where $d=4+2\epsilon$, $\mu$ represents the dimension of the coupling
constant, $p$ is the electron incoming momentum and $k$ is the photon one.
Thus,

$\displaystyle-i\Sigma(\not{p})=-e^{2}\left(\frac{(d-2)\not{p}}{2p^{2}}\left[A(m^{2})-A(\lambda^{2})+\left(\lambda^{2}-p^{2}-m^{2}\right)B(\lambda^{2},m^{2},p^{2})\right]+d\hskip
2.84526ptmB(\lambda^{2},m^{2},p^{2})\right)~{},$ (A.2)

where the integrals $A$ and $B$ are defined in Sec.A.5, together with their
divergent ($A_{d}$ and $B_{d}$) and finite ($A_{f}$ and $B_{f}$) parts.

In order to have the pole in the electron propagator at the physical mass with
residue equal to one, the renormalization proceeds by requiring that

$\displaystyle\delta Z_{m}=-\frac{\Sigma(\not{p}=m)}{m}\quad\&\quad\delta
Z_{2}=\frac{d\Sigma(\not{p})}{d\not{p}}\Big{|}_{\not{p}=m}~{}.$ (A.3)

Therefore,

$\displaystyle\delta Z_{m}$
$\displaystyle=\frac{e^{2}}{(4\pi)^{2}}\left(\frac{3}{\hat{\epsilon}}+3\ln\frac{m^{2}}{\mu^{2}}-4\right)~{},$
(A.4)

where here $\lambda^{2}=0$ has been set, and

$\displaystyle\frac{1}{\hat{\epsilon}}$
$\displaystyle=\frac{1}{\epsilon}+\gamma_{E}-\ln 4\pi~{},$ (A.5)

with $\gamma_{E}$ the Euler constant. In turn,

$\displaystyle\delta Z_{2}$
$\displaystyle=\frac{e^{2}}{(4\pi)^{2}}\left(\frac{1}{\hat{\epsilon}}+\ln\frac{m^{2}}{\mu^{2}}+2\ln\frac{m^{2}}{\lambda^{2}}-4\right)~{},$
(A.6)

where the infrared divergence as well as the ultraviolet one appear in this
renormalization factor.

###  A.2 Vacuum polarization

We compute here the correction coming from diagram (d) of Fig. 1, properly
renormalized with diagram (e). The expression for the one-loop amplitude given
by the diagram (d) reads

$\bar{u}_{\sigma^{\prime}}(\boldsymbol{p^{\prime}})(-ie\gamma^{\mu})u_{\sigma}(\boldsymbol{p})\left(\frac{-ig_{\mu\nu}}{q^{2}-\lambda^{2}}\right)i\Pi^{\nu\rho}(q)\mathcal{A}_{\rho}(q)~{},$
(A.7)

where

$\displaystyle i\Pi^{\nu\rho}(q)$
$\displaystyle=-de^{2}\mu^{-2\epsilon}\int\frac{d^{d}k}{(2\pi)^{d}}\frac{(-k^{2}+k\cdot
q+m^{2})g^{\nu\rho}+2k^{\nu}k^{\rho}-k^{\nu}q^{\rho}-q^{\nu}k^{\rho}}{(k^{2}-m^{2}+i\varepsilon)((k-q)^{2}-m^{2}+i\varepsilon)}~{},$
(A.8)

which, after the reduction to master integrals, has the form
$\Pi^{\nu\rho}(q)=\Pi(q^{2})\left(q^{2}g^{\nu\rho}-q^{\nu}q^{\rho}\right)$,
with

$\displaystyle\Pi(q^{2})=ie^{2}\frac{d}{(d-1)q^{2}}\left((2-d)A(m^{2})+\left(2m^{2}+(d-2)\frac{q^{2}}{2}\right)B(m^{2},m^{2},q^{2})\right).$
(A.9)

The effect of taking into account the counterterm reduces to subtract the
expression evaluated at vanishing momentum, so the counterterm is

$\delta
Z_{3}=\Pi(0)=\frac{e^{2}}{24\pi^{2}}\left(\frac{2}{\hat{\epsilon}}+2\ln\frac{m^{2}}{\mu^{2}}+1\right)~{}.$
(A.10)

It leads us to the amplitude

$\displaystyle\bar{u}_{\sigma^{\prime}}(\boldsymbol{p^{\prime}})(-ie\gamma^{\mu}\mathcal{A}_{\mu}(q))u_{\sigma}(\boldsymbol{p})\frac{q^{2}}{q^{2}-\lambda^{2}}\Pi_{R}(q^{2}),$
(A.11)

with $\Pi_{R}(q^{2})=\Pi(q^{2})-\Pi(0)$, corresponding to Eq.(2.5) .

###  A.3 Vertex modification of the electron-photon coupling

The expression for the sum of the diagrams (b) and (c) of Fig. 1 can be
explicitly written as
$\displaystyle{\bar{u}_{\sigma^{\prime}}\left(\boldsymbol{p}^{\prime}\right)(-ie\Gamma^{\mu})u_{\sigma}(\boldsymbol{p})A_{\mu}(q)}$,
with

$\displaystyle\Gamma^{\mu}$
$\displaystyle=-ie^{2}\int\frac{d^{d}k}{(2\pi)^{d}}\frac{\gamma^{\alpha}(\not{p^{\prime}}-\not{k}+m)\gamma^{\mu}(\not{p}-\not{k}+m)\gamma_{\alpha}}{((p^{\prime}-k)^{2}-m^{2}+i\varepsilon)((p-k)^{2}-m^{2}+i\varepsilon)(k^{2}-\lambda^{2}+i\varepsilon)}-ie^{2}\delta
Z_{1}~{}.$ (A.12)

The former diagrams are often written in terms of two form factors:

$\displaystyle\Gamma^{\mu}(p^{\prime},p)=F_{1}(q)\gamma^{\mu}+\frac{i\sigma^{\mu\nu}}{2m}q_{\nu}F_{2}(q)~{},$
(A.13)

with $\sigma^{\mu\nu}=\frac{i}{2}[\gamma^{\mu},\gamma^{\nu}]$.

The analysis is performed using a general $\mu$-index for the potential that
we should consider equal to $0$ for the Coulomb potential. The final result of
the renormalized vertex is:

$\displaystyle F_{1}(q)$
$\displaystyle=ie^{2}\left(2(t-2m^{2})\mathcal{C}(\lambda^{2},m^{2},t)+\frac{4A_{f}(m^{2})}{t-4m^{2}}+\frac{B_{f}(m^{2},m^{2},t)(3t-8m^{2})}{t-4m^{2}}+\frac{4B_{f}(0,m^{2},m^{2})(2m^{2}-t)}{t-4m^{2}}+\frac{2m^{2}-t}{8\pi^{2}(t-4m^{2})}\right)$
$\displaystyle+\delta Z_{2}^{f}~{},$ (A.14)

with $\delta Z_{2}^{f}=\delta Z_{2}-e^{2}/(16\pi^{2}\hat{\epsilon})$, and

$\displaystyle F_{2}(q)$
$\displaystyle=-\frac{4ie^{2}}{t-4m^{2}}\left(A_{f}(m^{2})+m^{2}B_{f}(m^{2},m^{2},t)-2m^{2}B_{f}(0,m^{2},m^{2})+\frac{im^{2}}{16\pi^{2}}\right)~{}.$
(A.15)

Let us notice that the infrared divergences in Eq. (A.3) stem only from the
triangle integral $\mathcal{C}(\lambda^{2},m^{2},t)$ and $\delta Z_{2}^{f}$,
while $F_{2}(q)$ has no infrared divergences. We also note that $F_{1}(q)$
vanishes linearly in $q^{2}$ for $q^{2}\to 0$. Renormalization takes place by
imposing the following condition on the total vertex function,
$\Lambda^{\mu}(p^{\prime},p)=\gamma^{\mu}+\Gamma^{\mu}(p^{\prime},p)$,

$\displaystyle\bar{u}(\boldsymbol{p}^{\prime})\Lambda^{\mu}(p,p^{\prime})u(\boldsymbol{p})\Big{|}_{\not{p}_{1}=\not{p}_{2}=m}=\bar{u}(\boldsymbol{p}^{\prime})\gamma^{\mu}u(\boldsymbol{p})\Big{|}_{\not{p}_{1}=\not{p}_{2}=m}$
(A.16)

The latter implies that

$\displaystyle\bar{u}(\boldsymbol{p})\Gamma^{\mu}(p,p)u(\boldsymbol{p})\Big{|}_{\not{p}=m}=(\delta
Z_{e}+\frac{1}{2}\delta
Z_{3})\bar{u}(\boldsymbol{p})\gamma^{\mu}u(\boldsymbol{p})=0~{},$ (A.17)

in agreement with the Ward identity. Equation (2.3) results by substituting in
the expressions given for $F_{1}(q)$ and $F_{2}(q)$ the master integrals in
the appendix A.5 and $\delta Z_{2}^{f}$.

###  A.4 Iterated diagram

The amplitude for the iteration of the Coulomb potential corresponding to
diagram (f) in Fig.1 is

$\displaystyle-i(ee^{\prime})^{2}(2\pi)\delta(E^{\prime}-E)\bar{u}_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\int\frac{d^{3}q}{(2\pi)^{3}}\frac{E\gamma^{0}+\vec{q}.\vec{\gamma}+m}{(({\boldsymbol{p}}^{\prime}-\boldsymbol{q})^{2}+\lambda^{2}+i\varepsilon)(({\boldsymbol{p}}-\boldsymbol{q})^{2}+\lambda^{2}+i\varepsilon)({\boldsymbol{p}}^{2}-\boldsymbol{q}^{2}+i\varepsilon)}u_{\sigma}({\boldsymbol{p}})$
(A.18)
$\displaystyle\equiv-i(ee^{\prime})^{2}(2\pi)\delta(E^{\prime}-E)\bar{u}_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\left((E\gamma^{0}+m)I+\frac{\vec{\gamma}.({\boldsymbol{p}}+{\boldsymbol{p}}^{\prime})}{2}J\right)u_{\sigma}({\boldsymbol{p}})~{},$

where

$I=\int\frac{d^{3}q}{(2\pi)^{3}}\frac{1}{(({\boldsymbol{p}}^{\prime}-\boldsymbol{q})^{2}+\lambda^{2}+i\varepsilon)(({\boldsymbol{p}}-\boldsymbol{q})^{2}+\lambda^{2}+i\varepsilon)({\boldsymbol{p}}^{2}-\boldsymbol{q}^{2}+i\varepsilon)}~{},$
(A.19)

and

$\frac{{\boldsymbol{p}}^{\prime}+{\boldsymbol{p}}}{2}J=\int\frac{d^{3}q}{(2\pi)^{3}}\frac{{\boldsymbol{q}}}{(({\boldsymbol{p}}^{\prime}-\boldsymbol{q})^{2}+\lambda^{2}+i\varepsilon)(({\boldsymbol{p}}-\boldsymbol{q})^{2}+\lambda^{2}+i\varepsilon)({\boldsymbol{p}}^{2}-\boldsymbol{q}^{2}+i\varepsilon)}~{},$
(A.20)

which are the same loop functions as introduced by Dalitz in Ref. [3] except
for a global factor $(2\pi)^{-3}$. However, we have now evaluated them using
dimensional regularization. After some manipulations we obtain:

$\displaystyle I$ $\displaystyle=\frac{-1}{8\pi
p\sin\frac{\theta}{2}\sqrt{4p^{4}\sin^{2}\frac{\theta}{2}+4p^{2}\lambda^{2}+\lambda^{4}}}\Big{(}\arctan\frac{\lambda
p\sin\frac{\theta}{2}}{\sqrt{4p^{4}\sin^{2}\frac{\theta}{2}+4p^{2}\lambda^{2}+\lambda^{4}}}+\frac{i}{2}\ln\frac{\sqrt{4p^{4}\sin^{2}\frac{\theta}{2}+4p^{2}\lambda^{2}+\lambda^{4}}+2p^{2}\sin\frac{\theta}{2}}{\sqrt{4p^{4}\sin^{2}\frac{\theta}{2}+4p^{2}\lambda^{2}+\lambda^{4}}-2p^{2}\sin\frac{\theta}{2}}\Big{)}~{},$
(A.21) $\displaystyle J$
$\displaystyle=-\frac{1}{2{\boldsymbol{p}}^{2}\cos^{2}\frac{\theta}{2}}\Big{(}{\cal
A}_{3}(\lambda^{2},{\boldsymbol{q}}^{2})+{\cal
B}_{3}(\lambda^{2},{\boldsymbol{p}}^{2})-(2{\boldsymbol{p}}^{2}+\lambda^{2})I\Big{)}$
$\displaystyle=\frac{(2{\boldsymbol{p}}^{2}+\lambda^{2})I}{2{\boldsymbol{p}}^{2}\cos^{2}\frac{\theta}{2}}-\frac{\arctan\frac{p\sin\theta/2}{\lambda}}{16\pi|{\boldsymbol{p}}|^{3}\cos^{2}\frac{\theta}{2}\sin\frac{\theta}{2}}+\frac{i}{16\pi|{\boldsymbol{p}}|^{3}\cos^{2}\frac{\theta}{2}}\ln\left(1-\frac{2i|{\boldsymbol{p}}|}{\lambda}\right)~{}.$
(A.22)

The master integrals in three dimensions ${\cal A}_{3}$ and ${\cal B}_{3}$ are
evaluated in Sec. A.5.1. Both functions $I$ and $J$ calculated here are in
agreement with their computation in Ref. [3]. By using the Dirac equation in
Eq.(A.18) we can rewrite it as

$\displaystyle-
ie^{4}(2\pi)\delta(E^{\prime}-E)\bar{u}_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\left((E\gamma^{0}+m)I+(E\gamma^{0}-m))J\right)u_{\sigma}({\boldsymbol{p}})$
$\displaystyle=-ie^{4}(2\pi)\delta(E^{\prime}-E)\bar{u}_{\sigma^{\prime}}({\boldsymbol{p}}^{\prime})\left(E\gamma^{0}(I+J)+m(I-J)\right)u_{\sigma}({\boldsymbol{p}})~{}.$
(A.23)

As discussed in Secs. 2.4 and 4.2, one needs to carefully take the limit
$\lambda\to 0$ for the iterated diagram $A_{3}$ in order to properly handle
the infrared divergences and its cancellation within partial-wave amplitudes.
This is a point obviated in Ref. [3] since it was interested in the Born
series associated to the Coulomb potential. The limit $\lambda\to 0$ can be
taken straightforwardly for the combination $I-J$ because of the factor
$\displaystyle{1-\frac{1}{\cos^{2}\theta/2}}=-\tan^{2}\frac{\theta}{2}$ that
vanishes for $\theta\to 0$. The resulting expression is the term between
square brackets multiplied by $m$ in Eq. (2.4).

The infrared divergent part for the sum $I+J$ is more delicate. In the term
$\sqrt{4p^{4}\sin^{2}\frac{\theta}{2}+4p^{2}\lambda^{2}+\lambda^{4}}$, which
repeatedly appears in the expression for $I$ in Eq. (A.21), we can neglect
$\lambda^{4}$ in front of $p^{2}\lambda^{2}$, but the latter cannot always be
neglected compared with $4p^{4}\sin^{2}\theta/2$, as the latter vanishes for
$\theta\to 0$. However, $\displaystyle{\lambda
p\sin\frac{\theta}{2}/\sqrt{4p^{4}\sin^{2}\frac{\theta}{2}+4p^{2}\lambda^{2}+\lambda^{4}}\to
0}$ for $\lambda\to 0$ because this ratio vanishes in any case, either when
$\lambda^{2}\ll p^{2}\theta^{2}$ or $p^{2}\theta^{2}\ll\lambda^{2}$. For the
last two terms in the expression for $J$ in Eq. (A.4) its asymptotic behavior
for $\lambda\to 0$ is straightforward, and it reads

$\displaystyle-\frac{\arctan\frac{p\sin\theta/2}{\lambda}}{16\pi|{\boldsymbol{p}}|^{3}\cos^{2}\frac{\theta}{2}\sin\frac{\theta}{2}}+\frac{i}{16\pi|{\boldsymbol{p}}|^{3}\cos^{2}\frac{\theta}{2}}\ln\left(1-\frac{2i|{\boldsymbol{p}}|}{\lambda}\right)\xrightarrow[\lambda\to
0]{}-\frac{1}{32|{\boldsymbol{p}}|^{3}\cos^{2}\frac{\theta}{2}\sin\frac{\theta}{2}}+\frac{i}{16\pi|{\boldsymbol{p}}|^{3}\cos^{2}\frac{\theta}{2}}\left(\ln\frac{2|{\boldsymbol{p}}|}{\lambda}-i\frac{\pi}{2}\right)$
(A.24)

Then, the term between square brackets multiplied by $E\gamma^{0}$ in Eq.
(2.4) results.

We have detected some typos in Ref. [3], which we list here with reference to
the numbered equations in [3]: $m$ should be changed to $-m$ in the numerator
of the second line of the equation before (2.3) and in the first line of this
equation. In addition, a global minus sign is missing in the expression for
$I$ in Eq. (2.6) in the limit $\lambda\to 0$. We stress that the final
expression for $M$ in the second line of Eq. (2.3) is correct.

###  A.5 Master Integrals

In this part of the Appendix we summarize the definitions and analytical
results for master integrals in $d=4+2\epsilon$ dimensions. To begin with, the
tadpole integral that is represented with $A(m^{2})$ has the expression:

$\displaystyle
A(m^{2})=\mu^{-2\epsilon}\int\frac{d^{d}q}{(2\pi)^{d/2}}\frac{1}{q^{2}-m^{2}+i\varepsilon}=-\frac{im^{2}}{(4\pi)^{2}}\left(\frac{1}{\hat{\epsilon}}-1+\ln\left(\frac{m^{2}}{\mu^{2}}\right)\right)~{},$
(A.25)

from where we can separate the divergent and finite parts of it in order to
define
$\displaystyle{A(m^{2})=\frac{A_{d}(m^{2})}{\hat{\epsilon}}+A_{f}(m^{2})}$:

$\displaystyle A_{d}(m^{2})=-\frac{im^{2}}{(4\pi)^{2}}\hskip
28.45274pt\&\hskip
28.45274ptA_{f}(m^{2})=-\frac{im^{2}}{(4\pi)^{2}}\left(\ln\left(\frac{m^{2}}{\mu^{2}}\right)-1\right)~{}.$
(A.26)

A similar procedure is used to describe the one loop sunset diagram,

$\displaystyle
B(m_{1}^{2},m_{2}^{2},p^{2})=\mu^{-2\epsilon}\int\frac{d^{d}q}{(2\pi)^{d/2}}\frac{1}{q^{2}-m_{1}^{2}+i\varepsilon}\frac{1}{(q-p)^{2}-m_{2}^{2}+i\varepsilon}=\frac{B_{d}}{\hat{\epsilon}}+B_{f}(m_{1}^{2},m_{2}^{2},p^{2})~{},$
(A.27)

where $B_{f}$ is finite. It then results that

$\displaystyle B_{d}$
$\displaystyle=-\frac{i}{(4\pi)^{2}}=\frac{A_{d}(m^{2})}{m^{2}}~{},$ (A.28)
$\displaystyle B_{f}(\lambda^{2},m^{2},p^{2})$
$\displaystyle=\frac{i}{(4\pi)^{2}}\left(2-\ln\frac{\lambda
m}{\mu^{2}}+\frac{m^{2}-\lambda^{2}}{2p^{2}}\ln\frac{\lambda^{2}}{m^{2}}-\frac{h}{2p^{2}}\left[\ln\frac{h+p^{2}+m^{2}-\lambda^{2}}{h-p^{2}-m^{2}+\lambda^{2}}+\ln\frac{h+p^{2}-m^{2}+\lambda^{2}}{h-p^{2}+m^{2}-\lambda^{2}}\right]\right)~{},$
$\displaystyle h$
$\displaystyle=\sqrt{(p^{2})^{2}+\lambda^{4}+m^{4}-2m^{2}p^{2}-2m^{2}\lambda^{2}-2p^{2}\lambda^{2}}~{},$

with $h$ the square root of the Källen function. For the case of equal masses,
$B_{f}$ simplifies as

$\displaystyle
B_{f}(m^{2},m^{2},p^{2})=\frac{i}{16\pi^{2}}\left(2-\ln\left(\frac{m^{2}}{\mu^{2}}\right)+\sigma(p^{2})\ln\left(\frac{\sigma(p^{2})-1}{\sigma(p^{2})+1}\right)\right)~{}.$
(A.29)

Other particular results derived from Eq. (A.28) are

$\displaystyle B_{f}(0,m^{2},m^{2})$
$\displaystyle=\frac{i}{(4\pi)^{2}}\left(2-\ln\frac{m^{2}}{\mu^{2}}\right)~{},$