text:
# Saturation Physics in Heavy Ion Collisions

## 1 Introduction

Saturation/Color Glass Condensate physics is based on the observation that the small-$`x`$ wave functions of ultrarelativistic hadrons and nuclei are characterized by a hard scale $`Q_s`$, known as the saturation scale glr ; mq ; BM ; mv . The scale $`Q_s`$ arises due to saturation of partonic densities at small-$`x`$ and is an increasing function of energy and atomic number of the nucleus bk ; LT ; mv . This large scale makes the strong coupling constant small $`\alpha _s(Q_s)1`$ leading to dominance of the classical gluonic fields in all high energy processes mv ; kjklw . Gluon production in high energy collisions is given by the classical field of the scattering color charges kmw . Corresponding gluon production cross section was found for $`pA`$ collisions in KM  and the effects of quantum evolution bk  were included in it in KT . The gluon production cross section for heavy ion collisions ($`AA`$) at the classical level has been studied both numerically KV  and analytically yuriaa . Since it is quite not clear at present how to include the effects of nonlinear quantum evolution bk  in the results of KV ; yuriaa , one has to construct models to describe the actual rapidity-dependent data produced in heavy ion collisions. Below we are going to show how some of these models, based on rather basic properties of saturation physics, provide a reasonably good description of RHIC data.

## 2 Particle Multiplicity from Saturation Models

### 2.1 Multiplicity at Mid-Rapidity Versus Centrality

Classical field $`A_\mu 1/g`$ leads to produced gluon multiplicity

$`{\displaystyle \frac{dN}{d^2kd^2bdy}}A_\mu A_\mu {\displaystyle \frac{1}{\alpha _s}}.`$ (1)
Gluon transverse momentum spectrum described by a single scale $`Q_s`$ can be written as

$`{\displaystyle \frac{dN}{d^2kd^2bdy}}={\displaystyle \frac{1}{\alpha _s}}f(k_{}/Q_s)`$ (2)
with some unknown function $`f(k_{}/Q_s)`$ to be determined by actual calculations. Integrating over $`k`$ and $`b`$ yields

$`{\displaystyle \frac{dN}{dy}}=\text{const}{\displaystyle \frac{1}{\alpha _s}}\pi R^2Q_s^2.`$ (3)
where the value of the constant is determined from $`f(k_{}/Q_s)`$. Following KN ; KL  we assume that the scale for the coupling constant in Eq. (3) is set by $`Q_s`$. (This step, of course, goes beyond the classical limit and assumes that Eq. (3) is valid even when running coupling corrections are included.) Then, as $`R^2A^{2/3}`$ and if $`Q_s^2A^{1/3}`$ mv ; kjklw , together with running coupling

$`\alpha _s(Q_s)={\displaystyle \frac{1}{b\mathrm{ln}Q_s^2/\mathrm{\Lambda }^2}}{\displaystyle \frac{1}{\mathrm{ln}A}}`$ (4)
we conclude from Eq. (3) that

$`{\displaystyle \frac{1}{A}}{\displaystyle \frac{dN}{d\eta }}\mathrm{ln}A.`$ (5)
For heavy ion experiments at different collision centralities we substitute $`A`$ by the number of participants $`N_{part}`$ so that Eq. (5) becomes

$`{\displaystyle \frac{1}{N_{part}}}{\displaystyle \frac{dN}{d\eta }}\mathrm{ln}N_{part}.`$ (6)
Eq. (6) allowed the authors of KN  to correctly predict the particle multiplicity at mid-rapidity at RHIC as a function of centrality at $`\sqrt{s}=130`$ GeV. A fit of particle multiplicity at other values of rapidity at $`\sqrt{s}=130`$ GeV taken from KL  is shown in Fig. 1.

### 2.2 Multiplicity as a Function of Energy

Eq. (3) allows one to test whether the scaling of total particle multiplicity with energy is consistent with saturation/Color Glass predictions. Using the fact that $`Q_s^21/x_{Bj}^\lambda `$ bk ; LT  in Eq. (3) and, for the moment, dropping the slower $`Q_s`$-dependence in $`\alpha _s`$ leads to

$`{\displaystyle \frac{dN/d\eta (\sqrt{s_1})}{dN/d\eta (\sqrt{s_2})}}=\left({\displaystyle \frac{\sqrt{s_1}}{\sqrt{s_2}}}\right)^\lambda .`$ (7)
Using PHOBOS data for total charge multiplicity at $`\sqrt{s}=130`$ GeV for most central collisions

$`{\displaystyle \frac{dN}{d\eta }}(\sqrt{s}=130\text{GeV})=\mathrm{\hspace{0.17em}555}\pm 12(stat)\pm 35(syst)`$ (8)
together with $`\lambda =0.25÷0.3`$ obtained in GBW  by analyzing HERA data, Kharzeev and Levin KL  predicted the total charge multiplicity at $`\sqrt{s}=200`$ GeV to be

$`{\displaystyle \frac{dN}{d\eta }}(\sqrt{s}=200\text{GeV})=\mathrm{\hspace{0.17em}616}÷634`$ (9)
in agreement with the later measured PHOBOS result

$`{\displaystyle \frac{dN}{d\eta }}(\sqrt{s}=200\text{GeV})=\mathrm{\hspace{0.17em}650}\pm 35(stat).`$ (10)
### 2.3 $`\frac{dN}{d\eta }`$ Versus $`\eta `$ and $`N_{part}`$

Describing rapidity distribution of the produced particles requires a little more modeling. Assuming $`k_T`$-factorization for particle production cross section

$`{\displaystyle \frac{d\sigma ^{AA}}{d^2kdy}}={\displaystyle \frac{2\alpha _s}{C_F}}{\displaystyle \frac{1}{\underset{¯}{k}^2}}{\displaystyle d^2q\varphi _A(\underset{¯}{q})\varphi _A(\underset{¯}{k}\underset{¯}{q})}`$ (11)
along with saturation-inspired unintegrated gluon distribution functions ($`\varphi _A(\underset{¯}{k})\alpha _s/\underset{¯}{k}^2`$ if $`k_{}>Q_s`$ and $`\varphi _A(\underset{¯}{k})S_{}/\alpha _s`$ if $`k_{}<Q_s`$) the authors of KL  produced an impressive fit of the PHOBOS data on charged particle multiplicities as functions of rapidity and centrality at $`\sqrt{s}=130`$ GeV shown in Fig. 2. The predictions made in KL  for similar multiplicity data at $`\sqrt{s}=200`$ GeV were also in good agreement with the later published BRAHMS data.

Phenomenological success of the saturation models presented above does not contradict the possibility of strong final state interactions leading to formation of quark-gluon plasma. As was argued in BMSS , thermalization in the saturation framework would not introduce any fundamentally new scales leaving Eqs. (3) and (6) practically unchanged. Late stage interactions are also not very likely to significantly modify the rapidity distribution of Fig. 2 due to causality constraints.

## 3 Non-Flow Contribution to $`v_2`$

The contribution of non-flow two-particle azimuthal correlations from the early stages of heavy ion collisions to the elliptic flow observable $`v_2`$ has been estimated in KTf  using a saturation-inspired model of particle correlations. Here we are going to construct a model-independent lower bound on these non-flow effects using $`v_2`$ extracted from the analysis of $`pp`$ data. We start with the definition of $`v_2(p_T)`$ for $`pp`$ corrected for uncertainty in the “reaction plane” definition (or, equivalently, defined through the two-particle correlation functions, such that $`v_2(p_T)<v_2>=<\mathrm{cos}2(\varphi _p\varphi _k)>_k`$ )

$`v_2^{pp}(p_T)={\displaystyle \frac{d^2k\frac{dN_{corr}^{pp}}{d^3pd^3k}\mathrm{cos}2(\varphi _p\varphi _k)}{\frac{dN^{pp}}{d^3p}\frac{dN^{pp}}{dy_k}+d^2k\frac{dN_{corr}^{pp}}{d^3pd^3k}}}\sqrt{{\displaystyle \frac{\frac{dN^{pp}}{dy_p}\frac{dN^{pp}}{dy_k}+\frac{dN_{corr}^{pp}}{dy_pdy_k}}{d^2pd^2k\frac{dN_{corr}^{pp}}{d^3pd^3k}\mathrm{cos}2(\varphi _p\varphi _k)}}}.`$ (12)
If we assume that the relative magnitude of the correlated terms in Eq. (12) compared to uncorrelated ones is roughly the same for all $`p_T`$, we can drop the former compared to the latter finding that Eq. (12) is approximately bounded from above by

