Gravitational Wave Hotspots: Ranking Potential Locations of Single-Source Gravitational Wave Emission
The steadily improving sensitivity of pulsar timing arrays (PTAs) suggests that gravitational waves (GWs) from supermassive black hole binary (SMBHB) systems in the nearby universe will be detectable sometime during the next decade. Currently, PTAs assume an equal probability of detection from every sky position, but as evidence grows for a non-isotropic distribution of sources, is there a most likely sky position for a detectable single source of GWs? In this paper, a collection of galactic catalogs is used to calculate various metrics related to the detectability of a single GW source resolvable above a GW background, assuming that every galaxy has the same probability of containing a SMBHB. Our analyses of these data reveal small probabilities that one of these sources is currently in the PTA band, but as sensitivity is improved regions of consistent probability density are found in predictable locations, specifically around local galaxy clusters.
Pulsar Timing Arrays (PTAs) are collections of millisecond pulsars (Foster & Backer 1990) whose joint timings will show correlations that are a specific signature of gravitational waves (GWs) passing between the Earth and the pulsars (Sazhin 1978; Detweiler 1979). Currently, there are three PTAs working collaboratively to detect GWs in this way (Hobbs et al. 2010). These arrays will soon have the sensitivity to detect single extragalactic sources of GWs, (Yardley et al. 2010; Ellis, Siemens, & Creighton 2012) which are resolvable above a GW background (Sesana 2013). To increase the efficiency of PTA observations and potentially decrease the time until a detection is made, various groups have considered the optimization of PTA observations (Lee, Jenet, & Price 2008; Burt, Lommen, & Finn 2011; Lee et al. 2012), but all have assumed an equal probability of detection across the sky. While all directions are equally likely to contain a GW source, there should exist more probable locations for the brightest GW source. And while the probability is quite small for there to exist a source that stands out above the background, the analysis in this paper identifies the potential locations of that source using what is currently known about the distribution of galaxies in the local universe, and allows others to enhance discussions on the optimization of a PTA (Anella et al. 2013).
Supermassive black hole binary (SMBHB) systems with periods of months to years are thought to be the most important source of gravitational waves (Jaffe & Backer 2003). Binaries like these form when galaxies containing nuclear black holes merge (Begelman, Blandford, & Rees 1980; Volonteri, Haardt, & Madau 2003). Corresponsingly, the number of such binaries should be greatest where galaxy mergers are more frequent; i.e., in galaxy clusters. All things being equal, the nearest clusters will play host to the brightest sources. This sugests that, as PTA sensitivities increase, we look toward the largest or richest nearby clusters (Virgo, Fornax, Norma, Perseus and Coma) as the most likely location of the earliest detectable SMBHB source. Here we make this expectation quantitative.
Below in §2, we use a compilation of several galaxy surveys to identify the mass, distance, and location of all galaxies within Mpc. From this data we estimate the SMBHB mass and the lifetime of GW emission in a detectable PTA band for each galaxy in §3, and together with the distance, we identify the probability of the existence of a detectable source in a given direction. §4 contains a full explanation of the results with sky position maps. We investigate the probability of detection given an increasing PTA sensitivity in §5 and a summary of our findings is found in §6.
2. Extragalactic Data Base
We require knowledge of the mass and distance of each galaxy in order to estimate the amplitude of the GWs emitted by a SMBHB that may exist at the center of that galaxy. We start by searching the Extragalactic Distance Database (EDD), created by Tully et al. (2009). This database is a compilation of many extragalactic surveys with the intention of compiling all visible galaxies within Mpc (z=0.03), where km/s/Mpc. The EDD111http://edd.ifa.hawaii.edu recently updated its records with the 2M++ galaxy reshift catalog which reaches complete out to 200 Mpc and is conservatively complete to 60 Mpc (Lavaux & Hudson 2011). This gives our sample the same completeness. Additionally, the EDD was updated with an extensive survey of all galaxies within 11 Mpc (Karachentsev, Makarov, & Kaisina 2013), later in the paper local sources will be highlighted and this recent addition to the EDD gives us great confidence in our ability to talk about neighboring galaxies.
The Lyon-Meudon Extragalactic Database222http://leda.univ-lyon1.fr (LEDA) (Paturel et al. 2003) is the largest database that the EDD draws on, but while LEDA compiles over three million objects (Vauglin et al. 2006) the EDD only gathers information from LEDA for galaxies that are found in other surveys (Tully et al. 2009) which at the time of this paper was just under a hundred thousand galaxies. Both the EDD and LEDA take advantage of the labeling mechanism started in the Principal Catalogue of Galaxies (PGC) (Paturel et al. 1989). The PGC number of every known galaxy is used to access the same galaxy across several surveys and to prevent any single galaxy from being counted more than once.
To obtain an accurate distance, the databases use various methods based on the known parameters of each galaxy. There are two main distance moduli calculated in LEDA, “mod0”, which is calculated from a distance catalog using the Tully-Fisher relation or the Faber-Jackson relation, and “modz”, which is calculated using redshift (Paturel et al. 1997). The EDD primarily uses the Tully-Fisher relation to get distances, but compares the results with other distances to assure a common scale (Tully et al. 2009). Distance measurement using the Tully-Fisher relation is only available for about of galaxies, with the rest being derived from redshift measurements.
2.1. Estimating Black Hole Masses
We calculate the total central black hole mass of a galaxy using the M- relation when an accurate central velocity dispersion () is found. This is available for of the galaxies. In an attempt to expand the pool of usable galaxies, B-band luminosity is used as a surrogate for mass. In LEDA, most known galaxies have a calculated total B-band luminosity, which we converted to a total luminosity using a bolometric correction (Buzzoni 2005). The number of galaxies for this work increased from 4,258 to 83,816. Fig 1 compares the luminosity and black hole mass for galaxies with known . While the relationship is most likely non-linear in lower mass regions, since PTAs are only sensitive to a certain range of black hole masses, solar masses (Sesana & Vecchio 2010), we can ignore the low mass trend and only use the clear linear trend in the PTA sensitive region. We advise the reader that the line looks like a poor fit to the data because of the invisible density gradient in the gray region of the plot.
For consistency, all gathered luminosities were converted to black hole mass using the above trend line which gives the power law,
where L is the corrected luminosity of each galaxy.
2.2. Sample Size
In Fig 2, we plot all found galaxies with a recorded B-band luminosity to look for survey biases. This plot was created with equal-area pixels generated using MEALpix333Program developed by GWAstro Research Group and available from http://gwastro.org. The empty areas in the plot are caused by the Milky Way galaxy, whose plane renders these areas of the sky unobservable, also called the zone of avoidance. The largest value pixels in this plot, appearing as black and dark gray, contain galaxy clusters. Significant clusters are labeled on the plots throughout this paper. The Virgo cluster contains approximately 1300 galaxies (Binggeli, Sandage, & Tammann 1985); the Fornax cluster contains around 60 galaxies (Jordán et al. 2007); the Norma cluser contains around 600 galaxies (Woudt et al. 2008); the Perseus cluster contains about 500 galaxies (Brunzendorf & Meusinger 1999); and the Coma cluster contains more than a thousand galaxies (Hammer et al. 2010). Since all pixels not obstructed by the plane of the Milky Way contain a number of galaxies that are within the same order of magnitude, we deem the distribution across the sky to be reasonably non-bias, particularly with respect to nearby sources, which are most important to us. Fig 3 contains a plot of all galaxies with a recorded B-band luminosity in the combined databases as a histogram over distance. We label the mean distance of some galaxy clusters, and in general the larger spikes in galaxies are due to clusters. The number of galaxies grow out to approximately 150 Mpc, and then fall off inversely with distance. From the completeness of the EDD, we feel confident that in this paper we achieve a reasonably accurate representation of our local universe, i.e. the galaxies within 150 Mpc.
Now that we have a significant sample size of galaxies, we cut down to only the galaxies with luminosities great enough to potentially harbor a detectable source; galaxies with central black holes larger than . In an attempt to remove bias to our closest neighbors, we remove the Andromeda galaxy and our galaxy from the dataset as well. The effect of this mass cut can be seen on the gray histogram overlayed on Fig 3, which leaves us with 75,486 galaxies. The region of sky within 50 Mpc, is the most effected by this cut, which makes sense given that less luminous galaxies will only be observable at smaller distances. The Virgo cluster is the only galaxy cluster strongly effected by this cut for similar reasons. After the mass cut was made, we cut down to only the galaxies within 150 Mpc of the Milky Way, which leaves us with 40,560 galaxies, just under half of the starting number. The effect of this distance cut verses sky position can be seen in Fig 4.
To estimate the detectibility of GWs with PTAs, we use the SMBHB mass and the distance of each galaxy to calculate the GW signal strength from each potential source. We also calculate the lifetime of emission in a detectable PTA band, and estimate the number of detectable sources. This allows us to identify the probability of the existence of a source at a given amplitude in a given direction.
We start with the simplified assumption that all galaxies contain a binary with equal mass black holes. Since we are only trying to estimate the “hotness” of sky positions relative to each other, we work in proportionalities. The overall factors are irrelevant at this stage of our understanding. A more mathematical description of these proportionalities can be found in the appendix.
We use MEALpix to divide the sky into equal-area pixels, and associate each galaxy with a particular pixel based on its sky location.
If we assume a given galaxy contains or contained a SMBHB then the probability that the SMBHB exists now is the ratio of the binary’s lifetime, , to the age of the universe. If we further assume that the source is detectable as long as the SMBHB exists then the expectation value of the number density of sources in a certain pixel, , is directly proportional to the lifetime, , of all sources in that pixel;
where N is the total number of galaxies in a given pixel and represents a particular galaxy in that pixel. serves as one of the metrics we use to characterize the ”hotness” of the GW sky. However, says nothing about the relative strength of sources, just the number of them. We therefore also use the metric that is proportional to the expectation value of the power in GWs emitted from a particular pixel. Consider the GW power contained in the residual response from a single source, P; the response in pulsar timing to a gravitational wave of strain is proportional to times the period of the binary (Jenet et al. 2004). However, we assume that all the SMBHBs have the same period, and drop the length of the period from our calculations. Therefore, the residual response is proportional to . The power in the pulsar timing residuals is the square of this response:
The total energy emitted over the lifetime of the source is proportional to the power times the lifetime of the source, . Therefore
where the integral is over all time. When we add up all the sources in a particular pixel we get
where N is the total number of galaxies in a given pixel and represents a particular galaxy in that pixel. In essence, this is the total amount of energy a PTA can expect to receive from this pixel over all time, and after dividing by the age of the universe would be the average power in the pixel. Therefore is also proportional to the expectation value of the power from this pixel and we have an expression for our second metric, ,
The above equations only depend on the GW strain, , and the lifetime of each SMBHB, . We use the standard dipole approximation from Jenet et al. (2004) to estimate the strain ;
where is the chirp mass of the SMBHB. The lifetime of the source, , is given by:
The chirp mass of a binary system, , is proportional to the total mass of the binary, , since we have already assumed that the masses in the binary are roughly equal,
We can now rewrite both the number of detectable sources and the GW signal strength in terms of and , which are the observed quantities gathered in §2.
Maps of the above quantities can be found in the next section.
There are many factors that contribute to the detection of a GW signal by a PTA, and so measuring only one quantity is insufficient to assess the likelihood of detecting a GW source in a particular region of the sky. For example, when calculating the number of detectable sources we rank the galaxies in terms of their likelihood of being ‘on’ during an observation, while when calculating GW signal strength we rank the sources in terms of the overall power they are expected to contribute to the pulsar timing band over the course of their lifetimes. These plots give different pictures of our local universe and are both needed to accurately understand a PTAs probability of making a detection.
4.1. Number Density of Detectable Sources
One way of finding probable locations for PTA detection is to look at the expectation value of the number of sources in any given direction , which we estimate using Eqn 10. This value is proportional to the probability of a detectable source currently being in the PTA band in a given pixel. Using the equations in the appendix, the amount of time each potential source spends in the PTA band is calculated, which is converted to a probability by dividing by the age of the universe. We find the total probability of a single source currently being detectable to be , with the ‘brightest’ location on the sky having a probability of currently containing a single source that stands out about the background in the PTA band. The probability of a detectable source currently being in the PTA band is plotted verses sky position in Fig 5, and as a function of distance in Fig 6. These plots are dominated by galaxy clusters in the distance range of 50 to 100 Mpc. Specifically the Coma, Norma and Perseus clusters, which all contain numerous massive galaxies, while also being at a close enough distance for those galaxies to be resolvable by a PTA. These plots highlight a region of space between 50 and 75 Mpc where a large portion of galaxies are a part of the Centaurus supercluster. This region of space which stretches down towards the Norma cluster, is partially in the zone of avoidance and contains the Great Attractor (Kocevski et al. 2007). While optical observations will likely not reveal the Great Attractor, PTAs have the potential to discover the source of this attraction from GW observations.
4.2. GW Power
is estimated using the total GW power for each source in a particular pixel and at a particular distance, integrated over its lifetime, which effectively weights the power from each source with the probability of whether or not it will be caught ‘on’ during an observation. Estimated using Eqn 11, is plotted verses sky position in Fig 7 and as a function of distance in Fig 8. These plots highlight a handful of local sources that have the potential to dominate a PTAs detection of a single GW source. In an attempt to show just how dominating a single source can be, we removed the brightest source, (RA h, dec , distance Mpc) and overlayed the total power as a function of distance in Fig 8. With the largest source of GW power removed from the plots, a small number of galaxies, specifically those around the Virgo cluster, continue to dominate the picture. Unlike the plots which highlight the region of space between 50 and 100 Mpc, the galaxies with the largest values are within a distance of 20 Mpc.
5. Probability of Detection With A PTA
For a PTA to detect a source, that source must not only be emitting GW radiation in a detectable band, but also at a sufficient amplitude. While the above plots provide an accurate picture of the local landscape for PTA detectable sources, we have not yet taken into account the threshold of a PTA to detect these sources.
Any PTA will have a minimum detectability threshold defined by its sensitivity. This threshold is directly related to the strain amplitude, of a potential GW source. Over time, the PTAs will continually lower this threshold as improvements to timing, data analysis, and amount of data continue. Eventually, virtually all potential sources in Section 4 will be detectable, and the maps presented are the best guides of where to focus efforts. However, the first sources detected will be the sources with the largest .
In this section, we make an educated guess as to where those first sources might be located based on our collected data. We start by considering the ten brightest potential sources, which we label case A. To mimic improvements to the sensitivity, we also consider two more cases, B C, that reduce the initial cutoff value by subsequent factors of . Our goal is to predict which of these bright sources will be detected first. Therefore, once all the sources above a given threshold have been identified, a probability is assigned to each source by considering its lifetime over the total lifetime of the sample, or .
In case A we consider the 10 sources with the highest amplitude according to Eqn 7. The likelihood of one of those sources being detected by a PTA is plotted in Fig 9(a). The pixel with the most likelihood contains two sources. Each of the two sources in that pixel has a lifetime of about one quarter of the total lifetime of all the detectable sources. In Fig 9(d), we see that half of these sources are within 20 Mpc.
When the sensitivity is increased by a factor of in case B, there are thirty three potential sources. Fig 9(b) plots the likelihood that one of these sources is detected. As in case A, the plot is dominated by one very likely pixel, however in this case there is only one source in that pixel and it contains 20 of the total lifetime of all detectable sources. This source has the largest value, discussed in §4.2, and while it has a smaller strain than the initial ten, its lifetime is six times longer than any of the potential sources in case A. This source is at a distance of 2.2 Mpc, and the distance distribution of the other sources is plotted in Fig 9(e), where a third of the potential sources are within 20 Mpc.
In case C, the sensitivity is increased one order of magnitude from case A, and there are 237 potential sources. The likelihood of detecting a certain source is plotted in Fig 9(c). This plot is unlike the Case A and Case B plots since it has sources that are more distributed in clusters. The pixel with the largest total lifetime contains three sources, and contains five percent of the total lifetime of all detectable sources at this cutoff. Fig 9(f) shows that these sources are evenly distributed throughout distance with a small “hotspot” appearing around the distance of the Virgo cluster.
Overall, the first two cases reveal very specific sources that dominate the likelihood of detection. It is important to note that as sensitivity is lowered the new sources that will become detectable have the potential for significantly longer lifetimes and these sources are more likely to be found in galaxy clusters, specifically the Virgo cluster.
In this paper, data is gathered from a compilation of at least complete galactic surveys out to 200 Mpc. These galaxies are all assumed to contain a SMBHB with equal mass black holes, and the total central mass is calculated using B-band luminosity as a surrogate for mass. This data set is then cut down to only detectable sources (i. e. sources with a total central mass larger than and with in a distance of 150 Mpc). From this data, two metrics are used to estimate the “hotness” of the GW sky: , which is proportional to the number density of sources in a given direction; and , which is proportional to the expectation value of the power in GWs emitted from a particular direction. Using , we are able to calculate the specific probability that each pixel contains a currently radiating GW source in the PTA band. While this number is very small for any individual pixel, we find a total probability of that one of the galaxies we considered is detectable. The ‘brightest’ location on the sky has a probability of currently containing a single source that stands out about the background in the PTA band. Overall the distribution of single sources potentially detectable by a PTA has a larger number density around local galaxy clusters. While the GW signal strength is dominated by a handful of sources, with the region of sky around the Virgo cluster having a larger number density of these ‘bright’ sources.
Work by Anella et al. (2013) has been done to show that if a new powerful timing pulsar was discovered in the direction of the Virgo cluster, a PTA would have twice the sensitivity to a region about around the discovered pulsar. Given that the sensitivity of a single pulsar to a GW source falls off as , where is the angle between a pulsar and a GW source (Burt, Lommen, & Finn 2011), we recommend focusing the search for new pulsars in the vicinity of the Virgo cluster.
As stated earlier in this paper, this is only a broad estimate that assumed an equal probability for all galaxies to contain a binary with equal mass black holes. Future work is planned to incorporate work done by Rosado & Sesana (2013), which uses the Sloan Digital Sky Survey and the Millenium simulation data (Springel et al. 2005) to search for SMBHBs in the redshift range of . Combining these two data sets is the next step towards creating realistic population distributions for single source GWs detectable by PTAs. While the distribution of detectable single sources will most likely scale with any new estimate, is more affected by these specific probabilities, and a future paper will address these factors.
In the above, the integration should be performed over that region of , , and where GWs would be detectable by a given PTA configuration.
Given the amplitude of the gravitational wave strain, , emitted by a particular SMBH system, the strength of the stochastic GW signal may be written as:
We want to make estimates of the relative strength of both and using observational data from galaxy surveys. In order to do this, we will make the following assumptions: 1) is proportional to the number of observable galaxies per unit solid angle in a given direction, 2) the evolution of the SMBH binaries are dominated by the effects of GW emission, 3) the probability of a galaxy harboring a SMBH is the same for all galaxies, 4) the chirp mass of the binary system is proportional to the total luminosity of the galaxy, , 5) all galaxies of interest have .
Using the fact that is independent of frequency, equation 1 may be written as:
where is the total time a SMBH with chirp mass is detectable in a given PTA and is given by:
The frequencies and are the lowest and highest detectable frequencies given the sensitivity of a PTA and a particular SMBH binary system. Note that both and depend on and . These frequencies are calculated using the following expression for the GW strain amplitude:
together with an expression for the minimum detectable strain of a PTA.
Let be the total number of galaxies per unit total luminosity, per unit red shift, per unit solid angle. Assumptions one and three allow us to write where is a constant. The total number of galaxies per unit solid angle in a given direction is given by
With assumption 4, we can write where is a constant. Putting this all into equation 3, we have
Given that is the fractional number of galaxies in a given direction with total luminosity between and located at a red shift between and , we can write the detectable number of systems as
where the angle brackets represent an average over all galaxies in a particular direction. Since , we find that
where , defined in the above equation, is a quantity that may be estimated from galaxy catalogue data and is directly proportional to the number of detectable sources per unit solid angle in a given direction. We can also define an idealized statistic which is valid for the case where the PTA can detect all frequencies down to some cuttoff frequency. In this case, we can ignore the frequency dependent terms and use
Like , is proportional to but it only depends on the properties of the galaxy distribution and not the details of a particular PTA. An estimate for is simply the sum of over all galaxies in a particular small region of the sky.
Using the same assumptions as with equation 2, one can show that:
Here, the range of integration is over all at a fixed frequency , unlike the expression for . Using assumption 2), we know that
This, together with equation 5 gives us the following:
where as defined above may be estimated from existing data and is proportional to the strength of the stochastic GW emission in a particular direction. An estimate for is obtained by summing over all galaxies in a particular direction.
Maps of both the number of detectable sources, and the strength of the stochastic GW signal, , calculated using data from extragalactic databases are presented in §4.
- Anella et al. (2013) Anella, R., Christy, B., Lommen, A., & Finn, L. S. 2013, ApJ in preparation
- Begelman, Blandford, & Rees (1980) Begelman, M. C., Blandford, R. D., & Rees, M. J. 1980
- Binggeli, Sandage, & Tammann (1985) Binggeli, B., Sandage, A., & Tammann, G. A. 1985, AJ, 90, 1681–1759
- Brunzendorf & Meusinger (1999) Brunzendorf, J., & Meusinger, H. 1999, A&AS, 139, 141–161
- Burt, Lommen, & Finn (2011) Burt, B. J., Lommen, A. N., & Finn, L. S. 2011, ApJ, 730, 17
- Buzzoni (2005) Buzzoni, A. 2005, MNRAS, 361, 725–742
- Detweiler (1979) Detweiler, S. 1979, ApJ, 234, 1100–1104
- Ellis, Siemens, & Creighton (2012) Ellis, J. A., Siemens, X., & Creighton, J. D. E. 2012, ApJ, 756, 175
- Foster & Backer (1990) Foster, R. S., & Backer, D. C. 1990, ApJ, 361, 300–308
- Hammer et al. (2010) Hammer, D., Verdoes Kleijn, G., Hoyos, C., den Brok, M., Balcells, M., Ferguson, H. C., Goudfrooij, P., Carter, D., Guzmán, R., Peletier, R. F., Smith, R. J., Graham, A. W., Trentham, N., Peng, E., Puzia, T. H., Lucey, J. R., Jogee, S., Aguerri, A. L., Batcheldor, D., Bridges, T. J., Chiboucas, K., Davies, J. I., del Burgo, C., Erwin, P., Hornschemeier, A., Hudson, M. J., Huxor, A., Jenkins, L., Karick, A., Khosroshahi, H., Kourkchi, E., Komiyama, Y., Lotz, J., Marzke, R. O., Marinova, I., Matkovic, A., Merritt, D., Miller, B. W., Miller, N. A., Mobasher, B., Mouhcine, M., Okamura, S., Percival, S., Phillipps, S., Poggianti, B. M., Price, J., Sharples, R. M., Tully, R. B., & Valentijn, E. 2010, ApJS, 191, 143–159
- Hobbs et al. (2010) Hobbs, G., Archibald, A., Arzoumanian, Z., Backer, D., Bailes, M., Bhat, N. D. R., Burgay, M., Burke-Spolaor, S., Champion, D., Cognard, I., Coles, W., Cordes, J., Demorest, P., Desvignes, G., Ferdman, R. D., Finn, L., Freire, P., Gonzalez, M., Hessels, J., Hotan, A., Janssen, G., Jenet, F., Jessner, A., Jordan, C., Kaspi, V., Kramer, M., Kondratiev, V., Lazio, J., Lazaridis, K., Lee, K. J., Levin, Y., Lommen, A., Lorimer, D., Lynch, R., Lyne, A., Manchester, R., McLaughlin, M., Nice, D., Oslowski, S., Pilia, M., Possenti, A., Purver, M., Ransom, S., Reynolds, J., Sanidas, S., Sarkissian, J., Sesana, A., Shannon, R., Siemens, X., Stairs, I., Stappers, B., Stinebring, D., Theureau, G., van Haasteren, R., van Straten, W., Verbiest, J. P. W., Yardley, D. R. B., & You, X. P. 2010, Classical and Quantum Gravity, 27, 084013
- Jaffe & Backer (2003) Jaffe, A. H., & Backer, D. C. 2003, ApJ, 583, 616–631
- Jenet et al. (2004) Jenet, F. A., Lommen, A., Larson, S. L., & Wen, L. 2004, ApJ, 606, 799–803
- Jordán et al. (2007) Jordán, A., Blakeslee, J. P., Côté, P., Ferrarese, L., Infante, L., Mei, S., Merritt, D., Peng, E. W., Tonry, J. L., & West, M. J. 2007, ApJS, 169, 213–224
- Karachentsev, Makarov, & Kaisina (2013) Karachentsev, I. D., Makarov, D. I., & Kaisina, E. I. 2013, AJ, 145, 101
- Kocevski et al. (2007) Kocevski, D. D., Ebeling, H., Mullis, C. R., & Tully, R. B. 2007, ApJ, 662, 224–235
- Lavaux & Hudson (2011) Lavaux, G., & Hudson, M. J. 2011, MNRAS, 416, 2840–2856
- Lee, Jenet, & Price (2008) Lee, K. J., Jenet, F. A., & Price, R. H. 2008, ApJ, 685, 1304–1319
- Lee et al. (2012) Lee, K. J., Bassa, C. G., Janssen, G. H., Karuppusamy, R., Kramer, M., Smits, R., & Stappers, B. W. 2012, Monthly Notices of the Royal Astronomical Society, 423, 2642–2655
- Paturel et al. (1989) Paturel, G., Fouque, P., Bottinelli, L., & Gouguenheim, L. 1989, A&AS, 80, 299–315
- Paturel et al. (1997) Paturel, G., Andernach, H., Bottinelli, L., di Nella, H., Durand, N., Garnier, R., Gouguenheim, L., Lanoix, P., Marthinet, M. C., Petit, C., Rousseau, J., Theureau, G., & Vauglin, I. 1997, A&AS, 124, 109–122
- Paturel et al. (2003) Paturel, G., Petit, C., Prugniel, P., Theureau, G., Rousseau, J., Brouty, M., Dubois, P., & Cambrésy, L. 2003, A&A, 412, 45–55
- Rosado & Sesana (2013) Rosado, P. A., & Sesana, A. 2013, ArXiv e-prints
- Sazhin (1978) Sazhin, M. V. 1978, Soviet Ast., 22, 36–38
- Sesana (2013) Sesana, A. 2013, Classical and Quantum Gravity, 30, 244009
- Sesana & Vecchio (2010) Sesana, A., & Vecchio, A. 2010, Phys. Rev. D, 81, 104008
- Springel et al. (2005) Springel, V., White, S. D. M., Jenkins, A., Frenk, C. S., Yoshida, N., Gao, L., Navarro, J., Thacker, R., Croton, D., Helly, J., Peacock, J. A., Cole, S., Thomas, P., Couchman, H., Evrard, A., Colberg, J., & Pearce, F. 2005, Nature, 435, 629–636
- Tully et al. (2009) Tully, R. B., Rizzi, L., Shaya, E. J., Courtois, H. M., Makarov, D. I., & Jacobs, B. A. 2009, AJ, 138, 323–331
- Vauglin et al. (2006) Vauglin, I., Prugniel, P., Courtois, H., Makarov, D., Petit, C., Mamon, G., & Paturel, G. 2006, In SF2A-2006: Semaine de l’Astrophysique Francaise, D. Barret, F. Casoli, G. Lagache, A. Lecavelier, and L. Pagani, eds., p. 365
- Volonteri, Haardt, & Madau (2003) Volonteri, M., Haardt, F., & Madau, P. 2003, The Astrophysical Journal, 582, 559
- Woudt et al. (2008) Woudt, P. A., Kraan-Korteweg, R. C., Lucey, J., Fairall, A. P., & Moore, S. A. W. 2008, Monthly Notices of the Royal Astronomical Society, 383, 445–457
- Yardley et al. (2010) Yardley, D. R. B., Hobbs, G. B., Jenet, F. A., Verbiest, J. P. W., Wen, Z. L., Manchester, R. N., Coles, W. A., van Straten, W., Bailes, M., Bhat, N. D. R., Burke-Spolaor, S., Champion, D. J., Hotan, A. W., & Sarkissian, J. M. 2010, MNRAS, 407, 669