Approximate Hamiltonian for baryons in heavy-flavor QCD

Aiming at relativistic description of gluons in hadrons, the renormalization group procedure for effective particles (RGPEP) is applied to baryons in QCD of heavy quarks. The baryon eigenvalue problem is posed using the Fock-space Hamiltonian operator obtained by solving the RGPEP equations up to second order in powers of the coupling constant. The eigenstate components that contain three quarks and two or more gluons are heuristically removed at the price of inserting a gluon-mass term in the component with one gluon. The resulting problem is reduced to the equivalent one for the component of three quarks and no gluons. Each of the three quark-quark interaction terms thus obtained consists of a spin-dependent Coulomb term and a spin-independent harmonic oscillator term. Quark masses are chosen to fit the lightest spin-one quarkonia masses most accurately. The resulting estimates for bbb and ccc states match estimates obtained in lattice QCD and in quark models. Masses of ccb and bbc states are also estimated. The corresponding wave functions are invariant with respect to boosts. In the ccb states, charm quarks tend to form diquarks. The accuracy of our approximate Hamiltonian can be estimated through comparison by including components with two gluons within the same method.


INTRODUCTION
Quark model represented baryons as bound states of three quarks, e.g. see [1,2]. In QCD, baryons are instead superpositions of states of quanta of quark and gluon fields. A priori, the number of quanta varies from three to infinity, across an infinite set of components. The range of momenta these quanta may have is also infinite. Their interactions diverge with their momenta. This article contributes to a development of a Hamiltonian approach to QCD that appears capable of filling the gap between the complex quantum-field picture and simple quark-model picture for hadrons [3]. We report results for heavy baryons that are obtained using our approach in its pilot version, which includes severe simplifications.
We limit the theory to quarks that have masses much greater than Λ QCD and we consider the weak coupling limit. Creation of quark-antiquark pairs is neglected. Components with more than one gluon are eliminated, by assuming that their dominant effect in the component with one gluon is that the gluon has a mass, allowed to be a function of the gluon kinematic relative momentum with respect to the quarks. We use secondorder perturbation theory to derive the resulting effective Hamiltonian for baryons that only acts in the component with three quarks. We compare the quark interaction terms in this Hamiltonian to similar terms in the Hamiltonian that only acts in the quark-anti-quark component in quarkonia, previously derived using the same method [3]. Masses of heavy baryons are estimated by extrapolating the coupling constant at the quark-mass scale from a formally infinitesimal value of the weak-coupling limit to the value implied by the known coupling constant at the scale of Z-boson mass. Quark masses are adjusted to the known spectra of heavy quarkonia. Our estimates for baryon masses contain no new parameters. The scale parameter for hadrons built from different flavors is fixed by a linear interpolation between its one-flavor values.
The concept of effective-gluon mass that we use is explained in Sec. II, followed by a brief outline of our method in Sec. III. The method is called the renormalization group procedure for effective particles (RGPEP). The second-order baryon eigenvalue problem is outlined in Sec. IV. Details of the effective quark-quark interaction terms in ccc and bbb baryons, implied by the gluon mass, are described in Sec. V, including a comparison with the case of heavy quarkonia. Sec. VI extends the calculation to the ccb and bbc baryons. The resulting estimates for baryon masses are described in Sec. VII. Comments concerning the RGPEP calculation of effective Hamiltonians in QCD in orders higher than second and beyond perturbation theory conclude the paper in Sec. VIII. Details of our fit to the spectra of quarkonia are described in App. A. Values of the RGPEP scale parameter we use are listed in App. B. Appendix C discusses dependence of harmonic oscillator frequencies on the scale parameter. Appendix D provides a detailed description of the baryon wave functions that are used in our estimates and App. E presents explicit formulas for the associated heavy-baryon masses.

