A Dipole formulas

One-loop corrections to gaugino (co-)annihilation into quarks in the MSSM

Abstract

We present the full supersymmetric QCD corrections for gaugino annihilation and co-annihilation into light and heavy quarks in the Minimal Supersymmetric Standard Model (MSSM). We demonstrate that these channels are phenomenologically relevant within the so-called phenomenological MSSM. We discuss selected technical details such as the dipole subtraction method in the case of light quarks and the treatment of the bottom quark mass and Yukawa coupling. Numerical results for the (co-)annihilation cross sections and the predicted neutralino relic density are presented. We show that the impact of including the radiative corrections on the cosmologically preferred region of the parameter space is larger than the current experimental uncertainty from Planck data.

pacs:
12.38.Bx,12.60.Jv,95.30.Cq,95.35.+d
1

I Introduction

Today there is striking evidence for the existence of a Cold Dark Matter (CDM) component in the universe, coming from a large variety of astronomical observations such as the rotation curves of galaxies, the inner motion of galaxy clusters, and the Cosmic Microwave Background (CMB), to name just a few. The Planck mission (1) has measured the CMB with previously unparalleled precision. These measurements, combined with the information from WMAP polarization data at low multipoles (2), allow to determine the dark matter relic density of the universe to

 ΩCDMh2=0.1199±0.0027, (1)

where denotes the present Hubble expansion rate in units of 100 km s Mpc.

The identification of the nature of CDM represents one of the biggest challenges for modern physics. One popular hypothesis is the existence of a new weakly interacting and massive particle (WIMP), which constitutes (at least a part of) the CDM. Besides the lack of direct experimental evidence, the biggest problem of this hypothesis is the fact that the Standard Model of particle physics (SM) does not contain a WIMP, since neutrinos are too light and can only form hot dark matter. This is a strong hint for physics beyond the Standard Model.

A well motivated example for an extension of the SM is the Minimal Supersymmetric Standard Model (MSSM). Under the assumption that a new quantum number, the so-called -parity, is conserved, the lightest supersymmetric particle (LSP) is stable. In many cases the LSP is the lightest of the four neutralinos , which is a mixture of the bino, wino, and two higgsinos, according to

 ~χ01 = Z1~B~B+Z1~W~W+Z1~H1~H1+Z1~H2~H2, (2)

and is probably the most studied dark matter candidate.

The time evolution of the neutralino number density is governed by a nonlinear differential equation, the Boltzmann equation (3)

 Missing \left or extra \right (3)

where the first term on the right-hand side containing the Hubble parameter stands for the dilution of dark matter due to the expansion of the universe. The second and third term describe the creation and annihilation of neutralinos. Both of these terms are proportional to the thermally averaged annihilation cross section . The creation is also proportional to the number density in thermal equilibrium , which for temperatures , being the lightest neutralino mass, is exponentially suppressed via

 neqχ∼exp{−mχT}. (4)

Therefore the creation rate drops to zero when the universe cools down. At some later point, the expansion of the universe will finally dominate over the annihilation, and the neutralino freezes out asymptotically.

Taking into account the possibility of co-annihilations between the neutralino and the other MSSM particles, the thermally averaged annihilation cross section can be written as (4); (5)

 ⟨σannv⟩=∑i,jσijvijneqineqχneqjneqχ, (5)

where the sum runs over all MSSM particles and , ordered according to etc.

As can be seen from Eq. (5), co-annihilations can occur not only if the LSP is involved, but also among several of its possible co-annihilation partners. However, depending on the exact MSSM scenario under consideration, not all of these contributions are numerically relevant. Indeed, by generalizing Eq. (4), the ratios of the occurring equilibrium densities are Boltzmann suppressed according to

 neqineqχ ∼ exp{−mi−mχT}. (6)

Consequently, only particles whose masses are close to can give sizeable contributions. In the MSSM, relevant particles can be light sfermions, in particular staus or stops, or other gauginos.

Once the Boltzmann equation for the total number density is solved numerically, the relic density is obtained via

 Ωχh2 = mχnχρcrit. (7)

Here, is the current neutralino number density after the freeze-out, obtained by solving the Boltzmann equation, and is the critical density of the universe. The theoretical prediction calculated in this way can be compared with the experimental data, i.e. the limits given in Eq. (1). This allows to identify the cosmologically preferred regions of the MSSM parameter space. The obtained constraint is complementary to information from collider searches, precision measurements, direct and indirect searches for CDM.

