Probing axion-like particles in leptonic decays of heavy mesons

In 1977, Peccei and Quinn introduced a new global chiral symmetry [1], known as $ U(1)_{\rm PQ} $ symmetry, to address the strong CP problem in QCD, that is, the absence of CP violation in the strong interactions and the neutron's electric dipole moment (EDM) being forbidden. At a energy scale $ f_a $, such a symmetry is assumed to be spontaneously broken, which results in the appearance of a pseudo-Nambu-Goldstone boson, namely the axion [1−5], whose mass $ m_a $ is constrained by the relation $ m_a \sim \frac{m_\pi f_\pi}{f_a} $, where $ m_\pi $ and $ f_\pi $ are the mass and decay constant of pion, respectively. This means that if the energy scale $ f_a $ is very high, the axion should be extremely light. The condition can be relaxed if we are not just limited to such QCD axion, but consider more general Axion-Like Particles (ALPs) [6−8]. In such cases, both the PQ symmetry breaking scale $ f_a $ and the ALP mass $ m_a $ can be considered as independent parameters.

People are interested in ALPs in many aspects. Theoretically, a large number of ALPs are predicted by string theories [9, 10] and supersymmerty[11, 12]. These particles may play an important role in the evolution of cosmology, most importantly, as a candidate of dark matter [13]. The cosmological and astronomical observations can set strict constraints on the very light ALPs [14−20]. The phenomenology of ALPs have also been extensively studied at the Large Hadron Collider (LHC) [21−29]. For example, they can be produced by the decays of on-shell Higgs/Z boson [6, 29], or participate as off shell mediators in the scattering processes [30]. At high-luminosity electron-positron colliders, ALPs with mass in the MeV-GeV range can also be explored non-resonantly or resonantly [31−36], which means they can be produced directly by coupling to the charged leptons or the gauge bosons, or produced through the decays of final mesons. For example, in Ref. [35] the production of ALPs at B factories is considered, including the contributions of $ e^{+}e^{-} \to \gamma a $ and $ e^+e^-\to\Upsilon\to\gamma a $. Recently, the BESIII experiment [37] use a data sample of $ \psi(3686) $ to get the upper limits of the branching fraction of the decay $ J/\psi \to \gamma a $ and set constrains on the coupling $ g_{a\gamma\gamma} $ in the mass range of $ 0.165\; {\rm GeV} \leq m_a \leq 2.84\; {\rm GeV} $. The fixed target experiments are also promising methods in searching ALPs [38−44], as the detectors can extend tens of meters and suitable for the detection of long-lived particles. The ALP productions though K meson decays are usually investigated in such cases [38−41].

Except the methods mentioned above, the decays of heavy-flavored mesons also provide a way to probe ALPs. On the one hand, a large number of bottomed and charmed mesons have been produced at B factories and τ-charm factory, respectively, and on the other hand, a larger range of $ m_a $ can be looked into in the decays of heavy mesons compared with that of K meson. We will focus on the decay processes $ h^- \to \ell^- \not{E} $, where $ h^-=D^-,\; D_s^-, \; B^- $, or $ B_c^- $; $ \ell^-=e^-,\; \mu^- $, or $ \tau^- $; $ \not{E} $ represents the missing energy. In the Standard Model, $ \not{E} $ is carried away by the antineutrino (see Fig. 1), and the corresponding partial width is

Figure 1.  The diagram of two-body leptonic decays of charged heavy meson.

where $ G_F $ is the Fermi coupling constant, $ f_h $ is the decay constant of h meson, $ V_{Qq} $ is the CKM matrix element, M and $ m_\ell $ are the masses of the meson and the lepton, respectively, and the neutrino mass is assumed to be zero. If an ALP can also be produced and has long enough lifetime to escape the detector, we will have $ \Gamma(h^-\to \ell^-\not{E})=\Gamma(h^-\to \ell^-\bar{\nu}_{\ell})+\Gamma(h^-\to \ell^-\bar{\nu}_{\ell}a) $. Then the experimental results of $ h^-\to \ell^- \not{E} $, which are presented in Table 1, can be used to constrain the coupling between the ALP and SM particles.

