Heavy Quarkonium in a Holographic Basis

We study the heavy quarkonium within the basis light-front quantization approach. We implement the one-gluon exchange interaction and a confining potential inspired by light-front holography. We adopt the holographic light-front wavefunction (LFWF) as our basis function and solve the non-perturbative dynamics by diagonalizing the Hamiltonian matrix. We obtain the mass spectrum for charmonium and bottomonium. With the obtained LFWFs, we also compute the decay constants and the charge form factors for selected eigenstates. The results are compared with the experimental measurements and with other established methods.


Introduction
Describing hadrons from quantum chromodynamics (QCD) remains a fundamental challenge in nuclear physics. Inspired by the discovery of a remarkable gauge/string duality [1], holographic QCD models, most notably the AdS/QCD [2], have been proposed as analytic semi-classical approximations to QCD (for a recent review, see Ref. [3]). In light of these phenomenological successes, as well as the recent progress in the ab initio nuclear structure calculations [4][5][6][7], the basis light-front quantization (BLFQ) [8] has been developed as a non-perturbative approach to address QCD bound-state problems from first principles.
BLFQ is based on the Hamiltonian formalism in light-front dynamics (LFD, [9]) in Minkowski space. The central task of the Hamiltonian approach is to diagonalize the QCD Hamiltonian operator, Here P ± = P 0 ± P 3 is the longitudinal momentum and the light-front quantized Hamiltonian operator, respectively. The eigenvalues directly produce the invariant-mass spectrum. The eigenfunctions, known as the light-front wavefunctions (LFWFs), play a pivotal role in the study of the hadron structures in deep inelastic scattering (DIS) [10] and deeply virtually Compton scattering (DVCS) [11]. In the Fock space expansion, Eq. (1) becomes a relativistic quantum many-body problem and can be solved by constructing and diagonalizing the many-body Hamiltonian matrix (see, e.g., [12] for a review).
The advantages of LFD are made explicit by BLFQ which can employ an arbitrary single-particle basis subject to completeness and orthonormality. By adopting a single-particle AdS/QCD basis, BLFQ naturally extends the AdS/QCD LFWFs to the multi-particle Fock sectors [8]. Furthermore, this basis preserves all the kinematical symmetries of the full Hamiltonian [13,14]. Such choice is in parallel with the no-core shell model (NCSM) used in non-relativistic quantum many-body theory [5]. Stateof-the-art computational tools developed in the many-body theory can be used to address the QCD eigenvalue problem [15]. BLFQ has been applied successfully to a range of non-perturbative problems, where, x = p + q /P + is the longitudinal momentum fraction of the quark, k ⊥ = p q⊥ − xP ⊥ is the relative transverse momentum, and r ⊥ is the transverse separation of the partons. κ is the strength of the confining potential. Note that the "light-cone Hamiltonian" has mass squared dimension, whose eigenvalues are the squared invariant masses. Following Brodsky and de Téramond [46], it is convenient to introduce the holographic coordinate ζ ⊥ = x(1 − x)r ⊥ , and its conjugate q ⊥ = k ⊥ / x(1 − x) ≡ −i∇ ζ ⊥ . In light-front holography, ζ ⊥ is mapped to the fifth coordinate z of the AdS space. In these coordinates, H sw is a harmonic oscillator (HO), Its eigenvalues follow the Regge trajectory M 2 = 2κ 2 (2n + |m| + 1). Its eigenfunctions are 2D HO functions in the holographic variables, Here q ⊥ = |q ⊥ |, θ = arg q ⊥ , and L m n (z) is the associated Laguerre polynomial. We adopt these functions as our basis. This basis has the advantage that in the many-body sector, it allows the exact factorization of the center-of-mass motion in the single-particle coordinates. This is a very valuable property, because the boson/fermion symmetrization/anti-symmetrization in the relative coordinates quickly becomes intractable, as the number of identical particles increases [13,14]. For this work, however, we do not have identical particles in the qq sector and we will use the relative coordinate. In future extensions, as sea quarks and gluons are added, it may be more advantageous to adopt single-particle coordinates.
The soft-wall Hamiltonian Eq. (3) is designed for massless quarks, and it is inherently 2-dimensional. For the heavy quarkonium systems, it should be modified to incorporate the quark masses and the longitudinal dynamics, Here V L is a longitudinal confining potential. Several longitudinal confining potentials have been proposed [47][48][49]. Here we propose a new longitudinal confinement which shares features with others proposed, where ∂ x ≡ (∂/∂x) ζ ⊥ . This term combined with the mass term from the kinetic energy forms a Sturm-Liouville problem, where α = 2mq(m q + mq)/κ 2 , β = 2m q (m q + mq)/κ 2 . The solutions of Eq. (7) are analytically known in terms of the Jacobi polynomial P and form a complete orthogonal basis on the interval [0, 1]. The corresponding eigenvalues are λ (α,β) l = (l + 1 2 (α + β))(l + 1 2 (α + β) + 1). (l = 0, 1, 2, · · · ) Comparing to other forms of longitudinal confinement, our proposal implements several attractive features. First, the basis functions resemble the known asymptotic parton distribution ∼ x α (1 − x) β with α, β > 0 [50]. This is our primary motivation for adopting the longitudinal confinement Eq. (6). Second, the basis function is also analytically known, which brings numerical convenience. Third, in the non-relativistic limit min{m q , mq} κ, the longitudinal confinement sits on equal footing with the transverse confinement, where together, they form a 3D harmonic oscillator potential, and rotational symmetry is manifest. This non-relativistic reduction also provides us an order-ofmagnitude estimate of the model parameters for our heavy quarkonium application. Fourth, in the massless limit max{m q , mq} κ, the longitudinal mode stays in the ground state and the longitudinal wavefunction χ 0 (x) = const. Thus we restore the massless model of Brodsky and de Téramond 1 .

