Hadronic bound states in SU(2) from Dyson–Schwinger equations

By using the Dyson–Schwinger/Bethe–Salpeter formalism in Euclidean spacetime, we calculate the ground state spectrum of J≤1\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$J\le 1$$\end{document} hadrons in an SU(2) gauge theory with two fundamental fermions. We show that the rainbow-ladder truncation, commonly employed in QCD studies, is unsuitable for a description of an SU(2) theory. This we remedy by truncating at the level of the quark–gluon vertex Dyson–Schwinger equation in a diagrammatic expansion. Results obtained within this novel approach show good agreement with lattice studies. These findings emphasize the need to use techniques more sophisticated than rainbow-ladder when investigating generic strongly interacting gauge theories.


Introduction
Quantum chromodynamics (QCD) is a strongly interacting gauge theory whose study has proven to be one of the most formidable challenges of modern theoretical physics. While the high-energy regime of QCD is by now relatively well explored in terms of perturbation theory, the arguably more interesting (and intrinsically non-perturbative) phenomena such as dynamical chiral symmetry breaking and confinement are yet to be fully understood.
One of the strategies which might lead to our better understanding of QCD is to investigate theories which are QCDlike, but have certain properties that make them technically less challenging than QCD itself. A prime example is provided by studies of SU(2) gauge theories with an even number of fermion flavors. Lattice simulations of these theories at non-zero chemical potential do not suffer from the sign problem, and such models thus provide ideal conditions to study the phase diagram of strongly interacting matter [1][2][3][4][5][6][7][8][9][10][11]. a e-mail: milan.vujinovic@uni-graz.at b e-mail: richard.williams@theo.physik.uni-giessen. de Here we wish to concentrate on the situation with two fundamentally charged Dirac fermions [1,8,10,11]. Such a theory may also be interesting in the context of a unified description of cold asymmetric Dark Matter (DM) and dynamical electroweak (EW) symmetry breaking [12][13][14], wherein the ground state hadronic spectrum at T = 0, μ = 0 is of great importance. It is exactly this hadronic spectrum that will be the central focus of our study.
However, as we will show in this paper, the RL truncation performs unsatisfactorily when adapted to an SU(2) theory with two fundamental flavors, even though the theory is expected to have QCD-like dynamics. We discuss possible reasons for this in more detail in Sect. 2. Here we only comment that we strongly believe that (most) of the inadequacy of RL method comes from its weak connection to the underlying gauge sector. Remedying this requires the use of beyond rainbow-ladder (BRL) techniques, with our preference toward those based on the diagrammatic expansion of quark-gluon vertex DSE [36][37][38][39][40][41][42][43][44][45][46]. While there are other BRL methods available [47][48][49][50][51][52][53], we choose the diagrammatic approach as it makes it easier to study the influence of the gauge sector on hadronic observables. Our aim in this paper is thus not only to provide a continuum calculation complementary to the lattice investigations of [12,13], but also to explicitly demonstrate the importance of using BRL methods when studying generic strongly interacting theories.
This manuscript is organized as follows. In Sect. 2 we discuss the DSEs relevant for our calculation, and also describe in some detail the approximations and model inputs we employ. In Sect. 3 we describe the extrapolation procedures used to obtain hadron masses, and provide estimates for errors coming from extrapolation. The results are discussed and compared to relevant lattice data in Sect. 4. We conclude in Sect. 5.

Framework
In a theory with two colors, both mesons and baryons (diquarks) can be described in terms of a two-body Bethe-Salpeter equation. For the meson where k stands for d 4 k/(2π) 4 and M ( p, P) is the meson amplitude with appropriate J PC quantum numbers, relative momentum p and total momentum P, and the meson wavefunction is The quark propagators are S(k ± ), at momenta k + = k + η P and k − = k − (1 − η)P, with k the loop momentum and η ∈ [0, 1] the momentum partition factor. In a covariant study, the results are independent of η: for concreteness, we work with η = 1/2. The final ingredient in Eq. (1) is the quark-antiquark 4-point interaction kernel K ( p, k, P). A diagrammatic representation of Eq. (1) for mesons is given in Fig. 1.
In order to solve the BSE, one clearly needs as input the quark propagator S( p). This Green function is decomposed as with Z f ( p 2 ) the quark wavefunction and M( p 2 ) the dynamical quark mass. The tree-level form is given by where Z m is the quark mass renormalization constant. The quark propagator satisfies its own DSE, see Fig. 2, and is given by Here, ν ( p, k) and D μν (k) are the full quark-gluon vertex and gluon propagator, respectively. Renormalization constants of the quark field and quark-gluon vertex are Z 2 and Z 1 f . They are related through a Slavnov-Taylor identity which takes a simple form when employing a miniMOM  The truncated two-body kernel in rainbow-ladder approximation scheme [54] in Landau gauge, Z 1 f = Z 2 / Z 3 with Z 3 the renormalization of the ghost propagator. The 4-point interaction kernel K ( p, k, P) of Eq. (1) is connected to the self-energy part ( p) of quark propagator DSE through the axial-vector Ward-Takahashi identity (axWTI) This identity encodes the chiral properties of the theory, and severely constrains the form of the BSE interaction kernel once an approximation for the quark DSE has been chosen. A direct connection is provided through the action of 'cutting' internal quark lines [36,37].