Channel	Experiment values
$ {\mathcal B}(B^-\to e^- \not{E}) $	$ \le 9.8 \times 10^{-7} $ [45]
$ {\mathcal B}(B^-\to \mu^- \not{E}) $	$ (5.3 \pm 2.0 \pm 0.9) \times 10^{-7} $ [46]
$ {\mathcal B}(B^-\to \tau^- \not{E}) $	$ (7.2 \pm 2.7 \pm 1.1 ) \times 10^{-5} $ [47]
$ {\mathcal B}(D_{s}^-\to e^- \not{E}) $	$ \le 0.83 \times 10^{-4} $ [48]
$ {\mathcal B}(D_{s}^-\to \mu^- \not{E}) $	$ (0.5294 \pm 0.0108 \pm 0.0085) \times 10^{-2} $ [49]
$ {\mathcal B}(D_{s}^-\to \tau^- \not{E}) $	$ (5.44 \pm 0.17 \pm 0.13 ) \times 10^{-2} $ [50]
$ {\mathcal B}(D^-\to e^- \not{E}) $	$ \le 8.8 \times 10^{-6} $ [51]
$ {\mathcal B}(D^-\to \mu^- \not{E}) $	$ (3.71 \pm 2.7 \pm 1.1) \times 10^{-3} $ [52]
$ {\mathcal B}(D^-\to \tau^- \not{E}) $	$ (1.20 \pm 0.24 \pm 0.12 ) \times 10^{-5} $ [53]

Table 1.  Experiment results of $ {\mathcal B}(B^- \to \ell^-\not{E}) $, $ {\mathcal B}(D^-\to $$ \ell^-\not{E}) $, and $ {\mathcal B}(D_s^-\to \ell^-\not{E}) $.

What we will do in this article is inspired by Refs. [54−57] where the processes $ h^- \to \ell^- \bar{\nu}_{\ell}a $ with a being an ALP were studied. Beyond the considerations for the case of pseudoscalar, Ref. [54] also considered the possibility that a being a vector particle. Here we will also probe such processes, but with three main differences. First, the hadron transition matrix elements will be calculated by using the Bethe-Salpeter (BS) method. This method has been widely used to study the weak decays [58−60], strong decays [61] and electromagnetic decays [62] of heavy mesons, and the theoretical results are consistent with the experimental data. Second, for the coupling of ALP with quarks, we consider the interference effects of different diagrams, and extract the allowed range of parameters in a more general way. Third, we will explore the ALP production from the decays of double heavy meson $ B_c $, that is $ B_c^- \to \ell^- \bar{\nu}_{\ell}a $, which could give the constraint in a larger mass range of ALP.

The paper is organized as follows. In Section II, we write the effective Lagrangian which describes the interaction between ALPs and charged fermions. Then we calculate the hadronic transition matrix elements of the $ h^-\to \ell^-\nu_{\ell} a $ processes by using the BS method. In Section III, we apply the experimental results to extract the allowed regions of the parameters which are then used to calculate the upper limits of the branching fractions of the $ B_c $ decay processes. Finally, we present the summary in Section IV.

Generally, the ALPs can interact both with the fermions and gauge bosons. Here in the case of leptonic decays of heavy mesons, we only consider the tree level diagrams shown in Fig. 2 which represent the interaction of the ALP with SM fermions. The contributions of gauge bosons are of higher order and assumed to be neglected. The effective Lagrangian, which describes flavor conserving processes, consists of dimension-five operators [54−57]

Figure 2.  Feynman diagrams of $ h \rightarrow \ell\bar{\nu}_\ell a $. (a), (b), (c) and (d) represent the ALP couples with the quark, antiquark, charged lepton and antineutrino, respectively.

where c_i is the effective coupling constant with the index i extending over all the fermions. As Refs. [55−57] did, in the last step we have used the equation of motion to rewrite the expression. This corresponds to the use of another operator basis, which will lead to the different upper limits of c_i. However, when we calculate the decay channels of $ B_c $ meson, two methods give quite close results. The dependence of the fermion mass $ m_i $ indicates that the ALP does not couple directly with the neutrinos and Fig. 2(d) needs not to be taken into account, which will simplify the calculation.

