text:
# 1 Introduction

## 1 Introduction

The description of various aspects of $`B_c`$-meson physics111For a review see  has received recently a lot of attention both from theoretical and experimental sides. $`B_c`$-meson is the only meson within SM, which contains two heavy quarks of different flavor. For this reason one might expect, that its spectroscopy, production mechanisms and decays differ significantly from those of charmonium, bottomonium as well as hadrons with one heavy quark. Indeed, contrary to the case of $`J/\mathrm{\Psi }`$ and $`\mathrm{{\rm Y}}`$-mesons, $`B_c`$-meson is a long living system, decaying trough electroweak interactions. To produce $`\overline{b}c`$-system in flavor conserving collisions, an additional heavy quark pair should be created, what makes the consideration of production mechanism complex already at the leading order in the coupling constant. The bound state heavy quarks dynamics, being very similar to that of other quarkonia, at the same time has some differences. With respect to the hierarchy of dynamical scales, $`B_c`$-meson stands among the families of charmonium $`\overline{c}c`$ and bottomonium $`\overline{b}b`$ and thus could be used to study both quantitatively and conceptually existing effective low energy frameworks for the description of bound state heavy quark dynamics, like NRQCD , pNRQCD  and vNRQCD . So, the study of $`B_c`$-meson gives us an opportunity to test SM in a completely new dynamical regime, explore it at a different angle, what may give us an additional insight into our understanding of SM in general.

Today, we already have many theoretical predictions, concerning various $`B_c`$-meson properties: spectroscopy , production  and decay channels . In the present work we will be interested in radiative corrections to the leptonic decay constant of $`B_c`$-meson, entering expression for its purely leptonic width. This particular decay mode of $`B_c`$-meson is not really interesting for a direct measurement — the necessarily spin-flip makes its branching fraction very small. However, the same quantity appears in the formulation of two-point QCD sum rules . And the latter may be used to extract an additional information about parameters of SM, namely the heavy quark masses. A first step towards systematic formulation of NNLO $`B_c`$-meson sum rules, treating $`B_c`$-meson as a Coulomb system, was made in , where NLO analysis was presented. In order to describe $`B_c`$-meson QCD sum rules at NNLO, one needs both 2-loop matching Wilson coefficient of QCD current with $`B_c`$ quantum numbers on corresponding NRQCD current and two point NRQCD correlator, computed at NNLO. Below we will present the results of matching procedure. These results are also interesting from the point of view of further development of computation technics used in dimensionally regulated NRQCD . Here we are dealing with threshold expansion with two mass scales $`m_b`$ and $`m_c`$, what is an extension of previously obtained results for matching of vector current with two equal heavy quark masses . The presence of two mass scales of course makes problem more complicated, however as we will see the whole calculation could be performed within the general framework of threshold expansion in dimensionally regulated NRQCD . Here we would like to stress, that we treat $`B_c`$-meson as a Coulomb system. So, in our matching procedure $`c`$-quark is considered to be heavy with respect to factorization scale $`m_c>\mu _{\mathrm{fac}}`$. As a consequence, it is impossible to make connection of our result with the results on the matching of heavy-light currents . On the other side the limit $`m_b=m_c`$ is well defined in our approach and in the case of vector current matching with equal quark masses we reproduce already known results of .

This paper is organized as follows. In section 2 we introduce necessary notation and present one-loop results obtained previously . Section 3 contains two-loop results for $`B_c`$-meson leptonic decay constant. Here we also describe shortly various steps of the performed calculation. And finally, in section 4 we present our conclusion.

## 2 $`B_c`$ \- meson leptonic constant: one loop result

The $`B_c`$-meson leptonic decay proceeds through a virtual $`W`$-boson. The $`W`$-boson couples to the $`B_c`$ through the axial-vector part of the charged weak current. All QCD effects, both perturbative and nonperturbative, enter into the decay rate through the decay constant $`f_{B_c}`$, defined by the matrix element

$$0|\overline{b}\gamma ^\mu \gamma _5c|B_c(P)=if_{B_c}P^\mu ,$$
(1)
where $`|B_c(P)`$ is the state of $`B_c`$-meson with four-momentum $`P`$. It has the standard covariant normalization $`B_c(P^{})|B_c(P)=(2\pi )^32E\delta ^3(P^{}P)`$ and its phase has been chosen so that $`f_{B_c}`$ is real and positive.

To separate short-distance ($`1/m_Q`$) and long-distance ($`1/m_Q`$) QCD contributions to $`f_{B_c}`$ we use the formalism of NRQCD . As a result the decay amplitude is factored into short-distance coefficients multiplied by NRQCD matrix elements. The Wilson short-distance coefficients could be then calculated further order by order in perturbation theory by matching a perturbative calculation in the full QCD with the corresponding perturbative calculation in NRQCD. It is the subject of this paper to calculate Wilson coefficient in front of NRQCD current up to 2-loop order in perturbative QCD.

Let us now discuss the matching procedure in some detail. We start with NRQCD lagrangian describing $`B_c`$-meson bound state dynamics

$$_{\mathrm{NRQCD}}=_{\mathrm{light}}+\psi _c^{}\left(iD_0+𝐃^2/(2m_c)\right)\psi _c+\chi _b^{}\left(iD_0𝐃^2/(2m_b)\right)\chi _b+\mathrm{},$$
(2)
where $`_{\mathrm{light}}`$ is the usual relativistic lagrangian for gluons and light quarks. The 2-component field $`\psi _c`$ annihilates charm quarks, while $`\chi _b`$ creates bottom antiquarks. The typical relative velocity $`v`$ of heavy quarks inside the meson provides a small parameter that can be used as a nonperturbative expansion parameter.

To express the decay constant $`f_{B_c}`$ in terms of NRQCD matrix elements, we express the axial-vector current $`\overline{b}\gamma ^\mu \gamma _5c`$ in terms of NRQCD fields. Only the $`\mu =0`$ component contributes to the matrix element (1) in the rest frame of the $`B_c`$-meson. This component of the current has the following operator expansion in terms of NRQCD fields:

$$\overline{b}\gamma ^0\gamma _5c=C_0(m_b,m_c)\chi _b^{}\psi _c+C_2(m_b,m_c)(𝐃\chi _b)^{}𝐃\psi _c+\mathrm{},$$
(3)
where $`C_0`$ and $`C_2`$ are Wilson short-distance coefficients that depend on the quark masses $`m_b`$ and $`m_c`$. By dimensional analysis, the coefficient $`C_2`$ is proportional to $`1/m_Q^2`$. The contribution to the matrix element $`0|\overline{b}\gamma ^0\gamma _5c|B_c`$ from the operator $`(𝐃\chi _b)^{}𝐃\psi _c`$ is suppressed by $`v^2`$ relative to the operator $`\chi _b^{}\psi _c`$, where $`v`$ is the relative velocity of heavy quarks inside $`B_c`$-meson. The dots in (3) represent other operators whose contributions are suppressed by higher powers of $`v^2`$.

The short-distance coefficient $`C_0`$ and $`C_2`$ can be determined by matching perturbative calculations of the matrix elements in the full QCD and NRQCD. A convenient choice for matching is the matrix element between the vacuum and the state $`|c\overline{b}`$ consisting of a $`c`$ and a $`\overline{b}`$ on their perturbative mass shells with nonrelativistic four-momenta $`p`$ and $`p^{}`$ in the center of momentum frame: $`𝐩+𝐩^{}=0`$. The matching condition is

$$0|\overline{b}\gamma ^0\gamma _5c|c\overline{b}|_{\text{QCD}}=C_00|\chi _b^{}\psi _c|c\overline{b}|_{\text{NRQCD}}+C_20|(𝐃\chi _b)^{}𝐃\psi _c|c\overline{b}|_{\text{NRQCD}}+\mathrm{},$$
(4)
where QCD and NRQCD represent perturbative QCD and perturbative NRQCD, respectively. At leading order in $`\alpha _s`$, the matrix element on the left hand side of (4) is $`\overline{v}_b(𝐩)\gamma ^0\gamma _5u_c(𝐩)`$. The Dirac spinors are

$`u_c(\text{p})`$ $`=`$ $`\sqrt{{\displaystyle \frac{E_c+m_c}{2E_c}}}\left(\begin{array}{c}\xi \\ \frac{\text{p}𝝈}{E_c+m_c}\xi \end{array}\right),`$ (7)
$`v(\text{p})`$ $`=`$ $`\sqrt{{\displaystyle \frac{E_b+m_b}{2E_b}}}\left(\begin{array}{c}\frac{(\text{p})𝝈}{E_b+m_b}\eta \\ \eta \end{array}\right),`$ (10)
where $`\xi `$ and $`\eta `$ are 2-component spinors and $`E_Q=m_Q^2+𝐩^2`$. Making a nonrelativistic expansion of the spinors to second order in $`𝐩/m_Q`$, we find

$$\overline{v}_b(𝐩)\gamma ^0\gamma _5u_c(𝐩)\eta _b^{}\xi _c\left(1\frac{1}{8}\left(\frac{m_b+m_c}{m_bm_c}\right)^2𝐩^2+\mathrm{}\right).$$
(11)
At leading order in $`\alpha _s`$, the matrix elements on the right hand side of (4) are $`\eta _b^{}\xi _c`$ and $`𝐩^2\eta _b^{}\xi _c`$. The short distance coefficients are therefore $`C_0=1`$ and

$$C_2=\frac{1}{8m_{\mathrm{red}}^2},$$
(12)
where $`m_{\mathrm{red}}=m_bm_c/(m_b+m_c)`$ is the reduced mass of the system.

In order to determine the short distance coefficients up to order $`\alpha _s^2`$, we must calculate the matrix elements on both sides of (4) up to $`\alpha _s^2`$. Up to this order it is sufficient to calculate radiative corrections only for the coefficient $`C_0`$, since the contribution proportional to $`C_2`$ is suppressed by additional power of $`v^2`$. The coefficient $`C_0`$ can be isolated by taking the limit $`𝐩0`$, in which case the matrix element of $`(𝐃\chi _b)^{}𝐃\psi _c`$ vanishes. Such calculation at the one-loop order was performed in , where the following result was obtained for $`C_0`$:

$$C_0=\mathrm{\hspace{0.33em}1}+\frac{\alpha _s(m_{\mathrm{red}})}{\pi }\left[\frac{m_bm_c}{m_b+m_c}\mathrm{log}\frac{m_b}{m_c}2\right].$$
(13)
$`m_{\mathrm{red}}`$ is the scale of the running coupling constant. To the accuracy of this calculation, any scale of order $`m_b`$ or $`m_c`$ would be equally correct.

Thus, the final result for the decay constant of the $`B_c`$ is

$$if_{B_c}M_{B_c}=C_00|\chi _b^{}\psi _c|B_c+C_20|(𝐃\chi _b)^{}𝐃\psi _c|B_c.$$
(14)
with short-distance coefficient $`C_0`$ up to next-to-leading order in $`\alpha _s`$ as given in (13), while $`C_2`$ is given to leading order in (12). The uncertainties consist of the perturbative errors in the short distance coefficients and an error of relative order $`v^4`$ from the neglected matrix elements that are higher order in $`v^2`$. The matrix element $`0|\chi _b^{}\psi _c|B_c`$ can be estimated using wavefunctions at the origin from nonrelativistic potential models:

$$|0|\chi _b^{}\psi _c|B_c|^2\mathrm{\hspace{0.33em}2}M_{B_c}\frac{3}{2\pi }|R(0)|^2.$$
(15)
The factor of $`2M_{B_c}`$ takes into account the relativistic normalization of the state $`|B_c`$.

