Hybrid geometric-random template placement algorithm for gravitational wave searches from compact binary coalescences

TL;DR

This paper introduces a hybrid geometric-random template placement algorithm for gravitational wave searches from compact binary systems, improving efficiency and reducing computational costs compared to traditional stochastic methods.

Contribution

A novel hybrid template placement method that automatically handles curvature effects and is faster than stochastic approaches, enhancing gravitational wave data analysis.

Findings

Template banks are 25% larger with stochastic methods.

Hybrid banks are generated nearly ten times faster.

The method maintains comparable fitting factors to stochastic banks.

Abstract

Astrophysical compact binary systems consisting of neutron stars and blackholes are an important class of gravitational wave (GW) sources for advanced LIGO detectors. Accurate theoretical waveform models from the inspiral, merger and ringdown phases of such systems, are used to filter detector data under the template based matched filtering paradigm. An efficient grid over the parameter space at a fixed minimal match has a direct impact on the overall time taken by these searches. We present a new hybrid geometric-random template placement algorithm for signals described by two masses and one spin magnitude parameters. Such template banks could potentially be used in GW searches from binary neutron stars and neutron star-blackhole systems. The template placement is robust and is able to automatically accommodate curvature and boundary effects with no fine tuning. We also compare these banks against vanilla-stochastic template banks and show that while both are equally efficient in the fitting-factor sense, the bank sizes are $\sim 25 \%$ larger in the stochastic method. Further, we show that the generation of the proposed hybrid banks can be sped-up by nearly an order of magnitude over the stochastic bank. Generic issues related to optimal implementation are discussed in detail. These improvements are expected to directly reduce the computational cost of gravitational wave searches.

Tables

Table 1. Table 1: Parameters used to generate the geometric-random and vanilla stochastic banks. The results for different sizes of 𝒰 𝒰 \mathscr{U} are summarized in in Fig. 4 . In Set-I, the parameter space is chosen by satisfying the NS-BH boundary mass i.e. components with individual mass ≤ 3 ​ M ⊙ absent 3 subscript M direct-product \leq 3\mathrm{M}_{\odot} are identified as NS with dimensionless spin magnitude in the range ± 0.4 plus-or-minus 0.4 \pm 0.4 .

Bank Parameter	Set-I	Set-II
Waveform model	$\mathrm{TaylorF2RedSpin}$	$\mathrm{TaylorF2RedSpin}$
Noise model	$\mathrm{aLIGOZeroDetHighPower}$	$\mathrm{aLIGOZeroDetHighPower}$
Lower cut-off frequency $f_{low}$	$20Hz$	$30Hz$
Higher cut-off frequency $f_{high}$	$2048Hz$	$1024Hz$
Mass of first object $m_1$	$[1,20]M_{\odot }$	$[3,15]M_{\odot }$
Mass of second object $m_2$	$[1,3]M_{\odot }$	$[1,3]M_{\odot }$
Spin of first object ${\chi }_1$	$[-0.98,0.98]$	$[-0.6,0.6]$
Spin of second object ${\chi }_2$	$[-0.4,0.4]$	$[-0.05,0.05]$
Size of $𝒰$	$\left[1-8\right]\times {10}^7$	$1\times {10}^7$
Minimal Match ${ℳ}_{\min }$	0.97	0.97

Table 2. Table 2: Summary of various template banks mentioned in this paper. The semi-analytic metric g i ​ j subscript 𝑔 𝑖 𝑗 g_{ij} for TaylorF2RedSpin TaylorF2RedSpin \mathrm{TaylorF2RedSpin} signal model was used in all cases. The Vanilla Stochastic algorithm can also work by directly calculating matches (instead of using g i ​ j subscript 𝑔 𝑖 𝑗 g_{ij} ) but in this head-to-head comparison, we have used the semi-analytic metric which speeds up the Vanilla Stochastic template placement significantly. The usage of the metric is compulsory for Geometric-Random placement.

Bank Parameters	Placement Algorithm	Size of $𝒰$	Bank Size	Execution Time (min)	Comments
Set-I of Table 1		$1\times {10}^7$	$694,422$	$375$	$25\%$ more templates $\times \left(8-10\right)$ faster
		$2\times {10}^7$	$749,705$	$482$
	Hybrid construction	$3\times {10}^7$	$777,113$	$616$
	Geometric-Random	$4\times {10}^7$	$798,269$	$885$
		$5\times {10}^7$	$812,570$	$990$
		$6\times {10}^7$	$824,541$	$1191$
		$8\times {10}^7$	$843,177$	$1712$
	Vanilla-Stochastic (lalapps_cbc_sbank)	$--$	$939,787$	$3666$
Set-II of Table 1	Hybrid construction Geometric-Random	$1\times {10}^7$	$107,547$	$69$	$25\%$ more templates $\times 11$ faster
	Vanilla-Stochastic (lalapps_cbc_sbank)	$--$	$134,563$	$762$
	$A_3^{\ast }$ seeded Stochastic	$--$	$128,185$	$--$
	Vanilla-Stochastic                   (top_down)	$--$	134, 426	$320$	$\times 2.3$ faster

Table 3. Table 3: Geometrical properties of truncated octahedron and pyramid where the truncated octahedron constructed using two truncated pyramid.

Name	Value
Side length of Pyramid $\left(PQ\right)$	$3a$
Side length of TO $\left(BN\right)$	$a$
Height of Pyramid $\left(OP\right)$	$\frac{3}{\sqrt{2}}a$
Height of TO along square face $\left(OS\right)$	$\sqrt{\frac{3}{2}}a$
Height of TO along hexagonal face $\left(OH\right)$	$\sqrt{2}a$
Height of TO along vertices $\left(ON\right)$	$\sqrt{\frac{5}{2}}a$
Volume of TO $\left(V\right)$	$8\sqrt{2}{a}^3$
Surface Area of TO $\left(S\right)$	$\left(6+12\sqrt{3}\right)a^2$
Volumetric Quotient of TO $\left(Q_v\right)$	$\approx 0.683292042$

Table 4. Table 4: These coordinates are the neighborhood positions of a TO where a TO is inscribed in a sphere with radius R 𝑅 R , where the center of a TO is placed at origin and the z 𝑧 z -axis goes through the centre of one of the square faces.

Neighborhood Type	Position
S-Neighborhood	$(\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}}),(-\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}}),$
square faces	$(-\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}}),(\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}})$,
	$(-\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}}),(\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}}),$
	$(-\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}}),(\frac{2R}{\sqrt{5}},-\frac{2R}{\sqrt{5}},\frac{2R}{\sqrt{5}})$
H-Neighborhood	$(0,\frac{4R}{\sqrt{5}},0),(0,-\frac{4R}{\sqrt{5}},0),(\frac{4R}{\sqrt{5}},0,0),$
hexagonal faces	$(-\frac{4R}{\sqrt{5}},0,0),(0,0,\frac{4R}{\sqrt{5}}),(0,0,-\frac{4R}{\sqrt{5}})$

Equations

$$\Lambda\big{(}\vec{\lambda}\big{)}=\mathrm{exp}\Big{\{}\big{<}s\big{|}h(\vec{\lambda})\big{>}-\frac{1}{2}\big{<}h(\vec{\lambda})\big{|}h(\vec{\lambda})\big{>}\Big{\}}.$$

$$\big{<}s\big{|}h(\vec{\lambda})\big{>}_{\Delta t}=4\;\mathrm{Re}\int_{f_{\mathrm{low}}}^{f_{\mathrm{high}}}\frac{\tilde{s}^{*}(f)\;\tilde{h}(f;{\vec{\lambda}})}{S_{n}(f)}e^{-2\pi if\Delta t}\ df,$$

