Seedless clustering in allsky searches for gravitationalwave transients
Abstract
The problem of searching for unmodeled gravitationalwave bursts can be thought of as a pattern recognition problem: how to find statistically significant clusters in spectrograms of strain power when the precise signal morphology is unknown. In a previous publication, we showed how “seedless clustering” can be used to dramatically improve the sensitivity of searches for longlived gravitationalwave transients. In order to manage the computational costs, this initial analysis focused on externally triggered searches where the source location and emission time are both known to some degree of precision. In this paper, we show how the principle of seedless clustering can be extended to facilitate computationallyfeasible, allsky searches where the direction and emission time of the source are entirely unknown. We further demonstrate that it is possible to achieve a considerable reduction in computation time by using graphical processor units (GPUs), thereby facilitating more sensitive searches.
pacs:
95.75.z,04.30.w,07.05.BxI Introduction
Longlived gravitationalwave transients (lasting ) constitute an interesting class of signals for secondgeneration detectors such as Advanced LIGO Harry, G. M. for the LIGO Scientific Collaboration (2010) and Advanced Virgo The Virgo Collaboration (2009). After reaching design sensitivity, Advanced LIGO expects to observe binary neutron stars mergers and neutronstar blackhole coalescences per year of science data Abadie et al. (2010). The standard searches for compact binary coalescences rely on matched filter template banks; see, e.g., Abadie et al. (2012a, b). More exotic sources of longlived transients, including emission from rotational instabilities in protoneutron stars Piro and Thrane (2012); Piro and Ott (2011); Piro and Pfahl (2007); Corsi and Mészáros (2009) and blackhole accretion disk instabilities Kiuchi et al. (2011); van Putten (2001, 2008), cannot be accurately modeled owing to theoretical uncertainties. However, searches for longlived bursts Aasi et al. (2013); Thrane et al. (2011); Thrane and Coughlin (2013) can be employed when a matched filter search is not possible. (There is a rich literature on short, subsecond gravitationalwave bursts and the different detection strategies available to detect them, but we focus here on longlived transients.)
In a crosscorrelation search such as Aasi et al. (2013); Thrane et al. (2011), the detection of gravitational waves can be thought of as a pattern recognition problem. The goal is to find tracks of excess strain crosspower, which appear as brighterthanexpected pixels on a signaltonoise ratio spectrogram (map). In a previous work Thrane and Coughlin (2013), we described how “seedless clustering” can be used to significantly enhance the sensitivity of searches for longlived transients when a trusted matched filter template bank is not available. We review the details of seedless clustering in Sec. III, but the basic idea is to integrate along many different cleverly chosen paths in a signaltonoiseratio spectrogram. This is in contrast to seedbased clustering algorithms which form clusters from bright spectrogram pixels called “seeds.”
The advantage of seedless clustering is most pronounced for long and weak signals Thrane and Coughlin (2013). For the waveforms considered in Thrane and Coughlin (2013), we found that seedless clustering can detect a gravitationalwave signal (at a fixed falsealarm and false dismissal rate) at a distance between – further than a seedbased clustering algorithm. This corresponded to an increased detection volume of –.
One of the challenges associated with seedless clustering is that it is, as a rule of thumb, more computationally expensive than seedbased alternatives. In Thrane and Coughlin (2013), we focused on applications to targeted searches, in which the sky location is tightly constrained and the time of the event is known to exist in some “onsource” window, thereby saving the extra computational cost associated with searching many sky positions and emission times.
In this work, we show how the seedless clustering formalism from Thrane and Coughlin (2013) can be extended to a highsensitivity, computationallyefficient, allsky search for longlived gravitational waves from arbitrary sky locations. There are two innovations which make this possible. First, by introducing a new random phase factor, we show that it is possible to efficiently scan the entire sky with a seedless clustering algorithm. Second, we take advantage of recent advances in computing to carry out our computations on graphical processor units (GPUs). Seedless clustering algorithms are “embarrassingly parallel” Foster (1995), which allows them to exploit the highly parallel architecture of GPUs. We show that an allsky search with seedless clustering is both computationally feasible and more sensitive than a seedbased algorithm. We outline the computational requirements for a realistic search and demonstrate the advantage of carrying out computations on GPUs.
The remainder of the paper is organized as follows. In Sec. II, we describe some of the general features and challenges of an allsky transient search. In Sec. III, we describe allsky stochtrack, an allsky algorithm which employs seedless clustering. In Sec. IV, we present the results of a sensitivity study comparing allsky stochtrack to a seedbased algorithm. In Sec. V, we describe the computational resources required for realistic searches and compare the algorithms’ performance on CPUs and GPUs. In Sec. VI, we offer concluding remarks and suggest directions for future research.
Ii The challenges of allsky radiometry
In this section, we outline some of the general features of an allsky transient search built on the principle of radiometry—using the time delay between two detectors to search data associated with a specific direction in the sky. We begin with strain time series and from detectors and , which are separated by a displacement . The data are split into segments (typically with a duration of ) and Fouriertransformed to create complexvalued strain spectrograms: and . Note that refers to sampling time whereas refers to segment start time.
Following Thrane et al. (2011); Thrane and Coughlin (2013); Aasi et al. (2013), the signaltonoise ratio spectrogram can be written as:
(1) 
Here, is a phase factor, which takes into account the time delay between detectors and ; is the speed of light. The phase factor rotates the crosspower signal in the complex plane so as to be real and positive. The term is a normalization factor, which uses neighboring segments to estimate the background ^{1}^{1}1 The factor may also include a directiondependent phase factor taking into account the relationship between the and polarizations of an elliptically polarized source. For the sake of simplicity, we do not include this additional phase factor. We expect the sensitivity to improve marginally by adding this phase factor by incorporating additional polarization information, though, at an increased computational cost. . Precise definitions of and are provided in Appendix A. If the source direction is known, for example, from an electromagnetic trigger (see Thrane and Coughlin (2013)), then it is straightforward to apply to the appropriate phase factor. When no electromagnetic trigger is available, it is necessary to search over multiple directions.
Consider the case where the source is located at but the filter is chosen for the direction , which introduces a timing error of:
(2) 
On average, the timing error reduces the signaltonoise ratio by
(3) 
Inspecting Eq. 3, we can infer the qualitative features of a signal in a spectrogram characterized by a timing error . For small values of , the apparent signal will be weaker than it would in the absence of a timing error. This is because some of the crosspower in Eq. 1 leaks into the imaginary direction. As crosses , the signal vanishes entirely before reappearing as negative signaltonoise ratio.
Graphically, large timing errors produce characteristic stripes in spectrograms; see Fig. 1. The bandwidth of each stripe is given by . The minimum stripe size is , where is the travel time between detectors and . For the two LIGO detectors, .
We can define a tolerance for the maximum possible timing error by requiring that we observe no less than, say, of the signaltonoise ratio. It follows that
(4) 
As frequency increases, the tolerable timing error decreases. For signals in the most sensitive part of LIGO’s band , the timing tolerance is . For highfrequency signals near , it is .
To summarize, radiometry relies on the use of a phase factor to characterize the time delay between two detectors (a function of source sky position). If the assumed sky position is incorrect, gravitationalwave crosspower leaks from the positive real direction into imaginary and/or negative components, which, in turn, leads to a reduced signaltonoise ratio. It may therefore be necessary to search many directions (with many time delays) in order to observe the signal with an acceptable signaltonoise ratio.
Iii Algorithm
In this section, we describe the details of an algorithm which uses seedless clustering to search for transient signals from all directions on the sky. We call it allsky stochtrack. We begin with a brief review of the stochtrack algorithm Thrane and Coughlin (2013), which will serve as a foundation on which to build.
The goal is to find the most significant cluster as determined by the value of the detection statistic ^{2}^{2}2 We note that alternative definitions of (with a different weighting scheme) are also possible; see Thrane and Coughlin (2013). ,
(5) 
where is defined in Eq. 1 and is the number of pixels in .
In any seedless clustering algorithm, is determined a priori by some set of rules (as opposed to by the data itself). In the stochtrack algorithm Thrane and Coughlin (2013), is chosen randomly from the set of quadratic Bézier curves Farin (1996) subject to the constraint that the curve persists for a duration . (Other parameterizable curves such as spline can be used as well.) Each randomly selected Bézier curve is referred to as a “template.” Each Bézier template is described by three timefrequency control points: , , and . The control points form a curve parameterized by :
(6) 
In order for the algorithm to have a high probability of guessing a close approximation to the true signal, many templates must be used. Fortunately, each template can be quickly generated from just six random numbers. By working with arrays of Bézier curves, stochtrack is able to carry out the sum in Eq. 5 for a large number of templates in parallel. As we shall see below, the parallel nature of the calculation lends itself to the use of GPUs.
The number of templates is denoted . For practical applications, it is typically chosen to be . In Thrane and Coughlin (2013), we described a default search with and a deep search (denoted stochtrack ) with ^{3}^{3}3 In Thrane and Coughlin (2013), the default and deep search are said to use and respectively. We believe the correct numbers are actually and as stated here. .
The allsky stochtrack algorithm builds on the foundation of stochtrack. First, we introduce “complex signaltonoise ratio”:
(7) 
This is necessary in order to preserve the complex phase information that encodes the direction of the source. Note that, unlike , is not defined for a particular direction.
Next, in addition to the six random control points, we add an additional random variable corresponding to the time delay between the detectors, which, as we saw in Sec. II, is a proxy for sky location. If we assume that the sky location of each transient is drawn from an approximately isotropic distribution, then the probability density function for time delay is a simple uniform distribution between .
Finally, we modify Eq. 5 to be
(8) 
The new sum described in Eq. 8 is carried out for many randomly selected clusters , each with a randomly selected time delay . By including a random time delay, the algorithm tries to guess not only the spectrographic shape of the signal, but also the appropriate phase factor that will minimize the timing error stripes shown in Fig. 1. Minimizing the timing error maximizes .
The addition of a new random variable comes at a cost. First, the allsky stochtrack algorithm will converge less quickly than stochtrack due to its expanded parameter space. Second, even if we imagine setting , allsky stochtrack templates span a larger space than the templates used in stochtrack, and so allsky stochtrack must contend with comparatively higher background. That said, we find that the extra cost is small. In the next section, we show that, for several signal models, allsky stochtrack achieves a sensitivity which is only slightly less than stochtrack, while searching a significantly expanded signal space.
Iv Sensitivity Study
In this section, we describe a study to determine the sensitivity of allsky stochtrack to four different longlived test waveforms using Monte Carlo and recolored LIGO noise. In order to compare with the baseline sensitivity of stochtrack, we use the same four waveforms used in Thrane and Coughlin (2013): two fallback accretion signals Piro and Thrane (2012) abbreviated FA 1 and FA 2, and two accretiondisk instability waveforms Santamaría and Ott (2011), abbreviated ADI 1 and ADI 2. The waveforms span durations of – and range in frequency from –. Additional details about the waveforms and the parameters used to generate them are provided in Appendix B. Additional information about the models behind the waveforms is available in Piro and Thrane (2012) and Santamaría and Ott (2011); see also Piro and Ott (2011); Piro and Pfahl (2007); Corsi and Mészáros (2009); Kiuchi et al. (2011); van Putten (2001, 2008).
Each waveform is injected into either Monte Carlo or recolored noise (initial LIGO noise which has been recolored to match the design sensitivity of Advanced LIGO while preserving nonstationary noise artifacts) ^{4}^{4}4 The data are taken in between GPS times 822917487 and 847549782. . The data are processed to form a complex signaltonoise ratio spectrogram (see Eq. 7). Following Thrane and Coughlin (2013), the ADI waveforms are analyzed in a band between – while the FA waveforms are analyzed in a band between –. The spectrogram resolution is except for FA 1, for which we use . Each spectrogram corresponds to of data. Data segments are constructed with overlapping Hann windows.
We characterize the sensitivity in terms of a detection distance, defined as the distance to which allsky stochtrack can detect a source with a falsealarm probability and a falsedismissal probability . We perform two series of tests. First, in order to compare allsky stochtrack with stochtrack, we inject each signal with an optimal orientation (faceon) and in an optimal sky location (where the detectors are most sensitive). The true source location is provided as input to stochtrack (and to the seedbased clustering algorithm, burstegard Prestegard and Thrane (2012)), but the allsky stochtrack algorithm is not provided any information about the true location of the source.
We show stochtrack results (from Thrane and Coughlin (2013)) for the default search () and for the deep search (), which is labeled: stochtrack . We compare these to new results obtained with the default allsky stochtrack () and allsky stochtrack ().
In the second series of tests, we inject signals at random sky locations (chosen from an isotropic distribution) and with random inclination and polarization angles . We expect that the detection distance for signals recovered with random values of will be of what is achieved for optimal sources based on the antenna response of our twodetector network.
Our hypothetical network consists of the Advanced LIGO detectors in Hanford, WA (H1) and Livingston, LA (L1) Harry, G. M. for the LIGO Scientific Collaboration (2010). We assume both detectors are operating at design sensitivity.
The results of the study are summarized in Tables 1 and 2 for Monte Carlo and recolored noise respectively. The Monte Carlo and recolored noise are processed identically except we apply a glitch identification Prestegard et al. (2012) cut when analyzing recolored noise ^{5}^{5}5In order to apply the algorithm from Prestegard et al. (2012), we assume that the source is optimally oriented with an optimal sky position.. For optimally oriented sources injected into Monte Carlo noise, we find that the allsky stochtrack can see sources – further than the seedbased burstegard, even though the burstegard algorithm is given the known sky location whereas allsky stochtrack is not. This corresponds to increased detection volume of –. For recolored noise, the improvement is – in distance and – in volume.
Repeating the Monte Carlo analysis with the computationally cheaper default version of allsky stochtrack (), we obtain distances of – times the distances obtained using burstegard. For recolored noise, these distances are – times the values obtained using burstegard. Note that while burstegard can detect the FA 1 waveform in recolored noise at greater distances than the default version of allsky stochtrack, this is very likely because the burstegard algorithm is supplied with the true source location. In an applestoapples comparison, seedless clustering using the default stochtrack is more sensitive than burstegard Thrane and Coughlin (2013).
Thus, the fact that burstegard can detect FA 1 signals in recolored noise at greater distances than allsky stochtrack is telling us that it is very useful to know where in the sky to look when trying to find FA 1 waveforms in recolored noise. This is, perhaps, not surprising since the FA 1 waveform is shorter and spans a greater bandwidth than the other waveforms we consider. Shorter signals are more prone to resemble nonstationary noise. Signals with larger bandwidth are more prone to loss of signal from phase factor mismatch (see Eq. 8).
We also present detection distance for sources with random values of , which are between – of the values obtained for the case of an optimal source.
waveform  algorithm  distance  volume  

