The One-Loop Correction to Heavy Dark Matter Annihilation
Grigory Ovanesyan, Nicholas L. Rodd, Tracy R. Slatyer, Iain W. Stewart
MMIT-CTP 4852
The One-Loop Correction to Heavy Dark Matter Annihilation
Grigory Ovanesyan, Nicholas L. Rodd, Tracy R. Slatyer, and Iain W. Stewart Physics Department, University of Massachusetts Amherst, Amherst, MA Center for Theoretical Physics, Massachusetts Institute of Technology, Cambridge, MA
We calculate the one-loop corrections to TeV scale dark matter annihilation in a model where thedark matter is described by an SU(2) L triplet of Majorana fermions, such as the wino. We use thisframework to determine the high and low-scale MS matching coefficients at both the dark matterand weak boson mass scales at one loop. Part of this calculation has previously been performedin the literature numerically; we find our analytic result differs from the earlier work and discusspotential origins of this disagreement. Our result is used to extend the dark matter annihilationrate to NLL (cid:48) (NLL+ O ( α ) corrections) which enables a precise determination of indirect detectionsignatures in present and upcoming experiments. I. INTRODUCTION
It is now well established that if dark matter (DM)is composed of TeV scale Weakly Interacting MassiveParticles (WIMPs) then its present day annihilation rateto produce photons is poorly described by the tree-levelamplitude. Correcting this shortcoming is important fordetermining accurate theoretical predictions for existingand future indirect detection experiments focussing onthe TeV mass range, such as H.E.S.S [1, 2], HAWC [3–5],CTA [6], VERITAS [7–9], and MAGIC [10, 11].The origin of the breakdown in the lowest order ap-proximation can be traced to two independent effects.The first of these is the so called Sommerfeld enhance-ment: the large enhancement in the annihilation crosssection when the initial states are subject to a long-rangepotential. In the case of WIMPs this potential is due tothe exchange of electroweak gauge bosons and photons.This effect has been widely studied (see for example [12–16]) and can alter the cross section by as much as severalorders of magnitude. The Sommerfeld enhancement isparticularly important when the relative velocity of theannihilating DM particles is low, as it is thought to be inthe present day Milky Way halo.The second effect is due to large electroweak Sudakovlogarithms of the heavy DM mass, m χ , over the elec-troweak scale, which enhance loop-level diagrams andcause a breakdown in the usual perturbative expansion.The origin of these large corrections can be traced to thefact that the initial state in the annihilation is not anelectroweak gauge singlet, and that a particular γ or Z final state is selected, implying that the KLN theoremdoes not apply [17–20]. While the importance of this ef-fect for indirect detection has only been appreciated morerecently (see for example [21–26]), it must be accountedfor, as it can induce O (1) changes to the cross section.Hryczuk and Iengo [21] (hereafter HI) calculated the one-loop correction to the annihilation rate of heavy winos to γγ and γZ , and found large corrections to the tree-levelresult, even after including a prescription for the Som-merfeld enhancement. These large corrections are symp-tomatic of the presence of large logarithms ln(2 m χ /m Z )and ln(2 m χ /m W ), which can generally be resummed us- ing effective field theory (EFT) techniques. This obser-vation has been made by a number of authors who intro-duced EFTs to study a variety of models and final states.The list includes the case of exclusive annihilation into γ or Z final states for the standard fermionic wino [24]and also a scalar version of the wino [23], as well as semi-inclusive annihilation into γ + X for the wino [22, 25, 26]and higgsino [26].In principle the EFT calculations are systematicallyimprovable to higher order and in a manner where theperturbative expansion is now under control. In order tofully demonstrate perturbative control has been regained,however, it is important to extend these works to higherorder. To this end, in this paper we extend the calcu-lation of exclusive annihilation of the wino, which hasalready been calculated to next-to-leading logarithmic(NLL) accuracy [24]. Doing so includes determining theone-loop correction in the full theory, as already consid-ered in HI. Nonetheless the results in that reference werecalculated numerically and are not in the form neededto extend the EFT calculation to higher order. As such,here we revisit that calculation and analytically deter-mine the DM-scale (high-scale) one-loop matching coef-ficients. We further calculate the electroweak-scale (low-scale) matching at one loop, thereby including the effectsof finite gauge boson masses. Taken together these twoeffects extend the calculation to NLL (cid:48) = NLL + O ( α )one-loop corrections, where α = g / π and g is theSU(2) L coupling. We estimate that our result reduces theperturbative uncertainty from Sudakov effects to O (1%),improving on the NLL result where the uncertainty was O (5%). Our calculation is complementary to the NLL (cid:48) calculation for the scalar wino considered in [23], andwhere relevant we have cross checked our work againstthat reference. In Sec. II we outline the EFT setup andreview the NLL calculation. Then in Sec. III we statethe main results of this work, the one-loop high and low-scale matching, leaving the details of their calculation toApp. A and App. C respectively. Detailed cross checkson the results are provided in App. B and App. D, whilstlengthy formulae are delayed till App. E. We compare ouranalytic results to the numerical ones of HI in Sec. IV andthen conclude in Sec. V. a r X i v : . [ h e p - ph ] D ec II. THE EFT FRAMEWORK
We begin by outlining the EFT framework for our cal-culation, and in doing so review the calculation of heavyDM annihilation to NLL, focussing on the treatment ofthe large logarithms that were partly responsible for thebreakdown in the tree-level approximation. We choosethe concrete model of pure wino DM – the same as usedin HI and [24] – to study these effects. Nevertheless weemphasise the point that the central aim is to quantifythe effect of large logarithms which can occur in manymodels of heavy DM, rather than to better understandthis particular model. Ultimately it would be satisfyingto extend these results to DM with arbitrary charges un-der a general gauge group to make the analysis less modelspecific. This is possible for GeV scale DM indirect de-tection where the tree-level approximation is generallyaccurate (see for example [27, 28]). Understanding thefull range of effects first in a simple model is an importantstep towards this goal.The model considered takes the DM to be a wino: anSU(2) L triplet of Majorana fermions. As already high-lighted, this is a simple example where both the Som-merfeld enhancement and large logarithms are impor-tant. Furthermore this model is of interest in its ownright. Neutralino DM is generic in supersymmetric theo-ries [29, 30]; models of “split supersymmetry” naturallyaccommodate wino-like DM close to the weak scale, whilethe scalar superpartners can be much heavier [31–33].DM transforming as an SU(2) L triplet has been studiedextensively in the literature, both within split-SUSY sce-narios [34–36] and more generally [14, 37, 38]. The modelaugments the Standard Model (SM) Lagrangian with L DM = 12 Tr ¯ χ (cid:0) i /D − M χ (cid:1) χ . (1)We take M χ = m χ I , such that in the unbroken theory allthe DM fermions have the same mass. After electroweaksymmetry breaking, the three states χ , , break into aMajorana fermion χ and a Dirac fermion χ + . A smallmass difference, δm , between these states is then gener-ated radiatively, ensuring that χ makes up the observedstable DM. Note, however, that both the charged andneutral states will be included in the EFT.An effective field theory for this model, NRDM-SCET,was introduced in [24] and used to calculate the rates forthe annihilation processes χχ → ZZ, Zγ, γγ . Specifi-cally the EFT generalizes soft-collinear effective theory(SCET) [39–42] to include non-relativistic dark matter(NRDM) in the initial state. Schematically the calcula-tion involves several steps. Firstly the full theory has tobe matched onto the relevant NRDM-SCET EW operatorsat the high scale of µ (cid:39) m χ . The qualifier EW indicatesthat this is a theory where electroweak degrees of free-dom – the W and Z bosons, top quark, and the Higgs– are dynamical, as introduced in [20, 43–46]. Theseoperators then need to be run down to the electroweakscale, µ (cid:39) m Z . At this low scale, we then match NRDM- SCET EW onto a theory where the electroweak degrees offreedom are no longer dynamical, NRDM-SCET γ . Thismatching accounts for the effects of electroweak symme-try breaking, such as the finite gauge boson masses. Atthis stage we can now calculate the low-scale matrix ele-ments which provide the Sommerfeld enhancement. Wenow briefly review each of these steps.The first requirement is to match NRDM-SCET EW and the full theory at the high scale µ m χ . The relevantoperators in the EFT to describe DM annihilation havethe following form: O r = 12 (cid:0) χ aTv iσ χ bv (cid:1) (cid:16) S abcdr B icn ⊥ B jd ¯ n ⊥ (cid:17) i(cid:15) ijk ( n − ¯ n ) k , (2)which is written in terms of the basic building blocksof the effective theory, and in the centre of momentumframe we can define v = (1 , , , n = (1 , ˆ n ), and¯ n = (1 , − ˆ n ) where ˆ n is the direction of an outgoinggauge boson. In more detail χ av is a non-relativistic two-component fermionic field of gauge index a correspondingto the DM and B ¯ n,n contain the outgoing (anti-)collineargauge bosons A µ ¯ n,n , which can be seen as B µn ⊥ = A µn ⊥ − k µ ⊥ ¯ n · k ¯ n · A µn + . . . , (3)where the higher order terms in this expression involvetwo or more collinear gauge fields. For B µ ¯ n ⊥ we simplyinterchange n ↔ ¯ n . The full form of B µn ⊥ can be found in[47], and is collinear gauge invariant on its own. Finallythe gauge index connection is encoded in S abcdr : S abcd = δ ab ( S cen S de ¯ n ) ,S abcd = ( S aev S cen )( S bfv S df ¯ n ) . (4)These expressions are written in terms of adjoint Wilsonlines of soft gauge bosons along some direction n , ¯ n , or v ; in position space the incoming Wilson line is S v ( x ) = P exp (cid:20) ig (cid:90) −∞ dsv · A v ( x + ns ) (cid:21) , (5)where the matrix A bcv = − if abc A av and for outgoing Wil-son lines the integral runs from 0 to ∞ .The fact there are only two possible forms of S abcdr means there are only two relevant NRDM-SCET opera-tors. An important requirement of the operators is thatthe incoming DM fields must be in an s -wave config-uration. Then being a two-particle state of identicalfermions, the initial state must be a spin singlet. Ifthe annihilation was p -wave or higher, it would be sup-pressed by powers of the low DM velocity relative tothese operators. The Wilson coefficients associated withthese operators are determined by the matching. Cal-culating to NLL only requires the tree-level result where C ( µ m χ ) = − C ( µ m χ ) = − πα ( µ m χ ) /m χ as an initialcondition. Here again α is the SU(2) L fine structureconstant. We extend this result to one loop in Sec. III.After matching, the next step is to evolve these op-erators down to the low scale, effectively resummingthe large logarithms ln(2 m χ /m Z ) and ln(2 m χ /m W ) thatcaused a breakdown in the perturbative expansion of thecoupling. This is done using the anomalous dimensionmatrix ˆ γ of the two operators (a matrix as the operatorswill in general mix during the running). In general thematrix can be broken into a diagonal piece γ W T , and anon-diagonal soft contribution ˆ γ S , asˆ γ = 2 γ W T I + ˆ γ S . (6)To NLL these results are given by [24]: γ W T = α π Γ g ln 2 m χ µ − α π b + (cid:16) α π (cid:17) Γ g ln 2 m χ µ , ˆ γ S = α π (1 − iπ ) (cid:18) − (cid:19) − α π (cid:18) (cid:19) . (7)Here the diagonal anomalous dimension has been writtenin terms of the SU(2) L one-loop β -function, b = 19 / g = 8 andΓ g = 8 (cid:0) − π (cid:1) , and we use the full SM particle con-tent for this evolution. Renormalization group evo-lution with the anomalous dimension also requires thetwo-loop β -function, and for this we take b = − / µdα /dµ = − b α / (2 π ) − b α / (8 π ). Below the DM matching scale,the spin of the DM is no longer important. As such theanomalous dimension determined in [24] for the fermionicwino should resum the same logarithms as those that ap-pear in the scalar case considered in [23], and we haveconfirmed they agree.We can then explicitly use the full anomalous dimen-sion to evolve the operators as follows: (cid:20) C X ± ( { m i } ) C X ( { m i } ) (cid:21) = e ˆ D X ( µ Z , { m i } )) P exp (cid:32)(cid:90) µ Z µ mχ dµµ ˆ γ ( µ, m χ ) (cid:33) × (cid:20) C ( µ m χ , m χ ) C ( µ m χ , m χ ) (cid:21) , (8)Let us carefully explain the origin and dependence of eachof these terms. Starting from the right, C and C arethe high-scale Wilson coefficients of the operators statedin Eq. (2), resulting from a matching of the full theoryonto NRDM-SCET EW . These only depend on the highscales, specifically µ m χ and m χ . Next the anomalousdimension ˆ γ is also a high scale object, and so only de-pends on m χ and now µ as it runs between the relevantscales. ˆ D X is a factor accounting for the low-scale match-ing from NRDM-SCET EW onto NRDM-SCET γ – a the-ory where the electroweak modes have been integratedout, see [20, 43–46]. It is a matrix as soft gauge boson This means we take m t ∼ m H ∼ m W,Z and integrate out allthese particles at the same time at the electroweak scale. exchanges can mix the operators. Furthermore ˆ D X islabelled by X to denote its dependence on the specific fi-nal state considered, γγ , γZ or ZZ . This object dependson the low-scale physics and so depends on µ Z and allthe masses in the problem, which we denote as { m i } .It contains both a resummation of low-scale logarithms(which can be carried out directly as in [43, 44] or moreelegantly with the rapidity renormalization group [48],see also [49]) as well as the low scale matching coefficientwhich does not necessarily exponentiate. Finally on theleft we have our final coefficients C X ± and C X , which asexplained below can be associated with the charged andneutral annihilation processes. In an all orders calcula-tion of all terms in Eq. (8), the scale dependence wouldcompletely cancel on the right hand side, implying that C X ± and C X depend only on the mass scales in the prob-lem and not µ m χ or µ Z . Nevertheless at any finite per-turbative order, the scale dependence does not cancelcompletely and so a residual dependence is induced inthese coefficients. We will exploit this to estimate theuncertainty in our results associated with missing higherorder terms.As we are performing a resummed calculation, the or-der to which we calculate is defined in terms of the largeelectroweak logarithms we can resum. In general thestructure of the logarithms can be written schematicallyas:ln CC tree ∼ ∞ (cid:88) k =1 (cid:104) α k ln k +1 (cid:124) (cid:123)(cid:122) (cid:125) LL + α k ln k (cid:124) (cid:123)(cid:122) (cid:125) NLL + α k ln k − (cid:124) (cid:123)(cid:122) (cid:125) NNLL + . . . (cid:105) , (9)where since Sudakov logarithms exponentiate, we havedefined the counting in terms of the log of the result.Furthermore all corrections are defined with respect tothe tree level result C tree ∼ O ( α ), which is a conven-tion we will follow throughout. With this definition ofthe counting, to perform the running in Eq. (8) to NLLorder, there are three effects that must be accounted for:1. high-scale matching at tree level; 2. two-loop cuspand one-loop non-cusp anomalous dimensions; and 3. thelow-scale matching at tree level, together with the ra-pidity renormalization group at NLL. To extend this toNNLL all three of these need to be calculated to one or-der higher. In between these two is the NLL (cid:48) result wepresent here, which involves determining both the highand low-scale matching at one loop. In terms of Eq. (9),this amounts to determining the leading k = 1 piece ofthe NNLL result. To the extent that O ( α ) correctionsare larger than those at O ( α ln( µ m χ /µ Z )), the NLL (cid:48) re-sult is an improvement over NLL and more importantthan NNLL.Before presenting the result of that calculation, how-ever, it is worth emphasising another advantage gainedfrom the effective theory. In addition to allowing us toresum the Sudakov logarithms, the effective theory alsoallows this problem to be cleanly separated from the issueof low-velocity Sommerfeld enhancement in the ampli-tude – in NRDM-SCET there is a Sommerfeld-Sudakovfactorization. At leading power the relevant SCET La-grangian contains no interaction with the DM field. Onthe other hand NRDM does contain soft modes, whichare responsible for running the couplings, however thesemodes do not couple the Sommerfeld potential to thehard interaction at leading power. Consequently matrixelements for the DM factorize from the matrix elementsof the states annihilated into. This allows for an all ordersfactorized formula for the DM annihilation amplitude inthis theory: M χ χ → X = 4 √ m χ P X (cid:104) s (cid:0) Σ X − Σ X (cid:1) + √ s ± Σ X (cid:105) , M χ + χ − → X = 4 m χ P X (cid:104) s ± (cid:0) Σ X − Σ X (cid:1) + √ s ±± Σ X (cid:105) . (10)Here X can be γγ , γZ or ZZ and P γγ = − e (cid:15) in ⊥ (cid:15) j ¯ n ⊥ (cid:15) ijk ˆ n k / (2 m χ ), whilst P γZ = cot ¯ θ W P γγ and P ZZ = cot ¯ θ W P γγ , with ¯ θ W the MS Weinberg angle.The key physics in this equation is that the contributionfrom Sommerfeld enhancement is captured in the terms s ij , whilst the contribution from electroweak logarithmsis in Σ Xi ; the two are manifestly factorized and can becalculated independently.The focus of the present work is to extend the calcu-lation of the Sudakov effects. In terms of the factorizedresult stated in Eq. (10) this amounts to an improved cal-culation of Σ Xi . Explicitly, from there we can see that: (cid:12)(cid:12) Σ X (cid:12)(cid:12) = σ (cid:26)(cid:26) SE χ + χ − → X σ tree χ + χ − → X , (cid:12)(cid:12) Σ X − Σ X (cid:12)(cid:12) = σ (cid:26)(cid:26) SE χ χ → X σ tree χ + χ − → X , (11)where (cid:8)(cid:8) SE denotes a calculation where Sommerfeld En-hancement is intentionally left out. To be even moreexplicit, we can write these Sudakov effects in terms ofthe Wilson coefficients in Eq. (8). Specifically we have:Σ X = C X ± C tree1 , Σ X − Σ X = C X C tree1 , (12)where as stated above C tree1 = − πα /m χ . III. THE ONE-LOOP CORRECTION
In this section we discuss the main results of this work,which includes analytic expressions for both the high andlow scale matching coefficients in the language introducedin the previous section. We start with reporting the re-sult of the calculation of the high-scale Wilson coeffi-cients C r to one loop. The details have been eschewed toApp. A. In short this calculation involves enumeratingand evaluating the 25 one-loop diagrams that mediate χ a χ b → W c W d in the unbroken full theory and thenmatching this result onto the NRDM-SCET EW opera-tors. For example, we evaluate diagrams such asand provide the analytic expression graph by graph. Herethe solid lines are DM particles and wavy lines are elec-troweak gauge bosons in the full theory above the DMscale. In addition we account for the counter term contri-bution, the change in the running of the coupling throughthe matching, and also ensure that the calculation main-tains the Sudakov-Sommerfeld factorization. Combiningall of these we find C ( µ ) = − πα ( µ ) m χ + α ( µ ) m χ (cid:20) µ m χ +2 ln µ m χ + 2 iπ ln µ m χ + 8 − π (cid:21) ,C ( µ ) = πα ( µ ) m χ − α ( µ ) m χ (cid:20) ln µ m χ +3 ln µ m χ − iπ ln µ m χ − π (cid:21) . (13)Here and throughout this section α ( µ ) is the couplingdefined below the scale of the DM mass, m χ . We explainthis distinction carefully in App. A. For each coefficientin Eq. (13) the first term represents the tree-level con-tribution. A cross check on this result is provided inApp. B, where we check that the µ dependence of thisresult properly cancels with that of the NLL resumma-tion from [24] for the O ( α ) corrections. The cancellationoccurs between our result in Eq. (13) and the running in-duced by the anomalous dimension stated in Eqs. (6) and(7); this can be seen clearly in Eq. (8) as these are theonly two objects that depend on µ m χ . As the anomalousdimension is independent of the DM spin, the logarithmsappearing in our high-scale matching coefficients shouldalso be, and indeed ours match those in the scalar calcu-lation of [23]. Of course the finite terms should not be,and are not, the same.We next state the contribution from the low-scalematching. Unsurprisingly, as this effect accounts for elec-troweak symmetry breaking effects such as the gauge bo-son masses, it is in general dependent upon the identity ofthe final states. Again this is a matching calculation andinvolves evaluating diagrams that appear in SCET EW ,but not SCET γ , and we provide three examples below. W/Z
W/Z
W/Z
Here springs with a line through them are collinear gaugebosons with energy ∼ m χ in the DM center-of-massframe, and springs without the extra line are soft gaugebosons with energy ∼ m W,Z . A central difficulty in thecalculation is accounting for the effects of electroweaksymmetry breaking, see for example [50] for a recent dis-cussion. In order to simplify this we make use of the gen-eral formalism for electroweak SCET of [20, 43–46], whichwe have extended to include the case of non-relativisticexternal states. We postpone the details to App. C. Theapproach breaks the full low-scale matching into a “soft”and “collinear” component, which are the labels associ-ated with the non-diagonal and diagonal contributionsrespectively, rather than the effective theory modes thatgive rise to them. This distinction is discussed further inApp. C. In our case, ˆ D X ( µ ) in Eq. (8) can be specifiedthroughexp (cid:104) ˆ D X ( µ ) (cid:105) = (cid:104) ˆ D s ( µ ) (cid:105) (cid:104) D χc ( µ ) I (cid:105) exp (cid:34)(cid:88) i ∈ X D ic ( µ ) I (cid:35) , (14)where again X can be γγ , γZ or ZZ , ˆ D s ( µ ) is the non-diagonal soft contribution and a matrix as it mixes theoperators, whilst D χc ( µ ) and D ic ( µ ) are the initial and fi-nal state diagonal contributions respectively. Note bothˆ D s ( µ ) and the identity matrix I are 2 × O ( α ), whereas the final state diagonalcontribution has its logarithmically enhanced contribu-tion resummed to all orders. Using this definition we findthat the components of the soft matrix are (see App. C):[ ˆ D s ] = 1 + α ( µ )2 π (cid:20) ln m W µ (1 − iπ ) + c W ln m Z µ (cid:21) , [ ˆ D s ] = α ( µ )2 π ln m W µ (1 − iπ ) , (15)[ ˆ D s ] = 1 + α ( µ )2 π ln m W µ (2 − iπ ) , [ ˆ D s ] = 1 . Here and throughout we use the shorthand c W = cos ¯ θ W and s W = sin ¯ θ W . Further, the diagonal contributionscan be written as: D χc ( µ ) =1 − α ( µ )2 π (cid:20) ln m W µ + c W ln m Z µ (cid:21) ,D ic ( µ ) = α ( µ )2 π (cid:34) ln m W µ ln 4 m χ µ −
12 ln m W µ − ln m W µ + c i ln m Z µ + c i (cid:21) , (16) This calculation can also be performed using the rapidity renor-malization group [48], but in order to make best use of earlierSCET calculations in SCET EW we will not use that formalismhere. where i = Z or γ and we have: c Z = 5 − s W − s W c W ,c γ = 1 − s W , (17)and c Z = − . − . i ,c γ = − . . (18)Analytic expressions for these last results are provided inApp. C and App. E, and we give numerical values hereas the expressions are lengthy. Note that we have distin-guished between factors of m W and m Z in all logarithms.The µ dependence of the low-scale matching is demon-strated to cancel with that in our high-scale matchingresult when the running is turned off, the details beingshown in App. D. We emphasise that this cross checkinvolves not only the µ dependence of the objects inEq. (14), but also the µ dependence of the high-scale co-efficients stated in Eq. (13) and further the SM SU(2) L and U(1) Y β -functions. The full µ cancellation is non-trivial – it requires the interplay between each of theseobjects. This ultimately provides us with confidence inthe results as stated. As a further check, our low-scalematching result does not depend on the spin of the DM.As such we should be again able to compare our result tothe scalar case calculated in [23]. In that work they onlyconsidered the γγ final state, and also neglected the im-pact of SM fermions. Restricting our calculation to thesame assumptions, we confirm that the µ dependence inour result matches theirs.Taking our results in combination, we can extend theNLL calculation to NLL (cid:48) . Of course we cannot show fullNNLL results in the absence of the higher order anoma-lous dimension calculation, nevertheless the results westate here determine the cross section with perturbativeuncertainties on the Sudakov effects reduced to the per-cent level. At O ( α ), our calculation accounts for allterms of the form α ln ( µ m χ /µ Z ), α ln ( µ m χ /µ Z ), and α ln ( µ m χ /µ Z ). The first perturbative term we are miss-ing at this order is α ln( µ m χ /µ Z ). Taking µ Z = m Z and m χ anywhere from m Z to 20 TeV, we estimate the ab-sence of these terms induces an uncertainty that is lessthan 1%, demonstrating the claimed accuracy.To combine the various results stated above into thecross section we take the factorized results in Eq. (10),and note that as the higher order Wilson coefficients havenothing to do with the Sommerfeld enhancement, theircontribution is included in the Σ terms as given explicitlyin Eq. (12). We know that at tree level s = s ±± = 1 Again note that all counting here is relative to the lowest or-der contribution, which occurs at C tree ∼ O ( α ). As such theabsolute order of the terms in this sentence is O ( α ). m χ [ TeV ] | Σ ������� ������������ �� χ + χ - → ��� � γ � γγ �������� ′ m χ [ TeV ] | Σ - Σ ������� ������������ �� χ � χ � → ��� � γ � γγ �������� ′ FIG. 1:
Here we show our NLL (cid:48) result for the electroweak corrections to the charged (left) and neutral (right) DM annihilations obtainedby adding the one-loop high and low-scale corrections to the NLL result. The result is in good agreement with the known NLL calculation,but with smaller uncertainty since the scale uncertainties have been reduced. The bands here are derived by varying the high scale between m χ and 4 m χ . m χ [ TeV ] | Σ ������� ������������ �� χ + χ - → ��� � γ � γγ �������� ′ m χ [ TeV ] | Σ - Σ ������� ������������ �� χ � χ � → ��� � γ � γγ �������� ′ FIG. 2:
As for Fig. 1, but showing a variation in the low-scale matching between m Z / m Z , rather than a variation of the high-scalematching. As can be seen the NLL (cid:48) contribution has reduced the low scale dependence in both charged and neutral DM annihilationcases, and is again consistent with the NLL result. and s ± = s ± = 0, implying that when the Sommerfeldenhancement can be ignored we can associate | Σ | withthe Sudakov contribution to χ + χ − annihilation and | Σ − Σ | with χ χ .For this reason, in Fig. 1 and Fig. 2 we show the con-tributions to | Σ | and | Σ − Σ | for LL, NLL and NLL (cid:48) .In both cases we see the addition of the one-loop correc-tions is completely consistent with the NLL results, sug-gesting that this approach has the Sudakov logarithmsunder control. In these plots we take a central value of µ m χ = 2 m χ and µ Z = m Z . In Fig. 1 the bands are de-rived from varying the high-scale matching between m χ and 4 m χ . Recall that if we were able to calculate thesequantities to all orders, they would be independent of µ ,and so varying these scales estimates the impact of miss-ing higher order terms. For the | Σ | NLL result, taking µ m χ = 2 m χ is a minimum in the range varied over, sowe symmetrise the uncertainties in order to more conser-vatively estimate the range of uncertainty. Similarly inFig. 2 we show the equivalent plot, but here the bandsare derived by varying the low scale µ Z from m Z / m Z . Improving on the high and low-scale matching, aswe have done here, should lead to a reduction in the scaleuncertainty. In all four cases shown this is clearly visibleand furthermore all results are still consistent with theNLL result within the uncertainty bands.We can also take this result and determine the impacton the full DM annihilation cross section into line pho-tons from γγ and γZ in this model, as we show in Fig. 3.We take the uncertainty on our final result to include thehigh and low-scale variations added in quadrature. ForH.E.S.S. limits we use [2], whilst for the CTA projection - - - - m χ [ TeV ] σ γγ + Z γ / v [ c m / s ] �� + ����� + ����� ′ + ��� - ���� + �� ( �� ) �� ′ γ + � ( �� ) ���� ����� ( ��� ) ��� ���������� ( ��� ) FIG. 3:
The impact of the NLL (cid:48) result on the full cross section,which includes the Sommerfeld Enhancement (SE), is shown to beconsistent with the lower orders result, suggesting the electroweakcorrections are under control. Also shown is the rate for the semi-inclusive process γ + X calculated to LL (cid:48) in [26]. In addition onthis plot we show current bounds from H.E.S.S. and projected onesfrom CTA, determined assuming 5 hours of observation time. Seetext for details. we assume 5 hours of observation time and use [37, 51].For both we assume an NFW profile with a local DM den-sity of 0.4 GeV/cm . We see again that our partial NLL (cid:48) results are consistent with the NLL conclusions. In thisfigure we also include the LL (cid:48) result for the semi-inclusiveprocess γ + X taken from Fig. 7 of [26], denoted by (BV).The semi-inclusive result is above our line photon result,except at low DM masses. Note that this work does notshow scale uncertainties, so the precise difference is hardto quantify numerically. IV. COMPARISON TO EARLIER WORK
In addition to using our results from the previous sec-tion in conjunction with the running due to the anoma-lous dimension, we can also consider the case where wetake our one-loop result in isolation. In this sense weshould be able to reproduce the initial problem of largelogarithms seen by HI. We show this in Fig. 4, comparedto the LL and NLL result. For Σ our one-loop resultis consistent with that from NLL, indicating the impor-tance of the α ln ( µ m χ /µ Z ) and α ln( µ m χ /µ Z ) correc-tions to C tree . For Σ − Σ , which starts at NLL, ourone-loop result is only consistent with the NLL expres-sion in the small m χ region. A digitized version of our cross section is available with the arXivsubmission or upon request.
For the | Σ | case we also show on that plot the equiva-lent curve for HI as extracted from Fig. 11 of their paper.From here it is clear that the qualitative shape of our re-sults agrees with theirs but that there is disagreement inthe normalization. This disagreement is already clear inFig. 3 and is more evident in Fig. 4. In Fig. 5 we analyzethis difference in more detail. In the left panel we showthe difference between their result and ours, showing ourcalculation with and without the low-scale matching in-cluded. Given the low-scale matching accounts for theelectroweak masses, which were included in HI, we wouldexpect including it to improve the agreement. This isseen, but it does not substantially relieve the tension.To further explore the difference, in the right panel ofFig. 5 we take our results and shift them down by a con-stant: 0.175 for the high only result and 0.137 for thehigh and low combination. Such a constant offset couldoriginate from a difference in m χ independent terms be-tween our result and HI. Unfortunately, however, a dif-ference in such terms could originate from almost anyof the graphs contributing to the result. Comparing ouranalytic expressions to the numerical results of HI wehave been unable to pinpoint the exact location of thedisagreement, although it is clear that we agree on the m χ dependence of the higher order corrections.Despite the discrepancy between our result and that ofHI, we emphasise that we have confidence in our resultas stated. This confidence is derived from the non-trivialcross checks we have performed on our result. In detail,these are • The cancellation in the O ( α ) corrections of the µ m χ dependence in our high-scale matching coeffi-cients, stated in Eq. (13), with the high-scale de-pendence entering from the anomalous dimension,as stated in Eqs. (6) and (7). This cancellation isdemonstrated in App. B; • In the absence of running, the cancellation in the O ( α ) corrections of the µ dependence between ourhigh and low-scale results, where the latter is statedin Eqs. (14), (15), (16), (17), and (18). This cancel-lation also depends on the SM SU(2) L and U(1) Y β -functions and is shown in App. D; • We have confirmed that the µ dependence in ourlow-scale result matches that in [23], when we elim-inate parts of our calculation in order to make thesame assumptions used in that work; • The form of the dominant µ independent terms inthe low-scale matching are in agreement with theresults of [20, 43–46], as discussed in App. C; and • We have confirmed that the framework used tocalculate the low-scale matching for our non-relativistic initial state kinematics, reproduces theresults of [20, 43–46] when we instead considermassless initial states as used in those references. m χ [ TeV ] | Σ ������� ������������ �� χ + χ - → ��� � γ � γγ ������ - ���� ( ���� ��� ) � - ���� ( �� ) m χ [ TeV ] | Σ - Σ ������� ������������ �� χ � χ � → ��� � γ � γγ ������ - ���� ( ���� ��� ) FIG. 4:
Similar to Fig. 1, but instead of displaying NLL (cid:48) curves we show our high and low-scale one-loop results including no runningfrom the anomalous dimension. For the case of χ + χ − annihilation we further show the equivalent result of HI, taken from Fig. 11 of theirwork (which only extends up to 3 TeV). There is evidently some discrepancy between the results. Note that at low masses where theSudakov logarithms are not too large, our result is consistent with the NLL result as would be expected. See text for details. m χ [ TeV ] | Σ ��� - ���� ���������� � - ���� ���� ( ���� ��� ) � - ���� ���� + ��� ( ���� ��� ) � - ���� ( �� ) m χ [ TeV ] | Σ ��� - ���� ���������� � - ���� ���� ( ���� ��� ) � �������� - ���� ���� + ��� ( ���� ��� ) � �������� - ���� ( �� ) FIG. 5:
We show the result of HI for | Σ | compared to two variations of our result. Firstly in the left panel we show our result with thehigh only or high and low-scale calculations compared to the result of HI, taken from Fig. 11 of their paper, demonstrating that there isa disagreement. In the right panel we take our results and shift them each by a m χ independent constant. The shifted results show thatabove ∼ m χ dependence of our result is in good agreement with HI. V. CONCLUSION
In this work we provide analytic expressions for the fullone-loop corrections to heavy wino dark matter annihila-tion, allowing the systematic resummation of electroweakSudakov logarithms to NLL (cid:48) for the line cross section.We have compared our result to earlier numerical cal-culations of such effects, finding results similar in be-haviour but quantitatively different. Our result is statedin a manner that can be straightforwardly extended tohigher order, with our result already reducing the pertur-bative uncertainty from Sudakov effects on this processto O (1%). Acknowledgements
The authors thank Aneesh Manohar for helpful discus-sions and for providing explicit cross checks for aspects ofour low-scale matching. We also thank Matthew Baum-gart, Tim Cohen, Ian Moult and Hiren Patel for help-ful discussions and comments. Feynman diagrams weredrawn using [52] and NLR thanks Joshua Ellis for assis-tance with its use. This work is supported by the U.S.Department of Energy under grant Contract NumbersDE-SC00012567, DE-SC0013999, DE-SC0011090 and bythe Simons Foundation Investigator grant 327942. NLRis supported in part by the American Australian Associ-ation’s ConocoPhillips Fellowship.
Appendix A: One-loop Calculation of χ a χ b → W c W d in the Full Theory In this appendix we outline the details of the high-scalematching calculation, which gives rise to the Wilson co-efficients stated in Eq. (13). These coefficients are deter-mined solely by the ultraviolet (UV) physics, allowing usto simplify the calculation by working in the unbrokentheory with m W = m Z = δm = 0. Combining this withthe heavy Majorana fermion DM being non-relativistic,there are only two possible Dirac structures that can ap-pear in the result: M A = (cid:15) ∗ µ ( p ) (cid:15) ∗ ν ( p ) (cid:15) σµνα p α i ¯ v ( p ) γ σ γ u ( p ) , M B = (cid:15) ∗ µ ( p ) (cid:15) ∗ ν ( p ) g µν ¯ v ( p ) /p u ( p ) , (A1)where p and p are the momenta of the incomingfermions, whilst p and p correspond to the outgoingbosons. The symmetry properties of these structures un-der the interchange of initial and final state particles al-low us to write our full amplitude as: M abcd = 4 πα m χ { [ B δ ab δ cd + B ( δ ac δ bd + δ ad δ bc )] M A + B ( δ ac δ bd − δ ad δ bc ) M B } . (A2)The above equation serves to define the Wilson coeffi-cients B r in a convenient form. These coefficients arerelated to the EFT coefficients of the operators definedin Eq. (2) and (4) via: C = ( − πα /m χ ) B , C = ( − πα /m χ ) B . (A3)For NLL accuracy we only need the tree-level value ofthese coefficients, which receive a contribution from s , t and u -channel type graphs and were calculated in [24].For completeness we state their values here: B (0)1 = 1 , B (0)2 = − , B (0)3 = 0 . (A4)Combining these with Eq. (A3), we see that the firstterms in Eq. (13) are indeed the tree-level contributionsas claimed.The operator associated with B was not discussed inthe earlier work of [24] as it cannot contribute to thehigh-scale matching calculation at any order, as we willnow argue. Firstly note that the B operator is skewunder the interchange a ↔ b . Due to the mass splittingbetween the neutral and charged states, present day anni-hilation is initiated purely by χ χ = χ χ , a symmetricstate that cannot overlap with B . One may worry thatexchange of one or more weak bosons between the initialstates – the hallmark of the Sommerfeld enhancement –may nullify this argument. But it can be checked thatif the initial states to such an exchange have identicalgauge indices, then so will the final states. As such B isnot relevant for calculating high-scale matching. Diagrams where a soft gauge boson is exchanged between an
In spite of this there are several reasons to calculate B here. From a practical point of view B gives us anadditional handle on the consistency of our result, whichwe check in App. B. Given that many graphs that gen-erate B and B also contribute to B , the consistencyof B provides greater confidence in the results for theoperators we are interested in. Further, from a physicspoint of view, although B is not relevant for high-scalematching when considering present day indirect detec-tion experiments, it could be relevant for calculating theannihilation rate in the early universe, where all statesin the DM triplet were present, to the extent that thenon-relativistic approximation is still relevant. For thisreasons we state it in case it is of interest for future work,such as expanding on calculations of the relic density atone loop (see for example [53–55]). Determining Matching Coefficients
Let us briefly review how matching coefficients are cal-culated at one loop. To begin with we can write the gen-eral structure of the UV and infrared (IR) divergences ofthe bare one-loop result for annihilation diagrams in thefull theory as: M fullbare = K(cid:15) + L(cid:15) IR + M(cid:15) UV + N (cid:18) (cid:15) UV − (cid:15) IR (cid:19) + C , (A5)where N is the coefficient associated with the variousscaleless integrals, and C is the finite contribution. Nowthe full theory is a renormalizable gauge theory, sowe know the additional counter-term and wavefunctionrenormalization contributions must be of the form: δ full = − M + N(cid:15) UV + D + E(cid:15) + F(cid:15) IR , (A6)where the values of D , E and F are scheme depen-dent. Nonetheless when calculating matching coefficientsit is easiest to work in the on-shell scheme for the wave-function renormalization factors, so below to denote thiswe add an “os” subscript to D , E and F . The reason thisscheme is the most straightforward, is that in any otherscheme when we map our Feynman amplitude calcula-tion for M full onto the S -matrix elements we want viathe LSZ reduction, there will be non-trivial residues cor-responding to the external particles. When using the on-shell scheme for the wave-function renormalization fac-tors, however, these residues are just unity, which sim-plifies the calculation as we can then ignore them. We initial and final state particle would in principle allow B tocontribute. Such a contributions would however be to the low-scale matching, which we discuss in App. C. As discussed there, B contributions to present day DM annihilation are power sup-pressed, and therefore do not contribute at any order in the lead-ing power effective theory. δ full with thebare results we obtain a UV finite answer: M fullren . = K + E os (cid:15) + L − N + F os (cid:15) IR + C + D os . (A7)In our calculation we will use dimensional regularizationto regulate both UV and IR divergences, which effectivelysets (cid:15) UV = (cid:15) IR , causing all scaleless integrals to vanish.Naively this seems to change the above argument, but aslong as we still use the correct counter-term in Eq. (A6)we find: M fullren . = K(cid:15) + L(cid:15) + M(cid:15) + C − M + N(cid:15) + D os + E os (cid:15) + F os (cid:15) = K + E os (cid:15) + L − N + F os (cid:15) + C + D os . (A8)Comparing this with Eq. (A7), we see that if we interpretall of the divergences in the final result as IR, then thismethod is equivalent to carefully distinguishing (cid:15) UV and (cid:15) IR throughout.In the EFT, with the above choice of zero massesand working on-shell with dimensional regularization, allgraphs are scaleless. At one loop they have the generalform: M EFTbare = O (cid:18) (cid:15) − (cid:15) (cid:19) + P (cid:18) (cid:15) UV − (cid:15) IR (cid:19) . (A9)Importantly if we have the correct EFT description ofthe full theory, then the two theories must have the sameIR divergences. Comparing Eq. (A9) to Eq. (A7), we seethis requires O = − K − E os and P = N − L − F os . TheEFT is again a renormalizable theory, so we can cancelthe UV divergences using δ EFT = ( K + E os ) (cid:15) − + ( L + F os − N ) (cid:15) − . Note as all EFT graphs are scaleless thereare no finite contributions that could be absorbed intothe counter-term, so in any scheme there is no finite cor-rection to δ EFT . Using this counter-term, we conclude: M EFTren . = K + E os (cid:15) + L − N + F os (cid:15) IR . (A10)Again note that for a similar argument to that in the fulltheory, if we had set (cid:15) UV = (cid:15) IR at the outset, then aslong as we still used the correct counter-term we wouldarrive at the same result.The matching coefficient is then obtained from sub-tracting the renormalized EFT from the renormalized full One may worry there could also be scaleless integrals of the form (cid:16) (cid:15) − − (cid:15) − (cid:17) , but the use of the zero-bin subtraction [56] en-sures such contributions cannot appear. theory result, so taking the appropriate results above weconclude: M fullren . − M EFTren . = C + D os . (A11)Comparing this with Eq. (A7), we see that providedwe have the correct EFT, then the matching coeffi-cient is just the finite contribution to the renormal-ized full-theory amplitude in the on-shell scheme. Eventhough this result makes explicit reference to a scheme in D on − shell , it is in fact scheme independent. The reasonfor this is that if we worked in a different scheme, al-though D would change, we would also have to accountfor the now non-trivial external particle residues that en-ter via LSZ. Their contribution is what ensures Eq. (A11)is scheme independent. Results of the Calculation
As outlined above, in order to obtain the matchingcoefficients we need the finite contribution to the renor-malized full theory amplitude. Now to compute this inthe particular theory we consider in this paper, we needto calculate the 25 diagrams that contribute to the one-loop correction to χ a χ b → W c W d . The diagrams areidentical to those considered in [57], where they defineda numbering scheme for the diagrams, grouping them bytopology and labelling them as T i for various i . We fol-low that numbering scheme here, but cannot use theirresults as they considered massless initial state fermionswhilst ours are massive and non-relativistic. In generalwe calculate the diagrams using dimensional regulariza-tion with d = 4 − (cid:15) to regulate the UV and IR, andwork in ‘t Hooft-Feynman gauge. Loop integrals are de-termined using Passarino-Veltman reduction [58], and wefurther make use of the results in [59–62] as well as Feyn-Calc [63, 64] and Package-X [65].In the EFT description of the full theory outlined inSec. II, the factorization of the matrix elements ensured aseparation between the Sommerfeld and Sudakov contri-butions. Yet for the full theory no clear separation existsand there will be graphs that contribute to both effects– in particular the graph T c considered below. The pur-pose of the Wilson coefficients we are calculating hereis to provide corrections to the Sudakov contribution –we do not want to spoil the EFT distinction by includ-ing Sommerfeld effects in these coefficients. In order tocleanly separate the contributions we take the relativevelocity of our non-relativistic initial states to be zero.This ensures that any contribution of the form 1 /v , char-acteristic of Sommerfeld enhancement, become power di-vergences and therefore vanish in dimensional regulariza-tion. This is different to the treatment in HI, where theycalculated the diagram without sending v → v andsubtracting the NRDM-SCET EW Sommerfeld graphs.1In our calculation the DM is a Majorana fermion. Itturns out that for almost all the graphs below the resultis identical regardless of whether we think of the fermionas Majorana or Dirac – a result that is also true at tree-level. The additional symmetry factors in the Majoranacase are exactly cancelled by the factors of 1 / T d and T d below, as well as closed fermionloop contributions to the counter-terms.Using the approach outlined above we now state thecontribution to B r as defined in Eq. (A2) graph by graph.Throughout we define L ≡ ln µ/ m χ . T a The result for this graph and its cross term is: B [1 a ]1 = α π (cid:20) − (cid:15) − (cid:15) (4 L + 2 iπ + 2) − L − L − iπL − π (cid:21) ,B [1 a ]2 = 12 B [1 a ]1 ,B [1 a ]3 = α π (cid:20) (cid:15) + 14 (cid:15) (2 L − iπ −
2) + 12 L − L − iπL + 17 π −
16 (2 + 7 iπ − (cid:21) . (A12)In calculating this graph in the non-relativistic limit viaPassarino-Veltman reduction there are additional spuri-ous divergences that must be regulated. The origin ofthese divergences is that Passarino-Veltman assumes themomenta appearing in the integrals to be linearly inde-pendent. But in the center of momentum frame if wetake v = 0, then p and p are identical and this assump-tion breaks down, leading to the divergences of the form( s − m χ ) − , where s = ( p + p ) . A simple way toregulate them is to give the initial states a small relativevelocity. This does not lead to a violation of our sepa-ration of Sommerfeld and Sudakov effects as this graphdoes not contribute to the Sommerfeld enhancement. Assuch this procedure introduces no 1 /v contributions tothe final result and the regulator can be safely removedat the end. This is the only diagram where this issueappears – if it occurred in a graph that did contribute tothe Sommerfeld effect we would need to use a differentregulator, or explicitly subtract the corresponding EFTgraph at finite v . T b This graph has a single crossed term and combining thetwo yields: B [1 b ]1 = B [1 b ]3 = 0 ,B [1 b ]2 = α π (cid:20) (cid:15) + 4 L + 2 (cid:15) + 4 L ( L + 1) − π − (cid:21) . (A13) T c The combination of this graph and its crossed term is: B [1 c ]1 = α π (cid:20) (cid:15) − L + 4 ln 2 (cid:21) ,B [1 c ]2 = 12 B [1 c ]1 ,B [1 c ]3 = α π (cid:20) (cid:15) − L + π − (cid:21) . (A14)Formally this graph also gives a contribution to the Som-merfeld enhancement in the full theory. Nevertheless aswe take v = 0 at the outset, the contribution here ispurely to the Sudakov terms. T d The contribution from this diagram vanishes in the non-relativistic limit, i.e. B [1 d ]1 = B [1 d ]2 = B [1 d ]3 = 0 . (A15) T a B [2 a ]1 = B [2 a ]2 = 0 ,B [2 a ]3 = α π (cid:20) (cid:15) + 2 L + iπ
24 + 1172 (cid:21) . (A16) T b For a scalar Higgs in the loop, the graph and its crossterm contribute: B [2 b ]1 = B [2 b ]2 = 0 ,B [2 b ]3 = α π (cid:20) (cid:15) + 2 L + iπ
12 + 1136 (cid:21) . (A17) T c There is no crossed graph associated with the graphabove as the gauge bosons running in the loop are realfields. As such taking just this graph gives B [2 c ]1 = B [2 c ]2 = 0 ,B [2 c ]3 = α π (cid:20) (cid:15) + 1 (cid:15) (cid:18)
34 (2 L + iπ ) + 178 (cid:19) + 38 (2 L + iπ ) + 178 (2 L + iπ ) + 9524 − π (cid:21) . (A18) T d There are two types of fermions that can run in the loop:the Majorana triplet fermion that make up our DM orleft-handed SM doublets. As with the gauge bosons theseSM fermions are taken to be massless and for generality we say there are n D of them. For the SM doublets thereis a crossed graph, whilst for the Majorana DM field thereis not, so that: B [2 d ]1 = B [2 d ]2 = 0 ,B [2 d ]3 = α π (cid:20) − (cid:18) (cid:15) + 43 L + 43 ln 2 −
59 + π (cid:19) − n D (cid:18) (cid:15) + 16 (2 L + iπ ) + 736 (cid:19)(cid:21) . (A19)If the DM had been a Dirac field instead, there wouldhave been a crossed graph and the result would be mod-ified such that the first line of B [2 d ]3 gets multiplied by2. The factor of 7 /
36 we find in the last line of B [2 d ]3 isconsistent with the expression found for this graph, butwith different kinematics, in [57], but disagrees with [70]. T e − h The four graphs shown above do not contribute to ourone-loop result; the graphs on the top row vanish atleading order for non-relativistic initial states, whilst theloops on the second line are both scaleless and so areidentically zero in dimensional regularization. As suchwe have: B [2 e − f ]1 = B [2 e − f ]2 = B [2 e − f ]3 = 0 . (A20) T a and T a The two graphs shown above have identical amplitudes.For each graph independently, the sum of it and its For the SM well above the electroweak scale n D = 12. In detail,for each generation there are four doublets: the lepton doubletand due to color, three quark doublets. As such for three gener-ations we have twelve left-handed SM doublets. B [3 a/ a ]1 = α π (cid:20) − (cid:15) + 2 − L(cid:15) − L +4 L − π (cid:21) ,B [3 a/ a ]2 = − B [3 a/ a ]1 ,B [3 a/ a ]3 = 12 B [3 a/ a ]1 . (A21) T b and T b As for T a and T a , these two graphs also have equalamplitudes. Again we provide the combination of eachwith its crossed graph: B [3 b/ b ]1 = α π (cid:20) (cid:15) + 2 L − π (cid:21) ,B [3 b/ b ]2 = − B [3 b/ b ]1 ,B [3 b/ b ]3 = 12 B [3 b/ b ]1 . (A22) T a Whether the above graph has a crossed graph associatedwith interchanging the initial states depends on the iden-tity of the initial state fermions. For Majorana fermionsthere is such a crossing, whilst for Dirac there is not.Despite this, in either case the combination of the graphand its crossing (where it exists) is the same in both casesand is simply: B [5 a ]1 = B [5 a ]2 = 0 ,B [5 a ]3 = α π (cid:20) − (cid:15) − L −
133 ln 2 −
83 + 23 iπ (cid:21) . (A23) T b As for T a the existence of a crossed graph depends onthe nature of the DM. Regardless again the result is thesame if we take it to be Dirac or Majorana, which is: B [5 b ]1 = B [5 b ]2 = 0 ,B [5 b ]3 = α π (cid:20) (cid:15) + 3 L + 3 ln 2 − (cid:21) . (A24) T a For a gauge boson in the loop we have: B [6 a ]1 = B [6 a ]2 = 0 ,B [6 a ]3 = α π (cid:20) − (cid:15) − L − − iπ (cid:21) . (A25)Note this graph and the remaining T type topologieshave no crossed graphs. T b In the case of a ghost loop we have: B [6 b ]1 = B [6 b ]2 = 0 ,B [6 b ]3 = α π (cid:20) − (cid:15) − L − − iπ (cid:21) . (A26) T c For a scalar Higgs we have an identical contribution to T b : B [6 c ]1 = B [6 c ]2 = 0 ,B [6 c ]3 = α π (cid:20) − (cid:15) − L − − iπ (cid:21) . (A27)4 T d As for T d the fermion in the loop could again be eitherDM or SM. Allowing there to be n D left-handed SM dou-blets we have B [6 d ]1 = B [6 d ]2 = 0 ,B [6 d ]3 = α π (cid:20)(cid:18) (cid:15) + 43 L + 43 ln 2 + 169 (cid:19) + n D (cid:18) (cid:15) + 13 L + 518 + 16 iπ (cid:19)(cid:21) . (A28)Here there is a symmetry factor of 1 / B [6 d ]3 would getmultiplied by 2 as this symmetry factor would not bepresent. T e and T f Both of these integrals are scaleless and vanish in dimen-sional regularization, so: B [6 e − f ]1 = B [6 e − f ]2 = B [6 e − f ]3 = 0 . (A29) T For the final graph we again have a crossed contribution,and combining the two gives: B [7]1 = α π (cid:20) − (cid:15) − L − (cid:21) ,B [7]2 = − B [7]1 ,B [7]3 = 12 B [7]1 . (A30) Counter-terms
To begin with, as B vanishes at tree level there are nocounter-term corrections to its value at one loop. Instead we only need to consider graphs that would contribute to B and B , of which there are three:The graph on the left corresponds to the internal wave-function and mass renormalization of the DM – renor-malization factors denoted as Z χ and Z m – whilst theremaining two graphs account for the renormalization ofthe DM and electroweak gauge boson interaction vertex g ¯ χ /W χ – here Z (which includes coupling and externalline wavefunction renormalization). Now if we calculatethe above three graphs, we find a contribution propor-tional to the tree-level amplitude M tree , as well as a termthat would contribute to B . The contribution to B iscancelled by the additional s -channel type counter-termgraphs not drawn, so the full counter-term contributionleaves only: (2 δ − δ χ − δ m ) M tree , (A31)where we have used Z i = 1 + δ i .Next, when determining the δ i we need to pick ascheme. As explained above, when calculating matchingcoefficients it is easiest to work in the on-shell scheme forwavefunction renormalization to ensure we do not have toworry about residues from the LSZ reduction. The mean-ing of the on-shell values of δ χ and δ m is clear, whereasfor δ we must write this out more explicitly. By defi-nition we know δ = δ g + δ W + δ χ , where δ g and δ W are the counter-terms for the coupling and gauge bosonwave-functions respectively. For the gauge-boson wave-function we use the on-shell scheme as usual. For thecoupling counter-term, however, we define it in the MSscheme. Since our full theory is defined with the DM asa propagating degree of freedom, this coupling is definedabove the m χ . In the EFT the DM is integrated out, sothe appropriate coupling for the matching is one definedbelow m χ . We put this issue aside for now and return toit in the next section.The above choices then define our scheme for δ ina manner that ensures all residues are still 1. With thisscheme, we can then calculate the relevant counter-termsand find: δ χ = − α π (cid:20) (cid:15) UV + 4 L + 4 ln 2 + 4 (cid:21) , (A32) δ m = − α π (cid:20) (cid:15) UV + 12 L + 12 ln 2 + 8 (cid:21) ,δ W = − α π (cid:20) n D − (cid:15) UV + 19 − n f (cid:15) IR + 163 L + 163 ln 2 (cid:21) ,δ g = − α π (cid:20) − n D (cid:15) UV (cid:21) ,δ = − α π (cid:20) (cid:15) UV + 19 − n D (cid:15) IR + 203 L + 203 ln 2 + 4 (cid:21) , n D is again the number of left-handed SM dou-blets. Recall that in determining the counter-terms wecannot neglect scaleless integrals as we did for the maincalculation, so their contribution has been included hereand we explicitly distinguish (cid:15) UV from (cid:15) IR . Substitutingthese results into Eq. (A31), we find the crossed contri-bution is: B [CT]1 = α π (cid:20) n D − (cid:15) IR + 83 L + 83 ln 2 + 4 (cid:21) ,B [CT]2 = α π (cid:20) − n D (cid:15) IR − L −
86 ln 2 − (cid:21) ,B [CT]3 = 0 . (A33)Interestingly the counter-term contribution is UV finite.This implies that the sum of all one-loop graphs beforeadding in counter-terms must be UV finite. Given thatwe used dimensional regularization to regulate both UVand IR divergences this cannot be immediately read offfrom our results, but going back to the integrals and keep-ing track of the UV divergences we confirmed that thesum is indeed UV finite.Note if our DM field had instead been a Dirac fermion,there would be several modifications to the above. Firstlythe L and ln 2 dependence in δ W and δ would be modi-fied, whilst the (cid:15) UV dependence in δ W and δ g would alsochange. In the combination stated in Eq. (A33) this onlychanges the L and ln 2 dependence, but in a way thatis exactly cancelled when we account for the scale of thecoupling in the next section. Scale of the Coupling
Throughout the above calculation we have treated theDM as a propagating degree of freedom and included itseffects in loop diagrams. This implies that the couplingused so far above in this appendix implicitly depends on n D +1 flavors – n D left-handed SM doublets and one Ma-jorana DM fermion – i.e. we have used α = α ( n D +1)2 ( µ ).In the EFT however, the DM is no longer a propagat-ing field and so the appropriate coupling is α ( n D )2 ( µ ). Atorder α , which we are working to at one loop, the dis-tinction will lead to a finite contribution because of thematching at the scale µ = m χ , which we calculate in thissection.Let us start by reviewing the standard treatment ofa running coupling in the MS scheme. This running iscaptured by the β -function, which is defined by β ( α ) = µdα /dµ , where here α is the renormalized coupling; thebare coupling is independent of µ . The β -function canbe written as: β ( α ) = − (cid:15)α − b π α + . . . , (A34)where we have expanded it to the order needed for thisthreshold matching analysis. At this order the LL solu- tion for the running of the coupling is: α ( µ ) = α ( µ )1 + α ( µ ) b π ln µµ . (A35)In order to determine the threshold matching correctionat the one-loop order we are working it suffices to simplydemand that the coupling is continuous at the scale m χ ,and this is captured by a difference in b . For our problemwe define b ( n D +1)0 to be the value above m χ and b ( n D )0 thevalue below. Then using Eq. (A35) to define α ( n D +1)2 ( µ )and α ( n D )2 ( µ ), it suffices to demand they match at a scale m χ , which gives: α ( n D +1)2 ( µ ) = α ( n D )2 ( µ ) (cid:20) α ( n D )2 ( µ )2 π (cid:16) b ( n D +1)0 − b ( n D )0 (cid:17) ln µm χ + . . . (cid:35) . (A36)So now we just need to determine b ( n D +1)0 − b ( n D )0 . Ingeneral for a theory containing just gauge bosons, Weylfermions (WF), Majorana fermions (MF) and chargedscalars (CS), we can write: b = 113 C A − (cid:88) i ∈ WF C ( R i ) − (cid:88) i ∈ MF C ( R i ) − (cid:88) i ∈ CS C ( R i ) . (A37)Our calculation has all four of these ingredients: elec-troweak gauge bosons, the left-handed SM fermions(which are Weyl because only one chirality couples tothe gauge bosons), the Majorana DM fermion and theHiggs. Then using C A = 2, C ( R ) = 1 / C ( R ) = 2 for theadjoint Wino, we conclude: b ( n D )0 = 43 − n D ,b ( n D +1)0 = 35 − n D . (A38)From this Eq. (A36) tells us that to the order we areworking: α ( n D +1)2 ( µ ) = α ( n D )2 ( µ ) (cid:34) − α ( n D )2 ( µ )4 π (cid:18) L + 83 ln 2 (cid:19)(cid:35) . (A39)Now as there is only a difference between the couplingsat next to leading order, this only corrects the tree levelresult stated in Eq. (A4). As such the impact of changingto the coupling defined below m χ , which is relevant forthe matching, is to add the following contribution: B [Matching]1 = α π (cid:20) − L −
83 ln 2 (cid:21) ,B [Matching]2 = − B [Matching]1 ,B [Matching]3 = 0 , (A40)6where after adding this contribution now here and in allearlier one-loop results we can simply take α = α ( n D )2 .As alluded to above, this result is modified for a Dirac DM fermion, but in a way exactly compensated by achange in the counter-term contribution. Combination
Combining the 25 graphs above with the counter-terms and the matching contributions, we arrive at the followingresult: B (1)1 = α π (cid:20) − (cid:15) − L + 12 iπ + 31 − n D (cid:15) − L − L − iπL − π (cid:21) ,B (1)2 = α π (cid:20) (cid:15) + 48 L − iπ + 55 − n D (cid:15) + 4 L + 6 L − iπL − π (cid:21) ,B (1)3 = α π (cid:20) n D −
72 ln 2 −
71 + 3 π (cid:21) , (A41)where recall L = ln µ/ m χ , n D is the number of SM left-handed doublets and now all (cid:15) = (cid:15) IR .As explained in detail at the outset of the calculation,the one-loop contribution to the matching coefficient isjust the finite part of this result. Combining this withthe tree-level term in Eq. (A4) and mapping back to C r using Eq. (A3) then gives us the Wilson coefficients inEq. (13), which we set out to justify.If instead we had a Dirac DM triplet rather than aMajorana, then the only impact on the above would befor B (1)3 , and we would instead have B (1)3 = α π (cid:20) n D −
72 ln 2 − (cid:21) . (A42) Appendix B: Consistency Check on the High-ScaleMatching
As a non-trivial check on our high-scale calculation,we can calculate the ln µ , or L in our case, pieces ofEq. (A41) independently by expanding the NLL results.To begin with, if we define C ≡ ( C C C ) T , then fromthe definition of the anomalous dimension we have: µ ddµ C ( µ ) = ˆ γ ( µ ) C ( µ ) . (B1)Next we expand the coefficients as a series in α : C ( µ ) = C (0) ( µ ) + C (1) ( µ ) + ... , where C (0) ( µ ) is the tree-levelcontribution and C (1) ( µ ) the one-loop result. Now wewant a cross check on the one-loop contribution, so weevaluate Eq. (B1) at O ( α ), giving µ dα dµ ∂C (0) ∂α + µ ∂C (1) ( µ ) ∂µ = ˆ γ − loop ( µ ) C (0) ( µ ) , (B2)and rearranging we arrive at: µ ∂C (1) ( µ ) ∂µ = ˆ γ − loop ( µ ) C (0) ( µ ) − µ dα dµ ∂C (0) ∂α . (B3) This equation shows that we can derive the µ and hence L dependence of the one-loop Wilson coefficient from theone-loop anomalous dimension and tree-level Wilson co-efficient, both of which are known from the NLL result.To be more explicit, we can write the bare Wilson coef-ficient as C bare = µ (cid:15) (cid:18) a(cid:15) + b(cid:15) + µ − independent (cid:19) = a(cid:15) + b + 2 aL(cid:15) + 2 aL + 2 bL + µ − independent , (B4)where in the second equality we swapped frpm ln µ to L and absorbed the additional ln 2 factors into the µ -independent term. From here we can write the renor-malized Wilson coefficient as C ren . = 2 aL + 2 bL + µ − independent , (B5)which we can then substitute into the left-hand side ofEq. (B3) to derive a and b for each Wilson coefficient.Doing this and then mapping back to B r using Eq. (A3),we find B (1)1 = α π (cid:20) − L(cid:15) − L − L − iπL + µ − ind . (cid:21) ,B (1)2 = α π (cid:20) L(cid:15) + 4 L + 6 L − iπL + µ − ind . (cid:21) ,B (1)3 = α π [0 + µ − ind . ] , (B6)in exact agreement with Eq. (A41). In particular, as B (0)3 = 0, we needed B (1)3 to be independent of L , as wefound.7 Appendix C: Low-Scale Matching Calculation
The focus of this appendix is to derive the low-scalematching conditions stated in Eqs. (14), (15), (16), (17),and (18). At this scale, the matching is from an effectivetheory where the W , Z , top and Higgs are dynamicaldegrees of freedom – NRDM-SCET EW – onto a theorywhere these electroweak modes have been integrated out– NRDM-SCET γ .In order to perform the calculation we will make use ofthe formalism of electroweak SCET developed in [20, 43–46]. As we are working in SCET, there are both collinearand soft gauge boson diagrams that will appear in theone-loop matching. In [20] it was proven that at one-loopthe total low-scale matching contribution from these softand collinear SCET modes can always be decomposedinto a contribution that is diagonal, in that it leads tono operator mixing, and another that is non-diagonal, asit does induce mixing. In their works, they then refer tothe diagonal parts as collinear and non-diagonal ones as soft , however we shall always use the term “diagonal” torefer to the contributions that have contributions fromboth soft and collinear diagrams, although we do use asubscript “c” for the diagonal piece. At one loop thematching amounts to evaluating the diagrams that ap-pear in NRDM-SCET EW but not NRDM-SCET γ . Thesediagrams can be broken into three classes:1. Wave-function diagrams correcting our initial non-relativistic states;2. Diagrams where a soft gauge boson is exchangedbetween two different external states; and3. Final state collinear diagrams, which are now cor-rections to collinear states.Each class will be discussed separately below. Before do-ing so, however, we first define our operators and outlinehow the low-scale matching proceeds at tree level.Unlike for the high-scale matching, here we only con-sider the two operators that match onto M A in Eq. (A1),as opposed to the third operator coming from M B . Thereason for this is the additional operator does not con-tribute to the low-scale matching calculation for presentday DM annihilation at any order in leading powerNRDM-SCET. To understand this note that the opera-tors coming from M A and M B have different spin struc-tures. In order to mix these structures we need to transferangular momentum between the states. The only low-scale graphs we can write down to do this are soft gaugeboson exchanges. The spin structure of the coupling of asoft exchange to an n -collinear gauge boson is /n and thecorresponding coupling to our non-relativistic DM fieldis /v . Neither coupling allows for a transfer of angular mo-mentum, demonstrating that these operators cannot mix.Unlike for the high-scale matching, we will not make useof the operator corresponding to M B for our low-scaleconsistency check, so we drop it from consideration atthe outset. Operator Definition and Tree-level Matching
Prior to electroweak symmetry breaking, the tworelevant operators in NRDM-SCET EW can be writtenschematically as: O = 12 δ ab δ cd χ a χ b W c W d , O = 14 ( δ ac δ bd + δ ad δ bc ) χ a χ b W c W d . (C1)Our notation here is schematic in the sense that we havesuppressed the Lorentz structure and soft Wilson lines.The form of these is written out explicitly in Eq. (2) andis left out for convenience as it appears in every operatorwritten down in this appendix. Further, in this equa-tion the factor of 1 / χ is a Majorana field this factor ensures the Feynmanrule associated with these operators has no additionalnumerical factor. Note also that the gauge bosons arelabelled as they are associated with a collinear direction.At tree-level the low-scale matching is effected simply bymapping the fields in these operators onto their brokenform. Explicitly we have: χ = 1 √ (cid:0) χ + + χ − (cid:1) ,χ = i √ (cid:0) χ + − χ − (cid:1) ,χ = χ ,W = 1 √ (cid:0) W + + W − (cid:1) ,W = i √ (cid:0) W + − W − (cid:1) ,W = s W A + c W Z . (C2)Substituting these into Eq. (C1) yields 22 operators inthe broken theory. Of these, 14 involve a W ± in the finalstate, so we will not consider them further. We definethe remaining 8 as:ˆ O = 12 χ χ A A , ˆ O = 12 χ χ Z A , ˆ O = 12 χ χ A Z , ˆ O = 12 χ χ Z Z , ˆ O = χ + χ − A A , ˆ O = χ + χ − Z A , ˆ O = χ + χ − A Z , ˆ O = χ + χ − Z Z , (C3)where again we have used the schematic notation ofEq. (C1), as we will for all operators in this appendix.At tree level, the operators in Eq. (C1) and (C3) are re-lated simply by the change of variables in Eq. (C2). Thismapping is performed by a 22 × D (0) s, − = s W s W s W c W s W c W s W c W s W c W c W c W s W s W c W s W c W c W . (C4)In terms of the calculation presented in the main text,what we actually want is the mapping onto the Sudakovfactors Σ, defined in Eq. (10), not the broken opera-tors in Eq. (C3). As given there, the s W and c W fac-tors are absorbed into P X , and so will not contributeto the Σ factors. Then ˆ O − represent the contribu-tions to neutral annihilation χ χ → X , represented byΣ − Σ , and ˆ O − the contributions to charged annihi-lation χ + χ − → X , represented by Σ . Accordingly wehave: ˆ D (0) s = (cid:20) (cid:21) . (C5)This provides the tree-level result we should use inEq. (14). Next we turn to calculating this one-loop low-scale matching in full, considering the three classes ofdiagrams that can contribute in turn. Initial State Wave-function Graphs
There are two graphs that fall under the categoryof initial state wave-function corrections, and these areshown below.Note here we follow the standard SCET conventions ofdrawing collinear fields as gluons with a solid line throughthem, whereas soft fields are represented simply by gluonlines. In these graphs, the soft gauge field can be either a W or Z boson. In either case the integral to be calculatedis: − g (cid:90) ¯ d d k µ (cid:15) [ k − m ] v · ( k + p ) , (C6)where g is the coupling – g for a W boson, c W g fora Z boson, p is the external momentum, k is the loopmomentum, m the gauge boson mass, and v is the veloc-ity associated with the non-relativistic χ field. Given ourinitial state is heavy, this is unsurprisingly exactly theheavy quark effective theory wave-function renormaliza-tion graph. The analytic solution can be found in e.g. [71, 72], and using this we find:= − iv · p α π (cid:20) (cid:15) + ln µ m (cid:21) , (C7)where α = g / π . Now in addition to the one-loopgraphs we drew above, at this order there will also bea counter-term of the form iv · p ( Z χ − Z χ = 1 + α ( µ )2 π (cid:20) (cid:15) − ln m W µ − c W ln m Z µ (cid:21) . (C8)Now each of our initial states will contribute Z / χ , im-plying that the contribution to ˆ D ( µ ) given in Eq. (14)is D χc ( µ ) = 1 − α ( µ )2 π (cid:20) ln m W µ + c W ln m Z µ (cid:21) , (C9)and the subscript c indicates this is a diagonal contri-bution in the sense that it leads to no operator mixing.This is exactly as in Eq. (16) and justifies this part of thelow-scale matching. Soft Gauge Boson Exchange Graphs
In this section we calculate the contribution from theexchange of a soft W or Z gauge boson between differentexternal final states. As these gauge bosons carry SU(2) L gauge indices, unsurprisingly these graphs will lead tooperator mixing. Consequently, in terms of the notationintroduced above these graphs will lead to non-diagonalcontributions. They will also induce diagonal terms, andwe will carefully separate the two below.Once separated, we will group the diagonal contribu-tion with those we get from the final state wave-functiongraphs we consider in the next subsection. The reasonfor this is that these diagonal contributions for photonand Z final states, as we have, were already evaluated in[46], and we will not fully recompute them here. In thatwork, however, the diagonal contribution was only statedin full. The breakdown into the soft boson exchange andfinal state wave-function graphs was not provided. Thisraises a potential issue because in that work all externalstates were taken to be collinear, not non-relativistic. Assuch, in this section we will explicitly calculate the softgauge boson exchange graphs for both kinematics anddemonstrate that the diagonal contribution is identicalin the two cases.Before calculating the graphs, we first introduce someuseful notation. At one loop the gauge bosons will havetwo couplings to the four external states. Each of thesecouplings will have an associated gauge index structure,and in order to deal with this it is convenient to introducegauge index or color operators T . This notation wasfirst introduced in [73, 74], and it allows the gauge index9structure to be organized generally rather than case bycase. Examples can be found in the original papers andalso in the SCET literature e.g. [20, 46, 75]. An examplerelevant for our purposes is the action of T on an SU(2) L adjoint, which is the representation of both our initialand final states: T χ a = ( T cA ) aa (cid:48) χ a (cid:48) = − i(cid:15) caa (cid:48) χ a (cid:48) , T W a = ( T cA ) aa (cid:48) W a (cid:48) = − i(cid:15) caa (cid:48) W a (cid:48) . (C10)In terms of this notation then, we can write the gauge in-dex structure of all relevant one-loop low-scale matchinggraphs as T i · T j , where i, j label any of the four externallegs. Because of this we label the result from these softexchange diagrams as S ij for the case of our kinematics– non-relativistic initial states and collinear final states– and we use S (cid:48) ij to denote the kinematics of [46] – allexternal states collinear. Following [20, 46], we take allexternal momenta to be incoming and further rapiditydivergences will be regulated with the ∆-regulator [76].Now let us turn to the graphs one by one. S ( (cid:48) )12 In this graph the soft gauge boson exchanged betweenthe initial state can be a W or Z boson. In either case,the value of this graph is: S = α π T · T (cid:20) (cid:15) − ln m µ (cid:21) , (C11) S (cid:48) = α π T · T (cid:20) (cid:15) − (cid:15) (cid:18) ln δ δ µ + iπ (cid:19) −
12 ln m µ + iπ ln m µ + ln m µ ln δ δ µ − π (cid:21) , where as above α = g / π and the identity g and m depend on whether this is for a W or Z . In S (cid:48) , δ and δ are the ∆-regulators. Unsurprisingly these onlyappear for the collinear kinematics for the initial state in S (cid:48) and not for the nonrelativistic kinematics in S . S ( (cid:48) )13 , S ( (cid:48) )14 , S ( (cid:48) )23 , and S ( (cid:48) )24 Again the exchanged soft boson can be a W or Z . Thesefour graphs are grouped together as they have a commonform, for example: S = α π T · T (cid:20) (cid:15) − (cid:15) ln δ µ −
14 ln m µ (C12)+ 12 ln δ µ ln m µ − π (cid:21) ,S (cid:48) = α π T · T (cid:20) (cid:15) − (cid:15) ln (cid:18) − δ δ µ w (cid:19) −
12 ln m µ + ln m µ ln (cid:18) − δ δ µ w (cid:19) − π (cid:21) . Then S ( (cid:48) )14 is given by the same expressions but with 3 →
4, whilst S ( (cid:48) )23 and S ( (cid:48) )24 are given by similar replacements.For the all collinear case we have defined the followingfunctions of the kinematics: w = w ≡ n · n = 12 n · n = ts ,w = w ≡ n · n = 12 n · n = us , (C13)where s , t , and u are the Mandelstam variables relevantfor all incoming momenta. The signs inside the logs inEq. (C12) can be understood by noting that as t < u <
0, whilst s >
0, we have w ij < S ( (cid:48) )34 Finally we have the graph above, which yields: S = α π T · T (cid:20) (cid:15) − (cid:15) (cid:18) ln δ δ µ + iπ (cid:19) −
12 ln m µ + iπ ln m µ + ln m µ ln δ δ µ − π (cid:21) ,S (cid:48) = S . (C14)This completes the list of graphs to evaluate. As writ-ten it appears that all graphs are non-diagonal from theirgauge index structure. However as we will now show, thecombinations of all graphs can be reduced to a diago-nal and non-diagonal piece. Firstly for the case of allcollinear external states we have: S (cid:48) + S (cid:48) + S (cid:48) + S (cid:48) + S (cid:48) + S (cid:48) ≡ (cid:88) (cid:104) ij (cid:105) S (cid:48) ij , (C15)0which serves to define (cid:104) ij (cid:105) . The part of this sum thatinvolving the rapidity regulators can be written as α π ln m µ (cid:88) (cid:104) ij (cid:105) T i · T j (cid:18) ln δ i µ + ln δ j µ (cid:19) . (C16)This can be simplified using the following identity: (cid:88) (cid:104) ij (cid:105) ( f i + f j ) T i · T j = − (cid:88) i f i T i · T i . (C17)If we identify f i = ln δ i /µ , then Eq. (C16) becomes:= − α π ln m µ (cid:88) (cid:104) ij (cid:105) T i · T i ln δ i µ , (C18)which is now diagonal in the gauge indices. For the re-maining terms that are independent of δ , we organisethem as follows: (cid:88) (cid:104) ij (cid:105) S (cid:48) ij = 12 [ S (cid:48) + S (cid:48) + S (cid:48) ]+ 12 [ S (cid:48) + S (cid:48) + S (cid:48) ]+ 12 [ S (cid:48) + S (cid:48) + S (cid:48) ]+ 12 [ S (cid:48) + S (cid:48) + S (cid:48) ] , (C19)where we used the fact S (cid:48) ij = S (cid:48) ji . Each of these groupscan now be simplified. For example, the first group canbe written as: S (cid:48) + S (cid:48) + S (cid:48) = α π ( T · T + T · T + T · T ) × (cid:20) −
12 ln m µ − π (cid:21) + α π T · T (cid:20) iπ ln m µ (cid:21) (C20) − α π T · T (cid:20) ln (cid:18) − ts (cid:19) ln m µ (cid:21) − α π T · T (cid:20) ln (cid:16) − us (cid:17) ln m µ (cid:21) , If we then use (cid:88) j,j (cid:54) = i T i · T j = − T i · T i , (C21) This and the gauge index identity stated below in Eq. (C21)follow simply from the fact (cid:80) i T i = 0 when it acts on gaugeindex singlet operators, see for example [20]. Eq. (C20) can be rewritten as:= α π T · T (cid:20)
12 ln m µ + π (cid:21) + α π T · T (cid:20) iπ ln m µ (cid:21) (C22) − α π T · T (cid:20) ln (cid:18) − ts (cid:19) ln m µ (cid:21) − α π T · T (cid:20) ln (cid:16) − us (cid:17) ln m µ (cid:21) . Repeating this for the remaining three terms in Eq. (C19)and reinserting the δ contributions, we can rewrite thecombination of all terms as: (cid:88) (cid:104) ij (cid:105) S (cid:48) ij ≡ (cid:88) (cid:104) ij (cid:105) ˆ S (cid:48) ij + (cid:88) i C i , (C23)where we have defined:ˆ S (cid:48) ij ≡ − α π ln m µ T i · T j U (cid:48) ij , (C24) C i ≡ α π T i · T i (cid:20)
14 ln m µ + π −
12 ln m µ ln δ i µ (cid:21) , and from the above we can see that: U (cid:48) = U (cid:48) = − iπ ,U (cid:48) = U (cid:48) = ln (cid:18) − ts (cid:19) ,U (cid:48) = U (cid:48) = ln (cid:16) − us (cid:17) . (C25)Thus as claimed, we have reduced (cid:80) (cid:104) ij (cid:105) S (cid:48) ij in Eq. (C23)into a diagonal and non-diagonal piece. Importantly wehave explicitly isolated the diagonal contribution C i , andas we will now show we get exactly the same diagonalcontribution for the kinematics of interest in this work.Before doing so, however, note that the irreduciblynon-diagonal contribution given in Eq. (C24) andEq. (C25) agrees with Eq. (150) in [20], where theygave the general form of U (cid:48) ij for the case of all externalcollinear particles: U (cid:48) ij = ln − n i · n j − i + . (C26)Next we repeat this procedure for (cid:80) (cid:104) ij (cid:105) S ij , where wehave non-relativistic fields in the initial state. As beforewe consider the contribution from the rapidity regulatorsat the outset, which for δ yield: α π ( T · T + T · T + T · T ) (cid:20)
12 ln m µ ln δ µ (cid:21) = − α π T · T (cid:20)
12 ln m µ ln δ µ (cid:21) , (C27)1where we again used Eq. (C21). An identical relation willhold for δ , and this time there is no δ or δ as the non-relativistic fields do not lead to rapidity divergences. Forthe remaining terms, we now organise them as follows: (cid:88) (cid:104) ij (cid:105) S ij = S + (cid:20) S + S + 12 S (cid:21) + (cid:20) S + S + 12 S (cid:21) . (C28)Evaluating each of the terms in square brackets and sim-plifying the gauge index structure as before, we arrive atthe following: (cid:88) (cid:104) ij (cid:105) S ij ≡ (cid:88) (cid:104) ij (cid:105) ˆ S ij + C + C , (C29)where we again have:ˆ S ij ≡ − α π ln m µ T i · T j U ij , (C30) C i ≡ α π T i · T i (cid:20)
14 ln m µ + π −
12 ln m µ ln δ i µ (cid:21) , and now U = 1 ,U = − iπ ,U = U = U = U = 0 . (C31)Critically, although the non-diagonal contribution is dif-ferent to the case of all collinear kinematics, we see thatthe diagonal function defined in Eq. (C30) is identical tothat in Eq. (C24). This justifies the claim made earlierthat the diagonal part of this calculation is the same forboth kinematics. As such we put the C i terms aside forthe moment, and return to them when we consider thefinal state collinear graphs.What remains here then is to evaluate the irreduciblynon-diagonal contribution: (cid:80) (cid:104) ij (cid:105) ˆ S ij . This essentiallyamounts to calculating the gauge index structure, whichthe use of gauge index operators has allowed us to putoff until now. In addition we need to recall that wehave a contribution to each graph from a W and Z bo-son exchange. As above we closely follow the approachin [20, 46], except accounting for the differences in ourkinematics. To this end, we begin by observing that af-ter electroweak symmetry breaking the unbroken SU(2) L and U(1) Y generators, t and Y , become α t · t + α Y · Y → α W ( t + t − + t − t + )+ α Z t Z · t Z + α em Q · Q , (C32)where α = α em /s W , α = α em /c W , α W = α , α Z = α /c W , and t Z = t − s W Q . This implies that we can write the full contribution as:ˆ D (1) s = α W ( µ )2 π ln m W µ − (cid:88) (cid:104) ij (cid:105)
12 ( t + t − + t − t + ) U ij + α Z ( µ )2 π ln m Z µ − (cid:88) (cid:104) ij (cid:105) t Zi t Zj U ij . (C33)Now the contribution on the first line is more compli-cated, because ( t + t − + t − t + ) U ij is a non-diagonal 22 × t Zi t Zj U ij is diagonal. Nev-ertheless we can simplify the non-diagonal part by usingthe following relation:12 ( t + t − + t − t + ) = t · t − t · t . (C34)Here t · t is again diagonal, and whilst t · t is non-diagonal, it is written in terms of the unbroken operatorsso that we can calculate it in the unbroken theory wherewe only have 2 operators not 22. Thus it is now a 2 × D s = ˆ D (0) s + ˆ D (1) s,W + ˆ D (1) s,Z , ˆ D (1) s,W = α W ( µ )2 π ln m W µ (cid:104) ˆ D (0) s S + D W ˆ D (0) s (cid:105) , ˆ D (1) s,Z = α Z ( µ )2 π ln m Z µ (cid:104) D Z ˆ D (0) s (cid:105) , (C35)where ˆ D (0) s is given in Eq. (C4) and as we will now demon-strate ˆ D s is effectively the matrix given in Eq. (15) thatwe set out to justify. In order to do this we have toevaluate the remaining terms: S ≡ − (cid:88) (cid:104) ij (cid:105) t i · t j U ij , D W ≡ (cid:88) (cid:104) ij (cid:105) t i · t j U ij , D Z ≡ − (cid:88) (cid:104) ij (cid:105) t Zi · t Zj U ij . (C36)The form of each of these matrices can be evaluated byacting with them on the operators – the unbroken op-erators in Eq. (C1) for S and the broken operators inEq. (C3) for D W/Z – where the action of the gauge indexoperators is given by Eq. (C10). Doing this, we find: S = (cid:20) − iπ − iπ iπ − (cid:21) , (C37)whilst D W, − = diag (0 , , , , − , − , − , − , D Z = − c W D W . (C38)2Substituting these results into Eq. (C35), we find:ˆ D s, − = s W [1 + G ( µ )] s W s W c W [1 + G ( µ )] s W c W s W c W [1 + G ( µ )] s W c W c W [1 + G ( µ )] c W s W [1 + H ( µ )] s W I ( µ ) s W c W [1 + H ( µ )] s W c W I ( µ ) s W c W [1 + H ( µ )] s W c W I ( µ ) c W [1 + H ( µ )] c W I ( µ ) , (C39)where we have defined: G ( µ ) ≡ α W ( µ )2 π ln m W µ (2 − iπ ) ,H ( µ ) ≡ α W ( µ )2 π ln m W µ (1 − iπ )+ c W α Z ( µ )2 π ln m Z µ ,I ( µ ) ≡ α W π ln m W µ (1 − iπ ) . (C40)From the form of ˆ D s given in Eq. (C39), we can againreduce this to a 2 × andΣ − Σ , exactly as we did for the tree-level low-scalematching. Doing this, the 2 × Final State Graphs
Finally we have the last contribution, which is the com-bination of final state collinear graphs as well as C + C ,as defined in Eq. (C30). As mentioned in the previoussubsection, this calculation has already been performedin [46], and given that the form of C i is the same for ourkinematics as it is for theirs, we take the result from theirwork. In that paper they calculated this diagonal contri-bution for all possible weak bosons. For our calculationwe are only interested in a final state photon or Z , forwhich they give: D Zc = α π (cid:20) F W + f S (cid:18) m Z m W , (cid:19)(cid:21) + 12 δ R Z + tan ¯ θ W R γ → Z ,D γc = α π [ F W + f S (0 , δ R γ + cot ¯ θ W R Z → γ . (C41)The various terms in these equations are outlined below.Nonetheless, once the full expressions are written out theanalytic result for the terms in Eq. (18) can be extractedas the terms independent of ln µ . To begin with we have: F W ≡ ln m W µ ln sµ −
12 ln m W µ − ln m W µ − π
12 + 1 , (C42)where note for our calculation s = 4 m χ . Next f S ( w, z )is defined as: f S ( w, z ) ≡ (cid:90) dx (2 − x ) x ln 1 − x + zx − wx (1 − x )1 − x , (C43)such that an explicit calculation gives us f S (cid:18) m Z m W , (cid:19) = 1 . ,f S (0 ,
1) = π − . (C44)Finally the R contributions are defined by: δ R Z ≡ Π (cid:48) ZZ ( m Z ) ,δ R γ ≡ Π (cid:48) γγ (0) , R γ → Z ≡ m Z Π Zγ ( m Z ) , R Z → γ ≡ − m Z Π γZ (0) , (C45)where Π (cid:48) ≡ ∂ Π( k ) /∂k and the various Π functionsare defined via the inverse of the transverse gauge bo-son propagator − i (cid:18) g µν − k µ k ν k (cid:19) (cid:20) k − m Z − Π ZZ ( k ) − Π Zγ ( k ) − Π γZ ( k ) k − Π γγ ( k ) (cid:21) . (C46)The form of the Π functions is not given explicitly in [46],but can be determined from the results of e.g. [70, 77].When doing so, there are two factors that must be ac-counted for. Firstly the Π functions must be calculatedin MS. This is because [46] accounts for the residuesexplicitly in (C41). If we used the on-shell scheme forexternal particles, as we did for the high-scale match-ing, we would double count the contribution from theresidues. Secondly we need to respect that the low-scalematching is performed above and below the electroweakscale, which means the Π functions for the photon and Z must be treated differently. Above the matching scalethe W , Z , top and Higgs are dynamical degrees of free-dom, but below it they are not. Light degrees of freedom Note there is a typo in Eq. B2 of [46], where R γ → Z and R Z → γ involved Π (cid:48) rather than Π. The expressions stated here are thecorrect ones, and we thank Aneesh Manohar for confirming thisand for providing a numerical cross check on our results for theseterms. Z contributions,we need to include all degrees of freedom – heavy andlight – in the loops, as the Z itself does not propagatebelow the matching and the light fermions are offshellin these loops. For the photon contributions, however,only the heavy degrees of freedom should be included.Accounting for these factors, we arrive at the following: δ R Z = α π (cid:20) − s W + 46 s W c W ln m Z µ Z +1 . − . i (cid:21) ,δ R γ = α π (cid:20) − s W ln m Z µ Z + 0 . (cid:21) , R γ → Z = α π (cid:20) − s W + 34 s W c W tan ¯ θ W ln m Z µ Z +0 . − . i (cid:21) , R Z → γ = α π (cid:20) s W c W ln m Z µ Z − . (cid:21) . (C47)Analytic forms for the Π functions are provided inApp. E, we do not provide the full expressions here asthey are lengthy. In order to determine the numericalvalues above we have used the following: m Z = 91 . ,m W = 80 .
385 GeV ,m t = 173 .
21 GeV ,m H = 125 GeV ,m b = 4 .
18 GeV ,m c = 1 .
275 GeV ,m τ = 1 . ,m s = m d = m u = m µ = m e = 0 GeV ,c W = m W /m Z . (C48)This completes the list of ingredients for Eq. (C41). Sub-stituting them into that equation gives exactly the rel-evant terms in Eqs. (16), (17), and (18), justifying thediagonal part of the low-scale matching. Note that theresults are insensitive to the precise values used for the m b and m c masses.We have now justified each of the pieces making up thelow-scale one-loop matching. All that remains is to crosscheck this result, which we turn to in the next appendix. Appendix D: Consistency Check on the Low-ScaleMatching
In this appendix we provide a cross check on the low-scale one-loop matching calculation, much as we did forthe high-scale result in App. B. Given that we already cross checked the high-scale result, we here make useof that to determine whether the ln µ contributions atthe low scale are correct. In order to do this, we takeEq. (8) and turn off the running, which amounts to set-ting µ m χ = µ Z ≡ µ . In detail we obtain: (cid:20) C X ± C X (cid:21) = e ˆ D X ( µ ) (cid:20) C ( µ ) C ( µ ) (cid:21) . (D1)Now as we have the full one-loop result, the ln µ depen-dence between these two terms must cancel at O ( α ) forany X , which we will now demonstrate.Before doing this in general, we first consider the sim-pler case where electroweak symmetry remains unbrokenand we just have a W W final state. In this case, asin general, to capture all µ dependence at O ( α ) we alsoneed to account for the β -function. If SU(2) L remainsunbroken, however, this is just simply captured in: α ( µ ) = α ( m Z ) + α ( m Z ) b π ln m Z µ , (D2)where b = (43 − n D ) /
6, with n D the number of SMdoublets. This follows directly from Eq. (A35). In theunbroken theory we can simply set c W = 1 and s W = 0,so if we do this and substitute our results from Eqs. (13),(14), (15), (16), (17) into Eq. (D1), then we find: C W ± = 1 m χ (cid:18) b c W − (cid:19) ln µ + µ − ind . ,C W = µ − ind . , (D3)Now we can calculate that c W = (2 n D − /
24, whichtaking n D = 12 exactly agrees with c Z in Eq. (17) when c W = 1 and s W = 0 as it must. Then recalling b fromabove we see that both coefficients are then µ indepen-dent at this order, demonstrating the required consis-tency.We now consider the same cross check in the full bro-ken theory. The added complication here is that for ourdifferent final states, γγ , γZ , and ZZ , the coupling isactually s W α , s W c W α , and c W α respectively. As wework in MS, we need to account for the fact that s W and c W are functions µ . Explicit calculation demonstratesthat the running is only relevant for the consistency of C X ± – the cancellation in C X is independent of the β -function at this order – and in fact we find: C X ± = 1 m χ (cid:32) b ( X )0 (cid:88) i ∈ X c i − (cid:33) ln µ + µ − ind . . (D4)To derive this we simply used Eq. (D2), with b → b ( X )0 ,leaving us to derive the appropriate form of b ( X )0 for X = γγ , γZ , ZZ . Firstly note that s W ( µ ) = α ( µ ) α ( µ ) + α ( µ ) ,c W ( µ ) = α ( µ ) α ( µ ) + α ( µ ) , (D5)4where α is the U(1) Y coupling. We can write a similarexpression to Eq. (D2) for α , but this time we have b (1)0 = − /
6. To avoid confusion we also now refer tothe SU(2) L b as b (2)0 = 19 / Z bosons in the final state, theappropriate β -function is: β ZZ = µ ddµ (cid:2) c W α (cid:3) . (D6)Combining this with Eq. (D5), we conclude that: b ( ZZ )0 = (cid:0) s W + 1 (cid:1) b (2)0 − s W c W b (1)0 = 19 + 22 s W c W . (D7)There is an additional factor of c W in this expressionthan if we were just calculating the β -function for α Z .The reason for this is that b ( ZZ )0 is the appropriate re-placement for b in Eq. (D2), which represents the cor-rection to α = c W α Z not α Z . Substituting this intoEq. (D4) along with the definition of c Z from Eq. (17)demonstrates consistency for the ZZ case.The case of two final state photons has to be treateddifferently, because of the fact our low-scale matching in-tegrated out the electroweak degrees of freedom, whichdid not include the photon. This means we need to usea modified version of the SU(2) L and U(1) Y couplingsthat only include the running due to the modes beingremoved. This amounts to accounting for the runningfrom the Higgs, W and Z bosons, and the top quark,which we treat as an SU(2) L singlet Dirac fermion to en-sure it is entirely removed through the matching. Doingso, the SM β -functions now evaluate to b (2) (cid:48) = 43 / b (1) (cid:48) = − /
18. Repeating the same calculation as weused to determine b ( ZZ )0 , we find that: b ( γγ )0 = (cid:16) b (1) (cid:48) + b (2) (cid:48) (cid:17) s W = 479 s W . (D8)Again, substituting this into Eq. (D4) shows that the twophoton case is also consistent. The final case γZ , but itis straightforward to see that in this case Eq. (D4) breaksinto two conditions that are satisfied if the ZZ and γγ cases are, so this is not an independent cross check.As such, in the absence of running, all the µ depen-dence in our calculation vanishes at O ( α ), as it must.But we emphasise that this is a non-trivial cross check,that involves all aspects of the calculation in the full bro-ken theory. Appendix E: Analytic Form of Π Here we state the analytic expressions for the MS elec-troweak Π functions for photon and Z boson, appropriatefor the matching from SCET EW to SCET γ . These re-sults can be determined using standard references, suchas [70, 77]. As the photon is a dynamical degree of free-dom above and below the matching, we only need to con-sider loop diagrams involving electroweak modes that areintegrated out through the matching. This simplifies theevaluation, and we have the following two functions:Π (cid:48) γγ (0) = α s W π (cid:26) −
169 ln µ m t + 3 ln µ m W + 23 (cid:27) , Π γZ (0) = α s W π (cid:26) m W s W c W ln µ m W (cid:27) . (E1)As the Z itself is being integrated out, we need to includeall relevant loops when calculating Π Zγ and Π (cid:48) ZZ . Inorder to simplify the expressions, we firstly introduce thefollowing expressions: β ≡ (cid:114) m s − , ξ ≡ (cid:114) − m s , (E2) λ ± ≡ s (cid:18) s − m + m ± (cid:113) ( s − m + m ) − s ( m − i(cid:15) ) (cid:19) . In terms of these we then define: a ( m , m ) ≡ m m − m ln m m ,b ( s, m ) ≡ iβ ln (cid:18) β + iβ − i (cid:19) ,b ( s, m ) ≡ − ξ ln 1 + ξ − ξ + iπξ ,c ( s, m ) ≡ − m s β (cid:18) β β + i ln β + iβ − i (cid:19) ,c ( s, m ) ≡ m s ξ (cid:18) ξξ − − ln 1 + ξ − ξ (cid:19) ,d ( s, m , m ) ≡ λ + ln (cid:18) λ + − λ + (cid:19) − ln ( λ + − λ − ln (cid:18) λ − − λ − (cid:19) − ln ( λ − − ,e ( s, m , m ) ≡ − s + ln (cid:18) λ + − λ + (cid:19) ∂λ + ∂s + ln (cid:18) λ − − λ − (cid:19) ∂λ − ∂s . (E3)5We can now write out the full expressions:Π Zγ ( m Z ) = α s W π (cid:26) − s W c W s W (cid:20) m Z − m Z ln µ m t − ( m Z + 2 m t ) b ( m Z , m t ) (cid:21) + 3 − s W c W s W (cid:20) m Z − m Z ln µ m b − ( m Z + 2 m b ) b ( m Z , m b ) (cid:21) + 6 − s W c W s W (cid:20) m Z − m Z ln µ m c − ( m Z + 2 m c ) b ( m Z , m c ) (cid:21) + 1 − s W c W s W (cid:20) m Z − m Z ln µ m τ − ( m Z + 2 m τ ) b ( m Z , m τ ) (cid:21) + m Z s W − c W s W (cid:20)
53 + iπ + ln µ m Z (cid:21) + 13 s W c W (cid:26)(cid:20)(cid:18) c W + 12 (cid:19) m Z + (cid:0) c W + 4 (cid:1) m W (cid:21) (cid:18) ln µ m W + b ( m Z , m W ) (cid:19) − (12 c W − m W ln µ m W + 13 m Z (cid:27)(cid:27) , (E4)and finallyΠ (cid:48) ZZ ( m Z ) = α s W π (cid:26) (cid:26) − s W + 32 s W c W s W (cid:20) − ln µ m t − b ( m Z , m t ) − ( m Z + 2 m t ) c ( m Z , m t ) + 13 (cid:21) + 34 s W c W m t c ( m Z , m t ) (cid:27) + 2 (cid:26) − s W + 8 s W c W s W (cid:20) − ln µ m b − b ( m Z , m b ) − ( m Z + 2 m b ) c ( m Z , m b ) + 13 (cid:21) + 34 s W c W m b c ( m Z , m b ) (cid:27) + 2 (cid:26) − s W + 32 s W c W s W (cid:20) − ln µ m c − b ( m Z , m c ) − ( m Z + 2 m c ) c ( m Z , m c ) + 13 (cid:21) + 34 s W c W m c c ( m Z , m c ) (cid:27) + 23 (cid:26) − s W + 8 s W c W s W (cid:20) − ln µ m τ − b ( m Z , m τ ) − ( m Z + 2 m τ ) c ( m Z , m τ ) + 13 (cid:21) + 34 s W c W m τ c ( m Z , m τ ) (cid:27) + 7 − s W + 16 s W s W c W (cid:20) − − ln µ m Z − iπ (cid:21) + 16 s W c W (cid:26)(cid:18) c W + 2 c W − (cid:19) (cid:18) ln µ m W + b ( m Z , m W ) (cid:19) + 13 (cid:0) c W − (cid:1) + (cid:20)(cid:18) c W + 2 c W − (cid:19) m Z + (cid:0) c W + 16 c W − (cid:1) m W (cid:21) c ( m Z , m W ) (cid:27) + 112 s W c W (cid:26) − (cid:18) ln µ m Z + d ( m Z , m Z , m H ) (cid:19) + (cid:0) m H − m Z (cid:1) e ( m Z , m Z , m H ) − ( m Z − m H ) m Z e ( m Z , m Z , m H ) −
23+ ( m Z − m H ) m Z (cid:18) ln m H m Z + d ( m Z , m Z , m H ) − a ( m Z , m H ) (cid:19)(cid:27)(cid:27) . (E5)6 [1] J. A. Hinton (HESS), New Astron. Rev. , 331 (2004),astro-ph/0403052.[2] A. Abramowski et al. (HESS), Phys. Rev. Lett. ,041301 (2013), 1301.1173.[3] G. Sinnis, A. Smith, and J. E. McEnery, in On recent de-velopments in theoretical and experimental general rela-tivity, gravitation, and relativistic field theories. Proceed-ings, 10th Marcel Grossmann Meeting, MG10, Rio deJaneiro, Brazil, July 20-26, 2003. Pt. A-C (2004), pp.1068–1088, astro-ph/0403096.[4] J. P. Harding and B. Dingus (HAWC), in
Proceedings,34th International Cosmic Ray Conference (ICRC 2015) (2015), 1508.04352.[5] J. Pretz (HAWC), in
Proceedings, 34th Interna-tional Cosmic Ray Conference (ICRC 2015) (2015),1509.07851.[6] M. Actis et al. (CTA Consortium), Exper. Astron. ,193 (2011), 1008.3703.[7] T. C. Weekes et al., Astropart. Phys. , 221 (2002),astro-ph/0108478.[8] J. Holder et al. (VERITAS), Astropart. Phys. , 391(2006), astro-ph/0604119.[9] A. Geringer-Sameth (VERITAS), in (2013), 1303.1406.[10] J. Flix Molina (MAGIC), in Proceedings, 40th Rencon-tres de Moriond on Very High Energy Phenomena in theUniverse (2005), pp. 421–424, astro-ph/0505313.[11] M. L. Ahnen et al. (Fermi-LAT, MAGIC), JCAP ,039 (2016), 1601.06590.[12] J. Hisano, S. Matsumoto, and M. M. Nojiri, Phys. Rev.Lett. , 031303 (2004), hep-ph/0307216.[13] J. Hisano, S. Matsumoto, M. M. Nojiri, and O. Saito,Phys. Rev. D71 , 063528 (2005), hep-ph/0412403.[14] M. Cirelli, A. Strumia, and M. Tamburini, Nucl. Phys.
B787 , 152 (2007), 0706.4071.[15] N. Arkani-Hamed, D. P. Finkbeiner, T. R. Slatyer, andN. Weiner, Phys. Rev.
D79 , 015014 (2009), 0810.0713.[16] K. Blum, R. Sato, and T. R. Slatyer, JCAP , 021(2016), 1603.01383.[17] M. Ciafaloni, P. Ciafaloni, and D. Comelli, Phys. Rev.Lett. , 4810 (2000), hep-ph/0001142.[18] P. Ciafaloni and D. Comelli, Phys. Lett. B446 , 278(1999), hep-ph/9809321.[19] P. Ciafaloni and D. Comelli, Phys. Lett.
B476 , 49 (2000),hep-ph/9910278.[20] J.-y. Chiu, A. Fuhrer, R. Kelley, and A. V. Manohar,Phys. Rev.
D80 , 094013 (2009), 0909.0012.[21] A. Hryczuk and R. Iengo, JHEP , 163 (2012), [Erra-tum: JHEP06,137(2012)], 1111.2916.[22] M. Baumgart, I. Z. Rothstein, and V. Vaidya, Phys. Rev.Lett. , 211301 (2015), 1409.4415.[23] M. Bauer, T. Cohen, R. J. Hill, and M. P. Solon, JHEP , 099 (2015), 1409.7392.[24] G. Ovanesyan, T. R. Slatyer, and I. W. Stewart, Phys.Rev. Lett. , 211302 (2015), 1409.8294.[25] M. Baumgart, I. Z. Rothstein, and V. Vaidya, JHEP ,106 (2015), 1412.8698.[26] M. Baumgart and V. Vaidya, JHEP , 213 (2016),1510.02470.[27] G. Elor, N. L. Rodd, T. R. Slatyer, and W. Xue, JCAP , 024 (2016), 1511.08787.[28] G. Elor, N. L. Rodd, and T. R. Slatyer, Phys. Rev. D91 ,103531 (2015), 1503.01773.[29] G. F. Giudice, M. A. Luty, H. Murayama, and R. Rat-tazzi, JHEP , 027 (1998), hep-ph/9810442.[30] L. Randall and R. Sundrum, Nucl. Phys. B557 , 79(1999), hep-th/9810155.[31] J. D. Wells, Phys. Rev.
D71 , 015013 (2005), hep-ph/0411041.[32] N. Arkani-Hamed and S. Dimopoulos, JHEP , 073(2005), hep-th/0405159.[33] G. F. Giudice and A. Romanino, Nucl. Phys. B699 ,65 (2004), [Erratum: Nucl. Phys.B706,487(2005)], hep-ph/0406088.[34] A. Arvanitaki, N. Craig, S. Dimopoulos, and G. Vil-ladoro, JHEP , 126 (2013), 1210.0555.[35] N. Arkani-Hamed, A. Gupta, D. E. Kaplan, N. Weiner,and T. Zorawski (2012), 1212.6971.[36] L. J. Hall, Y. Nomura, and S. Shirai, JHEP , 036(2013), 1210.2395.[37] T. Cohen, M. Lisanti, A. Pierce, and T. R. Slatyer, JCAP , 061 (2013), 1307.4082.[38] P. Ciafaloni, D. Comelli, A. De Simone, A. Riotto, andA. Urbano, JCAP , 016 (2012), 1202.0692.[39] C. W. Bauer, S. Fleming, and M. E. Luke, Phys. Rev. D63 , 014006 (2000), hep-ph/0005275.[40] C. W. Bauer, S. Fleming, D. Pirjol, and I. W. Stewart,Phys. Rev.
D63 , 114020 (2001), hep-ph/0011336.[41] C. W. Bauer and I. W. Stewart, Phys. Lett.
B516 , 134(2001), hep-ph/0107001.[42] C. W. Bauer, D. Pirjol, and I. W. Stewart, Phys. Rev.
D65 , 054022 (2002), hep-ph/0109045.[43] J.-y. Chiu, F. Golf, R. Kelley, and A. V. Manohar, Phys.Rev. Lett. , 021802 (2008), 0709.2377.[44] J.-y. Chiu, F. Golf, R. Kelley, and A. V. Manohar, Phys.Rev.
D77 , 053004 (2008), 0712.0396.[45] J.-y. Chiu, R. Kelley, and A. V. Manohar, Phys. Rev.
D78 , 073006 (2008), 0806.1240.[46] J.-y. Chiu, A. Fuhrer, R. Kelley, and A. V. Manohar,Phys. Rev.
D81 , 014023 (2010), 0909.0947.[47] C. W. Bauer, S. Fleming, D. Pirjol, I. Z. Rothstein, andI. W. Stewart, Phys. Rev.
D66 , 014017 (2002), hep-ph/0202088.[48] J.-Y. Chiu, A. Jain, D. Neill, and I. Z. Rothstein, JHEP , 084 (2012), 1202.0814.[49] T. Becher and M. Neubert, Eur. Phys. J. C71 , 1665(2011), 1007.4005.[50] A. Denner, L. Jenniches, J.-N. Lang, and C. Sturm,JHEP , 115 (2016), 1607.07352.[51] L. Bergstrom, G. Bertone, J. Conrad, C. Farnier, andC. Weniger, JCAP , 025 (2012), 1207.6773.[52] J. Ellis (2016), 1601.05437.[53] F. Boudjema, A. Semenov, and D. Temes, Phys. Rev. D72 , 055024 (2005), hep-ph/0507127.[54] N. Baro, F. Boudjema, and A. Semenov, Phys. Lett.
B660 , 550 (2008), 0710.1821.[55] N. Baro, F. Boudjema, G. Chalons, and S. Hao, Phys.Rev.
D81 , 015005 (2010), 0910.3293.[56] A. V. Manohar and I. W. Stewart, Phys. Rev.
D76 ,074002 (2007), hep-ph/0605001.[57] A. Fuhrer, A. V. Manohar, J.-y. Chiu, and R. Kelley, Phys. Rev.
D81 , 093005 (2010), 1003.0025.[58] G. Passarino and M. J. G. Veltman, Nucl. Phys.
B160 ,151 (1979).[59] A. Denner, H. Eck, O. Hahn, and J. Kublbeck, Nucl.Phys.
B387 , 467 (1992).[60] R. K. Ellis and G. Zanderighi, JHEP , 002 (2008),0712.1851.[61] G. ’t Hooft and M. J. G. Veltman, Nucl. Phys. B153 ,365 (1979).[62] R. K. Ellis, Z. Kunszt, K. Melnikov, and G. Zanderighi,Phys. Rept. , 141 (2012), 1105.4319.[63] R. Mertig, M. Bohm, and A. Denner, Comput. Phys.Commun. , 345 (1991).[64] V. Shtabovenko, R. Mertig, and F. Orellana (2016),1601.01167.[65] H. H. Patel, Comput. Phys. Commun. , 276 (2015),1503.01469.[66] M. Beneke and V. A. Smirnov, Nucl. Phys. B522 , 321(1998), hep-ph/9711391.[67] M. E. Luke, A. V. Manohar, and I. Z. Rothstein, Phys.Rev.
D61 , 074025 (2000), hep-ph/9910209.[68] N. Brambilla, A. Pineda, J. Soto, and A. Vairo, Rev.Mod.Phys. , 1423 (2005), hep-ph/0410047.[69] A. V. Manohar and I. W. Stewart, Phys. Rev. D62 ,074015 (2000), hep-ph/0003032.[70] D. Yu. Bardin and G. Passarino,
The Standard Model inthe Making: Precision Study of the Electroweak Interac-tions (1999).[71] I. W. Stewart, Nucl. Phys.
B529 , 62 (1998), hep-ph/9803227.[72] A. V. Manohar and M. B. Wise, Camb. Monogr. Part.Phys. Nucl. Phys. Cosmol. , 1 (2000).[73] S. Catani and M. H. Seymour, Phys. Lett. B378 , 287(1996), hep-ph/9602277.[74] S. Catani and M. H. Seymour, Nucl. Phys.
B485 ,291 (1997), [Erratum: Nucl. Phys.B510,503(1998)], hep-ph/9605323.[75] I. Moult, I. W. Stewart, F. J. Tackmann, and W. J.Waalewijn, Phys. Rev.
D93 , 094003 (2016), 1508.02397.[76] J.-y. Chiu, A. Fuhrer, A. H. Hoang, R. Kelley, and A. V.Manohar, Phys. Rev.
D79 , 053007 (2009), 0901.1332.[77] A. Denner, Fortsch. Phys.41