The standard calculation of the relic density is often carried out by a public dark matter code, such as micrOMEGAs (6) or DarkSUSY (7). Both of these codes evaluate the (co-)annihilation cross section at an effective tree level, including in particular running coupling constants and quark masses, but no loop diagrams. However, it is well known that higher-order loop corrections may affect the cross section in a sizeable way.

In order to ensure an adequate comparison with the very precise cosmological data, the uncertainties in the theoretical predictions have to be minimized. For a given supersymmetric mass spectrum, the main uncertainty on the particle physics side resides in the calculation of the annihilation cross sections , defined in Eq. (5), which govern the annihilation cross section and thus the relic density . It is the aim of the present work to improve on this point in the context of gaugino2 (co-) annihilation in the MSSM.

The impact of loop corrections on the annihilation cross section and the resulting neutralino relic density has been discussed in several previous analyses. The supersymmetric QCD (SUSY-QCD) corrections to the annihilation of two neutralinos into third-generation quark-antiquark pairs have been studied in Refs. (8); (9); (10). The corresponding electroweak corrections have been investigated in Refs. (11); (12); (13). Further studies are based on effective coupling approaches (14); (15), including the co-annihilation of a neutralino with a stau. SUSY-QCD corrections to neutralino-stop co-annihilation can be found in Refs. (16); (17); (18).

These analyses led to the common conclusion that radiative corrections are non-negligible in the context of relic density calculations, as they may influence the resulting theoretical prediction in a sizeable way. In particular, the impact of the corrections is in general larger than the experimental uncertainty of the WMAP or Planck data.

The aim of the present Paper is to extend the calculation of Refs. (8); (9); (10) to all gauginos in the initial and all quarks in the final state. We present the full corrections in supersymmetric QCD to the following annihilation and co-annihilation processes of gauginos into quark-antiquark pairs:

 ~χ0i~χ0j → q¯q, (8) ~χ0i~χ±k → q¯q′, (9) ~χ±k~χ±l → q¯q (10)

for , , and . The quark in Eq. (9) is the down/up-type quark of the same generation3 as the up/down-type quark . The corresponding Feynman diagrams at tree level are shown in Fig. 1.

This Paper is organized as follows: In Sec. II we specify the model framework, introduce our reference scenarios and discuss the phenomenology of gaugino (co-) annihilation. Sec. III contains technical details about the actual cross section calculation. We will discuss the subtleties of the dipole subtraction method for light quarks and the treatment of the bottom quark mass and Yukawa coupling. Aspects concerning the regularization and renormalization are kept rather short, as they can be found in Ref. (17). In Sec. IV we present our numerical results to illustrate the impact of the one-loop corrections on the cross section and the relic density, respectively. Finally, our conclusions are given in Sec. V.

Ii Phenomenology of gaugino annihilation and co-annihilation

Throughout this analysis, we work within the phenomenological MSSM (pMSSM)4, where the soft-breaking parameters are fixed at the input scale TeV according to the SPA convention (20). We choose to work with eleven free parameters, which are detailed in the following: The Higgs sector is fixed by the pole mass of the pseudoscalar Higgs boson , the higgsino mass parameter , and the ratio of the vacuum expectation values of the two Higgs doublets . The first and second generation squarks have a common soft-mass parameter , while the third generation squarks are governed by , the soft-mass parameter for the sbottoms and left-handed stops, and for the right-handed stops. All trilinear couplings are set to zero except for in the stop sector. In contrast to the three independent mass parameters in the squark sector, we have a single parameter for all sleptons. Finally, the gaugino and gluino sector is defined by the bino mass parameter , the wino mass parameter , and the gluino mass parameter . In the context of our analysis, the most interesting parameters are , , and , since they determine the decomposition of the neutralinos and charginos.

Within this setup, with the help of a scan over the parameter space, we have chosen three reference scenarios, which will be used to illustrate the numerical impact of the presented corrections. The corresponding input parameters as discussed above are listed in Tab. 1, while Tab. 2 summarizes the most important particle masses, mixings, and related observables.

