Theory for quarkonium: from NRQCD factorization to soft gluon factorization

We demonstrate that the recently proposed soft gluon factorization (SGF) is equivalent to the nonrelativistic QCD (NRQCD) factorization for heavy quarkonium production or decay, which means that, for any given process, these two factorization theories are either both valid or both violated. We use two methods to arrive at this conclusion. In the first method, we apply the two factorization theories to the physical process . Our explicit calculation shows that both SGF and NRQCD can correctly reproduce the low energy physics of full QCD, and the two factorizations are thus equivalent. In the second method, by using equations of motion, we successfully deduce SGF from NRQCD effective field theory. By identifying SGF with NRQCD factorization, we establish relations between the two factorization theories and prove the generalized Gremm-Kapustin relation as a byproduct. Compared with the NRQCD factorization, the advantage of SGF is that it resums the series of relativistic corrections originating from kinematic effects to all powers, yielding better convergence of the relativistic expansion.


I. INTRODUCTION
The widely used nonrelativistic QCD (NRQCD) factorization [1] has encountered notable difficulties in describing heavy quarkonium data. As the NRQCD factorization is based on the NRQCD effective field theory [2], it is likely rigorous, although for inclusive quarkonium production, only two-loop verification is available at present [3][4][5]. The main known problem of NRQCD factorization is the bad convergence of its relativistic expansion [6], which may be responsible for its deficiencies. Recently, a new factorization approach called soft gluon factorization (SGF) [7,8] was proposed to describe quarkonium production and decay. The aim of SGF is to resum the series of relativistic corrections originating from the kinematic effects in NRQCD, which are the main cause of the poor convergence of the relativistic expansion.
Nevertheless, SGF has not been well-established. In SGF, the hadronization of intermediate quark-antiquark pairs to physical quarkonium is described by nonperturb-ative soft gluon distributions (SGDs), which are only formally defined by QCD fields in a small loop momentum region [7]. Without an explicit definition of a small region, it is hard to prove the validity of SGF for physical processes. Furthermore, the unclear relation between SGF and NRQCD factorization makes it impossible to verify whether the kinematic effects have been correctly resummed.
In this paper, with the help of a new regulator, we provide a rigorous definition of a small region in SGF. We then present two strategies for exploring the relationship between SGF and NRQCD factorization. In the first strategy, we apply the two factorization theories to the physical process of and show that both SGF and NRQCD factorization can correctly reproduce all the low energy physics of full QCD in this process. In the second strategy, we argue that the SGF formula can be deduced from NRQCD at the operator level by using equations of motion. Both strategies demonstrate that SGF and NRQCD factorization are equivalent, which means that, for any process, the two factorizations theor-ies are either both valid or both violated. By identifying the two theories, we generate complete relations between the nonperturbative matrix elements in SGF and NR-QCD; this proves the generalized Gremm-Kapustin relation [9] as a byproduct.
The rest of this paper is organized as follows. In Sec. II, we study the exclusive process in SGF. We give a rigorous definition of nonperturbative matrix elements in SGF and show that the low energy physics of full QCD can be reproduced by SGF. In Sec. III, we discuss the equivalence between SGF and NRQCD factorization and establish relations between the nonperturbative matrix elements in SGF and NRQCD. In Sec. IV, we show that SGF can be deduced from NRQCD at the operator level. We summarize our results in Sec. V. Some of the technical details of our calculations are provided in Appendix A.

IZATION APPROACH
A. Factorization formula According to Ref. [8], for the exclusive decay process, we have the following SGF formula at the amplitude level: Ψ K n n where stands for the Dirac field of a heavy quark, is the projection operator defining the intermediate state , and the subscript "S" indicates that, to evaluate the matrix elements, one only includes integration regions where the off-shellness of all particles is much smaller than the heavy-quark mass. From the point view of the method of regions [10,11], the effect of "S" is to keep only small regions, which include everything except the hard regions.
For the process , the symmetries of QCD tell us that only are relevant, where the spectroscopic notation, with the superscript "[1]", denotes a color singlet. We thus have with the projection operators defined explicitly as , is the momentum of , is a polarization vector with , is the heavy-quark mass, is the mass of , the color projector is defined as , and the spin projection operator is defined as The hard parts can be perturbatively calculated according to the matching procedure discussed in [7,8].
To this end, we first replace with an on-shell colorsinglet state with momenta 1) on both sides of Eq. (3). The on-shell conditions result in q 0 |q| P q S which fixes and in the rest frame of . The remaining degrees of freedom of are removed by partial wave expansion, i.e., -wave expansion in this case. After the replacement, the l.h.s. of Eq. (3) becomes , An-Ping Chen, Yan-Qing Ma Chin. Phys. C 45, 013118 (2021) cc P 1) Note that the total momentum of the pair is fixed to the momentum of the physical quarkonium during the matching procedure in SGF. This is significantly different from the matching procedure in NRQCD factorization, where the total momentum of the pair is a free parameter.
where is the solid angle of the relative momentum in the rest frame, is the leptonic current with A µ cc cc cc 1 cc and is the hadronic part of the decay amplitude with the spinors of removed. The pair is projected to state by replacing the spinors of the pair by Similarly, we have Based on these equations, one can calculate , , and perturbatively.
Denoting the perturbative expansion of any quantity as , we have the following orthogonal relations [7]: Based on this, Eq. (3) results in the following matching relations: , (13a) By replacing with , we can obtain similar relations for . These relations enable us to calculate and perturbatively. For simplicity, we concentrate on the S-wave contribution in the rest of the paper.

