A QCD Sum-Rules Analysis of Vector ($1^{--}$) Heavy Quarkonium Meson-Hybrid Mixing

We use QCD Laplace sum-rules to study meson-hybrid mixing in vector ($1^{--}$) heavy quarkonium. We compute the QCD cross-correlator between a heavy meson current and a heavy hybrid current within the operator product expansion. In addition to leading-order perturbation theory, we include four- and six-dimensional gluon condensate contributions as well as a six-dimensional quark condensate contribution. We construct several single and multi-resonance models that take known hadron masses as inputs. We investigate which resonances couple to both currents and so exhibit meson-hybrid mixing. Compared to single resonance models that include only the ground state, we find that models that also include excited states lead to significantly improved agreement between QCD and experiment. In the charmonium sector, we find that meson-hybrid mixing is consistent with a two-resonance model consisting of the $J/\psi$ and a 4.3~GeV resonance. In the bottomonium sector, we find evidence for meson-hybrid mixing in the $\Upsilon(1S)$, $\Upsilon(2S)$, $\Upsilon(3S)$, and $\Upsilon(4S)$.


Introduction
Hybrids are hadrons containing explicit gluon degrees of freedom in addition to a constituent quark and antiquark. They are colour singlets and so should be allowed within QCD. However, they have yet to be conclusively identified in experiment (see, e.g., ref. [1] for a comprehensive review).
Hybrids can be broadly classified as having quantum numbers (i.e., J P C ) that are exotic or nonexotic. Exotic quantum numbers (e.g., 0 −− , 0 +− , 1 −+ , 2 +− ) are those not accessible to conventional quark-antiquark (qq) mesons; the rest of the quantum numbers are non-exotic and are accessible to both qq-mesons and hybrids. Looking for resonances with an exotic J P C is a promising hybrid search strategy being used, for example, at GlueX. Furthermore, hybrids with exotic J P C would be unable to quantum mechanically mix with qq-mesons (as no conventional meson could have the J P C in question), and so could perhaps appear as pure, unmixed states. In contrast, hybrids with non-exotic J P C are expected to mix with qq-mesons resulting in hadrons that would be superpositions of both conventional meson and hybrid.
In this article, we consider meson-hybrid mixing in (non-exotic) vector (1 −− ) charmonium (cc) and bottomonium (bb). The heavy quarkonium sectors have received considerable attention lately due primarily to the discovery of the XYZ resonances (see [2,3] for reviews and [4] for some recent developments). These XYZ resonances are a collection of hadrons many of whose properties (e.g., masses, widths, and decay rates) do not agree with quark model predictions [5]. Unsurprisingly, the XYZ resonances have generated a lot of discussion concerning outside-the-quark-model hadrons such as hybrids. We focus on 1 −− rather than some other J P C because more is known about the spectra of 1 −− heavy quarkonium than is known about the spectra for the other quantum numbers [6].
We investigate meson-hybrid mixing with QCD Laplace sum-rules (LSRs) [7,8,9,10]. Using the operator product expansion (OPE) [11], we compute the cross-correlator between a qq-meson current and a hybrid current (see (3) and (4) respectively below). In the cross-correlator calculation, we include leading-order (LO) QCD contributions from perturbation theory and non-perturbative corrections due to the four-dimensional (4d) and 6d gluon condensates as well as the 6d quark condensate. We then analyze several single and multi-resonance models of the hadron mass spectra that take known resonance masses as inputs. We determine which resonances couple to both currents and so can be considered mixed. The QCD sum-rules methodology has been applied to hadron mixing problems in a number of systems including pseudoscalar meson-glueball mixing [12], scalar meson-glueball mixing [13], 1 ++ charmonium hybrid-DD * molecule mixing [14], and open-flavour heavy-light meson-hybrid mixing [15].
We find that multi-resonance models that include excited states in addition to the ground state lead to significantly improved agreement between QCD and experiment when compared to single resonance models that include only the ground state. In addition, we show explicitly that the higher mass excited states make numerically significant contributions to the LSRs despite the tendency of LSRs to suppress such resonances. Finally, we find that meson-hybrid mixing in the charmonium sector is described well by a two-resonance model consisting of the J/ψ and a 4.3 GeV state such as the X(4260). In the bottomonium sector, we find evidence for meson-hybrid mixing in all of the Υ(1S), Υ(2S), Υ(3S), and Υ(4S).

