**6. The Nature of** *Y***(4260)**

The distinctive characteristics of the *Y*(4260) are its mass, which does not fit any charmonium in the same mass region, the suppression of open charm decays with respect to the *J*/*ψπ*+*π*<sup>−</sup> final state, and the appearance of the exotic charmonium *Zc*(3900) among its decay products. This interesting mix of properties is addressed in quite a few theoretical works, and like in other cases, the molecular, tetraquark, and several other explanations are invoked.

A support for the molecular picture was provided by the QCD lattice computations in [185], by QCD sum rules in [186], by a meson exchange model in [187], and also by the authors of [188], which favored it over the hadro-charmonium interpretation. Further arguments for *Y*(4260) being a molecule were based on the line shape study in [189], and the authors of [190] proposed an unconventional state with a large, but not completely dominant molecular component. An interesting paper [191] came up with a baryonic molecule concept, and the molecular hypothesis was also analyzed in [192–195].

On the contrary, the molecular scenario is strongly disfavored in [196] because of reasons related to the heavy quark spin symmetry and the molecular scenario was rejected in [197] in favor of a charmonium hybrid one. Here, the crux of the argument lies in an important separation between *Y*(4260) mass and its decay threshold. Further arguments to support the charmonium or hybrid-charmonium picture were given in the publications [198–200].

One should also mention different quark models [201–204] with some of them favoring the tetraquark description of *Y*(4260). The tetraquark hypothesis was also analyzed in the QCD sum rules study [205], and the coupled channels approach combined with the three-particle Faddeev equations was used to describe *Y*(4260) in [206].

The analysis of *Y*(4260) is within the CCQM [207] done in a similar way to the *Zc* case: its decay modes are analyzed in both the molecular and tetraquark scenario. With quantitative measurements related to *Y*(4260) not being very numerous, one can analyze the partial decay widths to *J*/*ψπ*+*π*<sup>−</sup> and open charm final states and see whether the latter ones are suppressed. The Feynman diagrams describing the studied transitions are drawn in Figure 7. The considered open charm final states include *DD*¯ , *DD*¯ <sup>∗</sup>, *D*∗*D*¯ , and *D*∗*D*¯ <sup>∗</sup>. As follows from the previous section, *Zc*(3900) is described as a molecular state (67).

**Figure 7.** Feynman diagrams of the *Y*(4260) decay to open charm (**a**) and *Zcπ* (**b**).

The molecular-type non-local interpolating current for *Y*(4260) is written as:

$$J\_{\rm Y^{\rm mol}}^{\mu}(\mathbf{x}) = \int d\mathbf{x}\_1 \dots \int d\mathbf{x}\_4 \,\delta\left(\mathbf{x} - \sum\_{i=1}^4 w\_i \mathbf{x}\_i\right) \Phi\_Y\left(\sum\_{i 
$$J\_{\rm Y^{\rm mol}; \mathbf{4q}}^{\mu} = \quad \frac{1}{\sqrt{2}} \left\{ (\vec{q}(\mathbf{x}\_3)\gamma\_5 \mathbf{c}(\mathbf{x}\_1)) \cdot (\mathcal{E}(\mathbf{x}\_2)\gamma^\mu \gamma\_5 \mathbf{q}(\mathbf{x}\_4)) - (\gamma\_5 \leftrightarrow \gamma^\mu \gamma\_5) \right\}\_{} \qquad (q = u, d)$$
$$

with:

$$w\_1 = w\_2 = \frac{m\_\mathcal{E}}{2(m\_q + m\_c)'} \qquad w\_3 = w\_4 = \frac{m\_q}{2(m\_q + m\_c)}$$

The matrix element corresponding to the open charm production is given by:

$$\mathcal{M}\left(Y\_{\mathfrak{u}}(p,\varepsilon\_{p}^{\mu}) \to D\_{1}^{0}(p\_{1}) + D\_{2}^{0}(p\_{2})\right) = \frac{9}{\sqrt{2}} \operatorname{\mathcal{Y}} \operatorname{\mathcal{Y}} \mathcal{S}\_{D\_{1}} \mathcal{Y} D\_{2} \tag{76}$$

$$\times \quad \int \frac{d^{4}k\_{1}}{(2\pi)^{4}i} \int \frac{d^{4}k\_{2}}{(2\pi)^{4}i} \tilde{\Phi}\_{Y}\left(-\Omega\_{q}^{2}\right) \tilde{\Phi}\_{D\_{1}}\left(-\ell\_{1}^{2}\right) \tilde{\Phi}\_{D\_{2}}\left(-\ell\_{2}^{2}\right)$$

$$\times \quad \left\{ \operatorname{tr}\left(\gamma\_{5} S\_{\mathfrak{c}}(k\_{1}) \Gamma\_{2} S\_{\mathfrak{u}}(k\_{3})\right) \operatorname{tr}\left(\gamma^{\mu} \gamma\_{5} S\_{\mathfrak{u}}(k\_{2}) \Gamma\_{1} S\_{\mathfrak{c}}(k\_{4})\right) - \left(\gamma\_{5} \leftrightarrow \gamma^{\mu}\gamma\_{5}\right) \right\},$$

where:

$$
\Gamma\_1 \otimes \Gamma\_2 = \begin{cases}
\gamma\_5 \otimes \gamma\_5 & \text{for } D\vec{D} \\
\epsilon\_{\nu\_1}^\* \gamma^{\nu\_1} \otimes \gamma\_5 & \text{for } D^\*\vec{D} \\
\epsilon\_{\nu\_1}^\* \gamma^{\nu\_1} \otimes \epsilon\_{\nu\_2}^\* \gamma^{\nu\_2} & \text{for } D^\*\vec{D}^\*
\end{cases}
\tag{77}$$

.

and the momenta are defined as:

$$\begin{aligned} \Omega\_{\eta}^{2} &= \begin{array}{c} \frac{1}{2} \sum\_{i \le j} q\_i q\_j, & q\_1 = -k\_1 - w\_1^{\underline{Y}} p\_{\prime} & q\_2 = k\_4 - w\_2^{\underline{Y}} p\_{\prime} & q\_3 = k\_3 - w\_3^{\underline{Y}} p\_{\prime} \\\\ k\_1 &= -k\_2 + w\_{\mu}^{\underline{D}} p\_{1\prime} & \ell\_2 = -k\_1 - w\_c^{\underline{D}} p\_{2\prime} & k\_3 = k\_1 + p\_{2\prime} & k\_4 = k\_2 + p\_1. \end{aligned} \end{aligned}$$

The decay into *Zc* + *π* involves a three-loop diagram, and the corresponding matrix element is:

$$M\left(Y\_{\boldsymbol{u}}(p,\boldsymbol{\epsilon}^{\mu}) \rightarrow Z\_{\boldsymbol{c}}^{+}(p\_{1},\boldsymbol{\epsilon}^{\nu}) + \pi^{-}\right) = \frac{9}{2} \operatorname{\mathcal{H}\mathcal{S}} \mathcal{Z}\_{\boldsymbol{c}} \mathcal{G}\_{\pi} \tag{78}$$

$$\times \quad \prod\_{j=1}^{3} \left[ \int \frac{d^{4}k\_{j}}{(2\pi)^{4}i} \right] \operatorname{\mathcal{S}}\_{Y}\left(-\Omega\_{\boldsymbol{q}}^{2}\right) \operatorname{\mathcal{S}}\_{\mathcal{Z}\_{\boldsymbol{c}}}\left(-\Omega\_{\boldsymbol{r}}^{2}\right) \operatorname{\mathcal{S}}\_{\pi}\left(-\ell^{2}\right) \varepsilon\_{\mu}(p) \epsilon\_{\nu}^{\*}(p\_{1}) \right. $$

$$\times \quad \times \quad \sum\_{\Gamma} \operatorname{tr}\left(\Gamma\_{1} \operatorname{\mathcal{S}}\_{\boldsymbol{c}}(k\_{1}) \Gamma\_{2} \operatorname{\mathcal{S}}\_{\boldsymbol{u}}(k\_{2})\right) \operatorname{tr}\left(\Gamma\_{3} \operatorname{\mathcal{S}}\_{\boldsymbol{u}}(k\_{3}) \Gamma\_{4} \operatorname{\mathcal{S}}\_{\boldsymbol{d}}(k\_{4}) \Gamma\_{5} \operatorname{\mathcal{S}}\_{\boldsymbol{c}}(k\_{5})\right) \ .$$

where:

$$\sum\_{\Gamma} \left[ \Gamma\_1 \otimes \Gamma\_2 \right] \cdot \left[ \Gamma\_3 \otimes \Gamma\_4 \otimes \Gamma\_5 \right] = \left[ \gamma\_5 \otimes \gamma\_5 \right] \cdot \left[ \gamma^\mu \gamma\_5 \otimes \gamma\_5 \otimes \gamma^\nu \right],$$

$$= \left[ \gamma^\mu \gamma\_5 \otimes \gamma^\nu \right] \cdot \left[ \gamma\_5 \otimes \gamma\_5 \otimes \gamma\_5 \right] - \left[ \gamma^\mu \gamma\_5 \otimes \gamma\_5 \right] \cdot \left[ \gamma\_5 \otimes \gamma\_5 \otimes \gamma^\nu \right].$$

and the momenta are defined as:

$$\begin{aligned} \Omega\_q^2 &= \frac{1}{2} \sum\_{i \le j} q\_i q\_j; \qquad q\_1 = -k\_1 - w\_1^\chi p\_\prime \quad q\_2 = k\_5 - w\_2^\chi p\_\prime \quad q\_3 = k\_2 - w\_3^\chi p\_\prime\\ \Omega\_r^2 &= \frac{1}{2} \sum\_{i \le j} r\_i r\_j; \qquad r\_1 = -k\_5 + w\_1^Z p\_{1\prime} \quad r\_2 = k\_1 + w\_2^Z p\_{1\prime} \quad r\_3 = k\_4 - w\_3^Z p\_{1\prime} \\ \ell\_1 &= \quad k\_3 + w\_4^\pi p\_{2\prime} \quad k\_4 = k\_3 + p\_{2\prime} \quad k\_5 = k\_1 - k\_2 + k\_3 + p. \end{aligned}$$

In the tetraquark scenario, the non-local *Y*(4260) current takes the form:

$$J\_{\mathbf{Y}^{\text{int}}}^{\mu}(\mathbf{x}) = \int d\mathbf{x}\_1 \dots \int d\mathbf{x}\_4 \,\delta\left(\mathbf{x} - \sum\_{i=1}^4 w\_i^{\mathbf{Y}} \mathbf{x}\_i\right) \Phi\_{\mathbf{Y}}\left(\sum\_{i$$

$$J\_{"\gamma^{\rm act};4,q}^{\mu} = \frac{1}{\sqrt{2}} \mathfrak{c}\_{abc} \mathfrak{c}\_{\rm dcc} \left\{ (q\_a(\mathbf{x}\_4) \mathbf{C} \gamma\_5 \mathbf{c}\_b(\mathbf{x}\_1)) (\bar{q}\_d(\mathbf{x}\_3) \gamma^\mu \gamma\_5 \mathbf{C} \bar{\mathbf{c}}\_c(\mathbf{x}\_2)) - (\gamma\_5 \leftrightarrow \gamma^\mu \gamma\_5) \right\}.\tag{80}$$

The matrix element of the decay into *DD*¯ is expressed as:

$$M\left(Y\_{\
u}^{\rm det}(p,\varepsilon\_{p}^{\mu}) \to D\_{1}^{0}(p\_{1}) + D\_{2}^{0}(p\_{2})\right) = \frac{6}{\sqrt{2}} \operatorname{\mathcal{G}}\_{Y} \operatorname{\mathcal{G}}\_{D\_{1}} \operatorname{\mathcal{G}}\_{D\_{2}} \tag{81}$$

$$\times \quad \int \frac{d^{4}k\_{1}}{(2\pi)^{4}i} \int \frac{d^{4}k\_{2}}{(2\pi)^{4}i} \operatorname{\mathcal{\Phi}}\_{Y}\left(-\operatorname{\mathcal{D}}\_{q}^{2}\right) \operatorname{\mathcal{\Phi}}\_{D\_{1}}\left(-\ell\_{1}^{2}\right) \operatorname{\mathcal{\Phi}}\_{D\_{2}}\left(-\ell\_{2}^{2}\right)$$

$$\times \quad \left\{ \operatorname{tr}\left(\gamma\_{5} \operatorname{\mathcal{S}}\_{\operatorname{\mathcal{E}}}(k\_{1}) \operatorname{\mathcal{I}}\_{2}^{0} \operatorname{\mathcal{S}}\_{\operatorname{\mathcal{E}}}(k\_{3}) \gamma^{\mu} \gamma\_{5} \operatorname{\mathcal{S}}\_{\operatorname{\mathcal{E}}}(k\_{2}) \operatorname{\mathcal{I}}\_{1}^{0} \operatorname{\mathcal{S}}\_{\operatorname{\mathcal{U}}}(k\_{4})\right) - (\gamma\_{5} \leftrightarrow \gamma^{\mu} \gamma\_{5}) \right\},$$

with the momenta:

$$\begin{aligned} \Omega\_q^2 &= \frac{1}{2} \sum\_{i \le j} q\_i q\_j; & q\_1 &= -k\_1 - w\_1^\chi p\_\prime & q\_2 &= -k\_2 - w\_2^\chi p\_\prime & q\_3 &= k\_3 - w\_3^\chi p\_\prime\\ \ell\_1 &= -k\_2 - w\_\varepsilon^D p\_1, & \ell\_2 &= -k\_1 - w\_\varepsilon^D p\_2, & k\_3 &= k\_1 + p\_2, & k\_4 &= k\_2 + p\_1. \end{aligned}$$

The matrix element of the decay into *Zcπ* is given by:

$$\begin{split} &\quad \mathcal{M}\left(\boldsymbol{\chi}\_{\boldsymbol{u}}^{\text{att}}(\boldsymbol{p},\boldsymbol{\epsilon}^{\mu}) \to \mathcal{Z}\_{\varepsilon}^{+}(\boldsymbol{p}\_{1},\boldsymbol{\epsilon}^{\mu}) + \boldsymbol{\pi}^{-}(\boldsymbol{p}\_{2})\right) = \boldsymbol{3}\,\boldsymbol{\operatorname{\mathcal{G}}}\boldsymbol{\mathcal{G}}\boldsymbol{\mathcal{Z}}\_{\varepsilon}\boldsymbol{\mathcal{G}}\_{\pi} \\ &\quad \times \quad \prod\_{j=1}^{3} \left[ \int \frac{d^{4}k\_{j}}{(2\pi)^{4}i} \right] \boldsymbol{\widetilde{\Phi}}\_{Y}\left(-\Omega\_{\eta}^{2}\right) \boldsymbol{\widetilde{\Phi}}\_{Z\_{\varepsilon}}\left(-\Omega\_{\tau}^{2}\right) \boldsymbol{\widetilde{\Phi}}\_{\pi}\left(-\boldsymbol{\ell}^{2}\right) \\ &\quad \times \quad \boldsymbol{\varepsilon}\_{\mu}(\boldsymbol{p})\boldsymbol{\varepsilon}\_{\nu}^{+}(\boldsymbol{p}\_{1}) \sum\_{\Gamma} \text{tr}\left[\Gamma\_{1}^{Y}\boldsymbol{S}\_{\varepsilon}(k\_{1})\Gamma\_{2}^{Z}\boldsymbol{S}\_{\mu}(k\_{2})\Gamma\_{2}^{Y}\boldsymbol{S}\_{\varepsilon}(k\_{3})\Gamma\_{1}^{Z}\boldsymbol{S}\_{d}(k\_{4})\gamma\_{5}\boldsymbol{S}\_{\mu}(k\_{5})\right], \end{split} \tag{82}$$

*Symmetry* **2020**, *12*, 884

with the momenta:

$$\begin{aligned} \Omega\_q^2 &= \ \frac{1}{2} \sum\_{i \le j} q\_i q\_j; \qquad q\_1 = -k\_1 - w\_1^\chi p\_\prime \quad q\_2 = -k\_3 - w\_2^\chi p\_\prime \quad q\_3 = k\_2 - w\_3^\chi p\_\prime\\ \Omega\_r^2 &= \ \frac{1}{2} \sum\_{i \le j} r\_i r\_j; \qquad r\_1 = k\_3 + w\_1^Z p\_\prime \quad r\_2 = k\_1 + w\_2^Z p\_\prime \quad r\_3 = -k\_4 + w\_3^Z p\_\prime \end{aligned}$$

$$\ell\_1 = \ -k\_4 - w\_4^\chi p\_{2\prime} \quad k\_4 = k\_1 - k\_2 + k\_3 + p\_{1\prime} \quad k\_5 = k\_1 - k\_2 + k\_3 + p.$$

Here, the summation over Γ is defined by:

$$\sum\_{\Gamma} = \left[\gamma\_5 \otimes \gamma^\mu \gamma\_5 - \gamma^\mu \gamma\_5 \otimes \gamma\_5\right]^Y \otimes \left[\gamma\_5 \otimes \gamma^\nu - \gamma^\nu \otimes \gamma\_5\right]^Z.$$

The considered decays comprise different combinations of pseudoscalar, vector, and axial-vector particles in the final state. The relevant expressions for the matrix elements and decay widths are written down:

$$\begin{aligned} M(V(p) \to P(p\_1) + P(p\_2)) &= \epsilon\_V^\mu q\_\mu G\_{VPP}, \qquad q = p\_1 - p\_2, \\ \Gamma(V \to PP) &= \frac{\left| \mathbf{p} \mathbf{1} \right|^3}{6\pi m^2} G\_{VPP}^2, \end{aligned}$$

$$\begin{aligned} M(V(p) \to A(p\_1) + P(p\_2)) &= \epsilon\_V^\mu \epsilon\_A^{\*\nu} \left( \mathbb{g}\_{\mu\nu} A + p\_1 \Box p\_\nu B \right), \\ \Gamma(V \to AP) &= \frac{|\mathbf{p\_1}|}{24\pi m^2} \left\{ \left( 3 + \frac{|\mathbf{p\_1}|^2}{m\_1^2} \right) A^2 + \frac{m^2}{m\_1^2} |\mathbf{p\_1}|^4 B^2 + \frac{m^2 + m\_1^2 - m\_2^2}{m\_1^2} |\mathbf{p\_1}|^2 AB \right\} \end{aligned}$$

$$\begin{aligned} M(V(p) \to V(p\_1) + P(p\_2)) &= \epsilon\_V^\mu \epsilon\_V^{\* \text{ } \text{\textquotedblleft}} \varepsilon\_{\mu \nu\_1 a \beta} p^a p\_1^\beta G\_{VVP}, \\ \Gamma(V \to VP) &= \frac{\left| \mathbf{p} \mathbf{1} \right|^3}{12 \pi} G\_{VVP}^2 \end{aligned}$$

$$M(V(p)\to V(p\_1)+V(p\_2)) = \epsilon\_V^\mu \epsilon\_V^{\*\nu\_1} \epsilon\_V^{\*\nu\_2} \left\{ p\_{1\mu} p\_{1\nu\_2} p\_{2\cdot\nu\_1} A + g\_{\mu\nu\_1} p\_{1\cdot\nu\_2} B + g\_{\mu\nu\_2} p\_{2\cdot\nu\_1} C + g\_{\nu\_1\nu\_2} p\_{1\cdot\mu} D \right\},$$

$$\begin{split} \Gamma(V \to V\_{1}V\_{2}) &= \frac{|\mathbf{p}\_{1}|^{3}}{24\pi m\_{1}^{2}m\_{2}^{2}} \Big\{ m^{2}|\mathbf{p}\_{1}|^{4}A^{2} + [|\mathbf{p}\_{1}|^{2} - 3m\_{1}^{2}]B^{2} + [|\mathbf{p}\_{1}|^{2} + 3m\_{2}^{2}]\big{C}^{2} \\ &+ [|\mathbf{p}\_{1}|^{2} + 3\frac{m\_{1}^{2}m\_{2}^{2}}{m^{2}}]D^{2} + |\mathbf{p}\_{1}|^{2}[m^{2} + m\_{1}^{2} - m\_{2}^{2}]AB + |\mathbf{p}\_{1}|^{2}[-m^{2} + m\_{1}^{2} - m\_{2}^{2}]A\big{C} \\ &+ [\mathbf{p}\_{1}|^{2}[m^{2} - m\_{1}^{2} - m\_{2}^{2}]AD + [2|\mathbf{p}\_{1}|^{2} - m^{2} + m\_{1}^{2} + m\_{2}^{2}]BC + [2|\mathbf{p}\_{1}|^{2} + m\_{1}^{2} + \frac{m\_{1}^{2}}{m^{2}}(m\_{2}^{2} - m\_{1}^{2})]BD \\ &+ [-2|\mathbf{p}\_{1}|^{2} - m\_{2}^{2} + \frac{m\_{2}^{2}}{m^{2}}(m\_{2}^{2} - m\_{1}^{2})]CD \Big]. \end{split}$$

The value of Λ*Zc* is set to 3.3 GeV, and guided by our experience, we assume that Λ*Y*(4260) = 3.3 ± 0.1 GeV. The numerical evaluation leads to the results presented in Table 4.

In both scenarios, the open charm decays are suppressed with respect to the *J*/*ψπ* decay channel. The discrimination between them is provided by the total decay width Γ[*Y*(4260)] = 55 ± 19 MeV, which is in contradiction with the molecular description. Thus, one can conclude that the CCQM approach favors the tetraquark structure of *Y*(4260).


**Table 4.** Decay widths of the selected *Y*(4260) transition in MeV.

#### **7. Bottomonium-Like States** *Zb***(10610) and** *Z <sup>b</sup>***(10650)**

Exotic quarkonia states appear also in the bottomonium sector: *Zb*(10610) and *Z*- *<sup>b</sup>*(10650) are two examples. Even though the exotic bottomonia masses tend to be significantly higher than the charmonia ones, the underlying dynamics is similar, and one finds the molecular, tetraquark, and other hypotheses in theoretical approaches that describe them.

*Zb*(10610) and *Z*- *<sup>b</sup>*(10650) were seen as molecules in the boson exchange model of [208], and the molecular picture was also favored in [209], where the spin structure of these two particles was analyzed. Further support of the molecular scenario came from the quark model based on a phenomenological Lagrangian used by the authors of [210] and also from other analyses preformed in [211] (QCD multipole expansion), [212] (effective field theory), [213] (pion exchange model), [214] (QCD sum rules, only *Z*- *<sup>b</sup>*(10650) included), [215] (heavy quark spin symmetry and coupled channels analysis),and [216] (coupled channels approach with pion exchange model). A different set of works supports, with various intensity, the tetraquark structure of the two bottomonia states. In [217], the conclusion followed from an effective diquark-antidiquark Hamiltonian combined with meson-loop induced effects. The authors of [218] based their analysis on the QCD sum rules and interpreted *Zb* and *Z*- *<sup>b</sup>* as axial-vector tetraquarks. The two works [219,220] also drew their conclusions from the QCD sum rules and allowed the tetraquark and molecular scenario. The former work suggested that *Zb* and *Z*- *<sup>b</sup>* could have both the diquark-antidiquark and molecular components (following from a mixed interpolating current). The latter one excluded neither the tetraquark nor molecular the interpretation of *Zb*(10610), and the idea of a mixed current appeared also. The mentioned analyses could be supplemented by numerous other works [221–241] where further ideas and approaches were exploited.

The theoretical analysis of the *Zb*(10610) and *Z*- *<sup>b</sup>*(10650) states by the CCQM was performed in [81]. The work assumed a molecular-type interpolating current, which is favored by most theoretical approaches when interpreting the experimental results. It is a natural choice reflecting the proximity of the particle masses to the corresponding thresholds:

$$\begin{aligned} m(Z\_b^+) &= 10607.2 \pm 2.0 \,\mathrm{MeV}, \qquad m(B^\* \bar{B}) &= 10604 \,\mathrm{MeV},\\ m(Z\_b^{'+}) &= 10652.2 \pm 1.5 \,\mathrm{MeV}, \qquad m(B^\* \bar{B}^\*) &= 10649 \,\mathrm{MeV}.\end{aligned}$$

The quantum numbers of the two states *IG*(*JPC*) = 1+(1+−) lead to the choice of (local) interpolating currents:

$$J\_{Z\_b^+}^{\mu} = \begin{array}{ccl} \frac{1}{\sqrt{2}} \left[ (\bar{d}\gamma\_5 b)(\bar{b}\gamma^\mu u) + (\bar{d}\gamma^\mu b)(\bar{b}\gamma\_5 u) \right] \end{array} ,\tag{83}$$

$$J\_{Z\_b^{t+}}^{\mu\nu} = \varepsilon^{\mu\nu a\p} (\bar{d}\gamma\_a b)(\bar{b}\gamma\_\beta \mu),\tag{84}$$

which guarantees that, when considering the transitions into *B*(∗)*B*¯(∗), the *Zb* state can decay only to the [*B*¯ <sup>∗</sup>*B* + *c*.*c*.] pair, while the *Z*- *<sup>b</sup>* state can decay only to a *<sup>B</sup>*¯ <sup>∗</sup>*B*<sup>∗</sup> pair. Decays into the *BB* channels are not allowed.

Further decay channels include a bottomonium particle accompanied with a charged light meson. Taking into account the *G* parity, which is conserved in strong interactions and kinematic considerations, only three possible bottomonium-meson decay channels are available: *Z*<sup>+</sup> *<sup>b</sup>* <sup>→</sup> <sup>Υ</sup> <sup>+</sup> *<sup>π</sup>*+,

*Z*<sup>+</sup> *<sup>b</sup>* <sup>→</sup> *hb* <sup>+</sup> *<sup>π</sup>*<sup>+</sup> and *<sup>Z</sup>*<sup>+</sup> *<sup>b</sup>* <sup>→</sup> *<sup>η</sup><sup>b</sup>* <sup>+</sup> *<sup>ρ</sup>*+. All mentioned *<sup>Z</sup>*<sup>+</sup> *<sup>b</sup>* transition can be arranged into three groups with respect to the spin kinematics:

$$\begin{aligned} 1^+ &\to 1^- + 0^- : & Z\_b^+ &\to \mathcal{Y} + \pi^+ \,, \quad Z\_b^+ &\to [\mathcal{B}^{\*0}B^+ + \text{c.c.}] \,, \quad Z\_b^+ &\to \eta\_b + \rho^+ \,, \\ 1^+ &\to 1^+ + 0^- : & Z\_b^+ &\to h\_b + \pi^+ \,, \\ 1^+ &\to 1^- + 1^- : & Z\_b^+ &\to \mathcal{B}^{\*0}B^{\*+} \,. \end{aligned}$$

The classification of the bottomonia particles based on their quantum numbers is shown in Table 5.


**Table 5.** The bottomonium states <sup>2</sup>*S*+1*<sup>L</sup> <sup>J</sup>*. We use the notation <sup>↔</sup> *<sup>∂</sup>* <sup>=</sup><sup>→</sup> *<sup>∂</sup>* <sup>−</sup> <sup>←</sup> *∂* .

The expressions for matrix elements and decay widths depend on the spin structure and are for the three cases as follows.

• For 1<sup>+</sup> <sup>→</sup> <sup>1</sup><sup>−</sup> <sup>+</sup> <sup>0</sup><sup>−</sup> transitions, the matrix element can be parameterized with two Lorentz structures:

$$<\langle 1^-(q\_1;\delta), 0^-(q\_2)|\, T|1^+(p;\mu)\rangle = \left(A \, g^{\mu\delta} + B \, q\_1^{\mu} q\_2^{\delta}\right) \varepsilon\_{\mu} \, \varepsilon\_{1\delta}^\*.\tag{85}$$

The invariant amplitudes *A* and *B* can be combined into the helicity amplitudes:

$$H\_{00} = -\frac{E\_1}{M\_1}A - \frac{M}{M\_1}|\mathbf{q}\_1|^2 B\_{\prime} \quad H\_{+1+1} = H\_{-1-1} = -A\_{\prime}$$

which are practical to express the decay width. For the derivation of the latter, it is useful to work in the rest frame of the initial particle, where <sup>|</sup>**q1**<sup>|</sup> <sup>=</sup> *<sup>λ</sup>*1/2(*M*2, *<sup>M</sup>*<sup>2</sup> <sup>1</sup>, *<sup>M</sup>*<sup>2</sup> <sup>2</sup>)/2*M* is the three-momentum and *E*<sup>1</sup> = (*M*<sup>2</sup> + *M*<sup>2</sup> <sup>1</sup> <sup>−</sup> *<sup>M</sup>*<sup>2</sup> <sup>2</sup>)/2*M* is the energy of the final state vector. Furthermore, the on-mass-shell character of the initial and final state particles is taken into account by *p*<sup>2</sup> = *M*2, *q*<sup>2</sup> <sup>1</sup> = *<sup>M</sup>*<sup>2</sup> <sup>1</sup>, *<sup>q</sup>*<sup>2</sup> <sup>2</sup> = *<sup>M</sup>*<sup>2</sup> <sup>2</sup>, and *<sup>p</sup>μεμ* = 0. One arrives at:

$$\Gamma = \frac{|\mathbf{q}\_1|}{24\pi M^2} \left\{ |H\_{+1+1}|^2 + |H\_{-1-1}|^2 + |H\_{00}|^2 \right\}. \tag{86}$$

• The matrix element for the 1<sup>+</sup> <sup>→</sup> <sup>1</sup><sup>+</sup> <sup>+</sup> <sup>0</sup><sup>−</sup> transitions is expressed through one covariant term only:

$$\langle \langle 1^+(q\_1;\delta), 0^-(q\_2) \vert T \vert 1^+(p;\mu) \rangle = \mathbb{C} \, q\_{1a} q\_{2\rho} \, \varepsilon^{a\beta\mu\delta} \, \varepsilon\_\mu \, \varepsilon\_{1\delta}^\* \,. \tag{87}$$

The decay with the formula can be written as:

$$
\Gamma = \frac{|\mathbf{q}\_1|^3}{12\pi M^2} \,\mathrm{C}^2,\tag{88}
$$

where one can note the *p*-wave suppression factor |**q1**| **3**.

• As shown in [81], the matrix element for 1<sup>+</sup> <sup>→</sup> <sup>1</sup><sup>−</sup> <sup>+</sup> <sup>1</sup><sup>−</sup> decay can be parameterized using three amplitudes:

$$\langle 1^-(q\_1;\delta), 1^-(q\_2;\rho) \vert \vert T \vert 1^+(p;\mu) \rangle = \left( B\_1 \varepsilon^{q\_1 q\_2 \rho \delta} q\_1^\mu + B\_2 \varepsilon^{q\_1 \mu \rho \delta} + B\_3 \varepsilon^{q\_2 \mu \rho \delta} \right) \varepsilon\_\mu \varepsilon\_\delta \varepsilon\_\rho. \tag{89}$$

The relation between the helicity amplitudes *Hλ*;*λ*1*λ*<sup>2</sup> (*λ* = *λ*<sup>1</sup> − *λ*2) and the invariant amplitudes can be shown to be:

$$\begin{aligned} H\_{0;+1;1} &= \ -H\_{0;-1;-1} = \ -E\_1 A\_1 - E\_2 A\_2 - M |\mathbf{q}\_1|^2 A\_5, \\ H\_{+1;+1;0} &= \ -H\_{-1;-1;0} = \frac{(E\_1 M - M\_1^2)}{M\_2} A\_1 + M\_2 A\_2 - \frac{M^2}{M\_2} |\mathbf{q}\_1|^2 A\_4, \\ H\_{-1;0;+1} &= \ -H\_{+1;0;-1} = M\_1 A\_1 + \frac{(E\_1 M - M\_1^2)}{M\_1} A\_2 - \frac{M^2}{M\_1} |\mathbf{q}\_1|^2 A\_3. \end{aligned}$$

The rate of the decay 1+(*p*) <sup>→</sup> <sup>1</sup>−(*q*1) + <sup>1</sup>−(*q*2), finally, reads:

$$\Gamma = \frac{|\mathbf{q}\_1|}{24\pi M^2} \cdot 2 \left\{ |H\_{0;+1+1}|^2 + |H\_{+1;+1}|^2 + |H\_{-1;0+1}|^2 \right\}.\tag{91}$$

Coming back to the CCQM description, one can write the non-local versions of Equations (83) and (84) as follows:

$$J\_{Z\_b^\mu}^\mu(\mathbf{x}) = \int d\mathbf{x}\_1 \dots \int d\mathbf{x}\_4 \,\delta\left(\mathbf{x} - \sum\_{i=1}^4 w\_i \mathbf{x}\_i\right) \Phi\_{Z\_b;Z\_b;Z\_b} \left(\sum\_{i 
$$J\_{Z\_b;4q}^\mu = \frac{1}{\sqrt{2}} \left\{ (\bar{d}(\mathbf{x}\_3)\gamma\_5 b(\mathbf{x}\_1))(\bar{b}(\mathbf{x}\_2)\gamma^\mu u(\mathbf{x}\_4)) + (\bar{d}(\mathbf{x}\_3)\gamma^\mu b(\mathbf{x}\_1))(\bar{b}(\mathbf{x}\_2)\gamma\_5 u(\mathbf{x}\_4)) \right\},$$
$$

$$\begin{array}{rcl}J\_{Z\_{b}^{\mu\nu}}^{\mu\nu}(\mathbf{x})&=&\int d\mathbf{x}\_{1}\ldots\int d\mathbf{x}\_{4}\,\delta\left(\mathbf{x}-\sum\_{i=1}^{4}w\_{i}\mathbf{x}\_{i}\right)\Phi\_{Z\_{b}}\left(\sum\_{i$$

The interaction Lagrangian is constructed in the usual way for *Zb*; in the case of *Z*- *<sup>b</sup>*, the stress tensor of the field is introduced *Z*- *<sup>b</sup>*, *μν* = *∂μZ*- *<sup>b</sup>*, *<sup>ν</sup>* − *∂νZ*- *<sup>b</sup>*, *<sup>μ</sup>*:

$$\mathcal{L}\_{\text{int},\text{Z}\_{\flat}} = \mathcal{g}\_{\text{Z}\_{\flat}} Z\_{\flat,\mu}(\mathbf{x}) \cdot \mathcal{J}\_{\text{Z}\_{\flat}}^{\mu}(\mathbf{x}) + \text{H.c.},\tag{94}$$

$$\mathcal{L}\_{\text{int},\mathcal{Z}\_b'} = \begin{array}{c} \mathbb{S}\mathcal{Z}\_b' \\ \hline 2M\_{\mathcal{Z}\_b'} \end{array} \\ Z\_{b,\mu\nu}'(\mathbf{x}) \cdot \mathbb{J}\_{\mathcal{Z}\_b'}^{\mu\nu}(\mathbf{x}) + \text{H.c.} \tag{95}$$

The factor 2*MZ*- *<sup>b</sup>* is put into the denominator in order to preserve the same physical dimensions of the *gZb* and *gZ*- *<sup>b</sup>* couplings. The link between these couplings and the size parameters is done via the compositeness condition, which is based on the evaluation of hadronic mass operators. The latter are written in the momentum space as:

$$\begin{split} \hat{\Pi}\_{Z\_{b}}^{\mu\nu}(p) &= \ \frac{9}{2} \prod\_{i=1}^{3} \int \frac{d^{4}k\_{i}}{(2\pi)^{4}i} \,\hat{\Phi}\_{Z\_{b}}^{2} \left(-\vec{\omega}^{2}\right) \\ &\times \quad \left\{ \text{tr}\left[\gamma\_{5}\mathbb{S}\_{1}(\hat{k}\_{1})\gamma\_{5}\mathbb{S}\_{3}(\hat{k}\_{3})\right] \text{tr}\left[\gamma^{\mu}\mathbb{S}\_{4}(\hat{k}\_{4})\gamma^{\nu}\mathbb{S}\_{2}(\hat{k}\_{2})\right] \\ &\quad + \text{tr}\left[\gamma^{\mu}\mathbb{S}\_{1}(\hat{k}\_{1})\gamma^{\nu}\mathbb{S}\_{3}(\hat{k}\_{3})\right] \text{tr}\left[\gamma\_{5}\mathbb{S}\_{4}(\hat{k}\_{4})\gamma\_{5}\mathbb{S}\_{2}(\hat{k}\_{2})\right] \right\}, \end{split} \tag{96}$$

$$\begin{split} \|\hat{\Pi}^{\mu\nu}\_{Z'\_b}(p)\| &=& -9\frac{\varepsilon^{\mu\rho\alpha\beta}\varepsilon^{\nu\rho\rho\sigma}}{M\_{Z'\_b}^2} \prod\_{i=1}^3 \int \frac{d^4k\_i}{(2\pi)^4\hat{\mathbf{i}}} \,\hat{\Phi}^2\_{Z'\_b} \left(-\vec{\omega}^2\right) \\ &\times \quad \text{tr}\left[\gamma\_\beta \mathcal{S}\_1(\hat{k}\_1)\gamma\_\alpha \mathcal{S}\_3(\hat{k}\_3)\right] \text{tr}\left[\gamma\_\beta \mathcal{S}\_4(\hat{k}\_4)\gamma\_\sigma \mathcal{S}\_2(\hat{k}\_2)\right], \end{split} \tag{97}$$

*Symmetry* **2020**, *12*, 884

where *ω* <sup>2</sup> = 1/2 (*k*<sup>2</sup> <sup>1</sup> + *<sup>k</sup>*<sup>2</sup> <sup>2</sup> + *<sup>k</sup>*<sup>2</sup> <sup>3</sup> + *k*1*k*<sup>2</sup> − *k*1*k*<sup>3</sup> − *k*2*k*3) and: ˆ *<sup>k</sup>*<sup>1</sup> <sup>=</sup> *<sup>k</sup>*<sup>1</sup> <sup>−</sup> *<sup>w</sup>*<sup>1</sup> *<sup>p</sup>* , <sup>ˆ</sup> *<sup>k</sup>*<sup>2</sup> <sup>=</sup> *<sup>k</sup>*<sup>2</sup> <sup>−</sup> *<sup>w</sup>*<sup>2</sup> *<sup>p</sup>* , <sup>ˆ</sup> *k*<sup>3</sup> = *k*<sup>3</sup> + *w*<sup>3</sup> *p*, ˆ *k*<sup>4</sup> = *k*<sup>1</sup> + *k*<sup>2</sup> − *k*<sup>3</sup> + *w*<sup>4</sup> *p* , *ε <sup>μ</sup>pαβ* = *p<sup>ν</sup> ε μναβ* .

A list of matrix elements for different decay reactions as predicted by the CCQM is given in what follows. For each element, we provide, in the last line of the corresponding expression, the form factor parametrization of the matrix element to be compared with the appropriate expression from Equations (85), (87), and (89). Beforehand, let us also define the argument of <sup>Φ</sup> *Zb* (*<sup>η</sup>* <sup>2</sup>). One has:

$$\begin{aligned} \vec{\eta}^2 &= \sum\_{i=1}^3 \eta\_i^2 & \eta\_1 &= \begin{array}{rcl} +\frac{1}{2\sqrt{2}} \left( 2k\_1 + (1 + w\_1 - w\_2)q\_1 + (w\_1 - w\_2)q\_2 \right), \end{array} \\\\ \eta\_2 &=& +\frac{1}{2\sqrt{2}} \left( 2k\_2 - (w\_3 - w\_4)q\_1 + (1 - w\_3 + w\_4)q\_2 \right), \\\eta\_3 &=& +\frac{1}{2} \left( (1 - w\_1 - w\_2)q\_1 - (w\_1 + w\_2)q\_2 \right), \end{aligned}$$

where *wi* denotes four body reduced masses *wi* <sup>=</sup> *mi*/ <sup>4</sup> ∑ *j*=1 *mj* and quarks are indexed as *q*<sup>1</sup> = *q*<sup>2</sup> = *b*, *q*<sup>3</sup> = *q*<sup>4</sup> = *d* = *u*.

• <sup>1</sup><sup>+</sup> <sup>→</sup> <sup>1</sup><sup>−</sup> <sup>+</sup> <sup>0</sup><sup>−</sup> matrix elements parametrized as in Equation (85):

$$\mathcal{M}^{\mu\delta}\left(Z\_{b}(p,\mu)\rightarrow\mathcal{Y}(q\_{1},\delta)+\pi^{+}(q\_{2})\right)=\frac{3}{\sqrt{2}}\,\mathfrak{S}z\_{b}\mathfrak{S}\chi\mathfrak{S}\pi$$

$$\int \frac{d^{4}k\_{1}}{(2\pi)^{4}i}\int \frac{d^{4}k\_{2}}{(2\pi)^{4}i}\tilde{\Phi}\_{\mathcal{Z}\_{b}}\left(-\vec{\eta}^{2}\right)^{2}\tilde{\Phi}\_{\mathcal{Y}}\left(-(k\_{1}+\upsilon\_{1}q\_{1})^{2}\right)\tilde{\Phi}\_{\mathcal{R}}\left(-(k\_{2}+\mathsf{u}\_{4}q\_{2})^{2}\right)$$

$$\times\quad\left\{\text{tr}\left[\gamma\_{5}S\_{1}(k\_{1})\gamma^{\delta}S\_{2}(k\_{1}+q\_{1})\gamma^{\mu}S\_{4}(k\_{2})\gamma\_{5}S\_{3}(k\_{2}+q\_{2})\right]$$

$$+\text{tr}\left[\gamma^{\mu}S\_{1}(k\_{1})\gamma^{\delta}S\_{2}(k\_{1}+q\_{1})\gamma\_{5}S\_{4}(k\_{2})\gamma\_{5}S\_{3}(k\_{2}+q\_{2})\right]\right\}$$

$$=\left.A\_{Z\_b\chi\_\pi}g^{\mu\delta} + B\_{Z\_b\chi\_\pi}q\_1^{\mu}q\_2^{\delta}\right|$$

$$M^{\mu\delta} \left( Z\_b'(p\_\tau \mu) \to \mathcal{Y}(q\_1, \delta) + \pi^+(q\_2) \right) = \mathfrak{Z}\_{Z\_b'} \mathfrak{Y} \mathcal{Y} \mathcal{G}\_\pi \frac{i\varepsilon^{\mu\nu\mu\delta}}{M\_{Z\_b'}} \tag{99}$$

$$\begin{split} \int \frac{d^4 k\_1}{(2\pi)^4 i} \int \frac{d^4 k\_2}{(2\pi)^4 i} \tilde{\Phi}\_{Z\_b'} \left( -\vec{\eta}^2 \right) \tilde{\Phi}\_Y \left( -(k\_1 + v\_1 q\_1)^2 \right) \tilde{\Phi}\_\pi \left( -(k\_2 + u\_4 q\_2)^2 \right) \\ \times \quad \text{tr} \left[ \gamma\_a S\_1(k\_1) \gamma^\delta S\_2(k\_1 + q\_1) \gamma\_\beta S\_4(k\_2) \gamma\_5 S\_3(k\_2 + q\_2) \right] \end{split} \tag{90}$$

$$=\left.A\_{Z'\_{\flat}Y\pi}\operatorname{g}^{\mu\delta} + B\_{Z'\_{\flat}Y\pi}\operatorname{q}^{\mu}\_{1}\operatorname{q}^{\delta}\_{2}\right.$$

$$M^{\mu\rho} \left( Z\_b(p,\mu) \to \eta\_b(q\_1) + \rho(q\_2,\rho) \right) = \frac{3}{\sqrt{2}} \mathcal{G} z\_b \mathcal{g}\_{\eta\rho} \mathcal{g}\_{\rho} \tag{100}$$

$$\int \frac{d^4 k\_1}{(2\pi)^4 i} \int \frac{d^4 k\_2}{(2\pi)^4 i} \breve{\Phi}\_{\mathcal{Z}b} \left( -\vec{\eta}^2 \right) \breve{\Phi}\_{\eta\_b} \left( -(k\_1 + v\_1 q\_1)^2 \right) \breve{\Phi}\_{\rho} \left( -(k\_2 + u\_4 q\_2)^2 \right)$$

$$\times \quad \left\{ \text{tr} \left[ \gamma\_5 S\_1(k\_1) \gamma\_5 S\_2(k\_1 + q\_1) \gamma^\mu S\_4(k\_2) \gamma^\rho S\_3(k\_2 + q\_2) \right] \right.$$

$$\left. \quad + \text{tr} \left[ \gamma^\mu S\_1(k\_1) \gamma\_5 S\_2(k\_1 + q\_1) \gamma\_5 S\_4(k\_2) \gamma^\rho S\_3(k\_2 + q\_2) \right] \right\}$$

$$=\left.A\_{Z\_{\mathbb{H}}\eta\_{\mathbb{H}^{\rho}}}\mathcal{g}^{\mu\rho} - B\_{Z\_{\mathbb{H}^{\mu}\mathbb{H}^{\rho}}}q\_2^{\mu}q\_1^{\rho}\right|$$

$$\begin{split} \mathcal{M}^{\mu\rho} \left( Z\_b^{\ell}(p\_\nu \mu) \to \eta\_b(q\_1) + \rho(q\_2, \rho) \right) &= \mathfrak{Z} g\_{Z\_b^{\ell}} g\_{\eta\_b} g\_{\rho} \frac{i\varepsilon^{\mu\rho\eta\delta}}{\mathcal{M}\_{Z\_b^{\prime}}} \\ \quad \quad \quad \quad \int \frac{d^4 k\_1}{(2\pi)^4 i} \int \frac{d^4 k\_2}{(2\pi)^4 i} \widetilde{\Phi}\_{Z\_b^{\prime}} \left( -\widetilde{\eta}^2 \right) \widetilde{\Phi}\_{\eta\_b} \left( -(k\_1 + v\_1 q\_1)^2 \right)^2 \widetilde{\Phi}\_{\rho} \left( -(k\_2 + u\_4 q\_2)^2 \right) \\ \quad \times \quad \quad \quad \quad \quad \left[ \gamma\_a S\_1(k\_1) \gamma\_5 S\_2(k\_1 + q\_1) \gamma\_\beta s\_4(k\_2) \gamma^\rho S\_3(k\_2 + q\_2) \right] \\ \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \end{split}$$

• <sup>1</sup><sup>+</sup> <sup>→</sup> <sup>1</sup><sup>+</sup> <sup>+</sup> <sup>0</sup><sup>−</sup> matrix elements parametrized as in Equation (87):

$$M^{\mathfrak{gl}}\left(Z\_{\mathfrak{b}}^{+}(p\_{\prime},\mathfrak{p}) \to h\_{\mathfrak{b}}(q\_{1},\delta) + \pi^{+}(q\_{2})\right) = \frac{3}{\sqrt{2}} \operatorname{g} z\_{\mathfrak{b}} \mathcal{g}\_{h\mathfrak{b}} \mathcal{g}\_{\mathfrak{n}} \tag{102}$$

$$\begin{split} & \quad \int \frac{d^{4}k\_{1}}{(2\pi)^{4}i} \int \frac{d^{4}k\_{2}}{(2\pi)^{4}i} \tilde{\Phi}\_{\mathfrak{Z}\_{\mathfrak{b}}} \left(-\vec{\eta}^{2}\right) \tilde{\Phi}\_{\mathfrak{h}\_{\mathfrak{b}}} \left(-(k\_{1} + v\_{1}q\_{1})^{2}\right) \tilde{\Phi}\_{\mathfrak{n}} \left(-(k\_{2} + u\_{4}q\_{2})^{2}\right) \\ & \quad \times \left\{\operatorname{tr}\left[\gamma\mathfrak{s}S\_{1}(k\_{1})\gamma\mathfrak{s} \cdot (2k\_{1}^{\delta})S\_{2}(k\_{1} + q\_{1})\gamma^{\mu}S\_{4}(k\_{2})\gamma\mathfrak{s}S\_{3}(k\_{2} + q\_{2})\right] \\ & \quad \quad + \operatorname{tr}\left[\gamma^{\mu}S\_{1}(k\_{1})\gamma\_{\mathfrak{b}} \cdot (2k\_{1}^{\delta})S\_{2}(k\_{1} + q\_{1})\gamma\_{\mathfrak{t}}S\_{4}(k\_{2})\gamma\_{\mathfrak{t}}S\_{3}(k\_{2} + q\_{2})\right] \right\} \end{split}$$

$$=\,^\varepsilon \varepsilon^{\mu\delta q\_1 q\_2} A\_{Z\_b l\_b;\pi} \ast$$

$$M^{\mu\delta} \left( Z\_b'(p\_\nu \mu) \to h\_b(q\_1, \delta) + \pi^+(q\_2) \right) = 3 \, \mathcal{G} Z\_{\bar{b}}^\* \mathcal{G}\_b \mathcal{G}\_\pi \frac{i\varepsilon^{\mu\nu\alpha\beta}}{M\_{Z\_b'}} \tag{103}$$

$$\begin{split} \int \frac{d^4 k\_1}{(2\pi)^4 i} \int \frac{d^4 k\_2}{(2\pi)^4 i} \tilde{\Phi}\_{Z\_b'} \left( -\vec{\eta}^2 \right) \tilde{\Phi}\_{\hbar \flat} \left( -(k\_1 + v\_1 q\_1)^2 \right) \tilde{\Phi}\_\pi \left( -(k\_2 + u\_4 q\_2)^2 \right) \\ \times \quad \text{tr} \left[ \gamma\_\mu \mathcal{S}\_1(k\_1) \gamma\_5 \cdot (2k\_1^\sharp) \mathcal{S}\_2(k\_1 + q\_1) \gamma\_\beta \mathcal{S}\_4(k\_2) \gamma\_5 \mathcal{S}\_3(k\_2 + q\_2) \right] \\ = \varepsilon^{\mu\delta q\_1 q\_2} A\_{z\_b' b\_b \pi} \, . \end{split}$$

The matrix elements describing decays to a pair of *B* mesons can be also listed within two groups depending on the quantum numbers. The argument of the *Zb*-vertex function *δ*<sup>2</sup> is defined as:

$$\delta^2 = \sum\_{i=1}^3 \delta\_i^2; \qquad \delta\_1 = \begin{array}{rcl} -\frac{1}{2\sqrt{2}} \left( k\_1 + k\_2 + (w\_1 - w\_2)q\_1 + (1 + w\_1 - w\_2)q\_2 \right), \\ + \frac{1}{2\sqrt{2}} \left( k\_1 + k\_2 + (1 - w\_3 + w\_4)q\_1 - (w\_3 - w\_4)q\_2 \right), \\ + \frac{1}{2} \left( k\_1 - k\_2 + (w\_1 + w\_2)q\_1 - (1 - w\_1 - w\_2)q\_2 \right). \end{array} \tag{104}$$

The quark indices are similar to the previous case *q*<sup>1</sup> = *q*<sup>2</sup> = *b*, *q*<sup>3</sup> = *q*<sup>4</sup> = *d* = *u*, *v*ˆ2 = *m*2/(*m*<sup>2</sup> + *m*4), *v*ˆ4 = *m*4/(*m*<sup>2</sup> + *m*4), *u*ˆ1 = *m*1/(*m*<sup>1</sup> + *m*3), and *u*ˆ3 = *m*3/(*m*<sup>1</sup> + *m*3).

• <sup>1</sup><sup>+</sup> <sup>→</sup> <sup>1</sup><sup>−</sup> <sup>+</sup> <sup>0</sup><sup>−</sup> matrix elements parametrized as in Equation (85):

$$M^{\mu\rho}\left(Z\_b^+(p,\mu)\to\bar{B}^0(q\_1)+B^{\*+}(q\_2,\rho)\right)=\frac{9}{\sqrt{2}}\,\mathfrak{S}\mathbf{z}\_b\mathcal{S}\mathbf{z}\mathcal{S}^\*\!\!/ $$

$$\begin{split} &\int \frac{d^4k\_1}{(2\pi)^4i}\int \frac{d^4k\_2}{(2\pi)^4i}\tilde{\Phi}\_{\mathbf{Z}\_b}\left(-\tilde{\delta}^2\right)\tilde{\Phi}\_B\left(-(k\_2+v\_4q\_1)^2\right)\tilde{\Phi}\_{\mathbf{B}^\*}\left(-(k\_1+u\_1q\_2)^2\right) \\ &\times \quad \text{tr}\left[\gamma^\mu\mathcal{S}\_1(k\_1)\gamma^\rho\mathcal{S}\_3(k\_1+q\_2)\right]\text{tr}\left[\gamma\_5\mathcal{S}\_4(k\_2)\gamma\_5\mathcal{S}\_2(k\_2+q\_1)\right] \end{split} \tag{105}$$

$$=\left.A\_{Z\_bBB^\*}\mathcal{G}^{\mu\rho} - B\_{Z\_bBB^\*}q\_2^{\mu}q\_1^{\rho}\right|$$

$$M^{\mu\mu}\left(Z\_b^+(p,\mu)\to B^+(q\_1,\delta)+B^+(q\_2)\right)=\frac{9}{\sqrt{2}}\,\mathcal{G}\_b\mathcal{G}\_B\cdot\mathcal{G}\_B\tag{106}$$

$$\int \frac{d^4k\_1}{(2\pi)^4i}\int \frac{d^4k\_2}{(2\pi)^4i}\,\widetilde{\Phi}\_{\mathcal{Z}\_b}\left(-\widetilde{\delta}^2\right)\,\widetilde{\Phi}\_{B^+}\left(-(k\_1+\mathfrak{H}\_1q\_1)^2\right)\,\widetilde{\Phi}\_B\left(-(k\_2+\mathfrak{H}\_2q\_2)^2\right)$$

$$\times\quad\text{tr}\left[\gamma\_5S\_1(k\_1)\gamma\_5S\_3(k\_1+q\_2)\right]\text{tr}\left[\gamma^\mu S\_4(k\_2)\gamma^\delta S\_2(k\_2+q\_1)\right]$$

$$=\quad A\_{Z\_bB^+B}\,\mathcal{g}^{\mu\delta}+B\_{Z\_bB^+B}\,\mathcal{q}^\mu\_1\mathcal{q}^\delta\_2.$$

• <sup>1</sup><sup>+</sup> <sup>→</sup> <sup>1</sup><sup>−</sup> <sup>+</sup> <sup>1</sup><sup>−</sup> matrix elements parametrized as in Equation (89):

$$\begin{split} &M^{\mu\delta\rho} \left( Z\_{b}^{\prime +} (p,\mu) \to B^{\*0} (q\_{1},\delta) + \bar{B}^{\*+} (q\_{2},\rho) \right) = 9\,\, \mathcal{G}\_{Z\_{b}^{\prime}} \mathcal{G}\_{B^{\*}} \mathcal{G}\_{B^{\*}} \frac{\varepsilon^{\mu\nu\delta\rho}}{\mathcal{M}\_{Z\_{b}^{\prime}}} \\ &\quad \int \frac{d^{4}k\_{1}}{(2\pi)^{4}i} \int \frac{d^{4}k\_{2}}{(2\pi)^{4}i} \tilde{\Phi}\_{Z\_{b}^{\prime}} \left( -\tilde{\delta}^{2} \right) \tilde{\Phi}\_{B^{\*}} \left( -(k\_{1} + \mathfrak{d}\_{1} q\_{1})^{2} \right) \tilde{\Phi}\_{B^{\*}} \left( -(k\_{2} + \mathfrak{d}\_{4} q\_{2})^{2} \right) \\ &\quad \times \quad \text{tr} \left[ \gamma\_{\mu} S\_{1} (k\_{1}) \gamma^{\delta} S\_{3} (k\_{1} + q\_{1}) \right] \text{tr} \left[ \gamma\_{\beta} S\_{4} (k\_{2}) \gamma^{\rho} S\_{2} (k\_{2} + q\_{2}) \right] \\ &= \quad B\_{1} q\_{1}^{\mu} \varepsilon^{q\_{1} q\_{2} \rho \delta} + B\_{2} \varepsilon^{q\_{1} \mu \rho \delta} + B\_{3} \varepsilon^{q\_{2} \mu \rho \delta} \text{.} \end{split}$$

With all the above theoretical expressions, one can proceed to the numerical evaluation of the decay widths. The first step is the adjustment of the size parameters Λ*<sup>Z</sup>* and Λ- *<sup>Z</sup>*. They are tuned so as to respect the observables measured by the Belle collaboration [35]:

$$
\Gamma\_{Z\_b}(BB^\*\pi) \quad = \begin{pmatrix} 25 \pm 7 \end{pmatrix} \text{MeV} \, \prime \qquad \mathcal{B}(Z\_b^+ \to [B^+B^{\*0} + \mathcal{B}^0 B^{\*+} ]) = 85.6^{+1.5 + 1.5}\_{-2.0 - 2.1} \,\%, \qquad \Gamma\_{Z\_b^{'+} \to \mathcal{B}^{\*+} \,\%} \qquad \Gamma\_{Z\_b^{'-} \to \mathcal{B}^{\*+} \,\%} \qquad \Gamma\_{Z\_b^{'+} \to \mathcal{B}^{\*+} \,\%} \tag{108}
$$
 
$$
\Gamma\_{Z\_b^{'+} \to \mathcal{B}^\*} \, (B^\*B^\*\pi) \quad = \quad (23 \pm 8) \, \text{MeV} \, \prime \qquad \mathcal{B}(Z\_b^{'+} \to \mathcal{B}^{\*+} \, B^{\*0}) = \mathcal{T}3.7^{+3.4 + 2.7}\_{-4.4 - 3.5} \,\% \,\prime \,\tag{108}
$$

leading to:

$$
\Lambda\_{Z\_b} = 3.45 \pm 0.05 \text{ GeV} \qquad \Lambda\_{Z'\_b} = 3.00 \pm 0.05 \text{ GeV} \,. \tag{109}
$$

With the decays into *B* pairs dominating all other decay channels, we approximate the total decay width as the sum of all herein evaluated channels. The CCQM gives:

$$
\Gamma\_{Z\_b} = 30.9^{+2.3}\_{-2.1} \,\mathrm{MeV} \,, \qquad \Gamma\_{Z'\_b} = 34.1^{+2.8}\_{-2.5} \,\mathrm{MeV} \,\mathrm{V} \,\tag{110}
$$

which is in fair agreement with (108). The predicted partial decay widths of *Zb*(10610) and *Z*- *<sup>b</sup>*(10650) particles are summarized in Table 6.


**Table 6.** Particle decay widths for the *Z*<sup>+</sup> *<sup>b</sup>* (10610) and *<sup>Z</sup>*<sup>+</sup> *<sup>b</sup>* (10650).

The *Zb* and *Z*- *<sup>b</sup>* decays are dominated [13] by <sup>Γ</sup>*Zb* (*B*+*B*∗<sup>0</sup> +*B*∗+*B*<sup>0</sup> )=(85.6+2.1 <sup>−</sup>2.9)% and Γ*Z*- *b* (*B*∗+*B*∗<sup>0</sup> )=(74+<sup>4</sup> <sup>−</sup>6)%, respectively, meaning that the bottomonia modes should not exceed 15 and 25 percent. This is observed for the *hb*(1*P*)*π*<sup>+</sup> final state; the other bottomonia channels are suppressed, but not so much as seen in the data:

$$\begin{array}{llll} \frac{\Gamma\left(Z\_{b} \rightarrow \operatorname{Y}(1S)\pi\right)}{\Gamma\left(Z\_{b} \rightarrow B\overline{B}^{\*} + c.c.\right)} & \approx & 0.29 \,, & \frac{\Gamma\left(Z\_{b} \rightarrow \eta\_{b}\rho\right)}{\Gamma\left(Z\_{b} \rightarrow B\overline{B}^{\*} + c.c.\right)} \approx 0.21 \,, \\\frac{\Gamma\left(Z\_{b}' \rightarrow \operatorname{Y}(1S)\pi\right)}{\Gamma\left(Z\_{b}' \rightarrow B^{\*}\overline{B}^{\*}\right)} & \approx & 0.56 \,, & \frac{\Gamma\left(Z\_{b}' \rightarrow \eta\_{b}\rho\right)}{\Gamma\left(Z\_{b}' \rightarrow B^{\*}\overline{B}^{\*}\right)} \approx 0.44 \,. \end{array}$$

The model also allows us to make predictions:

$$R\_{\mathbf{Y}(1S)\pi} = \frac{\Gamma(Z\_b \to \mathbf{Y}(1S)\pi)}{\Gamma(Z\_b' \to \mathbf{Y}(1S)\pi)} = 0.62 \pm 0.06, \qquad R\_{\eta\_b \rho} = \frac{\Gamma(Z\_b \to \eta\_b \rho)}{\Gamma(Z\_b' \to \eta\_b \rho)} = 0.59 \pm 0.06. \tag{111}$$

One can conclude that the CCQM provides, within a molecular picture, a fair description of *Zb*(10610)/*Z*- *<sup>b</sup>*(10650) states and related decay observables and catches the tendencies seen in experimental data. Some deviations are observed when the fraction of bottomonium in final states is considered.

#### **8. Summary and Conclusions**

The confined covariant quark model is an approach based on a non-local interaction Lagrangian of quarks and hadrons. It has many appealing features: a full Lorentz invariance, confinement, large applicability range (from mesons to exotic hadrons), inclusion of the electromagnetic interaction, and a limited number of free parameters. As a practical tool, it allows overcoming the difficulties related to the non-applicability of the perturbative approach for bound states in QCD. In this text, we used it to describe four quark exotic states *X*(3872), *Zc*(3900), *Y*(4260), *Zb*(10610), and *Z*- *<sup>b</sup>*(10650). We demonstrated that the CCQM had enough predictive power to make the distinctions between various hypothesis, with respect to the exotic quarkonia mostly related to their structure (molecular versus tetraquark one). At the same time, the model provides a good description of experimental data without large deviations and predictions for future measurements. Concerning the structure of the studied particles, the molecular picture is favored for *Zc*(3900), *Zb*(10610), and *Z*- *<sup>b</sup>*(10650) and the tetraquark one for *X*(3872) and *Y*(4260). These conclusions follow from the measured decay characteristics of the considered exotic states and the related model description: with the expected increase in the number and quality of experimental data, one may hope the quarkonia-structure puzzle will be solved in the years to come.

**Author Contributions:** Conceptualization, M.A.I.; methodology, M.A.I., A.Z.D. and S.D.; software, M.A.I.; validation, A.Z.D., S.D., M.A.I. and A.L.; formal analysis, M.A.I. and A.L.; investigation, A.Z.D., S.D., M.A.I. and A.L.; resources, A.Z.D., S.D., M.A.I. and A.L.; data curation, M.A.I. and A.L.; writing—original draft preparation, A.L.; writing—review and editing, A.Z.D., S.D., M.A.I. and A.L.; visualization, M.A.I.; supervision, M.A.I.; project administration, M.A.I. and A.L.; funding acquisition, A.Z.D, S.D., M.A.I. and A.L. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Joint Research Project of the Institute of Physics, Slovak Academy of Sciences (SAS), and Bogoliubov Laboratory of Theoretical Physics, Joint Institute for Nuclear Research (JINR), Grant No. 01-3-1135-2019/2023. A.Z.D., S.D., and A.L. acknowledge the funding from the Slovak Grant Agency for Sciences (VEGA), Grant No. 2/0153/17.

**Acknowledgments:** We would like to thank Thomas Gutsche, Jürgen Körner, Valery Lyubovitskij, Pietro Santorelli, A. Issadykov, F. Goerke, K. Xu, and G. G. Saidullaeva for their collaboration, by which the results discussed in this review have been obtained.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


c 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
