text:
# Physical limits to biochemical signaling

## I Introduction

A striking fact about biological systems is that single molecular events can have macroscopic consequences. The most famous example is of course the storage of genetic information in a single molecule of DNA, so that changes in the structure of this single molecule (mutations) can have effects on animal behavior and body plan from generation to generation . But there also are examples where the dynamics of individual molecular interactions can influence behavior on much shorter time scales. Thus we (and other animals) can see when a single molecule of rhodopsin in the rod cells of the retina absorbs a photon , and some animals can smell a single molecule of airborne odorant . Even if a single molecular event does not generate a specific behavior, it may still be that the reliability of behavior is limited by inevitable fluctuations associated with counting random molecular events. Thus the visual system has a regime where perception is limited by photon shot noise , and the reliability with which bacteria can swim up a chemical gradient appears to be limited by noise in the measurement of the gradient itself . It is an open question whether biochemical signaling systems within cells operate close to the corresponding counting noise limits.

The classical analysis of bacterial chemotaxis by Berg and Purcell provided a simple estimate and a clear intuitive picture of the noise in ‘measuring’ chemical concentrations. Their argument was that if we have a sensor with linear dimensions $`a`$, we expect to count an average of $`\overline{N}\overline{c}a^3`$ molecules when the mean concentration is $`\overline{c}`$. Each such measurement, however, is associated with a noise $`\delta N_1\sqrt{\overline{N}}`$. A volume with linear dimension $`a`$ can be cleared by diffusion in a time $`\tau _Da^2/D`$, so if we are willing to integrate over a time $`\tau `$ we should be able to make $`N_{\mathrm{meas}}\tau /\tau _D`$ independent measurements, reducing the noise in our estimate of $`N`$ by a factor of $`\sqrt{N_{\mathrm{meas}}}`$. The result is that our fractional accuracy in measuring $`N`$, and hence in measuring the concentration $`c`$ itself, is given by

$$\frac{\delta c}{\overline{c}}=\frac{\delta N}{\overline{N}}=\frac{1}{\sqrt{\overline{N}N_{\mathrm{meas}}}}=\frac{1}{\sqrt{Da\overline{c}\tau }}.$$
(1)
A crucial claim of Berg and Purcell is that this result applies when the sensor is a single receptor molecule, so that $`a`$ is of molecular dimensions, as well as when the sensor is the whole cell, so that $`a1\mu `$m. In particular, imagine a cell of radius $`R`$ which has receptor molecules of size $`a`$ on its surface. With just one receptor the limiting concentration resolution must be $`\delta c/c(Da\overline{c}\tau )^{1/2}`$, and if $`N_r`$ receptors are distributed sparsely over the cell surface we expect that they provide independent measurements, improving the resolution to $`\delta c/c(DN_ra\overline{c}\tau )^{1/2}`$. But as $`N_r`$ increases to the point where $`N_raR`$, this must saturate at $`\delta c/c(DR\overline{c}\tau )^{1/2}`$, presumably because of correlations among the concentration signals sensed by the different receptors.

The discussion by Berg and Purcell makes use of several special assumptions which we suspect are not required, and this leads to some clear questions:

* For interactions of a substrate with a single receptor, does Eq. (1) provide a general limit to sensitivity, independent of molecular details?
* Can we understand explicitly how correlations among nearby receptors result in a limit like Eq. (1), but with $`a`$ reflecting the size of the receptor cluster?
* Do the spatial correlations among nearby receptors have an analog in the time domain, so that there is a minimum averaging time required for noise reduction to be effective?

Finally, if we can establish Eq. (1) or its generalizations as a real limit on sensitivity for any signaling process, we would like to know if cells actually operate near this limit.

In most cases that we know about, biochemical signaling molecules are thought to interact with their receptors through some kinetic process which leads to the establishment of equilibrium between bound and unbound states. If this is the case, we can view the fluctuations in occupancy of a binding site as an example of thermal noise, and we can use the fluctuation–dissipation theorem rather than tracing through the consequences of different microscopic hypotheses about the nature of the interaction between signaling molecules and their targets. We begin with a simple example, to show that we can recover conventional results.

## II Binding to a single receptor

Consider a binding site for signaling molecules, and let the fractional occupancy of the site be $`n`$. If we do not worry about the discreteness of this one site, or about the fluctuations in concentration $`c`$ of the signaling molecule, we can write a kinetic equation

$$\frac{dn(t)}{dt}=k_+c[1n(t)]k_{}n(t).$$
(2)
This describes the kinetics whereby the system comes to equilibrium, and the free energy $`F`$ associated with binding is determined by detailed balance,

$$\frac{k_+c}{k_{}}=\mathrm{exp}\left(\frac{F}{k_BT}\right).$$
(3)
If we imagine that thermal fluctuations can lead to small changes in the rate constants, we can linearize Eq. (2) to obtain

$$\frac{d\delta n}{dt}=(k_+c+k_{})\delta n+c(1\overline{n})\delta k_+\overline{n}\delta k_{}.$$
(4)
But from Eq. (3) we have

$$\frac{\delta k_+}{k_+}\frac{\delta k_{}}{k_{}}=\frac{\delta F}{k_BT}.$$
(5)
Applying this constraint to Eq. (4) we find that the individual rate constant fluctuations cancel and all that remains is the fluctuation in the thermodynamic binding energy $`\delta F`$:

$$\frac{d\delta n}{dt}=(k_+c+k_{})\delta n+k_+c(1\overline{n})\frac{\delta F}{k_BT}.$$
(6)
Fourier transforming,

$$\delta n(t)=\frac{d\omega }{2\pi }\mathrm{exp}(i\omega t)\delta \stackrel{~}{n}(\omega ),$$
(7)
we can solve Eq. (6) to find the frequency dependent susceptibility of the coordinate $`n`$ to its conjugate force $`F`$,

$$\frac{\delta \stackrel{~}{n}(\omega )}{\delta \stackrel{~}{F}(\omega )}=\frac{1}{k_BT}\frac{k_+c(1\overline{n})}{i\omega +(k_+c+k_{})}$$
(8)
Now we can compute the power spectrum of fluctuations in the occupancy $`n`$ using the fluctuation–dissipation theorem:

$`\delta n(t)\delta n(t^{})`$ $`=`$ $`{\displaystyle \frac{d\omega }{2\pi }\mathrm{exp}[i\omega (tt^{})]S_n(\omega )}`$ (9)
$`S_n(\omega )`$ $`=`$ $`{\displaystyle \frac{2k_BT}{\omega }}\mathrm{}\left[{\displaystyle \frac{\delta \stackrel{~}{n}(\omega )}{\delta \stackrel{~}{F}(\omega )}}\right]`$ (10)
$`=`$ $`{\displaystyle \frac{2k_+c(1\overline{n})}{\omega ^2+(k_+c+k_{})^2}}.`$ (11)
It is convenient to rewrite this as

$$S_n(\omega )=(\delta n)^2\frac{2\tau _c}{1+(\omega \tau _c)^2},$$
(12)
where the total variance is

$`(\delta n)^2`$ $`=`$ $`{\displaystyle \frac{d\omega }{2\pi }S_n(\omega )}=k_BT{\displaystyle \frac{\delta \stackrel{~}{n}(\omega )}{\delta \stackrel{~}{F}(\omega )}}|_{\omega =0}`$ (13)
$`=`$ $`{\displaystyle \frac{k_+c(1\overline{n})}{k_+c+k_{}}}`$ (14)
$`=`$ $`\overline{n}(1\overline{n}),`$ (15)
and the correlation time is given by

$$\tau _c=\frac{1}{k_+c+k_{}}.$$
(16)
This is the usual result for switching in a Markovian way between two states; here it follows from the ‘macroscopic’ kinetic equations, plus the fact that binding is an equilibrium process.

The same methods can be used in the more general case where the concentration is allowed to fluctuate. Now we write

$$\frac{dn(t)}{dt}=k_+c(\stackrel{}{x}_0,t)[1n(t)]k_{}n(t),$$
(17)
where the receptor is located at $`\stackrel{}{x}_0`$, and

$$\frac{c(\stackrel{}{x},t)}{t}=D^2c(\stackrel{}{x},t)\delta (\stackrel{}{x}\stackrel{}{x}_0)\frac{dn(t)}{dt}.$$
(18)
Following the same steps as above, we find the linear response function

$`{\displaystyle \frac{\delta \stackrel{~}{n}(\omega )}{\delta \stackrel{~}{F}(\omega )}}`$ $`=`$ $`{\displaystyle \frac{k_+c(1\overline{n})}{k_BT}}{\displaystyle \frac{1}{i\omega [1+\mathrm{\Sigma }(\omega )]+(k_+\overline{c}+k_{})}}`$ (19)
$`\mathrm{\Sigma }(\omega )`$ $`=`$ $`k_+(1\overline{n}){\displaystyle \frac{d^3k}{(2\pi )^3}\frac{1}{i\omega +Dk^2}}`$ (20)
The “self–energy” $`\mathrm{\Sigma }(\omega )`$ is ultraviolet divergent, which can be traced to the delta function in Eq. (18); we have assumed that the receptor is infinitely small. A more realistic treatment would give the receptor a finite size, which is equivalent to cutting off the $`k`$ integrals at some (large) $`\mathrm{\Lambda }\pi /a`$, with $`a`$ the linear dimension of the receptor. If we imagine mechanisms which read out the receptor occupancy and average over a time $`\tau `$ long compared to the correlation time $`\tau _c`$ of the noise, then the relevant quantity is the low frequency limit of the noise spectrum. Hence,

$$\mathrm{\Sigma }(\omega D/a^2)\mathrm{\Sigma }(0)=\frac{k_+(1\overline{n})}{2\pi Da},$$
(21)
and

$$\frac{\delta \stackrel{~}{n}(\omega )}{\delta \stackrel{~}{F}(\omega )}=\frac{k_+\overline{c}(1\overline{n})}{k_BT}\left[i\omega \left(1+\frac{k_+(1\overline{n})}{2\pi Da}\right)+(k_+\overline{c}+k_{})\right]^1,$$
(22)
where $`\overline{c}`$ is the mean concentration. Applying the fluctuation–dissipation theorem once again we find the spectral density of occupancy fluctuations,

$$S_n(\omega )2k_+\overline{c}(1\overline{n})\frac{1+\mathrm{\Sigma }(0)}{\omega ^2(1+\mathrm{\Sigma }(0))^2+(k_+\overline{c}+k_{})^2}.$$
(23)
We note that the total variance in occupancy is unchanged since this is an equilibrium property of the system while coupling to concentration fluctuations serves only to change the kinetics.

Coupling to concentration fluctuations does serve to renormalize the correlation time of the noise,

$$\tau _c\tau _c[1+\mathrm{\Sigma }(0)].$$
(24)
The new $`\tau _c`$ can be written as

$$\tau _c=\frac{1\overline{n}}{k_{}}+\frac{\overline{n}(1\overline{n})}{2\pi Da\overline{c}},$$
(25)
so there is a lower bound on $`\tau _c`$, independent of the kinetic parameters $`k_\pm `$,

$$\tau _c>\frac{\overline{n}(1\overline{n})}{2\pi Da\overline{c}}.$$
(26)
Again, the relevant quantity is the low frequency limit of the noise spectrum,