We have used SPheno 3.2.3 (21) to obtain the physical mass spectrum from the given input parameters. The neutralino relic density and the numerical value of the branching fraction have been obtained using micrOMEGAs 2.4.1 (6) with the standard CalcHEP 2.4.4 (22) implementation of the MSSM. The only changes we introduced are that we have set as well as included a lower limit on the squark-width, which both do not influence the results concerning dark matter presented here, but will be relevant later in the discussion of the dipole subtraction method in Sec. III.2.

Our scenarios have been selected such that they fulfill the following constraints: In order to work with scenarios which are realistic with respect to the recent Planck measurements, we require the neutralino relic density to be in the vicinity of the limits given in Eq. (1). Let us note that we assume that the neutralino accounts for the whole amount of dark matter that is present in our universe. Moreover, we expect the relic density to be modified by our corrections to the (co-)annihilation cross section of the neutralino, so that we apply rather loose bounds at this stage.

Second, we require the mass of the lightest (“SM-like”) -even Higgs boson to agree with the observation at LHC,

 122 GeV≤mh0≤128 GeV, (11)

where we allow for a theoretical uncertainty of about 3 GeV on the value computed by SPheno. This uncertainty is motivated by higher-order corrections, which are at present not included in SPheno, see, e.g., Ref. (23). Finally, we impose the interval

 2.77⋅10−4≤BR(b→sγ)≤4.07⋅10−4 (12)

on the inclusive branching ratio of the decay . This corresponds to the latest HFAG value (24) at the confidence level.

As can be seen in Tab. 2, the selected scenarios fulfill the mentioned constraints within the required uncertainties. All channels with quark final states contributing to at least 0.1% of the total annihilation cross section are listed in Tab. 3, while in Tab. 4 we show the contributions of the different sub-channels, i.e. the different diagram classes shown in Fig. 1. We have grouped the contributions from -channel scalar exchange (contribution denoted ), the -channel contribution from vector boson exchange (), and the squark exchange in the - and -channels (). The contributions from the corresponding squared matrix elements are denoted , , and , while the interference terms are denoted by , , and . Note that negative numbers in Tab. 4 refer to destructive interferences.

In our scenario I, the dominant contribution to the total annihilation cross section is the co-annihilation between the LSP and the lighter chargino. The second most important channel is the co-annihilation between the two lightest neutralinos, while the pair-annihilation of the LSP is only the third most important channel. This hierarchy is explained as follows: First, as can be seen in Tab. 4, the dominant subchannels for this scenario are the exchange of a scalar in the -channel. More precisely, the value of is already large enough to favor bottom quarks in the final states due to the -enhanced bottom Yukawa coupling.

In the case of co-annihilation of the LSP with the second lightest neutralino, this process is mediated by the pseudoscalar Higgs-boson , whose mass GeV is relatively close to the total mass in the initial state, GeV. The same argument holds for the co-annihilation with the lighter chargino, which proceeds via the exchange of a charged Higgs boson ( GeV and GeV). Although these two processes are Boltzmann-suppressed, see Eq. (6), they are numerically more important than the LSP pair-annihilation, which is kinematically disfavored. Indeed, with GeV, the configuration is further away from the -resonance.

Finally, although they are kinematically even closer to the -resonance ( GeV), the pair annihilation of the lighter chargino or of the second lightest neutralino are highly suppressed by the Boltzmann factor of Eq. (6) and therefore numerically not relevant.

The main difference in our scenario II is the different setup in the Higgs sector. More precisely, the lower value of together with the higher wino mass modifies the composition of the lightest neutralino such that it has a larger higgsino component. Moreover, due to the smaller value of , down-type Yukawa couplings are less important such that the contribution of final states with top-quarks is larger as can be seen in Tab. 3. Co-annihilation of the LSP and the lighter chargino (proceeding through exchange in the -channel) remains the dominant contribution, while co-annihilation with the second lightest neutralino is less relevant in this scenario. This is explained by the different mass spectrum: The mass gap between and is larger than for scenario I, and the mass difference between and is even more important, leading to a stronger Boltzmann suppression.

Another difference with respect to scenario I arises in the LSP pair-annihilation. As can be seen in Tab. 4, the dominant subchannel in scenario II is the -channel exchange of a stop. This is due to the fact that the latter is much lighter than in scenario I ( GeV against GeV in scenario I). For the channel, the stop exchange is less relevant, since the kinematical configuration is very close to the resonance of the heavier -even Higgs boson, GeV and GeV.

