Variational analysis of low-lying states in supersymmetric Yang-Mills theory

We have calculated the masses of bound states numerically in N = 1 supersymmetric Yang-Mills theory with gauge group SU(2). Using the suitably optimised variational method with an operator basis consisting of smeared Wilson loops and mesonic operators, we are able to obtain the masses of the ground states and first excited states in the scalar, pseudoscalar and spin-1/2 sectors. Extrapolated to the continuum limit, the corresponding particles appear to be approximately mass degenerate and to fit into the predicted chiral supermultiplets. The extended operator basis including both glueball-like and mesonic operators leads to improved results compared to earlier studies, and moreover allows us to investigate the mixing content of the physical states, which we compare to predictions in the literature.


Introduction
The N = 1 supersymmetric Yang-Mills (SYM) theory is the minimal supersymmetric extension of the pure gauge sector of QCD, describing the strong interactions of gluons and their fermionic superpartners, the gluinos. Unbroken supersymmetry requires the physical bound states of N = 1 SYM to form supermultiplets of degenerate masses. Based on supersymmetric effective actions, the lowest-lying supermultiplet has originally been proposed in ref. [1] to be a chiral supermultiplet formed out of a scalar and a pseudoscalar meson and a spin-½ bound state of gluons and gluinos, called gluino-glue. The original analysis has been extended in refs. [2,3] to include also glueball operators which can mix with the meson operators carrying the same quantum numbers. The scalar and pseudoscalar bound states are therefore created by linear combinations of meson-like and glueball-like operators. The mass hierarchy and the amount of mixing are difficult to predict theoretically from the effective Lagrangian and are in general unconstrained by supersymmetry.

JHEP04(2019)150
Numerical simulations of N = 1 SYM on a space-time lattice are required to verify and extend the analytical predictions about the low energy effective theory. The study of the lightest bound states of N = 1 SYM with SU(2) gauge symmetry has been the main subject of several publications by our collaboration [4][5][6][7][8]. The results reveal the expected degeneracy of the members of the chiral multiplet in the continuum limit. In particular, the scalar glueball and the scalar meson masses, extrapolated to the continuum limit, have compatible values, which hints at mixing of the states in this channel.
In the present work we have optimised our techniques to extract the masses with respect to our previous work. This allows for the first time to investigate the mixing properties of the two lowest supermultiplets from the lattice gauge ensembles we have generated. In this article we explain these methods and their optimisations. We show, in particular, that the construction of an enlarged correlation matrix including both, the meson and glueball operators, is crucial to extract the physical states in the scalar channel. Using these techniques systematic errors stemming from excited state contributions can be controlled much better. We reanalyse the gauge configurations of previous investigations yielding more reliable and more precise results for the lightest masses. In addition, we are able to extract the masses of excited states in the scalar, pseudoscalar, and fermionic channel. In this way we provide first results concerning the possible formation of a supermultiplet of excited states on the lattice.
The formation of supermultiplets at the level of the excited states is strong evidence for the fact that the unavoidable breaking of supersymmetry on the lattice can be kept under control. The states with higher masses are affected stronger by lattice artefacts, and hence the predicted degeneracy of these states is more sensitive to the breaking of supersymmetry on the lattice. The present work represents the first step towards an investigation of the spectrum of excited states of N = 1 SYM. An alternative determination of supersymmetry violations beyond the ground state level would be the weighted sum of all energy differences in terms of the derivative of the Witten index, see [9].
In addition to our results on excited states, we confirm our earlier results concerning the mixing of glueball and meson operators in the ground state, and we are able to determine the glueball and the meson contents. Eventually this might lead to a better understanding of the conjectures concerning the low energy multiplets of the theory.
In this work we have optimised our methods for SYM with gauge group SU (2). The techniques are, however, also suitable for other theories. Currently we are applying the same methods successfully to SYM with gauge group SU(3), for which the results will be published soon. Our techniques are also applicable to lattice QCD, where similar investigations are being done.
It should be noted that our studies, and corresponding investigations in QCD as well, are quite challenging due to the noisy signals from glueball operators and the disconnected meson contributions. This paper is organised as follows: N = 1 SYM in the continuum and on the lattice is introduced in the next section. The current paper is particularly focused on the technical aspects, which are explained in sections 3, 4, and 5. The physical results, the improved signal for the ground states, the mixing, and the excited states, are finally presented in sections 6 and 7.