## 3 Two-loop result

In this section we consider the two-loop contribution to $`B_c`$-meson leptonic constant. An analogous calculation, but for the case of vector current and equal quark masses was done previously in . Here we are following the approach of , extending it, where it is necessary, to the case of two different mass scales involved in our problem. In the case of $`B_c`$-meson we consider the matching of QCD axial current on the NRQCD one.

To compute Wilson short distance coefficient $`C_0`$ up to $`\alpha _s^2`$ order, we replace $`B_c`$-meson state by a free quark-antiquark pair of on-shell heavy quarks with small relative velocity. In terms of the on-shell matrix elements for heavy quarks with small relative velocity: $`p=qm_1/(m_1+m_2)`$ and $`p^{}=qm_2/(m_1+m_2)`$ ($`q`$ is the momentum of the gauge $`W`$-boson) 222This and previous definitions of external heavy quark momenta lead to the same value of the matching coefficient, but latter is more convenient in two-loop calculations.. The matching equation has the following form

$$Z_{2,QCD}\mathrm{\Gamma }_{QCD}=C_0Z_{2,\mathrm{NRQCD}}Z_{J,\mathrm{NRQCD}}^1\mathrm{\Gamma }_{\mathrm{NRQCD}}+O(v^2),$$
(16)
where $`\mathrm{\Gamma }_{\mathrm{QCD}}`$ and $`\mathrm{\Gamma }_{\mathrm{NRQCD}}`$ are amputated bare axial-vector vertices in QCD and NRQCD, correspondingly. Here, we would like to note, that since terms of order $`v^2`$ are not needed to determine $`C_0`$ (see Eq. (4)), one can set relative momentum to zero and compute Green functions directly at the threshold. Note, that here we deal with partially conserved axial-vector current and its QCD renormalization constant is given by $`Z_{J,\mathrm{QCD}}=1`$, as in the case of vector current.

In order to obtain the coefficient $`C_0`$ at two loops one has to perform the following steps:

* sum over all contributions to $`\mathrm{\Gamma }`$’s, including the one-loop term, multiplied by the two-loop QCD on-shell wave function renormalization constant , $`Z_{J,\mathrm{QCD}}^1`$ and one-loop NRQCD-current renormalization constant (at one-loop it is equal to unity),
* perform the one-loop renormalization of the coupling and masses.

After these manipulations the remaining poles in $`ϵ`$ belong to the anomalous dimension of the corresponding NRQCD current and finite part correspond to short-distance constant $`C_0`$.

In step 1) calculating two-loop contribution to vector-axial vertex, it is convenient to write it in the following form

$$A_\mu =\overline{b}(p^{})\left[\gamma _\mu \gamma _5F_1(q^2)\frac{1}{m_1+m_2}q_\mu \gamma _5F_4(q^2)\right]c(p).$$
(17)
The formfactors $`F_1(q^2)`$ and $`F_2(q^2)`$ can be singled out by imposing appropriate projections on them and thus our calculation will reduce to the calculation of diagrams with nonvanishing numerators.

To simplify matching calculation at 2-loop order we use threshold expansion in relative velocity $`v`$ of . The latter is obtained by writing down contributions corresponding to hard $`(km_Q)`$, soft $`(km_Qv)`$, potential $`(k_0m_Qv^2,k_im_Qv)`$ and ultrasoft $`(km_Qv^2)`$ regions333Here $`k`$ is the integration momentum.. The contributions from the last three regions can be easily identified with NRQCD diagrams, that appear in the calculation of $`\mathrm{\Gamma }_{\mathrm{NRQCD}}`$. Thus, they simply drop out from the matching relation (16) and it will suffice to compute only contribution to threshold expansion of $`\mathrm{\Gamma }_{\mathrm{QCD}}`$ of the region, where all momenta are hard. Diagrams, whose hard contribution need to be computed, are depicted in Fig. 1. Both form factors, introduced above have Coulomb singularities at the threshold and diverge as $`1/v^2`$. However, these singularities appear only in the soft, potential and ultrasoft contributions, while the hard contribution is well defined in dimensional regularization directly at the threshold. Once the relative velocity is set to zero, the corresponding loop integrals are considerably simplified and depend only on two scales $`M`$ and $`m`$, where $`M`$ is the mass of the heaviest quark and $`m`$ is the mass of the lighter quark. All loop integrals with irreducible numerators can then be further reduced to the scalar integrals with increased space-time dimension . For the scalar integrals we derived recurrence relations using integration by parts method of , which allowed us to reduce all scalar integrals to sunsets. The sunset integrals obtained were further reduced to master integrals with the use of recurrence relations of . We would like to note, that our reduction procedure is valid for arbitrary values of heavy quark masses, while to calculate some master integrals an expansion in parameter $`r=m/M`$ was performed. All master integrals are of the sunset type discussed in the paper . In particular there we discuss in full details the procedure of expansion in parameter $`r`$. Here, we describe only shortly the calculation procedure and more detailed exposition will be given elsewhere. The results, obtained for each of diagrams in Fig. 1 could be found in Appendix A.

After performing the steps of renormalization procedure described above, we have the following $`\overline{\text{MS}}`$ expression for the anomalous dimension of NRQCR current, which first arises here at the two-loop order ($`r=m/M`$)

$`\gamma _{J,\mathrm{NRQCD}}={\displaystyle \frac{d\mathrm{ln}Z_{J,\mathrm{NRQCD}}}{d\mathrm{ln}\mu }}=\left\{C_F^2\left(2{\displaystyle \frac{(1r)^2}{(1+r)^2}}\right)+C_FC_A\right\}{\displaystyle \frac{\pi ^2}{2}}\left({\displaystyle \frac{\alpha _s}{\pi }}\right)^2+𝒪(\alpha _s^3).`$ (18)
This scale dependence is compensated by the scale dependence of the two-loop matching coefficient $`c_2(M/\mu )`$, which is defined as

$`C_0(\alpha _s,{\displaystyle \frac{M}{\mu }})=1+c_1(M/\mu ){\displaystyle \frac{\alpha _s(M)}{\pi }}+c_2(M/\mu )\left({\displaystyle \frac{\alpha _s(M)}{\pi }}\right)^2+\mathrm{},`$ (19)
In order to obtain $`c_2`$ in analytical form, we perform the expansion of master integrals, considering the parameter $`r=m/M`$ as small (in the practice $`r=m_c/m_b0.3`$). Below we give the result of expansion up to $`O(r^2)`$ inclusive. Separating the different colour group factors, in $`c_2(M/\mu )`$ in the $`\overline{\text{MS}}`$ we have:

$`c_2(M/\mu )`$ $`=`$ $`C_F^2c_{2,A}+C_FC_Ac_{2,NA}+C_FT_Fn_lc_{2,L}+C_FT_Fc_{2,H},`$ (20)
$`c_{2,A}`$ $`=`$ $`\left(2{\displaystyle \frac{(1r)^2}{(1+r)^2}}\right){\displaystyle \frac{\pi ^2}{4}}\mathrm{ln}\left({\displaystyle \frac{M^2}{\mu ^2}}\right)+{\displaystyle \frac{9}{32}}L^2+{\displaystyle \frac{\pi ^2}{3}}L+{\displaystyle \frac{41}{32}}L+\pi ^2\mathrm{ln}(2){\displaystyle \frac{91\pi ^2}{48}}+{\displaystyle \frac{29}{16}}`$
$`+r\left[{\displaystyle \frac{17\pi ^2}{12}}L{\displaystyle \frac{43}{16}}L6\zeta _3+{\displaystyle \frac{23\pi ^2}{24}}{\displaystyle \frac{3}{8}}\right]+r^2\left[{\displaystyle \frac{21}{8}}L^23\pi ^2L+{\displaystyle \frac{165}{16}}L+12\zeta _3{\displaystyle \frac{25\pi ^2}{12}}{\displaystyle \frac{147}{16}}\right],`$
$`c_{2,NA}`$ $`=`$ $`{\displaystyle \frac{\pi ^2}{4}}\mathrm{ln}\left({\displaystyle \frac{M^2}{\mu ^2}}\right)+{\displaystyle \frac{11}{16}}L^2+{\displaystyle \frac{13\pi ^2}{24}}L+{\displaystyle \frac{127}{96}}L{\displaystyle \frac{9}{4}}\zeta _3{\displaystyle \frac{1}{2}}\pi ^2\mathrm{ln}(2)+{\displaystyle \frac{41\pi ^2}{96}}{\displaystyle \frac{23}{48}}`$
$`+r\left[{\displaystyle \frac{25}{16}}L^2{\displaystyle \frac{\pi ^2}{24}}L+{\displaystyle \frac{65}{48}}L{\displaystyle \frac{2\pi ^2}{3}}{\displaystyle \frac{49}{32}}\right]+r^2\left[{\displaystyle \frac{15}{8}}L^2{\displaystyle \frac{23}{16}}L+{\displaystyle \frac{19\pi ^2}{48}}+{\displaystyle \frac{31}{36}}\right],`$
$`c_{2,L}`$ $`=`$ $`{\displaystyle \frac{1}{4}}L^2{\displaystyle \frac{11}{24}}L+{\displaystyle \frac{5\pi ^2}{144}}+{\displaystyle \frac{1}{12}}+r\left[{\displaystyle \frac{L^2}{2}}{\displaystyle \frac{L}{12}}\right]+r^2\left[{\displaystyle \frac{L}{12}}{\displaystyle \frac{L^2}{2}}\right],`$ (23)
$`c_{2,H}`$ $`=`$ $`{\displaystyle \frac{L^2}{4}}{\displaystyle \frac{11L}{24}}{\displaystyle \frac{19\pi ^2}{96}}+{\displaystyle \frac{1147}{288}}+r\left[{\displaystyle \frac{L^2}{2}}{\displaystyle \frac{L}{12}}{\displaystyle \frac{23\pi ^2}{96}}{\displaystyle \frac{103}{144}}\right]`$ (24)
$`+r^2\left[{\displaystyle \frac{L^2}{2}}+{\displaystyle \frac{89L}{60}}+{\displaystyle \frac{23\pi ^2}{96}}+{\displaystyle \frac{11527}{3600}}\right].`$
Here $`\zeta _3=1.202\mathrm{}`$, $`M=m_b`$, $`r=m_c/m_b`$, $`L\mathrm{ln}(r)`$ and we have taken all fermions with masses less than $`m_b`$ and $`m_c`$ as massless. The heavy quark masses, appearing in our formulae are pole quark masses. All calculations were performed with the help of computer algebra system FORM . As a check for our calculations, we employed a general covariant gauge. To test our FORM programs we also have reproduced already known results for the sort distance corrections to the vector current with equal masses.

Now let us discuss the expression for 2-loop anomalous dimension of $`B_c`$-meson NRQCD current. It is an exact formula in a sense, that we resumed a whole series in $`r`$. In the limit of equal quark masses it differs from the result obtained in the case of vector current. This fact is explained by the contribution, coming from Breit–Fermi potential, namely spin-spin and spin-orbital heavy quark interaction, which in the case of $`B_c`$-meson has a different expression compared to that computed for $`J/\mathrm{\Psi }`$ and $`\mathrm{{\rm Y}}`$-mesons . This result for anomalous expression could be also obtained considering NNLO corrections to the Coulomb Green function, given by Breit–Fermi and non-abelian potentials.

