format: text

text:
# HU-EP-02/45 Properties of heavy quarkonia and 𝑩_𝒄 mesons in the relativistic quark model

## I Introduction

The investigation of the properties of mesons composed of a heavy quark and antiquark ($`b\overline{b}`$, $`c\overline{c}`$, $`c\overline{b}`$) gives very important insight into heavy quark dynamics. Heavy quarkonia have a rich spectroscopy with many narrow states lying under the threshold of open flavour production. Excited states experience different decays among which there are radiative transitions to lower levels. The theoretical analysis shows that many properties of heavy quarkonia, including mass spectra and radiative decay rates, are significantly influenced by relativistic effects. Thus their inclusion is necessary for the correct description of the spectroscopy and the determination of quarkonium wave functions. Radiative decays are the most sensitive to relativistic effects. Some of these decays, which are forbidden in the exact nonrelativistic limit (so-called hindered transitions) due to the orthogonality of initial and final meson wave functions, have decay rates of the same order as the allowed ones. In the relativistic description of mesons an important role is played by properties of the confining quark-antiquark interaction, in particular its Lorentz structure. Thus comparison of theoretical predictions with experimental data can provide valuable information on the form of the confining potential. Such information is of great practical interest, since at present it is not possible to obtain the $`Q\overline{Q}`$ potential in the whole range of distances from the basic principles of quantum chromodynamics (QCD). As it is well known, the growing of the strong coupling constant with distance makes perturbation theory inapplicable at large distances (in the infrared region). In this region it is necessary to account for nonperturbative effects connected with the complicated structure of the QCD vacuum. All this leads to a theoretical uncertainty in the $`Q\overline{Q}`$ potential at large and intermediate distances. It is just in this region of large and intermediate distances that most of the basic meson characteristics are formed.

At present, a vast set of experimental data is available on the masses and different decays of heavy quarkonia. However, not all states predicted by theory have been observed yet, while the others need confirmation. Such missing or unconfirmed states are present both in charmonium ($`2^1S_0`$, $`1^1P_1`$, $`1^1D_2`$, $`1^3D_{2,3}`$) and bottomonium (spin-singlet $`{}_{}{}^{1}S_{0}^{}`$ and $`{}_{}{}^{1}P_{1}^{}`$ states, $`D`$ states). The different possibilities for their experimental observation are proposed and widely discussed in the literature elq ; suz ; kuang ; grnb ; gr ; grD . Radiative transitions from the spin-triplet levels with $`J^{PC}=1^{}`$ to these states as well as their subsequent radiative decays play an important role in these proposals. The missing charmonium states can also be searched in $`B`$ meson decays and identified by their radiative transitions elq ; suz . For this purpose, reliable relativistic predictions for the masses of these states and for the rates of radiative transitions involving them are necessary.

The properties of the $`B_c`$ meson are of special interest, since it is the only heavy meson consisting of two heavy quarks with different flavour. This difference of quark flavours forbids annihilation into gluons. As a result, the excited $`B_c`$ meson states lying below the $`BD`$ production threshold undergo pionic or radiative transitions to the ground pseudoscalar state which then decays weakly. There should be a rather rich set of such narrow states which are considerably more stable than corresponding charmonium or bottomonium states. The CDF collaboration cdf  reported the discovery of the $`B_c`$ ground state in $`p\overline{p}`$ collisions at Fermilab. More experimental data are expected to come in the near future from new hadronic colliders.

The purpose of this paper is to give a detailed analysis of mass spectra and radiative transitions in charmonium, bottomonium and $`B_c`$ mesons with the comprehensive account of the relativistic effects. This will allow one to get valuable information about the Lorentz structure of confining quark interactions from the comparison of obtained predictions with available experimental data. On the other hand, it will indicate the processes in which the missing states can be searched for.

The paper is organized as follows. In Sec. II we describe our relativistic quark model. The expression for the heavy quark-antiquark quasipotential with the account of relativistic (including retardation effects) and one loop radiative corrections is given in Sec. III. There it is applied to the calculation of the charmonium, bottomonium and $`B_c`$ meson mass spectra. In Sec. IV pseudoscalar and vector decay constants of the $`B_c`$ meson are calculated with the account of relativistic corrections and compared with other theoretical predictions. In Sec. V the relativistic expressions for the radiative transition matrix elements in the quasipotential approach are given. They are used for the calculation of the decay rates of radiative M1 and E1 transitions in Secs. VI and VII, respectively. The role of relativistic effects in these transitions is investigated. Special attention is payed to the influence of the Lorentz structure of the quark potential on the relativistic corrections to decay rates. Pure vector and scalar potentials as well as their mixture are considered. The obtained results are compared with available experimental data, and the possibilities for searching the missing states in bottomonium are discussed. Finally, our conclusions are given in Sec. VIII.

## II Relativistic quark model

In the quasipotential approach a meson is described by the wave function of the bound quark-antiquark state, which satisfies the quasipotential equation 3  of the Schrödinger type 4

$$\left(\frac{b^2(M)}{2\mu _R}\frac{𝐩^2}{2\mu _R}\right)\mathrm{\Psi }_M(𝐩)=\frac{d^3q}{(2\pi )^3}V(𝐩,𝐪;M)\mathrm{\Psi }_M(𝐪),$$
(1)
where the relativistic reduced mass is

$$\mu _R=\frac{E_1E_2}{E_1+E_2}=\frac{M^4(m_1^2m_2^2)^2}{4M^3},$$
(2)
and $`E_1`$, $`E_2`$ are given by

$$E_1=\frac{M^2m_2^2+m_1^2}{2M},E_2=\frac{M^2m_1^2+m_2^2}{2M}.$$
(3)
Here $`M=E_1+E_2`$ is the meson mass, $`m_{1,2}`$ are the quark masses, and $`𝐩`$ is their relative momentum. In the center of mass system the relative momentum squared on mass shell reads

$$b^2(M)=\frac{[M^2(m_1+m_2)^2][M^2(m_1m_2)^2]}{4M^2}.$$
(4)
The kernel $`V(𝐩,𝐪;M)`$ in Eq. (1) is the quasipotential operator of the quark-antiquark interaction. It is constructed with the help of the off-mass-shell scattering amplitude, projected onto the positive energy states. Constructing the quasipotential of the quark-antiquark interaction, we have assumed that the effective interaction is the sum of the usual one-gluon exchange term with the mixture of long-range vector and scalar linear confining potentials, where the vector confining potential contains the Pauli interaction. The quasipotential is then defined by mass1

$$V(𝐩,𝐪;M)=\overline{u}_1(p)\overline{u}_2(p)𝒱(𝐩,𝐪;M)u_1(q)u_2(q),$$
(5)
with

$$𝒱(𝐩,𝐪;M)=\frac{4}{3}\alpha _sD_{\mu \nu }(𝐤)\gamma _1^\mu \gamma _2^\nu +V_{\mathrm{conf}}^V(𝐤)\mathrm{\Gamma }_1^\mu \mathrm{\Gamma }_{2;\mu }+V_{\mathrm{conf}}^S(𝐤),$$

where $`\alpha _S`$ is the QCD coupling constant, $`D_{\mu \nu }`$ is the gluon propagator in the Coulomb gauge

$$D^{00}(𝐤)=\frac{4\pi }{𝐤^2},D^{ij}(𝐤)=\frac{4\pi }{k^2}\left(\delta ^{ij}\frac{k^ik^j}{𝐤^2}\right),D^{0i}=D^{i0}=0,$$
(6)
and $`𝐤=𝐩𝐪`$; $`\gamma _\mu `$ and $`u(p)`$ are the Dirac matrices and spinors

$$u^\lambda (p)=\sqrt{\frac{ϵ(p)+m}{2ϵ(p)}}\left(\begin{array}{c}1\\ \frac{𝝈𝐩}{ϵ(p)+m}\end{array}\right)\chi ^\lambda ,$$
(7)
with $`ϵ(p)=\sqrt{p^2+m^2}`$. The effective long-range vector vertex is given by

$$\mathrm{\Gamma }_\mu (𝐤)=\gamma _\mu +\frac{i\kappa }{2m}\sigma _{\mu \nu }k^\nu ,$$
(8)
where $`\kappa `$ is the Pauli interaction constant characterizing the anomalous chromomagnetic moment of quarks. Vector and scalar confining potentials in the nonrelativistic limit reduce to

$`V_V(r)`$ $`=`$ $`(1\epsilon )Ar+B,`$
$`V_S(r)`$ $`=`$ $`\epsilon Ar,`$ (9)
reproducing

$$V_{\mathrm{conf}}(r)=V_S(r)+V_V(r)=Ar+B,$$
(10)
where $`\epsilon `$ is the mixing coefficient.

The expression for the quasipotential of the heavy quarkonia, expanded in $`v^2/c^2`$ without and with retardation corrections to the confining potential, can be found in Refs. mass1  and mass , respectively. The structure of the spin-dependent interaction is in agreement with the parameterization of Eichten and Feinberg ef . The quasipotential for the heavy quark interaction with light antiquark without employing the expansion in inverse powers of the light quark mass is given in Ref. hlm . All the parameters of our model such as quark masses, parameters of the linear confining potential $`A`$ and $`B`$, mixing coefficient $`\epsilon `$ and anomalous chromomagnetic quark moment $`\kappa `$ are fixed from the analysis of heavy quarkonium masses mass  (see Sec. III) and radiative decays gf  (see Secs. VVII). The quark masses $`m_b=4.88`$ GeV, $`m_c=1.55`$ GeV and the parameters of the linear potential $`A=0.18`$ GeV2 and $`B=0.16`$ GeV have usual values of quark models. The value of the mixing coefficient of vector and scalar confining potentials $`\epsilon =1`$ has been determined from the consideration of the heavy quark expansion for the semileptonic $`BD`$ decays fg  and charmonium radiative decays gf . Finally, the universal Pauli interaction constant $`\kappa =1`$ has been fixed from the analysis of the fine splitting of heavy quarkonia $`{}_{}{}^{3}P_{J}^{}`$\- states mass1 . Note that the long-range magnetic contribution to the potential in our model is proportional to $`(1+\kappa )`$ and thus vanishes for the chosen value of $`\kappa =1`$. In the present paper we will take into consideration the retardation corrections as well as one-loop radiative corrections.

## III Heavy quarkonium and $`𝑩_𝒄`$ meson mass spectra

The heavy quark-antiquark potential with the account of retardation effects and one loop radiative corrections can be presented in the form of a sum of spin-independent and spin-dependent parts. For the spin-independent part we have mass

$`V_{\mathrm{SI}}(r)`$ $`=`$ $`{\displaystyle \frac{4}{3}}{\displaystyle \frac{\overline{\alpha }_V(\mu ^2)}{r}}+Ar+B{\displaystyle \frac{4}{3}}{\displaystyle \frac{\beta _0\alpha _s^2(\mu ^2)}{2\pi }}{\displaystyle \frac{\mathrm{ln}(\mu r)}{r}}`$ (16)
$`+{\displaystyle \frac{1}{8}}\left({\displaystyle \frac{1}{m_1^2}}+{\displaystyle \frac{1}{m_2^2}}\right)\mathrm{\Delta }\left[{\displaystyle \frac{4}{3}}{\displaystyle \frac{\overline{\alpha }_V(\mu ^2)}{r}}{\displaystyle \frac{4}{3}}{\displaystyle \frac{\beta _0\alpha _s^2(\mu ^2)}{2\pi }}{\displaystyle \frac{\mathrm{ln}(\mu r)}{r}}+(1\epsilon )(1+2\kappa )Ar\right]`$
$`+{\displaystyle \frac{1}{2m_1m_2}}(\{{\displaystyle \frac{4}{3}}{\displaystyle \frac{\overline{\alpha }_V}{r}}[𝐩^2+{\displaystyle \frac{(𝐩𝐫)^2}{r^2}}]\}_W`$
$`{\displaystyle \frac{4}{3}}{\displaystyle \frac{\beta _0\alpha _s^2(\mu ^2)}{2\pi }}\{𝐩^2{\displaystyle \frac{\mathrm{ln}(\mu r)}{r}}+{\displaystyle \frac{(𝐩𝐫)^2}{r^2}}({\displaystyle \frac{\mathrm{ln}(\mu r)}{r}}{\displaystyle \frac{1}{r}})\}_W)`$
$`+\left[{\displaystyle \frac{1\epsilon }{2m_1m_2}}{\displaystyle \frac{\epsilon }{4}}\left({\displaystyle \frac{1}{m_1^2}}+{\displaystyle \frac{1}{m_2^2}}\right)\right]\left\{Ar\left[𝐩^2{\displaystyle \frac{(𝐩𝐫)^2}{r^2}}\right]\right\}_W`$
$`+\left[{\displaystyle \frac{1}{4}}\left({\displaystyle \frac{1}{m_1^2}}+{\displaystyle \frac{1}{m_2^2}}\right)+{\displaystyle \frac{1}{m_1m_2}}\right]B𝐩^2,`$
where

$`\overline{\alpha }_V(\mu ^2)`$ $`=`$ $`\alpha _s(\mu ^2)\left[1+\left({\displaystyle \frac{a_1}{4}}+{\displaystyle \frac{\gamma _E\beta _0}{2}}\right){\displaystyle \frac{\alpha _s(\mu ^2)}{\pi }}\right],`$ (17)
$`a_1`$ $`=`$ $`{\displaystyle \frac{31}{3}}{\displaystyle \frac{10}{9}}n_f,`$ (18)
$`\beta _0`$ $`=`$ $`11{\displaystyle \frac{2}{3}}n_f.`$
Here $`n_f`$ is the number of flavours and $`\mu `$ is a renormalization scale. Note that for a quantity quadratic in the momenta we use the Weyl prescription bcp :

$$\{f(r)p^ip^j\}_W=\frac{1}{4}\{\{f(r),p^i\},p^j\}.$$

For the dependence of the QCD coupling constant $`\alpha _s(\mu ^2)`$ on the renormalization point $`\mu ^2`$ we use the leading order result

$$\alpha _s(\mu ^2)=\frac{4\pi }{\beta _0\mathrm{ln}(\mu ^2/\mathrm{\Lambda }^2)}.$$
(19)
In our numerical calculations we set the renormalization scale $`\mu =2m_1m_2/(m_1+m_2)`$ and $`\mathrm{\Lambda }=0.168`$ GeV, which gives $`\alpha _s=0.314`$ for $`m_1=m_2=m_c`$ (charmonium); $`\alpha _s=0.223`$ for $`m_1=m_2=m_b`$ (bottomonium); and $`\alpha _s=0.265`$ for $`m_1=m_c`$, $`m_2=m_b`$ ($`B_c`$ meson).

The spin-dependent part of the quark-antiquark potential for unequal quark masses with the inclusion of radiative corrections gupta ; ptn  can be presented in our model as follows:

$$V_{\mathrm{SD}}(r)=a𝐋𝐒+b\left[\frac{3}{r^2}(𝐒_1𝐫)(𝐒_2𝐫)(𝐒_1𝐒_2)\right]+c𝐒_1𝐒_2+d𝐋(𝐒_1𝐒_2),$$
(20)
$`a`$ $`=`$ $`{\displaystyle \frac{1}{4}}({\displaystyle \frac{1}{m_1^2}}+{\displaystyle \frac{1}{m_2^2}})\{{\displaystyle \frac{4}{3}}{\displaystyle \frac{\alpha _s(\mu ^2)}{r^3}}(1+{\displaystyle \frac{\alpha _s(\mu ^2)}{\pi }}[{\displaystyle \frac{7}{3}}{\displaystyle \frac{\beta _0}{12}}+\gamma _E({\displaystyle \frac{\beta _0}{2}}3)+{\displaystyle \frac{\beta _0}{2}}\mathrm{ln}(\mu r)`$ (24)
$`3\mathrm{ln}(\sqrt{m_1m_2}r)]){\displaystyle \frac{A}{r}}\}+{\displaystyle \frac{1}{m_1m_2}}{\displaystyle \frac{4}{3}}{\displaystyle \frac{\alpha _s(\mu ^2)}{r^3}}(1+{\displaystyle \frac{\alpha _s(\mu ^2)}{\pi }}[{\displaystyle \frac{1}{6}}{\displaystyle \frac{\beta _0}{12}}`$
$`+\gamma _E({\displaystyle \frac{\beta _0}{2}}{\displaystyle \frac{3}{2}})+{\displaystyle \frac{\beta _0}{2}}\mathrm{ln}(\mu r){\displaystyle \frac{3}{2}}\mathrm{ln}(\sqrt{m_1m_2}r)])`$
$`+\left({\displaystyle \frac{1}{m_1^2}}{\displaystyle \frac{1}{m_2^2}}\right){\displaystyle \frac{\alpha _s^2(\mu ^2)}{2\pi r^3}}\mathrm{ln}{\displaystyle \frac{m_2}{m_1}}+{\displaystyle \frac{1}{2}}\left({\displaystyle \frac{1}{m_1}}+{\displaystyle \frac{1}{m_2}}\right)^2(1+\kappa )(1\epsilon ){\displaystyle \frac{A}{r}},`$
$`b`$ $`=`$ $`{\displaystyle \frac{1}{3m_1m_2}}\{{\displaystyle \frac{4\alpha _s(\mu ^2)}{r^3}}(1+{\displaystyle \frac{\alpha _s(\mu ^2)}{\pi }}[{\displaystyle \frac{29}{6}}{\displaystyle \frac{1}{4}}\beta _0+\gamma _E({\displaystyle \frac{\beta _0}{2}}3)+{\displaystyle \frac{\beta _0}{2}}\mathrm{ln}(\mu r),`$ (27)
$`3\mathrm{ln}(\sqrt{m_1m_2}r)])+(1+\kappa )^2(1\epsilon ){\displaystyle \frac{A}{r}}\},`$
$`c`$ $`=`$ $`{\displaystyle \frac{4}{3m_1m_2}}\{{\displaystyle \frac{8\pi \alpha _s(\mu ^2)}{3}}([1+{\displaystyle \frac{\alpha _s(\mu ^2)}{\pi }}({\displaystyle \frac{5}{12}}\beta _0{\displaystyle \frac{11}{3}}[{\displaystyle \frac{m_1m_2}{m_1+m_2}}+{\displaystyle \frac{1}{8}}{\displaystyle \frac{m_1+m_2}{m_1m_2}}]\mathrm{ln}{\displaystyle \frac{m_2}{m_1}})]\delta ^3(r)`$ (31)
$`+{\displaystyle \frac{\alpha _s(\mu ^2)}{\pi }}[{\displaystyle \frac{\beta _0}{8\pi }}^2\left({\displaystyle \frac{\mathrm{ln}(\mu r)+\gamma _E}{r}}\right)+{\displaystyle \frac{21}{16\pi }}^2\left({\displaystyle \frac{\mathrm{ln}(\sqrt{m_1m_2}r)+\gamma _E}{r}}\right)])`$
$`+(1+\kappa )^2(1\epsilon ){\displaystyle \frac{A}{r}}\},`$
$`d`$ $`=`$ $`{\displaystyle \frac{1}{4}}({\displaystyle \frac{1}{m_1^2}}{\displaystyle \frac{1}{m_2^2}})\{{\displaystyle \frac{4}{3}}{\displaystyle \frac{\alpha _s(\mu ^2)}{r^3}}(1+{\displaystyle \frac{\alpha _s(\mu ^2)}{\pi }}[{\displaystyle \frac{7}{3}}{\displaystyle \frac{\beta _0}{12}}+\gamma _E({\displaystyle \frac{\beta _0}{2}}3)+{\displaystyle \frac{\beta _0}{2}}\mathrm{ln}(\mu r)`$ (34)
$`3\mathrm{ln}(\sqrt{m_1m_2}r)]){\displaystyle \frac{A}{r}}2(1+\kappa )(1\epsilon ){\displaystyle \frac{A}{r}}\}+({\displaystyle \frac{1}{m_1}}+{\displaystyle \frac{1}{m_2}})^2{\displaystyle \frac{\alpha _s^2(\mu ^2)}{2\pi r^3}}\mathrm{ln}{\displaystyle \frac{m_2}{m_1}},`$
where $`𝐋`$ is the orbital momentum and $`𝐒_{1,2}`$, $`𝐒=𝐒_1+𝐒_2`$ are the spin momenta. For the equal mass case ($`m_1=m_2=m`$) the second order in $`\alpha _s`$ contribution of the annihilation diagrams

$$\delta c=\frac{8\alpha _s^2(\mu ^2)}{3m^2}\left(1\mathrm{ln}2\right)\delta ^3(r)$$
(35)
must be added to the spin-spin interaction coefficient $`c`$ in Eq. (31).

The correct description of the fine structure of the heavy quarkonium mass spectrum requires the vanishing of the vector confinement contribution. This can be achieved by setting $`1+\kappa =0`$, i.e., the total long-range quark chromomagnetic moment equals zero, which is in accord with the flux tube buch  and minimal area bcp ; bv  models. One can see from Eq. (20) that for the spin-dependent part of the potential this conjecture is equivalent to the assumption about the scalar structure of confining interaction scal .

To calculate the heavy meson mass spectra with the account of all relativistic corrections (including retardation effects) of order $`v^2/c^2`$ and one-loop radiative corrections we substitute the quasipotential which is a sum of the spin-independent (16) and spin-dependent (20) parts into the quasipotential equation (1). Then we multiply the resulting expression from the left by the quasipotential wave function of a bound state and integrate with respect to the relative momentum. Taking into account the accuracy of the calculations, we can use for the resulting matrix elements the wave functions of Eq. (1) with the static potential 111This static potential includes also some radiative corrections. The remaining radiative correction term with logarithm in Eq. (16), also not vanishing in the static limit, is treated perturbatively.

$$V_{\mathrm{NR}}(r)=\frac{4}{3}\frac{\overline{\alpha }_V(\mu ^2)}{r}+Ar+B.$$
(36)
As a result we obtain the mass formula

$$\frac{b^2(M)}{2\mu _R}=W+a𝐋𝐒+b\left[\frac{3}{r^2}(𝐒_1𝐫)(𝐒_2𝐫)(𝐒_1𝐒_2)\right]+c𝐒_1𝐒_2+d𝐋(𝐒_1𝐒_2),$$
(37)
where

$$W=V_{\mathrm{SI}}+\frac{𝐩^2}{2\mu _R}.$$
(38)
The first term on the right-hand side of the mass formula (37) contains all spin-independent contributions, the second and the last terms describe the spin-orbit interaction, the third term is responsible for the tensor interaction, while the forth term gives the spin-spin interaction. The last term is not zero only for the unequal mass case $`m_1m_2`$ and leads to the mixing of triplet and singlet meson states with the total angular momentum $`J`$ equal to the orbital momentum $`L`$.

In Table 1 the calculated charmonium mass spectrum is compared with experimental data. For meson states we use the notation $`n{}_{}{}^{2S+1}L_{J}^{}`$, where $`n=n_r+1`$ and $`n_r`$ is the radial quantum number. Our predictions agree with PDG pdg  data within a few MeV. Our model correctly reproduces both the position of the levels and their fine and hyperfine splittings. In this Table we also give the recent Belle Collaboration data belle  on pseudoscalar $`\eta _c(1{}_{}{}^{1}S_{0}^{})`$ and $`\eta _c^{}(2{}_{}{}^{1}S_{0}^{})`$ states observed in $`B`$ decays. The measured mass of the ground spin-singlet state $`\eta _c(1{}_{}{}^{1}S_{0}^{})`$ is in good agreement with world averages and predictions of our model, while the radially excited $`\eta _c^{}(2{}_{}{}^{1}S_{0}^{})`$ state lies considerably higher than previous experimental indications and most of the theoretical predictions. If these data are confirmed, it will be difficult to accommodate such a small hyperfine splitting $`32`$ MeV (almost four times smaller than 117 MeV splitting for the ground state) in the framework of the quark model. 222The position of the $`\mathrm{\Psi }^{}(2{}_{}{}^{3}S_{1}^{})`$ can in principle be influenced by the nearby threshold of the open charm production.

Our prediction for the mass of $`h_c(1{}_{}{}^{1}P_{1}^{})`$ is consistent with the data from the Fermilab Experiment E760 E760  on $`p\overline{p}h_c(3526)\pi ^0J/\mathrm{\Psi }`$ which, however, need confirmation. The same is true for the indication of a $`1{}_{}{}^{3}D_{2}^{}`$ state with mass $`3836\pm 13`$ MeV in $`\pi ^\pm NJ/\mathrm{\Psi }\pi ^+\pi ^{}+`$anything E705 .

For the calculation of the bottomonium mass spectrum it is also necessary to take into account additional one-loop corrections due to the finite mass of the charm quark hm ; h ; m ; bsv . We considered these corrections within our model in Ref. efgmc  and found that they give contributions of a few MeV and are weakly dependent on the quantum numbers of the bottomonium states. The one-loop correction to the static $`Q\overline{Q}`$ potential in QCD due to the finite $`c`$ quark mass is given by m ; efgmc

$$\mathrm{\Delta }V(r,m_c)=\frac{4}{9}\frac{\alpha _s^2(\mu )}{\pi r}\left[\mathrm{ln}(\sqrt{a_0}m_cr)+\gamma _E+E_1(\sqrt{a_0}m_cr)\right],$$
(39)
where

$$E_1(x)=_x^{\mathrm{}}e^t\frac{dt}{t}=\gamma _E\mathrm{ln}x\underset{n=1}{\overset{\mathrm{}}{}}\frac{(x)^n}{nn!},$$

$`\gamma _E0.5772`$ is the Euler constant and $`a_0=5.2`$. The resulting bottomonium mass spectrum with the account of this correction is given in Table 2. We found that the small shift of the QCD parameter $`\mathrm{\Lambda }`$ from our previous mass  value 0.178 to 0.168 GeV (with all other parameters remaining fixed) allows us to get a good fit to the bottomonium mass spectrum with the account of finite charm quark mass corrections. The difference between the theoretical and experimental data is less than 3 MeV. Very recently CLEO Collaboration presented cleo  the first evidence for the production of the triplet $`\mathrm{{\rm Y}}(1D)`$ state in the four photon cascades starting from $`\mathrm{{\rm Y}}(3S)`$. In Table 2 we give their preliminary result for the mass of $`\mathrm{{\rm Y}}(1D_2)`$ state which is consistent with our prediction.

For the mesons consisting of quarks with different flavours ($`m_1m_2`$), such as the $`B_c`$ meson, the coefficient $`d`$ (34) in the spin-dependent part of the quark potential (20) is not equal to zero. This results in the mixing of singlet and triplet $`P`$ states with $`J=1`$,

$`nP1^{}`$ $`=`$ $`n{}_{}{}^{1}P_{1}^{}\mathrm{cos}\theta _{nP}+n{}_{}{}^{3}P_{1}^{}\mathrm{sin}\theta _{nP},`$ (40)
$`nP1`$ $`=`$ $`n{}_{}{}^{1}P_{1}^{}\mathrm{sin}\theta _{nP}+n{}_{}{}^{3}P_{1}^{}\mathrm{cos}\theta _{nP},`$ (42)
and of $`D`$ states with $`J=2`$,

$`nD2^{}`$ $`=`$ $`n{}_{}{}^{1}D_{2}^{}\mathrm{cos}\theta _{nD}+n{}_{}{}^{3}D_{2}^{}\mathrm{sin}\theta _{nD},`$ (43)
$`nD2`$ $`=`$ $`n{}_{}{}^{1}D_{2}^{}\mathrm{sin}\theta _{nD}+n{}_{}{}^{3}D_{2}^{}\mathrm{cos}\theta _{nD}.`$ (45)
For the $`B_c`$ meson the values of the mixing angles in our model are

$$\theta _{1P}=0.357,\theta _{2P}=0.405,\theta _{1D}=0.627.$$
(46)
In Table 3 we compare our model predictions for the mass spectrum of the $`B_c`$ meson with other quark model results eq ; gklt ; ful ; nl . We see that the differences between the predictions in most cases do not exceed 30 MeV. The only exceptions are masses of $`1D`$ states, which are $`5070`$ MeV heavier in our model. The fine and hyperfine splittings are also consistent with each other. All these predictions for the ground state pseudoscalar $`B_c`$ and vector $`B_c^{}`$ meson masses satisfy the bounds found by Kwong and Rosner kr :

$$6.194\mathrm{GeV}<M_{B_c}<6.292\mathrm{GeV}$$

and

$$6.284\mathrm{GeV}<M_{B_c^{}}<6.357\mathrm{GeV}.$$

In Ref. bvbc  the ground state $`B_c`$ mass was evaluated in perturbative QCD. Experimental data pdg  at present are available only for the $`B_c`$ ground state and have large error bars $`M_{B_c}=6.4\pm 0.4`$ GeV.

In the following sections we apply the masses and wave functions of $`\mathrm{\Psi }`$, $`\mathrm{{\rm Y}}`$ and $`B_c`$ mesons for the calculation of their decay constants and decay rates.

## IV Pseudoscalar and vector decay constants of the $`𝑩_𝒄`$ meson

The $`B_c`$ meson and its first excitations which lie below the $`BD`$ threshold are stable against strong decays, since they cannot annihilate into gluons. They can decay via electromagnetic and pionic transitions into the lightest pseudoscalar ground state $`B_c`$. The significant contribution to the $`B_c`$ total decay rate comes from the annihilation of the $`c`$ quark and $`\overline{b}`$ antiquark into the vector boson $`W^+`$ which decays into a lepton and a neutrino or a quark-antiquark pair. The weak annihilation decay rate is determined by the pseudoscalar constant of the $`B_c`$ meson.

The decay constants $`f_P`$ and $`f_V`$ of the pseudoscalar ($`P`$) and vector ($`V`$) mesons parameterize the matrix elements of the weak current between the corresponding meson and the vacuum. In the case of the $`B_c`$ meson they are defined by

$`0|\overline{b}\gamma ^\mu \gamma _5c|P(𝐊)`$ $`=`$ $`if_PK^\mu ,`$ (47)
$`0|\overline{b}\gamma ^\mu c|V(𝐊,\epsilon )`$ $`=`$ $`f_VM_V\epsilon ^\mu ,`$ (48)
where $`𝐊`$ is the meson momentum, $`P`$ corresponds to the pseudoscalar $`B_c`$ and $`V`$ to the vector $`B_c^{}`$ mesons, $`\epsilon ^\mu `$ and $`M_V`$ are the polarisation vector and mass of the vector $`B_c^{}`$ meson.

In the relativistic quark model the decay constants can be expressed through the meson wave function $`\mathrm{\Phi }_{P,V}(p)`$ in the momentum space and are given by gfm

$`f_{P,V}`$ $`=`$ $`\sqrt{{\displaystyle \frac{12}{M_{P,V}}}}{\displaystyle \frac{d^3p}{(2\pi )^3}\left(\frac{ϵ_c(p)+m_c}{2ϵ_c(p)}\right)^{1/2}\left(\frac{ϵ_b(p)+m_b}{2ϵ_b(p)}\right)^{1/2}}`$ (50)
$`\times \left\{1+\lambda _{P,V}{\displaystyle \frac{𝐩^2}{[ϵ_c(p)+m_c][ϵ_b(p)+m_b]}}\right\}\mathrm{\Phi }_{P,V}(p),`$
with $`\lambda _P=1`$ and $`\lambda _V=1/3`$. In the nonrelativistic limit $`p^2/m^20`$ these expressions for decay constants give the well-known formula

$$f_P^{\mathrm{NR}}=f_V^{\mathrm{NR}}=\sqrt{\frac{12}{M_{P,V}}}\left|\mathrm{\Psi }_{P,V}(0)\right|,$$
(51)
where $`\mathrm{\Psi }_{P,V}(0)`$ is the meson wave function at the origin $`r=0`$.

The calculated values of the pseudoscalar and vector decay constants of the $`B_c`$ meson in our model using the relativistic formula (50) are displayed in Table 4. They are compared with the ones calculated using the nonrelativistic expression (51) and other predictions of the nonrelativistic quark models eq ; ful , QCD sum rules gklt  and lattice NRQCD jw . We see that inclusion of relativistic corrections reduces the pseudoscalar decay constant $`f_{B_c}`$ by 20% and produces the difference between vector and pseudoscalar decay constants of approximately 70 MeV. The calculated values of these decay constants are consistent with lattice jw  and QCD sum rule gklt  predictions.

## V Radiative transitions in heavy quarkonia and $`𝑩_𝒄`$ mesons

To determine the rates of radiative decays ($`BA+\gamma `$) it is necessary to calculate the matrix element of the electromagnetic current $`J_\mu `$ between the initial ($`B`$) and final ($`A`$) meson states. In the quasipotential approach such matrix element has the form f

$$A|J_\mu (0)|B=\frac{d^3pd^3q}{(2\pi )^6}\overline{\mathrm{\Psi }}_{A𝐏}(𝐩)\mathrm{\Gamma }_\mu (𝐩,𝐪)\mathrm{\Psi }_{B𝐐}(𝐪),$$
(52)
where $`\mathrm{\Gamma }_\mu (𝐩,𝐪)`$ is the two-particle vertex function and $`\mathrm{\Psi }_{A,B}`$ are the meson wave functions projected onto the positive energy states of quarks and boosted to the moving reference frame. The contributions to $`\mathrm{\Gamma }`$ come from Figs. 1 and 2. The contribution $`\mathrm{\Gamma }^{(2)}`$ is the consequence of the projection onto the positive-energy states. Note that the form of the relativistic corrections resulting from the vertex function $`\mathrm{\Gamma }^{(2)}`$ explicitly depends on the Lorentz structure of the $`Q\overline{Q}`$-interaction. Thus the vertex function is given by

$$\mathrm{\Gamma }_\mu (𝐩,𝐪)=\mathrm{\Gamma }_\mu ^{(1)}(𝐩,𝐪)+\mathrm{\Gamma }_\mu ^{(2)}(𝐩,𝐪)+\mathrm{},$$
(53)
where

$$\mathrm{\Gamma }_\mu ^{(1)}(𝐩,𝐪)=e_1\overline{u}_1(p_1)\gamma _\mu u_1(q_1)(2\pi )^3\delta (𝐩_2𝐪_2)+(12),$$
(54)
and

$`\mathrm{\Gamma }_\mu ^{(2)}(𝐩,𝐪)`$ $`=`$ $`e_1\overline{u}_1(p_1)\overline{u}_2(p_2)\{𝒱(𝐩_2𝐪_2){\displaystyle \frac{\mathrm{\Lambda }_1^{()}(k_1^{})}{ϵ_1(k_1^{})+ϵ_1(q_1)}}\gamma _1^0\gamma _{1\mu }`$ (55)
$`+\gamma _{1\mu }{\displaystyle \frac{\mathrm{\Lambda }_1^{()}(k_1)}{ϵ_1(k_1)+ϵ_1(p_1)}}\gamma _1^0𝒱(𝐩_2𝐪_2)\}u_1(q_1)u_2(q_2)+(12).`$
Here $`e_{1,2}`$ are the quark charges, $`𝐤_1=𝐩_1𝚫;𝐤_1^{}=𝐪_1+𝚫;𝚫=𝐏𝐐`$;

$$\mathrm{\Lambda }^{()}(p)=\frac{ϵ(p)\left(m\gamma ^0+\gamma ^0(𝜸𝐩)\right)}{2ϵ(p)},ϵ(p)=\sqrt{p^2+m^2},$$

and

$`p_{1,2}`$ $`=`$ $`ϵ_{1,2}(p){\displaystyle \frac{p_A}{M_A}}\pm {\displaystyle \underset{i=1}{\overset{3}{}}}n^{(i)}(p_A)p^i,`$
$`q_{1,2}`$ $`=`$ $`ϵ_{1,2}(q){\displaystyle \frac{p_B}{M_B}}\pm {\displaystyle \underset{i=1}{\overset{3}{}}}n^{(i)}(p_B)q^i,`$
where $`n^{(i)}`$ are three four-vectors given by

$$n^{(i)\mu }(p)=\{\frac{p^i}{M},\delta _{ij}+\frac{p^ip^j}{M(E+M)}\},E=\sqrt{𝐩^2+M^2},i,j=1,2,3,$$

$`p_B=(E_B,𝐐)`$ and $`p_A=(E_A,𝐏)`$ are four-momenta of initial and final mesons.

It is important to note that the wave functions entering the current matrix element (52) cannot be both in the rest frame. In the initial $`B`$ meson rest frame, the final $`A`$ meson is moving with the recoil momentum $`𝚫`$. The wave function of the moving $`A`$ meson $`\mathrm{\Psi }_{A𝚫}`$ is connected with the wave function in the rest frame $`\mathrm{\Psi }_{A\mathbf{\hspace{0.17em}0}}\mathrm{\Psi }_A`$ by the transformation f

$$\mathrm{\Psi }_{A𝚫}(𝐩)=D_1^{1/2}(R_{L_𝚫}^W)D_2^{1/2}(R_{L_𝚫}^W)\mathrm{\Psi }_{A\mathbf{\hspace{0.17em}0}}(𝐩),$$
(56)
where $`R^W`$ is the Wigner rotation, $`L_𝚫`$ is the Lorentz boost from the rest frame to a moving one, and the rotation matrix $`D^{1/2}(R)`$ in the spinor representation is given by

$$\left(\genfrac{}{}{0pt}{}{\mathrm{1\; \hspace{0.17em}0}}{\mathrm{0\; \hspace{0.17em}1}}\right)D_{1,2}^{1/2}(R_{L_𝚫}^W)=S^1(𝐩_{1,2})S(𝚫)S(𝐩),$$
(57)
where

$$S(𝐩)=\sqrt{\frac{ϵ(p)+m}{2m}}\left(1+\frac{𝜶𝐩}{ϵ(p)+m}\right)$$

is the usual Lorentz transformation matrix of the four-spinor.

To calculate the radiative transition matrix element we adopt the following procedure. We substitute the vertex functions $`\mathrm{\Gamma }^{(1)}`$ and $`\mathrm{\Gamma }^{(2)}`$ given by Eqs. (54) and (55) in the decay matrix element (52) and take into account the wave function transformation (56). The resulting structure of this matrix element is rather complicated, because it is necessary to integrate both over $`d^3p`$ and $`d^3q`$. The $`\delta `$ function in expression (54) permits us to perform one of these integrations and thus this contribution can be easily calculated. The calculation of the vertex function $`\mathrm{\Gamma }^{(2)}`$ contribution is more difficult. Here, instead of a $`\delta `$ function, we have a complicated structure, containing the $`Q\overline{Q}`$ interaction potential in the meson. However, we can expand this contribution in powers of the heavy quark velocities $`v^2/c^2`$ and then use the quasipotential equation in order to perform one of the integrations in the current matrix element. It is easy to see that the vertex function $`\mathrm{\Gamma }^{(2)}`$ contributes already at the first order of the $`v^2/c^2`$ expansion.

We consider two main types of radiative transitions:
a) Magnetic dipole (M1) transitions which go with the spin flip of the quark ($`\mathrm{\Delta }S=1`$, $`\mathrm{\Delta }L=0`$) and thus the initial and final states belong to the same orbital excitation but have different spins. Examples of such transitions are vector to pseudoscalar ($`n{}_{}{}^{3}S_{1}^{}n^{}{}_{}{}^{1}S_{0}^{}+\gamma `$, $`nn^{}`$) and pseudoscalar to vector ($`n{}_{}{}^{1}S_{0}^{}n^{}{}_{}{}^{3}S_{1}^{}+\gamma `$, $`n>n^{}`$) meson decays.
b) Electric dipole (E1) transitions in which the orbital quantum number is changed ($`\mathrm{\Delta }L=1`$, $`\mathrm{\Delta }S=0`$) and thus the initial and final states belong to different orbital excitations but have the same spin. Examples of such transitions are $`n{}_{}{}^{3}S_{1}^{}n^{}{}_{}{}^{3}P_{J}^{}\gamma `$ ($`n>n^{}`$) and $`n{}_{}{}^{3}P_{J}^{}n^{}{}_{}{}^{3}S_{1}^{}\gamma `$ ($`nn^{}`$) decays.

## VI Radiative M1 transitions

### VI.1 M1 decay rates

The radiative M1 transition rate is given by gf

$$\mathrm{\Gamma }(BA+\gamma )=\frac{\omega ^3}{3\pi }(2J^{}+1)\left|𝓜_{BA}\right|^2,\mathrm{where}\omega =\frac{M_B^2M_A^2}{2M_B},$$
(58)
$`M_B`$ and $`M_A`$ are the initial and final meson masses, $`J^{}`$ is the total angular momentum of the final meson. The matrix element of the magnetic moment $`𝓜`$ is defined by

$$𝓜_{BA}=\frac{i}{2}\left[\frac{}{𝚫}\times A\left|𝐉(0)\right|B\right]_{𝚫=0},𝚫=𝐏𝐐,$$
(59)
where $`A\left|J_\mu (0)\right|B`$ is the matrix element of the electromagnetic current between initial ($`B`$) and final ($`A`$) meson states with momenta $`𝐐`$ and $`𝐏`$, respectively.

After inserting the vertex functions $`\mathrm{\Gamma }^{(1)}`$ and $`\mathrm{\Gamma }^{(2)}`$ from Eqs. (54) and (55) in the decay matrix element (52) with the account of the wave function transformation (56), we carry out the expansion in inverse powers of the heavy meson masses $`M_{B,A}`$, which are large due to the presence of two heavy quarks $`M_{B,A}m_Q+m_Q^{}`$. Then we calculate the matrix element of the magnetic moment operator (59) and get
(a) for the vector potential

$`𝓜_V`$ $`=`$ $`{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}\overline{\mathrm{\Psi }}_A(𝐩){\displaystyle \frac{e_1}{2ϵ_1(p)}}\{𝝈_1+{\displaystyle \frac{(1\epsilon )(1+2\kappa )[𝐩\times [𝝈_1\times 𝐩]]}{2ϵ_1(p)[ϵ_1(p)+m_1]}}`$ (63)
$`+{\displaystyle \frac{(1\epsilon )(1+\kappa )[𝐩\times [𝝈_2\times 𝐩]]}{ϵ_1(p)[ϵ_2(p)+m_2]}}`$
$`{\displaystyle \frac{ϵ_2(p)}{M_B}}\left(1+(1\epsilon ){\displaystyle \frac{M_Bϵ_1(p)ϵ_2(p)}{ϵ_1(p)}}\right)i\left[𝐩\times {\displaystyle \frac{}{𝐩}}\right]`$
$`+{\displaystyle \frac{1}{2M_B}}[𝐩\times [𝐩\times ({\displaystyle \frac{𝝈_1}{ϵ_1(p)+m_1}}{\displaystyle \frac{𝝈_2}{ϵ_2(p)+m_2}})]]\}\mathrm{\Psi }_B(𝐩)+(12),`$
(b) for the scalar potential

$`𝓜_S`$ $`=`$ $`{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}\overline{\mathrm{\Psi }}_A(𝐩){\displaystyle \frac{e_1}{2ϵ_1(p)}}\{(1+\epsilon {\displaystyle \frac{ϵ_1(p)+ϵ_2(p)M_B}{ϵ_1(p)}})`$ (66)
$`\times \left(𝝈_1{\displaystyle \frac{ϵ_2(p)}{M_B}}i\left[𝐩\times {\displaystyle \frac{}{𝐩}}\right]\right){\displaystyle \frac{\epsilon [𝐩\times [𝝈_1\times 𝐩]]}{2ϵ_1(p)[ϵ_1(p)+m_1]}}`$
$`+{\displaystyle \frac{1}{2M_B}}[𝐩\times [𝐩\times ({\displaystyle \frac{𝝈_1}{ϵ_1(p)+m_1}}{\displaystyle \frac{𝝈_2}{ϵ_2(p)+m_2}})]]\}\mathrm{\Psi }_B(𝐩)+(12).`$
Note that the last terms in Eqs. (63), (66) result from the wave function transformation (56) from the moving reference frame to the rest one. It is easy to see that in the limit $`p/m0`$ the usual nonrelativistic expression for the magnetic moment follows.

Since we are interested in radiative transitions between $`S`$ state (vector and pseudoscalar) mesons, it is possible to evaluate spin matrix elements using the relation $`𝝈_1=𝝈_2`$. Then, taking into account that both quarks are heavy ($`Q`$ and $`Q^{}`$), we further expand Eqs. (63), (66) in inverse powers of the heavy quark mass $`m_Q`$ up to the second order corrections to the leading contribution and get
(a) for the purely vector potential ($`\epsilon =0`$)

$$_V=\frac{e_Q}{2m_Q}\{1\frac{2𝐩^2}{3m_Q^2}+\frac{1+\kappa }{3}\frac{𝐩^2}{m_Q}(\frac{1}{m_Q}\frac{1}{m_Q^{}})\frac{𝐩^2}{6M_B}(\frac{1}{m_Q}+\frac{1}{m_Q^{}})\}(QQ^{}),$$
(67)
(b) for the purely scalar potential ($`\epsilon =1`$)

$$_S=\frac{e_Q}{2m_Q}\{2\frac{M_Bm_Q^{}}{m_Q}+\frac{𝐩^2}{2m_Q}(\frac{1}{m_Q^{}}\frac{1}{3m_Q})\frac{𝐩^2}{6M_B}(\frac{1}{m_Q}+\frac{1}{m_Q^{}})\}(QQ^{}).$$
(68)
Here $`\mathrm{}`$ denotes the matrix element between radial meson wave functions. For these matrix element calculations we use the meson wave functions obtained calculating their mass spectra.

For the quarks of the same flavour ($`m_Q=m_Q^{}`$, $`e_Q=e_{\overline{Q}^{}}`$) and $`\kappa =1`$ these expressions reduce further gf
(a) for the purely vector potential ($`\epsilon =0`$)

$$_V=2\frac{e_Q}{2m_Q}\left\{1\frac{2𝐩^2}{3m_Q^2}\frac{𝐩^2}{3M_Bm_Q}\right\},$$
(69)
(b) for the purely scalar potential ($`\epsilon =1`$)

$$_S=2\frac{e_Q}{2m_Q}\left\{3\frac{M_B}{m_Q}+\frac{𝐩^2}{3m_Q}\frac{𝐩^2}{3M_Bm_Q}\right\}.$$
(70)
### VI.2 Results and discussion

The resulting $`M1`$ radiative decay rates of charmonium, bottomonium and $`B_c`$ are presented in Tables 5-7. In these Tables we give predictions both for allowed ($`n{}_{}{}^{3}S_{1}^{}n^{}{}_{}{}^{1}S_{0}^{}+\gamma `$, $`n^{}=n`$) and hindered ($`n>n^{}`$) decays. For the calculation of allowed decay rates we use expanded expressions (67)–(70). For the hindered transitions, which are strongly suppressed in the nonrelativistic limit due to orthogonality of the initial and final state wave functions, relativistic effects are decisive. Thus for their calculations we use unexpanded expressions (63) and (66). In Tables 5-7 we present the photon energy $`\omega `$, the decay rates calculated discarding all relativistic corrections $`\mathrm{\Gamma }^{\mathrm{NR}}`$, as well as using relativistic expressions for purely vector $`\mathrm{\Gamma }^V`$, for purely scalar $`\mathrm{\Gamma }^S`$ and for the mixture (II) of vector and scalar potentials $`\mathrm{\Gamma }`$ with $`\epsilon =1`$. Note that in all these calculations we use the relativistic wave functions found calculating the meson mass spectra in Sec. III.

The M1 radiative decay rates are very sensitive to relativistic effects. Even for allowed transitions relativistic and nonrelativistic results differ significantly. An important example is the decay $`J/\mathrm{\Psi }\eta _c\gamma `$. It is well known that the nonrelativistic predictions for its rate are more than two times larger than the experimental data. As we see from Table 5, the inclusion of the relativistic corrections for purely scalar or purely vector potentials do not bring theoretical results in agreement with experiment. For the purely scalar potential the decay rate even increases by 15%. On the other hand for the purely vector potential relativistic effects decrease the decay rate by 25%, but such decrease is not enough: the theoretical result still deviates from experimental data by more than $`2\sigma `$555This is compatible with the estimate that relativistic effects can give contributions of order of 20–30% in charmonium. Only for the mixture of vector and scalar potentials (II) we get the necessary decrease of the decay rate which brings theory in agreement with experimental data for the $`J/\mathrm{\Psi }\eta _c\gamma `$ decay rate. For the hindered decay $`\mathrm{\Psi }^{}\eta _c\gamma `$ the decay rate calculated for the mixture of vector and scalar potentials is also in good agreement with experiment while the rates for pure potentials (especially the scalar one) are lower than the experimental value. In Table 5 we give predictions for decays involving the first radial excitation of the pseudoscalar state $`\eta _c^{}(2{}_{}{}^{1}S_{0}^{})`$ as well. Since there are two contradicting experimental measurements of its mass we calculated the rates using both values. The results obtained using recent Belle value belle  are given in parenthesis.

In Table 6 predictions for M1 decay rates and branching fractions of bottomonium are given. Since the hyperfine splitting in bottomonium is predicted to be small (around 60 MeV, see Table 2) the photon energies and hence decay rates of allowed M1 transitions are very small. This is one of the main reasons why no spin-singlet $`S`$-wave levels $`\eta _b(n{}_{}{}^{1}S_{0}^{})`$ have been observed yet. Our results show that relativistic effects for the favored mixture of vector and scalar potentials further decrease the allowed M1 decay rates.

Recently it was argued by Godfrey and Rosner grnb  that hindered transitions could be more favorable for discovering the $`\eta _b`$. Their analysis of different quark model predictions showed that most nonrelativistic models favor the $`\eta _b(1S)`$ production from $`\mathrm{{\rm Y}}(2S)`$ decays, while the account of relativistic corrections makes prospects for discovering $`\eta _b(1S)`$ in $`\mathrm{{\rm Y}}(3S)`$ radiative decays comparable to those in $`\mathrm{{\rm Y}}(2S)`$ decays. Our present relativistic consideration of these decays supports this observation. Indeed we see from Table 6 that relativistic effects significantly increase rates of hindered transitions. The hindered decay $`\mathrm{{\rm Y}}^{\prime \prime }(3S)\eta _b(1S)\gamma `$ has the largest branching fraction $`4.0\times 10^4`$, which is almost 2.7 times larger than the $`\mathrm{{\rm Y}}^{}(2S)\eta _b(1S)\gamma `$ decay branching fraction. Very recently CLEO Collaboration cleonb  searched for $`\eta _b(1S)`$ in such hindered M1 transition from $`\mathrm{{\rm Y}}^{\prime \prime }(3S)`$. No evidence of $`\mathrm{{\rm Y}}^{\prime \prime }(3S)\eta _b(1S)\gamma `$ transitions was found and rather strict upper limits on the branching fraction were set: $`(\mathrm{{\rm Y}}^{\prime \prime }(3S)\eta _b(1S)\gamma )<6\times 10^4`$, which rule out many previous phenomenological predictions reviewed in Ref. grnb . Our model result for the branching fraction of this decay is below but rather close to this experimental upper limit.

In Table 7 we give predictions for decay rates of M1 radiative transitions of the $`B_c`$ meson in our model in comparison with previous nonrelativistic quark model analysis eq ; gklt ; ful . We see that relativistic effects play an important role in $`B_c`$ meson M1 radiative decays. They reduce the rates of allowed decays and increase the rates of hindered transitions. The largest rates are predicted for the latter decays which are increased by relativistic effects almost by the factor of 3 and thus they are an order of magnitude larger than the rates of allowed M1 transitions.

## VII Radiative E1 transitions

### VII.1 E1 decay rates

The radiative E1 transition rate is given by gmf

$$\mathrm{\Gamma }(BA+\gamma )=\frac{\omega ^3}{3\pi }|𝐃_{BA}|^2,\mathrm{where}\omega =\frac{M_B^2M_A^2}{2M_B},$$
(71)
$`M_B`$ and $`M_A`$ are the initial and final meson masses. The matrix element of the electric dipole moment $`𝐃_{BA}`$ is defined by

$$𝐃_{BA}=i\frac{}{𝚫}A\left|J_0(0)\right|B|_{𝚫=\mathrm{𝟎}},𝚫=𝐏𝐐,$$
(72)
where $`A\left|J_\mu (0)\right|B`$ is the matrix element of the electromagnetic current between initial ($`B`$) and final ($`A`$) meson states with momenta $`𝐐`$ and $`𝐏`$, respectively.

We substitute expressions (52)–(55) in the definition of the electric dipole moment (72) and take into account the relativistic transformation of the wave function (56). Then, discarding some terms of order $`v^4/c^4`$ and higher, we get the following expressions for the electric dipole moment $`𝐃_{BA}`$ gmf  (indices 1,2 are changed to $`Q,Q^{}`$)

(a) for the purely vector potential

$`𝐃_V`$ $`=`$ $`{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}\overline{\mathrm{\Psi }}_A(𝐩)e_Q\{i{\displaystyle \frac{ϵ_Q^{}(p)}{M_B}}{\displaystyle \frac{}{𝐩}}{\displaystyle \frac{\left[𝝈_Q\times 𝐩\right]}{2ϵ_Q(p)[ϵ_Q(p)+m_Q]}}(1{\displaystyle \frac{ϵ_Q(p)}{M_B}}`$ (74)
$`+{\displaystyle \frac{2[M_Bϵ_Q(p)ϵ_Q^{}(p)]}{m_Q}}){\displaystyle \frac{\left[𝝈_Q^{}\times 𝐩\right]}{2M_B[ϵ_Q^{}(p)+m_Q^{}]}}\}\mathrm{\Psi }_B(𝐩)(QQ^{}),`$
(b) for the purely scalar potential

$`𝐃_S`$ $`=`$ $`{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}\overline{\mathrm{\Psi }}_A(𝐩)e_Q\{i{\displaystyle \frac{ϵ_Q^{}(p)}{M_B}}{\displaystyle \frac{}{𝐩}}{\displaystyle \frac{\left[𝝈_Q\times 𝐩\right]}{2ϵ_Q(p)[ϵ_Q(p)+m_Q]}}(1{\displaystyle \frac{ϵ_Q(p)}{M_B}}`$ (76)
$`{\displaystyle \frac{2[M_Bϵ_Q(p)ϵ_Q^{}(p)]}{m_Q}}){\displaystyle \frac{\left[𝝈_Q^{}\times 𝐩\right]}{2M_B[ϵ_Q^{}(p)+m_Q^{}]}}\}\mathrm{\Psi }_B(𝐩)(QQ^{}).`$
The operator $`i/𝐩`$ in Eqs. (74), (76) corresponds in the coordinate space to the operator $`𝐫`$. All other terms in these equations are relativistic corrections. Thus in the nonrelativistic limit the standard expression for the electric dipole moment is recovered.

It is easy to see that there are three different structures with respect to the orbital variables in Eqs. (74), (76): $`𝐫`$, $`[(𝝈_𝑸+𝝈_Q^{})\times 𝐩]`$ and $`[(𝝈_𝑸𝝈_Q^{})\times 𝐩]`$. Thus the matrix element of the electric dipole moment for the electromagnetic transition $`nJMLSn^{}J^{}M^{}L^{}S^{}+\gamma `$ can be presented in the form

$$𝐃_{V,S}=n^{}J^{}M^{}L^{}S^{}|𝒜(p^2)𝐫^{V,S}(p^2)[𝐒\times 𝐩]𝒞^{V,S}(p^2)[(𝐒_Q𝐒_Q^{})\times 𝐩]|nJMLS,$$
(77)
where functions $`𝒜(p^2)`$, $`^{V,S}(p^2)`$ and $`𝒞^{V,S}(p^2)`$ up to order $`𝐩^2/m_{Q(Q^{})}^2`$ are given by

$`𝒜(p^2)`$ $`=`$ $`{\displaystyle \frac{e_Qm_Q^{}e_Q^{}m_Q}{M_B}}+{\displaystyle \frac{e_Qm_Qe_Q^{}m_Q^{}}{M_B}}{\displaystyle \frac{𝐩^2}{2m_Qm_Q^{}}},`$ (78)
$`^V(p^2)`$ $`=`$ $`{\displaystyle \frac{e_Q}{4m_Q^2}}\left(1+{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q}}{\displaystyle \frac{7𝐩^2}{4m_Q^2}}{\displaystyle \frac{𝐩^2}{m_Q^{}^2}}\right)`$ (81)
$`{\displaystyle \frac{e_Q^{}}{4m_Q^{}^2}}\left(1+{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q^{}}}{\displaystyle \frac{7𝐩^2}{4m_Q^{}^2}}{\displaystyle \frac{𝐩^2}{m_Q^2}}\right)`$
$`{\displaystyle \frac{e_Q+e_Q^{}}{4M_B}}\left[{\displaystyle \frac{1}{m_Q}}{\displaystyle \frac{1}{m_Q^{}}}{\displaystyle \frac{𝐩^2}{4}}\left({\displaystyle \frac{1}{m_Q^3}}{\displaystyle \frac{1}{m_Q^{}^3}}\right)\right],`$
$`𝒞^V(p^2)`$ $`=`$ $`{\displaystyle \frac{e_Q}{4m_Q^2}}\left(1+{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q}}{\displaystyle \frac{7𝐩^2}{4m_Q^2}}{\displaystyle \frac{𝐩^2}{m_Q^{}^2}}\right)`$ (84)
$`+{\displaystyle \frac{e_Q^{}}{4m_Q^{}^2}}\left(1+{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q^{}}}{\displaystyle \frac{7𝐩^2}{4m_Q^{}^2}}{\displaystyle \frac{𝐩^2}{m_Q^2}}\right)`$
$`{\displaystyle \frac{e_Q+e_Q^{}}{4M_B}}\left[{\displaystyle \frac{1}{m_Q}}+{\displaystyle \frac{1}{m_Q^{}}}{\displaystyle \frac{𝐩^2}{4}}\left({\displaystyle \frac{1}{m_Q^3}}+{\displaystyle \frac{1}{m_Q^{}^3}}\right)\right],`$
$`^S(p^2)`$ $`=`$ $`{\displaystyle \frac{e_Q}{4m_Q^2}}\left(1{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q}}+{\displaystyle \frac{𝐩^2}{4m_Q^2}}+{\displaystyle \frac{𝐩^2}{m_Qm_Q^{}}}\right)`$ (87)
$`{\displaystyle \frac{e_Q^{}}{4m_Q^{}^2}}\left(1{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q^{}}}+{\displaystyle \frac{𝐩^2}{4m_Q^{}^2}}+{\displaystyle \frac{𝐩^2}{m_Qm_Q^{}}}\right)`$
$`{\displaystyle \frac{e_Q+e_Q^{}}{4M_B}}\left[{\displaystyle \frac{1}{m_Q}}{\displaystyle \frac{1}{m_Q^{}}}{\displaystyle \frac{𝐩^2}{4}}\left({\displaystyle \frac{1}{m_Q^3}}{\displaystyle \frac{1}{m_Q^{}^3}}\right)\right],`$
$`𝒞^S(p^2)`$ $`=`$ $`{\displaystyle \frac{e_Q}{4m_Q^2}}\left(1{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q}}+{\displaystyle \frac{𝐩^2}{4m_Q^2}}+{\displaystyle \frac{𝐩^2}{m_Qm_Q^{}}}\right)`$ (90)
$`+{\displaystyle \frac{e_Q^{}}{4m_Q^{}^2}}\left(1{\displaystyle \frac{2(Mm_Qm_Q^{})}{m_Q^{}}}+{\displaystyle \frac{𝐩^2}{4m_Q^{}^2}}+{\displaystyle \frac{𝐩^2}{m_Qm_Q^{}}}\right)`$
$`{\displaystyle \frac{e_Q+e_Q^{}}{4M_B}}\left[{\displaystyle \frac{1}{m_Q}}+{\displaystyle \frac{1}{m_Q^{}}}{\displaystyle \frac{𝐩^2}{4}}\left({\displaystyle \frac{1}{m_Q^3}}+{\displaystyle \frac{1}{m_Q^{}^3}}\right)\right].`$
The last structure in Eq. (77) proportional to $`[(𝐒_Q𝐒_Q^{})\times 𝐩]`$ leads to the spin-flip transitions. It vanishes for the $`c\overline{c}`$ and $`b\overline{b}`$ mesons, consisting of the quark and antiquark of the same flavour ($`m_Q=m_Q^{}`$, $`e_Q=e_Q^{}`$) since in that case $`𝒞^{V,S}(p^2)=0`$. The functions $`𝒜(p^2)`$, $`^{V,S}(p^2)`$ then simplify and coincide with the ones found previously gmf

$`𝒜(p^2)`$ $`=`$ $`e_Q\left({\displaystyle \frac{2m_Q}{M_B}}+{\displaystyle \frac{𝐩^2}{M_Bm_Q}}\right),`$
$`^V(p^2)`$ $`=`$ $`{\displaystyle \frac{e_Q}{2m_Q^2}}\left(1+{\displaystyle \frac{2(M_B2m_Q)}{m_Q}}{\displaystyle \frac{11𝐩^2}{4m_Q^2}}\right),`$
$`^S(p^2)`$ $`=`$ $`{\displaystyle \frac{e_Q}{2m_Q^2}}\left(1{\displaystyle \frac{2(M_B2m_Q)}{m_Q}}+{\displaystyle \frac{5𝐩^2}{4m_Q^2}}\right).`$
In the case of the $`B_c`$ meson $`𝒞^{V,S}(p^2)0`$ and thus relativistic corrections lead to spin-flip transitions ($`S^{}=S\pm 1`$) but only for decays involving mixed states $`nP1`$, $`nP1^{}`$ (40) or $`nD2`$, $`nD2^{}`$ (43). For all other transitions the spin-flip correction vanishes due to momenta relations (see Eq. (99) below).

Using the Wigner-Eckart theorem and relations for matrix elements of the tensor operator between coupled functions, one can rewrite Eq. (77) in the form

$`D_i`$ $`=`$ $`(1)^{J^{}+J+L+SM^{}}\left(\begin{array}{ccc}J^{}& 1& J\\ M^{}& i& M\end{array}\right)\sqrt{(2J+1)(2J^{}+1)}`$ (97)
$`\times \mathbf{[}\left\{\begin{array}{ccc}L^{}& J^{}& S\\ J& L& 1\end{array}\right\}\delta _{SS^{}}(n^{}L^{}||𝒜(p^2)r||nL\eta (J^{},L^{},J,L,S)n^{}L^{}||(p^2)p||nL)`$
$`\chi (J^{},L^{},J,L)\delta _{SS^{}\pm 1}n^{}L^{}||𝒞(p^2)p||nL\mathbf{]},`$
where

$$\eta (J^{},L^{},J,L,S)=(1)^{L+J+1}\sqrt{6}[S(S+1)(2S+1)]^{1/2}\left(\left\{\begin{array}{ccc}L^{}& J^{}& S\\ J& L& 1\end{array}\right\}\right)^1\left\{\begin{array}{cccc}J& S& L& \\ J^{}& S& L^{}& \\ 1& 1& 1& \end{array}\right\},$$
(98)
$$\chi (J^{},L^{},J,L)=(1)^{L+L^{}}\sqrt{\frac{2}{2J^{}+1}}\left\{\begin{array}{ccc}J& L& 1\\ 1& 1& L^{}\end{array}\right\}.$$
(99)
Here

$$\left(\begin{array}{ccc}J^{}& 1& J\\ M^{}& q& M\end{array}\right),\left\{\begin{array}{ccc}L^{}& J^{}& S\\ J& L& 1\end{array}\right\},\mathrm{and}\left\{\begin{array}{cccc}J& S& L& \\ J^{}& S& L^{}& \\ 1& 1& 1& \end{array}\right\}$$

are $`3j`$-, $`6j`$\- and $`9j`$-symbols, $`n^{}L^{}\mathrm{}nL`$ are reduced matrix elements.

The total E1 decay rate of the $`nJLS`$ state is obtained by summing the decay rates (71) over all possible values of $`M^{}`$ for a fixed value of $`M`$. The resulting expression is

$`\mathrm{\Gamma }^{V,S}(nJLSn^{}J^{}L^{}S^{}+\gamma )={\displaystyle \frac{4}{3}}\alpha \omega ^3`$ $`\{`$ $`\delta _{SS^{}}C^{1/2}(J^{},L^{},J,L,S)[(L^{}L)R_{n^{}L^{}}\left|𝒜(p^2)r\right|R_{nL}`$ (102)
$`{\displaystyle \frac{\eta (J^{},L^{},J,L,S)}{\sqrt{\mathrm{max}(L^{},L)}}}n^{}L^{}||^{V,S}(p^2)p||nL]`$
$`\delta _{SS^{}\pm 1}\chi (J^{},L^{},J,L)n^{}L^{}||𝒞^{V,S}(p^2)p||nL\}^2,`$
where

$$C(J^{},L^{},J,L,S)=\mathrm{max}(L^{},L)(2J^{}+1)\left\{\begin{array}{ccc}L^{}& J^{}& S\\ J& L& 1\end{array}\right\}^2.$$

The reduced matrix elements can be expressed through the usual matrix elements over radial wave functions $`R_{nL}(r)`$:
(a) for the transitions between $`P`$ and $`S`$ states

$`n^{}0pn1`$ $`=`$ $`R_{n^{}S}^{}|R_{nP},`$ (103)
$`n^{}0𝐩^2pn1`$ $`=`$ $`R_{n^{}S}^{\prime \prime \prime }|R_{nP}+2R_{n^{}S}^{\prime \prime }\left|{\displaystyle \frac{1}{r}}\right|R_{nP}2R_{n^{}S}^{}\left|{\displaystyle \frac{1}{r^2}}\right|R_{nP};`$ (104)
(b) for the transitions between $`D`$ and $`P`$ states

$`n^{}1pn2`$ $`=`$ $`\sqrt{2}\left(R_{n^{}P}\left|{\displaystyle \frac{1}{r}}\right|R_{nD}R_{n^{}P}^{}|R_{nD}\right),`$ (105)
$`n^{}1𝐩^2pn2`$ $`=`$ $`\sqrt{2}(R_{n^{}P}^{\prime \prime \prime }|R_{nD}+R_{n^{}P}^{\prime \prime }\left|{\displaystyle \frac{1}{r}}\right|R_{nD}`$ (107)
$`6R_{n^{}P}^{}\left|{\displaystyle \frac{1}{r^2}}\right|R_{nD}+6R_{n^{}P}\left|{\displaystyle \frac{1}{r^3}}\right|R_{nD}),`$
where the prime means differentiation of $`R`$ with respect to $`r`$.

### VII.2 Results and discussion

The results of numerical calculations of charmonium E1 radiative decay rates using Eqs. (98)–(105) are presented in Table 8. For calculations of photon energies $`\omega `$ we used the experimentally measured masses of charmonium $`S`$ and $`P`$ states. 666For decays involving $`\eta _c^{}(2S)`$, as in the case of M1 transitions, we use both experimental values of its mass, giving a prediction for the recent Belle value belle  in parenthesis. For masses of $`D`$ states we used our model predictions from Table 1. We give predictions for decay rates calculated in the nonrelativistic limit $`\mathrm{\Gamma }^{\mathrm{NR}}`$, for relativistic decay rates with pure vector $`\mathrm{\Gamma }^V`$ and scalar $`\mathrm{\Gamma }^S`$ potentials as well as for the mixture (II) of vector and scalar potentials $`\mathrm{\Gamma }`$ with $`\epsilon =1`$. As in the case of M1 decay rates calculations, we use the relativistic wave functions in our numerical analysis.

The results presented in Table 8 show that relativistic effects play an important role in E1 decays of charmonium. The most sensitive to the relativistic corrections are decays $`\mathrm{\Psi }^{}(2S)\chi _{cJ}(1P)+\gamma `$. Their account leads to the considerable reduction of the decay rates. The rates for the vector potential are reduced more significantly than for the scalar one. As a result, there arises an approximately twofold reduction of decay rates for the mixture of vector and scalar potentials with the value of mixing parameter $`\epsilon =1`$, bringing theoretical predictions in good agreement with experimental data. The large influence of relativistic corrections originates from the fact that the zero of the $`2S`$ wave function is close to the maximum of the $`1P`$ wave function. This results in a reduction of the leading order decay matrix element $`1P|r|2S|`$. Therefore relatively small relativistic corrections produce such a large effect. This observation is confirmed by the calculations of the $`\chi _{cJ}(1P)J/\mathrm{\Psi }+\gamma `$ decay rates. Here both initial and final states do not have zeros and relativistic contributions have usual values and lead to an approximately $`25\%`$ reduction of the decay rate. All theoretical predictions are in nice agreement with data. In Table 8 we also give predictions for E1 decay rates of charmonium $`1D`$ states. At present only the $`1{}_{}{}^{3}D_{1}^{}`$ state is experimentally observed. This state is considerably broader, since it lies above the $`D\overline{D}`$ threshold. The observed state $`\mathrm{\Psi }(3770)`$ can also have a significant $`2{}_{}{}^{3}S_{1}^{}`$ state admixture r . If we consider it to be a pure $`D`$ state, then using its measured total decay rate, we get the following predictions for the E1 radiative decay branching fractions:

$$(1{}_{}{}^{3}D_{1}^{}1{}_{}{}^{3}P_{0}^{}\gamma )1.5\%;(1{}_{}{}^{3}D_{1}^{}1{}_{}{}^{3}P_{1}^{}\gamma )0.6\%;(1{}_{}{}^{3}D_{1}^{}1{}_{}{}^{3}P_{2}^{}\gamma )0.03\%.$$

On the other hand, the tensor $`1{}_{}{}^{3}D_{2}^{}`$ and $`1{}_{}{}^{1}D_{2}^{}`$ charmonium states are under the threshold of open charm production, since their masses are slightly below the $`D\overline{D}^{}`$ threshold, and the decay of these states to $`D\overline{D}`$ is forbidden by parity and angular momentum conservation. Thus E1 radiative transitions are the main decay channels of these states.

The calculated decay rates of E1 radiative transitions in bottomonium are presented in Tables 9 and 10. The influence of relativistic effects in bottomonium is considerably less than in charmonium. The contribution of relativistic corrections does not exceed 10% almost for all decays. The only exceptions are decays $`\mathrm{{\rm Y}}^{\prime \prime }(3S)\chi _{bJ}(1P)+\gamma `$, where the leading contribution is substantially reduced due to the significantly different number of zeros in initial $`3S`$ and final $`1P`$ wave functions. For all $`SP+\gamma `$ transitions we find good agreement of our model predictions with experimental data.

The comparison of the theoretical predictions for the radiative decays of $`P`$ states of bottomonium $`\chi _b(nP)`$ with the experimental data is complicated by the fact that the total decay rates of these states are not measured yet. Experiment gives only branching fractions $`\mathrm{\Gamma }[\chi _b(nP)\mathrm{{\rm Y}}(n^{}S)+\gamma ]/\mathrm{\Gamma }^{\mathrm{total}}`$ ($`nn^{}`$) for these decays. Thus for this comparison it is necessary to get theoretical predictions for the total decay rates of $`\chi _b(nP)`$. The main decay channels of the bottomonium $`P`$ states are inclusive strong decays to gluon and quark states and radiative decays. The strong decays were extensively studied in the literature had ; kmrr  including leading-order QCD corrections. The two-gluon annihilation rates of $`{}_{}{}^{3}P_{0}^{}`$ states with the account of relativistic corrections are given in our model fgq88  by

$`\mathrm{\Gamma }({}_{}{}^{3}P_{0}^{}gg)`$ $`=`$ $`{\displaystyle \frac{8\alpha _s^2}{3M^2}}|{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}{\displaystyle \frac{m_Q^2}{Mp}}[{\displaystyle \frac{p}{ϵ_Q(p)}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}}`$ (109)
$`+(1{\displaystyle \frac{M}{2ϵ_Q(p)}})(2{\displaystyle \frac{ϵ_Q(p)}{p}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}})]\varphi _P(p)|^2.`$
To include leading-order QCD corrections kmrr  this expression should be multiplied by $`(1+10.0\alpha _s/\pi )`$ for $`\chi _{b0}`$ and $`(1+10.2\alpha _s/\pi )`$ for $`\chi _{b0}^{}`$. The corresponding expression for two-gluon annihilation rates of $`{}_{}{}^{3}P_{2}^{}`$ states reads as fgq88

$`\mathrm{\Gamma }({}_{}{}^{3}P_{2}^{}gg)`$ $`=`$ $`{\displaystyle \frac{8\alpha _s^2}{5M^2}}\{|{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}{\displaystyle \frac{m_Qϵ_Q(p)}{Mp}}[(2+{\displaystyle \frac{p}{ϵ_Q(p)}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}}{\displaystyle \frac{ϵ_Q(p)}{p}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}})`$ (114)
$`\times (1+{\displaystyle \frac{m_Q}{2ϵ_Q(p)[ϵ_Q(p)+m_Q]}}){\displaystyle \frac{2p^2}{3ϵ_Q(p)[ϵ_Q(p)+m_Q]}}]\varphi _P(p)|^2`$
$`+{\displaystyle \frac{1}{6}}|{\displaystyle }{\displaystyle \frac{d^3p}{(2\pi )^3}}{\displaystyle \frac{m_Q^2}{Mp}}[{\displaystyle \frac{p}{ϵ_Q(p)}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}}+3(2{\displaystyle \frac{ϵ_Q(p)}{p}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}})`$
$`{\displaystyle \frac{ϵ_Q(p)}{m_Q}}\left(2{\displaystyle \frac{M}{m_Q}}\right)\left(2{\displaystyle \frac{ϵ_Q(p)}{p}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}}\right)+{\displaystyle \frac{p^2}{m_Q[ϵ_Q(p)+m_Q]}}\left(2{\displaystyle \frac{M}{m_Q}}\right)`$
$`\times \mathbf{(}{\displaystyle \frac{1}{2}}({\displaystyle \frac{ϵ_Q(p)}{p}}\mathrm{ln}{\displaystyle \frac{ϵ_Q(p)+p}{ϵ_Q(p)p}}2)(1{\displaystyle \frac{3ϵ_Q^2(p)}{p^2}})+1\mathbf{)}]\varphi _P(p)|^2\},`$
with leading-order QCD corrections kmrr  accounted by the factors $`(10.1\alpha _s/\pi )`$ for $`\chi _{b2}`$ and $`(1+1.0\alpha _s/\pi )`$ for $`\chi _{b2}^{}`$. Here $`\varphi _P(p)`$ is the $`P`$ state radial wave function in momentum space. In the nonrelativistic limit $`p/m_Q0`$ and $`M2m_Q`$, Eqs. (109) and (114) reduce to the known expressions kmrr :

$`\mathrm{\Gamma }({}_{}{}^{3}P_{0}^{}gg)`$ $`=`$ $`{\displaystyle \frac{6\alpha _s^2\left|R_{nP}^{}(0)\right|^2}{m_Q^4}},`$ (115)
$`\mathrm{\Gamma }({}_{}{}^{3}P_{2}^{}gg)`$ $`=`$ $`{\displaystyle \frac{8\alpha _s^2\left|R_{nP}^{}(0)\right|^2}{5m_Q^4}}.`$ (116)
The calculation of radiative and relativistic corrections to the annihilation decay rates of $`{}_{}{}^{3}P_{1}^{}`$ and $`{}_{}{}^{1}P_{1}^{}`$ states is a very complicated problem which has not been solved yet. Thus to estimate their hadronic decay rates we use the tree-level nonrelativistic expressions kmrr

$`\mathrm{\Gamma }({}_{}{}^{3}P_{1}^{}q\overline{q}+g)`$ $`=`$ $`{\displaystyle \frac{8\alpha _s^3n_f}{9\pi m_Q^4}}\left|R_{nP}^{}(0)\right|^2\mathrm{ln}(m_Qr),`$ (117)
$`\mathrm{\Gamma }({}_{}{}^{1}P_{1}^{}ggg)`$ $`=`$ $`{\displaystyle \frac{20\alpha _s^3}{9\pi m_Q^4}}\left|R_{nP}^{}(0)\right|^2\mathrm{ln}(m_Qr),`$ (118)
$`\mathrm{\Gamma }({}_{}{}^{1}P_{1}^{}gg+\gamma )`$ $`=`$ $`{\displaystyle \frac{36}{5}}e_q^2{\displaystyle \frac{\alpha }{\alpha _s}}\mathrm{\Gamma }({}_{}{}^{1}P_{1}^{}ggg).`$ (119)
From the relativistic consideration of the decays of $`{}_{}{}^{3}P_{0,2}^{}`$ states, which shows that relativistic effects give $`10\%`$ contributions to the $`b\overline{b}`$ decay rates, we can expect that these formulae give a reasonable estimate of the corresponding decay rates. For numerical calculations of hadronic decay rates of $`\chi _b`$ states we use $`\alpha _s=0.18`$ obtained from the experimental ratio of $`\mathrm{\Gamma }(\mathrm{{\rm Y}}gg\gamma )/\mathrm{\Gamma }(\mathrm{{\rm Y}}ggg)`$ gr . The calculated partial decay rates and branching fractions for $`1P`$ and $`2P`$ states of the bottomonium are compared with available experimental data in Table 11. There we give both PDG pdg  averages and very recent CLEO cleo  data. We see that our predictions for the branching fractions for radiative transitions of $`\chi _b(nP)`$ to $`\mathrm{{\rm Y}}(n^{}S)`$ states are in good agreement with experiment. The only discrepancy (of $`2.5\sigma `$) is the CLEO value for $`(\chi _{b1}(2P)\mathrm{{\rm Y}}(2S)+\gamma )`$ transition which is approximately two times larger than our model prediction and PDG value. The CLEO Collaboration cleo  measured also two photon decays of $`\mathrm{{\rm Y}}(3S)`$ via $`\chi _b(1P_J)`$ states. They report the branching fraction for $`\mathrm{{\rm Y}}(3S)\chi _b(1P_J)+\gamma \mathrm{{\rm Y}}(1S)+\gamma \gamma `$ transitions summed over all the $`J`$ states:

$$(\mathrm{{\rm Y}}(3S)\mathrm{{\rm Y}}(1S)+\gamma \gamma )=(2.14\pm 0.22\pm 0.21)\times 10^3.$$
(120)
Using our model results in Tables 911 for corresponding decay rates, we get

$`(\mathrm{{\rm Y}}(3S)\chi _{b0}(1P)+\gamma \mathrm{{\rm Y}}(1S)+\gamma \gamma )`$ $`=`$ $`1.4\times 10^5,`$
$`(\mathrm{{\rm Y}}(3S)\chi _{b1}(1P)+\gamma \mathrm{{\rm Y}}(1S)+\gamma \gamma )`$ $`=`$ $`9.97\times 10^4,`$
$`(\mathrm{{\rm Y}}(3S)\chi _{12}(2P)+\gamma \mathrm{{\rm Y}}(1S)+\gamma \gamma )`$ $`=`$ $`9.96\times 10^4,`$
and the sum over all the $`J`$ states is equal to

$$(\mathrm{{\rm Y}}(3S)\mathrm{{\rm Y}}(1S)+\gamma \gamma )=2.04\times 10^3$$

in accord with CLEO data (120).

The CLEO Collaboration cleo  presented recently the first evidence for the production of the triplet $`\mathrm{{\rm Y}}(1D)`$ states in the four photon transitions $`3S2P+\gamma 1D+\gamma \gamma 1P+\gamma \gamma \gamma 1S+\gamma \gamma \gamma \gamma e^+e^{}+\gamma \gamma \gamma \gamma `$. The measured product branching fraction for these five decays is equal to $`(3.3\pm 0.6\pm 0.5)\times 10^5`$. In Table 12 we give the theoretical predictions for the branching fractions of such four photon decays in our model and in the recent quark model analysis of Godfrey and Rosner grD . In general, both theoretical predictions are consistent. In the analysis of Ref. grD , the dominant decay is $`3{}_{}{}^{3}S_{1}^{}2{}_{}{}^{3}P_{1}^{}+\gamma 1{}_{}{}^{3}D_{2}^{}+\gamma \gamma 1{}_{}{}^{3}P_{1}^{}+\gamma \gamma \gamma 1{}_{}{}^{3}S_{1}^{}+\gamma \gamma \gamma \gamma e^+e^{}+\gamma \gamma \gamma \gamma `$, while in our model the above and $`3{}_{}{}^{3}S_{1}^{}2{}_{}{}^{3}P_{2}^{}+\gamma 1{}_{}{}^{3}D_{3}^{}+\gamma \gamma 1{}_{}{}^{3}P_{2}^{}+\gamma \gamma \gamma 1{}_{}{}^{3}S_{1}^{}+\gamma \gamma \gamma \gamma e^+e^{}+\gamma \gamma \gamma \gamma `$ transitions have almost the same rate and dominate. In the last line we give the sum of all these decay channels. Both theoretical predictions agree with CLEO measurement.

An important problem of quarkonium physics is the search for bottomonium spin singlet states $`\eta _b(n{}_{}{}^{1}S_{0}^{})`$ and $`h_b(n{}_{}{}^{1}P_{1}^{})`$. In the previous section we discussed the possibilities of finding $`\eta _b`$ in radiative M1 decays. From Table 11 we see that these two states can be discovered simultaneously. Indeed the radiative decay $`h_b(1{}_{}{}^{1}P_{1}^{})\eta _b(1{}_{}{}^{1}S_{0}^{})+\gamma `$ with the photon energy of 480 MeV is the main decay channel of $`h_b`$ (the branching fraction of this decay exceeds 50%). Thus production of a few $`h_b`$ states, e.g., through $`\mathrm{{\rm Y}}^{\prime \prime }(3S)h_b(1{}_{}{}^{1}P_{1}^{})\pi \pi `$ or $`\mathrm{{\rm Y}}^{\prime \prime }(3S)h_b(1{}_{}{}^{1}P_{1}^{})\pi `$ decays, which branching fractions are predicted to be about 0.1–1% kyv , will give a good possibility of finding $`\eta _b`$.

In Tables 13 and 14 we compare our results for the E1 radiative decay rates of the $`B_c`$ meson with other quark model predictions eq ; gklt ; ful . Comparison of the calculations using relativistic $`\mathrm{\Gamma }`$ and nonrelativistic $`\mathrm{\Gamma }^{\mathrm{NR}}`$ formulae for decay rates shows that relativistic corrections do not exceed 20% in $`B_c`$ meson E1 decays. Most of the theoretical predictions for E1 transitions between $`P`$ and $`S`$ states of $`B_c`$ mesons given in Table 13 are compatible with each other. The largest differences occur for decays involving $`P1`$ and $`P1^{}`$ states which are the mixtures of spin singlet and spin triplet states (40) due to different mixing angles used by the authors. Note that for such transitions there are additional relativistic spin-flip contributions to decay rates (102) proportional to $`\chi (J^{},L^{},J,L)`$ (99) which are specific only for $`B_c`$ mesons. In general, our predictions are closer to the ones of Ref. gklt .

In Table 14 we present the E1 radiative decay rates of $`B_c`$ mesons where either the initial or final state is a $`D`$ wave state. Here we find rather large variations in theoretical predictions. The main reason of these distinctions is the difference in values of $`D`$ state masses, which for some states reaches 70 MeV (see Table 3). Since $`2P`$ and $`1D`$ states of the $`B_c`$ are rather close, such a difference significantly influences the energy of the emitted photon and thus the decay rates. For decays involving the mixed spin singlet and spin triplet states $`P1`$, $`P1^{}`$ (40) and $`D2`$, $`D2^{}`$ (43) the additional relativistic spin-flip contributions (102) are important, especially for transitions where both initial and final states are mixed states.

## VIII Conclusions

In this paper the mass spectra and radiative M1 and E1 decay rates of charmonium, bottomonium and $`B_c`$ mesons were calculated in the relativistic quark model based on the quasipotential approach in quantum field theory. Special attention was devoted to the role of the relativistic effects in these processes. Since both quarks in the considered mesons are heavy, the $`v/c`$ expansion was applied. In the mass spectra calculations retardation as well as one-loop radiative corrections were taken into account. We also included the one-loop correction due to the finite charm quark mass to the bottomonium mass spectrum. It was found that this correction is rather small efgmc  and its inclusion allows one to obtain an even better fit of the bottomonium excited states with the slightly shifted value of QCD parameter $`\mathrm{\Lambda }`$. The calculated mass spectra of charmonium and bottomonium agree with the experimental data within a few MeV. Comparison of our results for the $`B_c`$ meson mass spectrum with previous calculations showed that different predictions for ground state and low excitations agree within 30 MeV.

The pseudoscalar and vector decay constants of $`B_c`$ meson were calculated using the relativistic wave functions obtained during the mass spectrum calculations. It was found that relativistic effects reduce these constants by approximately 20% and produce the splitting between them of about 70 MeV.

It was shown that relativistic effects play a significant role in radiative decays of mesons. Their form strongly depends on the Lorentz structure of the quark-antiquark interaction. The most sensitive are radiative M1 decays, where even for allowed transitions they significantly influence predictions for the rates. An important example is the $`J/\mathrm{\Psi }\eta _c\gamma `$ transition, which is overestimated by a factor of more than two if the nonrelativistic approximation is used. It is argued that the inclusion of relativistic corrections for a pure scalar or vector confining potential is not enough to bring theoretical predictions in accord with experiment. Only for the specific mixture of these potentials (II) with the mixing coefficient $`\epsilon =1`$, the agreement can be obtained. For other decay rates this mixing of scalar and vector potentials also gives the best results. The hindered M1 transition rates are dominated by relativistic contributions and are significantly enhanced by them. The comparison of the allowed and hindered M1 rates in bottomonium shows that the latter provide better opportunity for the search of the missing pseudoscalar $`\eta _b`$ state of the bottomonium.

The analysis of radiative E1 transitions showed that the form of relativistic corrections is less dependent on the Lorentz structure of the quark potential than in the case of M1 transitions. However, for some decays, e.g., $`\mathrm{\Psi }^{}\chi _{cJ}\gamma `$, the consideration of the mixed (II) vector and scalar potentials (with the same value of $`\epsilon =1`$) is important for bringing decay rates in accord with experimental data. In general, all our predictions for radiative decay rates and branching fractions of charmonium and bottomonium agree with measured values. In the case of the $`B_c`$ meson radiative E1 decays an important additional relativistic correction to decay rates which causes the flip of the quark spin was found. This contribution (99) to the radiative E1 decay rate (102) is specific only for transitions involving mixed states $`nP1`$, $`nP1^{}`$ (40) or $`nD2`$, $`nD2^{}`$ (43) of $`B_c`$ and is caused by the difference of the $`c`$ and $`b`$ quark masses. Finally, a comparison of various quark model predictions for the radiative M1 and E1 decay rates of $`B_c`$ has been performed. These radiative transitions along with pionic ones are the main decay channels of the low lying excitations in the $`B_c`$ meson.

###### Acknowledgements.

The authors express their gratitude to A. Martin, M. Müller-Preussker and V. Savrin for support and discussions. We are grateful to J. L. Rosner for stimulating correspondence and quoting some of our preliminary results on bottomonium radiative transitions in Ref. grnb . Two of us (R.N.F and V.O.G.) were supported in part by the Deutsche Forschungsgemeinschaft under contract Eb 139/2-1, Russian Foundation for Basic Research under Grant No. 00-02-17768 and Russian Ministry of Education under Grant No. E00-3.3-45.


text:
# Complex Probabilities on 𝑅^𝑁 as Real Probabilities on 𝐶^𝑁 and an Application to Path Integrals

Present address:\]Clinton Group, 55 Water Street, New York, NY 10041

## Abstract

We establish a necessary and sufficient condition for averages over complex valued weight functions on $`R^N`$ to be represented as statistical averages over real, non-negative probability weights on $`C^N`$. Using this result, we show that many path-integrals for time-ordered expectation values of bosonic degrees of freedom in real-valued time can be expressed as statistical averages over ensembles of paths with complex-valued coordinates, and then speculate on possible consequences of this result for the relation between quantum and classical mechanics.

For a wide range of quantum mechanical systems, if real-valued time is analytically continued to imaginary values, the weight assigned to each path in the path integral formula for time-ordered vacuum expectation values becomes real and non-negative, permitting interpretation as a probability density. Paths distributed according to this probability density can then be generated by the Metropolis method, one version of which makes use of the Langevin equation.

For the Langevin equation applied to paths for imaginary time path-integrals, however, the coefficients entering the equation are themselves analytic functions of the action occurring in the path weight. As a consequence of this analyticity, the Langevin equation can be analytically continued, at least formally, to complex-valued weights Klauder , among which are those of the original path-integral with real-valued time. For complex path weights, the coefficients in the analytically continued Langevin equation and the coordinates in the path ensembles the equation generates are complex.

Numerical tests and some analytic results concerning the Langevin equation with complex coefficients yield a complicated picture subsequent . For some systems the method gives accurate answers, for others values of the coordinates blow up as the equation’s running time is increased, and in still other applications averages over complex Langevin trajectories converge but to wrong answers. A question which follows from these results is whether the failure of the Langevin equation for many complex-valued path weights is simply a limitation of the Langevin equation or, alternatively, arises from the nonexistence of equivalent real, non-negative weights on complex-valued trajectories. More generally, with no reference to the Langevin equation, under what circumstances can averages over complex weights on real coordinates be realized as statistical averages over real, non-negative weights on complex coordinates? The main result addressing this problem so far is a proof Salcedo  that an average over real coordinates of a complex polynomial multiplied by a complex Gaussian can be realized as an average over a real, non-negative weight on complex coordinates.

In the present article, by an explicit construction unrelated to the Langevin equation, we establish a necessary and sufficient condition for averages on $`R^N`$ over complex weight functions, with no explicit restriction on functional form, to be represented as averages over real, non-negative probability densities on $`C^N`$. For the particular case of $`N=1`$, we show this condition is fulfilled for all complex weights normalized to a total integral of 1.

As an application of our result for general $`N`$, we show that real-time path integrals for many systems with bosonic degrees of freedom can be represented as averages over ordinary statistical ensembles of paths through complex space. Among the systems which can be represented in this way are scalar field theories, in any space-time dimension, with self-interaction $`\varphi ^p`$, $`p`$ even, field theories, in any dimension, with bosonic degress of freedom which range only over a bounded region, and therefore lattice gauge field theories over any compact gauge group. We then briefly speculate on possible consequences of this result for the relation between quantum and classical mechanics.

The systems we consider all have time reduced to a finite, discrete lattice of points. We do not examine the existence of limits as the time lattice spacing is taken to zero. We also do not discuss algorithms for actually generating averages over the probability densities on $`C^N`$ which we prove exist.

For $`c(x_1,\mathrm{},x_N)`$ a complex-valued weight function on $`(x_1,\mathrm{},x_N)R^N`$, $`t(z_1,\mathrm{},z_N)`$ a real, non-negative weight on $`(z_1,\mathrm{},z_N)C^N`$, and $`p(r_1,\mathrm{},r_N)`$ a real, non-negative weight on $`(r_1,\mathrm{},r_N)[0,\mathrm{})^N`$, define the expectation values

$$\begin{array}{c}<f>_c=Z_c^1\underset{j}{}\left(_{\mathrm{}}^{\mathrm{}}dx_j\right)\times \hfill \\ \hfill f(x_1,\mathrm{},x_N)c(x_1,\mathrm{},x_N)\end{array}$$
(1a)
$$\begin{array}{c}<f>_t=Z_t^1\underset{j}{}\left(_0^{2\pi }d\theta _j_0^{\mathrm{}}r_jdr_j\right)\times \hfill \\ \hfill f(z_1,\mathrm{},z_N)t(z_1,\mathrm{},z_N),\end{array}$$
(1b)
$$\begin{array}{c}<f>_p=Z_p^1\underset{j}{}\left(_0^{\mathrm{}}dr_j\right)\times \hfill \\ \hfill f(r_1,\mathrm{},r_N)p(r_1,\mathrm{},r_N)\end{array}$$
(1c)
with $`z_j`$ given by $`r_jexp(i\theta _j)`$ in Eq. 1b. We require $`c(x)`$, $`t(z)`$ and $`p(r)`$ to have non-zero total integrals, then choose $`Z_c`$, $`Z_t`$ and $`Z_p`$ to give $`<1>_c`$, $`<1>_t`$, $`<1>_p`$ all the value $`1`$.

We will show that for any $`c(x_1,\mathrm{},x_N)`$, a necessary and sufficient condition for the existence of a $`t(z_1,\mathrm{},z_N)`$ with the same moments

$`<x_1^{m_1}\mathrm{}x_N^{m_N}>_c`$ $`=`$ $`<z_1^{m_1}\mathrm{}z_N^{m_N}>_t`$ (2)
for all integer $`m_j0`$, is that there exist a $`p(r_1,\mathrm{},r_N)`$ satisfying the bound

$`|<x_1^{m_1}\mathrm{}x_N^{m_N}>_c|`$ $``$ $`<r_1^{m_1}\mathrm{}r_N^{m_N}>_p`$ (3)
for all integer $`m_j0`$.

For notational simplicity, we consider first a detailed proof for $`N=1`$, then briefly summarize the generalization to arbitrary $`N`$.

To show that Eq. 2 implies Eq. 3 for some $`p(r)`$, we take absolute values of both sides of Eq. 2 and obtain

$`|Z_c^1{\displaystyle _{\mathrm{}}^{\mathrm{}}}x^mc(x)𝑑x|`$ $``$ $`Z_t^1{\displaystyle _0^{\mathrm{}}}r^mp(r)𝑑r,`$ (4a)
$`p(r)`$ $`=`$ $`r{\displaystyle _0^{2\pi }}𝑑\theta t(z).`$ (4b)
But $`Z_p`$ is identical to $`Z_t`$. Thus Eqs. 4 implies Eq. 3.

We now assume Eq. 3 and construct a $`t(z)`$ which satisfies Eq. 2. For $`z`$ given by $`rexp(i\theta )`$ as before and a real parameter $`\lambda >1`$, Eq. 3 implies we can define the functions

$`q_\lambda (\theta )`$ $`=`$ $`{\displaystyle \underset{m0}{}}{\displaystyle \frac{exp(im\theta )<x^m>_c}{\lambda ^m<r^m>_p}}1.`$ (5a)
$`s_\lambda (\theta )`$ $`=`$ $`\{1+2Re[q_\lambda (\theta )]\}.`$ (5b)
Eqs. 5 imply for $`m0`$

$`{\displaystyle _0^{2\pi }}{\displaystyle \frac{d\theta }{2\pi }}exp(im\theta )s_\lambda (\theta )`$ $`=`$ $`{\displaystyle \frac{<x^m>_c}{\lambda ^m<r^m>_p}}.`$ (6)
With the definitions

$`<f>_{\lambda t}`$ $`=`$ $`Z_{\lambda t}^1{\displaystyle _0^{2\pi }}𝑑\theta {\displaystyle _0^{\mathrm{}}}𝑑rrf(z)t_\lambda (z),`$ (7a)
$`t_\lambda (z)`$ $`=`$ $`r^1p(\lambda ^1r)s_\lambda (\theta ),`$ (7b)
and $`Z_{\lambda t}`$ chosen to make $`<1>_{\lambda t}=1`$, Eqs. 6 and 7 then yield

$`<z^m>_{\lambda t}`$ $`=`$ $`<x^m>_c,`$ (8)
which is Eq. 2 for $`N=1`$.

By a choice of $`\lambda `$ we can insure that $`t_\lambda (z)`$ is non-negative. From Eq. 7b it follows that $`t_\lambda (z)`$ is non-negative if $`s_\lambda (\theta )`$ is non-negative. Eqs. 3 and 5 then yield

$`s_\lambda (\theta )`$ $``$ $`12|q_\lambda (\theta )|,`$ (9a)
$``$ $`12{\displaystyle \underset{m1}{}}{\displaystyle \frac{|<x^n>_c|}{\lambda ^m<r^n>_p}},`$ (9b)
$``$ $`{\displaystyle \frac{\lambda 3}{\lambda 1}}.`$ (9c)
Eq. 9c imples that for $`\lambda 3`$, $`t_\lambda (z)`$ is non-negative.

Consider now arbitrary $`N`$. The proof that Eq. 2 implies Eq. 3 for some $`p(r_1,\mathrm{},r_N)`$ for arbitrary $`N`$ is an immediate extension of the $`N=1`$ proof. The proof that Eq. 3 implies Eq. 2 requires very little more effort. For general $`N`$, $`q_\lambda (\theta )`$ of Eq. 5a becomes

$$\begin{array}{c}q_\lambda (\theta _1,\mathrm{},\theta _N)=\hfill \\ \hfill \underset{m_i0}{}\frac{exp(i_im_i\theta _i)<x_1^{m_1}\mathrm{}x_N^{m_N}>_c}{\lambda ^{_im_i}<r_1^{m_1}\mathrm{}r_N^{m_N}>_p}1.\end{array}$$
(10)
Eq. 9c, insuring the positivity of $`s_\lambda (z_1,\mathrm{},z_N)`$ and therefore of $`t_\lambda (z_1,\mathrm{},z_N)`$ becomes

$`s_\lambda (\theta _1,\mathrm{},\theta _N)`$ $``$ $`3{\displaystyle \frac{2\lambda ^N}{(\lambda 1)^N}},`$ (11)
the right side of which is positive for

$`\lambda `$ $``$ $`{\displaystyle \frac{1}{1(\frac{2}{3})^{\frac{1}{N}}}}.`$ (12)
Returning to the case of $`N=1`$, we now show that a $`p(r)`$ exists fulfilling Eq. 3 for any $`c(x)`$ with $`<1>_c`$ of $`1`$proof . To do this, we make use of a theorem giving conditions for the existence of a $`p(r)`$ such that for all non-negative integers $`m`$

$`<r^m>_p=s_m,`$ (13)
for a specified sequence of real numbers $`s_0,s_1,\mathrm{}`$. From the $`s_m`$, for any non-negative integer $`n`$, define the $`(n+1)\times (n+1)`$ matrices

$`H_{ij}^{2n}`$ $`=`$ $`s_{i+j},`$ (14a)
$`H_{ij}^{2n+1}`$ $`=`$ $`s_{i+j+1},`$ (14b)
with $`0i,jn`$. Then according to a result of Stieltjes Akhiezer  a $`p(r)`$ obeying Eq. 13 exists if for all non-negative $`n`$

$`det(H^n)>0.`$ (15)
For any $`c(x)`$, we will construct a sequence $`s_0,s_1,\mathrm{},`$ obeying both Eqs. 3 and 15. Choose $`s_0`$ to be $`1`$. Now suppose $`s_0,s_1,\mathrm{},s_k,`$ have been found obeying Eqs. 3 and 15 for all $`nk`$. By isolating the dependence of $`det(H^{k+1})`$ on $`s^{k+1}`$, it is not hard to show that a positive real bound $`b`$ exists, such that Eq. 15 is fulfulled for $`n=k+1`$ for any $`s_{k+1}>b`$. We choose $`s_{k+1}`$ to be the greater of $`|<x^{k+1}>_c|`$ and $`b`$. By induction, the resulting sequence $`s_0,s_1,\mathrm{}`$ obeys Eqs. 3 and 15 for all $`m`$. It follows that a $`p(r)`$ exists fulfilling Eq. 13 and thus also Eq. 3.

The result we have just proved can not be extended automatically to $`N>1`$ since, for $`N>1`$, the analog of Eq. 15 is not a sufficient condition for the existence of a $`p(r_1,\mathrm{},r_N)`$ fulfilling the analog of Eq. 13 Berg .

For each $`c(x_1,\mathrm{},x_N)`$ obeying Eq. 3 for some $`p(r_1,\mathrm{},r_N)`$, the non-negative, real $`t_\lambda (z_1,\mathrm{},z_N)`$ we construct satisfying Eq. 2 is only one representative of a large family of such functions fulfilling Eq. 2. At the very least, given any $`p(r_1,\mathrm{},r_N)`$ obeying Eq. 3, an infinite family of distinct $`p^{}(r_1,\mathrm{},r_N)`$ which also obey Eq. 3 can be constructed by shifting some part of the weight carried by $`p(r_1,\mathrm{},r_N)`$ out to larger values of $`r_1,\mathrm{},r_N`$. For example, for positive $`a_1,\mathrm{},a_N`$, define $`p^{}(r_1,\mathrm{},r_N)`$ to be zero except on the region $`r_1a_1,\mathrm{},r_Na_N`$, where it is given by $`p(r_1a_1,\mathrm{},r_Na_N)`$. Each moment of $`p^{}(r_1,\mathrm{},r_N)`$ is greater than the corresponding moment of $`p(r_1,\mathrm{},r_N)`$, and thus satisfies Eq. 3. By means of Eqs. 10 \- 12, $`p^{}(r_1,\mathrm{},r_N)`$ can therefore be turned into a new $`t_\lambda ^{}(z_1,\mathrm{},z_N)`$ with moments satisfying Eq. 2.

An expression for $`q_\lambda (\theta _1,\mathrm{},\theta _N)`$ somewhat simpler than Eq. 10 can be found if $`p(r_1,\mathrm{},r_N)`$ in Eq. 3 is a product $`\widehat{p}(r_1)\mathrm{}\widehat{p}(r_N)`$ with moments $`<r^m>_{\widehat{p}}`$ which, for any $`R>0`$ and integer $`m0`$, obey

$`<r^m>_{\widehat{p}}`$ $``$ $`b_RR^m`$ (16)
for some $`b_R>0`$ independent of $`m`$. Then with the definition

$`f_\lambda (y)`$ $`=`$ $`{\displaystyle \underset{m0}{}}{\displaystyle \frac{y^m}{\lambda ^m<r^m>_{\widehat{p}}}},`$ (17)
Eq. 10 takes the form

$$\begin{array}{c}q_\lambda (\theta _1,\mathrm{},\theta _N)=\hfill \\ \hfill <f_\lambda [x_1exp(i\theta _1)]\mathrm{}f_\lambda [x_Nexp(i\theta _N)]>_c1.\end{array}$$
(18)
For many quantum mechanical systems, a finite volume, lattice approximation to the path integral in real-valued time leads to an expectation value fulfilling Eq. 3 and therefore respresentable as an ordinary statistical average over an ensemble of random paths through complex space. Among the systems which obey Eq. 3 for some $`p(x_1,\mathrm{},x_N)`$ are scalar field theories, in any space-time dimension, with self-interaction $`\varphi ^p`$, $`p`$ even, field theories, in any dimension, with bosonic degrees of freedom which range only over a bounded region, and therefore lattice gauge field theories over any compact gauge group. In addition, of course, for any system for which the real-valued time path integral is correctly handled by the complex Langevin equation, the probability distribution $`t(z_1,\mathrm{},z_N)`$ given by the limit as running time becomes large of the ensemble of points in $`C^N`$ struck by the Langevin trajectory fulfills Eq. 2. Since Eq. 2 implies Eq. 3, such a system then has a $`p(r_1,\mathrm{},r_N)`$ fulfilling Eq. 3.

As an example of a proof of Eq. 3, consider the path integral for a single anharmonic oscillator with time reduced to a periodic lattice of points labeled by a positive integer $`tT`$. With $`x_t`$ the position at time $`t`$, the time-ordered expectation of a polynomial $`f(x_1,\mathrm{},x_T)`$ becomes $`<f>_c`$ of Eq. 1a with

$$\begin{array}{c}c(x_1,\mathrm{},x_T)=\hfill \\ \hfill exp\left\{i\underset{t}{}\left[\mu \frac{(x_tx_{t+1})^2}{2\delta }\kappa x_t^4\delta \right]\right\},\end{array}$$
(19)
where $`\mu `$ is the particle’s mass, $`\kappa `$ the anharmonic spring constant, and $`\delta `$ the lattice spacing. Rotating each integral over real $`x_t`$ to an integral over $`exp(i\pi /8)y_t`$ with $`y_t`$ real, taking the absolute value of both sides of the rotated version of Eq. 1a, and using the inequality

$`(y_ty_{t+1})^2`$ $``$ $`2(y_t^2+y_{t+1}^2),`$ (20)
we obtain

$`|<x_1^{m_1}\mathrm{}x_N^{m_N}>_c|`$ $``$ $`{\displaystyle \frac{Z_{\widehat{p}}^T}{|Z_c|}}{\displaystyle \underset{t}{}}<r_t^{m_t}>_{\widehat{p}},`$ (21a)
$`\widehat{p}(r)`$ $`=`$ $`exp\left({\displaystyle \frac{\sqrt{2}r^2}{\mu \delta }}\kappa r^4\delta \right),`$ (21b)
with $`<\mathrm{}>_{\widehat{p}}`$ given by Eq. 1c for $`N=1`$. Eq. 21a then gives Eq. 3 with

$`p(x_1,\mathrm{},x_T)`$ $`=`$ $`{\displaystyle \underset{t}{}}\widehat{p}({\displaystyle \frac{x_t}{\lambda }}),`$ (22a)
$`\lambda `$ $`=`$ $`Max[1,{\displaystyle \frac{Z_{\widehat{p}}}{|Z_c|^{\frac{1}{T}}}}].`$ (22b)
The proof of Eq. 3 for an anharmonic oscillator extends easily to a $`\varphi ^p`$ field theory, $`p`$ even, for any dimension of space-time and to bosonic field theories, in any dimension, with with degress of freedom bounded by some positive real $`B`$. In a proof of Eq. 3 for bounded degrees of freedom, $`\widehat{p}(r)`$ in Eq. 21b becomes $`\delta (rB)`$.

Finally, we speculate briefly on the possibility that the degrees of freedom of the real world might actually be complex-valued, distributed according to a probability which, as we have shown possible in many cases, yields averages agreeing with those of conventional quantum mechanics. A version of this idea suggested by the complex Langevin equation is considered in Ref. Migdal . A problem encountered immediately by this proposal is that the macroscopic world, governed nearly by classical mechanics, exhibits only real-valued positions. This problem would be solved if, within the statistical ensemble of paths through complex-space, the paths for macroscopic variables happened nearly to obey classical mechanics and to have only small imaginary parts. A formulation of quantum mechanics incorporating this feature would, potentially, be capable of resolving many of the puzzles in the interpretation of quantum mechanics bell ; nuovocim . A discussion of the proposed resolution of these problems by decoherence criteria nbi ; others  and of the difficulties which these encounter appears in Ref. goldstein .

The speculation we offer here, in effect, is that quantum mechanics is a version of ordinary statistical mechanics but for paths through a complex-valued coordinate space. Although this prospect may sound implausible, the harmonic oscillator provides a simple model of how such behavior might conceivably come about.

Assume periodic boundary conditions for time period $`T=2M+1`$ and Fourier transform the oscillator coordinate $`x_t`$ at integer time $`t`$ according to the convention

$$x_t=\underset{0kM}{}\left[a_kcos(\frac{2\pi kt}{T})+b_ksin(\frac{2\pi kt}{T})\right],$$
(23)
with real coefficients $`a_k`$ and $`b_k`$, $`b_0`$ identically 0. The time-ordered expectation of a polynomial $`f(a_0,\mathrm{},a_M,b_1,\mathrm{},b_M)`$ is $`<f>_c`$ of Eq. 1a for

$$\begin{array}{c}c(a_0,\mathrm{},a_M,b_1,\mathrm{},b_M)=\hfill \\ \hfill exp\left[i\underset{0kM}{}(a_k^2+b_k^2)s_k\right],\end{array}$$
(24)
where $`s_k`$ is

$$s_k=\frac{2\mu }{\delta }[1cos(k)]\kappa \delta ,$$
(25)
and $`\delta `$ is the time lattice spacing.

In Eq. 1a expressed as integrals over real $`a_k`$ and $`b_k`$, for $`s_k`$ positive the integrals can be rotated, respectively, to $`\sqrt{i}c_k`$ and $`\sqrt{i}d_k`$ with real $`c_k`$ and $`d_k`$. For $`s_k`$ negative, the $`a_k`$ and $`b_k`$ integrals can be rotated, respectively, to $`\sqrt{i}c_k`$ and $`\sqrt{i}d_k`$. If $`\mu `$, however, is given a microscopic positive imaginary part $`ϵ`$ and a correctly chosen real part, there will be a single $`k`$ for which $`s_k`$ has no real part and a microscopic positive imaginary part. For this $`k`$ the integrals over $`a_k`$ and $`b_k`$ can be left pure real. Inverting the Fourier transform of Eq. 23, a statistical ensemble of complex-valued trajectories $`x_t`$ results with real, non-negative probability weight. A typical trajectory in this ensemble consists of the sum a real part, with amplitude of order $`ϵ^{1/2}`$, oscillating at the frequency predicted by classical mechanics, and a complex part independent of $`ϵ`$.

The preceding construction is, effectively, an application to the harmonic oscillator path integral of the method of steepest descent. A similar construction can be carried out for a free scalar field theory in any dimension of space time. For more complicated field theories, the method of steepest descent can also be applied to obtain a $`t(x_1,\mathrm{},x_N)`$ which includes real-valued trajectories obeying the classical equations of motion, but the $`t(x_1,\mathrm{},x_N)`$ produced in this way will only approximately fulfill Eq. 2. Whether for interacting systems $`t(x_1,\mathrm{},x_N)`$ can be found which fulfill Eq. 2 and include real-valued trajectories obeying classical equations of motion is an open question.