text:
# On the Topology of the Space of Contact Structures on Torus Bundles

## 1 Introduction

Very little seems to be known about the topology of the space $`\mathrm{\Xi }(M)`$ of contact structures on a given manifold $`M`$. (There is, as yet, no standard notation for that space. Since $`\xi `$ is a customary notation for contact structures, the letter $`\mathrm{\Xi }`$ seems an apt choice.) Some general results about the structure of $`\mathrm{\Xi }`$ can be found in \[10, Section 2.4\] and  (the latter has to be read with a certain amount of caution; see Mathematical Reviews 83k:58018). For our purposes we only need to observe that $`\mathrm{\Xi }(M)`$ is an open subset (in the $`C^k`$-topology for $`k1`$) of the Grassmannian manifold of all codimension one subbundles of the tangent bundle $`TM`$.

There also appear to be few results concerning $`\mathrm{\Xi }(M)`$ for specific $`M`$. One such result is due to Eliashberg. Let $`\xi _0=\mathrm{ker}(xdyydx+zdttdz)`$ be the standard contact structure on $`S^3^4`$ and denote by $`\mathrm{\Xi }_0(S^3,p_0)`$ the subspace of $`\mathrm{\Xi }(S^3)`$ which consists of contact structures isotopic to $`\xi _0`$ and with $`\xi (p_0)=\xi _0(p_0)`$ for some fixed $`p_0S^3`$. In  Eliashberg has shown that $`\mathrm{\Xi }_0(S^3,p_0)`$ is contractible.

This result can be rephrased as follows. Write $`\text{Diff}_{p_0}^+(S^3)`$ for the group of orientation preserving diffeomorphisms of $`S^3`$ that fix the contact plane $`\xi _0(p_0)`$, and denote by $`\text{Cont}(S^3,\xi _0)`$ the subgroup of $`\text{Diff}_{p_0}^+(S^3)`$ consisting of diffeomorphisms which preserve $`\xi _0`$. Then Eliashberg’s theorem says that the inclusion of $`\text{Cont}(S^3,\xi _0)`$ in $`\text{Diff}_{p_0}^+(S^3)`$ is a homotopy equivalence, cf. Section 3. Giroux  gives results about the connected components of other such contactomorphism groups.

The present paper is concerned with the topology of the space of contact structures on $`T^2`$-bundles over $`S^1`$. Specifically, we provide an essentially elementary argument showing the fundamental group of these spaces (with base point specified below) to contain an infinite cyclic subgroup. The basic techniques of the proof are classical (Gray stability and contact Hamiltonians), but a key ingredient is the recent work of Giroux  and Kanda  on the classification of contact structures on the $`3`$-torus. Throughout we assume contact structures to be coorientable, so they can be defined by global contact forms.

For results about related homotopical questions in symplectic topology see for instance  and .

## 2 Contact geometry of the $`3`$-torus

Let $`(x,y,\theta )`$ be standard angular coordinates on $`T^3=(/2\pi )^3`$. Define, for $`n`$,

$$\alpha _s^n=\mathrm{cos}(n\theta 2\pi s)dx\mathrm{sin}(n\theta 2\pi s)dy,s[0,1],$$

and $`\zeta _s^n=\mathrm{ker}\alpha _s^n`$. We shall abbreviate $`\zeta _s^1`$ to $`\zeta _s`$.

Giroux and Kanda have shown that the $`\zeta _0^n`$ provide a complete list, up to diffeomorphism, of so-called tight contact structures on the 3-torus (see  for the definition of ‘tight’). Notice, however, that the isotopy classification is more subtle. According to Eliashberg and Polterovich , for $`\varphi \text{SL}(3,)\text{Diff}^+(T^3)`$ the contact structures $`\varphi _{}\zeta _0`$ and $`\zeta _0`$ are isotopic if and only if $`\varphi `$ fixes the subspace $`^20`$ of the first homology group. Taken together, these results yield complete information about the set of connected components of $`\mathrm{\Xi }(T^3)`$, it being known that $`\pi _0(\text{Diff}^+(T^3))=\text{SL}(3,)`$.

The following is our main result about the topology of $`\mathrm{\Xi }(T^3)`$. Below we shall formulate corresponding statements for all other torus bundles over $`S^1`$.

###### Proposition 1.

For each $`n`$ the fundamental group $`\pi _1(\mathrm{\Xi }(T^3),\zeta _0^n)`$ based at $`\zeta _0^n`$ contains an infinite cyclic subgroup, generated by the loop $`\{\zeta _s^n:\mathrm{\hspace{0.17em}0}s1\}`$.

Together with the long homotopy exact sequence of the next section and the known results about the homotopy type of $`\text{Diff}(T^3)`$ this can be related to the topology of contactomorphism groups.

Notice that the loop $`\{\xi _s\}`$ is defined by a linear circle of contact forms,

$$\xi _s=\mathrm{ker}(\mathrm{cos}(2\pi s)\alpha _0+\mathrm{sin}(2\pi s)\alpha _{1/4}).$$

A pair of contact forms $`(\alpha _0,\alpha _{1/4})`$ with this property that any non-trivial linear combination of these forms is again a contact form was called a contact circle in , and it is natural to ask what can be said about such contact circles as elements of $`\pi _1(\mathrm{\Xi }(M))`$. In  it is shown that contact circles exist on all closed, orientable 3-manifolds.

Our proof of Proposition 1 makes essential use of the diffeomorphism classification of contact structures on $`T^3`$. It is worth emphasising that the analogue of Proposition 1 for more general torus bundles over the circle, formulated in Section 5 below, rests equally on the classification of contact structures on $`T^3`$, and not on the classification of contact structures on these spaces themselves. Again the essential loop is a contact circle.

Clearly, if a contact circle $`(\beta _1,\beta _2)`$ extends to a contact sphere $`(\beta _1,\beta _2,\beta _3)`$, i.e. a triple of contact forms such that any non-trivial linear combination is a contact form, then $`(\beta _1,\beta _2)`$ defines the trivial element in $`\pi _1(\mathrm{\Xi }(M))`$.

On the other hand, evaluation at a point $`pM`$ and identification of cooriented radial 2-planes in $`^3`$ with $`S^2`$ defines a diffeomorphism

$$\begin{array}{ccc}\hfill ^3S^2& & S^2\hfill \\ \hfill (\lambda _1,\lambda _2,\lambda _3)& & \mathrm{ker}(\lambda _1\alpha _1+\lambda _2\alpha _2+\lambda _3\alpha _3)(p),\hfill \end{array}$$

so a contact sphere always maps to the generator of $`\pi _2(S^2)`$ under this evaluation map and hence defines an element of infinite order in $`\pi _2(\mathrm{\Xi }(M))`$. Contact spheres exist, for instance, on the connected sum of left-quotients of $`\text{SU}(2)`$ and copies of $`S^1\times S^2`$, see \[7, Prop. 5.7\]. This does not contradict Eliashberg’s result about the contractibility of $`\mathrm{\Xi }_0(S^3,p_0)`$, of course, since in $`\mathrm{\Xi }(S^3)`$ we do not fix the contact plane at a point.

## 3 A homotopy exact sequence

The following considerations rest on the well-known concepts of Gray stability and contact Hamiltonians, cf. , . We only recall briefly the facts that we shall use:

Gray stability: Given a smooth family of contact structures $`\xi _t`$, $`t[0,1]`$, on a closed manifold $`M`$, there is a canonically defined time-dependent vector field $`X_t`$ tangent to $`\xi _t`$ whose flow $`\psi _t`$ satisfies $`\psi _t\xi _0=\xi _t`$ (that is, the differential $`T\psi _t`$ maps $`\xi _0`$ to $`\xi _t`$). Given a family of contact forms $`\alpha _t`$ defining $`\xi _t`$, with Reeb vector field $`R_t`$, this vector field is determined by the equations $`\alpha _t(X_t)=0`$ and $`i_{X_t}d\alpha _t=\mu _t\alpha _t\dot{\alpha }_t`$, with the function $`\mu _t`$ determined by $`\mu _t=\dot{\alpha }_t(R_t)`$, and with $`\dot{\alpha }_t`$ denoting the time derivative of $`\alpha _t`$. This $`X_t`$ is independent of the choice of $`\alpha _t`$, even though $`\mu _t`$ does depend on that choice.

Contact Hamiltonians: Given a smooth family of functions $`H_t`$ on a contact manifold $`(M,\xi =\mathrm{ker}\alpha )`$, a time-dependent vector field $`X_t`$ whose flow preserves $`\xi `$ is defined by

$$X_t=H_tR+Y_t,$$

where $`R`$ is the Reeb vector field of $`\alpha `$ and $`Y_t`$ is determined by

$$\alpha (Y_t)=0$$

and

$$i_{Y_t}d\alpha =dH_t(R)\alpha dH_t.$$

Indeed, the infinitesimal condition $`L_{X_t}\alpha =\rho _t\alpha `$ is equivalent to the equations above if we require $`H_t=\alpha (X_t)`$. Put differently, a time-dependent vector field $`X_t`$ whose flow preserves $`\xi =\mathrm{ker}\alpha `$ is completely determined by $`H_t=\alpha (X_t)`$, the Hamiltonian function of $`X_t`$.

Let $`(M,\xi _0)`$ be a closed contact manifold. Write $`\text{Diff}_0(M)`$ for the identity component of the diffeomorphism group of $`M`$ and $`\text{Cont}_0(M,\xi _0)`$ for its subgroup of contactomorphisms, i.e.

$$\text{Cont}_0(M,\xi _0)=\{\varphi \text{Diff}_0(M):\varphi _{}\xi _0=\xi _0\}.$$

We have an obvious inclusion $`i:\text{Cont}_0(M,\xi _0)\text{Diff}_0(M)`$ and, by Gray stability, a surjection

$$\begin{array}{ccc}\hfill \sigma :\text{Diff}_0(M)& & \mathrm{\Xi }_0(M)\hfill \\ \hfill \varphi & & \varphi _{}\xi _0,\hfill \end{array}$$

where $`\mathrm{\Xi }_0(M)`$ denotes the component of $`\mathrm{\Xi }(M)`$ containing $`\xi _0`$.

The following is a well-known (folklore) result:

###### Proposition 2.