absolute  %  %  
ADI 1  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  —  
ADI 2  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  —  
FA 1  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  —  
FA 2  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  — 
waveform  algorithm  distance  volume  

absolute  %  %  
ADI 1  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  —  
ADI 2  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  —  
FA 1  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  —  
FA 2  targeted seedbased  
targeted seedless  
targeted seedless  
allsky seedless  
allsky seedless  
…w/ random  —  — 
V Computing
The results from Section IV were obtained using graphical processor units (GPUs) on the LIGO Data Grid. In this section, we document how GPUs provide an efficient architecture for carrying out stochtrack and allsky stochtrack calculations. We compare the performance of the algorithm using both GPUs and CPUs. We utilize Kepler GK104s GPUs, which are capable of peak single precision floating point performance of according to the manufacturer. Each GPU card has memory. We use Intel Xeon E54650 CPUs.
For our benchmark test, we analyze spectrograms consisting of pixels () using the same deepsearch settings used to analyze the ADI 1 waveforms in the previous section. The computation time includes inputoutput tasks and other calculations, which do not take advantage of the GPU architecture. However, these computations correspond to a tiny fraction () of the total computation time. The results are summarized in Table. 3. We find that allsky stochtrack calculations can be carried out faster on GPUs than CPUs.
hardware  computation time 

