Tri-meson bound state $BBB^{\ast}$ via delocalized $\pi$~Bond

TL;DR

This paper explores the possibility of a bound state formed by three B mesons ($BBB^*$) using the one-pion exchange potential, suggesting such a tri-meson molecule could exist if its two-body subsystems are bound.

Contribution

It introduces a novel investigation of the $BBB^*$ three-meson system as a potential tri-meson molecular state using three-body Schrödinger equation analysis.

Findings

Predicted the possible existence of a $BBB^*$ bound state.

Established a link between two-body $BB^*$ molecular states and the three-body $BBB^*$ state.

Provided a theoretical framework for exploring tri-meson molecules.

Abstract

During the last decades, numerous exotic states which cannot be explained by the conventional quark model have been observed in experiment. Some of them can be understood as two-body hadronic molecules, such as the famous $X(3872)$, analogous to deuteron in nuclear physics. Along the same line, the existence of the triton leaves an open question whether there is a bound state formed by three hadrons. Since, for a given potential, a system with large reduced masses is more easier to form a bound state, we study the $BBB^{\ast}$ system with the one-pion exchange potential as an exploratory step by solving the three-body Schr\"odinger Equation. We predict that a tri-meson molecular state for the $BBB^{\ast}$ system is probably existent as long as the molecular states of its two-body subsystem $BB^*$ exist.

Tables

Table 1. Table 1: The coupling constants and meson masses in our calculation. The meson masses are taken from the PDG Patrignani:2016xqp

mass(MeV)	coupling constants
$m_{\pi }=139$	$g=0.57$
$M_B=5279$	$f_{\pi }=132.00$ MeV
$M_{B^{\ast }}=5325$

Table 2. Table 2: Channel dependent coefficients. Here, c 𝑐 c denotes the C-parity of the two-body system.

channel	isospin	$C(i,j)$	channel	$C(i,j)$
$B{B}^{\ast }$	$I=1$	$1/2$	$B^{+}{B}^{\ast +}$	$1/2$
	$I=0$	$-3/2$	$B^{+}{B}^{\ast 0}$	$-1/2$
$B{\overline{B}}^{\ast }$	$I=1$	$c/2$	$B^0{B}^{\ast +}$	$-1/2$
	$I=0$	$-3c/2$	$B^0{B}^{\ast 0}$	$1/2$

Table 3. Table 3: Bound state solutions of the N ​ N ​ N 𝑁 𝑁 𝑁 NNN system with isospin I 3 = 1 / 2 subscript 𝐼 3 1 2 I_{3}=1/2 . E 2 subscript 𝐸 2 E_{2} is the energy eigenvalue of its subsystem. E 3 subscript 𝐸 3 E_{3} is the reduced three-body energy eigenvalue relative to the break-up state of the N ​ N ​ N 𝑁 𝑁 𝑁 NNN system. E T subscript 𝐸 𝑇 E_{T} is the total three-body energy eigenvalue relative to the N ​ N ​ N 𝑁 𝑁 𝑁 NNN threshold. V B ​ O ​ ( 0 ) subscript 𝑉 𝐵 𝑂 0 V_{BO}(0) is the minimum of the BO potential. r r ​ m ​ s subscript 𝑟 𝑟 𝑚 𝑠 r_{rms} represents the root-mean-square radius of any two N 𝑁 N in the N ​ N ​ N 𝑁 𝑁 𝑁 NNN system. The S 𝑆 S -wave and D 𝐷 D -wave represent the probabilities for S 𝑆 S -wave and D 𝐷 D -wave components in any two N 𝑁 N in the N ​ N ​ N 𝑁 𝑁 𝑁 NNN system.

$\Lambda$(MeV)	$E_2$(MeV)	$E_3$(MeV)	$E_T$(MeV)	$V_{BO}\left(0\right)$(MeV)	S wave(%)	D wave(%)	$r_{rms}$(fm)
830.00	-0.18	-1.93	-2.11	-4.54	94.01	5.99	4.21
850.00	-0.67	-2.71	-3.38	-5.36	93.36	6.64	4.00
870.00	-1.23	-3.65	-4.88	-6.32	92.68	7.32	3.78
890.00	-1.88	-4.77	-6.66	-7.42	91.99	8.01	3.54
899.60	-2.23	-5.38	-7.62	-8.00	91.66	8.34	3.42
900.00	-2.25	-5.41	-7.66	-8.03	91.64	8.36	3.42
920.00	-3.05	-6.85	-9.90	-9.35	90.97	9.03	3.18
940.00	-3.98	-8.51	-12.49	-10.83	90.35	9.65	2.95
960.00	-5.03	-10.42	-15.45	-12.46	89.76	10.24	2.74
980.00	-6.21	-12.57	-18.78	-14.23	89.23	10.77	2.54
1000.00	-7.55	-14.97	-22.51	-16.14	88.73	11.27	2.37
1020.00	-9.04	-17.61	-26.65	-18.19	88.27	11.73	2.23
1040.00	-10.69	-20.51	-31.20	-20.37	87.84	12.16	2.10

Table 4. Table 4: Bound state solutions of the B ​ B ∗ 𝐵 superscript 𝐵 ∗ BB^{\ast} system with the isospin I 2 = 1 subscript 𝐼 2 1 I_{2}=1 . Λ Λ \Lambda is the parameter in the form factor. E 2 subscript 𝐸 2 E_{2} is the energy eigenvalue. The binding energy is − E 2 subscript 𝐸 2 -E_{2} . r r ​ m ​ s subscript 𝑟 𝑟 𝑚 𝑠 r_{rms} is the root-mean-square radius. α 𝛼 \alpha and β 𝛽 \beta are the probabilities for the components B ​ B ∗ 𝐵 superscript 𝐵 ∗ BB^{\ast} and B ∗ ​ B ∗ superscript 𝐵 ∗ superscript 𝐵 ∗ B^{\ast}B^{\ast} , respectively.

		$B{B}^{\ast }$		$B^{\ast }{B}^{\ast }$			Proportion
$\Lambda$(MeV)	$E_2$(MeV)	S wave(%)	D wave(%)	S wave(%)	D wave(%)	$r_{rms1}$(fm)	$\alpha$(%)	$\beta$(%)
1380	-2.11	99.13	0.87	99.21	0.79	1.51	91.91	8.09
1400	-2.94	99.12	0.88	99.36	0.64	1.30	90.22	9.78
1420	-3.93	99.14	0.86	99.49	0.51	1.15	88.47	11.53
1440	-5.08	99.16	0.84	99.59	0.41	1.03	86.69	13.31
1460	-6.40	99.19	0.81	99.68	0.33	0.94	84.91	15.09
1480	-7.88	99.22	0.78	99.74	0.26	0.86	83.14	16.86
1500	-9.54	99.25	0.75	99.80	0.20	0.80	81.40	18.60
1520	-11.38	99.29	0.71	99.84	0.16	0.75	79.70	20.30
1540	-13.39	99.32	0.68	99.88	0.12	0.71	78.07	21.93
1560	-15.59	99.36	0.64	99.91	0.09	0.67	76.50	23.50

Table 5. Table 5: Bound state solutios of the B ​ B ∗ 𝐵 superscript 𝐵 ∗ BB^{\ast} with the isospin I 2 = 0 subscript 𝐼 2 0 I_{2}=0 . Λ Λ \Lambda is the parameter in the form factor. E 2 subscript 𝐸 2 E_{2} is the energy eigenvalue. The binding energy is − E 2 subscript 𝐸 2 -E_{2} . r r ​ m ​ s subscript 𝑟 𝑟 𝑚 𝑠 r_{rms} is the root-mean-square radius. α 𝛼 \alpha and β 𝛽 \beta are the probabilities for the component B ​ B ∗ 𝐵 superscript 𝐵 ∗ BB^{\ast} and B ∗ ​ B ∗ superscript 𝐵 ∗ superscript 𝐵 ∗ B^{\ast}B^{\ast} , respectively.

		$B{B}^{\ast }$		$B^{\ast }{B}^{\ast }$			Proportion
$\Lambda$(MeV)	$E_2$(MeV)	S wave(%)	D wave(%)	S wave(%)	D wave(%)	$r_{rms1}$(fm)	$\alpha$(%)	$\beta$(%)
1040	-1.88	86.09	13.91	68.84	31.16	2.03	92.05	7.95
1060	-2.63	84.59	15.41	69.41	30.59	1.79	90.18	9.82
1080	-3.54	83.26	16.74	69.90	30.10	1.60	88.22	11.78
1100	-4.62	82.07	17.93	70.33	29.67	1.45	86.22	13.78
1120	-5.87	81.01	18.99	70.69	29.31	1.33	84.22	15.78
1140	-7.29	80.07	19.93	70.99	29.01	1.23	82.25	17.75
1160	-8.89	79.23	20.77	71.25	28.75	1.14	80.34	19.66
1180	-10.67	78.49	21.51	71.46	28.54	1.07	78.50	21.50
1200	-12.63	77.81	22.19	71.63	28.37	1.01	76.75	23.25
1220	-14.78	77.21	22.79	71.77	28.23	0.95	75.09	24.91

Table 6. Table 6: Bound state solutions of the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} with the isospin I 3 = 3 / 2 subscript 𝐼 3 3 2 I_{3}=3/2 . E 2 subscript 𝐸 2 E_{2} is the energy eigenvalue of its subsystem B ​ B ∗ 𝐵 superscript 𝐵 ∗ BB^{\ast} with the isospin I 2 = 1 subscript 𝐼 2 1 I_{2}=1 . E 3 subscript 𝐸 3 E_{3} is the reduced three-body energy eigenvalue relative to the break-up state of the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} system. E T subscript 𝐸 𝑇 E_{T} is the total three-body energy eigenvalue relative to the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} threshold. V B ​ O ​ ( 0 ) subscript 𝑉 𝐵 𝑂 0 V_{BO}(0) is the minimum of the BO potential. r r ​ m ​ s subscript 𝑟 𝑟 𝑚 𝑠 r_{rms} represents the root-mean-square radius of any two B 𝐵 B in the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} system. The S 𝑆 S -wave and D 𝐷 D -wave represent the probabilities for S 𝑆 S -wave and D 𝐷 D -wave components in any two B 𝐵 B in the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} . The α 𝛼 \alpha and β 𝛽 \beta denote the probabilities for the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} and B ​ B ∗ ​ B ∗ 𝐵 superscript 𝐵 ∗ superscript 𝐵 ∗ BB^{\ast}B^{\ast} components, respectively.