The map $`\sigma `$ has the homotopy lifting property with respect to cubes $`I^n=[0,1]^n`$ and smooth homotopies. Hence there is an exact sequence of homotopy groups

$$\mathrm{}\stackrel{\mathrm{\Delta }}{}\pi _k(\text{Cont}_0)\stackrel{i_\mathrm{\#}}{}\pi _k(\text{Diff}_0)\stackrel{\sigma _\mathrm{\#}}{}\pi _k(\mathrm{\Xi }_0)\stackrel{\mathrm{\Delta }}{}\pi _{k1}(\text{Cont}_0)\stackrel{i_\mathrm{\#}}{}\mathrm{}$$

Proof. We only need to prove the first statement. The usual proof that a Serre fibration yields an exact sequence of homotopy groups then still applies to give the desired conclusion.

This means we are given a commutative diagram

$$\begin{array}{ccc}I^n\times 0& & \text{Diff}_0(M)\\ & & \sigma & & \\ I^n\times I& & \mathrm{\Xi }_0(M),\end{array}$$

with the bottom map smooth in the $`I`$-factor; more precisely, a family $`\xi _{s,t}`$ of contact structures on $`M`$, continuous in $`sI^n`$, smooth in $`tI`$, and with all $`t`$-derivatives continuous in $`s`$. We also have a family $`\varphi _{s,0}\text{Diff}_0(M)`$ with $`(\varphi _{s,0})_{}\xi _0=\xi _{s,0}`$, and the aim is to find a lifting $`I^n\times I\text{Diff}_0(M)`$.

By Gray stability (with an added parameter $`s`$) one finds a family of vector fields $`X_{s,t}`$, continuous in $`s`$ and smooth in $`t`$, such that the flow of the $`t`$-dependent vector field $`X_{s,t}`$ defines a family of diffeotopies $`\psi _{s,t}`$ which is continuous in $`s`$ and satisfies $`(\psi _{s,t})_{}\xi _{s,0}=\xi _{s,t}`$. Set $`\varphi _{s,t}=\psi _{s,t}\varphi _{s,0}`$. Then $`(\varphi _{s,t})_{}\xi _0=\xi _{s,t}`$, so $`\varphi _{s,t}`$ is the desired lift of $`\xi _{s,t}`$. $`\mathrm{}`$

## 4 Proof of Proposition 1

We only show that $`\{\zeta _s=\zeta _s^1\}`$ is non-zero in $`\pi _1`$; a completely analogous argument applies to any multiple of this loop and to general $`n`$.

Our proof is by contradiction. Assume $`\{\zeta _s\}`$ defines the trivial element in $`\pi _1(\mathrm{\Xi }(T^3),\zeta _0)`$. Write $`\pi `$ for the canonical submersion

$$T^2\times T^2\times /2\pi T^3$$

and set $`\stackrel{~}{\zeta }_s=\pi ^{}\zeta _s`$ for each $`s[0,1]`$. Then $`\{\stackrel{~}{\zeta }_s\}`$ would define the trivial element in $`\pi _1(\mathrm{\Xi }_0^{}(T^2\times ))`$, where $`\mathrm{\Xi }_0^{}(T^2\times )`$ denotes the connected component containing $`\stackrel{~}{\zeta }_0`$ of the space of contact structures on $`T^2\times `$ which are invariant under shifts by $`2\pi `$ in $``$-direction.

Proposition 2 applies to $`T^2\times `$ if we replace $`\mathrm{\Xi }_0`$ by $`\mathrm{\Xi }_0^{}(T^2\times )`$ and $`\text{Diff}_0`$ by $`\text{Diff}_0^{}(T^2\times )`$, the identity component of the $``$-equivariant diffeomorphisms of $`T^2\times `$, i.e.

$$\text{Diff}_0^{}(T^2\times )=\{\varphi \text{Diff}(T^2\times ):\varphi (x,y,\theta +2\pi )=\varphi (x,y,\theta )+(0,0,2\pi )\}_0,$$

for the $``$-invariance resp. -equivariance guarantees that the vector field $`X_{s,t}`$ used in the proof of Proposition 2 still integrates to a global flow.

We would then have a homotopy $`\{\stackrel{~}{\zeta }_{s,t}\}`$ rel $`\{0,1\}`$ of $`\{\stackrel{~}{\zeta }_s\}`$ to the constant loop $`\{\stackrel{~}{\zeta }_0\}`$. We may assume that $`\stackrel{~}{\zeta }_{s,t}`$ is smooth in $`t`$ (by smoothing a homotopy which is constant near $`t=0,1`$).

Define a diffeomorphism

$$\begin{array}{ccc}\hfill \stackrel{~}{\varphi }_s:T^2\times & & T^2\times \hfill \\ \hfill (x,y,\theta )& & (x,y,\theta +2\pi s),\hfill \end{array}$$

so that $`\stackrel{~}{\varphi }_s\stackrel{~}{\zeta }_0=\stackrel{~}{\zeta }_{s,0}=\stackrel{~}{\zeta }_0`$. By Proposition 2 we find a path $`\psi (s)`$ in $`\text{Cont}_0(T^2\times ,\stackrel{~}{\zeta }_0)`$ that joins $`\stackrel{~}{\varphi }_0`$ with $`\stackrel{~}{\varphi }_1`$.

Let $`X_s`$ be the time-dependent vector field that integrates to $`\psi (s)`$ at time $`s`$. Choose a contact form $`\stackrel{~}{\alpha }_0`$ representing $`\stackrel{~}{\zeta }_0`$ and set $`H_s=\stackrel{~}{\alpha }_0(X_s)`$. Since the flow of $`X_s`$ preserves $`\stackrel{~}{\zeta }_0`$, the Hamiltonian vector field of $`H_s`$ coincides with $`X_s`$. The time-one-map of this Hamiltonian flow is $`\stackrel{~}{\varphi }_1`$. During the diffeotopy, the image of $`T^2\times [0,2\pi ]`$ under $`\psi (s)`$ stays above $`T^2\times \{2\pi (k1)\}`$ for some $`k`$ sufficiently large.

Now let $`\chi (\theta )`$ be a smooth function with $`\chi (\theta )0`$ for $`\theta 2\pi k`$ and $`\chi (\theta )1`$ for $`\theta 2\pi (k1)`$, and set $`\overline{H}_s(x,y,\theta )=H_s(x,y,\theta )\chi (\theta )`$. The Hamiltonian vector field $`\overline{X}_s`$ of $`\overline{H}_s`$ still integrates to a global flow, and the time-one-map of this flow defines a $`\stackrel{~}{\zeta }_0`$-preserving diffeomorphism

$$T^2\times [2\pi k,0]T^2\times [2\pi k,2\pi ].$$

But this would imply that the contact structures $`\zeta ^k`$ and $`\zeta ^{k+1}`$ on $`T^3`$ are diffeomorphic, contradicting the work of Giroux and Kanda. This contradiction proves the proposition. $`\mathrm{}`$

Here is an alternative proof of Proposition 1, based on the work of Eliashberg-Hofer-Salamon . Again arguing by contradiction, we assume that $`\{\zeta _s\}`$ is trivial in $`\pi _1(\mathrm{\Xi }(T^3),\zeta _0)`$. Consider the double cover

$$T_{4\pi }^3T^2\times /4\pi T^2\times /2\pi T^3,$$

and write $`\zeta _s^{}`$ for the lift of $`\zeta _s`$. Then $`\{\zeta _s^{}\}`$ is also a contractible loop in the corresponding space of contact structures.

The same argument as before yields a contact isotopy of $`(T_{4\pi }^3,\zeta _0^{})`$, beginning at the identity and ending with a shift by $`2\pi `$ in $`\theta `$-direction. In particular, this contact isotopy would move the pre-Lagrangian torus (see ) $`T^2\times \{0\}`$ to $`T^2\times \{2\pi \}`$ and therefore separate it completely from itself. This would contradict \[3, Thm. 3.8.3\], which uses a Floer homology argument to provide a lower bound for the number of intersection points between a pre-Lagrangian submanifold $`\mathrm{\Lambda }_0`$ and the image $`\mathrm{\Lambda }_1^{}`$ of a Legendrian submanifold $`\mathrm{\Lambda }_1\mathrm{\Lambda }_0`$ under a contact isotopy (The conditions $`O_1`$ and $`O_2`$ of that theorem are satisfied in the present context).

## 5 Contact circles on torus bundles

We now extend Proposition 1 to all torus bundles over the circle. For $`A\text{SL}_2`$, write $`M_A`$ for the quotient of $`T^2\times `$ under the transformation

$$\stackrel{~}{A}:(\left(\begin{array}{cc}x& \\ y& \end{array}\right),\theta )(A\left(\begin{array}{cc}x& \\ y& \end{array}\right),\theta +2\pi ).$$

On each $`M_A`$ there is a family of contact structures $`\zeta ^n`$, $`n`$, characterised up to fibre preserving isotopy by the property that it descends from a contact structure on $`^2\times `$ of the form

$$\mathrm{cos}f(\theta )dx\mathrm{sin}f(\theta )dy=0,$$

invariant under the transformation $`\stackrel{~}{A}`$, and with

$$2(n1)\pi <\underset{\theta }{sup}(f(\theta +1)f(\theta ))2n\pi ,$$

see .

###### Proposition 3.

For each $`A\text{SL}_2`$ and $`n`$, the fundamental group of $`\mathrm{\Xi }(M_A)`$ based at $`\zeta ^n`$ contains an element of infinite order.

To prove this proposition it is convenient to use the Thurston geometries adapted to the different choices of $`A`$, as done previously in . Rather than giving a case by case proof, we illustrate the method by a particular example which contains all the ideas required for the general argument.

Assume that $`\text{trace}A3`$. Then $`A`$ is conjugate in $`\text{GL}_2`$ to a matrix of the form

$$A^{}=\left(\begin{array}{cc}e^\gamma & 0\\ 0& e^\gamma \end{array}\right)$$

and $`M_A`$ is modelled on the solvable Lie group $`\text{Sol}^3`$ (the inhomogeneous Lorentz group), cf. . This means that we can write $`M_A`$ as a quotient of $`^3`$ under a group of transformations $`\mathrm{\Gamma }`$ generated by

$`(x,y,\theta )`$ $``$ $`(x+a_1,y+b_1,\theta )`$
$`(x,y,\theta )`$ $``$ $`(x+a_2,y+b_2,\theta )`$
$`(x,y,\theta )`$ $``$ $`(e^\gamma x,e^\gamma y,\theta +\gamma ),`$
where the vectors

$$\left(\begin{array}{cc}a_1& \\ b_1& \end{array}\right),\left(\begin{array}{cc}a_2& \\ b_2& \end{array}\right)$$

generate a lattice $`\mathrm{\Lambda }`$ in $`^2`$ invariant under $`A^{}`$.

Now define, for $`n`$,

$$\alpha _s^n=\mathrm{cos}\left(\frac{2\pi n\theta }{\gamma }2\pi s\right)e^\theta dx\mathrm{sin}\left(\frac{2\pi n\theta }{\gamma }2\pi s\right)e^\theta dy,s[0,1],$$

and set $`\zeta _s^n=\mathrm{ker}\alpha _s^n`$. This defines a contact structure on $`^3`$ that descends to $`M_A`$. Up to fibre preserving isotopy, $`\zeta _0^n`$ is equal to the $`\zeta ^n`$ studied by Giroux. Write $`\stackrel{~}{\zeta }_0`$ for the lift of $`\zeta _0=\zeta _0^1`$ to $`T^2\times `$.

Denote by $`\mathrm{\Xi }_0^A^{}(T^2\times )`$ the connected component of $`\stackrel{~}{\zeta }_0`$ in the space of contact structures on $`T^2\times `$ invariant under the transformation

$$\stackrel{~}{A}^{}:(x,y,\theta )(e^\gamma x,e^\gamma y,\theta +\gamma ),$$

where we think of $`T^2`$ as the quotient of $`^2`$ under the lattice $`\mathrm{\Lambda }`$. Write $`\text{Diff}_0^A^{}(T^2\times )`$ for the identity component of the diffeomorphisms of $`T^2\times `$ that commute with $`\stackrel{~}{A}^{}`$.

Then the argument in the proof of Proposition 1 goes through as before, with $`\mathrm{\Xi }_0^{}(T^2\times )`$ replaced by $`\mathrm{\Xi }_0^A^{}(T^2\times )`$, and $`\text{Diff}_0^{}(T^2\times )`$ by $`\text{Diff}_0^A^{}(T^2\times )`$. Assuming that $`\{\zeta _s\}`$ was trivial in $`\pi _1(\mathrm{\Xi }(M_A),\zeta _0)`$, we obtain a $`\stackrel{~}{\zeta }_0`$-preserving diffeomorphism

$$T^2\times [k\gamma ,0]T^2\times [k\gamma ,\gamma ].$$

On $`T^2\times \{l\gamma \}`$, $`l`$, the contact form $`\alpha _0^1`$ restricts to $`e^{l\gamma }dx`$, so the characteristic foliation on these tori is always given by $`dx=0`$. Since the characteristic foliation on a surface determines the germ of the contact structure along that surface, we can use the identity map on $`T^2`$ to glue the ends of $`T^2\times [k\gamma ,\epsilon \gamma ]`$, with $`\epsilon \{0,1\}`$, and obtain a contact structure on $`T^3`$, diffeomorphic to the standard structure $`\zeta ^k`$ of $`T^3`$ for $`\epsilon =0`$, and $`\zeta ^{k+1}`$ for $`\epsilon =1`$. The diffeomorphism above would induce a diffeomorphism between $`\zeta ^k`$ and $`\zeta ^{k+1}`$ on $`T^3`$, contradicting once again the work of Giroux and Kanda.

Hansjörg Geiges

Mathematisches Institut

Universität zu Köln

Weyertal 86–90

50931 Köln

Germany

e-mail: geiges@math.uni-koeln.de

Jesús Gonzalo

Departamento de Matemáticas

Universidad Autónoma de Madrid

28049 Madrid

Spain

e-mail: jesus.gonzalo@uam.es


format: text

text:
# Statistical Study of the Blue Straggler Properties in Galactic Globular Clusters

## 1. The Relative Frequency of Blue Stragglers

Our blue straggler star (BSS) catalog includes almost 3000 objects extracted from 62 Galactic globular cluster (GC) CMDs which are part of the photometrically homogeneous HST snapshot database (Piotto et al. 2002).

In each cluster, the number of BSS has been normalized to the total population of stars. We defined three different specific frequencies ($`F_{\mathrm{BSS}}`$) as the ratios of the BSS counts to the number of HB, RGB, and HB plus RGB stars, respectively. The results described below turned out to be independent from the adopted $`F_{\mathrm{BSS}}`$.

Figure 1 (upper panel) shows that there is a significant anticorrelation between $`F_{\mathrm{BSS}}`$ and the absolute total magnitude $`M_V`$ of the cluster. In particular, the faintest clusters in our sample (NGC 6717 and NGC 6838) present a BSS specific frequency that is more than a factor of 20 larger than that in the brightest clusters. Interestingly enough, PCC clusters (open circles), with the exception of NGC 5946 that shows an anomalously low value of $`F_{\mathrm{BSS}}`$, share the same distribution as the normal clusters.

A less significant dependence is found with the central cluster density. There is no clear correlation with the concentration parameter.

King (2002) demonstrated that the rate of stellar collisions in a King model GC is $`5\times 10^{15}(\mathrm{\Sigma }_0^3r_c)^{1/2}`$, where $`\mathrm{\Sigma }_0`$ is the central surface brightness in units of $`L_V`$pc<sup>-2</sup> (equivalent to $`\mu _v=26.41`$), and $`r_c`$ is the core radius in pc. As shown in Fig. 1 (lower panel), $`F_{BSS}`$ anticorrelates with this quantity. The clusters with the highest collision rates have the smallest relative number of BSS. The dispersion of the points in the lower panel is larger than in the upper panel. We note a saturation in the correlation for the clusters with collision rates smaller than $`10^9`$. The two clusters with the smallest total absolute magnitude do not seem to partecipate to this general trend.

## 2. Conclusions

The dependence of the BSS relative frequency with the cluster total mass and with the expected collision rate is noteworthy. Clearly, the environment strongly affects the number of BSS that a GC can form at the present time. Stellar encounters can favour binary destruction: separation of the soft binaries, merger of hard binaries. The latter might have given origin to a large number of BSS in the past, but now high central density clusters have much less primordial binaries to form BSS than low density objects. A decrease in the number of binaries means lower probability for a BSS to form. Of course, high density can also favour encounters, and therefore BSS formation via collision. The anticorrelation between the binary frequency and the collision rate of Fig. 1 suggests that direct collision of stars is not the the dominant mechanism that triggers the formation of BSS, while the number of surviving (primordial?) binaries might be the relevant parameter.

## References

King, I. R. 2002, Introduction to Classical Stellar Dynamics (in Russian) (Moscow: URSS), p. 214

Piotto, G. et al. 2002, A&A, 391, 945


text:
# Correlated Fluctuations between Luminescence and Ionization in Liquid Xenon

## Abstract

The ionization of liquefied noble gases by radiation is known to be accompanied by fluctuations much larger than predicted by Poisson statistics. We have studied the fluctuations of both scintillation and ionization in liquid xenon and have measured, for the first time, a strong anti-correlation between the two at a microscopic level, with coefficient $`0.80<\rho _{ep}<0.60`$. This provides direct experimental evidence that electron-ion recombination is partially responsible for the anomalously large fluctuations and at the same time allows substantial improvement of calorimetric energy resolution.

(EXO Collaboration)

The measurement of ionizing radiation in a liquiefied noble gas such as Ar, Kr, or Xe can be characterized by two parameters: $`W_e`$, the mean energy required to create a free electron-ion pair and the Fano factor, $`F_e`$, which parameterizes the fluctuations in the number of ion pairs. The Fano factor is defined by

$$\sigma _e=\sqrt{F_eN_e}$$
(1)
where $`\sigma _e^2`$ is the variance of the charge expressed in units of the electron charge, $`e`$, and $`N_e`$ is the number of ion pairs. $`F_e=1`$ corresponds to Poisson statistics. Fano originally predicted that the charge fluctuations would be sub-Poissonian, $`F_e<1`$, because the individual ion-pair creation processes are not independent once the additional constraint $`E=N_eW_e`$ is included, where $`E`$ is the energy deposited by the incident particle. The Fano factor for gaseous noble elements is relatively well understood, where for example the theory for argon predicts $`F_e=0.16`$ while experiment yields $`F_e=0.20`$.

The Fano factor for the liquid phase has been far more difficult to understand, even though the values of $`W_e`$ are similar to those of the gas. The pioneering theoretical work of Doke  in 1976 predicted $`F_e0.05`$ for liquid Xe (LXe). Experiments, however, showed $`F_e>20`$, that is, a variance 20 times worse than Poissonian and 400 times worse than predicted by Fano’s original argument. This discrepancy has not only raised important issues for the understanding of the phenomena of energy loss and conversion, but it is also of interest for experimental physics since this large Fano factor limits the resolution of calorimeters used in nuclear and particle physics.

Luminescence light provides a second process, complementary to ionization, with which to study energy loss and conversion phenomena. The properties of luminescence light emitted by ionizing radiation (scintillation) in liquid Ar, Kr, and Xe detectors have been extensively studied and have been exploited in calorimetry. As described elsewhere, scintillation photons are produced by the relaxation of a xenon excimer, $`\mathrm{𝚇𝚎}_\mathrm{𝟸}^{}2\mathrm{𝚇𝚎}+\gamma `$. The excimer is formed in two ways; 1) direct excitation by an energetic particle, 2) recombination of an electron with a xenon molecular ion, $`\mathrm{𝚇𝚎}_\mathrm{𝟸}^++𝚎\mathrm{𝚇𝚎}_\mathrm{𝟸}^{}`$ (recombination luminescence). The electron-ion recombination rate depends upon the ionization density and the applied electric field. The quantities $`W_p`$ and $`F_p`$ can be defined in analogy with the ionization case, where $`W_p`$ is the mean energy absorbed per emitted photon and $`F_p`$ parameterizes their variance. Fano’s theory was not intended to treat photon fluctuations and in this case $`F_p`$ is treated as a phenomenological parameter.