B. Perturbative calculation in full QCD
We first calculate according to Eq. (8).
The amplitude in full QCD can be decomposed as α s and up to order , one has [9] G =1 + where is the Spence function: In the above results, we have dropped imaginary parts that are irrelevant for our purpose. By inserting Eqs. (14) and (10) into Eq. (8), one gets

M
Note that in Ref. [9], is a free parameter, but to use these expressions for SGF, it needs to be the quarkonium mass.

C. Perturbative calculation of matrix elements in SGF
We now describe our method for calculating . As pointed out in Refs. [7,8], this quantity has been defined to include only small loop momentum regions. In what follows, we provide an explicit definition and choose a UV renormalization scheme. straightforward; the result is Consider the vertex correction in Fig. 1(d) as an example to explain the calculation at the one-loop level. The amplitude reads T S k where is an operator that forces the loop momentum to be in a small region [7,8], which will be defined explicitly by using the method of regions [10,11]. T S Before continuing, we note that, although the full QCD integral in Eq. (20) is well regularized by dimensional regularization, the manipulation to define generates unregularized integrals. Therefore, another regulator is needed to make our manipulation mathematically rigorous. We thus propose a new regularization at the full 1 + η QCD level by multiplying the power of all propagator denominators by ; this can regularize all possible divergences encountered in the derivation using the method of regions, including both ultraviolet and non-ultraviolet divergences 1) . A similar regularization method has been used in the light-cone ordered perturbation theory in Ref. [12] to regularize rapidity divergences. For the integral of interest, we obtain where is introduced to compensate for the change in mass dimension caused by the new regularization, and we assume to ensure that the theory can eventually be regularized by dimensional regularization.
With the new regularization at hand, we decompose the loop momentum into three domains: where the relation " " is understood as the negation of " ," is the complete integration domain, and an implicit cutoff scale exists to rigorously separate the three domains. This division satisfies F Then, one can split any original integral into the three domains:

Fig. 1. Feynman diagrams for the matrix element
An-Ping Chen, Yan-Qing Ma Chin. Phys. C 45, 013118 (2021) η η 1) Non-ultravoilet divergences are necessarily caused by singularities in denominators, which are clearly regularized because we have a power of for all denominators. The power of also makes ultraviolet divergences caused by one or more components of loop momentum going to infinity well regularized. Because effective field theories and factorization theories can be derived from the method of regions, our regularization method is so general that it can regularize all possible divergences in factorization theories and effective field theories. this definition involves a hard cutoff, which makes high order perturbative calculations inconvenient.
To obtain a more convenient definition, we introduce the operators , which expand the integ rand to a convergent power series of small quantities in each domain. We also define . We then have where the property has been used [11]. It is clear that the l.h.s. and the first term on the r.h.s. of this equation are equivalent in the low energy domain and that their difference is an integration in the hard domain, which is infrared safe but may be ultraviolet divergent. The difference can be interpreted as a different choice of renormalization scheme. Therefore, we arrive at our final definition: MS with the UV divergences removed by an renormalization scheme. and are conventionally respectively called the soft region and potential region, while the overlap region compensates for the double counting between and .
For the soft region contribution , we apply the operator to expand the integrand in Eq. (21) for small quantities. This results in the following integrals: where the term comes from the numerator in Eq. (21), and where are non-negative integers. As scaleless and infrared-finite integrals can be set to zero (which is another choice of scheme), we find that only integrals with are relevant: These integrals have pinch poles around that can be regularized by the new regulator . As integrals in can be obtained by expanding the term in the denominator in Eq. (29), it is clear that they cancel exactly with the contribution of the pinch poles around in . This cancellation demonstrates that the overlap region is conceptually important even though it contains only scaleless integrals that are usually set to zero in dimensional regularization. In fact, such overlap contributions have also been introduced as "zero-bin subtraction" in effective-theory treatments [13,14]; this re-η → 0 solves a number of problems in effective theories (NR-QCD and soft-collinear effective theory) [13], including the pinch singularity problem discussed above. By combining the soft region and the overlap region, we only need to consider the contribution from the gluon pole, and we can safely take the limit , which results in .
We emphasize that, if one wants to calculate sep-Theory for quarkonium: from NRQCD factorization to soft gluon factorization Chin. Phys. C 45, 013118 (2021) k 0 k k 0 = 0 arately, the correct order is to first integrate over with fixed . Otherwise, the poles around will be regularized by different regulators between the soft region and overlap region, which necessarily breaks the symmetries of the theory and makes the cancellation between the two regions impossible.
For , the terms in the expansion are proportional to where we can take as all divergences are well regularized by dimensional regularization. Then because and can be expressed as linear combinations of the denominators, if any of , or is nonzero, the integral can be decomposed to either scaleless and infrared-finite integrals or to integrals with pure virtual value. By keeping only the real part, we get Similarly, for the self-energy diagrams, one can derive , Summing these contributions, we obtain the matrix element at NLO before renormalization:

MS
Ultraviolet divergences in the above result can be removed by the renormalization procedure, yielding the renormalized matrix element Similarly, because we find that the real part of is proportional to , we have (38b) We find that the matrix element defined in SGF reproduces all infrared and Coulomb divergences in full QCD and that the obtained hard part is free of divergences. We therefore conclude that the SGF factorization holds at least at the one-loop level.

E. Validity of SGF
The correctness of SGF at the one-loop order can be understood in the following way. The full QCD results in Eq. (18) can also be reproduced by the method of regions (see Appendix A for details), in which is expressed as An-Ping Chen, Yan-Qing Ma Chin. Phys. C 45, 013118 (2021) Here, the first term on the r.h.s has nonzero support only in the hard domain and is thus infrared safe. It is straightforward to check that the low energy part of has been correctly reproduced by the matrix element in SGF, i.e., which leaves the corresponding short-distance hard part defined by the high energy part of , and thus, it is infrared safe. More precisely, we havê

MS MS
where denotes the removal of UV divergences by the subtraction scheme 1) .
The above one-loop argument can be generalized to all orders. By definition, the low energy part (the "small" region of loop momenta) in full QCD can be reproduced by the matrix element in SGF at any order in the expansion, with a proper definition of at the multi-loop level as discussed in Ref. [11]. Therefore, the short-distance hard part is perturbatively infrared-safe, which means that the SGF formula for the decay width of holds for all orders in .

A. NRQCD result
Ignoring operators involving gauge fields, the NR-QCD factorization for the decay amplitude is given by [9,15,16] where and are the two-component heavy quark fields in NRQCD, is defined as , and and are short-distance hard parts that can be perturbatively calculated. The first two orders in the expansion for are given by [9,15,16] given in Eq. (38). In Ref. [9], the S-wave contributions in Eqs. (42) and (44) are further resummed by using the generalized Gremm-Kapustin relation [9,17,18] ⟨0|O An |J/ψ⟩ These relations are obtained by computing the matrix element in the potential-model [17].

B. Equivalence between SGF and NRQCD factorization
The basic justification for the validity of SGF in this process is that the SGF matrix elements reproduce the low energy part of full QCD. As the NRQCD matrix elements also correctly reproduce the low energy part of full QCD, SGF is equivalent to NRQCD factorization. This Theory for quarkonium: from NRQCD factorization to soft gluon factorization Chin. Phys. C 45, 013118 (2021)

MS MS
1) Because overlap regions are scaleless which can be set to zero in dimensional regularization, a simpler way to obtain short-distance hard part is to use . But keep in mind that means to remove IR divergences by subtraction scheme. 2) Note that here the normal derivative is used instead of gauge covariant derivative. They are equivalent for our discussion because operators involving gauge fields are ignored. means that for any process, the two factorization formulas are either both valid or both broken.
Because the amplitude of can be factorized in both SGF and NRQCD factorization, we have (wave contributions are suppressed) W which can generate relations between the SGF matrix elements and NRQCD matrix elements, with denoting contributions from operators with explicit gauge fields. In fact, we can generate more relations by applying the two factorization formulas to any well defined QCD quantity : For example, if we choose , we have . To determine the corresponding , we replace by a pair with invariant mass 1) . Using the nonrelativistic expansion formulas given in [19], we have u v Here, we use nonrelativistic normalization for the spinors and and thus obtain where the extra factor in the first line appears because the NRQCD matrix elements have nonrelativistic normalization, whereas SGF matrix elements have relativistic normalization. As both the SGF matrix elements and NRQCD matrix elements maintain the same low energy physics and are renormalized in the same way, the coefficients should vanish at higher or-α s ders in , i.e., Thus, we get the relation Similarly, by choosing , we can obtain which provides complete relations between the SGF matrix elements and NRQCD matrix elements. Using these relations, we obtain which agrees with Eq. (45). We have thus proved the generalized Gremm-Kapustin relation by using the equivalence between SGF and NRQCD factorization. Based on our proof, it is clear that the terms in the relations are contributions from operators with explicit gauge fields.