$`v_2^{pp}(p_T){\displaystyle \frac{d^2k\frac{dN_{corr}^{pp}}{d^3pd^3k}\mathrm{cos}2(\varphi _p\varphi _k)}{\frac{dN^{pp}}{d^3p}\frac{dN^{pp}}{dy_k}}}\sqrt{{\displaystyle \frac{\frac{dN^{pp}}{dy_p}\frac{dN^{pp}}{dy_k}}{d^2pd^2k\frac{dN_{corr}^{pp}}{d^3pd^3k}\mathrm{cos}2(\varphi _p\varphi _k)}}}.`$ (13)
We want to estimate the contribution to $`v_2^{AA}(p_T)`$ of the non-flow correlations of the same physical origin as the ones giving rise to $`v_2^{pp}(p_T)`$ in Eq. (12). The contribution is

$`v_2^{AA}(p_T)|_{nonflow}={\displaystyle \frac{d^2k\frac{dN_{corr}^{AA}}{d^3pd^3k}\mathrm{cos}2(\varphi _p\varphi _k)}{\frac{dN^{AA}}{d^3p}\frac{dN^{AA}}{dy_k}}}\sqrt{{\displaystyle \frac{\frac{dN^{AA}}{dy_p}\frac{dN^{AA}}{dy_k}}{d^2pd^2k\frac{dN_{corr}^{AA}}{d^3pd^3k}\mathrm{cos}2(\varphi _p\varphi _k)}}}.`$ (14)
Using the fact that, approximately, both in saturation models and in the data $`dN/dyN_{part}`$ KN ; KL , we rewrite the right hand side of Eq. (13) as

$`v_2^{pp}(p_T)v_2^{AA}(p_T)|_{nonflow}\sqrt{{\displaystyle \frac{N^{AA}}{N^{pp}}}},`$ (15)
with $`N^{AA}`$ and $`N^{pp}`$ total particle multiplicities in $`AA`$ and $`pp`$ collisions, proportional to the average number of participants involved. Inverting Eq. (15) we obtain a lower bound

$`v_2^{AA}(p_T)|_{nonflow}v_2^{pp}(p_T)\sqrt{{\displaystyle \frac{N^{pp}}{N^{AA}}}}.`$ (16)
Preliminary analysis of $`pp`$ data yields $`v_2^{pp}(p_T)1`$ at high-$`p_T`$ rs . To get a lower bound we use $`N_{part}^{AA}=394`$ and $`N_{part}^{pp}=2`$ in Eq. (16) obtaining $`v_2^{AA}(p_T)|_{nonflow}\genfrac{}{}{0pt}{}{_>}{}7\%`$.

This work was supported in part by the U.S. Department of Energy under Grant No. DE-FG03-97ER41014. The preprint number is NT@UW-03-020.


format: text

text:
# The Galactic Bar

## 1. Introduction

The recognition that some galaxies contain central rectangular bars dates back at least as far as Curtis’ (1918) classification of a number of such systems as “$`\varphi `$-type spirals.” This fundamental property of disk galaxies was codified in the form we know today in Hubble’s (1926) classic “tuning fork” diagram, where barred and unbarred spirals form the prongs of the fork. Over the succeeding century, observations have improved, and in particular most recently infrared data have allowed us to view the stellar components of galaxies almost free from the confusing effects of dust obscuration. This mass of data has shown that bars are very common features, found in over a half of all spiral galaxies (Whyte et al. 2002).

It should therefore come as no particular surprise to find that the Milky Way is also a barred system. However, the commonness of such galaxies also raises the question of why one should care whether the Milky Way is barred or not. There are, after all, plenty of other galaxies that are oriented much more favorably than the Milky Way: not only is our own galaxy edge-on, rendering any bar almost invisible, but different parts of the Milky Way are different distances away, greatly complicating the geometry. Why not learn all there is to know about bars by studying the plethora of more accessible external systems? In fact, there are several good reasons why we might want to know about the Milky Way’s bar. First, it is natural to be inquisitive about our own neighborhood, and want to know what our galaxy would look like to an outsider (see Figure 1). Second, we have a mass of detailed information about the Milky Way that we have only been able to obtain because of its proximity. In order to put these data about small-scale structure, star streams, stellar populations, star formation, etc, into the broader context of galaxy evolution, we need to quantify the global properties of the galaxy in which they were found. Finally, and perhaps most importantly, although our position within the Galaxy complicates the understanding of its structure, it also provides us with information on the three-dimensional properties of its bar that is not available in any other system. Since we know very little about the three-dimensional shapes of bars, the Milky Way offers a unique laboratory for studying the structure of these systems.

In fact, our lack of understanding of their three-dimensional structure even leads to some ambiguity as to what exactly is meant by the term “bar.” In particular, it is not clear whether one should draw a distinction between a spheroidal (possibly triaxial) central bulge and a bar component. The different nomenclature may just reflect the limitations of our two-dimensional views of these systems: a triaxial structure that would be identified as a bar in a face-on system would be classified as a bulge component if the same system were viewed edge on. Indeed, there is now strong indirect evidence that bulges and bars are at the very least intimately related (Bureau & Freeman 1999, Merrifield & Kuijken 1999). In order not to prejudge the issue, this review will not draw any absolute distinction between the various putative components of the Galaxy in its central few kiloparsecs, but rather will look at the total distribution of stars and gas to see how far they are from axisymmetry, and hence how strongly the Milky Way is barred in that global sense.

In addition to referring the interested reader to the sophisticated treatment of the wealth of data available for such studies, we will attempt to develop some simple models that illustrate the particular challenges of investigating the bar in our own galaxy and the unique possibilities that its geometry admits. Section 2 looks at the properties of the bar that can be derived from the Milky Way’s gaseous components, while Section 3 explores the stellar distribution. In Section 4, we place the Milky Way in the context of other barred galaxies, and in Section 5 we take a brief look at what the future may hold for studies of the Galactic bar.

## 2. The Bar in Gas

Although most people think about distortions in starlight when considering bars, the earliest evidence that the Milky Way is barred came from it gaseous component. This is not surprising, as the Milky Way is much like an external system viewed edge-on, and in such cases a bar is almost impossible to detect photometrically. The best clues to the existence of a bar in an edge-on system come from the kinematics of its gaseous components, where the non-axisymmetric bar potential induces non-circular orbits in the gas.

In the case of the Milky Way, as Figure 2 shows, a plot of line-of-sight velocity versus Galactic longitude (an “$`l`$-$`v`$ diagram”) for the HI gas in the plane of the Galaxy reveals clear signatures of a non-axisymmetric distribution. There is a significant asymmetry between the gas properties at positive and negative longitudes, which is inconsistent with an axisymmetric disk. Further evidence for non-circular motions comes from the non-zero line-of-sight velocities of some of the HI gas at zero longitude: if the gas were following circular orbits, then all of its motion should be transverse to the line of sight at this point.

The strongest feature that illustrates these properties is classically known as the “3 kpc expanding arm,” because of its approximate radial location in the Galaxy. This feature, highlighted in Figure 2, is asymmetric about the center of the Galaxy, lying between Galactic longitudes of $`l^{}20^{}`$ and $`l^++35^{}`$. Tracing it through the center of the Galaxy reveals that it crosses zero longitude at a velocity of $`v_{\mathrm{los}}^0=53\mathrm{km}\mathrm{s}^1`$; presumably this loop in the $`l`$$`v`$ diagram re crosses $`l=0`$ at positive velocity, but this part of the feature is too faint to trace reliably. This non-zero line-of-sight velocity is inconsistent with an axisymmetric picture of the Galaxy with gas on circular orbits, as all gas motion in this direction should be transverse to the line of sight. It was this observation that led to the idea that the feature might be an expanding arm of material flung out from the center of the Galaxy in our direction. As we will see below, however, it is now understood as arising from the non-circular motions that occur in bars.

The interpretation of this feature as a bar is usually attributed to de Vaucouleurs (1964), who first made such a connection. However, his description of “focalized vortices plunging toward the nucleus” is not one that would be recognized today. It is therefore perhaps unsurprising that the properties he inferred for the bar turn out to be erroneous: even the side of the bar that lies closest to us was incorrectly identified. The earliest accurate modeling of this feature seems to have been undertaken by Shane (1971)111This paper is seldom noticed because it has the misfortune to be in one of the few volumes that is missing from the ADS, so is passed over in many literature reviews., who pointed out that this feature could be explained by an elliptical ring of gas with its major axis oriented at an angle of $`\varphi 20^{}`$ to the line of sight, and a short-to-long axis ratio of $`0.6`$. Subsequent analysis (e.g. Peters 1975) confirmed this result, and even recent full gas-dynamical simulations (e.g. Fux 1999) all find that a feature of this kind is required to fit the observed $`l`$-$`v`$ diagram.