We have simultaneously studied both ionization and scintillation in LXe using a special geometry ionization chamber equipped with a VUV photomultiplier tube (PMT) capable of detecting the $`175`$ nm scintillation light from the xenon. The chamber, shown in Figure 1, has a pancake geometry with a total volume of $`1.5`$ liters. It is built out of ultra-high-vacuum compatible parts, so that it can be baked and treated with hot xenon before filling to obtain good LXe purity. The chamber is cooled by immersion in a bath of HFE7100  that is in turn cooled by liquid nitrogen. The xenon is purified by an Oxysorb  cartridge during transfer from the storage cylinder to the chamber. Ionization is collected in a fiducial volume with radius $``$ 10 mm and height $``$ 6 mm defined by a small grid and a large cathode plane. The grid and cathode are made of an electroformed nickel mesh with 90% optical transparency. Drifted electrons are collected on a low-capacitance anode surrounded by a field defining ring held at the same ground potential. An electron transparency of 100% is obtained by maintaining the anode-to-grid field, ($`E_{AG}`$), double in value with respect to the drift field, ($`E_{KG}`$), for all the measurements reported here. All data was collected at a LXe temperature of $`169.2\pm 0.5`$ K and pressure of $`990_{50}^{+30}`$ Torr. The temperature gradients in the chamber were measured to be less than $`\pm 0.2`$ K.