The Feynman amplitude corresponding to each diagram can be written by using the meson wave function. Here, as we are considering the heavy flavored meson, which can be seen as a two-body bound state including at least one heavy quark (or antiquark), the instantaneous BS wave function is appropriate for the calculation. This kind of wave functions have been extensively used to study the decay processes of heavy mesons. For the pseudoscalar meson, its wave function has the form [63]

where P is the meson momentum, and q is the relative momentum between the quark and antiquark; $ q^\mu_\perp \equiv q^\mu- \frac{P \cdot q }{M^2} P^\mu = (0,\; \vec q) $; $ f_i $ is the radial wave function which can be achieved by solving the eigenvalue equations numerically (e.g. the detailed results for $ B_c $ meson can be found in Ref. [64]).

The amplitude corresponding to Fig. 2(a) can be written as

where $ c_Q $ is coupling constant; $ m_Q $ is the mass of quark Q; $ N_c=3 $ is the color factor; $ k_a $, $ k_\ell $, and $ k_{\bar\nu} $ are the momenta of ALP, charged lepton, and antineutrino, respectively; $ p_Q $ is the momentum of quark Q, which is related to P and q as $ p_Q = \frac{m_Q}{m_Q+m_q} P+q $ with $ m_q $ being the mass of the antiquark. To get Eq. (4), we have made an approximation, namely $ p_Q \approx \frac{m_Q}{m_Q+m_q} P+q_\perp $, so that the denominator is independent of the time component of the relative momentum q. After integrating out $ q_\perp $, $ \mathcal{M}_a $ becomes

where we have used $ A_Q = \frac{c_Q m_Q G_F V_{qQ}}{f_a \sqrt{2}} $ for short; $ F_1 $ and $ F_2 $ are the form factors of the hadronic transition matrix element, which are functions of the squared momentum transition $ (P-k_a)^2 $.

Similarly, the Feynman amplitude corresponding to Fig. 2(b) can be written as

where $ c_q $ is the coupling constant; $ m_q $ and $ p_q $ are the mass and momentum of the antiquark, respectively; $ p_q $ has the form $ p_q\approx \frac{m_q}{m_Q+m_q} P-q_\perp $. After integrating out $ q_\perp $, we get

where $ A_q=\frac{c_q m_q G_FV_{qQ}}{f_a \sqrt{2}} $; $ G_1 $ and $ G_2 $ are form factors.

When calculating the integral, one must be careful, since as the relative momentum $ \vec q $ changes, the denominator of the propagator may be zero, which will lead to a nonvanishing imaginary part of the form factors [65]. Specifically, we can write the propagator as

where θ is the angle between $ \vec q $ and $ \vec k_a $. a and b are expressed as

where $ E_a $ is the energy of the ALP. We will first integrate out $ \cos\theta $ analytically, and then integrate out $ |\vec q| $ numerically. As $ \vec q $ changes, $ a+b $ and $ a-b $ may have opposite signs, which means with a specific value of θ the pole can exist.

The Feynman amplitude corresponding to Fig. 2(c) can be written as

where $ c_\ell $ is the coupling constant. The hadronic part of the amplitude is proportional to the meson momentum [66], with the proportionality being the decay constant $ f_h $

Using Eq. (10) and defining $ A_\ell = \frac{c_\ell m_\ell G_FV_{qQ}}{f_a \sqrt{2}} $ for short, $ \mathcal{M}_c $ is simplified as

The total Feynman amplitude is $ \mathcal M=\mathcal M_a +\mathcal M_b +\mathcal M_c $. And the partial width of such a decay channel is achieved by finishing the three-body phase space integral

where the upper and lower limits have the following forms

In this section, we will first use the formulas obtained above to calculate the limits of the coupling constants, and then the results will be applied to investigate the similar decays of the $ B_c $ meson. We will focus on three different scenarios: first, the ALP couples only with a single fermion, namely the quark, the antiquark, or the charged lepton; second, the ALP couples with quarks and antiquarks, but not with leptons; third, the ALP couples with all the charged fermions with the same coupling constant.