Numerically, the scale dependence of the two-loop matching coefficient (20) is large. We expect it to be compensated by the scale dependence of $`B_c`$-meson wave function computed within the lattice NRQCD framework with the scale dependence of the results introduced by imposing ultraviolet cutoff/factorization scale. The results of nonrelativistic potential models should be considered as being obtained at a scale typical for $`B_c`$-meson bound state dynamics $`\mu 1`$ GeV. For numerical estimates we consider $`\mu =1`$ GeV as an adequate bound state scale and use its value for factorization scale $`\mu _{\mathrm{fac}}`$. In the case of $`B_c`$-meson and pole quark masses ($`m_b=4.8`$ GeV, $`m_c=1.65`$ GeV) we have

$`f_{B_c}=\left[11.48\left({\displaystyle \frac{\alpha _s(m_b)}{\pi }}\right)24.24\left({\displaystyle \frac{\alpha _s(m_b)}{\pi }}\right)^2\right]f_{B_c}^{NR}`$ (25)
Here we see that 2-loop corrections are large and constitute $`100\%`$ of 1-loop correction. Even less favorable situation was detected earlier for $`J/\mathrm{\Psi }`$ states . The variation of factorization scale does not help to solve this problem. So, at this point one may conclude that application of NRQCD to $`B_c`$-meson as well as for $`J/\mathrm{\Psi }`$ is not well justified. However, while in the case of $`J/\mathrm{\Psi }`$-meson 2-loop corrections far exceed 1-loop, one may speculate that situation in $`B_c`$-meson is marginal, in a sense that one may still obtain sensible results using NRQCD framework as well as treating $`B_c`$-meson as a heavy-light system. This way one may study the change in heavy quarks dynamics between bottomonium and charmonium families. Another possible reason for such large correction is our use of particular renormalization scheme $`\overline{\text{MS}}`$, which may turn out not to work well for heavy quark bound state systems. To prove or disprove this statement one has to calculate 2-loop corrections to any relation between physical observables, from which quarkonium wave function at the origin is eliminated .

## 4 Conclusion

In this paper we presented the results for the two-loop radiative corrections to the short distance coefficient in matching between QCD and NRQCD for the $`B_c`$-meson currents. It turns out that the numerical value of these corrections are large. There are several conclusions, which could be made here. First NRQCD description is inadequate not only for charmonium family but also for $`B_c`$-meson. To see whether it is true, we suppose to perform a systematical study of NRQCD $`B_c`$-meson sum rules at NNLO in one of our next publications. Second, is it the use of $`\overline{\text{MS}}`$ renormalization scheme, which makes things look nasty?. To answer this question one needs to compute some other decay mode of either $`J/\mathrm{\Psi }`$ or $`B_c`$ -mesons at 2-loop level and analyze perturbative corrections to the ratio of these decay modes.

We thank K. Chetyrkin, A. Davydychev, A. Grozin, M. Kalmykov, V. Kiselev, O. Tarasov and A. Pivovarov for fruitful discussion of the topics discussed in this paper. We especially grateful to V. Smirnov for encouragement, constant interest to this work and help with numerous checks of our calculation. This work was supported by DFG Forschergruppe ”Quantenfeldtheorie, Computeralgebra and Monte-Carlo-Simulation” (contract FOR 264/2-1) and by BMBF under grant No 05HT9VKB0.

## Appendix A 2 loop bare hard contribution

Here we collected results for hard contributions of diagrams, shown in Fig. 1. Below, to save space, we present only coefficient of $`\left(\alpha _s/\pi \right)^2(e^{\gamma _E}M^2/(4\pi \mu ^2))^{2\epsilon }`$ for each of the digrams, presented above. Though our calculations were performed in general covariant gauge, we present here for brevity only the results in the Feynman gauge.

* diagram 1 (Fig. 1a)

$`C_F^2\text{{}`$ $`{\displaystyle \frac{9}{32\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{3}{64}}{\displaystyle \frac{\pi ^2}{4}}{\displaystyle \frac{9L}{8}}+r[{\displaystyle \frac{9L}{8}}{\displaystyle \frac{\pi ^2}{2}}]+r^2[\pi ^2{\displaystyle \frac{9L}{8}}]\right)`$ (26)
$`+{\displaystyle \frac{9L^2}{4}}+{\displaystyle \frac{11}{12}}\pi ^2L+{\displaystyle \frac{L}{2}}3\zeta _3{\displaystyle \frac{299\pi ^2}{192}}{\displaystyle \frac{41}{128}}+r\left[3L^2+{\displaystyle \frac{3}{2}}\pi ^2L{\displaystyle \frac{9L}{8}}6\zeta _3+{\displaystyle \frac{\pi ^2}{3}}{\displaystyle \frac{15}{32}}\right]`$
$`+r^2\left[{\displaystyle \frac{7L^2}{8}}{\displaystyle \frac{19}{6}}\pi ^2L+{\displaystyle \frac{491L}{48}}+12\zeta _3{\displaystyle \frac{49\pi ^2}{24}}{\displaystyle \frac{401}{36}}\right]\text{}},`$
* diagram 2 (Fig. 1d, with selfenergy for quark with mass $`M`$)

$`C_F^2\text{{}`$ $`{\displaystyle \frac{3}{32\epsilon ^2}}(1{\displaystyle \frac{2}{r}})+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{1}{r}}[{\displaystyle \frac{3L}{8}}+{\displaystyle \frac{3}{8}}]{\displaystyle \frac{3L}{8}}{\displaystyle \frac{3}{64}}r[{\displaystyle \frac{3L}{8}}+{\displaystyle \frac{3}{4}}]+r^2[{\displaystyle \frac{9L}{8}}+{\displaystyle \frac{3}{4}}]\right)`$ (27)
$`{\displaystyle \frac{1}{r}}\left[{\displaystyle \frac{3}{8}}L^2+{\displaystyle \frac{3}{4}}L+{\displaystyle \frac{\pi ^2}{32}}+{\displaystyle \frac{13}{12}}\right]+{\displaystyle \frac{9}{8}}L^2+{\displaystyle \frac{13}{8}}L+{\displaystyle \frac{45\pi ^2}{64}}+{\displaystyle \frac{2041}{384}}`$
$`+r\left[{\displaystyle \frac{3}{8}}L^2{\displaystyle \frac{17}{4}}L+{\displaystyle \frac{\pi ^2}{4}}{\displaystyle \frac{63}{32}}\right]+r^2\left[{\displaystyle \frac{5}{8}}L^2+{\displaystyle \frac{31}{6}}L+{\displaystyle \frac{\pi ^2}{4}}+{\displaystyle \frac{1001}{288}}\right]\text{}},`$
* diagram 3 (Fig. 1d, with selfenergy for quark with mass $`m`$)

$`C_F^2\text{{}`$ $`{\displaystyle \frac{3}{32\epsilon ^2}}(12r)+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{3}{8}}L{\displaystyle \frac{51}{64}}+r[{\displaystyle \frac{3}{2}}L+{\displaystyle \frac{9}{8}}]r^2[{\displaystyle \frac{15}{8}}L+{\displaystyle \frac{3}{4}}]\right)`$ (28)
$`+{\displaystyle \frac{3}{4}}L^2+3L+{\displaystyle \frac{247\pi ^2}{192}}+{\displaystyle \frac{1057}{384}}r\left[3L^2+{\displaystyle \frac{29}{16}}L+{\displaystyle \frac{11\pi ^2}{96}}+{\displaystyle \frac{37}{48}}\right]`$
$`+r^2\left[{\displaystyle \frac{15}{4}}L^2{\displaystyle \frac{11}{8}}L{\displaystyle \frac{5\pi ^2}{8}}{\displaystyle \frac{21}{32}}\right]\text{}},`$
* diagram 4 (Fig. 1e, with two gluons coupled to quark with mass $`M`$)

$`C_FC_A\text{{}`$ $`{\displaystyle \frac{15}{64\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{19}{128}}{\displaystyle \frac{9}{16}}L+r[{\displaystyle \frac{9}{16}}L{\displaystyle \frac{\pi ^2}{16}}]+r^2[{\displaystyle \frac{\pi ^2}{16}}{\displaystyle \frac{9}{16}}L]\right)`$ (29)
$`+{\displaystyle \frac{9}{16}}L^2{\displaystyle \frac{3}{16}}L{\displaystyle \frac{7\pi ^2}{128}}+{\displaystyle \frac{329}{256}}+r\left[{\displaystyle \frac{15}{16}}L^2+{\displaystyle \frac{\pi ^2}{6}}L+{\displaystyle \frac{11}{16}}L{\displaystyle \frac{3}{4}}\zeta _3{\displaystyle \frac{25\pi ^2}{48}}{\displaystyle \frac{21}{64}}\right]`$
$`+r^2\left[{\displaystyle \frac{13}{16}}L^2{\displaystyle \frac{\pi ^2}{6}}L{\displaystyle \frac{7}{96}}L+{\displaystyle \frac{3}{4}}\zeta _3+{\displaystyle \frac{7\pi ^2}{16}}{\displaystyle \frac{115}{144}}\right]\text{}},`$
* diagram 5 (Fig. 1e, with two gluons coupled to quark with mass $`m`$)

$`C_FC_A\text{{}`$ $`{\displaystyle \frac{15}{64\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{15}{16}}L{\displaystyle \frac{\pi ^2}{16}}+{\displaystyle \frac{19}{128}}+r[{\displaystyle \frac{9}{16}}L+{\displaystyle \frac{\pi ^2}{16}}]r^2[{\displaystyle \frac{9}{16}}+{\displaystyle \frac{\pi ^2}{16}}]\right)`$ (30)
$`+{\displaystyle \frac{15}{8}}L^2+{\displaystyle \frac{\pi ^2}{4}}L{\displaystyle \frac{5}{16}}L{\displaystyle \frac{3}{4}}\zeta _3{\displaystyle \frac{221\pi ^2}{384}}+{\displaystyle \frac{473}{256}}+r\left[{\displaystyle \frac{3}{2}}L^2{\displaystyle \frac{\pi ^2}{4}}L+{\displaystyle \frac{3}{16}}L+{\displaystyle \frac{3}{4}}\zeta _3+{\displaystyle \frac{\pi ^2}{12}}{\displaystyle \frac{21}{64}}\right]`$
$`+r^2\left[2L^2+{\displaystyle \frac{\pi ^2}{4}}L{\displaystyle \frac{83}{96}}L{\displaystyle \frac{3}{4}}\zeta _3+{\displaystyle \frac{\pi ^2}{8}}+{\displaystyle \frac{353}{288}}\right]\text{}},`$
* diagram 6 (Fig. 1c, with two gluons coupled to quark with mass $`M`$)

$`C_F(C_A2C_F)\text{{}`$ $`{\displaystyle \frac{9}{64\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{3}{16}}L{\displaystyle \frac{27}{128}}{\displaystyle \frac{3}{16}}rL+{\displaystyle \frac{3}{16}}r^2L\right)+{\displaystyle \frac{3}{16}}L^2+{\displaystyle \frac{5}{8}}L{\displaystyle \frac{3}{16}}\zeta _3`$ (31)
$`+{\displaystyle \frac{\pi ^2}{8}}\mathrm{ln}(2)+{\displaystyle \frac{33\pi ^2}{128}}+{\displaystyle \frac{287}{256}}+r\left[{\displaystyle \frac{3}{16}}L^2L{\displaystyle \frac{3}{8}}\zeta _3+{\displaystyle \frac{\pi ^2}{4}}\mathrm{ln}(2)+{\displaystyle \frac{\pi ^2}{12}}{\displaystyle \frac{9}{64}}\right]`$
$`+r^2\left[{\displaystyle \frac{1}{16}}L^2{\displaystyle \frac{11}{48}}L+{\displaystyle \frac{3}{8}}\zeta _3{\displaystyle \frac{\pi ^2}{4}}\mathrm{ln}(2)+{\displaystyle \frac{5\pi ^2}{48}}+{\displaystyle \frac{71}{36}}\right]\text{}},`$
* diagram 7 (Fig. 1c, with two gluons coupled to quark with mass $`m`$)