The electron signal is detected by a low noise charge sensitive preamplifier with a cold input FET placed on the feedthrough of the chamber. After subsequent amplification and filtering, the charge signal $`f(t)`$ is recorded on a transient digitizer and fit to a model assuming drift to the anode of a gaussian ionization distribution represented by the function

$$f(t)=A\left[1+erf\left(\frac{\mathrm{\Delta }t}{C}\right)\right]e^{\frac{t}{\tau }}+\alpha $$
(2)
where $`\mathrm{\Delta }t=tt_0`$ and the term $`\alpha =De^{\frac{t}{\tau }}+100+0.019t`$ accounts for pile-up, ADC offset and a small integral non-linearity. We simultaneously fit for the amplitude $`A`$, the effective drift-time $`t_0`$, the rise-time parameter $`C`$, the pile-up offset $`D`$ and the integrator fall-time constant $`\tau `$. The fit fall-time is $`50\mu `$s for all the data, as expected. The parameter $`C`$ includes information about the multiplicity of the ionization event and its values range from 0.78 $`\mu `$s for a field of 4 kV/cm to 1.04 $`\mu `$s for 0.2 kV/cm. The decrease of this parameter with the increase of the electric field is consistent with the expected change in drift velocity.

The readout is calibrated in units of electron charge by injecting a signal into the input FET through a calibrated capacitor. The noise of the preamplifier was measured with a test signal and found to follow a gaussian distribution with width $`\sigma _N=381e`$.

The scintillation light is detected through the cathode grid by the 3-inch VUV PMT  with a solid angle that varies from 21.5% to 26% of $`4\pi `$ for signals produced in the region between cathode and grid active for ionization collection. The fraction of photoelectrons extracted from the PMT photocathode per scintillation photon produced in the LXe varies between 3.0% and 3.6% considering mesh and quartz window transparencies and photocathode quantum efficiency (22% including transparency of the PMT window). The PMT signal is read out using a 12 bit charge ADC calibrated by measuring the single-photon peak, which is well resolved from noise.

Data acquisition (DAQ) is triggered when there is a delayed coincidence between PMT and charge signals. The delay is required to account for the electron drift-time and varies with electric field. Once a coincidence is established, DAQ is inhibited until the event has been processed and the instruments cleared. About 25% of the collected data is removed off-line by analysis cuts that reject events that don’t fit model (2). These “bad events” are mostly due to multiple ionization deposits from Compton scattering, deposits at the edge of the ionization fiducial region, and cosmic ray background.

The LXe is excited by the 570 and 1064 keV $`\gamma `$’s (and their associated internal conversion electrons and X-rays) emitted by a 207Bi source. The source is prepared as a monoatomic layer electro-plated on a $`3\mu `$m diameter carbon fiber that is then threaded through the cathode grid. The thin source minimizes the effects of self-shadowing of the scintillation light and is essential to obtain acceptable scintillation resolution. The data presented here were collected at drift fields $`E_{KG}`$ of 0.2kV/cm, 0.5kV/cm, 1kV/cm, 2kV/cm, and 4kV/cm.

A two-dimensional scatter plot of ionization versus scintillation signals obtained at 4 kV/cm is shown in Figure 2. The two peaks in the ionization channel can be fit to a gaussian function plus a first order polynomial function. The energy resolution obtained from the ionization data, once the electronics noise is subtracted in quadrature, is in good quantitative agreement with the results of other authors, as shown for the 570 keV peak in Figure 3.

As shown in Figure 2, the scintillation spectrum has substantially poorer resolution and a larger continuum than the ionization spectrum. Differences between the two channels arise from the fact that the VUV photon acceptance depends somewhat upon the position of the energy deposited, as mentioned above, and from the different fiducial volumes. The resolution in the scintillation channel, as extracted by a fit with a gaussian plus a first degree polynomial, is $`\sigma =(16.5\pm 2.5)\%`$ at 570 keV, where the error refers to the range of values observed at different fields.

The observed ionization and scintillation spectra are qualitatively well reproduced in both peak positions and strengths by a detailed GEANT 4  simulation which includes the processes of internal capture, X-ray and Auger electron emission. An empirical gaussian smearing is used to match the GEANT 4 output with the data. However, in the absence of a microscopic model for the production of scintillation photons in LXe, the simulation could not be used to quantitatively fit the two dimensional spectrum.

Previous works  have reported evidence for what we call here “macroscopic” anti-correlation between scintillation light and ionization as a function of the electric field applied to the LXe. This effect, illustrated for our data in Figure 4, is thought to be due to the fact that scintillation occurs as ionization recombines. Different values of the electric field modify such recombination, changing the proportion between light and ionization but leaving their sum constant. The constancy of the sum is well reproduced by data, as shown, in our case, by Figure 4. We note that the absolute value of the sum does not directly correspond to the total energy deposited, as part of the energy is lost in other channels. We refer to this anti-correlation as “macroscopic” because is it an average behavior of the LXe and it is obtained by changing an external, macroscopic, parameter (the electric field).

In Figure 2, clear “microscopic” anti-correlation is evident in each of the two “islands” corresponding to the 570 keV and 1064 keV $`\gamma `$ lines. This effect, simply referred to as anti-correlation below, indicates that energy sharing between scintillation and ionization fluctuates event by event even though the sum of the two remains constant, in analogy with the macroscopic correlations described above. The microscopic nature of the anti-correlation may help in understanding the large Fano factor mentioned above. It also shows the way to dramatic improvements in the energy resolution of particle detectors, since the fluctuations clearly affect scintillation and ionization in opposite ways. Anti-correlations similar to the ones shown in Figure 2 are obtained at all fields and were not previously observed in LXe, although they were shown, without a detailed discussion, for heavy ions in liquid Argon, by Ref. .

The anti-correlation can be formally extracted from the data by fitting each of the peaks to a 2-dimensional gaussian distribution

$$G(N_e,N_p)=He^{(\frac{1}{2(1\rho _{ep}^2)}[\frac{\mathrm{\Delta }N_e^2}{\sigma _e^2}+\frac{\mathrm{\Delta }N_p^2}{\sigma _p^2}\frac{2\rho _{ep}\mathrm{\Delta }N_e\mathrm{\Delta }N_p}{\sigma _e^{}\sigma _p}])}$$
(3)
where $`H=1/2\pi \sigma _e^{}\sigma _p\sqrt{1\rho _{ep}^2}`$ is a normalization constant, $`\rho _{ep}=V_{ep}/\sigma _e^{}\sigma _p`$ is the correlation coefficient between ionization and scintillation signals, $`V_{ep}=(N_e\overline{N_e})(N_p\overline{N_p})`$ is the covariance, $`\mathrm{\Delta }N_e=N_e\overline{N}_e`$ and $`\mathrm{\Delta }N_p=N_p\overline{N}_p`$. While the instrumental noise on the scintillation channel is negligible and hence ignored, the variance for the ionization channel, $`\sigma _e^2=\sigma _e^2+\sigma _N^2`$, includes the term, $`\sigma _N=381e`$, that accounts for the readout noise mentioned above. This is a relatively small correction (e.g. 12% at 4 kV/cm for the 570 keV peak). Unlike the one-dimensional fits used for ionization and scintillation data separately, the 2-dimensional fitting function does not include a background term. A two-step procedure is used to first find an elliptical region of full width $`2.35\sigma `$ around each of the two $`\gamma `$-ray peaks and then perform the final fit. The limited range of the fit suppresses the K-shell internal capture electrons and Compton scattering events. We found that the two-dimensional gaussian fits tend to produce larger variances with respect to the one-dimensional fits, providing a conservative estimate of the real resolutions.

The correlation coefficients found by the fit are plotted in Figure 5 as a function of drift field. A clear anti-correlation exists at all fields for both the 570 keV and the 1064 keV lines. The errors shown in the Figure are systematic and are estimated as the full range of the variation of $`\rho _{ep}`$ when the different analysis cuts and fit intervals are changed. Statistical errors are negligible. The stability of the correlation with respect to different experimental circumstances (such as position in the chamber, size of the chamber and type of radiation) can be tested by repeating the fits for different subset of the data, each presenting a different range of effective drift-times, as derived from Eq. 2. The effective drift-time distribution, shown for 4 kV/cm in the inset of Figure 6, allows for the selection of events originating at different spatial positions along the drift field. The peak at large drift-times primarily corresponds to electrons that lose their energy and stop in the vicinity of the cathode plane where the source is. As shown in Figure 6 the correlation coefficient does not depend significantly on different effective drift times.

By diagonalizing the covariance matrix via a unitary rotation through angle $`\theta =2\rho _{ep}\sigma _e\sigma _p/(\sigma _e^2\sigma _p^2)`$, the correlated two-dimensional gaussian can be reduced to two uncorrelated gaussians with variances $`\sigma _{min}^2`$ and $`\sigma _{max}^2`$. This rotation is mathematically equivalent to measuring the variances along the major and minor axes of the two-dimensional ellipses. Figure 7 shows the ionization resolutions from the 2-dimensional fit and the minimum resolutions in the frame rotated by $`\theta `$ for the 570 keV gamma ray peak at all drift fields.