During the calculation process, the following numerical values of three groups of physical quantities are adopted. (1) The masses of constituent quarks are from Ref. [67]: $ m_u = 0.305 $ GeV, $ m_d = 0.311 $ GeV, $ m_s = 0.5 $ GeV, $ m_c = 1.62 $ GeV, and $ m_b = 4.96 $ GeV. (2) The relevant CKM matrix elements are from Particle Data Group (PDG) [68]: $ |V_{ub}| = 3.82 \times 10^{-3} $, $ |V_{cd}| = 0.211 $, $ |V_{cs}| = 0.975 $, and $ |V_{cb}| = 41.1 \times 10^{-3} $. (3) The decay constant of D, $ D_s $, and B mesons are the Lattice QCD results [69]: $ f_D = 212.6 $ MeV, $ f_{D_s}=249.9 $ MeV, and $ f_B = 188 $ MeV.

We consider separately the coupling between the ALP and each charged fermion. For example, by setting $ c_q $ and $ c_\ell $ to zero, we get the branching fraction of $ h^-\to\ell^-\bar\nu a $ with $ c_Q/f_a $ and $ m_a $ being unknown parameters. Then using the experimental results in Table 1, we get the upper limits of $ c_Q/f_a $ as functions of $ m_a $. The same is true for the other two cases. The results are shown in Fig. 3. One notices that the curves raise rapidly when $ m_a $ is large enough, which is mainly because the decreasing of phase space. The most stringent upper limits for different coupling constants come from different decay channels. For example, Fig. 3(a) shows that the most strict constraints of $ c_{b}/f_a $ and $ c_{c}/f_a $, which are of the order of $ 10^2\; {\rm TeV}^{-1} $ when $ m_a $ is less than $ 1 $ GeV, come from $ B^-\to\mu^-\bar\nu_\mu a $ and $ D_s^+\to e^+\nu_e a $, respectively. One can also notice that the upper limit of $ c_e/f_a $ in Fig. 3(c) is much higher than those of $ c_{Q(q)}/f_a $ set up by the same channels, which is because of the small $ m_e $.

Figure 3.  (color online) The upper limits of the coupling constants (a) $ c_Q / f_a $, (b) $ c_q / f_a $ and (c) $ c_{\ell} / f_a $ derived from the leptonic decays of D, $ D_s $, and B mesons.

Using above results, we can impose restrictions on $ \mathcal B(B_c^-\to \ell^- \bar{\nu}_\ell a) $, which are presented in Fig. 4. We can see that when $ m_a $ is less than 1 GeV, the branching ratios upper limits for ALP-quark coupling cases are of the order of $ 10^{-4} $, while for the ALP-tau coupling case, the order of magnitudes is $ 10^{-2} $, which is much larger than those of the ALP-$ e/\mu $ coupling cases. We also notice that all the curves, except the one corresponding to $ B_c^-\to\tau^-\bar\nu_\tau a $ in Fig. 4(a), show an ascending trend. This is because the maximum value of $ m_a $ permitted in the $ B_c $ decay is larger than that in the same decays of B and $ D_{(s)} $ mesons. For example, when we use $ c_{b}/f_a $ from $ B^-\to\mu^-\bar\nu_\mu a $ to constrain the branching ratios of $ B_c^-\to\mu^-\bar\nu_\mu a $, the result blows up before $ m_a $ reaches its maximum. The appearance of kinks in Fig. 4(c) is because we have used $ c_{e/\mu}/f_a $ from different channels when $ m_a $ takes different values.

Figure 4.  (color online) The upper limits of $ \mathcal{B}(B_c^- \to \ell^- \bar{\nu}_\ell a) $ in scenario 1.

In this scenario, we consider the case that the ALP couples with quarks but not with leptons. Two assumptions are made. First, all the coupling constants are assumed to be real numbers, and no additional relative phase is introduced. Second, the coupling constants of the ALP and light quarks, denoted by $ c_q $, are assumed to be equal, while for those of heavy quarks, $ c_Q $ ($ Q=c,\; b $) may have different values. We let $ c_q $ be a free parameter, then use the experimental results of B and $ D_{(s)} $ mesons to provide constraints on $ c_b $ and $ c_c $. Specifically, we first write the following expression,