$$\rho(\Delta t;{\vec{\lambda}})=\frac{\big{<}s\big{|}h(\vec{\lambda})\big{>}_{\Delta t}}{\sqrt{\big{<}h(\vec{\lambda})\big{|}h(\vec{\lambda})\big{>}_{\Delta t=0}}}.$$

$$\big{<}h(\vec{\lambda})\big{|}h(\vec{\lambda}+\Delta\vec{\lambda})\big{>}={\mathcal{M}}_{\mathrm{min}}.$$

$$1-g_{\mu \nu }\Delta {\lambda }^{\mu }\Delta {\lambda }^{\nu }={\mathcal{M}}_{\min },$$

$$g_{\mu\nu}=-\frac{1}{2}\frac{\partial^{2}\big{<}h(\vec{\lambda})\big{|}h(\vec{\lambda}+\Delta\vec{\lambda})\big{>}}{\partial\Delta\lambda^{\mu}\;\partial\Delta\lambda^{\nu}}\Big{|}_{\Delta\vec{\lambda}=0}$$

$$\begin{split}\theta_{0}&=\frac{5}{2^{\frac{1}{3}}}\bigg{(}\frac{1}{16\pi f_{0}m\eta^{3/5}}\bigg{)}^{\frac{5}{3}}\\ \theta_{3}&=\bigg{(}\frac{16\pi^{5}}{25}\frac{\theta_{0}^{2}}{\eta^{3}}\bigg{)}^{\frac{1}{5}}\\ \theta_{3s}&=\frac{113\chi_{r}\theta_{3}}{48\pi}\end{split}$$

$$A_i=\sqrt{\frac{1-{\mathcal{M}}_{\min }}{e_i}},$$

$${\xi }_i=\underset{j,k=1}{\sum ^3}{\mathcal{S}}_{ij}{\mathcal{R}}_{jk}^T{\theta }_k$$

$${\bar{g}}_{ij}=\underset{p,q,r,t=1}{\sum ^3}{\mathcal{S}}_{ip}^{-1}{\mathcal{R}}_{pq}^T{g}_{qr}{\mathcal{R}}_{rt}{\mathcal{S}}_{tj}^{-1}.$$

$${\bar{N}}_i^p=\sum _{j,k=1}^3{\mathcal{R}}_{ij}^T{\mathcal{S}}_{jk}{N}_k^p$$

$$\mathcal{FF}(h_{a})=\underset{h\in\mathscr{T}}{\mathrm{max}}\text{\hskip 2.84526pt}\big{<}h_{a}\big{|}h\big{>}$$

$$\widetilde{h}(f)=\mathcal{A}f^{-7/6}\exp\Big{\{}-i\big{[}\Psi(f)-\frac{\pi}{4}\big{]}\Big{\}}$$

$$\begin{split}\Psi(f)=2\pi ft_{0}+\phi_{0}+\frac{3}{128\eta v_{f}^{5}}\Bigg{\{}1+v_{f}^{2}\Big{[}\frac{55\eta}{9}+\frac{3715}{756}\Big{]}+v_{f}^{3}\big{[}4\beta-16\pi\big{]}+v_{f}^{4}\Big{[}\frac{3085\eta^{2}}{72}+\frac{27145\eta}{504}\\ +\frac{15293365}{508032}-10\sigma_{0}\Big{]}+v_{f}^{5}\Big{[}\frac{38645\pi}{756}-\frac{65\pi}{\eta}-\gamma_{0}\Big{]}\big{(}3ln(v_{f})+1\big{)}+v_{f}^{6}\Big{[}-\frac{6848\gamma_{E}}{21}-\frac{127825\eta^{3}}{1296}+\\ \frac{76055\eta^{2}}{1728}+\big{(}\frac{2255\pi^{2}}{12}-\frac{15737765635}{3048192}\big{)}\eta-\frac{640\pi^{2}}{3}+\frac{11583231236531}{4694215680}-\frac{6848ln(4v_{f})}{21}\Big{]}+\\ v_{f}^{7}\Big{[}-\frac{74045\pi\eta^{2}}{756}+\frac{37815\pi\eta}{1512}+\frac{77096675\pi}{254016}\Big{]}\Bigg{\}}\end{split}$$

$$\begin{split}\beta&=\frac{113\chi_{r}}{12}\\ \sigma_{0}&=\Big{(}-\frac{12769(4\eta-81)}{16(76\eta-113)^{2}}\Big{)}\chi_{r}^{2}\\ \gamma_{0}&=\Big{(}\frac{565(1713\eta^{2}+135856\eta-146597)}{2268(76\eta-13)}\Big{)}\chi_{r}\end{split}$$

$$\chi_{r}=\frac{1}{2}\Big{(}1-\frac{76}{113\eta}\Big{)}(\chi_{1}+\chi_{2})+\frac{1}{2}\frac{m_{1}-m_{2}}{m_{1}+m_{2}}(\chi_{1}-\chi_{2})$$

$${\chi }_{1,2}=\frac{{\mathbf{S}}_{1,2}\cdot {\hat{\mathbf{L}}}_N}{2m_{1,2}^2}.$$

$$Q_v=\frac{V}{\frac{4}{3}\pi R^3},$$

Full text

I Introduction

The direct detection of gravitational waves (GW) was vigorously pursued by the LIGO Scientific Collaboration (LSC) consisting of several hundred scientists which culminated with the discovery of the first gravitational wave event GW150914 Abbott et al. (2016a) in the twin advanced-LIGO (advLIGO) Harry and the LIGO Scientific Collaboration (2010) detectors. This event was determined to have been caused by the inspiral and merger of a spinning binary blackhole system of component masses $(36,\;29)\mathrm{M}_{\odot}$ located nearly 1.3 billion light years away from Earth. It was also the first observational evidence for the existence of stellar mass blackholes. Subsequently, advLIGO detectors also detected a second event GW151226 Abbott et al. (2016b) from the inspiral and merger of lighter compact objects. These detections mark the transition to the era of gravitational wave astronomy.

Several other kilometer scale detectors are under upgradation or construction at present around the globe - these include the French-Italian advanced Virgo detector Acernese et al. (2015) and the Japanese cryogenic KAGRA detector Somiya (2012); Aso et al. (2013). In India, the proposal for the advanced LIGO-India detector Iyer et al. (2011); Unnikrishnan (2013) has been approved and is expected to be built over the next few years. The network of these advanced detectors is expected to improve their overall science potential and herald a new wave in astronomy with the potential to observe the very early Universe and compliment information gathered by electromagnetic observations.

Binary blackholes and neutron stars are considered one of the most promising sources for the advanced terrestrial detectors. Precise theoretical model waveforms for GW emitted from the inspiral, merger and ringdown phases of such compact binary coalescences (CBC) are now available. These models are parametrized by the system’s intrinsic properties such as the component masses, spin etc. and extrinsic properties such as sky location, distance to the source, time at coalescence etc. Accurate theoretical models allow the use of matched filtering technique to search for weak GW signals buried in detector noise. As the signal parameters are not known a priori, one filters the data using a set of expected signals spanning the deemed parameter space. Each one of these expected signal corresponds to a single point in the parameter space, and are collectively known as the template bank. Coverage of the full range of search parameters using a finite grid of discrete points leads to an inevitable loss in the signal-to-noise ratio (SNR) which can be controlled by fixing the minimal match of the bank. The latter is often decided by striking a balance between desired detection efficiency and computational cost of carrying out the search.

Advanced detectors are reaching unprecedented sensitivities at low frequencies. Over the last few years, the development of theoretical spinning waveform models have also reached a mature stage. The combined effect of these factors is that one now needs to search over a significantly larger volume compared to initial LIGO era, in a parameter space spanning three or more dimension - posing new challenges for data analysis, which include devising efficient grid placement strategies. The practical issue of optimal template placement in matched-filter searches is an important open optimization problem Manca and Vallisneri (2010); Messenger et al. (2009).