$`S_n(\omega =0)`$ $`=`$ $`2k_+\overline{c}(1\overline{n}){\displaystyle \frac{1+\mathrm{\Sigma }(0)}{(k_+\overline{c}+k_{})^2}}`$ (27)
$`=`$ $`{\displaystyle \frac{2\overline{n}(1\overline{n})}{k_+\overline{c}+k_{}}}+{\displaystyle \frac{[\overline{n}(1\overline{n})]^2}{\pi Da\overline{c}}}.`$ (28)
If we average for a time $`\tau `$, then the root-mean-square error in our estimate of $`n`$ will be

$$\delta n_{\mathrm{rms}}=\sqrt{S_n(0)\frac{1}{\tau }},$$
(29)
and we see that this noise level has a minimum value independent of the kinetic parameters $`k_\pm `$,

$$\delta n_{\mathrm{rms}}>\frac{\overline{n}(1\overline{n})}{\sqrt{\pi Da\overline{c}\tau }}.$$
(30)
To relate these results back to the discussion by Berg and Purcell, we note that an overall change in concentration is equivalent to a change in $`F`$ by an amount equal to the change in chemical potential, so that $`\mathrm{\Delta }c/\overline{c}\mathrm{\Delta }F/k_BT`$. This means that there is an effective spectral density of noise in measuring $`c`$ given by

$$S_c^{\mathrm{eff}}(\omega )=\left(\frac{\overline{c}}{k_BT}\right)^2S_F(\omega ),$$
(31)
where the ‘noise force’ spectrum $`S_F(\omega )`$ is given by the fluctuation–dissipation theorem as

$$S_F(\omega )=\left|\frac{\delta \stackrel{~}{n}(\omega )}{\delta \stackrel{~}{F}(\omega )}\right|^2S_n(\omega )=\frac{2k_BT}{\omega }\mathrm{}\left[\frac{\delta \stackrel{~}{F}(\omega )}{\delta \stackrel{~}{n}(\omega )}\right].$$
(32)
In the present case we find that

$$S_c^{\mathrm{eff}}(\omega )=\frac{2\overline{c}^2}{k_+\overline{c}(1\overline{n})}\left[1+\frac{k_+(1\overline{n})}{2\pi Da}\right].$$
(33)
As before, the accuracy of a measurement which integrates for a time $`\tau `$ is set by

$$\delta c_{\mathrm{rms}}=\sqrt{S_c^{\mathrm{eff}}(0)\frac{1}{\tau }},$$
(34)
and we find again a lower bound which is determined only by the physics of diffusion,

$$\frac{\delta c_{\mathrm{rms}}}{\overline{c}}>\frac{1}{\sqrt{\pi Da\overline{c}\tau }}.$$
(35)
Note that this is (up to a factor of $`\sqrt{\pi }`$) exactly the Berg–Purcell result in Eq. (1).

## III Binding to multiple receptors

To complete the derivation of Berg and Purcell’s original results, we consider a collection of $`m`$ receptor sites at positions $`\stackrel{}{x}_\mu `$:

$`{\displaystyle \frac{dn_\mu (t)}{dt}}`$ $`=`$ $`k_+c(\stackrel{}{x}_\mu ,t)[1n_\mu (t)]k_{}n_\mu (t)`$ (36)
$`{\displaystyle \frac{c(\stackrel{}{x},t)}{t}}`$ $`=`$ $`D^2c(\stackrel{}{x},t){\displaystyle \underset{\mathrm{i}=1}{\overset{N}{}}}\delta (\stackrel{}{x}\stackrel{}{x}_\mu ){\displaystyle \frac{dn_\mu (t)}{dt}}.`$ (37)
From Eq. 37, we can write

$$\delta c(\stackrel{}{x_\nu },\omega )=\frac{i\omega \mathrm{\Lambda }}{2\pi ^2D}\delta \stackrel{~}{n}_\nu (\omega )+\frac{i\omega }{2\pi ^2}\underset{\mu \nu }{\overset{m}{}}\frac{\delta \stackrel{~}{n}_\mu (\omega )}{\left|\stackrel{}{x}\stackrel{}{x}_\mu \right|}_0^{\mathrm{}}\frac{k\mathrm{sin}\left(k\left|\stackrel{}{x}\stackrel{}{x}_\mu \right|\right)}{i\omega +Dk^2}𝑑k,$$
(38)
where $`\mathrm{\Lambda }`$ is the cut-off wave number; as before, the cut-off arises to regulate the delta function in Eq. 37, and is related to the size of the individual receptor. In the limit $`\left(\omega /D\right)^{1/2}1`$, we have

$$\delta c(\stackrel{}{x_\nu },\omega )=\frac{i\omega \mathrm{\Lambda }}{2\pi ^2D}\delta \stackrel{~}{n}_\nu (\omega )+\frac{i\omega }{4\pi D}\underset{\mu \nu }{\overset{m}{}}\frac{\delta \stackrel{~}{n}_\mu (\omega )}{\left|\stackrel{}{x}\stackrel{}{x}_i\right|},$$
(39)
and combining with Eq. 36 in Fourier space, we obtain

$`i\omega \delta \stackrel{~}{N}`$ $`=`$ $`\left[\left(k_+\overline{c}+k_{}\right){\displaystyle \frac{i\omega \mathrm{\Lambda }k_+(1\overline{n})}{2\pi ^2D}}\right]\delta \stackrel{~}{N}`$ (41)
$`+{\displaystyle \frac{i\omega k_+(1\overline{n})}{4\pi D}}{\displaystyle \underset{\nu =1}{\overset{m}{}}}{\displaystyle \underset{\mu \nu }{}}\delta \stackrel{~}{n}_\mu {\displaystyle \frac{1}{\left|\stackrel{}{x_\mu }\stackrel{}{x}_\nu \right|}}+mk_+\overline{c}(1\overline{n})\left({\displaystyle \frac{\delta \stackrel{~}{F}}{k_BT}}\right).`$
where we have defined $`\delta \stackrel{~}{N}(\omega )=_{\mu =1}^m\delta n_\mu (\omega )`$, and assumed the steady state fractional occupancies to be independent of the receptor site, $`\overline{n}_\mu =\overline{n}=k_+\overline{c}/\left(k_+\overline{c}+k_{}\right)`$ .

If we consider receptor cluster geometries such that the innermost sum is independent of $`\stackrel{}{x}_\nu `$, we can write the sum as

$$\underset{\nu =1}{\overset{m}{}}\underset{\mu \nu }{}\delta \stackrel{~}{n}_\mu \frac{1}{\left|\stackrel{}{x_\mu }\stackrel{}{x}_\nu \right|}=\varphi (m)\delta \stackrel{~}{N}$$
(43)
where

$$\varphi (m)=\underset{\mu =2}{\overset{m}{}}\frac{1}{\left|\stackrel{}{x_\mu }\stackrel{}{x}_1\right|}.$$
(44)
From the fluctuation–dissipation theorem, we find the spectrum of $`\delta \stackrel{~}{F}`$ and convert that to an equivalent concentration error as in Eq. 33:

$$\frac{\delta c_{\mathrm{rms}}}{\overline{c}}>\frac{1}{\sqrt{\pi D\overline{c}\tau }}\left(\frac{\mathrm{\Lambda }}{m\pi }+\frac{\varphi (m)}{2m}\right)^{1/2}.$$
(45)
For example, for receptors of radius $`b`$ uniformly distributed around a ring of radius $`a>b`$, we have $`\varphi (m)=mg_0/a`$, where $`g_0`$ is a geometric factor of order unity, and

$$\frac{\delta c_{\mathrm{rms}}}{\overline{c}}>\frac{1}{\sqrt{\pi D\overline{c}\tau }}\left(\frac{1}{mb}+\frac{g_0}{2a}\right)^{1/2}.$$
(46)
In summary, we find that the simple formula in Eq. (1) really does provide a general limit on the precision of concentration measurements by sensors of linear dimension $`a`$, at least in those cases where the interactions between the receptor and its ligand are passive. Further, there is a minimum level of receptor occupancy noise from Eq. (30), and a minimum correlation time from Eq. (26). Let us look at two examples to see how these limits compare with the performance of real cellular signaling mechanisms.

## IV Physical examples

### A Regulation of gene expression in bacteria

Expression of genes is controlled in part by the occupancy of promoter sites adjacent to the regions of DNA which code for protein . We thus can view gene expression as a readout mechanism for sensing promoter site occupancy, or even as a sensor for the concentration of the transcription factor proteins which bind to the promoter site. In a bacterium like E. coli, transcription factors are present in $`N_{\mathrm{TF}}100`$ copies in a cell of volume of $`1\mu \mathrm{m}^3`$ . If the transcription factor is a repressor then gene expression levels are determined by $`1n`$, while if it is an activator then expression is related to $`n`$; because $`\delta n_{\mathrm{rms}}\overline{n}(1\overline{n})`$ \[Eq. (30)\], fractional fluctuations in either $`A=n`$ or $`A=1n`$ are determined by

$$\frac{\delta A}{\overline{A}}=(1\overline{A})\frac{1}{\sqrt{\pi Da\overline{c}\tau }}$$
(47)
Direct measurements of diffusion constants for small proteins in the E. Coli cytoplasm yield $`D3\mu \mathrm{m}^2/\mathrm{s}`$ . A promoter site itself has linear dimensions $`a3\mathrm{nm}`$, and putting these factors together we find the crucial combination of parameters $`\pi Da\overline{c}3\mathrm{s}^1`$. In particular this means that the fluctuations in occupancy of the promoter site, averaged over a time $`\tau `$, are given by

$$\frac{\delta A}{\overline{A}}>(0.1)(1\overline{A})\left(\frac{100}{N_{\mathrm{TF}}}\right)^{1/2}\left(\frac{30\mathrm{s}}{\tau }\right)^{1/2}$$
(48)
Recent experiments  indicate that E. Coli achieves $`10\%`$ precision in control of gene expression at small values of $`\overline{A}`$. For this performance to be consistent with the physical limits, the transcription machinery must therefore integrate the promoter site occupancy for times of order one minute, even assuming that the translation from occupancy to expression level itself is noiseless. This integration can be provided by the lifetime of the mRNA transcripts themselves, which is $`3\mathrm{min}`$ in prokaryotes .

### B Control of the flagellar motor by CheY

The output of bacterial chemotaxis is control of the flagellar rotary motor . The phosphorylated form of the signaling protein CheY (CheY–P) binds to the motor and modulates the probability of clockwise versus counterclockwise rotation . Recent measurements  show that the probability $`p`$ of clockwise rotation depends very steeply on the concentration $`c`$ of CheY–P,

$$p=\frac{c^h}{c^h+c_{1/2}^h},$$
(49)
with $`h10`$ and $`c_{1/2}3\mu \mathrm{M}`$. Motors switch between clockwise and counterclockwise rotation as a simple random telegraph process, and for $`cc_{1/2}`$ the switching frequency is $`f1.5\mathrm{s}^1`$. If we view the motor as a sensor for the internal messenger CheY, then the observed behavior of the motor determines an equivalent noise level of

$$\delta c_{rms}=\left(\frac{p}{c}\right)^1\sqrt{p(1p)}\left(\frac{\tau _0}{\tau }\right)^{1/2},$$
(50)
where $`\tau _0`$ is the correlation time of the motor state; for the simple telegraph model it can be shown that $`\tau _0=2p(1p)/f`$. Using Eq. (49) we find

$$\frac{\delta c_{rms}}{c}=\frac{1}{h}\sqrt{\frac{2}{f\tau }}.$$
(51)
Thus, for $`cc_{1/2}`$, a single motor provides a readout of CheY–P concentration accurate to $`10\%`$ within two seconds. Given the dimensions of the flagellar motor’s C ring, $`a45\mathrm{nm}`$, with $`m34`$ individual subunits to which the CheY-P molecules bind , from Eq. 46 we find