The resolution improvement observed when correlating the ionization and scintillation signals is significantly better than the uncorrelated addition of the two spectra, where $`\sigma _{p+e}=(\sigma _p^2+\sigma _e^2)^{1/2}`$. At 4 kV/cm $`\sigma _{p+e}=3.7\%`$. When using the correlated nature of the two channels the resolution $`\sigma _{min}=3.0\%`$, a factor of 1.24 improvement. The improvement is more dramatic at the 0.2 kV/cm field where the ratio is $`9.7/5.4=1.8`$. While we found that the absolute values of the resolution change somewhat for different fitting methods and at different stages of the iterative procedure, the improvement factor is very stable. We note here that the main results of this work are to show that a microscopic correlation exists and it can be consistently used to improve the resolution. The ultimate resolution will have to be measured at a later stage with a detector of larger volume, better VUV light collection and, possibly, a different source and it is expected to be even better than the resolution $`\sigma _{min}`$ reported in Figure 7.

While the nature of the large Fano factor in liquefied noble gases is still not well understood, the different models considered generally involve fluctuations in the density of the energy deposited . Such fluctuations can arise from density variations in the medium or from the small number of delta electrons with energies between 1 and 30 keV that are produced by Rutherford scattering of atomic electrons in the liquid. When these electrons come to a stop they produce a dense cluster of electron-ion pairs whose recombination rate is relatively high. In both cases more recombination occurs where the ionization density is higher. Hence the statistics relevant for the fluctuations in the liquid is that of the (small) number of high ionization density regions (either the number of delta-electrons or the number of high density regions). In a sense, this would mean that the Fano factors computed in Eq. 1 are too small because the “wrong” $`N`$ is used.

As first suggested in Ref.  the recombination hypothesis above implies that the fluctuations of ionization and scintillation should be anti-correlated, since the high ionization-density regions which produce less free charge should emit a correspondingly larger amount of scintillation light. A linear combination of ionization and scintillation should therefore have smaller fluctuations than either alone. Note that since the Fano factors $`F_e`$ and $`F_p`$ are typically far in excess of their Poissonian values, the improvement should in principle be substantial and is not limited by the factor of $`\sqrt{2}`$ one would expect from the combination of two uncorrelated signals of equal resolution. The data presented here provides direct experimental evidence of this fluctuation and the associated recombination by observing a clear microscopic anti-correlation between ionization and scintillation. The correlation, $`\rho _{ep}`$, and the magnitude of the variances, $`\sigma _{max}^2`$ and $`\sigma _{min}^2`$, provide input to improve the model in terms of observables that are more directly related to the microscopic physics. It is interesting to note that from Figure 5 it appears that the magnitude of the correlation $`\rho _{ep}`$ decreases for increasing electric field (at least at low fields). This behavior is possibly connected to the dependence of the intensity of the recombination luminescence with the electric field, as expected in all the models discussed.

In conclusion, we have observed a clear microscopic anti-correlation between ionization and scintillation in liquid xenon. We have also measured the correlation coefficient between these two channels to be $`0.8<\rho _{ep}<0.6`$. This observation has two important consequences. First, it provides direct experimental evidence that recombination is a substantial source of the anomalously large Fano factors in the high density liquid state. Second, it provides a method for improving the resolution of liquefied noble gas calorimeters.

One of us (G.G.) is indebted to P. Picchi, F. Pietropaolo and T. Ypsilantis for early discussions on the subject, M. Moe for a careful reading of the manuscript and to R.G.H. Robertson for having stimulated the interest on this topic. This work was supported by US-DoE grant DE-FG03-90ER40569-A019 and by a Terman Fellowship of the Stanford University.


text:
# New Orbits for the 𝑛-Body Problem

## 1 Least Action Principle

Given $`n`$ bodies, let $`m_j`$ denote the mass and $`z_j(t)`$ denote the position in $`^2=`$ of body $`j`$ at time $`t`$. The action functional is a mapping from the space of all trajectories, $`z_1(t),z_2(t),\mathrm{},z_n(t)`$, $`0t2\pi `$, into the reals. It is defined as the integral over one period of the kinetic minus the potential energy:

$$A=_0^{2\pi }\left(\underset{j}{}\frac{m_j}{2}\dot{z}_j^2+\underset{j,k:k<j}{}\frac{m_jm_k}{z_jz_k}\right)𝑑t.$$

Stationary points of the action function are trajectories that satisfy the equations of motions, i.e., Newton’s law gravity. To see this, we compute the first variation of the action functional,

$`\delta A`$ $`=`$ $`{\displaystyle _0^{2\pi }}{\displaystyle \underset{\alpha }{}}\left({\displaystyle \underset{j}{}}m_j\dot{z}_j^\alpha \dot{\delta z}_j^\alpha {\displaystyle \underset{j,k:k<j}{}}m_jm_k{\displaystyle \frac{(z_j^\alpha z_k^\alpha )(\delta z_j^\alpha \delta z_k^\alpha )}{z_jz_k^3}}\right)dt`$
$`=`$ $`{\displaystyle _0^{2\pi }}{\displaystyle \underset{j}{}}{\displaystyle \underset{\alpha }{}}\left(m_j\ddot{z}_j^\alpha +{\displaystyle \underset{k:kj}{}}m_jm_k{\displaystyle \frac{z_j^\alpha z_k^\alpha }{z_jz_k^3}}\right)\delta z_j^\alpha dt,`$
and set it to zero. We get that

$$m_j\ddot{z}_j^\alpha =\underset{k:kj}{}m_jm_k\frac{z_j^\alpha z_k^\alpha }{z_jz_k^3},j=1,2,\mathrm{},n,\alpha =1,2$$
(1)
Note that if $`m_j=0`$ for some $`j`$, then the first order optimality condition reduces to $`0=0`$, which is not the equation of motion for a massless body. Hence, we must assume that all bodies have strictly positive mass.

## 2 Periodic Solutions

Our goal is to use numerical optimization to minimize the action functional and thereby find periodic solutions to the $`n`$-body problem. Since we are interested only in periodic solutions, we express all trajectories in terms of their Fourier series:

$$z_j(t)=\underset{k=\mathrm{}}{\overset{\mathrm{}}{}}\gamma _ke^{ikt},\gamma _k.$$

Abandoning the efficiency of complex-variable notation, we can write the trajectories with components $`z_j(t)=(x_j(t),y_j(t))`$ and $`\gamma _k=(\alpha _k,\beta _k).`$ So doing, we get

$`x(t)`$ $`=`$ $`a_0+{\displaystyle \underset{k=1}{\overset{\mathrm{}}{}}}\left(a_k^c\mathrm{cos}(kt)+a_k^s\mathrm{sin}(kt)\right)`$
$`y(t)`$ $`=`$ $`b_0+{\displaystyle \underset{k=1}{\overset{\mathrm{}}{}}}\left(b_k^c\mathrm{cos}(kt)+b_k^s\mathrm{sin}(kt)\right)`$
where

$`a_0=\alpha _0,a_k^c`$ $`=`$ $`\alpha _k+\alpha _k,a_k^s=\beta _k\beta _k,`$
$`b_0=\beta _0,b_k^c`$ $`=`$ $`\beta _k+\beta _k,b_k^s=\alpha _k\alpha _k.`$
Since we plan to optimize over the space of trajectories, the parameters $`a_0`$, $`a_k^c`$, $`a_k^s`$, $`b_0`$, $`b_k^c`$, and $`b_k^s`$ are the decision variables in our optimization model. The objective is to minimize the action functional.

ampl is a small programming language designed for the efficient expression of optimization problems Fourer et al. (1993). Figure 1 shows the ampl program for minimizing the action functional.

Note that the action functional is a nonconvex nonlinear functional. Hence, it is expected to have many local extrema and saddle points. We use the author’s local optimization software called loqo (see Vanderbei (1999), Vanderbei and Shanno (1999)) to find local minima in a neighborhood of an arbitrary given starting trajectory. One can provide either specific initial trajectories or one can give random initial trajectories. The four lines just before the call to solve in Figure 1 show how to specify a random initial trajectory. Of course, ampl provides capabilities of printing answers in any format either on the standard output device or to a file. For the sake of brevity and clarity, the print statements are not shown in Figure 1. ampl also provides the capability to loop over sections of code. This is also not shown but the program we used has a loop around the four initialization statements, the call to solve the problem, and the associated print statements. In this way, the program can be run once to solve for a large number of periodic solutions.

### 2.1 Choreographies

Recently, Chenciner and Montgomery (2000) introduced a new family of solutions to the $`n`$-body problem called choreographies. A choreography is defined as a solution to the $`n`$-body problem in which all of the bodies share a common orbit and are uniformly spread out around this orbit. Such trajectories are even easier to find using the action principle. Rather than having a Fourier series for each orbit, it is only necessary to have one master Fourier series and to write the action functional in terms of it. Figure 2 shows the ampl model for finding choreographies.

## 3 Stable vs. Unstable Solutions

Figure 3 shows some simple choreographies found by minimizing the action functional using the ampl model in Figure 2. The famous $`3`$-body figure eight, first discoverd by Moore (1993) and later analyzed by Chenciner and Montgomery (2000), is the first one shown—labeled FigureEight3. It is easy to find choreographies of arbitrary complexity. In fact, it is not hard to rediscover most of the choreographies given in Chenciner et al. (2001), and more, simply by putting a loop in the ampl model and finding various local minima by using different starting points.

However, as we discuss in a later section, simulation makes it apparent that, with the sole exception of FigureEight3, all of the choreographies we found are unstable. And, the more intricate the choreography, the more unstable it is. Since the only choreographies that have a chance to occur in the real world are stable ones, many cpu hours were devoted to searching for other stable choreographies. So far, none have been found. The choreographies shown in Figure 3 represent the ones closest to being stable.

Given the difficulty of finding stable choreographies, it seems interesting to search for stable nonchoreographic solutions using, for example, the ampl model from Figure 1. The most interesting such solutions are shown in Figure 4. The one labeled Ducati3 is stable as are Hill3\_15 and the three DoubleDouble solutions. However, the more exotic solutions (OrthQuasiEllipse4, Rosette4, PlateSaucer4, and BorderCollie4) are all unstable.