$E_2$(MeV)	$E_3$(MeV)	$E_T$(MeV)	$V_{BO}\left(0\right)$(MeV)	S wave(%)	D wave(%)	$r_{rms}$(fm)	$\alpha$(%)	$\beta$(%)
-0.18	-0.19	-0.38	-3.43	99.76	0.24	3.98	97.50	2.50
-0.48	-0.45	-0.93	-4.88	99.68	0.32	3.34	96.39	3.61
-0.89	-0.85	-1.74	-6.62	99.59	0.41	2.67	95.02	4.98
-1.43	-1.42	-2.85	-8.56	99.49	0.51	2.11	93.78	6.22
-2.11	-2.20	-4.31	-10.65	99.41	0.59	1.71	91.91	8.09
-2.94	-3.17	-6.11	-12.87	99.34	0.66	1.43	90.22	9.78
-3.93	-4.33	-8.26	-15.21	99.29	0.71	1.24	88.47	11.53
-5.08	-5.67	-10.75	-17.65	99.25	0.75	1.09	86.69	13.31
-6.40	-7.18	-13.58	-20.19	99.22	0.78	0.98	84.91	15.09
-7.88	-8.83	-16.71	-22.83	99.20	0.80	0.90	83.14	16.86
-9.54	-10.61	-20.16	-25.55	99.18	0.82	0.83	81.40	18.60
-11.38	-12.51	-23.89	-28.36	99.17	0.83	0.77	79.70	20.30
-13.39	-14.62	-28.01	-31.24	99.16	0.84	0.72	78.07	21.93
-15.59	-16.75	-32.34	-34.20	99.15	0.85	0.68	76.50	23.50
-17.97	-18.99	-36.95	-37.24	99.14	0.86	0.65	75.01	24.99

Table 7. Table 7: Bound state solutions of the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} with isospin I = 1 / 2 𝐼 1 2 I=1/2 . E 2 subscript 𝐸 2 E_{2} is the energy eigenvalue of its subsystem B ​ B ∗ 𝐵 superscript 𝐵 ∗ BB^{\ast} with the isospin I = 0 𝐼 0 I=0 . E 3 subscript 𝐸 3 E_{3} is the reduced three-body energy eigenvalue relative to the break-up state of the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} system. E T subscript 𝐸 𝑇 E_{T} is the total three-body energy eigenvalue relative to the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} threshold. V B ​ O ​ ( 0 ) subscript 𝑉 𝐵 𝑂 0 V_{BO}(0) is the minimum of the BO potential. r r ​ m ​ s subscript 𝑟 𝑟 𝑚 𝑠 r_{rms} represents the root-mean-square radius of any two B 𝐵 B in the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} system. The S 𝑆 S -wave and D 𝐷 D -wave represent the probabilities for S 𝑆 S -wave and D 𝐷 D -wave components in any two B 𝐵 B in the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} . The α 𝛼 \alpha and β 𝛽 \beta denote the probabilities for the B ​ B ​ B ∗ 𝐵 𝐵 superscript 𝐵 ∗ BBB^{\ast} and B ​ B ∗ ​ B ∗ 𝐵 superscript 𝐵 ∗ superscript 𝐵 ∗ BB^{\ast}B^{\ast} components, respectively.

$E_2$(MeV)	$E_3$(MeV)	$E_T$(MeV)	$V_{BO}\left(0\right)$(MeV)	S wave(%)	D wave(%)	$r_{rms}$(fm)	$\alpha$(%)	$\beta$(%)
-0.19	-0.32	-0.51	-2.15	94.66	5.34	3.89	97.68	2.32
-0.44	-0.64	-1.08	-3.04	92.56	7.44	3.27	96.66	3.34
-0.80	-1.13	-1.93	-4.19	90.43	9.57	2.69	95.36	4.64
-1.27	-1.82	-3.09	-5.57	88.49	11.51	2.24	93.80	6.20
-1.88	-2.72	-4.60	-7.14	86.72	13.28	1.93	92.05	7.95
-2.63	-3.82	-6.45	-8.89	85.09	14.91	1.70	90.18	9.82
-3.54	-5.11	-8.65	-10.78	83.60	16.40	1.53	88.22	11.78
-4.62	-6.57	-11.20	-12.82	82.27	17.73	1.40	86.22	13.78
-5.04	-7.13	-12.17	-13.56	81.84	18.16	1.36	85.52	14.48
-7.29	-10.00	-17.29	-17.25	80.06	19.94	1.21	82.25	17.75
-8.89	-11.93	-20.83	-19.63	79.15	20.85	1.13	80.34	19.66
-10.67	-14.00	-24.68	-22.12	78.36	21.64	1.07	78.50	21.50
-12.63	-16.20	-28.84	-24.70	77.66	22.34	1.02	76.75	23.25
-14.78	-18.52	-33.30	-27.38	77.03	22.97	0.97	75.09	24.91
-17.10	-20.96	-38.06	-30.15	76.46	23.54	0.93	73.52	26.48

Equations

$$|1,\frac{3}{2},\frac{3}{2}⟩$$

$$|1,\frac{3}{2},-\frac{3}{2}⟩$$

$$|1,\frac{3}{2},\frac{1}{2}⟩$$

$$|1,\frac{3}{2},-\frac{1}{2}⟩$$

$$|1,\frac{1}{2},\frac{1}{2}⟩$$

$$|1,\frac{1}{2},-\frac{1}{2}⟩$$

$$|0,\frac{1}{2},\frac{1}{2}⟩$$

$$|0,\frac{1}{2},-\frac{1}{2}⟩$$

$${\mathcal{L}}_P$$

$$\widetilde{{\mathcal{L}}_P}$$

$$\displaystyle\phi=\left(\begin{array}[]{cc}\frac{\pi^{0}}{\sqrt{2}}&\pi^{+}\\ \pi^{-}&-\frac{\pi^{0}}{\sqrt{2}}\\ \end{array}\right).$$

$$\displaystyle V_{BBB^{\ast}}=\left(\begin{array}[]{cccccc}0&V_{1}(\vec{r}_{ab})&V_{1}(\vec{r}_{ac})&V_{2}(\vec{r}_{ab})&V_{2}(\vec{r}_{ac})&0\\ V_{1}(\vec{r}_{ba})&0&V_{1}(\vec{r}_{bc})&V^{\prime}_{2}(\vec{r}_{ba})&0&V_{2}(\vec{r}_{bc})\\ V_{1}(\vec{r}_{ca})&V_{1}(\vec{r}_{cb})&0&0&V^{\prime}_{2}(\vec{r}_{ac})&V^{\prime}_{2}(\vec{r}_{bc})\\ V_{2}(\vec{r}_{ba})&V^{\prime}_{2}(\vec{r}_{ab})&0&V_{3}(\vec{r}_{ab})&V_{1}(\vec{r}_{bc})&V_{1}(\vec{r}_{ac})\\ V_{2}(\vec{r}_{ca})&0&V^{\prime}_{2}(\vec{r}_{ca})&V_{1}(\vec{r}_{cb})&V_{3}(\vec{r}_{ac})&V_{1}(\vec{r}_{ab})\\ 0&V_{2}(\vec{r}_{cb})&V^{\prime}_{2}(\vec{r}_{cb})&V_{1}(\vec{r}_{ca})&V_{1}(\vec{r}_{ba})&V_{3}(\vec{r}_{bc})\\ \end{array}\right).$$

$$\displaystyle H_{ab}=\left(\begin{array}[]{cccccc}T_{*}&V_{1}(\vec{r}_{ab})&0&V_{2}(\vec{r}_{ab})&0&0\\ V_{1}(\vec{r}_{ba})&T_{*}&0&V^{\prime}_{2}(\vec{r}_{ba})&0&0\\ 0&0&T&0&0&0\\ V_{2}(\vec{r}_{ba})&V^{\prime}_{2}(\vec{r}_{ab})&0&T_{**}+V_{3}(\vec{r}_{ab})+\delta M&0&0\\ 0&0&0&0&T_{*}+\delta M&V_{1}(\vec{r}_{ab})\\ 0&0&0&0&V_{1}(\vec{r}_{ba})&T_{*}+\delta M\\ \end{array}\right),$$

$$|1,1⟩$$

$$|1,-1⟩$$

$$|1,0⟩$$

$$|0,0⟩$$

$$\displaystyle\left(\begin{array}[]{ccc}T_{*}&V_{1}(\vec{r}_{ab})&V_{2}(\vec{r}_{ab})\\ V_{1}(\vec{r}_{ba})&T_{*}&V^{\prime}_{2}(\vec{r}_{ba})\\ V_{2}(\vec{r}_{ba})&V^{\prime}_{2}(\vec{r}_{ab})&T_{**}+V_{3}(\vec{r}_{ab})+\delta M\end{array}\right)\left(\begin{array}[]{c}\frac{1}{\sqrt{2}}\psi(\vec{r}_{ab})\\ \frac{1}{\sqrt{2}}\psi(\vec{r}_{ab})\\ \psi^{\prime}(\vec{r}_{ab})\end{array}\right)=E_{b}\left(\begin{array}[]{c}\frac{1}{\sqrt{2}}\psi(\vec{r}_{ab})\\ \frac{1}{\sqrt{2}}\psi(\vec{r}_{ab})\\ \psi^{\prime}(\vec{r}_{ab})\end{array}\right).$$

$$⟨f|S|i⟩={\delta }_{fi}+(2\pi {)}^4{\delta }^4(p_f-p_i)i{M}_{fi}={\delta }_{fi}-2\pi \delta (E_f-E_i)i{V}_{fi},$$

$$V_{fi}=-\frac{M_{fi}}{\sqrt{\prod _{f}2{p_f}^0\prod _{i}2{p_i}^0}}\approx -\frac{M_{fi}}{\sqrt{\prod _{f}2m_f\prod _{i}2m_i}},$$

