1.School of Physics and State Key Laboratory of Nuclear Physics and Technology, Peking University, Beijing 100871, China 2.Center for High Energy physics, Peking University, Beijing 100871, China 3.College of Physics and Communication Electronics, Jiangxi Normal University, Nanchang 330022, China 4.Collaborative Innovation Center of Quantum Matter, Beijing 100871, China Received Date:2020-05-19 Accepted Date:2020-09-23 Available Online:2021-01-15 Abstract: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 $J/\psi \to e^+e^-$. 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.
HTML
--> --> --> -->
A.Factorization formula
According to Ref. [8], for the exclusive decay process, we have the following SGF formula at the amplitude level:
where $ n $ denotes intermediate states. In general, $ n $ can contain dynamical soft partons (gluons or light quarks) in addition to a $ Q\bar Q $ pair. However, for simplicity, in this work we only discuss intermediate states without dynamical soft partons, which can nevertheless be treated similarly. The nonperturbative matrix elements $ \overline{R}^{n*}_{{ {\cal Q}}} $ are then defined as
where $ \Psi $ stands for the Dirac field of a heavy quark, $ {\cal{K}}_{n} $ is the projection operator defining the intermediate state $ n $, 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 $ J/\psi \to e^+e^- $, the symmetries of QCD tell us that only $ n = {^{{3}}{S}_{{1}}^{[{1}]}}, {^{{3}}{D}_{{1}}^{[{1}]}} $ are relevant, where the spectroscopic notation, with the superscript "[1]", denotes a color singlet. We thus have
where $ d = 4-2\epsilon $ is the space-time dimension, $ D_\mu $ is the gauge covariant derivative with $ \overline\Psi \overleftrightarrow{D}_\mu \Psi = \overline\Psi (D_\mu \Psi) - $$(D_\mu \overline\Psi)\Psi $, $ P $ is the momentum of $ J/\psi $, $ {\epsilon}_{S_z}^{\mu} $ is a polarization vector with $ P\cdot {\epsilon}_{S_z} = 0 $, $ m $ is the heavy-quark mass, $ M $ is the mass of $ J/\psi $, the color projector is defined as $ {\cal C}^{[1]} = 1/\sqrt{N_c} $, and the spin projection operator $ \mathbb{P}_{\alpha\beta} $ is defined as
The hard parts $ \hat{{\cal{A}}}^{n} $ can be perturbatively calculated according to the matching procedure discussed in [7,8]. To this end, we first replace $ J/\psi $ with an on-shell color-singlet state $ c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}}) $ with momenta①
$p_c = P/2+q, \quad p_{\overline{c}} = P/2-q$
(6)
on both sides of Eq. (3). The on-shell conditions $ p_c^2 = p_{\overline{c}}^2 = m^2 $ result in
$P\cdot q = 0, \quad \quad q^2 = m^2-M^2/4,$
(7)
which fixes $ q_0 $ and $ |{{q}}| $ in the rest frame of $ P $. The remaining degrees of freedom of $ q $ are removed by partial wave expansion, i.e., $ S $-wave expansion in this case. After the replacement, the l.h.s. of Eq. (3) becomes
and $ {\cal{A}}^{\mu }_{c\bar c} $ is the hadronic part of the decay amplitude with the spinors of $ c\bar{c} $ removed. The $ c\bar{c} $ pair is projected to state $ {^{{3}}{S}_{{1}}^{[{1}]}} $ by replacing the spinors of the $ c\bar{c} $ pair by
Based on these equations, one can calculate $ {\cal{A}}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})\to e^+e^-} $, $ \overline{R}^{ {^{{3}}{S}_{{1}}^{[{1}]}} *}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $ , and $ \overline{R}^{ {^{{3}}{D}_{{1}}^{[{1}]}} *}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $ perturbatively. Denoting the perturbative expansion of any quantity $ W $ as $ W = W^{(0)}+\alpha_s W^{(1)}+\alpha_s^2 W^{(2)}+\cdots $, we have the following orthogonal relations [7]:
By replacing $ J/\psi $ with $ c\bar{c}( {^{{3}}{D}_{{1}}^{[{1}]}}) $, we can obtain similar relations for $ \hat{{\cal{A}}}^{ {^{{3}}{D}_{{1}}^{[{1}]}}} $. These relations enable us to calculate $ \hat{{\cal{A}}}^{ {^{{3}}{S}_{{1}}^{[{1}]}}} $ and $ \hat{{\cal{A}}}^{ {^{{3}}{D}_{{1}}^{[{1}]}}} $ perturbatively. For simplicity, we concentrate on the S-wave contribution $ \hat{{\cal{A}}}^{ {^{{3}}{S}_{{1}}^{[{1}]}}} $ in the rest of the paper. 2B.Perturbative calculation in full QCD -->
B.Perturbative calculation in full QCD
We first calculate $ {\cal{A}}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}}) \to e^+e^-} $ according to Eq. (8). The amplitude $ {\cal{A}}^\mu_{c\bar c} $ in full QCD can be decomposed as
Note that in Ref. [9], $ M $is a free parameter, but to use these expressions for SGF, it needs to be the quarkonium mass. 2C.Perturbative calculation of matrix elements in SGF -->
C.Perturbative calculation of matrix elements in SGF
We now describe our method for calculating $ \overline{R}^{ {^{{3}}{S}_{{1}}^{[{1}]}} *}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $. 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. Up to order $ \alpha_s $, the corresponding Feynman diagrams are shown in Fig. 1, where the solid circle represents the operator $ \overline\Psi {\cal{K}}_{ {^{{3}}{S}_{{1}}^{[{1}]}}}\Psi $. The calculation at the tree level is straightforward; the result is Figure1. Feynman diagrams for the $ {^{{3}}{S}_{{1}}^{[{1}]}} $ matrix element
where $ T_S $ is an operator that forces the loop momentum $ k $ to be in a small region [7,8], which will be defined explicitly by using the method of regions [10,11]. Before continuing, we note that, although the full QCD integral in Eq. (20) is well regularized by dimensional regularization, the manipulation to define $ T_S $ generates unregularized integrals. Therefore, another regulator is needed to make our manipulation mathematically rigorous. We thus propose a new regularization at the full QCD level by multiplying the power of all propagator denominators by $ 1+\eta $; this can regularize all possible divergences encountered in the derivation using the method of regions, including both ultraviolet and non-ultraviolet divergences②. 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 $ \nu $ is introduced to compensate for the change in mass dimension caused by the new regularization, and we assume $ \eta\ll \epsilon $ to ensure that the theory can eventually be regularized by dimensional regularization. With the new regularization at hand, we decompose the loop momentum $ k_\mu = (k_0,{{k}}) $ into three domains:
where the relation "$ \lesssim $" is understood as the negation of "$ \gg $," $ D = \mathbb{R}^{d} $ is the complete integration domain, and an implicit cutoff scale exists to rigorously separate the three domains. This division satisfies
A possible definition of $ T_S $ is $ \{k\in D_s\cup D_p\} $. However, this definition involves a hard cutoff, which makes high order perturbative calculations inconvenient. To obtain a more convenient definition, we introduce the operators $ T^{(i)}\; (i \in \{h,s,p \}) $ , which expand the integ rand to a convergent power series of small quantities in each domain. We also define $ T^{(i,j,\cdots)}\equiv T^{(i)}T^{(j,\cdots)} $. We then have
where the property $ T^{(i,j)} = T^{(j,i)} $ 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:
$T_S = T^{(s)}+ T^{(p)}- T^{(s,p)},$
(26)
with the UV divergences removed by an $ \overline{ {{\rm{MS}}}} $ renormalization scheme. $ T^{(s)} $ and $ T^{(p)} $ are conventionally respectively called the soft region and potential region, while the overlap region $ T^{(s,p)} $ compensates for the double counting between $ T^{(s)} $ and $ T^{(p)} $. For the soft region contribution $ \overline{R}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})*,{(s)}}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $, we apply the operator $ T^{(s)} $ to expand the integrand in Eq. (21) for small quantities. This results in the following integrals:
where $ j_i $ 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 $ m_1 = m_2 = j_1 = j_2 = j_4 = j_5 = 0 $ are relevant:
These integrals have pinch poles around $ k_0 = 0 $ that can be regularized by the new regulator $ \eta $. As integrals in $ \overline{R}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})*,{(s,p)}}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $ can be obtained by expanding the $ k_0^2 $ term in the denominator in Eq. (29), it is clear that they cancel exactly with the contribution of the pinch poles around $ k_0 = 0 $ in $ \overline{R}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})*,{(s)}}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $. 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 resolves a number of problems in effective theories (NRQCD 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 $ \eta\to 0 $ , which results in
We emphasize that, if one wants to calculate $ \overline{R}^{ {^{{3}}{S}_{{1}}^{[{1}]}}*,(s)}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $ separately, the correct order is to first integrate over $ k_0 $ with fixed $ {{k}} $. Otherwise, the poles around $ k_0 = 0 $ 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 $ \overline{R}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})*,{(p)}}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $, the terms in the expansion are proportional to
where we can take $ \eta\to 0 $ as all divergences are well regularized by dimensional regularization. Then because $ k_0 $ and $ {{k}}\cdot {{q}} $ can be expressed as linear combinations of the denominators, if any of $ m_1, m_2, j_1, j_2 $, or $ j_3 $ 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
Ultraviolet divergences in the above result can be removed by the $ \overline{ {{\rm{MS}}}} $ renormalization procedure, yielding the renormalized matrix element
Similarly, because we find that the real part of $ \overline{R}^{ {^{{3}}{D}_{{1}}^{[{1}]}} *,(1)}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $ is proportional to $ \overline{R}^{ {^{{3}}{D}_{{1}}^{[{1}]}} *,(0)}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $, we have
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. 2E.Validity of SGF -->
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 $ {\cal{A}}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})\to e^+e^-} $ is expressed as
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 $ {\cal{A}}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})\to e^+e^-} $ 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 $ {\cal{A}}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})\to e^+e^-} $, and thus, it is infrared safe. More precisely, we have
where $ \overline{ {\rm{MS}}} $ denotes the removal of UV divergences by the $ \overline{ {\rm{MS}}} $ subtraction scheme③. 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 $ \alpha_s $ expansion, with a proper definition of $ T_S $ 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 $ J/\psi \to e^+e^- $ holds for all orders in $ \alpha_s $. -->
A.NRQCD result
Ignoring operators involving gauge fields, the NRQCD factorization for the decay amplitude $ J/\psi \to e^+e^- $ is given by [9,15,16]
where $ \psi $ and $ \chi^\dagger $ are the two-component heavy quark fields in NRQCD, $ \overleftrightarrow{\nabla} $ is defined as $\chi^\dagger \overleftrightarrow{\nabla} \psi = \chi^\dagger (\nabla \psi) - $$ (\nabla\chi^\dagger)\psi $, and $ s_n( {^{{3}}{S}_{{1}}^{[{1}]}}) $ and $ b_n( {^{{3}}{D}_{{1}}^{[{1}]}}) $ are short-distance hard parts that can be perturbatively calculated. The first two orders in the $ \alpha_s $ expansion for $ s_n( {^{{3}}{S}_{{1}}^{[{1}]}}) $ are given by [9,15,16]
with $ {\cal{G}}^\prime $, $ {\cal{H}}^\prime $ 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]
These relations are obtained by computing the matrix element $ \langle 0 \vert{\cal{O}}_{An} \vert J/\psi \rangle $ in the potential-model [17]. 2B.Equivalence between SGF and NRQCD factorization -->
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 means that for any process, the two factorization formulas are either both valid or both broken. Because the amplitude of $ J/\psi \to e^+e^- $ can be factorized in both SGF and NRQCD factorization, we have ($ D $-wave contributions are suppressed)
which can generate relations between the SGF matrix elements and NRQCD matrix elements, with $ {\cal{O}}(v^2) $ 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 $ W $:
For example, if we choose $ W = \frac{1}{4}(M+2m) (M-2m)$$\overline{R}^{ {^{{3}}{S}_{{1}}^{[{1}]}} *}_{J/\psi} $, we have $ \hat{W} = \frac{1}{4}(M+2m)(M-2m) $. To determine the corresponding $ w_n $, we replace $ J/\psi $ by a $ c\bar c $ pair with invariant mass $ M $⑤. Using the nonrelativistic expansion formulas given in [19], we have
where the extra factor $ \sqrt{2M} $ 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 $ w_n( {^{{3}}{S}_{{1}}^{[{1}]}}) $ should vanish at higher orders in $ \alpha_s $, i.e.,
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 $ {\cal{O}}(v^2) $ terms in the relations are contributions from operators with explicit gauge fields. We see that the equivalence between SGF and NRQCD factorization gives rise to the generalized Gremm-Kapustin relation. Moreover, by assuming these relations, the introduction of nonperturbative quantities $ \langle 0 \vert{\cal{O}}_{An} \vert J/\psi \rangle $ with $ n\geq 1 $ is unnecessary for exclusive processes. The dominant contributions of these quantities are purely kinematic and can be accounted for by the coefficients of $ \langle 0 \vert{\cal{O}}_{A0} \vert J/\psi \rangle $. Thus, the NRQCD factorization formula can be resummed to obtain
which is exactly the SGF formula: note Eq. (53), which relates $ \langle 0 \vert{\cal{O}}_{A0} \vert J/\psi \rangle $ to $ \overline{R}^{ {^{{3}}{S}_{{1}}^{[{1}]}} *}_{J/\psi} $. -->
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 $ \chi $ fields. Because we are not interested in gluon fields, we can replace $ D_0 $ by $ \nabla_0 $ and $ {{D}} $ by $ {{\nabla}} $. The equations of motion are usually used to replace operators involving $ \nabla_0 $ by operators involving $ {{\nabla}}^2 $ in NRQCD. Then, only spatial components $ {{\nabla}} $ appear in the NRQCD operators, which can be decomposed to a relative derivative and total derivative when acting on quark-antiquark bilinear operators. For example, beginning from the bilinear operator $ \chi^\dagger\psi $, in NRQCD, one can construct more operators by adding the relative derivative to obtain $ \chi^\dagger \overleftrightarrow{{\nabla}}^2 \psi $, the total derivative to obtain $ {{\nabla}}^2(\chi^\dagger \psi) $, or a combination of the derivatives. However, the equations of motion can also be used to replace the relative derivatives $ \overleftrightarrow{{\nabla}}^2 $ and $ \overleftrightarrow{ \nabla}_0 $ by total derivatives⑥, which results in the following matrix elements for exclusive processes:
where $ n_1 $ and $ n_2 $ are non-negative integers. As we are working in the rest frame of $ {\cal{Q}} $, by using integration by parts, we find that $ \nabla_0 $ gives rise to quarkonium mass $ M $, and $ {{\nabla}} $ vanishes. Therefore, in the factorization formula, the matrix element with $ n_1 = n_2 = 0 $ is enough to handle all the contributions in this series of matrix elements, although short-distance hard parts depend on the heavy quark mass$ m $as well as the quarkonium mass $ M $. 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. 2B.Inclusive processes -->
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 $ n_1 = n_2 = 0 $, with the short-distance hard parts depending on $ P^2 $, $ P\cdot P_X $, and $ P_X^2 $, where $ P_X $ is the momentum of unobserved particles $ X $. 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. Because the short-distance hard parts depend on $ P_X $, nonperturbative matrix elements and soft gluon distributions must be also functions of $ P_X $ [7]. Note the difference between the soft gluon distributions and shape functions introduced in the endpoint region [20-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].
APPENDIX A. FULL QCD RESULTS CALCULATED BY REGIONIn this appendix, we use the method of regions to calculate the amplitude $ {\cal{A}}^{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}}) \to e^+e^-} $ at the one-loop order. According to Eqs. (24) and (25), for the original integral $ F $, we have
As the first term on the r.h.s can be defined by integrals in the hard domain, it is infrared safe. We now apply Eq. (A3) to calculate the QCD correction to $ {\cal{A}}^\mu_{c+\bar c} $. We first consider the vertex correction. As showed in Ref. [9], the vertex correction can be expressed in terms of elementary integrals $ I_{111}, I_{011}, I_{-111}, I_{010} $ , and $ I_{110} $, with $ I_{abc} $ defined as
$\tag{A4}I_{abc} \equiv \mu^{2\epsilon} \int\frac{{\rm d}^dk}{(2\pi)^d} \frac{1}{[k^2+{\rm i}0^+]^a[k^2+2p_{ c}\cdot k +{\rm i}0^+]^b[k^2-2p_{\bar c}\cdot k +{\rm i}0^+]^c}.$
According to Eq. (A3), $ I_{abc} $ can be expressed as
The overlap contributions $ I_{abc}^{(h,s)}, I_{abc}^{(h,p)}, I_{abc}^{(h,s,p)}, \text{ and } I_{abc}^{(s,p)} $ and the soft region contribution $ I_{abc}^{(s)} $ are scaleless integrals and can be set to zero if there is no infrared divergence. The calculation of $ I_{abc}^{(s)}-I_{abc}^{(s,p)} $ and $ I_{abc}^{(h,s)}-I_{abc}^{(h,s,p)} $ is similar to $ \overline{R}^{ {^{{3}}{S}_{{1}}^{[{1}]}}*,(s)}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} -\overline{R}^{ {^{{3}}{S}_{{1}}^{[{1}]}}*,(s,p)}_{c\bar{c}( {^{{3}}{S}_{{1}}^{[{1}]}})} $; we have
It should be noted that the ultraviolet and infrared divergences in $ I_{abc}^{(h)} $ are not well distinguished. However, when the overlap contributions $ I_{abc}^{(h,s)}, I_{abc}^{(h,p)}, \text{ and } I_{abc}^{(h,p,s)} $ are subtracted from $ I_{abc}^{(h)} $, the infrared divergent parts in the hard region can be removed. Then, the $ 1/\epsilon $ divergences in $ I_{abc}^{(h)}-I_{abc}^{(h,s)}- $$I_{abc}^{(h,p)} + I_{abc}^{(h,s,p)} $ are converted to $ 1/\epsilon_{UV} $ divergences. Inserting Eqs. (A6), (A8), and (A9) into (A5), we find that the results for the elementary integrals are consistent with those in Refs. [9,25]. We also need to calculate the heavy-quark wave-function renormalization $ Z_Q $, which is given by
This agrees with the result in Ref. [26]. Based on the above results, one can reproduce the full QCD result in Eq. (15) by summing the contributions from the different domains: