Mass shifts of $^3P_J$ heavy quarkonia due to the effect of two-gluon annihilation

In this work, we calculate the nonrelativistic asymptotic behavior of the amplitudes of $q\overline{q}\rightarrow 2g \rightarrow q\overline{q}$ in the leading order of $\alpha_s$ (LO-$\alpha_s$) with $q\overline{q}$ in the $^3P_J$ channels. In the practical calculation we take the momenta of quarks and antiquarks on-shell and expand the amplitudes on the three-momentum of the quarks and antiquarks to order 6 and get three nonzero terms. The imaginary parts of the first term and the second term are the old. The real parts of the results have IR divergence. When applying the results to the heavy quarkonia, the corresponding amplitude of $q\overline{q}\rightarrow1g\rightarrow q\overline{q}$ with $q\overline{q}$ in the color octet $^3S_1$ channel is considered to absorb the IR divergence in a unitary way in the leading order of $v$ (LO-$v$). The finial results can be used to estimate the mass shifts of the $^3P_J$ heavy quarkonia due to the effect of two-gluon annihilation. The numerical estimation shows that the contributions to the mass shifts of $\chi_{c0,c1,c2}$ are about $1.23\sim1.58$ MeV, $1.57\sim1.86$ MeV and $5.92\sim5.45$ MeV when taking $\alpha_s\approx 0.25\sim0.35$.


I. INTRODUCTION
The energy spectrum of an elemental system is a basic question after the breakthrough of quantum mechanism. Currently, the energy spectrum of hadrons are still an unsolved problem in QCD due to the complex nonperturbative property. Many phenomenological models have been used to studies the energy spectrum of hadrons in the quark level such as the quark model [1], QCD sum rules [2], Dyson-Schwinger equation and Bethe-Salpeter equation [3], etc. In these calculations the annihilation effect whose imaginary and real parts correspond to the decay width and the mass shift is usually neglected. For heavy quarkonia, their inclusive decays can be well described by the effective theory nonrelativistic QCD (NRQCD) [4]. In NRQCD, the imaginary part of the coefficients of four fermions interactions are matched from the imaginary parts of the on-shell scattering amplitudes qq → 1g or 2g or 3g → qq or the decay widths of qq → 1g or 2g or 3g in perturbative QCD order by order.
In previous paper [5], we calculated the real parts of these coefficients in the leading order of α s (LO-α s ) in the 1 S 0 channel with the momenta of quarks and antiquarks off-shell and find the results are gauge invariant, while the similar calculation can not be directly extended to the 3 P J channels due to the gauge invariance. In this paper, we follow the idea of NRQCD and calculate the amplitudes of qq → 2g → qq in the 3 P J channels and qq → 1g → qq in the 3 S 1 channel with qq in color single and color octet states, respectively.
We organize the paper as follows. In Sec. II we give an introduction on the basic formula, in Sec. III we describe our calculation and present the analytic results for the coefficients to order 6 after the nonrelativistic expansion, in Sev. IV we estimate the effects to the mass shifts numerically and discuss the interesting properties of the results. which can be expressed as (1) where H 1 and H 8 are some nonperturbative matrix elements, M(qq( 3 P J ) 1 → qq( 3 P J ) 1 ) and M(qq( 3 S 1 ) 8 → qq( 3 S 1 ) 8 ) are the amplitudes with the momenta of the quarks and antiquarks on shell and the indexes 1 and 8 refer to the color signal and color octet states, respectively.
The amplitudes at quark level can be calculated perturbatively. In the perturbation theory, the corresponding Feynman diagrams for the amplitudes of qq( 3 P J ) 1 → 2g → qq( 3 P J ) 1 are shown in Fig. 1, and the the transition qq( 3 S 1 ) 8 → g → qq( 3 S 1 ) 8 are shown in Fig. 2.
In the center mass frame, we choose the momenta as follows.
with p 2 1 = p 2 2 = p 2 3 = p 2 4 = m 2 and m the mass of heavy quark. We can define P (E, 0, 0, 0), p i (0, p i ), and p f (0, p f ). For the heavy quark and antiquark pairs we can take |p i,f |/m The Feynman diagrams for qq → 2g → qq in the 3 P J channels in the leading order of α s . as small variables, and then can expand the expressions on these small variables. In the onshell case, we have the relations |p i |= |p f | p and E = m 2 + p 2 . This relation means we can not distinguish |p i | from |p f | and a non-uniqueness may happen when applying the final expressions to the bound states. There is no such non-uniqueness in the direct calculation of the imaginary parts since the momenta p i and p f appear independently after cutting the gluon lines. Fortunately, we can see the first two non-zero orders can be gotten uniquely due to the symmetry which will be discussed in the following.
To project the quark and antiquark pairs to the 3 P J state and the 3 S 1 state, we use the project matrix in the on-shell case [6,7] and have the following.
In our calculation, we only calculate the amplitudes in the perturbative QCD and do not go to match them with the corresponding amplitudes in NRQCD in the on-shell region, so we directly try to include the structure of the heavy quarkonia in our calculation. To include the information of the H( 3 P J ) in the perturbative QCD, we assume the following structure for the H( 3 P J ). (7) where the color factors 1/ √ N c and 1/ (N 2 c − 1)/2 are used to normalized the color parts to 1, φ 0,1 (|p|) refer to the wave functions of the H( 3 P J ) in the momentum space in the color single and color octet states, respectively. The relations between the wave functions φ 0,1 with the wave functions in the coordinate space are defined as Using the structure of the H( 3 P J ) and the above project matrices, the amplitudes can be expressed as follows. where and with D = 4 − 2ǫ and the color factor c the hard kernel T i and In the real bound states the values of |p i | and |p f | are independent which is different from the on-shell case, so we label |p i | and |p f | independently in the above original expressions.
To calculate G (a,b,c) (s i , s f ), we use the package Feyncalc [8] to do the trace of Dirac matrices in D-dimension and then expand the expressions on the variable p to a special order.
After the expansion, we use the tensor decomposition to re-expressed the loop integrations and finally use the package FIESTA [9] to do sector decomposition and then use Mathematica to do the analysic integration.
After the the loop integrations, the form of G (a+b,c) (s i , s f ) can be expressed as follows.
After getting the coefficients C (a+b,c) in , usually the properties of the integrations of angle and the sums of the spins are independently used to simplify the expressions as in Ref. [6].
In our calculation, for simplification we directly calculate the sums of the spins and the integrations of angles together. We define P (J, X, n) where X are some functions dependent onp i ,p f , ǫ(s i ) and ǫ * (s f ) withp i,f p i,f /|p i,f |, n = 0, 1, 2 and 3, J = 1, 2 and 3, P (X, n) and Q(J, X, n) are not dependent on J z whose manifest form are directly listed in Appendix A. Using the expressions of P (J, X, n) and Q(X, n), M( 3 P J ) and M( 3 S 1 ) can be calculated easily.

