format: text

text:
# Scaling of the 𝐵 and 𝐷 meson spectrum in lattice QCD

## I Introduction

Mesonic bound states consisting of a single heavy quark, $`b`$ or $`c`$, and a light quark, $`u`$, $`d`$ or $`s`$, as well as gluons, provide an interesting laboratory to study strong interactions. The typical momentum within such states is much lower than the mass of the heavy quark. This leads to a situation where the heavy quark becomes non-relativistic and the properties of the bound state are essentially determined by the light quark and the glue. At leading order the splittings within the spectrum become independent of the properties of the heavy quark, such as its mass $`m_Q`$ and spin $`s_Q`$, so that orbital and radial excitation energies are expected to match between the $`B`$-system and the $`D`$-system. The resulting approximate SU($`2N_h`$) symmetry, with $`N_h`$ denoting the number of heavy flavours, is usually referred to as heavy quark symmetry, see  and the references therein. At the next order, $`1/m_Q`$ effects give rise to fine structure in the spectrum, several times larger in the $`D`$-system than for the $`B`$, see e.g.  for a review.

The spectrum of $`B`$ and $`D`$ states is not yet well established experimentally  although several new results have been reported recently . Here we study the spectrum theoretically and from first principles using lattice QCD. This will aid the experimental search for new states. In the case of well-established states it will provide a test for the theory and/or the systematic errors in our calculation. Of key interest are decay matrix elements for $`B`$-factory experiments. Knowing how well the spectrum has been obtained gives confidence that we understand how to simulate $`B`$ and $`D`$ mesons reliably. This is important for the analysis of systematic errors in matrix element determinations.

To formulate heavy $`b`$ and $`c`$ quarks on the lattice, a naïve discretisation is inappropriate since the lattice spacings currently available are not small compared to the Compton wave length of those quarks ($`m_Qa>1`$). Presently there are two different formulations available to simulate heavy quarks, non-relativistic QCD (NRQCD)  and the heavy Wilson approach . For the $`b`$-quark on present lattices both approaches become essentially the same. However, in this regime, NRQCD is to be preferred since the inclusion of higher order correction terms is easily implemented.

In this publication we report on our calculations of the $`B`$-meson spectrum for two different values of the lattice spacing $`a`$. Together with the results of , which were obtained with the same methods at another value of the lattice spacing, we can investigate the dependence on $`a`$ of our results. Physical results must be independent of $`a`$ and hence we can perform a test of systematic errors inherent in our calculation. We find no such errors at a significant level. In addition, on our coarsest lattice, we were able to simulate the $`D`$-meson spectrum and compare to results using heavy Wilson methods on finer lattices (where NRQCD does not work well since $`am_c<1`$). Early results on our coarse lattice have already been published in .

Section II gives details of the simulations we performed and section III gives details of our fit procedure. Section IV gives our determination of the bare $`b`$ and $`c`$ quark masses. Section V discusses the behaviour of the splittings in the spectrum that we obtain. This includes fits to the dependence of the splittings on the mass of the heavy quark. Section VI compares the results in physical units at different values of the lattice spacing and with previous results as well as with experiment. Readers interested in our results for the physical meson spectrum could jump directly to this section. Section VII contains our conclusions and our best estimate for the $`B`$-spectrum, based on the combined input from three different values of the lattice spacing.

## II Simulation details

### A Gauge field action

Our calculation was performed on two sets of gauge field configurations, which were generated using the Wilson gauge action

$$S_G=\beta \underset{x,\mu <\nu }{}[1\frac{1}{3}\mathrm{Re}\mathrm{Tr}(U_{x,\nu }U_{x+\widehat{\nu },\mu }U_{x+\widehat{\mu },\nu }^+U_{x,\mu }^+)].$$
(1)
This action has lattice artifacts of $`𝒪(a^2)`$. For the bare gauge coupling $`\beta `$, we used $`5.7`$ and $`6.2`$. The lattice volumes and the number of configurations are given in table I. We will refer to these configurations by their respective $`\beta `$-values.

### B Light quark propagators

The light quark propagators have been generated with the use of the Sheikholeslami-Wohlert action, also known as the clover action ,

$`S_L`$ $`=`$ $`a^4{\displaystyle \underset{x}{}}[\overline{\psi }_x\psi _x+\kappa {\displaystyle \underset{\mu }{}}[\overline{\psi }_{x\widehat{\mu }}(\gamma _\mu 1)U_{x\widehat{\mu },\mu }\psi _x\overline{\psi }_{x+\widehat{\mu }}(\gamma _\mu +1)U_{x,\mu }^+\psi _x]`$ (3)
$`a\frac{1}{2}\mathrm{i}c_{\mathrm{sw}}\kappa {\displaystyle \underset{\nu ,\rho }{}}\overline{\psi }_xF_{\nu \rho ,x}\sigma _{\nu \rho }\psi _x].`$
On the configuration set with $`\beta =5.7`$ the clover coefficient $`c_{\mathrm{sw}}`$ is set to its tadpole-improved tree level value $`c_{\mathrm{sw}}=1.5667`$, as determined from the 4th root of the plaquette . This reduces the lattice spacing artifacts in the light quark propagators to $`𝒪(\alpha _sa,a^2)`$. At $`\beta =6.2`$ we used the non-perturbative determined value of $`c_{sw}=1.6138`$, which removes the $`𝒪(\alpha _sa)`$ artifacts from the light quark propagator as well.

In reference  the light hadron spectrum at $`\beta =6.2`$ has been calculated using the non-perturbative as well as the tadpole-improved tree level value for $`c_{\mathrm{sw}}`$. No significant differences in the meson and baryon spectrum could be resolved between the two values of $`c_{\mathrm{sw}}`$. From this we expect the difference between tadpole and non-perturbatively improved light quarks at $`\beta =6.2`$ to be well covered by the size of the statistical errors in our case as well. This allows us to compare our $`\beta =6.2`$ results to the tadpole-improved results at $`\beta =5.7`$ and in reference .

For each value of $`\beta `$ we used 3 different values for the hopping parameter $`\kappa `$. The actual values are detailed in table II. The table also contains the values of $`\kappa _c`$ and $`\kappa _s`$ from the UKQCD collaboration  used in our calculation. The use of these values is appropriate for the analysis in terms of chiral extrapolations and scale setting that we have done. We also carefully include systematic errors from different chiral extrapolations and associated uncertainties in setting the scale. A recent re-analysis by UKQCD of their light hadron spectrum  gives somewhat different values for $`\kappa _c`$ and $`\kappa _s`$. Our errors encompass any changes this would produce in our physical results.

### C Heavy quark propagators

The typical momentum scale inside a heavy light meson such as a $`B`$ or $`D`$ meson is of the $`𝒪(\mathrm{\Lambda }_{\mathrm{QCD}})`$, which is small compared to the mass of the heavy quark. Therefore the mass of the heavy quark $`m_Q`$ represents an irrelevant scale for the dynamics of the mesonic bound state and it is possible to simulate these states on lattices with a lattice spacing larger than the Compton wavelength of the heavy quark.

In our simulation we use a non-relativistic expansion of the heavy quark Hamiltonian, which is known as NRQCD .

$`H`$ $`=`$ $`H_0+\delta H,`$ (5)
$`H_0`$ $`:=`$ $`{\displaystyle \frac{𝐃^2}{2m_Q}},`$ (6)
$`\delta H`$ $`:=`$ $`c_4{\displaystyle \frac{g}{2m_Q}}𝝈𝑩+c_2{\displaystyle \frac{\mathrm{i}g}{8m_Q^2}}(𝐃𝑬𝑬𝐃)c_3{\displaystyle \frac{g}{8m_Q^2}}𝝈(𝐃\times 𝑬𝑬\times 𝐃)`$ (8)
$`c_1{\displaystyle \frac{(𝐃^2)^2}{8m_Q^3}}+c_5a^2{\displaystyle \frac{𝐃^{(4)}}{24m_Q}}c_6a{\displaystyle \frac{(𝐃^2)^2}{16nm_Q^2}}.`$
Please note that the rest mass term of $`H`$ has been omitted, resulting in a shift of the Hamiltonian, which is discussed in section IV. In the case of a heavy-light meson the NRQCD expansion has to be organised in powers of $`\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q`$ . Here this expansion is used up to $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^2)`$. We also include the $`𝒑^4`$ term, which is believed to be the leading term in $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$. The last two terms correct for discretisation errors from finite lattice spacing in respectively the spatial and temporal derivatives. $`n`$ is a stability parameter used in the evolution equation (II C). The matching coefficients $`c_1`$, …, $`c_6`$ are set to their tadpole-improved tree level values .

With the Hamiltonian $`H`$ and $`\delta H`$ the propagator of the heavy quark can be obtained from a Schrödinger-type evolution equation

$`G_{t+1}`$ $`=`$ $`\left(1a\frac{1}{2}\delta H\right)\left(1a\frac{1}{2n}H_0\right)^nU_4^+\left(1a\frac{1}{2n}H_0\right)^n\left(1a\frac{1}{2}\delta H\right)G_t\text{ for }t>1,`$ (10)
$`G_1`$ $`=`$ $`\left(1a\frac{1}{2}\delta H\right)\left(1a\frac{1}{2n}H_0\right)^nU_4^+\left(1a\frac{1}{2n}H_0\right)^n\left(1a\frac{1}{2}\delta H\right)\varphi _x.`$ (11)
With $`\varphi _x`$ we denote the source smearing function used on the initial time slice. At $`\beta =5.7`$ we use 20 different values for $`m_Q`$ in the range $`0.6am_Q20.0`$ and at $`\beta =6.2`$ we use 10 values in the range $`1.1am_Q6.0`$. Details, including the $`n`$ values, are given in the tables III and IV. For each value of $`\beta `$ we performed 3 different runs. At $`\beta =5.7`$ we label them A, C and S; for $`\beta =6.2`$ they are labeled H, N and P.

For the $`S`$-wave mesons at $`\beta =5.7`$ we used up to three different smearing functions, $`\varphi _{\mathrm{G},0}`$, $`\varphi _{\mathrm{G},1}`$ and $`\varphi _{\mathrm{G},2}`$, in the different runs. These are convolutions of Gaussian functions for the light and the heavy quark with radii as detailed in table V. The configurations were fixed to Coulomb gauge. A local sink will be denoted with $`\varphi _\mathrm{L}`$. In most cases our final $`\beta =5.7`$ results were obtained with both sink and source smearing.

For $`\beta =6.2`$ we use smearing for the heavy quark propagators only. In run H at $`\beta =6.2`$ we applied a hybrid procedure of Jacobi smearing  and fuzzing . For runs N and P we fixed the configurations to Coulomb gauge. We used hydrogenic wave functions $`\varphi _{\mathrm{Hg},1}`$, $`\varphi _{\mathrm{Hg},2}`$ and $`\varphi _{\mathrm{He},1}`$ for run N. The indices ‘g’ and ‘e’ denote wave functions of the ground and first excited state. The details are given in table VI. In the P run we used Gaussian smearing with two different radii, $`ar_Q=2.5`$ and $`5.0`$.

The spin operators applied to construct mesonic states with the correct quantum numbers are detailed in table 1 of reference .

### D Lattice spacing

In the quenched approximation one obtains different values for the lattice spacing, depending on the quantity it is determined from. This is expected to be caused by the strong coupling $`\alpha _s`$ running differently in the real world and the quenched theory.

We use the physical mass of the $`\rho `$-meson  to fix the lattice spacing. This procedure is justified from the typical gluon momentum in a $`B`$ or $`D`$ meson being of similar size to the momentum in a light meson such as the $`\pi `$ and $`\rho `$. Since heavyonium states probe a higher physical scale these are not appropriate to fix the scale for a heavy-light system in the quenched approximation. Using the $`\rho `$-scale should take care of most of the quenching effects.

The determination of $`m_\rho `$ is complicated by the chiral extrapolation required, see reference  for a review. At $`\beta =5.7`$ we use the result of . The result of the linear extrapolation in the light quark mass $`m_q`$ is quoted as the central value and the deviation of the quadratic fit is treated as a systematic uncertainty. At $`\beta =6.2`$ a linear extrapolation is reported in reference . We treat the difference to the 3rd order extrapolation from  as a systematic uncertainty. The numbers are compiled in table VII. We use:

$`\beta =5.7:a^1=1.116(12)(_0^{+56})\text{ GeV},a=0.1768(19)(_{88}^{+0})\text{ fm},`$ (13)
$`\beta =6.2:a^1=2.59(_{10}^{+6})(_0^{+9})\text{ GeV},a=0.0762(_{18}^{+29})(_0^{26})\text{ fm}.`$ (14)
For comparison, table VII also shows the lattice spacing as obtained from the string tension $`\sigma `$ and the bottomonium splitting $`\overline{\chi }_b\mathrm{{\rm Y}}`$. As a physical value for $`\sigma `$ we choose a result obtained from a potential model fit to the charmonium spectrum . The lattice numbers originate from . These results are in agreement with the outcome of the $`m_\rho `$ analysis. As explained above, the bottomonium system probes a different scale and the values obtained using it do not agree with the result from light spectroscopy .

## III Fitting techniques

### A Parametrisations

At $`\beta =5.7`$ we used several different smearings at source and sink. For hadron correlators with a local sink, we applied simultaneous vector fits, requiring the fitted mass(es) $`m_k`$ to agree for all propagators:

$`\varphi _\mathrm{L}(t)|\varphi _i(0)`$ $`=`$ $`{\displaystyle \underset{k=1}{\overset{n}{}}}A_{i,k}\mathrm{exp}(m_kt),1im,`$ (15)
$`A_{i,k}`$ $`=`$ $`\varphi _\mathrm{L}|\psi _k\psi _k|\varphi _i.`$ (16)
In the case of sink and source smearing, we used simultaneous matrix fits. In matrix fits, the fitted amplitudes are constrained in their relationship with each other as well:

$`\varphi _j(t)|\varphi _i(0)`$ $`=`$ $`{\displaystyle \underset{k=1}{\overset{n}{}}}B_{j,k}^{}B_{i,k}\mathrm{exp}(m_kt),1j,im,`$ (17)
$`B_{i,k}`$ $`=`$ $`\psi _k|\varphi _i.`$ (18)
The fitting techniques are described in more detail in reference . We found matrix fits to be more precise with respect to statistical errors. Due to the omission of the rest mass in eq. (II C) the fitted mass is shifted with respect to the bound state mass. We denote the result of the fit as the simulation mass, $`m_{\mathrm{sim}}`$. The determination of the shift will be discussed in section IV.

To extract mass splittings we applied two different procedures. One is to fit the masses as above, take their difference, and then calculate the error from the bootstrap or jackknife samples of the difference. With this procedure one can easily take advantage of using different smearings. In the case of a single smearing function, a ratio-fit provides an alternative . For this one divides the bootstrap or jackknife samples of the two propagators and fits the outcome with an exponential ansatz. The mass shift cancels out of the difference in both procedures.

### B Pseudo-scalar and vector meson

On the $`\beta =5.7`$ configurations the simulation masses for pseudo-scalar and vector mesons have been determined most accurately in run A. In this run we only used the smearing functions $`\varphi _{\mathrm{G},1}`$ and $`\varphi _{\mathrm{G},2}`$. We found the double exponential matrix fit with sink and source smearing to deliver the most precise result. For the fit range we choose the initial time slice $`t_{\mathrm{min}}`$ two time slices larger than the first time slice delivering a reasonable $`\chi ^2`$. In general we choose the number of dropped time slices multiplied by the number of propagators used for the fit to be larger than or equal to the number of fit parameters. The reason for this procedure is as follows. The first reasonable value of $`\chi ^2`$ is observed once the residual excitations are just masked by the statistical uncertainties, which allows for them to be still of similar size. Each excited data point can eat up one fit parameter. Dropping as many data points as fit parameters delivers a fit which is entirely dominated by statistical fluctuations. The residual fit range dependence of those fits becomes negligible against the statistical uncertainties. We judge $`\chi ^2`$ values resulting in a $`Q0.1`$ as reasonable, where $`Q`$ denotes the probability of a fit having an even higher value of $`\chi ^2`$. The final result is given in table VIII.

In run H at $`\beta =6.2`$ we only had one smearing function available. We extracted the final results from single exponential fits to the propagators with source and sink smearing. Their fit results turned out to be more precise than the ones from using a local sink. The final fit range was determined such that we observed a reasonable $`\chi ^2`$ and achieved stability of the fitted result against variation of the fit range. The results are displayed in table IX.

In the run N we used hydrogenic wave functions of different radii. We generated smeared local and smeared smeared meson propagators. However no cross correlators, e.g. $`\varphi _{\mathrm{Hg},1}`$ at sink and $`\varphi _{\mathrm{He},1}`$ at source, were calculated. Hence eq. (17) was inapplicable and we had to use vector fits in the case of smearing at sink and source as well.

In double exponential vector fits to two smearing functions, we observed extremely low values of $`Q`$. We observed that this is connected to unfortunate statistical fluctuations on certain time slices. However the fit parameters turned out to be stable with respect to variations of the fit range. These fits will be discussed in subsection V B in more detail. To obtain a more precise result for the $`S`$-wave ground states, we resorted to single exponential fits to single propagators and compared the outcome for the different smearing functions. This is shown in figure 1 for the pseudo-scalar propagator at $`am_Q=2.5`$. The octagons indicate the final result for each propagator, as determined from the $`Q`$-value after dropping two time slices. Within statistical errors all results are in reasonable agreement with each other. For the final result, which is also included in table IX, we choose the smeared-smeared $`\varphi _{\mathrm{Hg},1}`$ propagator. In the end these deliver the more accurate hyperfine splitting, due to superior noise cancellation between the pseudo-scalar and the vector meson state.

In this context it is interesting to note that the propagators with local sink and $`\varphi _{\mathrm{Hg},2}`$ source smearing plateau much later than the others, but the results are in agreement with those from other propagators.

For the $`S`$ wave states in the run P we only had the Gaussian smearing at the source with radius $`ar_0=2.5`$ and local sink. Since these propagators plateau quite late, we used $`t_{\mathrm{min}}=16`$ for the final result, the error bars are not competitive with those above. Since they are needed for the later analysis of the $`P`$ states we include them as well in table IX.

To describe physical bound states involving light $`u`$ and $`d`$ quarks, the results of tables VIII and IX have to be extrapolated in the light quark hopping parameter. On both sets of configurations, the difference between the critical and normal hopping parameter is smaller than the uncertainty we assigned to $`\kappa _c`$ in table II  and we use $`\kappa _c`$ in our extrapolations. The normal hopping parameter is the one for which the extrapolations deliver the physical $`m_\pi /m_\rho `$ ratio.

Due to the high statistical accuracy we achieved at $`\beta =5.7`$ in the pseudo-scalar case, a linear ansatz in $`am_q:=\frac{1}{2}(\frac{1}{\kappa }\frac{1}{\kappa _c})`$ in a full covariant fit to all three data points, results in a fit with $`\chi ^2/\mathrm{d}.\mathrm{o}.\mathrm{f}.>8/1`$ for $`am_Q<10`$. This corresponds to $`Q<0.004`$. The resulting curves do not describe the data. We carefully checked whether this is caused by a residual fit range dependence and found all the fit parameters including the would-be strange to non-strange meson splitting to be stable against variation of the fit range. This was done for an initial time slice $`t_{\mathrm{min}}`$ in a range from $`3`$ to $`6`$.

We therefore extracted our final result from a linear spline to the points with highest and lowest $`m_q`$ and use the deviation of a quadratic spline as a systematic uncertainty of the chiral extrapolation. An example for the extrapolation is given in figure 2. From the figure it is obvious that interpolations to extract the heavy-strange meson mass are insensitive to the different ansätze and we do not assign an uncertainty due to the different interpolations. However, in the case of the heavy strange meson, we are faced with the problem that $`\kappa _s`$ is highly sensitive to the quantity it is determined from. Our central value is interpolated to the $`\kappa `$ as determined from $`m_K/m_\rho `$, and the difference to the outcome for $`\kappa `$ corresponding to $`m_\varphi /m_\rho `$ is treated as an uncertainty of the quenched approximation. The results are presented in table X.

For $`\beta =6.2`$ the statistical accuracy is not as high and our data are well described by linear extrapolations. The results are presented in table XI.

## IV Heavy quark masses

### A Mass shift from dispersion relation

The omission of the rest mass term $`m_Q`$ in the Hamiltonian eq. (II C) causes most of the shift of the simulation mass $`m_{\mathrm{sim}}`$ with respect to the physical meson mass. The mass, $`m_{\mathrm{rel}}`$, of the meson can be determined from the relativistic dispersion relation of the meson $`E(\stackrel{}{p})=\sqrt{m_{\mathrm{rel}}^2+\stackrel{}{p}^2}`$, which gives

$$m_{\mathrm{rel}}=\frac{\stackrel{}{p}^2[E(\stackrel{}{p})E(0)]^2}{2[E(\stackrel{}{p})E(0)]}.$$
(19)
Here $`E(\stackrel{}{p})`$ denotes the total energy of the meson. The mass shift $`\mathrm{\Delta }_{\mathrm{rel}}`$ is defined as the difference

$$\mathrm{\Delta }_{\mathrm{rel}}:=m_{\mathrm{rel}}m_{\mathrm{sim}}.$$
(20)
This shift per heavy quark should be universal for all hadronic states simulated at the bare heavy quark mass $`m_Q`$.

In our calculation at $`\beta =5.7`$ we determined the mass shift from the difference in energy of the pseudo-scalar meson propagators with $`a|\stackrel{}{p}|=0`$ and $`2\pi /12`$. This was done in run C at $`\kappa =0.1400`$ with source smearing $`\varphi _{\mathrm{G},1}`$ and a local sink. At large values of $`m_Q`$, we found a single exponential fit to the ratio of the correlators to plateau much later than the fits to the individual propagators. This is reflected in a large fit range dependence of the jackknife difference of the masses of the individual fits, for time slices in which no plateau was observed in the ratio-fit. For our final result we choose a minimal $`t`$-value two time slices larger than the first $`t`$-value for which we obtained a decent $`\chi ^2`$ in a fit to the ratio of propagators. The final result is presented in table XII and figure 3.

We also tried simultaneous vector fits according to eq. (15) with two exponents. We used propagators with source smearing $`\varphi _{\mathrm{G},1}`$ and $`\varphi _{\mathrm{G},2}`$. The jackknifed difference of the fitted ground state mass is in agreement with the above procedure; however, the statistical uncertainties, especially for large values of $`m_Q`$, are larger.

For $`\beta =6.2`$ we calculated the mass shift in heavy quarkonia, since the statistical precision for heavy-light correlators at finite momentum was not sufficient. In the following, mass shifts from heavy quarkonia will be denoted by $`\mathrm{\Delta }_H`$. We simulated the vector-meson for $`a|\stackrel{}{p}|2\frac{2\pi }{24}`$. The kinetic mass $`m_1`$ was obtained from fits to the dispersion relations:

$`E_{\mathrm{sim}}(\stackrel{}{p})`$ $`=`$ $`m_0+{\displaystyle \frac{\stackrel{}{p}^2}{2m_1}}{\displaystyle \frac{\stackrel{}{p}^4}{8m_2^3}},`$ (22)
$`E_{\mathrm{sim}}(\stackrel{}{p})`$ $`=`$ $`m_0+{\displaystyle \frac{\stackrel{}{p}^2}{2m_1}}{\displaystyle \frac{\stackrel{}{p}^4}{8m_1^3}},`$ (23)
$`E_{\mathrm{sim}}(\stackrel{}{p})`$ $`=`$ $`m_0+{\displaystyle \frac{\stackrel{}{p}^2}{2m_1}},`$ (24)
with parameters $`m_0`$, $`m_1`$ and $`m_2`$. $`E_{\mathrm{sim}}(\stackrel{}{p})`$ denotes the simulation energy as determined from the propagator falloff. In the case of $`am_Q1.3`$ we used the ansätze (22) and (23); for the three heavier $`m_Q`$-values, (23) and (24). All fits gave fit parameters which were consistent within half of the statistical error. To obtain the shifts required for heavy-light spectroscopy we subtracted the simulation mass of the quarkonium vector-meson and divided by two. The final results are displayed in table XIII and figure 3. It is interesting to compare to the result from reference $`a\mathrm{\Delta }_H=1.29(2)`$ obtained at $`am_Q=1.22`$. Due to higher statistics, this result is much more precise. This value is included as a square into figure 3 and agrees well with the newer results.

### B Mass shift in perturbation theory

The mass shift $`\mathrm{\Delta }`$ can also be calculated in lattice perturbation theory :

$$\mathrm{\Delta }_{\mathrm{pert}}=Z_mm_QE_0.$$
(25)
Here $`Z_m`$ denotes the renormalisation constant connecting the bare lattice mass $`m_Q`$ with the pole mass and $`E_0`$ denotes the heavy quark self energy constant. In the perturbative expansion the 1-loop contributions from $`Z_m`$ and $`E_0`$ cancel each other to a large extent and the direct perturbative expansion of $`\mathrm{\Delta }_{\mathrm{pert}}`$ is much better behaved than either perturbative series on its own. The Lepage Mackenzie scale $`aq^{}`$  has been determined separately for $`\mathrm{\Delta }`$ and it is larger than for $`Z_m`$ or for $`E_0`$. The coefficients for

$$\mathrm{\Delta }_{\mathrm{pert}}=m_Q[1+\alpha _s(aq^{})\mathrm{\Delta }^{(1)}]$$
(26)
can be found in table XIV. We use the $`\alpha _P(aq=3.4)`$ values as determined from the $`1\times 1`$ Wilson loop  with 2-loop running in order to evolve to the respective $`aq^{}`$. For the final mass shift we assign a relative uncertainty of $`\alpha _s^2(aq^{})`$. Since $`\mathrm{\Delta }^{(1)}`$ is small, this is more conservative than the squared 1-loop contribution. The final results are displayed in table XII for $`\beta =5.7`$ and table XIII for $`\beta =6.2`$. For values of $`m_Q`$ not included in table XIV we interpolated linearly between the nearby values, which is completely sufficient within the claimed accuracy. The results for $`a\mathrm{\Delta }`$ from perturbation theory and the lattice simulation are compared in figure 3. Apart from possibly the low $`m_Q`$ region at $`\beta =5.7`$, the figure shows excellent agreement between the two ways of calculating the mass shift.

For $`\beta =5.7`$ the stability parameter $`n`$ differs in some cases between the perturbative results and the simulation. However for $`am_Q=4`$, where perturbative results exist for $`n=1`$ and $`2`$, the effect of $`n`$ is completely negligible: we obtain $`\mathrm{\Delta }_{\mathrm{pert}}=3.88(22)`$ vs $`3.89(24)`$. From a comparison of the simulation result of the runs A and C at $`\beta =5.7`$ we can also obtain evidence of the effect of the different $`n`$ on the simulation mass $`m_{\mathrm{sim}}`$. For $`am_Q=1.0`$ and $`\kappa =0.1400`$ we measure $`am_{\mathrm{sim},\mathrm{ps}}=0.6265(21)`$ and $`0.6248(21)`$ for $`n=5`$ and $`6`$ respectively. This difference is again completely negligible against the uncertainty we assign to $`a\mathrm{\Delta }`$, even if we enlarge it by a factor of $`3`$ to allow for a larger effect between $`n=4`$ and $`5`$. The former $`n`$ was used in the perturbation theory. Note also that this difference tends to be in the opposite direction to that in $`a\mathrm{\Delta }`$ implying that the effect of $`n`$ on the physical mass is reduced when compared to the shift.

Here it is interesting to note that physical mass differences like the hyperfine splitting $`m_{\mathrm{hpf}}=m_{\mathrm{sim},\mathrm{v}}m_{\mathrm{sim},\mathrm{ps}}`$ are even less sensitive to $`n`$. At the above mass parameter of $`am_Q=1.0`$ we measure $`am_{\mathrm{hpf}}=0.0833(20)`$ for $`n=4`$ and $`0.0835(20)`$ for $`n=5`$.

In summary the differences in $`n`$ between the different runs as well as the perturbative shifts can be neglected safely even at the high level of accuracy we achieved here. This leaves us with a discrepancy between $`\mathrm{\Delta }_{\mathrm{pert}}`$ and $`\mathrm{\Delta }_{\mathrm{rel}}`$ for our lowest $`m_Q`$-values, which is roughly twice as large as the uncertainty we assign to the perturbative result.

On the other hand for $`\beta =6.2`$ we observe excellent agreement between the precise result of  with the perturbative calculation at the relatively low value $`am_Q=1.22`$.

Given a value for the shift $`\mathrm{\Delta }`$ and the simulation mass $`m_{\mathrm{sim}}`$ from tables VIII, IX, XII and XIII, we can now calculate absolute masses for all the states. We do this for the ground state vector and pseudo-scalar mesons, both to fix the quark mass, as described in the next subsection. Moreover, we use the meson mass rather than the quark mass to discuss the $`m_Q`$ dependence, since it is more directly comparable to experiment. We frequently plot results against $`1/m_{\mathrm{sav}}`$, where $`m_{\mathrm{sav}}`$ is the spin-average of the ground state vector and pseudo-scalar mesons

$$m_{\mathrm{sav}}=\frac{1}{4}(3m_\mathrm{v}+m_{\mathrm{ps}}).$$
(27)
This is preferable to $`m_{\mathrm{ps}}`$ alone since the spin-averaging reduces the dependence on sub-leading spin-dependent terms.

### C Bare heavy quark mass

To determine the bare quark mass $`m_Q`$ corresponding to the $`b`$ and $`c`$-quark, we compared the mass of the spin-averaged $`S`$-wave meson denoted with an overbar, with the experimental result. We used 5313 MeV for the $`\overline{B}`$, 5405 MeV for the $`\overline{B}_s`$, 1973 MeV for the $`\overline{D}`$ and 2076 MeV for the $`\overline{D}_s`$ . For the interpolations we used spline-fits to three neighbouring points. The fits were done quadratically in $`m_Q`$ and $`1/m_Q`$ and no significant difference was observed between the two. From the strange and non-strange meson we obtained identical results for the quark masses:

$`am_b`$ $`=`$ $`1.64(5)(_5^{+8})(_7^{+0}),\beta =6.2`$ (29)
$`am_b`$ $`=`$ $`4.20(25)(5)(_{24}^{+0}),\beta =5.7`$ (30)
$`am_c`$ $`=`$ $`0.87(6)(3)(_{13}^{+0}),\beta =5.7`$ (31)
The errors as indicated in the parentheses give the uncertainty arising from the mass shift and the statistical and systematic uncertainty of the $`a`$ determination. The uncertainties associated with the simulation mass $`m_{\mathrm{sim}}`$ are completely negligible here. For $`m_b`$ we used the perturbative shifts $`\mathrm{\Delta }_{\mathrm{pert}}`$. Using $`\mathrm{\Delta }_H`$ at $`\beta =6.2`$ delivers $`am_b=1.59(_5^{+14})(_3^{+6})(_5^{+0})`$ and using $`\mathrm{\Delta }_{\mathrm{rel}}`$ at $`\beta =5.7`$ gives $`am_b=4.16(_{31}^{+53})(7)(_{30}^{+0})`$, which is agreement with $`\mathrm{\Delta }_{\mathrm{pert}}`$ but with larger error. For $`m_c`$ we used the simulation result $`\mathrm{\Delta }_{\mathrm{rel}}`$. Here $`\mathrm{\Delta }_{\mathrm{pert}}`$ would give $`am_c=1.02(8)(2)(_{10}^{+0})`$. The deviation from the result eq. (31) reflects the difference between the $`\mathrm{\Delta }`$-values at low $`am_Q`$ discussed above.

This careful analysis to fix the bare heavy quark mass is particularly necessary for fine structure splittings in the spectrum to be discussed in the next section. These are very sensitive to the heavy quark mass, generally as $`1/m_Q`$. In addition, any errors in the heavy quark mass must be fed into errors in the fine structure splittings in order to avoid underestimating those errors.

The bare masses do not scale with the lattice spacing as expected, because they are unphysical. A better quantity to consider would be the mass in the $`\overline{\mathrm{MS}}`$-scheme. This will be discussed in a future publication .

## V Mass dependence of level splittings

In this section we describe how the results for the level splittings are extracted from our data. The dependence of the different splittings on the light and heavy quark mass is also discussed.

### A Flavour dependent splittings

The mass difference between heavy-light states distinguished only by their strangeness survives into the static limit. Based on the ideas of heavy quark symmetry such splittings are expected to depend weakly on the mass of the heavy quark. If a spin-averaged combination is taken, the leading heavy quark mass dependence arises purely from the kinetic term in eq. (6). The size of the slope in $`1/m_{\mathrm{sav}}`$ then gives information on the difference in $`p_b^2`$ for the $`b`$-quark in the strange and non-strange states.

We calculated this splitting from the ground state S-wave results for $`\kappa _s`$ and $`\kappa _c`$. The result is highly sensitive to the reported uncertainties in the chiral extrapolation and the determination of the strange hopping parameter $`\kappa _s`$. At $`\beta =5.7`$ we determine the statistical uncertainties in a jackknife procedure applied to the difference of the individual results and at $`\beta =6.2`$ we use the bootstrap.

For $`\beta =5.7`$ our statistical errors are very small and we consider additional systematic uncertainties for the chiral extrapolation and $`\kappa _s`$. For $`\beta =6.2`$ the quality of our data is not as good and we give statistical errors only. Using the $`\kappa _s`$-value determined from the $`\varphi `$ would lead to a 9% increase of the result, which is small compared to our statistical errors. The results are displayed in tables XV and XVI. In figure 4 we plot the result for the spin-averaged splitting at $`\beta =5.7`$ versus the inverse of $`m_{\mathrm{sav}}`$, in order to display its heavy quark dependence.

The figure displays a clear increase of the splitting with decreasing heavy quark mass. To quantify the slope of this dependence we perform a linear fit of the splitting result versus $`1/m_{\mathrm{sav}}`$. The result, converted into physical units, is detailed in table XXV. The slope corresponds to a $`p_b^2`$ difference of $`[0.25(3)\text{ GeV}]^2`$, which is of the size of $`\mathrm{\Lambda }_{\mathrm{QCD}}^2`$, as expected.

Because of the larger statistical uncertainties, we do not observe a significant slope at $`\beta =6.2`$. The data can be fitted nicely to a constant.

A comparison of our results with the ones of  for the pseudo-scalar case is plotted in figure 5. In this plot we show the result for the strange quark as determined from the $`K/\rho `$ mass ratio only. Due to the large error bars we do not include the results obtained at the larger values of $`m_Q`$ for $`\beta =6.2`$. Within the accuracy of around 12% in the case of $`\beta =6.2`$ or better, no sign of scaling violations shows up in the plot. We also observe excellent agreement with the experimental result.

### B Radial excitations

In order to obtain a reasonably stable and long plateau for the radially excited $`S`$-states on the coarse lattice at $`\beta =5.7`$ we applied triple exponential matrix fits to the three smearing functions $`\varphi _{\mathrm{G},0}`$, $`\varphi _{\mathrm{G},1}`$ and $`\varphi _{\mathrm{G},2}`$. This was done for the run S for a single $`\kappa `$ of $`0.1400`$ only, which is approximately equal to the strange as determined from the $`K`$-meson. Since in reference  the dependence of the $`2^1S_01^1S_0`$ splitting on the light quark mass was found to be very small, a variation of less than 2% when fixing $`\kappa _s`$ from the $`K`$ or $`K^{}`$-meson, we can ignore any mismatch in our $`\kappa `$ vs $`\kappa _s`$ compared to the statistical uncertainties. We therefore treat our result as the answer for this splitting with $`\kappa _s`$ as determined from the $`K`$.

In figure 6 we show a typical example for the excellent stability of the simulation masses $`am_{\mathrm{sim}}`$ against variation of the starting point $`t_{\mathrm{min}}`$ of the fit range. The extent in $`t_{\mathrm{min}}`$ for which we can resolve the excited state is 5 time slices or $`0.28`$ GeV<sup>-1</sup>. The rate of its disappearance is set by the $`2S1S`$ splitting of $`600`$ MeV.

In the figure the first good value of $`Q`$ is observed for $`t_{\mathrm{min}}=2`$. To be safe with respect to residual excitations we quote the final result for a fit range starting at time slice $`4`$, which is the procedure described in subsection III B. The peak in $`Q`$ at $`t_{\mathrm{min}}=5`$ results from the fit becoming insensitive to the third exponential at this point.

The results for all 6 heavy quark masses are given in table XVII. In figure 7 we plot the heavy quark mass dependence of the spin-averaged splitting. The result shows a clear increase as the heavy quark mass is reduced. In table XXV we detail fit results for the offset and slope of this splitting with respect to $`1/m_{\mathrm{sav}}`$. From the assumption that the increase of the splitting with $`1/m_{\mathrm{sav}}`$ is caused by the difference in the kinetic energy $`p^2/2m_Q`$ between ground and radially excited states, the fitted slope gives a $`p_b^2`$ difference of $`[0.95(15)\text{ GeV}]^2`$. This is of the size of a few times $`\mathrm{\Lambda }_{\mathrm{QCD}}^2`$ as would be expected.

On our fine lattice, since no cross-correlators between the different smearings had been calculated, we used simultaneous vector fits in all cases. The differences between the smearings $`\varphi _{\mathrm{Hg},1}`$, $`\varphi _{\mathrm{He},1}`$ and $`\varphi _{\mathrm{Hg},1}`$ turned out to be too small for simultaneous fits with three exponentials and we had to resort to fits using 2 exponentials. Again we choose the starting point $`t_{\mathrm{min}}`$ of the fit range as described in subsection III B. With this procedure it is possible to extract reliable information on the excited state, as can be verified from the tables of reference  for the case of $`\mathrm{{\rm Y}}`$-spectroscopy. Using propagators with sink and source smearing in vector fits, leads to a suppression of the excited state contamination, which made it impossible to extract a signal for the excited state. Therefore we used propagators with smearing at the source and local sinks to extract the radially excited states.

As denoted earlier, the fits to these propagators are plagued by statistical fluctuations, which lead to quite large $`\chi ^2`$ and low $`Q`$ values. However, the fits describe the data reasonably and the fitted parameters are stable against variation of $`t_{\mathrm{min}}`$. This is shown in figure 8. For a fit range starting point $`t_{\mathrm{min}}5`$ we obtain $`Q>1\%`$, which is low compared to what we obtained in the other fits. However it is not that low, that the fit could be ruled out on statistical grounds. Together with the good stability of the fitted masses against variations of $`t_{\mathrm{min}}`$ as displayed on the right hand side of the figure, we believe that our signal is genuine. In this example we extract the final result from $`t_{\mathrm{min}}=8`$. The results for the radially excited $`S`$-wave for this lattice spacing are compiled in table XVIII.

### C Orbital excitations

Orbitally excited $`P`$-state mesons have been investigated at both of our lattice spacings. The possible states consist of four non-degenerate energy levels; total angular momentum $`J=0`$ and $`2`$ as well as two $`J=1`$ states. In the heavy quark symmetry picture a $`jj`$ coupled basis is appropriate. In the vicinity of the static limit, the $`J=2`$ and the higher of the $`J=1`$ states are close and separated from the $`J=0`$ and lower $`J=1`$. The former correspond to a light quark total angular momentum of $`j_l=\frac{3}{2}`$, the latter of $`j_l=\frac{1}{2}`$. We use an $`LS`$ coupled basis to study the states, but expect our $`{}_{}{}^{1}P_{1}^{}`$ and $`{}_{}{}^{3}P_{1}^{}`$ channels to mix, leading to the observation of the lower $`J^P=1^+`$ state with both operators. We will denote the state corresponding to $`j_l=\frac{3}{2}`$ with a prime.

At $`\beta =5.7`$ again we use one light hopping parameter, $`\kappa =0.1400`$. As in the case of the radial excitations, we treat this as the value corresponding to the strange quark as determined from the $`K`$-meson and the simulations have been performed in the run S. We used the derivatives of the smearings $`\varphi _{\mathrm{G},0}`$ and $`\varphi _{\mathrm{G},1}`$ at source and sink. The final results were obtained from double exponential simultaneous matrix fits to both smearings and are listed in table XIX.

The selection of the fit range proved to be very delicate for this $`a`$-value. The statistical error grows rapidly when increasing $`t_{\mathrm{min}}`$, since the signal to noise is exponentially related to the $`PS`$ splitting . We give an example in section V F, where the fine structure is discussed. For the lightest values of $`am_Q`$ and correlators with $`{}_{}{}^{3}P_{2}^{}`$ wave operators in the $`T`$-representation we obtained small values of $`Q`$ of a few permille. We include the corresponding mass values into the table for the sake of completeness and mark them in italics. However we disregard them in the further evaluation. The results are always in agreement with the ones obtained in the $`E`$-representation and we do not believe there is a serious problem with this, simply statistical fluctuations.

In the second to last column of table XIX we give the spin-averaged $`P`$-state result which we calculate as

$$m(P_{\mathrm{sav}})=\frac{1}{12}[1m(^3P_0)+3m(^1P_1)+3m(^3P_1)+5m(^3P_2E)].$$
(32)
The result is also shown in figure 9. For comparison we include the experimental result for the $`B_{sJ}^{}(5850)`$ resonance and the spin-average of the $`D_{s1}`$ and $`D_{s2}^{}`$ . The figure displays at most a mild heavy quark mass dependence. To quantify this, we report in table XXV on the offset and slope of this splitting in physical units. The slope is consistent with a $`p_b^2`$ difference of $`𝒪(\mathrm{\Lambda }_{\mathrm{QCD}}^2)`$, but is also consistent with zero.

$`P`$-states were also investigated in the run P for $`\beta =6.2`$. We choose $`am_Q=1.6`$, directly corresponding to the $`b`$-quark in eq. (29). For the light quarks we use the strange hopping parameter $`\kappa =0.1346`$. Since the results on our coarse lattice depend only weakly on $`m_Q`$ we take the outcome as the final answer for $`B_s`$.

In the simulations we used two different smeared sources together with local sinks. Again we used derivatives of Gaussian smearing functions. The masses were extracted from double exponential vector-fits to both propagators. We observe reasonable $`Q`$-values for all applied operators and include all channels into the spin-average

$$m(P_{\mathrm{sav}})=\frac{1}{12}[1m(^3P_0)+3m(^1P_1)+3m(^3P_1)+2m(^3P_2E)+3m(^3P_2T)].$$
(33)
The results for the fitted masses are displayed in table XX. The splitting between the spin-averaged $`P`$ and $`S`$ waves is given in tables XXI and XXII.

### D Radially excited $`P`$-states

Having available 2 different smearing functions at both of our lattice spacings, it is possible to obtain information on the radially excited $`P`$-states as well. In figures 10 we show the dependence of the fitted masses of the spin-averages of the $`1P`$ and the $`2P`$ on the starting point $`t_{\mathrm{min}}`$ of the fit range for two different values of $`m_Q`$. The figure displays a clear signal for an excited state and reasonable stability with respect to variations of $`t_{\mathrm{min}}`$. However the error grows rapidly with increasing $`t_{\mathrm{min}}`$. The $`Q`$-values of the $`{}_{}{}^{3}P_{2}^{}E`$-fit, which is the last of the individual states included in the spin-average to reach a plateau, are included in figure 17. Decent $`Q`$-values are observed for $`t_{\mathrm{min}}=3`$. Since this is a 6 parameter matrix fit to four propagators, we take our final result from $`t_{\mathrm{min}}=5`$. The results for the spin-averaged $`2P`$-state are summarised in table XIX. We give the splitting to the spin-averaged $`1S`$ and $`1P`$-states in table XXI. We do not observe a significant slope for the splitting with respect to $`1/m_{\mathrm{sav}}`$.

For $`\beta =6.2`$ we show the plateau in figure 11. Due to the finer lattice, the growth in the error with increasing $`t_{\mathrm{min}}`$ is smaller than before. For the mass parameters used here, an example of a $`Q`$-plot will be given in figure 18 below. Here the first decent $`Q`$ is observed for $`t_{\mathrm{min}}=2`$. Since this is a vector fit to 2 propagators, we take our final results from $`t_{\mathrm{min}}=6`$. The result and splittings are included in tables XX and XXII.

As noted when discussing radially excited $`S`$-states, due to our conservative selection of the fit range, we expect residual excitations to be negligible within the quoted statistical errors.

### E Hyperfine splittings

The mass difference between a pseudo-scalar and a vector $`S`$-wave meson is caused by the spin of the heavy quark. This hyperfine splitting is expected to vanish in the limit of infinitely heavy quark mass.

On our $`\beta =5.7`$ lattice we determined the hyperfine splitting $`m_{\mathrm{hpf}}`$ from the difference of the results in table VIII. A crucial ingredient in obtaining a small statistical error is to choose identical fitting ranges for both correlators. If one of them has a plateau at a larger value of $`t_{\mathrm{min}}`$ than the other, we took this larger value to obtain $`m_{\mathrm{hpf}}`$ from the difference of the fitted masses. The error in this procedure is estimated with a jackknife. The results are displayed in table XXIII.

The chiral extrapolation of the hyperfine splitting turned out to be less difficult than that for the simulation mass of the pseudo-scalar and vector mesons. The curvature seems to cancel out between them and we have been able to perform linear fits to extrapolate to the chiral limit. However, in order to be consistent, we assign a systematic uncertainty to the result. This uncertainty is obtained from the difference to the outcome of first extrapolating the individual mesons to the chiral limit and then determining the hyperfine splitting. In this case we use the quadratic extrapolations to the chiral limit to take a possible curvature into account.

In our chiral extrapolations of the hyperfine splitting we observe a negative slope with respect to the mass $`am_q`$ of the light quark, which is illustrated in figure 12. The left hand side gives an example of our chiral extrapolations and the right hand side shows the slopes measured at each value of $`am_Q`$. In order to construct a physically meaningful quantity, the latter has been multiplied by the strange quark mass, such that we can compare with experimental results for the difference between the strange and non-strange hyperfine splittings.

For $`B`$-mesons the light quark dependence of the hyperfine splitting is not well resolved experimentally, because of large uncertainties in the $`B_s`$ hyperfine splitting. For $`D`$-mesons the situation is much clearer and one observes an increase with the light quark mass. However the magnitude of the slope is largely dependent on whether you compare the $`D_s`$ hyperfine splitting with the hyperfine splitting of the $`D^+`$ or the $`D^0`$. We expect this difference in the experimental results to be mainly due to QED-effects, since these come in with opposite signs in the $`D^+`$ and the $`D^0`$. Since the $`D_s`$ is positively charged as well, QED effects should largely cancel when comparing the hyperfine splittings of the $`D_s`$ and the $`D^+`$ and one obtains a positive slope for the $`D`$-meson from the experiment.

Comparing our data to the experimental results, one observes our hyperfine splittings to be too small. This will be discussed in more detail in section VI. With respect to the slope, the result at the $`D`$ has clearly the wrong sign and its magnitude is approximately twice as large as that from the experiment. We did not observe this effect in our $`\beta =6.2`$ results, neither was it observed in . Both of these results did not achieve the high statistical accuracy we have at $`\beta =5.7`$ and also use values of the heavy quark mass at around the $`B`$ or heavier. For those values of $`m_Q`$ the light quark mass dependence at $`\beta =5.7`$ is also not that significant.

A slope of similar sign and size has been observed in the calculations of , although the authors did not comment on this. Reference  used a highly improved gluonic action with NRQCD heavy quarks on even coarser lattices and  a heavy clover action for the heavy quarks on a finer lattice $`\beta =6.0`$. A detailed comparison with these results will be given below in section VI. In this context it is interesting to note, the slope of the hyperfine splitting as a function of the quark mass turns out to be to small in light hadron spectroscopy as well .

The calculations listed above are performed in the quenched approximation, which could be a factor contributing to the wrong slope. In potential model language, which is not necessarily appropriate here, the hyperfine splitting is related to the square of the wave-function at the origin. This in turn depends on the light quark mass and is independent of $`m_Q`$ as $`m_Q\mathrm{}`$. The wrong slope could then reflect the fact that the wave-function at the origin is not increasing rapidly enough as the light quark mass increases. This is natural in the quenched approximation as the potential at the origin is weakened by the coupling constant running too quickly to zero at short distance.

An alternative scenario is one in which the coefficients of the relevant terms in the action, here $`c_{\mathrm{sw}}`$ in the light quark action, effectively carry some quark mass dependence that has not been included, leading to an underestimate of the hyperfine splitting at large $`am_q`$. In this case the effect would disappear as $`a`$ is reduced and this seems to be contradicted by results on finer lattices .

Another cause could be a problem in the chiral extrapolation itself. The experimental result for hyperfine splitting $`J/\psi \eta _c`$ in the charmonium system is smaller than the hyperfine splitting for the $`D`$ and $`D_s`$-mesons. If one considers charmonium as a $`D_c`$-meson, one has to conclude that there is a maximum of the hyperfine splitting as a function of the light quark mass for $`m_q<m_c`$. If this maximum is attained for $`m_q<m_s`$, our observation of a negative slope of the hyperfine for $`m_qm_s`$ would be in agreement with nature. In this case extrapolations from the strange region to lighter $`u`$ and $`d`$-quarks as well as the chiral limit would be impossible.

Our final results for the hyperfine splitting at $`\kappa _c`$ and $`\kappa _s`$ are given in table XXIII. In figure 13 we display the dependence of the hyperfine splitting on the spin-averaged heavy-light meson mass. A linear fit in $`m_{\mathrm{sav}}^1`$ for the five heaviest $`m_{\mathrm{sav}}`$ values gives reasonable values of Q. In table XXV we give the numerical outcome of this fit for the strange and non-strange hyperfine splittings. As expected from HQET the intercept always turns out to be zero within statistical errors.

In the H run for $`\beta =6.2`$ we determine the hyperfine splitting from ratio-fits. In order to determine this without excited state contamination we use a fit interval for which both of the individual correlators have reached a plateau. As noted above, no significant dependence on the light quark mass was observed and we were able to fit the results to a constant with reasonable values of $`Q`$. The result is given in table XXIV.

In case of the N and P run we determined the hyperfine splitting from the jackknife difference of masses obtained from the pseudo-scalar and vector meson propagator. In the N run we compared the outcome for the different smearings available, for different values of the starting point $`t_{\mathrm{min}}`$ of the fit range. An example is shown in figure 15. The different smearing functions lead to compatible answers for the hyperfine splitting. We use the outcome from the propagators with sink and source smearing $`\varphi _{\mathrm{Hg},1}`$ for our final result. In figure 16 we compare the outcome of the different runs at $`\beta =6.2`$. Clearly the outcome from the run N is the most precise. The result for the physical $`B_s`$ hyperfine splitting will be extracted from this results.

Having observed clear signals for the radially excited $`S`$-wave states on our coarse lattice, we also studied their hyperfine splittings. Unfortunately the statistical noise grows rapidly and we observe no clear signal for a non-zero splitting. Our results are given in figure 14, comparing the radially excited state hyperfine splitting to that of the ground state. Although we cannot give a value for the hyperfine splitting of the radially excited $`S`$-state, our results support the expectation that it should be equal to or smaller than the ground state splitting.

### F $`P`$-state fine structure

To extract the $`P`$-state fine structure we investigate the jackknife difference of the masses of the individual channels reported in tables XIX and XX. Because the statistical noise grows rapidly as $`t_{\mathrm{min}}`$ increased, this proved to be delicate. For $`\beta =5.7`$ this is illustrated in figure 17. For the matrix-fit to the $`{}_{}{}^{3}P_{2}^{}E`$ propagators we observe a jump in $`Q`$ for $`t_{\mathrm{min}}=3`$. However in the plot of the fit range dependence of the $`{}_{}{}^{3}P_{2}^{}E{}_{}{}^{3}P_{0}^{}`$ splitting the statistical uncertainty doubles between $`t_{\mathrm{min}}=3`$ and $`5`$. Therefore it is hard to tell whether there is a plateau or not.

We quote final results for $`t_{\mathrm{min}}=5`$, which corresponds to dropping 2 timeslices from the first reasonable $`Q`$-value. With this procedure we obtain a large statistical error and no significant splitting can be resolved. More aggressive fitting would have led to a result compatible with zero but with a statistical error of $`30`$ MeV. We give our final numbers in table XXI. For the splitting we used the same fit range for both channels, which leads to slight deviations from the direct difference of the results in table XIX.

For $`\beta =6.2`$ the situation is easier, as shown in figure 18. The noise on the splitting does not grow as fast as on the coarse lattice, because the $`PS`$ splitting is smaller in lattice units. We observe the first reasonable $`Q`$ values for $`t_{\mathrm{min}}=3`$. Since this is a 6 parameter fit to two propagators, we drop 3 time slices and quote the final result for $`t_{\mathrm{min}}=6`$. The results are given in table XXII. Here we also quote results for the splitting of the $`J=1`$ channels to the $`J=0`$ state. We reiterate that no significance should be attached to any difference in the results between the $`{}_{}{}^{3}P_{1}^{}`$ and the $`{}_{}{}^{1}P_{1}^{}`$ operators.

## VI The physical meson spectrum

In this section we determine the physical $`B`$ and $`D`$-meson spectrum and investigate scaling by comparing results at different values of the lattice spacing. We also compare with experimental results and other lattice calculations.

### A $`B`$-meson spectrum

At both of our lattice spacings we can simulate the $`b`$-quark directly. Here we discuss our results for the physical $`B`$-spectrum. Together with the findings of  we want to investigate the dependence of the individual splittings on the lattice spacing. The findings are compared both to the existing experimental results and lattice investigations performed by other groups within a similar framework using NRQCD .

The results of reference  are useful in that they work at a larger lattice spacing than we do here. There are a number of problems, however that make their results not directly comparable. For example, they do not use either smeared correlators or standard fitting techniques, and this will give rise to an unknown systematic error. In addition they do not see a difference between fixing the lattice spacing from $`m_\rho `$ using the clover action and from the charmonium $`1P1S`$ splitting. This is clearly seen on finer lattices . If this arises from overestimating $`a^1`$ from $`m_\rho `$ because of discretisation errors, then this is another source of systematic error. In particular, this feeds into the fixing of the bare $`b`$ or $`c`$ quark mass and into hyperfine splittings. Their final result for the splitting does not take into account the effect of any of the uncertainties in the bare quark mass determination. This is particularly important for the hyperfine splitting and causes their errors to be heavily underestimated.

The results of  overlap with, but are not as complete as ours.

Unfortunately there are no results for heavy clover fermions available that we can use. References  quote numbers for the $`B`$ spectrum. The first reference still uses extrapolations from the lighter quark masses into the $`b`$-region. Reference  determines the bare $`b`$-quark mass from heavyonium, which is not suitable for the heavy clover approach at the lattice spacings used . However we will later compare to their findings for the $`D`$-spectrum, since this problem is not so severe for charmonium at the lattice spacings used.

Results from taking the $`b`$-quark as a static source also exist for spin-independent and flavour splittings, which survive in the infinite mass limit, see for example . However we restrict the discussion here to a comparison with results simulated at the physical $`b`$-quark mass directly.

In the following we denote spin-averaged states by an overline.

We summarise our results for the $`B`$-spectrum in the tables XXVI and XXVII. As an example of the splitting between a strange and a non-strange meson we discuss the difference of the pseudo-scalar $`B_sB`$ splitting in figure 19. We observe no scaling violations between our results and the agreement with the experimental value is excellent.

The results from  given in figure 19 include the statistical errors only and are taken from their figure 15. In  additional uncertainties for the average of the results at the two finest lattices are mentioned in the text. The overall agreement with our results is good. They notice an upward jump, however, for their result on their finest lattice. Unfortunately our result from our finest lattice comes with large uncertainties so that we are unable to clarify whether there is any real effect here. Given the lack of scaling violations in the rest of the results, it seems unlikely to us. Table VII confirms that the scale from $`m_\rho `$ used by us and the scale from $`\sigma `$ used in  are very close. The mismatch of scales of $`3\%`$ can be neglected safely.

The results from  are also in agreement with ours. They use the $`K^{}/K`$ ratio to fix the strange quark mass. This reduces their results compared to that using the $`K/\rho `$ ratio. Assuming a shift of $`10`$ to $`20`$ MeV from this would increase the agreement. This is the size of the effect we observe on our coarse lattice from fixing the strange from $`\varphi /\rho `$ instead of the $`K`$.

Figure 20 shows the scaling of the radially excited $`B_s`$-meson. As discussed in subsection V B already, the extraction of a result for our finest value of $`a`$ turned out to be more problematic than anticipated, and we are left with quite large statistical uncertainties. However the final result is in good agreement with the result from our coarser lattice as well as the result of .

We also included a preliminary experimental result from the DELPHI collaboration for an admixture of the non-strange $`B^{}B`$ and $`B^{}{}_{}{}^{}B^{}`$ splitting . Assuming the hyperfine splitting of the two states to be of similar size, which our findings support, we observe reasonable agreement here. Table XXVI also contains results for the radial excitation energy of the vector state and the spin-averaged $`S`$-wave.

The orbital excitations are compared in figure 21. The lattice results for the spin-averaged strange $`P`$ state scale very well. The magnitude agrees nicely with the $`B_{sJ}^{}(5850)`$ resonance, which is expected to be an admixture of the two $`j_l=\frac{3}{2}`$ states.

Our results for the radially excited $`P`$-states are compared in figure 22. This is the first observation of a signal for these states in a lattice calculation. As the figure shows we get consistent results from the two different lattice spacings investigated. To the best of our knowledge radially excited $`P`$-states have not been observed yet experimentally.

The splittings discussed above are all essentially light quark quantities, which survive into the static limit. Their scaling or non-scaling says more about the light quark action than the heavy quark sector. The hyperfine splitting is of a different nature and from its scaling behaviour one can learn about how well the heavy quarks are being described on the lattice. We display this in figure 23. Our result for the strange hyperfine splitting together with the findings from  shows good scaling.

However the result is much smaller than the experimental value. Since the leading term in the hyperfine splitting arises from the $`𝝈𝑩`$ term in the action, eq. (8), the result for the splitting is sensitive to the coefficient $`c_4`$ and the inclusion of radiative corrections beyond tadpole improvement is required. Preliminary calculations  indicate that the inclusion of the 1-loop corrections would increase the hyperfine splittings on the order of $`10\%`$ for the lattice spacings used. The quenched approximation might also play a rôle here, since in light spectroscopy the hyperfine splittings turn out to be too small as well, see  for a review. This effect increases with increasing quark mass. Unfortunately reference , which investigates the effect of the inclusion of two flavours of dynamical quarks on the $`B`$-meson, does not give any evidence for an increase of the $`B^{}B`$ splitting due to sea-quark vacuum polarisation effects.

From the experience  in $`\mathrm{{\rm Y}}`$-spectroscopy using NRQCD, one could have expected to observe scaling violations in the hyperfine splitting. Using the same heavy quark Hamiltonian as we do,  reports an increase of $`50\%`$ for the $`\mathrm{{\rm Y}}\eta _b`$ splitting, within the range from $`\beta =5.7`$ up to $`\beta =6.2`$. The leading discretisation correction for the hyperfine splitting is $`𝒪((ap_{\mathrm{gluon}})^2)`$ . Typical gluon momenta for the $`\mathrm{{\rm Y}}`$-system are $`1`$ GeV, while for the $`B`$-system they are $`𝒪(\mathrm{\Lambda }_{\mathrm{QCD}})`$. From this one expects reduced scaling violations of $`10\%`$ in the $`B`$-system for our range of lattice spacings. This is the same size as our uncertainties on the hyperfine splitting and therefore consistent with the fact that no scaling violations show up in figure 23.

Results from  are for the chirally extrapolated splitting $`B^{}B`$ with statistical errors only, taken from their figure 17. In the text, the authors quote a result for the strange hyperfine splitting from the average of the two finest lattices, which is $`3`$ MeV higher than the same average for the non-strange hyperfine splitting. An upwards shift of $`3`$ MeV increases the already excellent agreement even further.

The results of , on the other hand, exhibit a clear disagreement to our findings as well as the findings of . We believe that this is because they have determined the bare $`b`$-quark mass $`am_b`$ using heavyonium.

The fine structure of the $`P`$-states is the last topic of this section. Unfortunately we have not been able to resolve this clearly on our coarsest lattice. The situation is displayed in figure 24, for the three sublevels which were resolved at $`\beta =6.0`$ and $`6.2`$. To investigate whether there is evidence for scaling violations in the fine structure we calculate the jackknifed difference of the highest and lowest state. This is shown in figure 25. The error bars turn out to be large and the figure is inconclusive. A more aggressive fit on the coarser lattice, as discussed in subsection V E would lead to the conclusion that scaling violations were seen, but we believe that further work is needed to resolve this question.

Our results for the $`B`$-meson spectrum, together with those of  do not show signs of residual lattice spacing dependence within the achieved accuracy. Therefore we can average the results for the different values of the lattice spacing $`a`$ to obtain our final results on the quenched $`B`$-meson excitation spectrum. The averages were obtained in the following way, for each value of $`a`$ we add the different uncertainties in quadrature to obtain a single value. Here we omitted those sources of uncertainty which are associated with the quenched approximation. These are the uncertainty arising from fixing the strange quark mass from different physical quantities and in the case of the results of  the additional uncertainty of the lattice spacing $`a`$ associated with the physical quantity used to fix $`a`$. At this step we also symmetrised with respect to unsymmetric uncertainties. The central values have been obtained from fitting the results to a constant in an uncorrelated fit with the above described uncertainties. This puts more weight on the more precise results than a simple average.

Our analysis at the individual values of $`a`$ does not include an uncertainty for the residual effect of the lattice spacing. Our actions are improved to $`𝒪(\alpha _sa,a^2)`$. Therefore for each value of $`a`$ we add the maximum of $`\alpha _sa\mathrm{\Lambda }_{\mathrm{QCD}}`$ and $`a^2\mathrm{\Lambda }_{\mathrm{QCD}}^2`$ in quadrature to the uncertainty used in the fit. We quote the smallest of these three so obtained uncertainties as our final uncertainty for the quenched $`B`$-meson spectrum. This way we quote an accuracy which is of the same size as the one we checked for scaling violations. Determining the final uncertainty from the $`\chi ^2`$ of the fit would reduce the uncertainty beyond this level. This procedure also ensures that residual lattice spacing artifacts are properly included if the achieved accuracy differs over the three individual results and the average is largely determined by the coarser lattices.

Our final result on the $`B`$-meson splitting spectrum in the quenched approximation is given in table XXVIII and figure 26.

The question of the effect of quenching on the spectrum goes beyond the scope of this paper. We refer the reader to . There the effects of the inclusion of 2 flavours of dynamical quarks on the spectrum in NRQCD have been investigated and compared to the findings of . No significant difference between the quenched and $`n_f=2`$ results for the $`1S`$, $`2S`$ and $`1P`$-states was found. In particular, as mentioned earlier, no significant sea quark effects were seen in the $`1S`$-hyperfine splitting. Since our investigations confirm scaling in the quenched heavy-light spectrum the conclusions of  are unchanged.

### B $`D`$-meson spectrum

In this section we discuss the $`D`$-meson spectrum and compare our result to existing lattice results as well as to experiment.

The convergence of the NRQCD expansion is particularly important in the $`D`$-range, where the expansion parameter $`\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q\frac{1}{4}`$. Useful results on the question of the convergence are contained in . There the authors study the complete NRQCD action to $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$. Here we include all terms up to $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^2)`$ and the relativistic correction to the kinetic energy in $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$. The authors of  calculate the heavy light kinetic masses using eq. (24) and show that the difference arising from $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$-terms is consistent with the expectation that they are sub-sub-leading in a $`\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q`$ expansion. The changes to the spin-averaged meson mass that we use to fix the quark mass are dominated by the $`\frac{p^4}{m_Q^3}`$ relativistic correction that we include. From this we conclude that remaining $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$ and higher order terms in the NRQCD expansion would only change the physical masses by at most a few percent. This allows us to use the results of  to estimate the changes in the hyperfine splitting which would be produced by these additional terms at fixed bare quark mass.

The authors of  find that the $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^2)`$ terms produce an effect somewhat smaller than a $`\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q`$ expansion might suggest, since they affect the hyperfine splitting indirectly. The only spin-dependent term at $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^2)`$ is a spin-orbit type interaction. At $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$ most terms produce a change of a few percent, but the one which is directly related to the $`𝝈𝑩`$-term: $`\{𝐃^2,𝝈𝑩\}`$ reduces the hyperfine splitting at the charm by $`20\%`$. However when including the other operators of $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$, the second largest effect comes from the $`𝝈(𝑬\times 𝑬+𝑩\times 𝑩)`$ operator, which is also spin-dependent and works in the opposite direction to the other one. The total effect of $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$ is below $`10\%`$. This is of the size of the naïve expectation for the suppression with respect to the leading term and not at all inconsistent with good convergence of the NRQCD expansion. Since we do not include these terms, we conclude, that we may be overestimating the quenched lattice hyperfine splitting of the $`D`$-meson by $`10\%`$.

Our results on the $`D`$-meson spectrum are summarised in table XXIX and in figure 27. The overall agreement to the experimentally observed spectrum is good. We will now discuss the individual splittings in more detail. We also compare our results to the lattice studies of . All of these results use the quenched approximation as well. The publications  apply the heavy clover approach , which has quite different systematic uncertainties from NRQCD for charm quarks.

The flavour dependent $`D_s^{()}D^{()}`$ splittings are in good agreement with the experimental results. Here it is interesting to note that our results reflect the increase of $`10`$ MeV from the $`B`$ to the $`D`$-meson system, which can already be expected from the good agreement of the slope with the experimental outcome in table XXV.

In figure 28 we compare our result for the strange to non-strange spin-averaged splitting to other lattice calculations. In order not to disguise other possible effects, we excluded the uncertainty of the strange quark mass from the plot. The results from  use the $`K^{}/K`$ ratio to define the strange quark. This should shift the results downwards, compared to fixing $`\kappa _s`$ from the $`K/\rho `$ ratio as used for the other results. The implications have already been discussed in the previous sub-section VI A. The uncertainties again allow for an upwards shift of these results by $`10`$ to $`20`$ MeV. It should be noted that we combined the results for the hyperfine splittings and the pseudo-scalar $`D_sD`$ splitting from  to obtain the spin-averaged splitting.

Because of different systematic uncertainties, it is particularly interesting to compare to the heavy clover results of . The results obtained with the use of the $`m_\rho `$-scale, which is the same as what we use, agree very well with ours. This agreement is expected since  uses the same light quark and gauge field action and this quantity is essentially determined by the light quarks and the gluon field. We conclude that the results for the flavour dependent $`D_sD`$ splitting agree well between the different approaches and, within the accuracy achieved, agree well with the experimental result.

For radially excited $`D_s^{()}^{}`$-mesons no experimental results are known to us. However the DELPHI collaboration reports on the observation of the non-strange $`D^{}^{}`$ . This result is still awaiting confirmation by the OPAL and the CLEO collaboration and its interpretation is disputed on the ground of its small experimentally observed width . The splitting between the DELPHI result and the $`D^{}`$ has a similar size to our $`D_s^{}{}_{}{}^{}D_s^{}`$ splitting. Reference  reports lattice results from the heavy-clover approach. From their plot we read $`D_s^{}D_s840(160)`$ MeV, which is in agreement with our findings. However this includes what the authors call “*continuum*” extrapolation out of a regime where the expansion parameter $`am_Q=𝒪(1)`$ is not small. We would prefer to compare to the unextrapolated results at the individual values of $`a`$.

Experimentally the only well established charmed $`P`$-states in the particle data book  are those which are expected to correspond to the states of total light angular momentum $`j_l=\frac{3}{2}`$. Recently the CLEO collaboration  claimed the observation of the $`D_1`$ state corresponding to $`j_l=\frac{1}{2}`$. CLEO gives a preliminary result of $`D_1=2461(_{34}^{+41})(10)(32)`$ MeV, which is slightly heavier but compatible in error bars to the $`D_1^{}=2425`$ MeV <sup>*</sup><sup>*</sup>*We quote the charge-average.. Our lattice calculation delivers the mass of the lighter of the two states. We did not observe a signal for an excited state slightly heavier than this.

In table XXIX our result for the $`\overline{D}_s(1P)\overline{D}_s`$ splitting is compared to the spin-average of the $`D_{s1}^{}`$ and $`D_{s2}^{}`$, the $`j_l=\frac{3}{2}`$ states. The agreement is reasonable. Reference  reports on the $`D_{s1}\overline{D}_s`$ splitting from a lattice study with the heavy clover approach. A comparison to our result for this splitting is given in figure 29. When using the same scale obtained from $`m_\rho `$ both lattice results agree very well with each other. The agreement with experiment is also good.

It should be noted however, that the experimental result included in figure 29 is not necessarily the same as ours. The experimental $`P`$-state corresponds to $`j_l=\frac{3}{2}`$. If the CLEO trend is confirmed and the $`D_1`$ is indeed heavier than the $`D_1^{}`$ and the same holds for the $`D_{s1}`$ states, then the lattice result also corresponds to $`j_l=\frac{3}{2}`$. If not the lattice result will correspond to $`j_l=\frac{1}{2}`$, but the two states will be so close, that any mismatch is well covered by the error bars.

Table XXIX contains our final result for the radially excited $`2P`$-state. This is the first result for this state from a lattice simulation.

As in the $`B`$-system the hyperfine splittings are too small when compared to the experimental result. Whereas in the $`B`$-system they were too small by $`40\%`$, here they are low by $`25\%`$. This could reflect a more severe quenching error for $`B`$-mesons. $`B`$-mesons are somewhat smaller states than $`D`$ mesons and probe slightly different scales. This is a sub-leading effect in a heavy quark symmetry picture, however. Alternatively, if the error comes from radiative corrections to the $`c_4`$ coefficient, that would need to increase with $`m_Q`$. That is seen by the authors of .