JHEP04(2019)150
2 Supersymmetric Yang-Mills theory on the lattice The Lagrangian L of the N = 1 SYM in Euclidean space-time, is similar to the one of one-flavour QCD. The main difference is that supersymmetry requires the gluino field λ to be a Majorana fermion field in the adjoint representation of the gauge group. Correspondingly the gauge covariant derivative is the adjoint one, given by (D µ λ) a = ∂ µ λ a + g f abc A b µ λ c . The gluino mass term breaks supersymmetry softly and the renormalised gluino mass must be tuned to zero to recover the full supersymmetry. In this work we focus on the theory with gauge group SU (2). The light mesons in the particle spectrum of SYM are similar to the mesons in QCD, apart from the difference that the constituent fermions are gluinos and not quarks. Due to this similarity, the scalar meson is called a-f 0 and the pseudoscalar meson is called a-η where the "a" indicates the adjoint representation. The mesons are expected to mix with the corresponding glueballs with the same quantum numbers. In addition to these particles, the low-energy spectrum contains a spin-½ bound state of gluons and gluinos, called gluino-glue, which has no analogy in QCD.
On the lattice, the fermion part of the action is discretised using the Wilson-Dirac operator D W , which effectively removes the doublers from the physical spectrum, but at the same time explicitly breaks chiral symmetry. As shown in ref. [10], chiral symmetry and supersymmetry can be recovered simultaneously in the continuum limit by tuning of the bare gluino mass. This tuning can be achieved in different ways. One way is to determine the point where the renormalised gluino mass vanishes from the supersymmetric Ward identities [11]. Another way is to employ the adjoint pion mass, which is defined in a partially quenched approach [12], and to extrapolate to the point where it vanishes. For non-zero lattice spacings, lattice artefacts introduce a small mismatch between these determinations. Because the adjoint pion mass can be measured quite precisely and effectively, we use this quantity to define our extrapolations to the chiral point. This corresponds to the point where supersymmetry should be restored in the continuum limit. N = 1 SYM is asymptotically free, and the running of the strong coupling has been investigated non-perturbatively in ref. [13]. The extrapolation to zero lattice spacing a → 0 corresponds to the weak coupling limit g → 0. The lattice spacing a used for the extrapolation to the continuum limit is measured in terms of the infrared scale w 0 defined from the gradient flow [8,14,15]. In order to suppress lattice discretisation effects, we use a tree level Symanzik improved action together with stout smearing on the gauge links in the Wilson-Dirac operator, see [6] for further details.
The production of the gauge field configurations is performed with the two-step polynomial hybrid Monte Carlo (PHMC) algorithm [4,16]. A mild sign problem of our lattice formulation arises from the integration of the Majorana fermions [17] and is corrected by reweighting.

JHEP04(2019)150
3 Techniques for the determination of light states In order to gain insights into the structure of the low energy spectrum of SYM, consisting of gluino-glue particles as well as mixtures of mesons and glueballs, we determine their masses as well as the glueball and meson content in the scalar and the pseudoscalar states. These quantities can be reliably extracted from the gauge field ensembles using the well-known variational method for which we have systematically optimised the required techniques and parameters.
Our approach is to use a set of basic interpolators of the physical states to which we then apply smearing techniques to create the operator basis used in the variational method. In sections 3 to 5 we explain in detail the techniques and the tunings used within this approach.
The mesons in SYM are all flavour-singlet ones, which require techniques for the measurement of all-to-all propagators. We have found that a preconditioned stochastic estimation is optimal for our lattice volumes. This is explained in section 3.3.

