NRQCD matrix elements for S-wave bottomonia and Gamma[eta_b(nS) ->gamma gamma] with relativistic corrections

We determine the leading-order nonrelativistic quantum chromodynamics (NRQCD) matrix element_Upsilon and the ratio_Upsilon, for Upsilon=Upsilon(nS) with n=1, 2, and 3 by comparing the measured values for Gamma[Upsilon ->e^+ e^-] with the NRQCD factorization formula in which relativistic corrections are resummed to all orders in the heavy-quark velocity v. The values for_Upsilon, which is the ratio of order-v^2 matrix element to_Upsilon, are new. They can be used for NRQCD predictions involving Upsilon(nS) and eta_b(nS) with relativistic corrections. As an application, we predict the two-photon decay rates for the spin-singlet states: Gamma[eta_b(1S) ->gamma gamma] = 0.512^{+0.096}_{-0.094} keV, Gamma[eta_b(2S) ->gamma gamma] = 0.235^{+0.043}_{-0.043} keV, and Gamma[eta_b(3S) ->gamma gamma] = 0.170^{+0.031}_{-0.031} keV.


I. INTRODUCTION
The pseudoscalar bottomonium η b (1S), which is the spin-singlet S-wave ground state, was first observed in the photon energy spectrum of the radiative Υ(3S) decay [1] and confirmed in the radiative Υ(2S) decay [2] by the BABAR Collaboration. The state was also confirmed by the CLEO Collaboration again in Υ(3S) → γη b (1S) [3]. So far, only the mass for the η b (1S) is known as m η b (1S) = 9390.9 ± 2.8 MeV [4], and any of its exclusive decay modes has not been observed, yet. Among its various decay modes, recent theoretical studies have been concentrated on relatively clean channels like η b → J/ψ J/ψ [5][6][7][8], η b → J/ψ γ [9], and others [10,11]. 1 On the other hand, the most elementary exclusive decay channel is η b → γγ, although it has a large background. With the decay mode Υ → e + e − of the spin-triplet partner, η b → γγ must be well described by the nonrelativistic quantum chromodynamics (NRQCD) factorization formulas for the electromagnetic decay of heavy quarkonia [12]. If one makes use of the heavy-quark spin symmetry, then one can make a rough estimate of the decay rate, whose branching fraction is ∼ 10 −5 , which is relatively greater than other channels listed above.
In the mean time, there has been a significant progress in the NRQCD calculations for Swave charmonium production and decay, in which relativistic corrections of all orders in the heavy-quark velocity v are resummed [28]. Precise determination of the wavefunction at the origin for the J/ψ was made based on this method [29][30][31][32]. This method has been applied to reconcile the large discrepancy between the theoretical prediction and the experimental 1 Throughout this Letter we use collective notations η b and Υ that indicate η b (nS) and Υ(nS), respectively, for n = 1, 2, and 3 unless a specific state is specified. results for the cross section σ[e + e − → J/ψ + η c ] at the B factories [33][34][35]. Therefore, it is worthwhile to improve the NRQCD prediction for Γ[η b → γγ] by taking into account the relativistic corrections to all orders in v. In order to carry out such an analysis, one needs to know the values for the color-singlet NRQCD matrix elements and those involving relativistic corrections. The NRQCD matrix element for the S-wave bottomonia can be determined by making use of the measured values for Γ[Υ → e + e − ] up to corrections of spin-symmetry breaking effects. Unfortunately, available order-v 2 NRQCD matrix element that has been fixed from lattice QCD simulations [36,37] suffers from large uncertainties originated from slow convergence of the cut-off regularization method.
In this Letter, we first determine the NRQCD matrix elements for the S-wave bottomonium states that are required to compute the relativistic corrections with considerably less uncertainties than available values, extending the method in Refs. [28,30]. As an application, we compute Γ[η b → γγ], in which corrections of order the strong coupling α s and relativistic corrections of all orders in α 0 s q 2n are included. Here, q is half the relative momentum of b andb in the bottomonium rest frame. The remainder of this Letter is organized as follows: In Section II, we estimate the NRQCD matrix elements for the S-wave bottomonium states by making use of the resummed NRQCD factorization formula against empirical data for the spin-triplet states. Our prediction for Γ[η b → γγ] is presented in Section III with the comparison with available predictions and we summarize in Section IV.