Rainbow-ladder
The 'rainbow' part of RL truncation refers to the replacement of the full quark-gluon vertex in Eq. (3) by i.e. its tree-level form augmented by a model dressing function, λ(k 2 ), that is, a function of the gluon momentum alone. The corresponding axWTI-preserving approximation for BSE kernel is that of one gluon exchange (the 'ladder'), which we show diagrammatically in Fig. 3.
In the RL approach, the model dressing function λ(k 2 ) of Eq. (5) is often combined with the dressing of the gluon propagator D μν (k 2 ) into a single model function, constructed to reproduce correctly some hadronic observables, usually m π and f π . While this method has shown considerable success in QCD phenomenology (see e.g. [55,56] for some of the limitations of the model), in an SU(2) theory the approach seems rather unsuitable, especially in the 1 ++ channel: see Table 2 for details.
There are two primary reasons why RL will not perform satisfactorily in a generic strongly interacting theory. One reason is with regards to its very limited interaction structure (γ ν × γ μ ) which offers no variation in interaction strength across different meson channels. The second is that the connection to the underlying gauge dynamics is typically lost in the construction of an effective quark-gluon interaction; this prevents adequate rescaling of parameters such as g 2 N c that cannot be translated into a re-parameterization of an effective model.

Beyond rainbow-ladder
A BRL approach which is well suited for studying the influence of underlying Yang-Mills sector on the hadronic observables is based on the quark-gluon vertex [41,44,46,[57][58][59][60][61][62]. Here, we focus on the truncated form of the DSE [46] shown in Fig. 4. Within this approximation, only the so-called non-Abelian (NA) diagram is kept in the quark-gluon vertex self-energy. The truncated kernel, consistent with constraints from chiral symmetry, is shown in Fig. 5.
So that the Bethe-Salpeter equation can be tackled, the evaluation of a fully self-consistent quark-gluon vertex is not performed. That is, the full calculated vertex (denoted by a red filled circle in Fig. 4) is not back-coupled into the non-Abelian diagram. Instead, the internal vertices (orange squares in Fig. 4) are modeled by Eq. (5) with λ(k 2 ) constructed such that it strongly resembles the tree-level projection of the full quark-gluon vertex at each iteration step; essentially, it depends upon a function (M 0 ) that encodes the interaction strength in terms of the dynamically generated quark mass. We used the parametrization Eq. (21) of Ref. [46], with modifications that account for the change N c = 2 and the rescaling of the gauge coupling g 2 (g 2 N c is left invariant). For (M 0 ) we use the functional form given in Eq. (22) of [46] with parameters a 2.44, b 1.79, c −0.20, d 0.30. The procedure described therein is used for solving the resultant coupled DSE system of a quark propagator and quark-gluon vertex.
We emphasize here that this model is, in a sense, highly constrained. Namely, once the input for the ghost and gluon propagators (which we will discuss shortly) and the trunca-  tion of the quark-gluon vertex DSE have been chosen, all other parts of the calculation are fixed. The BSE kernel follows from the axWTI, and the model dressing λ(k 2 ) of Eq. (5) follows from the tree-level projection of the full quark-gluon vertex. We will re-iterate this point in Sect. 4.1, when we provide an estimation of the model dependence.
The final ingredient which we need to specify in our calculation is the gluon propagator D μν (k). We work in Landau gauge, where this Green function takes the form with T μν (k) = δ μν − k μ k ν /k 2 the transverse projector with respect to momentum k. The gluon dressing function which we use is plotted in Fig. 6. The details of this function and its parametrization can be found in [63]. We point out that the gluon which we employ corresponds to a quenched DSE calculation. Ignoring the back-reaction of quarks onto the Yang-Mills sector is usually considered a good approximation for theories with QCD-like dynamics, as the corresponding effect on 'observables' like the chiral condensate, f π and others is quite small [64]. However, the quenched approximation should be reconsidered in theories which have (nearly) conformal, or 'walking' dynamics. Walking dynamics arises naturally in models with a relatively large number of light fundamentally charged fermions [65][66][67][68][69][70][71][72][73], or fermions belonging to higher-dimensional representations of the gauge group [74][75][76][77][78][79][80][81][82].
3 Bound states from space-like P 2 One of the consequences of working with Euclidean spacetime is that access to time-like quantities, such as masses of bound states, requires an analytic continuation of the component Green functions to complex momenta. While this is only a minor technicality thanks to many well-established techniques in the literature [22,49,[83][84][85], there are situations in which existing methods do not apply, or which are simply too complicated to implement. In this case, indirect methods can be employed that enable access to a limited number of time-like quantities [86,87].
In the next two sections we describe two techniques that have been widely used, and compare their performance in cases where direct analytic continuation is possible. This provides an estimate of the methods applicability to the study at hand.

Eigenvalue extrapolation
There are several means by which the mass spectrum of the BSE can be obtained. The most often used is through solution of Eq. (1), written as a matrix equation for simplicity This has solutions at discrete values of the bound state's total momentum P 2 = −M 2 i . By introducing the function λ(P 2 ) on the right, we obtain an eigenvalue equation whose boundstate solution correspond to λ(P 2 ) = 1.
Since λ(P 2 ) is a continuous function of P 2 , one can conceive that its continuation from space-like P 2 > 0 to timelike P 2 < 0 may be obtained through appropriate function fitting and extrapolation. The transformation of the eigenvalue g (λ) = 1 − 1/λ, see Ref. [26], removes a considerable degree of intrinsic curvature in the region close to the pole, rendering simple linear extrapolation viable provided the extrapolation is not far.
In the top panel of Fig. 7 we show the eigenvalue extrapolation of λ(P 2 ) for various J PC states. The data is first transformed via g (λ), before a linear fit f (P 2 ) = a + bx is performed. Finally, we plot the inverse function of g, λ fit = g −1 ( f (P 2 )) as solid lines. Exact results, obtained via calculation in the complex plane, are included as labeled points.

Inverse vertex function extrapolation
The second means to obtain the mass spectrum employs instead the inhomogeneous BSE for the vertex function i The obvious difference between this and the homogeneous BSE is the inhomogeneous term (0) i . Its introduction leads to several important changes to the solution. Setting the relative momentum p to zero, for convenience, we observe the appearance of poles as one approaches the bound state P 2 ∼ −M 2 . Then the determination of a bound-state mass is reduced to looking for zeros in 1/V i . Typically, the leading amplitude is used as point of reference, and one employs the method of biconjugate gradient (stabilized) for solution.
Restricting ourselves to space-like momenta P 2 requires once more the use of fit functions and extrapolation. Here, the most useful are rational polynomials Table 1 Results for vertex pole extrapolation for QCD rainbow-ladder in the chiral limit, compared with the result computed through direct analytic continuation. All units are in MeV. The points P 2 are taken from the region (0, L); the errors on extrapolated results come from the fitting procedure

J PC
Calc R (2,2) (L = 0.5) R (2,2) Note, that since the coefficients a i , b i are obtained through least-squares fitting, the resulting function is not a true Padé approximant. Regardless, the procedure appears quite reliable as can be seen in the bottom panel of Fig. 7. We summarize our results for vertex pole approximation in Table 1. The results obtained with eigenvalue extrapolation are not quoted as the method performs rather poorly, especially in the 1 ++ channel (see top panel of Fig. 7). In either of the extrapolation techniques there are two principal sources of uncertainty for the mass values. One comes from the fitting procedure, since the fit function coefficients (a i , b i of Eq. (10) for vertex pole method) come with their own error bars. These errors are straightforward to quantify, and the resulting uncertainties for the meson masses are quoted in parentheses in Table 1.
A second source of errors has to do with the applicability of the extrapolation procedure, as one would expect the whole method to become less reliable as one probes deeper into the P 2 < 0 region (i.e. the technique is less reliable for heavier mesons). Although it is very hard to quantify this, a comparison with exact results suggests that these effects are quite small for the inverse vertex approximation. In light of other systematic errors, present in both the continuum and lattice investigations of the SU(2) gauge theory, we will ignore this uncertainty in Sect. 4. As an additional check on the extrapolation method, we performed calculations with different fit ranges for P 2 , with the total momentum sampled in the region (0, L) (in GeV 2 ), and with L given in the table.
In the next section we employ the method with L = 0.5, which appears empirically to have the best performance.

Estimation of model dependence
As already highlighted, the majority of model dependence stems from the truncation of the quark-gluon vertex DSE. Other parts of the calculation are constrained either by the underlying gauge dynamics (i.e. the ghost and gluon propagator which are taken from appropriate lattice or continuum calculations) or by chiral symmetry (in the process of truncating the BSE kernel). Thus, we can test the sensitivity of the truncation by varying the solution of the quark-gluon vertex within the constraints imposed by chiral symmetry breaking and the axWTI.
The natural step is to dress the three-gluon vertex. This is motivated by both the 3PI formalism [88] and through the effective resummation of neglected diagrams in the full DSE for the quark-gluon vertex. This in turn enables us to give a rough estimate as to the impact of including additional corrections on our results. It is sufficient to describe the full three-gluon vertex in Landau gauge by its tree structure and one function of a symmetric variable The dressing function A(s 0 ) is obtained by solving the three-gluon vertex DSE in a 'ghost triangle' approximation, depicted in Fig. 9. The details of the calculation can be found in [63]. The resultant dressing function is shown in Fig. 8. Information available from continuum non-perturbative studies of the three-gluon vertex [63,[89][90][91] suggests that both the truncation of Fig. 9 and the restriction of possible tensor structures to the tree-level term, provide a reasonable phenomenological description of this Green function. The effect which the dressed three-gluon vertex has on the hadron masses can be seen in Table 2.
There is one further extension to our model that is possible, which is the inclusion of the so-called Abelian diagram in the quark-gluon vertex DSE, see Ref. [46]. This introduces no complications in the evaluation of the quark- Fig. 9 The truncated DSE for the three-gluon vertex. To ensure that bose-symmetry is maintained the right-hand side is averaged over all cyclic permutations Table 2 Chiral limit results for meson masses in rainbow-ladder (RL) and beyond rainbow-ladder (BRL) truncations, compared with lattice data for an SU(2) theory. All units are in TeV. Errors of the BRL results come from the extrapolation procedure. For the 0 ++ state, our continuum result is for an isoscalar; lattice results are forthcoming  (10) 3.08 (10) 3.3 ± 0.7 gluon vertex itself, and through the 'cutting' procedure it is straightforward to construct a solvable BSE kernel which is consistent with axWTI [37]. This BSE kernel would contain diagram with a new topology-the so-called crossed ladder diagram-which increases the algebraic and numerical effort considerably. However, in previous calculations the Abelian contribution has been shown to have a small effect on meson masses, typically less than 2% [92], which would similarly apply to our present investigation. For these reasons, and in light of other uncertainties of continuum and lattice investigations, we feel that it is justified to ignore this extension for now.

Discussion
Comparison of our results with the lattice [12,13] requires the scale to be set by equating the electroweak (EW) scale with the pseudoscalar meson ('pion') decay constant, i.e. v EW = f π = 246 GeV. This puts the theory under investigation in the context of dynamical EW symmetry breaking, otherwise known as Technicolor (TC) [93,94]. The drawbacks that the SU(2) model discussed here (and any other model with QCD-like dynamics) faces as a Technicolor template are by now well known. These include the problems with precision tests on flavor-changing neutral currents [95], and the composite 'Higgs boson' which is expected to be very heavy. This latter problem is seen here, whereupon we do not find an isoscalar scalar ('sigma') meson (a TC version of the Higgs boson) below 1.33 TeV. This situation, however, might change drastically if one con-siders explicitly the couplings to Standard Model particles [96], or more general EW embeddings [97]. Another promising approach to Technicolor phenomenology is to use nearly conformal theories as Technicolor templates [98][99][100][101]. As we are presently concerned with the QCD-like aspects of the model under investigation, we will not comment on its possible Technicolor applications further.
Ground state masses for various J PC mesons are shown for both the rainbow-ladder (RL) and beyond rainbow-ladder (BRL) truncations in Table 2, where they are additionally compared with the relevant lattice calculations. RL results were obtained by means of direct analytic continuation, while those of BRL were extrapolated from the region of space-like P 2 via the inverse vertex function. The pion decay constant, which is used to set the scale of the calculation, is evaluated via the relation [102]: where k ± = k ± P/2 and π is the pion BSE amplitude normalized according to the Nakanishi condition [103]. In QCD, the conventions employed in the above equation would correspond to the value f π = 93 MeV. When working in the region of space-like P 2 , the definition of Eq. (12) can be used without approximations only in the chiral limit, since the pion amplitude π can be obtained for the case P 2 → 0. For non-chiral quarks, and thus non-vanishing pion mass, the calculation would have to be set up for complex total momentum, which is a formidable task in a BRL setting [104]. For the discussion of results it would be useful to have an estimate on the mass of an isoscalar scalar meson, calculated in a method different from our DSE/BSE approach. Since the lattice results for this particle are yet to come, we will use the values obtained by means of group theory scaling, which for the model under investigation gives m σ ∈ [1, 1.5] TeV [96]. Taking this into consideration, it seems that the RL method fares well for the sigma meson, and to a lesser extent, the rho meson. In the 1 ++ channel, this truncation performs inadequately, with a result which deviates by about 30 % from the central lattice value. It is arguable whether or not one can modify the RL method so that it is better suited for an SU(2) theory, thus performing reasonably well for all considered mesons. Based on our current results, and given the limitations of the RL framework, we are skeptical toward this prospect.
On the other hand, the results of the BRL approach compare well with lattice, especially when employing the dressed three-gluon vertex. Since there are considerable error margins present in both the continuum and the lattice investigations, stronger statements about the agreement of our methods will have to wait for more refined calculations. Regarding the continuum calculation, dressing of the three-gluon vertex seems to lead to a better agreement with the discretized approach, but the overall impact of this modification is relatively mild, and all meson masses are rather robust in this respect. This leaves open the possibility that more elaborate modifications of our model (i.e. inclusion of additional diagrams and higher n-point Green functions in the quark-gluon vertex DSE) might not change the results appreciably. However, note that dynamical contributions that can collectively be termed 'pion cloud' effects are known to be important, and they are the focus of present and future investigations.
Apart from the chiral limit study, we also performed calculations for non-vanishing current quark masses. In Fig. 10 we demonstrate the validity of Gell-Mann-Oakes-Renner (GMOR) relation in the BRL approach, while in Fig. 11 we plot the masses of spin one mesons (in units of chiral limit f π ) as a function of current quark mass. Both plots correspond to a calculation with a bare three-gluon vertex. The results shown in Fig. 11 seem to compare well with the ones shown in Fig. 6 of [13]: however, a direct comparison is not possible since we don't have enough information to relate our m q to the ones employed in [13].
As a final remark, we note that the calculation of the baryonic spectrum in this theory does not require any additional effort. An SU(2) gauge theory possesses an enlarged (Pauli-Gürsey) flavor symmetry, which implies that chiral multiplets will contain both mesons and baryons (diquarks). In other words, a meson with J P quantum numbers will be degenerate with a J −P diquark. This degeneracy (which breaks down if the chemical potential is raised above some critical value μ c ) has been confirmed in numerous lattice investigations [2,4,5,12,13].

Conclusions and outlook
We presented a Dyson-Schwinger/Bethe-Salpeter calculation of ground state hadron masses in a theory with two colors and two fundamentally charged Dirac fermions. We employed a novel beyond rainbow-ladder method and obtained good agreement with lattice results for spin one mesons: however, improved calculations will be needed to reduce uncertainties in both lattice and continuum approaches.
For J = 0 mesons, we demonstrated that chiral dynamics is satisfied (i.e. the GMOR relation holds) and obtained the mass of the sigma meson to be in good agreement with the analysis based on group theory scaling. Additionally, we showed that the rainbow-ladder method performs unsatisfactorily in this strongly interacting template. This underlines the need to use more sophisticated techniques when studying generic non-Abelian gauge theories.
Besides masses, the beyond rainbow-ladder approach we outlined here can also be used to study hadronic decays and form factors. A first step toward accessing these quantities is to extend the calculation to complex Euclidean momenta. However, the technical complications which arise are considerable and are subject to future investigation.