One-Gluon exchange
Following Ref. [20] (cf. [22][23][24][25]), we introduce the one-gluon exchange in LFD. In the momentum space, this term reads, where S s,s,s , is the average momentum transfer. We have also introduced a gluon mass µ g to regularize the Coulomb singularity. This singularity is integrable and does not carry physical significance. The gluon mass is used to improve the numerics and will be taken small compared to other scales in this application. In Eq. (11), we have taken the total momentum P ⊥ = 0, P + = 1 by virtue of the boost invariance in LFD.
This one-gluon exchange introduces a logarithmic divergence [21,23,25,51]. This can be seen from the one-gluon exchange kernel by transverse power counting. In particular, the leading power comes from the spin non-flip spinor matrix elements S ↑↓↑↓ and S ↓↑↓↑ , which contain terms proportional to k 2 ⊥ or k 2 ⊥ . Such terms do not vanish in the large momentum limit 2 . In principle, this divergence can be handled by a proper renormalization [52][53][54][55][56]. However, such a procedure is usually rather challenging or even impractical. In this work, we adopt a simple counterterm technique [20,21,23,24,51] by replacing the spinor matrix elements S ↑↓↑↓ and S ↓↑↓↑ with, (s,s, s ,s = ↑↓↑↓ or ↓↑↓↑) (12) In its essence, this counterterm exactly removes the troubling k 2 ⊥ and k 2 ⊥ terms in the spinor matrix elements hence removing the UV divergence.