Finally, scenario II features a non-negligible contribution from co-annihilation of the LSP with the third neutralino, proceeding mainly through scalar exchange in the -channel. As before, this is due to the kinematical situation close to the -resonance ( GeV). In both scenarios I and II, the kinematical configuration is such that resonances with - or -exchange are not relevant.

The phenomenology of scenario III is rather different. Here, the lightest neutralino is mainly higgsino-like and the two lightest neutralinos together with the lighter chargino are almost mass-degenerate. Consequently, LSP pair annihilation is negligible, while co-annihilations with the lighter chargino and the second lightest neutralino as well as chargino pair-annihilation are the dominant processes. Moreover, the configuration at this parameter point is such that final states with first and second generation quarks are dominating (contributions from third-generation quarks only amount to 12%).

In contrast to the first two reference points, scenario III is characterized by important contributions from the - and -channel diagrams. This is explained by the fact that, here, the lightest stop is relatively light as compared to the annihilating gauginos ( GeV). Consequently, the squark propagator is numerically less suppressed. This configuration also leads to strong destructive interferences between the squark exchange and the -channel contributions. Finally, in this scenario, Higgs exchanges are negligible, since the kinematical configuration is above the corresponding resonances (e.g., GeV). However, the exchange of a - or -boson in the -channel gives sizeable contributions due to the higgsino nature of the lightest gauginos.

Iii Technical details

The radiative corrections at include the one-loop diagrams shown in Fig. 2 as well as the real gluon emission diagrams shown in Fig. 3. The loop contributions give rise to ultraviolet (UV) and infrared (IR) singularities. While the former are removed via renormalization, i.e. the introduction of appropriate counterterms, the latter cancel when including also the real emission of gluons (25). Altogether, the cross section at next-to-leading order (NLO) in is given by

 σNLO = ∫2dσV+∫3dσR, (13)

where the virtual part () and the real emission part () are integrated over the two- and three-particle phase space, respectively. We describe the different parts in the following.

iii.1 Calculation of loops

Here, we focus on the virtual part of the next-to-leading order cross section. The calculation of the loop diagrams is carried out in the dimensional reduction scheme (), and therefore all ultraviolet and infrared divergences are regulated dimensionally. The virtual part of the cross section is rendered UV-finite by redefining the original parameters of the theory. In our case we introduce counterterms for each of the parameters in the strong sector and specify conditions in order to determine the counterterms so that all UV divergences vanish.

In contrast to our previous analyses of neutrino annihilation (9); (10), here we use a hybrid on-shell/ renormalization scheme where the parameters , and are chosen as input parameters along with the heavy quark masses and . The trilinear couplings and the bottom quark mass are defined in the renormalization scheme, whereas all remaining input masses are defined on-shell. For details of the renormalization scheme see Ref. (17).

In this work we have improved on several aspects of the renormalization scheme. First, we have improved on the determination of the bottom quark mass in the -scheme. The bottom quark mass in the renormalization scheme is determined from a Standard Model analysis of sum rules (26). The transformation of the mass to the appropriate bottom quark mass in the renormalization scheme within the MSSM requires several steps as follows:

 mSM,¯¯¯¯¯¯¯MSb(mb) (1)⟶ mSM,¯¯¯¯¯¯¯MSb(Q) (2)⟶ mSM,¯¯¯¯¯¯¯DRb(Q) (3)⟶ mMSSM,¯DRb(Q). (14)

In the first step we use the three-loop renormalization group evolution to obtain the mass of the bottom quark at a scale (27) 5. Then in the second step at the final scale, while still considering only the Standard Model, we transform the bottom quark mass from the to the renormalization scheme using the two-loop relation (28)

 mSM,¯¯¯¯¯¯¯DRb(Q) = mSM,¯¯¯¯¯¯¯MSb(Q)[1−αeπ14CF+(α¯¯¯¯¯¯¯MSsπ)211192CACF−α¯¯¯¯¯¯¯MSsπαeπ(14C2F+332CACF) (15) +(αeπ)2(332C2F+132CFTnf)+…].

