text:
# Spin 1/2 as propagation on a lattice with symmetries modulo gauge transformations

## 1 Introduction

The mutual consistency of quantum mechanics and special relativity has remained a nontrivial issue, particularly with respect to locality  and quantum measurements , but also in describing a free particle. Quantum mechanics, on the one-particle level, makes a fundamental distinction between the roles of space and time which, it seems, can be overcome only in the framework of quantum field theory. But even in quantum field theory, in the usual line of argument, relativity is something that needs to be enforced. Moreover, if it is enforced by postulating unitary representations of the Poincaré group, there appear theoretical possibilities  which are never observed in Nature, such as continuous spin111Continuous spin arises if, in terms of induced representations, the little group of a light-like four-momentum is represented non-trivially in all of its components. or tachyons.

In this paper I suggest an ahistorical route to relativistic quantum mechanics as represented by the Dirac equation. I derive the equation, and with it the Lorentz invariance, from seemingly “non”-relativistic quantum theory. In fact, Dirac himself came to the conclusion in the 1950s  that the Michelson-Morley experiment, in view of symmetries present in quantum but not classical mechanics, had been overinterpreted as a support of special relativity. As for a derivation of Lorentz invariance from a “mechanism”, there is a well-known precedent: Maxwell , in a balance of working hypothesis and actual belief , utilized mechanistic ideas of electromagnetic fields which did not enforce relativity but got it right automatically. More recently, in constructing cellular automata, Bialynicki-Birula  noted that an automaton simulating the Weyl equation would require only very general conditions: a two-component wave function, an evolution that is linear and unitary, and (a vague remnant of relativity) that a wave function constant in space be also constant in time.

By the technicalities used (not by the direction of argument) the present paper is based on a discretization of the Dirac equation devised by Susskind in 1977 . The issue was to compensate for the doubling of the degrees of freedom encountered in replacing a derivative with an antihermitian difference

$$\frac{\mathrm{\Delta }f}{\mathrm{\Delta }x}=\frac{f(x+a)f(xa)}{2a}$$

A zero difference function, for example, is not only obtained from $`f=\mathrm{const}`$ but also from an alternating constant $`(1)^{x/a}`$ on the lattice sites. Susskind showed that spinorial degrees of freedom can be consistently assigned to different sites on a 3-dimensional lattice, thus thinning out the degeneracy of energy levels on a given lattice by a factor of four. The discretized Dirac equation resulting in this way is, in case of zero mass,

$$\begin{array}{ccc}\hfill i\dot{\psi }(x,y,z,t)& =& i\left(\psi (x+1,y,z,t)\psi (x1,y,z,t)\right)\hfill \\ & +& i\left(\psi (x,y+1,z,t)\psi (x,y1,z,t)\right)(1)^x\hfill \\ & +& i\left(\psi (x,y,z+1,t)\psi (x,y,z1,t)\right)(1)^{x+y}\hfill \end{array}$$
(1)
where $`\psi (x,y,z,t)`$ is a one-component wave function. Thus spin 1/2, usually thought of as “internal” to a point particle, can be encoded in a spatial arrangement of hopping amplitudes for a particle without internal structure.

In fact, equation (1) is the unique consequence of basic, not typically relativistic assumptions on the propagation of quantum particles living on the sites (as opposed to links or plaquettes) of a cubic spatial lattice:

* Locality: immediate propagation to nearest neighbours only
* Lattice symmetries are realised modulo gauge transformations

Also, time evolution will be assumed to be linear and unitary. There will remain two kinematical options, one of which will be discarded because it is infinitely slower than the other. Remarkably, the slow option is the one that would realise lattice symmetries in a strict sense, without accompanying gauge transformations.

In Section 2 the assumptions are specified; they include a general equation for linear, unitary propagation as it was already proposed by this author . In Section 3, the assumed invariances of the equation of propagation are evaluated, and equation (1) is derived. As for introducing particle mass, it is pointed out that an alternative to the standard term due to Susskind exists which avoids species doubling on infinite lattices. In Section 4, I present my Conclusions. In the Appendix, some omissions and simplifications anticipated in Sections 2 and 3 are justified.

## 2 Specifying the assumptions

### 2.1 Nearest-neighbour hopping

We assume, as in , that a quantum particle, initially in a position eigenstate, will “move” by gradually (differentiably in $`t`$) forming superpositions of nearest-neighbour eigenstates. Then a general state vector, given as a superposition of eigenstates with coefficients $`\psi (\stackrel{}{s},t)`$, will evolve according to

$$i\dot{\psi }(\stackrel{}{s},t)=\underset{\stackrel{}{n}}{}\kappa (\stackrel{}{s},\stackrel{}{n})\psi (\stackrel{}{s}+\stackrel{}{n},t)$$
(2)
where the sum runs over nearest neighbours, and where $`\kappa (\stackrel{}{s},\stackrel{}{n})`$ are complex hopping amplitudes whose properties are to be determined. The sum also includes an on-site hopping amplitude represented by $`\stackrel{}{n}=0`$.

Since the distances between nearest neighbours are all the same, we assume that all hopping amplitudes are of the same magnitude. With a suitable rescaling of time we thus assume

$$|\kappa (\stackrel{}{s},\stackrel{}{n})|=1\text{for all }\stackrel{}{s}\text{ and all }\stackrel{}{n}0$$
(3)
Unitarity of time evolution and hermiticity of the Hamiltonian will be taken for granted. Using the standard scalar product of wave functions,

$$\psi |\phi =\underset{\stackrel{}{s}}{}\overline{\psi (\stackrel{}{s})}\phi (\stackrel{}{s})$$

the linear operator acting on the rhs of (2) is hermitian if and only if

$$\kappa (\stackrel{}{s},\stackrel{}{n})=\overline{\kappa (\stackrel{}{s}\stackrel{}{n},\stackrel{}{n})}$$
(4)
### 2.2 Invariances modulo gauge transformations

For a free particle, the equation of motion should be “the same” at all times and locations, as well as after a rotation. Quantum mechanically, the arbitrariness of the phases of position eigenstates allows to interpret “the same” as “gauge equivalent”.

In a local gauge transformation, the wave function at each space-time point is multiplied by a phase factor. Thus

$$\psi (\stackrel{}{s},t)_{\mathrm{old}}=g(\stackrel{}{s},t)\psi (\stackrel{}{s},t)_{\mathrm{new}}|g(\stackrel{}{s},t)|=1$$

In terms of the new wave function, equation (2) involves the hopping amplitudes

$`\kappa (\stackrel{}{s},\stackrel{}{n})_{\mathrm{new}}`$ $`=`$ $`g(\stackrel{}{s}+\stackrel{}{n},t)\kappa (\stackrel{}{s},\stackrel{}{n})_{\mathrm{old}}g(\stackrel{}{s},t)^1\stackrel{}{n}0`$ (5)
$`\kappa (\stackrel{}{s},0)_{\mathrm{new}}`$ $`=`$ $`\kappa (\stackrel{}{s},0)_{\mathrm{old}}i\dot{g}(\stackrel{}{s},t)g(\stackrel{}{s},t)^1`$ (6)
If $`S`$ is a symmetry operation on the lattice (translation or rotation) the hopping amplitudes of equation (2) would in general change according to

$$\kappa (\stackrel{}{s},\stackrel{}{n})_{\mathrm{new}}=\kappa (S^1\stackrel{}{s},S^1\stackrel{}{n})_{\mathrm{old}}$$

Our assumption is that $`\kappa _{\mathrm{new}}`$ is a local gauge transform of $`\kappa _{\mathrm{old}}`$:

$$\kappa (\stackrel{}{s},\stackrel{}{n})_{\mathrm{new}}=g(\stackrel{}{s}+\stackrel{}{n})\kappa (\stackrel{}{s},\stackrel{}{n})_{\mathrm{old}}g(\stackrel{}{s})^1$$

Expressing this entirely in terms of $`\kappa _{\mathrm{old}}`$, and dropping the index, we have

$$g(\stackrel{}{s}+\stackrel{}{n})\kappa (\stackrel{}{s},\stackrel{}{n})g(\stackrel{}{s})^1=\kappa (S^1\stackrel{}{s},S^1\stackrel{}{n})$$
(7)
Since $`g`$ depends on $`S`$ we eventually write $`g(\stackrel{}{s},S)`$.

### 2.3 Maximal gauge fixing

Working in a particular gauge will greatly faciliate the evaluation of symmetries up to gauge transformations. Following the procedure of maximal gauge fixing as devised in  for Hamiltonian lattice gauge theories, let us choose to have

$$\kappa (x,y,z,\widehat{1})=1\text{ for all }x,y,z$$
(8)
This is accomplished, using (3), by a gauge transformation with a suitable behaviour in the $`\widehat{1}`$ direction:

$$g(x+1,y,z,t)=g(x,y,z,t)\kappa (x,y,z,\widehat{1})_{\mathrm{old}}^1$$