Computation of the masses using the variational method
A set of interpolating operators O i must be chosen in order to extract information about physical states within the variational method. In principle the only requirement is that the operators O i have the same quantum numbers under the lattice symmetry group as the physical state of interest. In practice the choice of this variational basis has a crucial influence on the precision of the results, and therefore we explain in detail the optimisation of these operators. To extract the masses of physical states using the variational method, the correlation matrix is built from the time-slice correlators of the set The solution of the generalised eigenvalue problem (GEVP) associated with the correlation matrix C(t), provides generalised eigenvalues λ (n) and their corresponding eigenvectors v (n) . In ref. [18] it has been shown that the generalised eigenvalues λ (n) (t, t 0 ) satisfy where m n are the masses of the physical states. In order to suppress contributions of higher excitations, it has been suggested in ref. [19] to use t 0 ≥ t/2. This, however, leads to very noisy eigenvalues λ (n) for the correlators considered here and we therefore use t 0 = 0 or t 0 = 1 when necessary. A first estimate can be obtained from the effective mass . A more precise mass estimation is obtained by fitting the eigenvalues λ (n) (t, t 0 ) to an exponential function of t.
The set of interpolating fields ideally should have large overlaps with the physical states of interest and small overlap with higher excited states, otherwise terms coming from higher order excitations might create large corrections to the expected exponential behaviour of the eigenvalues.

Interpolating operators
The basic interpolators for glueballs are built from gauge link loops that represent the spin and parity quantum number of the respective state. For the scalar glueball we use a sum of gauge plaquettes where P ij denotes a plaquette in the i-j plane. For the pseudoscalar glueball we use where the sum is over all rotations of the cubic group O h , and P C is the parity conjugate of the loop C, which is depicted in figure 1. The basic interpolating fields for the scalar mesons are When inserted into the correlation matrix, Wick contractions of these fields lead to connected and disconnected pieces where D −1 W (x, y) denotes the propagator from x to y (spin and group indices suppressed) and Γ represents 1 or γ 5 .
The interpolating field of the gluino-glue state is given by (3.10)

JHEP04(2019)150
The gluino-glue correlation matrix consists of an odd and an even part under time reversal Projections to the odd part C 1 (t) and to the even part C γ 4 (t) both provide valid correlators to be used in the GEVP. For the tuning of our methods we use C 1 (t) and for the final results we use a weighted sum of the results from both C 1 (t) and C γ 4 (t).

Stochastic estimator technique (SET)
The calculation of correlators including more than one gluino field in the interpolating operators requires the estimation of the fermion propagator from all to all lattice points. Since a complete determination of the inverse of the Dirac-Wilson operator D W is prohibitively expensive, we approximate D −1 W by using the stochastic estimator technique (SET) [20], which is improved by a truncated eigenmode approximation and even-odd preconditioning, see [21] for further details. The preconditioned approximation converges faster to the exact result and the numerical computation of the eigenspace of the preconditioned matrix is much faster than the one of the full matrix. For the inversion of the preconditioned Dirac matrix D p we use its Hermitean version, obtained by multiplication with γ 5 , since the singular value decomposition provides a better approximation.
The idea of SET is to solve the Dirac equation on a set of source noise vectors η i fulfilling the relation The entries of the vectors η i are Z 4 complex numbers of the form (±1 ± i)/ √ 2. The propagator D −1 p is then approximated by and D −1 p η i is calculated using a conjugate gradient solver. For the truncated eigenmode approximation the N E lowest eigenvalues λ i and eigenvectors |v i of the Hermitean operator γ 5 D p are computed. The truncated eigenvector approximation of the inverse matrix is The two approximations are easily combined: the noise vectors of SET are just projected to the subspace orthogonal to the one spanned by the lowest eigenvectors. In the end both contributions are summed, gaining also a speedup of the inversions due to the better condition number of the projected operator.

JHEP04(2019)150
For the case of plain Wilson fermions, investigated in the current study, the even-even part M ee of the Wilson-Dirac matrix is just the identity. More generally, the inverse of the block diagonal matrix M ee can be computed exactly. However, in order to simplify the smearing procedure, we have used a large number of stochastic sources to approximate the inverse of this matrix. These can be computed very efficiently.
We have optimised the parameters of this approximation in order to reduce the noise and to speed up the computations. The tuning of the number of noise vectors required for a reliable estimation of the disconnected piece is discussed in section 4.4. The eigenspace of the preconditioned matrix computed in the measurement of disconnected contributions is also used for a deflation of the inversions in other measurements, like the connected meson or the gluino-glue correlators. Concerning speedup, the optimisation is also machine dependent, e. g. our runs on KNL based machines require quite different parameters.

