Charmonium resonances from 2+1 flavor CLS lattices

TL;DR

This paper advances lattice QCD methods to study charmonium resonances, focusing on coupled-channel scattering and spectrum analysis using CLS ensembles, addressing the complexities of open decay channels.

Contribution

It introduces new techniques for calculating single-channel and coupled-channel scattering matrices for charmonium states on lattice ensembles.

Findings

Progress in determining scattering matrices in scalar and vector channels.

Updated results on charmonium spectrum in moving frames.

Enhanced understanding of open decay channels in lattice calculations.

Abstract

Many exotic charmonium resonances have been identified recently in experiment, however their nature and properties are mostly unknown. Algorithmic and theoretical progress in lattice calculations has enabled reliable numerical investigation of the spectrum below the strong decay threshold, while the study of charmonium resonances remains an open challenge. The main difficulty to overcome is the presence of many open decay channels which are coupled together, resulting in a complex finite volume quantization condition. We report on our recent progress towards the determination of single-channel and coupled-channel scattering matrices in the scalar and vector channels on CLS ensembles. We also present an update concerning the study of the charmonium spectrum in moving frames.

Equations

$$\frac{k^3\cot (\delta )}{\sqrt{s}}=A+Bs+C{s}^2.$$

$$k^3\cot (\delta )=-k^2\sqrt{k^2}$$

$${\tilde{K}}^{-1}=\left[\begin{array}{cc}{c}_{11}s+b_{11}& b_{12}\\ b_{12}& b_{22}\end{array}\right]$$

$$\Omega (\mu ,A(s))=\frac{\det (A)}{\det (({\mu }^2+A{A}^{†}{)}^{\frac{1}{2}})},$$

Full text

1 Introduction

During the last decade collider experiments have been able to identify many exotic QCD resonances that do not fit in the conventional quark model. These are possible candidates for pentaquarks, tetraquarks or hadronic molecules. The discovery of new exotic states has led to an increasing interest in understanding their properties from a phenomenological point of view. Lattice Monte Carlo simulations can provide an ab-initio approach to the study of QCD resonances, however a full consideration of all aspects and systematics of resonances on the lattice is still an open challenge, especially in the heavy flavor sector.

Here we focus on the charmonium spectrum and we present the recent progress of our collaboration towards the study of excited charmonium states and exotic resonances. We first determine the charmonium spectrum on the lattice in the single hadron approach, i.e. neglecting the effects of any strong decay modes. Here we extend the study of the spectrum from the rest frame [2] to two moving frames [6]. We are also exploring the Lüscher’s finite volume method to determine scattering amplitudes and extract scalar and vector charmonium resonances, extending the first results of Ref. [1].

2 Charmonium states in the single-hadron approach

The starting point of any lattice hadron spectroscopy calculation is the determination of the two-point correlation functions $C_{ij}(t)=\langle O_{i}(t)O_{j}^{\dagger}(0)\rangle$. The operators $O_{i}$ are classified according to their transformation properties on the lattice. Single- and multi- hadron operators are considered. The energy levels extracted from the solution of the Generalized EigenValue Problem (GEVP) using only single-hardon operators already provide information about the bound states of the charmonium spectrum and some indication about the resonant states above threshold.

We have measured the charmonium correlation functions on the $N_{f}=2+1$ ensembles labeled U101 and H105 generated by the CLS collaboration [4, 5]. The lattice spacing is $a=0.08636(98)(40)$ fm, and the pion and the kaon masses are equal to $m_{\pi}\simeq 280$ MeV and $m_{K}\simeq 467$ MeV. The lattice volume is $96\times 32^{3}$ for H105 and $128\times 24^{3}$ for U101. The gauge and the fermion fields fulfill open boundary conditions in the time direction. The bare charm quark mass $m_{c}$ is a parameter to be tuned. We have chosen two different values of $m_{c}$ corresponding to the $D$ meson mass being 80 MeV below and above its physical value to estimate the systematic errors and to explore how the charmonium spectrum and the exotic resonances are influenced by the mass of the charm quark.

