format: text

text:
# Gravitational Waves Probe the Coalescence Rate of Massive Black Hole Binaries

## 1 Introduction

Over the past decade, evidence has accumulated that Massive Black Holes111The literature has variable usage of the terms “Massive” and “Supermassive” to describe the $`10^{6\text{}9}`$ M black holes in the centers of galaxies (e.g., Hug, ????). For concision, we favor the former in this work. (MBHs) are ubiquitous in the spheroidal components of low-redshift galaxies (e.g., Magorrian et al., 1998), and that their masses are tightly correlated with the properties of the host galaxies. These initial studies noted that the MBH mass scaled with the mass of the spheroidal component of the parent galaxy as inferred from the luminosity although there was considerable scatter. More recently, a much tighter correlation has been observed between the MBH mass and the velocity dispersion in the spheroid (Ferrarese & Merritt, 2000; Gebhardt et al., 2000).

At higher redshift, evidence for the ubiquity of MBHs is more circumstantial. However, we do know that many distant galaxies harbor active nuclei for some portion of their life and, in turn, that the engine driving this activity is almost certainly gravitational accretion onto a central black hole. There are many outstanding questions regarding the era and duty cycle of AGN activity, the initial mass and mass evolution of MBHs, and the mass growth of MBHs during AGN episodes. Nonetheless, the existing evidence strongly favors the hypothesis that the present-day MBH population is the dormant remnant of the AGN population at high redshift.

During the last decade astrophysicists have also begun to understand the important role that mergers play in the history of galaxies. In the hierarchical clustering picture, mergers are by definition responsible for building up the mass of today’s galaxies. Numerical simulations have also shown that the merger history also determines the morphology of present-day galaxies: e.g., major mergers tend to produce elliptical galaxies. The strongest evidence for the mergers is the relatively flat density scaling laws in “some” spheroids (Milosavljević & Merritt, 2001). Simulations that reproduce this scaling actually calculate the merger of galaxy halos where the majority of the mass resides. The subsequent merger of the central galaxies within the halos is nearly guaranteed owing to the rapid process of dynamical friction for these dynamically soft stellar systems.

The ubiquity of MBHs in galaxies and the frequent growth of galaxies by merger leads to the obvious question about the fate of the pair of MBHs in the subsequent evolution of merged systems. The same dynamical friction process that brings together parent galaxies within merged halos will also lead to the sinking of the MBHs to the center of the merger daughter galaxy. Unfortunately, in the simplest case dynamical friction will “turn off” well before the MBH binary coalesces222We use the term merger to denote the joining of halo and galaxy pairs and the term coalescence for the joining of the members of a MBH binary to distinguish implicitly the two nominally sequential processes., or even before it reaches a regime where gravitational radiation losses will drive the binary’s evolution (Begelman et al., 1980). Various mechanisms have been suggested to evolve the MBH pair through this period.

If the pair can indeed evolve to the gravitational wave (GW) regime, it will over time lose angular momentum to gravitational radiation and eventually coalesce. The final inspirals of MBH pairs take about $`10^4`$ seconds, and are the most luminous gravitational wave events in the Universe. These events are one of the chief targets of the Laser Interferometer Space Antenna (LISA) GW interferometer satellite. In this paper, we concentrate on the quiescent evolution leading up to this state. Before the final inspiral, the GW amplitude is

$$h4.4\times 10^{17}M_8^{5/3}P_{\mathrm{yr}}^{2/3}D_{\mathrm{Gpc}}^1\frac{q}{(1+q)^2},$$
(1)
where $`h`$ is the strain or metric perturbation, $`M_8`$ is the total mass of the binary in units of $`10^8M_{}`$, $`q`$ is the mass ratio ($`q<1`$), $`P_{\mathrm{yr}}`$ is the observed GW period in years (the orbital period divided by 2 for the expected nearly-circular orbits), and $`D_{\mathrm{Gpc}}`$ is the distance in Gpc. The lifetime of the system is

$$t_{\mathrm{GW}}=1.1\times 10^6\mathrm{yr}M_8^{5/3}P_{\mathrm{yr}}^{8/3}\frac{(1+q)^2}{q}.$$
(2)
The precision of the rotation periods of millisecond period pulsars which is established via pulse arrival time measurements allows detection of the stochastic MBH-MBH background spectrum at nHz frequencies (Sazhin, 1978; Detweiler, 1979; Rajagopal & Romani, 1995). The detection process is very similar to that of the laser interferometers—passage of an electromagnetic signal through distorted space-time—except that in the case of pulsars no mirrors are needed as we have a distant precision clock, the rotating neutron star, sending us a periodic pulse train. Measurements of a single pulsar can place an upper limit on the presence of a spectrum of gravitational radiation (Kaspi et al., 1994; Lommen, 2001). In short, the limit on strain is given by the ratio of the timing precision and the measurement duration; e.g. 1 $`\mu `$s/10 yr $`3\times 10^{14}`$. While pulsar timing cannot detect individual events of the small amplitude given in equation (1), these measurements have a good prospect for detection of the stochastic background of GWs from the ensemble of these MBH-MBH coalescence events throughout the Universe333We refer to this as background radiation, although those interested in relic radiation from the first inflationary moments might choose to call this foreground radiation.. The GW perturbs pulse arrival times of a spatial array of pulsars in a correlated manner that is distinct from other known perturbations such as atomic time and ephemeris errors. A Pulsar Timing Array then acts as a nHz gravitational wave telescope capable of direct detection of the stochastic background radiation.

Questions remain. Are there dynamical processes that drive the MBH pair into the GW regime? Does this happen frequently? If so, can we observe the gravitational radiation via pulsar timing? If not, do we then see evidence of close binary MBHs “hung up” in the centers of many, or most, galaxies? If we observe neither the GW signal via precision pulsar timing nor inactive close pairs via high angular resolution studies, does this imply that some part of the paradigm—the ubiquity of black holes at all redshifts and the importance of mergers in galaxy evolution—is flawed? Some of these questions have been taken up by other groups over the years. In their important paper, Rajagopal & Romani (1995, hereafter RR), discussed several mechanisms driving the MBH pair into the GW regime, and calculated the GW spectrum under various assumptions. More recently, Gould & Rix (2000) follow up on the gas-dynamical mechanism for driving the MBH coalescence first mentioned by Begelman et al. (1980). Armitage & Natarajan (2002) explore the detailed interaction of the MBH binary with the accretion disk. Meanwhile, Milosavljević & Merritt (2001) and Yu (2001) have investigated the stellar-dynamical schemes. Finally, Menou et al. (2001) have investigated the effect that a time-dependent change in MBH “demographics” might have on the merger history of galaxies, and Phinney (2001) has shown an alternative method for calculating the GW spectrum for generic sources.

The advances in our understanding of MBHs along with steady progress in precision pulsar timing have led us to this paper. In §2 we first assemble the ingredients needed to model the stochastic GW background from MBH binaries, and in §3 we present the results of our calculations. In §4 we discuss the detection of gravitational radiation with a Pulsar Timing Array and the status of current experiments. In the concluding section we discuss the theoretical and observational prospects for a better understanding of the various ingredients to our calculations and measurements.

## 2 The Binary MBH Gravity Wave Spectrum

First, we need to define the terminology for a stochastic gravitational wave background spectrum that has been presented in the literature in a variety of ways (e.g., Burke (1975); Allen & Romano (1999); Maggiore (2000)). A single GW is denoted by the transverse-traceless part of the metric perturbation, $`h_{ab}(𝐱,t)=h_+(t)e_{ab}^+(𝐤)+h_\times (t)e_{ab}^\times (𝐤)`$, where plus signs and crosses refer to the two polarization degrees of freedom of a gravitational wave, and the $`e_{ab}^i`$ are basis tensors for the polarizations which are functions of the direction of propagation of the gravitational wave, $`𝐤`$. The two amplitudes combine with the two coordinate directions and the polarization position angle to yield five independent parameters of the single-wave metric perturbation. The stochastic background will excite both wave components equally and randomly $`|h_+|^2=|h_\times |^2`$ and $`h_+^{}h_\times =0`$, and all directions and position angles will be equally likely.

The total power spectral density of gravitational waves, $`S_h(f)`$, is defined by

$$\underset{P=+,\times }{}\stackrel{~}{h}_P(f)\stackrel{~}{h}_P^{}(f^{})=\frac{1}{2}\underset{a,b}{}\stackrel{~}{h}_{ab}(f)\stackrel{~}{h}_{ab}^{}(f^{})=\frac{1}{2}\delta (ff^{})S_h(f),$$
(3)
where $`\stackrel{~}{h}_P(f)`$ is the Fourier transform of the $`P=+,\times `$ polarization component of the metric strain tensor at frequency $`f`$, and $`\delta `$ is the Dirac delta function (which comes from our Fourier conventions and the definition of $`S_h(f)`$ as a spectral density with units of inverse frequency). There are two spectrum conventions in the literature: “two-sided” (defined for $`|f|<\mathrm{}`$) and “one-sided” ($`f0`$). They differ in amplitude by a factor of two: $`S_{\mathrm{one}}=2S_{\mathrm{two}}`$. We use the one-sided version in this paper. From this, we can define other useful quantities, such as the fractional contribution to the energy density of the Universe from a logarithmic interval of frequency,

$$\mathrm{\Omega }_{\mathrm{GW}}(f)=\frac{1}{\rho _c}\frac{d\rho _{\mathrm{GW}}}{d\mathrm{ln}f}=\frac{2\pi ^2}{3H_0^2}f^3S_h(f),$$
(4)
where $`H_0=100h_0`$km s-1 Mpc-1 is Hubble’s constant, and $`\rho _c=3H_0^2/8\pi G`$ is the critical density in an FRW (Freedman-Robertson-Walker) universe. Finally, we can also write down an expression relating the power spectrum to a “characteristic strain” spectrum, $`h_c(f)`$,

$$h_c(f)=\sqrt{fS_h(f)}.$$
(5)
We return in §4 to a discussion of the spatial spectrum of the stochastic background as we formulate the response of the Pulsar Timing Array to these waves.

In the following section, we will calculate the number density of sources in a given interval of strain and frequency, $`N(h,f)dhdf`$.444Here and throughout, we use the notation $`N()`$ to refer to the general concept of a number density function, and not to some specific function of the arguments. Thus, $`N(z,f)dzdf`$ refers to the number density of binaries in a redshift and frequency interval, whereas $`N(h,f)dhdf`$ refers to the number density of binaries in a strain and frequency interval. This is related to the strain power spectrum by

$$S_h(f)=_0^{\mathrm{}}h_{\mathrm{rms}}^2N(h,f)𝑑h.$$
(6)
Here, $`h_{\mathrm{rms}}`$ is defined so that the integrated power from a single event counted by $`N(h,f)`$ is $`𝑑fS_h(f)=h_{\mathrm{rms}}^2`$.

Thus, we need to combine various ingredients to calculate the observed gravitational wave spectrum due to binary MBHs: the galaxy merger rate, the black hole population demographics amongst galaxies, MBH binary dynamics and MBH binary gravitational wave emission, all of which enter $`N(h,f)`$. In order to calculate $`N(h,f)`$, we first consider a related quantity, $`\dot{N}(z)dz`$, the *rate* of coalescence events happening in redshift interval $`dz`$. Since this is a rate, the *number* of events of a given type in a time bin $`dt`$ is $`\dot{N}(z)dzdt=N(h,f)dhdf`$. We can relate this number to what we want by

$$N(h,f)=\dot{N}(z)\frac{dz}{dh}\frac{dt}{df}=\dot{N}(z)\frac{dz}{dh}\frac{dt}{dt_p}\frac{dt_p}{df_p}\frac{df_p}{df},$$
(7)
where the subscript $`p`$ denotes proper values in the rest frame at $`r(z)`$. Hence,

$$N(h,f)=\dot{N}[z(h)]\frac{dz}{dh}\left(f_p\frac{dt_p}{df_p}\right)\frac{(1+z)}{f}\dot{N}[z(h)]\frac{dz}{dh}\tau _{\mathrm{GW}}\frac{(1+z)}{f},$$
(8)
where we have introduced the gravitational wave timescale measured in the rest-frame:

$$\tau _{\mathrm{GW}}\left(f_p\frac{dt_p}{df_p}\right).$$
(9)
In this derivation we have assumed that the formation rate of the MBH binaries changes slowly on cosmological timescales, and that the binary enters the gravitational wave regime rapidly, so the rate is only a function of redshift and not, say, of the initial time of the halo merger. This should be an excellent approximation, as argued in Phinney (2001).

### 2.1 Galaxy Merger/Black Hole Coalescence Rate

We would therefore like to calculate $`\dot{N}(z,M_1,M_2)dzdM_1dM_2`$, the rate of black hole binary coalescence events observed at $`z=0`$ occurring in a given redshift interval $`dz`$ at $`z`$, due to MBHs in mass intervals $`dM_1`$ and $`dM_2`$. Because the formation mechanism for these objects in galactic halos is at best poorly understood, we will make a simplifying ansatz:

$$\dot{N}(z,M_1,M_2)=\nu (z)\frac{\varphi (M_1,z)\varphi (M_2,z)}{\varphi _{}^2},$$
(10)
where $`\nu (z)`$ is the total merger rate of the spheroids in which the MBHs reside in a notation similar to that of RR and $`\varphi (M,z)`$ is the mass function of black holes at redshift $`z`$ (normalized so $`𝑑M\varphi (M,z)=\varphi _{}(z)`$, the number density of MBHs at $`z`$). That is, we assume that the coalescence rate of black hole binaries is both independent of their masses and the properties of the halos in which they are present, and rapid on cosmological timescales. We will see that we are primarily interested in integrals over the full distribution, in which case the individual factors here can be seen as appropriately weighted averages over the multivariate distribution. We return to the determination of the black hole mass function in §2.2, and proceed next with formulation of the galaxy merger rate.

We start by considering a shell of some proper depth $`dl_p=cdt_p`$ at redshift $`z`$. This shell has a proper volume $`dV_p`$, or comoving volume $`dV_c`$ of

$$dV_p=4\pi d_A^2dl_p(1+z)^3dV_c,$$
(11)
where $`d_A`$ is the angular diameter distance. This distance is related to the FRW coordinate distance, $`r`$, by $`d_A=H_0^1[(a_0H_0r)/(1+z)]`$; note that, for example, Peebles (1993) refers to $`r(z)`$ itself as the “angular size distance”. The FRW distance is a function of $`z`$ and cosmological parameters as given by the dimensionless integral

$$a_0H_0r(z)=\frac{a_0H_0}{c}𝒮_k\left[\frac{c}{a_0H_0}_0^z\frac{dz^{}}{E(z^{})}\right],$$
(12)
where $`H_0=100h_0\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$ is the Hubble constant, $`a_0`$ is the FRW scale factor today, $`(\dot{a}/a)^2=H_0^2E^2(z)=H_0^2\left[\mathrm{\Omega }_m(1+z)^3+\mathrm{\Omega }_\mathrm{\Lambda }+(1\mathrm{\Omega }_m\mathrm{\Omega }_\mathrm{\Lambda })(1+z)^2\right]`$, $`\dot{a}=da/dt_p`$, $`𝒮_k(x)`$ is $`(\mathrm{sin}x,x,\mathrm{sinh}x)`$ for (closed, flat, open) geometries, and $``$ gives the curvature of space. In a flat Universe to which will specialize hereafter,

$$a_0H_0r(z)=_0^z\frac{dz^{}}{E(z^{})}(\mathrm{flat}).$$
(13)
We define the merger event rate per comoving unit volume per proper unit time at coordinate time $`t`$ as $`R(t)`$. The total event rate for the shell is $`R(t)dV_c`$. The observed rate from the shell requires correction of $`R`$ for the redshift of the time interval, $`dt_p/dt=1/(1+z)`$, and then substitution of earlier results:

$$g(t)dt=R(t)dV_cdt_p=R(t)dt_p\mathrm{\hspace{0.33em}4}\pi c(1+z)^2d_A^2=R(t)dt_p\mathrm{\hspace{0.33em}4}\pi c^3H_0^2(a_0H_0r)^2.$$
(14)
Converting the differential from time to redshift, we finally have

$$\nu (z)dz=g(t)dt=4\pi c^3R(z)H_0^3\frac{\left[a_0H_0r(z)\right]^2}{(1+z)E(z)}dz.$$
(15)
That is, $`g(t)𝑑t=\nu (z)𝑑z`$ gives the rate observed at $`z=0`$, with contributions from all along the past light cone. This is just the same as equation (12) of RR, except that we write our event rate as $`R`$ per unit proper time rather than their $`F_m`$ per unit redshift, although we calculate it as a function of $`z`$ (so $`F_mdz=Rdt_p`$). Also, they specialize to a flat universe with no cosmological constant. In particular, their equation (14) is equivalent to our expression with $`R(z)`$ equal to a constant.

Carlberg et al. (2000), following the method of Patton et al. (2000), estimate the galaxy merger rate out to $`z1`$ in the CNOC2 and CFGRS redshift surveys. They count the number of “kinematically close pairs”, those likely to be gravitationally infalling, based on their positions and redshifts. They find that roughly half of “close” pairs observed in projection at $`20h^1`$ kpc are physically close based on morphological disturbances of the pairs. The number of pairs detected is close to what is expected from the two-point correlation function. Because the total number of pairs in their sample at this separation is small, they determine the overall fraction at this separation by extrapolating the two-point correlation function of galaxies from measurements out to $`100h^1`$ kpc. Using dynamical arguments and results of simulations, they find that the average time for a merger starting at this separation of $`20h^1`$ kpc is $`T=0.3`$ Gyr, and that a fraction $`0.3`$ of such pairs will merge in this time. Approximately, then, $`/T1\mathrm{Gyr}^1`$ gives the fractional merger rate for these “close pairs,” which should hold roughly constant as long as the galaxy population is similar to that of the low-redshift universe.

Using these arguments along with observations of the number of close pairs from the CNOC2 and CFGRS surveys over $`0z1`$, Carlberg et al. (2000) find

$$R(z)\varphi _{}(z)R_g(z)\varphi _{}(z)(0.09\pm 0.05)(1+z)^{0.1\pm 0.5}\left(\frac{}{0.3}\right)\left(\frac{0.3\mathrm{Gyr}}{T}\right)\mathrm{Gyr}^1,$$
(16)
where $`R_g`$ is the merger rate for a single object and $`\varphi _{}(z)`$ is the comoving volume density of merging objects which for us is spheroids containing MBHs. We will from here on assume the merger-rate parameters to be equal to these canonical values given in equation (16). Crucially, we also assume that the measured merger rate per galaxy is equal to the desired merger rate. Inserting this merger rate into equation (15), we find

$$\nu (z)dz=0.03\mathrm{yr}^1\left[\frac{\varphi _{}(z)}{10^3h_0^3\mathrm{Mpc}^3}\right]\left[\frac{R_0}{0.09\mathrm{Gyr}^1}\right]\left[\frac{R_g(z)}{R_0}\right]\left[\frac{(a_0H_0r)^2}{(1+z)E(z)}\right]dz.$$
(17)
We have taken $`\varphi _{}=1.0\times 10^3h_0^3\mathrm{Mpc}^3`$ as our value for the number density of merging objects; see §2.2. We have left the evolution of the merger rate, $`R_g(z)/R_0`$ \[with normalization $`R_0=R_g(0)`$\], unspecified. Whereas Carlberg et al. (2000) find $`R_g(z)`$ approximately constant, Patton et al. (2002), a very similar group of authors using somewhat different data, find $`R_g(z)(1+z)^{1.3\text{}2.3}`$. To parameterize this difference (as well as other substantial observational disagreement on the merger rate as a function of time), we consider several other models, generically parameterizing the merger rate per unit time, $`R_g(z)`$, as a power-law in the expansion factor, $`R_g(z)=R_0(1+z)^\gamma `$. We will also allow the merger rate per unit redshift to be a power-law, which means that $`R_g(z)=R_0E(z)(1+z)(1+z)^{\gamma ^{}5/2}`$, where we use $`\gamma ^{}5/2`$ since that has $`\gamma ^{}=\gamma `$ when $`\mathrm{\Omega }_m=1`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }=0`$. In Figure 1 we show the merger rate with $`\gamma =0`$; increasing this exponent would increase the merger rate at higher redshift. Below, we will find that our limited knowledge of the merger rate is a main contributor to the uncertainty in the final GW spectrum. RR consider merger rates with both higher present-day normalizations and stronger redshift evolution (see Figure 1).

We note as an aside an alternative to these phenomenological and observational merger rates. Menou et al. (2001) have instead considered theoretical calculations of the merger rate of galaxy halos. They perform Monte Carlo simulations of the merger history of halos using semi-analytic methods (so-called “merger trees”) based on numerical simulations and the Press-Schechter formalism and its extensions. This allows them to vary the black hole mass function and its relationship to the underlying halo mass function as a function of time. In this work, however, they do not consider the detailed gravitational wave spectrum that results from their models.

Our model is useful beyond the nHz-$`\mu `$Hz regime. The LISA satellite will be sensitive to the total event rate of MBH binary inspirals as they go through their final coalescence. This event rate is just $`\nu (z)𝑑z`$, with $`\nu (z)`$ given by equation (15) or equation (17). We show this quantity in Figure 2. For the models we have been considering, the event rate is 0.01–1 binary coalescence events per year. However, only a subset of these events will be detected with sufficiently high signal-to-noise for a sufficient length of time over the different phases of the MBH binary merger to enable a measurement of their individual masses and distances (e.g., Hughes, 2002).

### 2.2 Black Hole Population Demographics

Magorrian et al. (1998) showed that essentially all “spheroids” (elliptical galaxies or the bulges of spirals) of mass $`M_{\mathrm{sph}}`$ have central black of mass $`M_{}0.006M_{\mathrm{sph}}`$. More recently Gebhardt et al. (2000), Ferrarese & Merritt (2000), Merritt & Ferrarese (2001a) and Tremaine et al. (2002) have shown that the data are consistent with an even tighter correspondence between the black hole mass and the velocity dispersion $`\sigma _v`$ of the galaxy raised to a power between about 4 and 5. These authors find a significantly smaller ratio of $`M_{}`$ to $`M_{\mathrm{sph}}`$ than that in Magorrian et al. (1998) with more careful treatment of the data. The $`\sigma _vM_{}`$ relation further implies that MBHs may grow primarily by other process such as accretion, rather than by the coalescence process discussed here (Merritt & Ferrarese, 2001b): if the spheroids grow simply by mergers, and the MBH mass initially tracks the spheroid mass, we might expect that the MBH mass would grow linearly with the mass of the spheroid, whereas the velocity dispersion raised to some power would not. There is some controversy, which is discussed by Merritt & Ferrarese (2001a), of the implications of these results for the detailed $`M_{}`$$`M_{\mathrm{sph}}`$ relationship.

Determination of either the mass function, the velocity dispersion function, or the luminosity function and the mass-to-light ratio for the population of galaxies with spheroids will provide us with the mass function of central black holes. Unfortunately, there have been no good measurements of the velocity dispersion distribution for galaxy spheroids. Instead, we use the original, less tight, correlation between black hole mass and spheroid mass. We will take the differential probability of black hole mass at a given spheroid mass, $`\varphi (M_{})dM_{}`$, to be lognormal with mean and dispersion given by

$$\mathrm{log}_{10}(M_{}/M_{\mathrm{sph}})=\mathrm{log}_{10}(0.0012)\pm 0.45$$
(18)
(Merritt & Ferrarese, 2001b). We then combine this with a theoretical model for the spheroid luminosity function (Fukugita et al., 1998; Magorrian & Tremaine, 1999, Maggorian, private communication) which is also lognormal with

$$\mathrm{log}_{10}(L_{\mathrm{sph}}/L_{})=\mathrm{log}_{10}(2\times 10^9)\pm 0.6.$$
(19)
Because we want an actual spatial density, we must normalize the luminosity function with the spheroid density, $`\varphi _{}=10^3h_0^3\mathrm{Mpc}^3`$. We also need the spheroid mass-to-light ratio, $`M/L=7.0(L/10^{10}L_{})^{0.2}M_{}/L_{}`$. Putting all this together by integrating over the two lognormal distributions gives us the MBH mass function, also in lognormal form, with

$$\mathrm{log}_{10}(M_{}/M_{})=\mathrm{log}_{10}(1.2\times 10^7)\pm 0.6$$
(20)
This expression is appropriate for elliptical and S0 galaxies, and for the spheroidal component of disks and spirals. There is also evidence (e.g., McLure & Dunlop, 2001, 2002) that the same or very similar relations hold for low-redshift ($`z<0.2`$) AGNs (Seyferts and QSOs), although there is some dispute regarding the constants of proportionality.

In Figure 3 we compare this mass function with the earlier model of Small & Blandford (1992), which was used in RR. Although the shapes of the distributions are similar, the more recent results on demographics predict that MBHs are more common: there are roughly an order of magnitude more at masses below $`10^9M_{\mathrm{}}`$.

This hard data on central black hole populations is of necessity at low redshift. In this work we are also concerned with the evolution of this population. The very existence of active galaxy populations at high redshift of course implies that some fraction of galaxies have had central black holes for a long time, but the detailed demographics are less certain. In their recent work, Menou et al. (2001) consider a scenario in which the fraction of spheroids containing an MBH is a function of redshift. In the following, such evolution is completely degenerate with changes in the merger rate, but for definiteness we will consider two scenarios, one in which the black hole mass function is constant with time, and another in which the normalizing number density follows a power law in $`(1+z)`$; this will allow us to combine a possible evolution in the merger rate with evolution in MBH populations as a single power law as in the previous subsection. In general we would certainly expect the peak of the mass function to also evolve over cosmological timescales, an issue we defer to later study.