II. ASSUMPTION OF GLUON MASS
Theoretically, a baryon state in heavy-flavor QCD is a superposition of states of virtual, point-like quarks and gluons, Components that include quark-anti-quark pairs are considered very small because quarks are heavy. In contrast, components with gluons are included because in canonical QCD gluons are massless. However, using the massless gluons in the expansion and limiting their number in a computation one expects to obtain the spectrum of excited baryons that gets dense toward the free quark threshold. The same feature is expected to occur in such computations of spectrum of quarkonia. Physically, the latter is observed to be not dense [4,5]. For example, the s-wave cc or bb mass splittings are on the order of half GeV and they do not decrease with the excitation number as they would if the interaction was purely of the Coulomb type, like in QED with massless photons. The mass splittings in the heavy baryon spectrum are also not expected to rapidly decrease with the excitation number.
To describe physical splittings, excitations of the gluon field must involve considerable energy. This requirement can be addressed using the concept of a gluon mass [3]. Introduction of a mass term for gluons in the canonical Hamiltonian of QCD would spoil its gauge-theory structure. Instead, we introduce a gluon mass in solving the eigenvalue problem of a Hamiltonian H t that is derived using the RGPEP [6], see below. The parameter t is the renormalization group parameter. It is useful to think about it as t = s 4 , where s has an interpretation of the size of effective particles. Canonical gluons are considered point-like. Instead of using canonical gluons, the goal is to represent a heavy baryon by a superposition where the quarks and gluons are the effective particles of size s. We introduce the gluon mass in the effective QCD eigenvalue problem for heavy baryons within the same computational scheme that we previously applied to heavy quarkonia [3]. Our leading principle is that the gluon mass µ t is the minimal price we have to pay for limiting the expansion in Eq. (2) to the first two terms. Such limitation makes sense because interactions in the Hamiltonian H t contain vertex form factors. The form factors are obtained by solving the RGPEP Eq. (16), with the initial condition provided by the canonical QCD Hamiltonian with regularization and counterterms. These form factors cause that interactions cannot change invariant masses of component states by amounts exceeding 1/s. Therefore, the effective gluons of size s cannot be as copiously produced as the point-like gluons can in the canonical representation of QCD. It is plausible that inclusion of a few effective components is sufficient to accurately describe a heavy-baryon solution.
Although the number of gluons G t that need to be included in the effective representation of a low-mass solution to the QCD eigenvalue problem is expected to be limited, direct inspection shows that inclusion of even a few components still leads to a mathematically difficult equations for their coupled-channel dynamics. These equations are hard to understand. The results presented in this paper follow from the effective eigenvalue problem that is obtained by using the hypothesis that the contribution of all components other than |3Q t and |3Q t G t may be approximated by inclusion of a gluon mass, µ t for the gluon in component |3Q t G t .
The mass assumption is falsifiable by extending the calculation to explicitly include more components and relegating the gluon mass ansatz to states with more gluons than one. The purpose would be to verify if the finite value of the RGPEP parameter s on the order of quark Compton wave-length is sufficient to prevent the spill of probability to states with many gluons, especially when the coupling constant is small. The latter situation is expected to occur for the quark masses that are much greater than Λ QCD . This is precisely the reason for us to study the dynamics of gluons using the RGPEP first in the context of heavy-flavor QCD. In order to simplify the problem and thus increase a chance of understanding the dynamics of effective gluons whose masses are likely to be much larger than Λ QCD , we exclude from the theory the quarks that have masses much smaller than Λ QCD . If the latter were included in the theory, they could appear in large numbers in the effective Fock-space basis and complicate the dynamics, as massless gluons do.

III. RGPEP FOR HADRONS
The RGPEP provides equations for calculating the renormalized Hamiltonian H t from the canonical one that includes regularization and counterterms. It is also used to calculate the counterterms. Eigenstates of H t define hadrons in terms of effective particle basis in the Fock space. We first consider QCD of only one flavor of heavy quarks, useful in discussing dynamics in baryons made of quarks of one flavor. The case of baryons made of two types of heavy quarks is discussed in Sec. VI. We calculate H t using expansion in powers of a formally infinitesimal coupling constant, up to terms of second-order. Results of our second-order calculations are later compared with results obtained in quark models and in lattice approach to QCD.

A. Canonical Hamiltonian
The Lagrangian for one-flavor QCD is We use the front form (FF) of Hamiltonian dynamics [7] and employ canonical quantization to derive the corresponding HamiltonianĤ can QCD in the gauge A + = 0, We adopt the FF notation of Refs. [8,9]. The Hamiltonian operator density, :Ĥ x + =0 : integrated over the front x + = 0, is expressed in terms of the quantum fieldŝ , u pσ and v pσ are the Dirac spinors, ε µ pσ is the transverse-gluon polarization vector, χ c and T c denote three-component color vector for quarks and eight-component color matrix vector for gluons, while σ and c stand for their spins and colors, respectively. We omit the hats and normal ordering symbols in further formulas. In second-order calculation, only two interaction terms count, the quark-gluon interaction, and the instantaneous quark-quark interaction, where

B. Regularization
H can QCD is regularized by inserting cutoff functions r 21.3 and r C 12.1 2 in the interaction vertices, as shown below and further explained in App. A of Ref. [3]. The terms that contribute to the baryon problem are The free, or kinetic term is where E q and E g are the FF quark and gluon energies, respectively. In the quark-gluon vertex, the first term corresponds to emission and the second to absorption of a gluon by a quark. Numbers 1, 2, 3 stand for sets of quantum numbers of particles 1, 2 and 3, and 123 includes integration over momenta and summation over spins and colors of particles 1, 2 and 3. In the factor The tilde over δ indicates an implicit factor 2(2π) 3 , multiplying the Dirac δ-function of momentum conservation. The function r 21.3 cuts off the large relative transverse momenta and small fractions of plus momenta for the particles involved in the vertex, cf. App. A in [3]. The instantaneous interaction of Eq. (8) yields C.

