text:
# An efficient and accurate quantum algorithm for the Dirac equation

(March 25, 2002)

## Abstract

An efficient quantum algorithm for the many-body three-dimensional Dirac equation is presented. Its computational complexity is dominantly linear in the number of qubits used to spatially resolve the 4-spinor wave function.

Dirac equation, quantum computing, quantum lattice gas, many-body relativistic quantum mechanics  preprint: Draft Revision 1.2

The first quantum algorithm to compute a path integral was introduced by Feynman in 1965. At that time he did not refer to it as such since it was not until in 1982 that he proposed the idea of using a quantum computer to efficiently simulate quantum mechanical many-body dynamics feynman-collection . In the second chapter of his manuscript on path integrals published with Hibbs feynman-65-1st-qlga , the problem is given to prove that the one-dimensional (1D) Dirac equation can be modeled by summing over all the possible zigzag paths of a single-speed particle traveling at the speed of light in a discrete two-dimensional space-time hopping from lattice site to lattice site and flipping left or right according to a deterministic local interaction rule. The amplitude a particular path contributes to the kernel is proportional to the number of its “collisions” or corners. In this way, the time evolution of the 2-spinor field of a single quantum particle is modeled by a “gas” of particles computing all paths simultaneously in a time-explicit fashion. This discrete and parallelized process of local collisions and the streaming along the lattice is described by a continuous effective field theory, the 1D Dirac equation, in the limit of the infinite lattice resolution. A solution to Feynman’s “quantum lattice gas” problem was published in 1984 by Jacobson and Schulman jacobson-jpamg84 . Here we give a three-dimensional (3D) solution.

In 1994, Bialynicki-Birula proposed a discrete model of the 3D Dirac equation implemented on a body-centered cubic lattice bialynicki-birula-prd94 . However, this model is 1st-order convergent (doubling the grid resolution merely doubles the numerical accuracy), problematic when modeling particle dynamics in an external potential. Although the model is unitary, it is specified using non-unitary matrices and requires ad hoc lattice partitioning if implemented in parallel. Furthermore, Bialynicki-Birula addresses only the one-body problem. Meyer published a subequent series of papers on the 1D quantum lattice gas algorithm, equivalent to Bialynicki-Birula’s algorithm, and cast in the form of Feynman’s original model meyer-collection . Meyer contributed 1D one-body numerical simulations and addressed the non-interacting lattice or checker-board problem using an additional rest particle. Yet he too did not address the many-body case nor the low-order numerical convergence issue.

Contemporaneously with Meyer, Succi published a series of papers on this subject emphasing the analogy between quantum mechanics and fluid mechanics: the connection between the Dirac equation and the Schroedinger equation to that between the kinetic Boltzmann equation and the Navier-Stokes equation of hydrodynamics succi-collection . Succi’s quantum lattice gas model on a cubic lattice for the 3D Dirac equation has, at the “kinetic” level, the particles undergoing mixing during free propagation and is again similar to Bialynicki-Birula’s model. Succi discusses the many-body case, but his algorithm runs into an “exponential complexity wall” succi-cpc02 .

Our quantum lattice-gas algorithm for the 3D Dirac equation is suited to direct implementation on a quantum computer using only two-qubit quantum gates and efficiently handles the many-body problem. For pedagoical purposes, first we state the simplest quantum lattice-gas algorithm on a cubic lattice. Then, we introduce an improved version that remedies two difficiencies: the checkboard problem of non-interacting sublattices and the low-order convergence. Finally, we recast our quantum algorithm to handle the many-body case in a second-quantized representation.

The relativistic quantum mechanical wave equation for a free particle is the linear Dirac equation

$$_t\psi =c\underset{i}{}\alpha _i_i\psi +i\frac{mc^2}{\mathrm{}}\beta \psi ,$$
(1)
where $`\psi `$ is a 4-spinor and the matrices $`\alpha _i`$ and $`\beta `$ satisfy the constraints $`\alpha _i^2=1`$, $`\beta ^2=1`$, $`\{\alpha _i,\alpha _j\}=0`$, and $`\{\beta ,\alpha _i\}=0`$ so that (1) is equivalent to the Klein-Gordon wave equation. Since the $`2\times 2`$ Pauli matrices $`\sigma _i`$ for $`i=x,y,z`$ satisfy these constraints, we can express $`\alpha _i`$ and $`\beta `$ as tensor products $`\alpha _i=a\sigma _i`$ and $`\beta =b\mathrm{𝟏}`$, where $`a`$, $`b`$, can be any two different Pauli matrices and where $`\mathrm{𝟏}`$ is the $`2\times 2`$ identity matrix. We choose $`a=\sigma _z`$ and $`b=\sigma _x`$. Then, the Dirac equation is

$$_t\psi =c\underset{i}{}\sigma _z\sigma _i_i\psi +i\sigma _x1\frac{mc^2}{\mathrm{}}\psi .$$
(2)
With the wave function defined on an infinite resolution cubical lattice at times separated by an infinitesimal duration $`\delta t`$ with the grid cell size the infinitesimal length $`\delta rc\delta t`$, the Heisenberg evolution

$$\psi ^{}=\psi +\delta \psi =e^{\mathrm{\Sigma }_i\sigma _z\sigma _i\delta r_ii\frac{mc^2}{\mathrm{}}\delta t\sigma _x\mathrm{𝟏}}\psi $$
(3)
corresponds exactly to (1) in the relativistic limit where $`\mathrm{}\omega mc^2`$ and $`\mathrm{}kmc`$.

The $`4\times 4`$ matrix

$$\sigma _z\sigma _z=\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right),$$
(4)
operating with the $`z`$ spatial derivative in (3) is diagonal whereas the matrices $`\sigma _z\sigma _x`$ and $`\sigma _z\sigma _y`$ for the $`x`$ and $`y`$ partial derivatives, respectively, are not diagonal.

We would like to transform (3) in such a way that all the matrices operating with the spatial partial derivatives are diagonal (and hence correspond to infinitesimal shifting along the orthogonal lattice directions). To do this, we need the two identities:

$$e^{i\frac{\pi }{4}\sigma _x}e^{\epsilon \sigma _z}e^{i\frac{\pi }{4}\sigma _x}=e^{\epsilon \sigma _y}e^{i\frac{\pi }{4}\sigma _y}e^{\epsilon \sigma _z}e^{i\frac{\pi }{4}\sigma _x}=e^{\epsilon \sigma _x},$$
(5)
that follow from $`e^{i\frac{\pi }{4}\sigma _i}=\frac{1}{\sqrt{2}}\left(1+i\sigma _i\right)`$ provided $`\epsilon `$ is infinitesimal. Then, using the identity $`\mathrm{𝟏}e^{i\theta a}=e^{i\theta \mathrm{𝟏}a}`$, the 2-spinor similarity transformations (5) can be generalized to 4-spinor transformatons

$`\left(\mathrm{𝟏}e^{i\frac{\pi }{4}\sigma _x}\right)e^{\epsilon \sigma _z\sigma _z}\left(\mathrm{𝟏}e^{i\frac{\pi }{4}\sigma _x}\right)`$ $`=`$ $`e^{\epsilon \sigma _z\sigma _y}`$ (6)
$`\left(\mathrm{𝟏}e^{i\frac{\pi }{4}\sigma _y}\right)e^{\epsilon \sigma _z\sigma _z}\left(\mathrm{𝟏}e^{i\frac{\pi }{4}\sigma _y}\right)`$ $`=`$ $`e^{\epsilon \sigma _z\sigma _x},`$
which we will use to diagonalize the $`x`$ and $`y`$ spatial derivative operators in (3). Using (6) and defining