We have first determined the charmonium spectrum in the single hadron approach on the ensemble U101 at rest and in moving frames. In both cases, many different lattice irreducible representations need to be considered in order to be able to identify the underlying continuum quantum numbers $J^{P}$ of the lattice energy levels. For further details about our results and our methodology, see Ref. [6].

3 Scattering in a finite box: Lüscher formalism

The study of the physical properties of a resonance on a Euclidean space-time requires a different formalism that is beyond the simple estimation of energy levels from the exponential decay of correlation functions of $\bar{c}c$-operators. The Lüscher formalism relates the phase shift to the difference between the energy of the non-interacting system to the measured energy levels in a finite box [3]. The most challenging aspect of such a method is that a very precise determination of the lattice spectrum is required possibly on many different volumes or in many different momentum frames to reliably determine the scattering amplitudes.

In the elastic case, assuming that only a single scattering channel is coupled to the resonance under investigation, the Lüscher formalism provides for each energy level a corresponding phase shift $p^{2l+1}\cot(\delta)$. In the case of two channels coupled together, the properties of the scattering matrix in the complex plane are accessible by fitting a suitable parameterization of the $K$-matrix to all energy levels [7]. We explore the elastic case in the vector channel and the coupled channel scenario for the scalar channel.

3.1 Vector channel

The lowest states beyond the $J/\psi$ in the charm vector meson channel are $\Psi_{2S}$, $\psi(3770)$ and $\psi(4040)$. $\psi(3770)$ is the first resonance just above the open charm $\bar{D}D$ decay threshold. Therefore we include in our basis both at rest and in moving frames the corresponding two-particle $\bar{D}D$ interpolators, neglecting the effects from a $3^{--}$ resonance and charm annihilation diagrams.

The extracted phase shift as function of $k^{2}=s/4-m_{D}^{2}$ is presented for U101 in Fig. 1 and for H105 in Fig. 1 at the charm quark mass corresponding to $\kappa_{c}=0.12522$. The first natural assumption to interpolate our data and parametrize the phase shift is a quadratic fit of the form

[TABLE]

We use the package TwoHadronsInBox for fitting our energy levels to the Lüscher formula [7]. The lowest energy level, corresponding to $J/\psi$, is excluded from our fit. The quadratic ansatz works reasonably well for H105 with a $\chi^{2}/\textrm{d.o.f.}=2.2$, while a higher $\chi^{2}/\textrm{d.o.f.}$ is obtained for U101, indicating that higher orders should possibly be included in Eq. 1.

A bound state appears as a pole of the $T$ matrix on the real axis below the scattering threshold in the first Riemann sheet. Such a pole appears if the condition

[TABLE]

is fulfilled. The above condition is depicted in Fig. 2. For the ensemble U101 there are two crossings of the red and blue curves corresponding to two bound states, while for H105 there is a single crossing corresponding to the $\Psi(2S)$. Given that both ensembles have been generated with the same bare lattice parameters but different lattices sizes, our findings lead to the conclusion that the subleading exponentially suppressed finite volume effects are relevant for U101 in this particularly delicate situation when $\psi(3770)$ lies almost on threshold. It results from the fact that the points corresponding the naive energy levels for the $\psi(3770)$ lie at negative $k^{2}$ in Fig 1, so this result is independent of the parametrization chosen. Therefore we will focus on the H105 ensemble for the scattering analysis of the near-threshold resonance $\psi(3770)$.

The study of a scattering channel is performed by looking to the structure of the $T$ matrix in the complex $s$-plane. There is a single pole on the real axis of the $T$-matrix in the first Riemann sheet, corresponding to the bound state $\Psi(2S)$, and there are two complex conjugated pairs of poles in the second sheet corresponding to the $\psi(3770)$ resonance, see Fig. 3. We have also explored different parametrization of our vector channel data to estimate the systematic errors and to find a functional form that is able to fulfill all physical requirements, such as the sanity check of Ref. [8]. The full results of our analysis will be presented in a forthcoming publication.