This transformation involves an evanescent coupling , which is identical to the strong coupling constant in a supersymmetric theory, whereas in QCD there is a subtle difference. The dots indicate higher order terms, which are not relevant here. The constants and are the usual color factors of QCD, and , respectively. The difference between the couplings to one-loop order is (29)

 αe = α¯¯¯¯¯¯¯DRs[1+α¯¯¯¯¯¯¯DRsπ{−TFLt2+CA4(2+L~g+∑i=1,2(L~g−L~qi)m2~qim2~g−m2~q1) (16) +CF4(∑i=1,2⎛⎝−1−2L~g+2L~qi+(−L~g+L~qi)m2~qim2~g−m2~qi⎞⎠m2~qim2~g−m2~qi+(−3−2L~g))}],

where stands for the strong coupling constant in the renormalization scheme in the MSSM and we have used the shorthand notation .

The last step requires adding threshold corrections from supersymmetric particles in the loop. Using the results of Ref. (29), we can write the final transformation as

 mSM,¯¯¯¯¯¯¯DRb(Q) = ζmbmMSSM,¯DRb(Q), (17)

where the coefficient can be expanded in the strong coupling constant,

 ζmb = 1+(α¯¯¯¯¯¯¯DRsπ)ζ(1)b+(α¯¯¯¯¯¯¯DRsπ)2ζ(2)b+O(α3s). (18)

Using the results for the coefficients and given in Ref. (29) and inverting Eq. (17) yields the final bottom quark mass in the renormalization scheme in the MSSM which is subsequently used in our analysis.

The bottom quark mass deserves our special attention as the Higgs exchange is the leading contribution to many of the (co-)annihilation cross sections in this analysis. Therefore the second improvement of this study with respect to our previous works involves an improvement to the Yukawa coupling of the bottom quark. Similar to our earlier analyses, we improve our full one-loop SUSY-QCD by including higher-order QCD and SUSY-QCD corrections to the Yukawa coupling.

Leading QCD and top-quark induced corrections to the Yukawa coupling of Higgs bosons to bottom quarks were calculated up to (30) and can be used to define an effective Yukawa coupling which includes these corrections as

 [h¯¯¯¯¯MS,QCD,Φb(Q)]2=[h¯¯¯¯¯MS,Φb(Q)]2[1+ΔQCD+ΔΦt], (19)

for each Higgs boson . The QCD corrections are explicitly given by

 ΔQCD =αs(Q)πCF174+α2s(Q)π2[35.94−1.359nf] +α3s(Q)π3[164.14−25.77nf+0.259n2f] (20) +α4s(Q)π4[39.34−220.9nf+9.685n2f−0.0205n3f],

and the top-quark induced corrections for each Higgs boson read

 Δht=ch(Q) [1.57−23logQ2m2t+19log2m2b(Q)Q2], (21) ΔHt=cH(Q) [1.57−23logQ2m2t+19log2m2b(Q)Q2], (22) ΔAt=cA(Q) [236−logQ2m2t+16log2m2b(Q)Q2], (23)

with

 Missing or unrecognized delimiter for \Biggr (24)

We exclude the one-loop part of these corrections as it is provided consistently through our own calculation.

In the MSSM, the Yukawa coupling to bottom quarks can be enhanced for large or large , and then effects even beyond the next-to-leading order should be included. Therefore, in addition, we incorporate these corrections that can be resummed to all orders in perturbation theory (31); (32). Denoting the resummable part by , we redefine the bottom quark Yukawa couplings as

 hMSSM,hb(Q) = h¯¯¯¯¯MS,QCD,hb(Q)1+Δb[1−Δbtanαtanβ], (25) hMSSM,Hb(Q) = h¯¯¯¯¯MS,QCD,Hb(Q)1+Δb[1+Δbtanαtanβ], (26) hMSSM,Ab(Q) = h¯¯¯¯¯MS,QCD,Ab(Q)1+Δb[1−Δbtan2β]. (27)

As a further improvement, we include also the leading NNLO contributions to as given in Ref. (33). These corrections are now also available for general, i.e. non-minimal sources of flavour violation (34). The electroweak one-loop corrections have also been computed and resummed analytically to all orders (35). We leave their implementation to future work.