Renormalized Hamiltonian
We call the regularized canonical Hamiltonian for quanta of size s = 0 the initial Hamiltonian, since it provides an initial condition for solving the RGPEP equation. Its solution defines a family of renormalized Hamiltonians H t , which are written in terms of the operators q t that create or annihilate particles of size s > 0, t = s 4 . The latter operators are defined by means of a unitary transformation U t , The idea is that nonzero size eliminates divergent integrals. Hence, the effective Hamiltonians cannot be sensitive to the cutoff parameters in the cutoff functions. This implies that the regularized canonical Hamiltonian needs to be supplemented with counter-terms, which ensures that the renormalized Hamiltonians do not depend on the regularization. The problem is to define U t that generates Hamiltonians H t ≡ H t (q t ) in a suitable operator basis. Instead of directly defining U t , we define H t = H t (q 0 ), in which the products of operators from the bare theory have the coefficients from renormalized theory. By definition, it obeys where prime denotes derivative with respect to the parameter t. G t is called a generator of the RGPEP. It is set to where H free is the free part of H t , and is identical with H free . The tilde above H t means that coefficients in front of interaction terms are multiplied by the square of total +-momentum entering the vertex. The RGPEP design guarantees that the interaction vertices that change invariant mass of interacting particles by more than 1/s are exponentially suppressed [6], cf. Eq. (20).
In the present work, we solve Eq. (16) using expansion in powers of renormalized coupling constant g t , up to second order, The zero-order term, H t0 , corresponds to terms without running coupling. The only difference between H t0 and the bare expression H free is the presence of effective creation and annihilation operators in place of bare ones, The solution for quark-gluon interaction term of order g t , apart from substitution of bare operators by effective ones, differs from the bare theory term of Eq. (12) by the presence of form factors f t 21.3 . Namely, with where M 2 21 is the square of free invariant mass of particles 1 and 2. The second-order terms that matter are the quark-quark interaction and quark self-interaction terms. The quark-quark interaction, H t2 QQ is a sum of two parts. The one that stems from the instantaneous interaction, differs from H QQ inst by effective particle operators and a form factor [3] f t 12.
The other one that stems from exchange of transverse gluons, where j µ ij =ū i γ µ u j , involves factors with S ij = M 2 ij − (m i + m j ) 2 and m i denoting the mass of particle i. In particular, 4 labels the gluon and m 4 = 0. Moreover, the gluon momentum is where (z) denotes the sign of z.
The differences between Eq. (24) and the corresponding equations for quarkonia [3], are: quark current j 22 instead of anti-quark currentj 2 2 , transposition of matrix t 4 22 , overall sign difference, color factor −2/3 instead of 4/3 (when acting on a color singlet state) and symmetrization factor 1/2 that is absent in quarkonia.
Renormalized quark self-interaction mass correction is where D.

Bound-state eigenvalue problem
Thanks to asymptotic freedom [10], g t = g is small in Hamiltonians with small t. This holds for λ = 1/s much larger than the scale of Λ QCD in the RGPEP scheme. We formally consider which allows us to simplify the eigenvalue problem In the baryon eigenstates represented using Eq. (2), we can neglect Fock sectors with more than three quarks, because of the first inequality in Eq. (32). In a matrix form, the eigenvalue problem reads where H t2 = H t2 QQ +H t2 δm and dots stand for the Fock components with more than one effective gluon and for the Hamiltonian terms that involve those components.

E. Gluon mass ansatz in the effective eigenvalue problem
Similarly to the case of heavy quarkonia [3], we remove the Fock components and Hamiltonian matrix elements that involve more than one gluon. The price for this removal is a gluon-mass ansatz for the component |3Q t G t . Indeed, a gluon mass term appears when one uses Gaussian elimination to express the component |3Q t 2G t in terms of component |3Q t G t . Our working hypothesis is that other terms can also be dropped once a mass ansatz is introduced. The reduced eigenvalue problem with the ansatz mass term µ 2 t is This two-component problem is reduced to an equation for the component |3Q t , using a transformation R succinctly described in [11]. Up to terms order g 2 and using notation r and l for right and left states in the matrix elements, one obtains

IV. 3Q EIGENVALUE PROBLEM
The three-quark component of a baryon satisfies the FF eigenvalue equation in which the state |3Q t is defined by FIG. 1. Self-interaction of effective quarks.

A.
Small-x dynamics The interaction kernel U t eff (12; 1 2 ) in Eq. (39) can be written in terms of relative momentum variables x 1/12 , κ ⊥ 1/12 and x 1 /12 , κ ⊥ 1 /12 (that is momenta relative to pair 12). The resulting expression has the same structure as in our quarkonium analysis [3], with the color factor 2/3 instead of 4/3. The gluon mass ansatz in baryons can be different from the one in quarkonia, since it is a function of gluon momentum relative to the three-quark subsystem instead of quark-antiquark subsystem. The small-x singular factors in the interaction do not produce divergences for the same reason as in quarkonia. The gluonexchange integral is finite because we assume that the gluon mass ansatz vanishes when x 5 → 0. An example of required behavior in quarkonia is µ 2 ∼ x δµ 5 κ 2 5 . In the notation used for baryons, the same behavior is described by µ 2 ∼ x δµ 4/12 κ 2 4/12 when x 4/12 → 0. Let us assume that µ 2 ∼ x δµ 4 κ 2 4 in the limit x 4 → 0. When x 4 goes to zero, so does x 4/12 = x 4 /(x 1 +x 2 ), and κ ⊥ 4/12 ≈ κ ⊥ 4 . Therefore, The same reasoning applies to every pair of the quarks that exchange a gluon. Similarly, in the quark self-interaction terms, the integration variables are x 4/i and κ ⊥ 4/i . In the small- Therefore, mass terms are finite when the gluon mass ansatz vanishes properly when x 4 → 0.

