Nucleon-$J/\psi$ and nucleon-$\eta_{c}$ scattering in $P_{c}$ pentaquark channels from LQCD

The lattice QCD simulation of $NJ/\psi$ and $N\eta_c$ scattering is performed at $m_\pi\simeq 266~$MeV in channels with all possible $J^P$. This includes $J^P=3/2^\pm$ and $5/2^\pm$ where LHCb discovered $P_c(4380)$ and $P_c(4450)$ pentaquark states in proton$-J/\psi$ decay. This is the first lattice simulation that reaches the energies 4.3-4.5 GeV where pentaquarks reside. Several decay channels are open in this energy region and we explore the fate of $P_c$ in the one-channel approximation in this work. Energies of eigenstates are extracted for the nucleon-charmonium system at zero total momentum for all quantum numbers, i.e. six lattice irreducible representations. No significant energy shifts are observed. The number of the observed lattice eigenstates agrees with the number expected for non-interacting charmonium and nucleon. Thus, we do not find any strong indication for a resonance or a bound state in these exotic channels within one-channel approximation. This possibly indicates that the coupling of $NJ/\psi $ channel with other two-hadron channels might be responsible for $P_c$ resonances in experiment. One of the challenges of this study is that up to six degenerate $J/\psi(p) N(-p)$ eigenstates are expected in the non-interacting limit due non-zero spins of $J/\psi$ and N, and we establish all of them in the spectra.