where the three terms on the right side of the last equation represent the contributions of Fig. 2(a), 2(b) and their interference, respectively; $ \tilde\Gamma $ denotes the decay width devided by the coupling constant sqaured, which can be calculated theoretically.

The experimental allowed region for $ c_Q $ and $ c_q $ is determined by Eq. (13). We investigate all the decay channels, and find that the decays of $ B^- $ and $ D^- $ mesons give the most stringent constraint. The area of the parameter space depends on $ m_a $. In Fig. 5, as two examples, we show the results when $ m_a $=0, and 1.6 GeV, respectively. One can see that when $ m=0 $ GeV, the two areas are comparable, while when $ m=1.6 $ GeV, the area for $ Q=c $ (just partially plotted) is more larger than that of $ Q=b $.

Figure 5.  (color online) The experimental allowed region (inside the ellipse) of the coupling constants. The green and blue areas correspond to $ Q=c $ and $ Q=b $, respectively.

Next, we scan the experimental allowed parameter area in Fig. 5, and calculate the limits of $ \mathcal{B}(B_c^- \to \ell^- \bar{\nu}_\ell a) $. Because both $ c_Q $ and $ c_q $ changing the sign does not affect the results, we only need to consider the region with $ c_q\ge 0 $. Our strategy is as follows. First, by solving Eq. (13), we obtain the boundary values of $ c_Q $ for a fixed $ c_q $,

where we have defined $ C_{Q(q)}=c_{Q(q)}/f_a $. The condition $ C_q^2 \le \frac{4 \tilde\Gamma_Q \delta\Gamma}{4 \tilde\Gamma_Q \tilde\Gamma_{q}-\tilde\Gamma_{Qq}^2} $ should be satisfied to make sure the quantity under the square root nonnegative. Second, we scan $ C_b $ and $ C_c $ together for a fixed $ C_q $ whose allowed region is determined by the small ellipse in Fig. 5.

Finally, we calculate the upper limits of $ \mathcal{B}(B_c^- \to \ell^- \bar{\nu}_\ell a) $ as functions of $ C_q $ with a specific $ m_a $. The results are shown in Fig. 6(a) and (c). Here we only present the result for the positive $ C_q $, and that of the negative $ C_q $ is symmetrical with it about the vertical axis. It can be seen that as $ m_a $ increases, the experimental allowed region of $ C_q $ shrinks, and the upper limit of $ \mathcal{B}(B_c^- \to \ell^- \bar{\nu}_\ell a) $ roughly increases with $ m_a $. This can be understood from Fig. (5), where as $ m_a $ changes from 0 GeV to 1.6 GeV, the allowed range of $ c_c/f_a $ becomes quite large, which leads to a larger upper limit. The interesting thing is that for some range of $ c_q $, there is also a nonvanishing lower limit of the branching fraction (see Fig. 6(b) and (d)). From Fig. 5 we can see the two ellipses are tilted, so at some $ c_q $, $ c_b $ and $ c_c $ cannot be zero at the same time, which leads to the nonvanishing $ \mathcal{B}(B_c^- \to \ell^- \bar{\nu}_\ell a) $.

Figure 6.  (color online) The upper and lower limits of $ \mathcal B(B_c \rightarrow \ell \bar{\nu}_\ell a) $ in scenario 2.

In this scenario, we assume that the ALP couples with all the charged fermions with the same coupling constant $ c_q=c_Q=c_l=c $, the upper limits of which as functions of $ m_a $ are shown in Fig. 7(a). We can see the $ B^-\to\mu\bar\nu_\mu a $ channel gives the most stringent restriction. We also notice that the curves have undulation around 1 GeV. To see why this happened, we take the $ B^-\to\mu\bar\nu_\mu a $ channel as an example, and plot the contributions of three Feynman diagrams in Fig. 2 and the corresponding interference terms. It can be seen that the effect of the lepton can be neglected because of its small mass, while the coupling between the ALP and quarks provides the main contribution, especially the interference term, which has sizable negative value when $ m_a<2 $ GeV. This causes the total result to be undulant around $ m_a=1 $ GeV, which is also transferred to the coupling constant.