In order to obtain some insight into the reason for the robustness of this result, it is instructive to see if it can be understood on the basis of a simple calculation. Specifically, can one go from the three basic observed quantities, $`l^+`$, $`l^{}`$ and $`v_{\mathrm{los}}^0`$ to obtain the three parameters of the ellipse, $`b`$, $`a`$ and $`\varphi `$? The values of the spatial coordinates, $`l^+`$ and $`l^{}`$, clearly do not provide enough information to solve for all three ellipse parameters, but obviously provide some constraint: the fact that they are not equal rules out an axisymmetric feature. One can therefore narrow down the range of possible ellipses using $`l^+`$ and $`l^{}`$, and, for example, solve for the axis ratio $`b/a`$ as a function of the feature’s angle $`\varphi `$ to the line of sight. This constraint is illustrated in Figure 3.

A further constraint comes from the non-zero velocity observed at zero longitude. This observation can be interpreted rather crudely in the epicyclic approximation, in which the elliptical orbit of the gas is constructed from a parent circular orbit of angular frequency $`\mathrm{\Omega }`$ and radius $`R_{\mathrm{bar}}`$ upon which an elliptical perturbation of radial amplitude $`X`$ and angular frequency $`\kappa `$ is superimposed (Binney & Tremaine 1987). At an arbitrary longitude, the line-of-sight velocity of gas following such motions is an ugly mixture of components of the circular motion plus the radial and tangential perturbations. However, at $`l=0`$ the only line-of-sight component is the radial perturbation, which must produce the observed value of $`v_{\mathrm{los}}^0`$. Writing the radial coordinate of gas as $`R=R_{\mathrm{bar}}+X\mathrm{cos}\kappa t`$, we find that $`\dot{R}=X\kappa \mathrm{sin}\kappa t`$. The situation is further simplified if it is assumed that the speed of circular orbits in the Milky Way, $`v_c`$, is approximately constant, as observed in most external spiral galaxies (although such an assumption must break down at small radii where the rotation speed must drop toward zero). For such a flat rotation curve, $`\kappa =\sqrt{2}\mathrm{\Omega }=\sqrt{2}v_c/R_{\mathrm{bar}}`$ (Binney & Tremaine 1987). Putting all this information together, we find that

$$\frac{X}{R_{\mathrm{bar}}}=\frac{v_{\mathrm{los}}^0}{\sqrt{2}v_c}\times \frac{1}{\mathrm{sin}(\sqrt{2}\varphi )},$$
(1)
and hence the ellipse’s axis ratio

$$\frac{b}{a}=\frac{R_{\mathrm{bar}}X}{R_{\mathrm{bar}}+X}=\frac{1X/R_{\mathrm{bar}}}{1+X/R_{\mathrm{bar}}},$$
(2)
which is just a function of $`\varphi `$. This function is also plotted in Figure 3.

The spatial and the kinematic constraints are both met where the lines intersect. It is therefore apparent that a self-consistent model requires an ellipse of gas with an axis ratio of $`b/a0.55`$ at an angle to the line-of-sight of $`\varphi 25^{}`$. A number of somewhat dubious approximations have gone into this argument – the epicyclic approximation has been pushed well beyond the perturbative regime where it is strictly valid, for example. However, full hydrodynamical modeling of the Galactic $`l`$$`v`$ diagrams in HI and CO come up with essentially identical parameters for this gas ring (Fux 1999). If this feature is interpreted as a signature of the Galactic bar, then we have a measure of its angle to the line-of-sight; as we shall see below, however, the axis ratio requires further interpretation.

## 3. The Bar in Stars

Seeing the stellar bar in the Milky Way presents particular difficulties. As has already been mentioned, it is viewed from a most unfavorable direction, but in addition it is seen through many magnitudes of optical extinction, and it covers many square degrees on the sky. These latter challenges were finally overcome when the COBE map of the sky was published, and Galactic astronomers began analyzing the parts of the data that the cosmologists discarded (see Figure 4). The COBE satellite’s infrared coverage spanned the optimal point in the spectrum where emission is still dominated by stars yet extinction is minimized, and the all-sky coverage of the survey meant that the large extent of the structure was no longer a problem, so finally some photometric insight into the bar’s properties could be obtained.

In external edge-on galaxies, bars are completely undetectable from photometric data. The only clues to their presence come from indirect indicators such as a plateau in the light distribution at small radii (e.g. de Carvalho & da Costa 1987). However, such features could equally well be axisymmetric disk features, so no firm conclusions can be drawn. One can do rather better if one measures optical spectra along the major axes of such galaxies, as both the stellar kinematics (as determined from absorption-line spectra) and gas kinematics (as seen in emission-line spectra) show complex structure in the line-of-sight velocity distribution as a function of position (Kuijken & Merrifield 1995), directly analogous to the structure seen in $`l`$$`v`$ diagram of the Milky Way (see Figure 2). It was through such studies that a link was established between the presence of a bar and boxy, or even double-lobed, isophotes in the central bulges of these systems (Bureau & Freeman 1999, Merrifield & Kuijken 1999). Such a link was predicted by numerical simulations, which showed that initially-thin bar structures should buckle perpendicular to their plane, creating triaxial objects with a double-lobed shape when viewed edge-on (e.g. Combes & Sanders 1981). The double-lobed structure apparent in some of the isophotes shown in Figure 4 therefore provides indirect evidence that the Milky Way is barred.

In the case of the Milky Way, however, we can obtain rather more direct evidence for the presence of a bar. Generally speaking, the complex geometry involved in observing a galaxy from the inside greatly hinders the interpretation of the data. On this occasion, though, it actually helps. The finite distance to the bar means that there is a parallax effect: the two ends of the bar are effectively viewed from different directions, causing the projected scale-length of the light distribution to be different on the two sides of the Galactic center. This asymmetry is clearly visible in the map of the light distribution shown in Figure 4.

Perhaps unsurprisingly, there is not enough information in such a two-dimensional map to derive a unique model for the three-dimensional distribution of starlight, even if quite strong symmetries are imposed on the bar \[as highlighted by Zhao (2000)\]. Nonetheless, a fairly consistent view has emerged as to the basic structural parameters of the bar using the observed asymmetry in the photometry, and once again a simple model provides a useful illustration as to why these parameters are fairly robustly measured.