V. HARMONIC OSCILLATOR CORRECTION TO COULOMB TERMS
Given Eq. (32), the effective Hamiltonian H eff t can be approximated by its non-relativistic (NR) limit. To define this limit, we introduce a set of convenient mo-mentum variables, cf. Refs. [12,13], where . The second equality in these equations holds only for equal masses. The nonrelativistic limit is defined as K/m → 0, Q/m → 0. It is valid because the relative momentum regions that significantly exceed the RGPEP scale λ m are suppressed by the exponentially-fast vanishing form factors in the interaction vertices of effective particles. In the leading NR approximation, the momenta K 12 and Q 3 are related to the Jacobi momenta: K 12 is the relative momentum of quark 1 with respect to 2 and Q 3 is the relative momentum of quark 3 with respect to the pair of quarks 1 and 2, see Ref. [12] for more details. Generically, we denote by K the relative momentum of a quark with respect to another quark with which it is involved in an interaction term, and we denote by Q the relative momentum of a spectator with respect to the pair in interaction. We introduce three sets of such relative momentum variables, arranged using the cyclic permutation of indices 123: K jk and Q i . In the NR limit, For equal quark masses, we write the baryon mass as M = 3m + B, divide Eq. (39) by 6m, take the NR limit and obtain where V ij C,BF = V C,BF ( K ij , K ij ) and W ij = W ( K ij − K ij ) are, respectively, the Coulomb term with Breit-Fermi (BF) corrections and the additional interaction resulting from the gluon mass ansatz. µ 12 = m/2, µ 3(12) = 2m/3 are the reduced masses. Both V and W are similar to the ones in the quarkonium case [3].
where ∆ K = K − K and the RGPEP form factor is The mass terms can be written similarly to the interaction terms because δm 2 1 t /(2m) = −(2π) −3 d 3K W ( K ). Using the Taylor expansion for the wave functions under the integrals, e.g., one can see that the first term cancels with a mass term. The second term is linear in momentum ∆ K 12 and gives zero after integration. The first non-vanishing term is the third one, quadratic in ∆ K 12 . This term provides a harmonic oscillator potential. Only terms with m = n are non-zero, and As in quarkonia, we assume that the ansatz µ 2 dominates ∆ K 2 in the relevant integration range. In this case, in Eq. (56), µ 2 /(µ 2 + ∆ K 2 ) ≈ 1, which further leads to the conclusion that w n for n = 1, 2, 3, corresponding to different directions in space, are the same. Thus, the effective oscillator interaction respects rotational symmetry in the Jacobi variables. However, there are only two independent relative momenta for three quarks. We distinguish one pair of quarks, e.g., 12, and rewrite the oscillators in terms of K 12 and Q 3 , Thus, we obtain an oscillator force between quarks 1 and 2 and an oscillator force between quark 3 and the pair 12.
Their strengths are in ratio 3/2 : 2. Since the ratio of corresponding reduced masses µ 3 (12) and µ 12 is 4/3, the frequencies of these oscillators are the same and equal This expression differs from the result for quarkonia by a factor √ 3/2, rendering ω 2 baryon /ω 2 meson = 3/4, assuming that m, λ and α are the same for mesons and baryons built from one flavor of heavy quarks. This result is very close to the ratio 5/8 suggested by models that employ the concept of gluon condensate in vacuum [14] or only inside hadrons [12].
The oscillator interaction may appear to be in contradiction with the linear confinement picture in QCD. However, the eigenvalues of the FF Hamiltonian are the baryon masses squared, in distinction from the instant form (IF) Hamiltonian eigenvalues that are the baryon energies, reducing to the baryon masses only for bound states at rest. At large distances between quarks, the quadratic potential in the FF corresponds to the linear potential in the IF of Hamiltonian dynamics [13].