### 2.3 Black Hole Binary Dynamics

Once the merger of two galaxy halos occurs, we must then follow the evolution of the central black holes. This is the subject of some debate. The paucity of obvious close binary nuclei imply that mergers between black-hole possessing-galaxies must be rare or short-lived. We know that the host galaxies themselves have often experienced at least one “major merger” in their lifetime. This is especially true of the most massive ellipticals harboring the most massive black holes. The recent study of the dynamics of MBHs in merging galaxies given by Yu (2001) builds on the seminal work of Begelman et al. (1980).

The evolution of the MBH binary proceeds in three stages. First, the pair spiral into what becomes the common single core of the merger remnant driven by dynamical friction, on a timescale given by

$$t_{\mathrm{DF}}\left(\frac{4\times 10^6\text{yr}}{\mathrm{log}N_{}}\right)\left(\frac{\sigma _c}{200\text{km}\text{s}^1}\right)\left(\frac{r_c}{100\text{pc}}\right)^2\left(\frac{M_{}}{10^8M_{}}\right)^1,$$
(21)
where the core of the galaxy has radius $`r_c`$, velocity dispersion $`\sigma _c`$ and contains $`N_{}`$ stars, and $`M_{}`$ is the mass of the smaller black hole.

The binary continues to harden (inspiraling and gaining energy) first through dynamical friction, subsequently by three-body interactions, and finally as a hard binary, when the semimajor axis is given by

$$a\frac{GM_{}}{4\sigma _c^2}=2.8\text{pc}\left(\frac{\sigma _c}{200\text{km}\text{s}^1}\right)^2\left(\frac{M_{}}{10^8M_{}}\right),$$
(22)
and the Keplerian period is given by

$$P4.4\times 10^4\text{yr}\left(\frac{M_{}}{10^8M_{}}\right)\left(\frac{\sigma _c}{200\text{km}\text{s}^1}\right)^3.$$
(23)
At this separation, however, the timescale for evolution via gravitational radiation, given by equation (2), is much longer than a Hubble time. So the typical pair must somehow evolve to a much tighter orbit with $`a0.02\text{pc}`$ in order to ever complete the merger. At this separation the Keplerian period is $`P30\text{yr}`$ and orbital speed is $``$4000 km s-1. Various dynamical mechanisms have been proposed for this process. Begelman et al. (1980) first discussed the dynamics in detail and noted the difficulty of reaching the rapid GW inspiral regime, as the MBH binary empties the “loss cone” in the stellar distribution function of stars with which the MBH binary can interact. RR surveyed the field at the time, and decided that while “unassisted stellar dynamics” would bring only a small fraction to the GW regime, “external perturbations may, however, cause efficient inspiral.” Examining these possible perturbations further, Gould & Rix (2000) and Armitage & Natarajan (2002) revisited the earlier idea (Begelman et al., 1980) that gas dynamics—the same physics as planetary migrations—could be responsible for this evolution, and indeed for fueling some or all quasar activity. Recently, Milosavljević & Merritt (2001) have done detailed N-body calculations of a MBH-MBH binary falling into a spherical potential of stars, and discuss the issue of potential stalling of the binary coalescence. Yu (2001) has examined in yet greater detail the stellar-dynamical schemes in tri-axial galaxies. Yu concludes that indeed there should be a sizable population of binary MBHs remaining in massive post-merger galaxies, in the absence of such gas-dynamical effects, but notes the difficulties of detecting the remnant binaries. On the other hand, Zhao et al. (2002) argue that the central stellar distribution of galaxies is considerably more complicated (cuspy and asymmetric) keeping the loss cone populated, possibly rescuing the stellar-dynamic scenarios.

Another possibility is that the binary will be driven to the inspiral regime by the presence of a third MBH from a second merger event. Tidal forces from the third body—or equivalent perturbations from a non-axisymmetric galactic potential—can decrease the merger time by an order of magnitude and significantly increase the eccentricity of the binary via the “Kozai mechanism” (Blaes et al., 2002).

In the following, we will assume that some mechanism like gas physics is indeed successful. Thus, all relevant pairs enter the gravitational wave regime and eventually shed their orbital energy in the form of gravitational radiation. The final ten thousand seconds or so of this evolution are thought to be the most energetic GW events in the (present-day) Universe, and will be easily detected by the LISA satellite over a wide range of MBH masses and redshifts (Hughes, 2002). Here, though, we are concerned with the quiescent quasi-Keplerian evolution that precedes the infall events.

### 2.4 Gravity Waves from Binary MBHs

Now, we will briefly review the emission of gravitational radiation from binary black holes. First, we will need the amplitude of the gravitational radiation emitted by a binary at distance $`D`$. This is given by (Peters & Matthews, 1963; Thorne, 1989)

$$h_{\mathrm{rms}}(r)=4\sqrt{\frac{2}{5}}\left(\frac{G}{c^3}\right)^{5/3}\left(\frac{2\pi }{P_p}\right)^{2/3}\left(\frac{c}{D}\right),$$
(24)
where $`=[M_1M_2(M_1+M_2)^{1/3}]^{3/5}`$ is the “chirp mass” of the system, and $`P_p`$ is the proper rest-frame period of the binary. The observed frequency of the radiation in the $`n`$th harmonic will be $`f=nf_p/(1+z)`$, where $`f_p1/P_p`$. Only terms with $`n2`$ are nonzero, and $`n=2`$ for the circular orbits to which we will restrict ourselves hereafter (Quinlan, 1996; Quinlan & Hernquist, 1997). However, if the binaries are driven together by the presence of a third MBH (Blaes et al., 2002), the eccentricity may be significant and populate higher-frequency harmonics. This would only serve to increase the amplitude of gravitational radiation and populate the higher harmonics (by an eccentricity-dependent factor). This amplification would offset, or exceed, the loss resulting from expulsion of the third body which is assumed to coalescence in our model. Equation 24 differs by a factor of $`\sqrt{3/4}`$ from that given in RR. Instead of the “characteristic strain” from Thorne (1989), which includes this factor to recover the signal-to-noise for an earth-bound detector, we just use the angle-averaged mean-square strain (i.e., eq.  is $`1/\sqrt{2}`$ times the maximum rms strain). We point out for ease of interpretation of this and later formulae that $`(GM/c^3)`$ has units of time. In a cosmological setting, we replace $`D`$ with $`ca_0H_0r(z)/H_0`$.

We will also need the characteristic timescale of the emission, defined in equation (9). Using the Kepler formula, $`f^2a^3=(2\pi )^2G(M_1+M_2)`$ and the formula for energy loss from gravitational waves, the gravitational wave timescale is (e.g., Peters & Matthews, 1963; Shapiro & Teukolsky, 1983)

$$\tau _{\mathrm{GW}}^1\left(\frac{1}{f_p}\frac{df_p}{dt_p}\right)=\frac{96}{5}\left(\frac{G}{c^3}\right)^{5/3}\left[\pi f(1+z)\right]^{8/3}.$$
(25)
where we have used $`f_p=f(1+z)/n`$ and set $`n=2`$ for circular orbits.

## 3 Strain Spectrum

Now we can combine all of the ingredients into an observable strain power spectrum. We have

* $`\nu (z)dz`$, the volume merger rate observed at $`z=0`$ of galaxies with MBHs located between $`z`$ and $`z+dz`$;
* $`\varphi (M_{})dM_{}`$ the distribution of MBH masses, normalized to the number density of galaxies whose merger rate is given by $`\nu (z)`$;
* $`h_{\mathrm{rms}}(z,f,M_1,M_2)`$, the strain observed at $`z=0`$ for a binary of given masses with frequency $`f`$ at redshift $`z`$; and
* $`\tau _{\mathrm{GW}}`$, the gravitational radiation timescale.

Putting all this together using equation (7), the strain distribution is

$$N(h,f,M_1,M_2)dhdfdM_1dM_2=N(z,f,M_1,M_2)\frac{dz}{dh}dhdfdM_1dM_2$$
(26)
where the distribution in redshift is

$`N(z,f,M_1,M_2)`$ $`=`$ $`\nu \left(z\right){\displaystyle \frac{\varphi \left(M_1\right)\varphi \left(M_2\right)}{\varphi _{}^2}}\left(1+z\right)\tau _{GW}(M_1,M_2,z,f_p){\displaystyle \frac{1}{f}}`$ (27)
$`=`$ $`{\displaystyle \frac{5}{96}}{\displaystyle \frac{\varphi \left(M_1\right)\varphi \left(M_2\right)}{\varphi _{}^2}}\left({\displaystyle \frac{G}{c^3}}\right)^{5/3}\nu \left(z\right)\left(1+z\right)^{5/3}\left(\pi f\right)^{8/3}{\displaystyle \frac{1}{f}}`$
$`=`$ $`{\displaystyle \frac{5}{96}}{\displaystyle \frac{4\pi c^3}{H_0^3}}{\displaystyle \frac{\varphi \left(M_1\right)\varphi \left(M_2\right)}{\varphi _{}^2}}\left({\displaystyle \frac{G}{c^3}}\right)^{5/3}R\left(z\right){\displaystyle \frac{\left[a_0H_0r\left(z\right)\right]^2}{E\left(z\right)\left(1+z\right)^{8/3}}}\left(\pi f\right)^{8/3}{\displaystyle \frac{1}{f}}.`$
First, note that this expression differs from equation (20) of RR; they neglected to include the factor of $`(1+z)`$ needed to convert from proper rest-frame time interval $`t_p`$ to observed interval $`t`$. Second, note that equation (27) factors into functions of mass, redshift and frequency. This means that marginal densities (integrals over one or more of the variables) are easy to calculate and have the same dependence on the remaining parameters as the unmarginalized distribution. This will make the following calculations particularly simple. However, it is certain that a fuller treatment of the model, incorporating the likely dependence of the halo merger and MBH coalescence rate simultaneously upon mass and redshift, including a complete accounting for the probability of a given MBH binary of successfully inspiraling to reach the GW regime, will no longer display this simplicity.

We now formulate the strain spectrum using the definitions of equations (56):

$$h_c(f)^2=f𝑑h𝑑M_1𝑑M_2h_{\mathrm{rms}}^2(z,f,M_1,M_2)N(h,f,M_1,M_2).$$
(28)
We can simplify this expression considerably by using $`N(h,\mathrm{})dh=N(z,\mathrm{})dz`$, which allows us to remove the $`dz/dh`$ factor in equation (26). Thus

$$h_c(f)^2=f𝑑z𝑑M_1𝑑M_2h_{\mathrm{rms}}^2(z,f,M_1,M_2)N(z,f,M_1,M_2).$$
(29)
Using the results of equations (2429) and the developments of the previous sections, the strain spectrum is given by

$$h_c(f)^2=𝑑z\frac{4\pi c^3}{3}\left[\pi f(1+z)\right]^{4/3}\left(G/c^3\right)^{5/3}_{}\frac{R(z)}{H_0E(z)},$$
(30)
where we use $`\mathrm{}_{}`$ to denote averages over the MBH population. In equation (30)

$$^{5/3}_{}=𝑑M_1𝑑M_2\frac{\varphi (M_1)\varphi (M_2)}{\varphi _{}^2}^{5/3}=\left(2.3\times 10^7M_{}\right)^{5/3},$$
(31)
where we have used the MBH mass function of §2.2, and have assumed that the number density of MBHs and the halos whose merger rate is described by equation (17) are the same.

We can now also separate out the $`z`$ integration to see the simple form of the result:

$`h_c\left(f\right)^2`$ $`=`$ $`{\displaystyle \frac{4\pi c^3}{3}}\left(G/c^3\right)^{5/3}_{}\left(\pi f\right)^{4/3}R_0{\displaystyle \frac{R\left(z\right)}{R_0}\frac{dz}{H_0E\left(z\right)\left(1+z\right)^{4/3}}}`$ (32)
$`=`$ $`\left(1.8h_0\times 10^{15}\right)^2\left({\displaystyle \frac{}{2.3\times 10^7M_{}}}\right)^{5/3}_{}\left({\displaystyle \frac{f}{\text{yr}^1}}\right)^{4/3}{\displaystyle \frac{R\left(z\right)}{R_0}\frac{dz}{E\left(z\right)\left(1+z\right)^{4/3}}}.`$
In the last line we have set the number density of merging galaxies and their present day merger rate to

$`\varphi _{}`$ $`=`$ $`1.0\times 10^3h_0^3\mathrm{Mpc}^3,\text{and}`$
$`R_g`$ $`=`$ $`0.09\mathrm{Gyr}^1.`$ (33)
Phinney (2001) derived an expression equivalent to equation (32), and noted the weak dependence of $`h_c(f)`$ on cosmology owing to a cancellation of $`a_0H_0r(z)`$ factors. Only the $`E(z)`$ factor remains in the integral. However, the details of the distribution $`N(h,f)`$ do depend more strongly on the cosmology.