The values of $`g`$ on a plane with a constant $`x`$ coordinate are still free, and can be used to fix $`\kappa (\stackrel{}{s},\widehat{2})`$ on that plane. Let us choose to have

$$\kappa (0,y,z,\widehat{2})=1\text{ for all }y,z$$
(9)
which requires

$$g(0,y+1,z,t)=g(0,y,z,t)\kappa (0,y,z,\widehat{2})_{\mathrm{old}}^1$$

Finally, the values of $`g`$ along the line $`x=y=0`$ can be chosen so that

$$\kappa (0,0,z,\widehat{3})=1$$
(10)
Any further gauge transformation $`g`$ that is not constant throughout the lattice will destroy at least one of the conditions (8)-(10).

## 3 Propagation on a simple cubic lattice

### 3.1 Symmetries used

For a translation by a vector $`\stackrel{}{a}`$ we have

$$S^1\stackrel{}{n}=\stackrel{}{n}S^1\stackrel{}{s}=\stackrel{}{s}\stackrel{}{a}$$
(11)
For a rotation by $`90^{}`$ about the $`\widehat{1}`$ axis,

$$\begin{array}{ccc}\hfill S^1\widehat{1}& =& \hfill \widehat{1}\\ \hfill S^1\widehat{2}& =& \hfill \widehat{3}\\ \hfill S^1\widehat{3}& =& \hfill \widehat{2}\end{array}S^1(x,y,z)=(x,z,y)$$
(12)
For a rotation by $`90^{}`$ about the $`\widehat{3}`$ axis,

$$\begin{array}{ccc}\hfill S^1\widehat{1}& =& \hfill \widehat{2}\\ \hfill S^1\widehat{2}& =& \hfill \widehat{1}\\ \hfill S^1\widehat{3}& =& \hfill \widehat{3}\end{array}S^1(x,y,z)=(y,x,z)$$
(13)
### 3.2 Determining the hopping amplitudes

#### 3.2.1 Evaluating translations

Let $`S`$ in equation (7) be a translation as specified in (11). Thus

$$g(\stackrel{}{s}+\stackrel{}{n},\stackrel{}{a})\kappa (\stackrel{}{s},\stackrel{}{n})g(\stackrel{}{s},\stackrel{}{a})^1=\kappa (\stackrel{}{s}\stackrel{}{a},\stackrel{}{n})\stackrel{}{n}0$$
(14)
Putting $`\stackrel{}{n}=\widehat{1}`$ in (14) and using (8) we see that $`g(\stackrel{}{s},\stackrel{}{a})`$ must be independent of the coordinate $`x`$,

$$g(x,y,z,\stackrel{}{a})=g(y,z,\stackrel{}{a})\text{for all }\stackrel{}{a}$$

Now putting $`\stackrel{}{n}=\widehat{2}`$ and $`\stackrel{}{a}=\widehat{1}`$ in (14), we find

$$\kappa (x1,y,z,\widehat{2})=\kappa (x,y,z,\widehat{2})\left(g(y+1,z,\widehat{1})g(y,z,\widehat{1})^1\right)$$

Solving the recursion in $`x`$ and using gauge condition (9), we have

$$\kappa (x,y,z,\widehat{2})=e^{ix\alpha (y,z)}\text{ where }e^{i\alpha (y,z)}=g(y,z,\widehat{1})g(y+1,z,\widehat{1})^1$$

In fact, $`\alpha `$ must be independent of $`y`$ and $`z`$ since by re-inserting the last equation in (14) (with $`\stackrel{}{n}=\widehat{2}`$) and considering $`\stackrel{}{a}=\widehat{2},\widehat{3}`$ we encounter an $`x`$ dependence of $`\mathrm{exp}(ix\alpha (y,z))`$ on the lhs and $`\mathrm{exp}(ix\alpha (y1,z))`$ or $`\mathrm{exp}(ix\alpha (y,z1))`$, respectively, on the rhs. Hence, $`\alpha (y1,z)=\alpha (y,z1)=\alpha (y,z)`$ so that

$$\kappa (x,y,z,\widehat{2})=e^{i\alpha x}$$
(15)
Using (15), equation (14) with $`\stackrel{}{n}=\widehat{2}`$, $`\stackrel{}{a}=\widehat{1},\widehat{2}`$ can now be read as a recursion relation determining the gauge transformations $`g(\stackrel{}{s},\widehat{1})`$ and $`g(\stackrel{}{s},\widehat{2})`$ up to their values at $`x=y=0`$. We find

$`g(y,z,\widehat{1})`$ $`=`$ $`g(0,z,\widehat{1})e^{i\alpha y}`$ (16)
$`g(y,z,\widehat{2})`$ $`=`$ $`g(0,z,\widehat{2})`$ (17)
To obtain restrictions on the hopping amplitude in the $`\widehat{3}`$ direction, we now insert (16) and (17) into (14), using $`\stackrel{}{n}=\widehat{3}`$ and $`\stackrel{}{a}=\widehat{1},\widehat{2}`$. Thus

$`g(0,z+1,\widehat{1})\kappa (x,y,z,\widehat{3})g(0,z,\widehat{1})^1`$ $`=`$ $`\kappa (x1,y,z,\widehat{3})`$
$`g(0,z+1,\widehat{2})\kappa (x,y,z,\widehat{3})g(0,z,\widehat{2})^1`$ $`=`$ $`\kappa (x,y1,z,\widehat{3})`$
Taking into account the gauge condition (10) the recursions are readily resolved, yielding

$$\kappa (x,y,z,\widehat{3})=e^{ix\beta (z)}e^{iy\gamma (z)}$$
(18)
where

$`e^{i\beta (z)}`$ $`=`$ $`g(0,z+1,\widehat{1})^1g(0,z,\widehat{1})`$
$`e^{i\gamma (z)}`$ $`=`$ $`g(0,z+1,\widehat{2})^1g(0,z,\widehat{2})`$
#### 3.2.2 Evaluating unitarity

Applying (4) to (8), (15), (18) we obtain

$`\kappa (x,y,z,\widehat{1})`$ $`=`$ $`1`$
$`\kappa (x,y,z,\widehat{2})`$ $`=`$ $`e^{i\alpha x}`$ (19)
$`\kappa (x,y,z,\widehat{3})`$ $`=`$ $`e^{ix\beta (z1)}e^{iy\gamma (z1)}`$
#### 3.2.3 Evaluating $`\mathrm{𝟗𝟎}^{\mathbf{}}`$ rotations about the $`𝒙`$ axis

Let $`S`$ in equation (7) be the rotation specified by (12). Putting $`\stackrel{}{n}=\widehat{1},\widehat{2},\widehat{3}`$,

$`g(x+1,y,z)\kappa (x,y,z,\widehat{1})g(x,y,z)^1`$ $`=`$ $`\kappa (x,z,y,\widehat{1})`$
$`g(x,y+1,z)\kappa (x,y,z,\widehat{2})g(x,y,z)^1`$ $`=`$ $`\kappa (x,z,y,\widehat{3})`$ (20)
$`g(x,y,z+1)\kappa (x,y,z,\widehat{3})g(x,y,z)^1`$ $`=`$ $`\kappa (x,z,y,\widehat{2})`$
whence, using (8),(15), (18) and (19),

$`g(x+1,y,z)g(x,y,z)^1`$ $`=`$ $`1`$
$`g(x,y+1,z)e^{i\alpha x}g(x,y,z)^1`$ $`=`$ $`e^{ix\beta (y1)}e^{iz\gamma (y1)}`$
$`g(x,y,z+1)e^{ix\beta (z)}e^{iy\gamma (z)}g(x,y,z)^1`$ $`=`$ $`e^{i\alpha x}`$
By the first of these equations, $`g`$ must not depend on $`x`$. Thus, in the third and second equation, the only $`x`$ dependence occurs in the exponentials, implying

$$e^{i\beta (z)}=e^{i\alpha }=e^{i\beta (y1)}\text{for all }y,z$$

Thus $`\beta =\mathrm{const}`$, and there remain two possibilities,

$$e^{i\alpha }=e^{i\beta }=\pm 1$$
(21)
#### 3.2.4 Evaluating $`\mathrm{𝟗𝟎}^{\mathbf{}}`$ rotations about the $`𝒛`$ axis

Now let $`S`$ in equation (7) be the rotation specified by (13). Putting $`\stackrel{}{n}=\widehat{1},\widehat{2}`$,

$`g(x+1,y,z)\kappa (x,y,z,\widehat{1})g(x,y,z)^1`$ $`=`$ $`\kappa (y,x,z,\widehat{2})`$
$`g(x,y+1,z)\kappa (x,y,z,\widehat{2})g(x,y,z)^1`$ $`=`$ $`\kappa (y,x,z,\widehat{1})`$
Using (8),(15),(19) we obtain