Smearing techniques
The different interpolating fields used in the variational method are constructed by applying smearing techniques to the basic interpolating fields defined in 3.2.
We use APE-smearing [22] on the gauge links in order to create smeared gluino-glue and smeared glueball operators. We use a smearing parameter APE = 0.4 for smearing the gluino-glue interpolators and APE = 0.5 for smearing the glueball interpolators.
For the construction of the fermion source we use gauge invariant Jacobi smearing [23] λ →λ = F λ with the smearing operator Here N J is the integer Jacobi smearing level, κ J is a Jacobi smearing coefficient, C J is a normalisation constant andŨ i (x) are APE smeared adjoint gauge links. The tuning of the parameters N J , κ J and C J and the smearing of the gauge links is explained in section 4. Smearing the fermions fields is equivalent to replacing the propagator D −1 W with the smeared propagator For the disconnected piece, this translates to smearing the source and sink vectors in equation (3.15) The connected piece is calculated using standard delta sources on which F † , D −1 W and F are subsequently applied:D

JHEP04(2019)150
4 Optimising the methods Significant sources of possible systematic errors are non-optimal parameters of the smearing techniques, the operator basis in the variational method, the stochastic estimation of disconnected propagator pieces, and the fitting methods and parameters. In this section we discuss their influence and optimisation.

Jacobi smearing
Aiming to achieve a good signal-to-noise ratio for the meson measurements, we systematically optimise the Jacobi smearing parameters κ J , C J , the smearing level N APE of the smeared gauge fieldsŨ and the Jacobi smearing levels N J . As explained in [24] there is a critical value of the parameter κ c J . For values smaller than κ c J , the smearing operator F converges in the limit N J → ∞, while for values larger than κ c J it diverges. In order to smear efficiently and at the same time to avoid large numerical errors we choose a value which is just above κ c J . With this choice the smearing radius is less than 6 lattices spacings up to smearing levels N J = 100 [24]. The normalisation factor C J is chosen such that the values of the correlation functions stay at the same order of magnitude for large and small smearing levels. We use In principle, the final results do not depend on the values for κ J and C J . The efficiency could, however, depend on the choice of the values of the smearing parameters and one could re-tune them for every ensemble. In our experience from different ensembles of SYM with gauge group SU(2) and SU(3), the dependency on these parameters is very mild, and it is not necessary to always find the optimal values. We have therefore kept fixed the values stated here.
For the tuning of the smearing levels we have used SU(2) gauge ensembles at β = 1.75 and κ = 0.14925. There are in total 6800 thermalised configurations of which we have measured every 64th for our tests. We use 40 stochastic estimators for the disconnected pieces.

Presmearing of the gauge fields
The application of Jacobi smearing defined with unsmeared gauge links introduces additional noise to the signal of the disconnected contribution to the a-f 0 meson, see figure 2. We have observed that this noise can be suppressed by using smeared gauge fields in the Jacobi smearing. We have tested APE and Stout smearing together with different smearing levels for the preparation of the gauge field. Our results indicate that N APE = 20, APE = 0.5 leads to a noise reduction for Jacobi smearing levels up to N J = 80. The results are not very sensitive to N APE , and values of N APE between 10 and 80 all feature sufficient noise suppression. Larger values for N APE lead again to a degradation of the signal, see N APE = 160 in figures 2, 3. Using Stout smearing instead of APE smearing didn't improve the noise suppression in our tests and we therefore stay with APE smearing.    Figure 3. Left: jackknife-error of the a-η correlator averaged over the interval t ∈ [3,12] plotted against the Jacobi smearing level, using different smearing levels for the gauge-fieldŨ in the smearing kernel (3.17). Right: effective mass of the a-η correlator using different levels of APE smeared gauge fields in the Jacobi smearing. The Jacobi smearing level is fixed to N J = 20. Using smeared gauge fields more effectively suppresses excited state contributions than Jacobi smearing with unsmeared gauge fields (blue). The correlator has been normalised to 1 at t = t 0 .
In the case of the a-η correlator, using smeared gauge fields in the Jacobi smearing has a different effect than on the a-f 0 correlator. It does not suppress the noise that is introduced by Jacobi smearing (see figure 3). However, when smeared gauge fields are used instead of unsmeared ones, Jacobi smearing more effectively suppresses excited state contributions to the correlator. Again, the results are not very sensitive to the exact value of N APE as long as it is in a suitable range 10 ≤ N APE ≤ 80.
Considering both a-f 0 and a-η correlators, the smearing level N APE = 20 is chosen, as it appears to be a good choice for suppressing noise and excited state contributions.