The simplest plausible model for the bar is a triaxial Gaussian ellipsoid, with principal axes of scalelengths $`a`$, $`b`$, and $`c`$ in the $`x`$, $`y`$ and $`z`$ directions, respectively (where $`z`$ is the distance perpendicular to the Galactic plane. If viewed edge-on from a large distance at an angle $`\varphi `$ to the $`x`$-axis, the projected distribution of light would also be Gaussian with a scale length $`c`$ perpendicular to the plane and a scale length

$$A(\varphi )=\left[\frac{\mathrm{sin}^2\varphi }{a^2}+\frac{\mathrm{cos}^2\varphi }{b^2}\frac{\mathrm{sin}^2\varphi \mathrm{cos}^2\varphi \left(\frac{1}{a^2}\frac{1}{b^2}\right)}{\frac{\mathrm{cos}^2\varphi }{a^2}+\frac{\mathrm{sin}^2\varphi }{b^2}}\right]^{1/2}$$
(3)
in the projected plane. However, since we are not an infinite distance away, observations at $`\pm l`$ observe the bar from somewhat different directions, which will yield different projected scalelengths, $`A^+=A(\varphi +l)`$ and $`A^{}=A(\varphi l)`$. Thus, we essentially have two observables at any given longitude, $`\overline{A}(l)=(A^++A^{})/2`$ and $`f(l)=A^+/A^{}`$; at $`|l|15^{}`$, for example, we can see from Figure 4 that $`\overline{A}6.0^{}`$ and $`f1.25`$. There are three in-plane physical quantities that we would like to derive – $`a`$, $`b`$ and $`\varphi `$ – so clearly the problem is under-determined, illustrating the degeneracy in deprojecting from two dimensions to three. However, we can solve for $`a`$ and $`b/a`$ as functions of $`\varphi `$, as we did above for the gaseous component; the resulting values for $`b/a`$ as a function of $`\varphi `$ are also shown in Figure 3.

If we adopt the angle of $`\varphi 25^{}`$, as derived from the gaseous component, we see that the stellar component must have an axis ratio of $`b/a0.6`$. Folding in the value of the vertical scaleheight, $`c`$, as measured in Figure 4, we find that the bar is a triaxial structure with axis ratios $`a:b:c1:0.6:0.4`$. Again, although this result has been derived on the basis of a simplified model, it reproduces the results of much more sophisticated calculations: Fux’s (1999) matching of N-body simulations to the COBE map and Binney, Gerhard & Spergel’s (1997) direct deprojection of the same data both predict aspect ratios of $`a:b:c1:0.6:0.4`$. This good agreement must occur because the amount of information on the structure of the bar is rather limited in the two-dimensional COBE map, and all the various modeling techniques latch onto essentially the same constraint arising from the asymmetry in the distribution of starlight.

Note that this analysis uses all the light from the central region, and does not attempt to separate out distinct bar and bulge components. The fact that it is possible to produce a coherent picture of the bar region without a separate axisymmetric bulge component adds further weight to the argument that these features should be viewed as different projections of the same physical object rather than as distinct entities. It is also worth noting that this measurement illustrates why it is worth persisting with the challenging task of studying the Milky Way’s bar: such a fundamental measurement of the three-dimensional shape of a bar has not been possible in any other system.

One use that we can immediately make of this unique measurement is to compare the shape of the bar in gas to that in the stellar component. Since the gas disk is essentially confined to the Galactic plane, the stellar bar is clearly rounder perpendicular to the plane than the gaseous bar. It is also apparent from Figure 3 that the stellar component is also rounder in the plane, with an axis ratio of $`b/a0.6`$ as compared to the value of $`0.5`$ for the gaseous component. This difference in shape can be straightforwardly understood if the stars do not lie on the closed orbits that the gas follows, but instead oscillate about these orbits, adding a more random component to the mean streaming, thus puffing out the structure of the stellar bar.

It is not yet clear where the stars that currently make up this triaxial bar formed. Some have argued that the stars predate the bar, with little star formation once the bar has formed (e.g. Cole & Weinberg 2002). However, there is also strong evidence for continuing star formation in the bar. As Figure 5 shows, $`l`$$`v`$ diagrams for OH/IR stars measured by Sevenster (1999) reveal that those with massive progenitors (which implies they must be young) pick out the same 3 kpc arm feature seen in the gas. Since this feature is a clear signature of material in an elliptical bar-like structure, it is apparent that star formation in the bar is an on-going process.

Such continuing star formation suggests that one should view the puffed-up nature of the stellar bar as arising from a continuing “heating” process, in which the stars are progressively scattered away from the planar gas distribution in which they formed, making their distribution rounder over time. Indeed, the lower mass stars in Figure 5, which are on average older than the more massive stars, show the somewhat larger scatter in velocity that one would expect from a heated population. Another good tracer of intermediate and old populations is provided by planetary nebulae (PNe), which are relatively easy to identify from their emission lines, and whose kinematics can be easily measured from these lines; an analysis of PNe in the region of the Galactic center (Beaulieu et al. 2000) reveals a distribution in the $`l`$$`v`$ diagram very similar to the older OH/IR stars in Figure 5, indicating that they have also had extra random motion injected into their velocities.

A simple model for this extra random component involves adding an amplitude $`\delta `$ in quadrature to the gaseous highly-flattened ellipsoidal distribution, in order to represent the excursions that the heated stars make from the closed gaseous orbits. This heuristic model makes no attempt to explain the details of the scattering process and the resulting complex distribution function of the stars, but such an isotropic swelling of the bar is the simplest manner in which the heating of the stellar component might occur.

Figure 6 shows the effect on the shape of adding progressively larger random components of amplitude $`\delta `$ to an initial distribution shaped like the Milky Way’s gaseous component with an aspect ratio $`a:b:c1:0.5:0`$. Given the simple-minded nature of this calculation, it is perhaps surprising that one can simultaneously reproduce the values of both $`b/a0.6`$ and $`c/a0.4`$ with a single value for the amplitude of the random motions, $`\delta 0.4a`$. This agreement may just be a fortuitous accident, but perhaps it sheds some light on the mechanism by which galactic bars become heated.

## 4. The Milky Way as a Barred Galaxy

We now look beyond the Milky Way to try and place it in the broader context of spiral galaxies in general. Again, this is partly just to satisfy a natural curiosity about our home galaxy. However, it also has a broader significance if we are to use the results derived for the Milky Way to infer anything about galaxies in general: it is always somewhat risky to draw general inferences from a single object, but it is downright foolhardy to generalize in this way if we already know that the Milky Way is in some way unusual.

The first rigorous attempt to classify the Milky Way as if it were an external galaxy was made by de Vaucouleurs & Pence (1978), who looked at a variety of qualitative and quantitative observations of the Milky Way, and concluded that it should be designated an SABbc(rs) system. One has to suspect, however, that this large string of characters was more a reflection of uncertainty than of a detailed understanding of the Galaxy’s structure. In the case of bar strength, for example, the paper has very little to say on how the classification was made. This reticence is not particularly surprising, as almost nothing was known about the Galactic bar at that time. The SAB classification therefore just reflects the authors’ desire to hedge their bets against future measurements.

Twenty-five years on, we have much more data and information, including the results reviewed in this paper, so we are now in a position to revisit this question of the Milky Way’s classification in comparison to other galaxies. We can also now exploit the greater access to computers and quantitative data in order to move away from subjective classifications by eye to objective measures of quantities such as bar strength. However, even today quantifying bar strength turns out to be a non-trivial issue, since the term “strength” is not well defined: is a very small, highly elliptical feature at the center of a galaxy a stronger or weaker bar than a much larger feature that is closer to axisymmetric? A variety of measures have been deployed, ranging from simple observational quantities like deprojected isophotal flattenings (Martin 1995) to complex physically-motivated algorithms involving calculating the non-radial gravitational field of the galaxy (Buta & Block 2001).

In the current analysis, we recognize the intrinsically two-dimensional nature of the classification by specifying both an isophotal ratio for the bar as viewed face on, and the fraction of the total galaxy luminosity that lies within this isophote. Figure 7 shows these parameters for the galaxies in the Ohio State Bright Spiral Galaxy Survey as well as the Milky Way. It is apparent from this figure that the conventional SA/SAB/SB classifications of bar strength are primarily driven by the isophote flattening rather than the luminosity of the bar, although there is significant scatter in the quantitative measurements for each of the qualitative classes. Also, somewhat frustratingly for anyone who would like to classify the Milky Way definitively into one of these traditional categories, our galaxy seems to lie in a region of this parameter space where one finds SA, SAB and SB galaxies. It might appear that we have not come very far from de Vaucouleurs and Pence’s tentative SAB classification after all! However, the uncertainty in the classical classification of the Milky Way bar now reflects the limitations of these qualitative schemes rather than any paucity of data on the Milky Way’s bar.

## 5. The Future

At some level, the Galactic bar is a solved problem. As outlined in this review, there is now a mass of data that provides a coherent picture of the Milky Way as a barred galaxy. There are certainly subtleties that remain to be discovered, but it might reasonably be argued that bars are not subtle things, and we already know the salient properties of the Galactic bar to the accuracy that is likely to be of any interest in studies of galactic structure.

Although there is some truth in this statement, it also reflects our ignorance of the fine-scale properties of bars, which arises from the difficulty of studying these objects in detail at the large distances of external galaxies. There is also a good prospect that important clues to the formation and evolution of bar structures, which are lacking in the broad-brush picture, can be gleaned from this detailed information. The Milky Way therefore offers a unique laboratory for the detailed study of a bar (albeit one viewed from a rather disadvantageous direction).

It is already known that bars can have quite complex structures, with smaller nuclear bars sometimes found within bigger bars at a variety of orientations \[as first noted by de Vaucouleurs (1974)\]. These nested bars are quite common, arising in around a quarter of all barred galaxies (Laine et al. 2002), so it would not be too surprising to find similar complexity in the arrangement of stars in the Milky Way’s bar. Indeed, preliminary analysis of the 2MASS star counts by Alard (2001) has shown that the asymmetry seen in the COBE data that reveals the presence of a bar (see Section 3) flips in sign within the central couple of degrees, indicative of a secondary bar at a radically different orientation. Our closeness to the Galactic bar means that there is the prospect of probing much closer to the center of the Galaxy than is possible in external barred systems, so we might hope to go down at least two more orders of magnitude in scale to determine whether this bar-within-a-bar is just the beginning of a hierarchy of structures nested like Russian dolls.

Huge data sets like the 2MASS star counts are certainly taking the study of the Milky Way’s structure to a new level. The analysis by Cole & Weinberg (2002) showed that the bar can be detected not only in the asymmetry of the counts, but also more directly by using infrared carbon stars as standard candles to get a crude measure of the third dimension. However, the real revolution will arrive when we unlock the other three dimensions of phase space by measuring the velocities of large samples of stars. There is already some kinematic data available from studies of line-of-sight velocities and proper motions of stars in low-extinction windows toward the Galactic center, and even this limited dynamical information is quite a challenge to reconcile with detailed models of the bar (Häfner et al. 2000), although line-of-sight contamination by non-bar stars may well compromise the analysis. When the next generation astrometric satellite GAIA is launched (ESA 2000), we will have measures of both line-of-sight velocities and proper motions for stars all across the Galaxy, including those in the bar region, as well as obtaining parallax distances to these objects to an accuracy of $`10\%`$ which will almost eliminate any problems from line-of-sight contamination. By assigning ages through position on the color–magnitude diagram, one will be able to detect the effects of any secular evolution in the orbit structure \[as was done so successfully to measure evolution in the Galactic disk using data from GAIA’s predecessor Hipparcos (Dehnen & Binney 1998)\]. These data will push dynamical analysis into a new regime where very sophisticated modeling will be required to study the details of the orbit families that are populated in the bar. The challenges to modeling and interpreting such huge data sets are far from trivial, but they will provide a completely new perspective on the formation and evolution of galactic bars.

The history of astronomy provides a long list of discoveries which have shown that we are nowhere particularly special in the Universe: Copernicus moved us out of the center of the Solar System; Shapley moved the Solar System to a dull suburb of the Milky Way; Hubble proved that the Milky Way is just one of billions of similar stellar metropolises. It should therefore come as no particular surprise to discover that the Milky Way is like the majority of spiral galaxies in that it contains a bar. Although the ordinariness of our situation is in some ways rather a disappointment, it does have one major up side: when in the future we want to learn about generic features of the Universe like bars in disk galaxies, we will not have very far to look.

## References

Alard, C. 2001, A&A, 379, L44

Binney, J., Gerhard, O. & Spergel, D. 1997, MNRAS, 288, 365

Binney, J. & Merrifield, M. 1998, Galactic Astronomy (New Jersey: Princeton University Press)

Binney, J. & Tremaine, S. 1987, Galactic Dynamics (New Jersey: Princeton University Press), 120

Bureau, M. & Freeman, K. C. 1999, AJ, 118, 126

Buta, R. & Block, D. 2001, ApJ, 550, 243

Cole, A. A. & Weinberg, M. D. 2002, ApJ, 574, L43

Combes, F. & Sanders, R. H. 1981, A&A, 96, 164

Curtis, H. 1918, Pub. Lick Obs., 13, 12

de Carvalho, R. R. & da Costa, L. N. 1987, A&A, 171, 66

Dehnen, W. & Binney, J. J. 1998, MNRAS, 298, 387

de Vaucouleurs, G. 1964, in IAU Symp. 20, The Galaxy and the Magellanic Clouds, ed. F. J. Kerr (Canberra: Australian Academy of Science), 195

de Vaucouleurs, G. 1974, in IAU Symp. 58, The Formation and Dynamics of Galaxies, ed. J. R. Shakeshaft (Dordrecht: Reidel), 335

de Vaucouleurs, G. & Pence, W. D. 1978, AJ, 83, 1163

ESA 2000, GAIA: Composition, Formation and Evolution of the Galaxy, Technical Report ESA-SCI(2000)4. (The scientific case is available on-line at http://astro.estec.esa.nl/GAIA).

Fux, R. 1997, A&A, 327, 983

Fux, R. 1999, A&A, 345, 787

Häfner, R., Evans, N. W., Dehnen, W. & Binney, J. 2000, MNRAS, 314, 433

Hubble, E. P. 1926, ApJ, 64, 321

Kuijken, K. & Merrifield, M. R. 1995, ApJ, 443, L13

Laine, S., Shlosman, I., Knapen, J. H. & Peletier, R. F. 2002, ApJ, 567, L97

Martin, P. 1995, AJ, 109, 2428

Merrifield, M. R. & Kuijken, K. 1999, A&A, 345, L47

Peters, W. L. III 1975, ApJ, 195, 617

Shane, W. W. 1971, A&A, 16, 118

Sevenster, M. N. 1999, MNRAS, 310, 629

Whyte, L., et al. 2002, MNRAS, 336, 1281

Whyte, L., et al. 2003, in preparation

Zhao, H. S. 2000, MNRAS, 316, 418


text:
# Discrete Hamiltonian evolution and quantum gravity

## I Introduction

In all approaches to quantum gravity the idea of a discrete structure of spacetime at the Planck scale is fundamental. In lattice approaches, such as Regge calculus, dynamical triangulations, and other path integral methods, discreteness in space and time is built in right from the startloll . In canonical methods, such as loop quantum gravity, discreteness of space arises as a rather natural consequence of Dirac quantization av1 ; av2 ; av3 ; length , and it is possible that a form of discrete time evolution may emerge from a regulated action of the Hamiltonian constraint. This is already manifested to some extent in spin foam modelsfoam , which may be viewed as ”covariantizations” of loop quantum gravity models.

In attempts to construct a path integral for quantum gravity, a standard starting point is a classical non-Hamiltonian discretization using triangulations. It is interesting to ask whether a phase space path integral for quantum gravity can be defined starting from a classical discretization that is Hamiltonian. If possible, this would provide a regularization that is more closely tied to the classical theory than triangulations or spin foam methods. There are two main approaches to the problem: (i) fix all coordinate gauge conditions classically and discretize the reduced Hamiltonian system, or (ii) attempt a discretization of the constrained system and study its consequences for the Hamiltonian structure. The primary requirement for canonical quantization of a discrete system is to ensure that the fundamental Poisson bracket is preserved under discrete classical evolution. In the second approach, the consequences of discretization on the constraint algebra also require scrutiny.

These questions have been studied to some extent in the literature. Discrete evolution schemes for unconstrained Hamiltonian dynamics appear in bander , with further refinements in moncrief , where schemes preserving the Poisson bracket are presented. More recently, the applicability of discrete evolution schemes to constrained Hamiltonian theories have been discussed gp , with a view to applicability to quantum gravity.

From elementary numerical methods, it is known that there are implicit and explicit approaches to discretization of differential equations. Although either approach may be useful for studying classical evolution, explicit methods are more natural for quantization. This is because discrete Hamiltonians that are functions of canonical variables at more than one time instant are effectively non-local in time, so it is unclear whether canonical quantization makes sense. It is therefore necessary to study classical discretizations for which the discrete Hamiltonian does not depend on more than one time point. This may restrict somewhat the potential applicability of implicit schemes.

For the purpose of developing a fundamental lattice approach to quantum gravity starting from a discretization of Hamilton’s equations, it is necessary to ensure that canonical structures are not affected by the procedures used. This is because classical preservation of Poisson brackets is a requirement for a clean approach to canonical quantization. The main purpose of this paper is to investigate this aspect. We study certain issues associated with discrete time evolution. We do not address space discretization, since our main concern is the extent to which time discretization can maintain Hamiltonian structures for quantization. Specifically we address the following questions: (i) What are the general classes of discretization? (ii) What are the conditions required to maintain the fundamental Poisson brackets under evolution for a given Hamiltonian system? (iii) Are there restrictions on allowed Hamiltonians? (iv) How does the approach work in simple systems, and (v) How viable is the approach for quantum gravity?

Our starting points are general forms of explicit and implicit discretizations containing arbitrary functions of phase space variables. In Section II we derive the necessary conditions on these functions to preserve the fundamental Poisson bracket for both types of schemes. We also study the consequences of discretization for constraint evolution for the parametrized particle. In Section III we apply the approaches to the parametrized scalar field, FRW cosmological models, and BF theory. In the final section we dicsuss some consequences of our results for applicability to quantum gravity.

## II Discrete Hamiltonian evolution

As for differential equations, there are two general approaches that may be utilised to obtain a discretization of Hamilton’s equations that preserve the canonical Poisson bracket. The explicit schemes, as the name suggests, give explicit functions of the form

$$p_{n+1}=f(q_n,p_n),q_{n+1}=g(q_n,p_n)$$
(1)
for evolution. This is in contrast to implicit schemes, where the general equations are of the form

$$F(p_{n+1},q_{n+1},p_n,q_n)=0,G(p_{n+1},q_{n+1},p_n,q_n)=0.$$
(2)
We consider both types of discretization, and outline the conditions under which the schemes preserve Poisson brackets, and thus have a utility for quantization at the discrete level.

### II.1 Explicit schemes

Consider a one particle system with canonical phase space variables $`(q,p)`$, and arbitrary Hamiltonian $`H(q,p)`$. A first attempt at time discretization of Hamilton’s equations with step $`ϵ`$ might be the simplest finite differencing scheme

$`q_{n+1}`$ $`=`$ $`q_n+ϵ{\displaystyle \frac{H(q_n,p_n)}{p_n}},`$
$`p_{n+1}`$ $`=`$ $`p_nϵ{\displaystyle \frac{H(q_n,p_n)}{q_n}}`$ (3)
However there is an immediate problem with this: it does not preserve the Poisson bracket under evolution for even the standard Hamiltonian $`H(q,p)=p^2/2m+V(q)`$; the order $`ϵ^2`$ term does not vanish unless $`V(q)=0`$. To rectify this it is necessary to introduce higher order terms in $`ϵ`$ in the defining equations of the discretization. This might be done using explicit or implicit finite differencing schemes. We choose a general explicit approach and study its consequences. The reasons are simplicity and potential utility for canonical quantization.

The defining equations for the discrete evolution we study are the following:

$`q_{n+1}`$ $`=`$ $`q_n+ϵ{\displaystyle \frac{H(q_n,p_n)}{p_n}}+ϵ^2\alpha (q_n,p_n)`$ (4)
$`p_{n+1}`$ $`=`$ $`p_nϵ{\displaystyle \frac{H(q_n,p_n)}{q_n}}+ϵ^2\beta (q_n,p_n).`$ (5)
$`\alpha ,\beta `$ are at this stage arbitrary functions. The Poisson bracket of the evolved variables is

$`\{q_{n+1},p_{n+1}\}`$ $`=`$ $`1ϵ\left[\{q_n,{\displaystyle \frac{H}{q_n}}\}+\{p_n,{\displaystyle \frac{H}{p_n}}\}\right]`$ (6)
$`+ϵ^2\left[\{\alpha ,p_n\}+\{q_n,\beta \}\{{\displaystyle \frac{H}{p_n}},{\displaystyle \frac{H}{q_n}}\}\right]`$
$`+ϵ^3\left[\{{\displaystyle \frac{H}{p_n}},\beta \}\{\alpha ,{\displaystyle \frac{H}{q_n}}\}\right]`$
$`+ϵ^4\{\alpha ,\beta \}`$
This equation displays the conditions required for preserving the canonical Poisson bracket. The first order term in $`ϵ`$ vanishes identically because the Hamiltonian is a function of the discrete phase space variables at a single time point. (See below for discussion of a discretization where this is not the case.) The remaining task is to see what functions $`H`$, $`\alpha `$, and $`\beta `$ lead to vanishing of the three higher order terms in $`ϵ`$. It is obviously better to fix $`H`$ by looking at systems of interest, and then attempt solutions for $`\alpha `$ and $`\beta `$.

In order to obtain a general procedure for solving these conditions, it is useful to consider first the Hamiltonian independent $`ϵ^4`$ term. Its solutions can be classified into one of the following four cases: (i) $`\alpha =0`$, (ii) $`\beta =0`$, (iii) $`\alpha (q_n),\beta (q_n)`$ or $`\alpha (p_n),\beta (p_n)`$, and (iv) $`\alpha =\beta 0`$. For a given Hamiltonian, the other conditions can be written out explicitly.

The restrictions that arise may be seen by considering the standard Hamiltonian

$$H(q,p)=\frac{p^2}{2m}+V(q)$$
(7)
The conditions have a unique solution: Case (ii) gives

$$\alpha (q)=\frac{1}{m}\frac{V(q)}{q}$$
(8)
The remaining cases (i),(iii) and (iv) do not give a solution. This result is surprising given that there are numerous possibilities for discretizing classical differential equations. The lesson is that preserving canonical Poisson brackets imposes rather strong conditions on allowable discretizations.

As a further illustration consider the parametrized particle. This case is also useful for seeing what happens to the constraints under discrete evolution. The canonical variables are $`(t,p^t)`$ and $`(q,p)`$, which depend on a time parameter $`\lambda `$. The Hamiltonian constraint is

$$H:=p^t+\frac{p^2}{2m}+V(q)=0$$
(9)
Following (5), discrete evolution with lapse function $`N`$ is defined by

$`q_{n+1}`$ $`=`$ $`q_n+Nϵ{\displaystyle \frac{p_n}{m}}+ϵ^2\alpha (q_n,p_n)`$
$`p_{n+1}`$ $`=`$ $`p_nNϵ{\displaystyle \frac{V(q_n)}{q_n}}+ϵ^2\beta (q_n,p_n).`$
$`p_{n+1}^t`$ $`=`$ $`p_n^t`$
$`t_{n+1}`$ $`=`$ $`t_n+Nϵ.`$ (10)
From these equations it is clear that preservation of the Poisson bracket goes through in the same way as for the unparametrized case, with $`\beta =0`$ and $`\alpha `$ given by (8) up to a factor of $`N`$.

The second check is the evolution of constraints using (10). Is $`H_{n+1}H(t_{n+1},p_{n+1}^t,q_{n+1},p_{n+1})=H_n`$? A straightforward computation using the evolution equations gives

$`H_{n+1}`$ $`=`$ $`p_n^t+{\displaystyle \frac{1}{2m}}\left[p_nϵNV+ϵ^2\beta \right]^2+V(q_{n+1})`$ (11)
$`=`$ $`H_n+{\displaystyle \frac{ϵ^2}{2m}}\left[N^2(V)^2+\beta ^2+2p_n\beta +2m\alpha \right]+𝒪(ϵ^3),`$
where $`V=V(q_n)/q_n`$. From this it is clear that $`H_{n+1}`$ does not equal $`H_n`$: the Taylor expansion of the potential contains terms that cannot be cancelled by any choice of the functions $`\alpha ,\beta `$ (which are in any case effectively fixed by the requirement of preserving the canonical Poisson bracket). The discretization has made the constraint second class.

At this stage there are two possible ways to proceed: (i) impose the Hamiltonian constraint strongly since it is not preserved in time, or (ii) fix the free function $`N`$ such that the additonal terms on the right hand side of (11) vanish. The latter amounts to $`N`$ becoming a fixed function of $`n`$, and hence the time step $`ϵN`$ also becomes a function of $`n`$. Thus, making the evolution consistent via (ii) requires a variable time step in this sense, which concomitantly means fixing the lagrange multiplier.

From the theory of constrained Hamiltonian systems it is possible to see that these two choices are in fact equivalent: Gauge fixing leads to first class constraints becoming second class pairs with the gauge conditions. Since we started with a first class system (the parametrized particle), which became second class as a result of time discretization, we can ask what time gauge condition has been implicitly imposed by the discretization that led to the Hamiltonian constraint not being preserved in time.

To answer this question let us recall the procedure for the continuous case in the canonical theory, where time gauge fixing is accomplished by setting a function of the phase space variables to equal the time parameter. Consider a phase space function $`f(t,p^t,q,p)`$ and impose the gauge condition

$$\chi :=\lambda f=0.$$
(12)
This fixes parameter time $`\lambda `$. Preserving the gauge condition under evolution requires

$$\frac{d\chi }{d\lambda }=\frac{\chi }{\lambda }+\{\chi ,NH\}=0.$$
(13)
This fixes the lapse as a function of $`f`$

$$N=\frac{1}{\{f,H\}}$$
(14)
Conversely, fixing $`N`$ first leads via the above equations to a connection with the canonical gauge fixing function $`f`$. Thus a discretization may also be determined by fixing the lapse, which may be viewed as imposing a variable time step. This of course, is what is routinely done in numerical relativity, with the difference that in that case it is not necessary to preserve the canonical Poisson bracket for quantization.

For the discrete evolution we have defined, this procedure may be used to extract the gauge fixing condition. This is done by first imposing $`H_{n+1}=0`$ strongly since it is second class. This gives an equation for $`N`$. The gauge fixing function $`f`$ is then determined by solving (14), which in general is rather messy. For example, for the above discretization of the parametrized particle, it is not possible to solve for $`N`$ to all orders in $`ϵ`$ due to the presence of all orders of $`N`$ in $`H_{n+1}`$ (11). (This is because the expansion of $`V(q_{n+1})`$ in powers of $`\mathrm{\Delta }q_n=q_{n+1}q_n`$ contains powers of $`N`$.) Thus only an approximate solution for $`N`$ is possible, up to a specified order in $`ϵ`$.

The lessons from this section are two-fold: (i) The conditions on allowable explicit time discretizations are rather strong if canonical structures are to be preserved. (ii) For constrained theories with time reparametrization invariance, time discretization is equivalent to rather complicated gauge fixing conditions.

### II.2 Implicit schemes

There are implicit schemes for time discretizations discussed in the literature. We mention two examples below before introducing some generalizations. One of these uses a leap-frog numerical scheme to construct a consistent Hamiltonian evolution moncrief . The other method is implicit and uses a ”Hamiltonian” that depends on phase space variables at neighbouring time steps gp .

The leap-frog discretization scheme for the standard Hamiltonian $`p^2/2m+V(q)`$ is moncrief

$`q_{n+1}`$ $`=`$ $`q_n+ϵ{\displaystyle \frac{p_n}{m}}{\displaystyle \frac{ϵ^2}{2m}}{\displaystyle \frac{V(q_n)}{q_n}}`$
$`p_{n+1}`$ $`=`$ $`p_n{\displaystyle \frac{ϵ}{2}}\left({\displaystyle \frac{V(q_n)}{q_n}}+{\displaystyle \frac{V(q_{n+1})}{q_{n+1}}}\right).`$ (15)
This is an implicit scheme since the $`V(q_{n+1})`$ term in the $`p_{n+1}`$ equation is determined by first computing $`q_{n+1}`$. Without this term it is, up to factors, the same as the one we derived above from the general form (5). It preserves the canonical Poisson bracket under evolution.

The non-local Hamiltonian scheme gp  is rather unusual. Evolution equations are derived from a discrete action, which may be viewed as a function of coordinates $`q_n`$ and velocities $`v_n`$:

$$S[q_n,v_n]=\underset{n}{}L_n(q_n,v_n)=\underset{n}{}\left[v_n(q_{n+1}q_n)ϵH(v_n,q_n)\right].$$
(16)
This resembles a Hamiltonian form. The conjugate momentum is defined by $`p_{n+1}:=L_n/q_{n+1}=v_n`$. This leads to the following discrete equations for the canonically conjugate pair $`(q_n,p_n)`$ (Eqns. (37,43-44) of Ref. gp ):

$`q_{n+1}`$ $`=`$ $`q_n+ϵ{\displaystyle \frac{H(q_n,p_{n+1})}{p_{n+1}}}`$
$`p_{n+1}`$ $`=`$ $`p_nϵ{\displaystyle \frac{H(q_n,p_{n+1})}{q_n}}`$ (17)
It is in the sense manifested by these equations that the discretised phase space Hamiltonian is non-local in time.

The requirement of preserving the Poisson bracket for more general implicit schemes is at first sight more complicated than the procedure we followed for the explicit case. The computation of Poisson brackets requires implicit differentials assuming the validity of the implicit function theorem for the discretization. For a given discretization this leads to restrictions on the form of the Hamiltonian. It therefore appears that a general result may be difficult to obtain, and that one is stuck with a case by case analysis.

There is however a general approach for finding a large class of implicit schemes that are automatically consistent. This involves (i) viewing discrete evolution as a canonical transformation from the variables $`(q_n,p_n)`$ to the variables $`(q_{n+1},p_{n+1})`$, and (ii) requiring that the canonical transformation resemble a discretization of Hamilton’s equations. This approach gives two general classes of implicit discretization starting from a given Hamiltonian $`H(q,p)`$.

To see how this approach works it is useful to consider each of the four classes of generating functions for canonical transformations. As a first example consider a generating function of the form

$$I(q_{n+1},p_n)=q_{n+1}p_n+ϵF(q_{n+1},p_n)+ϵ^2G(q_{n+1},p_n)+\mathrm{}$$
(18)
The standard canonical transformation rules give the equations

$`q_{n+1}`$ $`=`$ $`q_n+ϵ{\displaystyle \frac{F}{p_n}}+ϵ^2{\displaystyle \frac{G}{p_n}}+\mathrm{}`$
$`p_{n+1}`$ $`=`$ $`p_nϵ{\displaystyle \frac{F}{q_{n+1}}}ϵ^2{\displaystyle \frac{G}{q_{n+1}}}+\mathrm{}.`$ (19)
The discretization prescription is to fix the functional form of $`F`$ to be that of the given Hamiltonian $`H(q_n,p_n)`$ with the replacement $`q_nq_{n+1}`$, ie set

$$F=H(q_{n+1},p_n).$$
(20)
The result is an automatically consistent class of implicit discretizations where the function $`G`$ is arbitrary. The limit $`ϵ0`$ gives the continuous time Hamilton equations.

The second example of this type is provided by the generating function

$$I(q_n,p_{n+1})=q_np_{n+1}+ϵF(q_n,p_{n+1})+ϵ^2G(q_n,p_{n+1})+\mathrm{}$$
(21)
The transformation equations are now

$`q_{n+1}`$ $`=`$ $`q_n+ϵ{\displaystyle \frac{F}{p_{n+1}}}+ϵ^2{\displaystyle \frac{G}{p_{n+1}}}+\mathrm{}`$
$`p_{n+1}`$ $`=`$ $`p_nϵ{\displaystyle \frac{F}{q_n}}ϵ^2{\displaystyle \frac{G}{q_n}}+\mathrm{}`$ (22)
These again give the continuum Hamilton equations in the $`ϵ0`$ limit with $`F=H(q_n,p_{n+1})`$ for arbitrary $`G`$. (A special case of this scheme arises from the action viewpoint discussed in gp .)

The remaining two types of transformations with generating functions depending on the pairs $`(q_n,q_{n+1})`$ and $`(p_n,p_{n+1})`$ do not give a form resembling discretized Hamilton equations, so the two cases discussed above appear to exhaust this approach for generating consistent implicit schemes.

As for the explicit scheme of the last section, it is readily verified that constraints are not preserved generically with these two types of implicit discretization. For the parametrized particle (9) for example, the evolution equations with the first method (19) (without the $`ϵ^2`$ terms, and with lapse $`N`$) are

$`q_{n+1}`$ $`=`$ $`q_n+Nϵ{\displaystyle \frac{p_n}{m}}.`$
$`p_{n+1}`$ $`=`$ $`p_nNϵ{\displaystyle \frac{V(q_{n+1})}{q_{n+1}}}.`$
$`p_{n+1}^t`$ $`=`$ $`p_n^t`$
$`t_{n+1}`$ $`=`$ $`t_n+Nϵ.`$ (23)
Substitution of these into $`H_n=p_n^t+p_n^2/2m+V(q_n)`$ shows that $`H_{n+1}H_n`$.

In summary, we have seen that both the explicit and implicit schemes are quite restrictive in their respective settings, although with the latter applicable to any Hamiltonian. In particular, the schemes we have considered are two-step, involving only variables at times $`n`$ and $`n+1`$. Multi-step general schemes cannot be represented as canonical transformations in any simple way for the implicit case, and checking of Poisson bracket preservation appears to be much more tedious for the explicit case. It is also apparent from the simple example of the parametrized particle that constraints are not preserved under discrete time evolution.

## III Applications

It is useful to see whether there are viable implementations of the explicit and implicit discretization discussed in the last section. In the following we consider a few examples of models with Hamiltonian constraints to see if there are solutions for the functions $`\alpha ,\beta `$ in (5) that preserve the canonical Poisson bracket. To do this we follow the procedure described above. We also examine the time gauge fixing conditions introduced by the implicit discretizations (19-22).

### III.1 The parametrized scalar field

As a first application we look at the parametrized scalar field in a curved spacetime. Its Hamiltonian is

$$H=N\left[p^t+\frac{1}{2\sqrt{q}}\mathrm{\Pi }^2+\frac{1}{2}\sqrt{q}q^{ab}_a\varphi _b\varphi +\sqrt{q}V(\varphi )\right]$$
(24)
where the canonical coordinates are $`(t,p^t)`$ and $`(\varphi ,\mathrm{\Pi })`$, and $`q_{ab}`$ is the spatial background metric. The discrete evolution equations are

$`\varphi _{n+1}`$ $`=`$ $`\varphi _n+ϵN{\displaystyle \frac{\mathrm{\Pi }}{\sqrt{q}}}+ϵ^2\alpha (\varphi _n,\mathrm{\Pi }_n)`$
$`\mathrm{\Pi }_{n+1}`$ $`=`$ $`\mathrm{\Pi }_n+ϵ\left[_a\left(\sqrt{q}q^{ab}_b\varphi N\right)V^{}N\sqrt{q}\right]+ϵ^2\beta (\varphi _n,\mathrm{\Pi }_n)`$
$`t_{n+1}`$ $`=`$ $`t_n+Nϵ`$
$`p_{n+1}^t`$ $`=`$ $`p_n^t.`$ (25)
From this one calculates $`\{\varphi _{n+1},\mathrm{\Pi }_{n+1}\}`$ to try to find functions $`\alpha `$ and $`\beta `$ such that $`\{\varphi _{n+1},\mathrm{\Pi }_{n+1}\}=1`$. It turns out that this is not possible except for the special case $`\varphi =\varphi (t)`$, where it reduces to the case of the parametrized particle (9).

### III.2 Cosmological models

We consider two models here, FRW coupled to a scalar field and deSitter space. These serve to illustrate the main problem with discretization arising from the Hamiltonian constraint of general relativity: the kinetic term contains a mixture of gravitational configuration and momentum variables, which makes it difficult to find a discretization that preserves the canonical commutation rules under evolution for the explicit scheme.

For both models we use geometrodynamics phase space variables $`(q_{ab},\stackrel{~}{\pi }^{ab})`$, and the parametrization

$$q_{ab}=a(t)e_{ab}\stackrel{~}{\pi }^{ab}=\frac{1}{3}P(t)e^{ab}$$
(26)
where $`e_{ab}`$ is the flat Euclidean $`3`$metric. (The main features of the analysis are unchanged for the triad-connection phase space variables.)

The Hamiltonian constraint for flat FRW minimally coupled to a scalar field is

$$H(a,p)=\frac{1}{6}a^{1/2}p^2+\frac{1}{2}a^{3/2}\mathrm{\Pi }^2a^{3/2}V(\varphi ),$$
(27)
where $`V(\varphi )`$ is a potential for the scalar field $`\varphi `$, and $`\mathrm{\Pi }`$ is its canonical momentum. For De Sitter space, the potential is replaced by the cosmological constant $`\mathrm{\Lambda }`$.

The discrete evolution equations obtained using (5) are

$$a_{n+1}=a_nϵN\frac{1}{3}a_n^{1/2}p_n+ϵ^2\alpha (a_n,p_n)$$
(28)
for both the FRW and deSitter case. The momentum evolution equation is

$$p_{n+1}=p_n+ϵN\left[\frac{1}{12}a_n^{1/2}p_n^2+\frac{3}{2}a_n^{5/2}\mathrm{\Pi }^2\frac{3}{2}a_n^{1/2}V(\varphi )\right]+ϵ^2\beta (a_n,p_n)$$
(29)
for the FRW case, and

$$p_{n+1}=p_n+ϵN\left[\frac{1}{12}a_n^{1/2}p_n^2+\frac{3}{2}a_n^{5/2}\mathrm{\Lambda }\right]+ϵ^2\beta (a_n,p_n)$$
(30)
for deSitter. (We ignore the evolution equations for the scalar field, since these are not relevant for the point we wish to make.) The main observation again is that computing $`\{a_{n+1},p_{n+1}\}`$ shows that there is no choice for $`\alpha (a_n,p_n)`$ and $`\beta (a_n,p_n)`$ which preserves the Poisson bracket to all non-vanishing orders in $`ϵ`$.

Thus we conclude that this general form of discretization is not applicable to models derived from general relativity.

### III.3 $`BF`$ theory

This example of a fully constrained dynamics is of particular interest in the context of spin foam models. Once again we focus on time discretization, which is what is relevant for preserving the Poisson bracket under evolution.

Abelian case: the continuum Hamiltonian is given by

$$H=\left(\lambda ϵ^{ab}F_{ab}+\mu _aE^a\right)d^3x$$
(31)
with the two constraints $`F_{ab}=0`$ and $`_aE^a=0`$ and Lagrange multipliers $`\lambda `$ and $`\mu `$. The discrete evolution equations for the canonical variables $`A_a`$ and $`E^a`$ are

$`A_{a,n+1}`$ $`=`$ $`A_{a,n}ϵ_a\mu _n,`$ (32)
$`E_{n+1}^a`$ $`=`$ $`E_n^a+ϵ2ϵ^{ab}_b\lambda _n.`$ (33)
From these it is obvious that the Poisson brackets are trivially preserved, even without considering terms of order $`ϵ^2`$.

Non-Abelian case: The only difference is that the canonical variables are now Lie-algebra valued. Thus the Hamiltonian is

$$H=(\lambda ^iϵ^{ab}F_{ab}^i+\mu ^iD_aE^{ai})d^3x$$
(34)
where $`D_a=_a+A_a`$ is the covariant derivative. As a result, the discrete evolution now becomes non-trivial:

$`A_{a,n+1}^i`$ $`=`$ $`A_{a,n}^iϵD_a\mu _n^i+ϵ^2\alpha (A_{a,n}^i,E_n^{ai})`$ (35)
$`E_{n+1}^{ai}`$ $`=`$ $`E_n^{ai}+ϵ\left[2ϵ^{ab}D_b\lambda _n^i+ϵ^{ijk}E_n^{aj}\mu _n^k\right]+ϵ^2\beta (A_{a,n}^i,E_n^{ai}).`$ (36)
As a consequence of the additional terms compared to the Abelian case, it is straightforward to check that there is no solution for $`\alpha `$ and $`\beta `$ that leads to a preservation of the Poisson bracket.

For the implicit schemes the Poisson bracket is of course preserved by the general argument given in the last section. The constraint associated with time evolution in BF theory is a combination of the two ”internal symmetry” constraints, as occurs for example in 2+1 gravity bengt . Thus it is expected that there is some partial breaking of these symmetries due to time discretization. Given this, it is useful to see what are the (partial) gauge fixing conditions induced by the time discretization.

Denoting the continuum gauge conditions by $`f^i(A,E)=0`$ and $`g^i(A,E)=0`$, the relation between these and the lagrange multipliers $`\lambda ^i`$ and $`\mu ^i`$ are given by the requirement that the gauges are preserved in time. Schematically this is

$$\{f,F[\lambda ]+G[\mu ]\}=0,\{g,F[\lambda ]+G[\mu ]\}=0.$$
(37)
For the discretization, one can in principle extract the functions $`\lambda ^i`$ and $`\mu ^i`$ by setting the time evolved constraints to zero, ie. $`G_{n+1}^iG^i(A_{n+1},E_{n+1})=0`$ and $`F_{n+1}^iF^i(A_{n+1},E_{n+1})=0`$. One can then attempt to solve Eqn. (37) for the gauge conditions $`f`$ and $`g`$. To get an idea of what the resulting equations look like, the first task is to extract $`\lambda ^i`$ and $`\mu ^i`$ from the evolved constraints. For the first implicit scheme (19) for example, where $`H=H[A_{n+1},E_n]`$, these equations are

$`G_{n+1}^i`$ $`=`$ $`G^i+ϵ\left[{\displaystyle \frac{\delta H}{\delta A_{n+1}^i}}+ϵ^{ijk}\left(E_n^{ak}{\displaystyle \frac{\delta H}{\delta E_n^j}}A_{n+1}^j{\displaystyle \frac{\delta H}{\delta A_{n+1}^k}}\right)\right]ϵ^2ϵ^{ijk}{\displaystyle \frac{\delta H}{\delta E_n^j}}{\displaystyle \frac{\delta H}{\delta A_{n+1}^k}}=0`$ (38)
$`F_{n+1}^i`$ $`=`$ $`F_n^i+ϵ\left[{\displaystyle \frac{\delta H}{\delta E_n^i}}+2ϵ^{ijk}A_n^j{\displaystyle \frac{\delta H}{\delta E_n^k}}\right]+ϵ^2ϵ^{ijk}{\displaystyle \frac{\delta H}{\delta E_n^j}}{\displaystyle \frac{\delta H}{\delta E_n^k}}=0`$ (39)
up to order $`ϵ`$, where we have suppressed $`ϵ^{ab}`$ and the spatial indices ($`a,b\mathrm{}`$). Inserting

$$\frac{\delta H}{\delta E_n^i}=ϵ^{ijk}A_{n+1}^j\mu _n^k$$
(40)
and

$$\frac{\delta H}{\delta A_{n+1}^i}=2\left(\lambda _n^i+ϵ^{ijk}A_{n+1}^j\lambda _n^k\right)$$
(41)
into the equations above, one can in principle solve for $`\lambda `$ and $`\mu `$, and substitute in (37) to get the functions $`f`$ and $`g`$. These equations give an indication of the rather complex gauge fixing induced by this discretization.

Note that although we are dealing here with time discretization only, this may still affect the first class nature of the purely internal symmetry constraints: the first class property is preserved automatically under evolution if $`G_{n+1}^i=0`$ follows from $`G_n^i=0`$, without fixing any Lagrange multipliers. We have seen above that this is not the case. Therefore, although time discertization does not affect the internal symmetry algebra at fixed time ($`n`$), it nevertheless makes the constraints second class with respect to evolution in time.

From the above examples, the pattern of calculation is clear for application to quantum gravity. It is unlikely that explicit schemes will work due to non-conservation of the Poisson bracket. The two classes of implicit schemes can be pursued to some extent at least, in the connection-triad canonical variables. A first attempt might be to keep space continuous and discretize time for the general form of the constraints with the Barbero-Immirzi parameter. It is possible there are choices of this parameter that simplify the form of the constraints and gauge fixing conditions.

## IV Discussion

We have studied general classes of explicit and implicit time discretizations of Hamilton’s equations containing arbitrary functions. The general forms were motivated by simplicity, and relative ease in establishing the primary criterion that the Poisson bracket is preserved under discrete evolution.

The explicit class of schemes we studied preserves the Poisson bracket for only a limited class of Hamiltonians. Therefore its use for quantum gravity is rather limited. The two classes of implicit schemes have more potential utility in that these apply to any Hamiltonian. However the resulting gauge conditions are quite complicated.

Our approach appears to exhaust explicit and implicit two-step discretization schemes that preserve the Poisson bracket. In both cases, continuum first class constraints are not preserved in time. Therefore gauge fixing conditions are implicit in the discretization, and can in principle be explicitly recovered. Given this one can ask what is the utility of this lattice approach to gauge fixing, versus the alternative path of gauge fixing in the continuum theory followed by discretization (where preservation of constraints is obviously not an issue). For example, the problem of time is ”solved” by gauge fixing in both approaches. Furthermore, at the level of quantization, the old question of whether different gauge fixings lead to the same quantum gravity model would reemerge in a lattice disguise.

Acknowledgements We thank Thomas Thiemann for discussions. This work was supported in part by the Natural Science and Engineering Research Council of Canada.