$`g(x+1,y,z)g(x,y,z)^1`$ $`=`$ $`e^{i\alpha y}`$
$`g(x,y+1,z)g(x,y,z)^1`$ $`=`$ $`e^{i\alpha x}`$
The solution to these recursion relations is

$$g(x,y,z)=e^{i\alpha xy}g(0,0,z)$$

The $`xy`$ dependent factor drops out when inserted in (7) with $`\stackrel{}{n}=\widehat{3}`$, leaving

$$g(0,0,z+1)\kappa (x,y,z,\widehat{3})g(0,0,z)^1=\kappa (y,x,z,\widehat{3})$$

Using (18) we arrive at

$$g(0,0,z+1)e^{i\alpha x}e^{iy\gamma (z)}g(0,0,z)^1=e^{i\alpha y}e^{ix\gamma (z)}$$

Considering the $`x`$ and $`y`$ dependences we obtain

$$e^{i\gamma }=e^{i\alpha }$$
(22)
Thus, using (21), the hopping amplitudes are determined up to the choice of $`\alpha =0`$,

$$\kappa (x,y,z,\widehat{n})1$$
(23)
or $`\alpha =\pi `$,

$$\begin{array}{ccc}\hfill \kappa (x,y,z,\widehat{1})=& \kappa (x,y,z,\widehat{1})& =1\hfill \\ \hfill \kappa (x,y,z,\widehat{2})=& \kappa (x,y,z,\widehat{2})& =(1)^x\hfill \\ \hfill \kappa (x,y,z,\widehat{3})=& \kappa (x,y,z,\widehat{3})& =(1)^{x+y}\hfill \end{array}$$
(24)
### 3.3 Staticity of the scalar solution

Option (23) for the hopping amplitudes would also result from postulating strict invariance under the lattice symmetries, as already studied in . Its continuum limit was found to be the nonrelativistic Schrödinger equation of a scalar particle. However, the time scale222$`\kappa `$ is an inverse time by equation (2). The lattice spacing $`a`$ emerges from Taylor expansions of next-neighbour terms; hence, its dimension is always cancelled by that of a spatial derivative. The derivatives act on wave functions in the continuum limit, so they are independent of $`a`$. on which the wave functions would evolve was found to be $`\lambda ^2/\kappa a^2`$, where $`\kappa `$ denotes the nearest-neighbour hopping amplitude, and $`\lambda `$ is a length scale of the wave function. $`\lambda ^2/\kappa a^2`$ is also the time scale of unitary cellular automata simulating scalar particles . In contrast, the time scale of option (24) can be seen to be $`\lambda /\kappa a`$ from the initial choice of scale in (3) and its modification in (27); again, $`\kappa `$ denotes the nearest-neighbour hopping amplitude. If $`a`$ is very small (like the Planck length), we obviously have

$$(\lambda /a)^2\lambda /a$$

so option (23) tends to a static (non-kinetic and, in this sense, non-particle) scenario relative to (24).

### 3.4 Recovering the massless Dirac equation

This section reviews standard procedure with lattice fermions. Using (24) the hopping equation (2) reads

$$\begin{array}{ccc}\hfill i\dot{\psi }(x,y,z,t)& =& \left(\psi (x+1,y,z,t)+\psi (x1,y,z,t)\right)\hfill \\ & +& \left(\psi (x,y+1,z,t)+\psi (x,y1,z,t)\right)(1)^x\hfill \\ & +& \left(\psi (x,y,z+1,t)+\psi (x,y,z1,t)\right)(1)^{x+y}\hfill \end{array}$$

Equation (1) is recovered by the gauge transformation $`\psi _{\mathrm{old}}=i^{x+y+z}\psi _{\mathrm{new}}`$. For the new $`\psi `$,

$$\begin{array}{ccc}\hfill i\dot{\psi }(x,y,z,t)& =& i\left(\psi (x+1,y,z,t)\psi (x1,y,z,t)\right)\hfill \\ & +& i\left(\psi (x,y+1,z,t)\psi (x,y1,z,t)\right)(1)^x\hfill \\ & +& i\left(\psi (x,y,z+1,t)\psi (x,y,z1,t)\right)(1)^{x+y}\hfill \end{array}$$

Since the alternating sign factors are strongly fluctuating when viewed on a length scale much larger than the lattice spacing $`a`$, there can be no smooth solution to equation (1). However, the equation is solved by a superposition of wave functions of the form

$$\begin{array}{cc}\psi _{00}(x,y,z,t)\hfill & \psi _{01}(x,y,z,t)(1)^y\hfill \\ \psi _{10}(x,y,z,t)(1)^x\hfill & \psi _{11}(x,y,z,t)(1)^{y+x}\hfill \end{array}$$
(25)
where $`\psi _{AB}(x,y,z,t)`$ is assumed to be smooth in the sense that it varies from a lattice site to the next in $`𝒪(a)`$ at most. As suggested by the double index, the space of solutions is a tensor product. Multiplication by $`(1)^x`$, for example, interchanges the presence/absence of that factor in the wave function, hence it is represented by the Pauli matrix $`\sigma _1`$ acting on the first index, and by $`\sigma _1\mathrm{𝟏}`$ acting on both indices. Similarly, differentiating along the $`x`$ direction gives an extra minus sign depending on whether the factor $`(1)^x`$ is present or absent; this corresponds to the action of $`\sigma _3\mathrm{𝟏}`$. Thus the right-hand side of (1) combines matrix factors and spatial differences into

$$\sigma _3\mathrm{𝟏}i\mathrm{\Delta }_x+\sigma _1\sigma _3i\mathrm{\Delta }_y+\sigma _1\sigma _1i\mathrm{\Delta }_z$$
(26)
The tensor products are readily seen to satisfy the algebraic relations of the Dirac $`\alpha `$ matrices. The difference operations asymptotically tend to $`2a/x`$, $`2a/y`$, $`2a/z`$ in the continuum limit $`a0`$. We may absorb the factor of $`2a`$ in a redefinition of the time parameter, thus recovering the massless Dirac equation

$$i\frac{\psi }{t^{}}=i\alpha _k\frac{\psi }{x_k}t^{}=2at$$
(27)
### 3.5 Mass terms

In spatial continuum, massless Dirac particles have chiral symmetry. This symmetry gets broken if a mass term is introduced. Therefore it is not unsatisfactory to find that mass terms on a lattice may require the breaking of a lattice symmetry. Susskind’s mass term is an on-site hopping amplitude $`\mu (1)^{x+y+z}`$ so that equation (1) becomes

$$i\dot{\psi }(x,y,z,t)=\left(i\mathrm{\Delta }_x+i(1)^x\mathrm{\Delta }_y+i(1)^{x+y}\mathrm{\Delta }_z+\mu (1)^{x+y+z}\right)\psi (x,y,z,t)$$

The last term breaks the invariance (modulo gauge transformations) under translations by one lattice unit, while invariance under translations by two units is preserved. Susskind’s mass term requires a doubling of the dimension of the space of solutions, since the functions (25) need to be complemented by analogous functions with an extra factor of $`(1)^z`$. This expands the tensor products (26) of the Hamiltonian to

$$\sigma _3\mathrm{𝟏}\mathrm{𝟏}i\mathrm{\Delta }_x+\sigma _1\sigma _3\mathrm{𝟏}i\mathrm{\Delta }_y+\sigma _1\sigma _1\sigma _3i\mathrm{\Delta }_z+\mu \sigma _1\sigma _1\sigma _1$$

The enlarged space of states is also recovered by acting on (25) with a symmetry of equation (1), Susskind’s version of the parity operation

$$\psi (x,y,z,t)(1)^{x+y+z}\psi (x,y,z,t)$$

In Monte Carlo simulations, which can only use lattices with a finite number of sites, there is no natural distinction between smooth and strongly fluctuating wave functions. Thus the degeneracy related to the above parity operation can only be suppressed at the expense of some arbitrariness.

On infinite lattices, however, an extra factor of $`(1)^z`$ to the functions (25) does make a difference. It may therefore be of interest to note that Susskind’s term is not the only possibility of introducing mass. For example, we may allow for a variation of the magnitude of the hopping amplitude in the $`x`$ direction so as to violate (3) while keeping (4),

$$\kappa (x,y,z,\widehat{1})=i+i\mu (1)^x\kappa (x,y,z,\widehat{1})=i+i\mu (1)^x$$

Physically this would correspond to an alternating variation of the lattice spacing. The additional term in hopping equation (1) is

$$i\mu (1)^x\left(\psi (x+1,y,z,t)+\psi (x1,y,z,t)\right)$$

Repeating the arguments that lead to (26) we obtain the operator

$$\sigma _3\mathrm{𝟏}i\mathrm{\Delta }_x+\sigma _1\sigma _3i\mathrm{\Delta }_y+\sigma _1\sigma _1i\mathrm{\Delta }_z+2i\mu (\sigma _1\sigma _3)\mathrm{𝟏}$$

