Our understanding of the internal structure of mesons from strong interaction physics can be tested by the decay rates of mesons annihilating to photons or through radiative transitions with the emission of a photon. In this work we use lattice QCD to calculate widths of the charm-anticharm meson decay processes 𝜂 𝑐 → 𝛾𝛾 and 𝐽/𝜓 → 𝛾𝜂 𝑐 .
The literature for 𝜂 𝑐 → 𝛾𝛾 is very unclear. The PDG [1] fits multiple sets of products of branching fractions to obtain a result with 7% uncertainty, but individually the experimental results are typically much less accurate and have a significant spread of central values. Lattice QCD does not make a clear case either: a result in the quenched approximation [2] is in tension with the PDG fit value; as are results with 𝑢/𝑑 quarks in the sea [3,4]. A further result with 𝑢/𝑑 quarks in the sea [5] is in agreement with the fit to experiment but with a large uncertainty. The work presented in these proceedings includes the effect of 𝑢/𝑑, 𝑠 and 𝑐 sea quarks and provides the first accurate lattice QCD calculation of this process.
Past results for the process 𝐽/𝜓 → 𝛾𝜂 𝑐 are in better shape, with several experimental [6][7][8][9] and lattice QCD [10][11][12][13] results in existence. Still, the lattice QCD results consistently give partial decay widths somewhat higher than the PDG average [1], which uses only the Crystal Ball [6] and CLEO results [7]. Our calculation of this process, including the realistic sea, improves the lattice QCD accuracy.
Our calculations are performed on Highly Improved Staggered Quark (HISQ) gauge ensembles generated by the MILC collaboration [14,15] that include the effect of 2 + 1 + 1 quarks in the sea. Lattice spacings range from 𝑎 ≈ 0.15 fm down to 𝑎 ≈ 0.06 fm. The light (𝑢 and 𝑑) quark masses 𝑚 𝑙 -taken to be the same -are one-fifth the strange quark mass or set to their physical value. Valence 𝑐 quarks also use the HISQ action. Ensemble details are given in Table 1 1: Details of the MILC gluon field configurations. Column 2 gives the lattice spacings. The number of lattice points in the spatial, 𝑁 𝑥 , and temporal, 𝑁 𝑡 , directions are given in column 3. The next three columns give the masses of the sea quarks in lattice units, while the final column supplies the mass of the valence charm. Set 3A is the same as set 3 except the valence charm quark has been deliberately mistuned.
Our 𝜂 𝑐 → 𝛾𝛾 calculation involves extracting the matrix element between the 𝜂 𝑐 and two on-shell photons, i.e., with squared 4-momentum 𝑞 2 = 0. To do this we construct a 3-point function from 𝑐 quark propagators
with
where 𝑗 𝜇 and 𝑗 𝜈 are c𝛾 𝜇 𝑐 and c𝛾 𝜈 𝑐, and 𝑂 𝜂 𝑐 is a pseudoscalar or temporal axial current that couples to pseudoscalar charmonium states. A diagram of this setup is shown in Fig. 1 depicting the position of the operator insertions for the 𝜂 𝑐 and vector currents 𝛾 1 and 𝛾 2 .
< l a t e x i t s h a 1 _ b a s e 6 4 = " By calculating a weighted integral over 𝑡 𝛾 1 of our 3-point correlation function to set the external photon on-shell [16,17], we get a resulting 2-point correlation function
We can then fit C𝜇𝜈 (simultaneously with the 𝜂 𝑐 2-point correlation function 𝐶 𝜂 𝑐 ) to the form
with 𝑓 (𝐸, 𝑡) = 𝑒 -𝐸𝑡 + 𝑒 -𝐸 ( 𝑁 𝑡 -𝑡 ) .
(
We have chosen 𝜔 1 = 𝑀 phys,latt 𝜂 𝑐 /2 in our calculation, where this mass is the value obtained from connected 𝜂 𝑐 correlation functions in the physical continuum limit. The sequential propagator (that which is between 𝑡 𝛾 1 and 𝑡 𝛾 2 in Fig. 1) is given spatial momentum purely in the 𝑦 direction -since the vector polarisation demands that q 1 has a component in that direction -using twisted boundary conditions [18,19]. The twist angle 𝜃 and momentum 𝑞 𝑦 1 are related by
where 𝑁 𝑥 is the number of lattice points in each spatial direction. In practice, our choice of 𝜔 1 results in photon 1 being exactly on-shell while photon 2 will be slightly off-shell. We account for this by extracting
where 𝑏 0 is the amplitude from the fit in eq. ( 4) and 𝑎𝑀 sim 𝜂 𝑐 is the mass of the 𝜂 𝑐 on each ensemble, both of which are in lattice units.
A fit to 𝐹 latt (0, 𝑞 2 ) allows us to extract 𝐹 (0, 0), which relates to the decay rate for this process by
where 𝑄 𝑐 = 2/3 is the electric charge of the 𝑐 quark and 𝛼 = 1/137 is the fine structure constant.
We fit 𝐹 latt (0, 𝑞 2 2 ) to the following form to determine 𝐹 (0, 0) in the continuum:
+𝜅 (0) sea,𝑢𝑑𝑠 𝛿 sea,𝑢𝑑𝑠 1 + 𝜅 (1) sea,𝑢𝑑𝑠 (𝑎 Λ) 2 + 𝜅 (2) sea,𝑢𝑑𝑠 (𝑎 Λ) 4 .
The 1 -𝑞 2 2 /𝑀 2 pole term addresses the fact that one of the two photons is off-shell by the inclusion of a pole [3,10]. We take the pole mass be around the 𝐽/𝜓 mass by setting a prior of 3.0(0.3) GeV for this term in our fit. The 𝛿 terms allow for the mistuning of valence and sea quark masses. The scale Λ is tuned using the Empirical Bayes criterion such that we accept the value that maximises the Bayes factor [20]. The final term allows for discretisation effects -up to (𝑎 Λ) 4 -coming from the sea, with Λ = 1 GeV. Fig. 2 (left) shows the result of the fit, with the value of 𝐹 (0, 0) in the continuum limit depicted by the black star. At leading order (LO) in nonrelativistic QCD (NRQCD) we can write
and so by taking our results for 𝐹 (0, 𝑞 2 2 ) and those of 𝑀 𝐽/𝜓 and 𝑓 𝐽/𝜓 from [21] we can assess how well LO NRQCD determines this ratio despite missing subleading terms in the velocity expansion of the heavy quark. We perform a fit with a form analogous to that in eq. ( 9) and show the result in Fig. 2 (right). This result demonstrates that LO NRQCD is already a (surprisingly) good approximation to this ratio.
We now compare our lattice QCD result to the decay widths as determined in experiment. We plot a comparison in Fig. 3 (left). The red star and band is the result from this work. The blue points are determined from experimental measurement of Γ(𝜂 𝑐 → 𝑖)Γ(𝜂 𝑐 → 𝛾𝛾)/Γ total (𝜂 𝑐 ), where the decay channel 𝑖 is given on the y-axis, combined with the experimental measurement of the branching fraction B (𝜂 𝑐 → 𝑖), with uncertainties added in quadrature. The blue square is from the product of branching fractions from 𝐽/𝜓 → 𝛾𝜂 𝑐 and 𝜂 𝑐 → 𝛾𝛾 from BESIII [22] combined with the PDG average for the branching fraction B (𝐽/𝜓 → 𝛾𝜂 𝑐 ). The PDG fit is shown by the blue dashed line and band, but it should be noted that the quality of the fit is poor (𝜒 2 = 118 with dof= 81). Our result disagrees with this PDG fit value by over 4𝜎.
We also compare our result with results using NRQCD in Fig 3 (right). The LO NRQCD result is depicted by the green band and dashed line, while the green circles show NRQCD results that include higher order corrections. The red star is our result, which has a significantly smaller uncertainty. (blue circles). The blue band is the PDG fit. Right: comparison with other theory calculations that used NRQCD. The green band is the LO NRQCD result, while the green circles show results from higher-order calculations. Our HISQ result is given by the red star.
We now turn to the decay process 𝐽/𝜓 → 𝛾𝜂 𝑐 . For an electromagnetic current 𝑗 𝜇 𝑐 = c𝛾 𝜇 𝑐 we can relate the matrix element between the 𝐽/𝜓 and 𝜂 𝑐 and form factor 𝑉 (𝑞 2 ) through
where 𝜖 𝐽/𝜓 𝜎 is the polarisation of the 𝐽/𝜓. At 𝑞 2 = 0 we can calculate the decay width to a real photon, relating it to the form factor by
where the spatial momentum | ì 𝑞| is given by
For this decay we calculate the 𝑞 2 dependence of 𝑉 -which will allow us to access to the width of the Dalitz decay Γ(𝐽/𝜓 → 𝜂 𝑐 𝑒 + 𝑒 -) at a later stage -and use this spread of data to access 𝑉 (𝑞 2 = 0). On the lattice it is useful to define a form factor V,
since we only calculate one of the two possible diagrams. These are identical as the photon can be emitted by either the 𝑐 or c. Figure 4 depicts the setup of our calculation with the relative positions of the 𝐽/𝜓, 𝜂 𝑐 and photon current insertion shown. We apply a twist to the propagator between the 𝜂 𝑐 and current 𝑗 𝜇 to give it spatial momentum. By using multiple twists we obtain form factors across the entire 𝑞 2 range.
Our final results for V (0) and Γ(𝐽/𝜓 → 𝛾𝜂 𝑐 ) will come from a simultaneous fit of the V (𝑞 2 ) data from all ensembles. In these proceedings, however, we provide preliminary results by fitting the form factor data to the form 𝐴𝑞 2 + 𝐵 and simply interpolating to 𝑞 2 = 0 for each individual ensemble. An example of one such fit is shown for Set 6 in Fig. 5 (left). The red circles are the lattice data, the blue band is the result of the fit, and the black star is its value at 𝑞 2 = 0. We then perform a chiral-continuum fit to this interpolated data for each ensemble using the form:
+ 𝜅 (0) sea,𝑢𝑑𝑠 𝛿 sea,𝑢𝑑𝑠 1 + 𝜅 (1) sea,𝑢𝑑𝑠 (Λ𝑎) 2 + 𝜅 (2) sea,𝑢𝑑𝑠 (Λ𝑎) 4 ,
the results of which are shown in Fig. 5 (right). The red points are the 𝑉 (0) data from the lattice (with the grey square denoting the mistuned 𝑚 𝑐 that is not included in this preliminary fit), the blue band the fit result and the black star the physical result in the continuum.
-0.010 -0.005 0.000 0.005 0.010 0.015 q 2 GeV 2 By using eqs. ( 12) and ( 14), we can convert V (0) from the lattice to the decay width Γ(𝐽/𝜓 → 𝛾𝜂 𝑐 ). This allows a comparison of the 𝐽/𝜓 → 𝛾𝜂 𝑐 process at 𝑞 2 = 0 with experiment, which we show in Fig. 6 (left). The green star and band is our result, the blue circles those from Crystal Ball [6] and CLEO [7], while the black cross and grey band are their average as given in the PDG data tables [1]. Finally, in Fig 6 (right) we compare our result for V (0) with those from earlier lattice calculations (blue points) where we restrict the comparison to those that used multiple lattice spacings. The ETM result includes the effect of 2 quark flavours in the sea [12], while the previous HPQCD number comes from a 2 + 1 flavour calculation [13]. The green star and band is once again our new result. with those from Crystal Ball, CLEO (blue circles) and their average as reported in the PDG data tables (black cross and band). Right: lattice results for V (0). Our result is shown by the green star and band, while a previous HPQCD result and a value from the ETM Collaboration are given by the blue circles.
Our result for the decay width for 𝜂 𝑐 → 𝛾𝛾 is the first accurate calculation from lattice QCD and gives a significantly clearer picture of this process. This work should also motivate more accurate results from experiment, which would allow a test of the Standard Model, similar to that in [21] for 𝐽/𝜓 → 𝛾 → 𝑒 + 𝑒 -.
We have also provided preliminary results for the decay width for the process 𝐽/𝜓 → 𝛾𝜂 𝑐 . The accuracy is improved over previous results and our measurement includes a realistic sea. Our full study will use results from the entire 𝑞 2 range to give an accurate determination of Γ(𝐽/𝜓 → 𝛾𝜂 𝑐 ) as well as allowing us to calculate the width of the Dalitz decay Γ(𝐽/𝜓 → 𝜂 𝑐 𝑒 + 𝑒 -).
We plan work on several follow-on studies. We can use the same approach as that for 𝜂 𝑐 → 𝛾𝛾 to determine the rate for 𝜋 0 → 𝛾𝛾 for comparison to experimental results. We can also tackle 𝜂 𝑏 → 𝛾𝛾, which has not yet been measured experimentally, to provide a robust prediction. Additionally, the presence of at least one off-shell photon in the processes 𝜋 → 𝛾 * 𝛾 ( * ) and 𝜂 𝑐 → 𝛾 * 𝛾 ( * ) gives access to associated form factors in these decays that can be used to further understand meson structure.