# Calculation of mass and width of unstable molecular state using the
developed Bethe-Salpeter theory

Xiaozhao Chen<EMAIL_ADDRESS>Department of Fundamental Courses, Shandong
University of Science and Technology, Taian, 271019, China  Xiaofu Lü
Department of Physics, Sichuan University, Chengdu, 610064, China  Institute
of Theoretical Physics, The Chinese Academy of Sciences, Beijing 100080, China
CCAST (World Laboratory), P.O. Box 8730, Beijing 100080, China  Xiurong Guo
Department of Fundamental Courses, Shandong University of Science and
Technology, Taian, 271019, China  Zonghua Shi  Department of Fundamental
Courses, Shandong University of Science and Technology, Taian, 271019, China
Qingbiao Wang  Department of Fundamental Courses, Shandong University of
Science and Technology, Taian, 271019, China

###### Abstract

Applying the developed Bethe-Salpeter theory for dealing with resonance, we
investigate the time evolution of molecular state composed of two vector
mesons as determined by the total Hamiltonian. Then exotic meson resonance
$\chi_{c0}(3915)$ is considered as a mixed state of two unstable molecular
states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$, and the mass and width for
physical resonance $\chi_{c0}(3915)$ are calculated in the framework of
relativistic quantum field theory. In this actual calculation, we minutely
show how to obtain the correction for energy level of resonance and to exhibit
the key features of dispersion relation in an extended Feynman diagram. The
numerical results are consistent with the experimental values.

###### pacs:

12.40.Yx, 14.40.Rt, 12.39.Ki

##  I introduction

Hadronic molecule structure has been proposed to interpret the internal
structure of exotic meson resonance for many years Swanson (2004); Törnqvist
(2004). In the previous works Swanson (2004); Törnqvist (2004); Liu and Zhu
(2009); Chen and Lü (2015); Zhao et al. (2022), molecular states were
considered as meson-meson bound states and homogeneous Bethe-Salpeter (BS)
equation was frequently used to investigate molecular states. Solving
homogeneous BS equations for meson-meson bound states, the authors of these
works obtained the masses and BS wave functions. The mass of meson-meson bound
state was regarded as mass of exotic meson resonance. However, all decay
channels of resonance should contribute to its physical mass and the
correction for energy level of molecular state due to decay channels has
seldom been considered Swanson (2004); Törnqvist (2004); Liu and Zhu (2009);
Chen and Lü (2015); Zhao et al. (2022); Chen and Lü (2018); Chen et al.
(2020). Fortunately, recent fundamental research Chen and Lü (2023) has
noticed that hadron resonance should be regarded as an unstable two-body
system, and developed BS theory for dealing with the dynamics of coupled
channels in the framework of relativistic quantum field theory. Though Ref.
Chen and Lü (2023) illuminated the physical meaning of the developed Bethe-
Salpeter theory for dealing with resonance, many details in computational
process have not been presented. In this paper, we will comprehensively and
systematically show the theoretical approach about unstable molecular state
composed of two heavy vector mesons, and this approach is applied to
investigate exotic meson resonance $\chi_{c0}(3915)$ Esposito et al. (2017),
once named $\emph{X}(3915)$, which is considered as a mixed state of two
unstable molecular states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$.

Since resonance is an unstable state which decays spontaneously into other
particles, the molecular state composed of two heavy vector mesons should not
be a stationary vector-vector bound state. To investigate this unstable two-
body system, we suppose that at some given time this unstable state has been
prepared to decay and then study the time evolution of this system as
determined by the total Hamiltonian. This prepared state can be described by
the ground-state BS wave function for vector-vector bound state at the times
$t_{1}=0$ and $t_{2}=0$. In our previous works Chen and Lü (2015); Chen et al.
(2020), the most general form of BS wave functions for the bound states
created by two vector fields with arbitrary spin and definite parity has been
given. According to the effective theory at low energy QCD, we have
investigated the light meson interaction with quarks in heavy vector mesons
and obtained the interaction kernel between two heavy vector mesons derived
from one light meson ($\sigma$, $\omega$, $\rho$, $\phi$) exchange Chen et al.
(2013); Chen and Lü (2015). Solving BS equation with this interaction kernel,
we have obtained the mass and BS wave function for bound state composed of two
vector mesons Chen and Lü (2015); Chen et al. (2016). After providing the
description for the prepared state, we can study the time evolution of BS wave
function and obtain the pole corresponding to resonance through the scattering
matrix element.

The crucial point of our resonance theory is that the scattering matrix
element between bound states is calculated in the framework of relativistic
quantum field theory. According to dispersion relation, the total matrix
element between a final state and an initial bound state should be calculated
with respect to arbitrary value of the final state energy Chen and Lü (2023).
It is necessary to note that the total energy of the final state extends over
the real interval while the initial state energy is specified. For the initial
bound state composed of two heavy vector mesons, we have given the generalized
Bethe-Salpeter (GBS) amplitude for four-quark state describing this meson-
meson structure Chen and Lü (2018); Chen et al. (2020), which should be
specified. Because the value of the final state energy is an arbitrary real
number over the real interval, we may obtain several closed channels derived
from the interaction Lagrangian and all open and closed channels should
contribute to the mass of physical resonance. For exotic resonance
$\chi_{c0}(3915)$, we consider three open decay channels $J/\psi\omega$,
$D^{+}D^{-}$, $D^{0}\bar{D}^{0}$ and one closed channel $D^{*}\bar{D}^{*}$
from the effective interaction Lagrangian at low energy QCD. Mandelstam’s
approach is applied to calculate the matrix element between bound states with
respect to arbitrary value of the final state energy, which are exhibited by
extended Feynman diagrams. Finally, we obtain the correction for energy level
of resonance $\chi_{c0}(3915)$ and the physical mass is used to calculate the
decay width of physical resonance $\chi_{c0}(3915)$.

The structure of this article is as follows. In Sec. II we give the revised
general form of GBS wave functions for meson-meson bound states as four-quark
states. The mass and GBS wave function for the mixed state of two bound states
$D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ is obtained in instantaneous
approximation. Section III gives the traditional technique to calculate the
matrix element with mass of meson-meson bound state, which is applied to
investigate the decay modes $\chi_{c0}(3915)\rightarrow J/\psi\omega$,
$\chi_{c0}(3915)\rightarrow D^{+}D^{-}$ and $\chi_{c0}(3915)\rightarrow
D^{0}\bar{D}^{0}$. Section IV gives the developed Bethe-Salpeter theory. In
Sec. V we emphatically introduce the matrix element between bound states with
respect to arbitrary value of the final state energy. Three open decay
channels $J/\psi\omega$, $D^{+}D^{-}$, $D^{0}\bar{D}^{0}$ and one closed
channel $D^{*}\bar{D}^{*}$ are considered. In Sec. VI we obtain the physical
mass and width for unstable molecular state. Our numerical results are
presented in Sec. VII and we make some concluding remarks in Sec. VIII.

##  II GBS wave function of meson-meson bound state as a four-quark state

In this paper, we investigate the light meson interaction with the light
quarks in heavy mesons. As in effective theory at low energy QCD, the
interaction Lagrangian for the coupling of light quark fields to light meson
fields is Chen and Lü (2018)

$\begin{split}&\mathscr{L}^{\text{eff}}_{I}=ig_{0}\left(\begin{array}[]{ccc}\bar{u}&\bar{d}&\bar{s}\end{array}\right)\gamma_{5}\left(\begin{array}[]{ccc}\pi^{0}+\frac{1}{\sqrt{3}}\eta&\sqrt{2}\pi^{+}&\sqrt{2}K^{+}\\\
\sqrt{2}\pi^{-}&-\pi^{0}+\frac{1}{\sqrt{3}}\eta&\sqrt{2}K^{0}\\\
\sqrt{2}K^{-}&\sqrt{2}\bar{K}^{0}&-\frac{2}{\sqrt{3}}\eta\end{array}\right)\left(\begin{array}[]{c}u\\\
d\\\ s\end{array}\right)\\\
&+ig^{\prime}_{0}\left(\begin{array}[]{ccc}\bar{u}&\bar{d}&\bar{s}\end{array}\right)\gamma_{\mu}\left(\begin{array}[]{ccc}\rho^{0}+\omega&\sqrt{2}\rho^{+}&\sqrt{2}K^{*+}\\\
\sqrt{2}\rho^{-}&-\rho^{0}+\omega&\sqrt{2}K^{*0}\\\
\sqrt{2}K^{*-}&\sqrt{2}\bar{K}^{*0}&\sqrt{2}\phi\end{array}\right)_{\mu}\left(\begin{array}[]{c}u\\\
d\\\
s\end{array}\right)+g_{\sigma}\left(\begin{array}[]{cc}\bar{u}&\bar{d}\end{array}\right)\left(\begin{array}[]{c}u\\\
d\end{array}\right)\sigma.\end{split}$ (1)

From this effective interaction Lagrangian at low energy QCD, we have to
consider that the heavy meson is a bound state composed of a quark and an
antiquark and investigate the interaction of light meson with quarks in heavy
meson. The quark current $J_{\mu}$ coupling with light vector meson and the
quark scalar density $J$ coupling with $\sigma$ meson can be obtained. In this
section, our attention is only focused on the bound state composed of two
vector mesons and some errors in previous works are revised.

###  II.1 BS wave function for bound state composed of two vector mesons

If a bound state with spin $j$ and parity $\eta_{P}$ is created by two
Heisenberg vector fields with masses $M_{1}$ and $M_{2}$, respectively, its BS
wave function is defined as

$\chi_{P(\lambda\tau)}^{j}(x_{1}^{\prime},x_{2}^{\prime})=\langle
0|TA_{\lambda}(x_{1}^{\prime})A^{\dagger}_{\tau}(x_{2}^{\prime})|P,j\rangle=\frac{1}{(2\pi)^{3/2}}\frac{1}{\sqrt{2E(P)}}e^{iP\cdot
X}\chi_{P(\lambda\tau)}^{j}(X^{\prime}),$ (2)

where $P$ is the momentum of the bound state,
$E(p)=\sqrt{\mathbf{p}^{2}+m^{2}}$,
$x_{1}^{\prime}=(\mathbf{x}_{1}^{\prime},it_{1})$,
$x_{2}^{\prime}=(\mathbf{x}_{2}^{\prime},it_{2})$,
$X=\eta_{1}x_{1}^{\prime}+\eta_{2}x_{2}^{\prime}$,
$X^{\prime}=x_{1}^{\prime}-x_{2}^{\prime}$ and
$\eta_{1,2}=M_{1,2}/(M_{1}+M_{2})$. Making the Fourier transformation, we
obtain the BS wave function in the momentum representation

$\begin{split}\chi^{j}_{P}(p_{1}^{\prime},p_{2}^{\prime})_{\lambda\tau}=\frac{1}{(2\pi)^{3/2}}\frac{1}{\sqrt{2E(P)}}(2\pi)^{4}\delta^{(4)}(P-p_{1}^{\prime}+p_{2}^{\prime})\chi^{j}_{\lambda\tau}(P,p),\end{split}$
(3)

where $p$ is the relative momentum of two vector fields and we have
$P=p_{1}^{\prime}-p_{2}^{\prime}$,
$p=\eta_{2}p_{1}^{\prime}+\eta_{1}p_{2}^{\prime}$, $p_{1}^{\prime}$ and
$p_{2}^{\prime}$ are the momenta carried by two vector fields, respectively.
The polarization tensor of the bound state
$\eta_{\mu_{1}\mu_{2}\cdots\mu_{j}}$ can be separated,

$\chi_{\lambda\tau}^{j}(P,p)=\eta_{\mu_{1}\mu_{2}\cdots\mu_{j}}\chi_{\mu_{1}\mu_{2}\cdots\mu_{j}\lambda\tau}(P,p),$
(4)

where the subscripts $\lambda$ and $\tau$ are derived from these two vector
fields. The polarization tensor $\eta_{\mu_{1}\mu_{2}\cdots\mu_{j}}$ describes
the spin of the bound state, which is totally symmetric, transverse and
traceless:

$\eta_{\mu_{1}\mu_{2}\cdots}=\eta_{\mu_{2}\mu_{1}\cdots},\quad
P_{\mu_{1}}\eta_{\mu_{1}\mu_{2}\cdots}=0,\quad\eta_{\mu_{1}\mu_{2}\cdots}^{\mu_{1}}=0.$
(5)

From Lorentz covariance, we have

$\begin{split}\chi_{\mu_{1}\cdots\mu_{j}\lambda\tau}=&p_{\mu_{1}}\cdots
p_{\mu_{j}}[g_{\lambda\tau}f_{1}+(P_{\lambda}p_{\tau}+P_{\tau}p_{\lambda})f_{2}+(P_{\lambda}p_{\tau}-P_{\tau}p_{\lambda})f_{3}+P_{\lambda}P_{\tau}f_{4}+p_{\lambda}p_{\tau}f_{5}]\\\
&+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\lambda}p_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\tau}p_{\lambda})f_{6}\\\ &+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\lambda}p_{\tau}-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\tau}p_{\lambda})f_{7}\\\ &+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\lambda}P_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\tau}P_{\lambda})f_{8}\\\ &+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\lambda}P_{\tau}-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\tau}P_{\lambda})f_{9}\\\ &+p_{\mu_{1}}\cdots
p_{\mu_{j}}\epsilon_{\lambda\tau\xi\zeta}p_{\xi}P_{\zeta}f_{10}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}p_{\xi}f_{11}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}P_{\xi}f_{12}\\\
&+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}p_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}p_{\lambda})f_{13}\\\
&+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}p_{\tau}-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}p_{\lambda})f_{14}\\\
&+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}P_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}P_{\lambda})f_{15}\\\
&+(p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}P_{\tau}-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}P_{\lambda})f_{16}\\\
&+p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\lambda}g_{\mu_{2}\\}\tau}f_{17}+p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\lambda\xi\zeta}p_{\xi}P_{\zeta}\epsilon_{\mu_{2}\\}\tau\xi^{\prime}\zeta^{\prime}}p_{\xi^{\prime}}P_{\zeta^{\prime}}f_{18}\\\
&+(p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\lambda}\epsilon_{\mu_{2}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}+p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\tau}\epsilon_{\mu_{2}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta})f_{19}\\\
&+(p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\lambda}\epsilon_{\mu_{2}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}-p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\tau}\epsilon_{\mu_{2}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta})f_{20},\end{split}$
(6)

where $\\{\mu_{1},\cdots,\mu_{j}\\}$ represents symmetrization of the indices
$\mu_{1},\cdots,\mu_{j}$. In fact, the relative momenta
$p_{\mu_{1}},\cdots,p_{\mu_{j}},p_{\lambda},p_{\tau}$ represent the orbital
angular momenta. There should be 20 scalar functions $f_{i}(P\cdot
p,p^{2})(i=1,\cdots,20)$ in Eq. (6). In Ref. Chen et al. (2020), three tensor
structures are omitted. In this paper, these missing terms are added as the
last three terms in Eq. (6). Using the transversality condition Chen et al.
(2013); Chen and Lü (2015)

$p^{\prime}_{1\lambda}\chi^{j}_{\lambda\tau}(P,p)=p^{\prime}_{2\tau}\chi^{j}_{\lambda\tau}(P,p)=0$
(7)

and considering the properties of BS wave function under space reflection, we
obtain the revised general form of BS wave functions for the bound states
created by two massive vector fields with arbitrary spin and definite parity
(see details in Chen and Lü (2015)), for $\eta_{P}=(-1)^{j}$,

$\begin{split}\chi_{\lambda\tau}^{j}(P,p)=&\frac{1}{\mathcal{N}^{j}}\eta_{\mu_{1}\cdots\mu_{j}}[p_{\mu_{1}}\cdots
p_{\mu_{j}}(\mathcal{T}^{1}_{\lambda\tau}\Phi_{1}+\mathcal{T}^{2}_{\lambda\tau}\Phi_{2})+\mathcal{T}^{3}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi_{3}+\mathcal{T}^{4}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi_{4}\\\
&+\mathcal{T}^{5}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi_{5}+\mathcal{T}^{6}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi_{6}],\end{split}$
(8)

for $\eta_{P}=(-1)^{j+1}$,

$\begin{split}\chi_{\lambda\tau}^{j}(P,p)=&\frac{1}{\mathcal{N}^{j}}\eta_{\mu_{1}\cdots\mu_{j}}(p_{\mu_{1}}\cdots
p_{\mu_{j}}\epsilon_{\lambda\tau\xi\zeta}p_{\xi}P_{\zeta}\Phi^{\prime}_{1}+\mathcal{T}^{7}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi^{\prime}_{2}+\mathcal{T}^{8}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi^{\prime}_{3}+\mathcal{T}^{9}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi^{\prime}_{4}\\\
&+\mathcal{T}^{10}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi^{\prime}_{5}+\mathcal{T}^{11}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi^{\prime}_{6}+\mathcal{T}^{12}_{\mu_{1}\cdots\mu_{j}\lambda\tau}\Phi^{\prime}_{7}),\end{split}$
(9)

where $\mathcal{N}^{j}$ is normalization, the independent tensor structures
$\mathcal{T}^{i}_{\lambda\tau}$ are given in Appendix A, $\Phi_{i}(P\cdot
p,p^{2})$ and $\Phi^{\prime}_{i}(P\cdot p,p^{2})$ are independent scalar
functions. Scalar functions $f_{i}$ in Eq. (6) are the linear combinations of
$\Phi_{i}$ and $\Phi^{\prime}_{i}$.

###  II.2 Kernel between two heavy vector mesons

In this paper, we assume that the isoscalar $\chi_{c0}(3915)$ is a mixed state
of two unstable molecular states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ with
spin-parity quantum numbers $0^{+}$. In this section, we only investigate the
mixed state of two stable bound states $D^{*0}\bar{D}^{*0}$ and
$D^{*+}D^{*-}$, and BS wave function for this system is a linear combination
of two components as

$\begin{split}\chi^{D^{*}\bar{D}^{*},j}_{\lambda\tau}(P,p)=&\frac{1}{\sqrt{2}}\chi^{D^{*0}\bar{D}^{*0},j}_{\lambda\tau}(P,p)+\frac{1}{\sqrt{2}}\chi^{D^{*+}D^{*-},j}_{\lambda\tau}(P,p),\end{split}$
(10)

where

$\begin{split}\chi^{D^{*0}\bar{D}^{*0},j}_{\lambda\tau}(P,p)=&\chi^{j}_{\lambda\tau}(P,p)\bigg{(}-\bigg{|}\frac{1}{2},-\frac{1}{2}\bigg{\rangle}\bigg{)}^{D^{*0}}\otimes\bigg{|}\frac{1}{2},\frac{1}{2}\bigg{\rangle}^{\bar{D}^{*0}},\\\
\chi^{D^{*+}D^{*-},j}_{\lambda\tau}(P,p)=&\chi^{j}_{\lambda\tau}(P,p)\bigg{(}-\bigg{|}\frac{1}{2},\frac{1}{2}\bigg{\rangle}\bigg{)}^{D^{*+}}\otimes\bigg{|}\frac{1}{2},-\frac{1}{2}\bigg{\rangle}^{D^{*-}},\end{split}$
(11)

and $P$ becomes the total momentum for the mixed state of two meson-meson
bound states, $\chi^{j}_{\lambda\tau}(P,p)$ is the component wave function in
the momentum representation;
$(-|\frac{1}{2},-\frac{1}{2}\rangle)\otimes|\frac{1}{2},\frac{1}{2}\rangle$
and
$(-|\frac{1}{2},\frac{1}{2}\rangle)\otimes|\frac{1}{2},-\frac{1}{2}\rangle$
are the isospin wave functions of pure bound states $D^{*0}\bar{D}^{*0}$ and
$D^{*+}D^{*-}$, respectively; $\chi^{D^{*0}\bar{D}^{*0},j}_{\lambda\tau}$ and
$\chi^{D^{*+}D^{*-},j}_{\lambda\tau}$ represent BS wave functions for the
bound states of two vector mesons, which are the eigenstates of Hamiltonian
without considering the coupled-channel terms. These eigenstates have the same
quantum numbers. The error in Ref. Chen et al. (2020) has been revised. As
usual the momentum for the mixed state of two bound states is set as
$P=(0,0,0,iM_{0})$ in the rest frame.

Let $D^{*}_{l}$ denote one of $D^{*0}$ and $D^{*+}$, and $l=u,d$ represents
the $u$ or $d$ antiquark in heavy vector meson $D^{*0}$ or $D^{*+}$,
respectively; $\bar{D}_{l}^{*}$ denotes the antiparticle of $D_{l}^{*}$. From
Eq. (8), we can obtain BS wave function describing pure bound state
$D^{*}_{l}\bar{D}^{*}_{l}$