In Figure 4 we show the dimensionless integral,

$$I_{h^2}\frac{R(z)}{R_0}\frac{dz}{E(z)(1+z)^{4/3}},$$
(34)
which is essentially the same as $`(1+z)^{1/3}`$ in Phinney (2001). We plot $`I_{h^2}`$ for a variety of cosmologies, merger rates and black-hole mass-function histories. We parameterize the latter two by the single exponent $`\gamma `$, defined such that $`\varphi (M_1,z)\varphi (M_2,z)R(z)(1+z)^\gamma `$ or $`E(z)(1+z)^{\gamma ^{}3/2}`$ as in §§2.12.2. The maximum redshifts of 3 and 5 used in Figure 4 are motivated by the distribution of objects with redshift which is discussed below in §3.2. For $`\gamma 2`$, we are most sensitive to $`z<3`$ and the cutoff is irrelevant. For stronger merger-rate evolution, we are sensitive to more distant events. Figure 4 shows the weak dependence upon both the cosmology ($`\mathrm{\Omega }_m`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }`$) and the maximum redshift of MBH binary formation. The results also show the weak dependence on whether we consider the merger rate as a power law in redshift (with exponent $`\gamma `$) or time ($`\gamma ^{}`$). Similarly, if we allow this exponent to vary by as much as one or two units, we again get very little change. The dependence on the maximum redshift is obviously somewhat stronger if we increase the power-law exponent $`\gamma `$, thereby putting more events at higher redshift. There remains a prefactor of $`R_0`$ in $`h_c^2(f)`$—the overall normalization of the merger rate—which can have a strong effect. Of course, each of these effects can change things by a factor of two. In summary, the uncertainty in $`h_c^2(f)`$ is at least an order of magnitude.

### 3.1 The Strain Distribution

We would like to estimate the entire probability distribution of $`h_c^2(f)`$ about this mean. The simplest way to approach this is to use the method of *characteristic functions*, which are simply the Fourier transforms of probability densities. In order to do this, first consider the quantity

$$h_c^2(f)\frac{df}{f}=𝑑z𝑑M_1𝑑M_2h_{\mathrm{rms}}^2(z,f,M_1,M_2)\widehat{N}(z,f,M_1,M_2)𝑑f.$$
(35)
In this equation, the quantity $`\widehat{N}(z,f,M_1,M_2)dzdfdM_1dM_2`$ gives the actual number in a cell of size $`dz\times df\times dM_1\times dM_2`$ at a given frequency and bandwidth $`f`$ and $`df`$. It can be shown (see Appendix A) that the characteristic function of the distribution of $`h_c^2(f)df/f`$ is given by

$$\mathrm{ln}\phi (t)=𝑑z𝑑M_1𝑑M_2N(z,f,M_1,M_2)\left[e^{ith_{\mathrm{rms}}^2(z,f,M_1,M_2)}1\right]𝑑f.$$
(36)
While one can show that this integral exists for the cases we are considering, it cannot be expressed in closed form. For this reason, the distribution itself—the inverse Fourier transform of $`\varphi (t)`$—is quite difficult to calculate.

To start, we can use the fact that the coefficients of the Taylor expansion of $`\mathrm{ln}\varphi (t)`$ around $`t=0`$ give $`\kappa _k`$, the ‘semi-invariant moments’ or ‘cumulants’ of the distribution (i.e., the mean, variance, skewness, kurtosis, etc. for $`k=1,2,3,4,\mathrm{}`$), which we can read off as

$$\kappa _k(f)=𝑑z𝑑M_1𝑑M_2\left[h_{\mathrm{rms}}^2(z,f,M_1,M_2)\right]^kN(z,f,M_1,M_2)𝑑f.$$
(37)
Note that these are the moments of $`h_c^2(f)df/f`$, not $`h_c^2(f)`$ itself, which can only be measured with a finite bandwidth. For example, the variance of $`h_c^2(f)df`$ is

$$\mathrm{var}\left[h_c^2(f)df\right]=f^2\mathrm{var}\left[h_c^2(f)df/f\right]=f^2𝑑z𝑑M_1𝑑M_2h_{\mathrm{rms}}^4N(z,f,M_1,M_2)𝑑f.$$
(38)
We measure $`h_c^2(f)`$ itself with some finite bandwidth, $`\mathrm{\Delta }f`$, as $`h_c^2(f)_{\mathrm{\Delta }f}h_c^2(f)𝑑f/\mathrm{\Delta }f`$. The variance on this quantity decreases as the bandwidth increases.

The variance of the MBH density $`N(h,f)`$ under consideration is

$$\mathrm{var}\left[h_c^2(f)\right]=4\pi c^3H_0\frac{32}{5}\left(\frac{\mathrm{\Delta }f}{f}\right)^1\left(G/c^3\right)^5_{}𝑑z\left[a_0H_0r(z)\right]^2R(z)/E(z).$$
(39)
As before, this separates into an average over the MBH mass function and an integral over redshift. Considering the latter, we see that the integral in question diverges as $`z0`$, with the integrand behaving as $`dz/z^2`$. Indeed, all of the integrals with $`k>1`$ diverge: the moments do not exist. Nonetheless, the characteristic function and the distribution itself remain well-defined. This divergence occurs because the strain from a single event is proportional to the inverse of the distance; although the mean-square accumulated strain (the power spectrum) converges, the dispersion from small distances diverges.

We have performed numerical experiments with the simpler problem of a uniform MBH binary distribution in a non-expanding, Euclidean universe, so $`N(z)z^2`$ and $`h_{\mathrm{rms}}^2(z)1/z^2`$ for all $`z<z_{\mathrm{max}}`$—an approximation to the situation at $`z1`$. We find that the distribution of strains is highly skewed, with a power-law tail toward larger amplitudes. This tail is due to the possible contribution of low-redshift pairs which give a contribution $`h^2(z)1/z^2`$.

### 3.2 Other Quantities

Our derivation also allows us to calculate other quantities related to the distribution of strains. First, we can ask what is the number of binaries in a frequency interval,

$$N(f)df=𝑑z𝑑M_1𝑑M_2N(z,f,M_1,M_2)𝑑f.$$
(40)
Using equation (27), and again under the assumptions we have made so far that the $`\varphi _{}`$ is independent of redshift, the $`z`$ integration separates out, and we see first that $`N(f)df/ff^{8/3}df/f`$, a steeply-falling function of frequency. In full,

$$N(f)df=5.6\times 10^5df/f\left(\frac{}{4.1\times 10^6M_{}}\right)^{5/3}_{}\left(\frac{f}{\text{yr}^1}\right)^{8/3}𝑑z\frac{R(z)}{R_0}\frac{[a_0H_0r(z)]^2}{E(z)(1+z)^{8/3}}.$$
(41)
$`^{5/3}_{}`$ is $`\left(4.1\times 10^6M_{}\right)^{5/3}`$ for our MBH mass function which is different from the value above owing to a different weighting function. In Figure 5 we show the dimensionless integral occurring in this expression for a variety of cosmologies and merger histories. The stochastic GW background at nHz (yr-1) frequencies is therefore the result of nearly a million simultaneous binaries across the Universe.

Similarly, we can ask, what is the range in redshift that we are probing? This is just

$$N(z)\frac{R_g(z)}{R_0}\frac{[a_0H_0r(z)]^2}{E(z)(1+z)^{8/3}}.$$
(42)
which we show in Figure 6. For larger values of the exponent ($`\gamma 2`$), relatively more events occur at higher redshift, and the dependence of the observed strain on the maximum redshift will be stronger. More specifically, we can calculate the mean redshift,

$$z(f)=\frac{𝑑zzN(z,f)}{𝑑zN(z,f)}.$$
(43)
Again, the integrals separate and we find that this average redshift is independent of frequency:

$$z=\frac{𝑑zR(z)\left[(a_0H_0r)^2/E(z)\right](1+z)^{5/3}}{𝑑zR(z)\left[(a_0H_0r)^2/E(z)\right](1+z)^{8/3}}1.$$
(44)
This is independent of the details of the MBH population and the normalization of the merger rate, but dependent on the merger rate evolution and the cosmology. In Figure 7 we show the average redshift for a variety of cosmologies and merger histories.

In the previous sections we calculated the strain spectrum, $`h_c^2(f)`$, and variants. Using the above developments, we can calculate a slightly different quantity: the average strain amplitude that contributes to $`h_c^2(f)`$. This is just

$`h^2(f)`$ $`=`$ $`{\displaystyle \frac{h_c^2(f)}{N(f)}}f^{4/3},`$ (45)
which rises owing to the more rapid decrease of number density with frequency relative to the spectrum amplitude decrease. That is, higher frequency signals are from small numbers of rare, powerful events in comparison with lower frequency signals.

Finally, we can also calculate quantities related to the mass distribution that contribute to the strain spectrum. Again due to the factorization of $`N(M_1,M_2,z,f)`$ these are independent of redshift and frequency. For any function of the masses, we see that

$`F(M_1,M_2)={\displaystyle \frac{F(M_1,M_2)^{5/3}_{}}{^{5/3}_{}}},`$ (46)
where as above $`\mathrm{}_{}`$ refers to averages over the MBH distribution. One particular quantity of interest is simply the average chirp mass:

$$=2.9\times 10^6M_{}.$$
(47)
which, as would be expected, is close to the average of the MBH distribution itself (Fig. 3). Another is the mass ratio,

$$qM_1/M_2=11.$$
(48)
Because the MBH mass distribution is itself so wide (Fig. 3), the gravitational waves sample a fairly wide range of black-hole mass ratios. In particular, we are sensitive to minor mergers as well as major mergers. However, these quantities change to $`=5\times 10^6M_{}`$ and $`q=1.7`$ if a minimum black hole mass of $`10^6M_{}`$ is assumed. Such a minimum hole mass is warranted based on considerations of the strong dependence of dynamical friction time scale on mass (Xu & Ostriker, 1994). The values in equations (4748) can be expressed as individual masses $`M_1,M_2(1.3\times 10^7,1.1\times 10^6)M_{}`$, or $`(7.4\times 10^6,4.4\times 10^6)M_{}`$ if the distribution is cut off at $`10^6M_{}`$.

In a similar vein, we can weight these functions by the signal amplitude they produce, rather than the number of binaries contributing. Since $`h_{\mathrm{rms}}^2^{10/3}`$ (eq. ), this changes the exponents in equation 46 from $`5/3`$ to $`+5/3`$. Such a weighting gives a signal-weighted average chirp mass of $`9\times 10^7M_{}`$ and mass ratio of $`6`$. Again, we see that a small number of high-mass black holes are responsible for the bulk of the signal. Similarly, the signal-weighted average redshift is somewhat lower than the number-weighted redshift of Figure 7.

### 3.3 Results

In Figure 8 we consider two models of galaxy mergers and the black hole population, and show the resultant spectrum of characteristic strain along with several realizations of a Monte-Carlo simulation of the spectrum. We perform the Monte Carlo by simply Poisson-sampling from $`N(z,f,M_1,M_2)`$ on a grid of the parameters and summing the contributions of $`h_{\mathrm{rms}}^2(z,f,M_1,M_2)`$. We have chosen $`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, the constant halo merger rate given in equation (17), the MBH mass function from equation (20), and a maximum redshift of $`z_{\mathrm{max}}=3`$. The upper set of curves have a relatively strong power-law index $`\gamma =2`$ describing the galaxy merger rate; the lower curves have $`\gamma =0`$. We have assumed a bandwidth of $`\delta f/f0.12`$. We note that at higher frequencies the simulation is somewhat below the expected mean spectrum. This is due to the high-amplitude tails mentioned in the discussion of the full strain distribution, §3.1 above: small number statistics and numerical inaccuracies mean that we miss the large spikes.

We expect the actual dispersion to be greater at higher frequencies—smaller numbers of objects are contributing more due to the steeply-falling function $`N(f)`$ (eq. ). However, just as the mean of our simulation falls below the expectation value, so does the dispersion around the mean: we miss the large excursions. This is a larger effect for the relatively high value $`\gamma =2`$ than at $`\gamma =0`$, representing comparatively fewer merger events at low redshift.

However, a central result of this paper is that the *slope* of the spectrum is independent of all of this. As in Phinney (2001), we find that $`h_c^2(f)f^{4/3}`$, or $`S_h(f)f^{7/3}`$, or $`\mathrm{\Omega }_{GW}(f)f^{2/3}`$. The overall amplitude, however, will vary up and down by the factors calculated above and shown in Figure 4. As mentioned above, this means that, on the one hand, we cannot use the information from any observed slope to limit the physics. On the other hand, it potentially provides a check that we are observing the background from a coalescing population (although it is not sensitive to the details of that population). Indeed, other proposed sources of stochastic gravitational waves often have different slopes. Early universe cosmological backgrounds, for example, generically have $`\mathrm{\Omega }(f)\text{const}`$ and hence $`h_c^2(f)f^2`$.