VI. TWO FLAVORS OF HEAVY QUARKS
Several new elements appear when one of the three quarks, say quark 3, is of different flavor than the other two. Besides smaller particle-exchange symmetry and the fact that β 1 = β 3 = 1/3, a new feature emerges that the NR effective quark masses are modified.
The reason for NR mass modification is that, when we deal with two different flavors of quarks, the constant term that cancels completely in Eq. (59) for the same flavor no longer does so for different flavors. A finite function of x and κ ⊥ is left and it multiplies ψ (1, 2, 3). This effect is small, but in principle ought to be considered. The correction shifts the minimal invariant mass squared value around which the NR approximation is obtained. Namely, the optimal values of β i around which one expands are slightly altered, cf. Eqs. (48) to (51).
The shifts spoil the rotational symmetry of the secondorder Coulomb and harmonic oscillator potentials. For b and c quarks, the deviation from spherical symmetry appears to be on the order of a few percent. It depends on the gluon mass.
This effect is certainly going to change in calculations of higher order than second, because it depends on the gluon mass ansatz and the ansatz will be replaced by theory. Since this effect is relatively small, we neglect it in what follows. Apart form the neglected effect, the Coulomb interactions between quarks are not altered when flavors differ.
The harmonic oscillator forces between pairs of quarks depend on the quark masses. Instead of Eq. (60), in which a common coefficient w is omitted, one obtains where and w 31 = w 23 . Using the reduced masses µ 12 and µ 3 (12) , one can write the frequencies squared for the Jacobi oscillation modes 12 and 3(12) as The frequencies strongly depend on the quark masses and the RGPEP scale parameter λ. The latter needs to be adjusted to justify the assumption that the component |3Q t in the baryon eigenstate is dominant. This is unlikely in the case when λ is much greater than the heavy quark mass and additional components are significantly involved in the dynamics. One can solve our effective eigenvalue problem for the component |3Q t and thus obtain approximate baryon masses and wave functions when λ is optimized by demanding agreement with data for quarkonia.
The effective eigenvalue equation for heavy baryons in QCD with two heavy flavors, implied by our gluon mass hypothesis, is where ∆ denotes Laplacian, reduced masses are µ 12 = m 1 /2, µ 3(12) = 2m 1 m 3 /(2m 1 + m 3 ) and The baryon mass eigenvalue is obtained from the eigenvalue E, We omitted the BF spin-dependent terms and the RG-PEP form factors whose numerical inclusion requires the fourth-order RGPEP calculation. So, Eq. (66) only accounts for interactions of order α. The associated quarkonium eigenvalue equation is [3] where ω 2 12 is given in Eq. (A8).

VII. SKETCH OF TRIPLY HEAVY BARYON SPECTRA
Our sketch of the heavy baryon spectra is drawn using parameters α, m c and m b . To estimate these quantities we use data for heavy quarkonia. Our premise that m Q λ k B ∼ α µ, where k B is the strong Bohr momentum and µ is the quark reduced mass, is formally satisfied in the limit of infinitesimal strong coupling constant by setting λ ∼ √ α m Q , where m Q is a suitable quark mass parameter. Results of our weak coupling expansion are extrapolated to the values of α between about one quarter and one half, depending on the system under consideration and assuming that the functional form of effective Hamiltonian as a function of α, does not change significantly [15].
The numbers we obtain are listed including four or even five significant digits only because the data we approximately reproduce [4] provide that many digits. We ignore data error bars and use our analytic expressions. The Coulomb effects are estimated in first-order perturbation theory around the relevant oscillator solutions. Spins are ignored, but we do include effects due to the Pauli exclusion principle for fermions.

A.
Adjustment of parameters α, m b , mc and λ The coupling constant dependence on λ is set to the well-known approximate function, cf. Ref. [10], The quark masses are assumed to be independent of λ because their dependence is not known yet in the RG-PEP. Confinement poses a conceptual difficulty concerning the definition of quark mass [4]. Quantitative estimates of quark masses would need the RGPEP calculation to at least fourth order while we consider only second. Formulas for running masses of quarks in other approaches, like in Eq. (9.6) in [4], do not concern mass terms in the FF Hamiltonian H t . At the current level of crude approximation and not knowing the masses precisely, we assume that m b and m c can be treated as constants in the range of values of λ that we use in fitting data, see below.
If our calculations of H t and its eigenvalues were exact, observables we obtain would be independent of λ, which hence could be chosen arbitrarily. Since we solve the RGPEP equation only up to order α and we introduce a gluon mass ansatz to reduce the eigenvalue problem to the hadron dominant Fock component, the effective dynamics we obtain may provide a reasonable approximation only in a certain window of values of λ [15].
Note that keeping λ proportional to the square root of α secures proportionality of the resulting hadron binding energies to α 2 [16], which resembles analogous scaling in QED. This scaling is maintained with our oscillator terms because their frequencies emerge proportional to α 2 . Knowing that the observed quarkonium mass spectra can be characterized as intermediate between the Coulomb and oscillator spectra [5], we expect that the harmonic oscillator frequencies obtained from QCD may be comparable in size with the strong-interaction Rydberg-like constant R = µ(4α/3) 2 /2.
In the case of quarkonia, we set wherem QQ is the average mass of quark and antiquark that form a heavy meson, such as J/ψ, Υ or B + c . The quark masses and unknown values of a QQ and b QQ are fitted to the spectra of heavy quarkonia. Separate fits for a set of bb states and a set of cc states give us most suitable λ bb and λ cc , and quark masses m b and m c . This set of numbers allows us to fix values of a QQ and b QQ in the linear formula of Eq. (72). With a QQ and b QQ fixed, we test Eq. (72) by comparing our theoretical spectrum of B c particles, computed for λ bc given by Eq. (72), with experimental data. The agreement is satisfactory, cf. Sec. VII B. The adjustment of constants a QQ and b QQ reflects the current lack of knowledge of the values of λ at which one can most accurately approximate different hadron eigenvalue problems using merely their lowest Fock components and gluon mass ansatz. Details of our fits of two quark masses, m c and m b at most suitable values of λ bb and λ cc , are described in App. A.
In the case of baryons, we set wherem 3Q is the average mass of the three quarks that form a lowest Fock component of a baryon at scale λ 3Q . We set a 3Q = a QQ , b 3Q = b QQ because it is the simplest possible choice that secures λ bbb = λ bb and λ ccc = λ cc , which is a natural choice. It turns out to yield estimates for bbb and ccc spectra that resemble results of other approaches, see below.
Our estimates are quite crude. We ask two questions. One is if the oscillator terms that follow from the assumption of gluon mass are capable of providing a reasonable first approximation to heavy hadrons. Provided that in the case of heavy quarkonia the answer is yes, the other question is what character of the heavy baryons spectrum one expects using the assumption that effective gluons develop a mass. To address these qualitative questions, we ignore the BF spin-dependent terms and we estimate strong-Coulomb effects by evaluating expectation values of the corresponding interaction terms in the oscillator eigenstates. Details of our baryon wave functions are described in App. D. Comparison with other approaches, including lattice estimates, suggests that our extremely simple oscillator picture and thus possibly also the gluon mass hypothesis, appear reasonable. Reliable estimates of better accuracy require fourth-order solution to the RGPEP Eq. (16). and these values are associated with α(λ cc ) = 0.3926 and ω cc = 321.6 MeV. The resulting charmonium masses are illustrated in the right panel of Fig. 3.
Values of λ bb , λ cc , and m b , m c allow us to fix These coefficients imply, according to Eq. (72), The middle panel of Fig. 3 shows the comparison of experimental masses of B c and B c (2S) and an average of different predictions for a mass of B * c with our theoretical levels. Note that because we fit the masses of spin-one quarkonia while we neglect spin-dependent interactions, we present our mass estimates as for 1 − . Because there are no experimental data to compare for spin-one B * c , we provide an average of various theoretical predictions [17]. The agreement is satisfactory, given that we do not expect our estimates to be precise. Equation (78) gives C.