For the interested reader, a java applet can be found at Vanderbei (2001) that allows one to watch the dynamics of each of the systems presented in this paper (and others). This applet actually integrates the equations of motion. If the orbit is unstable it becomes very obvious as the bodies deviate from their predicted paths.

### 3.1 Ducati3 and its Relatives

The Ducati3 orbit first appeared in Moore (1993) and has been independently rediscovered by this author, Broucke Broucke (2003), and perhaps others. Simulation reveals it to be a stable system. The java applet at Vanderbei (2001) allows one to rotate the reference frame as desired. By setting the rotation to counter the outer body in Ducati3, one discovers that the other two bodies are orbiting each other in nearly circular orbits. In other words, the first body in Ducati3 is executing approximately a circular orbit, $`z_1(t)=e^{it}`$, the second body is oscillating back and forth roughly along the $`x`$-axis, $`z_2(t)=\mathrm{cos}(t)`$, and the third body is oscillating up and down the $`y`$-axis, $`z_3(t)=i\mathrm{sin}(t)`$. Rotating so as to fix the first body means multiplying by $`e^{it}`$:

$`\overline{z}_1(t)`$ $`=`$ $`e^{it}(e^{it})=1`$
$`\overline{z}_2(t)`$ $`=`$ $`e^{it}\mathrm{cos}(t)=(1+e^{2it})/2`$
$`\overline{z}_2(t)`$ $`=`$ $`e^{it}i\mathrm{sin}(t)=(1e^{2it})/2.`$
Now it is clear that bodies 2 and 3 are orbiting each other at half the distance of body 1. So, this system can be described as a Sun, Earth, Moon system in which all three bodies have equal mass and in which one (sidereal) month equals one year. The synodic month is shorter—half a year.

This analysis of Ducati3 suggests looking for other stable solutions of the same type but with different resonances between the length of a month and a year. Hill3\_15 is one of many such examples we found. In Hill3\_15, there are 15 sidereal months per year. Let Hill3\_$`n`$ denote the system in which there are $`n`$ months in a year. All of these orbits are easy to calculate and they all appear to be stable. This success suggests going in the other direction. Let Hill3\_$`\frac{1}{n}`$ denote the system in which there are $`n`$ years per month. We computed Hill3\_$`\frac{1}{2}`$ and found it to be unstable. It is shown in Figure 6.

In the preceding discussion, we decomposed these Hill-type systems into two $`2`$-body problems: the Earth and Moon orbit each other while their center of mass orbits the Sun. This suggests that we can find stable orbits for the $`4`$-body problem by splitting the Sun into a binary star. This works. The orbits labeled DoubleDouble$`n`$ are of this type. As already mentioned, these orbits are stable.

Given the existence and stability of FigureEight3, one often is asked if there is any chance to observe such a system among the stars. The answer is that it is very unlikely since its existence depends crucially on the masses being equal. The Ducati and Hill type orbits, however, are not constrained to have their masses be equal. Figure 5 shows several Ducati-type orbits in which the masses are not all equal. All of these orbits are stable. This suggests that stability is common for Ducati and Hill type orbits. Perhaps such orbits can be observed.

## 4 Limitations of the Model

The are certain limitations to the approach articulated above. First, the Fourier series is an infinite sum that gets truncated to a finite sum in the computer model. Hence, the trajectory space from which solutions are found is finite dimensional.

Second, the integration is replaced with a Riemann sum. If the discretization is too coarse, the solution found might not correspond to a real solution to the $`n`$-body problem. The only way to be sure is to run a simulator.

Third, as mentioned before, all masses must be positive. If there is a zero mass, then the stationary points for the action function, which satisfy (1), don’t necessarily satisfy the equations of motion given by Newton’s law.

Lastly, the model, as given in Figure 1, can’t solve 2-body problems with eccentricity. We address this issue in the next section.

## 5 Elliptic Solutions

An ellipse with semimajor axis $`a`$, semiminor axis $`b`$, and having its left focus at the origin of the coordinate system is given parametrically by:

$$x(t)=f+a\mathrm{cos}t,y(t)=b\mathrm{sin}t,$$

where $`f=\sqrt{a^2b^2}`$ is the distance from the focus to the center of the ellipse.

However, this is not the trajectory of a mass in the $`2`$-body problem. Such a mass will travel faster around one focus than around the other. To accomodate this, we need to introduce a time-change function $`\theta (t)`$:

$$x(t)=f+a\mathrm{cos}\theta (t),y(t)=b\mathrm{sin}\theta (t).$$

This function $`\theta `$ must be increasing and must satisfy $`\theta (0)=0`$ and $`\theta (2\pi )=2\pi `$.

The optimization model can be used to find (a discretization of) $`\theta (t)`$ automatically by changing param theta to var theta and adding appropriate monotonicity and boundary constraints. In this manner, more realistic orbits can be found that could be useful in real space missions.

In particular, using an eccentricity $`e=f/a=0.0167`$ and appropriate Sun and Earth masses, we can find a periodic Hill-Type satellite trajectory in which the satellite orbits the Earth once per year.

## 6 Sensitivity Analysis

The determination of stability vs. instability mentioned the previous sections was done empirically by simulating the orbits with a integrator and very small step sizes. Two integrators were used: a midpoint integrator and a $`4`$-th order Runge-Kutta integrator. Orbits that are claimed to be stable were run for several hours of cpu time (which corresponds to many thousands of orbits) without falling apart. Orbits that are claimed to be unstable generally became obviously so in just a few seconds of cpu time, which corresponds to only a few full orbits. In this section, we describe a Floquet analysis of stability and present this measure of stability for the various orbits found.

For simplicity, in this section we assume that all masses are equal to one. Let $`\xi ^{}(t)=(z^{}(t),\dot{z}^{}(t))`$ be a particular solution to

$$\dot{\xi }=A(\xi )$$

where

$$A(z(t),\dot{z}(t))=(\dot{z}(t),a(z(t)))$$

and

$$a(z)=(a_1(z),\mathrm{},a_n(z))$$

and

$$a_j(z)=\underset{k:kj}{}\frac{z_jz_k}{z_jz_k^2},j=1,2,\mathrm{},n.$$

Consider a nearby solution $`\xi (t)`$:

$`\dot{\xi }(t)`$ $`=`$ $`A(\xi (t))`$
$``$ $`A(\xi ^{}(t))+A^{}(\xi ^{}(t))(\xi (t)\xi ^{}(t))`$
$`=`$ $`\dot{\xi }^{}(t)+A^{}(\xi ^{}(t))(\xi (t)\xi ^{}(t)).`$
Put $`\mathrm{\Delta }\xi =\xi \xi ^{}`$. Then $`\dot{\mathrm{\Delta }\xi }=A^{}(\xi ^{}(t))\mathrm{\Delta }\xi .`$ A finite difference approximation yields

$`\mathrm{\Delta }\xi (t+h)`$ $`=`$ $`\mathrm{\Delta }\xi (t)+hA^{}(\xi ^{}(t))\mathrm{\Delta }\xi (t)`$
$`=`$ $`\left(I+hA^{}(\xi ^{}(t))\right)\mathrm{\Delta }\xi (t).`$
Iterating around one period, we get:

$$\mathrm{\Delta }\xi (T)=\left(\underset{i=0}{\overset{n1}{}}\left(I+hA^{}(\xi ^{}(t_i))\right)\right)\mathrm{\Delta }\xi (0),$$

where $`h=T/n`$ and $`t_i=iT/n`$.

The following perturbations, which are associated with invariants of the physical laws, are unimportant in calculating $`\mathrm{\Delta }\xi (T)`$:

$$\left[\begin{array}{c}\mathrm{\Delta }z\\ \mathrm{\Delta }\dot{z}\end{array}\right]=\left[\begin{array}{c}e_1\\ e_1\\ e_1\\ 0\\ 0\\ 0\end{array}\right],\left[\begin{array}{c}e_2\\ e_2\\ e_2\\ 0\\ 0\\ 0\end{array}\right],\left[\begin{array}{c}0\\ 0\\ 0\\ e_1\\ e_1\\ e_1\end{array}\right],\left[\begin{array}{c}0\\ 0\\ 0\\ e_2\\ e_2\\ e_2\end{array}\right],\left[\begin{array}{c}Rz_1\\ Rz_2\\ Rz_3\\ R\dot{z}_1\\ R\dot{z}_2\\ R\dot{z}_3\end{array}\right],\frac{1}{2}\left[\begin{array}{c}3\dot{z}_1+2z_1\\ 3\dot{z}_1+2z_1\\ 3\dot{z}_1+2z_1\\ 3a_1\dot{z}_1\\ 3a_2\dot{z}_2\\ 3a_3\dot{z}_3\end{array}\right],$$

where $`R`$ denotes rotation by $`90^{}`$. The first two of these perturbations correspond to translation. The next two correspond to moving frame of reference and the last two correspond to rotation, and dilation. Dilation is explained below. Of course, all positions and velocities are evaluated at $`t=0`$. Vector $`e_i`$ denotes the $`i`$-th unit vector in $`^2`$.

Consider spatial dilation by $`\rho `$ together with a temporal dilation by $`\theta `$:

$$Z_j(t)=\rho z_j(t/\theta ).$$

Given that the $`z_j`$’s are a solution, it is easy to check that

$$\ddot{Z}_j(t)=\frac{\rho ^3}{\theta ^2}\underset{kj}{}\frac{Z_j(t)Z_k(t)}{Z_j(t)Z_k(t)^2}.$$

Hence, if mass is to remain fixed, we must have that $`\rho ^3=\theta ^2`$:

$$Z_j(t)=\rho z_j(t/\rho ^{3/2})\dot{Z}_j(t)=\rho ^{1/2}\dot{z}_j(t/\rho ^{3/2}).$$

To find the perturbation direction corresponding to this dilation, we differentiate with respect to $`\rho `$ at $`\rho =1`$:

$$\frac{d}{d\rho }\left[\begin{array}{c}\rho z_j(t/\rho ^{3/2})\\ \rho ^{1/2}\dot{z}_j(t/\rho ^{3/2})\end{array}\right]|_{\rho =1}=\left[\begin{array}{c}\frac{3}{2}\dot{z}_j+z_j\\ \frac{3}{2}a_j\frac{1}{2}\dot{z}_j\end{array}\right].$$

For checking stability, we project any initial perturbation onto the null space of $`P^T`$, where

$$P=\left[\begin{array}{cccccc}e_1& e_2& 0& 0& Rz_1& (3\dot{z}_1+2z_1)/2\\ e_1& e_2& 0& 0& Rz_2& (3\dot{z}_2+2z_2)/2\\ e_1& e_2& 0& 0& Rz_3& (3\dot{z}_3+2z_3)/2\\ 0& 0& e_1& e_2& R\dot{z}_1& (3a_1\dot{z}_1)/2\\ 0& 0& e_1& e_2& R\dot{z}_2& (3a_2\dot{z}_2)/2\\ 0& 0& e_1& e_2& R\dot{z}_3& (3a_3\dot{z}_3)/2\end{array}\right].$$

The projection matrix is given by

$$\mathrm{\Pi }=IP(P^TP)^1P^T.$$

From the fact that $`z_1+z_2+z_3=0`$ and $`\dot{z}_1+\dot{z}_2+\dot{z}_3=0`$, it follows that all columns of $`P`$ are mutually orthogonal except for the 5-th and 6-th columns. Hence, $`P^TP`$ is not a purely diagonal matrix.

Let

$$\mathrm{\Lambda }_n=\left(\underset{i=0}{\overset{n1}{}}\left(I+hA^{}(\xi ^{}(t_i))\right)\right).$$

We say that an orbit is stable if all eigenvalues of

$$\underset{n\mathrm{}}{lim}\mathrm{\Lambda }_n\mathrm{\Pi }$$

are at most one in magnitude.

### 6.1 Stable Orbits

We computed $`\mathrm{\Lambda }_n`$ for $`n=10^6`$. Table 1 shows maximum eigenvalues for those orbits that seemed stable from simulation. Table 2 shows maximum eigenvalues for those orbits that appeared unstable when simulated.

Acknowledgements. The author received support from the NSF (CCR-0098040) and the ONR (N00014-98-1-0036).


text:
# A 3×10⁹ solar mass black hole in the quasar SDSS J1148+5251 at 𝑧=6.41

## 1 Introduction

The vast energy requirements of the most luminous quasars can be met by a model invoking extraction of gravitational potential energy from matter falling toward a supermassive ($``$ a billion solar mass) black hole. The correlation between black-hole mass and galaxy mass observed in the local Universe strongly suggests that the most luminous quasars will reside in the most massive galaxies in the most massive dark matter halos. The luminous quasars at redshifts $`z>6`$ found in the Sloan Digital Sky Survey (SDSS) therefore pinpoint the earliest massive objects to form (Fan et al. 2001, 2003).