Variational operator basis
A variational basis for the GEVP can be built from different smearing levels of the interpolating operators. There is a trade-off between the cost required to build a large variational basis from many different smearing levels and the gain of new information from such operators. Each new smearing level requires additional inversions for the gluino-glue correlator and for the connected piece of the meson correlator. Therefore it is important to find those smearing levels which are most relevant for the extraction of the masses. For this purpose we consider the set of operators constructed from the smearing levels N Smearing ∈ {0, 5, 15, . . . , 95}, where Jacobi smearing is used for the meson interpolators and APE smearing is used for the glueball and gluino-glue interpolators. From these sets we have systematically chosen different subsets to analyse which smearing levels allow the most efficient estimation of the lowest two states in each sector.
For each number n of operators we have picked the following sets: 1. Small smearing levels: only the smallest n smearing levels of the full set are taken into account.
2. Large smearing levels: only the largest n smearing levels of the full set are taken into account.
3. Uniformly distributed: out of the full set we chose the smallest and the largest one and n − 2 additional smearing levels uniformly distributed in between.
4. Mid-large, uniformly distributed: out of the full set we chose a medium level (N Smearing = 35) and the largest one, and n − 2 additional smearing levels uniformly distributed between them.
To judge the quality of a mass determination with a chosen correlation matrix, the effective masses of the lowest two states at fixed t are used as estimators. Note that we use a rather small value of t here to demonstrate the suppression of excited states by choosing relevant operators. The fitting intervals for the final extraction of the masses are chosen more carefully (cf. section 6). The results from the described procedure are shown in figures 4 and 6 together with the best estimate for the mass obtained by fitting the eigenvalues using the full set of interpolating operators.
The results for the dependency of the effective masses on the smearing levels is similar for all three particles: smallest smearing levels obviously lead to unwanted contributions of higher excitations so that the effective masses at small t are significantly higher than the best estimate. Using only the largest smearing levels leads to the strongest suppression of these excited state contributions. The uniformly distributed subsets also suppress excited state contributions, but not quite as much as the large smearing levels. Note that the error of the effective masses is almost constant for the different combinations of the smearing levels. We conclude that it is optimal to use a set of rather large smearing levels for the mass estimations.

Number of stochastic estimators
With the optimised Jacobi smearing parameters we also tested the influence of the number of SET estimators on the error of the disconnected pieces. For this estimation we used every fourth configuration of our test ensemble. The aim is to find a value for the number of stochastic estimators where the error from the stochastic estimators is much smaller than the gauge noise. We find that higher smearing levels require less stochastic estimators, see figure 5, and the scalar correlator requires more stochastic estimators than the pseudoscalar one. At around 20 stochastic estimators the error is dominated by the gauge noise (the stochastic error is smaller than 15% of the gauge noise) for all tested smearing levels except for the unsmeared scalar correlator.