The Correlator
We consider the following cross-correlator between quarkonium meson current j (m) and quarkonium hybrid current [16] where is the dual gluon field strength tensor and νρωζ is the totally antisymmetric Levi-Civita symbol. The function Π in (2) probes 1 −− states. We calculate the correlator (1) within the OPE in which perturbation theory is supplemented by nonperturbative terms, each of which is the product of a perturbatively computed Wilson coefficient and a nonzero vacuum expectation value, i.e., a QCD condensate. In addition to perturbation theory, we include OPE terms proportional to the 4d and 6d gluon condensates and the 6d quark condensate defined respectively by αG 2 = α s : G a ωφ G a ωφ : J 2 = Tr : J ν J ν : where  (13).
In (9), the sum on the right-hand side is over quark flavours. We use the vacuum saturation hypothesis [7] to express J 2 in terms of the 3d quark condensate qq = : q α i q α i : resulting in where κ quantifies deviations from exact vacuum saturation. Throughout, we set κ = 2 (see, e.g., ref. [10] and references cited therein). The diagrams that contribute to (1) at LO in the coupling g s are shown in Figure 1 where we have suppressed a second set of similar diagrams in which the quark line runs clockwise. Wilson coefficients are computed using the fixed-point gauge method (see [17,18], for example), and divergent integrals are handled using dimensional regularization in D = 4 + 2 dimensions at MS renormalization scale µ. As in [19], we use the following convention for a dimensionally regularized γ 5 : We employ TARCER [20], a Mathematica package that implements the recurrence algorithm of [21,22], to express dimensionally regularized integrals in terms of a small set of master integrals. An exact calculation of each needed master integral is either found in [23,24] or is a well-known one-loop result. We denote the OPE computation of Π from (2) as Π (OPE) which we then decompose as where the superscripts in (13) correspond to the labels of the diagrams in Figure 1. For Π (I) , we find an exact -dependent result where m is a heavy quark mass (i.e., m c or m b ), Γ is the gamma function, and p F q (· · · ; · · · ; z) are generalized hypergeometric functions (see [25], for example). Expanding (14) in and dropping terms polynomial in z as they will not contribute to the LSR, we find For the sake of brevity, we do not include an explicit expression for the derivative term on the right-hand side of (16). (Note that (16) is ultimately superseded by (27), and we provide a complete expression for the latter.) Expanding the remaining terms on the right hand side of (13) in , we find Perturbation theory (16) contains a nonlocal divergence. Following [14,15], this divergence is eliminated via operator mixing under renormalization. The meson current (3) is renormalization-group (RG) invariant and so we only need consider the operator mixing of the hybrid current (4) which induces operator mixing with (3) and with where C 1 and C 2 are as-yet-undetermined renormalization constants. Substituting (23) into (1) (in D rather than four dimensions) gives Diagram RI Diagram RII Figure 2: Renormalization-induced Feynman diagrams that provide a LO perturbative contribution to the mixed correlator. The square insertion denotes the current (22).
The last two terms on the right-hand side of (24) each generate a new renormalization-induced Feynman diagram, the pair of which are shown in Figure 2. Note that a square insertion represents the current (22). Evaluating these two diagrams and choosing C 1 and C 2 such that the right-hand side of (24) is free of nonlocal divergences, we find as well as an updated expression for Π (I) from (13) that is free of nonlocal divergences where, again, we have omitted polynomials in z as they will not contribute to the LSR. In summary, taking operator mixing into account, the LO QCD expression Π (OPE) can be decomposed as in (13) with the terms on the right-hand side given by (27) and (17)-(21).