The recent discovery of SDSS J114816.64+525150.3 (hereafter SDSS J1148+5251) at a redshift of $`z=6.43`$ by Fan et al. (2003) makes it the most distant known quasar, observed only 840 million years after the beginning of the universe (cosmological parameters of $`H_0=70\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$, $`\mathrm{\Omega }_\mathrm{M}=0.3`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ are assumed throughout). The very high luminosity of this quasar ($`M_{1450}=27.8`$) suggests that it contains a black-hole mass of at least several billion solar masses.

We have obtained near-infrared spectroscopy of SDSS J1148+5251 to search for the Mg ii emission line (rest wavelength 2799 Å) redshifted into the $`K`$-band. This line is particularly useful since it is usually found close to the systemic redshift, unlike higher ionization lines such as C iv which typically show large blueshifts $`1000`$ km s-1 (Richards et al. 2002). The redshift for this quasar given by Fan et al. was based on the position of a heavily absorbed Ly $`\alpha `$ line and the onset of the Lyman break and is therefore quite uncertain. An accurate redshift is important for determining the strength of the proximity effect and for molecular line searches. Additionally, it has recently been shown that the width of the Mg ii emission line can be used to derive the mass of the quasars central black hole in the same way that the H $`\beta `$ line has traditionally been used (McLure & Jarvis 2002).

We describe the spectroscopic observations and measurements made from them in §2. The importance of the new redshift determination is discussed in §3. We determine an estimate for the mass of the black hole in §4. Our conclusions are presented in §5.

## 2 Near-infrared spectroscopy

### 2.1 Observations

We observed SDSS J1148+5251 with the United Kingdom Infra-Red Telescope (UKIRT) equipped with the UKIRT Imager Spectrometer (UIST) on 2003 January 10. UIST is a 1-5 $`\mu `$m imager-spectrometer with a 1024x1024 InSb array. Conditions were photometric with seeing at $`K`$-band approximately 0.6 arcsec. The $`HK`$ grism was employed with a 0.48 arcsec wide slit giving continuous wavelength coverage from 1.4 to 2.5 microns with a resolving power of $`500`$. The pixels have size 10.9 Å spectrally and 0.12 arcsec spatially. The target was nodded along the slit a distance of 12 arcsec in between exposures to enable accurate sky subtraction. The total integration time of the observation was 6480 seconds.

The data frames were flat-fielded and then adjacent pairs were subtracted from each other to subtract the sky emission and produce positive and negative spectra of the target. These data were combined and any pixels deviating from the mean at the $`3\sigma `$ level were excluded since this indicates contamination by cosmic rays. The residual background level was subtracted from the combined image. The positive and negative quasar spectra were extracted from this image using apertures of 1.0 arcsec and combined.

Flux calibration and atmospheric extinction corrections were made by dividing the spectrum by that of a star of spectral type F6 which was observed immediately preceding the quasar. An accurate wavelength calibration was made by comparison to an Argon arc lamp. This solution was then checked against the positions of about 20 sky emission lines giving an rms uncertainty in the wavelength calibration of 1 Å. The absolute flux scale of the spectrum was set by comparison with the magnitude measured in a snapshot $`K`$-band image. The image gave $`K=16.86\pm 0.05`$ consistent with the $`K^{}=16.91`$ quoted by Fan et al. (2003).

### 2.2 Emission line fitting

The spectrum of SDSS J1148+5251 is shown in Fig. 1. It shows a broad emission line at approximately the position expected for redshifted Mg ii. To determine the wavelength and width of the Mg ii emission line, a line-fitting process was applied which involves the determination of the underlying continuum and Fe ii emission surrounding the Mg ii line. Due to the blending of the Mg ii line emission and the Fe ii features in the 2700 Å – 2900 Å region, the fitting of the continuum and Fe ii emission is confined to two bands either side of the Mg ii line (rest-frame 2250 Å – 2700 Å and 2900 Å – 3100 Å).

The Fe ii emission in this spectral region is modeled using an Fe ii template based on an archival HST FOS spectrum of the narrow line Seyfert galaxy I Zw 1, which is notable for its extremely strong Fe ii emission. Due to the fact that the Fe ii emission in I Zw 1 is relatively narrow (FWHM $`900`$ kms-1), this template can be smoothed in order to match the much broader Fe ii features typical of powerful quasars. The fitting of the Fe ii emission template and AGN continuum are performed simultaneously, with the amplitude and FWHM of the Fe ii iron template left as free parameters. The AGN continuum is modeled as a power law which is normalized with the continuum magnitude at rest-frame 1280 Å of $`AB_{1280}=19.10`$ given by Fan et al. (2003). The spectral index of the power law is a free parameter and the fitting produces a value of $`\alpha =0.2`$ (defining $`f_\nu \nu ^\alpha `$, where $`f_\nu `$ is the flux-density at frequency $`\nu `$). This continuum slope is typical of quasars (e.g. Vanden Berk et al. 2001). The goodness of fit for the Fe ii +continuum fit is determined from a chi-square test using the noise array from the observation.

Once the minimum chi-square fit has been determined the best-fitting combination of continuum and Fe ii emission is then subtracted from the spectrum, leaving the isolated Mg ii emission line in the region 2700 Å – 2900 Å. The Mg ii emission line itself is then modeled using two identical gaussians which represent the line doublet ($`\lambda \lambda 2796,2802`$). During the fitting of the line profile, the amplitude, FWHM and central wavelength of the gaussians are treated as free parameters, while the doublet separation is held fixed at its laboratory value. The best fit to the Mg ii emission was determined with a chi-square test. Although more complex deconvolutions of the line profile are obviously possible, due to the signal-to-noise constraints of the spectrum the above procedure was adopted because it was able to effectively describe the Mg ii line profile with the minimum level of complexity.

We also performed an alternative fitting procedure in which the continuum was fixed at two points with low expected Fe ii emission and then the Mg ii and Fe ii features fit simultaneously. The results from this independent fitting method are consistent with those using the former method. The ratio of Fe ii /Mg ii emission lines in high redshift quasars is a powerful diagnostic of early element production inside supernovae (Yoshii, Tsujimoto & Kawara 1998). Unfortunately we do not have sufficient signal to noise to place meaningful constraints on this ratio using either of our two fitting methods.

The best fit Mg ii line doublet in SDSS J1148+5251 is at a redshift of $`z=6.41\pm 0.01`$. The deconvolved FWHM of each component of the doublet is $`6000_{600}^{+1100}`$ km s-1. The rest-frame equivalent width of the Mg ii line is 35 Å. The typical equivalent width of Mg ii lines in lower redshift SDSS quasars is 32 Å (Vanden Berk et al. 2001). The fact that the equivalent width of Mg ii in SDSS J1148+5251 is similar to that of other quasars indicates that the UV continuum emission from this quasar is not strongly beamed. Fan et al. (2003) noted that the Ly $`\alpha `$ equivalent width of this and other $`z6`$ quasars are smaller than lower redshift quasars. Pentericci et al. (2002) observed the C iv line in the $`z=6.28`$ quasar SDSS J1030+0524 and also found an equivalent width typical of lower redshift quasars. These results suggest that the low Ly $`\alpha `$ equivalent widths are primarily due to Ly $`\alpha `$ absorption and not a fundamental difference between the emission properties of these quasars.

## 3 An accurate redshift determination

The redshift of $`z=6.43\pm 0.05`$ for this quasar given in Fan et al. (2003) is based on the position of the heavily absorbed Ly $`\alpha `$ line and the onset of the Lyman break. The lack of distinguishable N v emission or any other lines in the optical spectrum mean that a more accurate redshift could not be obtained. In our spectrum we are able to measure the location of the Mg ii line and hence obtain a more accurate redshift for this quasar. We find a redshift of $`z=6.41\pm 0.01`$. The Mg ii line is particularly useful for redshift estimation compared to, for example, C iv which could be observed in the $`J`$-band. The low-ionization Mg ii line is usually found close to the systemic redshift, whereas the high-ionization C iv line often shows blueshifts of several thousand km s-1 (Richards et al. 2002). Another important fact is that there is no systematic velocity offset for the Ly $`\alpha `$ and Mg ii lines in a large sample of SDSS quasars (Vanden Berk et al. 2001), so we can predict the intrinsic, unabsorbed Ly $`\alpha `$ line peak.

Accurate redshifts for the most distant quasars are important for several reasons. The rest-frame UV spectra of $`z6`$ quasars show a deficit of flux shortward of Ly $`\alpha `$ due to absorption by neutral hydrogen in the intergalactic medium (Fan et al. 2001; Becker et al. 2001). These observations have led to the possibility that we are now witnessing the epoch of reionization. Whilst accurate quasar redshifts are not essential for determining the ionization state of the universe as a function of redshift, they are very important when dealing with the state of the IGM close to the quasar.

Assuming that the IGM in general is neutral at the redshift of the quasar, the strong ionizing field of the quasar will ionize a local region around the quasar allowing Ly $`\alpha `$ photons to be transmitted. This proximity effect has been observed as a lack of absorption lines in the spectra of quasars close to the peak of Ly $`\alpha `$ (Bajtlik, Duncan & Ostriker 1988). Cen & Haiman (2000) showed that the size of this ionized region and the transmitted flux blueward of the peak of Ly $`\alpha `$ depend upon the number of ionizing photons which have been emitted from the quasar, i.e. the product of its lifetime and luminosity. Haiman & Cen (2002) apply this to the spectrum of the $`z=6.28`$ quasar SDSS J1030+0524 to show that it is not strongly gravitationally lensed and to derive a lower limit on how long the quasar has been active. Such an analysis depends critically on knowing where the peak of the Ly $`\alpha `$ line should be in the spectrum.

The spectrum of SDSS J1148+5251 presented by Fan et al. (2003) also shows zero flux at Ly $`\alpha `$ redshifts of 5.7 to 6.33, suggestive of a neutral universe. However, there is some flux apparent on the blue wing of Ly $`\alpha `$ which is likely due to the proximity effect. The marked location of Ly $`\alpha `$ at $`z=6.43`$ given by Fan et al. is at the peak of the Ly $`\alpha `$ emission. However, adopting the new redshift for the quasar of $`z=6.41`$ shifts the location of Ly $`\alpha `$ to 9010 Å, which is co-incident with a steep decline in flux to less than half of the peak value. Hence the revised redshift indicates a substantially more neutral IGM close to the quasar. The amount and extent of transmitted flux on the blue side of Ly $`\alpha `$ is similar to that in SDSS J1030+0524 and the conclusions about this quasar made by Haiman & Cen (2002) and Pentericci et al. (2002) probably apply equally to SDSS J1148+5251.

The evolutionary state of the host galaxies of high redshift quasars can be probed by (sub)millimeter observations of dust and molecular gas. Particularly important tracers are the emission lines from transitions of CO molecules and carbon atoms (Blain et al. 2000). However current millimeter receivers have a fairly small bandwidth corresponding to $`2000`$ km s-1. Hence systemic redshifts must be known to an accuracy of $`0.01`$ to ensure that the emission line will fall within the frequency range of the observation. Our observations of SDSS J1148+5251 provide a redshift of sufficient accuracy for millimeter line searches.

## 4 Estimating the black-hole mass

Reverberation mapping studies of low redshift quasars have been successfully employed to estimate the masses of their black holes (Netzer & Peterson 1997). This method determines the radial distance of the broad line emitting region (BLR) from the nucleus ($`R_{\mathrm{BLR}}`$), by observing the time lag ($`1`$ year) between variations in the UV/optical continuum and the line strength. Combining this distance with the velocity of the emission line gas ($`V_{\mathrm{BLR}}`$, which is closely related to the observed FWHM of the line) gives a virial estimate of the black-hole mass ($`M_{\mathrm{bh}}=G^1R_{\mathrm{BLR}}V_{\mathrm{BLR}}^2`$). This estimate assumes that the dynamics of the gas are dominated by gravitational forces (see Peterson & Wandel 2000 for a discussion). These reverberation studies have revealed a correlation between the continuum luminosity and the radius of the BLR which can be used to estimate black-hole masses when only the line FWHM and luminosity are available (Kaspi et al. 2000).

Traditionally, the H $`\beta `$ emission line has been used for this work, since it is readily observable in the optical spectra of low redshift quasars where the reverberation mapping calibration of luminosity and BLR has been derived. Recently, McLure & Jarvis (2002) have cross-calibrated this relationship from the H $`\beta `$ line to the Mg ii line using a sample of AGN with reverberation mapping measures of $`R_{\mathrm{BLR}}`$, so that black-hole mass estimates can be derived to higher redshifts. They note that the Mg ii line (which has an ionization potential close to that of H $`\beta `$) is a much better choice than higher ionization lines like C iv whose dynamics may be strongly influenced by outflows. McLure & Jarvis derived a calibration between $`R_{\mathrm{BLR}}`$ and $`L_{3000}`$ (the monochromatic luminosity at 3000 Å) and performed a fit to the H $`\beta `$ and Mg ii FWHMs, finding a linear relationship with virtually identical normalization. The resulting black-hole mass estimator has the form

$$\frac{M_{\mathrm{bh}}}{\text{ }\mathrm{M}_{}}=3.37\left(\frac{\lambda L_{3000}}{10^{37}\mathrm{W}}\right)^{0.47}\left(\frac{\mathrm{FWHM}(\mathrm{Mg}\mathrm{II})}{\mathrm{kms}^1}\right)^2.$$
(1)
Our spectrum covers the rest wavelength region 3000 Å and we have determined the Mg ii FWHM. Therefore we can apply this formula to derive an estimate of the black-hole mass in SDSS J1148+5251. Using $`\lambda L_{3000}=6.2\times 10^{39}`$ W and FWHM(Mg ii) = 6000 km s-1 we get $`M_{\mathrm{bh}}=3\times 10^9`$ $`\mathrm{M}_{}`$. Assuming that $`z6`$ quasars are not fundamentally different to low redshift quasars (evidence supporting this is that their spectra are similar), the accuracy of this black-hole mass estimator with respect to those from reverberation mapping is a factor of 2.5 ($`1\sigma `$).

This is an extremely massive black hole to be in existence such a short time (840 Myr) after the beginning of the universe. Assuming that this black hole is residing in a proportionately massive halo ($`10^{13}`$ $`\mathrm{M}_{}`$) then standard theories of structure formation predict that this is one of the rarest, most massive halos to have collapsed by this time (e.g. Haiman & Loeb 2001; Fan et al. 2001). There are two caveats to the derived black-hole mass. Since the estimator depends on $`L_{3000}^{0.47}`$, if we have over-estimated the intrinsic, isotropic luminosity we will have over-estimated the black hole mass. If the luminosity function of $`z6`$ quasars has a steep slope then there are expectations that a significant fraction of the SDSS $`z6`$ quasars will have their luminosities enhanced by either gravitational lensing or beaming (Turner 1991; Wyithe & Loeb 2002a). Fan et al. (2003) discussed $`K`$’-band images of the three new $`z6`$ quasars and found that none of them show evidence for gravitational lensing. Our spectrum found that the equivalent width of the Mg ii line is typical of quasars which shows that the continuum luminosity is not strongly enhanced by beaming. Therefore, the available evidence suggests that there is minimal amplification of the intrinsic quasar luminosity.

It is interesting to consider the maximum luminosity that can be attained by an accreting black hole of this mass. This is usually considered to be the Eddington limit at which the outward radiation pressure equals the inward gravitational attraction. The Eddington limit for a black hole of mass $`3\times 10^9`$ $`\mathrm{M}_{}`$ is $`4\times 10^{40}`$ W. Using a bolometric correction of 7 to get the bolometric luminosity from $`\lambda L_{3000}`$ (Wandel, Peterson & Malkan 1999) we determine a bolometric luminosity for SDSS J1148+5251 of $`L=4\times 10^{40}`$ W. Therefore we find that the black hole in this quasar is radiating the maximal amount possible. This supports the often made assumption when relating quasar luminosities to black hole and halo masses that the highest redshift quasars are accreting at the Eddington limit (e.g. Fan et al. 2001; Wyithe & Loeb 2002b). It is not too surprising to find that SDSS J1148+5251 has the minimum size black hole given its luminosity. As mentioned above, the short cosmic time available for the growth of supermassive black holes by this redshift means black holes with such masses will only have arisen within rare collapsed peaks. In such conditions, there is likely to be a plentiful gas supply to fuel the quasar and the black hole is expected to be accreting close to its Eddington limit.

## 5 Conclusions

We have presented $`H`$ and $`K`$ spectra of the quasar SDSS J1148+5251 which currently has the highest known redshift. We clearly detect the broad Mg ii emission line. Fitting the location of this emission line we derive a redshift of $`z=6.41\pm 0.01`$. The equivalent width of this line is similar to that of lower redshift SDSS quasars, suggesting that the optical continuum is not strongly beamed. We obtain a virial estimate of the black-hole mass of $`M_{\mathrm{bh}}=3\times 10^9`$ $`\mathrm{M}_{}`$. The quasar is radiating at the Eddington luminosity, in line with the expectation that it is easier to fuel a supermassive black hole than create one at such an early epoch. Similar observations of other $`z>6`$ quasars would reveal if they are all radiating at the Eddington luminosity.

We thank the staff at the UKIRT for their excellent technical support. The United Kingdom Infrared Telescope is operated by the Joint Astronomy Centre on behalf of the U.K. Particle Physics and Astronomy Research Council. We thank the referee for comments. CJW thanks the National Research Council of Canada for support.