Basis light-front quantization
In Sect. 2, we derived the effective light-cone Hamiltonian operator for quarkonium, which reads, In this section, we will focus on solving the eigenvalue equation non-perturbatively, to obtain the mass spectrum and the LFWFs ψ J m J (k ⊥ , x, s,s) ≡ k ⊥ , x, s,s|ψ J P C m J . Here J, P , C and m J are the total spin, parity, charge parity, and the magnetic projection of the state, respectively. Our strategy is to construct a finite-dimensional matrix from the effective Hamiltonian operator and then diagonalize it numerically. To do this, we first adopt a basis expansion, Here ψ J m J (n, m, l, s,s) ≡ n, m, l, s,s|ψ J P C m J are the LFWFs in the BLFQ basis. As mentioned above, the basis functions Φ nml (q ⊥ , x) ≡ φ nm (q ⊥ )χ l (x) are generated by the confining part of the Hamiltonian -the kinematic energy plus the confinement. As such, the light-front holographic QCD results serve as our first approximation. The matrix elements of the Hamiltonian within this basis are, n , m , l , s ,s |H eff |n, m, l, s,s = n , m , l , s ,s |V g |n, m, l, s,s We take advantage of the conservation of the angular momentum in the transverse plane, and choose a particular magnetic projection m J : m + s +s = m J . In order to carry out practical calculations, we truncate the infinite basis space to a finite size by restricting the quantum numbers according to 2n + |m| + 1 ≤ N max , l ≤ L max . N max controls the range of momenta covered by the harmonic oscillator basis as it is related to a transverse IR regulator λ ir ∼ b/ √ N max and a UV regulator Λ uv ∼ √ N max b [57], where b is the HO basis scale and we have taken b = κ, unless otherwise stated. L max controls the basis resolution in the longitudinal direction.
The relativistic bound states are identified by three quantum numbers J P C . However, P is broken by the finite basis truncation, because the parity transformation is dynamical in LFD. Instead, one can exploit the so-called mirror parityP x =R x (π)P that is related to P [58], Then P can be extracted from the LFWFs by, Similarly, C can be obtained from the LFWFs by, Rotation is also a dynamical symmetry in LFD. As a result, J is not conserved and the mass degeneracy for different magnetic projections m J is lifted by the Fock-space truncation. Nevertheless, in a nonrelativistic system, such as the heavy quarkonium, the discrepancy is small and we can still extract J by counting the multiplicity of the nearly-degenerate mass eigenstates. It is also instructive to assign states the non-relativistic quantum numbers n 2S+1 L J , where n, S and L are the radial quantum number, the spin and the orbital angular momentum, respectively. L, S are related to parity and charge conjugation through P = (−1) L+1 , C = (−1) L+S . The radial quantum number n can be deduced from the mass hierarchy of the spectrum. The total spin S, though not an exact quantum number, can be obtained from its expectation value, ψ J P C m J | S 2 |ψ J P C m J = S(S + 1). L, S, J are constrained by the angular momentum addition |L − S| ≤ J ≤ L + S. Reconstructing these quantum numbers allows us to identify the states and to compare with experimental data and with other methods.

Numerical Results
In this work, we focus on charmonium and bottomonium, where the fermion masses are equal (m q = mq) and heavy (m q κ). The Hamiltonian matrix element (Eq. (16)) involves a four-dimensional integral (two in the radial direction, two in the longitudinal direction). They are evaluated using Gauss quadratures. The number of the quadrature points N rad (N lfx ) is taken to be at least twice N max (L max ). Then the obtained matrix is diagonalized using LAPACK software [59]. We fix the bottomonium coupling α s (M bb ) = 0.25 and obtain the charmonium coupling α s (M cc ) = 0.3595 using the leading-order pQCD evolution of strong coupling (N f = 5, M cc 3.5 GeV, M bb 9.5 GeV). We then fit the parameters κ, m q (m c , m b ) by minimizing the root-mean-squared (r.m.s.) deviation from the experimental masses [60] of the states below the DD or BB threshold. The details of the model parameters are summarized in Table 1. In these calculations, we take the regulators µ g = 0.02 GeV κ, N max = L max = 8, 16, 24 (N max = 8 is available in the supplemental materials), and the number of quadrature points N rad = N lfx = 64.
Previous BLFQ study of positronium [20,21] shows that the continuum limit N max → ∞, L max → ∞, µ g → 0 can be reached through extensive calculation and successive extrapolations. However, we shall not investigate the N max , L max and µ g extrapolations in this paper because of the numerical efforts involved and also the presence of several fitted parameters. Comparing the results with N max = L max = 16 and N max = L max = 24 in Table 1, where the model parameters κ, m q were fitted separately and turned out to be very close ( 1% change). The r.m.s. deviations from the measured spectra are also comparable. We also studied the trend of the mass eigenvalues with decreasing µ g in the range of 10 −4 GeV ≤ µ g ≤ 3 × 10 −1 GeV at fixed N max = L max , and found that the mass eigenvalues are well converged with respect to µ g . Table 1: Summary of the model parameters. The coupling α s and the gluon mass µ g are fixed (see the text). The confining strength κ and the quark mass m q are fitted with m J = 0 using the experimental data below the DD or BB threshold. The r.m.s. deviations for the m J = 0 spectrum, δM m J =0 , and the r.m.s. average-m J spectrum, δM , are computed for states below the threshold.

