The study of light invisible particles in B c decays

,


Introduction
The Standard Model (SM) is extremely successful. However, it is considered to be an effective field theory which is valid only up to certain energy scale. For example, it will be invalid at the Planck scale, with gravity giving large contribution. Far below that, there are many arguments supporting that new physics (NP) will appear at the TeV scale. The NP can show itself as the missing energy in the collision at the pp or e + e − colliders. For example, CODEX-b at the LHCb experiment is proposed to probe for GeV-scale longlived particles [1]. If we assume the possible new particle to be the candidate for the dark matter (DM), the high energy collision will provide a powerful way to detect such particles. Among the DM candidates, the weakly interacting massive particle (WIMP), which appears in many theoretical models, has attracted extensive attention (see [2] for reviews). The WIMP annihilation cross section is constrained by the observed dark matter density, which sets the lower bound of the WIMP mass to a few GeV (the so-called Lee-Winberg limit [3]). However, this result is model-dependent. If the DM is nonfermionic and the weak mass scales or weak interactions are not assumed [4], this constraint can be relaxed, and more lower mass, such as a few keV, will be possible. Theoretically, this kind of light dark matter (LDM) can have different spins [5], for example, it can be a scalar particle [6], sterile neutrino [7], or hidden vector particle [8]. The MeV-scale LDM is proposed [9,10] to explain the unexpected emission of 511 keV photons from the galaxy center. Experimentally, the parameter space for the WIMP with mass larger than several GeV has been severely constrained by the recent experiment [11], which also provides a strong motivation for the study of the sub-GeV dark matter.
The LDM emission from the heavy meson decays is an interesting approach for such studies. Phenomenologically, the LDM of some hidden sector can weakly interact with -1 -

JHEP03(2019)028
the SM fermions through different ways. For example, it can couple directly to the Higgs boson [12,13]. Or there are some connectors with quantum numbers of both SM and hidden sectors. Such connector can be a chiral fermion [14] or a dark gauge boson [15]. At the energy level of heavy mesons, these processes will be greatly suppressed by the large mass in the propagator of the connector or by the small coupling constant between the connector and the SM fermions. By a model-independent way, we can introduce an effective Lagrangian to describe phenomenologically the interaction between the invisible particles and SM fermions. This method has been extensively used in refs. [16][17][18][19][20][21][22] to study the flavor-changing neutral current (FCNC) processes of K, D, and B mesons. The SM background comes from the decays with νν in the final states, which has small branching fraction and makes the detection of NP possible. The difference between the experimental results for M → M f / E and the theoretical predictions for M → M f νν in the SM will set the constrains for the LDM emission channels, where M and M f are the masses of the initial and final mesons, respectively.
The same analysis can also be applied to the B c meson. As consisting of a heavy quark and a heavy antiquark with different flavors, this meson is unique. It can only decay through weak interaction, and either the b quark or thec antiquark can be a spectator. Therefore more possible decay modes involving the invisible particles are allowed. Experimentally, there are abundant B c samples are collected at the LHC [23][24][25], which gives us the chance to study its various decay channels precisely, especially the rare decays. Until now there is no experimental data for such decays of the B c meson available, so we expect detections in the near future. Theoretically, many methods have been applied to study the semileptonic, nonleptonic, pure leptonic, and the FCNC processes of the B c meson [26][27][28][29][30]. In this work, we will apply the instantaneous Bethe-Salpeter method to calculate the hadronic transition amplitude when both the initial and final mesons are heavy. This method has been used extensively to study the weak decays of B q mesons [31,32], and gotten consistent results with experiments. In the LDM emission processes of the B c meson, this method is still valid and the calculation steps are the similar to those in the SM.
The rest of the paper is organized as follows: in section 2, we first construct the effective Lagrangian which describes the coupling between quarks and light invisible particles. Then by comparing the theoretical and experimental results, we extract the upper and lower limits of the Wilson coefficients. In Set. III, these limits are used to constrain the branching fractions of the decay channel B c → hχχ with h and χ being the final meson and the invisible particle, respectively. Finally, we give the summary and perspective in section 4.

χ is a scalar
At the quark level, the χ emission processes of the heavy meson can be described by the effective Lagrangian [18],