Extended variational basis
In general, an operator transforming according to a given irreducible representation of the lattice symmetry group has a non-zero overlap with all the eigenstates of the Hamiltonian with same quantum numbers. Mixing occurs if two operators share the same transformation properties, independently of their fermion or gluon field-content. In the case of the 0 ++ or 0 −+ channels the relevant operators can be constructed from glueball-like combinations of Wilson loops and meson-like operators of the formλΓλ. Therefore the variational basis (3.2) can be enlarged to include the most general mixing between glueball and meson operators. Taking both kinds of operators into account, the full correlation matrix has the following form Each entry of C(t) is a submatrix consisting of the correlators among different interpolating fields, where O gb stands for glueball-like operators and O meson stands for meson-like operators. Since the meson and the glueball operators are constructed from quite different components, it is expected that their mutual overlap is small. Using this larger correlation matrix we expect to obtain significantly improved signals.
Scalar 0 ++ channel. The results of our calculations show that the enlarged variational basis leads to a more effective suppression of excited states contributions than using meson or glueball operators alone. Even the minimal choice of using only one meson and one glueball operator appears to be sufficient to extract the masses of the lowest two states in the scalar channel, see figure 6. Using more than one glueball or more than one meson operator does not improve this estimation, but provides better access to the excited states. We therefore conclude that it is crucial to include both meson and glueball operators in the variational basis to reliably extract the ground state in the scalar channel.
Pseudoscalar 0 −+ channel. In this channel the estimation of masses does not improve when the enlarged basis also includes glueball operators. The lowest two states can be analysed sufficiently by using only meson operators, see figure 7, which means that the glueball operators do not mix significantly with these states. This is in accord with the observation that the entries in the off-diagonal blocks of the correlation matrix (5.1) are zero (within rather large errors).

Supermultiplets
Using the methods explained above we have reanalysed the low-lying spectrum using the gauge configurations of our previous investigations [5]. The parameters of the configurations, on which the present investigation is based, are compiled in table 1.
The main improvements compared to the previous analysis are the use of the variational method in all three channels including glueball and meson operators in the 0 ++ channel (in previous investigations the variational method was only applied to the glueball analyses), the inclusion of the additional gluino-glue correlator C γ4 (t), and the use of smeared gauge fields in the definition of Jacobi smearing. Due to these improvements, the contributions from higher states to the correlators and noise stemming from the Jacobi smearing are suppressed so that the effective mass plateaus can be estimated more reliably and smaller values of t min can be used as the lower boundary of the fitting intervals. Consequently, we not only obtained more reliable and more precise estimates of the ground state masses, but we were also able to extract the masses of the first excited states, see figure 8.
For this investigation the full set of smearing levels under consideration consists of the gluino-glue smearing levels N APE = (5, 15, . . . , 95) and the glueball smearing levels N APE = (4, 8, . . . , 64). Note, that in contrast to the analysis above, we used a more conservative choice of smearing levels for the meson operators, namely N J = (0, 3, 20, 40). It is preferable to use not more than four smearing levels in the 0 ++ channel since this assures that the eigenvalues of the GEVP are clearly separated in the region of the fit and additional operators do not reduce the contamination by excited states further. We therefore use combinations of two meson operators, e.g. N J = (20, 40) and two glueball operators N AP E = (24, 48) in this channel. As explained above, including the glueball operator in the GEVP for the pseudoscalar channel does not improve the results, thus we have analysed the pseudoscalar meson a-η and the glueball gb −+ separately. The fit intervals were determined at each ensemble individually based on the formation of effective mass plateaus and the stabilisation of the χ 2 /d.o.f of the fit. Using the mixed basis of glueball and meson operators allowed to include smaller time separations t min for the scalar channel compared to the a-η . For the ensemble β = 1.9, κ = 0.14415 e. g. we used t min = 4 for the ground state and t min = 3 for the excited state in the scalar channel whereas we used t min = 7 for ground and excited states of the a-η . The errors were obtained statistically using the standard Jackknife procedure together with binning in order to properly take into account autocorrelations.
Due to the increased precision, we observed that in the chiral extrapolation the gluinoglue mass deviated strongly from the assumed linear behaviour at the ensemble with the largest adjoint pion mass at β = 1.9, κ = 0.1433. This resulted in a χ 2 /d.o.f ≈ 7.5 of the linear fit to the chiral limit. We therefore neglected this data point for the gluino-glue. This also contributes to the fact that the updated result for the gluino-glue mass is smaller than the one from our previous investigation. We also tested to use a second order polynomial instead, which gave a consistent result, but a slightly larger error.
The results for the masses, extrapolated to the chiral and to the continuum limit, are collected in table 2. The ground state masses, extrapolated to the continuum limit, agree quite well with each other within errors. Similarly to the ground states, the respective excited states also seem to form a mass degenerate supermultiplet of approximately three times the ground state mass. Interestingly the 0 −+ -glueball ground state and the excited state of the a-η both have similar masses, but they appear as two independent states in the variational method which do not mix.
In the continuum limit, the two-particle decay threshold is well below the mass of the excited supermultiplet, meaning that decays into particles of the ground state are allowed. Since the operators presently used in the GEVP very likely have only small overlap with two-particle states, the full consideration of the effects of the strong decay channels would require to include in our variational basis also multi-particle operators, such as two glueballs with opposite momenta or a glueball and a gluino-glue. The resulting spectrum would be enlarged and would include also energy levels corresponding to twoparticles interacting on a periodic box, fulfilling a quantisation condition that determines the resonance properties, such as the phase shift, of the excited supermultiplet [25]. Such calculations require, however, a very good precision that is out of reach in the singlet sectors within the current statistics. The excited energy levels within the single-hadron approximation are already a good estimate of the mass of the excited supermultiplet.
Compared to our previous investigations [5] the extracted ground state masses turn out to be slightly smaller. The differences to the old results can be viewed as systematic errors of the old results, which are now much better under control. This includes the correct estimation of the fitting intervals and the exclusion of the gluino-glue mass at β = 1.9, κ = 0.1433 in the chiral extrapolation.   Table 1. Summary of the ensembles, masses of the lowest (superscript 0) and excited (superscript 1) states in lattice units, and mixing coefficient of the scalar channel. The ensembles have been presented in [4][5][6]. For consistency we use only the ensembles with one level of stout smearing. In some cases the data could not be extracted, because the fit intervals could not be determined reliably.  Table 2. Ground and excited state masses in lattice units extrapolated to the chiral limit (κ = κ c ), and the extrapolations to the continuum limit in units of the scale w 0 . The subscript 0 ++ denotes the result in the mixed channel, gg the gluino-glue, a-η the pseudoscalar meson and gb −+ the pseudoscalar glueball. The masses in the continuum limit are given as w 0 m. The quoted errors are statistical.