Estimates of masses of heavy baryons
The fit to quarkonia described in Sec. VII B establishes optimal values of λ for all baryons. Values of the coupling constant are obtained from Eq. (71). The optimal values we obtain for these parameters are listed in App. B. The resulting masses of heavy baryons are shown in Fig. 4. Labels of states describe internal orbital motion of quarks, where the first part of a label corresponds to the motion of quark 1 with respect to quark 2 and the second part corresponds to the motion of quark 3 with respect to the pair of quarks 1 and 2. For example, in the state 1P 1S, the pair 12 is in a p-wave without radial excitation, while the quark 3 in its motion with respect to the pair 12 is in an s-wave state without radial excitation. In 1S2S, both 1 with respect to 2 and 3 with respect to 12 are in an s-wave state but the latter is radially excited. States A, B, C and D in ccc and bbb correspond to the second excitation of harmonic oscillator with excitation energy 2ω above the ground state. These states have spin-momentum wave functions that are symmetrized in a way due for fermions in colorless states. Details of the harmonic oscillator basis wave functions are described in App. D. Analytical formulas for masses of baryons are given in App. E.
The values of masses we obtain for bbb and ccc baryons agree well with model calculations [18][19][20][21][22][23][24] including quark-diquark [25] and hypercentral approxima- tions [26,27], bag models [28][29][30], Regge phenomenology [31,32], sum rules [33][34][35][36], pNRQCD [37], Dyson-Schwinger approach [38,39] and lattice studies [40][41][42][43][44][45][46], where comparison is available. As an example of comparison, we note that the ground state of ccc is assigned masses from 4733 MeV to 4796 MeV, by different lattice calculations, with an average of 4768 MeV. Our result is 4797 MeV, differing by 29 MeV, or 0.6 % from the average. For bbb, the average of two lattice results we have identified is 14369 MeV, and our result is 14346 MeV, which is 23 MeV difference, or 0.2 %. These comparisons refer to Table I in Ref. [27] that summarizes results of calculations of masses of Ω ccc and Ω bbb reported in twenty different articles. Ground state of bbc is also very close to the lattice result [40]. In contrast, our ccb differs by about 300 MeV from the lattice. We comment on this feature below. Comparison with lattice calculations reported in Ref. [41] shows that our splittings in bbb differ only by about 10 %. In case of ccc [43], the difference of splittings does not exceed 20 %. This degree of agreement is surprising in view of the complexity of lattice calculations in comparison with the simplicity of our effective Hamiltonian calculation.
The prominent feature visible in Fig. 4 is the extraordinary magnitude of splittings in the ccb baryons. It is a consequence of large oscillator frequency in cc subsystem in Eq. (64), due to large ratio of λ ccb /m c , in which the scale parameter λ ccb is large in comparison with m c due to m b . This separation of scales may make precise calcu-lations of masses of excited ccb baryons difficult. Since the harmonic excitation is so high, it is likely that components with gluons of mass on the order of 1 GeV have to be included in a nonperturbative way.
The surprising feature that the very crude, first approximation based on the RGPEP, with no free parameters left after adjusting quark masses and scale to bb and cc data, produces in an elementary analytic way similar splittings to the ones resulting from advanced calculations, is further illustrated in Fig. 5. It presents splittings in a second band of harmonic oscillator caused by Coulomb interactions. Splittings m D − m C , m C − m B , m B − m A are in relation 2:1:5, which is the general result in the first order of perturbation theory for harmonic oscillator perturbed by any potential [47].
Interestingly, analogous lattice QCD splittings with spin dependent interactions turned off [41], also appear in the ratios 2:1:5. These results suggest that the RG-PEP constituent picture with a gluon mass ansatz may be grasping the physics of lowest-mass heavy baryons. Since experimentally triply heavy baryons are difficult to produce and detect [48], their theoretical understanding using standard techniques is weakly motivated and hence also limited [49,50]. Therefore, the ease with which our method yields results for heavy baryons in agreement with complex approaches suggests that application of the RGPEP in fourth order and including components with one or more effective gluons in the eigenvalue problem beyond perturbation theory, are worth attempting. Qualitative picture of triply heavy baryon mass spectrum implied by the second-order RGPEP in heavy-flavor QCD and our gluon mass ansatz. The figure shows excitations above the ground states 1S1S, whose absolute masses are written at the bottom of each column. The ccb spectrum displays extraordinarily large mass excitation for states 1P 1S, 2S1S and 1D1S in ccb. Such high excitations are associated with formation of cc-diquarks, bound by a harmonic force that is strong because the charmed quarks are much lighter than the bottom quarks. Much less pronounced splittings appear in the bbc baryons. See the text for further discussion.