$$\frac{\delta c_{\mathrm{rms}}}{c}\frac{1}{22}\left(\frac{2\mathrm{s}}{\tau }\right)^{1/2},$$
(52)
where we have taken the size of the individual receptor binding site to be $`b1\mathrm{nm}`$, and $`D3\mu \mathrm{m}^2/\mathrm{s}`$ as above. Hence, the collection of receptors comprising the motor are able to measure the CheY-P concentration with $`5\%`$ precision within two seconds. This is in agreement with our earlier result obtained by observing the switching statistics of the motor to within a factor of two.

## V Concluding remarks

In conclusion, we have derived from statistical mechanics considerations the physical limits to which biological sensors that rely on the binding of a diffusing substrate can measure its concentration. Our approach complements and extends the classic work by Berg and Purcell. For a single receptor, we arrive at their earlier intuitive result, which states that the accuracy in measurement of concentration is limited by the noise associated with the arrival of discrete substrate molecules at the receptors. Our approach extends in a straightforward way to multiple receptors without relying on additional considerations; for this case, our result demonstrates more transparently the role of multiple receptors in improving the measurement accuracy, as well as the saturating limit in this improvement set by the receptor cluster size. Relevant internal or external signaling molecules are often present in low copy numbers, and their concentration in turn regulates downstream biochemical networks crucial to the cell’s functions. For two experimentally well-studied examples, we show that the cell’s performance is consistent with the counting noise limits in measuring the concentration of these signaling molecules.


format: text

text:
# Single and Binary Black Holes and their Influence on Nuclear Structure

## 1. Introduction

The effect of a supermassive black hole on its stellar surroundings depends, as so often in stellar dynamics, on how one imagines the system evolved to its present state. Collisional relaxation times are too long to have affected the stellar distribution in all but the densest nuclei (Faber et al. 1997), hence the structure and kinematics of nuclei are fossil relics of the interactions between stars and black holes. In the simplest scenario, the black hole grows by accreting gas on a time scale long compared with the orbital periods of the surrounding stars (Peebles 1972; Young 1980). This “adiabatic growth” model makes fairly definite predictions about the distribution of stars near the black hole, predictions which are consistent with the steep density cusps observed in faint galaxies but which can not explain the flatter profiles at the centers of bright galaxies. But galaxies merge, implying the formation of binary black holes (Begelman, Blandford & Rees 1980) which are efficient at displacing matter as they spiral together. The dynamics of black hole binaries in galactic nuclei are complex; among the unanswered questions are the long-term efficiency of stellar dynamical processes at extracting energy from a binary, and whether decay of black hole binaries ever stalls at separations too great for the efficient emission of gravitational waves. But the binary black hole model is at least crudely consistent with the observed dependence of nuclear structure on galaxy luminosity. This article summarizes theoretical work on the single and binary black hole models and suggests avenues for future progress.

## 2. Preliminaries

A black hole of mass $`M_{}`$ embedded in a galactic nucleus will strongly affect the motion of stars within a distance $`r=r_h`$, the “radius of influence.” A standard definition for $`r_h`$ is

$$r_h\frac{GM_{}}{\sigma ^2}10.8\mathrm{pc}\left(\frac{M_{}}{10^8M_{}}\right)\left(\frac{\sigma }{200\mathrm{km}\mathrm{s}^1}\right)^2$$
(1)
with $`\sigma `$ the 1D velocity dispersion of the stars at $`rr_h`$. This definition had its origin in the isothermal sphere model for galactic nuclei; $`\sigma `$ is independent of $`r`$ in such a model and $`r_h`$ is the radius at which the circular velocity around the black hole equals $`\sigma `$. We now know that nuclei are power laws in the stellar density, $`\rho r^\gamma `$, and that $`\gamma `$ can lie anywhere between $`0`$ and $`2.5`$ (Lauer, this volume). The velocity dispersion in a power-law nucleus is only constant if $`\gamma =0`$ or $`2`$, hence a definition like equation (1) is problematic. One alternative would be to define $`r_h`$ as the root of $`\sigma ^2(r)GM_{}/r=0`$. A simpler definition, which will be adopted in this article, is the radius at which the enclosed mass in stars is twice the black hole mass:

$$M_{}(r<r_h)=2M_{}.$$
(2)
This definition is exactly equivalent to equation (1) when $`\rho (r)=\sigma ^2/2\pi Gr^2`$, the singular isothermal sphere. For an arbitrary power law,

$$\rho (r)=\rho _0\left(\frac{r}{r_0}\right)^\gamma ,$$
(3)
equation (2) implies

$$r_h=r_0\left(\frac{3\gamma }{2\pi }\frac{M_{}}{\rho _0r_0^3}\right)^{1/(3\gamma )}.$$
(4)
Note that the “theorist’s” convention is adopted here, in which $`\gamma `$ is the power-law index of the space (not projected) density.

## 3. The Adiabatic Growth Model

If a black hole grows at the center of a stellar system through the accretion of gas, the stellar density in the core will also grow as the black hole’s gravity pulls in nearby stars (Peebles 1972; Young 1980). The change in the stellar density can be computed straightforwardly if it is assumed that the black hole grows on a time scale long compared with stellar orbital periods. This is reasonable, since even Eddington-limited accretion requires $`10^8`$ yr to double the black hole mass, and orbital periods throughout the region dominated by the black hole are $`\begin{array}{c}<\hfill \\ \hfill \end{array}10^6`$ yr. Under these assumptions, the adiabatic invariants $`𝐉`$ associated with the stellar orbits are conserved as the black hole grows and the phase-space density $`f`$ remains fixed when expressed in terms of the $`𝐉`$. Computing the final $`f`$ becomes a simple matter of expressing the final orbital integrals in terms of their initial values under the constraint that the adiabatic invariants remain fixed (Young 1980).

In spherical potentials, the adiabatic invariants are the angular momentum $`L`$ and the radial action $`I=2_r_{}^{r_+}\sqrt{2\left[E\mathrm{\Phi }(r)\right]L^2/r^2}`$, where $`\mathrm{\Phi }(r)`$ is the gravitational potential and $`r_\pm `$ are pericenter and apocenter radii. It may be shown (e.g. Lynden-Bell 1963) that orbital shapes remain nearly unchanged when $`L`$ and $`I`$ are conserved, implying that an initially isotropic velocity distribution $`f(E)`$ remains nearly isotropic after the black hole grows (though not exactly isotropic – see below). The final $`f`$ corresponding to an initially isotropic $`f`$ is then simply

$`f_f(E_f,L)`$ $`=`$ $`f_i(E_i,L)`$ (5)
$``$ $`f_f(E_f)`$
where $`E_f`$ is related to $`E_i`$ through the condition $`I_f(E_f,L)=I_i(E_i,L)`$.

While the form of the nuclear density profile before the black hole appeared is not known, $`N`$-body studies of structure formation suggest that power laws are generic (Power et al. 2003 and references therein) and observed nuclear density profiles are often well described as power laws even on scales $`rr_h`$. Setting $`\rho _ir^{\gamma _0}`$, $`\mathrm{\Phi }_ir^{2\gamma _0}`$ ($`0<\gamma _0<2`$), the initial distribution function becomes

$$f_i(E_i)E_i^\beta ,\beta =\frac{6\gamma _0}{2(2\gamma _0)}(0<\gamma _0<2).$$
(6)
To compute $`f_f(E_f)`$ we need a relation between $`E_f`$ and $`E_i`$; we restrict attention to the region within the black hole’s sphere of influence by setting $`E_f=v^2/2GM_{}/r`$. The radial action $`I(E,L)`$ in the power-law model can not be computed analytically for every $`(E,L)`$, but for certain orbits $`E_f(E_i)`$ has a simple form. For instance, circular orbits have $`I=0`$, and conservation of angular momentum implies $`r_iM_i(r_i)=r_fM_{}`$ or $`r_fr_i^{4\gamma _0}`$. Thus $`E_fr_f^1r_i^{\gamma _04}E_i^{(\gamma _04)/(2\gamma _0)}`$, or

$$E_i(E_f)^{(2\gamma _0)/(4\gamma _0)}.$$
(7)
The same relation turns out to be precisely correct for radial orbits as well and is nearly correct at intermediate eccentricities (Gondolo & Silk 1999). Thus we can write

$$f_f(E_f)=f_i(E_i)E_i^\beta (E_f)^\delta ,\delta =\frac{6\gamma _0}{2(4\gamma _0)}$$
(8)
and the final density profile within the sphere of influence of the black hole is

$`\rho _f(r)`$ $`=`$ $`{\displaystyle f_f(v)d^3v}{\displaystyle _{\mathrm{\Phi }(r)}^0}(E)^\delta \sqrt{E\mathrm{\Phi }(r)}𝑑E`$ (9)
$``$ $`r^\gamma ,\gamma =2+{\displaystyle \frac{1}{4\gamma _0}}.`$
For $`0<\gamma _0<2`$, $`\gamma `$ varies only between $`2.25`$ and $`2.5`$; the slope of the final density profile within $`r_h`$ is almost independent of $`\gamma _0`$.

The form of $`\rho _f(r)`$ at $`rr_h`$ must be computed numerically (e. g. Young 1980; Cipollina & Bertin 1994; Cipollina 1995; Quinlan, Hernquist & Sigurdsson 1995; Gondolo & Silk 1999). Figure 1 shows $`\rho _f(r)`$ when $`\rho _i(r)r^{\gamma _0}`$. Defining $`r_{cusp}`$ to be the radius at which the inner and outer power laws intersect, one finds

$$r_{cusp}=\alpha r_h,0.19\begin{array}{c}<\hfill \\ \hfill \end{array}\alpha \begin{array}{c}<\hfill \\ \hfill \end{array}0.22,0.5\gamma _01.5.$$
(10)
In early treatments of the adiabatic growth model (Peebles 1972; Young 1980), the black hole was assumed to grow inside of a constant-density isothermal core. The index of the power-law cusp that forms from this initial state is $`\gamma =1.5`$, compared with the limiting valure $`\gamma =2.25`$ as $`\gamma _00`$ in the power-law models. This difference can be traced to differences in the central density profile:

$$\rho _i(r)=\rho _0\times \left(1+C_1r+C_2r^2+\mathrm{}\right).$$
(11)
The isothermal model has $`C_1=0`$ (an “analytic core”) implying a phase space density that tends to a constant value at low energies. Other sorts of cores have $`C_10`$ and $`f`$ diverges at low energies; for instance, the core produced by setting $`\gamma =0`$ in $`\rho (r)=r^\gamma (1+r)^4`$ has $`f(E)[E\mathrm{\Phi }(0)]^1`$. In fact models with finite central $`\rho `$’s can be found that generate final cusp slopes anywhere in the range $`1.5\gamma 2.25`$ (Quinlan, Hernquist & Sigurdsson 1995). There is probably no way of ruling out an analytic core in the progenitor galaxy on the very small scales that are relevant to the later formation of a cusp, hence the adiabatic growth model is compatible with any final slope in the range $`1.5\begin{array}{c}<\hfill \\ \hfill \end{array}\gamma \begin{array}{c}<\hfill \\ \hfill \end{array}2.5`$. The upper limit could even be extended beyond $`2.5`$ if $`\gamma _0>2`$.