$\begin{split}\chi_{\lambda\tau}^{0^{+}}(P^{D\bar{D}},p)=\frac{1}{\mathcal{N}_{D\bar{D}}^{0^{+}}}[\mathcal{T}^{1}_{\lambda\tau}\mathcal{F}_{1}(P^{D\bar{D}}\cdot
p,p^{2})+\mathcal{T}^{2}_{\lambda\tau}\mathcal{F}_{2}(P^{D\bar{D}}\cdot
p,p^{2})].\end{split}$ (12)

$P^{D\bar{D}}$ represents the momentum of pure bound state in the rest frame,
whose fourth component is different from the one of $P$. This BS wave function
should satisfy the equation

$\chi^{0^{+}}_{\lambda\tau}(P^{D\bar{D}},p)=-\int\frac{d^{4}p^{\prime}}{(2\pi)^{4}}\Delta_{F\lambda\theta}(p_{1}^{\prime})\mathcal{V}_{\theta\theta^{\prime},\kappa^{\prime}\kappa}(p,p^{\prime};P^{D\bar{D}})\chi^{0^{+}}_{\theta^{\prime}\kappa^{\prime}}(P^{D\bar{D}},p^{\prime})\Delta_{F\kappa\tau}(p_{2}^{\prime}),$
(13)

where $\mathcal{V}_{\theta\theta^{\prime},\kappa^{\prime}\kappa}$ is the
interaction kernel, $P^{D\bar{D}}=(0,0,0,iM_{D\bar{D}})$,
$p_{1}^{\prime}=p+P^{D\bar{D}}/2$, $p_{2}^{\prime}=p-P^{D\bar{D}}/2$,
$\Delta_{F\lambda\theta}(p_{1}^{\prime})$ and
$\Delta_{F\kappa\tau}(p_{2}^{\prime})$ are the propagators for the spin 1
fields,
$\Delta_{F\lambda\theta}(p_{1}^{\prime})=(\delta_{\lambda\theta}+\frac{p^{\prime}_{1\lambda}p^{\prime}_{1\theta}}{M_{1}^{2}})\frac{-i}{p_{1}^{\prime
2}+M_{1}^{2}-i\varepsilon}$,
$\Delta_{F\kappa\tau}(p_{2}^{\prime})=(\delta_{\kappa\tau}+\frac{p^{\prime}_{2\kappa}p^{\prime}_{2\tau}}{M_{2}^{2}})\frac{-i}{p_{2}^{\prime
2}+M_{2}^{2}-i\varepsilon}$, $M_{1}=M_{D^{*}_{l}}$ and
$M_{2}=M_{\bar{D}_{l}^{*}}$. We emphasize that the kernel $\mathcal{V}$ is
defined in two-body channel so $\mathcal{V}$ is not complete interaction. The
kernel in homogeneous BS equation (13) plays a central role for making two-
body system to be a stable bound state, and the solution of homogeneous BS
equation (13) should only describe bound state.

To construct the interaction kernel between $D^{*}_{l}$ and $\bar{D}^{*}_{l}$,
we consider that the effective interaction is derived from one light meson
($\sigma$, $\rho^{0}$, $V_{1}$ and $V_{8}$) exchange Chen et al. (2013); Chen
and Lü (2015); Chen et al. (2016), shown as Fig. 1. The flavor-SU(3) singlet
$V_{1}$ and octet $V_{8}$ states of vector mesons mix to form the physical
$\omega$ and $\phi$ mesons as

$\phi=-V_{8}\cos\theta+V_{1}\sin\theta,\quad\omega=V_{8}\sin\theta+V_{1}\cos\theta,$
(14)

where the mixing angle $\theta=38.58^{\circ}$ was obtained by KLOE Ambrosino
et al. (2009). Then the exchanged mesons should be the octet $V_{8}$ and
singlet $V_{1}$ states, and the relation of the octet-quark coupling constant
$g_{8}$ and the singlet-quark coupling constant $g_{1}$ has the form

$g_{\phi}=-g_{8}\cos\theta+g_{1}\sin\theta,\quad
g_{\omega}=g_{8}\sin\theta+g_{1}\cos\theta,$ (15)

where the meson-quark coupling constants $g_{\omega}^{2}=2.42/2$ and
$g_{\phi}^{2}=13.0$ were determined by QCD sum rules approach Reinders et al.
(1985).

Figure 1:  The light meson exchange between two heavy vector mesons.

In Fig. 1, $VM$ represents $D_{l}^{*}$ and $\overline{VM^{\prime}}$ represents
$\bar{D}_{l}^{*}$. From the Lorentz-structure, the matrix elements of quark
scalar density $J$ and quark current $J_{\alpha}$ can be expressed as

$\begin{split}\langle&VM^{\varrho}(p_{1}^{\prime})|J|VM^{\vartheta}(q_{1}^{\prime})\rangle=\frac{1}{2\sqrt{E_{D^{*}_{l}}(p_{1}^{\prime})E_{D^{*}_{l}}(q_{1}^{\prime})}}\\\
&\times\bigg{\\{}[\varepsilon^{\varrho}(p_{1}^{\prime})\cdot\varepsilon^{\vartheta}(q_{1}^{\prime})]h^{(\text{s})}_{1}(w^{2})-h^{(\text{s})}_{2}(w^{2})\frac{1}{M_{1}^{2}}[\varepsilon^{\varrho}(p_{1}^{\prime})\cdot
q_{1}^{\prime}][\varepsilon^{\vartheta}(q_{1}^{\prime})\cdot
p_{1}^{\prime}]\bigg{\\}},\end{split}$ (16a)
$\begin{split}\langle&\overline{VM^{\prime}}^{\varrho^{\prime}}(-p_{2}^{\prime})|J|\overline{VM^{\prime}}^{\vartheta^{\prime}}(-q_{2}^{\prime})\rangle=\frac{1}{2\sqrt{E_{\bar{D}^{*}_{l}}(-p_{2}^{\prime})E_{\bar{D}^{*}_{l}}(-q_{2}^{\prime})}}\\\
&\times\bigg{\\{}[\varepsilon^{\varrho^{\prime}}(-p_{2}^{\prime})\cdot\varepsilon^{\vartheta^{\prime}}(-q_{2}^{\prime})]\bar{h}^{(\text{s})}_{1}(w^{2})-\bar{h}^{(\text{s})}_{2}(w^{2})\frac{1}{M_{2}^{2}}[\varepsilon^{\varrho^{\prime}}(-p_{2}^{\prime})\cdot(-q_{2}^{\prime})][\varepsilon^{\vartheta^{\prime}}(-q_{2}^{\prime})\cdot(-p_{2}^{\prime})]\bigg{\\}},\end{split}$
(16b)
$\begin{split}\langle&VM^{\varrho}(p_{1}^{\prime})|J_{\alpha}|VM^{\vartheta}(q_{1}^{\prime})\rangle=\frac{1}{2\sqrt{E_{D^{*}_{l}}(p_{1}^{\prime})E_{D^{*}_{l}}(q_{1}^{\prime})}}\\\
&\times\bigg{\\{}[\varepsilon^{\varrho}(p_{1}^{\prime})\cdot\varepsilon^{\vartheta}(q_{1}^{\prime})]h_{1}^{(\text{v})}(w^{2})(p_{1}^{\prime}+q_{1}^{\prime})_{\alpha}-h^{(\text{v})}_{2}(w^{2})\\{[\varepsilon^{\varrho}(p_{1}^{\prime})\cdot
q_{1}^{\prime}]\varepsilon^{\vartheta}_{\alpha}(q_{1}^{\prime})\\\
&+[\varepsilon^{\vartheta}(q_{1}^{\prime})\cdot
p_{1}^{\prime}]\varepsilon^{\varrho}_{\alpha}(p_{1}^{\prime})\\}-h_{3}^{(\text{v})}(w^{2})\frac{1}{M_{1}^{2}}[\varepsilon^{\varrho}(p_{1}^{\prime})\cdot
q_{1}^{\prime}][\varepsilon^{\vartheta}(q_{1}^{\prime})\cdot
p_{1}^{\prime}](p_{1}^{\prime}+q_{1}^{\prime})_{\alpha}\bigg{\\}},\end{split}$
(16c)
$\begin{split}\langle&\overline{VM^{\prime}}^{\varrho^{\prime}}(-p_{2}^{\prime})|J_{\beta}|\overline{VM^{\prime}}^{\vartheta^{\prime}}(-q_{2}^{\prime})\rangle=\frac{1}{2\sqrt{E_{\bar{D}^{*}_{l}}(-p_{2}^{\prime})E_{\bar{D}^{*}_{l}}(-q_{2}^{\prime})}}\\\
&\times\bigg{\\{}[\varepsilon^{\varrho^{\prime}}(-p_{2}^{\prime})\cdot\varepsilon^{\vartheta^{\prime}}(-q_{2}^{\prime})]\bar{h}_{1}^{(\text{v})}(w^{2})(-p_{2}^{\prime}-q_{2}^{\prime})_{\beta}\\\
&-\bar{h}^{(\text{v})}_{2}(w^{2})\\{[\varepsilon^{\varrho^{\prime}}(-p_{2}^{\prime})\cdot(-q_{2}^{\prime})]\varepsilon^{\vartheta^{\prime}}_{\beta}(-q_{2}^{\prime})+[\varepsilon^{\vartheta^{\prime}}(-q_{2}^{\prime})\cdot(-p_{2}^{\prime})]\varepsilon^{\varrho^{\prime}}_{\beta}(-p_{2}^{\prime})\\}\\\
&-\bar{h}_{3}^{(\text{v})}(w^{2})\frac{1}{M_{2}^{2}}[\varepsilon^{\varrho^{\prime}}(-p_{2}^{\prime})\cdot(-q_{2}^{\prime})][\varepsilon^{\vartheta^{\prime}}(-q_{2}^{\prime})\cdot(-p_{2}^{\prime})](-p_{2}^{\prime}-q_{2}^{\prime})_{\beta}\bigg{\\}},\end{split}$
(16d)

where $p^{\prime}_{1}=(\mathbf{p},ip_{10}^{\prime})$,
$p^{\prime}_{2}=(\mathbf{p},ip_{20}^{\prime})$,
$q^{\prime}_{1}=(\mathbf{p}^{\prime},iq_{10}^{\prime})$,
$q^{\prime}_{2}=(\mathbf{p}^{\prime},iq_{20}^{\prime})$,
$w=q_{1}^{\prime}-p_{1}^{\prime}=q_{2}^{\prime}-p_{2}^{\prime}$ is the
momentum of the light meson and $\mathbf{w}=\mathbf{p}^{\prime}-\mathbf{p}$;
$h(w^{2})$ and $\bar{h}(w^{2})$ are scalar functions, the four-vector
$\varepsilon(p)$ is the polarization vector of heavy vector meson with
momentum $p$ and $E_{D_{l}^{*}}(p)=\sqrt{\mathbf{p}^{2}+M_{D_{l}^{*}}^{2}}$.
Taking away the external lines including normalizations and polarization
vectors $\varepsilon^{\varrho}_{\theta}(p_{1}^{\prime})$,
$\varepsilon^{\vartheta}_{\theta^{\prime}}(q_{1}^{\prime})$,
$\varepsilon^{\varrho^{\prime}}_{\kappa}(-p_{2}^{\prime})$,
$\varepsilon^{\vartheta^{\prime}}_{\kappa^{\prime}}(-q_{2}^{\prime})$, we
obtain the interaction kernel from one light meson ($\sigma$, $\rho^{0}$,
$V_{1}$ and $V_{8}$) exchange Chen et al. (2013); Chen and Lü (2015)

$\begin{split}\mathcal{V}&{}_{\theta\theta^{\prime},\kappa^{\prime}\kappa}(p,p^{\prime};P^{D\bar{D}})\\\
=&h_{1}^{(\text{s})}(w^{2})\frac{-ig_{\sigma}^{2}}{w^{2}+M_{\sigma}^{2}}\bar{h}_{1}^{(\text{s})}(w^{2})\delta_{\theta\theta^{\prime}}\delta_{\kappa^{\prime}\kappa}+\left(\frac{-ig_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}+\frac{-ig_{1}^{2}}{w^{2}+M_{\omega}^{2}}+\frac{-ig_{8}^{2}}{w^{2}+M_{\phi}^{2}}\right)\\{h_{1}^{(\text{v})}(w^{2})\bar{h}_{1}^{(\text{v})}(w^{2})\\\
&\times(p_{1}^{\prime}+q_{1}^{\prime})\cdot(-p_{2}^{\prime}-q_{2}^{\prime})\delta_{\theta\theta^{\prime}}\delta_{\kappa^{\prime}\kappa}-h_{1}^{(\text{v})}(w^{2})\bar{h}_{2}^{(\text{v})}(w^{2})\delta_{\theta\theta^{\prime}}[-(p_{1}^{\prime}+q_{1}^{\prime})_{\kappa^{\prime}}q_{2\kappa}^{\prime}-p_{2\kappa^{\prime}}^{\prime}(p_{1}^{\prime}+q_{1}^{\prime})_{\kappa}]\\\
&-h_{2}^{(\text{v})}(w^{2})\bar{h}_{1}^{(\text{v})}(w^{2})[q_{1\theta}^{\prime}(-p_{2}^{\prime}-q_{2}^{\prime})_{\theta^{\prime}}+(-p_{2}^{\prime}-q_{2}^{\prime})_{\theta}p_{1\theta^{\prime}}^{\prime}]\delta_{\kappa^{\prime}\kappa}+h_{2}^{(\text{v})}(w^{2})\bar{h}_{2}^{(\text{v})}(w^{2})[-q_{1\theta}^{\prime}\delta_{\theta^{\prime}\kappa^{\prime}}q_{2\kappa}^{\prime}\\\
&+q_{1\theta}^{\prime}\delta_{\theta^{\prime}\kappa}(-p_{2\kappa^{\prime}}^{\prime})-\delta_{\theta\kappa^{\prime}}p_{1\theta^{\prime}}^{\prime}q_{2\kappa}^{\prime}+\delta_{\theta\kappa}p_{1\theta^{\prime}}^{\prime}(-p_{2\kappa^{\prime}}^{\prime})]\\},\end{split}$
(17)

where $g_{\sigma}=\frac{B(M_{\sigma})}{f_{\sigma}}=\frac{299}{60}$ Lü et al.
(1996, 1997), $g_{\rho}^{2}=2.42$ Reinders et al. (1985), and these terms
containing $M_{1,2}$ are neglected because the masses of heavy mesons are
large. Using the method above, we can obtain the interaction kernels from
one-$\rho^{\pm}$ exchange Chen et al. (2016).

###  II.3 Instantaneous approximation

####  II.3.1 The extended Bethe-Salpeter equation

Substituting BS wave function given by Eq. (12) and the kernel (17) into BS
equation (13), we find that the integral of one term on the right-hand side of
(12) has contribution to the one of itself and the other term. Ignoring the
cross terms, one can obtain two individual equations:

$\mathcal{F}^{1}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})=-\int\frac{d^{4}p^{\prime}}{(2\pi)^{4}}\Delta_{F\lambda\theta}(p_{1}^{\prime})\mathcal{V}_{\theta\theta^{\prime},\kappa^{\prime}\kappa}(p,p^{\prime};P^{D\bar{D}})\mathcal{F}^{1}_{\theta^{\prime}\kappa^{\prime}}(P^{D\bar{D}}\cdot
p^{\prime},p^{\prime 2})\Delta_{F\kappa\tau}(p_{2}^{\prime}),$ (18)
$\mathcal{F}^{2}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})=-\int\frac{d^{4}p^{\prime}}{(2\pi)^{4}}\Delta_{F\lambda\theta}(p_{1}^{\prime})\mathcal{V}_{\theta\theta^{\prime},\kappa^{\prime}\kappa}(p,p^{\prime};P^{D\bar{D}})\mathcal{F}^{2}_{\theta^{\prime}\kappa^{\prime}}(P^{D\bar{D}}\cdot
p^{\prime},p^{\prime 2})\Delta_{F\kappa\tau}(p_{2}^{\prime}),$ (19)

where $\mathcal{F}^{1}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})=\mathcal{T}^{1}_{\lambda\tau}\mathcal{F}_{1}(P^{D\bar{D}}\cdot
p,p^{2})$ and $\mathcal{F}^{2}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})=\mathcal{T}^{2}_{\lambda\tau}\mathcal{F}_{2}(P^{D\bar{D}}\cdot
p,p^{2})$. Comparing the tensor structures in both sides of Eqs. (18) and
(19), respectively, we obtain

$\begin{split}\mathcal{F}_{1}(P^{D\bar{D}}\cdot
p,p^{2})=\frac{1}{p_{1}^{\prime
2}+M_{1}^{2}-i\varepsilon}\frac{1}{p_{2}^{\prime
2}+M_{2}^{2}-i\varepsilon}\int\frac{d^{4}p^{\prime}}{(2\pi)^{4}}V^{0^{+}}_{1}(p,p^{\prime};P^{D\bar{D}})\mathcal{F}_{1}(P^{D\bar{D}}\cdot
p^{\prime},p^{\prime 2}),\end{split}$ (20) $\begin{split}p_{2}^{\prime
2}\mathcal{F}_{2}(P^{D\bar{D}}\cdot p,p^{2})=\frac{1}{p_{1}^{\prime
2}+M_{1}^{2}-i\varepsilon}\frac{1}{p_{2}^{\prime
2}+M_{2}^{2}-i\varepsilon}\int\frac{d^{4}p^{\prime}}{(2\pi)^{4}}V^{0^{+}}_{2}(p,p^{\prime};P^{D\bar{D}})q_{2}^{\prime
2}\mathcal{F}_{2}(P^{D\bar{D}}\cdot p^{\prime},p^{\prime 2}),\end{split}$ (21)

where $V^{0^{+}}_{1}(p,p^{\prime};P^{D\bar{D}})$ and
$V^{0^{+}}_{2}(p,p^{\prime};P^{D\bar{D}})$ are derived from the interaction
kernel between $D_{l}^{*}$ and $\bar{D}_{l}^{*}$. In instantaneous
approximation, we set the momentum of exchanged meson as $w=(\mathbf{w},0)$.
Then Eqs. (20) and (21) become two relativistic Schr$\ddot{o}$dinger-like
equations (see details in Refs. Chen and Lü (2015); Chen et al. (2016))

$\left(\frac{b_{1}^{2}(M_{D\bar{D}})}{2\mu_{R}}-\frac{\mathbf{p}^{2}}{2\mu_{R}}\right)\Psi_{1}^{0^{+}}(\mathbf{p})=\int\frac{d^{3}w}{(2\pi)^{3}}V_{1}^{0^{+}}(\mathbf{p},\mathbf{w})\Psi_{1}^{0^{+}}(\mathbf{p},\mathbf{w}),$
(22)
$\left(\frac{b_{2}^{2}(M_{D\bar{D}})}{2\mu_{R}}-\frac{\mathbf{p}^{2}}{2\mu_{R}}\right)\Psi_{2}^{0^{+}}(\mathbf{p})=\int\frac{d^{3}w}{(2\pi)^{3}}V_{2}^{0^{+}}(\mathbf{p},\mathbf{w})\Psi_{2}^{0^{+}}(\mathbf{p},\mathbf{w}),$
(23)

where $\Psi_{1}^{0^{+}}(\mathbf{p})=\int
dp_{0}\mathcal{F}_{1}(P^{D\bar{D}}\cdot p,p^{2})$,
$\Psi_{2}^{0^{+}}(\mathbf{p})=\int dp_{0}p_{2}^{\prime
2}\mathcal{F}_{2}(P^{D\bar{D}}\cdot p,p^{2})$,
$\mu_{R}=E_{1}E_{2}/(E_{1}+E_{2})=[M_{D\bar{D}}^{4}-(M_{1}^{2}-M_{2}^{2})^{2}]/(4M_{D\bar{D}}^{3})$,
$b^{2}(M_{D\bar{D}})=[M_{D\bar{D}}^{2}-(M_{1}+M_{2})^{2}][M_{D\bar{D}}^{2}-(M_{1}-M_{2})^{2}]/(4M_{D\bar{D}}^{2})$,
$E_{1}=(M_{D\bar{D}}^{2}-M_{2}^{2}+M_{1}^{2})/(2M_{D\bar{D}})$ and
$E_{2}=(M_{D\bar{D}}^{2}-M_{1}^{2}+M_{2}^{2})/(2M_{D\bar{D}})$. The potentials
between $D_{l}^{*}$ and $\bar{D}_{l}^{*}$ up to the second order of the
$p/M_{D^{*}_{l}}$ expansion are