$$V(\vec{r})=\frac{1}{(2\pi {)}^3}\int d^3\vec{q}{e}^{i\vec{q}\cdot \vec{r}}V(\vec{q}{F}^2(\vec{q}).$$

$$F_i(q)=\frac{{\Lambda }^2-m_{\pi }^2}{{\Lambda }^2-q_i^2}=\frac{{\Lambda }^2-m_{\pi }^2}{{\tilde{\Lambda }}_i^2+{\vec{q}}_i^2},$$

$${\tilde{\Lambda }}^{(\prime )2}={\Lambda }^2-(\Delta M^{(\prime )}{)}^2,$$

$$V_{B{B}^{\ast }\to B{B}^{\ast }}(\vec{r})$$

$$V_{B{B}^{\ast }\to B^{\ast }{B}^{\ast }}(\vec{r})$$

$$V_{B^{\ast }{B}^{\ast }\to B^{\ast }{B}^{\ast }}(\vec{r})$$

$$\Psi (\vec{r})={\psi }_S(\vec{r})+{\psi }_D(\vec{r}),$$

$${\chi }_{nl}(r)=\sqrt{\frac{(2\lambda {)}^{2l+3}n!}{\Gamma (2l+3+n)}}{r}^l{e}^{-\lambda r}{L}_n^{2l+2}(2\lambda r),n=1,2,\mathrm{3...}$$

$${\int }_0^{\infty }{\chi }_{im}(r){\chi }_{in}(r)r^{2}dr={\delta }_{ij}{\delta }_{mn}.$$

$$\psi (\vec{r})$$

$${\psi }^{\prime }(\vec{r})$$

$$\Psi (\Lambda ,{\vec{r}}_{ab})=\frac{1}{\sqrt{2}}\psi (\Lambda ,{\vec{r}}_{ab})|B_a{B}_b^{\ast }⟩+\frac{1}{\sqrt{2}}\psi (\Lambda ,{\vec{r}}_{ab})|B_a^{\ast }{B}_b⟩+{\psi }^{\prime }(\Lambda ,{\vec{r}}_{ab})|B^{\ast }{B}^{\ast }⟩.$$

$$\displaystyle V_{a}=\left(\begin{array}[]{cccccc}0&V_{1}(\vec{r}_{ab})&V_{1}(\vec{r}_{ac})&V_{2}(\vec{r}_{ab})&V_{2}(\vec{r}_{ac})&0\\ V_{1}(\vec{r}_{ba})&0&0&V^{\prime}_{2}(\vec{r}_{ba})&0&0\\ V_{1}(\vec{r}_{ca})&0&0&0&V^{\prime}_{2}(\vec{r}_{ac})&0\\ V_{2}(\vec{r}_{ba})&V^{\prime}_{2}(\vec{r}_{ab})&0&V_{3}(\vec{r}_{ab})&0&V_{1}(\vec{r}_{ac})\\ V_{2}(\vec{r}_{ca})&0&V^{\prime}_{2}(\vec{r}_{ca})&0&V_{3}(\vec{r}_{ac})&V_{1}(\vec{r}_{ab})\\ 0&0&0&V_{1}(\vec{r}_{ca})&V_{1}(\vec{r}_{ba})&0\\ \end{array}\right).$$

$$\displaystyle V_{bc}=\left(\begin{array}[]{cccccc}0&0&0&0&0&0\\ 0&0&V_{1}(\vec{r}_{bc})&0&0&V_{2}(\vec{r}_{bc})\\ 0&V_{1}(\vec{r}_{cb})&0&0&0&V^{\prime}_{2}(\vec{r}_{bc})\\ 0&0&0&0&V_{1}(\vec{r}_{bc})&0\\ 0&0&0&V_{1}(\vec{r}_{cb})&0&0\\ 0&V_{2}(\vec{r}_{cb})&V^{\prime}_{2}(\vec{r}_{cb})&0&0&V_{3}(\vec{r}_{bc})\\ \end{array}\right)$$

$$E_2$$

$$=$$

$$\psi ({\vec{r}}_{ab},{\vec{r}}_{ac})$$

Full text

I Introduction

In the last few decades, numerous exotic states named as “XYZ” as well as charm-strange mesons beyond the conventional quark model have been reported by many experimental collaborations. For a review of these exotic states, we recommend Refs. Klempt:2007cp ; Klempt:2009pi ; Brambilla:2010cs ; Olsen:2014qna ; Oset:2016lyh ; Chen:2016qju ; Chen:2016spr ; Esposito:2016noz ; Lebed:2016hpi ; Hosaka:2016pey ; Dong:2017gaw ; Guo:2017jvc ; Olsen:2017bmm ; Francis:2016hui ; Karliner:2017qjm ; Eichten:2017ffp . Some of them can be understood in the hadronic molecular picture Guo:2017jvc which is an analog of the deuteron as a loosely bound state of a proton and a neutron. In their formation, the one-pion exchange potential (OPEP) plays an important role, e.g. in the formation of deuteron and the $X(3872)$ Machleidt:1987hj ; Tornqvist:2003na , due to its long-range property. Analogously, the existence of the triton arouses interest in the study of the three-hadron system, especially after the large accumulation of experimental data, which might give some hints about the existence of this kind of bound states. In general, to solve the three-body problem, one should solve the Faddeev equations rigorously Malfliet:1968tj ; Eichmann:2009qa ; Ishii:1995bu ; Eichmann:2011vu ; Huang:1993yd ; Ishii:1993np ; Ishikawa:2002ti ; SanchisAlepuz:2011jn ; Elster:2008hn ; Eichmann:2008ef ; Popovici:2010ph ; Fujiwara:2003wr ; Stadler:1991zz . However, for a specific system, one can do some approximation to simplify the problem, such as the Fixed Center Approximation (FCA) in the study of the $X(2175)$ as a resonance of the $\phi K\bar{K}$ system MartinezTorres:2008gy and the approximation on unitary chiral dynamics on the $\pi K\bar{K}$ and $\pi\pi\eta$ systems MartinezTorres:2011vh . The FCA has also been widely applied to other systems, such as the $KK\bar{K}$ Torres:2011jt , the $D\bar{D}^{\ast}K$ and $\bar{D}D^{\ast}K$ Ren:2018pcd , the $J/\psi K\bar{K}$ MartinezTorres:2009xb , the $NDK$, $\bar{K}DN$ and $ND\bar{D}$ Xiao:2011rc , the $N\bar{K}K$ Xie:2010ig , the $BD\bar{D}$ and $BDD$ Dias:2017miz ; Jido:2008kp ; MartinezTorres:2008kh systems. There are many other studies that employ FCA method discussed in Refs. Bayar:2012rk ; Oset:2012gi ; MartinezTorres:2010ax ; Liang:2013yta ; Bayar:2015oea ; YamagataSekihara:2010qk ; Roca:2010tf ; Xie:2011uw ; Debastiani:2017vhv . The isobar familism has also been applied to discuss three-body systems, such as the strange dibaryon resonance in the $\bar{K}NN-\pi\sigma N$ system Ikeda:2007nz , the effect of the $\Delta(1236)$ isobar on the three nucleon bound states Hajduk:1979yn , and other systems MartinezTorres:2008kh ; Ikeda:2008ub ; Gal:2013dca ; Dreissigacker:1981az . The dimer familism is another approximation method for a three-body system where a composite field is intoduced to describe its two-body subsystem when rescattering with a third particle, which has been applied to the three-hadron systems Konig:2015aka ; Konig:2016yka ; Wilbring:2017fwy ; Schmidt:2018vvl . Recently, a series of studies Hammer:2017uqm ; Hammer:2017kms ; Meng:2017jgx of a three-particle system in a finite volume via the dimer field have be proposed to gain a insights about a three-body system in a discretized space such as used in lattice QCD. It is worth mentioning that the system $BB^{*}B^{*}-B^{*}B^{*}B^{*}$ has been studied recently in Ref. Garcilazo:2018rwu based on the colored interaction for its subsytems. By solving exactly the Faddeev equations for the tri-meson system, the authors find a bound state about 90 $\mathrm{MeV}$ below three B mesons threshold. Similar discussions on the $\Omega NN$ and $\Omega\Omega N$ systems can be found in Ref. Garcilazo:2019igo .

As discussed above, the OPEP plays an important role in binding the two-hadron system. From another point of view, one can view it as a pion shared by the two constituents and form a bound state. It can be regarded as a bond similar to the $\sigma~{}bond$ in hydrogen molecules. There is another kind of bond called delocalized $\pi~{}bond$ universally existing in benzene molecules, which is a pair of electrons shared by the six carbon atoms. A simple extension is replacing the carbon atoms by hadrons. We have studied the role of the delocalized $\pi~{}bond$ in forming the three-body bound state for the double heavy tri-meson systems, i.e. $DD^{\ast}K$, $D\bar{D}^{\ast}K$, $BB^{\ast}\bar{K}$ and $B\bar{B}^{\ast}\bar{K}$ Ma:2017ery , based on the sufficient information of it sub two-body system and the Born-Oppenheimer Approximation (BOA) which works well for a system with several heavy and light particles Moroz:2014eba . The crucial idea is to use Born-Oppenheimer (BO) potential for considering the influence of the light part on the dynamics of the heavy part. Therefore, it is a fascinating idea whether the delocalized $\pi~{}bond$ and the BOA could be applied to a three-heavy system, such as the $BBB^{\ast}$ system with large reduced mass.

The same three bottomed meson system has been studied in Ref.  Wilbring:2017fwy by calculating the scattering amplitudes between the $Z_{b}(10610)$ or the $Z_{b}(10650)$ and the bottomed meson. The universal bound states of three bottomed mesons from the Efimov effect has been ruled out. As in their calculation, only the contact interaction is included which might be the reason why they do not find a bound state. After including the long-range OPEP, the case might be different. Thus we solve the three-body Schrödinger equation to discuss the $BBB^{\ast}$ system by considering the OPEP. Without an assumption about its two-body subsystem, i.e. the molecular nature of the $Z_{b}(10610)$ or $Z_{b}(10650)$, we focus on the three-body bound state as a function of the binding energy of its subsystem. Hopefully, the present extensive investigations will be useful to deepen our the understanding of a system made of three heavy particles.

This paper is organized as follows. After the introduction, the formalism and the inputs for the $BBB^{\ast}$ system are presented in Sec. II. The dynamics of the two-body subsystem and the corresponding BO potential are given in Sec. III and Sec. IV, respectively. By constructing a proper interpolating wave functions for the $BBB^{\ast}$ in Sec. V, we solve the three-body Schrödinger equation in Sec. VI. Numerical results and discussions are given in the following section. The summary is presented in the last section. Some technicalities are relegated to the appendix.

II Formalism and the inputs

The BOA has been successfully used in few-body system with several heavy and light particles Braaten:2014qka ; Moroz:2014eba . For a three-body system with one light and two heavy mesons, such as the $DD^{*}K$ Ma:2017ery system, the three-body Schrödinger equation is divided into two sub equations, one is the motion of light meson with two static sources and the other one is the equation for the two heavy mesons with the BO potential Ma:2017ery , which reflects the influence of the light meson on the dynamics of the two heavy mesons. If the interaction between the light meson and the heavy meson is attractive, it would make the two heavy mesons come closer, thus facilitating the formation of a bound state of the whole system. However, for a three heavy meson system, although the application of BOA is not straightforward, one can employ the underlying permutation symmetry, which means the corresponding dynamics is invariant under the interchange of any two constituents, for the system and continue to use the BOA for the calculation.

The OPEP indicates that there is only one virtual pion exchanged by any two constituents as shown in Fig. 1.

One can use $a$, $b$ and $c$ to label the three mesons in the original channel, i.e. $B_{a}^{*}B_{b}B_{c}$. It changes into $B_{a}B_{b}^{*}B_{c}$ via one-pion exchange (OPE) between $a$ and $b$, and the channel $B_{a}B_{b}^{*}B_{c}$ changes into $B_{a}B_{b}B_{c}^{*}$ through the OPE between $b$ and $c$. When the virtual pion arises between $a$ and $c$, it returns back into the original channel $B_{a}^{*}B_{b}B_{c}$. Within this scenario, the virtual pion is not localized between any two constituents but rather shared by the whole system. It is very similar with the benzene molecule which has a pair of electrons shared by the six carbon atoms, which is called $delocalized~{}\pi~{}bond$ in molecular physics. Since the three constituents have the same probability to be the $B$ and $B^{\ast}$ mesons, one can write the system as $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$. Furthermore, the order of the $a$, $b$ and $c$ labels of the three mesons is artificial, as the system is invariant under the interchange of the $a$, $b$ and $c$. This interchange symmetry will help to simplify our calculations. The point is that one can count the influence of each heavy meson on the dynamics of the other two mesons one by one. In other words, we can divide the system into three two-body subsystems $a~{}b$, $b~{}c$ and $a~{}c$. In each subsystem, one should add the BO potential from the remaining one. The existence of a negative common eigenvalue for the three subsystems may partly answer whether there is three-body bound state for the three heavy system. For simplicity, we call this method as $Born-Oppenheimer~{}potential~{}method$ ($BO~{}potential~{}method$).

Before performing the calculation, we define the isospin wave functions of the $BBB^{\ast}$ systems as $|I_{2},I_{3},I_{3z}\rangle$ with $I_{2}$ the isospin of the sub-$BB^{\ast}$ system. $I_{3}$ and $I_{3z}$ represent the total isospin of the three-body system and its $z$ direction, respectively. One thus obtains the isospin wave functions of the $BBB^{\ast}$ system,

[TABLE]

Since $BB^{\ast}$ can couple with $B^{\ast}B^{\ast}$ via OPE, the coupled channel effect is not negligible. We only consider the next close $BB^{\ast}B^{\ast}$ channel in our calculation. If we distinguish the specific locations of the constituents as $a$, $b$ and $c$, there are six channels in total, i.e. $B_{a}^{\ast}B_{b}B_{c}$, $B_{a}B_{b}^{\ast}B_{c}$, $B_{a}B_{b}B_{c}^{\ast}$, $B_{a}^{\ast}B_{b}^{\ast}B_{c}$, $B_{a}^{\ast}B_{b}B_{c}^{\ast}$ and $B_{a}B_{b}^{\ast}B_{c}^{\ast}$.

The Lagrangians with SU(2) chiral symmetry (we only consider OPE) and C-parity conservation read

[TABLE]

where the heavy flavor meson fields $P$ and $P^{*}$ represent $P=(B^{-},\bar{B}^{0})$, $P^{*}=(B^{*-},\bar{B}^{*0})$, respectively. Its corresponding heavy anti-meson fields $\widetilde{P}$ and $\widetilde{P}^{*}$ represent $\widetilde{P}=(B^{+},B^{0})$, $\widetilde{P}^{*}=(B^{*+},B^{*0})$. $\phi$ is the pion matrix

[TABLE]

We use the pion decay constant $f_{\pi}=132~{}\mathrm{MeV}$ Zhao:2014gqa . The pionic coupling constant $g\!=\!0.57$ is extracted from the width of $D^{*+}$ by assuming heavy quark flavor symmetryPatrignani:2016xqp . All the parameters and input datas are listed in Table 1. Here, we neglect isospin breaking effect and use the masses of their charged particles.

Under SU(2) chiral symmetry, the OPE interaction is of order $\mathcal{O}(p^{0})$ for the three-body system. In this paper, we only take into account the OPEP to the $\mathcal{O}(p^{0})$ order. Thus, there are four kinds of effective potentials. We use $V_{1}$ to denote the effective potential for the interaction $BB^{\ast}\to B^{\ast}B$. $V_{2}$ and $V^{\prime}_{2}$ denote the process $BB^{\ast}\to B^{\ast}B^{\ast}$ and its reverse, respectively. $V_{3}$ represents diagonal process $B^{\ast}B^{\ast}\to B^{\ast}B^{\ast}$. Since the interactions are physical, the effective potentials should be unitary which gives $V^{\prime}_{2}=(V_{2})^{{\dagger}}$. $\vec{r}_{ij}$ is used to denote the relative displacement between the $i$-th and $j$-th particles. Thus, the effective potentials of the three-body system in the channel space $|\mathcal{BBB}\rangle\vcentcolon=\{B_{a}^{\ast}B_{b}B_{c},B_{a}B_{b}^{\ast}B_{c},B_{a}B_{b}B_{c}^{\ast},B_{a}^{\ast}B_{b}^{\ast}B_{c},B_{a}^{\ast}B_{b}B_{c}^{\ast},B_{a}B_{b}^{\ast}B_{c}^{\ast}\}$ take the following form,

[TABLE]

Its graphical illustration is shown in Fig. 2.

III The break-up state and two-body subsystem

For the three bottomed meson system, suppose one of the constituent is infinitely away from the remaining two mesons. The system can be divided into a two-body subsystem plus a free meson. A bound state solution of the two-body subsystem indicates a break-up state for the three-body system, i.e. a two-body bound state plus a free meson. In the OPE model, as there is no direct interaction between two $B$ mesons, one could expect a break-up state with the subsystem $BB^{*}$ with quantum number $J^{P}=1^{+}$ and a free meson $B$. We can detach the subsystem $B_{a}B_{b}^{*}$ first, and explore its binding solution. The Hamiltonian of the subsystem in the channel space $|\mathcal{BBB}\rangle\vcentcolon=\{B_{a}^{\ast}B_{b}B_{c},B_{a}B_{b}^{\ast}B_{c},B_{a}B_{b}B_{c}^{\ast},B_{a}^{\ast}B_{b}^{\ast}B_{c},B_{a}^{\ast}B_{b}B_{c}^{\ast},B_{a}B_{b}^{\ast}B_{c}^{\ast}\}$ reads

[TABLE]

where the $T_{*}=-({1}/{2\mu_{*}})\nabla_{ab}^{2}$ and $T_{**}=-({1}/{2\mu_{**}})\nabla_{ab}^{2}$ are the relative kinetic energy for the $BB^{*}$ and $B^{*}B^{*}$ in their center-of-mass frame, respectively, with $\mu_{*}=({M_{B}M_{B^{*}}})/({M_{B}+M_{B^{*}}})$, $\mu_{**}=({M_{B^{*}}})/{2}$, $\nabla_{ab}^{2}=({1}/{r_{ab}})({d^{2}}/{dr_{ab}^{2}})r_{ab}-({\overrightarrow{L_{ab}}^{2}})/({r_{ab}^{2}})$. Here $\overrightarrow{L_{ab}}$ is the angular momentum operator between meson $a$ and $b$. We also have the mass gap $\delta M=M_{B^{*}}-M_{B}$. The effective potentials $V_{1}$, $V_{2}$, $V^{\prime}_{2}$ and $V_{3}$ depend on the isospin of the specific channels, thus we rewrite the $BB^{*}$ wave functions with fixed isospin

[TABLE]

For the specific channels $B_{a}^{*}B_{b}$, $B_{a}B_{b}^{*}$ and $B_{a}^{*}B_{b}^{*}$, the Schrödinger equation in the channel space $\{B_{a}^{*}B_{b},~{}B_{a}B_{b}^{*},~{}B_{a}^{*}B_{b}^{*}\}$ takes the form

[TABLE]

Based on this, we can derive the scattering amplitude at the tree level

[TABLE]

where the $T$-matrix is the interaction part of the $S$-matrix and the $M$ is defined as the invariant matrix element. In the second equation we have applied the first order of Born series expansion on the Lippmann-Schwinger equation with $V_{fi}$ being the effective potential. The relation between the scattering amplitude $M_{fi}$ and the potential $V_{fi}$ is

[TABLE]

where $p_{f(i)}$ and $m_{f(i)}$ denote the four-momentum and the mass of the final (initial) state.

In the calculation, $p_{1}(E_{1},\vec{p})$ and $p_{2}(E_{2},-\vec{p})$ denote the four-momenta of the initial state particles in the center-of-mass system shown in Fig. 3, while $p^{\prime}_{1}(E^{\prime}_{1},\vec{p^{\prime}})$ and $p^{\prime}_{2}(E^{\prime}_{2},-\vec{p^{\prime}})$ denote the four-momenta of the final state particles, respectively. $q=p^{\prime}_{1}-p_{1}=(E^{\prime}_{1}-E_{1},\vec{p^{\prime}}-\vec{p})=(E_{2}-E^{\prime}_{2},\vec{q})$ is the transferred four-momentum. For convenience, we always use $\vec{q}=\vec{p^{\prime}}_{1}-\vec{p}_{1}$ and $\vec{k}=(\vec{p^{\prime}}_{1}+\vec{p}_{1})/2$ instead of $\vec{p^{\prime}}$ and $\vec{p}$ in the calculations. The effective potential in coordinate space can be derived by Fourier transformation

[TABLE]

To take into account in a rough way the substructure of each vertex, a monopole form factor

[TABLE]

with $m_{\pi}$ the pion mass and

[TABLE]

is used to suppress the contribution from UV energies. Here, $\Delta M=M_{B}^{\ast}-M_{B}$ and $\Delta M^{\prime}=(M_{B}^{\ast}-M_{B})/2$. As the parameter $\Lambda$ is related to non-perturbative QCD, it cannot be well determined. Here we only explore its effect on the binding energy of the $BB^{*}$ with the quantum number $J^{P}=1^{+}$ system. To solve the time-independent Schrödinger equation in coordinate space, the potential $V(\vec{q},\vec{k})$ in momentum space can be transformed in to that in coordinate space as shown in the Appendix.

The isosinglet and isotriplet $BB^{*}$ potentials $V_{BB^{*}\rightarrow BB^{*}}(\vec{r})$ in coordinate space are shown in Fig. 4 (a) and (b), respectively, with $\Lambda=1440~{}\mathrm{MeV}$. The isosinglet potential $V_{BB^{*}\rightarrow BB^{*}}(\vec{r})$ is repulsive which does not indicate a bound solution. Nevertheless, the potential $V_{BB^{*}\rightarrow B^{*}B^{*}}(\vec{r})$ for the isosinglet is attractive as shown in Fig. 5. So there is still the possibility of a binding solution. On the contrary, the isotriplet potential $V_{BB^{*}\rightarrow BB^{*}}(\vec{r})$ is attractive, while its potential $V_{BB^{*}\rightarrow B^{*}B^{*}}(\vec{r})$ is repulsive. These potentials in coordinate space can be expressed as

[TABLE]

[TABLE]

[TABLE]

with $\tilde{m}^{(\prime)2}_{\pi}=m_{\pi}^{2}-\Delta M^{(\prime)2}$. The tensor operator $\hat{S}_{T}$ has the form $\hat{S}_{T}=3(\vec{r}\cdot\hat{\vec{\epsilon_{b}}})(\vec{r}\cdot\hat{\vec{\epsilon_{a}}}^{{\dagger}})-\hat{\vec{\epsilon_{b}}}\cdot\hat{\vec{\epsilon_{a}}}^{{\dagger}}$ with $\epsilon$ the polarization vector of $B^{\ast}$. The $C_{\pi}(i,j)$ are channel dependent coefficients, summarized in Table 2. The $c$ in Table 2 represents the C-parity of the corresponding channel.

Since the tensor operator $\hat{S}_{T}$ leads to S-D wave mixing, the contributions from D-wave should be taken into account. Thus the wave function $\Psi(\vec{r})$ has two parts

[TABLE]

with $\psi_{S}(\vec{r})$ and $\psi_{D}(\vec{r})$ the $S$-wave and $D$-wave functions, respectively. In the matrix method, we use Laguerre polynomials

[TABLE]

as a set of orthogonal basis with the normalization condition

[TABLE]

Thus the total wave function can be expanded as

[TABLE]

where $\phi_{S}$ and $\phi_{D}$ are the angular part of the spin and orbital wave function for the $S$-wave (${}^{3}S_{1}$) and $D$-wave (${}^{3}D_{1}$) states, respectively. $a^{(\prime)}_{i}$ and $b^{(\prime)}_{i}$ are the corresponding expansion coefficients of S-wave and D-wave, respectively. After solving the coupled-channel Schrödinger equation with the S-D wave mixing, we obtain the binding energy $E_{b}$ and its corresponding wave function $\Psi(\Lambda,\vec{r}_{ab})$ for a given parameter $\Lambda$. Thus the wave function has the form

[TABLE]

Here, the wave function $\Psi(\Lambda,\vec{r}_{ab})$ is normalized. If we choose the value of the parameter $\Lambda=1440$ MeV for instance, one finds a loosely bound state for the isospin triplet system with a binding energy of 5.08 MeV, when the quantum number is $J^{P}=1^{+}$. There is also a loosely bound state for the isospin singlet system, when the quantum number is $J^{P}=1^{+}$. If the value of the parameter is chosen at $\Lambda=1107.7$ MeV, the isospin singlet and triplet systems have the same binding energy of 5.08 MeV. The dependence of the binding energy on the parameter $\Lambda$ will be given in Tables 4-5 and discussed in Sec. VII.

IV Born-Oppenheimer potential

As discussed in Sec. II, the BO potential reflects the influence of one of the mesons on the dynamics of the other two. For the $BBB^{\ast}$ (labeled as $a$, $b$ and $c$) system, one can derive the BO potential from $a$ for the $bc$ system. The procedure is divided into the following three steps:

• •

Considering that the particle $b$ and $c$ are static with the separation $r_{bc}$, one can separate the degree of freedom of $a$ from the three-body system.

• •

We assume the distance $r_{bc}$ is a parameter. The mesons $b$ and $c$ are static, and have one-pion interactions with meson $a$, which can be viewed as two static sources.

• •

We explore the dynamics for the meson $a$ in the limit $r_{bc}\to\infty$, and subtract the binding energy for the break-up state which is trivial for the three-body bound state.

Within this scheme, we divide the motion of the system into two parts, one is the motion of the meson $a$ relative to the mesons $b$ and $c$. The other one is the relative motion between mesons $b$ and $c$ in the presence of the BO potential from $a$.

As illustrated in Figs. 6, we use $\vec{r}_{bc}$ to denote the relative displacement between $b$ and $c$. Further, $\vec{r}_{ab}$ and $\vec{r}_{ac}$ represent the displacement of the meson $a$ relative to the meson $b$ and $c$, respectively. One can separate the effective potentials for the meson $a$

[TABLE]

from the Eq. (12). The remaining part

[TABLE]

in Eq. (12) is the potential between $bc$. As discussed in the previous section, one can obtain the two-body binding energy

[TABLE]

The $\psi(\vec{r}_{ab})$ and $\psi^{\prime}(\vec{r}_{ab})$ in the above equation are the eigenstate wave functions in Eq. (31).

In OPE model, as the virtual pion can only be exchanged between two of the $BBB^{*}$ subsystems, the wave function of $a$ can be either $\frac{1}{\sqrt{2}}\psi(\vec{r}_{ab})|B_{a}^{*}B_{b}B_{c}\rangle+\frac{1}{\sqrt{2}}\psi(\vec{r}_{ab})|B_{a}B_{b}^{*}B_{c}\rangle+\psi^{\prime}(\vec{r}_{ab})|B_{a}^{*}B_{b}^{*}B_{c}\rangle$ with pion exchanged between $a$ and $b$ or $\frac{1}{\sqrt{2}}\psi(\vec{r}_{ac})|B_{a}^{*}B_{b}B_{c}\rangle+\frac{1}{\sqrt{2}}\psi(\vec{r}_{ac})|B_{a}B_{b}B_{c}^{*}\rangle+\psi^{\prime}(\vec{r}_{ac})|B_{a}^{*}B_{b}B_{c}^{*}\rangle$ with pion exchanged between $a$ and $c$. The final wave function for the meson $a$ should be the superposition of these two components

[TABLE]

For simplicity, we neglect the mass difference for the $BB^{*}$ and $B^{*}B^{*}$ in the kinetic operator, i.e. $T_{**}\approx T_{*}$. Then the Hamiltonian of the meson $a$ is

[TABLE]

Accordingly, one can obtain the energy eigenvalue of the meson $a$

[TABLE]

where in the second step Eq. (53) and the symmetry between $b$ and $c$ are used. Since both the two-body energy eigenvalue $E_{2}$ and the wave functions $\psi_{b}$, $\psi_{c}$ depend on the parameter $\Lambda$, $E_{a}$ is also a function of $\Lambda$.

We take the parameter $\Lambda=1440$ MeV as an example and plot $E_{a}$ for the isospin triplet of the $BB^{*}$ system in Figs. 7(a). As shown in the figure, the energy of the meson $a$ has a minimum $-17.65$ MeV when $r_{bc}=0$, which corresponds to the limit that the mesons $b$ and $c$ are on top of each other and the system is reduced to the $bc$-$a$ quasi-two-body system. When $r_{bc}\to\infty$, then $E_{a}$ tends to the two-body energy eigenvalue $E_{2}$, i.e.$-5.08$ MeV. This corresponds to the situation that the meson $b$ is infinitely far away from the meson $c$. Then the meson $a$ can only form a two-body bound state with either $b$ or $c$. It is not a three-body bound state anymore, but rather a two-body bound state plus a free meson state. In fact, this is nothing but the break-up state that we have discussed in the earlier sections. We also plot $E_{a}$ for the isospin triplet of the $BB^{*}$ system in Fig. 8(a), taking the parameter $\Lambda=1107.7$ MeV. Similarly to the above, $E_{a}$ tends to the two-body energy eigenvalue $-5.08$ MeV. Therefore, we should subtract the limiting value $E_{2}$ when investigating the three-body bound state for the $BBB^{*}$ system. We define the BO potential as

[TABLE]

In other words, the BO potential between $b$ and $c$ is the energy eigenvalue of the meson $a$ relative to that of the break-up state.

V The configurations of the three-body systems

In the OPE model, there is only one pion exchanged between any two constituents in the $BBB^{\ast}$ system. The constituents will change themselves from vector mesons into pseudo-scalar mesons or vice verse when they exchange one pion. Each constituent has the same probability to be a vector meson or a pseudo-scalar meson. Thus, the symbol $\ast$ is shared among them. Since only one virtual pion occurs in the $BBB^{\ast}$ molecule, the virtual pion also be shared by the three mesons. We can thus write the $BBB^{\ast}$ as $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$.

The BO potential can describe the contribution for the one meson on the dynamics of the two remaining mesons as we have discussed in the last section. Assuming that the meson $b$ and $c$ are much heavier than the meson $a$, then we can use the Born-Oppenheimer approximation to separate the degree of freedom of $a$ from the three-body system. In other words, it is a kind of an adiabatic approximation that we divide the degrees of freedom of the three-body system into a light one and a heavy one. The motion of the light degree of freedom is the motion of meson $a$ relative to the three-body centre of mass. The motion of the heavy degree of freedom is the relative motion between meson $b$ and $c$. When exploring the dynamics for the meson $a$, we can assume the meson $b$ and $c$ are static with the distance $r_{bc}$. Then the three-body system can be simplified as a two-body system consisting of mesons $b$ and $c$ but with an additional BO potential generated by the meson $a$. Overall, only the meson $a$ can be separated from the system due to the fact that this meson is much lighter than the other ones. A separation in this way can be a good approximation for this system. With the same procedure that we derived Eq. (54), we obtain the wave functions $\psi(\vec{r_{ab}},\vec{r_{ac}})$ for the meson $a$. The remaining degree of freedom is the relative motion between meson $b$ and $c$ that can be described by a wave function assumed as $\Phi(\vec{r_{be}})$, to be determined from the Schrödinger equation. Then the total wave function of the system then has the form

[TABLE]

Nevertheless, the true system $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ is that the three mesons have little mass difference. Every meson can be considered to be a lighter one and separated from the three-body system. Thus, the system has the three basic simplification schemes. That is we can divide the system $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ into three kinds of two-body subsystems, i.e., $B_{a}^{(\ast)}B_{b}^{(\ast)}$ with the BO potential created by the meson $c$, $B_{b}^{(\ast)}B_{c}^{(\ast)}$ with the BO potential created by the meson $a$ and $B_{a}^{(\ast)}B_{c}^{(\ast)}$ with the BO potential created by the meson $b$ as shown in Figs. 9. These three simplification schemes can be regarded as three kinds of basic configurations. The eigenstates of the three-body system should be combinations of them. As the most simplest combination, one might expect the three-body eigenstate should be the superposition of the three kinds of basic configurations. We use the $\psi_{\not{a}}$, $\psi_{\not{b}}$, $\psi_{\not{c}}$ to denote these three configurations. The configuration wave function $\psi_{\not{a}}$ represent the configuration that we omit the meson $B_{a}^{(\ast)}$ and add the corresponding BO potential instead. Similarly, $\psi_{\not{b}}$, $\psi_{\not{c}}$ denote the configurations with the BO potentials provided by the mesons $B_{b}^{(\ast)}$ and $B_{c}^{(\ast)}$, respectively.

Taking the configuration function $\psi_{\not{a}}$ as an example, we separate the motion of the $B_{a}^{(\ast)}$ relative to the other mesons $B_{b}^{(\ast)}$ and $B_{c}^{(\ast)}$ where their relative displacement $r_{bc}$ is regarded as a parameter as shown in Fig. 9(a). The wave function of the $B_{a}^{(\ast)}$ has been discussed in the last section and can be written as $\psi(\vec{r}_{ab},\vec{r}_{ac})$. The remaining degree of freedom is the relative motion between $B_{b}^{(\ast)}$ and $B_{c}^{(\ast)}$, which can be taken as $\Phi(\vec{r}_{bc})$. Thus we have the configuration function $\psi_{\not{a}}=\Phi(\vec{r}_{bc})\psi(\vec{r}_{ab},~{}\vec{r}_{ac})$. The other two wave functions $\psi_{\not{b}}$ and $\psi_{\not{c}}$ can be obtained analogously, i.e. $\psi_{\not{b}}=\Phi(\vec{r}_{ac})\psi(\vec{r}_{ab},~{}\vec{r}_{bc})$, $\psi_{\not{c}}=\Phi(\vec{r}_{ab})\psi(\vec{r}_{bc},~{}\vec{r}_{ac})$ which correspond to the Fig. 9(b) and Fig. 9(c), respectively. If we regard the three configuration functions as a set of basis states, then the basis constitute a configuration space $\{\psi_{\not{a}},~{}\psi_{\not{b}},~{}\psi_{\not{c}}\}$. The three-body eigenstate expressed as a superposition of the three kinds of basic configurations can be described as a state vector in this configuration space. Thus, as an interpolating wave function, the three-body wave functions can be written as

[TABLE]

where the $\Phi(\vec{r}_{bc})$, $\Phi(\vec{r}_{ac})$ and $\Phi(\vec{r}_{ab})$ are undetermined functions that need to be solved. The $\alpha$, $\beta$ and $\gamma$ are the expansion coefficients. According to Eq. (54), we rewrite the three basic configuration functions in the channel space $\{B_{a}^{\ast}B_{b}B_{c},B_{a}B_{b}^{\ast}B_{c},B_{a}B_{b}B_{c}^{\ast},B_{a}^{\ast}B_{b}^{\ast}B_{c},B_{a}^{\ast}B_{b}B_{c}^{\ast},B_{a}B_{b}^{\ast}B_{c}^{\ast}\}$ as

[TABLE]

which can be expanded as a set of Laguerre polynomials

[TABLE]

Here the subscript $i$ is the order of Laguerre polynomials. We define the $i^{\rm th}$ order of the configuration functions as $\psi_{\not{a}}^{i}=\phi_{i}(\vec{r}_{bc})\psi(\vec{r}_{ab},~{}\vec{r}_{ac})$, $\psi_{\not{b}}^{i}=\phi_{i}(\vec{r}_{ac})\psi(\vec{r}_{ab},~{}\vec{r}_{bc})$ and $\psi_{\not{c}}^{i}=\phi_{i}(\vec{r}_{ab})\psi(\vec{r}_{bc},~{}\vec{r}_{ac})$. Further, $C$ is a normalisation constant.

We expect the three-body bound state that we seek for can be expressed as a state vector in the configuration space $\{\psi_{\not{a}},~{}\psi_{\not{b}},~{}\psi_{\not{c}}\}$. However, the configuration functions in Eq. (66) are not an orthogonal basis. Thus we orthonormalize the $\{\psi_{\not{a}},~{}\psi_{\not{b}},~{}\psi_{\not{c}}\}$ into a new basis $\{\tilde{\psi}_{\not{a}},~{}\tilde{\psi}_{\not{b}},~{}\tilde{\psi}_{\not{c}}\}$. We use $\tilde{\psi}_{\not{a}}^{i}$, $\tilde{\psi}_{\not{b}}^{i}$ and $\tilde{\psi}_{\not{c}}^{i}$ to denote the $i^{\rm th}$th order of the new configuration functions $\tilde{\psi}_{\not{a}}$, $\tilde{\psi}_{\not{b}}$ and $\tilde{\psi}_{\not{c}}$, respectively. Then we have

[TABLE]

where the $x_{ij}$ is a parameter matrix which will be determined later. The $N_{i}$ are normalization coefficients. The parameter matrix $x_{ij}$ in the three configuration functions are the same due to the interchange symmetry for the $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ system.

Since the $i^{\rm th}$ order configuration function $\tilde{\psi}_{\not{a}}^{i}$ should be orthogonal with the any order of the other configuration function $\tilde{\psi}_{\not{a}}^{j}$, one can gets the orthogonalization condition

[TABLE]

which gives

[TABLE]

This equation will determine the parameter matrix $x_{ij}$. Considering the normalization of the $i^{\rm th}$ order configuration function $\tilde{\psi}_{\not{a}}^{i}$

[TABLE]

one can obtain the normalization equation for the $N_{i}$ as

[TABLE]

After solving the equations for $x_{ij}$ and $N_{i}$, we obtain an orthonormalized configuration basis. This basis constitutes a orthonormalized configuration space. Then the eigenvector for the three-body system $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ can be written as a vector in the configuration space $\{\tilde{\psi}_{\not{a}},~{}\tilde{\psi}_{\not{b}},~{}\tilde{\psi}_{\not{c}}\}$. Therefore, we have

[TABLE]

where the $\tilde{\alpha}_{i}$, $\tilde{\beta}_{i}$ and $\tilde{\gamma}_{i}$ are the $i^{\rm th}$ order expansion coefficients.

VI Three-body Schrödinger equation

As we have discussed in previous sections, if the three-body binding energy is below the break-up threshold, the three-body system will disintegrate into a two-body system and a free meson. Since we only focus on the three-body bound state, we could make an energy shift and remove the energy eigenvalue $E_{2}$ for the break-up state and define a reduced Hamiltonian for the three-body system as

[TABLE]

The explicit form of $H$ is

[TABLE]

where $T_{*}=-({1}/{2\mu_{*}})\nabla_{ab}^{2}$, $T=-({1}/{2\mu})\nabla_{ab}^{2}$, $T_{**}=-({1}/{2\mu_{**}})\nabla_{ab}^{2}$, $T^{\prime}_{*}=-({1}/{2\mu^{\prime}_{*}})\nabla_{\xi}^{2}$, $T^{\prime}=-({1}/{2\mu^{\prime}})\nabla_{\xi}^{2}$, $T^{\prime}_{**}=-({1}/{2\mu^{\prime}_{**}})\nabla_{\xi}^{2}$ are the kinetic energy operators and the corresponding reduced masses are $\mu_{*}=({M_{B}M_{B^{*}}})/({M_{B}+M_{B^{*}}})$, $\mu={M_{B}}/{2}$, $\mu_{**}={M_{B^{*}}}/{2}$, $\mu^{\prime}_{*}=({(M_{B}+M_{B^{*}})M_{B}})/({2M_{B}+M_{B^{*}}})$, $\mu^{\prime}=({2M_{B}M_{B^{*}}})/({2M_{B}+M_{B^{*}}})$, $\mu^{\prime}_{**}=({2M_{B^{*}}M_{B}})/({2M_{B^{*}}+M_{B}})$. Here $\nabla_{ab}^{2}=({1}/{r_{ab}})({d^{2}}/{dr_{ab}^{2}})r_{ab}-({\overrightarrow{L}_{ab}^{2}}/{r_{ab}^{2}})$ and $\nabla_{\xi}^{2}=({1}/{\xi})({d^{2}}/{d\xi^{2}})\xi-({\overrightarrow{L}_{\xi}^{2}}/{\xi^{2}})$ with $\vec{\xi}={\vec{r}_{ab}}/{2}-\vec{r}_{bc}$. $\vec{r_{bc}}$ is the direction of the meson $b$ relative to the meson $c$. $\overrightarrow{L}_{ab}$ is the angular momentum operator between mesons $a$ and $b$. $\overrightarrow{L}_{\xi}$ is the relative angular momentum operator between two-body centre of mass for the meson $a~{}b$ and the meson $c$. The mass gap is $\delta M=M_{B^{*}}-M_{B}$.

The total Hamiltonian for the three-body system in the configuration space $\{\tilde{\psi}_{\not{a}},~{}\tilde{\psi}_{\not{b}},~{}\tilde{\psi}_{\not{c}}\}$ can be written as

[TABLE]

with $H_{\not{m}\not{n}}=\langle\tilde{\psi}_{\not{m}}|H|\tilde{\psi}_{\not{m}}\rangle$ $(m,~{}n=a,~{}b,~{}c)$.

The total reduced Hamiltonian for the three-body system $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ in the configuration space $\{\tilde{\psi}_{\not{a}},~{}\tilde{\psi}_{\not{b}},~{}\tilde{\psi}_{\not{c}}\}$ can be expressed as

[TABLE]

with $\mathcal{H}_{\not{m}\not{n}}=\langle\tilde{\psi}_{\not{m}}|\mathcal{H}|\tilde{\psi}_{\not{m}}\rangle$ $(m,~{}n=a,~{}b,~{}c)$. Thus we have

[TABLE]

The matrix element of the $\mathcal{H}_{\not{a}\not{a}}$ can be written as

[TABLE]

where, in the last step, the interchange symmetry in the $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ system is used. Similarly, we also have

[TABLE]

There are two independent matrices

[TABLE]

[TABLE]

to be determined, where we have used the abbreviations $\phi_{ab}^{i}$, $\phi_{bc}^{i}$, $\phi_{ac}^{i}$, $\psi_{ab}$, $\psi_{bc}$, $\psi_{ac}$ for $\phi(\vec{r}_{ab})^{i}$, $\phi(\vec{r}_{bc})^{i}$, $\phi(\vec{r}_{ac})^{i}$, $\psi(\vec{r}_{ab})$, $\psi(\vec{r}_{bc})$, $\psi(\vec{r}_{ac})$ respectively. The expression for the $\mathcal{H}_{\not{c}\not{a}}$ can be easily obtained by the replacement $c\rightarrow b,~{}b\rightarrow c$ on the expression for the $\mathcal{H}_{\not{b}\not{a}}$. Similarly, the expression for the $\mathcal{H}_{\not{c}\not{b}}$ is obtained by the replacement $c\rightarrow b,~{}b\rightarrow a,~{}a\rightarrow c$ on the expression for the $\mathcal{H}_{\not{b}\not{a}}$. In fact, interchange invariance for the $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ system can simplify the calculation, i.e. $\mathcal{H}_{\not{c}\not{a}}=\mathcal{H}_{\not{c}\not{b}}=\mathcal{H}_{\not{b}\not{a}}$ and $\mathcal{H}_{\not{a}\not{a}}=\mathcal{H}_{\not{b}\not{b}}=\mathcal{H}_{\not{c}\not{c}}$.

Based on the above discussion, the three-body Schrödinger equation can finally be written as

[TABLE]

where the energy eigenvalue $E_{3}$ is the reduced three-body energy eigenvalue. The total energy eigenvalue relative to the $BBB^{\ast}$ mass threshold is $E_{T}=E_{3}+E_{2}$. Solving the three-body Schrödinger equation may partly answer whether the three-body system has a loosely bound state or not.

VII Application to the $NNN$ system

In order to verify the feasibility of the Born-Oppenheimer potential method for the three heavy system, we apply it to the three nucleon system. Since there is sufficient experimental data for this system, we can apply the formalism introduced above to investigate its binding energy and illustrate the feasibility of our formalism. As we know, the triton and the helium-3 (3He) nucleus are the two possible bound states of the $NNN$ system, both of them have the quantum numbers $I(J^{P})=\frac{1}{2}({\frac{1}{2}}^{+})$ but have different isospin on its $z$ direction. They have the same structure and the binding energy if isospin symmetry breaking effect is neglected. The calculation on the three-nucleon system is much more straightforward than the $BBB^{*}$ system, as there are no other coupled channels. For simplicity, we only write down the isospin wave functions of the triton and helium-3 nuclei, which are

[TABLE]

The Lagrangian reads

[TABLE]

where the $g_{N}=14.70$ is the coupling constants (we use here the pseudoscalar coupling, which is fine to the order we are working, see e.g. Ref. Bernard:1989fe ). $N=(\psi_{p},\psi_{n})$ is the nucleon doublet. Further, $\vec{\tau}=\{\tau_{1},\tau_{2},\tau_{3}\}$ are the Pauli matrices, and $\vec{\pi}=\{\frac{1}{\sqrt{2}}(\pi^{+}+\pi^{-}),\frac{i}{\sqrt{2}}(\pi^{+}+\pi^{-}),\pi^{0}\}$ are the $\pi$ fields. With a procedure similar to the one discussed in the Sec. II-VI, we can investigate the properties of the break-up state formed by a deuteron and a free nucleon as well as the three-body bound states. As discussed in the above sections, there is only one free parameter $\Lambda$ in the monopole form factor introduced in Sec. III, which reflects, in a rough way, the internal structure of the interacting hadrons. In other words, the size of hadron is proportional to $1/\Lambda$, which is still unknown from the fundamental theory. Thus the parameter $\Lambda=899.60~{}\mathrm{MeV}$ is fixed by the binding energy $E_{2}=2.23$ MeV of deuteron in our calculations. With the so determined parameter $\Lambda$, we can obtain the BO potentials for the $NNN$ system using the formalism in the Sec. II-VI, with just the replacement of the effective potential $V_{BB*}$ by $V_{NN}$ in the calculations. This potential reads

[TABLE]

where the $C_{\pi}^{NN}(i,j)$ is the channel-dependent coefficient for the two-nucleon system, $M_{N}$ is the mass of nucleon, $g_{N}$ is the pion-nucleon coupling constant, and $\sigma_{1}$ and $\sigma_{2}$ are the spin Pauli matrix for the nucleon 1 and 2 in the scattering process $1+2\rightarrow 3+4$.

After solving the three-body Schrödinger equation, i.e.Eq. (123), one can obtain the dependence of the binding of the three-nucleon system on the parameter $\Lambda$ (Table 3). As shown in the table, there is a three-body bound state with the reduced binding energy and the total three-body bound energy in the range of $1.93$-$20.51$ $\mathrm{MeV}$ and $2.11$-$31.20$ $\mathrm{MeV}$, respectively, when the parameter $\Lambda$ varies from $830~{}\mathrm{MeV}$ to $1040~{}\mathrm{MeV}$. The corresponding isospin singlet two-body subsystem $NN$ has the binding enrgy in the range of $0.18$-$10.69$ $\mathrm{MeV}$. The root-mean-square of the system decreases from 4.21 $\mathrm{fm}$ to 2.10 $\mathrm{fm}$ when the parameter increases. Once the parameter $\Lambda=899.60~{}\mathrm{MeV}$ fixed by the deuteron binding energy, the reduced three-body binding energy and the total binding energy relative to the three free nucleons are 5.38 $\mathrm{MeV}$ and 7.62 $\mathrm{MeV}$, respectively. The later is comparable with the empirical binding energies of the triton (8.48 $\mathrm{MeV}$) and helium-3 (7.80 $\mathrm{MeV}$) nuclei. Noet again that there is no numerical difference between the binding energies of triton and helium-3 in the calculation, as the isospin breaking has not be considered.

For a better illustration of the binding property, we plot the dependence of the reduced three-body binding energy on the two-body binding energy of its deuteron subsystem. As shown in Fig. 10, the binding energy of the three-nucleon system becomes larger when the binding energy of its subsystem $NN$ increases. There are two red points in the figure, the left red point is the critical point which indicates the lower limit of the required binding energy of the deuteron to form a three-body bound state. It is very interesting that even though the binding energy of its subsystem is zero, there is a small binding energy of the three-nucleon system, which is 1.71 $\mathrm{MeV}$. This is reminsicent of a Borromean state, where a three-body system may have a bound state despite the fact that none of its subsystems forms a bound state. The other red point is our numerical result of the binding energy of triton or helium-3. It is a little below the experimental values since in our calculations we use the BOP method to construct our interpolating wave functions, which can be regarded as a version of the variational principle. As we know, this always give an upper limit of the energy of a system.

VIII Numerical results on the $BBB^{*}$ system

The application of the BOP method to the three-nucleon system has verified its feasibility to some extent. Now we return to the system mainly discussed in this paper, i.e. the three $B$ mesons system $BBB^{*}$. There is only one free parameter $\Lambda$ in the monopole form factor which is undetermined in our calculations. For the deuteron case, the parameter $\Lambda$ is within the range $0.8\sim 1.5~{}\mathrm{GeV}$. One would expect that the size of heavier bottom system is smaller than the size of deuteron, leading to a larger $\Lambda$. Thus, we vary the parameter $\Lambda$ from $0.9~{}\mathrm{GeV}$ to $1.6~{}\mathrm{GeV}$ to study whether the $BBB^{\ast}$ system is bound or not.

In order to show the properties of the two-body interactions for the $BB^{\ast}$, we first present the numerical results for the break-up state in Tables 4-5. We plot the effective potentials for the $BB^{\ast}$ in Fig. 4 and Fig. 5, where the regularization parameter is fixed at 1440 MeV. In these figures, (a) and (b) correspond to the isospin $I=0$ and $I=1$ cases, respectively. After carefully solving the coupled-channel Schrödinger equation with the treatment of the S-D wave mixing, we find loosely bound states for both cases.

For the isospin triplet case, i.e. $I_{2}=1$, the dependence of the binding energy of the two-body $BB^{*}$ system on the regularization parameter $\Lambda$ is shown in Table 4. The energy threshold of the break-up state for the $BBB^{\ast}$ is just the two-body energy eigenvalue of the $BB^{*}$ plus the mass of the three static free meson. We use $E_{2}$ denote the energy eigenvalue of the $BB^{*}$. When the parameter $\Lambda$ varies from $1380~{}\mathrm{MeV}$ to $1560~{}\mathrm{MeV}$, there is a bound state solution with the binding energy $2.11\sim 15.59~{}\mathrm{MeV}$ and the root-mean-square radius $1.51\sim 0.67~{}\mathrm{fm}$. The $S$-wave component takes over $99.13\%\sim 99.36\%$ comparing to the value $0.87\%\sim 0.64\%$ for $D$-wave. The $BB^{*}$ and $B^{*}B^{*}$ channels have probabilities $91.91\%\sim 76.50\%$ and $8.09\%\sim 23.50\%$, respectively. The proportion of the $B^{*}B^{*}$ channel and D-wave component are relatively small. However, as the value of the regularization parameter $\Lambda$ increases, the $B$ and $B^{*}$ interacting with pion is more like a point particle, the proportion of the $B^{*}B^{*}$ channel increases greatly, while the D-wave component decreases. We plot the radial wave functions of the $S$-wave and $D$-wave Fig. 11(a) for the system $BB^{\ast}$, obviously, the bound state we have found is ground state.

For the isospin singlet case, i.e. $I_{2}=0$, the dependence of the binding energy of the two-body $BB^{*}$ system on the regularization parameter $\Lambda$ is shown in Table 5. We also use $E_{2}$ to denote the energy eigenvalue of the $BB^{*}$. When the parameter varies from $1040~{}\mathrm{MeV}$ to $1220~{}\mathrm{MeV}$, there is a bound state solution with binding energy $1.88\sim 14.78~{}\mathrm{MeV}$ and the root-mean-square radius $2.03\sim 0.95~{}\mathrm{fm}$. The $S$-wave component is $86.09\%\sim 77.21\%$ compared to the value $13.91\%\sim 22.79\%$ for $D$-wave. The $BB^{*}$ and $B^{*}B^{*}$ channels have probabilities $92.05\%\sim 75.09\%$ and $7.95\%\sim 24.91\%$, respectively. The proportion of the $B^{*}B^{*}$ channel and D-wave component are relatively small which is similar with the case of isospin triplet. As the value of the regularization parameter $\Lambda$ increases, the proportion of the $B^{*}B^{*}$ channel increases greatly. Different from the case of isospin triplet the S-wave component decrease and D-wave increase as $\Lambda$ increases. As the parameter $\Lambda$ increases, all of the effective potentials become stronger. The S-wave potential increases faster than the D-wave potential for the isospin triplet case, while it is reverse for the isospin singlet case. In order to check whether the bound state we have found is the ground state, we also plot the radial wave functions of the $S$-wave and $D$-wave Fig. 11(b) for the system $BB^{\ast}$.

In Sec. II, we have listed the isospin wave functions of the $BBB^{*}$ which are expressed as $|I_{2},I_{3},I_{3z}\rangle$. After solving the three-body Schrödinger equation via the method of Sec. V, we find that all of these isospin eigenstates have bound state solutions. As long as the two-body system $BB^{*}$ has a loosely bound state, the three-body system $BBB^{*}$ is most likely to have a loosely bound state, too. We have collected the dependence of the three-body bound state solutions on the two-body binding energy in Tables 6-7.

The bound state solutions for the state $|1,\frac{3}{2},\pm\frac{1}{2}(\pm\frac{3}{2})\rangle$ are shown in Table 6. The three-body binding energy relative to their break-up states is 5.67 MeV, when the parameter $\Lambda$ is chosen at 1440 MeV and the two-body binding energy of their subsystems $BB^{\ast}$ is 5.08 MeV. To search for the dependence on the binding energy of the two-body system $E_{2}$, we change the parameter $\Lambda$. It turns out that if the value of the $E_{2}$ varies from $-0.18~{}\mathrm{MeV}$ to $-17.97~{}\mathrm{MeV}$, then the reduced three-body energy eigenvalue $E_{3}$ decreases from $-0.19~{}\mathrm{MeV}$ to $-18.99~{}\mathrm{MeV}$ and the total three-body energy eigenvalue $E_{T}$ decreases from $-0.38~{}\mathrm{MeV}$ to $-36.95~{}\mathrm{MeV}$. The structure of the three-body bound state is a regular triangle with the root-mean-square length of one side decreasing from 3.98 fm to 0.65 fm. In order to illustrate the strength of the BO potential, we also collect its minimum $V_{BO}(0)$ in the table within the range of -3.43$\sim$-37.24 MeV as the $E_{2}$ increases. As $E_{2}$ increases, the effective attraction between $B$ and $B^{*}$ becomes stronger, the BO potential becomed deeper, so then the three-body system becomes tighter and has a larger binding energy. From the results in the table, we can also see that the dominant wave between any two $B^{(*)}$ in the $BBB^{\ast}$ is $S$ wave and the dominant channel is the $BBB^{\ast}$ instead of the $BB^{\ast}B^{\ast}$ channel. For comparison, we plot the wave functions for any two $B^{(*)}$ in the $BBB^{\ast}$ system and that for the two-body $BB^{*}$ system in Fig. 11(a) with $\Lambda=1440$ MeV. The shapes of these exhibit little difference. From another perspective, one more $B$ meson has little effect on the size of the system but greatly increases the binding energy.

For the state $|0,\frac{1}{2},\pm\frac{1}{2}\rangle$, we also find a loosely bound solution, which is shown in Table 7. The three-body binding energy relative to their break-up states is 7.18 MeV, when the parameter $\Lambda$ is chosen at 1107.7 MeV and the two-body binding energy of their subsystems $BB^{\ast}$ is 5.08 MeV. In order to show the dependence on the binding energy of the two-body system $E_{2}$, we also change the parameter $\Lambda$. We find that if the value of the $E_{2}$ varies from $-0.19~{}\mathrm{MeV}$ to $-17.10~{}\mathrm{MeV}$, then the reduced three-body energy eigenvalue $E_{3}$ decreases from $-0.32~{}\mathrm{MeV}$ to $-20.96~{}\mathrm{MeV}$ and the total three-body energy eigenvalue $E_{T}$ decreases from $-0.51~{}\mathrm{MeV}$ to $-38.06~{}\mathrm{MeV}$. The structure of the three-body bound state is a regular triangle with the root-mean-square length of one side decreasing from 3.89 fm to 0.93 fm. As an illustration for the strength of the BO potential, we also list its minimum $V_{BO}(0)$ in the table within the range of $-$2.15$\sim-30.15$ MeV as $E_{2}$ increases. Similar to the $I_{3}=\frac{3}{2}$ case, the dominant wave between any two $B^{(*)}$ in the $BBB^{\ast}$ is $S$ wave and the dominant channel is the $BBB^{\ast}$ one. In order to show that one more $B$ meson has little effect on the size of the system, we also plot the wave functions for any two $B^{(*)}$ in the $BBB^{\ast}$ system and that for the two-body $BB^{*}$ system in Fig. 11(b) with $\Lambda=1107.7$ MeV. Here we chose the parameter $\Lambda=1107.70$ MeV for a better comparison with the $|1,\frac{3}{2},\pm\frac{1}{2}(\pm\frac{3}{2})\rangle$ case, since both cases have the same two-body binding energy 5.08 MeV.

The state $|1,\frac{1}{2},\pm\frac{1}{2}\rangle$ also has a loosely bound solution. However, in our calculation the states $|1,\frac{1}{2},\pm\frac{1}{2}\rangle$ and $|1,\frac{3}{2},\pm\frac{1}{2}(\pm\frac{3}{2})\rangle$ are degenerate. This is due to the fact that in the OPE model only two-body interactions are considered, and the two-body interaction is only depends on the total isospin of the two interacting mesons. The states $|1,\frac{1}{2},\pm\frac{1}{2}\rangle$ and $|1,\frac{3}{2},\pm\frac{1}{2}(\pm\frac{3}{2})\rangle$ have the same two-body interaction but may have different three-body interaction. If we further consider the calculation to the next-to-next leading order, this degeneracy may disappear. The calculation that contains three-body interactions via pion exchange is quite complicated, which is left for the further work.

The numerical results show that the three-body binding energy $|E_{3}|$ increases as the two-body binding energy $|E_{2}|$ increases. One may wonder whether there is a critical value of the $|E_{2}|$, below which the three-body system has no bound state solution. After lots of calculations, it turns out that all of the isospin eigenstates have no such critical value. That is to say, no matter what a small value for the two-body binding energy $|E_{2}|$, as long as the two-body system $BB^{*}$ has a loosely bound state, the three-body system $BBB^{*}$ probably has a loosely bound state. To show this conclusion explicitly, we plot the dependence of the three-body binding energy on the variety of two-body binding energy in Fig. 12. The isospin eigenstate $|1,\frac{3}{2},\pm\frac{1}{2}(\pm\frac{3}{2})\rangle$ and $|1,\frac{1}{2},\pm\frac{1}{2}\rangle$ cases are shown in Fig. 12(a), where $E_{I=3/2}^{BBB^{*}}$ and $E_{I=1}^{BB^{*}}$ denote the reduced three-body binding energy and two-body binding energy, respectively. When the two-body binding energy $E_{I=1}^{BB^{*}}$ approaches 0 MeV, the reduced three-body binding energy $E_{I=3/2}^{BBB^{*}}$ approaches a small value of about 0.06 MeV. Similarly, we also plot the dependence curve for the $|0,\frac{1}{2},\pm\frac{1}{2}\rangle$ case in Fig. 12(b), where $E_{I=1/2}^{BBB^{*}}$ and $E_{I=0}^{BB^{*}}$ denote the reduced three-body binding energy and two-body binding energy in this case, respectively. As shown in the figure, the reduced three-body binding energy $E_{I=1/2}^{BBB^{*}}$ also has a small value of about 0.12 MeV, when the two-body binding energy $E_{I=0}^{BB^{*}}$ approaches zero.

IX Summary and Discussion

In the present paper, we have performed an extensive study on the possibility of the $BBB^{\ast}$ system to form tri-meson molecules. Based on the Born-Oppenheimer potential method as well as the OPE scheme, we derived the three-body Schrödinger equation for the system $BBB^{\ast}$. Since the regularization parameter $\Lambda$ is difficult to be pinned down, we choose the parameter in the range of 0.9$\sim$1.6 GeV and show various bound state solutions of the $BBB^{\ast}$ system. After careful treatments of the S-D wave mixing and the coupled-channel $BB^{\ast}B^{\ast}$ effects, we found that all of the isospin eigenstates expressed by the $|I_{2},I_{3},I_{3z}\rangle$ have bound state solutions. For instance, in the states $|1,\frac{3}{2},\pm\frac{1}{2}(\pm\frac{3}{2})\rangle$ and $|1,\frac{1}{2},\pm\frac{1}{2}\rangle$, the three-body binding energy relative to their break-up states is 5.67 MeV, when the parameter $\Lambda$ is chosen at 1440 MeV and the two-body binding energy of their subsystems $BB^{\ast}$ is 5.08 MeV. In the state $|0,\frac{1}{2},\pm\frac{1}{2}\rangle$, the three-body binding energy relative to their break-up states is 7.18 MeV, when the parameter $\Lambda$ is chosen at 1107.7 MeV and the two-body binding energy of their subsystems $BB^{\ast}$ is also 5.08 MeV. After careful calculations, we find no critical value for the two-body binding energy, which indicates the lower limit of the required binding energy of their subsystem $BB^{\ast}$ to form a three-body bound state. That is to say, no matter how small the two-body binding energy is, as long as the two-body subsystem $BB^{*}$ has a loosely bound state, the three-body system $BBB^{*}$ is most likely to have a loosely bound state, too.

The BO potential method we have used in this paper is an adiabatic approximation that divide the degrees of freedom of the motion for the three-body system into a light one and a heavy one. Then we simplify the three-body system into a two-body system only with heavy degrees of freedom but with an additional BO potential generated by the relative light meson. Since the system $B_{a}^{(\ast)}B_{b}^{(\ast)}B_{c}^{(\ast)}$ has little mass difference on its constituents, the motion of every constituent can be regarded as a light degree of freedom. Therefore, the eigenstates of the three-body system should be the combinations of all of the possible cases. As the most simplest combination, one might expect the three-body eigenstate should be a superposition all of the possible cases. In other words, the three-body bound state solutions we have listed in the last section are approximate solutions. It may be that the strict solutions will be a more complicated combinations. To answer this question requires further study.

Our calculations are based on the OPE scheme, which is leading order in the chiral power counting (neglecting contact interactions). Since only one virtual pion occurs in the $BBB^{\ast}$ molecule, the virtual pion is also shared by the three mesons. Therefore, the three constituents in the $BBB^{\ast}$ system share one virtual pion which corresponds to a delocalized pion bond. It is attractive and strong enough to make them form a three-body molecular state.

To summarize briefly, with the delicate efforts of the long-range one-pion exchange, the S-D wave mixing and coupled-channel effects, we have investigated the existence of the loosely bound tri-meson molecules $BBB^{\ast}$ and find that it is very easy to form a tri-meson molecular state as long as its two-body subsystem $BB^{*}$ has a molecular state. Hopefully, the present extensive investigations will be useful to the understanding of the few-body hadronic systems and the future well-developed experiments on hadron collisions will provide us with a platform to seek out the tri-meson molecules.

Acknowledgements

This work is supported in part by the DFG (Grant No. TRR110) and the NSFC (Grant No. 11621131001) through funds provided to the Sino-German CRC 110 “Symmetries and the Emergence of Structure in QCD”. The work of UGM was also supported by the Chinese Academy of Sciences (CAS) President’s International Fellowship Initiative (PIFI) (Grant No. 2018DM0034) and by VolkswagenStiftung (Grant No. 93562). QW is also supported by the Thousand Talents Plan for Young Professionals and research startup funding at SCNU.

Appendix: Some helpful functions and Fourier transformations

The functions $Y(\tilde{m}_{\pi}r)$ and $Z(\tilde{m}_{\pi}r)$ in Eqs. (25)-(27) are defined as

[TABLE]

[TABLE]

where

[TABLE]

Fourier transformation formulas read

[TABLE]

[TABLE]

[TABLE]

where $\hat{S}_{T}=3(\vec{r}\cdot\hat{\vec{\epsilon_{b}}})(\vec{r}\cdot\hat{\vec{\epsilon_{a}}}^{{\dagger}})-\hat{\vec{\epsilon_{b}}}\cdot\hat{\vec{\epsilon_{a}}}^{{\dagger}}$.

The polarization vector is the S-D wave space have the following substitution

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]