Since $`i\sigma _1\sigma _3=\sigma _2`$ and since $`\sigma _2\mathrm{𝟏}`$ anticommutes with the first three tensor products we recover the Dirac matrices in the form

$$\sigma _3\mathrm{𝟏}=\alpha _1\sigma _1\sigma _3=\alpha _2\sigma _1\sigma _1=\alpha _3\sigma _2\mathrm{𝟏}=\beta $$

## 4 Conclusions

We have derived Susskind’s discretization of the Dirac equation from assumptions which apparently do not anticipate special relativity. While time was assumed to run continuously, spatial coordinates were confined to a lattice (reminiscent of a stack of particle detectors). Locality, too, was imposed in an unrelativistic sense, assuming that propagation from some position will, within a short interval of time, reach the nearest neighbours only. The intrinsically quantum-mechanical assumption was that the amplitudes of propagation will respect the symmetries of the lattice to the extent they have to in quantum theory, namely up to phase shifts of position eigenstates.

On the basis of these assumptions, a non-relativistic and a relativistic option appeared at the same stage in section 3.3. As it happened, this was simultaneously the alternative between strict symmetry and symmetry modulo gauge transformations. It was a matter of kinematical speed, rather than principle, that the non-relativistic option was discarded.

What insight do we gain by this route to the Dirac equation? I think it explains the preferred role of spin 1/2 in the Standard Model, since no internal structure of a particle (of a kind living on lattice sites) was assumed, and yet the Dirac equation resulted. In particular, in the case of zero mass where continuous spin is a possibility consistent with Poincaré invariance, it was just the massless version of the Dirac equation which emerged. More generally, I think the unity rather than mere consistency of special relativity and quantum theory—even in a “non”-relativistic formulation of the latter—is emerging here. Finally, gauge transformations turn out to be as fundamental to the propagation of free particles as to particle interactions.

It would be interesting to determine the “internal” degrees of freedom on lattices with other than simple cubic structure, especially with some of the infinitely many close-pack structures .

## Appendix

### Time-dependence ruled out

To simplify the notation, hopping amplitudes were so far considered as functions of the spatial coordinates only. Here we show that in the gauge we were using there is, in fact, no other possibility consistent with time-translation invariance modulo gauge transformations.

Let $`g(\stackrel{}{s},t)`$ be the gauge transformation accomplishing the shift $`S`$ of the hopping amplitudes by a time $`\delta t`$. Equation (7) then reads

$$g(\stackrel{}{s}+\stackrel{}{n},t)\kappa (\stackrel{}{s},\stackrel{}{n},t)g(\stackrel{}{s},t)^1=\kappa (\stackrel{}{s},\stackrel{}{n},t\delta t)$$
(28)
Putting $`\stackrel{}{n}=\widehat{1}`$ and using our gauge condition (8), which holds at all times, we obtain $`g(x+1,y,z,t)=g(x,y,z,t)`$ and hence

$$g(x,y,z,t)=g(y,z,t)$$

Now putting $`\stackrel{}{n}=\widehat{2}`$ and $`x=0`$ in (28) and using gauge condition (9) we moreover obtain

$$g(y,z,t)=g(z,t)$$

Finally, putting $`\stackrel{}{n}=\widehat{3}`$ and $`x=y=0`$ and using gauge condition (10) we see that the gauge factor can only be a function of the time parameter:

$$g(z,t)=g(t)$$

But then $`g`$ has no effect at all in equation (28), and the nearest-neighbour hopping amplitudes must be strictly invariant under a time shift:

$$\kappa (\stackrel{}{s},\stackrel{}{n},t)=\kappa (\stackrel{}{s},\stackrel{}{n},t\delta t)$$

### On-site hopping gauged away

On-site hopping amplitudes consistent with the spatial symmetries (modulo gauge transformations) can always be gauged away, as we now show. Thus it was justified to omit them in the previous sections (excluding section 3.5 where translational symmetry was partially broken).

The analogue of equation (7), using gauge transformation (6), would be

$$\kappa (\stackrel{}{s},0,t)i\dot{g}(\stackrel{}{s},t)g(\stackrel{}{s},t)^1=\kappa (S^1\stackrel{}{s},0,t)$$
(29)
Let us reconsider the gauge factors which accomplished spatial translations of the nearest-neighbour amplitudes (section 3.2.1), taking into account the final expressions of the amplitudes as given by (24). Equation (14) with $`\stackrel{}{n}`$ put equal to $`\widehat{1},\widehat{2},\widehat{3}`$ can then be read as a recursion relation determining the $`x,y,z`$ dependence of $`g(x,y,z,\stackrel{}{a})`$ for a given shift vector $`\stackrel{}{a}`$. The arguments leading to those expressions did anticipate that $`\kappa (x,y,z,\widehat{n})`$ would be independent of time, but this was justified in the previous section.

Those recursion relations implicit in (14) determine $`g(x,y,z,\stackrel{}{a})`$, at each instant of time, up to a global phase factor. That factor could depend on $`t`$, allowing for an expression of the form

$$g(x,y,z,t,\stackrel{}{a})=g_0(x,y,z,\stackrel{}{a})h(t,\stackrel{}{a})$$

Inserting this in (29) the time-independent $`g_0`$ drops out, leaving

$$\kappa (\stackrel{}{s},0,t)i\dot{h}(t,\stackrel{}{a})h(t,\stackrel{}{a})^1=\kappa (\stackrel{}{s}\stackrel{}{a},0,t)$$

In particular, putting $`\stackrel{}{a}=\widehat{1},\widehat{2},\widehat{3}`$ and defining

$$\stackrel{}{c}=(c_1,c_2,c_3)c_n(t)=i\dot{h}(t,\widehat{n})h(t,\widehat{n})^1$$

we see that on-site hopping amplitudes can at most take the form

$$\kappa (\stackrel{}{s},0,t)=\kappa (\stackrel{}{0},0,t)+\stackrel{}{s}\stackrel{}{c}(t)$$
(30)
But $`\stackrel{}{c}(t)`$ must vanish due to rotational symmetry (modulo gauge transformations) by an argument similar to the above for translations. Reconsidering section 3.2.3 we find that equations (20) in conjunction with (24) determine the gauge factor up to a time-dependent global factor,

$$g(x,y,z,t,R_x)=g_0(x,y,z,R_x)h(t,R_x)$$

In equation (29), again, the time-independent $`g_0`$ drops out, leaving

$$\kappa (\stackrel{}{s},0,t)i\dot{h}(t,R_x)h(t,R_x)^1=\kappa (R_x^1\stackrel{}{s},0,t)$$

Taken at the origin $`\stackrel{}{s}=0`$ the equation implies that $`\dot{h}`$ must vanish. Hence, $`\kappa (\stackrel{}{s},0,t)`$ must be strictly invariant under the $`90^{}`$ rotation about the $`x`$ axis, which implies that $`\stackrel{}{c}(t)`$ can at most have an $`x`$ component. This latter possibility can finally be ruled out by reconsidering the $`90^{}`$ rotation about the $`z`$ axis as in section 3.2.4.

The remaining term of (30) can be removed by a gauge transformation only dependent on time, satisfying $`\dot{g}(t)=i\kappa (\stackrel{}{0},0,t)g(t)`$.


format: text

text:
# Evolution of binaries producing Type Ia supernovae, luminous supersoft X-ray sources, and recurrent novae

## 1. Introduction

Type Ia supernovae (SNe Ia) are one of the most luminous explosive events of stars. Recently, SNe Ia have been used as good distance indicators that provide a promising tool for determining cosmological parameters (Perlmutter et al. 1999; Riess et al. 1998) because of their almost uniform maximum luminosities. However, one of the most important unresolved problems on SNe Ia is that we do not know the exact progenitor systems of SNe Ia.

It is widely accepted that the exploding star itself is an accreting white dwarf (WD) in a binary (Nomoto 1982; Nomoto, Thielemann, & Yokoi 1984). However, the companion star (and thus the observed binary system) is not known. Several objects have ever been considered, which include merging double white dwarfs (Iben & Tutukov 1984; Webbink 1984), recurrent novae (Starrfield, Sparks, & Truran, 1985), symbiotic stars (Munari & Renzini 1992) etc (Livio 2000, for recent summary).

Among these models, the systems of double white dwarfs (the double degenerate model) are not theoretically supported; when two WDs merge, carbon ignites not at the center but near the surface layer of the C+O WD; so that the C+O WD is converted into an O+Ne+Mg WD (Saio & Nomoto 1998). If the total mass of the merger exceeds the Chandrasekhar mass, the merged WD does not explode as an SN Ia but collapses to a neutron star (Nomoto & Kondo 1991).

On the other hand, the single degenerate systems has become more promising, in which a C+O WD gains mass from a non-degenerate companion star and grows its mass up to the Chandrasekhar mass. Such promising progenitors include luminous supersoft X-ray sources (SSS), recurrent novae (RNe), and symbiotic stars. However, the most important problem in these single degenerate systems is that we did not know the mechanism how to grow the mass of the WDs. Instead, there are effective mass-losing mechanisms such as nova ejection and common envelope evolutions that get rid of the mass accumulated on the WD. Here, I review and summarize our recent developments of the mechanisms that enable to grow the mass of WDs in the single degenerate systems and give a unified picture of evolutions of SSSs, recurrent novae, and SNe Ia.