III. THE ANALYTIC RESULTS FOR THE ASYMPTOTIC BEHAVIOR
In the practical calculation, we expand the expressions on p to order 6. Since the calculation is taken with the momenta on-shell, the gauge invariance is manifest. The final result can be expressed as where the subindexes p means to expand the expressions on p. For the real bound states, the terms p 2 and p 4 only receive contributions from the terms |p i ||p f | and 1 2 (|p i ||p f | 3 +|p i | 3 |p f |), respectively, since only they are nonzero. The term p 6 receives the contributions both from the terms 1 2 (|p i ||p f | 5 +|p i | 5 |p f |) and |p i | 3 |p f | 3 , which results in that one can not distinguish them in a unitary from.
The imaginary parts of c J,i and d i are expressed as follows. with In the following discussion, we define Re[c f in J,i ] as the finite part of Re[c J,i ] with the C IR related parts being subtracted.
Using the above expressions and the quasi potential method, one has the following relation for the corresponding effective potential in the LO-α s .
In the LO-v, the corresponding decay widths of H( 3 P J ) to the light hadrons (l.h) from the above diagrams which is labeled as Γ( 3 P J → l.h) are expresses as and the corresponding mass shifts labeled as ∆M( 3 P J ) are expressed as where we have used the relation Eq.(27) shows when one goes to discuss the mass shifts of the H( 3 P J ) states due to the two-gluon annihilation effects in the LO-α s , the color octet contribution should also be considered. The IR divergence in M( 3 P J ) can be absorbed in a unitary way by M( 3 S 1 ) in the LO-v. The absorbed form is unique and same with the case of the one-loop radiative corrections to the decay width Γ( 3 P J → l.h). This is nature if one goes to match the above results with the corresponding NRQCD coefficients.
In the literature, the decay widths Γ( 3 P J → l.h) in the LO-v and the NLO-α s can be expressed as follows [10].
where n f is the number of light quarks, n f = 3 for charmonium and n f = 4 for bottomonium states, C 0 and C 2 are expressed as and H 1 is related to the derivative of the wave function through the relation: After the following replacement similarly with that in [11]: we can see that the IR divergences in Eq.(27) and Eq. (29) are absorbed in the same form.
In the literature, the relation between H 8 and the wave function of the high Fock state is not given. From the comparison between Eq. (27) and Eq. (29) we can get the following relation: Combing Eq. (27) with Eq. (29), finally one can get with Eq. (34) can be used to estimate the mass shifts of H( 3 P J ) in the LO-α s and the LO-v due to the two-gluon annihilation effect.

IV. NUMERICAL RESULT AND CONCLUSION
The main results of our calculation are the expressions of the coefficients c J,i and the mass shifts ∆M( 3 P J ). In the LO-v and the LO-α s , if one assumes that the contribution from the H 8 related term is small, then the ratios between the mass shifts and the decay widths can be expressed as where the similar result for the 1 S 0 state is also presented. We can see that the ratio for the 3 P 0 state is much larger than the ratios for the 3 P 2 and 1 S 0 states and the ratios are positive for the 3 P J states and negative for the 1 S 0 state. taken from Ref. [13]. The similar estimation can be applied to the bottomonium. Comparing these numerical results with the corresponding results of η c [5], we can find that the mass shifts of χ cJ are very different. These properties mean the corrections to different states can not be subtracted or hidden in a unified way. Combing the numerical results, one can get ∆M( 3 P 0 ) − ∆M( 1 S 0 ) ≈ 9.0 ∼ 8.6 MeV with α s ≈ 0.25 ∼ 0.35, correspondingly. This numerical result suggests that the annihilation effects should be considered seriously when try to understand the spectrum of heavy quarkonia precisely, especially when some decay channels with large decay widths are opened.  Another interesting property is that although the decay width of the 3 P 1 states to twogluon intermediated state is zero, the corresponding mass shift is nonzero.
In summary, the real part of the nonrelativistic asymptotic behavior of the amplitudes of qq → 2g → qq in the 3 P J channels is discussed in the LO-α s . By expanding the expressions on the three-momentum of quarks and antiquarks, the expressions are calculated to order 6. The imaginary part of the first 2 terms of our results are the same with those given in the references. The real part of our results can be used to estimate the mass shifts of the 3 P J heavy quarkonia due to the two-gluon annihilation effect. In the LO-α s and the LO-v, we get the following properties: (1)