## 4 Gravity Wave Detection with a Pulsar Timing Array

We turn now to the experimental approaches to detection of the stochastic background from MBH binaries. We begin with a statement of the response of pulsar timing to a plane GW. In the following subsection we discuss the angular and temporal basis functions required for detection of a stochastic background. In the final subsection the current limits and future prospects are discussed in relation to the expected signal level presented in Figure 8.

### 4.1 Plane Wave Response

Precision timing of the arrival of pulses from a pulsar allows the direct detection of gravitational radiation traversing the sightline (Sazhin, 1978; Detweiler, 1979). The spacetime metric strain perturbation $`h`$ can be viewed as a perturbation from unity of the index of refraction for electromagnetic radiation in flat space. For a plane GW traversing the entire sightline to a pulsar the pulse arrival time is advanced, or retarded, in proportion to the incomplete cycles of the wave traversed by the pulse at the pulsar and in the solar system. Observable effects are possible only when the GW amplitude, or state of the spacetime metric, changes locally or at the pulsar or both in ways distinguishable from pulsar spin and other parameters. A useful analytic approach is to formulate the apparent redshift $`Z`$ of a pulsar spin frequency which arise from gravitational radiation. Observationally $`Z`$ is estimated by taking the derivative of the residuals that result after fitting pulse arrival times to a model for the pulsar’s spin properties and astrometric coordinates, and binary parameters if relevant. For a single plane wave the observable may be written as

$$Z(\alpha ,\beta ,\gamma ;t_\mathrm{r})=\frac{(\alpha ^2\beta ^2)}{2(1+\gamma )}[h_+(ct_\mathrm{r})h_+(ct_\mathrm{e}\gamma l)]+\frac{\alpha \beta }{(1+\gamma )}[h_\times (ct_\mathrm{r})h_\times (ct_\mathrm{e}\gamma l)],$$
(49)
where $`(\alpha ,\beta ,\gamma )`$ are the direction cosines of the pulsar with respect to the GW, $`(h_+,h_\times )`$ are the polarization components of the GW, and $`(t_\mathrm{e},t_\mathrm{r})`$ are the times of emission (e) and reception (r) of the pulses from the pulsar which is at a distance $`l`$.

The $`t_\mathrm{r}`$ terms in equation (49) will be correlated between different pulsars. Therefore timing measurements of an array of pulsars across the sky, the Pulsar Timing Array, creates a gravitational wave telescope that is sensitive to a spectrum of waves with periods less than the duration of the measurements which is years, or $``$ 10 nHz (e.g. Backer, 1993). The five degrees of freedom in equation (49) correspond to the five degrees of freedom in the trace-free space-time metric. The $`t_\mathrm{e}`$ terms in equation (49) will be uncorrelated between different pulsars, and yet will have comparable amplitude. These terms create an irreducible noise background in addition to any limitations from the measurement errors and properties of the pulsars themselves or the intervening turbulent plasma through which the signals propagate. Lommen & Backer (2001) apply these ideas to a search for plane waves from both our Galactic Center and nearby galaxies which host MBHs based on the hypothesis that the central objects are binary.

### 4.2 Fitting Algorithms for a Stochastic Background

In the preceding sections we considered the production of a stochastic background of gravitational radiation from coalescences of MBH binaries. This is summarized in the form of a characteristic strain spectrum, $`h_c(f)`$, Figure 8. For the purposes of analysis of timing data from a spatial array of pulsars we need to consider the full tensor perturbation field $`h_{ij}(𝐱,t)`$. All elements are statistically independent owing to the random superposition of many sources from many directions. This superposition destroys the simple angular plane wave pattern described in equation (49). However, Burke (1975) has shown that 5/8 of the variance in the stochastic fluctuations can be extracted using the 5 quadrupole spherical harmonic $`Y_m^l`$ functions as the angular basis.

Temporal modulations of each spherical harmonic term can be described using either Fourier frequency coefficients as the basis, or polynomials (Foster & Backer, 1990), or orthogonal polynomials (Stinebring et al., 1990). The polynomial approaches are particularly suited to this analysis owing to both the steep spectrum expected for the MBH-MBH GW emission (eq. ) and the need to fit for the spin parameters of the stars. Given any choice of the temporal basis function, the amplitude $`h_{\mathrm{rms}}(f)`$ for any temporal term is formed by the square root of the sum of the squares of the 5 spherical harmonic coefficients. Conversion to $`h_c(f)`$ requires consideration of the spectral window function which will be addressed in a later work.

A second approach to the angular basis has been suggested by Hellings (1990). The tensor field stated above can be transformed into a spectrum of waves $`\stackrel{~}{h}_{ij}(𝐪,f)`$. These can be decomposed into waves traveling along the three cardinal directions, each with independent polarizations:

$$\left[\begin{array}{ccc}h_{11}& h_{12}& h_{13}\\ h_{12}& h_{22}& h_{23}\\ h_{13}& h_{23}& h_{33}\end{array}\right]=\left[\begin{array}{ccc}h_{x,+}& h_{x,\times }& 0\\ h_{x,\times }& h_{x,+}& 0\\ 0& 0& 0\end{array}\right]+\left[\begin{array}{ccc}h_{y,+}& 0& h_{y,\times }\\ 0& 0& 0\\ h_{y,\times }& 0& h_{y,+}\end{array}\right]+\left[\begin{array}{ccc}0& 0& 0\\ 0& h_{z,+}& h_{z,\times }\\ 0& h_{z,\times }& h_{z,+}\end{array}\right]$$
(50)
In this case one of the six terms in the decomposed format is not independent owing to the trace-free property. The Doppler patterns of the three waves can then be summed to create Doppler patterns for each tensor element. Thus the full variance of the stochastic background can be sampled. As in the spherical harmonic case each term will be temporally modulated, and the quantity $`h_{\mathrm{rms}}`$ is formed by the square root of the sum of the squares of 5 independent terms.

### 4.3 Current Limits and Future Prospects

While the preceding subsection outlines approaches to detection of the GW background using an array of pulsars, measurements from a single object can establish an upper limit on the background. Kaspi et al. (1994) has stated the best limit at nHz frequencies using precision timing of PSR B1855+09 relative to the world’s best atomic time scale. A second pulsar in their study, PSR B1937+21, was shown to be unstable on long time scales, an effect they attribute to internal structure of that neutron star. Thorsett & Dewey (1996) improved the limit using further statistical considerations. At a typical frequency of 5 nHz (8 yr period) they conclude that the energy density in gravitational radiation per logarithmic frequency interval is less than $`\mathrm{\Omega }_gh^26\times 10^8`$, where $`h`$ in this expression is the ratio of Hubble’s constant to 100 km s-1.

Energy density scales quadratically with the time derivative of the metric strain. Therefore the level of $`\mathrm{\Omega }_gh^2`$ that can be detected by pulsar timing scales quadratically with the rms timing residual $`R`$, and as $`T^5`$ with the duration of the experiment $`T`$, which is the inverse of the minimum frequency sampled. Recently Lommen & Backer (2002) have extended the Kaspi study of PSR B1855+09 by more than doubling the experiment duration to 17 y. The new upper limit on energy density is roughly an order of magnitude lower at frequencies near 5 nHz. A full report on this work is in preparation. An approximate statement of the Kaspi-Lommen limit on characteristic strain is given in Figure 8, which shows that the measurements are approaching a level of significance given our current model of the Universe. An alternate way of stating this is that the measurements are placing useful constraints on the uncertain parameters in our model. The ratio of the pulsar detection level to the characteristic strain spectrum level scales as $`T^{13/6}`$.

The Pulsar Timing Array experiment, which will use both existing and new data sets, can improve on the Kaspi et al. (1994) result in three obvious ways: (a) smaller timing residuals owing to combined data sets, new and upgraded telescopes and better data acquisition techniques; (b) more objects which both provide the capability of actual detection as opposed to upper limits and, if sufficiently numerous, provide a “root N” advantage; and (c) longer experiment duration. In Figure 8 we show a reasonable goal for the future: 200 ns timing precision over eight years. Current measurement series are already achieving this level of precision for a few objects.

## 5 Summary

We have calculated the spectrum of the stochastic background of gravitational radiation from the Universe of coalescing binary black holes in the centers of galaxies with simple parameterizations of the current uncertainties. This is followed by a discussion of the influence of the stochastic background on precision timing measurements of pulsars which includes upper limits on the background and use of an array of pulsars for direct detection. This work was motivated by improvements both in our knowledge of the present-day massive black hole population and the rate of galaxy mergers as well as new pulsar measurements.

Following the approach of past authors, we construct the spectrum from coalescing binary MBHs from: the galaxy merger rate, the black hole population demographics amongst galaxies, MBH binary dynamics and MBH binary gravitational-wave emission. Our galaxy merger rate is based on observations of close pairs and an estimate of their dynamical friction time scale. We convert from the galaxy mass function to the black hole mass function using recent determinations of the correlation between black hole mass and spheroid mass.

Our “fiducial” model has rapid evolution in the merger rate per unit time per galaxy as a function of redshift, and a constant mass function of MBHs out to $`z=3`$. The characteristic strain spectrum predicted for this model, $`h_c(f)10^{16}(f/\mathrm{yr}^1)^{2/3}`$, is just below the latest observational limits at $`f0.2\mathrm{yr}^1`$. Whereas the slope and general character of this prediction agree with the earlier work of Rajagopal & Romani (1995), our predicted amplitude is somewhat higher. This increase is the result of the order of magnitude increase in the local MBH number density, which is tempered by a somewhat lower merger rate at high redshift.

With our formulation we can calculate other quantities. The number of binaries contributing for unit bandwidth ($`\delta f/f1`$) at nHz frequencies is approximately $`10^6`$. It is evident from our simulations that the variance in the strain spectrum in the same frequency band is roughly 50%, considerably larger than the $`10^3`$ that would be expected if all events were weighted equally. This indicates that the spread around the mean is due to the possibility of a small number of high-amplitude, nearby events. The median redshift for the binaries contributing to the spectrum is $`1`$ if the merger rate evolves relatively slowly with redshift ($`\gamma 2`$), but can be large for strong evolution of the merger rate. The mean “chirp” mass is $`3\times 10^6M_{}`$ and the mean mass ratio is $`11`$; the individual masses are then $`10^7`$ and $`10^6M_{}`$.

Our calculations allow us to see directly how the remaining lack of understanding of some of these physical processes impacts our predictions. These are areas for future research. The most important factors are:

* The value of the present day galaxy merger rate is a matter of considerable debate. The value we have quoted, $`0.09\mathrm{Gyr}^1`$ per galaxy from Carlberg et al. (2000) is known to no better than 50%. Moreover, this value is measured for massive $`L_{}`$ galaxies; does it apply to the wider range harboring MBHs?
* The evolution of the merger rate for moderate redshifts ($`z1`$) is even less well-known. Looking back still further, the evolution of the merger rate at $`z>1`$ is yet more difficult to pin down. The high redshift surveys, in progress and planned, are crucial for progress.
* We have not taken into account the details of MBH binaries dynamics: under what circumstances, if any, do MBH binaries make it to the GW regime? We have assumed that all galaxy mergers promptly lead to coalescence of their central black holes. A more complex model needs to consider both a delay between merger and coalescence which may depend on mass and the possibility of expulsion in a triple system.
* Throughout, we have simplified the MBH coalescence rate by splitting the MBH population from the galaxy merger rate, as in equation (10). Although this ansatz is useful for examining the dependences of the results upon the various ingredients, a full treatment of the time- and mass-dependent merger rate of galaxies (or the black holes directly) needs to be included. Appropriate results could come from considerably improved observations, N-body simulations, or formalisms such as the Press-Schechter calculation of the halo mass function.
* Finally, we have assumed a simple model for the black hole population, combining the spheroid mass function and MBH demographics, as in §2.2. This does not take into account more recent results relating the MBH mass with the gravitational potential of its parent spheroid via the velocity dispersion (Gebhardt et al., 2000; Ferrarese & Merritt, 2000). Moreover, we assume that the evolution in the MBH mass function can be modeled simply as a scaling of the number density, rather than the perhaps more realistic scenario where the mean and spread of the mass also evolves.