CPU  
GPU 
Using our benchmark tests, we estimate the computational requirements for fullfledged gravitational searches running stochtrack and allsky stochtrack on GPUs. (Interestingly, stochtrack and allsky stochtrack take about the same time to run given identical parameters.) We consider two analyses: one targeted (using stochtrack) and one allsky (using allsky stochtrack). For both analyses, we assume an analysis band of (following Aasi et al. (2013)). For the targeted analysis, we assume that search analyzes external triggers, e.g., from gammaray bursts; see Aasi et al. (2013). Following Aasi et al. (2013), we assume that the search is carried out in a wide onsource window. For the targeted analysis, we further assume that timeshift analyses are carried out in order to evaluate the significance of candidate events; see, e.g., Abbott et al. (2004). For the allsky analysis, we assume .
Before we present estimates of computational cost, it will be useful to define a new variable: , the number of templates per of bandwidth. This variable is useful since, all else equal, bigger bands must be analyzed with more templates than smaller bands due to the increased size of the template parameter space. We chose in order to facilitate comparisons with the ADI 1 and ADI 2 results given in Tables 1 and 2. However, we note that waveforms FA 1 and FA 2 are analyzed in a wide band, six times wider than the ADI analysis band. Thus, corresponds to (the default search) in the FA 1 and FA 2 analysis analysis band. corresponds to (more sensitive than stochtrack ) in the FA 1 and FA 2 analysis analysis band.
Given our assumptions, the estimated computational time for a triggered stochtrack search with GPUs is:
(9) 
Here is the number of GPUs. The estimated computational time for an allsky search with GPUs is:
(10) 
From Eq. 9, we conclude that GPUs can facilitate a deepsearch sensitivity with stochtrack using modest computational resources. From Eq. 10, we conclude that a yearlong allsky analysis with defaultsensitivity allsky stochtrack can also be carried out using reasonable computational resources.
Since we know that allsky stochtrack sensitivity can improve significantly with added templates, it would be advisable to follow up on of the loudest events identified by the allsky analysis, with a deeper search. This would add only a marginal increase to the computational burden while ensuring that a marginal detection is promoted to a strong detection (or revealed to be a noise fluctuation).
The sensitivity of an allsky search with allsky stochtrack can be increased after the analysis has commenced (supposing, for example, that more GPUs become available) through the use of intermediate data files. Namely, we recommend recording for each spectrogram. If multiple runs of the analysis are carried out, one can choose the largest value of among each run and for every spectrogram in a simple postprocessing step ^{6}^{6}6Do not forget to use a different random seed for each run, good reader.. In other words, it is easy to combine the results from three runs with in order to obtain results identical to a single search. This parallelizability can be exploited to plan for a computationally conservative analysis, while being ready for a more aggressive analysis, should the resources be available.
Vi Conclusions and Future Work
In previous work, we proposed a new seedless clustering algorithm called stochtrack and demonstrated how it could significantly improve the sensitivity of searches for longlived, unmodeled gravitationalwave transients. Here we extend the principle of stochtrack to the case of an allsky search, when there is no external trigger telling us where on the sky to look. We compare the sensitivity of allsky stochtrack to that of a seedbased algorithm (which takes the true sky direction as input), and find that, for the most part, allsky stochtrack is significantly more sensitive, even though it is searching for the signal in a much larger parameter space.
We point out that stochtrack and allsky stochtrack are “embarrassingly parallel” algorithms and we perform benchmark tests using CPUs and GPUs. We find that GPUs can carry out stochtrack and allsky stochtrack calculations ten times faster than CPUs. We estimate the computational cost of realistic analyses, and show that interesting investigations can be carried out in a reasonable amount of time with a modest number of GPUs.
While we present allsky stochtrack as a tool for allsky analyses, it should also be very helpful in targeted analyses in which the sky localization of the external trigger is large compared to the pointspread function of the gravitationalwave detector network. Instead of drawing the time delay variable (see Eq. 8) from a distribution derived from an isotropic prior, it is straightforward to draw it from a distribution corresponding to a particular patch of sky. This hybrid solution provides an efficient alternative to running stochtrack for many different directions.
We previously mentioned in Thrane and Coughlin (2013) the possibility of using stochtrack to search for compact binaries. In general, compact binaries can be wellmodeled, and so it is expected that matched filtering is the optimal search strategy. However, there are good reasons to explore alternative methods:

Improved robustness and redundancy with an alternative method.

Investigate potentially challenging corners of parameter space, e.g., systems with nonnegligible spin and/or eccentricity.

Detect exotic systems and/or new physics which are not included in matched filter template banks.
In order to place this discussion in context and to motivate future work, we close by reporting the results of a sensitivity study for detecting the coalescence of two neutron stars with stochtrack. We consider the case of an optimally oriented system at an optimal sky location. We assume the signal is confined to a onsource region as in previous searches triggered by gammaray bursts Abadie et al. (2012a). We find that such a binary neutron star coalescence can be detected in Advanced LIGO Monte Carlo noise using stochtrack with and at a distance of . For comparison, the best upper limit from initial LIGO and Virgo on binary neutron star coalescence coincident with gamma ray bursts is Abadie et al. (2012a). It is probable that the sensitivity of stochtrack to binary neutron stars can be enhanced with additional tuning; see Fig. 2. Thus, the application of stochtrack and allsky stochtrack to compact binary coalescence signals appears promising and worthy of future work.
Acknowledgements.
We thank Stuart Anderson, Juan Barayoga, and Fred Donovan for assistance with GPUs. We thank Anthony Piro for sharing the fallback accretion waveforms used in this analysis. We thank Tanner Prestegard for helpful comments on a draft of this paper. ET is a member of the LIGO Laboratory, supported by funding from United States National Science Foundation. LIGO was constructed by the California Institute of Technology and Massachusetts Institute of Technology with funding from the National Science Foundation and operates under cooperative agreement PHY0757058. MC is supported by National Science Foundation Graduate Research Fellowship Program, under NSF grant number DGE 1144152. This paper has been assigned document number LIGOP1400010.Appendix A Formalism
The signaltonoise ratio spectrogram can be written as
(11) 
Here is an estimator for crosspower:
(12) 
and is an estimator for its variance
(13) 
Here is the normalization from a discrete Fourier transform and is a filter function, which accounts for the time delay between detectors and as well as the detector responses. Typically is defined such that is an unbiased estimator for gravitationalwave power Thrane et al. (2011). The variables and are the autopower spectral densities for detectors and in the segments neighboring .
Appendix B Model Parameters
This section reproduces details about the test waveforms from Thrane and Coughlin (2013). The FA waveforms Piro and Ott (2011); Piro and Thrane (2012) are described by the following parameters: initial protoneutron star mass , maximum neutron star mass , a dimensionless factor related to the supernovae explosion energy –, and the radius of the protoneutron star . The values of these parameters for FA 1 and FA 2 are given in Table 5. The ADI waveforms Santamaría and Ott (2011) are parameterized by black hole mass , dimensionless spin parameter , the fraction of the accretion disk mask that forms clumps –, and the torus mass . The values of these parameters for ADI 1 and ADI 2 are given in Table 6.
waveform  duration (s)  – (Hz)  