QCD Laplace Sum-Rules
The function Π from (2) satisfies a dispersion relation where Π on the left-hand side is to be identified with the QCD prediction Π (OPE) ; ImΠ(t) is the hadronic spectral function; t 0 is the hadron threshold parameter; and · · · represents subtraction constants, collectively a quadratic polynomial in Q 2 . To eliminate these subtraction constants as well as local divergences in Π (OPE) and to accentuate the resonance contributions of the hadronic spectral function to the integral on the right-hand side of (28), we apply the Borel transform with Borel parameter τ to formulate the 0 th -order LSR [7] On the right-hand side of (30), we use a "resonance(s) plus continuum" model where ρ (had) represents the resonance content of the spectral function (to be discussed further in Section 4), θ is the Heaviside step function, and s 0 is the continuum threshold. Then, we define the continuumsubtracted 0 th -order LSR To compute R 0 (τ, s 0 ), we use the following identity relating the Borel transform to the inverse Laplace transformL −1 [7]: where c is any real number for which f (Q 2 ) is analytic for Re(Q 2 ) > c. Generalized hypergeometric functions of the form p F p−1 have a branch cut along the positive real semi-axis originating at the branch point z = 1. As such, in the complex Q 2 -plane, Π (OPE) (Q 2 ) is analytic except for a branch cut along the negative real semi-axis originating at a branch point and deform the integration contour on the right-hand side to that shown in Figure 3. Then, we apply definitions (30) and (32) to find where and, from (27) and (17)- (21), we get For both Π (I) and Π (II) , the first integral on the right-hand side of (34) converges and the second vanishes for η → 0 + . For Π (III) -Π (VI) , however, each integral diverges although their sum is finite. To isolate this finite contribution, we first expand the imaginary parts (38)-(41) near t = 4m 2 : where is analytic in a neighbourhood about t = 4m 2 . When (42) is inserted into the first integral on the righthand side of (34), the part of the result stemming from the p(t) term converges whereas the parts stemming from the (t − 4m 2 ) −2 and (t − 4m 2 ) −1 terms diverge. Focusing on these divergent parts, we have for η → 0 + . In (44), erf is the error function and E n is the exponential integral function On the right-hand side of (44), note that the terms proportional to η −3/2 and η −1/2 diverge whereas the remaining terms are finite. Next, we consider the contributions of Π (III) -Π (VI) to the second integral on the right-hand side of (34). Parameterizing for θ i = −π + to θ f = π − , we find that for η → 0 + . When (44) and (48) are added together, the divergent terms which go like η −3/2 and η −1/2 cancel leaving a finite result. Finally, collecting together (34), (35), (42), (44), and (48), we have where, again, p(t) is given in (43), and the imaginary parts ImΠ (I) and ImΠ (II) are given in (36) and (37). The integral on the right-hand side of (49) can be evaluated analytically; however, the result is long and so we omit it for the sake of brevity. Renormalization-group improvement [26] implies that the strong coupling and quark mass in the simplified (34) get replaced by corresponding running quantities evaluated at renormalization scale µ, i.e., α s → α s (µ) and m → m c,b (µ). At one-loop in the MS renormalization scheme, we have for charmonium For charmonium, we set µ to m c ; for bottomonium, we set µ to m b . Finally, we use the following values for the gluon and quark condensates [27,28,29]:

Analysis and Results
To extract hadron properties from the LSR (49) we must first select an acceptable range of τ values, i.e., a Borel interval (τ min , τ max ). To do so, we follow the same methodology as in [14,15,30,31]. To choose τ max , we demand that the LSR converge in the following sense: the magnitude of the 4d gluon condensate contribution (stemming from Π (II) ) must be less than one-third that of the perturbative contribution (stemming from Π (I) ), and the magnitude of the sum of the 6d gluon and quark condensate contributions (stemming from Π (III) -Π (VI) ) must be less than one-third that of 4d gluon condensate contribution. For charmonium, we find τ max = 0.6 GeV −2 ; for bottomonium, we find τ max = 0.2 GeV −2 . To choose τ min , we consider the pole contribution i.e., the ratio of the LSR's hadron contribution to its hadron plus continuum contribution, and demand that it be at least 10%. In both the charmonium and bottomonium analyses, the value of τ min selected using this prescription depends weakly on s 0 , a parameter not known at the outset. Hence, we first choose reasonable seed values for s 0 : s 0 = 25 GeV 2 for charmonium and s 0 = 130 GeV 2 for bottomonium. When input into (61), these two seed values correspond to τ min = 0.1 GeV −2 for charmonium and τ min = 0.01 GeV −2 for bottomonium. After making predictions for s 0 through the optimization procedure explained below, we then update τ min using the new, predicted value of s 0 . In all cases considered, the effect on τ min was insignificant.
Next, we turn our attention to ρ (had) from (31). As ρ (had) represents the resonance(s) portion of the hadronic spectral function, it contains those hadrons which couple to both the meson current (3) and the hybrid current (4). Such hadrons can be thought of as mixtures that have a qq-meson and a hybrid component. Our analysis approach is to input known vector heavy quarkonium resonances into ρ (had) in order to test them for meson-hybrid mixing. In Table 1, we list all vector charmonium resonances that have a Particle Data Group entry in [6], and in Table 2, we do the same for bottomonium. (Note that, in Table 1, states named with a ψ or J/ψ have I G = 0 − whereas those named with an X have unknown I G .) All resonances listed in the two tables have widths 100 MeV. In general, LSRs are insensitive to resonance widths of up to several hundred MeV, and so, we ignore the widths of individual resonances. But, for a cluster of resonances for which the mass difference between successively heavier states is 250 MeV, we amalgamate the cluster into a single resonance with nonzero effective width. And so, we consider a variety of ρ (had) of the form where n is the number of distinct resonances (or clusters of resonances) and where each ρ or, for a resonance cluster, a rectangular pulse with effective width Γ i = 0 in which the resonance strength is uniformly distributed over are mixing parameters related to the combined effect of coupling to the hybrid and qq-meson currents. A state with both qq-meson and hybrid components has ξ i = 0; a pure qq-meson or pure hybrid state has ξ i = 0. The specific models for which we present results are defined for the charmonium and bottomonium sectors in Tables 3 and 4 respectively.
Substituting (62) into (32) gives   Table 3: A representative collection of hadron models analyzed in the charmonium sector.
Model  Table 4: A representative collection of hadron models analyzed in the bottomonium sector. Model As a specific example, consider a ρ (had) that has three resonances with masses {m 1 , m 2 , m 3 }. If the first two resonances are narrow (i.e., Γ 1 = Γ 2 = 0) and the third has Γ 3 = 0, then and For particular choices of {m i } n i=1 and {Γ i } n i=1 , the quantities {ξ i } n i=1 and s 0 are extracted as best-fit parameters to (65). More precisely, we partition the Borel interval (τ min , τ max ) into N = 20 equal length subintervals with {τ j } N j=0 and define With the specific ρ (had) (t) given in (67), for example, eqn. (69) becomes Minimizing (69) gives predictions for {ξ i } n i=1 and s 0 corresponding to the best fit agreement between QCD and the hadronic model in question. For the models defined in Tables 3 and 4, our results are shown in  Tables 5 and 6 respectively. Rather than present each ξ i , we instead present ζ and ξ i ζ where The errors included are associated with the strong coupling reference values (54)-(55), the quark mass parameters (56)-(57), the condensates (58)-(60), and an allowed ±0.1 GeV variability in the renormalization scale [32]. We also allow for the end points of the Borel interval to vary by half the value of τ min , i.e., 0.05 GeV −2 in the charmonium sector and 0.005 GeV −2 in the bottomonium sector. We don't vary κ from (11) as the numerical contribution to the LSR (49) stemming from the 6d quark condensate diagram is negligible. Our results are most sensitive to varying the quark mass parameters.