In future work, we plan on addressing these issues. Many of the observational uncertainties in the model are at high redshift. Modern ‘semi-analytic’ models (e.g., Cattaneo et al., 1999; Haehnelt & Kauffmann, 2000) go a long way toward elucidating these issues within the context of particular scenarios for structure formation. Observationally, the quasar luminosity function (Boyle et al., 2000; Fan et al., 2001), OH Megamasers (Darling & Giovanelli, 2002) and “X-type” radio sources (Merritt & Ekers, 2002) are other tools for understanding the high-redshift MBH population. For this work, a detailed high-resolution search for more nearby MBH binaries themselves, perhaps in the form of double-nucleus QSOs (but see Kochanek et al., 1999) and Ultra-Luminous Infrared Galaxies (Murphy et al., 2001), will be a crucial tool. Analytical and N-body studies of the dynamics and energetics of MBH binaries within galaxies will provide population information on MBH-MBH mass ratios and orbit sizes crucial to future work. These will allow us to also constrain the dynamics of the MBH binary population: if and when do they stall prior to the BW regime?

We can, of course, use our current model to at least see the effect of these uncertainties. We can find a lower limit to the effect of MBH mass-function evolution by using a suitable initial mass function for the whole redshift range. If the average MBH mass has increased by an order of magnitude (say), then this would lower the characteristic strain (eq. ) by no more than $`10^{5/6}7`$.

The current experimental limit on the low-frequency gravitational wave background has improved over the pioneering work of Kaspi et al. (1994), and will be accurately stated by Lommen & Backer (in preparation). The new limit is lower at lower frequencies while the expected spectrum is rising. This new result will constrain parameters of our model spectrum on the high side of the current estimate. Ongoing work with an array of millisecond pulsars has the prospect of significant improvement of detection capability in the coming decade. Techniques for optimal fitting of Pulsar Timing Array data need to be further developed to meet the demands of the new measurements.

We thank Jon Arons, Ray Carlberg, Ron Hellings, Yuri Levin, Andrea Lommen, John Maggorian, Sterl Phinney, Philip Stark and members of the Center for Particle Astrophysics for helpful conversations. AHJ acknowledges support from NSF KDI grant 9872979 and NASA LTSA grant NAG5-6552, and PPARC in the UK. DCB acknowledges support from NSF grant AST-9731106 which partially supported the Pulsar Timing Array experiment.

## Appendix A Appendix

In this appendix, we calculate the probability distribution of the quantity $`h_c^2(f)`$ as given by equation (28). For generality, consider some quantity

$$y=𝑑zg(z)N(z)$$
(A1)
where $`N(z)dz`$ gives the probability of some event happening in $`(z,z+dz)`$ (i.e., a Poisson rate), and $`g(z)`$ is any function. We want to know $`P(y)`$, the distribution function of $`y`$. We start with a few important facts from probability theory.

First, the characteristic function (or moment generating function) of a distribution is defined as the Fourier transform of the distribution, i.e.,

$$\phi (t)=e^{ity}=𝑑ye^{ity}P(y).$$
(A2)
This has the property that the moments of the distribution are then given by

$$y^k=i^k\frac{d^k\phi (t)}{dt^k}|_{t=0}$$
(A3)
and the “cumulants” or “connected moments” or “semi-invariants” ($`\kappa _k`$, the mean, variance, skewness, kurtosis, …, for $`k=1,2,3,4,\mathrm{}`$) are given by the same expression, substituting $`\mathrm{ln}\phi `$ for $`\phi `$.

Next, we use the fact (trivial from Fourier theory) that the characteristic function of the sum of two variables is the product of the individual functions, and that the characteristic function of $`ay`$, where $`a`$ is a scalar and $`y`$ is the random variable is $`\phi (at)`$, where $`\phi (t)`$ is the characteristic function for $`y`$.

Finally, we need the characteristic function for a Poisson distribution, $`P(n)=r^ne^r/n!`$ where $`r`$ is the rate. This is just $`\phi (t)=\mathrm{exp}\left[r(e^{it}1)\right]`$.

Putting all of these facts together gives the characteristic function for $`y`$, which is just the sum (integral) of Poisson variables with rates $`N(z)dz`$ times scalars $`g(z)`$. Thus

$$\mathrm{ln}\phi (t)=𝑑zN(z)\left[e^{itg(z)}1\right].$$
(A4)
The Taylor series of $`\mathrm{ln}\phi `$ around $`t=0`$ then gives the cumulants of the distribution. They are just

$$\kappa _k=𝑑z\left[g(z)\right]^kN(z)$$
(A5)
(Note that the mean, $`k=1`$, is what would be trivially expected.)

In our case, we have $`N(z)dzN(z,f,M_1,M_2)dzdM_1dM_2`$ and $`g(z)h_{\mathrm{rms}}^2(z,f,M_1,M_2)`$ for $`y=h_c^2(f)df/f`$, giving equation (36) above.


text:
# 1 Introduction

## 1 Introduction

Supergravity solutions of intersecting D-branes are relatively easy to find as long as they are sufficiently smeared . Supergravity solutions of the localized intersections are far more difficult to find. Starting with the work of  there has been a steady enterprise of attempts to construct such supergravity solutions . However, to date, there are no known techniques for determining the gravitational back reaction due to a general intersecting brane configurations that arise in string theory.

For brane intersections involving a D6-brane, there are special techniques which allow the explicit construction of certain localized intersections . Here, one takes advantage of the fact that the M-theory lift of the geometry of the D6-brane near its core is an ALE space which is essentially the flat $`^4`$. This allowed the construction of the completely localized supergravity solution of D2 parallel to D6, as well as D4 ending on the D6. This method however was limited in its applicability to the region near the core of the D6.

The aim of this paper is to construct the fully localized supergravity solution of D2 parallel to D6 without restricting to the near core region of the D6. The construction of the supergravity solution turns out to be possible for this case due to the fact that there exists a simple ansatz which reduces the problem to a single linear differential equation  which is separable  and admits a regular boundary condition . Therefore, the problem can be solved using elementary methods. The solutions we obtain are nonetheless very interesting. Instead of taking the near horizon limit of the D6-brane as was done in , one can consider the near horizon limit of the D2-brane. This gives rise to a geometry where a D6-brane is slicing through the near horizon geometry of the D2-branes. From the open string point of view, this corresponds to taking the decoupling limit which keeps only the gauge fields on the D2-brane and the charged fundamental matter arising from strings stretching between the D2 and the D6 branes. It is therefore the 2+1 dimensional version of the holographic dual of gauge theory with flavors considered recently by Katz and Karch in  for the 3+1 dimensional case. There are several advantages for considering the case of 2+1 dimensions over 3+1. The field theory in 2+1 dimensions is superrenormalizable even after adding the fundamental matter, in contrast to 3+1 dimensional theory which looses the asymptotic freedom.111There are, however, constructions involving orientifolds which maintain conformal invariance . The supergravity solution takes the full gravitational back reaction of the D6-branes into account in contrast to the 3+1 example where the effect of D7-brane is treated in the probe approximation. One can therefore think of the localized D2/D6 solution as the “cleaner” version of the AdS/CFT-like correspondence with fundamental flavors.

The near horizon limit of the D2/D6 solution captures the full RG flow of the weakly coupled supersymmetric Yang-Mills theory with fundamental matter in the UV to a superconformal fixed point in the IR along the lines of . In fact, certain qualitative features of precisely this RG flow was anticipated in . Our explicit supergravity solution confirms the expectation of . As an added bonus, these supergravity solutions can be used to illustrate the mirror symmetry of Intriligator and Seiberg  in the language of AdS/CFT correspondence.

## 2 The Solution

Let us begin by describing the explicit construction of the supergravity solution. The idea is to start with a lift of D6-brane to the Taub-NUT geometry in M-theory and to consider the effect of placing large number of M2-branes in this background. To describe this background, one employs the ansatz

$`ds^2`$ $`=`$ $`H^{2/3}(dt^2+dx_1^2+dx_2^2)+H^{1/3}(dy^2+y^2d\mathrm{\Omega }_3^2+ds_{TN}^2)`$
$`F`$ $`=`$ $`dtdx_1dx_2dH^1`$ (2.1)
where $`ds_{TN}^2`$ is the metric of the Taub-NUT space. We find it convenient to parameterize the coordinates of the Taub-NUT space so that the metric takes the form

$$ds_{TN}^2=\left(1+\frac{2m}{r}\right)(dr^2+r^2(d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2))+\left(\frac{1}{1+\frac{2m}{r}}\right)(4m)^2\left(d\psi +\frac{1}{2}\mathrm{cos}\theta d\varphi \right)^2$$
(2.2)
where the coordinate take on values with range $`0r`$, $`0\theta \pi `$, $`0\varphi 2\pi `$, and $`0\psi 2\pi `$. The parameter $`m`$ is related to the radius $`R`$ of the circle of the Taub-NUT metric at infinity by the formula

$$R=4m.$$
(2.3)
For small radius

$$r=\frac{z^2}{8m}m,$$
(2.4)
the Taub-NUT metric simplifies to

$$ds_{TN}^2=dz^2+z^2d\mathrm{\Omega }_3^2.$$
(2.5)
The ansatz (2.1) is a solution to the equation of motion of 11 dimensional supergravity if $`H`$ solves the harmonic equation in the background of $`^4\times \text{Taub-NUT}`$ space

$$^2H=0,$$
(2.6)
except at the location of the M2-brane source. In order to maximize the symmetry of the problem to simplify the analysis, let us consider the case where the M2-brane source is placed at the origin $`y=r=0`$. One can than take $`H(y,r)`$ to be a function of two variables, and the harmonic equation (2.6) becomes

$$\left(\frac{1}{1+\frac{2m}{r}}\right)\left(\frac{^2}{r^2}+\frac{2}{r}\frac{}{r}\right)H(y,r)+_y^2H(y,r)=0.$$
(2.7)
Our task is simply to solve this differential equation. To this end, it is convenient to separate variables

$$H(y,r)=1+Q_{M2}\frac{d^4p}{(2\pi )^4}e^{ipy}H_p(r),$$
(2.8)
where $`Q_{M2}`$ is the membrane charge in the standard normalization

$$Q_{M2}=32\pi ^2N_2l_p^6.$$
(2.9)
Then, $`H_p`$ satisfies

$$\left(\frac{1}{1+\frac{2m}{r}}\right)\left(\frac{^2}{r^2}+\frac{2}{r}\frac{}{r}\right)H_p(r)p^2H_p(r)=0.$$
(2.10)
Since this is a second order differential equation, there are two solutions. The one which decays at large $`r`$ can be written in a closed analytic form

$$H_p(r)=c_pe^{pr}𝒰(1+pm,2,2pr),$$
(2.11)
where $`𝒰(a,b,z)`$ is the confluent hypergeometric function . The normalization factor $`c_p`$ is fixed to

$$c_p=\frac{\pi ^2}{8}\frac{1}{m^2}(pm)^2\mathrm{\Gamma }(pm)$$
(2.12)
by requiring that in the $`m\mathrm{}`$ limit keeping $`z^2=8mr`$ fixed (which is equivalent to looking at $`rm`$), $`H_p(r)`$ becomes (using 13.3.3 of )

$$H_p(z)=\frac{\pi ^2}{2z^2}pzK_1(pz)$$
(2.13)
whose Fourier transform is

$$\frac{d^4p}{(2\pi )^4}e^{ipy}H_p(z)=\frac{1}{(y^2+z^2)^3}.$$
(2.14)
We finally arrive at the statement

$$H(y,r)=1+Q_{M2}𝑑p\frac{(py)^2J_1(py)}{4\pi ^2y^3}H_p(r)$$
(2.15)
where we have reduced (2.8) to an integral over a single variable by exploiting the spherical symmetry in $`p`$-space.

The solution (2.15) combined with the ansatz (2.1) is the main result of this paper. Dimensional reduction of this solution along the $`\psi `$ coordinate of the Taub-NUT geometry (2.2) will give rise to the solution type IIA supergravity describing D2 localized along the world volume of D6. The metric part of the solution is given by

$`ds^2`$ $`=`$ $`H(y,r)^{1/2}\left(1+{\displaystyle \frac{2m}{r}}\right)^{1/2}(dt^2+dx_1^2+dx_2^2)`$
$`+H(y,r)^{1/2}\left(1+{\displaystyle \frac{2m}{r}}\right)^{1/2}(dy^2+y^2d\mathrm{\Omega }_3^2)+H(y,r)^{1/2}\left(1+{\displaystyle \frac{2m}{r}}\right)^{1/2}(dr^2+r^2d\mathrm{\Omega }_2^2),`$
where $`4m=R=g_sl_s`$. It is clear that when we set $`Q_{M2}=0`$, the solution reduces to the supergravity solution containing only the D6-branes. Similarly, in the $`m0`$ limit,

$$H(y,r)=1+\frac{Q_{D2}}{(y^2+r^2)^{5/2}},$$
(2.17)
where using $`l_p=g_s^{1/3}l_s`$,

$$Q_{D2}=\frac{3}{64m}Q_{M2}=6\pi ^2g_sN_2l_s^5,$$
(2.18)
which agrees with the supergravity solution of the D2 by itself including the numerical factors. Although (2.15) is left in an integral form, the expression is completely explicit and the final integration can be done numerically if desired. To demonstrate this point, we have computed