$`C_F(C_A2C_F)\text{{}`$ $`{\displaystyle \frac{9}{64\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{9}{16}}L{\displaystyle \frac{27}{128}}{\displaystyle \frac{3}{16}}rL+{\displaystyle \frac{3}{16}}r^2L\right){\displaystyle \frac{9}{8}}L^2+{\displaystyle \frac{3}{4}}L{\displaystyle \frac{3}{16}}\zeta _3`$ (32)
$`+{\displaystyle \frac{5}{8}}\pi ^2\mathrm{ln}(2)+{\displaystyle \frac{77\pi ^2}{128}}+{\displaystyle \frac{215}{256}}+r\left[{\displaystyle \frac{3}{8}}L^2+{\displaystyle \frac{25}{32}}L{\displaystyle \frac{3}{8}}\zeta _3{\displaystyle \frac{3}{4}}\pi ^2\mathrm{ln}(2){\displaystyle \frac{9\pi ^2}{32}}{\displaystyle \frac{53}{32}}\right]`$
$`+r^2\left[{\displaystyle \frac{3}{16}}L^2{\displaystyle \frac{5}{4}}L+{\displaystyle \frac{3}{8}}\zeta _3+{\displaystyle \frac{3}{4}}\pi ^2\mathrm{ln}(2)+{\displaystyle \frac{11\pi ^2}{24}}+{\displaystyle \frac{137}{64}}\right]\text{}},`$
* diagram 8 (Fig. 1b)

$`C_F(C_A2C_F)\text{{}`$ $`{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{3}{16}}{\displaystyle \frac{\pi ^2}{16}}\right)+{\displaystyle \frac{7}{24}}\pi ^2L{\displaystyle \frac{5}{8}}L{\displaystyle \frac{15}{8}}\zeta _3{\displaystyle \frac{3}{4}}\pi ^2\mathrm{ln}(2){\displaystyle \frac{3\pi ^2}{16}}{\displaystyle \frac{35}{12}}`$ (33)
$`+r\left[{\displaystyle \frac{3}{16}}L^2+{\displaystyle \frac{\pi ^2}{24}}L{\displaystyle \frac{5}{32}}L+{\displaystyle \frac{3}{4}}\zeta _3+{\displaystyle \frac{\pi ^2}{2}}\mathrm{ln}(2){\displaystyle \frac{\pi ^2}{32}}+{\displaystyle \frac{59}{64}}\right]`$
$`+r^2\left[{\displaystyle \frac{11}{16}}L^2{\displaystyle \frac{\pi ^2}{12}}L+{\displaystyle \frac{11}{6}}L{\displaystyle \frac{3}{4}}\zeta _3{\displaystyle \frac{\pi ^2}{2}}\mathrm{ln}(2){\displaystyle \frac{35\pi ^2}{48}}{\displaystyle \frac{2119}{576}}\right]\text{}},`$
* diagram 9 (Fig. 1f)

$`C_FC_A\text{{}`$ $`{\displaystyle \frac{19}{128\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{33}{256}}{\displaystyle \frac{19}{32}}L+{\displaystyle \frac{19}{32}}rL{\displaystyle \frac{19}{32}}r^2L\right)+{\displaystyle \frac{19}{16}}L^2{\displaystyle \frac{7}{32}}L+{\displaystyle \frac{95\pi ^2}{768}}+{\displaystyle \frac{1315}{512}}`$ (34)
$`r\left[{\displaystyle \frac{19}{16}}L^2+{\displaystyle \frac{5}{64}}L\right]+r^2\left[{\displaystyle \frac{19}{16}}L^2+{\displaystyle \frac{5}{64}}L\right]\text{}},`$
* diagram 10 (Fig. 1g)

$`C_FC_A\text{{}`$ $`{\displaystyle \frac{1}{128\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{3}{256}}{\displaystyle \frac{1}{32}}L+{\displaystyle \frac{1}{32}}rL{\displaystyle \frac{1}{32}}r^2L\right)+{\displaystyle \frac{1}{16}}L^2{\displaystyle \frac{1}{32}}L+{\displaystyle \frac{5\pi ^2}{768}}+{\displaystyle \frac{73}{512}}`$ (35)
$`+r\left[{\displaystyle \frac{1}{64}}L{\displaystyle \frac{1}{16}}L^2\right]+r^2\left[{\displaystyle \frac{1}{16}}L^2{\displaystyle \frac{1}{64}}L\right]\text{}},`$
* diagram 11 (Fig. 1h, with light quark loop)

$`C_FT_Fn_l\text{{}`$ $`{\displaystyle \frac{1}{8\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{1}{2}}L{\displaystyle \frac{1}{16}}{\displaystyle \frac{1}{2}}r(1r)L\right)L^2{\displaystyle \frac{5\pi ^2}{48}}{\displaystyle \frac{67}{32}}+r(1r)\left[L^2+{\displaystyle \frac{1}{4}}L\right]\text{}},`$
* diagram 12 (Fig. 1h, with quark of mass $`M`$ in the loop)

$`C_FT_F\text{{}`$ $`{\displaystyle \frac{1}{4\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left({\displaystyle \frac{1}{2}}L+{\displaystyle \frac{5}{48}}{\displaystyle \frac{1}{2}}r(1r)L\right){\displaystyle \frac{1}{2}}L^2{\displaystyle \frac{1}{3}}L{\displaystyle \frac{\pi ^2}{8}}{\displaystyle \frac{3}{32}}`$ (37)
$`+r\left[{\displaystyle \frac{1}{2}}L^2+{\displaystyle \frac{1}{3}}L+{\displaystyle \frac{\pi ^2}{6}}{\displaystyle \frac{211}{144}}\right]+r^2\left[{\displaystyle \frac{1}{2}}L^2{\displaystyle \frac{13}{30}}L{\displaystyle \frac{\pi ^2}{6}}+{\displaystyle \frac{6007}{3600}}\right]\text{}},`$
* diagram 13 (Fig. 1h, with quark of mass $`m`$ in the loop)

$`C_FT_F\text{{}`$ $`{\displaystyle \frac{1}{4\epsilon ^2}}+{\displaystyle \frac{1}{\epsilon }}\left(L+{\displaystyle \frac{5}{48}}{\displaystyle \frac{1}{2}}r(1r)\right)2L^2{\displaystyle \frac{2}{3}}L+{\displaystyle \frac{13\pi ^2}{96}}{\displaystyle \frac{563}{288}}`$ (38)
$`+r\left[L^2+{\displaystyle \frac{1}{4}}L{\displaystyle \frac{7\pi ^2}{32}}+{\displaystyle \frac{3}{4}}\right]+r^2\left[L^2+{\displaystyle \frac{5}{4}}L+{\displaystyle \frac{13\pi ^2}{32}}\right]\text{}},`$

where $`L\mathrm{ln}r`$ and $`n_l`$ is the number of light fermions.


format: text

text:
# The Vela Pulsar’s proper motion and parallax derived from VLBI observations

## 1 Introduction

The Vela pulsar is an archetypal young pulsar. It was the first to be observed to glitch (Radhakrishnan & Manchester, 1969) and was always associated with the large Vela supernova remnant (SNR). Since Vela is a young pulsar and displays a significant amount of glitching behavior and timing noise (Cordes et al., 1988), no accurate proper motion measurements from timing observations have been possible; it is only through direct proper motion observations that we can trace its path on the sky. This is true of nearly all young pulsars, which are usually the most interesting. These have, for example, the possibility of discovering the associated birth SNR. The first accurate measurements of Vela’s proper motion were produced using the Parkes-Tidbinbilla Interferometer (PTI) (Bailes et al., 1989). This 275 km baseline at 1.6 GHz gave a resolution of 140 mas. This provided the first accurate measurement of the proper motion and thus, by extrapolation, the birth site. We have extended this work to allow the third dimension to be derived via the parallax. Measurement of relative motion on milli-arcsecond scales is a challenging task and several attempts have been made, with various success. Previous observations have been made with ground based optical observations, with connected radio interferometers (PTI and the VLA) and the Hubble Space Telescope (HST). We now present results from the Australian Long baseline array (LBA).

We have 18 observations of the position and can provide the most accurate measurements of parallax and proper motion. Recent measurements of the position of Vela against a number of background sources has been reported by De Luca et al. (2000); Caraveo et al. (2001). Our more accurate values, using completely independent methods, are compared with theirs, confirming the parallax values. Differences in the derivation of the proper motion are discussed. We find good agreement if we account for the probable galactic rotation of the reference stars.

The distance to the Vela pulsar was originally estimated from the SNR distance by comparison to the Cygnus loop and IC443. This lead to an estimate of 500 pc (Milne, 1968), however it has long been argued that this is an overestimate. Our results have been foreshadowed by the predictions in several papers that have made the case for the distance to the Vela pulsar to be drastically reduced. Analysis of X-ray observations of the pulsar made with ROSAT (Page et al., 1996) and Chandra (Pavlov et al., 2001b), of the SNR also with ROSAT (Bocchino et al., 1999) and with optical absorption (Jenkins & Wallerstein, 1995; Cha et al., 1999) have all suggested that the original distance estimate should perhaps be halved. All of these results, however, are to some extent model dependent and therefore doubt has existed over their accuracy.

Furthermore, our results allow the testing of the prediction of Spruit & Phinney (1998) on the alignment of the spin axes of pulsars with their proper motion vector. The spin axis can be derived from the high resolution images from Chandra. These allow the symmetry, and thus presumably the spin, axes to be directly discerned. Only two cases have been tested so far, the Crab (Caraveo & Mignani, 1999) and Vela. The extra accuracy we can provide refines the conclusion reached by Pavlov et al. (2000); Helfand et al. (2001) for the Vela alignment. A note of caution should be raised in that an alternate explanation for the X-ray nebula around Vela has been put forward (Radhakrishnan & Deshpande, 2001) which explains the structure in terms of the particle beam from the polar cap. This description requires that only the projected spin axis lies along the proper motion axis. If the model of Spruit & Phinney (1998) is correct any misalignment of the axes allows the us to estimate the timescale of the impulse in terms of the initial rotation period, a value of great importance to understanding the core collapse. The accuracy, therefore, of the alignment of the axes is of great theoretical importance.

## 2 Observations and Data reduction

The LBA array is an Australian national facility and is usually made up of six telescopes, three operated by the CSIRO (Australian Telescope Compact Array, Parkes and Mopra), two by the University of Tasmania (Hobart and Ceduna) and one in South Africa operated by the Hartebeesthoek Radio Astronomical observatory. In addition DSN telescopes at Tidbinbilla are often included. The observations reported here are between the Tidbinbilla 70m and the Hobart 26m antennae, a baseline of 832 km and at two 16-MHz wide frequencies; 2.29 and 8.425-GHz. In addition there was one observation with the Australian array, 32-MHz centered on 8.417-GHz, although only the Hobart to Tidbinbilla data was directly used for the results in this paper. Details of all the telescopes and their parameters are available via the Australian Telescope National Facility (ATNF) website111http://www.atnf.csiro.au/vlbi. All data were recorded using beam-switching on a $``$5 minute cycle, using the extra galactic phase-reference source Vela-G ($`\alpha _{J2000}=08^h33^m22^s.31563`$, $`\delta _{J2000}=44^{}41^{}38^{\prime \prime }.71463`$). Vela-G is part of the Radio Reference Frame and was an International Celestial Reference Frame (ICRF) candidate source (Ma et al., 1998).