⟨0|O
An |J/ψ⟩ n ≥ 1 We see that the equivalence between SGF and NR-QCD factorization gives rise to the generalized Gremm-Kapustin relation. Moreover, by assuming these relations, the introduction of nonperturbative quantities with is unnecessary for exclusive processes. The dominant contributions of these quantities are purely kinematic and can be accounted for by the coefficients of . Thus, the NRQCD factorization formula can be resummed to obtain which is exactly the SGF formula: note Eq. (53), which relates to .

FACTORIZATION
A. Exclusive processes We now show that, in general, SGF can be deduced from NRQCD at the operator level. The equations of motion of heavy quark fields in NRQCD are given by [2] ( with a similar equation for fields. Because we are not interested in gluon fields, we can replace by and by . The equations of motion are usually used to replace operators involving by operators involving in NRQCD. Then, only spatial components appear in the NRQCD operators, which can be decomposed to a relative derivative and total derivative when acting on quarkantiquark bilinear operators. For example, beginning from the bilinear operator , in NRQCD, one can construct more operators by adding the relative derivative to obtain , the total derivative to obtain , or a combination of the derivatives. However, the equations of motion can also be used to replace the relative derivatives and by total derivatives 1) , which results in the following matrix elements for exclusive processes: where and are non-negative integers. As we are working in the rest frame of , by using integration by parts, we find that gives rise to quarkonium mass , and vanishes. Therefore, in the factorization formula, the matrix element with is enough to handle all the contributions in this series of matrix elements, although short-distance hard parts depend on the heavy quark mass as well as the quarkonium mass . This is simply the SGF formula. It is also clear that the SGF resums a series of power corrections originating from the kinematic effects in NRQCD.

B. Inclusive processes
For inclusive quarkonium processes, we can also use equations of motion to decompose the NRQCD matrix elements by Using integration by parts, we can eliminate all the matrix elements except , with the short-distance P 2 P · P X P 2 X P X X hard parts depending on , , and , where is the momentum of unobserved particles . The final result is the SGF formula for inclusive quarkonium processes. It again resums a series of power corrections originating from the kinematic effects in NRQCD. P X P X Because the short-distance hard parts depend on , nonperturbative matrix elements and soft gluon distributions must be also functions of [7]. Note the difference between the soft gluon distributions and shape functions introduced in the endpoint region [20][21][22][23]. The purpose of the shape functions is to resum large logarithms in the endpoint region in the NRQCD factorization framework, which are defined at a fixed power in the relativistic expansion (usually leading power). The SGF with soft gluon distributions is designed to resum a power series of the relativistic expansion and can be applied at both the endpoint region and other regions. If SGF is applied at the endpoint region, the large logarithms can be naturally resummed by the renormalization group equations of the soft gluon distribution. A detailed discussion of this topic will be presented in a separate work [24].

V. SUMMARY
In summary, taking as an example, we have demonstrated that the SGF is equivalent to the NR-QCD factorization for heavy quarkonium production or decay. We have also shown that SGF can be deduced from NRQCD effective field theory at the operator level by using equations of motion. To arrive at this conclusion, we introduced a new regulator and defined the SGF matrix elements rigorously. Based on the equivalence between the two factorizations, we derived explicit relations between the SGF matrix elements and NRQCD matrix elements and proved the generalized Gremm-Kapustin relation.
The results obtained in this paper imply that, for any given process, SGF and NRQCD factorization are either both valid or both violated. Therefore, SGF is valid for all orders in perturbation theories for many processes in which NRQCD factorization has been proved. This provides a solid theoretical foundation for SGF. Compared with NRQCD factorization, SGF effectively resums a subset of relativistic correction terms originating from kinematic effects; this can reduce theoretical uncertainties and may thus provide a better description of experimental data.

ACKNOWLEDGMENTS
We would like to thank K.T. Chao and C. Meng for useful discussions.
Theory for quarkonium: from NRQCD factorization to soft gluon factorization Chin. Phys. C 45, 013118 (2021) 1) During the matching procedure in SGF, We indeed using onshell conditions to fix relative momentum. See discussion after Eq. (7).

LATED BY REGION
In this appendix, we use the method of regions to cal-A cc( 3 S [1] 1 )→e + e − F culate the amplitude at the one-loop order.
According to Eqs. (24) and (25) Defining we have As the first term on the r.h.s can be defined by integrals in the hard domain, it is infrared safe. The remaining contributions and are given by An-Ping Chen, Yan-Qing Ma Chin. Phys. C 45, 013118 (2021)