$$f(r)=\frac{512m^3r^3}{Q_{M2}}(H(0,r)1)$$
(2.19)
numerically. The normalization was chosen so that $`f(r)=1`$ for $`r0`$. The result of this computation is illustrated in figure 1. The result clearly illustrates the cross-over between asymptotics (2.14) and (2.17) for small and large $`r`$, respectively.

## 3 The Decoupling Limit

Now that we have constructed the supergravity solution of the localized intersection of D2 and D6 in type IIA supergravity, let us consider taking its near horizon limit which gives rise to a holographic dual of the gauge theory on the D2-branes. We will scale $`l_s`$ to zero keeping the two dimensional gauge coupling

$$g_{YM2}^2=g_sl_s^1$$
(3.1)
fixed. In this limit, the gauge coupling on the six brane

$$g_{YM6}^2=(2\pi )^4g_sl_s^3=(2\pi l_s)^4g_{YM2}^2$$
(3.2)
goes to zero so that the dynamics on the six-brane decouples. In order to identify the corresponding near horizon geometry on the supergravity side, we also scale the radial coordinates so that

$$Y=\frac{y}{l_s^2},U=\frac{r}{l_s^2}$$
(3.3)
is fixed. In this limit, the harmonic function due to the D6-brane source scales as

$$1+\frac{2m}{r}=1+\frac{g_sl_s}{2r}=1+\frac{g_{YM2}^2}{2U},$$
(3.4)
whereas the harmonic function due to the D2-brane source (2.15) scales as

$$H(Y,U)=\frac{1}{l_s^4}h(Y,U),$$
(3.5)
where

$`h(Y,U)=`$
$`\pi ^2g_{YM2}^4N_2{\displaystyle 𝑑P\frac{(PY)^2J_1(PY)}{Y^3}P^2\mathrm{\Gamma }\left(\frac{g_{YM2}^2P}{4}\right)e^{PU}𝒰(1+\frac{g_{YM2}^2P}{4},2,2PU)}.`$
We have also scaled the integration variable

$$p=\frac{P}{l_s^2}$$
(3.7)
so that the string length $`l_s`$ does not appear anywhere in the definition of $`h(Y,U)`$.

We now have all the ingredients to explicitly write down the supergravity solution for the decoupled D2/D6 system. Generalizing slightly to the case with multiple coincident D6-branes, the solution takes the form
$`{\displaystyle \frac{ds^2}{l_s^2}}=h(Y,U)^{1/2}\left(1+{\displaystyle \frac{g_{YM2}^2N_6}{2U}}\right)^{1/2}(dt^2+dx_1^2+dx_2^2)`$      $`+h(Y,U)^{1/2}\left(1+{\displaystyle \frac{g_{YM2}^2N_6}{2U}}\right)^{1/2}(dY^2+Y^2d\mathrm{\Omega }_3^2)`$      $`+h(Y,U)^{1/2}\left(1+{\displaystyle \frac{g_{YM2}^2N_6}{2U}}\right)^{1/2}(dU^2+U^2d\mathrm{\Omega }_2^2)`$    $`h(Y,U)=`$      $`\pi ^2g_{YM2}^4N_2{\displaystyle 𝑑P\frac{(PY)^2J_1(PY)}{Y^3}P^2\mathrm{\Gamma }\left(\frac{g_{YM2}^2N_6P}{4}\right)e^{PU}𝒰(1+\frac{g_{YM2}^2N_6P}{4},2,2PU)}.`$

The only dependence of this metric on the string length $`l_s`$ is in the overall normalization which is what one expects for the supergravity dual of a quantum field theory. Note also that although we no longer have the “1” in the harmonic function of the D2-brane in taking the decoupling limit, we are left with the “1” in the harmonic function of the D6-brane. The effect of the D2-brane is therefore to “warp” not only the ALE region but also the asymptotically flat region of the D6-brane geometry.

## 4 RG Flow and Mirror Symmetry

Let us now interpret the various features of the supergravity solution (3) from the point of view of the field theory dual. By construction, (3) is dual to maximally supersymmetric SYM in 2+1 dimensions with gauge group $`SU(N_2)`$ further coupled to $`k=N_6`$ flavors of massless hypermultiplets in the fundamental representation. The coupling to the hypermultiplets reduces the number of unbroken supersymmetries from 16 to 8. The supergravity solution (3) asymptotes to the geometry of the near horizon D2-brane in the large $`U`$ limit. This suggests that the dynamics of this theory is dominated by the free gluons in the UV as one expects for a superrenormalizable theory. In the small $`U`$ region, the geometry of (3) asymptotes to $`AdS_4\times S_7/_k`$. This is the superconformal field theory one expects to find on M2-brane probing $`^4\times (^4/_k)`$.

Now, the superconformal theory on the M2-brane on an orbifold does not have a simple Lagrangian formulation. One way to define such a theory without relying on string theory is to define it as an IR fixed point of a different theory which has a Lagrangian formulation. The 2+1 SYM with $`k`$ flavors is one concrete example of a UV theory which flows to this superconformal field theory in the IR. Roughly speaking, by compactifying one of the directions of $`^4/_k`$ into a Taub-NUT space, one has embedded the dynamics of M2 into the dynamics of D2, and it is the latter which has a good Lagrangian description. By making the D2 warp the $`^4\times \text{Taub-NUT}`$ geometry due to its gravitational backreaction, one constructs the supergravity dual to the decoupled theory on the D2-brane. The supergravity solution (3) encodes this full renormalization group flow in the language of AdS/CFT correspondence. Similar observations can be found in the earlier work of .

It turns out that there is a different way to embed the superconformal field theory on the M2 on $`^4\times (^4/_k)`$ as the IR fixed point of a field theory with a Lagrangian description. This is the $`_k`$ quiver theory of the 2+1 SYM. From the point of view of branes, this amounts to considering the decoupling limit of M2-branes on $`^3\times S_1\times (^4/_k)`$. The supergravity solution for such a brane configuration is easy to find. They are simply

$$ds^2=H^{2/3}(dt^2+dx_1^2+dx_2^2)+H^{1/3}(dy^2+y^2d\mathrm{\Omega }_2^2+dz^2+dr^2+r^2ds_{Lens}^2)$$
(4.1)
where $`ds_{Lens}^2`$ is the metric on the Lens space which is the base of the $`^4/_k`$ viewed as a cone, and

$$H(y,z,r)=1+\underset{n=\mathrm{}}{\overset{\mathrm{}}{}}\frac{Q_{M2}}{(r^2+y^2+(z2\pi nR)^2)^3}.$$
(4.2)
It is possible to take the decoupling limit of this solution keeping

$$\frac{R}{l_s^2}=\frac{g_s}{l_s}=g_{YM2}^2=\text{fixed}$$
(4.3)
which will give rise to a different supergravity background describing the renormalization group flow of the 2+1 dimensional $`_k`$ quiver theory flowing to the same superconformal field theory.

What we have here is a pair of supergravity solutions, both of which asymptotes to the same $`AdS_4\times S_7/_k`$. It is therefore a holographic realization of two different RG flows which flow to the same conformal field theory in the far IR. This is mirror symmetry. Although the metric on the supergravity side asymptotes to the same thing near the core, the geometry away from the core of the two solutions are clearly different from each other. This illustrates quite explicitly in the AdS/CFT language the basic fact that mirror symmetry is an equivalence only for the far IR of a pair of field theories.

The basic idea behind the embedding of the superconformal field theory into a Lagrangian field theory was to compactify one of the dimensions either in the $`^4`$ or the $`^4/^k`$. The freedom to choose between the two was the basis for mirror symmetry. Let us now consider what happens if one compactifies both so that we have $`(^3\times S^1)\times \text{Taub-NUT}`$. Now there are two ways to reduce the same geometry from M-theory to type IIA. Let us for the sake of the argument reduce on the circle in the Taub-NUT. This will give rise to a D2/D6 system with one of the direction transverse to the D2-brane compactified. In the decoupling limit, compactness of the directions transverse to the world volume of the brane is an indication that the underlying gauge theory is a $`U(\mathrm{})/`$ theory because of the presence of images. According to the argument of , it is better to view this as a theory with one extra dimension. From the point of view of the supergravity, the same picture manifests itself in the fact that in the near horizon limit, where the backreaction of the D2-brane dominates, the proper size of the $`S^1`$ transverse to the D2 shrinks as one approaches the boundary. At the point where this proper size becomes smaller than the string length, the supergravity description of this geometry become unreliable, and following the argument of , one is instructed to go to the T-dual picture, where the D2 becomes a D3.

Unfortunately, the same T-duality maps the D6 to a D5. T-duality in supergravity is not capable of handling this map except for the case where the D5 is completely smeared. We are therefore unable to provide a purely supergravity description of the 3+1 dimensional UV fixed point for the decoupled theory on M2 in $`(^3\times S^1)\times \text{Taub-NUT}`$. The problem of finding this supergravity solution was attempted most recently in . Let us note in passing that at least the solution for the case of the smeared D5 considered in section 4.1 of  can be obtained from (2) by applying the T-duality rules for supergravity.

Another interesting aspect of the decoupled theory of M2 on $`(^3\times S^1)\times \text{Taub-NUT}`$ is the fact that its UV description on the field theory side is precisely the defect field theory introduced in . In fact, one of the main motivations of  was to find the purely gravitational holographic description of the defect field theory. By reducing from M-theory to IIA on the circle of the Taub-NUT and T-dualizing to IIB on the circle in $`^3\times S^1`$, we arrive at a defect theory consisting only of the D5 defects. The Lagrangian of this theory was worked out in . Alternatively, one could have reduced from M-theory to IIA on the circle in $`^3\times S^1`$ and T-dualizing along the circle of the Taub-NUT arriving at the theory with NS defects. The relation between the two ways of going from M-theory to type IIB is the natural extension of mirror symmetry. Clearly, from the point of view of the type IIB theory, this equivalence is S-duality. What one learns here is that the origin of mirror symmetry in 2+1 dimensional gauge theory is the S-duality of the defect field theory in 3+1 dimensions.

This idea of embedding a 2+1 dimensional theory in some UV structure to make mirror symmetry manifest is not a new idea. Embedding to string theory was exploited for this goal some time ago in . Embedding of the 2+1 dimensional theory into the 3+1 dimensional defect field theory amounts to taking the decoupling limit of . Although the formulation of defect field theories was strongly motivated by string theory, their existence is independent of string theory. The relation between S-duality of the defect field theory and the mirror symmetry of its dimensionally reduced theory can be studied by exploiting their holographic duality, instead of their embedding, to string theory. The former is physically economical.

By embedding a pair of 2+1 dimensional mirror theories into a pair of S-dual defect field theories, one obtains a mirror duality which applies at all scales, not just in the far IR. So the embedding into the defect field theories can be interpreted as an intricate UV modification of the mirror pair theories so as to extend their range of validity beyond the far IR. For the Abelian case, this issue was addressed in . For the non-Abelian case, the natural UV modification appears to involve a theory with one extra dimension and defects.

## 5 Some Generalizations

In the previous sections, we discussed mainly the localized D2/D6 supergravity configuration where all of the D2 and D6 are coincident. Let us now consider some generalizations.

### 5.1 Separating D2 from D6

One simple generalization one can consider is separating the D2-brane from the D6-brane. For the sake of concreteness, let us consider the case where there is one of each of D2 and D6.

From the point of view of M-theory, this is a configuration of a single M2-brane at a generic point on the Taub-NUT background (2.2). Therefore, the supergravity solution is given by the same ansatz (2.1) where $`H`$ is a solution of the harmonic equation (2.6) but with a source located at a generic point in the Taub-NUT. From the point of view of finding a localized D2/D6 solution of the type IIA supergravity equation of motion, however, one is only interested in sources which are smeared along the 11-th coordinate $`\psi `$. Let us therefore take the M2-brane source to be smeared evenly along $`\psi `$ as this would also simpilfy the analysis. The harmonic function $`H`$ is then a solution of

$$(_{TN}^2+_y^2)H=2\pi ^4Q_{M2}\delta ^4(\stackrel{}{y}\stackrel{}{y}_0)\delta ^3(\stackrel{}{r}\stackrel{}{r}_0)\delta (\psi \psi _0).$$
(5.1)
The factor of $`2\pi ^4`$ arises from the fact that

$$_r^2\left(\frac{1}{r^6}\right)=2\pi ^4\delta ^8(r)$$
(5.2)
in 8 dimensions. Smearing along $`\psi `$ and separating the variables as was done in (2.8), one finds that $`H_p(\stackrel{}{r})`$ is a solution of

$$\left(_r^2\frac{2mp^2}{|\stackrel{}{r}|}p^2\right)H_p(\stackrel{}{r})=\frac{2\pi ^4Q_{M2}}{2\pi R}\delta ^3(\stackrel{}{r}\stackrel{}{r}_0),R=4m$$
(5.3)
with appropriate boundary conditions. This equation is of the form