How do these predictions compare with the data? Observed luminosity profiles are well described as power laws at the smallest resolvable radii, and in the case of faint ellipticals, $`M_V\begin{array}{c}>\hfill \\ \hfill \end{array}20`$, the observed range of slopes is $`1.5\begin{array}{c}<\hfill \\ \hfill \end{array}\gamma \begin{array}{c}<\hfill \\ \hfill \end{array}2.5`$ (Lauer, this volume). This is precisely the range in $`\gamma `$ predicted by the adiabatic growth model. However in bright galaxies, $`\gamma `$ extends down to $`0`$. A natural interpretation is that the steep cusps in faint galaxies are a result of adiabatic black hole growth, while some additional mechanism, like mergers, has acted to modify the profiles in the brighter galaxies.

Some fine tuning is still required to reproduce the luminosity profiles of galaxies with $`1.5\begin{array}{c}<\hfill \\ \hfill \end{array}\gamma \begin{array}{c}<\hfill \\ \hfill \end{array}2`$. These intermediate slopes require a shallow, core-like initial profile, and if the initial core radius exceeded $`r_h`$, the final profile will exhibit an upward inflection at $`rr_h`$ (e.g. Figure 2 of Young 1980; Figure 1). Such inflections are rarely if ever seen; observed profiles have slopes that decrease smoothly inward. A way out is to require that the black hole mass exceeds the initial core mass so that its growth obliterates the core; alternatively, all nuclei with $`\gamma \begin{array}{c}<\hfill \\ \hfill \end{array}2`$ may have been the products of mergers.

An ingenious attempt to reconcile the adiabatic growth model with both steep and shallow cusps was made by van der Marel (1999). He postulated the existence of isothermal cores in the progenitor galaxies with core masses scaling as $`L^{1.5}`$, with $`L`$ the total galaxy luminosity. Since $`M_{}L`$, $`M_{}/M_{core}L^{0.5}`$ and the black holes in faint galaxies would grow to dominate their cores, producing a cusp profile that approximates the featureless power laws of faint galaxies. In bright galaxies, van der Marel argued that the upward inflection at $`r\begin{array}{c}<\hfill \\ \hfill \end{array}r_h`$ would be difficult to resolve; hence these galaxies would exhibit nearly unperturbed core profiles. Van der Marel’s model is intriguing, although the assumed relation between core size and galaxy luminosity is ad hoc, and the assumption of large, pre-existing cores does not fit naturally into any current model of galaxy formation.

The predictions of the adiabatic growth model can be altered if the initial galaxy is anistropic, non-spherical or rotating. The simplest such case to treat is a spherical nucleus containing only circular orbits, with or without rotation (Young 1980; Quinlan, Hernquist & Sigurdsson 1995; Ullio, Zhao & Kamionkowski 2001). The derivation given above applies without approximations to this case: if the initial density profile is a power law, $`\rho _ir^{\gamma _0}`$, the final profile is also a power law with $`\gamma =2+1/(4\gamma _0)`$, for any $`\gamma _0`$ in the range $`[0,3]`$. As shown above, the isotropic model yields approximately the same final slope for $`0<\gamma _0<2`$, suggesting that even extreme tangential anisotropies have little effect on the final profile. The effect of radial anistropies has apparently never been investigated; this is worth doing since real galaxies often show evidence for significant radial anisotropies (e.g. Kronawitter et al. 2000).

Black hole growth in axisymmetric nuclei with analytic cores was considered by Leeuwin & Athanassoula (2000). They found little dependence of the final cusp slope on the degree of flattening. Merritt & Quinlan (1998) grew black holes of various masses in a triaxial $`N`$-body model that was formed via gravitational collapse, and found $`\rho r^2`$ at $`r\begin{array}{c}<\hfill \\ \hfill \end{array}r_h`$. Their models evolved to axisymmetry after growth of the black hole, but it is now known that triaxiality can be maintained throughout the black hole’s sphere of influence (Poon & Merritt 2002, 2003). Triaxial potentials support a wide range of different orbit families and the effects of black hole growth on such models have not been examined in detail.

In the spherical geometry, adiabatic growth of a black hole induces a mild anisotropy in the stellar motions at $`r\begin{array}{c}<\hfill \\ \hfill \end{array}r_h`$ due to the slightly different ways that circular and eccentric orbits respond to the changing potential; the net effect is a decrease in the average orbital eccentricity (Young 1980; Goodman & Binney 1984; Quinlan, Hernquist & Sigurdsson 1995; Figure 1). If the progenitor galaxy is rotating, growth of the black hole tends to increase $`V_{rot}`$ more rapidly than $`\sigma `$ (Lee & Goodman 1989; Leeuwin & Athanassoula 2000), although again the effect is slight.

To summarize: the adiabatic growth model is limited in its ability to reproduce the full range of luminosity profiles observed in galactic nuclei. The model predicts power-law profiles at $`r\begin{array}{c}<\hfill \\ \hfill \end{array}r_h`$ with logarithmic slopes $`1.5\begin{array}{c}<\hfill \\ \hfill \end{array}\gamma \begin{array}{c}<\hfill \\ \hfill \end{array}2.5`$. This nicely brackets the range of slopes observed in the nuclei of galaxies fainter than $`M_V20`$. However slopes less than $`\gamma 1.5`$ are not naturally produced by the adiabatic growth model, and some fine tuning is required to avoid inflections in the profile at $`rr_h`$.

## 4. The Binary Black Hole Model

The adiabatic growth model was proposed (Peebles 1972) before the importance of galaxy interactions and mergers was appreciated. We now know that supermassive black holes have been present in at least some spheroids since redshifts of $`z6`$ (e.g. Fan et al. 2001), and we believe that most galaxies have experienced at least one major merger since that time; indeed the era of peak quasar activity may coincide with the era of galaxy assembly via mergers (e.g. Caviliere & Vittorini 2000). If a nucleus forms via the merger of two galaxies containing pre-existing black holes, the net effect on the nuclear density profile is roughly the opposite of what the adiabatic growth model predicts: the black holes displace matter as they spiral into the center (Figure 2). This is a natural way to account for the shallow nuclear profiles in bright galaxies. The process may be understood as a sort of dynamical friction, with the “heavy particles” (the black holes) transferring their kinetic energy to the “light particles” (the stars). However most of the energy transfer takes place after the two black holes have come within each other’s sphere of influence, and in this regime, the interaction with the background is dominated by another mechanism, the gravitational slingshot (Saslaw, Valtonen & Aarseth 1974). The massive binary ejects passing stars at high velocities, removing them from the nucleus and simultaneously increasing its binding energy.

We begin by reviewing the dynamics of massive binaries in fixed, homogeneous backgrounds, then consider the more difficult problem of a binary located at the center of an inhomogeneous and evolving galaxy.

### 4.1. Dynamics of Massive Binaries

Consider a binary system consisting of two black holes with masses $`M_1`$ and $`M_2`$, where $`pM_2/M_11`$ and $`M_{12}M_1+M_2`$. Let $`a`$ be the semi-major axis of the Keplerian orbit and $`e`$ the orbital eccentricity. The binding energy of the binary is

$$|E|=\frac{GM_1M_2}{2a}=\frac{G\mu M_{12}}{2a},$$
(12)
with $`\mu =M_1M_2/M_{12}`$ the reduced mass. The relative velocity of the two black holes, assuming a circular orbit, is

$$V_{bin}=\sqrt{\frac{GM_{12}}{a}}=658\mathrm{km}\mathrm{s}^1\left(\frac{M_{12}}{10^8M_{}}\right)^{1/2}\left(\frac{a}{1\mathrm{pc}}\right)^{1/2}.$$
(13)
A binary is called “hard” when its binding energy per unit mass, $`|E|/M_{12}=G\mu /2a`$, exceeds $`\sigma ^2`$. For concreteness, a binary will here be called hard if its separation falls below $`a_h`$, where

$$a_h\frac{G\mu }{4\sigma ^2}0.27\mathrm{pc}\left(1+p\right)^1\left(\frac{M_2}{10^7M_{}}\right)\left(\frac{\sigma }{200\mathrm{km}\mathrm{s}^1}\right)^2.$$
(14)
Other definitions of $`a_h`$ are possible (e. g. Hills 1983; Quinlan 1996).

Stars passing within a distance $`3a`$ of the center of mass of a hard binary undergo a complex interaction with the two black holes, followed almost always by ejection at velocity $`\sqrt{\mu /M_{12}}V_{bin}`$ (Saslaw, Valtonen & Aarseth 1974). Each ejected star carries away energy and angular momentum, causing the semi-major axis, eccentricity and orientation of the binary to change and the local density of stars to drop. If the stellar distribution is assumed fixed far from the binary and if the contribution to the potential from the stars is ignored, the rate at which these changes occur can be computed by carrying out scattering experiments of massless stars against a binary whose orbital elements remains fixed during each interaction.

The results of the scattering experiments can be summarized via a set of dimensionless coefficients $`H,J,K,L,\mathrm{}`$ which define the mean rates of change of the parameters characterizing the binary and the stellar background (Hills & Fullerton 1980; Roos 1981; Hills 1983, 1992; Baranov 1984; Mikkola & Valtonen 1992; Quinlan 1996; Merritt 2001, 2002). These coefficients are functions of the binary mass ratio, eccentricity and hardness but are typically independent of $`a`$ in the limit that the binary is very hard. The hardening rate of the binary is given by

$$\frac{d}{dt}\left(\frac{1}{a}\right)=H\frac{G\rho }{\sigma }$$
(15)
with $`\rho `$ the density of stars at infinity. The mass ejection rate is

$$\frac{dM_{ej}}{d\mathrm{ln}(1/a)}=JM_{12}$$
(16)
with $`M_{ej}`$ the mass in stars that escape the binary. The rate of change of the binary’s orbital eccentricity is

$$\frac{de}{d\mathrm{ln}(1/a)}=K.$$
(17)
The diffusion coefficient describing changes in the binary’s orientation is

$$\mathrm{\Delta }\vartheta ^2=L\frac{m_{}}{M_{12}}\frac{G\rho a}{\sigma }$$
(18)
with $`m_{}`$ the stellar mass. Additional coefficients describe the rate of diffusion of the binary’s center of mass, or “Brownian motion” (Merritt 2001).

The binary hardening coefficient $`H`$ reaches a constant value of $`16`$ in the limit $`aa_h`$, with a weak dependence on $`M_2/M_1`$ (Hills 1983; Mikkola & Valtonen 1992; Quinlan 1996). In a fixed background, equation (15) therefore implies that a hard binary hardens at a constant rate:

$$\frac{1}{a(t)}\frac{1}{a_h}H\frac{G\rho }{\sigma }\left(tt_h\right),tt_h,a(t_h)=a_h.$$
(19)
If the supply of stars remains steady, hardening continues at this rate until the components of the binary come close enough together that the emission of gravitational radiation is important. In this regime, gravity wave coalescence takes place in a time:

$`t_{gr}`$ $`=`$ $`{\displaystyle \frac{5}{256F(e)}}{\displaystyle \frac{c^5}{G^3\mu M_{12}^2}}a^4,`$
$`F(e)`$ $`=`$ $`(1e^2)^{7/2}\left(1+{\displaystyle \frac{73}{24}}e^2+{\displaystyle \frac{37}{96}}e^4\right)`$ (20)
(Peters 1964). Coalescence in a time $`t_{gr}`$ occurs when $`a=a_{gr}`$, where

$$\frac{a_h}{a_{gr}}75F^{1/4}\frac{p^{3/4}}{\left(1+p\right)^{3/2}}\left(\frac{\sigma }{200\mathrm{km}\mathrm{s}^1}\right)^{7/8}\left(\frac{t_{gr}}{10^9\mathrm{yr}}\right)^{1/4}.$$
(21)
The $`M_{}\sigma `$ relation has been used to express $`M_{12}`$ in terms of $`\sigma `$. For mass ratios $`p`$ of order unity and $`e0`$, the binary must decay by a factor of $`10^2`$ in order for gravitational radiation to induce coalescence in a time shorter than $`10^9`$ yr. Less decay is required if the binary is eccentric or if $`M_2M_1`$.