JHEP03(2019)028
Experimental bound [41][42][43] SM prediction [37][38][39][40] Invisible particles bound  where q and q f are the Dirac spinor fields of the initial and final quarks, respectively; g s1 and g s2 are the phenomenological coupling constants. This Lagrangian is model-independent. And for specific models, the four-particle vertex may be generated at the tree or loop level [16][17][18][19][20] by introducing other new particles. In this work, we will not focus on any specific model, but consider the FCNC processes of B c meson induced by such effective operators. Theoretically, there are many studies [33][34][35][36] of the FCNC processes of the B c meson, while the corresponding detection is still missing. So we cannot use the experimental data of B c meson to set constraints on the coupling constants. Our strategy is in the opposite direction. That is, the allowed-region of the coupling constants from other processes are used to constraint the branching ratios of the B c decays. Experimentally, there are data for such decays of B and D mesons. The corresponding channels are The experimental bounds for their branching ratios are listed in table 1. Within the SM, the missing energy / E represents the νν pair, and the branching fractions are calculated in refs. [37][38][39][40]. The difference between theoretical predictions and experimental bound allows the existence of NP. Here the NP processes are described by the Feynman diagrams in figure 1.
For the B − → π − (K − )χχ processes, only the scalar current gives contribution to the transition amplitude, which can be written as  Table 2. Parameters in the form factors of the B → π(K) processes [44].
where P and P f are the momenta of the initial or final mesons, respectively; m q and m q f are the masses of quarks; s is defined as (P − P f ) 2 . In the first step, the equation of motion is used. The hadronic transition matrix is parameterized as the form factors f + and f 0 .
Here we adopt the results of the QCD light-cone sum rules (LCSR) [44], where the form factors are constructed as The corresponding parameters are presented in table 2.
For the B − → ρ − (K * − )χχ processes, only the pseudoscalar current gives contribution to the transition amplitude, which has the form, where is the polarization vector of the final meson; M and M f are the masses of the initial and final mesons, respectively; A 0 , A 1 , A 2 , and A 3 are form factors.  Table 3. Parameters in the form factors of the B → ρ(K * ) processes with k max = 2 [45].
In ref. [45], Bharucha et al. also used the LCSR, but with a different parameterization method, to write the form factors as where F 1 , F 2 , F 3 , and F 4 represent A 0 , A 1 , A 12 , and V , respectively; P i (s) = (1 − s/m 2 R,i ) −1 represents the pole structure. And z(s) is defined as The related parameters are listed in table 3. The form factors A 12 and A 3 are related to A 1 and A 2 by By finishing the three-body phase space integral, we get the branching ratios where λ(x, y, z) = x 2 + y 2 + z 2 − 2xy − 2xz − 2yz is the Källen function; m χ is the mass of the invisible particle; θ is the angel between the three-dimensional momenta P χ and P f   in the momentum center frame of invisible particles; Γ B − is the total width of B − meson; Ω = 2 originates from the final two invisible particles being identical. For the annihilation processes of B 0 , D 0 , B 0 s mesons, that is M → χχ, only the pseudoscalar current contributes to the decay amplitude, which has the form, where f M is the decay constant of the initial meson, which has the values: 216 GeV and f D 0 = 0.230 GeV [46]. By finishing the two-body phase space integral we get the partial width, In figure 2, we plotΓ = Γ/|g s2 | 2 as a function of m χ . One can see that they all have same trend that decrease when m χ gets larger, because the phase space gets smaller. By comparing the theoretical predictions and the experimental upper limits (the third column of table 1) of the branching ratios for these decays, we can set the upper bounds for the effective coupling constants g s1 and g s2 with specific mass of the invisible particle. The results are shown in figure 3 (represented by the solid lines). One can see that as m χ increasing, the upper limits of the effective coupling constants get more and more larger. The reason is simple: larger m χ means more suppression from the phase space. So from these decay channels we can set more stringent upper limits for the effective coupling constants when 2m χ is not close to the threshold. The B → π(ρ)χχ channel gives smaller bound of |g s1 | 2 (|g s2 | 2 ) compared with the B → K(K * )χχ channel. The B → χχ mode gives the most stringent upper bound of |g s2 | 2 , because the two-body phase space is larger than the three-body case. We also present the result from D 0 decay, which is larger due to its smaller mass. For B 0 s , the experimental results are still missing. Once the experimental data for the annihilation channel are available, they can also be used to set the upper limit of |g s2 | 2 .
The discussions above are model-independent except the calculations of the hadronic transition matrix. If we introduce some specific models, more information can be extracted. For example, in ref. [14], McKeen introduced an effective Lagrangian, where Y is a heavy fermion and Λ is the coupling constant. With the B 0 s −B 0 s mixing (see figure 4), a mass difference can be estimated to be to m q g s1 . Experimentally, the latest value of ∆M Bs−Bs is 17.63±0.11 (stat)±0.02 (syst) ps −1 , which comes from LHCb [47]. This sets an upper limit for the contribution of the light invisible particle. And by assuming m Y 400 GeV [14] we can estimate the upper limit of g s1 which is shown by the green solid line in figure 3a.
One can see it's a very loose restriction compared with other decay modes. So we will not use this result in the B c case.
The lower bound for the effective couplings can be set by considering the relic density [48,49] where χ is the candidate of dark matter which assumed to be a scalar and SM singlet. It can annihilate into the SM particles by introducing the effective Lagrangian [12,16] -7 -

