The studies on $Z \to \Upsilon(1S)+g+g$ at the next-to-leading-order QCD accuracy

In this paper, we carry out the next-to-leading-order (NLO) studies on $Z \to \Upsilon(1S)+g+g$ via the color-singlet (CS) $b\bar{b}$ state. We find the newly calculated NLO QCD corrections to this process can significantly influence its leading-order (LO) results, and greatly improve the dependence on the renormalization scale. By including the considerable feeddown contributions, the branching ratio $\mathcal{B}_{Z \to \Upsilon(1S)+g+g}$ is predicted to be $(0.56 \sim 0.95)\times 10^{-6}$, which can reach up to $19\% \sim 31\%$ of the LO predictions given by the CS dominant process $Z \to \Upsilon(1S)+b+\bar{b}$. Moreover, $Z \to \Upsilon(1S)+g+g$ also seriously affect the CS predictions on the $\Upsilon(1S)$ energy distributions, especially when $z$ is relatively small. In summary, for the inclusive $\Upsilon(1S)$ productions in $Z$ decay, besides $Z \to \Upsilon(1S)+b+\bar{b}$, the gluon radiation process $Z \to \Upsilon(1S)+g+g$ can provide indispensable contributions as well.

The leading-order (LO) color-singlet (CS) predictions obtained by calculating the CS dominant process Z → Υ(1S) + b +b are only at the 10 −6 order [2][3][4]. Subsequently Li et al. [4] accomplish the next-to-leading-order (NLO) QCD corrections to this bb pair associated channel, pointing out the higher-order terms can give rise to a 24% enhancement to the total decay width.
Reviewing the inclusive J/ψ productions via Z decay [2][3][4][5][6][7][8][9][10][11][12][13][14][15], besides the αα 2 s -order process Z → J/ψ + c +c that serves as the leading role in the CS predictions, the electromagnetic processes Z → ff γ * with γ * → J/ψ (f = l, u, d, s, c, b) and the gluon fragmentation processes Z → f qfq g * with g * → J/ψgg (f q = u, d, s, c, b) can also provide nonnegligible contributions. While, due to the suppression by m 2 c m 2 Z [3], the total decay width of the other αα 2 s -order process Z → J/ψ +g +g is only about two orders of magnitudes smaller than that of Z → J/ψ + c +c. As for the Υ productions, the situations become just the opposite. The relative significances of the electromagnetic and gluon fragmentation processes are much less important than the J/ψ case, since the larger value of m b than m c highly suppress the denominator of the propagators γ * (→ Υ) and g * (→ Υgg). However, for Z → Υ + g + g, the stated above suppresion effect ( m 2 c m 2 Z ) will be largely weaken by replacing m c with m b , subsequently making this process to be indispensable in comparison with Z → Υ + b +b [2].
Moreover, the Υ energy distributions in Z → Υ+g +g and Z → Υ+b+b may be thoroughly different. This can be understood by that the former process is seriously suppressed by the [11,16,17], so the value of z concerning the maximum dΓ dz should be small; however, as a result of the b quark fragmentation, the dominant contributions of Z → Υ + b +b exist in the region of large z. From these points of view, the process Z → Υ + g + g would be crucial for Z decaying to the inclusive Υ, deserving a separate investigation.
Seeing that all the existing studies on Z → Υ(1S) + g + g are only accurate to the first order in α s , to investigate the effects of the higher-order terms, in this paper we will for the first time carry out the NLO QCD corrections to this process. In general, for the The superscripts "q" and "u g " denote the light quarks (u, d, s) and the ghost particles, respectively.
The rest of the paper is organized as follows: In Sec. II, we give a description on the calculation formalism. In Sec. III, the phenomenological results and discussions are presented. Section IV is reserved as a summary.