The first three epochs observations of the Vela pulsar (Table 1) were made using MkIII/MkIIIA VLBI recording systems, using recording mode B, which provides 14 contiguous frequency channels each 2 MHz wide in right-hand circular polarization at 2.3 GHz. Sampling was one-bit (two-level) for all data. None of the Mark III/IIIA observations used the pulsar gating nor binning described later. Only the MkIII positions were available for this analysis. The MkIII/MkIIIA VLBI data were processed at the Washington Mark IIIA correlator, located at the U.S. Naval Observatory (USNO). Phase models were applied with the CALC 6.0 package, supplying the best possible values of station locations, clock models and earth orientation parameters. The raw data were then exported in ”FRNGX” format and further analyzed following the phase-referencing technique described by Lestrade et al. (1990), and used the software SPRINT, developed by J.-F. Lestrade.

The S2 recordings system and correlator for LBA were commissioned in the mid-1990s and observations were transferred to this system. In all observations, two separate frequency channels of 16 MHz bandwidth were recorded with two-bit (four-level) sampling. Most of the observations of the Vela pulsar using the S2 system were made with one 16 MHz channel at 2.3 GHz (13 cm) and one at 8.4 GHz (3 cm), taking advantage of the dual 13/3 cm (right-circular polarization) receivers at both Hobart and Tidbinbilla. In the final epoch, using the whole LBA at 8.4 GHz, two 16 MHz channels were used to form one contiguous 32 MHz band.

The S2 correlator is operated by ATNF at its Marsfield headquarters in Sydney. This wonderfully flexible device can pulsar bin (as opposed to pulsar gate) with thirty-two bins across the Vela pulsar’s period, thereby limiting the dispersion to individual channels rather than the entire band. There is, however, no significant dispersion at these observing frequencies for Vela. We used only those pulsar bins containing significant flux, usually one or two out of thirty two, gaining up to a factor of 5.6 in signal to noise.

Post correlation the data were fringe fitted at the nominal pulsar position using the binned flux, with ten minute solution interval in AIPS222AIPS, Astronomical Image Processing System, developed and maintained by the NRAO., and then exported. Further processing was done with DIFMAP (Shepherd, 1997). The data were flagged and averaged to generate a statistical weight. The offsets of Vela-G from its nominal reference position were found by fitting a point source to the visibilities. The position of Vela was shifted (phase rotated) to correct for these offsets. We then selfcalibrated the Vela phases in AIPS at 5 minutes intervals, and copied those corrections to the Vela-G visibilities. The final offsets of Vela-G from the reference position, all very small, were used to calculate the final Vela position.

Lobe ambiguity was not an issue even though this was a single baseline experiment, as observing with two widely separated frequencies broke any degeneracy. The final observation, which was with two 8.4-GHz bandpasses, had multiple baselines with which to identify the correct lobe. Only the Tidbinbilla to Hobart baseline was used for measurements of the positions.

## 3 Results

We find a pulsar position of $`\alpha _{J2000}=08^h35^m20^s.61153\pm 0.00002`$, $`\delta _{J2000}=45^{}10^{}34^{\prime \prime }.8755\pm 0.0003`$ for a reference epoch of 2000.0. The proper motion is $`\mu _{\alpha \mathrm{cos}\delta }=49.61\pm 0.06,\mu _\delta =29.8\pm 0.1`$ $`\text{mas yr}^1`$ and a parallax of $`3.4\pm 0.2`$ mas, equivalent to a distance of $`293_{17}^{+19}`$ pc.

### 3.1 Sources of error in phase referenced observations

#### 3.1.1 Source structure of the phase reference source

Whilst strictly the pulsar itself is the phase-reference in this experiment (with a binned flux of 5 and 0.6 Jy at 2.3 and 8.4 GHz) the effect of significant structure is the same in either case. If the source is non-symmetrical the incoming wavefront arriving at the baseline will not have the assumed intrinsic zero phase, degrading the source fitting. Therefore for our final observation we used the full array to allow imaging of both Vela-G and the pulsar. As we were looking for local structure we were able to use selfcalibration.

Following the Chandra observations which unveiled the fine X-ray structure around the pulsar we carefully searched for radio structure at VLBI resolutions in the off pulse data (calibrated with the on pulse phases). We found no discernible structure on baselines between 2.6 and 15.6 M$`\lambda `$ (equivalent to a full width half maximum of 16 – 95 mas) within 10′′ of the pulsar. The peak point source flux was 0.8 mJy/beam in line with expectations given the RMS of 0.1 mJy. This is of no great surprise as the radio structure around the pulsar has been reported to resolve out at resolutions finer than a few arcseconds (Dodson et al., 2003 Accepted; Lewis et al., 2002).

As pulsar radiospheres have angular sizes much less than 1 mas, they may be treated as a point source. (Scattering caused by the ISM may cause observable angular broadening of the pulsar image but is unlikely to bias the observed position significantly). The phase-reference source has not been optically identified but is most likely a distant AGN, and might well have source structure that produces systematic errors in the positions determined at each epoch for the pulsar. The fact that in most observations we used the Vela pulsar as the reference source for the weaker Vela G does makes no essential difference in this respect. With this in mind, the final epoch of our observations was used to image the reference source with the full LBA array. No structure brighter than 1 mJy was found around Vela-G (18 mJy) at 8.4 GHz, with an image RMS of 0.3 mJy. The phase residuals had a RMS of 4. We note here that Vela-G is a well known extra-galactic source and considered for the ICRF (candidate source ICRF J083322.3-444138) (Ma et al., 1998).

#### 3.1.2 Source structure of the ionospheric delay

A considerable phase error can be caused by the different ionospheric delay encountered over the different antennae. Several approaches exist to address these problems, the best being using multiple frequencies in the VLBI observations and solving for the ionosphere as part of the reduction (Brisken et al., 2000). Unfortunately the S2 system used in our experiments does not have enough spanned bandwidth to allow this approach to be used. We therefore we have looked into using the measured Total Electron Count (TEC) from GPS observations (Walker & Chatterjee, 1999). However the best data is on a 5 grid. The separation of our phase-reference and the target source is 0.7 and therefore GPS data cannot provide a useful correction.

As no correction was possible we modeled the data quality using DIFWRAP (Lovell, 2000) which allows an error estimate that includes all possible contributions without attempting to identify them. This approach involves exploring a range of model parameters to identifying the range of acceptable fits, which should be the bounds of the multi-dimensional one sigma contour. The errors found are indeed approximately equivalent to one sigma, as confirmed by fact that the proper motion fit has a reduced $`\chi ^2`$ of 1.1.

### 3.2 Space velocity

Two corrections need to be applied to our results to the true local motion of the Vela pulsar in its local environment. Our observations are directly tied to the ICRF, therefore the solar peculiar motion and the galactic rotation contribute to the observed proper motion, and need to be removed. We have used the solar constants from Dehnen & Binney (1998) of $`10\pm 0.36,5.25\pm 0.62,7.17\pm 0.38`$ $`\text{km s}^1`$ in galactic coordinates. We have used a flat rotation curve ($`\mathrm{\Omega }_0=220`$ $`\text{km s}^1`$, $`R_0=8.5`$ kpc; Fich et al. (1989)), which produces a local proper motion of $`5.4\mu _l`$ $`\text{mas yr}^1`$. This compares with $`5.7\mu _l`$ $`\text{mas yr}^1`$ from the local values for the Oort constants, as found by Feast & Whitelock (1997). The corrections are shown in the Table 2. The dominant source of error is from the uncertainties in the solar peculiar motion parameters. We have used the measured uncertainties in our observations and combined those with the models. We have ignored the possibilities of systematic errors or alternate models. This gives us an angular motion, at Vela’s local standard of rest, of $`\mu _\alpha =38.5\pm 1.2`$ $`\text{mas yr}^1`$, $`\mu _\delta =23\pm 1.5`$ $`\text{mas yr}^1`$ or $`\mu _{}=45\pm 1.3`$ $`\text{mas yr}^1`$ at a position angle of $`301^{}\pm 1.8`$. The pulsars’ transverse space velocity is therefore $`62\pm 2`$ $`\text{km s}^1`$. The PA no longer lies quite along the spin axis (e.g. $`310^{}\pm 1.5`$ Helfand et al. (2001), or $`307^{}\pm 2`$ Pavlov et al. (2001a)), which may strengthen the case that the impulse timescale was not quite long enough to average the off-axis component to zero.

## 4 Discussion

### 4.1 Fitting methods

We used the publicly available proper motion fitting routines, PMPAR333http://nacho.princeton.edu/w̃alterfb/pmpar/pmpar.html, created by W. Brisken. As each epoch had quite different observation spans and phase stability great care was required to ensure that the error estimates were accurate. The traditional approach has been to assume that the errors are a fraction of the beam size, but these ignore variations in the observing conditions during the experiment. We used DIFWRAP (Lovell, 2000) to measure the complete range of errors and derived error estimates which are realistic. Where we only had the archival positions, and not the data, we have used the median value of the difference between the formal errors and the errors found with DIFWRAP. It is particularly important to get the errors correct, as the effect to be measured is small and the variation between the data quality is large.

De Luca et al. (2000) point out that, as they calculate the proper motion using data collected on nearly identical day-numbers, their proper motion measurement is not contaminated by the parallax. They refer to this as a ‘pure’ proper motion. This is a concern with very short time baselines, which would blend the parallax and proper motion. We, however, have data spanning seven years and the correlation between the parameters is low.

### 4.2 Comparison of proper motion with other studies

Our calculated proper motion and parallax can also be compared with recent values obtained by optical proper motion studies and phase-referencing VLBI. We ignore the historic observations, which were blighted by poor resolution and low elevation, and concentrate on the two Radio VLBI observations, this one and Bailes et al. (1989) and the four optical observations, one purely ground based (Nasuti et al., 1997), one ground and HST (Markwardt & Ogelman, 1994) and two reports from the HST data set (De Luca et al., 2000; Caraveo et al., 2001) for which we take only the latest results. The results found by each of these studies is shown in Table 3. The table shows the proper motions in right ascension and declination for the pulsar along with the reference for this work. The proper motions listed do not account for the rotation of the Galaxy nor the peculiar motion of the Sun in the local galactic potential being only with reference to the calibrators used. The radio reference is an extra-galactic source, and the optical references are a significant number of field stars. Both styles of observations are internally consistent to within two standard deviations (see Figure 2), with improvement of errors over time. With the reduction of errors, however, the radio and optical results are steadily becoming less compatible.

The major advantage our observation have over the HST ones (other than the formal resolution being approximately one hundred times better) is that the reference is tied to a well defined reference frame. For the optical observations only field stars were available. The distance to the reference stars in the HST must be significantly greater than that of the pulsar, otherwise they would have had observable parallaxes themselves. Not knowing a distance we have assumed that they lie between two and ten kpc and calculated what the apparent proper motion would be for these limits using the standard flat rotation curve. We find that the galactic rotation contributes between -5 – -3 $`\text{mas yr}^1`$ to the proper motion in $`l`$. As all the sources would have proper motions within this range the scatter is within the HST errors of 1 mas. This contribution was not included in the calculations of Caraveo et al. (2001), and when it is the two sets of observations are consistent. Figure 2 includes a line representing the range that having the reference source at two to ten kpc would have contributed to our result. This means that the radio results are now consistent with all the optical results. This contribution, of course, needs to be removed to produce the correct space velocity and position angle for the pulsar.

## 5 Conclusions