$$X_\theta ^{(1)}e^{i\theta \sigma _x}\mathrm{𝟏}=\left(\begin{array}{cccc}\mathrm{cos}\theta & 0& i\mathrm{sin}\theta & 0\\ 0& \mathrm{cos}\theta & 0& i\mathrm{sin}\theta \\ i\mathrm{sin}\theta & 0& \mathrm{cos}\theta & 0\\ 0& i\mathrm{sin}\theta & 0& \mathrm{cos}\theta \end{array}\right),$$
(7)
$$X_\theta ^{(2)}\mathrm{𝟏}e^{i\theta \sigma _x}=\left(\begin{array}{cccc}\mathrm{cos}\theta & i\mathrm{sin}\theta & 0& 0\\ i\mathrm{sin}\theta & \mathrm{cos}\theta & 0& 0\\ 0& 0& \mathrm{cos}\theta & i\mathrm{sin}\theta \\ 0& 0& i\mathrm{sin}\theta & \mathrm{cos}\theta \end{array}\right),$$
(8)
and

$$Y_\theta ^{(2)}\mathrm{𝟏}e^{i\theta \sigma _y}=\left(\begin{array}{cccc}\mathrm{cos}\theta & \mathrm{sin}\theta & 0& 0\\ \mathrm{sin}\theta & \mathrm{cos}\theta & 0& 0\\ 0& 0& \mathrm{cos}\theta & \mathrm{sin}\theta \\ 0& 0& \mathrm{sin}\theta & \mathrm{cos}\theta \end{array}\right),$$
(9)
and $`S_ie^{\sigma _z\sigma _z\delta r_i}`$, the spatial displacement operators in the Heisenberg representation of the evolution equation (3) can be written

$`e^{\sigma _z\sigma _x\delta r_x}`$ $`=`$ $`Y_{\frac{\pi }{4}}^{(2)}S_xY_{\frac{\pi }{4}}^{(2)}e^{\sigma _z\sigma _y\delta r_y}=X_{\frac{\pi }{4}}^{(2)}S_yX_{\frac{\pi }{4}}^{(2)}`$ (10)
$`e^{\sigma _z\sigma _z\delta r_x}=S_z,`$
so the evolution equation itself can be rewritten as

$$\psi ^{}=Y_{\frac{\pi }{4}}^{(2)}S_xY_{\frac{\pi }{4}}^{(2)}X_{\frac{\pi }{4}}^{(2)}S_yX_{\frac{\pi }{4}}^{(2)}S_zX_{\frac{mc^2\delta t}{\mathrm{}}}^{(1)}\psi .$$
(11)
This has the form of a quantum lattice-gas algorithm with local interaction (“collision”) operators $`X_{\frac{\pi }{4}}^{(2)}`$, $`Y_{\frac{\pi }{4}}^{(2)}`$ and $`X_{\frac{mc^2\delta t}{\mathrm{}}}^{(1)}`$, as well as lattice-directed displacement (“streaming”) operators $`S_x`$, $`S_y`$, and $`S_z`$.

For numerical purposes, we would like to represent the wave function on a finite resolution grid with cell size $`\delta r\mathrm{\Delta }r`$ and update time $`\delta t\mathrm{\Delta }t`$. In this approximation, the wave function becomes a discrete field that exists only at the spacetime grid points $`\stackrel{}{x}_{\mathrm{}}`$ for $`\mathrm{}=1\mathrm{}L^3`$ and $`t_n`$ for $`n=0,1,2,\mathrm{}`$

$$\psi (\stackrel{}{x}_{\mathrm{}},t_n)=\left(\begin{array}{c}\alpha (\stackrel{}{x}_{\mathrm{}},t_n)\\ \beta (\stackrel{}{x}_{\mathrm{}},t_n)\\ \mu (\stackrel{}{x}_{\mathrm{}},t_n)\\ \nu (\stackrel{}{x}_{\mathrm{}},t_n)\end{array}\right)$$
(12)
and the operators $`S_i`$ for $`i=x,y,`$ or $`z`$ induce a finite displacement $`S_i\psi (\stackrel{}{x})\psi (\stackrel{}{x}_{\mathrm{}}+\sigma _z\sigma _z\mathrm{\Delta }r\widehat{x}_i)`$ of the components of the 4-spinor only along lattice directions:

$$S_x\psi (x_{\mathrm{}},y_{\mathrm{}},z_{\mathrm{}})=\left(\begin{array}{c}\alpha (x_{\mathrm{}}+\mathrm{\Delta }r,y_{\mathrm{}},z_{\mathrm{}})\\ \beta (x_{\mathrm{}}\mathrm{\Delta }r,y_{\mathrm{}},z_{\mathrm{}})\\ \mu (x_{\mathrm{}}\mathrm{\Delta }r,y_{\mathrm{}},z_{\mathrm{}})\\ \nu (x_{\mathrm{}}+\mathrm{\Delta }r,y_{\mathrm{}},z_{\mathrm{}})\end{array}\right),$$
(13)
$$S_y\psi (x_{\mathrm{}},y_{\mathrm{}},z_{\mathrm{}})=\left(\begin{array}{c}\alpha (x_{\mathrm{}},y_{\mathrm{}}+\mathrm{\Delta }r,z_{\mathrm{}})\\ \beta (x_{\mathrm{}},y_{\mathrm{}}\mathrm{\Delta }r,z_{\mathrm{}})\\ \mu (x_{\mathrm{}},y_{\mathrm{}}\mathrm{\Delta }r,z_{\mathrm{}})\\ \nu (x_{\mathrm{}},y_{\mathrm{}}+\mathrm{\Delta }r,z_{\mathrm{}})\end{array}\right),$$
(14)
and

$$S_z\psi (x_{\mathrm{}},y_{\mathrm{}},z_{\mathrm{}})=\left(\begin{array}{c}\alpha (x_{\mathrm{}},y_{\mathrm{}},z_{\mathrm{}}+\mathrm{\Delta }r)\\ \beta (x_{\mathrm{}},y_{\mathrm{}},z_{\mathrm{}}\mathrm{\Delta }r)\\ \mu (x_{\mathrm{}},y_{\mathrm{}},z_{\mathrm{}}\mathrm{\Delta }r)\\ \nu (x_{\mathrm{}},y_{\mathrm{}},z_{\mathrm{}}+\mathrm{\Delta }r)\end{array}\right).$$
(15)
These streaming operators are classical data shifting operators causing global permutations of the components of the 4-spinor wave function across the lattice and on a quantum computer can be implemented by 2-qubit local swap operators yepez-pre99 . The collision operators act independently on each node of the lattice and cause local quantum entanglement between component pairs of the 4-spinor. The streaming operators in turn propagate this local on-site entanglement to next nearest neighbors so that eventually quantum entanglement covers the entire lattice.

It is possible to rewrite (11) as a finite difference equation on a body-centered cubical lattice. The resulting set of coupled finite difference equations are similar to the finite difference representation of the 3D Dirac equation given by Bialynicki-Birula in 1994 bialynicki-birula-prd94 . A drawback of expressing the algorithm as a finite-difference equation is its unsuitability for a quantum computer implementation using two-qubit quantum gates whereas our manifestly unitary expression (11) is suitable.

A continuous effective field theory for $`\psi =(\alpha ,\beta ,\mu ,\nu )`$ follows in the continuum limit of the emergent finite-difference equations by Taylor expanding in $`\mathrm{\Delta }rx_{\mathrm{}+1}x_{\mathrm{}}`$ and in $`\mathrm{\Delta }tt_{n+1}t_n`$. We obtain