Figure 7.  (color online) (a) shows the upper limits of the coupling constants in scenario 3. (b) shows the contributions to $ \mathcal{B}(B^- \to e^- \bar{\nu}_e a) f_a^2 /c_{i}^{2} $ by each Feynman diagram and the interference terms.

Using the upper limit of the coupling constant given by the $ B^-\to\mu\bar\nu_\mu a $ channel, we obtain the constraints of the branching fraction of $ B_c^- \to \ell^- \bar{\nu}_{\ell} a $, which are shown in Fig. 8(a). We can see, for $ \mathcal{B}(B_c^- \to e^- \bar{\nu}_e a) $ and $ \mathcal{B}(B_c^- \to \mu^- \bar{\nu}_{\mu} a) $, the upper limit is about $ 10^{-6} $ when $ m_a $ is less than 1 GeV, and then it keeps increasing until being larger than one. The reason for this is that the phase space of the B decay is small compared with that of the $ B_c^- $ decay when $ \ell=e $ or μ. For $ \mathcal{B}(B_c^- \to \tau^- \bar{\nu}_e a) $, the upper limit is around $ 10^{-5} $ when $ m_a<1 $ GeV, and then there is a fluctuation with the peak value of $ 5.4\times 10^{-4} $ at $ m_a=3.4 $ GeV, and finally the branching fraction goes to zero because of the vanishing phase space. It is worth to compare these results with the branching fractions of $ B_c^-\to\ell^-\bar\nu_\ell $ in SM. Using Eq. (1) and (10) we get $ \mathcal B(B_c^-\to \ell\bar\nu_\ell)= 2.2\times 10^{-9},\; 9.5\times 10^{-5} $, and $ 2.2\times 10^{-2} $ for $ \ell=e,\; \mu $, and τ, respectively. One can see that $ \mathcal B(B_c^-\to e^-\bar\nu_e) $ is much smaller than the upper limit of $ \mathcal B(B_c^-\to e^-\bar\nu_e a) $ for the chiral suppression; on the contrary, $ \mathcal B(B_c^-\to \tau^-\bar\nu_\tau) $ is much larger than the upper limit of $ \mathcal B(B_c^-\to \tau^-\bar\nu_\tau a) $.

Figure 8.  (color online) (a) shows the upper limit of $ \mathcal{B}(B_c^- \to \ell^- \bar{\nu}_{\ell} a) $ in scenario 3. (b), (c), and (d) shows the upper limits of the differential branching fraction with $ \ell=e $, μ, and τ, respectively.

In Fig. 8(b), (c), and (d), we present the upper limits of the differential branching fractions as functions of the lepton energy. For $ \ell $ being e or μ, the results are similar because of the small lepton mass, while for $ \ell=\tau $, the distribution shape is quite different. We can see that as $ m_a $ gets larger, the peak moves left. Experimentally, the detection of the unmonoenergetic charged lepton in the decay $ h^-\to\ell^-\not{E} $ may indicate the existence of ALP. Around the peak, it has the largest experimental allowed probability to find the charged lepton.

To conclude, we have studied the ALP production through the processes $ h^- \to \ell^- \bar{\nu}_\ell a $. The instantaneous BS wave functions of the heavy mesons are applied to compute the branching fractions of such decay channels. We adopt three scenarios, that is the ALP coupling only to one charged fermion, the ALP coupling only to quarks, and the ALP coupling to all the charged fermions with the same coupling constant. In each scenario, by comparing the theoretical and experimental results, we get the upper limits of the coupling constants, which are then used to calculate the upper limits of the branching fractions of the $ B_c^- \to \ell^- \bar{\nu}_{\ell} a $ channel. For the second scenario, we also get the nonzero lower limit of the branching fraction at some range of $ c_q $. We hope this work could be helpful for the future detection of the ALP through heavy meson decays.