Mixing between glueballs and mesons
The numerical results of the extended variational approach indicate mixing between glueball and mesonic states. Such a mixing has been predicted in the literature on the structure of the lowest chiral supermultiplets [2,3]. From the calculated correlation matrices we can get insights about the nature of the physical states with respect to their mesonic and glueball content. From the eigenvectors v n of the GEVP (3.2) the corresponding physical states |n can be reconstructed and decomposed into a glueball contribution |φ (g) and a meson contribution |φ (m) : where v , respectively, and |Ω is the vacuum state. Note that |n , |φ (g) and |φ (m) are not normalised here.
ki c * (m) kj (7.5) Now we define the glueball and the meson contents of the physical state |n as the overlap of this state with the glueball and meson contributions We would like to point out that the glueball contribution |φ (g) and the meson contribution |φ (m) are not necessarily orthogonal to each other and therefore c (g)2 n and c (m)2 n in general do not add up to 1. The formulae above are independent of the time separation t. Therefore, we calculate the glueball and meson contents at each t and perform a fit to a constant over the region where the values form a plateau to extract the final results.
We have determined the glueball and meson contents c (g) n and c (m) n of the ground state and excited state in the scalar channel using again two meson and two glueball operators as in the analysis of the masses, see figure 9 and table 3. Extrapolated to the continuum limit, we find that the ground state has a glueball content of c (g) = 0.82(3) and a meson content of c (m) = 0.62 (6), meaning that the ground state in the scalar channel shows significant mixing of glueball and meson contents. The excited state appears to have the opposite glueball and meson contents than the ground state, namely c (g) = 0.6(1) and c (m) = 0.91 (6). While the ground state is more glueball like, the excited state is more meson like. This is especially visible in the chiral extrapolation at our smallest lattice spacing, figure 9. Interestingly the squares c (g)2 and c (m)2 do add up to 1 within errors as in the orthogonal case.

