Next-to-leading-log renormalization-group running in heavy-quarkonium creation and annihilation

In the framework of potential NRQCD, we obtain the next-to-leading-log renormalization-group running of the matching coefficients for the heavy quarkonium production currents near threshold. This allows to obtain S-wave heavy-quarkonium production/annihilation observables with next-to-leading-log accuracy within perturbative QCD. In particular, we give expressions for the decays of heavy quarkonium to e^+e^- and to two photons. We also compute the O(m\alpha^8\ln^3\alpha) corrections to the Hydrogen spectrum.

Heavy quark-antiquark systems near threshold are characterized by the small relative velocity v of the heavy quarks in their center of mass frame. This small parameter produces a hierarchy of widely separated scales: m (hard), mv (soft), mv 2 (ultrasoft), ... . The factorization between them is efficiently achieved by using effective field theories, where one can organize the calculation as various perturbative expansions on the ratio of the different scales effectively producing an expansion in v. The terms in these series get multiplied by parametrically large logs: ln v, which can also be understood as the ratio of the different scales appearing in the physical system. Again, effective field theories are very efficient in the resummation of these large logs once a renormalization group (RG) analysis of them has been performed. This will be the aim of this paper for annihilation and production processes near threshold.
We will restrict ourselves, in this paper, to the situation where Λ QCD ≪ mα 2 s (to be implicit in what follows), which is likely to be relevant, at least, for t-t production near threshold.
NRQCD [1] has an ultraviolet (UV) cut-off ν N R = {ν p , ν s } satisfying mv ≪ ν N R ≪ m. At this stage ν p ∼ ν s . ν p is the UV cut-off of the relative three-momentum of the heavy quark and antiquark, p. ν s is the UV cut-off of the three-momentum of the gluons and light quarks.
Potential NRQCD (pNRQCD) [2] (see [3,4] for details) is defined by its particle content and cut-off ν pN R = {ν p , ν us }, where ν p is the cut-off of the relative three-momentum of the heavy quarks and ν us is the cut-off of the three-momentum of the gluons and light quarks. They satisfy the following inequalities: |p| ≪ ν p ≪ m and p 2 /m ≪ ν us ≪ |p|. Note that no gluons or light quarks with momentum of O(|p|) are kept dynamical in pNRQCD. The motivation to integrate out these degrees of freedom is that they do not appear as physical (on-shell) states near threshold. Nevertheless, they can appear off-shell and, since their momentum is of the order of the relative three-momentum of the heavy quarks, integrating them out produces non-local terms (potentials) in three-momentum space. Indeed, these potentials encode the non-analytical behavior in the transfer momentum of the heavy quark, k = p − p ′ , of the order of the relative three-momentum of the heavy quarks.
The matching process, which basically means the computation of the potentials, is carried out for a given external incoming (outcoming) momentum p (p ′ ). Therefore, one has to sum over all of them in the pNRQCD Lagrangian, since they are still physical degrees of freedom as far as their momentum is below ν p . In position space, this means that an integral over x, the relative distance between the heavy quarks, appears in the Lagrangian when written in terms of the heavy quark-antiquark bilinear.
Within pNRQCD, integrals over p (or x) appear when solving the Schrödinger equation that dictates the dynamics of the heavy quarkonium near threshold. At lower orders, these integrals are finite effectively replacing p by ∼ mα s . Nevertheless, at higher orders in quantum mechanics perturbation theory and/or if some singular enough operators are introduced (as it will be the case of the heavy quarkonium production currents) singularities proportional to ln ν p appear. These must be absorbed by the potentials or by the matching coefficients of the currents. We will describe how to resum the logarithms associated to this cutoff within pNRQCD.
A RG analysis for non-relativistic systems have been addressed before in Refs. [5][6][7], where they match to an effective theory called vNRQCD. On physical terms, this theory should be equivalent to the previously defined pNRQCD once the RG evolution has been performed and the soft degrees of freedom have been integrated out, as only ultrasoft gluons and light fermions and potential quarks are left as dynamical degrees of freedom. We will compare with their results. In some cases disagreement will be found.
Let us now describe the matching between QCD and pNRQCD within an RG framework. For the case where no divergences proportional to ln ν p appear, the procedure reduces to the results of Ref. [8] to which we refer for the notation and background material necessary to follow this paper.
We first address the procedure that gives the running of the potentials. One first does the matching from QCD to NRQCD. The latter depends on some matching coefficients: c(ν s ) and d(ν p , ν s ), which can be obtained order by order in α s (with ν p = ν s ) following the procedure described in Ref. [9]. The c(ν s ) stand for the coefficients of the operators that already exist in the theory with only one heavy quark (ie. HQET) and the d(ν p , ν s ) stand for the coefficients of the four heavy fermion operators. The starting point of the renormalization group equation can be obtained from these calculations by setting ν p = ν s = m (up to a constant of order one). In principle, we should now compute the running of ν p and ν s . The running of the c(ν s ) can be obtained using HQET techniques [10]. The running of the d(ν p , ν s ) is more complicated. At one-loop, ν p does not appear and we effectively have d(ν p , ν s ) ≃ d(ν s ), whose running can also be obtained using HQET-like techniques [8]. At higher orders, the dependence on ν p appears and the running of the d(ν p , ν s ) becomes more complicated. Fortunately, we need not compute the running of d in this more general case because, as we will see, the relevant running of the d for near threshold observables can be obtained within pNRQCD. The next step is the matching from NRQCD to pNRQCD. The latter depends on some matching coefficients (potentials). They typically have the following structure:Ṽ (c(ν s ), d(ν p , ν s ), ν s , ν us , r). After matching, any dependence on ν s disappears since the potentials have to be independent of ν s . Therefore, they could be formally written asṼ (c(1/r), d(ν p , 1/r), 1/r, ν us , r). These potentials can be obtained order by order in α s following the procedure of Refs. [2][3][4]. The integrals in the matching calculation would depend on a factorization scale µ, which should correspond either to ν s or to ν us . In the explicit calculation, they could be distinguished by knowing the UV and infrared (IR) behavior of the diagrams: UV divergences are proportional to ln ν s , which should be such as to cancel the ν s scale dependence inherited from the NRQCD matching coefficients, and IR divergences to ν us . In practice, however, as far as we only want to perform a matching calculation at some given scale µ = ν s = ν us , it is not necessary to distinguish between ν s and ν us (or if working order by order in α s without attempting any log resummation).
Before going into the rigorous procedure to obtain the RG equations of the potentials, let us first discuss their structure on physical grounds. As we have mentioned, the potential is independent of ν s . The independence of the potential with respect ν s allows us to fix the latter to 1/r that, in a way, could be understood as the matching scale for ν s 1 . Therefore, 1/r, the point where the multipole expansion starts, would also provide with the starting point of the renormalization group evolution of ν us (up to a constant of order one). The running of ν us can then be obtained following the procedure described in Refs. [11,8]. At the end of the day, we would haveṼ (c(1/r), d(ν p , 1/r), 1/r, ν us , r), where the running on ν us is known and also the running in 1/r if the d is ν p -independent. So far, the only explicit dependence of the potential on ν p appears in the d. Nevertheless, the potential is also implicitly dependent on the three-momentum of the heavy quarks through the requirement 1/r ∼ p ≪ ν p , and also through ν us , since ν us needs to fulfill p 2 /m ≪ ν us ≪ |p|. This latter requirement holds if we fix ν us = ν 2 p /m (this constraint tells you how much you can run down ν us in the potential before finding the cutoff ν 2 p /m caused by the cutoff of p).
Within pNRQCD, the potentials should be introduced in the Schrödinger equation. This means that integrals over the relative three-momentum of the heavy quarks take place. When these integrals are finite one has p ∼ 1/r ∼ mα s and p 2 /m ∼ mα 2 s . Therefore, one can lower ν us down to ∼ mα 2 s reproducing the results obtained in Ref. [8]. In some cases, in particular in heavy quarkonium creation, the integrals over p are divergent, and the log structure is dictated by the ultraviolet behavior of p and 1/r. This means that we can not replace 1/r and ν us by their physical expectation values but rather by their cutoffs within the integral over p. Therefore, for the RG equation of ν p , the anomalous dimensions will depend (at leading order) onṼ (c(ν p ), d(ν p , ν p ), ν p , ν 2 p /m, ν p ) 2 and the running will go from ν p ∼ m down to ν p ∼ mα s . Note that, at this stage, a single cutoff, ν p , exists and the correlation of cutoffs can be seen. The importance of the idea that the cutoffs of the non-relativistic effective theory should be correlated was first realized by Luke, Manohar and Rothstein in Ref. [5] (for an application to QED see [12]). Note also that at the matching scale ν p ∼ m, what it would be the ultrasoft cutoff is also of order m. In this sense it should be understood the statement in Ref. [5] that ultrasoft gluons appear at the scale m, a point that becomes relevant within a RG approach.
With the above discussion in mind, the matching between NRQCD and pNRQCD could be thought as follows. One does the matching by computing the potentials order by order in α s at 1 In practice, the potential is often first obtained in momentum space so that one could then set ν s = k. Note, however, that this is not equivalent to fix ν s = 1/r, since finite pieces will appear after performing the Fourier transform.
2 Roughly speaking, this result can be thought as expanding ln r around ln ν p in the potential ie.
The ln(ν p r) terms may give subleading contributions to the anomalous dimension when introduced in divergent integrals over p. The discussion at this stage is not very rigorous and a more precise discussion would require a full detailed study within dimensional regularization, which goes beyond the aim of this work. Nevertheless, we do not expect it to change the underlying idea, although it deserves further investigations.
the matching scale ν p = ν s = ν us following the procedure of Refs. [2][3][4] (by doing the matching at a generic ν p some of the running is trivially obtained). The structure of the potential at this stage then readsṼ (c(ν p ), d(ν p , ν p ), ν p , ν p , ν p ) (and similarly for the derivatives with respect ln r of the potential). This provides the starting point of the renormalization group evolution of ν us (up to a constant of order one). The running of ν us can then be obtained following the procedure described in Refs. [11,8]. For the final point of the evolution of ν us , we choose ν us = ν 2 p /m. At the end of the day, we obtainṼ (c(ν p ), d(ν p , ν p ), ν p , ν 2 p /m, ν p ) ≡Ṽ (ν p ).
The running of ν p goes from ν p = m (this was fixed when the matching between QCD and NRQCD was done) up to the physical scale of the problem ν p ∼ mα s . If the running of the NRQCD matching coefficients is known, the above result gives the complete running of the potentials. The procedure to get the running of the c is known at any finite order. For the d it is just known at one-loop order, since, at this order, it is only ν s -dependent. Nevertheless, at higher orders, dependence on ν p appears. Therefore, the above method is not complete unless an equation for the running of ν p is provided. This is naturally given within pNRQCD. It appears through the iteration of potentials. Let us consider this situation more in detail.
We first remind what the Hamiltonian in pNRQCD for the singlet sector is (see Ref. [8] for notation and further details): where C f = (N 2 c − 1)/(2N c ) and we will set c k = c 4 = 1 (we will only eventually use c 4 for tracking of the contribution due to this term). The propagator of the singlet is (formally) At leading order (within an strict expansion in α s ) the propagator of the singlet reads 0000000000000000000000000000 0000000000000000000000000000 1111111111111111111111111111 1111111111111111111111111111 If we were interested in computing the spectrum at O(mα 6 s ), one should consider the iteration of subleading potentials (δh s ) in the propagator as follows: In general, if these potentials are singular enough, these contributions will produce logarithmic divergences due to potential loops. These divergences can be absorbed in the matching coefficients, D d,s and D S 2 ,s , of the local potentials (those proportional to the δ (3) (r)) providing with the renormalization group equations of these matching coefficients in terms of ν p . Let us explain how it works in detail. Since the singular behavior of the potential loops appears for |p| ≫ α s /r, a perturbative expansion in α s is licit in G c (E), which can be approximated by Therefore, with the accuracy we aim at in this paper, a practical simplification follows from the fact that the Coulomb potential, −C f αs r , can be considered to be a perturbation as far as the computation of the ln ν p ultraviolet divergences is concerned. This means that the computation of the anomalous dimension can be organized within an expansion in α s and using the free propagators G (0) c . Moreover, each G (0) c produces a potential loop and one extra power of m in the numerator, which kills the powers in 1/m of the different potentials. This allows the mixing of potentials with different powers in 1/m. One typical example would be the diagram in Fig. 1. 3 The computation of this diagram would go as follows: Using δ (3) (r) = |r = 0 r = 0|, we can see that the relevant computation reads (instead of α Vs one could use α s since the non-trivial running of α Vs is a subleading effect. Nevertheless, we keep α Vs since it allows to keep track of the contributions due to the Coulomb potentials) where D = 4 + 2ǫ and q = p − p ′ . This divergence is absorbed in D d,s contributing to its running at next-to-leading-log (NLL) order as follows Therefore, even without knowing the running of the d (which need to be known at NLL order in this case), we can obtain the running of the potential (one can also think of trading Eq. (7) into an equation for d, which is the only unknown parameter within the potential). This is so because D d,s is only needed with LL accuracy in the right-hand side of Eq. (7).
The above method deals with the resummation of logs due to the hard, soft and ultrasoft scales. Nevertheless, for some specific kinematical situations even smaller scales could appear. Their study, however, goes beyond the aim of this paper. In any case, pNRQCD can be considered to be the right starting point to study these kinematical situations.
The matching scale between QCD and NRQCD is ν p ∼ ν s ∼ m. On the other hand, the matching scale between NRQCD and pNRQCD is also the hard scale: ν p ∼ ν 2 p /m ∼ m. Therefore, one could wonder about the necessity of using the intermediate theory NRQCD. This is indeed the attitude in Refs. [5,12,6,7], where they directly perform the matching between QCD to an effective field theory: vNRQCD that, once the RG evolution has been performed and the soft degrees of freedom have been integrated out, should be physically equivalent to pNRQCD with ν p ∼ mα s . One motivation for going through NRQCD is that it allows to perform the factorization of the hard scale within an effective field theory framework. In fact, a full factorization of the different regions of momentum that ought to be integrated out is achieved within pNRQCD. This extremely simplifies the matching process since one deals with only one scale (region of momentum) in the loops at each step. In the matching between QCD and NRQCD only hard loops need to be considered, whereas in the matching between NRQCD and pNRQCD only soft loops need to be considered. Moreover, the structure of the UV cutoffs of the theory is better understood in this way. For instance, one can see that all the explicit dependence of the potentials on ν p is inherited from the d matching coefficients. Within a diagrammatic approach the factorization of the different regions of momentum have been achieved using the threshold expansion [13].
Let us now consider the case of the electromagnetic current, which will provide an example where to apply the above discussion. The procedure is analogous to the potentials. We first do the matching from QCD to NRQCD: We will just concentrate in the coefficient b 1,NR . Within NRQCD, it should be understood as a function of ν p and ν s , ie. b 1,NR (ν p , ν s ). One should first obtain the matching conditions at the hard scale. This has been computed up to two loops [14] but we will only need the one-loop expression [15]: since we only aim to a NLL resummation in this paper. If we compare with the previous discussion of the potentials, the matching coefficients d play the role of b 1 . Therefore, within pNRQCD, we will need b 1,NR (ν p , ν p ) ≡ b 1,NR (ν p ). We first have to consider b 1,NR (ν p , ν s ). In this case, unlike for the d's, there is no running due to ν s at the order of interest. This can be easily seen in the Coulomb gauge. Moreover, the matching from NRQCD to pNRQCD creates the potentials but let b 1 unchanged since soft loops or HQET-like calculations give zero correction to b 1 at the order of interest. Formally, or, in other words, the matching condition reads B 1,pNR (b 1,NR (ν p ), ν us = ν p ) = b 1,NR (ν p ). The running of ν us is also trivial as there is none at the order of interest (this has to do with the fact that we are dealing with an annihilation process). Therefore, we finally have We can see that we are in the analogous situation to the running of D (2) d,s (ν p ) versus the running of d(ν p , ν p ). We now need the RG equation for B 1 (ν p ). This demands to obtain the ultraviolet corrections to the current within pNRQCD keeping track of the contributions due to the different potentials. Fortunately, this calculation has already been done and we can extract the relevant information from Ref. [16]. The computation goes along the same lines than in the example of Fig. 1. The explicit diagrams to be computed for the RG running of B s (ν p ) are given in Fig. 2 (where s  denotes the spin). From this figure, we can clearly illustrate the structure of the computation. O(1/m) corrections to h (0) s only need one potential loop to kill the 1/m coefficient. O(1/m 2 ) corrections to h (0) s need two potential loops to kill the 1/m 2 coefficient and so on. In the situation with more than one potential loop, the additional potential loops can be produced without additional 1/m factors coming from the potential only if Coulomb potentials are introduced. This explains why the 1/m potential needs zero Coulomb potential insertions, the 1/m 2 potentials need one Coulomb potential insertions and the 1/m 3 term needs two Coulomb potential insertions (for the running of D S 2 ,s we expect a similar structure). In principle, this would be be a never ending story unless there is an small parameter that tells us how far we have to go in the calculation in order to achieve some given accuracy. This is indeed so. The 1/m potential is a NLL effect [8] and therefore higher powers in D (1) s produce NNLL effects or beyond. On the other hand, the introduction of Coulomb potentials brings powers in α s , which suppresses the order of the calculation. In our case, for a NLL calculation, the maximum power of the anomalous dimension should be α 2 s . This means that with zero α Vs insertions (O(1/m) potentials) the matching coefficient (D (1) s ) has to be known with NLL accuracy, with one α Vs insertion (O(1/m 2 ) potentials) the matching coefficients (D (2) ) have to be known with LL accuracy and with two α Vs insertions (O(1/m 3 ) potentials) the matching coefficients must have no running (this explains why only c 4 is considered at this order).
From the above discussion, the RG equation reads 4 where C A = N c , and the RG-improved matching coefficients of the potentials can be read from Ref. [11,8] with the assignment 1/r → ν p and ν us → ν 2 p /m (see also [17,4] for calculations of the potentials at finite orders in α s ). We have kept the spin explicit so our results will also be valid for the pseudoscalar current: with the matching condition [18]: Eq. (11) gives subleading effects within an strict expansion in α s . Therefore, it can be approximated to and the solution reads where Our evaluation can be compared with the result obtained using the vNRQCD formalism [7]. We agree for the spin-dependent terms but differ for the spin-independent ones. The disagreement still holds if we consider QED with light fermions (C f → 1, C A → 0, T F → 1). Agreement is found if we consider QED without light fermions (C f → 1, C A → 0, n f → 0, T F → 1). If we expand our results in α s , we can compare with earlier results in the literature. By following the discussion in Ref. [7], we can relate our results with the correction to the wave-function at the origin as defined in Ref. [16]. We obtain where we have expanded up to second order in ln(ν p ) = ln(mα s ) with α s ≡ α s (ν p ). The first term reproduces the leading log term [14,19] (see also [20]), the β 0 -independent O(α 3 s ln 2 α s ) terms reproduce Kniehl and Penin results [16] and we agree with the complete O(α 3 s ln 2 α s ) term computed by Manohar and Stewart [7] (the sign of difference for the β 0 -dependent terms displayed in Ref. [7] is due to the fact that in Ref. [7] the expansion was made with α s (m) whereas here we have chosen α s (mα s )). Nevertheless, disagreement with this last evaluation appears at higher orders in the expansion in α s (we have explicitly checked this for the O(α 4 s ln 3 α s ) terms). As far as we can see, the disagreement seems to be due to the fact that they have different expressions for the RG improved potentials [6,7]. 5 By setting ν p ∼ mα s , B s (ν p ) includes all the large logs at NLL order in any (inclusive enough) S-wave heavy-quarkonium production observable we can think of. For instance, the decays to e + e − and to two photons at NLL order read where V and P stand for the vector and pseudoscalar heavy quarkonium, we have fixed ν p = m Q C f α s /n, α s = α s (ν p ), and (Ψ n (z) = d n ln Γ(z) dz n and Γ(z) is the Euler Γ-function) which has been read from Ref. [19] (see also [23]). Working along similar lines one could 5 The running of the Coulomb potential is not needed for the precision of the above calculation. Nevertheless, the running obtained in vNRQCD [21] also disagrees with the one obtained in pN-RQCD [11]. At this respect, we would like to report on a recent computation [22] of the 4-loop double log term of the Coulomb potential proportional to C 3 A β 0 that agrees with the pNRQCD result and disagrees with the vNRQCD one. easily obtain NLL expressions for other heavy quarkonium observables in the study of t-t production near threshold or in sum rules of bottomonium. Note that for t-t production near threshold there already exists a (partial) NNLL RG improved evaluation within the vNRQCD formalism [24]. Since we disagree for the RG improved expression for the electromagnetic current matching coefficient, this discrepancy would also propagate to that evaluation.
In conclusion, by using the method of Ref. [8] and incorporating the idea [5] of correlating the cut-offs of the effective theory, we have given the first steps towards the creation of a comprehensive system of RG equations in pNRQCD once the scale ν p enters into the game. We have used this formalism to compute the running of the matching coefficients of the vector and pseudoscalar currents and disagreement with the results obtained using the vNRQCD framework [7] has been found. Our results allow to obtain S-wave heavy quarkonium production observables with NLL accuracy. We have explicitly illustrated this point for heavy-quarkonium decays to e + e − and to two photons. We have also computed the O(mα 8 ln 3 α) corrections to the Hydrogen spectrum in the Appendix.
With the above discussion, we may also try to see whether we are able to obtain the O(mα 8 ln 3 α) contributions to the Hydrogen spectrum. It goes beyond the scope of this paper to perform a detailed analysis. Here, we will just see that under some assumptions from Ref. [12], we are able to obtain the O(mα 8 ln 3 α) correction to the Hydrogen spectrum. We will use the notation of Ref. [25]. In that paper, Hydrogen-like atoms with light fermions were considered. Here, we will use their results on the strict Hydrogen limit (no light fermions: n f = 0).
According to Ref. [12], the O(mα 8 ln 3 α) correction to the Hydrogen energy can be obtained from the anomalous dimension due to diagrams of the type of Fig 2.a in Ref. [12] or, in our case, to the diagram in Fig. 1. The argument that led to this conclusion was that the O(mα 8 ln 3 α) terms had the highest possible log power that could appear from a NNNLL evaluation of the energy and that, in order to achieve such power, it was necessary to mix with NNLL logs. The latter only appear in the LL evaluation of D (2) d , which, indeed, only produces a single log [12] (see also [8]). The other point was that the NLL evaluation of the potentials would only produce single logs unless mixed with LL running. Therefore, the diagrams with the highest possible power of D d,s at NLL. The first picture represents the calculation in terms of the free quark-antiquark propagator G (0) c and the small rectangles the potentials. The picture below is the representation within a more standard diagrammatic interpretation in terms of quarks and antiquarks. The delta potentials are displayed as local interactions and the Coulomb potential as an extended in space (but not in time) object.