There are at least two approaches to solving the problem of efficient template placement - (a) via deterministically placed points through the tessalations of a geometrical lattice Cokelaer (2007); Babak et al. (2006) or (b) via stochastically placed points by choosing them at random Harry et al. (2009) over the parameter space. The geometric method requires the metric over the signal manifold and has been used extensively for CBC searches over 2D mass parameters by arranging the grid points in a hexagonal lattice. A variant of the geometrical placement method for aligned spin CBC systems has been explored Brown et al. (2012) for the TaylorF2 Buonanno et al. (2009); Poisson and Will (1995) signal model. Lack of availability of metric and also the intricate fine tuning required to avoid uncovered regions arising from variations in curvature across the parameter space also make it difficult to generalize the geometric placement methods to higher dimension.

The stochastic template bank is constructed from random proposals drawn from a uniform distribution over the deemed parameter space that are accepted as a new template point only if the new proposal is far (in minimal match sense) from existing templates in the bank. Such banks are easy to implement, robust and can work even without the explicit knowledge of the metric on the signal manifold by using brute force match calculations. This approach has been demonstrated to be more space efficient than a square lattice in 2D but less efficient than a hexagonal lattice. Such banks can be calculated for higher dimensions as well Ajith et al. (2014).

An earlier attempt to combine the geometric and stochastic methods Capano et al. (2016) by seeding a stochastic bank with a pre-fabricated geometric bank has also been tried and has been demonstrated to improve the efficiency marginally.

We present a new hybrid algorithm for template placement in 3D (two mass components and one reduced spin magnitude) parameter space for gravitational wave searches from compact binary coalescence by combining the efficiency of optimal geometrical placement and the robustness and ease of stochastic placement algorithm. This geometric-random bank placement method uses a local truncated octahedral lattice to place the templates and requires the metric over the signal manifold.

Plan of the Paper:

The paper is organized as follows:

In Section II we recapitulate the fundamentals of templated matched filtering technique used for gravitational wave searches from CBC and review the definition of the metric on the signal manifold with a view to set the notation used in the paper and also define key terms. Numerous papers on template placement strategies have been written over the last few years and many new ideas have emerged Brown et al. (2012); Capano et al. (2016); Harry et al. (2009); Ajith et al. (2014): we provide a concise review of these efforts for the benefit of the reader. In particular, we elucidate the stochastic algorithm by casting it in two different ways which are algorithmically equivalent. This sets the stage for Section III, where we present a new hybrid geometric-random template placement algorithm and give a detailed explanation of the issues concerning its optimal implementation. We argue that by construction this new method can not be less space efficient than the vanilla-stochastic placement method. We also construct an explicit template bank for neutron star-blackhole searches using this new algorithm and compare it against a vanilla-stochastic template bank, noting the improvement in overall bank size as well as the time taken to generate the bank. The $\mathrm{TaylorF2RedSpin}$ signal model Ajith (2011) was used for this purpose. The template banks presented in this work have been calculated using the $\mathrm{aLIGOZeroDetHighPower}$ Shoemaker (2010) noise sensitivity curve for advLIGO.

In Section IV, the template banks generated by the hybrid method are tested and validated against the vanilla stochastic bank. We present the fitting factor results using $50,000$ signal injections using $\mathrm{TaylorF2RedSpin}$ and TaylorF2 aligned spin waveform models using software implemented in the LALApps package of the LIGO Algorithm Library LIGO Data Analysis Software Working Group and show that the two banks are nearly identical in efficiency. We also compare the hybrid bank against a $A_{3}^{*}$ lattice seeded vanilla-stochastic bank and show that the former is more efficient.

Finally, in Section V, we summarize the main results and make some comments related to several key issues related to the new method, indicating a possible way for extending it to higher dimensions.

II Metric on the signal manifold

In this section we shall quickly summarize the basics of template based matched filtering technique used in GW searches from CBC. A basic assumption in matched filtering based searches is that the astrophysical GW signal buried in the detector noise is faithfully represented by the signal model used in the search over the range of search parameters.

The signal manifold ${\bm{\lambda}}$ is the set of all possible GW signals $h({\vec{\lambda}})\equiv h(t;{\vec{\lambda}})$ characterized by the parameter vector $\vec{\lambda}$. It is customary to represent the corresponding frequency domain signal as $\tilde{h}(\vec{\lambda})\equiv\tilde{h}(f;\vec{\lambda})$. The detector output $s(t)$ consists of detector noise $n(t)$ and a possible gravitational wave signal of unknown parameters. The additive noise throws the signal out from this manifold to the space of all possible functions. In order to find the point in signal manifold closest to the detector output $s(t)$, the latter is projected over the signal manifold by calculating the maximum likelihood over ${\vec{\lambda}}$ which serves as the detection statistic. For additive Gaussian noise, the likelihood $\Lambda$ is given by Finn (1992)

[TABLE]

The inner product \big{<}s\big{|}h\big{>} is the complex cross-correlation Sathyaprakash and Dhurandhar (1991) between the detector output $s(t)$ and the gravitational wave signal $h(\vec{\lambda})$ weighed inversely by the noise power spectral density of the detector Babak et al. (2013):

[TABLE]

where $S_{n}(f)$ is the one-sided noise power spectral density defined by \big{<}\tilde{n}(f)\big{|}\tilde{n}^{*}(f^{{}^{\prime}})\big{>}=\frac{1}{2}S_{n}(|f|)\;\delta(f-f^{{}^{\prime}}), asterisk ∗ denotes the complex conjugation operator, the frequency range $f_{\mathrm{low}}\leq f\leq f_{\mathrm{high}}$ marks the effective bandwidth of the detector and $\Delta t$ is the time delay between these two signals. The signal-to-noise ratio (SNR) $\rho$ after filtering $s(t)$ is defined as:

[TABLE]

Without any loss of generality, we assume that all template waveforms are normalized such that \big{<}h(\vec{\lambda})\big{|}h(\vec{\lambda})\big{>}_{\Delta t=0}=1. From Eq. (1) it is clear that this allows us to use the log-likelihood function (or equivalently, the SNR) maximized over the parameters, as the detection statistic.

The log-likelihood function can be maximized over all time lags ($\Delta t$) by using Fast Fourier Transform (FFT) based convolution as shown in Eq. (2). It can be maximized over other extrinsic parameters analytically Sathyaprakash and Dhurandhar (1991). On the other hand, as one cannot maximize the log-likelihood function over intrinsic parameters analytically, a brute force approach is needed. In this case, a discrete grid of points is placed to cover the intrinsic parameter space. One evaluates the log-likelihood surface at each of these points and the maximum is suitably determined.

The template bank consisting of these discrete set of points on ${\bm{\lambda}}$ is constructed using a control parameter ${\mathcal{M}}_{\mathrm{min}}$; commonly known as minimal match in GW literature. This parameter is chosen such that the minimum overlap of an arbitrary vector in the signal manifold (within the deemed parameter space) and at least one template in the bank never drops below this value. The art of template placement lies in maximizing the inter-template separation without violating this constraint with the aim of achieving the smallest bank size. In this way, we can map the template placement to the sphere covering problem Conway and Sloane (1999); Prix (2007) with spherical cells of radius equal to $\sqrt{1-{\mathcal{M}}_{\mathrm{min}}}$. One wants the smallest number of overlapping templates (i.e. metric spheres) to fully cover the space (i.e. leave no holes). The equation of such a spherical cell, centered at $\vec{\lambda}\in\bm{\lambda}$, is given by

[TABLE]