If the binary manages to shrink by such a large factor, the damage done to its stellar surroundings will be considerable. The mass ejected by the binary in decaying from $`a_h`$ to $`a_{gr}`$ is given by the integral of equation (16):

$$M_{ej}=M_{12}_{a_{gr}}^{a_h}\frac{J(a)}{a}𝑑a.$$
(22)
Figure 1.2 shows $`M_{ej}`$ as a function of the mass ratio $`M_2/M_1`$ for $`\sigma =200`$ km s<sup>-1</sup> and various values of $`t_{gr}`$. The mass ejected in reaching coalescence is of order $`M_{12}`$ for equal-mass binaries, and several times $`M_2`$ when $`M_2M_1`$. A black hole that grew to its current size through a succession of mergers should therefore have displaced a few times its own mass in stars. If this mass came mostly from stars that were originally in the nucleus, the density within $`r_h`$ would drop drastically and the hardening would slow. Without some way of replenishing the supply of stars, decay could stall at a separation much greater than $`a_{gr}`$. This is the “final parsec problem”: how to avoid stalling and bring the black holes from their separation of $`1`$ pc, when they first form a hard binary, to $`10^2`$ pc, where gravity wave emission is efficient.

Changes in the binary’s orbital eccentricity (equation 17) are potentially important because the gravity wave coalescence time drops rapidly as $`e1`$ (equation 20). For a hard binary, scattering experiments give $`K(e)K_0e(1e^2)`$, with $`K_00.5`$ for an equal-mass binary (Mikkola & Valtonen 1992; Quinlan 1996). The dependence of $`K`$ on $`M_2/M_1`$ is not well understood. The implied changes in $`e`$ as a binary decays from $`a=a_h`$ to $`a_{gr}`$ are modest, $`\mathrm{\Delta }e\begin{array}{c}<\hfill \\ \hfill \end{array}0.2`$, for all initial eccentricities. The rms change in the orientation of the binary’s spin axis (equation 18) is $`\delta \theta 2(m_{}/M_{12})^{1/2}\mathrm{log}^{1/2}(a_h/a)`$, which is of order one degree or less if $`m_{}M_{}`$. Reorientations of the binary’s angular momentum affect the spin direction of the coalesced black hole and the direction of any associated jet (Merritt 2002).

### 4.2. Binary Black Holes in Galaxies

The scattering experiments summarized above treat the binary’s environment as fixed and homogeneous. In reality, the binary is embedded at the center of an inhomogeneous and evolving galaxy, and the supply of stars that can interact with it is limited.

In a fixed spherical galaxy, stars can interact with the binary only if their pericenters lie within $`\times a`$, where $``$ is of order unity. Let $`L_{lc}=a\sqrt{2\left[E\mathrm{\Phi }(a)\right]}\sqrt{2GM_{12}a}`$, the angular momentum of a star with pericenter $`a`$. The “loss cone” is the region in phase space defined by $`LL_{lc}`$. The mass of stars in the loss cone is

$`M_{lc}(a)`$ $`=`$ $`m_{}{\displaystyle 𝑑E_0^{L_{lc}}𝑑LN(E,L^2)}`$ (23)
$`=`$ $`m_{}{\displaystyle 𝑑E_0^{L_{lc}^2}𝑑L^24\pi ^2f(E,L^2)P(E,L^2)}`$
$``$ $`8\pi ^2GM_{12}m_{}a{\displaystyle 𝑑Ef(E)P_{rad}(E)}.`$
Here $`P`$ is the orbital period; in the final line, $`f`$ is assumed isotropic and $`P`$ has been approximated by the period of a radial orbit of energy $`E`$. An upper limit to the mass that is available to interact with the binary is $`M_{lc}(a_h)`$, the mass within the loss cone when the binary first becomes hard; this is an upper limit since some stars that are initially within the loss cone will “fall out” as the binary shrinks. Assuming a singular isothermal sphere for the stellar distribution, $`\rho r^2`$, and taking the lower limit of the energy integral to be $`\mathrm{\Phi }(a_h)`$, equation (23) implies

$$M_{lc}(a_h)3\mu .$$
(24)
We can compute the change in $`a`$ that would result if the binary interacted with this entire mass, by using the fact the mean energy change of a star interacting with a hard binary is $`3G\mu /2a`$ (Quinlan 1996). Equating the energy carried away by stars with the change in the binary’s binding energy gives

$$\frac{3}{2}\frac{G\mu }{a}dM\frac{GM_1M_2}{2}d\left(\frac{1}{a}\right)$$
(25)
or

$$\mathrm{ln}\left(\frac{a_h}{a}\right)3\frac{\mathrm{\Delta }M}{M_{12}}\frac{9\mu }{M_{12}}9\frac{p}{(1+p)^2},pM_2/M_1$$
(26)
if $`\mathrm{\Delta }M`$ is equated with $`M_{lc}`$. Only for very low mass ratios ($`M_2\begin{array}{c}<\hfill \\ \hfill \end{array}10^3M_1`$) is this decay factor large enough to satisfy equation (21), but the time required for such a small black hole to reach the nucleus is likely to exceed a Hubble time (Merritt 2000). Hence even under the most favorable assumptions, the binary would not be able to interact with enough mass to reach gravity-wave coalescence.

But the situation is even worse than this, since not all of the mass in the loss cone will find its way into the binary. The time scale for the binary to shrink is comparable with stellar orbital periods, and some of the stars with $`r_{peri}a_h`$ will only reach the binary after $`a`$ has fallen below $`a_h`$. We can account for the changing size of the loss cone by writing

$`{\displaystyle \frac{dM}{dt}}`$ $`=`$ $`{\displaystyle _{E_0(t)}^{\mathrm{}}}{\displaystyle \frac{1}{P(E)}}{\displaystyle \frac{dM_{lc}}{dE}}𝑑E`$ (27)
$`=`$ $`8\pi ^2GM_{12}m_{}a(t){\displaystyle _{E_0(t)}^{\mathrm{}}}f_i(E)𝑑E,`$
where $`M(t)`$ is the mass in stars interacting with the binary and $`f_i(E)`$ is the initial distribution function; setting $`P(E_0)=t`$ reflects the fact that stars on orbits with periods less than $`t`$ have already interacted with the binary and been ejected. Combining equations (25) and (27),

$$\frac{d}{dt}\left(\frac{1}{a}\right)24\pi ^2Gm_{}_{E_0(t)}^{\mathrm{}}f_i(E)𝑑E.$$
(28)
Solutions to equation (28) show that a binary in a singular isothermal sphere galaxy stalls at $`a_h/a2.5`$ for $`M_2=M_1`$, compared with $`a_h/a10`$ if the full loss cone were depleted (equation 26). In galaxies with shallower central cusps, decay of the binary would stall at even greater separations.

It is therefore entirely possible that uncoalesced binaries exist at the centers of many galaxies, particularly galaxies like large ellipticals with low central densities. But there is circumstantial evidence that long-lived black hole binaries are rare. Some radio galaxies show clear evidence of a recent flip in the black hole’s spin orientation, as would occur when two black holes coalesce (Dennett-Thorpe et al. 2002), and their numbers are consistent with a coalescence rate that is roughly equal to the galaxy merger rate (Merritt & Ekers 2002). The almost complete lack of correlation between jet directions in Seyfert galaxies and the angular momenta of their disks (Ulvestad & Wilson 1984; Kinney et al. 2000) also suggests that black hole coalescences were common in the past. There are few if any “smoking gun” detections of binary black holes among AGN (e.g. Halpern & Eracleous 2000); the best, but still controversial, case is OJ 287 (Pursimo et al. 2000).

Below are discussed some additional mechanisms that have been proposed for extracting energy from binary black holes and hastening their decay. When the agent interacting with the binary is stars, continued decay of the binary implies continued destruction of the stellar cusp. However for most of the mechanisms discussed below, the detailed effects on the stellar distribution have yet to be worked out.

#### 4.2.1 Scattering of stars into the loss cone

Destruction of a pre-existing stellar cusp generates strong gradients in the phase-space density at $`LL_{lc}`$, the angular momentum of an orbit that lies just outside the loss cone. A small perturbation can deflect a star on such an orbit into the loss cone. This process has been studied in detail in the context of gravitational scattering of stars into the tidal disruption sphere of a single black hole (Lightman & Shapiro 1977; Cohn & Kulsrud 1978). Once a steady-state flow of stars into the loss cone has been set up, the distribution function near $`L_{lc}`$ has the form:

$$f(E,L)\frac{1}{\mathrm{ln}(1/R_{lc})}\overline{f}(E)\mathrm{ln}\left(\frac{R}{R_{lc}}\right)$$
(29)
(Lightman & Shapiro 1977), where $`RL^2/L_c^2(E)`$, $`L_c(E)`$ is the angular momentum of a circular orbit of energy $`E`$, and $`\overline{f}`$ is the distribution function far from the loss cone, assumed to be isotropic. The mass flow into the central object is $`m_{}(E)𝑑E`$, where

$$(E)dE=4\pi ^2L_c^2(E)\left\{\frac{dr}{v_r}\underset{R0}{lim}\frac{(\mathrm{\Delta }R)^2}{2R}\right\}\frac{\overline{f}}{\mathrm{ln}(1/R_{lc})}dE.$$
(30)
The quantity in brackets is the orbit-averaged diffusion coefficient in $`R`$. Yu (2002) evaluated the contribution of two-body scattering to the decay rate of a massive binary and found that it was usually too small to overcome the stalling that occurs when the loss cone is first emptied.

Standard loss cone theory as applied by Yu (2002) assumes a quasi-steady-state distribution of stars in phase space near $`L_{lc}`$. This assumption is appropriate at the center of a globular cluster, where relaxation times are much shorter than the age of the universe, but is less appropriate for a galactic nucleus, where relaxation times almost always greatly exceed a Hubble time (e.g. Faber et al. 1997). The distribution function $`f(E,L)`$ immediately following the formation of a hard binary is approximately a step function,