3.2 Scalar channel

The $0^{++}$ charmonium channel is not yet experimentally well understood. Beyond the ground state, the $\chi_{c_{0}}$, two possible scalar resonances have been identified: the recently discovered $\chi_{c0}(3860)$ [9] and the $X(3915)$ [10], for which an alternative interpretation as a $2^{++}$ state has been suggested in Ref. [11]. Both states are quite above the $\bar{D}D$ threshold and the $X(3915)$ is slightly below the $\bar{D}_{s}D_{s}$ threshold. In our setup we expect a significant influence of both $\bar{D}D$ and $\bar{D}_{s}D_{s}$ on our energy spectrum and we include their respective interpolators in our basis for the GEVP. For this first analysis, we however exclude $J/\Psi-\omega$ and $\bar{D}^{*}D^{*}$ interpolators and we do not consider a possible $2^{++}$ resonance. We assume that these scattering channels have negligible influence.

If two decay channels are coupled, the $K^{-1}$-matrix can be parametrized in terms of three functions of $s$, the simplest ansatz for them is a polynomial in powers of $s$. We have tried various possibilities and we find the form

[TABLE]

is the best compromise between constraining our fitting parameters while keeping a good $\chi^{2}/\textrm{d.o.f}\approx 1$ (notations for $\tilde{K}^{-1}$ are as in Ref. [7]). In the coupled-channel case, the quality of fit can be judged by looking to the difference between our extracted energy levels and the zeros of the $\Omega$ functions in the various lattice irreducible representations. The $\Omega$ function, as introduced in Ref. [7], is defined as

[TABLE]

with $A(s)$ related to the box-quantization matrix and $\mu$ a regulator parameter. Our analysis is performed combining the U101 and H105 data from two different frames together. The quality of our fits is presented in Fig. 4, where each point on the horizontal axis represents one of our lattice energy levels and the crossing of the blue curve to the real axis is the predicted energy level from our parametrization.

In the case of two coupled channels, there are four Riemann sheets in the complex plane. These four sheets are labeled according to ($sign(Im[k_{m_{D}}]),~{}sign(Im[k_{m_{D_{s}}}])$), where $k_{m_{i}}^{2}=s/4-m_{i}^{2}$, in $\bar{D}D$ and $D_{s}D_{s}$ scattering. In Fig. 5 we show the preliminary results for the $T$-matrix in the complex plane for each Riemann sheet. We plan to investigate the pole structures emerging in all complex planes and to extract the physical parameters of the scalar resonances in the near future.

4 Conclusions

The vector resonance $\psi(3770)$ requires a delicate tuning of the lattice volume and of the charm mass such that the strong decay in $\bar{D}D$ is allowed. The study of scalar resonances is challenging given the presence of multiple decay thresholds in the region of experimental interest, however we have been already able to identify a good parameterization of the $K$ matrix. There are still many systematic uncertainties and assumptions that have to be addressed to provide quantitative results that can be directly matched with experiments.

5 Acknowledgements

We thank Ben Hörz for his help and support with the package TwoHadronsInBox. Our code that implements the distillation method is written within the framework of the Chroma software package [12]. The inversion of the Dirac operator is performed used the multigrid solver of Ref. [13, 14]. The simulations were performed on the Regensburg iDataCool cluster, and the SFB/TRR 55 QPACE 2 [15] and QPACE 3 machines. We thank our colleagues in CLS for the joint effort in the generation of the gauge field ensembles which form a basis for the here described computation. The Regensburg group was supported by the Deutsche Forschungs-gemeinschaft Grant No. SFB/TRR 55. M. P. acknowledges support from the EU under grant no. MSCA-IF-EF-ST-744659 (XQCDBaryons). S. Prelovsek was supported by Slovenian Research Agency ARRS (research core funding No. P1-0035 and No. J1-8137).