$\begin{split}V_{1}^{0^{+}}(\mathbf{p},\mathbf{w})=&\frac{h_{1}^{(\text{s})}(w^{2})}{2E_{1}}\frac{g_{\sigma}^{2}}{w^{2}+M_{\sigma}^{2}}\frac{\bar{h}_{1}^{(\text{s})}(w^{2})}{2E_{2}}+h_{1}^{(\text{v})}(w^{2})\bar{h}_{1}^{(\text{v})}(w^{2})\\\
&\times\bigg{(}\frac{g_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}+\frac{g_{1}^{2}}{w^{2}+M_{\omega}^{2}}+\frac{g_{8}^{2}}{w^{2}+M_{\phi}^{2}}\bigg{)}\left(-1-\frac{4\mathbf{p}^{2}+5\mathbf{w}^{2}}{4E_{1}E_{2}}\right),\end{split}$
(24)
$\begin{split}V_{2}^{0^{+}}(\mathbf{p},\mathbf{w})=&\frac{h_{1}^{(\text{s})}(w^{2})}{2E_{1}}\frac{g_{\sigma}^{2}}{w^{2}+M_{\sigma}^{2}}\frac{\bar{h}_{1}^{(\text{s})}(w^{2})}{2E_{2}}\left(1-\frac{\mathbf{w}^{2}}{M_{1}^{2}}\right)+h_{1}^{(\text{v})}(w^{2})\bar{h}_{1}^{(\text{v})}(w^{2})\bigg{(}\frac{g_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}\\\
&+\frac{g_{1}^{2}}{w^{2}+M_{\omega}^{2}}+\frac{g_{8}^{2}}{w^{2}+M_{\phi}^{2}}\bigg{)}\left(-1-\frac{2\mathbf{p}^{2}+2\mathbf{w}^{2}}{4M_{1}^{2}}-\frac{2\mathbf{p}^{2}+2\mathbf{w}^{2}}{4E_{1}E_{2}}\right).\end{split}$
(25)

Solving these two equations (22) and (23), respectively, one can obtain the
eigenvalues $b_{1}^{2}(M_{D\bar{D}})$ and $b_{2}^{2}(M_{D\bar{D}})$ and the
corresponding eigenfunctions $\Psi_{1}^{0^{+}}(\mathbf{p})$ and
$\Psi_{2}^{0^{+}}(\mathbf{p})$. From $\Psi_{1}^{0^{+}}$ and
$\Psi_{2}^{0^{+}}$, it is easy to obtain $\mathcal{F}_{1}$ and
$\mathcal{F}_{2}$, respectively.

Because the cross terms are small, we can take the ground state BS wave
function to be a linear combination of two eigenstates
$\mathcal{F}^{10}_{\lambda\tau}(P^{D\bar{D}}\cdot p,p^{2})$ and
$\mathcal{F}^{20}_{\lambda\tau}(P^{D\bar{D}}\cdot p,p^{2})$ corresponding to
lowest energy in Eqs. (18) and (19). Then in the basis provided by
$\mathcal{F}^{10}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})=\mathcal{T}^{1}_{\lambda\tau}\mathcal{F}_{10}(P^{D\bar{D}}\cdot
p,p^{2})$ and $\mathcal{F}^{20}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})=\mathcal{T}^{2}_{\lambda\tau}\mathcal{F}_{20}(P^{D\bar{D}}\cdot
p,p^{2})$, BS wave function $\chi^{0^{+}}_{\lambda\tau}$ is considered as

$\chi^{0^{+}}_{\lambda\tau}(P^{D\bar{D}},p)=\frac{1}{\mathcal{N}_{D\bar{D}}^{0^{+}}}[\mathcal{C}_{1}\mathcal{F}^{10}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})+\mathcal{C}_{2}\mathcal{F}^{20}_{\lambda\tau}(P^{D\bar{D}}\cdot
p,p^{2})].$ (26)

Substituting (26) into BS equation (13) and comparing the tensor structures in
both sides, we obtain an eigenvalue equation in instantaneous approximation
Chen and Lü (2015)

$\left(\begin{array}[]{cc}\frac{b_{10}^{2}(M_{D\bar{D}})}{2\mu_{R}}-\lambda&H_{12}\\\
H_{21}&\frac{b_{20}^{2}(M_{D\bar{D}})}{2\mu_{R}}-\lambda\end{array}\right)\left(\begin{array}[]{c}\mathcal{C}^{\prime}_{1}\\\
\mathcal{C}^{\prime}_{2}\end{array}\right)=0,$ (27)

where we have the matrix elements

$\begin{split}H_{12}=H_{21}=&\int\frac{d^{3}p}{(2\pi)^{3}}\Psi_{10}^{0^{+}}(\mathbf{p})^{*}\int\frac{d^{3}w}{(2\pi)^{3}}h_{1}^{(\text{v})}(w^{2})\\\
&\times\left(\frac{g_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}+\frac{g_{1}^{2}}{w^{2}+M_{\omega}^{2}}+\frac{g_{8}^{2}}{w^{2}+M_{\phi}^{2}}\right)\bar{h}_{1}^{(\text{v})}(w^{2})\frac{\mathbf{w}^{2}}{E_{1}E_{2}}\Psi_{20}^{0^{+}}(\mathbf{p},\mathbf{w}),\end{split}$
(28)

and $b_{10}^{2}(M_{D\bar{D}})/(2\mu_{R})$ and
$b_{20}^{2}(M_{D\bar{D}})/(2\mu_{R})$ are the eigenvalues corresponding to
lowest energy in Eqs. (22) and (23), respectively; $\Psi_{10}^{0^{+}}$ and
$\Psi_{20}^{0^{+}}$ are the corresponding eigenfunctions. From this equation,
we can obtain the eigenvalues and eigenfunctions which contain the
contribution from the cross terms. Some errors in our previous works have been
revised.

####  II.3.2 Form factors of heavy meson

To calculate these heavy vector meson form factors $h(w^{2})$ describing the
heavy meson structure, we have to know the wave function of heavy vector meson
$D_{l}^{*}$ in instantaneous approximation. For heavy vector mesons, the
authors of Refs. Burden et al. (1997); Maris et al. (1998); Ivanov et al.
(1999, 2007) have obtained their BS amplitudes in Euclidean space:

$\begin{split}\Gamma_{\lambda}^{V}(K,k)=\frac{1}{\mathcal{N}^{V}}\bigg{(}\gamma_{\lambda}+K_{\lambda}\frac{\gamma\cdot
K}{M_{V}^{2}}\bigg{)}\varphi_{V}(k^{2}),\end{split}$ (29)

where $K$ is the momentum of heavy meson, $k$ denotes the relative momentum
between quark and antiquark in heavy meson, $M_{V}$ is heavy vector meson
mass, $\Gamma_{\lambda}^{V}(K,k)$ is transverse
($K_{\lambda}\Gamma_{\lambda}^{V}(K,k)=0$), $\mathcal{N}^{V}$ is
normalization, and $\varphi_{V}(k^{2})$ is scalar function fixed by providing
fits to observables. The charmed meson $D_{l}^{*}$ is composed of $c$-quark
and $l$-antiquark. As in heavy-quark effective theory (HQET) Neubert (1994),
we consider that the heaviest quark carries all the heavy-meson momentum and
obtain BS wave function of $D_{l}^{*}$ meson

$\begin{split}\chi_{\lambda}(K,k)=\frac{-1}{\gamma\cdot(k+K)-im_{c}}\frac{1}{\mathcal{N}^{D_{l}^{*}}}\bigg{(}\gamma_{\lambda}+K_{\lambda}\frac{\gamma\cdot
K}{M_{D_{l}^{*}}^{2}}\bigg{)}\varphi_{D_{l}^{*}}(k^{2})\frac{-1}{\gamma\cdot
k-im_{l}},\end{split}$ (30)

where $K$ is set as the momentum of heavy meson in the rest frame, $k$ becomes
the relative momentum between $c$-quark and $l$-antiquark, $m_{c,l}$ are the
constituent quark masses,
$\varphi_{D_{l}^{*}}(k^{2})=\varphi_{\bar{D}_{l}^{*}}(k^{2})=\exp(-k^{2}/\omega_{D^{*}}^{2})$
and $\omega_{D^{*}}$=1.50 GeV Ivanov et al. (2007). The components of this BS
wave function are $4\times 4$ matrices, which can be written as Chen et al.
(2009)

$\begin{split}\chi_{\lambda}(K,k)=\Psi_{\lambda}^{S}+\Psi^{V}_{\lambda,\mu}\gamma_{\mu}+\Psi^{T}_{\lambda,\mu\nu}\sigma_{\mu\nu}+\Psi^{AV}_{\lambda,\mu}\gamma_{\mu}\gamma_{5}+\Psi^{Pse}_{\lambda}\gamma_{5},\end{split}$
(31)

and the coefficient corresponding to $\gamma_{\mu}$ is

$\begin{split}\Psi_{\lambda,\mu}^{V}=\frac{1}{4}\text{Tr}[\gamma_{\mu}\chi_{\lambda}(K,k)].\end{split}$
(32)

Substituting Eq. (30) into (32), we can obtain the heavy vector meson wave
function in instantaneous approximation

$\begin{split}\Psi_{ij}^{D_{l}^{*}}(\mathbf{k})&=\int
dk_{4}\frac{1}{N^{D_{l}^{*}}}\exp\bigg{(}\frac{-\mathbf{k}^{2}-k_{4}^{2}}{\omega_{D_{l}^{*}}^{2}}\bigg{)}\frac{\mathbf{k}^{2}/3+k_{4}^{2}+m_{c}m_{l}}{(\mathbf{k}^{2}+k_{4}^{2}+m_{c}^{2})(\mathbf{k}^{2}+k_{4}^{2}+m_{l}^{2})}\delta_{ij}~{}~{}i,j=1,2,3.\end{split}$
(33)

In the previous works Chen et al. (2013); Chen and Lü (2015); Chen et al.
(2016), we have obtained the form factors for the vertices of heavy vector
meson $D_{l}^{*}$ coupling to scalar meson ($\sigma$)

$\begin{split}-\frac{h_{1}^{(\text{s})}(w^{2})}{2E_{1}}=&\frac{\bar{h}_{1}^{(\text{s})}(w^{2})}{2E_{2}}=F_{1}({\mathbf{w}^{2}}),~{}~{}~{}h_{2}^{(\text{s})}(w^{2})=\bar{h}_{2}^{(\text{s})}(w^{2})=0,\\\
F_{1}(\mathbf{w}^{2})=&\int\frac{d^{3}k}{(2\pi)^{3}}\bar{\Psi}^{D_{l}^{*}}\bigg{(}\mathbf{k}+\frac{2E_{c}(k)}{E_{D_{l}^{*}}+M_{D_{l}^{*}}}\mathbf{w}\bigg{)}\sqrt{\frac{E_{l}(k)+m_{l}}{E_{l}(k+w)+m_{l}}}\\\
&\times\left\\{\frac{E_{l}(k+w)-E_{l}(k)+2m_{l}}{2\sqrt{E_{l}(k+w)E_{l}(k)}}-\frac{\mathbf{k\cdot
w}}{2\sqrt{E_{l}(k+w)E_{l}(k)}[E_{l}(k)+m_{l}]}\right\\}\Psi^{D_{l}^{*}}(\mathbf{k}),\end{split}$
(34)

and to vector meson ($\rho$, $V_{1}$ and $V_{8}$)

$\begin{split}h_{1}^{(\text{v})}(w^{2})=&h_{2}^{(\text{v})}(w^{2})=\bar{h}_{1}^{(\text{v})}(w^{2})=\bar{h}_{2}^{(\text{v})}(w^{2})=F_{2}({\mathbf{w}^{2}}),~{}~{}~{}h_{3}^{(\text{v})}(w^{2})=\bar{h}_{3}^{(\text{v})}(w^{2})=0,\\\
F_{2}(\mathbf{w}^{2})=&\frac{2\sqrt{E_{D_{l}^{*}}M_{D_{l}^{*}}}}{E_{D_{l}^{*}}+M_{D_{l}^{*}}}\int\frac{d^{3}k}{(2\pi)^{3}}\bar{\Psi}^{D_{l}^{*}}\bigg{(}\mathbf{k}+\frac{2E_{c}(k)}{E_{D_{l}^{*}}+M_{D_{l}^{*}}}\mathbf{w}\bigg{)}\sqrt{\frac{E_{l}(k)+m_{l}}{E_{l}(k+w)+m_{l}}}\\\
&\times\left\\{\frac{E_{l}(k+w)+E_{l}(k)}{2\sqrt{E_{l}(k+w)E_{l}(k)}}+\frac{\mathbf{k\cdot
w}}{2\sqrt{E_{l}(k+w)E_{l}(k)}[E_{l}(k)+m_{l}]}\right\\}\Psi^{D_{l}^{*}}(\mathbf{k}),\end{split}$
(35)

where $E_{c,l}(p)=\sqrt{\mathbf{p}^{2}+m_{c,l}^{2}}$ and $\Psi^{D_{l}^{*}}$ is
the wave function of heavy vector meson expressed as Eq. (33). In this paper,
some errors in our previous works have been revised. Equations (22) and (23)
can be solved numerically with these form factors, and then the eigenvalue
equation (27) can be solved. The masses $M_{D\bar{D}}$ and wave functions of
pure bound states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ with spin-parity
quantum numbers $0^{+}$ can be obtained.

Considering the interaction kernels from one-$\rho^{\pm}$ exchange and using
the coupled-channel approach (see details in Ref. Chen et al. (2016)), we can
calculate the mass $M_{0}$ of the mixed state of two pure bound states
$D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ with $0^{+}$. Since the mixing of
component wave functions causes the change of energy, the fourth component of
$P^{D\bar{D}}$ in the original BS wave function becomes the total energy of
mixed state, and $\chi_{\lambda\tau}^{0^{+}}(P^{D\bar{D}},p)$ in Eq. (26)
becomes

$\begin{split}\chi_{\lambda\tau}^{0^{+}}(P,p)=&\frac{1}{\mathcal{N}^{0^{+}}}\\{[(p^{\prime}_{1}\cdot
p^{\prime}_{2})g_{\lambda\tau}-p^{\prime}_{2\lambda}p^{\prime}_{1\tau}]\mathcal{C}_{1}\mathcal{F}_{10}(P\cdot
p,p^{2})\\\ &+[p_{1}^{\prime 2}p_{2}^{\prime
2}g_{\lambda\tau}+(p_{1}^{\prime}\cdot
p_{2}^{\prime})p^{\prime}_{1\lambda}p^{\prime}_{2\tau}-p_{2}^{\prime
2}p^{\prime}_{1\lambda}p^{\prime}_{1\tau}-p_{1}^{\prime
2}p^{\prime}_{2\lambda}p^{\prime}_{2\tau}]\mathcal{C}_{2}\mathcal{F}_{20}(P\cdot
p,p^{2})\\}.\end{split}$ (36)

We emphasize that the mass $M_{0}$ of meson-meson bound state should not be
the mass of physical resonance. Substituting Eq. (36) into (10), we obtain BS
wave function $\chi^{D^{*}\bar{D}^{*},0^{+}}_{\lambda\tau}(P,p)$ for the mixed
state of two bound states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ with
$0^{+}$.

###  II.4 GBS wave function for four-quark state

The heavy meson is a bound state consisting of a quark and an antiquark and
the meson-meson bound state is actually composed of four quarks. We have to
give GBS wave function of meson-meson bound state as a four-quark state. If a
bound state with spin $j$ and parity $\eta_{P}$ is composed of four quarks,
its GBS wave function can be defined as Chen and Lü (2018)

$\chi^{j}_{P}(x_{1},x_{3},x_{4},x_{2})=\langle
0|T\mathcal{Q}^{\mathcal{C}}(x_{1})\bar{\mathcal{Q}}^{\mathcal{A}}(x_{3})\mathcal{Q}^{\mathcal{B}}(x_{4})\bar{\mathcal{Q}}^{\mathcal{D}}(x_{2})|P,j\rangle,$
(37)

where $P$ is the momentum of the four-quark bound state, $\mathcal{Q}$ is the
quark operator and its superscript is a flavor label. From translational
invariance, this GBS wave function can be written as

$\chi^{j}_{P}(x_{1},x_{3},x_{4},x_{2})=\frac{1}{(2\pi)^{3/2}}\frac{1}{\sqrt{2E(P)}}e^{iP\cdot
X}\chi^{j}_{P}(X^{\prime},x,x^{\prime}),$ (38)

where
$X=\eta_{1}(\eta_{1}^{\prime\prime}x_{1}+\eta_{3}^{\prime\prime}x_{3})+\eta_{2}(\eta_{4}^{\prime\prime}x_{4}+\eta_{2}^{\prime\prime}x_{2})$,
$X^{\prime}=(\eta_{1}^{\prime\prime}x_{1}+\eta_{3}^{\prime\prime}x_{3})-(\eta_{4}^{\prime\prime}x_{4}+\eta_{2}^{\prime\prime}x_{2})$,
$x=x_{1}-x_{3}$, $x^{\prime}=x_{2}-x_{4}$, $\eta_{1}+\eta_{2}=1$,
$\eta^{\prime\prime}_{1,3}=m_{\mathcal{C},\mathcal{A}}/(m_{\mathcal{C}}+m_{\mathcal{A}})$,
$\eta^{\prime\prime}_{2,4}=m_{\mathcal{D},\mathcal{B}}/(m_{\mathcal{D}}+m_{\mathcal{B}})$
and $m_{\mathcal{A},\mathcal{B},\mathcal{C},\mathcal{D}}$ are the quark
masses. In the momentum representation, GBS wave function of four-quark bound
state becomes

$\begin{split}\chi^{j}_{P}(p_{1},p_{3},p_{4},p_{2})=\frac{1}{(2\pi)^{3/2}}\frac{1}{\sqrt{2E(P)}}(2\pi)^{4}\delta^{(4)}(P-p_{1}+p_{3}-p_{4}+p_{2})\chi^{j}(P,p,k,k^{\prime}),\end{split}$
(39)

where $p_{1},p_{3},p_{4},p_{2}$ are the momenta carried by the fields
$\mathcal{Q}^{\mathcal{C}}$, $\bar{\mathcal{Q}}^{\mathcal{A}}$,
$\mathcal{Q}^{\mathcal{B}}$, $\bar{\mathcal{Q}}^{\mathcal{D}}$; $p$, $k$,
$k^{\prime}$ are the conjugate variables to $X^{\prime}$, $x$, $x^{\prime}$,
respectively; and $p=\eta_{2}(p_{1}-p_{3})+\eta_{1}(p_{2}-p_{4})$,
$k=\eta_{3}^{\prime\prime}p_{1}+\eta_{1}^{\prime\prime}p_{3}$,
$k^{\prime}=\eta_{4}^{\prime\prime}p_{2}+\eta_{2}^{\prime\prime}p_{4}$. In the
hadronic molecule structure, $p$ is the relative momentum between two mesons
in molecular state, $k$ and $k^{\prime}$ are the relative momenta between
quark and antiquark in these two mesons, respectively, shown as Fig. 2. This
work is aimed to investigate the bound state composed of two vector mesons. In
Fig. 2, $VM$ represents the vector meson with mass $M_{1}$,
$\overline{VM^{\prime}}$ represents the anti-particle of vector meson
$VM^{\prime}$ with mass $M_{2}$, and $MS$ represents the meson-meson bound
state.

Figure 2:  Generalized Bethe-Salpeter wave function for four-quark state in
the momentum representation. The solid lines denote quark propagators, and the
unfilled ellipses represent Bethe-Salpeter amplitudes.

In Fig. 2, there are three two-body systems: a meson-meson bound state and two
quark-antiquark bound states. We define BS wave functions of these two-body
systems as $\chi^{j}_{P}(p_{1}^{\prime},p_{2}^{\prime})$,
$\chi_{p_{1}^{\prime}}(p_{1},p_{3})$, $\chi_{p_{2}^{\prime}}(p_{4},p_{2})$,
respectively. BS wave function for the bound state of two vector mesons has
been given by Eq. (3) and BS wave functions of two vector mesons are

$\begin{split}\chi_{p_{1}^{\prime}}(p_{1},p_{3})_{\lambda}=\frac{1}{(2\pi)^{3/2}}\frac{1}{\sqrt{2E(p^{\prime}_{1})}}(2\pi)^{4}\delta^{(4)}(p_{1}^{\prime}-p_{1}+p_{3})\chi_{\lambda}(p_{1}^{\prime},k),\end{split}$
(40)
$\begin{split}\chi_{p_{2}^{\prime}}(p_{4},p_{2})_{\tau}=\frac{1}{(2\pi)^{3/2}}\frac{1}{\sqrt{2E(p^{\prime}_{2})}}(2\pi)^{4}\delta^{(4)}(p_{2}^{\prime}+p_{4}-p_{2})\chi_{\tau}(p_{2}^{\prime},k^{\prime}),\end{split}$
(41)