Let us finally briefly comment on possible Sommerfeld enhancement effects of gaugino (co-)annihilation in the MSSM. As is well known, potentially large loop corrections arise from long-range interactions of WIMPs before their annihilation, which are mediated by bosons with masses well below the WIMP mass . More precisely, these corrections become relevant when the Bohr radius is smaller than the interaction range or . This condition is almost never realized in the MSSM with an LSP mass of TeV, and weak gauge and Higgs boson masses GeV (36). In particular, it is never realized in our scenarios I–III. Note, however, that in scenario III the pair annihilation of charginos, which are almost mass degenerate with the LSP, contributes about 10% to the total cross section, so that this (subleading) channel would indeed be enhanced through multiple exchanges of massless photons if resummed to all orders (37); (38). These calculations, which involve either numerical solutions of coupled Schrödinger equations or analytic calculations in an effective non-relativistic theory, are postponed to future work.

iii.2 Dipole subtraction method

As mentioned above, the real gluon emission shown in Fig. 3 needs to be included in order to cancel the remaining infrared (IR) singularities in the virtual part of the cross section (25). However, this is not as straightforward as in the ultraviolet case, since the two contributions reside in the differential cross sections and , which are integrated over different phase spaces. Moreover, working in dimensions, the soft and collinear divergencies appearing in the virtual contribution can be explicitly isolated and appear as single and double poles, and , while the divergencies in the real corrections arise from the phase space integration over the gluon phase space. In addition, quasi-collinear divergencies can appear in including large logarithmic corrections of the form , which cancel against logarithms of the same form in .

For these reasons, and generally speaking, a separate numerical evaluation of the two phase-space integrations in Eq. (13) cannot lead to numerically stable results. There are two approaches to render both of these terms separately infrared and collinear safe and therefore numerically evaluable: The so-called phase-space slicing method (39) and the dipole subtraction method (40); (41); (42). In the present work, we shall use the latter, which we will describe in the following.

The dipole subtraction method renders the integrands in Eq. (13) seperately finite by adding and subtracting an auxiliary cross section . Using dimensional regularization, this is done according to

 Missing or unrecognized delimiter for \Big (28)

where in the last term on the right-hand side the three-particle phase-space integral is factorized into the two-particle phase-space integral of and the integration over the one-particle phase-space of the radiated gluon. The auxiliary cross section , acting as a local counterterm for , has to possess the same pointwise singular behavior as and has to be analytically integrable over the gluon phase space in dimensions. Then, on the one hand, reproduces the potentially soft or collinear singular terms in the real corrections, such that one ends up with a convenient form for numerically performing the three-particle phase-space integration in Eq. (28). On the other hand, cancels all single and double poles appearing in in a way that the sum is rendered finite even in the limit . In addition, can be written in such a way, that it also cancels all quasi-collinear divergencies.

The dipole contributions to the matrix elements of real corrections in the case of final state radiation can be written in the general form

 ∣∣MR∣∣2=∑i,j∑k≠i,jDijk+⋯=Dgq,¯q+Dg¯q,q+…. (29)

This expression encodes the singular structure of the real radiation matrix element as a summation over so-called emitter-spectator pairs, singled out over the two Born-level external particles in all possible ways, and the dots stand for further infrared and collinear finite terms. Here, and run over the final state particles connected to the emitter through a splitting process as depicted in Fig. 4, and stands for the spectator particle, which is needed to maintain conservation of gauge-group charges and total momentum.

The general structure of the associated matrix element of can then be rewritten as

 ∣∣MA∣∣2 = ∑i,j∑k≠i,jDijk (30) = ∑i,j∑k≠i,jVij,k(pi,~pij,~pk)⊗∣∣MB(~pij,~pk)∣∣2.

The universal product form on the right hand side mimics the factorization of in the soft and collinear limit. It encodes the two-step process of the Born-level production of an emitter-spectator pair with momenta and followed by the decay of the emitter described by as represented by the box in Fig. 5. The are matrices in color and helicity-space of the emitter and the symbol stands for phase space convolution and possible helicity and color sums between and the exclusive Born-level matrix element . They become proportional to the Altarelli-Parisi splitting functions in the collinear region and to eikonal factors in the soft region (42); (43).

In addition, Eq. (30) allows for a factorizable mapping of the three-particle phase-space spanned by , , and onto the two-particle phase-space represented by the emitter and spectator momenta and times the single gluon phase space spanned by as implied by Eq. (28). It takes the form

 ~pμk = λ1/2(Q2,m2ij,m2k)λ1/2(Q2,(pi+pj)2,m2k)(pμk−Q⋅pkQ2Qμ) (31) +Q2+m2k−m2ij2Q2Qμ