4.2
Mass    [51]. We compare our results with the experimental data from the particle data group (PDG, [60], cf. [61]). The r.m.s. deviations are computed for charmonium (bottomonium) below the DD (BB) threshold (see Table. 1). Comparing with the recent results in Ref. [43], our approach improves the charmonium and bottomonium mass spectra from light-front holography.
The LFWFs can be used to calculate the transition amplitudes. Here we consider the decay constants. These quantities are useful for computing the decay widths and constraining the Standard Model parameters. In the non-relativistic limit, they are proportional to the wavefunctions at the origin and, therefore, test the short-distance physics of the model. The decay constants for a scalar S, a pseudo-scalar P , a axial-vector A, and a vector V states are defined as, respectively, where e µ λ (p) is the spin vector for the vector boson, Here ⊥ ± = (1, ±i)/ √ 2. Note that due to the charge conjugation symmetry, the decay constants of 0 ++ and 1 +− vanish. In LFD, the decay constants can be computed from the "+" component of the current decay constants [GeV ] Figure 2: The decay constants of scalar and vector quarkonia as compared with PDG data [60] as well as Lattice QCD (Lattice) [62][63][64][65] and Dyson-Schwinger (DSE) [32] approaches. We extrapolate N max using second-order polynomials in N −1 max and adopt the difference between the extrapolated value (N max extr.) and the N max = 24 value as the uncertainty (not including systematic errors). Note that Υ(nD) are referring to the vector bottomonia Υ(n 3 D 1 ). matrix elements, which read [10], The decay constants from this work are plotted in Fig. 2. We also list the PDG data [60], as well as results from Lattice QCD (Lattice, [62][63][64][65]), and the Dyson-Schwinger equation (DSE, [32]). The PDG data are extracted from the dilepton decay widths Γ ee (for vectors) and diphoton decay widths Γ γγ (for pseudo-scalars). We obtain results for three successive sets of basis regulators, N max = L max = 8, 16, 24, where the parameters κ and m q are fitted to the mass spectrum separately, and α s , µ g are kept fixed. We then extrapolate N max using simply polynomials in N −1 max and estimate the uncertainty associated with N max from the difference between the extrapolated and the N max = 24 results. While the resultant masses are close as we mentioned, the decay constants show noticeable residual regulator dependence. This may not be a surprise as the decay constant probes the short-distance physics, whereas the basis is chosen to emulate confinement. We expect slower convergence for the decay constants than the masses. From Fig. 2, our calculated decay constants are in reasonable agreement with the known experimental measurements as well as Lattice and DSE results. This is encouraging since the measured decay constants are not used in our fits. However, compared to Lattice and DSE results, most of our results are systematically larger than the PDG data. This is likely due to the systematic errors of our model and can be improved by more realistic models for the LF Hamiltonian and/or by including higher Fock spaces. Fig. 2 also includes decay constants for the D-wave states. In the non-relativistic limit, these quantities should vanish. The small but non-vanishing D-wave decay constants in our results indicate the mixing of the S-wave component, as expected from the relativistic treatment.
The LFWFs also provide direct access to other hadronic observables. Here we study the fictitious charge form factor with the photon coupling only to the quark (but not the anti-quark). Of course, this quantity is not a physical observable. Nevertheless, it provides important insight to the system. In particular, this charge form factor at small momentum transfer yields the r.m.s. radius of the hadron, In LFD within the impulse approximation, the form factors can be obtained from the Drell-Yan-West