$`_t\left(\begin{array}{c}\alpha \\ \beta \\ \mu \\ \nu \end{array}\right)`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Delta }r}{\mathrm{\Delta }t}}_x\left(\begin{array}{c}\beta \\ \alpha \\ \nu \\ \mu \end{array}\right)+i{\displaystyle \frac{\mathrm{\Delta }r}{\mathrm{\Delta }t}}_y\left(\begin{array}{c}\beta \\ \alpha \\ \nu \\ \mu \end{array}\right)`$
$`+{\displaystyle \frac{\mathrm{\Delta }r}{\mathrm{\Delta }t}}_z\left(\begin{array}{c}\alpha \\ \beta \\ \mu \\ \nu \end{array}\right)+i{\displaystyle \frac{mc^2}{\mathrm{}}}\left(\begin{array}{c}\mu \\ \nu \\ \alpha \\ \beta \end{array}\right)+𝒪(c\mathrm{\Delta }r,\mathrm{\Delta }t),`$
which is exactly the Dirac equation (1) when $`\mathrm{\Delta }t\mathrm{\Delta }r\epsilon `$ are infinitesimal and when the partial derivative with respect to time is defined as $`_t\psi \frac{\psi ^{}\psi }{\mathrm{\Delta }t}`$. (11) gives rise to perfectly unitary evolution of the discretized wave function and, therefore, is an unconditionally stable numerical algorithm. The effective field theory (An efficient and accurate quantum algorithm for the Dirac equation) is 1st-order convergent in space.

Our basic approach to improve the accuracy of the quantum algorithm is to set the grid size $`\mathrm{\Delta }r`$ to be smaller than the Compton wavelength $`\lambda =\frac{h}{mc}`$ of the modeled particle

$$\mathrm{\Delta }r\epsilon \frac{h}{mc},$$
(17)
and to introduce a small temporal scale that is much smaller than $`\frac{\lambda }{c}`$

$$\mathrm{\Delta }t\epsilon ^2\frac{h}{mc^2}.$$
(18)
The diffusive ordering condition of spatial and temporal fluctuations typical of random walk processes, $`\mathrm{\Delta }r^2=\nu \mathrm{\Delta }t`$, provides a context to understand the scaling behavior of the small parameter $`\epsilon `$. According to (17) and (18), the diffusive transport coefficient is $`\nu =\frac{h}{\epsilon m}`$ and the particle velocity is $`\frac{\mathrm{\Delta }r}{\mathrm{\Delta }t}=\frac{c}{\epsilon }`$, which approaches infinity as $`\epsilon 0`$. In this limit, the velocity of the modeled quantum particle is relatively small, hence the resulting effective field theory should correspond to the non-relativistic limit of the Dirac equation as $`\epsilon 0`$.

To diagonalize the streaming operators in (10), we used a fixed and finite rotation angle $`\frac{\pi }{4}`$ independent of the grid resolution. We will now diagonalize the streaming operators using a small rotation angle proportional to $`\mathrm{\Delta }t`$. By (18), the rotation angle is $`\theta =\frac{mc^2\mathrm{\Delta }t}{h}=\epsilon ^2`$, which is dependent on the grid resolution. The displacement operators in the Dirac equation (1) can be represented by interleaving streaming and collision operators on a cubical lattice as follows:

$$e^{\sigma _z\sigma _x\delta r_x}E_xS_x^{2,4}Y_{\frac{\epsilon }{2}}^{(2)}S_x^{2,4}Y_{\frac{\epsilon }{2}}^{(2)}S_x^{1,3}Y_{\frac{\epsilon }{2}}^{(2)}S_x^{1,3}Y_{\frac{\epsilon }{2}}^{(2)}$$
(19)
and

$$e^{\sigma _z\sigma _y\delta r_y}E_yS_y^{2,4}X_{\frac{\epsilon }{2}}^{(2)}S_y^{2,4}X_{\frac{\epsilon }{2}}^{(2)}S_y^{1,3}X_{\frac{\epsilon }{2}}^{(2)}S_y^{1,3}X_{\frac{\epsilon }{2}}^{(2)},$$
(20)
where the superscripts on the streaming operators refer to individual components of the 4-spinor. The streaming operators $`S_i=S_i^{2,3}S_i^{1,4}`$ in (13) and (14) are now separated by collision operators. This kind of interleaving of streaming and collision operators removes the spurious check-board effect of independent sublattice dynamics that otherwise occurs yepez-ijmpc00b ; yepez-cpc01 . So far we treated the non-diagonal operators $`e^{\sigma _z\sigma _x\delta r_x}`$ and $`e^{\sigma _z\sigma _y\delta r_y}`$, but not the displacement operator $`e^{\sigma _z\sigma _z\delta r_z}`$ because no such improvement exists since it is diagonal. However, if instead of using the Dirac matrix $`\sigma _z\sigma _z`$, we use an alternative non-diagonal representation for the $`z`$-direction partial derivative, then we can again employ interleaving. Therefore, we consider this alternate form of the Dirac equation

$$_t\psi =c\sigma _z\sigma _x_x\psi +c\sigma _z\sigma _y_y\psi +c\sigma _y\mathrm{𝟏}_z\psi +i\sigma _x1\frac{mc^2}{\mathrm{}}\psi .$$
(21)
Now the displacement operator in (21) for the $`z`$-direction can be re-expressed in a fashion similar to (19) and (20) as

$$e^{\sigma _y\mathrm{𝟏}\delta r_z}E_zS_z^{2,3}X_{\frac{\epsilon }{2}}^{(1)}S_z^{2,3}X_{\frac{\epsilon }{2}}^{(1)}S_z^{1,4}X_{\frac{\epsilon }{2}}^{(1)}S_z^{1,4}X_{\frac{\epsilon }{2}}^{(1)}.$$
(22)
Then instead of (10), we use (19), (20), and (22) for an improved quantum algorithm

$$\psi (t+\mathrm{\Delta }t)=E_xE_yE_zX_\epsilon ^{(1)}\psi (t).$$
(23)
In (23) we have appended a collision operator $`X_\epsilon ^{(1)}`$ to produce the mass term in the Dirac equation.

It is possible to derive a finite-difference equation representation of the quantum lattice-gas algorithm (23) by carrying out all the collision and streaming operations symbolically. The result is no longer expressible on the body-centered cubic lattice. Nevertheless, once again, a continuous effective field theory for $`\psi =(\alpha ,\beta ,\mu ,\nu )`$ follows in the continuum limit and Taylor expanding in $`\mathrm{\Delta }r`$ and in $`\mathrm{\Delta }t`$:

$`_t\left(\begin{array}{c}\alpha \\ \beta \\ \mu \\ \nu \end{array}\right)`$ $`=`$ $`{\displaystyle \frac{mc^2}{\mathrm{}}}[\mathrm{\Delta }r_x\left(\begin{array}{c}\beta \\ \alpha \\ \nu \\ \mu \end{array}\right)+i\mathrm{\Delta }r_y\left(\begin{array}{c}\beta \\ \alpha \\ \nu \\ \mu \end{array}\right)`$
$`+`$ $`i\mathrm{\Delta }r_z\left(\begin{array}{c}\mu \\ \nu \\ \alpha \\ \beta \end{array}\right)+i\left(\begin{array}{c}\mu \\ \nu \\ \alpha \\ \beta \end{array}\right)]+𝒪(\mathrm{\Delta }r^2,\mathrm{\Delta }t),`$
which approximates the Dirac equation (21) when $`\mathrm{\Delta }t`$ is small. The effective field theory (An efficient and accurate quantum algorithm for the Dirac equation) is 1st-order convergent because of the error term $`𝒪(\mathrm{\Delta }t)`$. With the evolution operator $`E=E_xE_yE_z`$, we define the dual operator $`\stackrel{~}{E}E_x^{}E_y^{}E_z^{}`$, by taking the adjoint of the collision operators and reversing the streaming directions. Then, it is possible use a symmetrized evolution operator yepez-cpc01

$$\psi (t+\mathrm{\Delta }t)=\stackrel{~}{E}Ee^{\mathrm{\Delta }t^2}\psi (t),$$
(25)
which is better than 2nd-order accurate in space, as demonstrated in Fig. 1.

Our quantum algorithm for the many-body Dirac equation can be expressed in terms of 2-qubit gates that conserve particle number acting on an initial ket with 4 qubits per lattice node, $`|\mathrm{\Psi }=_{r=1}^{L^3}|q_1(r)|q_2(r)|q_3(r)|q_4(r)`$. With $`\widehat{a}_\alpha ^{}`$, $`\widehat{a}_\alpha `$, and $`\widehat{n}=\widehat{a}_\alpha ^{}\widehat{a}_\alpha `$ denoting the creation, annihilation, and number operator, respectively, of the $`\alpha `$th qubit ($`1\alpha 4L^3`$), the collision operators are

$`\widehat{X}_{\alpha \beta }`$ $`=`$ $`\mathrm{𝟏}i\mathrm{sin}\theta (\widehat{a}_\alpha ^{}\widehat{a}_\beta +\widehat{a}_\beta ^{}\widehat{a}_\alpha )+(\mathrm{cos}\theta 1)(\widehat{n}_\alpha +\widehat{n}_\beta )`$ (26)
$`2\mathrm{cos}\theta \widehat{n}_\alpha \widehat{n}_\beta `$
$`\widehat{Y}_{\alpha \beta }`$ $`=`$ $`\mathrm{𝟏}+\mathrm{sin}\theta (\widehat{a}_\alpha ^{}\widehat{a}_\beta \widehat{a}_\beta ^{}\widehat{a}_\alpha )+(\mathrm{cos}\theta 1)(\widehat{n}_\alpha +\widehat{n}_\beta )`$ (27)
$`2\mathrm{cos}\theta \widehat{n}_\alpha \widehat{n}_\beta ,`$
where $`\alpha `$ and $`\beta `$ index different qubits at the same site. Then (7-9) are rewritten as $`X_\theta ^{(1)}\widehat{X}_{13}\widehat{X}_{24}`$, $`X_\theta ^{(2)}\widehat{X}_{12}\widehat{X}_{34}`$, and $`Y_\theta ^{(2)}\widehat{Y}_{12}\widehat{Y}_{34}`$. Hence, $`2L^3`$ applications of either $`\widehat{X}_{\alpha \beta }`$ or $`\widehat{Y}_{\alpha \beta }`$ are required for a single collision step. Streaming occurs by successive application of the interchange operator $`\widehat{S}_{\mu \nu }=\mathrm{𝟏}+\widehat{a}_\mu ^{}\widehat{a}_\nu +\widehat{a}_\nu ^{}\widehat{a}_\mu \widehat{n}_\mu +\widehat{n}_\nu `$ yepez-pre99 . $`(L1)^3`$ number of applications of $`\widehat{S}_{\mu \nu }`$ ($`\mu `$ refers to one qubit-component at some site and $`\nu `$ to the same component at its neighboring site) are required to stream one qubit-component along a cubic lattice direction. The total evolution operator $`\widehat{E}`$ is the product of collision operators $`\widehat{X}`$ and $`\widehat{Y}`$ and streaming operators $`\widehat{S}`$ corresponding to algorithm (11) or some variant of (23) depending on the desired degree of numerical accuracy. With the new ket $`|\mathrm{\Psi }^{}(t+\mathrm{\Delta }t)=\widehat{E}|\mathrm{\Psi }(t)`$, the resulting probability of finding a particle at site $`\stackrel{}{x}`$ is $`P(\stackrel{}{x})=_{i=1}^4\mathrm{\Psi }^{}|\widehat{n}_{\alpha +i}|\mathrm{\Psi }^{}`$, where $`\alpha `$ the index of the 1st qubit at $`\stackrel{}{x}`$.

The computational complexity of one time step scales as $`C=\rho _c2L^3+\rho _s(L1)^3`$, where $`\rho _c`$ and $`\rho _s`$ are the number of operations per node for collisions and streaming. For the simplest algorithm (11), $`\rho _c=5`$ and $`\rho _s=12`$, and for the improved algorithm (23), $`\rho _c=13`$ and $`\rho _s=24`$. Both $`\rho _c`$ and $`\rho _s`$ double when we use a symmetrized rule like (25) but are a fixed-cost overhead. With $`Q=4L^3`$ qubits, the size of the Hilbert space is exponential $`2^Q`$, whereas the complexity $`C=\frac{\rho _c}{2}Q+\rho _s[Q\frac{3}{4}(2Q)^{\frac{2}{3}}+\frac{3}{2}(2Q)^{\frac{1}{3}}1]`$ for all the versions of our quantum algorithm is dominantly linear in $`Q`$.

I would like to thank the Air Force Office of Scientific Research for supporting this work.


format: text

text:
# The Sub-Luminal Parsec-scale Jet of Mkn 501

## 1 Introduction

The BL Lac object Mkn 501 (1652+398, J1653+3945, DA426) has been well-studied at radio wavelengths (e.g. Mufson et al., 1984; Gabuzda et al., 1992), but interest in the source was rejuvenated by the discovery of TeV gamma-ray emission (Quinn et al., 1996). A prolonged high-state at TeV energies in 1997 included activity on timescales of several hours, implying the TeV gamma-rays originate in a relatively compact area (see, e.g., the review of Catanese & Weekes, 1999). The detection of correlated X-ray and TeV gamma-ray variability in the other well-studied TeV source, Mkn 421, is strong evidence in favor of the X-ray emission being the high end of the synchrotron component of the spectral energy distribution, with the TeV emission arising from inverse Compton scattering of photons by the synchrotron emitting electrons, and such a model is also widely accepted for Mkn 501.

On the parsec scale, VLBI observations have revealed that a jet emerges from the core at a position angle of $``$180, and bends by $``$90 within the first $``$2 milli-arcseconds (mas). The jet extends to the east until $``$20 mas from the core, when it bends further, finally reaching the position angle of $``$45 seen on the kilo-parsec scale (Conway & Wrobel, 1995; Giovannini et al., 1999). The parsec-scale jet of Mkn 501 is one-sided, and it is assumed that this jet is relativistically doppler boosted, while the counter-jet is doppler de-boosted to such an extent it is invisible at the current sensitivity of VLBI observations. The jets are believed to originate in the accretion disk surrounding a central super-massive black hole, which for Mkn 501 has been suggested to have a mass of 108.93±0.21 M (Falomo, Kotilainen & Treves 2002).

A number of different motions have been reported for components in the Mkn 501 jet, ranging from 0.27$`\pm `$0.02 mas yr-1 (Gabuzda et al., 1994) to 2.4 mas yr-1 (Giovannini et al., 1999). These studies have relied on only a small number of epochs, typically less than four, and, in hindsight, have probably under-estimated the errors in locating jet components. As illustrated by Piner et al. (1999) for Mkn 421, reliable determination of component motions generally requires larger, multi-epoch data sets.

Here we study the parsec-scale jet of Mkn 501 from twelve VLBI observations spanning 4.28 years, at frequencies of 5, 8 and 15 GHz. Mkn 501 lies at a redshift of 0.034 (Wills & Wills, 1974) which, for the value of H$`{}_{0}{}^{}=65`$ km s-1 Mpc-1 adopted throughout his paper, corresponds to a distance of 155 Mpc. At this distance, an angular separation of 1 mas corresponds to a projected linear distance of 0.72 pc.

## 2 VLBI observations

We have compiled data from three programs for this study: the Very Long Baseline Array (VLBA) 2 cm survey Kellermann et al. (1998, see also http://www.cv.nrao.edu/2cmsurvey), the Radio Reference Frame Image Database (Fey & Charlot, 1997, see also http://rorf.usno.navy.mil/rrfid.shtml), and a single VLBI Space Observatory Programme (VSOP; Hirabayashi et al., 1998, 2000) observation. The data used in this study are summarized in Table 1.

The VLBA 2 cm survey is being undertaken at multiple epochs to study the properties and evolution of over 100 active galactic nuclei. The Mkn 501 observations consisted of eight scans at one hour intervals of typically 5 minutes duration. Data were recorded with a bandwidth of 64 MHz using 1 bit samples and left-circular polarization. An image derived from the March 1997 data used in this paper was published by Kellermann et al. (1998). The full ten-station array was used at all epochs except the last, for which the Saint Croix telescope was unavailable.

The Radio Reference Frame Image Database (RRFID) of the U.S. Naval Observatory (USNO) is a program to regularly image the radio sources used for precise astrometry. The VLBA was used at all epochs, although the October 1995 observation was made without the Mauna Kea and North Liberty telescopes. The June 1998 observation was made with the addition of the Fairbanks 26 m (Alaska), Green Bank 20 m (West Virginia), Kokee Park 20 m (Hawaii), Medicina 32 m (Italy), Ny Alesund 20 m (Norway), Onsala 20 m (Sweden), and Westford 18 m (Massachusetts) telescopes. Typically, four scans of $``$3 minutes were made, with bandwidths of 16 MHz for the first two epochs, 32 MHz for the 8 GHz observations of the last two epochs, and 64 MHz for the 15 GHz observation in April 1996. Right circular polarization is recorded for all RRFID observations. An image from the 1995 April 12 epoch at 8.4 GHz was published by Fey & Charlot (1997).

The 5 GHz VSOP observation, in April 1998, was made over a 13 hour period with the HALCA satellite, the VLBA and the Effelsberg 100 m (Germany) telescope. Interferometric fringes to the satellite were detected from tracking passes totaling 7 hours. In the standard VSOP observing mode, 32 MHz of two-bit sampled, left circular polarization data is recorded. Although made at the lowest frequency considered here, the long baselines to the orbiting telescope result in the synthesized beam-size for this observation being the smallest of these data. The VSOP data considered here were combined with a 1.6 GHz VSOP observation one day later to derive a spectral index map of the source (Edwards et al., 2000a).

## 3 Analysis

The data were fringe-fit in AIPS and we have imaged all data ourselves using the Difmap package. Beams were calculated using natural weighting (uvweight=0,$``$1 in Difmap), with the exception of the VSOP observation for which uniform weighting (uvweight=2,0 in Difmap) is more appropriate (Hirabayashi et al., 2000). An image from the 15 GHz RRFID observation in April 1996 is shown in Figure 1.

Model-fitting of the images was carried out in Difmap. As inspection of Figure 1 reveals, at most epochs four jet components were required in addition to the core to provide a good representation of the data. We have labeled these C1 through C4, with C1 being the component farthest from the core. Circular gaussian components were fitted at all epochs. Full details of the model-fits are given in Table 2. Reduced $`\chi ^2`$ values for the fits are not given in the Table, as they are dependent on the way the data from the different programs was reduced. Thus, while model fits represent a minimum in $`\chi ^2`$ for the given number of components, the comparison of $`\chi ^2`$ values between data from different programs is potentially misleading (see also Piner et al., 1999).

The positions of all model-fit jet components are plotted as a function of time in Figure 2. In order to determine the uncertainty in the component location for the purposes of determining component motions, we have projected the beam major-axis onto the line joining the component and core and then taken a fraction of this projected length as the error in position. For the extended, outermost component, C1, we conservatively adopted half a projected beamwidth for the uncertainty. For C2 and C3 we used one quarter of the projected beamwidth, and for the innermost component, C4, we used one-eighth of the beamwidth. Motions were determined by weighted linear fits to the data. As shown in Figure 2, both C1 and C4 show little evidence of motion over the 4.28 year period. Both C2 and C3 show clear evidence of motion, with an apparent component speeds for C2 of 0.6$`\pm `$0.1 $`c`$, and for C3 of 0.3$`\pm `$0.1 $`c`$. The apparent speeds along the jet, as opposed to radial separations, are, within errors, the same. These speeds supersede the preliminary values reported in Piner et al. (2002).

## 4 Discussion

Our component locations agree well with those reported from a contemporaneous 5 GHz observation made in June 1996 as part of the VLBA Pre-launch Survey (Fomalont et al., 2000). We can also extrapolate our derived motions and compare them with model-fits to the 5 GHz observations at epochs 1987.4 (Gabuzda et al., 1992), and 1989.3 (Gabuzda et al., 1994). The extrapolated motion of C3 is consistent with the positions of the K2 of Gabuzda et al. (1992) and Gabuzda et al. (1994), assuming uncertainties $``$3 times larger than the $`\pm `$0.1 mas adopted by these workers. The K1 of Gabuzda et al. (1992) and Gabuzda et al. (1994), with similarly increased uncertainties, is consistent with the extrapolated motion of C2, particularly if the speed of C2 lies at the lower end of the range determined in §3 (assuming a constant motion). These identifications were qualitatively suggested by Edwards et al. (2000b), but are quantitatively borne out by the fuller analysis presented here.

The speeds of C4 and C1 are formally consistent with zero, i.e. they appear to be stationary components. Before considering this further, we reconcile this result for C4 with the observations of Marscher (1999), who reported the detection at 22 GHz of a resolved component between 0.5 and 1 mas from the core, with an apparent motion derived from three epochs between April and August 1997 of 0.96$`\pm `$0.1 mas yr-1, corresponding to 2.3$`\pm `$0.2 $`c`$. This location is consistent with our C4, however we do not see such rapid motion over the four year period. The motion reported by Marscher corresponds to 0.27 mas in the 0.29 years the observations spanned. In our data, C4 ranges between 0.59 mas from the core (1998 April 7) and 0.83 mas from the core (1999 July 19), a range of 0.21 mas, similar in magnitude to that of Marscher (1999). Any attempt at further interpretation is complicated by the fact that there are likely to be frequency-dependent offsets in the separation of components from the core (see, e.g., Lobanov, 1998) in our data, which would be most important for C4.

Stationary components have been reported for a number of sources in the past, with a detailed study being made as part of the multi-epoch monitoring program of Jorstad et al. (2001a). The monitoring revealed that the super-luminal speeds detected for these of EGRET-detected blazars were much faster than for the general population of bright compact radio sources, however evidence was also found for at least one stationary component in 27 of the 42 sources (Jorstad et al., 2001a). Jorstad et al. suggested that the stationary components within several parsecs of the core were associated with standing recollimation shocks caused by pressure imbalances at the boundary between the jet and the surrounding medium. In contrast, the stationary components further from the core tended to be associated with bends in the parsec-scale jet. There is support for this scenario in our data. C4 is located at a projected distance of $``$0.5 pc from the core, and is quite plausibly associated with a recollimation shock. C1, on the other hand, is an extended component, which 1.6 GHz VLBI imaging has revealed is associated with a significant change in the jet from a bright “spine” to a limb-brightened morphology (Giovannini et al., 1999).

If we assume our fastest observed pattern speed (0.6$`c`$ for C2) reflects the bulk apparent speed of the jet, then we can solve for the intrinsic speed and angle to the line-of-sight, provided we also have an estimate of the Doppler beaming factor. A Doppler factor $`\delta 10`$ is inferred from the TeV observations of this source (e.g., Tavecchio, Maraschi, & Ghisellini 1998). The VSOP observations yield our best measurement of the radio-core brightness temperature, 4$`\times 10^{11}`$ K. This is consistent with a Doppler factor of $``$10 if the source is in equipartition (Readhead, 1994), but it is also consistent with lower Doppler factors if equipartition is violated (Kellermann, 2002). If we accept the values of 0.6$`c`$ and 10 for the apparent bulk speed and Doppler factor, then the Lorentz factor of the Mkn 501 jet is $`\gamma =5`$ ($`v=0.98c`$) and its angle to the line-of-sight is $`\theta =0.7\mathrm{°}`$. Such a small angle to the line-of-sight may be expected of a $`\gamma `$-ray blazar, although subluminal apparent speeds are in general not expected (see the Monte Carlo simulations of Lister, 1999).

Alternative kinematics that do not place such tight constraints on the angle to the line-of-sight assign the Doppler factor measurement to the TeV-emitting region (on the light-day size scale) and the apparent bulk speed to the VLBI jet (on the light-year size scale), and allow a change in the bulk Lorentz factor or angle to the line-of-sight in the intermediate region. If the jet in the TeV-emitting region has, e.g., $`\theta =5\mathrm{°}`$ and $`\gamma =7`$ (enforcing $`\delta =10`$), then a decrease in the Lorentz factor to $`\gamma =2`$ would reproduce the observed apparent speed in the VLBI jet. Such a deceleration of electron-positron jets close to the core is proposed by Marscher (1999) for the TeV blazars. A change in angle to the line-of-sight, perhaps accompanying the large bend in the jet seen $``$2 mas from the core, cannot by itself reproduce the observed values; a jet with $`\delta =10`$ has a minimum Lorentz factor of 5, and a jet with $`\gamma =5`$ can only have an apparent speed of 0.6$`c`$ in the large-angle solution for $`\theta >90\mathrm{°}`$. Any set of kinematic parameters must also be constrained by the one-sided appearance of the source; the example above with $`\theta =5\mathrm{°}`$ and $`\gamma =2`$ would have a jet-to-counterjet brightness ratio greater than $``$200, somewhat higher than the limit that can be placed from our observations.

Similar values for the Doppler factor and apparent jet speed apply to the other well-studied TeV blazar, Mkn 421 (Piner et al., 1999). From these two sources, it appears that TeV blazars as a class may either have very small angles to the line-of-sight ($`\theta <1\mathrm{°}`$), or may decelerate significantly between the TeV-emitting region and the parsec scale.

It is notable that no new component has emerged from the core after the prolonged TeV high-state in 1997. A component with a speed similar to that of C2 or C3 would now be $``$0.5 mas from the core and would have been detected at the latter epochs. This suggests that events which give rise to extended TeV (and associated X-ray) activity are different in nature to those which result in the production of new VLBI components (see also Marscher, 1999). Mkn 421 and Mkn 501, for which the inverse-Compton component of the spectral energy distribution (SED) peaks at TeV energies, have sub-luminal component speeds and apparently no new component emerging after epochs of TeV activity. In contrast, sources with the inverse-Compton component of the SED peaking at GeV energies tend to have the emergence of new, super-luminal, VLBI components associated with GeV flaring states (Jorstad et al., 2001b). The detection of more TeV gamma-ray sources by the next generation of air Cerenkov telescopes will enable these apparent trends to be investigated more quantitatively.

This research has made use of the United States Naval Observatory Radio Reference Frame Image Database (RRFID), and the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. Ken Kellermann and Alan Fey are particularly thanked for the provision of fringe-fit data from the VLBA 2cm survey and RRFID observations, respectively. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. We gratefully acknowledge the VSOP Project, which is led by the Institute of Space and Astronautical Science in cooperation with many organizations and radio telescopes around the world. BGP acknowledges support from Whittier College’s Newsom Endowment.


text:
# Cosmic Distances: Current Odds and Future Perspectives

## 1. Introduction

Dating back to Baade (1944, 1958) the interplay between investigations on stellar populations and on radial variables as primary distance indicators is a well-established route to improve the knowledge on the stellar content of simple and complex stellar systems. During the last half century paramount observational and theoretical efforts have been devoted to single out the evolutionary and pulsation properties of population I and population II stars. In particular, spectroscopic measurements disclosed that the metallicity was not the leading physical parameter to split Baade’s populations. Actually, Galactic RR Lyrae stars present a metallicity distribution that ranges from $`[Fe/H]2.2`$ up to solar metallicities (Preston 1964; Lub 1977). At the same time, Classical Cepheids have been discovered in dwarf galaxies such as Phoenix and IC 1613 that are 1 dex more metal-poor ($`[Fe/H]1.9`$; $`[Fe/H]1.3`$) than the Large and the Small Magellanic Cloud (LMC, $`[Fe/H]0.3`$; SMC, $`[Fe/H]0.7`$). This is the reason why we subdived the two most popular groups of variables located in the Cepheid instability strip in low and intermediate-mass radial variables.

The microlensing experiments (EROS, OGLE, MACHO) provided an unprecedented amount of photometric data. The number of variables for which pulsation parameters (periods, mean magnitude and colors, amplitudes) are available have doubled and in some cases increased by one order of magnitude. In this context, HST also played a crucial role, and indeed the two fundamental projects on stellar distances succeeded in the measurement of Cepheids in several spirals of the Virgo cluster (Freedman et al. 2001; Saha et al. 2001). This unique opportunity allowed a substantial improvement in the calibration of secondary distance indicators, and in the evaluation of the Hubble constant.

This is the coarse-grained scenario and no doubts that the perspectives for the near future are quite promising. However, the fine-grained scenario still presents several unsettled problems. In the following we discuss the impact that they might have on cosmic distance determinations and on stellar populations.

## 2. Classical Cepheids and intermediate-mass stars

Classical Cepheids are the most popular distance indicators, since they are bright stars and they can be easily identified in the Galaxy, in Local Group (LG) galaxies, and in spiral galaxies of the Virgo and Fornax clusters. Pros and cons of the Cepheid distance scale have been discussed countless in the recent literature (Sandage et al. 1999; Tanvir 1999; Feast 1999; Freedman et al. 2001). Plain physical arguments suggest that individual Cepheid distances can be obtained using the Period-Luminosity-Color (PLC) relation. The use of the PLC relation is hampered by two problems: i) two accurate mean magnitudes are required; ii) the use of optical bands does require an accurate knowledge of the reddening. To overcome the first problem the two $`H_0`$ projects decided to properly sample the V band and then to adopt a constant amplitude ratio to improve the accuracy of the mean I magnitude. This empirical trick is supported by theoretical predictions. Two different routes can be followed to overcome the latter problem: i) near-infrared (NIR) magnitudes are marginally affected by reddening. However, we still lack accurate mean J and K magnitudes for Cepheids in external galaxies. Ground-based NIR observations with 8m class telescopes and/or H-band NICMOS observations will become available in the near future, but certainly NGST will play a major role to improve both the accuracy and the sample size. ii) the Wesenheit magnitudes (W=V-2.45(V-I)) are reddening free. This approach mimics the use of a PLC relation in the V, V-I band, but the color coefficient attain slightly different values. The strengths and weaknesses mentioned before are the main reasons why Cepheid distances are estimated using the Period-Luminosity (PL) relation. The main problem in using the PL relation is its intrinsic width in effective temperature. This thorny problem is severe in the optical but negligible in the NIR bands (Madore et al. 1987).

On top of these problem there is the long-standing debate concerning the universality of both PL and PLC relations, i.e. the fact that the zero-point and the slopes of these relations might depend on the metal content. During the last few years this crucial point has been addressed in several theoretical and empirical investigations (Sasselov et al. 1997; Kennicutt et al. 1998; Macri et al. 2001). These investigations do not agree on the magnitude of the effect, but they agree on the sign. In particular, they found that metal-rich Cepheids, at fixed period, are brighter than metal-poor ones. Theoretical models constructed by adopting different physical assumptions concerning the coupling between radial displacements and convection (Chiosi, Wood, & Capitanio 1993; Bono et al. 1999a,b; Alibert et al. 1999) also do agree on the sign among themselves but do not agree on the magnitude of the metallicity dependence. What is really puzzling is that both linear and nonlinear Cepheid models are at odds concerning the sign, since they predict that metal-rich Cepheids, at fixed period, are fainter than metal-poor ones. Two pieces of evidence somehow support theoretical predictions: i) theory and observations suggest that metal-rich Cepheids are redder than metal-poor ones; ii) recent detailed analysis of V, I, K-band data of a sizable sample of MC Cepheids support the theoretical sign (Groenewegen & Oudmaijer 2000; Groenewegen 2000; Storm et al. 2000). This problem requires a firm solution, since it affects the calibration of both zero-point and slope. The approach adopted by Pietrzynski et al. (2002) to measure Cepheids on a substantial portion of NGC 300 is quite promising to settle this problem.

However, recent theoretical and empirical results do suggest that the metal content marginally affects the PL and the PLC relation of First Overtone (FO) Cepheids. As a matter of fact, Bono et al. (2002) found that predicted and empirical $`PL_K`$ relations (K band data from the Two Micron All Sky Survey) and Wesenheit functions (V,I band data from the Optical Gravitational Lensing Experiment) supply mean distances to the Magellanic Clouds (MCs) that agree very well with each other. In particular they estimated a distance for the LMC of $`18.53\pm 0.08`$ (theory) and $`18.48\pm 0.13`$ (obser.), as well as for the SMC of $`19.04\pm 0.11`$ (theory) and $`19.01\pm 0.13`$ (obser.). The reasons for the agreement are manifolds: i) K-band data and Wesenheit function are marginally affected by uncertainties on reddening corrections; ii) The $`PL_K`$ relation and the Wesenheit function presents a mild dependence on metal content (Bono et al. 1999b). iii) The previous intrinsic features do apply to fundamental (F) Cepheids, however, they are magnified for FO Cepheids, since the width in temperature of the latter pulsators is systematically smaller than for the former ones. Therefore, distances based on these variables are less affected by the typical spread of PL relations.

The improvement in using FOs to constrain the accuracy of the Cepheid distance scale is undoubtful. The simultaneous detection of both F and FO Cepheids in nearby stellar systems should allow us to estimate on a quantitative basis the systematic errors affecting optical and NIR PL/PLC relations. Unfortunately, current photometric data for FO Cepheids are scanty, and indeed we still lack detailed information for FOs in nearby dwarfs such as Phoenix (Caldwell et al. 1988) but also for spiral galaxies such as M 31. Fig. 1 shows the period distribution of F and FO Cepheids in the Galaxy, M31, and in the MCs. Data plotted in this figure clearly show that the mean metallicity affects the period distribution. In fact, the peak decreases from $`\mathrm{log}P=0.70.8`$ in the Galaxy and in M 31 ($`Z0.02`$) to $`\mathrm{log}P=0.2`$ in the SMC ($`Z0.004`$). The same outcome applies to FO periods, and indeed the period distribution peaks around $`\mathrm{log}P=0.3`$ (LMC) and around $`\mathrm{log}P=0.00.1`$ (SMC). The sample of Galactic FOs is too small to draw any firm conclusion. The main drawback that limits the use of FOs is that they present short periods (see Fig. 1), and, therefore, they are on average 1-2 mag fainter than F variables (see Fig. 2).

The detection of a sizable sample of FOs in nearby galaxies is important to improve the accuracy of Cepheid distances as well as to improve our knowledge on He-burning phases and on the occurrence of blue loops. The lower limit in the period distribution of FOs depends on the minimum mass whose blue loop crosses the instability strip. Current evolutionary predictions (Bertelli, Bressan, & Chiosi 1985; Castellani, Chieffi, & Straniero 1990; Stothers & Chin 1993; Bono et al. 2000; Limongi et al. 2000) suggest that the occurrence and the extent in temperature of blue loops strongly depend on the chemical compositions, on the efficiency of the mass loss as well as on the input physics and on the physical assumptions adopted to handle mixing processes in the convective layers.

Up to now a detailed comparison between theory and observations has been performed for field stars in the MCs and for a few LMC clusters and in particular for NGC 1866. Although several theoretical and observational (HST) studies have been devoted to this cluster we still lack firm quantitative constraints on the size of the convective core around the Turn-Off region (Testa et al. 1999; Barmina et al. 2002; Walker et al. 2002). This cluster might be crucial to address some problems concerning the Main Sequence fitting and the Cepheid distance scale, since it hosts more than 20 Cepheids (Walker et al. 2001). However, the pulsation parameters are only available for half of them and we still lack homogeneous photometric data for bright and faint cluster stars.

## 3. RR Lyrae variables and low-mass stars

RR Lyrae stars are very good tracers of the old stellar component, and play a key role in the absolute age determination of Globular Clusters. Although, RR Lyrae stars present several undoubtful advantages, the difference between distance estimates based on different calibrations (Baade-Wesselink method, HB models) of the $`M_V`$ vs \[Fe/H\] relation is systematically larger than the empirical uncertainties, and cover a range of $`0.3`$ mag on the distance modulus. This indicates that either the reddening accuracy is poor or RR Lyrae distances are still affected by systematic uncertainties. Moreover, recent theoretical (Caputo et al. 2000) and empirical (Layden 2000) studies confirm that the $`M_V`$ vs \[Fe/H\] relation is not linear when moving from metal-poor to metal-rich RR Lyrae.

Some of the problems affecting the RR Lyrae distance scale can be overcome using the K-band Period-Luminosity ($`PL_K`$) relation. In a seminal empirical investigation Longmore et al. (1990) demonstrate that RR Lyrae do obey a $`PL_K`$ relation in this band. This finding was further strengthened by a recent theoretical investigation (Bono et al. 2001) suggesting that RR Lyrae obey to a very tight $`PLZ_K`$ relation connecting the period, the luminosity, the K-band absolute magnitude, and the metallicity. Current models also disclosed that the $`PLZ_K`$ relation is, in contrast with the $`M_V`$ vs \[Fe/H\] relation, marginally affected by off-ZAHB evolution and by a spread in stellar masses (see Fig. 3). The intrinsic accuracy of previous predictions was tested using the prototype RR Lyr. Recently Benedict et al. (2002) on the basis of new astrometric data collected with FGS 3@HST, provided an accurate estimate of the absolute trigonometric parallax of RR Lyr. The uncertainty of the new parallax ($`\pi _{abs}=3.82\pm 0.20`$ mas) is $`3`$ times smaller than the uncertainty on the Hipparcos parallax ($`\pi _{abs}=4.38\pm 0.59`$ mas). Interestingly, Bono et al. (2002) using the predicted $`PLZ_K`$ relation found a pulsation parallax $`\pi _{puls}=3.858+/0.131`$ mas which agrees quite well with the HST parallax and presents a smaller formal error.

It is noteworthy that HB stars are the crossroad of several astrophysical problems such as the second-parameter problem, the UV-upturn in ellipticals (Ferguson, this volume), as well as the dependence of $`M_V(RR)`$ on metallicity. Detailed extragalactic samples of HB and RR Lyrae stars are only available for MCs and a few LG dwarfs (Mateo 1998; Bersier & Wood 2002). New homogeneous data of HB stars in stellar systems with different chemical and dynamical histories might play a crucial role to constrain their evolution properties.

## 4. Final remarks

Some problems affecting the evolutionary and pulsation properties of classical Cepheids and RR Lyrae stars can be settled during the next ten years using the unprecedented sensitivity and spatial resolution of the CCD camera that are available (ACS) or will become available (WFPC3) on board of HST. Sizable samples of FO Cepheids might be detected and measured across the disk of M 31. The number of orbits necessary to accomplish this experiment is rather modest, since this galaxy hosts regions rich of Cepheids and also because the apparent magnitudes around the minimum luminosity, assuming a distance modulus of 24.5 and a mean reddening of E(B-V)=0.08, roughly range from $`V22`$ ($`\mathrm{log}P=0.6`$) to $`V24`$ ($`\mathrm{log}P=0.1`$). The same conclusions apply to RR Lyrae stars, and indeed their apparent magnitudes range from $`V25`$ to $`V26`$ (Clementini et al. 2001). Note that up to now we still lack complete sample of HB and RR Lyrae stars in elliptical galaxies such as M 32.

A substantial improvement concerning the intrinsic accuracy of the Cepheid distance scale will also be provided by NGST. The superior sensitivity in the NIR bands will allow us to provide accurate measurements of Cepheids mean K magnitudes up to the Virgo and the Fornax cluster. The new data will supply robust distance determinations based on both the $`PL_K`$ and the $`PLC(V,VK)`$ relations. At the same time, NGST will supply the unprecedent opportunity to estimate the distance of several galaxies up to the outskirts of the LG using simultaneously NIR relations for RR Lyrae ($`PLZ_K`$) and Classical Cepheids. Future astrometric missions such as the SIM, and GAIA may also supply accurate empirical constraints concerning the accuracy of absolute (zero-point) and relative (slope) distance estimates as well as on their dependence on metallicity.

The previous improvements would supply tight constraints on the systematic uncertainties affecting primary distance indicators, and hopefully a determination of the Hubble constant with a global accuracy of the order of 5-7%. The next step that can allow us to nail down the systematic uncertainty affecting the evaluation of $`H_0`$ is to only use Cepheids and to by pass secondary distance indicators. Data available in the literature clearly show that distance determinations to the Coma cluster based on different zero-points and secondary indicators range from $`34.64\pm 0.25`$ (SBF in the K-band, Jensen et al. 1999) to $`35.29\pm 0.11`$ (various methods, Tammann et al. 1999), while the $`H_0`$ values range from $`85\pm 10`$ to $`60\pm 6`$ $`kms^1Mpc^1`$.

The main advantage in using Classical Cepheids is that we know the physical mechanisms that drive the pulsation instability in these objects. On the other hand, the experiment is challenging from an observational point of view, since it is necessary to detect and measure Cepheids in giant spirals of the Coma cluster. This means that we should be able to perform accurate photometry down to $`V3132`$ mag with a spatial resolution that is at least a factor of five better than current HST capabilities. We also note that this limit magnitude would supply the unique opportunity to estimate the age of old stellar component in LG galaxies up to the outer fringes (NGC 3109, Antlia), and in turn to supply a robust lower limit to the age of the Universe. Even though we expect that these two parameters are tightly correlated we could agree with Baade’s statement: Although I am quite certain that I do not mistake pink elephants for pink mice an unprejudiced check is always reassuring (as quoted by Feast 2000).

It is a pleasure to thank N. Panagia for a detailed reading of a draft of this paper. I am indebted to my collaborators for helpful discussions and suggestions.

## References

Alibert, Y., Baraffe, I., Hauschildt, P., & Allard, F. 1999, A&A, 344, 551

Baade, W. 1944, ApJ, 100, 79

Baade, W. 1958, in Stellar Populations, ed. D.J.K. O’Connell (Amsterdam: North Holland), 3

Baade, W. & Swope, H. H. 1965, AJ, 70, 212

Barmina, R., Girardi, L., & Chiosi, C. 2002, A&A, 385, 847

Benedict, G. F., et al. 2002, AJ, 123, 473

Bersier, D., & Wood, P. R. 2002, AJ, 123, 840

Bertelli, G., Bressan, A., & Chiosi, C. 1985, A&A, 150, 33

Bono, G., Caputo, F., Castellani, V., & Marconi, M. 1999b, ApJ, 512, 711

Bono, G., Caputo, F., Castellani, V., Marconi, M., & Storm, J. 2001, MNRAS, 326, 1183

Bono, G., Caputo, F., Castellani, V., Marconi, M., & Storm, J. 2002, MNRAS, 332, L78

Bono, G., Castellani, V., Marconi, M., 2002, ApJ, 565, L83

Bono, G., Groenewegen, M. A. T., Marconi, M., Caputo, F. 2002, ApJ, 574, L33

Bono, G., et al. 2000, ApJ, 543, 955

Bono, G., Marconi, M., & Stellingwerf, R. F. 1999a, ApJS, 122, 167

Caldwell, N., & Schommer, R. 1988, in The extragalactic distance scale, (San Francisco: ASP), 77

Caputo, F., Castellani, V., Marconi, M., & Ripepi, V. 2000, MNRAS, 316, 819

Castellani, V., Chieffi, A., & Straniero, O. 1990, ApJS, 74, 463

Chiosi, C., Wood, P. R., & Capitanio, N. 1993, ApJS, 86, 541

Clementini, G., Federici, L., Corsi, C., Cacciari, C., Bellazzini, M., & Smith, H. A. 2001, 2001, ApJ, 559, L109

Feast, M. 1999, PASP, 111, 775

Feast, M. 2000, JHA, 31, 29

Freedman, W. L., HST Key Project 2001, ApJ, 553, 47

Groenewegen, M. A. T. 2000, A&A, 363, 901

Groenewegen, M. A. T., & Oudmaijer, R. D. 2000, A&A, 356, 849

Jensen, J. B., Tonry, J. L., & Luppino, G. A. 1999, ApJ, 510, 71

Kennicutt, R. C. Jr., HST Key Project 1998, ApJ, 498, 181

Layden, A. 2002, in IAU Colloq. 185, Radial and Nonradial Pulsations as Probes of Stellar Physics, ed. C. Aerts, T. R. Bedding & J. Christensen-Dalsgaard (San Francisco: ASP), 118

Limongi, M., Straniero, O., & Chieffi, A. 2000, ApJS, 129, 625

Longmore, A. J., Dixon, R., Skillen, I., Jameson, R. F., Fernley, J. A., 1990, MNRAS, 247, 684

Lub, J. 1977, A&AS, 29, 345

Macri, L. M. et al. 2001, ApJ, 559, 243

Madore, B. F., Welch, D. L., McAlary, C. W., & McLaren, R. A. 1987, ApJ, 320, 26

Mateo, M. 1998, ARA&A, 36, 435

Mochejska et al. 2001, AJ, 122, 2477

Pietrzynski, G., Gieren, W., Fouqué, P., & Pont, F. 2002, AJ, 123, 789

Preston, G. W. 1964, ARA&A, 2, 23

Saha, A., Sandage, A., Tammann, G. A., Dolphin, A. E., Christensen, J., Panagia, N., Macchetto, F. D. 2001, ApJ, 562, 314

Sandage, A., Bell, R. A., & Tripicco, M. J. 1999, ApJ, 522, 250

Sasselov, D. D., et al. 1997, A&A, 324, 471

Storm, J., Carney, B. W., Gieren, W. P., Fouqué, P., & Fry, A. M. 2000, in IAU Coll. 176, The Impact of Large-Scale Surveys on Pulsating Star Research, eds. L. Szabados & D. Kurtz (San Francisco: ASP), 145

Stothers, R. B., & Chin, C. W. 1993, ApJ, 412, 294

Tammann, G. A., Sandage, A., Reindl, B. 1999, in 19th Texas Symposium on Relativistic Astrophysics and Cosmology, astro-ph/9904360

Tanvir, N. R. 1999, in Post-Hipparcos Cosmic Candles, eds. A. Heck & F. Caputo (Dordrecht: Kluwer), 17

Testa, V. et al. 1999, AJ, 118, 2839

Walker, A. R., et al. 2002, to be submitted

Walker, A. R., Raimondo, G., Di Carlo, E., Brocato, E., Castellani, V., & Hill, V. 2001, ApJ, 560, L139