where $p_{1}^{\prime}$ and $p_{2}^{\prime}$ are the momenta of two vector
mesons, respectively, $p_{1}^{\prime}=p+\eta_{1}P$,
$p_{2}^{\prime}=p-\eta_{2}P$ and $\eta_{1,2}=M_{1,2}/(M_{1}+M_{2})$. Applying
the Feynman rules and comparing with Eq. (39), we obtain the revised GBS wave
function for four-quark state describing the bound state composed of two
vector mesons with arbitrary spin and definite parity Chen and Lü (2018); Chen
et al. (2020)

$\begin{split}\chi^{j}(P,p,k,k^{\prime})=\chi_{\lambda}(p_{1}^{\prime},k)\chi^{j}_{\lambda\tau}(P,p)\chi_{\tau}(p_{2}^{\prime},k^{\prime}).\end{split}$
(42)

From Eq. (30), we obtain BS wave functions of vector mesons

$\begin{split}&\chi_{\lambda}(p_{1}^{\prime},k)=\frac{-1}{\gamma^{\mathcal{C}}\cdot
p_{1}-im_{\mathcal{C}}}\frac{1}{\mathcal{N}^{V}}\bigg{(}\gamma_{\lambda}+p_{1\lambda}^{\prime}\frac{\gamma\cdot
p_{1}^{\prime}}{M_{V}^{2}}\bigg{)}\varphi_{V}(k^{2})\frac{-1}{\gamma^{\mathcal{A}}\cdot
p_{3}-im_{\mathcal{A}}},\\\
&\chi_{\tau}(p_{2}^{\prime},k^{\prime})=\frac{-1}{\gamma^{\mathcal{B}}\cdot
p_{4}-im_{\mathcal{B}}}\frac{1}{\mathcal{N}^{\bar{V}^{\prime}}}\bigg{(}\gamma_{\tau}+p_{2\tau}^{\prime}\frac{\gamma\cdot
p_{2}^{\prime}}{M_{\bar{V}^{\prime}}^{2}}\bigg{)}\varphi_{\bar{V}^{\prime}}(k^{\prime
2})\frac{-1}{\gamma^{\mathcal{D}}\cdot p_{2}-im_{\mathcal{D}}}.\end{split}$
(43)

In this section, we consider a mixed state of two bound states
$D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ with spin-parity quantum numbers
$0^{+}$. In Fig. 2, $VM$ and $\overline{VM^{\prime}}$ become $D_{l}^{*}$ and
$\bar{D}_{l}^{*}$, respectively, and in Eq. (37) the flavor labels
$\mathcal{C}=\mathcal{D}$ and $\mathcal{A}=\mathcal{B}$ represent $c$-quark
and $l$-quark, respectively. From Eqs. (10), (42) and (43), we obtain the GBS
wave function for meson-meson bound state as a four-quark state

$\begin{split}\chi^{D^{*}\bar{D}^{*},0^{+}}(P,p,k,k^{\prime})=&\frac{1}{\sqrt{2}}\chi^{D^{*0}}_{\lambda}(p_{1}^{\prime},k)\chi^{D^{*0}\bar{D}^{*0},0^{+}}_{\lambda\tau}(P,p)\chi^{\bar{D}^{*0}}_{\tau}(p_{2}^{\prime},k^{\prime})\\\
&+\frac{1}{\sqrt{2}}\chi^{D^{*+}}_{\lambda}(p_{1}^{\prime},k)\chi^{D^{*+}D^{*-},0^{+}}_{\lambda\tau}(P,p)\chi^{D^{*-}}_{\tau}(p_{2}^{\prime},k^{\prime}),\end{split}$
(44)

where

$\begin{split}&\chi^{D_{l}^{*}}_{\lambda}(p_{1}^{\prime},k)=\frac{-1}{\gamma\cdot
p_{1}-im_{c}}\frac{1}{\mathcal{N}^{D_{l}^{*}}}\bigg{(}\gamma_{\lambda}+p_{1\lambda}^{\prime}\frac{\gamma\cdot
p_{1}^{\prime}}{M_{D_{l}^{*}}^{2}}\bigg{)}\varphi_{D_{l}^{*}}(k^{2})\frac{-1}{\gamma\cdot
p_{3}-im_{l}},\\\
&\chi^{\bar{D}_{l}^{*}}_{\tau}(p_{2}^{\prime},k^{\prime})=\frac{-1}{\gamma\cdot
p_{4}-im_{l}}\frac{1}{\mathcal{N}^{\bar{D}_{l}^{*}}}\bigg{(}\gamma_{\tau}+p_{2\tau}^{\prime}\frac{\gamma\cdot
p_{2}^{\prime}}{M_{\bar{D}_{l}^{*}}^{2}}\bigg{)}\varphi_{\bar{D}_{l}^{*}}(k^{\prime
2})\frac{-1}{\gamma\cdot p_{2}-im_{c}}.\end{split}$ (45)

###  II.5 Normalization of BS wave function

####  II.5.1 Heavy vector meson

Here, we determine normalizations $\mathcal{N}^{D^{*0}}$ and
$\mathcal{N}^{D^{*+}}$. The authors of Refs. Ivanov et al. (1999, 2007)
employed the ladder approximation to solve the BS equation for quark-antiquark
state, and the reduced normalization condition for the BS wave function of
$D_{l}^{*}$ meson given by Eq. (30) is

$\begin{split}\frac{-i}{(2\pi)^{4}}\frac{1}{3}\int
d^{4}k\bar{\chi}_{\lambda}(K,k)\frac{\partial}{\partial
K_{0}}[S_{F}(k+K)^{-1}]S_{F}(k)^{-1}\chi_{\lambda}(K,k)=(2K_{0})^{2},\end{split}$
(46)

where $S_{F}(p)^{-1}$ is the inverse propagator for quark field and the factor
$1/3$ appears because of the sum of three transverse directions.

####  II.5.2 Molecular state

The reduced normalization condition for $\chi_{\lambda\tau}^{0^{+}}(P,p)$
expressed as Eq. (36) is

$\begin{split}&\frac{-i}{(2\pi)^{4}}\int
d^{4}p\bar{\chi}_{\tau^{\prime}\lambda^{\prime}}(P,p)\frac{\partial}{\partial
P_{0}}[\Delta_{F\lambda^{\prime}\lambda}(p+P/2)^{-1}\Delta_{F\tau\tau^{\prime}}(p-P/2)^{-1}]\chi_{\lambda\tau}(P,p)=(2P_{0})^{2},\end{split}$
(47)

where $\Delta_{F\beta\alpha^{\prime}}(p)^{-1}$ is the inverse propagator for
the vector field with mass $m$,
$\Delta_{F\beta\alpha^{\prime}}(p)^{-1}=i(\delta_{\beta\alpha^{\prime}}-\frac{p_{\beta}p_{\alpha^{\prime}}}{p^{2}+m^{2}})(p^{2}+m^{2})$
Chen and Lü (2018). After determining normalization $\mathcal{N}^{0^{+}}$, we
automatically obtain the normalized BS wave function for the mixed state of
two components $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ given by Eq. (10).
Immediately, the normalized GBS wave function for meson-meson bound state as a
four-quark state expressed as Eq. (44) is obtained.

##  III scattering matrix element from four-quark state to final state

In experiments two strong decay modes of $\chi_{c0}(3915)$ have been observed:
$J/\psi\omega$ and $D^{+}D^{-}$. The narrow state $\chi_{c0}(3915)$ was
discovered in 2005 Choi et al. (2005) by the Belle collaboration and for a
long time a series of experiments Aubert et al. (2008); Uehara et al. (2010);
Lees et al. (2012); Vinokurova et al. (2015) only observed one strong decay
mode of $\chi_{c0}(3915)$: $J/\psi\omega$ denoted as $c^{\prime}_{1}$. In 2020
the LHCb Collaboration observed another decay channel $D^{+}D^{-}$ Aaij et al.
(2020) denoted as $c^{\prime}_{2}$. Though the neutral channel
$D^{0}\bar{D}^{0}$ still has not been observed, this neutral channel should
exist for the isospin conservation, which is denoted as $c^{\prime}_{3}$. In
this section, we present the traditional technique to calculate decay width
for these processes and revise some errors in previous works Chen and Lü
(2018); Chen et al. (2020).

###  III.1 Decay channel $J/\psi\omega$ with respect to mass of bound state

Mandelstam’s approach is a technique based on BS wave function for evaluating
the general matrix element between bound states Lurié (1968). Applying
Mandelstam’s approach, we have obtained the scattering matrix element from a
four-quark state to a heavy meson plus a light meson Chen and Lü (2018) in the
momentum representation, as shown in Fig. 3. In this work, we retain only the
lowest order term of the two-particle irreducible Green’s function. In Fig. 3,
$VM$ and $\overline{VM^{\prime}}$ still represent $D_{l}^{*}$ and
$\bar{D}_{l}^{*}$, respectively; $HM$ represents $J/\psi$ with momentum $Q$
($Q^{2}=-M_{J/\psi}^{2}$) and $LM$ represents $\omega$ with momentum
$Q^{\prime}$ ($Q^{\prime 2}=-M_{\omega}^{2}$). The momentum of the initial
state is set as $P=(0,0,0,iM_{0})$ in the rest frame, and $M_{0}$ is the mass
of the mixed state of two pure bound states $D^{*0}\bar{D}^{*0}$ and
$D^{*+}D^{*-}$, which should not be the physical mass of resonance. It is
necessary to emphasize that the momenta in the final state satisfy
$Q+Q^{\prime}=P$ in this section. Here, we consider that in the final state
the light vector meson $\omega$ is an elementary particle and the heavy vector
meson $J/\psi$ is a bound state of $c\bar{c}$. From Eq. (43), we obtain the BS
wave function of heavy vector meson $J/\psi$

$\begin{split}\chi_{\nu}(Q,q)=\frac{-1}{\gamma\cdot(q+Q/2)-im_{c}}\frac{1}{\mathcal{N}^{J/\psi}}\bigg{(}\gamma_{\nu}+Q_{\nu}\frac{\gamma\cdot
Q}{M_{J/\psi}^{2}}\bigg{)}\varphi_{J/\psi}(q^{2})\frac{-1}{\gamma\cdot(q-Q/2)-im_{c}},\end{split}$
(48)

where $\varphi_{J/\psi}(q^{2})=\exp(-q^{2}/\omega_{J/\psi}^{2})$ and
$\omega_{J/\psi}$=0.826 GeV was obtained from lattice QCD (see details in Ref.
Chen and Lü (2018)). The reduced normalization condition for BS wave function
of $J/\psi$ meson expressed as Eq. (48) is

$\begin{split}\frac{-i}{(2\pi)^{4}}\frac{1}{3}\int
d^{4}q\bar{\chi}_{\nu}(Q,q)\frac{\partial}{\partial
Q_{0}}[S_{F}(q+Q/2)^{-1}S_{F}(q-Q/2)^{-1}]\chi_{\nu}(Q,q)=(2Q_{0})^{2},\end{split}$
(49)

where the factor $1/3$ appears for the three transverse directions are summed.
Normalization $\mathcal{N}^{J/\psi}$ can be determined. These momenta in Fig.
3 become

$\begin{split}&p_{1}=(Q+Q^{\prime})/2+p+k,~{}~{}p_{2}=(Q+Q^{\prime})/2-Q+p+k,~{}~{}p_{3}=k,~{}~{}p_{4}=Q^{\prime}+k,\\\
&q=Q^{\prime}/2+p+k,~{}~{}k^{\prime}=Q^{\prime}+k,~{}~{}p_{1}^{\prime}=p+P/2,~{}~{}p_{2}^{\prime}=p-P/2,~{}~{}Q+Q^{\prime}=P.\end{split}$
(50)

Figure 3:  The lowest order matrix element between bound states in the
momentum representation.

Using the Heisenberg picture, we obtain the total matrix element from the
initial state $|P~{}\text{in}\rangle$ to a final state $\langle
Q,Q^{\prime}~{}\text{out}|$

$\begin{split}-iR_{(c^{\prime}_{1};b)a}(M_{0})=\langle
Q,Q^{\prime}~{}\text{out}|P~{}\text{in}\rangle=-i(2\pi)^{4}\delta^{(4)}(Q+Q^{\prime}-P)T_{(c^{\prime}_{1};b)a}(M_{0}),\end{split}$
(51)

where $T_{(c^{\prime}_{1};b)a}(M_{0})$ is the $T$-matrix element with mass
$M_{0}$ for channel $c^{\prime}_{1}$. According to Mandelstam’s approach, we
obtain

$\begin{split}T_{(c^{\prime}_{1};b)a}(M_{0})=\frac{ig_{\omega}\varepsilon_{\mu}^{\varrho^{\prime}}(Q^{\prime})\varepsilon_{\nu}^{\varrho}(Q)}{(2\pi)^{9/2}\sqrt{2E_{J/\psi}(Q)}\sqrt{2E_{\omega}(Q^{\prime})}\sqrt{2E(P)}}\bigg{(}\frac{1}{\sqrt{2}}\mathcal{M}^{c^{\prime}_{1},D^{*0}\bar{D}^{*0}}_{\nu\mu}+\frac{1}{\sqrt{2}}\mathcal{M}^{c^{\prime}_{1},D^{*+}D^{*-}}_{\nu\mu}\bigg{)},\end{split}$
(52)

where $\varepsilon_{\nu}^{\varrho=1,2,3}(Q)$ and
$\varepsilon_{\mu}^{\varrho^{\prime}=1,2,3}(Q^{\prime})$ are the polarization
vectors of $J/\psi$ and $\omega$, respectively, $\varepsilon^{\varrho}(Q)\cdot
Q=\varepsilon^{\varrho^{\prime}}(Q^{\prime})\cdot Q^{\prime}=0$, and

$\begin{split}\mathcal{M}^{c^{\prime}_{1},D^{*}_{l}\bar{D}^{*}_{l}}_{\nu\mu}=&\int\frac{d^{4}kd^{4}p}{(2\pi)^{8}}\frac{1}{N^{J/\psi}}\frac{\varphi_{J/\psi}(q^{2})}{p_{2}^{2}+m^{2}_{c}}\frac{1}{p_{1}^{2}+m^{2}_{c}}\frac{1}{N^{D_{l}^{*}}}\frac{\varphi_{D_{l}^{*}}(k^{2})}{p_{3}^{2}+m^{2}_{l}}\frac{1}{N^{\bar{D}_{l}^{*}}}\frac{\varphi_{\bar{D}_{l}^{*}}(k^{\prime
2})}{p_{4}^{2}+m^{2}_{l}}\\\ &\times\text{Tr}[(\gamma\cdot
p_{2}+im_{c})\gamma_{\nu}(\gamma\cdot
p_{1}+im_{c})\gamma_{\lambda}(\gamma\cdot
p_{3}+im_{l})\gamma_{\mu}(\gamma\cdot
p_{4}+im_{l})\gamma_{\tau}\chi_{\lambda\tau}^{0^{+}}(P,p)].\end{split}$ (53)

Here $\chi_{\lambda\tau}^{0^{+}}(P,p)$ is expressed as Eq. (36). In Eq. (53)
the trace of the product of 8 $\gamma$-matrices contains 105 terms and the
resulting expression has been given in Appendix B of Ref. Chen and Lü (2018).
In our approach, the $p$ integral is computed in instantaneous approximation.
To calculate this tensor
$\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}$, we have
given a simple method in Ref. Chen and Lü (2018). It is obvious that the
tensor $\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}$ only
depends on $Q$ and $Q^{\prime}$, so in Minkowski space
$\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}$ can be
expressed as

$\begin{split}&\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}\\\
&=g_{\nu\mu}U_{1}(Q^{\prime},Q)+Q^{\prime}_{\nu}Q_{\mu}U_{2}(Q^{\prime},Q)+Q^{\prime}_{\nu}Q^{\prime}_{\mu}U_{3}(Q^{\prime},Q)+Q_{\nu}Q^{\prime}_{\mu}U_{4}(Q^{\prime},Q)+Q_{\nu}Q_{\mu}U_{5}(Q^{\prime},Q),\end{split}$
(54)

where $U_{i}(Q^{\prime},Q)(i=1,\cdots,5)$ are scalar functions. The above
expression is multiplied by these tensor structures $g_{\nu\mu}$,
$Q^{\prime}_{\nu}Q_{\mu}$, $Q^{\prime}_{\nu}Q^{\prime}_{\mu}$,
$Q_{\nu}Q^{\prime}_{\mu}$, $Q_{\nu}Q_{\mu}$, respectively; and a set of
equations is obtained

$\begin{split}g_{\nu\mu}\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}&=U_{1}^{\prime}=4U_{1}+(Q^{\prime}\cdot
Q)U_{2}+Q^{\prime 2}U_{3}+(Q^{\prime}\cdot Q)U_{4}+Q^{2}U_{5},\\\
Q^{\prime}_{\nu}Q_{\mu}\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}&=U_{2}^{\prime}=(Q^{\prime}\cdot
Q)U_{1}+Q^{\prime 2}Q^{2}U_{2}+Q^{\prime 2}(Q^{\prime}\cdot
Q)U_{3}+(Q^{\prime}\cdot Q)^{2}U_{4}+Q^{2}(Q^{\prime}\cdot Q)U_{5},\\\
Q^{\prime}_{\nu}Q^{\prime}_{\mu}\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}&=U_{3}^{\prime}=Q^{\prime
2}U_{1}+Q^{\prime 2}(Q^{\prime}\cdot Q)U_{2}+Q^{\prime 2}Q^{\prime
2}U_{3}+Q^{\prime 2}(Q^{\prime}\cdot Q)U_{4}+(Q^{\prime}\cdot Q)^{2}U_{5},\\\
Q_{\nu}Q^{\prime}_{\mu}\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}&=U_{4}^{\prime}=(Q^{\prime}\cdot
Q)U_{1}+(Q^{\prime}\cdot Q)^{2}U_{2}+Q^{\prime 2}(Q^{\prime}\cdot
Q)U_{3}+Q^{2}Q^{\prime 2}U_{4}+Q^{2}(Q^{\prime}\cdot Q)U_{5},\\\
Q_{\nu}Q_{\mu}\mathcal{M}^{c^{\prime}_{1},D_{l}^{*}\bar{D}_{l}^{*}}_{\nu\mu}&=U_{5}^{\prime}=Q^{2}U_{1}+Q^{2}(Q^{\prime}\cdot
Q)U_{2}+(Q^{\prime}\cdot Q)^{2}U_{3}+Q^{2}(Q^{\prime}\cdot
Q)U_{4}+Q^{2}Q^{2}U_{5},\end{split}$ (55)

where $U^{\prime}_{i}$ are numbers. Subsequently, we numerically calculate
$U^{\prime}_{i}$ and solve this set of equations. The values of $U_{i}$ can be
obtained.

Then we can obtain the decay width with mass of meson-meson bound state for
channel $J/\psi\omega$

$\begin{split}\Gamma_{1}(M_{0})=&\int
d^{3}Qd^{3}Q^{\prime}(2\pi)^{4}\delta^{(4)}(Q+Q^{\prime}-P)\sum_{\varrho^{\prime}=1}^{3}\sum_{\varrho=1}^{3}(2\pi)^{3}|T_{(c^{\prime}_{1};b)a}(M_{0})|^{2},\end{split}$
(56)

where $P=(0,0,0,iM_{0})$,
$Q=(\mathbf{Q}(M_{0}),i\sqrt{\mathbf{Q}^{2}(M_{0})+M^{2}_{J/\psi}})$,
$Q^{\prime}=(-\mathbf{Q}(M_{0}),i\sqrt{\mathbf{Q}^{2}(M_{0})+M^{2}_{\omega}})$
and
$\mathbf{Q}^{2}(M_{0})=[M_{0}^{2}-(M_{J/\psi}+M_{\omega})^{2}][M_{0}^{2}-(M_{J/\psi}-M_{\omega})^{2}]/(4M_{0}^{2})$.
To calculate the decay width $\Gamma_{1}(M_{0})$, we use the transverse
condition $\varepsilon^{\varrho}(Q)\cdot
Q=\varepsilon^{\varrho^{\prime}}(Q^{\prime})\cdot Q^{\prime}=0$ and the
completeness relation. It is necessary to emphasize that the decay width
$\Gamma_{1}(M_{0})$ expressed as Eq. (56) is not the decay width of physical
resonance.