## 2. Accretion wind evolution

First, we briefly explain how the mass of WDs in the single degenerate systems grow toward the Chandrasekhar mass limit. The reason that the standard binary evolution theory fails to explain the evolutionary path to SNe Ia and recurrent novae lies mainly in the theoretical prediction of a common envelope formation and the ensuing spiral-in process, because these processes inhibit the growth of white dwarfs in binary systems. As discussed in many previous papers on binary evolution, it has been widely accepted that a hydrogen-rich envelope on a mass-accreting white dwarf expands to a red giant size when its mass accretion rate, $`\dot{M}_{\mathrm{acc}}`$, exceeds a critical limit, i.e.,

$$\dot{M}_{\mathrm{cr}}=8.3\times 10^7(\frac{M_{\mathrm{WD}}}{M_{}}0.4)M_{}\text{ yr}^1,$$
(1)
(see, e.g., Fig. 9 of Nomoto 1982) and easily forms a common envelope. Once a common envelope is formed, two stars begin to spiral-in each other as a result of viscous drag, thus producing a double degenerate system (Iben & Tutukov 1984; Webbink 1984).

However, Hachisu, Kato, & Nomoto (1996) found that white dwarfs begin to blow strong winds ($`v_{\mathrm{wind}}1000`$ km s-1 and $`\dot{M}_{\mathrm{wind}}>1\times 10^6M_{}`$ yr-1) when $`\dot{M}_{\mathrm{acc}}>\dot{M}_{\mathrm{cr}}`$, as illustrated in Fig. 1, thus preventing the binary from collapsing. The binary does not spiral-in but keeps its separation. Hydrogen burns steadily, and therefore the helium layer of the white dwarf can grow at a rate of $`\dot{M}_{\mathrm{He}}\dot{M}_{\mathrm{cr}}`$. The other transferred matter is blown in the wind ($`\dot{M}_{\mathrm{wind}}\dot{M}_{\mathrm{acc}}\dot{M}_{\mathrm{cr}}`$), where $`\dot{M}_{\mathrm{wind}}`$ is the wind mass loss rate. We call this accretion wind because it begins to blow when the accretion rate exceeds the critical rate.

Therefore, we have to revise Fig. 9 of Nomoto (1982), which has been used widely in binary evolution scenarios to show the basic processes of mass-accreting white dwarfs. The revision to Nomoto’s Fig. 9 is shown in Fig. 1. The most important difference of Fig. 1 from Nomoto’s Fig. 9 is that white dwarfs increase their masses in a much wider parameter region. In our new diagram, the status of a mass-accreting white dwarf is specified by the following three phases: (1) accretion wind phase ($`\dot{M}_{\mathrm{acc}}>\dot{M}_{\mathrm{cr}}`$); (2) steady shell-burning phase ($`\dot{M}_{\mathrm{std}}<\dot{M}_{\mathrm{acc}}<\dot{M}_{\mathrm{cr}}`$); and (3) intermittent shell-flash phase ($`\dot{M}_{\mathrm{acc}}<\dot{M}_{\mathrm{std}}`$), where $`\dot{M}_{\mathrm{std}}`$ means the lowest limit of the mass accretion rate for steady hydrogen shell-burning. The new growing region of white dwarfs, $`\dot{M}_{\mathrm{acc}}>\dot{M}_{\mathrm{std}}`$, is much wider than old Nomoto’s (1982) narrow strip, $`\dot{M}_{\mathrm{std}}<\dot{M}_{\mathrm{acc}}<\dot{M}_{\mathrm{cr}}`$.

### 2.1. RX J0513.9-6951 as an example of accretion wind phase

As an example of the accretion wind phase, we have modeled long-term optical light curves of the LMC supersoft X-ray source RX J0513.9-6951 based on our wind model of mass-accreting white dwarfs. The source shows a quasi-periodic variability, the period of which is $`170`$ d and the amplitude of which is $`1`$ mag (Alcock et al. 1996). The source has a relatively longer ($`130`$ d) optical high-state ($`V16.2`$) and a shorter ($`40`$ d) optical low-state ($`V17.0`$) with the transition time of several days. X-rays were observed only in the optical low-states. This indicates that a strong wind blows only in the optical high-state because X-rays are self-absorbed by the wind. To examine the idea of the strong wind, we have calculated envelope solutions of mass-accreting WDs and found that optically thick winds really occur when the mass accretion rate exceeds a critical rate of $`\dot{M}_{\mathrm{cr}}`$ for $`M_{\mathrm{WD}}>0.7M_{}`$ even in a low metallicity of $`Z=0.004`$ of LMC.

We attribute the $`100150`$ d bright phase to our strong wind phase with an expanded white dwarf photosphere, whereas the $`40`$ d less bright phase to our static (no-wind) phase with a shrinking envelope of the WD. The duration and the brightness of the optical high state is explained by the model of a WD mass, $`M_{\mathrm{WD}}=1.35M_{}`$, with an extending accretion disk edge over the Roche lobe size during the wind phase as illustrated in Fig. 2. Irradiation of the accretion disk by the WD photosphere plays an essential role to reproduce the optical light curve because a large part of the optical light comes from the accretion disk.

The switching mechanism from the optical high to low state is understood as follows: Strong winds from the WD themselves make the separation increase, because the high velocity wind do not carry away large angular momentum. The winds further strip off the very surface layer of the companion star. As a result, the mass transfer rate decreases and it finally stops the strong winds. After the wind stops, the surface of the companion recovers in a thermal timescale of the envelope and the mass transfer rate increases up and strong winds begin to blow again. Thus, the semi-periodic cycles repeat.

Light curve fittings indicate that the optical high state is initiated by a sudden increase in the mass transfer rate of $`\dot{M}_{\mathrm{acc}}1\times 10^5M_{}`$. Then, the envelope of the white dwarf expands to $`R_{\mathrm{ph}}0.1R_{}`$, not exceeding over the accretion disk size of $`2R_{}`$. Our light curve fittings are shown in Fig. 3. In our models, about 40% of the transferred matter can be accumulated on the white dwarf but the rest is blown off in the wind. This suggests that it takes $`1\times 10^5`$ yr for the white dwarf of $`1.35M_{}`$ to grow to $`1.38M_{}`$ and to explode as a Type Ia supernova.

### 2.2. Two evolutionary paths to RNe and SNe Ia

Based on the new mechanism of the accretion wind evolution as shown in Fig. 1, Hachisu, Kato, & Nomoto (1999a) and Hachisu et al. (1999b) have proposed two paths to SNe Ia (see also, Li & van den Heuvel 1997). In these two paths, WDs can accrete mass continuously from the companion and grow at a rate of $`\dot{M}_{\mathrm{cr}}`$ of eq. (1). Therefore, candidates for SN Ia progenitors are systems at the final stages on these paths; i.e., one is a supersoft X-ray source (SSS) consisting of a WD and a lobe-filling, more massive, main-sequence (MS) companion as shown in the late stage of Fig. 4 (the WD+MS systems; see also Li & van den Heuvel 1997) and the other is a symbiotic star consisting of a WD and a mass-losing red giant (RG) as shown in the late stage of Fig. 5 (the WD+RG systems, see also Hachisu et al. 1996). Both the systems contribute to main parts of the SN Ia birth rate (Hachisu et al. 1999a, 1999b).

Hachisu et al. (1999a, 1999b) have followed many binary evolutions until SN Ia explosions. The final outcome is summarized in Fig. 6 for the initial orbital period, $`P_0`$ (in units of days), and the initial donor mass, $`M_{d,0}`$, which is representing $`M_{\mathrm{RG},0}`$ (the initial mass of the red-giant component), or $`M_{\mathrm{MS},0}`$ (the initial mass of the slightly evolved main-sequence component). In this figure, the initial mass of the white dwarf is assumed to be $`M_{\mathrm{WD},0}=1.0M_{}`$. The right-hand region, longer orbital periods of $`P_0100`$—1000 d, represents our results for the WD+RG system. In the left-hand region, shorter orbital periods of $`P_00.3`$—10 d, are shown the results for the WD+MS system.

We plot the initial parameter regions of binaries (thin solid curves in Fig. 6) that finally produce an SN Ia, which are corresponding to stage (F) in Fig. 4 or stage (D) in Fig.5. Left diagram in Fig. 6 shows such two regions of the WD+MS systems (SSS channel) and the WD+RG systems (symbiotic channel). We have added in this figure the final parameter regions (thick solid curves) just before an SN Ia explosion occurs. Starting from the initial region encircled by the thin solid curves, binary systems evolve and can explode as an SN Ia in the regions enclosed by the thick solid curves. Hachisu & Kato (2001) showed that some recurrent novae lie on the boarder of their parameter region that can produce an SN Ia (see also Kato 2001 in this volume). Among six recurrent novae with known orbital periods, five fall in the final regions of SN Ia progenitors just before the explosion. Only T Pyx is far out of the regions of SNe Ia.