From the normalization of waveforms, it is clear that ${\mathcal{M}}_{\mathrm{min}}\leq 1$. For high values of this parameter, the LHS can be Taylor series expanded upto leading order in the small parameter $\Delta\vec{\lambda}$:

[TABLE]

where $\Delta\lambda^{\mu}$ are the components of the vector $\Delta\vec{\lambda}$ and

[TABLE]

is the metric over the signal manifold which is essentially the Fischer information matrix projected on the intrinsic parameter space and calculated using standard covariance matrix method Poisson and Will (1995); Owen (1996); Owen and Sathyaprakash (1999). The match between two nearby points in the parameter space can be calculated easily using the metric. Eq. (5) can be re-arranged as $g_{\mu\nu}\Delta\lambda^{\mu}\Delta\lambda^{\nu}=(1-{\mathcal{M}}_{\mathrm{min}})$, and identified to be the equation of an ellipsoid in 3D centered at a point $\vec{\lambda}$. For higher dimensions it represents a hyper-ellipsoid. We shall refer to this as the minimal-match ellipsoid elsewhere in the paper. In GW searches from CBC sources, template banks are usually constructed at ${\mathcal{M}}_{\mathrm{min}}\sim 0.97$ Babak et al. (2013). This corresponds to a loss in detection rate of $\sim 9\%$ assuming uniform distribution of such sources.

II.1 State of the art in template placement

As mentioned earlier, the geometric and stochastic template placement algorithms are two broad class of methods used in searches for gravitational wave signals from compact binary coalescences.

Previous searches for GW signals from nonspinning compact binaries in initial-LIGO and Virgo Abbott et al. (2009a, b); Abadie et al. (2010) data have used the metric based geometrical hexagonal template placement in two dimensions Cokelaer (2007); Babak et al. (2006); Owen (1996); Owen and Sathyaprakash (1999); Sathyaprakash and Dhurandhar (1991). Templates are placed in chirp time coordinates $\{\tau_{0},\tau_{3}\}$ instead of component masses $\{m_{1},m_{2}\}$, since the templates are almost uniformly spaced in the former. This process starts by initializing a template point and then finding neighboring points in a $A_{2}^{*}$ hexagonal lattice, where the centre of the hexagons represent the position of individual templates. Hexagonal tiling offers the most efficient space filling in two dimensions, which optimizes the number of template points and in turn reduces the total computational cost of the search. To construct this geometric bank one require the semi-analytic metric on the signal manifold assumed to be slowly varying over the parameter space. The curvature effects leads to some loss of efficiency in this strategy. This geometric placement was initially demonstrated Cokelaer (2007) for 2PN SPA family of waveforms. At present, the metric for 3.5PN SPA waveforms Keppel et al. (2013) is available which allows template placement for such waveforms as well.

Geometric template placement in higher dimensional intrinsic parameter space (e.g. component masses and spins) has several problems. First, the metric may not be available for such signal manifolds and secondly, (unlike the hexagonal packing in 2D), optimal geometrical placement in higher dimensions with curvature are not well known. It is further complicated from curvature and boundary effects. To mitigate the curvature issues which leads to rapidly changing metric components over the parameter space, new coordinates are being explored in which the metric is slowly varying Ajith et al. (2014), but they do not solve the problem completely. Recent studies Brown et al. (2012); Harry et al. (2014) have explored geometrical placement in higher dimensions for aligned-spin BNS and NSBH systems for some specific waveform families. In this method, one constructs a metric on the parameterized coordinates - taken to be the coefficients of the 3.5PN TaylorF2 expansion of the orbital phase Pai and Arun (2013), instead of the usual chirp time coordinates. Since the metric is globally flat in these new coordinates, one can globally transform it into a Euclidean coordinate system. Finally, a principal coordinate analysis facilitates the projection to an effective lower dimensional parameter space, which can be covered by a grid placed in a hexagonal lattice. For NSBH Harry et al. (2014) systems, the parameter space is covered by stacking several two dimensional hexagonal lattices along the direction of the minor axis. This method is available in the PyCBC software package Nitz et al. (2016); Usman et al. (2016); Dal Canton et al. (2014); Allen et al. (2012).

An alternative approach of template placement is the so-called stochastic method Harry et al. (2009) where one starts with an empty template bank (or a set of seed points) to which random points, drawn from a uniform distribution over the deemed parameter space, are appended in an iterative fashion. At every iterative step, a new random point is proposed to be included in the template bank: this proposal is rejected if it happens to lie too close to the points that are already in the existing list otherwise it is accepted. For each accepted proposal, the rejection rate is determined and the process terminates when this rate exceeds a certain threshold averaged over the last few acceptances. We shall call this the bottom-up approach in building the list of templates by considering random proposals one by one, retaining only the valid ones. This has been encoded in the LSC Algorithm Library (LAL) software suite via the program lalapps_cbc_sbank.

The stochastic bank placement algorithm can also be cast in a top-down fashion. In this alternate implementation, one starts with an empty template list $\mathscr{T}$ and a list of very large number of proposals $\mathscr{U}$ distributed uniformly over the deemed parameter space. One picks a random point $r\in\mathscr{U}$ and appends it to $\mathscr{T}$, following which all points from $\mathscr{U}$ that lie within the minimal match distance from $r$ are removed. The process continues in an iterative manner until all points from $\mathscr{U}$ are exhausted. We call this the top-down approach as the template-bank is microfabricated out of a large block of random proposals by paring it down to the desired shape and order.

Both these methods are algorithmically equivalent - however, the top-down approach is more useful in projecting a geometric structure over the stochastic template bank leading to the hybrid geometric-random placement algorithm described in this paper. In Table 2 we demonstrate that the top-down approach is also faster by a factor $\gtrsim 2$ over the traditional bottom-up approach as one is able to eliminate many proposals for a single accepted proposal using efficient computational data structures such as binary search trees (BST Hibbard (1962)). The number of templates generated by both these methods is nearly identical.

The stochastic method is relatively easier to implement and does not require the metric over the signal manifold per se. The distance between the proposed point and the elements of the current template bank can be directly calculated by evaluating the match inner product Eq. (2). This brute-force approach for match calculation allows it to be extended easily to higher dimensional parameter spaces and overcome irregular boundary effects. The disadvantages of stochastic method include the requirement of high computational time as several million proposals have to be processed to guarantee adequate coverage and the fact that it generates substantially more templates than the geometric bank. The computational time can be reduced, if we use the metric (if available) for match calculation. But the intrinsic stochastic nature of the algorithm leads to the inefficiency in grid placement.

Another instance of template placement developed for aligned spin binary black hole (BBH) searches has explored a combination of geometric and stochastic approaches Capano et al. (2016); Privitera et al. (2014). In this method, at first one generates an aligned-spin geometric hexagonal lattice template bank upto some valid range of parameters which is then used as a “seed” bank for the stochastic placement - thereby accelerating the placement. This method generates $\sim 5.5\%$ fewer template points than the stochastic method and has been used as part of the uber template bank used in CBC searches in the data from the first observational run of aLIGO Abbott et al. (2016c).

III A new geometric random algorithm for template placement

In this section we present the metric based hybrid geometric-random template placement algorithm in three dimensions using a truncated octahedral lattice. Such lattices are the Dirichlet-Voronoi polytope of body-centred cubic $A_{3}^{*}$ lattice Schürmann and Vallentin (2006). The latter provide optimal coverage for conformally flat spaces where the metric-coefficients are constant Prix (2007). It is interesting to note that this is in line with Lord Kelvin’s conjecture Thomson (1887) according to which, truncated octahedron based space filling is optimal in flat 3D space. A truncated octahedron is a 14-sided space filling polyhedron and has the highest volumetric quotient (Eq. 18) which makes it suitable for the lattice structure. Other geometric properties are tabulated in Appendix B.