###  III.2 Decay channel $D^{+}D^{-}$ with respect to mass of bound state

Considering the lowest order term of the two-particle irreducible Green’s
function, we obtain the interaction between two heavy vector mesons derived
from a light meson exchange. Applying Mandelstam’s approach, we can obtain the
$T$-matrix element with mass $M_{0}$ for channel $c^{\prime}_{2}$, which can
be represented graphically by Fig. 4. In Fig. 4, $PM$ and
$\overline{PM^{\prime}}$ represent pseudoscalar mesons $D^{+}$ and $D^{-}$,
respectively; $Q_{1}$ and $Q_{2}$ represent the momenta of final particles,
$Q_{1}^{2}=-M_{D^{+}}^{2}$, $Q_{2}^{2}=-M_{D^{-}}^{2}$ and in this section
$Q_{1}+Q_{2}=P$.

Figure 4:  Matrix element for decay channel $D^{+}D^{-}$. The momenta in the
final state satisfy $Q_{1}+Q_{2}=P$. $w$ represents the momentum of the
exchanged light meson.

To simplify the computational process, we use the vertex function for the
exchanged light meson, heavy pseudoscalar and vector mesons; and then Fig. 4
can be reduced to Fig. 5. From the Lorentz-structure, we obtain the matrix
elements of quark scalar density $J$ and quark current $J_{\alpha}$ between
heavy pseudoscalar and vector mesons

$\begin{split}\langle
PM(Q_{1})|J|VM^{\vartheta}(p_{1}^{\prime})\rangle=&\frac{1}{\sqrt{2E_{D^{+}}(Q_{1})}\sqrt{2E_{D^{*}_{l}}(p_{1}^{\prime})}}[Q_{1}\cdot\varepsilon^{\vartheta}(p_{1}^{\prime})]h^{(s)}_{4}(w^{2}),\end{split}$
(57a)
$\begin{split}\langle\overline{PM^{\prime}}(Q_{2})|J|\overline{VM^{\prime}}^{\vartheta^{\prime}}(-p_{2}^{\prime})\rangle=&\frac{1}{2\sqrt{E_{D^{-}}(Q_{2})E_{\bar{D}^{*}_{l}}(-p_{2}^{\prime})}}[Q_{2}\cdot\varepsilon^{\vartheta^{\prime}}(-p_{2}^{\prime})]\bar{h}^{(s)}_{4}(w^{2}),\end{split}$
(57b)
$\begin{split}\langle&PM(Q_{1})|J_{\alpha}|VM^{\vartheta}(p_{1}^{\prime})\rangle=\frac{1}{2\sqrt{E_{D^{+}}(Q_{1})E_{D^{*}_{l}}(p_{1}^{\prime})}}\\\
&\times\\{h_{4}^{(\text{v})}(w^{2})\\{[Q_{1}\cdot\varepsilon^{\vartheta}(p_{1}^{\prime})](Q_{1}+p_{1}^{\prime})_{\alpha}-[(Q_{1}+p_{1}^{\prime})\cdot(Q_{1}-p_{1}^{\prime})]\varepsilon_{\alpha}^{\vartheta}(p_{1}^{\prime})\\}\\\
&-h^{(\text{v})}_{5}(w^{2})\\{[Q_{1}\cdot\varepsilon^{\vartheta}(p_{1}^{\prime})](Q_{1}-p_{1}^{\prime})_{\alpha}-(Q_{1}-p_{1}^{\prime})^{2}\varepsilon_{\alpha}^{\vartheta}(p_{1}^{\prime})\\}\\},\end{split}$
(57c)
$\begin{split}\langle&\overline{PM^{\prime}}(Q_{2})|J_{\beta}|\overline{VM^{\prime}}^{\vartheta^{\prime}}(-p_{2}^{\prime})\rangle=\frac{1}{2\sqrt{E_{D^{-}}(Q_{2})E_{\bar{D}^{*}_{l}}(-p_{2}^{\prime})}}\\\
&\times\\{\bar{h}_{4}^{(\text{v})}(w^{2})\\{[Q_{2}\cdot\varepsilon^{\vartheta^{\prime}}(-p_{2}^{\prime})](Q_{2}-p_{2}^{\prime})_{\beta}-[(Q_{2}-p_{2}^{\prime})\cdot(Q_{2}+p_{2}^{\prime})]\varepsilon^{\vartheta^{\prime}}_{\beta}(-p_{2}^{\prime})\\}\\\
&-\bar{h}^{(\text{v})}_{5}(w^{2})\\{[Q_{2}\cdot\varepsilon^{\vartheta^{\prime}}(-p_{2}^{\prime})](Q_{2}+p_{2}^{\prime})_{\beta}-(Q_{2}+p_{2}^{\prime})^{2}\varepsilon^{\vartheta^{\prime}}_{\beta}(-p_{2}^{\prime})\\}\\},\end{split}$
(57d)

where $p^{\prime}_{1}=(\mathbf{p},ip_{10}^{\prime})$,
$p^{\prime}_{2}=(\mathbf{p},ip_{20}^{\prime})$,
$Q_{1}=(\mathbf{Q}_{D},iQ_{10})$, $Q_{2}=(-\mathbf{Q}_{D},iQ_{20})$,
$w=p_{1}^{\prime}-Q_{1}=p_{2}^{\prime}+Q_{2}=p-(Q_{1}-Q_{2})/2$ is the
momentum of the light meson and $\mathbf{w}=\mathbf{p}-\mathbf{Q}_{D}$;
$E_{D^{\pm}}(p)=\sqrt{\mathbf{p}^{2}+M_{D^{\pm}}^{2}}$, $h(w^{2})$ and
$\bar{h}(w^{2})$ become vertex functions.

Figure 5:  Reduced matrix element for decay channel $D^{+}D^{-}$.

Now, we introduce the vertex function for the exchanged light meson, heavy
pseudoscalar and vector mesons, shown as Fig. 6. The charmed meson $D^{+}$ is
composed of $c$-quark and $d$-antiquark. For heavy pseudoscalar mesons, the
authors of Refs. Burden et al. (1997); Maris et al. (1998); Ivanov et al.
(1999, 2007) also gave their BS amplitudes in Euclidean space:

$\begin{split}\Gamma^{P}(K,k)=\frac{1}{\mathcal{N}^{P}}i\gamma_{5}\varphi_{P}(k^{2}),\end{split}$
(58)

where $K$ is the momentum of heavy meson, $k$ denotes the relative momentum
between quark and antiquark in heavy meson, $\mathcal{N}^{P}$ is
normalization, and $\varphi_{P}(k^{2})$ is scalar function fixed by providing
fits to observables. Using the approach introduced in Sec. II.3.2, we can
obtain the heavy pseudoscalar meson wave function in instantaneous
approximation

$\begin{split}\Psi^{D^{+}}(\textbf{k})&=\int
dk_{4}\frac{1}{4}\text{Tr}\bigg{\\{}\gamma_{5}\frac{-1}{\gamma\cdot(k+K)-im_{c}}\frac{1}{\mathcal{N}^{D^{+}}}i\gamma_{5}\varphi_{D^{+}}(k^{2})\frac{-1}{\gamma\cdot
k-im_{d}}\bigg{\\}}\\\ &=\int
dk_{4}\frac{-i}{N^{D^{+}}}\exp\bigg{(}\frac{-\textbf{k}^{2}-k_{4}^{2}}{\omega_{D}^{2}}\bigg{)}\frac{\textbf{k}^{2}+k_{4}^{2}+m_{c}m_{d}}{(\textbf{k}^{2}+k_{4}^{2}+m_{c}^{2})(\textbf{k}^{2}+k_{4}^{2}+m_{d}^{2})},\end{split}$
(59)

where $m_{c,d}$ are the constituent quark masses,
$\varphi_{D^{+}}(k^{2})=\varphi_{D^{-}}(k^{2})=\exp(-k^{2}/\omega_{D}^{2})$
and $\omega_{D}$=1.50 GeV Ivanov et al. (2007). Then we can apply the method
given in Refs. Chen et al. (2013); Chen and Lü (2015); Chen et al. (2016) to
obtain the explicit forms of the vertex function for heavy pseudoscalar meson
$D$ and vector meson $D_{l}^{*}$ coupling to scalar meson ($\sigma$)

$\begin{split}-\frac{h^{(s)}_{4}(w^{2})}{2E_{1}}=&\frac{\bar{h}_{4}^{(s)}(w^{2})}{2E_{2}}=F_{4}(\textbf{w}^{2}),\\\
F_{4}(\textbf{w}^{2})=&\frac{1}{\sqrt{M_{D^{*+}}^{2}-M_{D^{+}}^{2}}}\int\frac{d^{3}k}{(2\pi)^{3}}\bar{\Psi}^{D^{+}}\bigg{(}\textbf{k}+\frac{E_{c}(k)}{E_{D^{+}}}\textbf{w}\bigg{)}\sqrt{\frac{E_{d}(k)+m_{d}}{E_{d}(k+w)+m_{d}}}\\\
&\times\left\\{\frac{E_{d}(k+w)-E_{d}(k)+2m_{d}}{2\sqrt{E_{d}(k+w)E_{d}(k)}}-\frac{\textbf{k}\cdot\textbf{w}}{2\sqrt{E_{d}(k+w)E_{d}(k)}[E_{d}(k)+m_{d}]}\right\\}\Psi^{D^{*+}}(\textbf{k}),\end{split}$
(60)

and to vector meson ($\rho$, $V_{1}$ and $V_{8}$)

$\begin{split}h_{4}^{(\text{v})}(w^{2})=&h_{5}^{(\text{v})}(w^{2})=\bar{h}_{4}^{(\text{v})}(w^{2})=\bar{h}_{5}^{(\text{v})}(w^{2})=F_{5}({\textbf{w}^{2}}),\\\
F_{5}(\textbf{w}^{2})=&\frac{1}{\sqrt{M_{D_{l}^{*}}^{2}-M_{D^{+}}^{2}}}\frac{2\sqrt{E_{D_{l}^{*}}E_{D^{+}}}}{E_{D_{l}^{*}}+E_{D^{+}}}\int\frac{d^{3}k}{(2\pi)^{3}}\bar{\Psi}^{D^{+}}\bigg{(}\textbf{k}+\frac{E_{c}(k)}{E_{D^{+}}}\textbf{w}\bigg{)}\sqrt{\frac{E_{l}(k)+m_{l}}{E_{d}(k+w)+m_{d}}}\\\
&\times\left\\{\frac{E_{d}(k+w)+E_{l}(k)}{2\sqrt{E_{d}(k+w)E_{l}(k)}}+\frac{\textbf{k}\cdot\textbf{w}}{2\sqrt{E_{d}(k+w)E_{l}(k)}[E_{l}(k)+m_{l}]}\right\\}\Psi^{D_{l}^{*}}({\textbf{k}}),\end{split}$
(61)

Figure 6:  Vertex function for the exchanged light meson, heavy pseudoscalar
and vector mesons.

where $E_{d}(p)=\sqrt{\mathbf{p}^{2}+m_{d}^{2}}$, $\Psi^{D^{+}}$ and
$\Psi^{D_{l}^{*}}$ are the wave functions of heavy pseudoscalar and vector
mesons expressed as Eqs. (59) and (33), respectively.

Taking away the external lines including normalizations and polarization
vectors $\varepsilon^{\vartheta}_{\lambda}(p_{1}^{\prime})$,
$\varepsilon^{\vartheta^{\prime}}_{\tau}(-p_{2}^{\prime})$ in Eq. (57), we
obtain the interaction from one light meson ($\sigma$, $\rho^{0}$, $V_{1}$ and
$V_{8}$) exchange

$\begin{split}\mathcal{V}_{\lambda\tau}(Q_{1},Q_{2},\textbf{p})=&-2E_{1}F_{4}(\textbf{w}^{2})\frac{-ig_{\sigma}^{2}}{w^{2}+M_{\sigma}^{2}}2E_{2}F_{4}(\textbf{w}^{2})Q_{1\lambda}Q_{2\tau}\\\
&-4F_{5}(\textbf{w}^{2})\left(\frac{-ig_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}+\frac{-ig_{1}^{2}}{w^{2}+M_{\omega}^{2}}+\frac{-ig_{8}^{2}}{w^{2}+M_{\phi}^{2}}\right)F_{5}(\textbf{w}^{2})(p_{1}^{\prime}\cdot
p_{2}^{\prime})Q_{1\lambda}Q_{2\tau},\\\ \end{split}$ (62)

where $E_{1}=E_{2}=M_{0}/2$, and $w=(\mathbf{w},0)$. The interaction from
one-$\rho^{\pm}$ exchange becomes

$\begin{split}\mathcal{V}^{\prime}_{\lambda\tau}(Q_{1},Q_{2},\textbf{p})=-4F_{5}(\textbf{w}^{2})\frac{-i2g_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}F_{5}(\textbf{w}^{2})(p_{1}^{\prime}\cdot
p_{2}^{\prime})Q_{1\lambda}Q_{2\tau}.\\\ \end{split}$ (63)

These momenta in Fig. 5 become

$\begin{split}w=(\mathbf{p}-\mathbf{Q}_{D}(M_{0}),0),~{}Q_{1}+Q_{2}=P,~{}p_{1}^{\prime}=p+P/2,~{}p_{2}^{\prime}=p-P/2,\end{split}$
(64)

where $P=(0,0,0,iM_{0})$, $Q_{1}=(\mathbf{Q}_{D}(M_{0}),iM_{0}/2)$,
$Q_{2}=(-\mathbf{Q}_{D}(M_{0}),iM_{0}/2)$ and
$\mathbf{Q}_{D}^{2}(M_{0})=[M_{0}^{2}-(M_{D^{+}}+M_{D^{-}})^{2}]/4$. For decay
channel $D^{+}D^{-}$, we obtain the total matrix element

$\begin{split}-iR_{(c^{\prime}_{2};b)a}(M_{0})=\langle
Q_{1},Q_{2}~{}\text{out}|P~{}\text{in}\rangle=-i(2\pi)^{4}\delta^{(4)}(Q_{1}+Q_{2}-P)T_{(c^{\prime}_{2};b)a}(M_{0}),\end{split}$
(65)

where $T_{(c^{\prime}_{2};b)a}(M_{0})$ is the $T$-matrix element with mass
$M_{0}$ for channel $c^{\prime}_{2}$. From Fig. 5, we obtain

$\begin{split}T_{(c^{\prime}_{2};b)a}(M_{0})=\frac{-i}{(2\pi)^{9/2}\sqrt{2E_{D^{+}}(Q_{1})}\sqrt{2E_{D^{-}}(Q_{2})}\sqrt{2E(P)}}\bigg{(}\frac{1}{\sqrt{2}}\mathcal{M}^{c^{\prime}_{2}}+\frac{1}{\sqrt{2}}\mathcal{M}^{\prime
c^{\prime}_{2}}\bigg{)},\end{split}$ (66)

where

$\begin{split}\mathcal{M}^{c^{\prime}_{2}}=&\int\frac{d^{4}p}{(2\pi)^{4}}\mathcal{V}_{\lambda\tau}(Q_{1},Q_{2},\mathbf{p})\chi_{\lambda\tau}^{0^{+}}(P,p),\end{split}$
(67a) $\begin{split}\mathcal{M}^{\prime
c^{\prime}_{2}}=&\int\frac{d^{4}p}{(2\pi)^{4}}\mathcal{V}^{\prime}_{\lambda\tau}(Q_{1},Q_{2},\mathbf{p})\chi_{\lambda\tau}^{0^{+}}(P,p).\end{split}$
(67b)

Here $\chi_{\lambda\tau}^{0^{+}}(P,p)$ is expressed as Eq. (36). The $p$
integral is also computed in instantaneous approximation. Then the decay width
with mass of meson-meson bound state for channel $D^{+}D^{-}$ becomes

$\begin{split}\Gamma_{2}(M_{0})=&\int
d^{3}Q_{1}d^{3}Q_{2}(2\pi)^{4}\delta^{(4)}(Q_{1}+Q_{2}-P)(2\pi)^{3}|T_{(c^{\prime}_{2};b)a}(M_{0})|^{2}.\end{split}$
(68)

The decay width $\Gamma_{2}(M_{0})$ also is not the width of physical
resonance.

###  III.3 Decay channel $D^{0}\bar{D}^{0}$ with respect to mass of bound
state

Since the $\chi_{c0}(3915)$ state is a isoscalar, there should exist the
neutral channel $D^{0}\bar{D}^{0}$. In Figs. 4, 5 and 6, $PM$ and
$\overline{PM^{\prime}}$ represent pseudoscalar mesons $D^{0}$ and
$\bar{D}^{0}$, respectively. Following the same procedure as for charged
channel $D^{+}D^{-}$, we can obtain the $T$-matrix element
$T_{(c^{\prime}_{3};b)a}(M_{0})$ and the decay width $\Gamma_{3}(M_{0})$ with
mass $M_{0}$ for neutral channel $c^{\prime}_{3}$. The decay width
$\Gamma_{3}(M_{0})$ should not be the width of physical resonance.

##  IV the developed Bethe-Salpeter theory

Secs. II and III give the traditional technique to deal with molecular state
in present particle physics. These masses of meson-meson bound states were
regarded as masses of resonances Swanson (2004); Törnqvist (2004); Liu and Zhu
(2009); Chen and Lü (2015); Zhao et al. (2022) and used to calculate decay
widths of resonances Chen and Lü (2018); Chen et al. (2020), which should not
be impeccable. To deal with resonance in the framework of relativistic quantum
field theory, we considered the time evolution of molecular state as
determined by the total Hamiltonian and provided the developed Bethe-Salpeter
theory in Ref. Chen and Lü (2023).

Because the time evolution of molecular state is determined by the total
Hamiltonian, exotic meson resonance should be considered as an unstable meson-
meson molecular state. According to the developed Bethe-Salpeter theory for
dealing with resonance Chen and Lü (2023), this unstable state has been
prepared to decay at given time, and the prepared state can be regarded as a
bound state with ground-state energy. Solving BS equation for arbitrary meson-
meson bound state, one can obtain the mass $M_{0}$ and BS wave function
$\chi_{P}(x^{\prime}_{1},x^{\prime}_{2})$ for this bound state with momentum
$P=(\mathbf{P},i\sqrt{\mathbf{P}^{2}+M_{0}^{2}})$. Setting $t_{1}=0$ and
$t_{2}=0$ in the ground-state BS wave function, we obtain a description for
the prepared state (ps)

$\begin{split}\mathscr{X}^{\text{ps}}_{a}=\chi_{P}(\mathbf{x}_{1}^{\prime},t_{1}=0,\mathbf{x}_{2}^{\prime},t_{2}=0)=\frac{1}{(2\pi)^{3/2}}\frac{1}{\sqrt{2E(P)}}e^{i\mathbf{P}\cdot\mathbf{X}}\chi_{P}(\mathbf{X}^{\prime}).\end{split}$
(69)

Now it is necessary to consider the total Hamiltonian

$\begin{split}H=K_{I}+V_{I},\end{split}$ (70)

where $K_{I}$ represents the interaction responsible for the formation of
stationary bound state and $V_{I}$ stands for the interaction responsible for
the decay of resonance. Then the time evolution of this system determined by
the total Hamiltonian $H$ has the explicit form

$\begin{split}\mathscr{X}(t)=e^{-iHt}\mathscr{X}^{\text{ps}}_{a}=\frac{1}{2\pi
i}\int_{C_{2}}d\epsilon e^{-i\epsilon
t}\frac{1}{\epsilon-H}\mathscr{X}^{\text{ps}}_{a},\end{split}$ (71)

where $(\epsilon-H)^{-1}$ is the Green’s function and the contour $C_{2}$ runs
from $ic_{r}+\infty$ to $ic_{r}-\infty$ in energy-plane. The positive constant
$c_{r}$ is sufficiently large that no singularity of $(\epsilon-H)^{-1}$ lies
above $C_{2}$. The time-dependent wave function $\mathscr{X}(t)$ provides a
complete description of the system for $t>0$. Since $H\neq K_{I}$, this system
should not remain in the prepared state $\mathscr{X}^{\text{ps}}_{a}$. Then at
arbitrary time $t$ the probability amplitude of finding the system in the
state $\mathscr{X}^{\text{ps}}_{a}$ is

$\begin{split}\mathscr{A}_{a}=(\mathscr{X}^{\text{ps}}_{a},\mathscr{X}(t))=\frac{1}{2\pi
i}\int_{C_{2}}d\epsilon\frac{e^{-i\epsilon t}}{\epsilon-
M_{0}-(2\pi)^{3}T_{aa}(\epsilon)}.\end{split}$ (72)