## 3. Case BB mass transfer

It has been shown that the ejecta of the recurrent nova U Sco are helium-rich as much as He/H$`2`$ by number (Williams et al. 1981). The accretion disk of U Sco is also suggested to be helium-rich (Hanes 1985). Many of the SSSs have a stronger He II line than Balmer lines, thus also suggesting that the companion of the WD+MS systems (SSS channel) is helium-rich. In this section, therefore, we describe the reason that the cool component of the WD+MS systems has a helium-rich envelope.

Previous works consider only the initial systems consisting of a more massive AGB star with a C+O core and a less massive main-sequence star. This system undergoes a common envelope evolution and it finally yields a binary system of a mass-accreting C+O WD and a lobe-filling MS/subgiant star. Thus, they do not include another important evolutionary path, in which a more massive component fills up its inner critical Roche lobe when it develops a helium core of $`0.82.0M_{}`$ in its red giant phase as shown in the early stage of Fig. 4.

We consider, for example, a close binary with the primary mass of $`M_{1,i}=7M_{}`$ and the secondary mass of $`M_{2,i}=2M_{}`$ with the initial separation of $`a_i80600R_{}`$. When the primary has evolved to a red giant developing a helium core of $`M_{1,\mathrm{He}}1.01.6M_{}`$, it fills up its inner critical Roche lobe. Then, the binary undergoes a common envelope (CE) evolution, which yields a much more compact close binary of a naked helium star of $`M_{1,\mathrm{He}}1.01.6M_{}`$ and a main-sequence star of $`M_2=2M_{}`$ with the final separation of $`a_{\mathrm{f},\mathrm{CE}}450R_{}`$. Here, we assume a relation of

$$\frac{a_{\mathrm{f},\mathrm{CE}}}{a_i}\alpha _{\mathrm{CE}}\left(\frac{M_{1,\mathrm{core}}}{M_{1,i}}\right)\left(\frac{M_2}{M_{1,i}M_{1,\mathrm{core}}}\right),$$
(2)
with the efficiency factor of $`\alpha _{\mathrm{CE}}=1.0`$ for the common envelope evolutions.

After the common envelope evolution, the radii of the inner critical Roche lobes of the primary and the secondary become $`R_1^{}0.36a_{\mathrm{f},\mathrm{CE}}`$ and $`R_2^{}0.4a_{\mathrm{f},\mathrm{CE}}`$, respectively. Since the secondary radius should be smaller than its inner critical Roche lobe, i.e., $`1.7R_{}R_2<R_2^{}`$ ($`1.7R_{}`$ is the radius of a $`2M_{}`$ ZAMS star), the initial separation $`a_i`$ should exceeds $`80R_{}`$. The upper bound of the initial separation is obtained from the maximum radius of the $`7M_{}`$ star which has formed a helium core, i.e., $`a_i2R_{1,\mathrm{max}}<2\times 300R_{}`$. Thus, the allowable range of the initial separations is $`80R_{}<a_i<600R_{}`$ for a pair of $`M_{1,i}=7M_{}`$ and $`M_{2,i}=2M_{}`$. Then, we have $`R_1^{}1.518R_{}`$ and $`R_2^{}1.720R_{}`$ after the first common envelope evolution.

After the hydrogen-rich envelope is stripped away and hydrogen shell burning vanishes, the naked helium core contracts to ignite a central helium burning and becomes a helium main sequence star. When the helium star forms a C+O core, its helium envelope expands to fills its inner critical Roche lobe again. The helium is transferred stably to the secondary on an evolutionary time scale of $`\tau _{\mathrm{EV}}10^5`$ yr because the mass ratio is smaller than $`0.79`$ ($`q=M_1/M_2<0.79`$). The resultant mass transfer rate is $`\dot{M}_110^5M_{}`$ yr-1. No common envelope is formed for such a low rate. After the helium envelope is lost, the primary becomes a C+O WD of $`M_{\mathrm{C}+\mathrm{O}}0.91.1M_{}`$ and the separation increases by 10%$``$40%, i.e., $`a_{\mathrm{f},\mathrm{He}}(1.11.4)a_{\mathrm{f},\mathrm{CE}}(470)R_{}`$. Here, we assume the conservations of the total mass and angular momentum during the helium mass transfer to obtain the separation after the helium mass transfer, $`a_{\mathrm{f},\mathrm{He}}`$. The secondary receives almost pure helium matter of $`\mathrm{\Delta }M_{\mathrm{He}}0.10.5M_{}`$ to form a helium-enriched envelope. The secondary’s hydrogen content decreases to $`X0.6`$ if helium is completely mixed into the star. For the $`9+2.5M_{}`$ case, much more helium is transferred, but much less helium for the $`6+2M_{}`$ case.

## 4. Shrink of orbital separation during superwind phase

Based on the population synthesis analysis, Yungelson & Livio (1998) claimed that almost no realization frequency is derived for the original WD+RG model (Hachisu et al. 1996). In this section, we point out that very wide binaries with initial separations of $`a_i>1500R_{}`$, which were not included in their analysis, are essentially important in our SN Ia modeling.

The more massive component (the initial mass of $`M_{1,i}`$) of a binary first evolves to a red giant (AGB stage) and fills its inner critical Roche lobe. After a common envelope phase, the more massive component leaves a C+O WD, and the separation of the binary decreases by a factor of $`a_f/a_i1/401/50`$ from eq.(2) for $`M_{\mathrm{WD}}1M_{}`$ and $`M_21M_{}`$ because a $`1M_{}`$ WD descends from a main-sequence star of $`M_{1,i}78M_{}`$. Yungelson & Livio (1998) assumed that the separation of interacting binaries is smaller than $`a_i<1500R_{}`$. Then, the widest binaries have separations of $`a_f<3040R_{}`$ after common envelope evolution. The orbital period is $`P_0<20`$ days for $`M_{\mathrm{WD},0}1M_{}`$ and $`M_{\mathrm{RG},0}1M_{}`$. There is no SN Ia region of the WD+RG systems for $`P_0<20`$ days as seen from Fig. 6. Thus, Yungelson & Livio concluded that we cannot expect any SN Ia explosions from the right-hand SN Ia region (WD+RG system) in Fig. 6.

The reason that Yungelson & Livio’s modeling fails to reproduce our SN Ia progenitors is their assumption of $`a_i<1500R_{}`$. In what follows, we show that WD+RG binaries wide enough to have $`P_01001000`$ d are born from initially very wide binaries with $`a_i150040000R_{}`$.

A star with a zero-age main-sequence mass of $`M_{1,i}<8M_{}`$ ends its life by ejecting its envelope in a wind of relatively slow velocities ($`v1040`$ km s-1) as shown in Fig. 5. These wind velocities are as low as the orbital velocities of binaries with separations of $`a_i150040000R_{}`$ for a pair of $`M_{1,i}7M_{}`$ and $`M_{2,i}1M_{}`$. When the wind velocity is as low as the orbital velocity, we have the separation shrinking due to effective loss of the angular momentum by the outflowing gas, i.e., $`\dot{a}/a2\dot{M}_1/M_2<0`$ (see Hachisu et al. 1999a).

Once the binary system begins to shrink, its evolution becomes similar to a common envelope evolution (see stages B-D in Fig. 5). Thus, the separation is reduced by a factor of $`1/401/50`$, i.e., $`a_f301000R_{}`$ for $`M_{1,i}7M_{}`$ and $`M_{2,i}1M_{}`$. The orbital period becomes $`P_0153000`$ days for a pair of $`M_{\mathrm{WD},0}1M_{}`$ and $`M_21M_{}`$ (stage D in Fig. 5). These initial sets of the parameters are very consistent with the initial conditions of our WD+RG progenitor systems.

## 5. Type Ia supernova rates

We estimate the SN Ia rate from the WD+RG and WD+MS systems in our Galaxy by

$$\nu =0.2𝑑qd\mathrm{log}A\frac{dM}{M^{2.5}}\mathrm{yr}^1,$$
(3)
where $`q`$, $`A`$, and $`M`$ are the mass ratio ($`q=M_2/M_11`$), the initial separation in units of solar radius, and the primary initial mass in solar mass units (from eq.  of Iben & Tutukov 1984). The SN Ia rate is essentially the same as the volume of the initial phase space consisting of $`q`$, $`A`$, and $`M`$ as shown in Fig. 6.