II. NRQCD MATRIX ELEMENTS FOR THE Υ(nS)
In this section, we briefly review the method to determine the NRQCD matrix elements for Υ at leading and subleading order in q 2 based on the strategy for the charmonium counterpart in Ref. [30]. The results are compared with those of lattice QCD calculations.
The NRQCD factorization formula for the electromagnetic decay of the S-wave quarkonium H is a linear combination of nonperturbative NRQCD matrix elements O n H that are classified in powers of v, where O n is the NRQCD operator. The factorization is achieved at the amplitude level and the ratio q 2n H of the order-q 2n matrix element to the leading one are all, in general, independent. In addition, the ratios q 2n H have power-ultraviolet divergences that must be regulated and, therefore, the values can even be negative under subtraction. In lattice QCD calculations, this subtraction is made by making use of the hard-cut-off regularization whose convergence is slow, resulting in large uncertainties [28].
However, in an electromagnetic decay, in which the color-singlet contributions dominate, one can calculate the quarkonium wavefunction of the leading heavy-quark-antiquark (QQ) Fock state up to corrections of relative order v 2 if one knows the static, spin-independent QQ potential exactly. The authors of Refs. [28,30] have constructed the generalized ver- H ] n , of the Gremm-Kapustin relation [38] to resum a class of relativistic corrections. The method has been devised to be consistent with dimensional regularization of these power-ultraviolet divergent matrix elements.
The resultant formula for the decay rate of Υ → e + e − , in which relativistic corrections of all color-singlet QQ operator matrix elements are resummed, is given by [29,30] where m Υ is the Υ mass, O 1 Υ is the color-singlet NRQCD matrix element for the electromagnetic decay of the Υ at leading order in v, and v 2 The order-α 2 s corrections to Γ[Υ → e + e − ] (Refs. [39,40]) contain a strong dependence on the NRQCD factorization scale. If one were to include those corrections in Eq. (1) and use it to determine O 1 Υ , then O 1 Υ would also contain a strong dependence on the NRQCD factorization scale, which would cancel in other quarkonium decay and production processes only if the short-distance coefficients were calculated through relative order α 2 s . Generally, short-distance coefficients for quarkonium processes have not been calculated beyond relative order α s . For this reason, we omit the order-α 2 s corrections to the leptonic width in Eq. (1). Nevertheless, if one includes the order-α 2 s corrections and take the factorization scale to be m b , the resultant NRQCD matrix elements are increased by about a factor of 40%.
We briefly discuss the method employed in this Letter to compute O 1 Υ and q 2 Υ . We follow the method given in Ref. [30] and make use of the Cornell potential model [41].
By using the Schrödinger equation we can express O 1 Υ and q 2 Υ as functions of the parameters of the Cornell potential model, which are the mass parameter in the Schrödinger equation, the string tension, and the Coulomb strength of the Cornell potential. The mass parameter can be expressed in terms of the 1S-2S mass splitting [30], which we compute from the masses of Υ(1S) and Υ(2S). The value of the string tension, which is universal, is taken from lattice measurements as 0.1682 ± 0.0053 GeV 2 [42][43][44][45]. Finally, the Coulomb strength parameter is determined by constraining the rate (1) to be consistent with the experimental value [4] and solving the resulting nonlinear equation numerically. Because of this, the value of the Coulomb strength parameter is chosen differently for each quarkonium. between this analysis and that for the S-wave charmonium in Ref. [30] is that there are no measured data for Γ[η b → γγ]. Therefore, we use the spin-triplet data only.
By carrying out these calculations, the Coulomb strength parameter is fixed as 9.955 for Υ(1S), 10.960 for Υ(2S), and 11.127 for Υ(3S), respectively. From these we obtain our results for O 1 Υ and q 2 Υ , which are tabulated in Table I Table I reflect the uncertainties arising from m b , Γ[Υ → e + e − ], string tension, α s , and the ignorance of the spin-dependent interactions of the potential in the Schrödinger equation [30], all of which are added in quadrature. The values for the leading-order matrix elements O 1 Υ in Table I have been used to predict the inclusive charm production in Υ decays [46] and the heavy quarkonium production associated with a photon in e + e − annihilation [47]. The values for q 2 Υ in Table I are new. The value for Υ(1S) has errors significantly less than those of the available lattice QCD calculations [36,37]. One can reduce theoretical uncertainties   by considering the dependence on m b that also appears in the short-distance coefficients of factorization formulas. Therefore, we present the sources of errors in Table I. Unlike the S-wave charmonium case in Ref. [30], the uncertainties of O 1 Υ and q 2 Υ due to the errors of the heavy-quark mass are insignificant.