Such a template bank would be applicable for gravitational wave searches from the compact binaries described adequately by their component masses and a single “reduced spin” parameter $\chi_{r}$ using $\mathrm{TaylorF2RedSpin}$ signal model. This signal model is constructed using post-Newtonian (PN) template family of gravitational waveforms from inspiralling compact binaries with non-precessing spins. Here, the spin effects are captured by the parameter $\chi_{r}$ defined as a mass-weighted, linear combination of individual dimensionless spin magnitudes. The details of the signal model can be looked up in Appendix A.

In order to construct the three-dimensional geometric random template bank, we require the metric on the parameter space. The metric $g_{ij}$ for $\mathrm{TaylorF2RedSpin}$ approximant varies rapidly in $\{m,\eta,\chi_{r}\}$ coordinate system where $m$ denotes the total mass and $\eta$ denotes the symmetric mass ratio. In order to enhance the efficiency of the algorithm, one places the templates in dimensionless chirp time coordinates $\{\theta_{0},\theta_{3},\theta_{3s}\}$, in which the metric components are slowly varying over the parameters. These new coordinates are defined as:

[TABLE]

However it is to be noted that the curvature effects do not vanish completely due to the above coordinate transformations. As such, Kelvin’s conjecture does not hold directly in this non-flat space. In the following sections, we show that the truncated-octahedral design can still be used for spaces with slowly varying curvature. This is achieved by merging the stochastic bank placement algorithm along with a local $A_{3}^{*}$ lattice.

The method outlined in Algorithm (1) proceeds by initializing three lists as follows:

1. (a)

$\mathscr{U}$: a list of uniformly distributed random points sprayed over the deemed parameters space 2. (b)

$\mathscr{T}$: an empty list for template points 3. (c)

$\mathscr{K}$: an empty temporary list

At first $\mathscr{K}$ is initialized or seeded with a random point chosen from $\mathscr{U}$ which is immediately appended to the list of templates $\mathscr{T}$ as its first element. The possible 14 TO neighbours of this initial point are then calculated and appended to $\mathscr{K}$ followed by the removal of the first element. Only those TO neighbors which fall within the parameter space and are farther than $\mathcal{M}_{min}$ from existing points in the $\mathscr{T}$ and $\mathscr{K}$ can be appended. This completes one iteration of the inner loop of the algorithm. We continue these steps until $\mathscr{K}$ becomes empty. The latter happens when no new valid TO neighbors are appended in successive iterations of the inner loop of the algorithm which ultimately exhausts all the points in $\mathscr{K}$. At this stage, all points from $\mathscr{U}$ that are within the minimal match distance from the elements in $\mathscr{T}$ are deleted; following which $\mathscr{K}$ is re-seeded with a new random point from $\mathscr{U}$. Termination of template placement algorithm occurs when $\mathscr{U}$ is exhausted.

The advantage of casting the stochastic template placement algorithm (see Sec II.1) in a top-down fashion is clear when we consider the extreme case where no TO lattice neighbours can be found for any proposal point in $\mathscr{U}$. In this case, it is clear from Algorithm 1 that the geometric-random algorithm naturally falls back to the vanilla-stochastic placement algorithm. This brings in the robustness of the latter to the proposed new algorithm. This contingency arises for a small fraction ($\sim 5-10\%$) of points in $\mathscr{U}$ - due to boundary edge effects and small uncovered patches arising from curvature effects. By construction the geometric-random placement method presented here is more efficient than vanilla-stochastic algorithm.

Also note that the TO, whose neighbors are added to $\mathscr{K}$ in Lines 8-9 of Algorithm (1) refers to the one inscribed inside the minimal match ellipsoid given in Eq. (5) above. The gradual change in curvature (changing metric components) leads to changing orientation and size of these ellipsoids which affects the size and orientation of the inscribed truncated-octahedron. This in turn, affects the coordinates of its 14 neighbors. In this way the placement algorithm automatically responds to the curvature effects. This is equivalent to assuming flat local patches of the signal manifold which can be optimally covered by $A_{3}^{*}$ lattice.

III.1 Implementation

We now highlight some salient features for efficient implementation of the algorithm.

III.1.1 Initializing $\mathscr{U}$ with uniform random points

In order to generate a random list of points over the parameter space, at first we calculate the minimum and maximum possible values of dimensionless chirp times $\theta_{0},\theta_{3}$ and $\theta_{3s}$ for the range of parameters over which the bank is to be placed. A set of random points are generated in $\{m_{1}-m_{2}\}$ space along the boundaries marking the constraints for component masses, chirp masses and mass ratios. Afterwards, a coordinate transformation is taken to the dimensionless chirp time coordinate system $\{\theta_{0}-\theta_{3}\}$ which allows the calculation of the extremities in these coordinates. The extreme values of $\theta_{3s}$ are evaluated using $|\chi_{r}|\leq 1$ in Eq. (7). Once these extremities are known, uniform random points are proposed within these values - out of which only those which lie within the specified mass and spin ranges are retained.

The final banksize depends strongly on the initial density of points in $\mathscr{U}$ upto a critical value, beyond which it does not change very much. Empirically, we find this number to be about $\sim 1/10$ of the total volume calculated in the dimensionless coordinates. The latter depends on the metric and varies across the parameter space due to which we estimate an average cell volume by sampling different points across this space. The total number of cells needed to cover the space is given by the ratio of the total volume divided by the average volume of a cell. Monte Carlo integration method can be used to estimate the total volume.

III.1.2 Deletion of points within minimal-match ellipsoid

We discuss an efficient way to delete points from $\mathscr{U}$ that lie within minimal match distance from the template bank points. For a fixed minimal match value, Eq. (5) represents the surface of an ellipsoid with semi axes length:

[TABLE]

where $e_{i}$’s are eigenvalues of the metric $g_{ij}$ evaluated at the ellipsoid centre. The orientation of the ellipsoid is defined by matrix $\mathcal{R}$ composed of eigenvectors of the metric $g_{ij}$, where $\mathcal{R}_{ij}$ is $j^{th}$ component of the $i^{th}$ eigenvector. To visualize the deletion process, let us consider a cell centred around a template. The points in $\mathscr{U}$ lying inside the cell can be determined by calculating the match between this template with the all other points in $\mathscr{U}$ but such a brute force approach will not be computationally efficient. We need to find a smaller set of points against which we can identify these neighbors more efficiently. To this end, we can first use binary search tree (BST) Hibbard (1962); Douglas (1959); Windley (1960) data structure to identify a smaller list of points from $\mathscr{U}$ that lie within a sphere of radius equal to the largest semi axes of the cell. From this smaller list of points, we then proceed to delete the ones that satisfy $g_{ij}\Delta\lambda^{i}\Delta\lambda^{j}\leq 1-{\mathcal{M}}_{\mathrm{min}}$.

This can be further refined by binning the template points in $\theta_{0}$ such that each bin contains about $\sim 1000$ templates. For each such subset, the removal of points from $\mathscr{U}$ proceeds by identifying those points that lie in the same bin (within some acceptable margin in $\theta_{0}$) and then applying the above strategy to eliminate the points. This refinement is made possible due to the fact that match values decrease monotonically with the difference in parameters and is most sensitive to changes in $\theta_{0}$. Additional binning in other two coordinates may possibly improve the computational efficiency of removing the points from $\mathscr{U}$ even more.

III.1.3 Global co-ordinate transformation