In field theory the operator $T(\epsilon)$ is just the scattering matrix with
energy $\epsilon$, and $T_{aa}(\epsilon)$ is the $T$-matrix element between
two bound states, which is defined as

$\begin{split}\langle a~{}\text{out}|a~{}\text{in}\rangle=\langle
a~{}\text{in}|a~{}\text{in}\rangle-i(2\pi)^{4}\delta^{(4)}(P-P)T_{aa}(\epsilon).\end{split}$
(73)

Because of the analyticity of $T_{aa}(\epsilon)$, we define

$\begin{split}T_{aa}(\epsilon)=\mathbb{D}(\epsilon)-i\mathbb{I}(\epsilon),\end{split}$
(74)

where $\epsilon$ approaches the real axis from above, $\mathbb{D}$ and
$\mathbb{I}$ are the real and imaginary parts, respectively. In experiments,
many exotic particles are narrow states and their decay widths are very small
compared with their energy levels, i.e., $(2\pi)^{3}\mathbb{I}(M_{0})\ll
M_{0}$. This situation is ordinarily interpreted as implying that both
$(2\pi)^{3}|\mathbb{D}(\epsilon)|$ and $(2\pi)^{3}\mathbb{I}(\epsilon)$ are
also very small quantities, as compared to $M_{0}$. Therefore, we can expect
that $[\epsilon-M_{0}-(2\pi)^{3}T_{aa}(\epsilon)]^{-1}$ has a pole on the
second Riemann sheet

$\begin{split}\epsilon_{\text{pole}}\cong
M_{0}+(2\pi)^{3}[\mathbb{D}(M_{0})-i\mathbb{I}(M_{0})]=M-i\frac{\Gamma(M_{0})}{2},\end{split}$
(75)

where $\Delta M=(2\pi)^{3}\mathbb{D}(M_{0})$ is the correction for energy
level of resonance and $M=M_{0}+(2\pi)^{3}\mathbb{D}(M_{0})$ is the physical
mass for resonance. This pole at $\epsilon_{\text{pole}}$ describes the
resonance. The mass $M_{0}$ of two-body bound state is obtained by solving
homogeneous BS equation, which should not be the mass of physical resonance.
$\Gamma(M_{0})$ with mass $M_{0}$ also should not be the width of physical
resonance, which should depend on its physical mass $M$. We will minutely show
the computational process of $T$-matrix element between two bound states
$T_{aa}(\epsilon)$ in the next section.

##  V $T$-matrix element $T_{aa}(\epsilon)$

When there is only one decay channel, we can use the unitarity of
$T_{aa}(\epsilon)$ to obtain Goldberger and Watson (1964)

$\begin{split}2\mathbb{I}(\epsilon)=\sum_{b}(2\pi)^{4}\delta^{(3)}(\mathbf{P}_{b}-\mathbf{P})\delta(E_{b}-\epsilon)|T_{ba}(\epsilon)|^{2},\end{split}$
(76)

where $P_{b}=(\mathbf{P}_{b},iE_{b})$ is the total energy-momentum vector of
all particles in the final state and the $T$-matrix element $T_{ba}(\epsilon)$
is defined as $\langle
b~{}\text{out}|a~{}\text{in}\rangle=-i(2\pi)^{4}\delta^{(3)}(\mathbf{P}_{b}-\mathbf{P})\delta(E_{b}-\epsilon)T_{ba}(\epsilon)$.
The delta-function in Eq. (76) means that the energy $\epsilon$ in scattering
matrix is equal to the total energy $E_{b}$ of the final state, and $\sum_{b}$
represents summing over momenta and spins of all particles in the final state.
For $E_{b}=\epsilon$, we also denote the total energy of the final state by
$\epsilon$ and $\mathbb{I}(\epsilon)$ becomes a function of the final state
energy. Using dispersion relation for the function $T_{aa}(\epsilon)$, we
obtain

$\begin{split}\mathbb{D}(\epsilon)=-\frac{\mathcal{P}}{\pi}\int_{\epsilon_{M}}^{\infty}\frac{\mathbb{I}(\epsilon^{\prime})}{\epsilon^{\prime}-\epsilon}d\epsilon^{\prime},\end{split}$
(77)

where the symbol $\mathcal{P}$ means that this integral is a principal value
integral and the variable of integration is the total energy
$\epsilon^{\prime}$ of the final state. To calculate the real part, we need
calculate the function $\mathbb{I}(\epsilon^{\prime})$ of value of the final
state energy $\epsilon^{\prime}$, which is an arbitrary real number over the
real interval $\epsilon_{M}<\epsilon^{\prime}<\infty$. As usual the momentum
of initial bound state $a$ is set as $P=(0,0,0,iM_{0})$ in the rest frame and
$\epsilon_{M}$ denotes the sum of all particle masses in the final state. We
suppose that the final state $b$ may contain $n$ composite particles and
$n^{\prime}$ elementary particles in decay channel $c^{\prime}$. From Eq.
(76), we have

$\begin{split}\mathbb{I}(\epsilon^{\prime})=&\frac{1}{2}\int
d^{3}Q^{\prime}_{1}\cdots d^{3}Q^{\prime}_{n^{\prime}}d^{3}Q_{1}\cdots
d^{3}Q_{n}(2\pi)^{4}\delta^{(4)}(Q^{\prime}_{1}+\cdots+Q_{n}-P^{\epsilon^{\prime}})\sum_{\text{spins}}|T_{(c^{\prime};b)a}(\epsilon^{\prime})|^{2},\end{split}$
(78)

where $Q^{\prime}_{1}\cdots Q^{\prime}_{n^{\prime}}$ and $Q_{1}\cdots Q_{n}$
are the momenta of final elementary and composite particles, respectively;
$P^{\epsilon^{\prime}}=(0,0,0,i\epsilon^{\prime})$,
$T_{(c^{\prime};b)a}(\epsilon^{\prime})$ is the $T$-matrix element with
respect to $\epsilon^{\prime}$, and $\sum_{\text{spins}}$ represents summing
over spins of all particles in the final state. In Eq. (78) the energy in
scattering matrix is equal to the total energy $\epsilon^{\prime}$ of the
final state $b$, which is an arbitrary real number over the real interval
$\epsilon_{M}<\epsilon^{\prime}<\infty$. The mass $M_{0}$ and BS amplitude of
initial bound state $a$ have been specified and the value of the initial state
energy in the rest frame is a specified value $M_{0}$. From Eq. (78), we have
$\mathbb{I}(\epsilon^{\prime})>0$ for $\epsilon^{\prime}>\epsilon_{M}$ and
$\mathbb{I}(\epsilon^{\prime})=0$ for
$\epsilon^{\prime}\leqslant\epsilon_{M}$, which is the reason that the
integration in dispersion relation (77) ranges from $\epsilon_{M}$ to
$+\infty$.

If there are several decay channels, we should write instead

$\begin{split}\mathbb{I}(\epsilon^{\prime})=&\frac{1}{2}\sum_{c^{\prime}}\int
d^{3}Q^{\prime}_{1}\cdots d^{3}Q^{\prime}_{n^{\prime}}d^{3}Q_{1}\cdots
d^{3}Q_{n}(2\pi)^{4}\delta^{(4)}(Q^{\prime}_{1}+\cdots+Q_{n}-P^{\epsilon^{\prime}})\sum_{\text{spins}}|T_{(c^{\prime};b)a}(\epsilon^{\prime})|^{2},\end{split}$
(79)

where $\sum_{c^{\prime}}$ represents summing over all open and closed
channels. Because the total energy $\epsilon^{\prime}$ of the final state
extends from $\epsilon_{M}$ to $+\infty$, we may obtain several closed
channels derived from the interaction Lagrangian. Assuming that resonance
$\chi_{c0}(3915)$ is a mixed state of two components $D^{*0}\bar{D}^{*0}$ and
$D^{*+}D^{*-}$, we obtain one closed channel $D^{*}\bar{D}^{*}$ derived from
the interaction Lagrangian (1), denoted as $c^{\prime}_{4}$. Since bound state
lies below the threshold, i.e., $M_{0}<M_{D^{*}}+M_{\bar{D}^{*}}$, the closed
channel $c^{\prime}_{4}$ can not occur inside the physical world.

###  V.1 Channel $J/\psi\omega$ with respect to arbitrary value of the final
state energy

From Eq. (51), we obtain the total matrix element between the final state
$\langle J/\psi(Q),\omega(Q^{\prime})~{}\text{out}|$ and the specified initial
four-quark state $|P~{}\text{in}\rangle$

$\begin{split}-iR_{(c^{\prime}_{1};b)a}(\epsilon^{\prime})=\langle
Q,Q^{\prime}~{}\text{out}|P~{}\text{in}\rangle=-i(2\pi)^{4}\delta^{(4)}(Q+Q^{\prime}-P^{\epsilon^{\prime}})T_{(c^{\prime}_{1};b)a}(\epsilon^{\prime}),\end{split}$
(80)

where the total energy $\epsilon^{\prime}$ of the final state extends from
$\epsilon_{c^{\prime}_{1},M}$ to $+\infty$, i.e.,
$\epsilon_{c^{\prime}_{1},M}<\epsilon^{\prime}<\infty$ and
$\epsilon_{c^{\prime}_{1},M}=M_{J/\psi}+M_{\omega}$.
$T_{(c^{\prime}_{1};b)a}(\epsilon^{\prime})$ is the bound state matrix element
with respect to $\epsilon^{\prime}$ for channel $c^{\prime}_{1}$, shown as
Fig. 7. It is necessary to emphasize that the energy in the two-particle
irreducible Green’s function is equal to the final state energy
$\epsilon^{\prime}$ while the mass $M_{0}$ and BS amplitude of initial bound
state is specified. We have introduced _extended Feynman diagram_ in Ref. Chen
and Lü (2023) to represent arbitrary value of the final state energy. In Fig.
7, the quark momenta in left-hand side of crosses depend on the final state
energy and the momenta in right-hand side depend on the initial state energy,
i.e., $p_{1}-p_{2}-p_{3}+p_{4}=Q+Q^{\prime}=P^{\epsilon^{\prime}}$ and
$p_{1}^{\prime}-p_{2}^{\prime}=P$. When $\epsilon^{\prime}=M_{0}$, the crosses
in Fig. 7 disappear and then Fig. 7 becomes Fig. 3;
$T_{(c^{\prime}_{1};b)a}(\epsilon^{\prime}=M_{0})$ is the $T$-matrix element
with mass $M_{0}$ for channel $c^{\prime}_{1}$ expressed as Eq. (52). Though
the $T$-matrix element $T_{(c^{\prime}_{1};b)a}(\epsilon^{\prime})$ has the
same form expressed as Eq. (52), these momenta should become

$\begin{split}&p_{1}=(Q+Q^{\prime})/2+p+k,~{}~{}p_{2}=(Q+Q^{\prime})/2-Q+p+k,~{}~{}p_{3}=k,~{}~{}p_{4}=Q^{\prime}+k,\\\
&q=Q^{\prime}/2+p+k,~{}k^{\prime}=Q^{\prime}(M_{0})+k,~{}~{}p_{1}^{\prime}=p+P/2,~{}~{}p_{2}^{\prime}=p-P/2,~{}~{}Q+Q^{\prime}=P^{\epsilon^{\prime}},\end{split}$
(81)

where $P=(0,0,0,iM_{0})$, $P^{\epsilon^{\prime}}=(0,0,0,i\epsilon^{\prime})$,
$Q^{\prime}(M_{0})=(-\mathbf{Q}(M_{0}),i\sqrt{\mathbf{Q}^{2}(M_{0})+M_{\omega}^{2}})$,
$Q=(\mathbf{Q}(\epsilon^{\prime}),i\sqrt{\mathbf{Q}^{2}(\epsilon^{\prime})+M_{J/\psi}^{2}})$,
$Q^{\prime}=(-\mathbf{Q}(\epsilon^{\prime}),i\sqrt{\mathbf{Q}^{2}(\epsilon^{\prime})+M_{\omega}^{2}})$
and $\mathbf{Q}^{2}(\epsilon^{\prime})=[\epsilon^{\prime
2}-(M_{J/\psi}+M_{\omega})^{2}][\epsilon^{\prime
2}-(M_{J/\psi}-M_{\omega})^{2}]/(4\epsilon^{\prime 2})$. The initial state is
considered as a four-quark state, so the specified GBS amplitude of initial
state should be

$\begin{split}\Gamma^{D_{l}^{*}}_{\lambda}(p_{1}^{\prime},k)\chi^{0^{+}}_{\lambda\tau}(P,p)\Gamma^{\bar{D}_{l}^{*}}_{\tau}(p_{2}^{\prime},k^{\prime}),\end{split}$
(82)

where $k^{\prime}$ depends on $P$. Then we obtain the function
$\mathbb{I}_{1}(\epsilon^{\prime})$ for channel $J/\psi\omega$

$\begin{split}\mathbb{I}_{1}(\epsilon^{\prime})&=\frac{1}{2}\int
d^{3}Qd^{3}Q^{\prime}(2\pi)^{4}\delta^{(4)}(Q+Q^{\prime}-P^{\epsilon^{\prime}})\sum_{\varrho^{\prime}=1}^{3}\sum_{\varrho=1}^{3}|T_{(c^{\prime}_{1};b)a}(\epsilon^{\prime})|^{2}.\end{split}$
(83)

Figure 7:  Matrix element with respect to $\epsilon^{\prime}$ for channel
$J/\psi\omega$. The momenta in the final state satisfy
$Q+Q^{\prime}=P^{\epsilon^{\prime}}$ and the momentum of the initial state is
$P$. The final state energy extends from $\epsilon_{M}$ to $+\infty$ while the
initial state energy is specified, and the crosses mean that the momenta of
quark propagators depend on the final state energy $\epsilon^{\prime}$.

###  V.2 Channel $D^{+}D^{-}$ with respect to arbitrary value of the final
state energy

The $T$-matrix element with respect to $\epsilon^{\prime}$ for channel
$c^{\prime}_{2}$ can be represented graphically by Fig. 8. The total energy
$\epsilon^{\prime}$ of the final state extends from
$\epsilon_{c^{\prime}_{2},M}$ to $+\infty$, i.e.,
$\epsilon_{c^{\prime}_{2},M}<\epsilon^{\prime}<\infty$ and
$\epsilon_{c^{\prime}_{2},M}=M_{D^{+}}+M_{D^{-}}$. In Fig. 8, the crosses mean
that the momenta of quark propagators and the momentum $w$ of the exchanged
light meson depend on $Q_{1}$ and $Q_{2}$, i.e.,
$p_{1}-p_{2}-p_{3}+p_{4}=p_{1}-p_{2}-q_{3}+q_{4}=Q_{1}+Q_{2}=P^{\epsilon^{\prime}}$
and $p_{1}^{\prime}-p_{2}^{\prime}=P$.

Figure 8:  Matrix element with respect to $\epsilon^{\prime}$ for channel
$D^{+}D^{-}$. The momenta in the final state satisfy
$Q_{1}+Q_{2}=P^{\epsilon^{\prime}}$ and the momentum of the initial state is
$P$. $w$ represents the momentum of the exchanged light meson. The crosses
mean that the momenta of quark propagators and the momentum $w$ of the
exchanged light meson depend on the final state energy $\epsilon^{\prime}$.

We still use the vertex function to calculate the $T$-matrix element with
respect to $\epsilon^{\prime}$ for channel $c^{\prime}_{2}$. However,
different from the ordinary vertex function, we should introduce the vertex
function with respect to $\epsilon^{\prime}$, which is shown as Fig. 9. In
Fig. 9, $Q_{1}$ depends on $P^{\epsilon^{\prime}}$, $p_{1}^{\prime}$ depends
on $P$ and the crosses mean that the momenta of quark propagators and the
momentum $w$ of the exchanged light meson depend on the final state energy
$\epsilon^{\prime}$. Using the approach introduced in Sec. III.2, we can
obtain the explicit forms for the vertex functions with respect to
$\epsilon^{\prime}$, and then Fig. 8 can be reduced to Fig. 10. In Fig. 10, we
have $Q_{1}+Q_{2}=P^{\epsilon^{\prime}}$, $p_{1}^{\prime}-p_{2}^{\prime}=P$
and the crosses lie on the right-hand side of light meson propagator, which
implies that the momentum $w$ of the exchanged light meson depends on $Q_{1}$
and $Q_{2}$.

Figure 9:  Vertex function for the exchanged light meson, heavy pseudoscalar
and vector mesons with respect to $\epsilon^{\prime}$. $Q_{1}$ depends on
$P^{\epsilon^{\prime}}$ and $p_{1}^{\prime}$ depends on $P$. The crosses mean
that the momenta of quark propagators and the momentum $w$ of the exchanged
light meson depend on the final state energy $\epsilon^{\prime}$. Figure 10:
Reduced matrix element with respect to $\epsilon^{\prime}$ for channel
$D^{+}D^{-}$. The crosses mean that the momentum $w$ of the exchanged light
meson depends on $Q_{1}$ and $Q_{2}$.

From Eq. (65), we obtain the total matrix element between the final state
$\langle D^{+}(Q_{1}),D^{-}(Q_{2})~{}\text{out}|$ and the mixed state of two
pure bound states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$

$\begin{split}-iR_{(c^{\prime}_{2};b)a}(\epsilon^{\prime})=\langle
Q_{1},Q_{2}~{}\text{out}|P~{}\text{in}\rangle=-i(2\pi)^{4}\delta^{(4)}(Q_{1}+Q_{2}-P^{\epsilon^{\prime}})T_{(c^{\prime}_{2};b)a}(\epsilon^{\prime}),\end{split}$
(84)

where $T_{(c^{\prime}_{2};b)a}(\epsilon^{\prime})$ is the bound state matrix
element with respect to $\epsilon^{\prime}$ for channel $c^{\prime}_{2}$,
shown as Fig. 10. When $\epsilon^{\prime}=M_{0}$, the crosses in Figs. 8, 9
and 10 disappear and then these three extended Feynman diagrams become Figs.
4, 6 and 5, respectively; $T_{(c^{\prime}_{2};b)a}(\epsilon^{\prime}=M_{0})$
is the $T$-matrix element with mass $M_{0}$ for channel $c^{\prime}_{2}$
expressed as Eq. (66). Though the $T$-matrix element
$T_{(c^{\prime}_{2};b)a}(\epsilon^{\prime})$ has the same form expressed as
Eq. (66), these momenta should become

$\begin{split}&w=(\mathbf{p}-\mathbf{Q}_{D}(\epsilon^{\prime}),0),~{}Q_{1}+Q_{2}=P^{\epsilon^{\prime}},~{}p_{1}^{\prime}=p+P/2,~{}p_{2}^{\prime}=p-P/2,\end{split}$
(85)

where $P=(0,0,0,iM_{0})$, $P^{\epsilon^{\prime}}=(0,0,0,i\epsilon^{\prime})$,
$Q_{1}=(\mathbf{Q}_{D}(\epsilon^{\prime}),i\epsilon^{\prime}/2)$,
$Q_{2}=(-\mathbf{Q}_{D}(\epsilon^{\prime}),i\epsilon^{\prime}/2)$ and
$\mathbf{Q}_{D}^{2}(\epsilon^{\prime})=[\epsilon^{\prime
2}-(M_{D^{+}}+M_{D^{-}})^{2}]/4$. The coefficients $E_{1}$ and $E_{2}$ in
interaction $\mathcal{V}_{\lambda\tau}(Q_{1},Q_{2},\textbf{p})$ given by Eq.
(62) should become
$E_{1}(\epsilon^{\prime})=E_{2}(\epsilon^{\prime})=\sqrt{\epsilon^{\prime}M_{0}}/2$.
Then we obtain the function $\mathbb{I}_{2}(\epsilon^{\prime})$ for channel
$D^{+}D^{-}$

$\begin{split}\mathbb{I}_{2}(\epsilon^{\prime})&=\frac{1}{2}\int
d^{3}Q_{1}d^{3}Q_{2}(2\pi)^{4}\delta^{(4)}(Q_{1}+Q_{2}-P^{\epsilon^{\prime}})|T_{(c^{\prime}_{2};b)a}(\epsilon^{\prime})|^{2}.\end{split}$
(86)

###  V.3 Channel $D^{0}\bar{D}^{0}$ with respect to arbitrary value of the
final state energy