with

 ~pμij=Qμ−~pμk,Qμ=pμi+pμj+pμk, (32)

and the emitter mass . For further definitions see App. A. Since the only describe the emitter decay and are essentially independent of the Born-level cross section part, they need to be calculated only once. In the case of (anti)quark-gluon splitting, the are given by

 ⟨s∣∣Vij,k∣∣s′⟩ = 2g2sμ2ϵ CF{21−~zj(1−yij,k)−~vij,kvij,k[1+~zj+m2qpipj+ϵ(1−~zj)]}δss′ = ⟨Vij,k⟩δss′ (33)

with , (or ) and (or ). and are the emitter-spins, is the strong coupling, and for . The dimensional regularization scale , which drops out of the final result , is set to be equal to the renormalization scale. For further definitions see App. A. The auxiliary cross section can then be constructed using Eq. (33) together with Eq. (30).

The virtual dipole contributions of Eq. (28) can be rewritten as

 ∫2[dσV+∫1dσA]ϵ=0=∫2[dσV+dσB⊗I2(ϵ,μ2;{pa,ma})]ϵ=0, (34)

where the index runs over the two Born-level final states. Following Eq. (30), due to the phase-space mapping denoted in Eqs. (31) and (32), explicitly depends on the gluon phase-space spanned by only via the universal factor . Therefore this integration can be performed analytically once and for all as given in Eq. (34). Following again Ref. (41) for final state radiation, the function can then be written as

 I2(ϵ,μ2;{pa,ma}) = −g2s8π2(4π)ϵΓ(1−ϵ)∑j1T2j∑k≠jTj⋅Tk[T2j(μ2sjk)ϵ(Vj(sjk,mj,mk;ϵ)−π23) (35) +Γj(μ,mj;ϵ)+γjlog(μ2sjk)+γj+Kj+O(ϵ)],

where , the are the color charges in the representation of the associated particle, and and run over all possible emitter-spectator combinations as in Eq. (29). For the expressions of , and see App. A. The function can be further decomposed into a ()-symmetric and singular (S) and a non-symmetric and non-singular (NS) part

 Vj(sjk,mj,mk;ϵ) = V(S)j(sjk,mj,mk;ϵ) (36) +V(NS)j(sjk,mj,mk).

Dealing with light quarks, numerical instabilities of the relevant dipoles are encountered if the mass of a final state quark drops far below the initial state energy, e.g., for first generation quarks. Therefore we take the light quarks (, and ) as massless. That makes it necessary to also take into account the associated and , which are given in App. A for all possible combination of emitter-spectator masses.

These are all the parts needed for calculating the relevant dipole contributions. The concrete calculation of the various elements has been performed in a semi-automatic way using FeynCalc (44) for simplifying the Dirac algebra. Furthermore, since we are dealing with initial-state neutralinos which are Majorana fermions, one needs to take care of handling the fermion flow correctly. This has been done following Ref. (45).

In addition, we have introduced a minimal width for squark propagators to avoid numerical instabilities due to quasi-collinear singularities in the - and -channels. This minimal value is set to 10 GeV. For our typical scenarios we have explicitly checked that on the one hand this renders the three-particle phase-space integration stable, but on the other hand has no relevant impact on the relic density.

Iv Numerical results

iv.1 Impact on the cross section

In this Section we discuss the impact of our full corrections on the cross sections of the processes in Eqs. (8) – (10). In Fig. 6 we present the cross sections of the most relevant gaugino (co-)annihilation channels of Tab. 3 against , the momentum in the center-of-mass frame. We show the cross sections at tree-level (black dashed line), at one loop (blue solid line), and the corresponding value obtained with micrOMEGAs/CalcHEP (orange solid line). Also shown as gray shaded areas are the Boltzmann velocity distributions of the dark matter particles (in arbitrary units). Their maxima may coincide with the maximal cross sections (induced, e.g., by Higgs resonances as in the upper two plots), but depending on the particle masses the maximal cross section often also sits on the shoulder of the Boltzmann distribution (cf. the two central plots). The lower parts of each plot show the different ratios between the three cross sections (second item in the legends).

The upper left plot of Fig. 6 shows the cross sections for the channel for scenario I. The most striking feature is the large