We have measured the proper motion and parallax of the Vela pulsar to an unprecedented accuracy ($`\mu _{\alpha \mathrm{cos}\delta }=49.61\pm 0.06,\mu _\delta =29.8\pm 0.1`$ $`\text{mas yr}^1`$ , $`\pi =3.4\pm 0.2`$ mas), and have been able to convert these back to the space velocity and position angle of the pulsar in its local environment with greater precision that previously possible ($`62\pm 2`$ $`\text{km s}^1`$at $`301^{}\pm 1.8`$), because of the unambiguity in the radio reference frame. This allows the precise comparison of the Vela X-ray nebula symmetry axis and the proper motion of the Vela pulsar, opening insights into the timescale of the core collapse processes.

## 6 Acknowledgments

The Long Baseline Array is part of the Australia Telescope, funded by the Commonwealth of Australia and operated by ATNF and the University of Tasmania as a National Facility. This research has made use of the SIMBAD database operated at CDS, Strasbourg, France. This research has made use of NASA’s Astrophysics Data System Abstract Service. The author would like to expressly thank Warwick Wilson, head of Engineering for the ATNF, for always being prepared to listen to even the most outlandish correlator configuration request and often implementing them. Dr Brisken provided extremely helpful remarks and comments, as did the referee, Dr Caraveo, and we thank both.


text:
# Direct measurement of spatial Wigner function with area-integrated detection

## Abstract

We demonstrate experimentally a novel technique for characterizing transverse spatial coherence using the Wigner distribution function. The presented method is based on measuring interference between a pair of rotated and displaced replicas of the input beam with an area-integrating detector, and it can be superior in regimes when array detectors are not available. We analyze the quantum optical picture of the presented measurement for single-photon signals and discuss possible applications in quantum information processing.

The problem of characterizing spatial coherence of optical fields appears in a variety of applications. In biomedical optics, spatial coherence is being studied as a means to reveal the internal structure of scattering media, thus providing novel modes of tomographic imaging.WaxYangJOSAA02  In optical communications, spatial coherence of light sources is necessary to achieve efficient coupling into fiber systems. This aspect becomes even more critical in quantum communications,MichKiraSCI00 ; KurtMayePRL00 ; BeveBrouPRL02  where non-classical sources of radiation are used, and detection is often based on subtle interference effects that require good mode matching.

In this paper we present a novel technique for characterizing spatial coherence of light in the Wigner representation.WignPhR32 ; BastOpC78 ; Dragoman  Our technique has several features that make it distinct from previously proposed schemes. First, it is based on a series of measurements performed with an area integrating detector that collects all available light power, rather than on spatially resolved detection employed is some of the previous proposals.BrenLohmOpC82 ; McAlBeckOpL95  Therefore, it can be used in regimes where array detectors are not available, for example in mid- and far infrared, or when single photon signals are involved. In such cases, our technique provides the optimal signal-to-noise ratio, analogously to the case of Fourier transform spectroscopy which is more efficient than a spectrometer with a single scanning detector. Compared to the recently proposed heterodyne schemes for measuring the Wigner function,WaxThomOpL96 ; LeeReilOpL99  our method is self-referencing and it does not require any auxiliary sources of radiation.

Our technique is based on the following observation. For a complex quasi-monochromatic field $`E(𝝃)`$ varying in the transverse plane parameterized with the two-dimensional vector $`𝝃`$, the Wigner function is defined according toBastOpC78 ; Dragoman

$$W(𝒙,𝒌)=\frac{1}{\pi ^2}d^2𝝃e^{2ik𝝃}E^{}(𝒙𝝃)E(𝒙+𝝃).$$
(1)
Here $`𝒌`$ is the transverse spatial frequency vector and the angular brackets denote a statistical average. The above definition can be easily transformed to the form:

$$W(𝒙,𝒌)=\frac{1}{\pi ^2}d^2𝝃[e^{i𝒌𝝃}E(𝒙𝝃)]^{}e^{i𝒌𝝃}E(𝒙+𝝃).$$
(2)
This expression shows that the Wigner function at a point $`(𝒙,𝒌)`$ is given by a spatially integrated overlap of the field $`e^{i𝒌𝝃}E(𝒙+𝝃)`$ with a complex conjugate of its replica rotated by $`180^0`$, which corresponds to the transformation $`𝝃𝝃`$. The field $`e^{ik𝝃}E(𝒙+𝝃)`$ can be obtained from the original input by two simple transformation: displacement in space by $`𝒙`$ and changing the direction of propagation by $`𝒌`$, which in the paraxial approximation results in multiplying the field by a phase factor $`e^{ik𝝃}`$.

The above idea can be implemented to measure the Wigner function in a setup depicted in Fig. 1(a). The input beam is steered using the mirror M1 into a three-mirror Sagnac interferometer. The displacement and the tilt of the mirror control the point at which the Wigner function is measured. The 50:50 beam splitter BS generates a pair of replicas of the displaced and tilted input beam travelling in the opposite directions. The interferometer contains a Dove prism whose base forms $`45^0`$ with the plane of the interferometer. Transmission through the Dove prism rotates each of the counterpropagating beams by $`90^0`$, which adds up to the $`180^0`$ relative rotation required to measure the Wigner function. The operation of the setup can be most easily understood by analyzing propagation of the replicas of an exemplary image through the interferometer, as shown in Fig. 1(b). The counterpropagating replicas of the input field are recombined at the output port of the interferometer, and the emerging beam is focused with the help of a lens on an area-integrating detector.

The intensity $`I`$ recorded by the detector can be decomposed into a sum of three terms, $`I=I_1+I_2+I_{12}`$. The first two of them are given by the space-integrated input intensity $`|E(𝝃)|^2`$ and they remain constant as long as the aperture of the interferometer does not clip off any of the input field. The third term $`I_{12}`$, originating from the interference between the counterpropagating replicas is proportional to the right hand side of Eq. (2):

$$I_{12}=\frac{\pi ^2\mathrm{cos}\phi }{2}W(𝒙,𝒌)$$
(3)
thus giving the value of the Wigner function. Here $`\phi `$ is the relative phase between the interfering replicas. The complete Wigner function can be therefore scanned by measuring the detector photocurrent as a function of the position and the tilt of the steering mirror, and subtracting the constant pedestal. We note that the factor relating the interference term to the Wigner function depends on the relative phase $`\phi `$ between the interfering fields, and it can assume in principle either a positive or negative sign. We also note that if we extract the beam leaving through the second port of the interferometer, for example with the help of an optical circulator, then the constant pedestal can be removed by subtracting photocurrents measured at both the output ports of the interferometer.

We used the presented technique to measure the Wigner function for several beams derived from a He-Ne laser. For simplicity, we restrict ourselves only to one spatial dimension, parallel to the plane of the interferometer. The transverse distributions of the electric field for the prepared beams were factorizable in the directions parallel and perpendicular to the interferometer plane, so that the input field can be effectively considered as one dimensional. The steering mirror was mounted on a pair of translation and rotary stages (Physik Instrumente, models M–014.D0 and M–036.D0) controlled by a computer, which also recorded the value of the photodiode current. The interferometer was constructed using 1 inch diameter mirrors and a beam splitter, and a Dove prism 15 mm in cross-section. The beam emerging from the interferometer was focused with an $`f=25`$ mm lens on a detector with active area 3.6 mm $`\times `$ 3.6 mm.

In Fig. 2 we depict the measured Wigner functions. The first two plots were obtained for a laser beam filtered through a single-mode fiber and (a) collimated or (b) made lightly divergent. These data show Gaussian Wigner functions with diffraction-limited widths along the principal axes. The tilt of principle axis in Fig. 2(b) reflects positive correlations between the position and the direction of propagation characteristic for a diverging beam. Fig. 2(c) shows the Wigner function for an input obtained by a reflection from the front and the back surfaces of glass wedge. The measured Wigner distribution function shows two slightly converging Gaussian beams with different central transverse wave vectors. The oscillating pattern located between the Gaussians is a signature of their mutual coherence. We fitted our experimental data to a theoretical model of a Wigner function describing a pair of partially coherent Gaussian beams, which yielded the degree of coherence equal to $`1.02\pm 0.04`$. An additional feature seen in Fig. 2(c) is the presence of ripples on one of the Gaussian peaks, at twice the spatial frequency of the fringe pattern at the center of the distribution. These ripples can be attributed to interference with a third, much weaker beam generated by a secondary reflection from the glass wedge.

In the theoretical analysis presented above, we assumed that the range of transverse spatial frequencies in the input beam is narrow enough to make diffraction effects across the setup negligible. In a geneneral case, the measured Wigner function describes the distribution of the electric field in the plane in which the two replicas experience relative rotation. It is therefore critical that the replicas counterpropagating in the Sagnac interferomer are inverted at the same distance from the beam splitter BS that generated them. Diffraction after the Dove prism is irrelevant, as it does not change the total intensity which is the only quantity measured on the output. Also, when diffraction cannot be neglected, the parameterization of the Wigner distribution function obtained from the displacement and the tilt of the steering mirror needs to be compensated for the propagation distance from the steering mirror to the inversion plane. The compensation is given by a simple linear transformation of the Wigner function variables.BastOpC78 ; Dragoman

One important regime of optical measurements where array detectors are not easily available are single-photon signals. When the output of the interferometer is measured using a detector with single-photon sensitivity, for example an avalanche photodiode, the presented setup can be used to characterize transverse coherence of the output produced by single photon sources.MichKiraSCI00 ; KurtMayePRL00  Development of such sources is presently a subject of extensive research effort owing to prospective applications in quantum cryptography.BeveBrouPRL02  As secure quantum communication is most likely to be based on single-mode fiber links, efficient matching of the modal structure of the source to the fiber is necessary to achieve high transmission rates and close options for eavesdropping attacks.

The Wigner distribution function defined in Eq. (1) characterizes second-order statistical properties of fluctuating classical electromagnetic fields. However, it admits also a quantum mechanical interpretation which has a number of interesting consequences for our work. If the input signal consists of a single quasimonochromatic photon, then $`E(𝝃)`$, when properly normalized, is the wave function of that photon for the transverse spatial degrees of freedom.IBB  This correspondence has been recently used to propose encoding of information into spatial modes of a single photon characterized by different orbital angular momenta.LeacPadgPRL02  Separation of these modes was performed by interference between mutually rotated replicas in a manner similar to that underlying our scheme.

Following the quantum mechanical picture, the spatial degree of freedom of a photon is another optical realisation of a continuous-variable quantum system, alternative to the quadrature of a single light mode. There exists an analogy between the scheme proposed here and a direct measurement of the quadrature Wigner function of a light mode.WallVogePRA96 ; BanaRadzPRA99  In the latter scheme, the photon number parity operator is the equivalent of the spatial rotation between the interfering replicas that yields the Wigner function in the present case. This analogy has interesting consequences when generalized to composite systems. Consider a pair of photons with perfectly correlated positions, described by a pure entangled wave function. Such a pair exemplifies the Einstein-Podolsky-Rosen paradox, and it forms a source of entanglement. For example, it can be used to test Bell’s inequalities by sending the two photons to a pair of Sagnac interferometers and measuring the joint spatial Wigner function at selected points of the phase space. Such a measurement reveals quantum nonlocality, analogously to the joint measurement of displaced parity operators on a pair of quadrature entangled beams that has been discussed previously.BanaWodkPRA98

Concluding, we presented a method for measuring transverse spatial coherence of optical beams in the Wigner representation. The scheme is particularly beneficial in situations when array detectors are not available. It has a number of implications in quantum optics and quantum information processing: it can be used as a diagnostic tool for novel sources of nonclassical radiation, and also probe spatial entanglement in optical systems. This research was supported by ARDA Quantum Computing Program.