In Figs. 8, 9 and 10, $PM$ and $\overline{PM^{\prime}}$ represent pseudoscalar
mesons $D^{0}$ and $\bar{D}^{0}$, respectively. Following the same procedure
as for charged channel $D^{+}D^{-}$, we can obtain the $T$-matrix element
$T_{(c^{\prime}_{3};b)a}(\epsilon^{\prime})$ with respect to
$\epsilon^{\prime}$ and the function $\mathbb{I}_{3}(\epsilon^{\prime})$ for
neutral channel $c^{\prime}_{3}$. Here, the total energy $\epsilon^{\prime}$
of the final state extends from $\epsilon_{c^{\prime}_{3},M}$ to $+\infty$,
i.e., $\epsilon_{c^{\prime}_{3},M}<\epsilon^{\prime}<\infty$ and
$\epsilon_{c^{\prime}_{3},M}=M_{D^{0}}+M_{\bar{D}^{0}}$.

###  V.4 Closed channel $D^{*}\bar{D}^{*}$

The final state $\langle D^{*},\bar{D}^{*}~{}\text{out}|$ can be written as

$\begin{split}\langle
D^{*},\bar{D}^{*}~{}\text{out}|=\frac{1}{\sqrt{2}}\langle
D^{*0},\bar{D}^{*0}~{}\text{out}|+\frac{1}{\sqrt{2}}\langle
D^{*+},D^{*-}~{}\text{out}|.\end{split}$ (87)

The total energy $\epsilon^{\prime}$ of the final state extends from
$\epsilon_{c^{\prime}_{4},M}$ to $+\infty$, i.e.,
$\epsilon_{c^{\prime}_{4},M}<\epsilon^{\prime}<\infty$ and
$\epsilon_{c^{\prime}_{4},M}=M_{D^{*}_{l}}+M_{\bar{D}^{*}_{l}}$. Considering
the lowest order term of the two-particle irreducible Green’s function, we can
obtain the $T$-matrix element between the final state $\langle
D^{*}_{l},\bar{D}^{*}_{l}~{}\text{out}|$ and the initial four-quark state,
which can be represented graphically by Fig. 11. In Fig. 11, $VM$ and
$\overline{VM^{\prime}}$ still represent $D_{l}^{*}$ and $\bar{D}_{l}^{*}$,
respectively; $Q_{1}$ and $Q_{2}$ still represent the momenta of final
particles, but $Q_{1}^{2}=-M_{D_{l}^{*}}^{2}$ and
$Q_{2}^{2}=-M_{\bar{D}_{l}^{*}}^{2}$; the crosses mean that the momenta of
quark propagators and the momentum $w$ of the exchanged light meson depend on
$Q_{1}$ and $Q_{2}$, i.e.,
$p_{1}-p_{2}-p_{3}+p_{4}=p_{1}-p_{2}-q_{3}+q_{4}=Q_{1}+Q_{2}=P^{\epsilon^{\prime}}$,
and $p_{1}^{\prime}-p_{2}^{\prime}=P$.

Figure 11:  Matrix element for closed channel $D_{l}^{*}\bar{D}_{l}^{*}$. The
momenta in the final state satisfy $Q_{1}+Q_{2}=P^{\epsilon^{\prime}}$ and the
momentum of the initial state is $P$. $w$ represents the momentum of the
exchanged light meson. The crosses mean that the momenta of quark propagators
and the momentum $w$ of the exchanged light meson depend on the final state
energy $\epsilon^{\prime}$.

To calculate the $T$-matrix element with respect to $\epsilon^{\prime}$ for
channel $c^{\prime}_{4}$, we also introduce the form factor of heavy meson
with respect to $\epsilon^{\prime}$, which is shown as Fig. 12. Using the
approach introduced in Sec. II.3.2, we can obtain the explicit forms for the
heavy meson form factors $h(w^{2})$ with respect to $\epsilon^{\prime}$, and
then Fig. 11 can be reduced to Fig. 13. In Fig. 13, we have
$Q_{1}+Q_{2}=P^{\epsilon^{\prime}}$, $p_{1}^{\prime}-p_{2}^{\prime}=P$, and
the crosses lie on the right-hand side of light meson propagator, which
implies that the momentum $w$ of the exchanged light meson depends on $Q_{1}$
and $Q_{2}$.

Figure 12:  Heavy meson form factor with respect to $\epsilon^{\prime}$.
$Q_{1}$ depends on $P^{\epsilon^{\prime}}$ and $p_{1}^{\prime}$ depends on
$P$. The crosses mean that the momenta of quark propagators and the momentum
$w$ of the exchanged light meson depend on the final state energy
$\epsilon^{\prime}$. Figure 13:  Reduced matrix element for closed channel
$D_{l}^{*}\bar{D}_{l}^{*}$. The crosses mean that the momentum $w$ of the
exchanged light meson depends on $Q_{1}$ and $Q_{2}$.

Using the Heisenberg picture, we obtain the total matrix element between the
final state $\langle D^{*}_{l}(Q_{1}),\bar{D}^{*}_{l}(Q_{2})~{}\text{out}|$
and the mixed state of two pure bound states $D^{*0}\bar{D}^{*0}$ and
$D^{*+}D^{*-}$

$\begin{split}-iR_{(c^{\prime}_{4};b)a}(\epsilon^{\prime})=\langle
Q_{1},Q_{2}~{}\text{out}|P~{}\text{in}\rangle=-i(2\pi)^{4}\delta^{(4)}(Q_{1}+Q_{2}-P^{\epsilon^{\prime}})T_{(c^{\prime}_{4};b)a}(\epsilon^{\prime}).\end{split}$
(88)

According to Mandelstam’s approach, the $T$-matrix element becomes

$\begin{split}T_{(c^{\prime}_{4};b)a}(\epsilon^{\prime})=\frac{1}{2}\sum_{l=u,d}\frac{-i\varepsilon_{\mu}^{\varrho^{\prime}}(Q_{2})\varepsilon_{\nu}^{\varrho}(Q_{1})}{(2\pi)^{9/2}\sqrt{2E_{D^{*}_{l}}(Q_{1})}\sqrt{2E_{\bar{D}^{*}_{l}}(Q_{2})}\sqrt{2E(P)}}(\mathcal{M}^{c^{\prime}_{4}}_{\nu\mu}+\mathcal{M}^{\prime
c^{\prime}_{4}}_{\nu\mu}),\end{split}$ (89)

where $\varepsilon_{\nu}^{\varrho=1,2,3}(Q_{1})$ and
$\varepsilon_{\mu}^{\varrho^{\prime}=1,2,3}(Q_{2})$ become the polarization
vectors of $D^{*}_{l}$ and $\bar{D}^{*}_{l}$, respectively, and

$\begin{split}\mathcal{M}^{c^{\prime}_{4}}_{\nu\mu}=&\int\frac{d^{4}p}{(2\pi)^{4}}\mathcal{V}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})\chi_{\lambda\tau}^{0^{+}}(P,p),\end{split}$
(90a) $\begin{split}\mathcal{M}^{\prime
c^{\prime}_{4}}_{\nu\mu}=&\int\frac{d^{4}p}{(2\pi)^{4}}\mathcal{V}^{\prime}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})\chi_{\lambda\tau}^{0^{+}}(P,p).\end{split}$
(90b)

Here $\chi_{\lambda\tau}^{0^{+}}(P,p)$ is expressed as Eq. (36),
$\mathcal{V}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})$ and
$\mathcal{V}^{\prime}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})$ represent
the interactions derived from one light meson ($\sigma$, $\rho^{0}$, $V_{1}$
and $V_{8}$) exchange and one-$\rho^{\pm}$ exchange, respectively.

Now, we determine the interactions
$\mathcal{V}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})$ and
$\mathcal{V}^{\prime}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})$.
Structurally similar to Eq. (16), we obtain the vertices of the heavy vector
meson interaction with light meson

$\begin{split}\langle
VM^{\varrho}(Q_{1})|J|VM^{\vartheta}(p_{1}^{\prime})\rangle=&\frac{1}{2\sqrt{E_{D^{*}_{l}}(Q_{1})E_{D^{*}_{l}}(p_{1}^{\prime})}}[\varepsilon^{\varrho}(Q_{1})\cdot\varepsilon^{\vartheta}(p_{1}^{\prime})]h^{(\text{s})}_{1}(w^{2}),\end{split}$
(91a)
$\begin{split}\langle\overline{VM^{\prime}}^{\varrho^{\prime}}(Q_{2})|J|\overline{VM^{\prime}}^{\vartheta^{\prime}}(-p_{2}^{\prime})\rangle=&\frac{1}{2\sqrt{E_{\bar{D}^{*}_{l}}(Q_{2})E_{\bar{D}^{*}_{l}}(-p_{2}^{\prime})}}[\varepsilon^{\varrho^{\prime}}(Q_{2})\cdot\varepsilon^{\vartheta^{\prime}}(-p_{2}^{\prime})]\bar{h}^{(\text{s})}_{1}(w^{2}),\end{split}$
(91b) $\begin{split}\langle
VM^{\varrho}(Q_{1})|J_{\alpha}|VM^{\vartheta}(p_{1}^{\prime})\rangle=&\frac{1}{2\sqrt{E_{D^{*}_{l}}(Q_{1})E_{D^{*}_{l}}(p_{1}^{\prime})}}\\{[\varepsilon^{\varrho}(Q_{1})\cdot\varepsilon^{\vartheta}(p_{1}^{\prime})]h_{1}^{(\text{v})}(w^{2})(Q_{1}+p_{1}^{\prime})_{\alpha}\\\
&-h^{(\text{v})}_{2}(w^{2})\\{[\varepsilon^{\varrho}(Q_{1})\cdot
p_{1}^{\prime}]\varepsilon^{\vartheta}_{\alpha}(p_{1}^{\prime})+[\varepsilon^{\vartheta}(p_{1}^{\prime})\cdot
Q_{1}]\varepsilon^{\varrho}_{\alpha}(Q_{1})\\}\\},\end{split}$ (91c)
$\begin{split}\langle\overline{VM^{\prime}}&{}^{\varrho^{\prime}}(Q_{2})|J_{\beta}|\overline{VM^{\prime}}^{\vartheta^{\prime}}(-p_{2}^{\prime})\rangle\\\
=&\frac{1}{2\sqrt{E_{\bar{D}^{*}_{l}}(Q_{2})E_{\bar{D}^{*}_{l}}(-p_{2}^{\prime})}}\\{[\varepsilon^{\varrho^{\prime}}(Q_{2})\cdot\varepsilon^{\vartheta^{\prime}}(-p_{2}^{\prime})]\bar{h}_{1}^{(\text{v})}(w^{2})(Q_{2}-p_{2}^{\prime})_{\beta}\\\
&-\bar{h}^{(\text{v})}_{2}(w^{2})\\{[\varepsilon^{\varrho^{\prime}}(Q_{2})\cdot(-p_{2}^{\prime})]\varepsilon^{\vartheta^{\prime}}_{\beta}(-p_{2}^{\prime})+[\varepsilon^{\vartheta^{\prime}}(-p_{2}^{\prime})\cdot
Q_{2}]\varepsilon^{\varrho^{\prime}}_{\beta}(Q_{2})\\}\\},\end{split}$ (91d)

where $w=p-(Q_{1}-Q_{2})/2$ is the momentum of the light meson, $h(w^{2})$ and
$\bar{h}(w^{2})$ are the heavy meson form factors with respect to
$\epsilon^{\prime}$. Similarly, taking away the external lines including
normalizations and polarization vectors $\varepsilon^{\varrho}_{\nu}(Q_{1})$,
$\varepsilon^{\vartheta}_{\lambda}(p_{1}^{\prime})$,
$\varepsilon^{\varrho^{\prime}}_{\mu}(Q_{2})$,
$\varepsilon^{\vartheta^{\prime}}_{\tau}(-p_{2}^{\prime})$, we obtain the
interaction from one light meson ($\sigma$, $\rho^{0}$, $V_{1}$ and $V_{8}$)
exchange

$\begin{split}\mathcal{V}&{}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})\\\
=&-2E_{1}(\epsilon^{\prime})F_{1}(\mathbf{w}^{2})\frac{-ig_{\sigma}^{2}}{w^{2}+M_{\sigma}^{2}}2E_{2}(\epsilon^{\prime})F_{1}(\mathbf{w}^{2})\delta_{\nu\lambda}\delta_{\tau\mu}+F_{2}(\mathbf{w}^{2})\left(\frac{-ig_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}+\frac{-ig_{1}^{2}}{w^{2}+M_{\omega}^{2}}+\frac{-ig_{8}^{2}}{w^{2}+M_{\phi}^{2}}\right)\\\
&\times
F_{2}(\mathbf{w}^{2})\\{(Q_{1}+p_{1}^{\prime})\cdot(Q_{2}-p_{2}^{\prime})\delta_{\nu\lambda}\delta_{\tau\mu}-\delta_{\nu\lambda}[-(Q_{1}+p_{1}^{\prime})_{\tau}p_{2\mu}^{\prime}+Q_{2\tau}(Q_{1}+p_{1}^{\prime})_{\mu}]\\\
&-[p_{1\nu}^{\prime}(Q_{2}-p_{2}^{\prime})_{\lambda}+(Q_{2}-p_{2}^{\prime})_{\nu}Q_{1\lambda}]\delta_{\tau\mu}-p_{1\nu}^{\prime}\delta_{\lambda\tau}p_{2\mu}^{\prime}+p_{1\nu}^{\prime}\delta_{\lambda\mu}Q_{2\tau}-\delta_{\nu\tau}Q_{1\lambda}p_{2\mu}^{\prime}+\delta_{\nu\mu}Q_{1\lambda}Q_{2\tau}\\},\end{split}$
(92)

where
$E_{1}(\epsilon^{\prime})=E_{2}(\epsilon^{\prime})=\sqrt{\epsilon^{\prime}M_{0}}/2$,
and $w=(\mathbf{w},0)$. The interaction from one-$\rho^{\pm}$ exchange becomes

$\begin{split}\mathcal{V}&{}^{\prime}_{\nu\lambda,\tau\mu}(Q_{1},Q_{2},\mathbf{p})\\\
=&F_{2}(\mathbf{w}^{2})\frac{-i2g_{\rho}^{2}}{w^{2}+M_{\rho}^{2}}F_{2}(\mathbf{w}^{2})\\{(Q_{1}+p_{1}^{\prime})\cdot(Q_{2}-p_{2}^{\prime})\delta_{\nu\lambda}\delta_{\tau\mu}-\delta_{\nu\lambda}[-(Q_{1}+p_{1}^{\prime})_{\tau}p_{2\mu}^{\prime}+Q_{2\tau}(Q_{1}+p_{1}^{\prime})_{\mu}]\\\
&-[p_{1\nu}^{\prime}(Q_{2}-p_{2}^{\prime})_{\lambda}+(Q_{2}-p_{2}^{\prime})_{\nu}Q_{1\lambda}]\delta_{\tau\mu}-p_{1\nu}^{\prime}\delta_{\lambda\tau}p_{2\mu}^{\prime}+p_{1\nu}^{\prime}\delta_{\lambda\mu}Q_{2\tau}-\delta_{\nu\tau}Q_{1\lambda}p_{2\mu}^{\prime}+\delta_{\nu\mu}Q_{1\lambda}Q_{2\tau}\\}.\end{split}$
(93)

These momenta in Fig. 13 become

$\begin{split}w=(\mathbf{p}-\mathbf{Q}_{D^{*}}(\epsilon^{\prime}),0),~{}~{}Q_{1}+Q_{2}=P^{\epsilon^{\prime}},~{}~{}p_{1}^{\prime}=p+P/2,~{}~{}p_{2}^{\prime}=p-P/2,\end{split}$
(94)

where $P=(0,0,0,iM_{0})$, $P^{\epsilon^{\prime}}=(0,0,0,i\epsilon^{\prime})$,
$Q_{1}=(\mathbf{Q}_{D^{*}}(\epsilon^{\prime}),i\epsilon^{\prime}/2)$,
$Q_{2}=(-\mathbf{Q}_{D^{*}}(\epsilon^{\prime}),i\epsilon^{\prime}/2)$ and
$\mathbf{Q}_{D^{*}}^{2}(\epsilon^{\prime})=[\epsilon^{\prime
2}-(M_{D^{*}_{l}}+M_{\bar{D}^{*}_{l}})^{2}]/4$.

Substituting Eqs. (92) and (93) into (90), we obtain the explicit forms for
tensors $\mathcal{M}^{c^{\prime}_{4}}_{\nu\mu}$ and $\mathcal{M}^{\prime
c^{\prime}_{4}}_{\nu\mu}$. The $p$ integral is also computed in instantaneous
approximation. $\mathcal{M}^{c^{\prime}_{4}}_{\nu\mu}$ and
$\mathcal{M}^{\prime c^{\prime}_{4}}_{\nu\mu}$ only depend on $Q_{1}$ and
$Q_{2}$, which can be calculated by means of the method given in Sec. III.1.
Applying Eq. (78), we obtain the function $\mathbb{I}_{4}(\epsilon^{\prime})$
for closed channel $D^{*}\bar{D}^{*}$

$\begin{split}\mathbb{I}_{4}(\epsilon^{\prime})&=\frac{1}{2}\int
d^{3}Q_{1}d^{3}Q_{2}(2\pi)^{4}\delta^{(4)}(Q_{1}+Q_{2}-P^{\epsilon^{\prime}})\sum_{\varrho^{\prime}=1}^{3}\sum_{\varrho=1}^{3}|T_{(c^{\prime}_{4};b)a}(\epsilon^{\prime})|^{2}.\end{split}$
(95)

##  VI physical mass and width of resonance

For resonance $\chi_{c0}(3915)$, the dispersion relation (77) becomes

$\begin{split}\mathbb{D}(M_{0})=&-\frac{\mathcal{P}}{\pi}\int_{\epsilon_{c^{\prime}_{1},M}}^{\infty}\frac{\mathbb{I}_{1}(\epsilon^{\prime})}{\epsilon^{\prime}-M_{0}}d\epsilon^{\prime}-\frac{\mathcal{P}}{\pi}\int_{\epsilon_{c^{\prime}_{2},M}}^{\infty}\frac{\mathbb{I}_{2}(\epsilon^{\prime})}{\epsilon^{\prime}-M_{0}}d\epsilon^{\prime}-\frac{\mathcal{P}}{\pi}\int_{\epsilon_{c^{\prime}_{3},M}}^{\infty}\frac{\mathbb{I}_{3}(\epsilon^{\prime})}{\epsilon^{\prime}-M_{0}}d\epsilon^{\prime}\\\
&-\frac{1}{\pi}\int_{\epsilon_{c^{\prime}_{4},M}}^{\infty}\frac{\mathbb{I}_{4}(\epsilon^{\prime})}{\epsilon^{\prime}-M_{0}}d\epsilon^{\prime},\end{split}$
(96)

where $\epsilon_{c^{\prime}_{1},M}=M_{J/\psi}+M_{\omega}$,
$\epsilon_{c^{\prime}_{2},M}=M_{D^{+}}+M_{D^{-}}$,
$\epsilon_{c^{\prime}_{3},M}=M_{D^{0}}+M_{\bar{D}^{0}}$, and
$\epsilon_{c^{\prime}_{4},M}=M_{D_{l}^{*}}+M_{\bar{D}_{l}^{*}}$. From Eq.
(75), we obtain that the physical mass of resonance $\chi_{c0}(3915)$ is
$M=M_{0}+(2\pi)^{3}\mathbb{D}(M_{0})$. Replacing $M_{0}$ by $M$ in Eqs. (52)
and (56), we recalculate the matrix element $T_{(c^{\prime}_{1};b)a}(M)$ and
obtain the width $\Gamma_{1}$ for physical decay model
$\chi_{c0}(3915)\rightarrow J/\psi\omega$. Replacing $M_{0}$ by $M$ in Eqs.
(66) and (68), we recalculate the matrix element $T_{(c^{\prime}_{2};b)a}(M)$
and obtain the width $\Gamma_{2}$ for physical decay model
$\chi_{c0}(3915)\rightarrow D^{+}D^{-}$. For the isospin conservation, it is
easy to obtain the width $\Gamma_{3}$ for physical decay model
$\chi_{c0}(3915)\rightarrow D^{0}\bar{D}^{0}$.

##  VII numerical result