VIII. CONCLUSION
The effective Hamiltonians we finesse for heavy quarkonia and baryons from QCD of charm and beauty quarks using our gluon mass ansatz, lead to the baryon mass spectra in the ball park of expectations from other approaches to physics of ccc and bbb systems. In addition, the Hamiltonians suggest that quarks c form tight diquarks in ccb baryons. Diquarks are less likely in bbc baryons. Other approaches do not foresee tight diquarks in ccb baryons. This feature may thus distinguish a physically proper approach in future. However, such tight diquarks are hard to excite and mass splittings due their excitation are comparable or even exceed values of the gluon mass one may expect in theory. In that case, the highly excited baryon component with a heavy gluon may be large and our approximation to the three-quark compo-  Fig. 4. The harmonic oscillator and Coulomb potential are sufficient to qualitatively reproduce the pattern of splittings presented in Fig. 11 in Ref. [41]. However, quantitatively our pilotstudy splittings are about twice smaller. Precise agreement is not expected because our approximation is very crude and does not include spin interactions. Hence, it cannot match splittings between the C states and between the D states obtained in other calculations that include spin effects.
nent as dominant may be invalid. Calculations that treat the highly excited baryons as having significant components with one heavy effective gluon may yield smaller masses than our approximation based on the dominance of the three-quark component. If it were the case, the RGPEP approach would still apply, but in the domain of hadron physics in which gluons appear as constituents in competition with quarks for probability of appearance. Taking into account that the method of RGPEP that we use is invariant under boosts and that it is a priori capable of providing a relativistic theory of hadrons in terms of a limited number of their effective constituents with suitably adjusted size, an extension of the RGPEP calculation to fourth order appears worth undertaking. It is certainly needed for verifying if the gluon mass ansatz we introduced provides an adequate representation of dynamics of gluons in the presence of heavy color sources. Fourth-order Hamiltonian is also needed for control on the spin splittings and rotational symmetry.
The ratio 8/6 of harmonic oscillator frequencies in heavy quarkonia and triply heavy baryons is close to the ratio 8/5 obtained for u and d constituent quarks in models using the concept of gluon condensate. If this is not accidental, one may hope that the RGPEP formalism shall apply also to light hadrons as built from constituent quarks and massive gluons, the latter nearly decoupled after generating effective interactions for quarks on the way down in λ toward 1/fm [15]. But even for heavy baryons alone, the effective oscillator picture provides simple wave functions that can be used in description of relativistic processes that involve heavy hadrons.

ACKNOWLEDGMENTS
María Gómez-Rocha acknowledges financial support from European Centre for Theoretical Studies in Nuclear Physics and Related Areas, and from MINECO FPA2016-75654-C2-1-P. Stan G lazek acknowledges support of ECT* during his visit in February 2017.