ADI 1  –  
ADI 2  –  
FA 1  –  
FA 2  – 
References
 Harry, G. M. for the LIGO Scientific Collaboration (2010) Harry, G. M. for the LIGO Scientific Collaboration, Classical Quantum Gravity 27, 084006 (2010).
 The Virgo Collaboration (2009) The Virgo Collaboration, Advanced Virgo Baseline Design (2009), URL https://tds.egogw.it/itf/tds/file.php?callFile=VIR0027A09.%pdf.
 Abadie et al. (2010) J. Abadie et al., Class. Quant. Grav 27, 173001 (2010).
 Abadie et al. (2012a) J. Abadie et al., Astrophys. J. 760, 12 (2012a).
 Abadie et al. (2012b) J. Abadie et al., Phys. Rev. D 85, 082002 (2012b).
 Piro and Thrane (2012) A. L. Piro and E. Thrane, Astrophys. J. 761, 63 (2012).
 Piro and Ott (2011) A. L. Piro and C. D. Ott, Astrophys. J. 736, 108 (2011).
 Piro and Pfahl (2007) A. L. Piro and E. Pfahl, Astrophys. J. 658, 1173 (2007).
 Corsi and Mészáros (2009) A. Corsi and P. Mészáros, Astrophys. J. 702, 1171 (2009).
 Kiuchi et al. (2011) K. Kiuchi, M. Shibata, P. J. Montero, and J. A. Font, Phys. Rev. Lett. 106, 251102 (2011).
 van Putten (2001) M. H. P. M. van Putten, Phys. Rev. Lett. 87, 091101 (2001).
 van Putten (2008) M. H. P. M. van Putten, Astrophys. J. Lett. 684, 91 (2008).
 Aasi et al. (2013) J. Aasi et al., Phys. Rev. D 88, 122004 (2013).
 Thrane et al. (2011) E. Thrane, S. Kandhasamy, C. D. Ott, et al., Phys. Rev. D 83, 083004 (2011).
 Thrane and Coughlin (2013) E. Thrane and M. Coughlin, Phys. Rev. D 88, 083010 (2013).
 Foster (1995) I. Foster, Designing and Building Parallel Programs: Concepts and Tools for Parallel Software Engineering (Addison Wesley, 1995).
 Farin (1996) G. Farin, Curves and Surfaces for CAGD, Fourth Edition: A Practical Guide (Academic Press, 1996).
 Santamaría and Ott (2011) L. Santamaría and C. D. Ott, LIGO DCC p. T1100093 (2011), https://dcc.ligo.org/LIGOT1100093v2/public.
 Prestegard and Thrane (2012) T. Prestegard and E. Thrane, LIGO DCC p. L1200204 (2012), https://dcc.ligo.org/cgibin/DocDB/ShowDocument?docid=93146.
 Prestegard et al. (2012) T. Prestegard, E. Thrane, et al., Classical Quantum Gravity 29, 095018 (2012).
 Abbott et al. (2004) B. Abbott et al., Phys. Rev. D 69, 122001 (2004).