$$f(E,L)\{\begin{array}{cc}\overline{f}(E),\hfill & L>L_{lc}\hfill \\ 0,\hfill & L<L_{lc}\text{,}\hfill \end{array}$$
(31)
much steeper than the $`\mathrm{ln}L`$ dependence in a collisonally relaxed nucleus (equation 29). Since the transport rate in phase space is proportional to the gradient of $`f`$ with respect to $`L`$, steep gradients imply an enhanced flux into the loss cone. The total mass in stars consumed by the binary can exceed the predictions of the standard model by factors of a few, implying greater cusp destruction and more rapid decay of the binary (Milosavljevic & Merritt 2003b). This time-dependent loss cone refilling might be particularly effective in a nucleus that that continues to experience mergers or accretion events, in such a way that the loss cone repeatedly returns to an unrelaxed state with its associated steep gradients.

#### 4.2.2 Re-ejection

Unlike the case of tidal disruption of stars by a single black hole, a star that interacts with a massive binary remains inside the galaxy and is available for further interactions. In principle, a single star can interact many times with the binary before being ejected from the galaxy or falling outside the loss cone; each interaction takes additional energy from the binary and hastens its decay. Consider a simple model in which a group of $`N`$ stars in a spherical galaxy interact with the binary and receive a mean energy increment of $`\mathrm{\Delta }E`$. Let the original energy of the stars be $`E_0`$. Averaged over a single orbital period $`P(E)`$, the binary hardens at a rate

$$\frac{d}{dt}\left(\frac{GM_1M_2}{2a}\right)=m_{}\frac{N\mathrm{\Delta }E}{P(E)}.$$
(32)
In subsequent interactions, the number of stars that remain inside the loss cone scales as $`L_{lc}^2a`$ while the ejection energy scales as $`a^1`$. Hence $`N\mathrm{\Delta }Ea^1a^1a^0`$. Assuming the singular isothermal sphere potential for the galaxy, one finds

$$\frac{a_h}{a(t)}1+\frac{\mu }{M_{12}}\mathrm{ln}\left[1+\frac{m_{}N\mathrm{\Delta }E}{2\mu \sigma ^2}\frac{tt_h}{P(E_0)}\right]$$
(33)
(Milosavljevic & Merritt 2003b). Hence the binary’s binding energy increases as the logarithm of the time, even after all the stars in the loss cone have interacted at least once with the binary. Re-ejection would occur differently in nonspherical galaxies where angular momentum is not conserved and ejected stars could miss the binary on their second passage. However there will generally exist a subset of orbits defined by a maximum pericenter distance $`\begin{array}{c}<\hfill \\ \hfill \end{array}a`$ and stars scattered onto such orbits can continue to interact with the binary.

#### 4.2.3 Chaotic loss cones

Loss cone dynamics are qualitatively different in non-axisymmetric (triaxial or bar-like) potentials, since a much greater number of stars may be on “centrophilic” orbits which take them near to the black hole(s) (Gerhard & Binney 1985). Triaxial models need to be taken seriously given recent demonstrations (Poon & Merritt 2002, 2003) that black hole nuclei can be stably triaxial and that the fraction of mass on centrophilic – typically chaotic – orbits can be large. The frequency of pericenter passages, $`r_{\mathrm{peri}}<d`$, for a chaotic orbit of energy $`E`$ in a triaxial black-hole nuclus is roughly linear in $`d`$, $`N(r_{\mathrm{peri}}<d)d\times A(E)`$ (Merritt & Poon 2003). The total rate at which stars pass within a distance $`a`$ of the massive binary is therefore

$$\frac{dM}{dt}aA(E)M_c(E)𝑑E$$
(34)
where $`M_c(E)`$ is the mass on centrophilic orbits. The implied feeding rate can be comparable to that in a spherical nucleus with a constantly-refilled loss cone, and orders of magnitude greater than in a loss cone that is re-supplied via star-star interactions (Merritt & Poon 2003). Even transient departures from axisymmetry, for instance during mergers, might result in substantial loss cone refilling due to this mechanism.

#### 4.2.4 Multiple black holes

If binary decay stalls, an uncoalesced binary may be present in a nucleus when a third black hole, or a second binary, is deposited there following a subsequent merger. The multiple black hole system that forms will engage in its own gravitational slingshot interactions, eventually ejecting one or more of the black holes from the nucleus (though probably not from the galaxy). This process has been extensively modelled assuming a fixed potential for the galaxy (e.g. Valtaoja, Valtonen & Byrd 1989; Mikkola & Valtonen 1990; Valtonen et al. 1994). The effect on the stellar distribution of $`N\begin{array}{c}>\hfill \\ \hfill \end{array}2`$, interacting black holes is not well understood, though $`N`$-body simulations with $`1020`$ massive particles and a “live” background show that the black holes displace $`10`$ times their own mass before being ejected (Merritt & Milosavljevic 2003). The separation and eccentricity of the dominant binary can change dramatically during each interaction and this may be an effective way to shorten the gravity wave coalescence time. In a wide, hierarchical triple, the eccentricity of the dominant binary oscillates through a maximum value of $`\sqrt{15\mathrm{cos}^2i/3}`$, $`|\mathrm{cos}i|<\sqrt{3/5}`$, with $`i`$ the mutual inclination angle (Kozai 1962). Blaes, Lee & Socrates (2002) estimate that the coalescence times in equal-mass, hierarchical triples can be reduced by factors of $`10`$ via the Kozai mechanism.

#### 4.2.5 Gas dynamics

If the inner $`1`$ pc of the nucleus contains a mass in gas comparable to $`M_2`$, torques from the gas will cause the orbit of the smaller black hole to decay in a time of order the gas accretion time (Syer & Clarke 1995; Ivanov, Papaloizou & Polnarev 1999). Given standard assumptions about accretion disk viscosities, the gas-dynamical decay rate would exceed that from gravity wave emission for $`a>a_{acc}`$, where

$$\frac{a_{acc}}{a_h}1\times 10^3\left(\frac{p}{0.1}\right)^{2/5}\left(\frac{\sigma }{200\mathrm{km}s^1}\right)^2$$
(35)
(Armitage & Natarajan 2002). If the orbit of the secondary is strongly inclined with respect to the accretion disk around the larger black hole, its passages through the disk could generate periodic outbursts, and this has been suggested as a model for the $`12`$ yr cycle of optical flaring observed in the blazar OJ 287 (Lehto & Valtonen 1996). Gas deposition of the required magnitude almost certainly occurred during the quasar epoch, although it is less clear that this mechanism is effective for galaxies in the current universe.

## 5. $`N`$-Body Studies

Unless great care is taken, $`N`$-body studies of binary black hole dynamics are unlikely to give an accurate picture of the evolution expected in real galaxies. This follows from the result (§ 4.2) that time scales for two-body scattering of stars into the binary’s loss cone are of order the Hubble time or somewhat longer in real galaxies. In $`N`$-body simulations, relaxation times are shorter by factors of $`N/10^{11}`$ than in real galaxies, hence the long-term evolution of the binary is likely to be dominated by spurious loss cone refilling, Brownian motion of the black holes and other noise-driven effects (Milosavljevic & Merritt 2003b). The stalling that is predicted in the absence of loss cone refilling in real galaxies (Figure 4) can only be reproduced via $`N`$-body codes if the mean field is artificially smoothed and the black hole binary is “nailed down” (e.g. Quinlan & Hernquist 1997).

$`N`$-body studies are most useful at characterizing the early stages of binary formation and decay, or simulating the disruptive effects of a single black hole on an infalling galaxy; indeed scattering experiments (§ 4.1) are almost useless in these regimes. Due to algorithmic limitations, most $`N`$-body studies (e.g. Ebisuzaki, Makino & Okumura 1991; Makino et al. 1993; Governato, Colpi & Maraschi 1994; Makino & Ebisuzaki 1996; Makino 1997; Nakano & Makino 1999a,b; Hemsendorf, Sigurdsson & Spurzem 2002) have been based on galaxy models with unrealistically large cores. $`N`$-body merger simulations using realistically dense initial conditions (Holley-Bockelmann & Richstone 2000; Merritt & Cruz 2001; Merritt et al. 2002) show that the black hole in the larger galaxy is efficient at tidally disrupting the steep cusp in the infalling galaxy, producing a remnant with only slightly higher central density than that of the giant galaxy initially. This result helps explain the absence of dense cusps in bright galaxies (Forbes, Franx & Illingworth 1995) and the persistence of the “core fundamental plane” in the face of mergers (Holley-Bockelmann & Richstone 1999).

Quinlan & Hernquist (1997) studied the evolution of a black hole binary inside cuspy models with $`\rho r^1`$ and $`\rho r^2`$ and a range of black hole masses and particle numbers, $`N2\times 10^5`$. Their $`N`$-body code was unable to simulate an actual merger and all of their detailed results were derived from initial conditions consisting of a single galaxy into which two “naked” black holes were symmetrically dropped. This artificial starting configuration produced substantial evolution of the cusp before the formation of the binary. The late evolution of the binary was found to be strongly dependent on $`N`$, due in part to spurious Brownian motion of the black hole particles. The cores that formed were characterized by strong velocity anisotropies, $`\sigma _t\sigma _r`$, due to the ejection of stars on eccentric orbits.

Milosavljevic & Merritt (2001) followed the evolution of cuspy ($`\rho r^2`$) galaxy models containing black holes, starting from pre-merger initial conditions and continuing until the binary separation had decayed a factor of $`10`$ below $`a_h`$. The initially steep nuclear cusps were converted to shallower, $`\rho r^1`$ profiles shortly after the black holes had formed a hard binary; thereafter the nuclear profile evolved slowly toward shallower slopes. The decay rate of the binary was found not to be strongly dependent on $`N`$, probably due to the fact that the loss cone was continuously refilled by two-body scattering in their collisional simulation (Milosavljevic & Merritt 2003b). The velocity anisotropies created during formation of the core were much milder than in the simulations of Quinlan & Hernquist (1997), similar in magnitude to the anisotropies predicted by the adiabatic growth model (Figure 1).

An $`N`$-body test of the time-dependent loss cone dynamics described in §4.2 would be possible using a collisional code with $`N\begin{array}{c}>\hfill \\ \hfill \end{array}10^7`$ (Milosavljevic & Merritt 2003b). Such large particle numbers demand a novel $`N`$-body algorithm coupled with special-purpose hardware like the GRAPE-6.

## 6. Observational Evidence for the Binary Black Hole Model

Figure 3 shows that a massive binary must eject of order its own mass in reaching a separation of $`a_{gr}`$ if $`M_2M_1`$, or several times $`M_2`$ if $`M_2M_1`$. These numbers should be interpreted with caution since: (1) binaries might not decay as far as $`a_{gr}`$, or the final stages of decay might be driven by some process other than energy exchange with stars; (2) the definition of “ejection” used in Figure 3 is escape of a star from an isolated binary, and does not take into account the confining effect of the nuclear potential; (3) the effect of repeated mergers on nuclear density profiles, particularly mergers involving very unequal-mass binaries, is poorly understood. Nevertheless, a clear prediction of the binary black hole model is that galactic mergers should result in the removal of a mass of order $`M_{12}`$ from the nucleus. As in the adiabatic growth model, we are handicapped in testing the theory by lack of knowledge of the primordial nuclear profiles. A reasonable guess is that all galaxies originally had steep, $`\rho r^2`$ density cusps, since these are generic in the low-luminosity ellipticals which are least likely to have been influenced by mergers, and since the adiabatic growth model predicts $`\rho r^2`$ (Figure 1).

The “mass deficit” is defined as the difference in integrated mass between the observed density profile, and a $`\rho r^{\gamma _0}`$ profile extrapolated inward from the break radius $`r_b`$:

$$M_{def}4\pi _0^{r_b}\left[\rho (r_b)\left(\frac{r}{r_b}\right)^{\gamma _0}\rho (r)\right]r^2𝑑r.$$
(36)
Milosavljevic et al. (2002) and Ravindranath, Ho & Filippenko (2002) computed $`M_{def}`$ in samples of “core”-profile elliptical galaxies. The former authors found a good proportionality between $`M_{def}`$ and $`M_{}`$, with $`M_{def}/M_{}1`$ for $`\gamma _0=1.5`$ and $`M_{def}/M_{}10`$ for $`\gamma _0=2`$ (Figure 4). These numbers are within the range predicted by the binary black hole model, given the uncertainties associated with the effects of multiple mergers. Ravindranath, Ho & Filippenko (2001) computed black hole masses using a shallower assumed slope for the $`M_{}\sigma `$ relation, $`M_{}\sigma ^{3.75}`$, and found a steeper, nonlinear dependence of $`M_{def}`$ on $`M_{}`$.

More rigorous tests of the binary black hole model will require a better understanding of the expected form of $`\rho (r)`$. As discussed above, while the best current $`N`$-body simulations suggest $`\rho r^1`$ following binary formation (Milosavljevic & Merritt 2001), the simulations are dominated by noise over the long term. If the decay stalls, the predicted density profile can be very different: a “hole” forms inside of $`3a_{stall}`$ (e.g. Figure 1 of Zier & Biermann 2001). Central minima may in fact have been seen in the luminosity profiles of a few galaxies (Lauer et al. 2002).

Other processes could result in energy exchange between binary black holes and stars. If the binary eventually coalesces, the gravitational radiation carries a linear momentum leading to a recoil of the coalesced black hole at a velocity $`v_{\mathrm{recoil}}10^210^3`$ km s<sup>-1</sup> (Bekenstein 1973; Fitchett 1983; Eardley 1983), and possibly even higher if the black holes were rapidly spinning prior to coalescence (Redmount & Rees 1989). Recoil velocities of this order would eject the black hole from the center of the nucleus and its subsequent infall would displace stars. Quantitative evaluation of this effect will require more accurate estimates of $`v_{\mathrm{recoil}}`$ based on fully general-relativistic calculations of black hole coalescence.

Equating $`M_{def}`$ with the mass removed by the binary black hole is only justified if galaxy mergers leave nuclear density profiles nearly unchanged in the absence of the gravitational slingshot. This is known to be the case in equal-mass mergers between galaxies with power-law cusps (Barnes 1999; Milosavljevic & Merritt 2001), though less is known about the effects of unequal-mass mergers. It is difficult to see how the break radii in galaxies with $`r_br_h`$ (e.g. NGC 3640, 4168) can be accounted for by the binary black hole model unless mergers sometimes induce changes in $`\rho (r)`$ well beyond the radius of influence of the black holes.

Another unsolved problem is the persistence of steep, $`\rho r^2`$ nuclear profiles in fainter galaxies; even these galaxies should have undergone at least one merger since their formation. One possibility is that faint ellipticals experienced their last major merger at a time when spheroids contained a much larger gas fraction than at present (e.g. Kauffmann & Haehnelt 2000), and that binary decay was driven by gas dynamics during this last merger.

A major focus of future work should be to calculate the evolution of $`\rho (r)`$ as predicted by the various scenarios for binary decay discussed in this article.

An abridged version of this article will appear in Carnegie Observatories Astrophysics Series, Vol. 1: “Coevolution of Black Holes and Galaxies,” edited by L. C. Ho (Cambridge: Cambridge Univ. Press). I thank M. Milosavljevic for his detailed reading of the manuscript. This work was supported by NSF grants AST 00-71099 and AST 02-06031, by NASA grants NAG5-6037 and NAG5-9046, and by grant HST-AR-08759 from STScI.


text:
# Experimental Phase-Sensitive Cloning

## Abstract

This paper has been withdrawn.


text:
# CNO abundances in the Damped Ly𝛼 clouds

## 1. Introduction

In the previous ESO workshop dedicated to the Production and Distribution of C,N,O elements held in Garching in May 1985 the information on CNO elements coming from absorption systems along distant QSO were just at the starting point and the first tentative results were reviewed by Max Pettini. We had to wait for one decade for the first N measurement in a damped Ly$`\alpha `$ system (DLA), even longer for O, and the time for C has probably not arrived yet.

DLAs are intervening absorption systems with large hydrogen column densities particularly appropriate for the determination of accurate chemical abundances. Typical abundances are $``$2.5$`<`$ \[Fe/H\]$`<`$ $``$1 with a moderate chemical evolution. The increase in metallicity is of the order of 0.3 dex per unit of decreasing redshift and suggests that DLAs do not trace the population of galaxies responsible for the bulk of star formation in the universe. A floor of abundances at about -2.5 is observed at all redshifts up to z=5.28 (Songaila et al 2002) quite suggestive of a prompt enrichment at very early phases. At redshif 5 the age of the universe is about one billion years and there is hope to find young systems where we can better recognize the imprints of the first generation of stars.

The observed abundances are gas phase abundances and in general suffer from possible dust and ionization corrections. The large neutral hydrogen column densities of the DLAs shield the clouds from the IGM field effectively. New evidence comes from the recently measured Ar which is particularly sensitive to ionization effects. The \[Ar/Si\] is often very close to solar or mildly deficient, which rules out significant ionization (Molaro et al 2001, Vladilo et al 2003). Dust is definitely present in the DLAs as revealed by a number of indications such as the systematic deficiency of Fe compared with Zn, the anti-correlation between refractory elements and the elemental column density, the correlations between \[Fe/Zn\] and H2 and the correlation between \[Si/Fe\] and \[Zn/Fe\] (Molaro 2001 and references therein). CNO elements are all mildly refractory and are virtually free from this complication (Savage & Sembach 1996).

The determination of the abundances of the CNO elements have proven to be very challenging. This because the most accessible lines for CNO are either saturated, despite the low metallicities of the DLAs, or they fall in spectral regions strongly contaminated by the neutral H absorptions of the Ly$`\alpha `$ forest. The OI 1302 Å line is generally nicely placed redwards the Ly $`\alpha `$ forest but is always found strongly saturated. Conversely the OI 1355 Å has never been detected and provides, together with the previous line, high and low bounds to the oxygen abundance. A similar situation holds for C where the only suitable lines are the CII resonance lines at 1334.53 Å and CII 1036.33 Å which are heavily saturated in all DLAs. A more favourable circumstance is found for NI which has several multiplets offering a relevant dynamical range in the line strenghts. However, the more accessible NI multiplets are $`\lambda \lambda `$ 1134.16, 1134.42, 1134.98 Å and 1199.55, 1200.22, 1200.71 Å which fall in the Ly$`\alpha `$ forest.

## 2. Oxygen

First measurements of O were performed at last through observations of OI 925 Å and OI 950 Å in the Ly-$`\alpha `$ forest, associated with the zabs=3.39 DLA towards QSO 0000-2621 (Molaro et al 2000). Partial contamination by hydrogen clouds is taken into account using a model for the absorber constructed on lines falling redwards the Ly $`\alpha `$ forest. This approach has been followed in the DLA at zabs=4.49 towards the J 0307-4945 (Dessaugues-Zavadsky et al (2001), and towards QSO 0347-381 through the measurement of a dozen OI lines (Levshakov et al 2002). Probably the best available measure due to the simplicity of the absorber and to the quality of Keck data has been realized by Prochaska et al (2001) in the DLA at z=2.844 towards QSO 1946+76 were several OI lines have been detected allowing an OI measure of log(OI)=14.819 $`\pm `$ 0.007. Pettini et al (2002) targeted few systems with low elemental abundances and low HI column density and succeded in measuring O from an almost unsaturated OI 1302 Å in a couple of DLAs. OI measurements have been provided for other DLAs towards QSO 1202-0725 (D’Odorico et al 2003) QSO 2059-360 (Dessaugues-Zavadsky et al 2003), QSO 2243-6031, QSO 1104-1805 (Lopez et al 1999, Lopez et al 2002). A total of a dozen of systems with some information on the O abundance is presently available.

The \[O/Zn\] and \[O/Fe\] are shown in Fig 1 with filled and empty circles respectively, and the solar O abundance taken from the recent revision by Holweger (2001). O and Zn do not require any dust correction and two \[O/Zn\] reveal a moderate enhancement in $`\alpha `$-elements at the level of 0.2 dex or so, while in one case the ratio is of $``$ 0.6 dex. There are also few \[O/Fe\]$``$ 0.2 which give stringent values. Most likely these systems are dust-free, as QSO 0000-2621 where \[O/Zn\]$``$\[O/Fe\]. If dust is present the effect is to increase the intrinsic Fe abundance and further lower the \[O/Fe\] ratio. Overall there is a low overabundance of O relative to Zn, which is consistent with previous indications obtained from \[S/Zn\] measurements or from \[Si/Fe\] once expected dust depletion is accounted for (Centurión et al 2000, Vladilo 2001).

In the sample there are two DLAs with O and S measured. The \[O/S\] ratio is $``$ 0 showing that there is no particular behaviour of O with respect to the other $`\alpha `$ elements. Si is generally available but it is a refractory element and may be partially depleted onto dust grains. In all but two DLAs the \[O/Si\] is found $``$ 0, which implies that Si is only moderately locked up into dust grains in DLAs. If this is the case it is rather interesting to consider the behaviour of Si with comparison of the undepleted element Zn. The average value of the 26 DLAs with both elements measured is $`<`$\[Si/Zn\]$`>`$=-0.07 $`\pm `$0.2, providing additional evidence of a moderate, if any, enhancement of $`\alpha `$-elements. Overall the low $`\alpha `$ over iron-peak element ratio appears as a prominent characteristic of the DLA. Cases of genuinely $`\alpha `$ enhancement exist such as the DLA towards J0347-3819 but they look as an exception rather than the rule (Levshakov et al 2002).

## 3. Carbon

In the survey by Prochaska et al 2001 there are 14 lower limits of C abundances, all obtained from the saturated CII 1335 Å line. CI is observed in a few cases but it traces very cool gas and being strongly ionised is not very informative of the total C. Levshakov et al (2002) and Lopez et al (2002) through detailed modelling of all absorption features in the effort to account for the saturation of the CII lines, claimed measurements with a reasonable error. In these two DLAs Zn is also measured and for both \[C/Zn\]$``$ 0.25, which suggests a rather normal behaviour. Highly supersolar \[Si/C\] values at the level of 0.5-1.5 dex are predicted by PopIII SNe, hypernovae or pair-instability SNe of very massive stars (Umeda and Nomoto 2001) for which there is no evidence in DLAs where \[Si/C\] is solar. Better prospects to measure C rely on the observations of very high redshift DLAs if abundances finally decrease, or in the class of sub-DLAs (cfr Dessauges-Zavadski, these proceedings). In the zabs=4.383 DLA towards BR 1202-0725 the CII line is found only moderately saturated and \[C/Fe\]=-0.03 ($`\pm `$0.15) and \[Si/C\]=0.15 $`\pm `$0.11 are derived (D’Odorico et al 2003). In the zabs=5.8 DLA towards SDSS 1044-0125 identified by Songaila and Cowie (2002), in which $`\mathrm{log}`$ N(HI)= 20.5 cm-2 and \[Fe/H\]=$``$2.65, the CII 1334 transition is unsaturated providing a \[C/Fe\]=-0.6, which suggests a possible drop in the C abundance at very high redshift. At redshift about 6 we may be within 100 million years from the epoch of star formation, which puts tight constraints to the C build-up processes (Molaro et al 2003).

## 4. Nitrogen

Among the CNO triplet N was the first to be observed. In their first attempt Pettini et al (1995) failed in detecting the element in the DLA towards QSO 2348-146 but established a very low upper limit at \[N/H\]$`<`$$``$3.15, which has been slightly improved with Keck observations to \[N/H\]$`<`$$``$3.28 (Prochaska and Wolfe 1999). This failure was followed by the detection of N in the DLA at zabs=3.39 towards QSO 0000-2621 and towards QSO 1331+17 with the N abundances \[N/H\]=-2.77 $`\pm `$ 0.17, and \[N/H\]= -2.73, respectively (Molaro et al 1996, Kulkarni et al 1996). The three DLAs showed a difference by $``$ 0.6 dex in their \[N/Si\] ratios and introduced since the beginning the notion of significant scatter in the \[N/Si\] abundances. It is just amazing that this became evident with the first two systems and that the first one to be observed was a Nitrogen poor system.

The early findings were confirmed and extended on larger samples, of which the most important for the number of objects are those of Centurión et al (1998), Lu et al (1998), Pettini et al (2002) and Prochaska et al (2002). In the following we use the complete and updated compilation of Centurión et al (2002, Table 9) who contributed with some new N observations as well. The sample of N observations or limits now comprises 33 systems, out of which 28 have at least one of O, Si or S measured and provide informative ratios.

In Fig 2 the \[N/$`\alpha `$\] are plotted versus \[$`\alpha `$/H\]. For the $`\alpha `$ element either of O, S or Si are taken in the order. As discussed in the previous section they trace each other rather well so that it does not matter very much which $`\alpha `$ element is used. The data for the BCG are also shown in the figure with smaller symbols.

The data show a peculiar pattern. Most of the sample form a plateau at \[N/Si\]$``$ $``$0.82 ($`\pm `$ 0.13), while few points show lower N/O ratios. The possibility of a gap between low and high values rather than a pure scatter has been pointed out by Prochaska et al (2002) who suggested a bimodal distribution. To the two measurements and to one upper limit with low N/O considered originally by Prochaska et al we add the value towards J0307-4945.. Dessauges-Zavadsky (2001), rejected by Prochaska et al (2002) for a possible ionization effect, and the most recent one towards QSO 2059-360 from Dessauges-Zavadsky et al (2002). Two new upper limits quite separate from the high values and consistent with a bimodal distribution are provided by Pettini et al ( 2002). The average of the four determinations is \[N/Si\]$``$ $``$1.45($`\pm `$ 0.05), strongly suggestive of a second plateau even with a narrower dispersion. The low N/O values are mostly observed at the lowest metallicities, but with an overlap over the higher plateau in the range -2.0 $`<`$\[O/H\]$`<`$ -1.5. In Fig 3 the plot of \[N/H\] versus \[Si/H\] clearly shows the primary-like behaviour of N in both plateau.

The high N/O values look as an extension of the BCG towards lower metallicities and give definitive evidence of a plateau, which was not so clear from the BCG only. The sharing of the same N/O between DLAs and BCG suggest that the two populations experienced similar chemical evolution.

N synthesis is thought to occur in Intermediate Mass Stars (IMS) (4 to 7 M$`_{}`$) through Hot Bottom Burning in the AGB phase. The nucleosynthesis by IMS implies a delay in the N production when compared to that of O as well as a sensitivity of N production to the low mass regime of the IMF. The general interpretation of the observed data is that N is primary at low metallicities and secondary at high metallicity (Henry et al 2000). The Henry et al model is able to reproduce the plateau by reducing the SFR and thus by preventing the rapid O build up by massive stars. For high SFR the elemental build up is fast and the system enters very soon in the secondary production regime without undergoing a primary behaviour at all. Thus the presence of an extended plateau in the N/O suggests low SFR in the DLA, which may also explain the failure of imaging DLA galaxies. The lower plateau is a new feature of DLA with no other counterparts. Izotov et al (1999) argued that the low N/Si were resulting from an inflated Si column density by additional Si formed in the ionized region. However, the recent direct O measurements provided five genuine low N/O obtained by means of truly oxygen observations, shown in Fig 2, which rule out such possibility.

### 4.1. About the origin of the low N/O plateau

Low N/O values have been interpreted in a context of a delayed N released model as the result of a relatively young age ( Molaro 2001, Pettini 2002,) or as being the system observed close to a recent episode of star formation (Lu et al 1998, Centurión et al 1998). According to Henry et al the necessary scale time to attain a significant N production by IMS is of at least 250 Myrs. It is therefore possible that DLAs are catched in the early formation times when the IMS had not yet expelled their N. However, as Prochaska et al pointed out one should expect to find values everywhere below the plateau and it is difficult to account for the fact that the low N/O show a characteristic value at \[N/O\]=$``$1.4. In the low N/O plateau the scatter is so small that these systems are unlikely in a transient phase towards the upper plateau.

A second possibility to explain the low N/O values is a variation in the IMF. A top-heavy IMF or an IMF truncated below a mass threshold of $``$ 5 M$`_{}`$ which exclude low mass stars has been suggested by Prochaska et al (2002) . However, it remains unclear what produces the change in the IMF. As discussed by Henry at this meeting it is rather appealing to think that it is related to the formation of PopIII stars

A further possibility is that we are dealing with variable N yields as a function of the metallicity since the low N/O values occur preferably at the lowest metallicities. This may affect the efficiency of the 3rd dredge-up or of the mass-loss rate along the AGB phase somewhat altering the characteristic N yields of the IMS. However, it is not clear how to account for the fact that both high and low N/O values occur at equal metallicities.

Fresh clues for the interpretation of the low N/O plateau can be obtained plotting the \[N/Si\] values with respect to \[N/H\] as shown in Fig 4. The distribution does not show a degeneracy with the metallicity any longer but reveal two different regimes with a transition occurring at \[N/H\]$``$ $``$3. The low N/O are in the region with \[N/H\] $`<`$ -3 and the high N/O for \[N/H\] $`>`$ -3. Since there is no obvious link between N yields and the amount of N content in the systems we think that what we see is most likely a temporal effect. In fact, since the O abundance is more sensitive to the SFR history while N depends from IMS evolutionary time-scales, N could be a better indicator of the age of system. Thus the two regimes reflect probably two epochs.

We propose that the higher plateau is produced by IMS after they had time to evolve and release their N, while the lower plateau is produced by more massive, namely M$`>`$ 8 M$`_{}`$, and relatively younger stars. At later times (or \[N/H\]$``$ -3) the contribution from the massive stars is overwhelmed by that from the IMS. In fact a hint of connection between the two plateaux can be recognized in the two DLAs which lie just below the primary plateau at \[N/Si\] $``$ -1.1 and \[N/H\]$``$ -2.8 as shown in Fig 4. No objects are expected to be found below the low N/O plateau and therefore it would be rather crucial to improve the three upper limits to the levels at which they can probe the existence of the plateau. We note that for the DLA at zabs=2.456 towards QSO 1409+095 the upper limit is only 0.17 dex above the plateau and in fact there is a hint of the presence of the NI 1200 Ålines (Pettini et al 2002, their Fig 2).

Standard models do not predict N synthesis by massive stars but the situation is far from being settled (Heger and Woosley 2002). Chieffi et al (2001) found that also Pop III stars of intermediate mass do experience thermal pulses and the third dredge-up arguing that stars with primordial composition are sources of C and N. If we take the models of Limongi and Chieffi (2002) for zero metallicity, computed with FRANEC code CF85, with no mix of the H shell and a mass loss parameter of $`\eta `$=6, and we integrate with a Salpeter IMF over the whole mass spectrum 4-80 M$`_{}`$ we get a \[N/O\]=-0.85, namely very close to the higher plateau. On the other hand the integration in the mass range between 15 and 80 M$`_{}`$ gives \[N/O\]=-1.5, which is very close to the lower plateau. The models for low metallicity stars with rotational mixing or overshooting of Meynet and Maeder (2002) predict N production for massive stars and can also reproduce the ratio observed in the lower plateau. The observations of DLAs might actually offer crucial guidelines for N production in massive stars.

The fraction of systems on the two plateaux and their relative abundances can be envisaged as two important tests to probe the origin of the two plateaux. Presenly the fraction of the systems forming the low N/O plateau is of 25% of the total. We can compute the expected ratio assuming that we need roughly 250 Myrs for significant N production and computing the age of the systems from the formation epoch to the observed redshift by means of a look- back time $`t`$= 0.538\[(1+z)/10\]-3/2 Gyrs. The formation epoch is unknown but to get the observed statistics we need to place it at z$``$ 4, which is quite conceivable. Few systems with low N/O have zabs of about 2-3 suggesting a rather continuous formation of the DLA. The case of QSO 1202-0725 with high N/O and zabs=4.38 is rather remarkable and pushes back the epoch of star formation considerably. If the synthesis of N observed by IMS require about 0.25-0.5 Gyrs then the time of star formation is shifted to a redshift greater than six.

An $`\alpha `$\- element enhancement over iron-peak element is expected in the case that massive stars are responsible for the low N/O plateau. Fe is produced considerably by Type Ia SNae and should lag O similarly to N. In Fig. 5 the \[Si/Zn\] and \[Si/Fe\] are plotted versus \[N/Si\] Unfortunately there are no \[Si/Zn\] values for DLAs in the low N/O plateau and the \[Si/Fe\] are between 0.2 and 0.4 which is rather common among the DLAs. In absence of dust these ratios would imply an $`\alpha `$ enhancement, but unfortunately there is no Zn information, and the test is not presently conclusive.

## 5.  Acknowledgments

It is a pleasure to acknowledge my collaborators, and in particular Miriam Centurión for allowing me to present here preliminary results of common work. Marco Limongi and Alessandro Chieffi are acknowledged for making the computation of N/O ratios with their model.

## References

Centurión M., Bonifacio P., Molaro P., & Vladilo G. 1998, ApJ, 509, 620

Centurión M., Bonifacio P., Molaro P., & Vladilo G. 2000, ApJ, 536, 540

Dessauges-Zavadsky M., D’Odorico S., McMahon R. G. et al. 2001, A&A, 370, 426

Dessauges-Zavadsky M. et al. 2002, in preparation

D’Odorico, V. et al. 2003, in preparation

Henry R. B. C., Edmunds M. G., & Köppen J. 2000, ApJ, 541, 660

Holweger H. 2001, in Solar and Galactic Composition, ed. R. F. Wimmer-Schweingruber, AIP Conf. Proceed., 598, 23

Izotov Y. I., & Thuan T. X. 1999, ApJ, 511, 639

Kulkarni V.P., Huang K., Green R.F. et al. 1996, MNRAS 279, 197

Levshakov S. A., Dessauges-Zavadsky M., D’Odorico S., and Molaro, P. 2002, ApJ, 565, 696

Limongi M., Chieffi M. 2002 PASA 19, 246

Lopez S., Reimers D., Rauch M., et al 1999, ApJ 513, 598

Lopez, S., Reimers D., D’Odorico S., Prochaska J. X. 2002, A&A 385, 778

Lu L., Sargent W. L. W., & Barlow T. A. 1998, AJ, 115, 55

Meynet G. & Maeder A. 2002, A&A, 390, 561

Molaro, P. D’Odorico, S., Fontana, A., Savaglio S., Vladilo, G., 1996, A&A 308,1

Molaro P., Bonifacio P., Centurión M. et al. 2000, ApJ, 541, 54

Molaro P., Levshakov S. A, D’Odorico S., et al 2001, ApJ, 549, 90

Molaro, P., 2001 in 17th IAP Astrophysics Colloquium, 307

Molaro, P., et al in preparation, 2003

Pettini M., Lipman K., & Hunstead R. W. 1995, ApJ, 451, 100

Pettini M., Ellison S., Bergeron J., & Petitjean P. 2002, A&A, 391, 21

Prochaska J. X., Wolfe A. M., Tytler D. et al. 2001, ApJS, 137, 21

Prochaska J.X., Henry R.B.C., O’Meara J.M. et al. 2002, PASP, 114, 933

Songaila, A. Cowie L.L 2002, AJ 123, 2183

Savage B. D. Sembach K. R. 1996, ARA&A 34, 279

Umeda H., Nomoto K., & Nakamura T. 2000, in The First Stars, eds. A. Weiss, T. Abel, and V. Hill, (Springer, Heidelberg), 150

Vladilo G. 2002a, ApJ, 569, 295

Vladilo G., Centurión M., Molaro M., Peroux C., & D’Odorico, V., 2003, in prep.