10
. Figure 3: The charge form factors G 0 (Q 2 ) for η c , J/ψ, χ c0 , η b , Υ, and χ b0 as extrapolated from N max = L max = 8, 16, 24 using second-order polynomials in N −1 max . The difference between the extrapolated and the N max = 24 values are used to quantify the uncertainty, which does not include systematic errors. Note that, for the scales we are showing, the charge form factor of η b and Υ are on top of each other.   (8) Lattice [70] 0.063(1) 0.066(2) 0.095 (6) DSE [71] 0.048(4) 0.052 (3) formula within the Drell-Yan frame P + = P + [66], where q = P − P , and Q 2 = −q 2 = q 2 ⊥ . For (pseudo) scalars, it directly produces the charge form factor G 0 (Q 2 ) = I 0,0 (Q 2 ). For (axial) vector mesons, due to the violation of the rotational symmetry, there exists some ambiguity on finding the physical (Sachs) form factors G 0 , G 1 , and G 2 from I m J ,m J [67]. We adopt the prescription of Grach and Kondratyuk [67], which has been shown to be free of zero-mode contributions in some analytical models [68,69]. The charge form factor according to this prescription reads, where η = Q 2 /(4M 2 ), M is the mass of the hadron. Fig. 3 shows the charge form factors of η c , J/ψ, χ c0 , η b , Υ, and χ b0 . Table 2 lists the r.m.s. radii of the first few states. From our results, the radius of J/ψ is close to, but slightly larger than, that of η c . This is consistent with the Lattice and DSE results and can be understood from the non-relativistic point of view. Bottomonia are in general smaller than charmonia, a result, again, that can be drawn from the non-relativistic argument. The comparison of the first few available results shows that radii from our approach are in qualitative agreement with those of other approaches [70,71], though ours are systematically smaller than the Lattice and DSE results. This observation is consistent with the trend of the decay constants.

Summary and Outlook
We studied heavy quarkonium based on holographic QCD and a realistic one-gluon exchange on the light front. We proposed a longitudinal confining potential to incorporate quark mass and longitudinal dynamics. We solved the bound-state problem in the Basis Light-Front Quantization approach. We calculated the spectroscopy, decay constants and charge form factors for the charmonium and bottomonium. The comparison with the experimental data and results from Lattice QCD and Dyson-Schwinger equation shows reasonable agreement among the available observables. As such, we improve the light-front holography from the first approximation to QCD.
The LFWFs may readily be used in the study of, e.g., radiative transitions [26] and diffractive vector meson production in ultra-peripheral heavy ion collisions [72]. These observables will also serve as stringent tests of our model in different regimes.
This model can also be applied to heavy-light and light-light systems, where the light-front holography provides a good initial approximation [45]. The phenomenological confinement can be improved by a better understanding of the string/gauge duality as well as a more complete derivation of the interquark potentials from various first-principle approaches to QCD [73]. Ultimately, the phenomenological confining interaction should be replaced by the QCD Hamiltonian.
This work can be extended to higher Fock sectors to incorporate sea-quark and gluon degrees of freedom as well. The treatment of the many-body dynamics is essential for obtaining realistic predictions for states above the thresholds. One of the major challenges of the light-front Hamiltonian approach however is that explicitly including many gluons, as well as many quark-antiquark pairs, in the Fock space expansion quickly becomes numerically expensive. The ever increasing computational capacity and the progress in ab initio many-body calculations represent growing opportunities for understanding the strong interaction in the Basis Light-Front Quantization approach. For systems involving light quarks, however, the naïve Fock sector truncation may not suffice, as addressing dynamical chiral symmetry breaking is essential [74]. The coherent basis (see, e.g., [75]) and the light-front coupled-cluster method [76] are two promising approaches for dealing with collective modes and dynamical symmetry breaking in the light-front Hamiltonian formalism.
Supplementary data for "quarkonium in a holographic basis" Yang Li, Pieter Maris, Xingbo Zhao, and James P. Vary April 5, 2016 Table 1: Summary of the model parameters. The coupling α s and the gluon mass µ g are fixed. We choose the bottomonium α s (9.5 GeV) = 0.25 and then obtain the charmonium coupling α s (3.5 GeV) 0.36 from pQCD evolution of the strong coupling. The confining strength κ and the quark mass m q (m q = mq) are fitted with m J = 0 using the experimental data below the DD or BB threshold. The r.m.s. deviations for the m J = 0 spectrum, δM m J =0 , and the r.m.s. average-m J spectrum, δM , are computed for states below the threshold. The numbers for κ, m q , δM m J =0 and δM are rounded.