Conclusions
In this work we have presented in detail our improved techniques for the estimations of bound state masses in N = 1 supersymmetric Yang-Mills theory. We have used the variational method in all three channels, including a combined basis of glueball and mesonic . operators in the scalar channel. The enlarged basis allows for the first time an investigation of the mixing between these two classes of operators. Furthermore, we have optimised the Jacobi smearing, which now includes presmearing of the gauge fields and included an additional correlator C γ 4 (t) in the gluino-glue channel. The new techniques reduce excited state contributions and therefore allow more reliable estimations of the fit-intervals and smaller values of t min for the extraction of the ground state masses. We have reanalysed the gauge configurations of our previous investigations and improved our estimates of the ground state masses. Furthermore, we have obtained a first estimation of the masses of excited states, which have been out of reach in previous studies without the optimised variational basis and smearing techniques. We have combined the results for the masses of the ground states and first excited states from several different lattice spacings in an extrapolation to the continuum limit. The first interesting observation is the formation of supermultiplets. The ground state masses of the scalar, pseudoscalar, and fermion channel become approximately degenerate when extrapolated to the chiral and continuum limit. This is in line with our previous results and indicates the formation of a chiral supermultiplet of lightest states. The average mass of the lightest supermultiplet is around w 0 m (0) = 0.97(4).
Another very interesting result is the possible formation of a chiral supermultiplet of excited states at w 0 m (1) ≈ 3.1. The masses of these states in lattice units are around 0.5 a −1 , i. e. around half of the cutoff scale. It is quite unexpected that states at these high energies are not more affected by the supersymmetry breaking lattice artefacts. For a complete analysis a more detailed investigation of the bound state spectrum is required, including also other multiplets and a larger set of quantum numbers. We are currently adding, for instance, an investigation of the baryonic states of the theory [26].
We have investigated the mixing between glueball and mesonic operators in the scalar and pseudoscalar channels. In the scalar channel we have found significant mixing between JHEP04(2019)150 glueball and meson contents. In the pseudoscalar channel we found no hints of a significant mixing. The pseudoscalar ground state is clearly dominated by the mesonic contribution, whereas the scalar ground state has an apparent predominant glueball contribution.
Our results might help for a better understanding of the low energy effective action [1][2][3] and the conjectures presented in [27,28]. In [2,3], which refines and extends the analysis of [27], the nature of the lowest two chiral supermultiplets has been investigated on the basis of an effective action including mesonic and glueball-type degrees of freedom. The lightest states were conjectured to be of the glueball type, if mixing is not too strong. The argument is based on the perturbation of the effective theory by a small gluino mass, which leads to a splitting of the multiplets. In the mesonic multiplet the pseudoscalar meson becomes lighter than the scalar meson. Drawing on the proof [29] that the lightest scalar state, which has overlap with the 0 ++ glueball operator, is not heavier than the lightest pseudoscalar state, which has overlap with the 0 −+ glueball operator, they conclude that the multiplet of glueball states must be lighter than the multiplet of mesons.
On the other hand, in ref. [28] the lightest states are conjectured to be of mesonic type with a small mixing of the glueballs. Their analysis employs an effective Lagrangian with an arbitrary mixing angle between the glueball (R-charge 0) and meson (R-charge 2) multiplet. Generally, the effective Lagrangian allows either mesonic or glueball-type states to be the lightest ones, depending on an unknown parameter. The argument for the ordering of states is then based on the large-N c equivalence of SU(N c ) SYM and QCD-like theories. In QCD the mesonic η appears to be much lighter than the scalar glueball, which leads to the conjectures about the ordering of states and only a small mixing for SYM. A distinction between the scalar and pseudoscalar channels is, however, not being made.
Our numerical findings in the scalar channel, where the lightest state is dominantly of glueball type, are consistent with the predictions of [2,3,27]. On the other hand, it appears that the pseudoscalar ground state is dominated by the mesonic contribution with a negligible mixing of the glueball, which would correspond to the scenario advocated in [28]. Also, in contrast to the above arguments, we do not find the same mixing angle for the scalar and the pseudoscalar channel. A detailed consideration of the effects of a small gluino mass could shed more light on these questions.
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.