I. INTRODUCTION
Gell-Mann indicated in his 1964 paper [1] that along simple baryons (qqq) and mesons (qq) also states with more complex structure could exist. This could be particles with mesonic quantum numbers qqqq (tetraquarks) or states with baryonic quantum numbers qqqqq or qqqqqq (pentaquarks and dibaryons). In recent years a wide variety of those states were observed by several experiments. Until 2015, almost all of the experimentally confirmed exotic hadrons -so called XYZ states-carried mesonic quantum numbers. The confirmed exception with baryonic quantum numbers is dibaryon deuteron.
In 2015, two peaks in proton-J/ψ invariant mass with minimal flavor structure of uudcc were observed by LHCb [2]. Discovery was later confirmed by the model independent study in 2016 by the same collaboration [3]. Two resonances were observed: the lighter with the mass M 1 4380 MeV has broad width Γ 2 205 MeV, while the heavier state with M 2 4450 MeV is narrower with Γ 2 40 MeV. Resonances were identified as hidden charm pentaquarks P c . LHCb finds the best fit for spin-parity assignments (J Strong interaction allows P c uudcc to decay to a charmonium (cc) and a nucleon (uud), as well as to different combinations of a charmed meson (qc) and a charmed * ursa.skerbis@ijs.si † sasa.prelovsek@ijs.si charmed baryon (qqc). Some of the allowed decay channels are summarized in Table VII in appendix A. There were several attempts to explain the origin and the structure of two charmed pentaquarks phenomenologically. They were predicted in the coupled channel dynamics in [4,5], where the significant componentD * Σ c couples also to N J/ψ via the vector meson exchange. They were studied as hadronic molecules of a charmed meson and a charmed baryon, for example in [6][7][8][9][10][11], or charmonium and nucleon, for example in [12]. All listed studies report to find an indication for P c . A review of hadronic molecules including P c can be found in [13]. These studies are often based on the phenomenological meson-exchange models [14][15][16][17]. Several studies considered compact diquark-diquark-antiquark internal structure and find P c , for example [18][19][20][21][22]. The heavier P c was shown to be compatible with the kinematical effects of rescattering from N χ c1 to N J/ψ [22]. Other attempts to explain P c as a coupled channel effect consider rescattering with a charmed meson and a charmed baryon [23][24][25]. Enhancements corresponding to P c were also found as a consequence of the triangle singularity transitions in the kinematics of Λ b → K − pJ/ψ [22,26]. Further phenomenological studies can be found in the recent reviews [13,[27][28][29][30].
On the other hand, there is no knowledge on P c resonances based on the first-principle lattice QCD at present. The lattice simulations of systems with flavorccuud have never reached energies 4.3 − 4.5 GeV where the pentaquarks reside. They mostly considered interactions of a nucleon and a charmonium near threshold. Slightly attractive interaction was found in the preliminary study of N J/ψ and N η c scattering for m π = 293 − 598 MeV, where no bound state or reso-nances appear [31]. Similar conclusion was previously in the quenched simulation [32] at m π 197 MeV. Interestingly, the dynamical NPLQCD study of N η c scattering in s-wave [33] rendered a bound state about 20 MeV below threshold for very heavy m π 800 MeV. The simulation [34] presents preliminary results for N J/ψ and N η c potentials and phase shifts in s-wave using HALQCD method in one-channel approximation. These were extracted up to the energies 0.2 GeV above threshold at a very heavy m u/d . Attractive interaction was found in all channels explored, but not attractive enough to form bound states or resonances. Similar approach and conclusions were obtained in an earlier study [35,36]. The hadroquarkonium picture was considered in [37], where the staticcc potential V (r) was extracted for m c → ∞ as function of distance r in the presence of the nucleon. The potential is found shift down only by a few MeV due to the presence of the nucleon.
Effective field theory for s-wave quarkonium-nucleon system near threshold is developed in [38], where lowenergy constants are determined from lattice data [31,32,35,36]. This approach also does not feature bound states or near-threshold resonances.
The aim of this lattice study is to establish whether P c resonances could appear in the one-channel approximation for N J/ψ or N η c scattering. This is the first lattice study that reaches energies 4.3 − 4.5 MeV where P c resonances are found in experiment. We neglect contributions from other channels and coupled-channel effects. This study could shed light on whether the experimental P c resonances are crucially related to the coupledchannel effects or not. The purpose is to calculate the eigen-energies of the interacting N J/ψ system on the lattice in the rest frame for all possible J P . This includes previously studied partial wave = 0 and for the first time also > 0. The eigen-energies are then compared with (i) the energies expected in the limit when N and J/ψ do not interact and (ii) the energies expected based on the experimental P c resonances in the one-channel approximation. Analogous approach is followed for the N η c system. The results shed light on the fate of P c in the one-channel approximation.
Let us discuss the expectation for the spectra of N M system (M = J/ψ or η c ) in the limit when hadrons do not interact, as this will be an important reference case. The momenta p = n 2π L of each hadron are discretized due to the periodic boundary conditions on our lattice with L 2 fm ( 1 ). The total energies of non-interacting (n.i.) N (p)M (−p) system are (1) with n ∈ N 3 , while p with be denoted by p for simplicity from here on. The E H=N,M (p) are single-hadron energies measured for various momenta on our lattice ( Table   1 Symbol L denotes the lattice size, while denotes partial wave. 1); they would satisfy E H (p) = (m 2 H + p 2 ) 1/2 in continuum. The non-interacting energies (1) for both channels are presented in Fig. 1, together with the location of P c resonances from experiment. In order to capture the resonances region, we explore both channels up to the energy region which captures the states with relative momenta p 2 ≤ 2 (2π/L) 2 and slightly above. The N J/ψ scattering is investigated up to approximately E ≤ m sa + 1.6 GeV, while N η c is studied up to E ≤ m sa + 1.5 GeV, where reference energy is m sa = (m ηc + 3m J/ψ )/4 (16).
Dashed lines represent the non-interacting energies E n.i. = EN (p) + Ecc(−p) (1) for the nucleon-charmonium system with total momentum zero on our lattice: black for p 2 = 0, red for p 2 = 1 and blue for p 2 = 2, where p 2 is given in units of (2π/L) 2 . The experimental masses of both Pc are also given. Energies are presented with respect to spin-averaged charmonium mass (mη c + 3m J/ψ )/4 (16).
Furthermore, it is important to understand the effect of the non-zero spins of the scattering particles N and J/ψ in the non-interacting limit. Even in the continuum infinite-volume theory, several linearly-independent eigenstates have degenerate energy in a given channel J P . Those are eigenstates with different partial waves and total spins S that render certain J P . For N J/ψ channel with J P = 1 2 + , for example, there are two linearly independent physical states that have the same energy (1) in the non-interacting limit and exist for relative momenta p > 0 (p = 0 does not render partial waves > 0). Our aim is to find all linearlyindependent eigenstates and we will indeed confirm this pair of nearly degenerate eigenstates 2 . The effect of finite lattice implies that several channels J P contribute to a given lattice irreducible representation (irrep), which further enhances the number of degenerate eigenstates N (p)M (−p) at given p. Up to six degenerate eigenstates are expected at nonzero relative momenta for the N J/ψ channel and up to two are expected for N η c . One of the 2 They appear in irrep G + 1 that contains J P = 1 2 + . There are two levels with p = 1 in G + 1 of Fig. 4.
challenges is to extract all linearly independent eigenstates and we indeed establish all of them in our simulation, as will be evidenced from high degeneracies in the Figs. 4 and 7.
The existence of a pentaquark resonance P c would modify the eigen-energies of the nucleon-charmonium system with respect to the non-interacting case discussed above. We address the fate of P c in one-channel approximation, i.e. when a given pentaquark resonance couples only to N J/ψ or only to N η c , while it is decoupled from other two-hadron channels. If such a narrow resonance P c (4450) exist with a given spin-parity J P , the Lüscher formalism predicts one additional eigenstate with respect to the non-interacting case in the irreducible representations that contains this J P . This extra state is expected to have an energy approximately M Pc ±Γ Pc , while nearby states can get shifted in energy. Analogous expectation applies for the broader P c (4380), which is experimentally found in another spin-parity channel, while energy shifts would be larger in this case.
One of our goals is therefore to look for possible extra eigen-states in the spectrum, which could signal the presence of pentaquarks within one-channel approximation. This underlines the importance of finding a complete spectrum of eigenstates with all expected nearlydegenerate energy levels. Only then one could relate an extra level to a possible signature for P c . Spectra in all lattice irreducible representations are considered in our work and those include also J P = 3/2 ± , 5/2 ± that are particularly relevant for pentaquark searches.
The paper is organized as follows. The variational approach to extracting eigen-energies is reviewed II. The operators for single hadrons and for nucleon-meson systems are discussed in Sections III and IV, underlying complications arising from the fact that scattering particles N and J/ψ carry spin. Construction of two-hadron correlators is described in Sections V, followed by details of the lattice simulation in Section VI. Section VII presents the results and compares them to previous simulations. We end by Conclusions and Outlook.

II. EXTRACTING EIGEN-ENERGIES
The aim of this work is to extract energies of the eigenstates for N J/ψ and N η c systems, as well as for the relevant single hadrons N , J/ψ and η c . The eigen-energies of a certain system are obtained by computing N × N correlation matrices C on the lattice where O i=1,..,N denote interpolators that have quantum numbers of the desired system. In order to reliably extract the excited states, it is favorable to employ a large number of relevant interpolators, like for example in [39].
The information on the eigenstates is obtained by inserting complete sum of eigenstates |n . We employ the widely used variational approach by solving the generalized eigenvalue problem (GEVP) [40,41] where t 0 = 2 is employed. The eigen-energies E n are extracted from the resulting eigenvalues by one-exponential fits in the plateau region.

III. SINGLE-HADRON OPERATORS
In order to determine the non-interacting energies of the nucleon-charmonium system, we separately compute nucleon and charmonium correlators with momenta p 2 = 0, 1, 2 (in units of (2π/L) 2 ). Three standard nucleon operators and two standard operators for each charmonia are employed at each momentum p.

IV. TWO-HADRON OPERATORS AND EXPECTED DEGENERACIES IN THE NON-INTERACTING CASE
The operators for two-hadron system are needed to compute eigen-energies of these systems from the correlation functions. We consider only the system with total momentum zero since parity P is a good quantum number in this case, which simplifies the study. Operators are therefore of the form where each hadron is projected to a definite momentum p 2 = 0, 1, 2. All possible quantum numbers J P are considered since J P of P c pentaquarks are not reliably established from experiment yet. Therefore we consider all six irreducible representations Operators in partial-wave method The operators which have good total angular momentum J, its third component m J , angular momentum and total spin of particles S in continuum are referred as to partial-wave (PW) operators [42,43] The operator has good parity P = P 1 P 2 (−1) . The N η c system carries S = 1 2 , while N Jψ can have S = 1 2 or 3 2 . A given combination of J P can be obtained with multiple combinations of (S, ) and such channels can be coupled in the continuum infinite-volume (e.g. example in Eq. 2).
On the finite cubic lattice, the operators O J,m J form a reducible representation with respect to the lattice group O h . One has to employ operators which transform according to irreducible representations Γ P . Those are listed in Table I, where several J P contribute to a given irrep Γ P . Therefore partial wave operators (8) where the subduction coefficients S J,m J Γ,r are given in the Appendix of [44]. All these operators for nucleon-vector and nucleon-pseudoscalar systems with p 2 = 0, 1 are explicitly written in Appendix C of [43].
A simple example of the N J/ψ operator for p = 0 that transforms according to Let us explicitly present also a more non-trivial example, where N J/ψ with p 2 = 1 transform according to H − . Here the relations (8,9) lead to a number of linearlydependent operators , We find a basis of three linearly independent operators, which are chosen to be the operators in the first line.
Those are listed in Table IX and employed for the correlation matrices for N J/ψ system. The first operator, for example, has the form and these are not explicitly incorporated in the correlation matrix.
Along similar lines, we found linearly independent partial-wave operators for all irreducible representations and relative momenta. Our choices of linearlyindependent sets (J, , S) are given in appendix B 3 . The number of linearly independent operator-types for each |p| is provided in the Table II. For each of the operatortypes listed in Appendix B we in fact implement six operators: three choices for the nucleon (Eq. 6) times two choices for charmonium (Eq. 7). The final number of operators used in the computation of the correlation matrix for each irrep is given in Table III. irrep N J/ψ 6 × 6 5 × 6 4 × 6 4 × 6 10 × 6 9 × 6 TABLE III. Number N of interpolators used to compute the N × N correlation matrix C (3) for a given row of irrep. The study of N J/ψ scattering in H − , for example, uses N = 10 × 6 = 60 interpolators (10 operator types with 3 nucleon and 2 vector choices for each), so correlation matrix of size 60 × 60 is computed.
We note that the number of linearly independent operator-types in Table II agrees with the number of linearly-independent operators-types obtained using the projection method [43], helicity method [43] and Clebsch-Gordan decomposition in [45,46].

B. Expected number of degenerate states in the non-interacting limit
In the non-interacting limit, one expects several degenerate N (p)M (−p) eigenstates for most of J P (or irreps) and relative momenta p > 0. Let us try to understand this first in the continuum limit of an infinite volume. Different combinations of ( , S) lead to a given J P (| − S| ≤ J ≤ | + S|) due to the non-zero spins of the scattering particles. An example is given in Eq. (2), while more examples can be deduced from Tables VIII and IX. The linearly independent combinations ( , S) represent linearly independent eigenstates, so each of them should feature as an independent eigenstate in the spectrum.
This will remain true on the lattice on a finite volume and the corresponding discrete symmetry group. However, in this case also different spins J P can contribute to a given irrep Γ P as listed in Table I. Linearly independent combinations (J, , S) that subduce to given irrep Γ P now present linearly-independent eigenstates. These are listed in appendix B and their number is summarized in Table II for each Γ P and relative momentum p 2 . In the non-interacting limit, each of those linearly independent eigenstates should appear in the spectrum. Therefore the expected number of degenerate eigenstates for each row is given in Table II. These degeneracies can be lifted by the presence of the interactions.

V. CONSTRUCTION OF TWO-HADRON CORRELATORS
In the one-channel approximation there is no contraction connecting charmonium meson and nucleon interpolator, as shown in Figure 2. Therefore single hadron correlation functions can be simulated separately and then later combined to the two-hadron correlation function. All two-hadron correlators in our study can be expressed in terms of which are pre-computed for all combinations of p , p = 0, 1, 2, i, i = x, y, z and m s , m s = 1/2, −1/2, as well as for all configurations c. We omit the Wick contractions with the disconnected charm contributions in the study of pentaquark as well as charmonium systems. These contractions induce charmonium decays to the light hadrons and have been omitted in most of the previous lattice studies related to charmonium-like systems. Let us give an example of two-hadron correlator corresponding to operator O G − 1 ,r=1 (Eq. 10) at the sink and conjugate operator at the source. Given that there is no contractions connecting J/ψ and N interpolators (Fig. 2), the twohadron correlator can be expressed as The two-hadron correlation function is a sum of products of pre-calculated single-hadron correlators. Products of nucleon and vector correlators, that are calculated on individual configurations c, are averaged over configurations, as indicated by the large parentheses ... . All N × N correlation functions were constructed in a similar manner, where the number of interpolators n is given in the Table III. We calculated the correlation functions for all rows of irreps and then averaged them over the rows to gain better statistical accuracy.  The simulation is performed on the N f = 2 ensemble with parameters listed in table IV, that was generated in context of the work [47,48]. Wilson Clover action was used for light u/d quarks, while Fermilab approach [49] was employe for the charmed quarks. Strange quark is not dynamical in this simulation; we expect that dynamical strange quark is not crucial for this channel which does not involve strange valence quarks. Decay channels to hadrons with strange quarks also do not feature as important channels in the phenomenological studies of P c . Further details about the ensembles and treatment of the charm quark may be found in [47,48,50,51]. Periodic boundary conditions in space and anti-periodic boundary conditions in time are employed for fermions. The small spatial size (in comparison to other currently used lattices) L ≈ 2 fm of our lattice has a practical advantage since only levels up to relative momentum p 2 ≤ 2 (2π/L) 2 have to be extracted to cover the P c resonance region. A larger spatial size would imply denser energies (1) and require extracting a larger number of energy states, which would be more challenging in view of degeneracies appearing in the spectra.

VI. LATTICE SETUP
The quark fields are smeared according to the full distillation method [52], which allows the calculation of all the necessary correlation matrices. The smearing of the light quarks in the nucleon is based on N v = 48 lowest Laplace eigenvectors v k where k = 1, .., N v . The smearing of the charmed quarks in the charmonium is based on 96 eigenvectors. Perambulator τ k k (t , t) represents the quark propagator from the Laplace eigenvector k at time t to the Laplace eigenvector k at time t . All lightquark and charmed-quark perambulators are precomputed and saved for all distillation indices k, k = 1, .., N v and t, t = 1, .., N T . From those we compute charmonium two-point functions and nucleon two-point functions (13) needed for this study. Charmonium two-point function is given as a product to two perambulators in Eq. (11) of [52], while nucleon two-point function is a product of three perambulators in Eq. (19) of [52]. The numerical cost of computing these correlators from given perambulators comes from summing over the distillation indices (k, k , ..). The cost is considerable for the nucleon twopoint functions, which have to be evaluated for all combinations of momenta and polarizations at source and sink (13).
In the Fermilab approach to the discretization of charm quarks [49], one can attribute physical significance only to the quantities where the rest energy of the charm quarks cancels. This cancels, for example, in the difference between the energy of theccuud system and the charmoniumcc. We choose the mass of the spin-averaged 1S charmonium as the reference energy since this mass is least prone to the discretization errors. We will compare the lattice energies E lat − m lat cc,sa to the experimental masses m exp Pc − m exp cc,sa .

VII. RESULTS
The resulting eigen-energies of the single hadrons (N , J/ψ, η c ) and two hadron systems (N J/ψ, N η c ) are presented in this Section. They are obtained from the correlated one-exponential fits (5) and the errors are determined using jack-knife method.

A. Single hadron energies
The single-hadron energies E H (p) on the employed lattice are needed to determine the non-interacting energies of two hadrons (1). The fitted energies for various momenta are collected in Table V. They arise from 2 × 2 correlation matrices where both meson operators (Eq. 7) are used, while the first and the third operators (Eq. 6) are used for the nucleon. Results for E H on each resampled jack-knife will be used to determine E n.i and energy-shifts ∆E (1) in the next subsection.    (5) for the lowest four eigenstates of the N J/ψ system in G + 2 irrep. This gives the eigen-energy En in the plateau region. We observe all N (p)J/ψ(−p) eigenstates, expected in the non-interacting limit: this number is 0, 1 and 3 states for p 2 = 0, 1 and 2, respectively (Table  II). No additional eigenstate is found. The non-interacting energies (1) are indicated by the dashed lines. Color-coding of different eigenstates is arbitrary.

B. Two hadron energies
In the non-interacting limit, the two-hadron energies are equal to the sum of energies of the individual hadrons (Eq. 1) given in Table V. Resonances in LQCD manifest themselves by the non-zero energy shifts with respect to the non-interacting ones [41,[53][54][55][56]. Therefore we calculate the energy spectrum for all irreducible representations Γ P of the lattice group O h , which contain contributions from various J P as indicated in Table I. Channels with possible P c candidates are contained in the irreps G ± 2 for J P = 5 2 ± and in irreps H ± for J P = 5 2 ± or 3 2 ± .
Eigen energies are calculated from the correlation matrices as described in Section V. The size of the calculated correlation matrices for all irreps is displayed in Table III. These large correlation matrices render rather noisy eigenvalues, therefore we restricted our analysis to a somewhat smaller subset in Table VI, where each operator type (listed in Appendix B) is represented by two meson operators (7) and the first nucleon operator from (6). All resulting eigen-energies are obtained from the correlated one-exponential fits of the eigenvalues in the time range t = [7, 10] 4 . The energies of theccuud system 4 Two-exponential fits starting from earlier time-slices lead to com-will be provided with respect to the spin-averaged charmonium mass (m ηc + 3m J/ψ )/4, where the rest energies of the valence charm quark-pair cancels.

N J/ψ channel
The charmed P c resonances were observed as two peaks in the spectrum of proton−J/ψ invariant masses. Therefore, this channel could be the promising one in which to look for the charmed pentaquarks.
The final eigen-energies of N J/ψ system are presented in Figure 4 and in Table XI for all six irreducible representations. An example of the effective energies in irrep G + 2 is given in Figure 3. The energies of the states that are not plotted are higher and have typically (much) larger error-bars. Let us compare this spectrum separately with the non-interacting limit and with a scenario featuring P c .
• First we compare the resulting energies in Fig. 4 with the expectation from the non-interacting limit.
The non-interacting energies of N (p)J/ψ(−p) (Eq. 1) are given by the horizontal dashed lines: these energies are obtained as a sum separate lattice energies of N and J/ψ given in Table V. The observed energies in Fig. 4 agree with the non-interacting energies within errors. The number of the expected degenerate eigenstates in the non-interacting limit is listed in Table II and in Fig. 4. Such multiplicities of levels arise due to the non-zero spin of the scattering particles N and J/ψ. The energies in Fig. 4 indicate that we observe exactly the same pattern of degeneracy as in the non-interacting limit. The use of carefully-constructed interpolators [42,43] were crucial for this. We find no extra eigenstates in addition to those expected in the non-interacting limit. So, the observed lattice spectrum is rather close to the non-interacting case.
• Next we compare the energies with the analytic prediction based on the existence of P c (4450) or P c (4380) in Fig. 5. The aim is to explore the one-channel scenario where P c is coupled to N J/ψ and decoupled from other two-hadron decay channels. Such scenario renders an additional eigenstates near M Pc ± Γ Pc with respect to the noninteracting case (Fig. 5). The favored channels J P = 5/2 + for P c (4450) and 3/2 − for P c (4380) [2] are considered 5 . We assume that P c resides only in a single partial wave ( , S) and that there is no interaction in the other channels. The Breit-Wignertype dependence of the phase shift is employed, patible results. In few cases these are less stable and have larger errors. 5 The same conclusions apply for other possible J P listed in the introduction.
FIG. 4. Energies of N J/ψ eigenstates in the one-channel approximation for all lattice irreducible representations. The quantum numbers J P that contribute to each irrep are listed on the top. Each box represents one eigenstate. The centre of the box represents it energy En, while height represents to 2σE n . Ten lowest eigen-energies are, for example, shown in the irreducible representation H − . The number of the observed near-degenerate states agrees with the expected number of states in the non interacting limit. This number is given in Table II and indicated in the plot; it arises due to the non-zero spins of the nucleon and J/ψ. Dashed lines represent non-interacting energies EN (p) + E J/ψ (−p) for different value of relative momentum p: black for p 2 = 0, red for p 2 = 1 and blue for p 2 = 2 where p 2 is given in units of (2π/L) 2 . The dash-dotted (green and turquoise) lines correspond to experimental masses of Pc states. The observed spectrum shows no significant energy shifts or additional eigenstates. Energies are presented with a respect to spin-averaged charmonium mass (mη c + 3m J/ψ )/4 (16). Color-coding of the eigen-energies is arbitrary. with the resonance parameters M Pc and Γ Pc taken from experiment [2] The relation between the eigen-energies E and the phase shift δ for the scattering of particles with arbitrary spin was derived in [55]. This reduces to the well-known Lüschers relation [41,53,54,56] for a resonance that appears only in the channel (J P , , S) (17). The analytic predictions for the energies (orange circles) are obtained by solving Eq. (17) for E, assuming the continuum dispersion relation E(p) = (p 2 + m 2 N ) 1/2 + (p 2 + m 2 J/ψ ) 1/2 and employing m N,J/ψ determined from the lattice. Other levels in a given irrep remain intact and have noninteracting energies (brown circles). The predictions in Fig. 5 show all possible irreps and choices of partial waves ≤ 2 where the P c resonance with given J P could reside.
The scenario with a P c resonance in Fig. 5 features an additional eigenstate (with respect to the non-interacting case) in the energy range roughly M Pc ± Γ Pc , while some of the other levels near the resonance region get shifted. So, this scenario predicts one eigen-state more (with respect to Table II) in the explored energy region within the corresponding irreducible representation. We do not observe such an additional eigenstate, so such scenario featuring P c is not supported by our lattice data.
In summary, our lattice spectra in Fig. 4 are roughly in agreement with the predictions for almost noninteracting N and J/ψ. These spectra do not support the scenario where a P c resonance couples only to N/ψ decay channel and is decoupled from other channels. These results indicate that the existence of P c resonance within a one-channel N J/ψ scattering is not favored in QCD. This might suggest that the strong coupling between the N J/ψ with other channels might be responsible for the existence of the P c resonances in experiment. Future lattice simulations of the coupled-channel scattering will be needed to confirm or refute this hypothesis.

N ηc channel
The final eigenenergies of N η c system are presented in Figure 7 and Table X for all six irreducible representa- 6 The resonance mass in (17)  tions. An example of effective energies in irrep H − is given in Fig. 6. We observe all N (p)ηc(−p) eigenstates expected in the noninteracting limit: this number is 0, 1 and 2 states for p 2 = 0, 1 and 2, respectively (Table II). No additional eigenstate is found. The non-interacting energies (1) are indicated by the dashed lines. Color-coding of eigen-energies is arbitrary.
• First we compare the resulting energies in Fig. 7 with the expectation from the non-interacting limit.
The non-interacting energies of N (p)η c (−p) (Eq. 1) are indicated by the horizontal dashed lines. We find that the observed energies agree with the noninteracting N η c within sizable errors of our calculation. The number of expected degenerate eigenstates in the non-interacting limit is listed in Table II and in Fig. 7. The degree of degeneracy is smaller than in N J/ψ case since η c does not carry spin. We observe exactly the same pattern of degeneracy, while no additional eigenstates are found.
• The scenario featuring P c , coupled dominantly to N η c and largely decoupled from other channels 7 , would predict an additional eigenstate in the energy range roughly M Pc ± Γ Pc . This is based on an analogous argument, presented in more detail for N J/ψ system. Such an additional eigenstate is not observed in Fig. 7, so this scenario is not favored by our lattice results.
In summary, Figure 7 shows that no additional eigenstate or significant energy shift is observed with respect to non-interacting N and η c . We conclude that there is no strong indication for a P c resonance in one-channel approximation for N η c scattering.

C. Comparison with previous lattice simulations
Finally, we compare our conclusions with previous lattice simulations of N J/ψ and N η c systems.  Table II and indicated in the plot; it arises due to the non-zero spin of the nucleon. This number is smaller for N ηc than for N J/ψ since J/ψ carries spin one, while ηc is spinless. Dashed lines represent non-interacting energies EN (p) + Eη c (−p) (Eq. 1) for different values of the relative momenta p: black for p 2 = 0, red for p 2 = 1 and blue for p 2 = 2, where p 2 is given in units of (2π/L) 2 . The dash-dotted (green and turquoise) lines correspond to experimental masses of Pc states. The observed spectrum shows no significant energy shifts or additional eigenstates. Energies are presented with respect to spin-averaged charmonium mass (mη c + 3m J/ψ )/4 (16).
The s-wave interaction between a nucleon and J/ψ(η c ) was studied using the Lüscher formalism in dynamical [31] and quenched [32] QCD. All calculated scattering lengths a 0 , are consistent with zero within 1 or 2 sigma, implying very small attractive interaction. As a 0 ∝ ∆E, these results agree with our result ∆E ≈ 0 within errors.
The NPLQCD collaboration observed the energy shift ∆E ≈ −20 MeV of N (0)η c (0) ground state for very heavy m π 800 MeV in J = 1 2 − (G − 1 ) channel, which was almost independent of the volume L 3.4 − 6.7 fm [33]. We can neither confirm neither refute such an energyshift given the errors of our present calculation and different m π . This study concludes that N η c bound state 20 MeV below threshold exists at such a heavy pion mass. Further accurate studies are needed to explore possible existence of such a bound state at physical quark masses. If this state is theoretically confirmed, looking for its experimental signatures will be of prime interest.
Two studies considered potential between N and M (M = J/ψ or η c ) in s-wave as a function of distance withing the HALQCD method in a dynamical [34] and a quenched [35,36] simulation. The light-quark mass was larger than physical with the nucleon mass m N 1.8 GeV in [34] and m π = 640 − 870 MeV in [35]. They find weakly-attractive interaction near threshold in three channels explored: J P = 1 2 − , 3 2 − for N J/ψ and 1 2 − for N η c . The resulting interaction was not strong enough to form bound states or resonances, but the most interesting experimental region 4.3−4.5 GeV was not explored. The absence of a very pronounced interactions in this system agrees with our conclusions.
Finally, the study [37] of hadroquarkonium picture considered the staticcc (m c → ∞) as function of distance betweenc and c in presence of the nucleon. The N f = 2 + 1 CLS ensemble at m π 223 MeV was employed. The shift of the potential due to the presence of the nucleon was extracted by an impressive precision. This shifts was found to be down only by a few MeV. Such a shift is compatible with our results given the uncertainties on our eigen-energies.
Lattice results [31,32,35,36] were used for determination of parameters in effective field theory [38]. No signs of quarkonium-nucleon bound state is found within range of applicability for described EFT.

VIII. CONCLUSIONS AND OUTLOOK
We perform a N f = 2 lattice QCD simulation of N J/ψ and N η c scattering in the one-channel approximation, where N denotes a proton or a neutron. This is the first study that reaches the energies, where the charmed pentaquarks P c resonances were observed in N J/ψ decay by the LHCb experiment. The resulting energies of eigenstates in Figures 4 and 7 are compared to the analytic predictions of (i) a scenario with non-interacting nucleoncharmonium system and (ii) a scenario featuring a P c resonance coupled to a single channel. The non-interacting spectrum E = E N (p) + Ec c (−p) with p = n2π/L is discrete due to periodic boundary conditions on the lattice of finite size L. We find that the extracted lattice spectra is consistent with the prediction of an almost noninteracting nucleon-charmonium system within errors of our calculation. The scenario based on a Breit-Wignertype P c resonance, coupled solely to N J/ψ or to N η c , is not supported by our lattice data. The results indicate that the existence of P c resonance within a one-channel scattering is not favored in QCD. This might suggest that the strong coupling between the N J/ψ with other two-hadron channels might be responsible for the existence of the P c resonances in experiment. Future lattice simulations of coupled-channel scattering are needed to investigate this hypothesis.
One of the challenges in extracting the eigenstates of N J/ψ system in the current simulation is related to the high number of almost-degenerate eigenstates. This degeneracy arises due to the non-zero spins of N and J/ψ, since a number of different partial waves ( , S) can couple to a certain channel J P . Furthermore, several J P contribute to a certain lattice irreducible representation due to the reduced symmetry on the lattice. As a result, the non-interacting scenario predicts up to six degenerate linearly-independent N (p)J/ψ(−p) eigenstates with the relative momentum p 2 = 2( 2π L ) 2 for a given row and a given lattice irreducible representation. We establish all such eigenstates with p 2 ≤ 2( 2π L ) 2 , together with the pattern of degeneracies expected from the noninteracting case. The use of carefully-constructed interpolators [42,43] were vital for this.
Our work is only the first step towards exploring the dynamics of the charmed pentaquark channels by means of the lattice QCD simulations. Although P c uudcc resonances were experimentally observed so far only in the proton-J/ψ decay, they are allowed to strongly decay to a nucleon and a charmonium as well as to a charmed meson and a charmed baryon. The coupled-channel lattice study of the relevant channels is a challenging task left for the future simulations. These will contain also the Wick contractions connecting a meson and a baryon, which are absent for the case of a nucleon-charmonium system considered here. The number of eigenstates in a given irreducible representation will become even denser than in the present study. The extraction of the small non-zero energy shifts will require a particular effort in reducing the statistical error, especially for the underlying baryonic component at zero and non-zero relative momenta. The study of the systems with the non-zero total momentum will render additional information on the scattering matrices, but poses additional challenges since spin J and parity P are no longer good quantum numbers, even in the continuum. Therefore, understating the realistic P c resonances through the energies of eigenstates by means of a rigorous Lüscher-type approach brings considerable challenges. In light of this, it would be valuable to explore if there are any other ways to investigate these interesting systems by means of the first-principle lattice QCD.

Acknowledgments
We are grateful to M. Padmanath for valuable discussions and help concerning the cross-checks for the nucleon. We would like to kindly thank D. Mohler and C.B. Lang for allowing us to use the quark perambulators, generated during our previous joint projects. We acknowledge A. Hasenfratz for sharing with us the gauge configurations employed here. We thank L. Leskovec for sharing GEVP code for an independent cross-check. This work was supported by Research Agency ARRS (research core funding No. P1-0035 and No. J1-8137) and DFG grant No. SFB/TRR 55.
Note added. -Recently, LHCb reported a discovery of a new pentaquark state P + c (4312) that is also observed in the pJ/ψ invariant mass [57]. All the discovered narrow P c lie near Σ + cD 0( * ) threshold, indicating that coupling of pJ/ψ to this channel might be important for giving rise to P c in experiment. This might provide a possible explanation as to why P c resonances are not observed in our lattice study of pJ/ψ channel in the approximation where it is decoupled from all other channels. Recent LHCb results are in line with our conclusion that the coupling of pJ/ψ with other channels might be crucial for the existence of P c resonances.