Appendix A: Fits to masses of well-established quarkonia
Once the coupling constant α as a function of λ is set, the eigenvalues of Eq. (69) estimated by evaluating expectation values of the Coulomb terms in known eigenstates of the oscillator part of the effective Hamiltonian, with and m = m c or m = m b , are used to evaluate corresponding masses from the formula These are compared with data [4]. Thus, Υ(1S), Υ(2S) and χ b1 (1P ) are used to find best values m b and λ bb , using χ 2 . We obtain, Using α(λ bb ) = 0.2664, we find that ω bb = 268.8 MeV.
For charmonia, J/ψ, ψ(2S) and χ c1 (1P ) masses are used in the same way to find best values of m c and λ cc , which turn out to be Using α(λ cc ) = 0.3926, we obtain ω cc = 321.6 MeV. Two heavy mesons made of quarks b and c were observed, B c and B c (2S) [4]. Fits to quarkonia fix masses of quarks, hence, the only free parameter is λ bc , which we fix by assuming Eq. (72). Without any freedom left we plot in Fig. (3) the spectrum of B c using the same Eqs. (A1) and (A2), but with The physical quarkonia masses are read from On Fig. 6, there are also two green triangles at the bottom of the plot. They indicate the values of m b and m c . Furthermore, Figure 6 presents also α(λ) given in Eq. (71)  Harmonic oscillator basis for baryons is constructed from products of wave functions of two harmonic oscillators associated with relative motion of particles 1 and 2 (with momentum K 12 ), and relative motion of particle 3 of momentum Q 3 with respect to the pair 12.
Eigenfunctions for relative motion of two particles with reduced mass µ interacting with harmonic oscillator force characterized by frequency ω are where ν = 1/(2µω), k is the radial excitation number, l is the orbital angular momentum number of the state, m is the projection of angular momentum on the z axis, Y lm are spherical harmonics [4] and L a k (x) are generalized Laguerre polynomials. The first two polynomials are L a 0 (x) = 1, L a 1 (x) = 1 + a − x. The normalization factors are so that d 3 p (2π) 3 ψ * klm ψ k l m = δ kk δ ll δ mm . Finally, the energies are We use a convenient notation for products of wave functions of two harmonic oscillators, (k 12 + 1)(l 12 ) m12 (k 3(12) + 1)(l 3 (12) ) m 3 (12) = |ψ k12l12m12 ψ k 3(12) l 3(12) m 3 (12) , where index 12 corresponds to harmonic oscillator between 1 and 2, with ν 12 = 1/(2µ 12 ω 12 ), and index 3 (12) corresponds to harmonic oscillator between 3 and 12, with ν 3(12) = 1/(2µ 3(12) ω 3 (12) ). For example, the ground state is |1S 0 1S 0 ≡ |1S1S , while |1P 1 2S 0 ≡ |1P 1 2S is the state with harmonic oscillator between 1 and 2 excited to the first orbital excitation with angular momentum projection on z-axis equal 1 and the harmonic oscillator between 3 and 12 being radially excited. The quantum numbers m 12 and m 3 (12) are omitted below, unless they are relevant. We consider states whose excitation energies are at most 2ω 12 or 2ω 3(12) or ω 12 + ω 3 (12) . That is, we consider states 1S1S, 1P 1S, 1S1P , 2S1S, 1S2S, 1D1S, 1S1D and 1P 1P . The total angular momentum of a baryon is conserved, therefore, each state of the basis should have definite orbital angular momentum L. Since states 1P m12 1P m 3 (12) do not have definite angular momentum, we introduce instead the following states with angular momenta L = 2, 1 and 0, respectively, |2, +2 = |1P 1 1P 1 , |2, +1 = . . . (D5) where only the highest L z state is written explicitly. The construction of these states is done in accordance with the rules of adding angular momenta in quantum mechanics and we use convention defined in [4] in the tables of Clebsch-Gordan coefficients. Because quarks have spin, we also need to construct the spin wave functions. We define spin-3/2 quadruplet, which is fully symmetric with respect to exchange of any pair of quarks, and another spin-1/2 doublet, which we call (1/2) A and which is 12-antisymmetric, For completeness, we describe our construction of states of definite total angular momentum. First consider ccb and bbc systems, where quarks 1 and 2 are identical and 3 is different. Because quarks 1 and 2 are identical, the total spin-momentum wave function has to be 12symmetric; color-singlet wave function is antisymmetric. Therefore, one must add orbital angular momentum and spin respecting the Pauli exclusion principle for fermions. For example, the states |1P 1S are 12-antisymmetric and to obtain 12-symmetric spin-momentum wave function we can combine them only with spin (1/2) A , which is also 12-antisymmetric. States |1S1P are 12-symmetric and we can combine them only with spin 3/2 and (1/2) S . The list of possible states is summarized in Table I. To obtain the explicit formulas for the wave functions, we use Clebsch-Gordan tables, as in Eqs. (D5) to (D7).

Symmetric wave functions
In the case of three identical quarks, we need to use fully symmetric wave functions. One can symmetrize the wave functions given above. The ground state 1S1S wave function is fully symmetric in momentum, and we can combine it only with spin 3/2. By the way, symmetrization of 1S1S with (1/2) S gives zero. In this case the wave function is the same as in the case of only two quarks being identical, 0ω, 3 2 where J z = +3/2, +1/2, −1/2, −3/2 is the projection of baryon spin on z-axis.
After symmetrization of the oscillator once-excited states 1P 1S and 1S1P , one is left with only two linearly independent multiplets of states, whose wave functions are (we write only the highest J z state in each multiplet, more information is available in Table II) x = 4ν 3(12) ν 12 .
We list the Coulomb interaction expectation values for ccb and bbc states in Fig. 4. For baryons ccc and bbb, we also make use of Eqs. (E1), (E3) and (E4). Formulas for ground states and once orbitally excited states do not change. For identical quarks x = 3, andṼ Energies for states A, B, C and D need to be evaluated separately. We haveṼ V D = 43 20 . (E20)