Sources of errors
Our results are now compared with those for the ground-state S-wave bottomonium obtained from a lattice QCD simulation. The results from the quenched approximation are given in Ref. [36]. We quote the updated results of the unquenched analysis in Ref. [37]: The leading-order NRQCD matrix element is O 1 1S = 4.10(1)(9)(41) GeV 3 and the ratio 1S ranges from about −5 GeV 2 to about 2 GeV 2 [37]. Here, the subscript 1S indicates the average of η b (1S) and Υ(1S). Our central value for the O 1 Υ(1S) is about 25 % smaller than that of Ref. [37], which is greater than the quenched case [36] by about a factor of 2. In the case of q 2 Υ(1S) , our result is consistent with that in Ref. [37] but ours has uncertainties significantly smaller than that of the lattice result.

III. TWO-PHOTON WIDTHS FOR THE η b
In this section, we predict Γ[η b → γγ] by making use of the NRQCD matrix elements determined in Section II. In fact, the NRQCD matrix element O 1 η b that appears in the factorization formula for Γ[η b → γγ] might be different from O 1 Υ by a relative order v 2 , which breaks the approximate heavy-quark spin symmetry [12]. We recall that the effect of spin-symmetry breaking in the low-lying S-wave charmonia J/ψ and η c is not significant [30].
Therefore, the errors in the approximation O 1 η b (nS) ≈ O 1 Υ(nS) may be insignificant based on the fact that v 2 Υ ≪ v 2 J/ψ = 0.22 [30]. As in the leptonic decay of the Υ, we include the relativistic corrections to Γ[η b → γγ] to  all orders in v. The resultant factorization formula is given by [30,31] where the relativistic corrections are incorporated into the function g( The input parameters for the numerical calculations are chosen in a similar way in Ref. [30] for Γ[η c → γγ]. (1)], we omit the order-α 2 s corrections to the two-photon width of the η b , whose result is available in Ref. [48].
The resultant predictions for Γ[η b → γγ] are tabulated in Table II Table I. We also include the errors of using O 1 Υ and q 2 Υ , which are of relative order v 2 set to be 0.1. From the order-α 2 s corrections to the electromagnetic widths of the Υ and η b [39,40,48], we find that the order- Refs. [15][16][17][18][20][21][22][23] agree with our prediction within errors, while some models [13,14,19], which does not use the heavy-quark spin symmetry, apparently underestimate the rate in comparison with ours. Our result Γ[η b (1S) → γγ] = 0.512 +0.096 −0.094 agrees with the most recent potential-NRQCD prediction 0.54 ± 0.15 keV in Ref. [27] in which the leading relativistic corrections are included, while it is smaller than another potential-NRQCD prediction in Ref. [26]. Note that we have borrowed O 1 Υ and q 2 Υ for O 1 η b and q 2 η b in Eq. (2) after taking into account the errors of spin-symmetry breaking effect as v 2 ∼ 0.1 because → γγ] is consistent with a recent potential-NRQCD prediction in Ref. [27], in which the leading relativistic corrections are included.
These are not sufficient to observe the mode η b (1S) → γγ. However, we expect that this channel can be observed at the superKEKB or superB factory if more data are accumulated. The CERN Large Hadron Collider is expected to produce about 5 × 10 9 η b 's with the integrated luminosity ∼ 300 fb −1 [6], with which one can probe about 45000 events of η b (1S) → γγ. We anticipate such a stage against which our predictions can be tested.