The large uncertainty of $`20`$ MeV on our result arises from the chiral extrapolations used in the lattice spacing determination and the way in which this feeds into the fixing of the bare quark mass. Naively we expect a doubling of the relative error, because a larger value of $`a`$ requires a smaller value of $`am_Q`$ to deliver the same physical $`m_{\mathrm{sav}}`$. This smaller value of $`am_Q`$ gives a larger hyperfine splitting $`am_{hpf}`$. When converting to physical units it picks up the uncertainty of $`a`$ for the second time. In fact a factor of four is seen because of the flattening of the relation between the mass shift $`\mathrm{\Delta }`$ and bare heavy quark mass $`m_Q`$, see figure 3, as well as the steepening up of the hyperfine splitting curve for large values of $`m_{\mathrm{sav}}`$ in figure 13.

In figure 30 we compare our results to the results from  obtained in NRQCD. We choose their result in $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^2)`$ as most relevant for this comparison. The good agreement with our results is in fact misleading. They fix their $`c`$-quark mass from charmonium instead of the $`D`$. Fixing from the $`D`$ would lead to a larger value of $`am_c`$ and lower hyperfine splitting, see subsection IV C and reference .

It is interesting to compare the result for the hyperfine splitting between NRQCD and heavy clover quarks. This is also done in figure 30 for the $`D_s`$ hyperfine splitting. The heavy clover results of  appear to be higher than our result. However this is a result of their higher choice of scale coming from $`J/\psi `$ instead of from $`m_\rho `$. This is confirmed by the findings of . Using a scale from $`m_\rho `$ gives a result which agrees with ours, using a scale from $`J/\psi `$ agrees with reference . It should be noted that in  the bare quark mass is determined from the $`D`$, where as in  it is determined from charmonium. For the heavy clover approach at $`\beta =6.0`$ these differences are negligible within statistical errors .

As discussed at the beginning of this subsection, the inclusion of the terms $`𝒪((\mathrm{\Lambda }_{\mathrm{QCD}}/m_Q)^3)`$ contributing to the hyperfine splitting would decrease our result by $`10\%`$. However the heavy clover approach requires similar correction terms in the Hamiltonian to achieve this level of accuracy.

The agreement of the NRQCD and the heavy clover result for the spin-dependent hyperfine splitting is encouraging, since the systematic uncertainties are quite different. The light quark content plays only a minor rôle for the hyperfine splitting, which depends essentially only on the heavy quark Hamiltonian. In NRQCD the leading contribution to the hyperfine splitting comes from the $`𝝈𝑩`$ term in the action, whereas for heavy clover this is split between the kinetic hopping term and the clover term $`\sigma _{\nu \rho }F_{\nu \rho }`$, with the latter becoming more important as the lattice spacing becomes coarser. Both these actions give rise to systematic errors in the hyperfine splitting from mass-dependent radiative corrections to coefficients and neglected higher order terms, each at the $`10\%`$ level, so the differences could have been significantly larger than observed.

## VII Discussion

We present an extensive study of the $`B`$ and $`D`$-meson spectrum using NRQCD heavy quarks and clover light quarks in the quenched approximation.

Our results include the splitting between the strange and the non-strange meson, hyperfine splittings, radially and orbitally excited states. For the first time in a lattice calculation we obtained a result on radially excited $`P`$-wave states. For spin-independent splittings we observe good agreement with experimental results. However, our result for the spin-dependent hyperfine splitting turns out to be too low in comparison to experiment. This is a well known effect in quenched hadron spectroscopy. Furthermore, in the present calculation hyperfine splittings are also affected by the neglect of radiative corrections in the matching of lattice NRQCD to continuum QCD.

Using two different values of the lattice spacing in the $`B`$-spectrum together with the results of  allows for a detailed investigation of the residual lattice spacing dependence of our final results. No scaling violations are observed within the achieved accuracy. Of particular interest is the scaling of the $`B_s^{}B_s`$ splitting, which depends heavily on the properties of the heavy quark content of the theory. Here scaling violations could be ruled out with an accuracy of $`10\%`$. The $`P`$-fine structure has not been resolved for all values of the lattice spacings and further work is needed for this quantity.

Our results on the $`B`$-meson spectrum are summarised in table XXVIII and figure 26 together with the findings of . In addition to the uncertainties considered in the analysis at the individual values of $`a`$, the quoted uncertainties also contain an estimate for the residual lattice spacing artifacts of $`𝒪(\alpha _sa,a^2)`$. The table gives our final results for the $`B`$-meson spectrum in the quenched approximation.

Our final results on the $`D`$-meson spectrum are shown in table XXIX and figure 27 above. This is our final result for a lattice spacing of $`a^11.1`$ GeV and does not include an estimate of the residual lattice spacing artifacts of $`𝒪(\alpha _sa,a^2)`$. For the above value of $`a`$, this corresponds to $`13\%`$, which is of similar size to or smaller than the otherwise achieved accuracy.

We compared our results to lattice results of other collaborations obtained with NRQCD or in the heavy clover framework. In general we observe good agreement. Discrepancies which appear at first sight could be traced to underestimated errors in these other results or the use of different scales when converting the lattice results into physical units. The excellent agreement of our results with the results obtained in the heavy clover approach is noteworthy because of the different systematics of these approaches.

These results are the most complete lattice results on the $`B`$ and $`D`$-meson spectrum to date.

### Acknowledgements

We would like to acknowledge useful discussions with Peter Boyle and Peter Lepage.

J.H. was supported by a Marie Curie research fellowship by the European commission under ERB FMB ICT 961729, by PPARC and the National Science Foundation. S.C. acknowledges fellowships from the Royal Society of Edinburgh and the Alexander von Humboldt Stiftung. A.A.K. was supported by the Research for the Future Program of the Japanese Society for the Promotion of Science. J.Sh. would like to thank members of the theoretical physics group at the University of Glasgow for their hospitality during an extended visit. Support from an UK PPARC Visiting Fellowship PPA/V/S/1997/00666, is gratefully acknowledged.

This work was supported by the DOE under DE-FG02-91ER40690, PPARC under GR/L56343 and NATO under CRG/94259.

The gauge configurations at both values of $`\beta `$ and the light quark propagators at $`\beta =6.2`$ have been generously provided by the UKQCD-collaboration. The simulations at $`\beta =5.7`$ have been performed at NERSC supported by DOE and the ones at $`\beta =6.2`$ at EPCC in Edinburgh supported by PPARC.