$$\left(_r^2+\frac{2k\nu }{|\stackrel{}{r}|}+k^2\right)G(\stackrel{}{r},\stackrel{}{r}_0)=\delta ^3(\stackrel{}{r}\stackrel{}{r}_0)$$
(5.4)
if one identifies

$$H_p(\stackrel{}{r})=\frac{2\pi ^4}{2\pi R}G(\stackrel{}{r},\stackrel{}{r}_0),\nu =imp,k=ip.$$
(5.5)
Precisely this equation with the appropriate boundary conditions was considered in  and the solution was found to be

$$G(\stackrel{}{r},\stackrel{}{r}_0)=\frac{\mathrm{\Gamma }(1i\nu )}{4\pi |\stackrel{}{r}\stackrel{}{r}_0|}\frac{1}{ik}\left(\frac{}{y}+\frac{}{x}\right)𝒲_{i\nu ,1/2}(ikx)_{i\nu ,1/2}(iky)$$
(5.6)
where $`_{a,b}(z)`$ and $`𝒲_{a,b}(z)`$ are the Whittaker functions, and

$$x=r+r_0+|\stackrel{}{r}\stackrel{}{r}_0|,y=r+r_0|\stackrel{}{r}\stackrel{}{r}_0|.$$
(5.7)
If one sends the source $`\stackrel{}{r}_0`$ to zero, $`G(\stackrel{}{r},\stackrel{}{r}_0)`$ simplifies to

$$G(\stackrel{}{r})=\frac{1}{4\pi r}\mathrm{\Gamma }(1i\nu )𝒲_{i\nu ,1/2}(2ikr).$$
(5.8)
Using the identity

$$𝒲_{\kappa ,\mu }=z^{\mu +1/2}e^{z/2}𝒰(\frac{1}{2}+\mu \kappa ,1+2\mu ,z)$$
(5.9)
and the identification (5.5), one can show that (5.8) is equivalent to (2.11) including all the numerical factors.

Using (5.6) and the identification (5.5), one can write an explicit expression for the supergravity solution of the localized D2/D6 configuration. Since harmonic equations are linear, it is straightforward to generalize this to the case where the D2 is distributed arbitrarily in transverse coordinates. It is also straightforward to generalize this solution to the case where there are multiple D6’s as long as all of the D6’s are coincident. Simply set

$$\nu =iN_6mp.$$
(5.10)
By scaling $`\stackrel{}{r}_0=\alpha ^{}\stackrel{}{U}_0`$ keeping $`\stackrel{}{U}_0`$ fixed, one can take the decoupling limit of the D2/D6 as we did in section 3. From the point of view of the field theory dual, this corresponds to turning on a vacuum expectation value of some of the adjoint scalars, so that the matter fields in the fundamental acquire mass.

### 5.2 Separating the D6

Another possible generalization one might consider is to separate the D6’s from one another. The M-theory lift of this configuration is the multi-centered Taub-NUT geometry whose metric is given by

$$ds^2=Vd\stackrel{}{r}^2+\frac{(4m)^2}{V}(d\psi +\stackrel{}{\omega }d\stackrel{}{r})^2$$
(5.11)
where

$$V(\stackrel{}{r})=1+\underset{i=1}{\overset{N_6}{}}\frac{2m}{|\stackrel{}{r}\stackrel{}{r}_i|},V=4m\times \stackrel{}{\omega }.$$
(5.12)
Applying the same ansatz as (2.1) gives rise to a new harmonic equation (2.6). Although it is not absolutely necessary to do so, let us smear the M2 source along the $`\psi `$ coordinates. This will simplify the analysis. Then, the harmonic equation (2.6) can be written explicitly as

$$\left(_r^2p^2V(\stackrel{}{r})\right)H_p(\stackrel{}{r})=\frac{2\pi ^4Q_{M2}}{2\pi R}\delta ^3(\stackrel{}{r}\stackrel{}{r}_0).$$
(5.13)
For a geneneral multi-centered Taub-NUT background, this is still a difficult equation to solve. The case of all the D6 being coincident gives rise to a single centered Taub-NUT considered in the previous section.

It turns out that double centered Taub-NUT also admits natural coordinates in which the harmonic equation (5.13) separates . Let us consider this case in some detail.

Consider a double centered Taub-NUT with $`N_1`$ coincident centers at $`\stackrel{}{r}_1`$ and $`N_2`$ coincident centers at $`\stackrel{}{r}_2`$. Without loss of generality, we can set $`\stackrel{}{r}_2=\stackrel{}{r}_1`$. One can then introduce the so called prolate spheroidal coordinates

$$\xi =\frac{|\stackrel{}{r}\stackrel{}{r}_1|+|\stackrel{}{r}\stackrel{}{r}_2|}{2L},\eta =\frac{|\stackrel{}{r}\stackrel{}{r}_1||\stackrel{}{r}\stackrel{}{r}_2|}{2L}$$
(5.14)
where $`2L=|\stackrel{}{r}_1\stackrel{}{r}_2|`$ is the distance between the two centers. The ranges of these coordinates are $`1<\xi `$ and $`1<\eta <1`$. The countours of fixed $`\xi `$ and $`\eta `$ are illusterated in figure 2.

Let $`\varphi `$ denote the angular coordinate around a symmetry axis defined by $`\stackrel{}{r}_1\stackrel{}{r}_2`$. Then, the set of coordinates $`\xi `$, $`\eta `$, and $`\varphi `$ specifies a point $`\stackrel{}{r}`$. In these coordinates, the harmonic equation (5.13) becomes

$`\left(_\xi (\xi ^21)_\xi +_\eta (1\eta ^2)_\eta +\left({\displaystyle \frac{1}{\xi ^21}}+{\displaystyle \frac{1}{1\eta ^2}}\right)_\varphi ^2L^2(\xi ^2\eta ^2)Vp^2\right)H_p(\xi ,\eta ,\varphi )`$
$`={\displaystyle \frac{2\pi ^4Q_{M2}}{2\pi LR}}\delta (\xi \xi _0)\delta (\eta \eta _0)\delta (\varphi \varphi _0)`$ (5.15)
where

$$(\xi ^2\eta ^2)V=(\xi ^2\eta ^2)+\frac{2m(N_1+N_2)\xi }{L}\frac{2m(N_1N_2)\eta }{L}.$$
(5.16)
To further analyze this problem, it is convenient to consider the solution to the equation

$$\left(_\eta (1\eta ^2)_\eta +\frac{k^2}{1\eta ^2}p^2L^2\eta ^22p^2mL(N_1N_2)\eta \right)B_{\lambda ,k}(\eta )=\lambda B_{\lambda ,k}(\eta ).$$
(5.17)
This equation takes the form

$$\frac{1}{\mathrm{sin}\theta }\frac{d}{d\theta }\mathrm{sin}\theta \frac{d}{d\theta }B_{\lambda ,k}+\left(\lambda \frac{k^2}{\mathrm{sin}^2\theta }+2p^2mL(N_1N_2)\mathrm{cos}\theta +p^2L^2\mathrm{cos}^2\theta \right)B_{\lambda ,k}=0$$
(5.18)
after making the change of variables

$$\eta =\mathrm{cos}\theta ,0<\theta <\pi .$$
(5.19)
In the $`L0`$ limit, this equation becomes the Legendre equation. For finite $`L`$, one expects a discrete spectrum of eigenvalues $`\lambda _n`$ and its associated eigenfunction $`B_{\lambda _n,k}(\eta )`$. They can be determined either using numerical methods, perturbation theory, or by expanding

$$B_{\lambda ,k}(\eta )=e^{pL\eta }\underset{s=0}{\overset{\mathrm{}}{}}c_sP_{s+k}^k(\eta ),$$
(5.20)
and deriving a recursion relation for the coefficients $`c_s`$ along the lines of . We will consider these eigenfunctions to be orthonormalized so that

$$_0^\pi 𝑑\theta \mathrm{sin}(\theta )B_{\lambda _m,k}(\mathrm{cos}\theta )B_{\lambda _n,k}(\mathrm{cos}\theta )=_1^1𝑑\eta B_{\lambda _m,k}(\eta )B_{\lambda _n,k}(\eta )=\delta _{mn}.$$
(5.21)
One now sees that upon parameterizing

$$H_p(\xi ,\eta ,\psi )=\underset{n,k}{}A_{\lambda _n,k}(\xi )B_{\lambda _n,k}(\eta )e^{ik\varphi },$$
(5.22)
where $`B_{\lambda _n,k}`$ are orthonormal, $`A_{\lambda _n,k}(\xi )`$ satisfies

$`\left(_\xi (\xi ^21)_\xi {\displaystyle \frac{k^2}{\xi ^21}}p^2(L^2\xi ^2+2m(N_1+N_2)L\xi )\lambda _n\right)A_{\lambda _n,k}(\xi )`$
$`={\displaystyle \frac{2\pi ^4Q_{M2}}{2\pi LR}}B_{\lambda _n,k}(\eta _0)e^{k\varphi _0}\delta (\xi \xi _0).`$ (5.23)
Although somewhat complicated, this is a linear inhomogenious ordinary differential equation which can be solved numerically or using the method of . One can then evaluate

$$H(y,\xi ,\eta ,\varphi )=1+\frac{d^4p}{(2\pi )^4}e^{ipy}\underset{n,k}{}A_{\lambda _n,k}(\xi )B_{\lambda _n,k}(\eta )e^{ik\varphi }$$
(5.24)
which when substituted into (2.1) gives rise to the supergravity solution of the D2 parallel to two collections of D6-branes.

## 6 Concluding Remarks

The main goal of this paper was the construction of the localized D2/D6 supergravity solution. We have constructed this solution explicitly as an integral expression in (2.15). Using this solution, it was possible to construct a holographic dual to 2+1 dimensional Yang-Mill with matter in the fundamental representation, and describe mirror symmetry in the language of AdS/CFT correspondence.

It would be interesting to explore the standard holographic observables: entropy, Wilson loop, and correlation functions, for this supergravity solution. It would also be interesting to explore the mapping of observables between the mirror pairs from the point of view of holography.

The D2/D6 system appears to be unique in providing a conceptually clean setup to add flavors to AdS/CFT. The D3/D7 system suffers from the lack of asymptotic freedom, and the D1/D5 system suffers from the lack of moduli-space. One can still describe the decoupled theory on D1/D5 in Born-Oppenheimer approximation, but that does not appear to be compatible with the holographic duality, which has as a starting point a stationary supergravity solution with definite configuration of static sources. One manifestation of this difficulty is the fact that a localized supergravity solution of D1 coincident with D5 does not even appear to exist .

The key to the simplicity is the fact that the D6-brane lifts to a Taub-NUT geometry in M-theory which is purely geometrical. Furthermore, the geometry is sufficiently regular both near the core and at infinity. This is what made generalization of the D2/D6 intersection in the near core region considered in  to the full Taub-NUT geometry possible. It would be very interesting to see if the D4/D6 intersection considered in  can be extended in a similar manner by taking advantage of the simplicity of the Taub-NUT geometry.

A different simple generalzation is the supergravity dual of the non-commutative gauge theory with fundamental matter. Such a solution can be found by applying the same T-duality transformation considered in  or by following the twist operation for the dipole theories introduced in . The resulting supergravity background is similar to (3) but with the metric along the D2-brane worldvolume replaced by

$$\frac{ds^2}{l_s^2}=h(Y,U)^{1/2}\left(1+\frac{g_{YM2}^2N_6}{2U}\right)^{1/2}\left(dt^2+\frac{dx_1^2+dx_2^2}{1+\mathrm{\Delta }^4h(Y,U)^1\left(1+\frac{g_{YM2}^2N_6}{2U}\right)^1}\right)+\mathrm{}$$
(6.1)
where $`2\pi \mathrm{\Delta }^2=\theta ^{12}`$ is the non-commutativity parameter along the D2-brane world volume.

There are other generalizations one might consider. For example, one can separate the D6’s completely, so that the M-theory background becomes that of a multi-center Taub-NUT. One might also consider finding the supergravity solution for a decoupled theory on M2 probing a manifold of $`Sp(2)`$ holonomy222The theory on the brane probe was considered in . Mirror symmetry for these theories were first considered from the point of view of brane constructions in . considered in . For these cases, separation of variables do not appear to work and more sophisticated methods for solving for the the harmonic function must to be employed.

## Acknowledgements

We would like to thank K. Hori, N. Itzhaki, A. Kapustin, A. Karch, D. Marolf, and S. Mukhi for discussions. We also thank the Aspen Center for Physics and Institute des Hautes Études Scientifiques where part of this work was done. This work is supported in part by DOE grant DE-FG02-90ER40542 and by the Marvin L. Goldberger fellowship.