JHEP03(2019)028
where H is the SM Higgs doublet. As calculated in ref. [12], the annihilation cross section has the form where x f = m f /m χ and f refers to the SM fermions. Combining eqs. (13) and (15), we can get the lower limit forΛ, which is related to g s1 and g s2 by [12] where G F is the Fermi coupling constant; V q 1 q 2 is the Cabibbo-Kobayashi-Maskawa (CKM) matrix element. So the lower limits of g s1 and g s2 are also obtained which are represented by the dashed lines in figure 3. The shadow areas is allowed by both constraints (meson decays and relic density). When the parameters are in this region, the invisible scalar particle can be a candidate of the DM. For the regions below the dashed lines, the scalar invisible particle can also be possible to exist as a portal DM, which is a mediator between SM and hidden sectors [9, 50-52].
For the decays of B meson, when the final meson is a pseudoscalar, the second and the fourth terms in eq. (2.17) will not contribute to the decay. The FCNC process can be induced by the scalar or vector current, and the transition amplitude has the form where P 1 and P 2 are the four-dimensional momenta of χ and χ † , respectively; h − is π − or K − . The hadronic transition matrix element is parameterized the same as that in eq. (2.2) or eq. (2.4), and the form factors are expressed in eq. (2.3) or eq. (2.5).
The transition amplitude receives the contribution from two terms in the effective Lagrangian, and the partial width can be written as Here we have defined Γ 1(3) = dP S 3 |T 1(3) | 2 , which are independent of the effective coupling constants. The interference terms are proved to be zero.   Comparing the theoretical predictions and the experimental upper bound of these channels, we give the possible relations of the modulus square of the effective coupling constants, which are presented in figure 5. In this figure, the area below the colored line is allowed experimentally with a specific mass of χ. One can see that as m χ increasing, the allowed region gets larger and larger.
If the final meson is a vector, the situation is a little more complicated, because this time the decay processes can be induced by the second, the third, and the fourth operators in the effective Lagrangian. The transition amplitude is where h * − represents ρ − or K * − . Here we need to consider two kinds of hadronic transition matrix elements. h * − |q f γ µ γ 5 q|B − is parameterized the same as eq. (2.4).
where V (s) is expressed by eq. (2.5) and the parameters are given in table 3. The relationship between three effective couplings can be achieved by comparing the theoretical results and the experimental upper limits. Numerical calculation indicates the cross terms can also been neglected. In figure 6, we show that the experimentally allowed region is that under the colored plane which corresponding a specific mass of χ.

The decay modes of the B c meson
In the previous section, we have used the results of LCSR to study the FCNC processes of B meson. This method is suitable for the heavy-light state. For the B c meson, which consists of a heavy quark and a heavy antiquark, we choose the BS method to study its where P is the momentum of the meson; p 1 and the p 2 are the momenta of the quark and antiquark, respectively; m 1 and m 2 are the masses of the quark and antiquark, respectively; q is the relative momentum between quark and antiquark; χ P (q) is the BS wave function; V is the interaction kernel. For B c meson, we can safely make an instantaneous approximation for V , that is V (P, k, q) ≈ V (P, k ⊥ , q ⊥ ), where q ⊥ = q − P ·q √ P 2 P , and the same is for k ⊥ . By defining the Salpeter wave function ϕ(q ⊥ ) = i dq 0 2π χ P (q), we reduce eq. (3.1) to the three-dimensional form, which can be solved numerically. ϕ(q ⊥ ) is constructed from P , q ⊥ , Dirac gamma matrices, and some scalar function of q 2 ⊥ . We take the 0 − and 1 − states as examples, whose Salpeter wave functions are [54]