The accretion winds of WDs open new two channels to SNe Ia. The case BB mass transfer in the SSS channel extends the phase volume towards the low mass ratio region as shown in Fig. 6. On the other hand, the orbital shrinkage during the superwind phase in the symbiotic channel enlarges the phase space towards the large initial separation also shown in Fig. 6. Therefore, the estimated rate of WD+RG/WD+MS systems becomes close to the observed rate in our Galaxy, i.e., $`\nu 0.003`$ yr-1, which comes from the SSS channel of $`\nu _{\mathrm{MS}}0.001`$ yr-1 and the symbiotic channel of $`\nu _{\mathrm{RG}}0.002`$ yr-1.

## References

Alcock, C. et al. 1996, MNRAS, 280, L49

Belczyński, K., & Mikolajewska, J. 1998, MNRAS, 296, 77

Hachisu, I., & Kato, M. 2001, ApJ, 558, 323

Hachisu, I., Kato, M., Kato, T, & Matsumoto, K. 2000, ApJ, 528, L97

Hachisu, I., Kato, M., & Nomoto, K. 1996, ApJ, 470, L97

Hachisu, I., Kato, M., & Nomoto, K. 1999a, ApJ, 522, 487

Hachisu, I., Kato, M., Nomoto, K., & Umeda, H. 1999b, ApJ, 519, 314

Hanes, D. A. 1985, MNRAS, 213, 443

Iben, I. Jr., & Tutukov, A. V. 1984, ApJS, 54, 335

Kato, M. 2001, in The Physics of Cataclysmic Variables and Related Objects (this volume), in press

Livio, M. 2000, in Type Ia Supernovae: Theory and Cosmology, ed. J. C. Niemeyer & J. W. Truran (Cambridge: Cambridge Univ. Press), 33

Li, X.-D., & van den Heuvel, E. P. J. 1997, A&A, 322, L9

Munari, U., & Renzini, A. 1992, ApJ, 397, L87

Nomoto, K. 1982, ApJ, 253, 798

Nomoto, K., & Kondo, Y. 1991, ApJ, 367, L19

Nomoto, K., Thielemann, F., & Yokoi, K. 1984, ApJ, 286, 644

Patterson, J. et al. 1998, PASP, 110, 380

Perlmutter, S. et al. 1999, ApJ, 1999, 517, 565

Riess, A. G. et al. 1998, AJ, 116, 1009

Saio, H., & Nomoto, K. 1998, ApJ, 500, 388

Starrfield, S., Sparks, W. M., & Truran, J. W. 1985, ApJ, 291, 136

Webbink, R. F. 1984, ApJ, 277, 355

Williams, R. E., Sparks, W. M., Gallagher, J. S., Ney, E. P., Starrfield, S. G., & Truran, J. W. 1981, ApJ, 251, 221

Yungelson, L., & Livio, M. 1998, ApJ, 497, 168


text:
# The apeNEXT Project

## 1 OVERVIEW

Computer power requirements for lattice gauge theory (LGT) simulations have increased exponentially in the last 20 years. Estimates of these requirements, in terms of computing power and memory needs for large scale simulations, such as studies of the hadronic spectrum with dynamical fermions close to the physical limit have been published recently  and extensively debated at this conference . It is increasingly clear that unrestricted access to computing resources of several sustained Tflops will be needed in the next few years. For about 15 years, a large fraction of the compute cycles used by LGT simulations have been provided by LGT-optimized massively parallel processors. However, presently available systems (see  for a recent review) are unable to reach the required performance target. The development of yet another generation of dedicated LGT engines is still one open option, even if different approaches to the problem are currently discussed .

This paper describes in detail one such project, apeNEXT, carried out jointly by DESY, INFN and the University of Paris-Sud (a similar project is described in ). This paper is structured as follows: the next section introduces requirements and guidelines of the apeNEXT architecture. Section 3 outlines the proposed structure of the new system, at the hardware level, while section 4 discusses software issues. Section 5 is a brief status report.

## 2 PROJECT REQUIREMENTS AND GUIDELINES

Our main goal is the development of a massively parallel system for “compute intensive” and “memory intensive” LGT simulations. The former case refers to dynamical fermions on not too large lattices (typical sizes of $`48^3\times 96`$ corresponding to system of $`L=2\mathrm{}4fm`$ and $`a=0.1\mathrm{}0.05fm`$) and small dynamic quark masses, while the latter case corresponds to quenched simulations on very large lattices ($`100^3\times 100÷200`$) and large $`\beta `$ ($`L=1.5\mathrm{}2.0fm`$ and $`a=0.1\mathrm{}0.02fm`$) relevant for $`b`$ physics with limited extrapolation in the mass of the heavy quark.

These physical requirements translate into the capability to handle online data structures of the order of several hundred Gbyte to 1 Tbyte, allowing the use of optimized algorithms that trade memory for performance, and to transfer from disk to memory and viceversa the tens of Tbyte corresponding to the full set of propagators needed for weak interaction studies (see  for an accurate estimate of these figures). In terms of computing power, we aim for systems delivering several TFlops.

We try to minimize architectural changes with respect to APEmille and to leverage on technology improvements, but we introduce improvements in some key architectural areas:

* Floating point arithmetics is upgraded to double precision.
* Typical LGT kernels have a ratio $`R4\mathrm{}6`$ between the number of floating point operations and the required operands. APEmille was able to provide enough memory bandwidth without any intermediate storage hierarchy between registers and main memory. In apeNEXT longer memory latencies (in units of processor clock cycles) reduce effective bandwidths for data bursts of the size of gauge matrices or spin-color fields. We introduce additional features to handle this problem.
* APEmille is not able to overlap remote data transfer with computation, so the system must idle waiting for data from remote nodes to arrive. This loss of performance can be partially avoided, by allowing concurrent link activity and computations. We also support concurrent data transfers on more than one data link.

## 3 THE apeNEXT ARCHITECTURE

The apeNEXT architecture is based on a three-dimensional grid of processing nodes, connected by data links between first-neigbours with periodic boundary conditions. Each node is a completely independent unit, with a processor and its memory bank. It executes independently its own copy of the same program, so the whole system is a SPMD processor. Data transfers between nodes follow the message-passing paradigm: the program on the destination node must explicitely receive data sent by the source node. This is done however with small latency ($`100ns`$).

apeNEXT is based on one building block, namely the J&T processor. It controls program flow, performs arithmetic operations and handles data links to neighbour nodes. A block diagram is shown in figure 1.

J&T has several similarities with previous APE processors:

* An interface to main memory, through a data bus of 128 bits, plus error detection and correction bits. The memory interface supports Double Data Rate (DDR) Synchronous Dynamic Rams (SDRAM), that will be widely used in standard PCs in the next two years. A set of 8 plus 1 memory chips provides 256 Mbytes of main (program and data) memory for each node.
* a large Register File (RF), holding 256 registers of 64+64 bits. There are three read ports, one write port and one bi-directional ports. Three read ports feed data to the arithmetic box, results enter via the write port, and the bi-dirirectional port exchanges data with memory and the queue system (see later).
* The arithmetic box performs the “normal” floating-point operation $`a\times b+c`$ on complex (or pair of real) numbers, corresponding to 8 (or 4) standard floating point operations per clock cycle. We obtain 1.6 Gflops peak performance with a 200 MHz clock. Pipeline length is 10 clock cycles. The arithmetic box also performs arithmetic and logic operations on pairs of 64 bit integer values and conversions from/to floating point format.

The new features of J&T include:

* A software-controlled program cache, where heavily used compute kernels are pre-loaded.
* An address-generation unit, which operates independently of the arithmetic unit.
* A serial control interface for initialization, exception handling and debugging.
* A queue system for pre-fetching local and remote data.
* Seven data links with a bandwidth of 200 Mbytes/s each.

Here we discuss in some details the memory-queue system and the data links.

The apeNEXT memory interface has a peak bandwidth of 3.2 Gbyte/sec (one complex data word per clock cycle, able to sustain algorithms with $`R4`$). In practice, a startup latency of about 12 cycles affects sustained performances. We plan to solve this problem through a program controlled prefetch mechanisms. During execution of iteration $`i`$ of a critical kernel, data structures needed for iteration $`(i+1)`$ are moved from memory to an intermediate storage element, a data queue, close (in terms of bandwidth and latency) to the register file. The size of the queue is 1024 complex words. The user program starts data transfers from memory to the queue. Memory accesses are local (from local memory to local queue) or remote (from memory to the queue in a neighbour processor) and hence may have different transfer times. The queue system, however, ensures delivery of the data to the receiving register file in the same order in which memory accesses were scheduled.

Prefetch reduces substantially remote bandwith requirements. Consider the evaluation of the Wilson-Dirac operator, for a lattice point on the border between two processing nodes. A total of 12 complex data words (fermionic degrees of freedom) plus possibly 9 complex data words (gauge fields) must be moved on a data link in the time interval used for the evaluation of the operator, corresponding to $`320`$ complex normal operations. If we (unrealistically) assume 100% efficiency, we require a transfer rate of about 1 byte per clock cycle, i.e. 200 Mbytes/sec. Sites at the edges (corners) of the physical region mapped onto each processor need two (three) data transfers. If the links in different directions operate concurrently bandwidth requirements are not increased.