II. CALCULATION FORMALISM
Based on the nonrelativistic quantum chromodynamics [18], the decay width of Z → Υ(χ bJ ) + g + g 1 can be factorized as whereΓ 1 for Υ, and n = 3 P J (J = 0, 1, 2) for χ bJ . The procedures for dealing with the soft singularities involved in Z → bb[ 3 P [1] J ] + g + g (J = 0, 1, 2) have been described detailedly in our previours paper [19,20], so here we just give a brief presentation on the bb[ 3 S [1] 1 ] related calculations. The NLO short distance coefficients for n = 3 S 1 can be written aŝ whereΓ Γ Virtual is the virtual corrections, consisting of the contributions from the one-loop diagrams (Γ Loop ) and the counter terms (Γ CT ).Γ Real stands for the real corrections, which include the soft terms (Γ S ), hard-collinear terms (Γ HC ), and hard-noncollinear terms (Γ HC ). ForΓ Real , three processes are involved: There are 177 Feynman diagrams in total, including 6 diagrams forΓ Born , 111 diagrams forΓ Virtual (30 counter-terms, 12 self-energy, 30 3-points, 27 4-points, 12 5-points), and 60 diagrams forΓ Real (42 ggg v , 6 ggg av , 6 gqq, and 6 gu gūg ), as representatively shown in Fig.   1. ggg v and ggg av denote the vector and axial-vector parts of Z → bb[ 3 S To isolate the ultraviolet (UV) and infrared (IR) divergences, we adopt the dimensional regularization with D = 4 − 2 . The on-mass-shell (OS) scheme is employed to set the renormalization constants for the heavy quark mass (Z m ), heavy quark filed (Z 2 ), and gluon filed (Z 3 ). The modified minimal-subtraction (M S) scheme is used for the QCD gauge coupling (Z g ). The renormalization constants are [19][20][21] where γ E is the Euler's constant, T F n f is the one-loop coefficient of the β-function, and β 0 = 11 3 C A − 4 3 T F n lf . n f (= 5) and n lf (= n f − 2) are the number of active quark flavors and light quark flavors, respectively. In SU(3), the color factors are given by T F = 1 2 , C F = 4 3 , and C A = 3. The two-cutoff slicing strategy is utilized to subtract the IR divergences in Γ Real [22].
The package MALT@FDC that has been employed to preform the NLO QCD corrections to several heavy quarkonium related processes [19,20,[23][24][25] is used to deal withΓ Virtual ,Γ S , andΓ HC . The FDC package [26] serves as the agent for calculating the hard-noncollinear partΓ HC . Both the cancellation of the −2(−1) -order divergences and the independence on the cutoff parameters δ s,c have been checked carefully, as shown in Figs. 2 and 3.
As a crosscheck, taking the same input parameters, we have reproduced the NLO results of σ e + e − →J/ψ+g+g in Refs. [27,28].
The superscripts "DR" and "F D" denote the direct and feeddown contributions, respectively.  The total decay widths of Z → Υ(1S) + g + g are listed in Table. I. To demonstrate the dependence on the renormalization scale µ r , the results for µ r = 2m b and µ r = m Z are presented simultaneously. From the data in this table, one can observe i) For the direct productions, when µ r = 2m b , the QCD corrections diminish the LO results by about 5%, and cause a 60% enhancement for µ r = m Z . In addition, incorporating these higher-order terms notably weaken the dependence on µ r . As is illustrated in Fig. 4, the line referring to "DR N LO " decreases much more slowly than that for "DR LO " with the increase of µ r .
iii) The dependence on the mass of the b quark is mild, e.g., varying m b by ± 0.2 GeV from the central value of 4.9 GeV just results in a 10% variation of the total decay width.
Summing up the direct and feeddown contributions, we finally obtain where the theoretical uncertainty is induced by the choices of the values of µ r (2m b ∼ m Z ) and m b (4.9 ∼ 5.1 GeV).
In Fig. 5, the Υ(1S) energy distributions in Z → Υ(1S) + g + g are presented with z defined as m Z . For the direct productions, when µ r is equal to 2m b , the QCD corrections are positive for z < 0.52, and negative in the remaining scope of z; however, for µ r = m Z , these corrections keep always positive for all the available z values. With µ r being relatively small, the higher-order terms can significantly enhance the differential decay width. Taking where "N LO" represents the sum of the contribution of the LO terms and that of the QCD corrections. Including the feedown contributions can further increase dΓ dz by about 20% ∼ 30% for most values of z. Note that, phenomenological and theoretical arguments suggest that the NRQCD factorization holds only when the quarkonium is produced at relatively large momentum, thus in Fig Comparing to the results with m b = 4.9 GeV in Tab. I, one can obtain where the uncertainty arises from the variation of µ r in [2m b , m Z ]. This ratio suggests that, for the total decay width, the contributions via Z → Υ(1S) + g + g is comparable with the magnitude of the QCD corrections to Z → Υ(1S) + b +b. In Fig. 6, the comparison of the Υ(1S) energy distributions in Z → Υ(1S) + g + g and Z → Υ(1S) + b +b are presented, where "Υgg" denotes the contributions via the former process up to the NLO accuracy in α s , including the feeddown effects; "Υbb" stands for the direct contributions of the bb pair associated process at the LO QCD accuracy. One can find, for relatively small z, adding the Υgg contributions can increase the Υbb predictions to a surprisingly large extent. Such a remarkable enhancement on dΓ dz is almost the same in size as the QCD corrections to Z → Υ(1S) + b +b. All these points obviously reveal the phenomenological importance of Z → Υ(1S) + g + g for Z decaying to the inclusive Υ(1S).
Note that, within the NRQCD framework, for the Υ production in association with final gluon(s) via Z decay, in addition to the CS process Z → bb[ 3 S [1] 1 ]+g+g that we focus on in our present paper, the color-octet (CO) channels Z → bb[ 1 S [8] 0 , 3 S 1 ] + g + g, together with the potential enhancement via the NLO corrections to these CO channels 2 , may partly compensate for the suppression caused by LDMEs, subsequently making the CO contributions nonnegligible.
Of course, whether this is indeed the case depends on the future rigorous NLO calculations for these CO production channels.

IV. SUMMARY
In this manuscript, we for the first time perform the complete NLO studies on the process Z → Υ(1S) + g + g via the CS bb states. We find the impacts of the QCD corrections on the LO results are significant, including both the total and differential decay widths. In addition, these higher-order terms markedly weaken the dependence of the theoretical predictions on µ r . By incorporating the substantial feeddown effects via Υ(2, 3S) and χ bJ (1, 2, 3P ), B Z→Υ(1S)+g+g is scattered in the range (0.56 ∼ 0.95) × 10 −6 . This value is about 19% ∼ 31% of the LO predictions given by Z → Υ(1S) + b +b, which is responsible for the main 2 Such as the fragmentation process of Z → qqg * with g * → bb[ 3 S contributions in the CS predictions. Moreover, for small E Υ(1S) , Z → Υ(1S) + g + g has vital influence on the Υ(1S) energy distributions. In view of these points, for Z decaying to the inclusive Υ(1S), in addition to Z → Υ(1S) + b +b, the process Z → Υ(1S) + g + g is also phenomenologically important.