In dimensionless chirp time coordinates, most of the ellipsoidal cells have semi-axes ratio around $1000:10:1$. This implies that the number of points inside this flat and elongated ellipsoid is a small fraction of the total number of points contained in the sphere with radius equal to semi-major axis. This undermines the efficiency of deleting the points using the technique outlined in Section III.1.2. We use a conformal coordinate transformation in such a way that one of the cells transforms to a unit sphere. Note that due to curvature effects, the same transformation when applied to other cells does not guarantee them to change to unit spheres but mitigates the problem of highly asymmetric semi-axes ratio to a large extent. We use the eigenvalues and eigenvectors of the metric at a fiducial point to transform the coordinates ($\theta_{k}\rightarrow\xi_{k}$) as,

[TABLE]

where, $\mathcal{S}_{ij}=\sqrt{e_{i}}\delta_{ij}$ is a diagonal scaling matrix and the rotational matrix $\mathcal{R}_{ij}$ is constructed from the $i^{\text{th}}$ component of the $j^{\text{th}}$ eigenvector of the metric calculated at the fiducial point. This is carried out by first binning the templates in $\theta_{0}$. The transformation in Eq. (9) is carried out using the metric of the centre point in each bin.

The metric at all other points transform as:

[TABLE]

It is trivial to check that the metric at the fiducial points transforms to a $3\times 3$ identity matrix.

This global coordinate transformation Eq. (9) causes other nearby cells to become almost spherical with semi-axes in the ratio $\sim\ 1:1:1$. In Fig. 1 we show the effect of global coordinate transformation on nearby cells. This is doubly advantageous: not only does the ball-point query volume for the BST searches decrease (leading to fewer points to deal with); the points that lie within the inscribed minimal match ellipsoid now occupy a large fraction of its volume leading to increased efficiency in removal of points.

III.1.4 Locally placed truncated octahedral lattice

The template placement problem in 3D can be mapped to the tiling problem such that a minimum number of similar cells is used without leaving any region uncovered. As such, truncated octahedral lattice becomes a natural choice for this problem.

As stated earlier, the TO shares 14 faces with its neighbours in the lattice. Assuming a TO inscribed in a sphere of radius $R$, the coordinates of its 14 lattice neighbors $N_{k}^{p}$ are available in Table 4. Here $k=(1,2,3)$ is an index on the coordinates of the $p-th$ neighbor. When mapped to the template placement problem, we need to consider TO’s inscribed within the minimal match ellipsoids. In this case, the coordinates of the lattice neighbors $\bar{N}_{i}^{p}$ can be calculated using the rotation and scaling matrices:

[TABLE]

Fig. 2 shows the 14 neighbours of a TO inscribed in a elliptical cell.

Due to boundary effects, all 14 TO neighbors of a point need not necessarily be part of the template bank. The following conditions must be checked for:

1. (a)

The point is inside the deemed parameter space and also satisfies $\eta\leq 1/4$. The latter corresponds to the condition that $\theta_{k}$’s can be inverted to yield physical masses. 2. (b)

The point is not within the minimal match distance of existing points in $\mathscr{T}$ and $\mathscr{K}$.

Check (b) above ensures that we do not double count the neighbors.

As shown in Algorithm (1), we start from a random point in the parameter space (by seeding $\mathscr{K}$) and tessellate with local TO lattice. Because of curvature and boundary effects, it is not guaranteed that these tessellations cover the entire parameter space. This is marked by the exhaustion of $\mathscr{K}$ as the placement proceeds. At this stage, we need to re-seed $\mathscr{K}$ and continue the process iteratively until all points in $\mathscr{U}$ have been used up.

III.1.5 Choice of initial point and variations in banksize

The template placement algorithm starts from a randomly chosen point in the parameter space by seeding $\mathscr{K}$ which is copied over as the first element of the template bank list $\mathscr{T}$.

One can start from any point in parameter space: in the present work we have started from “mid-point” corresponding to component masses $m_{1,2}=(m_{1,2}^{max}+m_{1,2}^{min})/2$ with individual spins $\chi_{1,2}=(\chi_{1,2}^{max}+\chi_{1,2}^{min})/2$. We have checked that starting from other points (e.g. the extremities of the parameter space) also work quite well. This initial choice results in a minor fluctuation of the template bank size and it is recommended that we start from a point well inside the deemed parameter space where the local variations in curvature are less. This ensures maximum tiling before we need to re-seed $\mathscr{K}$.

We demonstrate these fluctuations in template bank size by constructing several template banks for compact binary systems whose bank parameters are given in Set II of Table 1 starting from different locations in the deemed parameter space. The random seed used to initialize $\mathscr{U}$ was kept the same to eliminate bias. The five different starting points were taken to be the centre and extremities of the parameters space respectively and as expected, different choices of the initial seed point resulted in slightly different number of templates in the final bank. We generated an average of $111,257$ templates with $\sim 3\%$ fluctuation which is quite insignificant.

IV Test of the Algorithm

In this section we construct a template bank using our algorithm and demonstrate its validity for CBC searches. We also compare its performance against a bank generated using the vanilla stochastic method available in the LSC Algorithm Library Suite (LALSuite) LIGO Data Analysis Software Working Group . We need a metric over the dimensionless chirp coordinates and at present two such models are available for use, namely; $\mathrm{TaylorF2RedSpin}$ and IMRPhenomB signal model Kalaghatgi et al. (2015). Our demonstration makes use of the $\mathrm{TaylorF2RedSpin}$ model.

IV.1 Construction of the template banks for reduced-spin binary system

We generate a template bank using $\mathrm{TaylorF2RedSpin}$ metric in $(\theta_{0},\theta_{3},\theta_{3s})$ parameter space for aligned spin compact binary system. The range of parameters is chosen such that the mass of the first object lies between $1-20\mathrm{M}_{\odot}$ with dimensionless spin magnitude in the range $\pm 0.98$. The mass of the second object is taken between $1-3\mathrm{M}_{\odot}$ with dimensionless spin magnitude in the range $\pm 0.4$. The NS-BH boundary mass condition is satisfied i.e any object with mass $\leq 3\mathrm{M}_{\odot}$ is considered to be a neutron star with corresponding limit on spin magnitude. This template bank can be used for binary neutron star (BNS) and neutron-star-blackhole (NSBH) searches.

We construct several template banks by varying the number of uniformly sprayed random points over the parameter space, $\mathscr{U}$. We consder sizes of $\mathscr{U}$ varying between $1-8\times 10^{7}$ points. In all the cases, the placement proceeds from a point corresponding to individual masses $m_{1,2}=(\frac{21}{2},\frac{4}{2})\mathrm{M}_{\odot}$ and individual spin magnitudes $\chi_{1}=\chi_{2}=0.0$ . The full specification of input arguments for template generation are given in Set-I of Table 1.

The corresponding final template bank size are listed in Table 2. In Fig. 3 we show the two dimensional projections of such a bank along $\hat{\theta}_{3s},\hat{\theta}_{3}$ and $\hat{\theta}_{0}$ directions respectively. The hat $(\hat{)}$ over the symbols denote unit normal vector along the specified direction. For comparison, we also generated a template bank using the vanilla stochastic method for the same parameters as listed in Set-I of Table 1. The vanilla stochastic bank placement was terminated at a point when there were 1000 rejected proposals for each accepted proposal averaged over the last 10 acceptances. This stochastic bank was found to contain $939,787$ templates which is $\sim 25\%$ larger than the geometric-random bank. The computational run time was also recorded and we observed that the geometric-random method took $482$ minutes while the vanilla stochastic method took $3666$ minutes to execute on a single, unloaded processor, which is $\sim 8$ times faster.

IV.2 Validation of the template banks

We investigate the performance of both the geometric-random and the vanilla stochastic template banks against a set of signal injections from the reduced spin $\mathrm{TaylorF2RedSpin}$ signal model. In this section we summarize the results of this comparison and demonstrate that the two banks are nearly identical in performance.