(3.2)
In Mandelstam formalism, the hadronic transition matrix element can be expressed as the overlap integral of the BS wave functions of the initial and final mesons. With the instantaneous approximation, it can be reduced to the overlap integral of Salpeter wave functions. To make the calculation simple, we just keep the positive energy parts of the wave functions which give the main contribution. The transition amplitude is [55]

The SM backgroud
In the Standard Model, the missing energy in the decay processes B − c → D ( * )− s + / E is carried by the (anti)neutrino. The corresponding Feynman diagrams are presented in figure 7. It can be described by an effective Lagrangian where G F is the Fermi coupling constant; α is the fine structure constant; θ W is the Weinberg angle; V q 1 q 2 is the CKM matrix element; X l (x t ) is the Inami-Lim function [56], which has the form The hadronic transition matrix element is calculated by eq. (3.3). The branching fraction is achieved by finishing the three-body phase space integral, which is presented in table 4 to compare with the results of other models. The errors come from varying the parameters in our model by ±5%.
There are also the B c → B u νν processes, which is induced by c → u at the quark level. The Feynman diagrams for such channels are given in figure 8. The corresponding effective Lagrangian is
In the above equation, we have defined X l (x q ) =D(x q , y l )/2, and the Inami-Lim function D(x q , y l ) is expressed as [56] D(x q , y l ) = 1 8 x q y l x q − y l y l − 4 y l − 1 2 ln y l + 1 8 where we have used x q = m 2 q /M 2 W and y l = m 2 l /M 2 W . The transition amplitude has the form and s , or B ( * ) , and the hadronic transition matrix elements are calculated with eq. (3.3). By finishing the three-body phase space integral, we get the decay widths expressed as the product of the squared effective coupling constant |g si | 2 and the quantity Γ i . Γ i is independent of the coupling constants, and can be calculated by taking a specific value of m χ . In figure 9, we plot them as functions of m χ . One can see that they all decrease when m χ gets larger, because the phase space gets smaller. With the same value of m χ , Γ 1 from the B c → P χχ channels is larger than Γ 2 from the B c → V χχ channels due to the different effective vertex in the amplitude. For Γ 1 , the B c → D s channel gives a larger result than that of the B c → D channel when m χ < 1.9 GeV. When m χ gets even larger, the phase space suppression will be important. For Γ 2 , this turning point is about 1.53 GeV. We also notice that Γ i of the c → u processes is two orders of magnitude less than that of the b → d(s) processes. This comes from both the smaller phase space and smaller m q for the former case.
The upper limits of the squared coupling constants |g s i | 2 with different values of m χ have been given in figure 3. Combining the results in figure 3 and figure 9, we can make predictions of the upper limits of the branching ratios of the B c → P (V )χχ channels. In figure 10 we present the results which are represented by the red solid lines. For the B c → Dχχ and B c → D s χχ processes, we use the upper limits extracted from the B → πχχ and B → Kχχ channels, respectively (see figure 3a). For the B c → D * s χχ channel, the result of B → K * χχ is used. For the B c → D * χχ case, there are two processes available -13 -

JHEP03(2019)028
to set the upper limits, namely B → ρχχ and B → χχ. The later one gives more stringent constraint, which is applied here. As there is no experimental result now for the D → π+ / E, we cannot give the constraints for the B c → Bχχ channel. The |g s2 | 2 extracted from the D 0 → χχ is used to set the upper limit for the branching ratio of the B c → B * χχ channel.
The upper limits of the branching ratios of B c → P (V )χχ channels are of the order of 10 −6 when m χ is not close its maximum value. For B c → Dχχ and B c → D ( * ) s χχ, a specific feature appears. When m χ is less than about 1.5 GeV, the branching ratios increases slowly with m χ ; after that, the branching ratios decreases rapidly to zero. It is the result of a combination of the increasing |g si | 2 and decreasing Γ i . One notices that in figure 9a the Γ of the D * (s) case is smaller than that of the D (s) case, however, the branching ratios of the former are several times larger than that of the later, because the experimental upper bound of B and D mesons in table 1  E modes, the upper limits of their branching ratios are very close to that of the corresponding χχ channel. This is because the SM background is small, for example, the branching ratios of B c → B ( * ) νν is of the order of 10 −14 . This provides a way to test our results. If the future experiments find a quit large branching ratio of such channels compared with the SM prediction, it definitely indicates the existence of some new physics. We also present the lower bounds of the branching ratio which come from the constraints of the relic density if χ is a dark matter. They are represented by the dashed lines in figure 10. The shadow areas are the allowed regions of the branching ratios. It should be pointed that these lower bounds are model-dependent.

The B c → P (V )χχ † process
If χ is a pseudoscalar, the Lagrangian L 2 is applied. When the final meson is a pseudoscalar, the transition amplitude has the form and for the vector meson case, the transition amplitude can be written as 14) The hadronic transition matrix elements are also calculated with eq. (3.3). But this situation is more complicated, because there are two or three operators contribute to the decay. Similar to the subsection 2.2, we will neglect the cross terms, and keep the ones proportional to |g pi | 2 . These terms are named as Γ i , which are independent of the effective coupling constants but depend on the mass of the light dark matter. To compare the contribution of different terms, we calculate them by giving a specific value of m χ . The results are presented in figure 12. We can see all the Γ i s are decreasing     with m χ , which is due to the suppression of phase space. When the final meson is D or D s , Γ 1 is large than Γ 3 , while for B, the situation is very different. When the final meson is a vector, Γ 2 and Γ 4 are close to each other, but both larger than Γ 4 . The approach we get the upper limits of the decay width is as follows. When m χ is given, there is an experimental allowed region for the effective coupling constants which are presented in figure 5 and figure 6. So we scan the parameter space to get the maximum value of the partial decay width. The results are given in figure 12. For the B c → B ( * ) χχ † channels, as there are no constraints for the effective coupling constants available now, so the upper limits of their branching ratios cannot be calculated.
One notices that, the figure 12a and figure 10a, 10b are almost exactly the same though the calculation processes are very different. For the later, χ is a scalar particle, and only the scalar coupling operator takes effect. For the former, χ is a pseudoscalar, and two effective operators make contribution to the branching ratios. When we scan the parameter space in figure 5, which are right triangular regions, we find that if the right endpoint of the hypotenuse is taken, the branching ratio will achieve the maximum. This also means only the scalar operator should be considered. So when we calculate the upper limit of the branching ratios, the operator contributes to the decay modes of figure 12a is just the same as that contributes to figure 10a, 10b, which makes their results are the same. But the condition in figure 12b is quite different. One can see the upper limits of the branching ratios are  larger than those in figure 10. When we scan the parameter space in figure 6, we find the axial vector operator provides most of the contribution, which is different with the case in figure 10c, 10d, where only the pseudoscalar operator takes effect. For the B − c → D * − χχ † channel, there is a kink when m χ is larger than 2 GeV. This is because the corresponding operator which has the most important contribution turns to the vector from the axial vector.
In figure 13 we present the differential distribution of the upper limits of the widths as a function of s. As examples, two cases with m χ = 0 GeV and 0.25(M −M f ) are considered both for χ being a scalar or a pseudoscalar particle. In figure 13a and figure 13c, the lines for the decay modes with χ being a scalar or a pseudoscalar coincide. The reason for this   is the same as that mentioned in the previous paragraph. The mass of χ determines the lower bound of s, namely, the left starting points of the curves. It is interesting to notice that for the B c → D ( * ) χχ † and B c → D ( * ) s χχ † channels, the peaks of the distribution curves are always in the position s = 16 GeV to 18 GeV, which is almost independent of m χ . The distribution curves for the B c → B * χχ mode are a little bit different. When m χ is very small, there is no peak. For comparison, we also plot the differential widths for the B c → D ( * ) s νν channels, which are smaller than those of invisible particles channels in most regions of s. For the channels B c → D s νν, B c → D * s νν, and B c → B ( * ) νν, their decay widths are too small to be shown in figure 13c, figure 13d, and figure 13e, respectively.

Conclusion
We have studied the light invisible particles (m χ is less than several GeV) through the rare decays of the B c meson. These particles can be the candidates of the light dark matter when the parameters taking specific values constrained by the B meson decays and the relic density. Both the scalar and pseudoscalar cases are considered. Effective Lagrangians which contain the dimension-six operators are constructed to generate such processes. The effective coupling constants are constrained by the experimental results for the B and D decays with missing energy. Then the upper limits of the branching fractions of the B c → P (V )χχ and B c → P (V )χχ † channels are calculated. For the former, when the final meson is D ( * ) s , the largest value of the upper limits is of the order of 10 −6 ; for the later, the largest value is of the order of 10 −5 when the final meson is D * s . Although the results change with m χ , their orders of magnitude almost have no change if m χ is not close to the threshold. Considering that the SM background is very small for some channels, we hope that the future experiments will find something new through such processes or set more stringent constraints for them.