text:
# Absolute Wavelength Shifts – A new diagnostic for rapidly rotating stars

## 1. Introduction

Much of what we know about rotating stars is deduced from their spectra. However, with increased rotational velocity, spectral lines become smeared-out, making the measurement and interpretation of stellar spectral features more difficult.

Recent developments in astrometry now permit the accurate determination of stellar radial motion without involving spectroscopy or the Doppler principle. Although astrometry has normally not been associated with studies of stellar motion in the radial direction, astrometric quantities (e.g., proper motions) are sensitive to also the radial-velocity component (although only as a second-order effect). The accuracies recently realized in space astrometry have made it practical to determine such astrometric radial velocities (Dravins, Lindegren, & Madsen 1999; Madsen, Dravins, & Lindegren 2002). The differences between these values (giving motions of the stellar centers of mass with even sub-km s-1 accuracies) and apparent spectroscopic velocities reveal wavelength shifts intrinsic to stellar atmospheres, such as gravitational redshifts or those caused by gas flows on stellar surfaces.

With current techniques, the most accurate results are obtained for stars in open clusters, moving through space with essentially the same velocity vector. Parallaxes give the distances, while proper-motion vectors show the fractional change with time of the cluster’s angular size. The latter equals the time derivative of distance, yielding the radial velocity. Such young open clusters contain a significant number of rapidly rotating stars (Fig. 1).

## 2. Lineshifts intrinsic to stellar atmospheres

On Hipparcos, an observing program was carried out for stars in moving clusters, yielding astrometric radial-velocity solutions for more than 1,000 stars (Madsen et al. 2002), including about one hundred stars in the Hyades and Ursa Major groups. The error budgets are somewhat complex, but accuracies can reach around 0.5 km s-1 (Madsen 2003).

Differences between astrometrically determined radial velocities and the apparent spectroscopic ones reveal lineshifts intrinsic to stellar atmospheres. Fig. 2 indicates that F-type spectra have enhanced blueshifts relative to cooler ones (as actually expected from hydrodynamic simulations), a trend that seems to continue to even hotter stars (where corresponding hydrodynamic models are not yet developed).

The likely origin of such wavelength shifts can be traced back to stellar surface structure. While this can be studied in detail on the Sun, corresponding phenomena in other stars are mainly seen as subtle asymmetries in photospheric spectral lines. Bright (hot) elements in solar (and stellar) granulation are rising, causing local blueshifts. A spectral line averaged over the stellar disk then obtains a convective blueshift since the bright (and blueshifted) elements statistically contribute a greater number of photons. Three-dimensional hydrodynamic models reproduce both the granulation patterns and the resulting wavelength shifts (Dravins & Nordlund 1990; Allende Prieto et al. 2002).

### 2.1.  Effects for rotating stars

Increased stellar rotation produces an increased centrifugal force and lowers the effective surface gravity. Analogous to modeled differences between hydrodynamic atmospheric models for dwarfs and lower-gravity subgiants, more vigorous granulation might then be expected to develop (since the convective energy flux must then be carried by a lower-density gas). The greater velocity amplitudes and temperature contrasts become visible in integrated starlight as enhanced convective blueshifts. Further effects could perhaps arise from Coriolis forces, which might influence the granulation structure in rapidly rotating stars.

A tantalizing possibility is the potential of observing meridional flows across stellar surfaces. A circulation pattern between the equator and poles will involve rising (or sinking) gases at the equator, with the opposite flows near the poles. Stars seen equator-on (with a large observed value of Vsin$`i`$) will then have their spectra influenced by Doppler shifts from the systematic upflow (or downflow) patterns there, while stars seen pole-on (with a small value of Vsin$`i`$) will show lineshift signatures of the opposite sign.

Since various effects of atmospheric dynamics will combine, a segregation of the various contributions will require a detailed physical modeling, in particular exploiting the fact that the exact wavelength displacements depend on the properties of various spectral lines. For example, convective blueshifts differ among lines of different strength (formed at different depths), of different excitation potential and ionization level (formed preferentially in differently hot elements), and of wavelength region (with different intensity contrasts between the hotter and cooler elements).

For stars in the Hyades, a Vsin$`i`$ dependence of the intrinsic wavelength shift is observed – spectra of rapidly rotating stars appear blueshifted relative to those of slower rotators (Fig. 3). While such a dependence is thus perhaps not unexpected, any interpretations must still be viewed with caution since this dependence also correlates with that on temperature (spectra of earlier-type stars are blueshifted relative to those from cooler ones; Fig. 2), and rapid rotation statistically dominates for earlier types.

#### Spectroscopic effects of stellar rotation

The accurate determination of spectroscopic wavelength shifts in rotating stars can involve some complications since the precise value may depend on precisely how the wavelength shifts are measured. For example, radial-velocity instruments often apply a cross correlation of the observed spectrum with a spectrum template taken from a slowly rotating star. However, for a rapidly rotating early-type star, a typical template may have a mismatch exceeding 1 km s-1, due mainly to the rotational broadening and the ensuing blending of spectral lines (Verschueren & David 1999; Griffin, David, & Verschueren 2000).

Even modest rotational velocities in sharp-lined late-type stars may cause wavelength displacements of the spectral-line bottoms used for radial-velocity determinations. For sufficiently rapid rotation, when asymmetric line components originating near the stellar limbs begin to affect the wings of the profile integrated over the stellar disk, the intrinsic line asymmetries may become enhanced (Gray & Toner 1985; Gray 1986; Smith, Livingston, & Huang 1987; Dravins & Nordlund 1990).

A conclusion to be drawn is (contrary to a sometimes expressed belief) that spectroscopic observations of rapidly rotating stars need to be made not only at a high photometric signal-to-noise ratio, but also at very high spectral resolution – much higher than that corresponding to the stellar rotational velocity itself.

## 3.  A new diagnostic for rapidly rotating stars

In the past, absolute lineshifts could be studied only for the Sun (since the relative Sun-Earth motion is known from planetary-system dynamics, and does not rely on spectral measurements). For other stars, this is now becoming possible thanks to three separate developments: (a) Astrometric measurements permit to accurately determine stellar radial motion without using any spectroscopy; (b) The availability of high-resolution spectrometers with accurate wavelength calibration (such as designed for exoplanet searches), and (c) Accurate laboratory wavelengths for several atomic species. These developments appear particularly significant for the study of stars with complex or smeared-out spectra, such as rapidly rotating ones, or such with large-scale atmospheric motions or mass loss. For a further discussion, see Dravins (2003).

##### Acknowledgments.

This work is supported by The Swedish Research Council, the Swedish National Space Board, and The Royal Physiographic Society in Lund.

## References

Allende Prieto, C., Asplund, M., García López, R. J., Lambert, D. L. 2002, ApJ, 567, 544

Dravins, D. 2003, in N.E. Piskunov, W.W. Weiss, & D.F. Gray (eds.), Modelling of Stellar Atmospheres, Proc. IAU Symp. No. 210 (San Francisco: ASP), in press

Dravins, D., Lindegren, L., Madsen, S. 1999, A&A, 348, 1040

Dravins, D., Lindegren, L., Madsen, S., Holmberg, J. 1997, in ESA SP-402, Proc. Hipparcos - Venice ’97, p. 733

Dravins, D., Nordlund, Å. 1990, A&A, 228, 203

Głȩbocki, R., Stawikowski, A. 2000, Acta Astron., 50, 509

Gray, D.F. 1986, PASP, 98, 319

Gray, D.F. Toner, C.G. 1985, PASP, 97, 543

Griffin, R.E.M., David, M., Verschueren, W. 2000, A&AS, 147, 299

Madsen, S. 2003, A&A, in press

Madsen, S., Dravins, D., Lindegren, L. 2002, A&A, 381, 446

Smith, M.A., Livingston, W., Huang, Y. 1987, PASP, 99, 297

Verschueren, W., David, M. 1999, A&AS, 136, 591


text:
# Binary evolution models with rotation

## 1. Introduction

The evolution of a single star can be strongly influenced by its rotation (e.g., Heger & Langer 2000; Meynet & Maeder 2000), and evolutionary models of rotating stars are now available for many masses and metallicities. While the treatment of the rotational processes in these models is not yet in a final stage (e.g., magnetic dynamo processes are just about to be included; Heger et al., this volume), they provide first ideas of what rotation can really do to a star.

Effects of rotation, as important they are in single stars, can be much stronger in the components of close binary systems: Estimates of the angular momentum gain of the accreting star in mass transferring binaries show that critical rotation may be reached quickly (Packet 1981; Langer et al. 2000). Therefore, we need binary evolution models which include a detailed treatment of rotation in the stellar interior, as in recent single star models. However, in binaries, tidal processes as well as angular momentum accretion need to be considered at the same time. Some first such models are now available and are discussed below.

These models provide evidence for rotational processes in binaries being essential for some of the most exciting cosmic phenomena, which may occur exclusively in binaries: Type Ia supernovae, the main producers of iron and cosmic yardstick to measure the accelerated expansion of the universe. Gamma-ray bursts from collapsars — which current stellar models with rotation preclude to occur in single stars (Heger et al., this volume; Woosley & Heger, this volume) — may provide the most powerful explosions in the universe and trace star formation to its edge. And black holes, even though also formed by single stars, need a companion to become “visible”.

## 2. How much matter can stars accrete from a binary companion?

As mentioned above, non-magnetic accretion, i.e. accretion via a viscous disk or via ballistic impact, transports angular momentum and can lead to a strong spin-up of the mass gaining star. For disk accretion, it appears plausible that the specific angular momentum of the accreted matter corresponds to Kepler-rotation at the stellar equator; this leads to a spin-up of the whole star to critical rotation when its initial mass is increased by about 20% (e.g., Packet 1981). Can accretion continue beyond this?

Theoretically, this appears possible, as viscous processes may transport angular momentum outwards through the star, the boundary layer, and the accretion disc (e.g., Paczynski 1991). However, as the star may be strongly rotating, its wind mass loss may be dramatically increased (Langer 1997, 1998), which may render the mass transfer process inefficient.

Observations of massive post-mass transfer binary systems constrain this effect. Table 1 lists parameters of four different kinds of massive close binary systems which give opposite answers. The two O stars in the Case A binary (mass transfer starts while both stars undergo core hydrogen burning) V 729 Cyg have a mass ratio of 3.5 but the same spectral type and visual flux. Clearly, an initial mass ratio close to 1 is required to get close to the observed current mass ratio. However, as during Case A the primary star (as we designate the initially more massive star in a binary) loses just about half of its mass, a mass ratio of at most 2 could be produced were the secondary (the initially less massive star in a binary) not allowed to accrete. Another system showing strong evidence for accretion is the massive X-ray binary Wray 977; it would require that stars of $`>40\mathrm{M}_{}`$ form neutron stars to explain this system without accretion (Wellstein & Langer 1999).

Several Galactic short period WNE+O binaries, on the other hand, can not be understood had the O star accreted substantial amounts from the WNE progenitor (Petrovic & Langer 2002). While those might have formed through common envelope evolution — for which little accretion is expected — the key X-ray binary 4U 1700-37 has such a short period that a major accretion phase can be excluded. However, as Case C evolution (mass transfer starts after core helium burning of the primary star) would lead to a compact object much more massive than 2.4$`\mathrm{M}_{}`$, Case B evolution (mass transfer starts just after core hydrogen burning of the primary star) is most likely here (Clark et al. 2002).