Following Apostolatos Apostolatos (1995), the “fitting factor” $\mathcal{FF}(h_{a})$ is defined as a measure of the maximum match over the template bank $\mathscr{T}$ for a putative injected signal $h_{a}$:

[TABLE]

The parameters of the injected signal (chosen from within the deemed parameter space over which the bank is placed) is chosen at random and may not coincide with that of a template point. The mismatch $(1-\mathcal{FF}(h_{a}))$ indicates the fractional loss of optimal SNR. Because differences between true waveforms and waveform models will always reduce the fitting factor, banks are usually tested to achieve fitting factors slightly larger than the minimal match used in their construction. The minimal match at which the template banks are constructed ($0.97$ in this case) is somewhat arbitrary - and is arrived by carefully balancing the computational cost of the search against the desired detection efficiency.

For the banks generated as per the parameters in Set-I of Table 1, we injected $50,000$ signals from reduced-spin $\mathrm{TaylorF2RedSpin}$ waveform model and calculated the fitting factor using Eq. (12). The lalapps_cbc_sbank_sim program as implemented in LALApps package of the LSC Algorithm Library Suite (LALSuite) LIGO Data Analysis Software Working Group was used for this calculation. The range of mass and spin parameters of the injected signals were chosen to be the same as that of the bank and were drawn from a uniform distribution in this range. Other extrinsic parameter were drawn from (i) a uniform random distribution over all possible sky locations (ii) a fixed inclination angle corresponding to edge-on orientation of the plane of the binary, and (iii) a fixed luminosity distance of 1 Mpc. For the banks (corresponding to different initial sizes of $\mathscr{U}$), we found that $\sim 92-94\%$ signals were recovered with fitting factor $\geq 0.97$ for both geometric-random bank as well as for the stochastic bank. The results are shown in Fig. 4. The fluctuation at higher significant digit was ignored due to the limited sample size of injections. These results are graphically depicted in Fig. 4 and demonstrate the equivalence of the two banks along with the computational speed up and corresponding efficiency of template placement using the new method.

Although the injection signal model is identical to the one used to construct the template banks, we still find that the fitting factors fall below the minimal match of the bank for a small fraction of the injections. This alludes to the fact that both template banks have holes in certain regions of the deemed parameter space where neighbouring templates do not provide adequate coverage. This arises from curvature effects.

The fitting factor depends on the parameters of the injected signals - e.g. masses, spin components, sky position and inclination angle of the binary’s plane to the line of sight etc. For non-precessing signals including only the dominant radiation mode, the fitting factor depends only on the intrinsic parameters. To understand the systematics, it is convenient to represent it as a function of two parameters while averaging out over the remaining ones.

In Fig. 5LABEL:sub@fig_ff_a–LABEL:sub@fig_ff_b, we show the histogram of “minimum” fitting factors Harry et al. (2014) (for geometric-stochastic template bank) over various combinations of intrinsic parameters of the compact binary system where both the signals and template waveforms are generated from $\mathrm{TaylorF2RedSpin}$ waveform model. The bank performs well to match the injected signals throughout the parameter space for NSBH and BNS systems except for those regions where both the total mass and reduced spins are high.

In order to carry out a high-precision test comparing the efficiency of the hybrid bank with that of the vanilla-stochastic bank, we need a large number of injections to calculate the fitting factors upto high significant figures. To this end, we construct both the template banks using $\mathrm{TaylorF2RedSpin}$ metric and parameters given in Set-II of Table 1. The geometric random bank was constructed by initializing $\mathscr{U}$ with $1\times 10^{7}$ uniform random points in dimensionless chirp-time coordinates whereas the stochastic bank code was set to terminate when the rejection rate reached a value $0.9996$ averaged over the last $10$ acceptances. The geometric random bank was found to contain $107,547$ templates whereas the stochastic bank contained ($25\%$ more) $134,563$ template points. In this case, the geometric-stochastic bank took $\gtrsim 11$ times less time than the vanilla-stochastic bank on a single unloaded processor.

We quantify the performance of these two template banks by computing fitting factors for two different injection families of aligned spin waveforms: $\mathrm{TaylorF2RedSpin}$ and TaylorF2. $50,000$ injections were made in both cases where the intrinsic parameters of the injected waveforms were chosen from Set-II of Table 1 and other parameters chosen as earlier. For the case of geometric-random template bank created using local TO lattice, we found $0.354\%$ signals were recovered below a fitting factor of 0.97 for $\mathrm{TaylorF2RedSpin}$ injections, and $2.892\%$ signals are recovered below this level for TaylorF2 injections. In the latter case, the injection model is different from the one used to construct the template bank - hence it is expected that the coverage for TaylorF2 will be less. The corresponding numbers for the vanilla-stochastic bank are found to be $0.514\%$ ($\mathrm{TaylorF2RedSpin}$ injections) and $2.408\%$ (TaylorF2 injections) respectively. From these numbers it is evident that the geometric-stochastic bank is equally efficient as the vanilla-stochastic bank. These results are summarized in Fig. 6.

As mentioned earlier, the $A_{3}^{*}$ lattice provides optimal coverage for conformally flat spaces in 3D. One may be tempted therefore to use such a lattice as a seed for stochastic template placement. We have investigated this approach and compared it with the geometric-random and vanilla-stochastic algorithms. At first, the deemed parameter space (corresponding to parameters in Set-II of Table 1) was covered by a $A_{3}^{*}$ lattice using the metric at a putative point to determine the dimensions of the unit cell. Using a point in this volume for which the unit cell had the smallest dimension, we placed $31,732$ $A_{3}^{*}$ lattice points to entirely cover the deemed parameter space. Using these as seed points for stochastic placement, the final bank size was found to have $128,185$ templates. This is marginally ($\sim 5\%$) smaller than the vanilla-stochastic bank which has $134,563$ templates as reported above. The geometric-random bank outperforms this $A_{3}^{*}$ seeded stochastic bank by a good margin of more than $16\%$.

As remarked earlier, the vanilla-stochastic algorithm can be cast in two different ways. The traditional bottom-up approach has been implementation in the LSC Algorithm Library and has been used in this work for comparison with the geometric-random bank. In order to compare it head-to-head with the top-down approach we implemented it in software and ran it for the exact same parameters as given in Set II of Table 1. As expected, the top-down implementations gives nearly identical bank sizes ($0.1\%$ difference in size) but takes less than half the time as the bottom-up LAL implementation. The computational speed comes from the fact that in this implementation, one can reject many proposals that lie within the minimal match ellipsoid centred around a single accepted proposal. Efficient computational data structures like binary search trees are readily available for such queries.

The summary of various templates banks referred to in the above discussion is available in Table 2.

V Discussion and outlook

Templated matched filtering is the mainstay of gravitational wave detection pipelines. With unprecedented improvement in low frequency sensitivity of advanced detectors, and the availability of theoretical spinning waveform models, it has become imperative to conduct these searches over increasingly larger volumes in higher dimensional parameter spaces. For such cases, the stochastic algorithm is used for template placement as it is easily scalable to higher dimensions. But it is computationally expensive and by design not the most space efficient.

This paper introduces a new template placement algorithm in 3D with an attempt to combine the space efficiency of $A_{3}^{*}$ lattice along with the robustness of stochastic placement algorithm. Such a template bank can be used in gravitational wave searches from binary neutron stars and neutron-star-blackhole compact binary systems where the waveform is described by two mass parameters and a mass-weighted spin magnitude parameter providing coverage for aligned-spin systems.