The apeNEXT links are designed for the above requirements. Each link delivers a block of data words to the destination queues. The transfer, proceeeds concurrently with local node operations. An interlock mechanisms stallss the receiving node if the latter tries to pop data not yet arrived at the destination queue.

Data communication provides a synchronization mechanism: the user program starts a data transfer to a remote node over one link. At the same time it enables data reception from another link. If the sender is ahead in time with respect to the receiver, the former will be forced to wait till synchronization is achieved. In most cases relevant for LGT simulations, data sends and receive are in opposite directions, but more general communication patterns, such as systolic paths, can be established.

Each node has 7 data links. Six links make up the regular three-dimensional network, while the “seventh” link is available for I/O to the host system (see later, for more details). Each link uses Cyclic Redundancy Check (CRC) to identify errors. Each data block is divided in packets of 128 bits (plus 16 CRC bits). Corrupted packets are re-trasmitted.

apeNEXT machines use Processing Boards (PB) as the basic hardware building block. Each PB has 16 processing elements (an array of $`4\times 2\times 2`$ nodes). A set of 16 PBs is assembled in a Crate, a system of $`4\times 8\times 8`$ nodes with communication links on the backplane. Crates are connected by cable links, building larger systems of $`4n\times 8\times 8`$ nodes. Each node is also connected to the slow control interface. Each PB can be connected to a host computer across the “seventh link”, which handles heavy data I/O. The number of hosts and of PBs using the “seventh link” is arbitrary. Therefore the I/O bandwidth of apeNEXT can be scaled according to the specific applications.

## 4 SOFTWARE

apeNEXT programs are written in the TAO programming language, also used in earlier APE machines. Our goal is almost complete compatibility, in the sense that old APEmille program will only need to be recompiled. Optimal efficiency on apeNEXT will of course involve some machine specific tuning of the codes. In particular programs will have to handle the queues, which are accessible at TAO level.

We are also developing a C compiler, based on open source compilers. Our target is to treat TAO and C on almost equal footing, also from the point of view of efficiency. From first experiences with a prototype version we hope that C codes reach performances comparable with TAO.

At the backend level of the compilation chain, we plan to use the instruction scheduling techniques used in APEmille to optimize pipeline usage. We are also in the process of adopting further optimization tools at the assembly level.

## 5 PROJECT STATUS

All apeNEXT hardware building blocks are in an advanced stage of development: the first PB prototype has just been delivered, while the J&T processor will start to be fabricated at a silicon foundry later this year. We expect all prototype buiding blocks to be available next spring. Our next important goal is a complete apeNEXT crate (400 GFlops peak processing power) running at the end of 2002. Subject to funding availability, we look forward to Tflops class machines in the year 2003.


text:
# ASTROPHYSICAL POLARIMETRY OF COSMOLOGICAL SOURCES

## 1 Introduction

Many astronomical objects are hidden from direct view by clouds of gas accompanied by small solid particles called dust by astronomers. If scattering particles are also present, a halo of scattered light is observable, with tangential position angles in a centro-symmetric pattern around the direction on the sky of the hidden object. For example, during a total solar eclipse, a beautiful outer atmosphere becomes visible (the “corona”) and much of this light is simply Thomson-scattered from the unseen bright solar photosphere. For Thomson scattering, spectral features are broadened somewhat because of the electrons’ thermal motions.

A radio galaxy is a galaxy with powerful radio emission, typically from two “lobes” on $`10^6`$ light year scales, symmetrically placed about the galaxy. The detected continuum, which is spatially resolved on the sky, is mostly starlight from the host galaxies. By contrast, a quasar is a pointlike optical continuum source with a similar radio manifestation. This powerful and variable continuum source may be thermal emission from gas heated as it falls toward the event horizon of a supermassive black hole.

The quasar pointlike optical sources are surrounded by fuzzy light, which is in fact mostly starlight from its host galaxy. It is now known that the powerful distant radio galaxies and quasars differ only in orientation with respect to the line of sight! In particular the powerful radio galaxies contain quasars hidden from direct view by dust clouds. We detect the hidden quasars via their halos of scattered light, just as in the case of the solar corona. The observed perfect tangential and centrosymmetric polarization patterns are evidentally unaffected by the long traverse through the universe.

Data on the material in this review are analyzed quantitatively by Wardle et al. (1997).

## 2 The Unified Model

Let’s look closely at the optical spectra of radio galaxies and quasars. Figure 1 shows a spectrum of the typical powerful radio galaxy Cygnus A, with a modest redshift of 0.07: the continuum light is mostly starlight. The emission lines however must come from clouds of ionized gas. The emission lines are due to recombination or collisional excitation of species photoionized by the quasar continuum. The spectrum is quite different from that of laboratory gases: the transitions which violate electric dipole selection rules are quite strong. This indicates a very low density gas. In such a gas, every time an excited state is populated, whether or not it’s metastable, a photon results from the subsequent spontaneous radiative de-excitation. The Einstein “A” coefficient, which is many orders of magnitude smaller for the forbidden transitions than for the permitted ones, makes no difference. Every excitation results in a “bankable” photon, whether the decay time is very long or very short.

The finite spectral width of both kinds of emission lines comes from bulk motions of the ionized gas clouds, and indicate internal motions within the emitting region of $`1,000`$ km/sec.

Figure 2 shows the spectrum of the quasar (or “broad line radio galaxy,” for any astronomers reading this) called 3C382. It consists of some starlight, the same type of low-density emission lines seen in the radio galaxy Cygnus A above, plus two more components. One is a mysterious continuum component, which is powerful and variable and may come from thermal radiation by optically thick matter accreting onto a supermassive black hole; and the other is the broad bases on the permitted lines only. These must come from a family of ionized gas clouds with a $`10,000`$ km/sec velocity dispersion, and a relatively high density so that the permitted lines are much much stronger than the forbidden lines. At these densities the gas is like that in the laboratory: the excited states are well populated for both permitted and forbidden transitions, so their ratio is determined by the Einstein A’s in this case.

During the 1980’s it was deduced from polarimetry and other data that the latter two components, the variable (and thus extremely compact) continuum source and the broad-emission-line region, are in fact present in powerful radio galaxies as well as in quasars, and can be detected via scattered polarized light. The electric vectors are always perpendicular to the radio source symmetry axes, and thus the photons’ last flights before being scattered into the line of sight were along the radio axes. From the point of view of the nuclear light, the other directions are blocked by a torus-shape collection of dusty (neutral) gas clouds. Figure 3 shows a cartoon.

The scattered light around a radio galaxy manifests as a bipolar reflection region (Fig. 4) and a polarized light spectrum identical to the total-flux spectra of quasars, complete with the nuclear continuum and the broad wings on the permitted emission lines.

The Unified Model for active galactic nuclei and quasars is reviewed in Antonucci 1993.

## 3 Implications for Particle Physicists

These scattering regions are visible at redshifts up to $`>2`$, though the data aren’t as pretty as for the nearby case shown here. The scatterers are sometimes small dust particles and sometimes free electrons, or a combination of the two. The point is, the electric vectors are all exactly perpendicular to the direction to a single point, where the hidden quasar lies. Any propagation effects that would spoil this perpendicular relationship are tightly constrained in size. One example is foreground gravitational lensing, which rotates the radius vectors due to shear, but leaves the electric vectors of the scattered light unchanged. Another is the hypothetical Lorentz-symmetry-violating effect discussed by Kostelecsky and Mewes. (Their constraint depends on the ratio of the wavelength of observation to the source distance, so will become several orders of magnitude more powerful when X-ray polarization can be used.)

Finally let me apprise you of a related effect which can sometimes be used for the same purpose. It is less accurate and robust. Consider the radio maps themselves. The radio photons are also polarized, but by a different effect: here the photons are intrinsically polarized because they derive from the synchrotron process. Since the double-lobed radio sources are axisymmetric to zeroth order, the net radio polarization integrated over the entire source tends to be either parallel or perpendicular to the structural axis. This just follows from the overall approximate axisymmetry of the radio morphology. Its polarization is not trivially understood from first principles as for the scattered optical light. And because it depends on the detailed “gastrophysics” or environmental “weather,” it isn’t precise; for some objects the polarization angle isn’t related to the symmetry axis at all. However, statistically we can compare the radio structural axes with the integrated radio polarization electric vector position angles, and limit deviations in a statistical sense: clearly if all the polarization angles (or all of the structural axes) were rotated by over a radian by some propagation effect, the statistical correlation with the radio axes would be destroyed.

## Acknowledgments

My polarimetry work is partially funded by NSF grants NSF AST96-17160 and NSF AST00-98719. Also, I thank Pat Ogle for Fig. 4 and for the Wardle et al. reference; Don Osterbrock for Figs. 1 and 2; and M. Urrey and Padovani for Fig. 3.

## References