Considering the isospin conservation, we employ the constituent quark masses
$m_{u}=m_{d}=0.33$ GeV, the heavy quark mass $m_{c}=1.55$ GeV Ebert et al.
(2006) and the meson masses $M_{\sigma}=0.45$ GeV, $M_{\omega}=0.782$ GeV,
$M_{\rho^{0}}=M_{\rho^{\pm}}=0.775$ GeV, $M_{\phi}=1.019$ GeV,
$M_{D^{*0}}=M_{D^{*+}}=2.007$ GeV, $M_{D^{0}}=M_{D^{+}}=1.865$ GeV,
$M_{J/\psi}=3.097$ GeV Workman et al. (2022). Without an adjustable parameter,
we numerically solve the eigenvalue equation (27) and obtain the masses and
wave functions of pure bound states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$
with spin-parity quantum numbers $0^{+}$. Considering the cross terms between
these two pure bound states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$ and using
the coupled-channel approach, we obtain the mass $M_{0}$ of the mixed state
with $0^{+}$. Then $M_{0}$ and GBS wave function
$\chi^{D^{*}\bar{D}^{*},0^{+}}(P,p,k,k^{\prime})$ given in Eq. (44) are used
to evaluate the matrix elements $T_{(c^{\prime}_{1};b)a}(M_{0})$ and
$T_{(c^{\prime}_{2};b)a}(M_{0})$ with mass of meson-meson bound state, and the
decay widths $\Gamma_{1}(M_{0})$ and $\Gamma_{2}(M_{0})$ with mass of meson-
meson bound state should not be the width of physical resonance. From Eqs.
(52), (53) and (81), we calculate the $T$-matrix element
$T_{(c^{\prime}_{1};b)a}(\epsilon^{\prime})$ with respect to
$\epsilon^{\prime}$ for channel $J/\psi\omega$. From Eqs. (62), (63), (66),
(67) and (85), we calculate the $T$-matrix element
$T_{(c^{\prime}_{2};b)a}(\epsilon^{\prime})$ with respect to
$\epsilon^{\prime}$ for channel $D^{+}D^{-}$. From Eqs. (89), (90), (92), (93)
and (94), we calculate the $T$-matrix element
$T_{(c^{\prime}_{4};b)a}(\epsilon^{\prime})$ with respect to
$\epsilon^{\prime}$ for closed channel $D^{*}\bar{D}^{*}$. From Eqs. (83),
(86) and (95), we calculate the functions $\mathbb{I}_{1}(\epsilon^{\prime})$
over $\epsilon_{c^{\prime}_{1},M}<\epsilon^{\prime}<\infty$,
$\mathbb{I}_{2}(\epsilon^{\prime})$ over
$\epsilon_{c^{\prime}_{2},M}<\epsilon^{\prime}<\infty$,
$\mathbb{I}_{3}(\epsilon^{\prime})$ over
$\epsilon_{c^{\prime}_{3},M}<\epsilon^{\prime}<\infty$, and
$\mathbb{I}_{4}(\epsilon^{\prime})$ over
$\epsilon_{c^{\prime}_{4},M}<\epsilon^{\prime}<\infty$, respectively. By doing
the numerical calculation, we obtain the mass correction $\Delta
M=(2\pi)^{3}\mathbb{D}(M_{0})$ and the physical mass $M$ for resonance
$\chi_{c0}(3915)$. Finally, the physical mass is used to recalculate these
strong decay widths $\Gamma_{1}(\chi_{c0}(3915)\rightarrow J/\psi\omega)$,
$\Gamma_{2}(\chi_{c0}(3915)\rightarrow D^{+}D^{-})$ and
$\Gamma_{3}(\chi_{c0}(3915)\rightarrow D^{0}\bar{D}^{0})$. Some errors in Ref.
Chen and Lü (2023) have been revised. $M$ and $\Gamma$ should be the observed
mass and full width in experiments. Our numerical results for resonance
$\chi_{c0}(3915)$ are in good agreement with the experimental data, which are
presented in Table 1.

Table 1: Mass $M$ and full width $\Gamma$ for physical resonance $\chi_{c0}(3915)$. $M_{0}$ is the mass of mixed state of two bound states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$, $\Delta M$ is the calculated correction due to all open and closed channels, and $\Gamma_{i}$ is the calculated width of $i$th decay channel. (Dimensioned quantities in MeV.) Quantity |  $M_{0}$ |  $\Delta M$ |  $M$ |  $\Gamma_{1}$ |  $\Gamma_{2}$ |  $\Gamma_{3}$ |  $\Gamma$
---|---|---|---|---|---|---|---
this work |  3953.7 |  $-$30.9 |  3922.8 |  22.5 |  1.5 |  1.5 |  25.5
PDGWorkman et al. (2022) |  |  |  3921.7$\pm$1.8 |  |  |  |  18.8$\pm$3.5

It is necessary to emphasize that there is not an adjustable parameter in our
approach. We require the meson-quark coupling constants $g$ and the parameters
$\omega_{H}$ in BS amplitudes of heavy mesons to calculate the mass and decay
width of physical resonance. The meson-quark coupling constants can be
determined by QCD sum rules approach Reinders et al. (1985), and these
parameters in BS amplitudes of heavy mesons are fixed by providing fits to
observables Ivanov et al. (1999, 2007); Kawanai and Sasaki (2011). Our
approach also involves the constituent quark masses $m_{u}$, $m_{d}$, and the
heavy quark mass $m_{c}$. According to the spontaneous breaking of chiral
symmetry, the light quarks ($u,d,s$) obtain their constituent masses because
the vacuum condensate is not equal to zero, and the heavy quark mass $m_{c}$
is irrelevant to vacuum condensate. Normally, the value slightly greater than
a third of nucleon mass is employed as the constituent mass of light quark.
The value of heavy quark mass $m_{c}$ can be determined by the experimental
mass of charmonium system $J/\psi$. Of course, the values of these parameters,
including $g$, $\omega_{H}$, $m_{u}$, $m_{d}$ and $m_{c}$, are values in
respective ranges. Simultaneously varying these parameters in respective
ranges, we find that the uncertainties of numerical results are at most 5%.
Despite the large uncertainty of meson mass $M_{\sigma}$, it has been found
that the uncertainties of numerical results from meson mass $M_{\sigma}$ are
also very small in our previous works Chen and Lü (2015); Chen et al. (2016);
Chen and Lü (2018); Chen et al. (2020) and Refs. Ding (2009); Zhao et al.
(2022). Therefore, in our approach the calculated mass and decay width are
uniquely determined.

Up to now, a theoretical approach from QCD to investigate resonance which is
regarded as an unstable two-body system has been established. In this paper,
we only explore exotic meson resonance which is considered as an unstable
molecular state composed of two heavy vector mesons. The extension of our
approach to more general resonances is straightforward, while the interaction
Lagrangian may be modified. More importantly, it is most reasonable and
fascinating to investigate resonance as far as possible from QCD. In the
framework of quantum field theory, the nonperturbative contribution from the
vacuum condensates can be introduced into BS wave function Chen et al. (2016)
and the two-particle irreducible Green’s function, and then the calculated
mass and decay width of resonance will contain more inspiration of QCD.

##  VIII conclusion

Exotic resonance $\chi_{c0}(3915)$ is considered as a mixed state of two
unstable molecular states $D^{*0}\bar{D}^{*0}$ and $D^{*+}D^{*-}$, and we
investigate the time evolution of meson-meson molecular state as determined by
the total Hamiltonian. According to the developed Bethe-Salpeter theory, the
total matrix elements for all decay channels should be calculated with respect
to arbitrary value of the final state energy. Because the total energy of the
final state extends from $\epsilon_{M}$ to $+\infty$, we consider three open
decay channels $J/\psi\omega$, $D^{+}D^{-}$, $D^{0}\bar{D}^{0}$ and one closed
channel $D^{*}\bar{D}^{*}$ from the effective interaction Lagrangian at low
energy QCD, which are exhibited by extended Feynman diagrams. Using the
developed Bethe-Salpeter theory, we calculate the mass $M$ and full width
$\Gamma$ of physical resonance $\chi_{c0}(3915)$, which are in good agreement
with the experimental data. Obviously, our work can be extended to more
general resonances.

###### Acknowledgements.

This work was supported by the National Natural Science Foundation of China
under Grants No. 11705104, 11801323 and 52174145; Shandong Provincial Natural
Science Foundation, China under Grants No. ZR2016AQ19 and ZR2016AM31; and
SDUST Research Fund under Grant No. 2018TDJH101.

##  Appendix A tensor structures in the general form of BS wave functions

The tensor structures in Eqs. (8) and (9) are given below Chen and Lü (2015);
Chen et al. (2020)

$\mathcal{T}_{\lambda\tau}^{1}=(p^{2}+\eta_{1}P\cdot p-\eta_{2}P\cdot
p-\eta_{1}\eta_{2}P^{2})g_{\lambda\tau}-(p_{\lambda}p_{\tau}+\eta_{1}P_{\tau}p_{\lambda}-\eta_{2}P_{\lambda}p_{\tau}-\eta_{1}\eta_{2}P_{\lambda}P_{\tau}),$
$\begin{split}\mathcal{T}_{\lambda\tau}^{2}=&(p^{2}+2\eta_{1}P\cdot
p+\eta_{1}^{2}P^{2})(p^{2}-2\eta_{2}P\cdot
p+\eta_{2}^{2}P^{2})g_{\lambda\tau}\\\ &+(p^{2}+\eta_{1}P\cdot
p-\eta_{2}P\cdot
p-\eta_{1}\eta_{2}P^{2})(p_{\lambda}p_{\tau}+\eta_{1}P_{\lambda}p_{\tau}-\eta_{2}P_{\tau}p_{\lambda}-\eta_{1}\eta_{2}P_{\lambda}P_{\tau})\\\
&-(p^{2}-2\eta_{2}P\cdot
p+\eta_{2}^{2}P^{2})(p_{\lambda}p_{\tau}+\eta_{1}P_{\lambda}p_{\tau}+\eta_{1}P_{\tau}p_{\lambda}+\eta_{1}^{2}P_{\lambda}P_{\tau})\\\
&-(p^{2}+2\eta_{1}P\cdot
p+\eta_{1}^{2}P^{2})(p_{\lambda}p_{\tau}-\eta_{2}P_{\lambda}p_{\tau}-\eta_{2}P_{\tau}p_{\lambda}+\eta_{2}^{2}P_{\lambda}P_{\tau}),\end{split}$
$\begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{3}=&\frac{1}{j!}p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\lambda}(p^{2}+2\eta_{1}P\cdot
p+\eta_{1}^{2}P^{2})[(p^{2}-2\eta_{2}P\cdot
p+\eta_{2}^{2}P^{2})(p+\eta_{1}P)_{\tau}\\\ &-(p^{2}+\eta_{1}P\cdot
p-\eta_{2}P\cdot p-\eta_{1}\eta_{2}P^{2})(p-\eta_{2}P)_{\tau}]\\\
&-p_{\mu_{1}}\cdots p_{\mu_{j}}[(p^{2}-2\eta_{2}P\cdot
p+\eta_{2}^{2}P^{2})(p_{\lambda}p_{\tau}+\eta_{1}P_{\lambda}p_{\tau}+\eta_{1}P_{\tau}p_{\lambda}+\eta_{1}^{2}P_{\lambda}P_{\tau})\\\
&-(p^{2}+\eta_{1}P\cdot p-\eta_{2}P\cdot
p-\eta_{1}\eta_{2}P^{2})(p_{\lambda}p_{\tau}+\eta_{1}P_{\lambda}p_{\tau}-\eta_{2}P_{\tau}p_{\lambda}-\eta_{1}\eta_{2}P_{\lambda}P_{\tau})],\end{split}$
$\begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{4}=&\frac{1}{j!}p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\tau}(p^{2}-2\eta_{2}P\cdot
p+\eta_{2}^{2}P^{2})[(p^{2}+\eta_{1}P\cdot p\\\ &-\eta_{2}P\cdot
p-\eta_{1}\eta_{2}P^{2})(p+\eta_{1}P)_{\lambda}-(p^{2}+2\eta_{1}P\cdot
p+\eta_{1}^{2}P^{2})(p-\eta_{2}P)_{\lambda}]\\\ &+p_{\mu_{1}}\cdots
p_{\mu_{j}}[(p^{2}+2\eta_{1}P\cdot
p+\eta_{1}^{2}P^{2})(p_{\lambda}p_{\tau}-\eta_{2}P_{\lambda}p_{\tau}-\eta_{2}P_{\tau}p_{\lambda}+\eta_{2}^{2}P_{\lambda}P_{\tau})\\\
&-(p^{2}+\eta_{1}P\cdot p-\eta_{2}P\cdot
p-\eta_{1}\eta_{2}P^{2})(p_{\lambda}p_{\tau}+\eta_{1}P_{\lambda}p_{\tau}-\eta_{2}P_{\tau}p_{\lambda}-\eta_{1}\eta_{2}P_{\lambda}P_{\tau})],\end{split}$
$\begin{split}\mathcal{T}^{5}_{\mu_{1}\cdots\mu_{j}\lambda\tau}=&\frac{1}{j!}(p^{2}+2\eta_{1}P\cdot
p+\eta_{1}^{2}P^{2})(p^{2}-2\eta_{2}P\cdot
p+\eta_{2}^{2}P^{2})p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\lambda}g_{\mu_{2}\\}\tau}\\\
&-\frac{1}{j!}p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\tau}(p^{2}-2\eta_{2}P\cdot
p+\eta_{2}^{2}P^{2})(p+\eta_{1}P)_{\lambda}\\\
&-\frac{1}{j!}p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}g_{\mu_{1}\\}\lambda}(p^{2}+2\eta_{1}P\cdot
p+\eta_{1}^{2}P^{2})(p-\eta_{2}P)_{\tau}\\\ &+p_{\mu_{1}}\cdots
p_{\mu_{j}}(p_{\lambda}p_{\tau}+\eta_{1}P_{\lambda}p_{\tau}-\eta_{2}P_{\tau}p_{\lambda}-\eta_{1}\eta_{2}P_{\lambda}P_{\tau}),\end{split}$
$\begin{split}\mathcal{T}^{6}_{\mu_{1}\cdots\mu_{j}\lambda\tau}=&p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\lambda\xi\zeta}p_{\xi}P_{\zeta}\epsilon_{\mu_{2}\\}\tau\xi^{\prime}\zeta^{\prime}}p_{\xi^{\prime}}P_{\zeta^{\prime}},\end{split}$
$\begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{7}=&-(2p^{2}+\eta_{1}P\cdot
p-\eta_{2}P\cdot p)p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}p_{\xi}\\\
&+(2\eta_{1}\eta_{2}P\cdot p+\eta_{2}p^{2}-\eta_{1}p^{2})p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}P_{\xi}\\\
&+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}p_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}p_{\lambda},\end{split}$
$\ \begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{8}=&-(P\cdot
p)p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}p_{\xi}+p^{2}p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}P_{\xi}\\\
&-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}p_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}p_{\lambda},\end{split}$
$\begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{9}=&-(2P\cdot
p+\eta_{1}P^{2}-\eta_{2}P^{2})p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}p_{\xi}\\\
&+P\cdot(\eta_{2}p-\eta_{1}P+2\eta_{1}\eta_{2}P)p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}P_{\xi}\\\
&+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}P_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}P_{\lambda},\end{split}$
$\begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{10}=&-P^{2}p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}p_{\xi}+(P\cdot
p)p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\tau\xi}P_{\xi}\\\
&-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}P_{\tau}+p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}P_{\lambda},\end{split}$
$\begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{11}=&(p^{2}+\eta_{1}P\cdot
p-\eta_{2}P\cdot p-\eta_{1}\eta_{2}P^{2})p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\lambda}\epsilon_{\mu_{2}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}\\\
&-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}\tau\xi\zeta}p_{\xi}P_{\zeta}(p-\eta_{2}P)_{\lambda},\end{split}$
$\begin{split}\mathcal{T}_{\mu_{1}\cdots\mu_{j}\lambda\tau}^{12}=&(p^{2}+\eta_{1}P\cdot
p-\eta_{2}P\cdot p-\eta_{1}\eta_{2}P^{2})p_{\\{\mu_{3}}\cdots
p_{\mu_{j}}g_{\mu_{1}\tau}\epsilon_{\mu_{2}\\}\lambda\xi\zeta}p_{\xi}P_{\zeta}\\\
&-p_{\\{\mu_{2}}\cdots
p_{\mu_{j}}\epsilon_{\mu_{1}\\}_{\lambda}\xi\zeta}p_{\xi}P_{\zeta}(p+\eta_{1}P)_{\tau}.\end{split}$

## References

  * Swanson (2004) E. S. Swanson, Phys. Lett. B 588, 189 (2004).
  * Törnqvist (2004) N. A. Törnqvist, Phys. Lett. B 590, 209 (2004).
  * Liu and Zhu (2009) X. Liu and S. L. Zhu, Phys. Rev. D 80, 017502 (2009).
  * Chen and Lü (2015) X. Chen and X. Lü, Eur. Phys. J. C 75, 98 (2015).
  * Zhao et al. (2022) M.-J. Zhao, Z.-Y. Wang, C. Wang, and X.-H. Guo, Phys. Rev. D 105, 096016 (2022).
  * Chen and Lü (2018) X. Chen and X. Lü, Phys. Rev. D 97, 114005 (2018).
  * Chen et al. (2020) X. Chen, X. Lü, R. Shi, X. Guo, and Q. Wang, Phys. Rev. D 101, 014009 (2020).
  * Chen and Lü (2023) X. Chen and X. Lü, Eur. Phys. J. C 83, 499 (2023).
  * Esposito et al. (2017) A. Esposito, A. Pilloni, and A. D. Polosa, Phys. Rep. 668, 1 (2017).
  * Chen et al. (2013) X. Chen, R. Liu, R. Shi, and X. Lü, Phys. Rev. D 87, 065013 (2013).
  * Chen et al. (2016) X. Chen, X. Lü, R. Shi, and X. Guo, Nucl. Phys. B 909, 243 (2016).
  * Ambrosino et al. (2009) F. Ambrosino et al. (KLOE Collaboration), JHEP 07, 105 (2009).
  * Reinders et al. (1985) L. Reinders, H. Rubinstein, and S. Yazaki, Phys. Rep. 127, 1 (1985).
  * Lü et al. (1996) X. Lü, Y. Liu, and E. Zhao, Chinese Physics Letters 13, 652 (1996).
  * Lü et al. (1997) X. Lü, Y. Liu, and E. Zhao, Science in China (Series A) 27, 361 (1997).
  * Burden et al. (1997) C. J. Burden, L. Qian, C. D. Roberts, P. C. Tandy, and M. J. Thomson, Phys. Rev. C 55, 2649 (1997).
  * Maris et al. (1998) P. Maris, C. D. Roberts, and P. C. Tandy, Phys. Lett. B 420, 267 (1998).
  * Ivanov et al. (1999) M. A. Ivanov, Y. L. Kalinovsky, and C. D. Roberts, Phys. Rev. D 60, 034018 (1999).
  * Ivanov et al. (2007) M. A. Ivanov, J. G. Körner, S. G. Kovalenko, and C. D. Roberts, Phys. Rev. D 76, 034018 (2007).
  * Neubert (1994) M. Neubert, Phys. Rep. 245, 259 (1994).
  * Chen et al. (2009) X. Chen, B. Wang, X. Li, X. Zeng, S. Yu, and X. Lü, Phys. Rev. D 79, 114006 (2009).
  * Choi et al. (2005) S.-K. Choi et al. (Belle Collaboration), Phys. Rev. Lett. 94, 182002 (2005).
  * Aubert et al. (2008) B. Aubert et al. (BABAR Collaboration), Phys. Rev. Lett. 101, 082001 (2008).
  * Uehara et al. (2010) S. Uehara et al. (Belle Collaboration), Phys. Rev. Lett. 104, 092001 (2010).
  * Lees et al. (2012) J. P. Lees et al. (BABAR Collaboration), Phys. Rev. D 86, 072002 (2012).
  * Vinokurova et al. (2015) A. Vinokurova et al. (Belle Collaboration), JHEP 06, 132 (2015).
  * Aaij et al. (2020) R. Aaij et al. (LHCb Collaboration), Phys. Rev. D 102, 112003 (2020).
  * Lurié (1968) D. Lurié, _Particles and Fields_ (Interscience Publishers, New York, 1968).
  * Goldberger and Watson (1964) M. L. Goldberger and K. M. Watson, _Collision Theory_ (Wiley, New York, 1964).
  * Ebert et al. (2006) D. Ebert, R. N. Faustov, and V. O. Galkin, Phys. Lett. B 634, 214 (2006).
  * Workman et al. (2022) R. L. Workman et al. (Particle Data Group), Prog. Theor. Exp. Phys. 2022, 083C01 (2022).
  * Kawanai and Sasaki (2011) T. Kawanai and S. Sasaki, Phys. Rev. Lett. 107, 091601 (2011).
  * Ding (2009) G.-J. Ding, Phys. Rev. D 79, 014001 (2009).