The truncated octahedron which is a Dirichlet-Voronoi polytope of the $A_{3}^{*}$ lattice, inscribed within the minimal match ellipsoid is used as a unit cell for the geometric placement. Such lattices are known to provide optimal coverage for conformally flat 3D spaces. While the signal manifold is not globally flat, one can assume local flat patches and use such cells to cover them optimally. The interface to the stochastic placement algorithm is made by spraying random points over the parameters space - which are removed if found within the minimal match ellipsoid of any template. We discuss how this merger of methods is able to handle the issues arising out of varying curvature and irregular boundaries. The nuances of its optimal implementation are discussed in detail. We show in a direct comparison with stochastic algorithm that the new method generates significantly fewer templates and is computationally more efficient in Table2. We now make a couple of comments related to the new template placement algorithm presented in this paper.

One of the key issues with the geometric template placement algorithms (e.g. geometric hexagonal placement Cokelaer (2007)) is related to the amount of fine tuning needed in the method to account for curvature effects. For example, in 2D the hexagonal cells used to place templates need to be whittled or pared to a slightly smaller size so that no holes are left due to the changes in relative orientation of neighboring hexagonal cells. This reduces the overall efficiency of placement leading to increase in template bank sizes. But more importantly, one loses the ability to make generic codes that are usable for different waveform models. In this regard, the geometric-random algorithm proposed here is robust against such fine-tuning by design. The template placement proceeds by first spraying a large number of random points over the parameters space which are later removed if found within the minimal match ellipsoid of the templates in the bank. Suppose a small portion of the deemed space is left uncovered due to curvature effects, it would lead to some residual points that are not removed from $\mathscr{U}$. As evident from Algorithm 1, these points are revisited in subsequent iterations where a random point (out of the residuals) is added to the template bank leading to complete coverage.

We would like to point out that the template banks shown in this paper make use of the $\mathrm{TaylorF2RedSpin}$ waveform family which model only the inspiral part of the evolution of the compact objects. They are not ideal for BBH searches directly, where a significant part of the SNR is expected to be contributed from the merger and ringdown phases. However the method developed in this paper is quite general and can be used for placing templates for IMR waveform families as well. If analytical metrics for IMR approximants (IMRPhenomDKhan et al. (2016); Husa et al. (2016), SEOBNRv2Taracchini et al. (2014); Pürrer (2014), etc.) were available, the proposed method would be an elegant and efficient solution for covering the entire combined BNS+NSBH+BBH space. This will thereby mesh well with modern LIGO searches which are carried out over a large range of parameters with combined BNS+NSBH+BBH banks using a mix of inspiral- only and IMR templates.

The assumption that the metric on the signal manifold is slowly varying and is locally flat is crucial for space efficiency. One can imagine a hypothetical case where this assumption does not hold true (e.g. metric coefficients are random at every point) in which case, the geometric-random algorithm will effectively fall back to the top-down version of the vanilla-stochastic template placement by design. In other words, the number of templates in the bank from the new method will not exceed the stochastic template bank in the limiting case. We have also shown in a direct comparison that the top-down stochastic bank implementation is computationally more efficient and should be used where the metric is available. Incorporating an intelligent way of spraying the random proposals (instead of drawing them from a uniform distribution) over the parameter space may lead to further optimization of this method.

On the other hand if the metric was perfectly flat (and given in coordinates where it was constant), the hybrid construction would fall back on a perfect $A_{3}^{*}$ lattice.

Finally, the geometric-random placement method presented here for 3D is generically scalable to signal manifolds in higher dimensions by using the appropriate $A_{n}^{*}\;(n>3)$ lattices.

Appendix A The $\mathrm{TaylorF2RedSpin}$ signal model for gravitational waves from inspiraling compact binary coalescence

The TaylorF2 “reduced spin” waveform model in frequency domain is given byAjith (2011)

[TABLE]

where, the amplitude $\mathcal{A}$ depends on the component masses, distance to the source, sky position, orientation of the binary’s plane, and $\Psi(f)$ is the instantaneous phase which can be explicitly written as:

[TABLE]

where $t_{0}$ is the time of arrival of the signal at the detector marking the epoch at which the instantaneous frequency takes a fiducial value, $\phi_{0}$ is the corresponding phase, $v_{f}\equiv(\pi mf)^{1/3}$ is the instantaneous velocity, $m=m_{1}+m_{2}$ is the total mass and $\eta=m_{1}m_{2}/(m_{1}+m_{2})^{2}$ is the symmetric mass ratio of the binary and $\gamma_{E}$ is the Euler gamma constant.

The spin effects are encoded through $\beta$, $\sigma_{0}$ and $\gamma_{0}$ which appear at 1.5PN, 2PN and 2.5PN phase terms respectively and are given by:

[TABLE]

where the reduced spin parameter $\chi_{r}$ is defined as the weighted sum of individual spins $\chi_{1}$ and $\chi_{2}$ of the component masses as:

[TABLE]

The individual spins of the components $\chi_{1,2}$ are the projections of their spin vectors $\bm{S}_{1,2}$ along the Newtonian orbital angular momentum vector $\bm{{L}}_{N}$ and defined as:

[TABLE]

Appendix B Space filling truncated octahedron

A polyhedron is a three dimensional solid that has a finite number of polygon faces. One can fill a 3D space completely without any overlap or gap through the tessellations of space-filling polyhedra. Examples of such space-filling polyhedron include cube, hexagonal prism etc. Solution to such space-filling problems find many practical applications like optimal placement of a network of communication towers (Alam and Haas, 2006). The template placement problem addressed in this paper can also be mapped to an optimal space-filling problem in curved space. The geometric properties of optimal space-filling polyhedra can be understood from the volumetric quotient $(Q_{v})$ defined as

[TABLE]

where $V$ is the volume of polyhedron, $R$ is maximum distance from its center to any vertex and $S$ is the surface area of the polyhedron. Note that $(Q_{v})$ is the inverse of the thickness, a common measure for the efficiency of a given covering. The polyhedra with the highest $Q_{v}$ is expected to provide the optimum coverage.

The TO is a 14-faced Archimedean solid, with 8 hexagonal faces, 6 square faces and has 24 vertices. It is generated by joining two regular pyramids upside down and cutting a pyramid from all six vertices in such a way that the length of all the sides generated are equal. Thus a truncated octahedron of side $a$ can be created by removing six pyramids of side $a$ from an octahedron of side $3a$. Fig. 7 shows the TO obtained from two pyramids. The geometric properties of a TO and pyramid are given in Table 3.

Suppose a three dimensional volume is covered by the tessellations of TO cells. Each face of such a cell can be shared by another neighboring TO cell. There are two kind of neighboring cells: ones that share the square faces which we shall call S-neighbors and the others that share the hexagonal faces which we shall refer to as H-neighbors. Each TO has a maximum of six S-neighbors and eight H-neighbors. The distance between S-neighbors is twice the height of TO along square face and for H-neighbors, twice of height of TO along the hexagonal face. When a TO is inscribed in a sphere of radius $R$ such that the $z$-axis goes through the centre of one of the square faces, then, the sides of the squares and hexagons $a$ is given by $a=\sqrt{\frac{2}{5}}R$, the distance from the centre to each of the S-neighbors is equal to $\sqrt{\frac{12}{5}}R$, and the distance to the H-neighbors is $\frac{4}{\sqrt{5}}R$. The coordinates of all the 14 neighbours are listed in Table 4.

Acknowledgements.

A.S. would like to thank LS for motivation. SR thanks IIT Gandhinagar for research fellowship and fellow graduate students (Amit Reza, Chakresh Singh, Md. Yousuf) for useful discussions and help with the manuscript. We thank Samarth Vaijanapurkar (B.Tech. 2016 student of IIT Gandhinagar) for useful discussions. The authors thank Duncan Brown for valuable feedback. Last but not least, the authors thank the anonymous referees whose careful reading of the manuscript, detailed queries and many suggestions have helped to improve the content and presentation of this work.