Discussion
As can be seen from Tables 5 and 6, in both the charmonium and bottomonium sectors, the inclusion of a third heavy resonance cluster in the analysis significantly improves the fit between QCD and experiment as measured by (69). The improvement is particularly dramatic for bottomonium. It is important to note that these third resonance clusters make large contributions to the LSR, i.e., the right-hand side of (65), despite the fact that high mass states are suppressed relative to low mass states due to the exponentially decaying kernel. As a quantitative measure of the excited state signal strength, consider   (5) the ratio of the third resonance's net contribution to the LSR to the sum (of the magnitudes) of the contributions made by all three resonances. In the charmonium sector, evaluating (72) for model 3 from Table 5 gives 0.43. In the bottomonium sector, evaluating (72) for model 3 from Table 6 gives 0.35. Thus the signal strength of the excited state is significant, as expected by its clear effect of reducing the χ 2 -values in Tables 5 and 6. Including one or more resonance widths in the analysis has almost no impact on the quality of fit between QCD and experiment as can be seen from the value of the minimized χ 2 of model 3-5 in Table 5 and models 3-4 in Table 6. This is unsurprising given the general insensitivity of LSRs to resonance width.
In both charmonium and bottomonium sectors, including a fourth resonance or resonance cluster in ρ (had) leads to a χ 2 that minimizes at s 0 ≈ m 2 4 , i.e., the heaviest resonance essentially merges with the continuum, contrary to the initial assumption articulated in (31) that there is a separation between resonance physics and the continuum. Furthermore, as can be seen from both Tables 5 and 6, the two-resonance scenario model 2 also suffers from this problem which gives us another reason to disfavour it compared to the three-resonance models.
Focusing on the three-resonance models in the charmonium sector (model 3-5 in Table 5), we find a nonzero mixing parameter for the J/ψ; essentially no evidence for mixing in the ψ(2S), ψ(3770) resonance cluster; and a large mixing parameter corresponding to a resonance (or resonance cluster) of mass (or average mass) 4.3 GeV. We investigated the effect of varying the mass of the third resonance, m 3 , from 4.0 GeV-4.6 GeV. We found that the minimum value of the χ 2 was indeed lowest for m 3 = 4.3 GeV, about one-third the value for either m 3 = 4.0 GeV or m 3 = 4.6 GeV.
Given the lack of evidence for meson-hybrid mixing in the ψ(2S), ψ(3770) resonance cluster, it is reasonable to exclude it from ρ (had) . As can be seen from models 6-7 in Table 5, doing so has a small effect on the fitted values of ξ 1 , ξ 3 , and s 0 as well as the minimum value of the χ 2 .
Focusing on the three-resonance models in the bottomonium sector (models 3-4 in Table 6), we find a nonzero mixing parameter for all three resonances, i.e., the Υ(1S), the Υ(2S), and the Υ(3S), Υ(4S) resonance cluster, indicating that all have qq-meson and hybrid components. In summary, the best agreement between our QCD predictions and experiment is achieved with threeresonance models in both the charmonium and the bottomonium sectors although, in the charmonium sector, omitting the second heaviest resonance cluster has minimal effect on the results. In fact, qq-mesonhybrid mixing in the charmonium sector is well-described by a two resonance model consisting of the J/ψ and a second state with mass 4.3 GeV. It has been hypothesized that the X(4260) might be a resonance with a significant hybrid component [33,34,35]. Our results are certainly consistent with this idea. In the bottomonium sector, our results indicate that there is nonzero qq-meson-hybrid mixing in the Υ(1S), the Υ(2S), and in the Υ(3S), Υ(4S) pair.