We conclude from these observed binary systems that in some massive close (Case A or B) binaries, the mass transfer process is nearly conservative, while in others it is strongly non-conservative. In the following, we show that including rotational physics into binary evolution calculations allows to recover both features within the same physical model.

## 3. Binary models with rotating components

We have constructed binary evolution models using the code of Wellstein et al. (2001), but including the physics of rotation as in the single star models of Heger et al. (2000) for both components. In addition, spin-orbit coupling according to Zahn (1977) has been added, and rotationally enhanced winds are implemented as in Langer (1998). The specific angular momentum of the accreted matter is assumed to be that of Kepler rotation at the stellar equator in the case of disk accretion, and determined by integrating the equation of motion of a test particle in the Roche potential in case the accretion stream impacts directly on the secondary star (Wellstein 2001).

Fig. 1 shows the evolution of the equatorial rotation velocity in a system starting out with a $`16\mathrm{M}_{}`$ and a $`15\mathrm{M}_{}`$ star in a 3 day orbit. The initial rotational velocity of both stars is unimportant since they evolve quickly into rotation which is synchronous with the orbital revolution, due to spin-orbit coupling. Each of the three thermal time scale mass transfer phases through which this system evolves (rapid Case A, Case AB, and Case ABB; see Wellstein et al. 2001) leads to a strong spin-up of the secondary star and an equally drastic spin-down of the primary (see Langer 1998, for the purely mechanical spin-down effect).

The accretion induced spin-up can bring the secondary close to critical rotation and thus strong rotationally enhanced mass loss sets in. The difference between the mass overflow rate $`\dot{M}_1`$ and the wind mass loss rate of the secondary is the net accretion rate $`\dot{M}_2`$ of the secondary star. Fig. 2 shows that during much of the rapid Case A mass transfer — during which the accreting star increases its mass from 15$`\mathrm{M}_{}`$ to $`23\mathrm{M}_{}`$ —, it is $`\dot{M}_1=\dot{M}_2`$. This is possible due to two factors: 1) During direct impact accretion, the specific angular momentum of the accreted material is only a fraction of the respective Kepler angular momentum. 2) As the secondary star fills a significant fraction of its Roche volume, tidal coupling removes part of the accreted angular momentum and feeds it into the orbit.

As the orbit widens significantly during the later evolution, the specific angular momentum of the accreted matter increases and tidal forces weaken. Thus, while the accretion efficiency is close to 1 during Case A, it is less than 0.1 later, resulting in a time average over the whole evolution of 0.67. For systems which start out with a wider orbit, we find a low accretion efficiency throughout (Fig. 3). For $`16\mathrm{M}_{}+15\mathrm{M}_{}`$ systems, we find a critical initial period of $`8`$d beyond which the secondary accretes only little.

These models imply that, despite the angular momentum problem for the accretion star (Packet 1981), quasi-conservative evolution of massive close binaries is possible. However, already in early Case B systems, the accretion efficiency may be strongly reduced compared to binary models without rotation.

## 4. Black hole formation in binaries

While single stars more massive than $`25\mathrm{M}_{}`$ are supposed to form black holes (Fryer 2002), it is a challenge to binary evolution modellers to find a way to make black hole binaries — binary systems consisting of a non-degenerate star and a black hole. The reason is that the initially more massive star, which is supposed to end up as black hole, undergoes tremendous mass loss (Fig. 2, and see Wellstein & Langer 1999).

Brown, Lee & Bethe (1999) suggested the Case C channel to form the Galactic low mass BH binaries — a low mass star in close orbit with a black hole of typically 6$`\mathrm{M}_{}`$ (Orosz 2002), i.e., starting with a wide system consisting of a (say) 25$`\mathrm{M}_{}`$ and a $`5\mathrm{M}_{}`$ star. This suggestion is strongly supported by the massive post-Case B X-ray binary 4U 1700-37 (Table 1): it shows that in Case B systems, even a 60$`\mathrm{M}_{}`$ primary forms a compact object of only 2.4$`\mathrm{M}_{}`$ (Clark et al. 2002).

The problem of how wide the Case C channel is for such massive stars is difficult. Current single star models suggest no significant radius increase late during or after core helium burning at solar metallicity. On the other hand, Smartt’s (2002) abundance analysis of Galactic supergiants suggests that they do not undergo “blue loop” evolution — which contradicts the current models. The problem is less severe at low metallicity, where many models predict a gradual radius increase during core helium burning.

It is even harder to produce Cyg X-1 type BH binaries, with black hole masses as high as 14$`\mathrm{M}_{}`$. Brown et al. (2001) argue that this requires Case C evolution and an upturn of the initial-final mass relation for single stars above $`60\mathrm{M}_{}`$, as predicted by Langer (1987).

## 5. Spins of compact objects and gamma-ray bursts

The inclusion of rotationally enhanced magnetic fields in the evolution of massive single stars (Heger et al., this volume) has improved the agreement between observed ($`15\mathrm{}150`$ms) and predicted ($`>5`$ms) spin periods of young neutron stars — with the consequence that obtaining gamma-ray bursts from black hole formation in single stars (collapsars; Woosley & Heger, this volume) seems difficult at present.

From Sections 3 and 4 above, we conclude that the initially more massive stars in massive close binaries are even less likely to produce a gamma-ray burst. First of all, they lose so much mass that even stars with a very large initial mass may not even form a black hole but rather a neutron star (see also Wellstein & Langer 1999). And secondly, Figure 1 shows how drastic these stars are spun down as a consequence of their heavy mass loss. The 16$`\mathrm{M}_{}`$ star in the computed binary system is expected to produce a neutron star with an initial spin period of more than one second!

The only way to avoid both drawbacks is to employ Case C evolution, which leads to a core evolution as in single stars. However, even then the CO-core of the star needs to be spun-up significantly to produce a collapsar and a gamma-ray burst — a possibility suggested by Brown et al. (1999, 2000) in the context of common envelope evolution and spiral-in. No detailed models for this scenario exist at present.

The initially less massive star in a massive binary, on the other hand, accretes large amounts of angular momentum and will thus acquire a larger core spin than a corresponding single star (see Fig. 1). It is thus conceivable that accretion stars are the progenitors of asymmetric supernova explosions and rapidly spinning compact objects. Those sufficiently massive to transform into a Wolf-Rayet star during core helium burning, or those which lose their envelope in a reverse Case C mass transfer, may possess all required ingredients to produce a gamma-ray burst within the collapsar model.

## 6. Progenitors of Type Ia supernovae

Supernovae of Type Ia are supposed to occur in close binary systems where the accreting component is a CO-white dwarf. The currently favoured scenario assumes that a white dwarf of initially less than $`1\mathrm{M}_{}`$ grows in mass due to quiescent nuclear shell burning of hydrogen and/or helium, which is accreted from a non-degenerate companion. In order to explode, it must reach the Chandrasekhar-mass, i.e., it needs to accrete at least $`0.4\mathrm{M}_{}`$. While single white dwarfs are slow rotators (Kawaler, this volume), the idea that Type Ia supernova progenitors are significantly spun-up due to accretion appears to be supported by recent measurements of the rotational velocity of the white dwarf in cataclysmic variables (Starrfield, this volume).

In order for the degenerate CO-core to grow in mass, accreting white dwarfs need to possess a helium burning shell source. However, like in AGB stars, which have a degenerate CO-core in their center, the helium shell source in accreting white dwarfs tends to turn thermally unstable, causing large oscillations in nuclear energy production, as well as in the white dwarf radius (Fig. 4). It is unclear whether the binary system can survive these oscillations intact.

However, Yoon et al. (2003) found that the consideration that angular momentum is accreted along with the matter by the white dwarf results in different shell source properties. The spin-up of the outer part of the white dwarf leads to a lower effective gravity and to a more relaxed helium shell source. Additionally, differential rotation leads to shear mixing at the location of the shell source, increasing its geometrical thickness and often turning it stable. This effect increases the likelihood of a given system to reach the Chandrasekhar-mass at all.

However, not only does the thickness of the helium shell source increase, but the chemical profiles through it become shallower (Fig. 5). As thus most of the helium burning occurs at a lower helium concentration, oxygen is produced at a much higher abundance than in corresponding non-rotating models, and carbon is produced correspondingly less. As the combustion of oxygen yields less energy than that of carbon, rotational effects influence the brightness of Type Ia supernovae. As the rotation rate depends on the amount of accreted material, and the average amount of accreted matter increases with the initial metallicity of the binary system (Langer et al. 2000), a drop of the brightness of Type Ia supernovae with higher metallicity — or an increase with redshift — might occur.

Finally, rotation also affects the evolution of white dwarfs in the so called sub-Chandrasekhar mass scenario, as outlined by Yoon & Langer (2003; also this volume).

## References

Bohannan B., Conti P.S., 1976, ApJ 204, 797

Brown G.E., Lee C.-H., Bethe H.A., 1999, New Astron. 4, 313

Brown G.E., Lee C.-H., Wijers R.A.M.J., Lee H.K., Israelian G., Bethe H.A., 2000, New Astron. 5, 191

Brown G.E., Heger A., Langer N., Lee C.-H., Wellstein S., Bethe H.A., 2001, New Astron. 6, 457

Clark J.S., Goodwin S.P., Crowther P.A., Kaper L., Fairbairn M., Langer N., Brocksopp C., 2002, A&A 392, 909

Fryer C.L., in A massive star odyssey: from main sequence to supernova, proc. IAU-Symp. 212, (San Francisco: ASP), K.A. van der Hucht, A. Herrero & C. Esteban, eds., in press

Heger A., Langer N., 2000, ApJ 544, 1016

Heger A., Langer N., Woosley S.E., 2000, ApJ 528, 368

Kaper L., Lamers H.J.G.L.M., Ruymaekers E., van den Heuvel E.P.J., Zuidervijk E.J., 1995, A&A 300, 446

Langer N., 1987, A&A 171, L1

Langer N., 1997, in Luminous Blue Variables: Massive Stars in Transition, A. Nota, H.J.G.L.M. Lamers, eds, ASP Conf. Ser., Vol. 120, p. 332

Langer N., 1998, A&A 329, 551

Langer N., Deutschmann A., Wellstein S., Höflich P., 2000, A&A 362, 1046

Meynet G., Maeder A., 2000, A&A 361, 101

Orosz J.A., 2002, in A massive star odyssey: from main sequence to supernova, proc. IAU-Symp. 212, (San Francisco: ASP), K.A. van der Hucht, A. Herrero & C. Esteban, eds., in press

Packet W., 1981, A&A 102, 17

Paczynski B., 1991, ApJ 370, 597

Petrovic J., Langer N., 2002, in A massive star odyssey: from main sequence to supernova, proc. IAU-Symp. 212, (San Francisco: ASP), K.A. van der Hucht, A. Herrero & C. Esteban, eds., in press

Scheithauer S., 2000, Diploma thesis, Potsdam University

Smartt S.J., 2002, in A massive star odyssey: from main sequence to supernova, proc. IAU-Symp. 212, (San Francisco: ASP), K.A. van der Hucht, A. Herrero & C. Esteban, eds., in press

Wellstein S., 2001, PhD thesis, Potsdam University

Wellstein S., Langer N., 1999, A&A 350, 148

Wellstein S., Langer N., Braun H., 2001, A&A 369, 939

Yoon S.-C., Langer N., 2003, A&A, in preparation

Yoon S.-C., Langer N., 2002, in: The Physics of Cataclysmic Variables and Related Objects, ASP Conf. Ser. Vol 261, 2002, B.T. Gaensicke et al., eds., p 79

Yoon S.-C., Langer N., Scheithauer S., 2003, A&A, in preparation

Zahn J.-P., 1977, A&A 57, 383