A hybrid geometricrandom template placement algorithm for gravitational wave searches from compact binary coalescences
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 geometricrandom 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 starblackhole 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 vanillastochastic template banks and show that while both are equally efficient in the fittingfactor sense, the bank sizes are larger in the stochastic method. Further, we show that the generation of the proposed hybrid banks can be spedup 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.
pacs:
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 advancedLIGO (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 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 FrenchItalian 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 LIGOIndia 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 signaltonoise 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 matchedfilter 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 prefabricated 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 geometricrandom 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 geometricrandom 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 vanillastochastic placement method. We also construct an explicit template bank for neutron starblackhole searches using this new algorithm and compare it against a vanillastochastic template bank, noting the improvement in overall bank size as well as the time taken to generate the bank. The signal model Ajith (2011) was used for this purpose. The template banks presented in this work have been calculated using the 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 signal injections using 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 lattice seeded vanillastochastic 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 is the set of all possible GW signals characterized by the parameter vector . It is customary to represent the corresponding frequency domain signal as . The detector output consists of detector noise 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 , the latter is projected over the signal manifold by calculating the maximum likelihood over which serves as the detection statistic. For additive Gaussian noise, the likelihood is given by Finn (1992)
(1) 
The inner product is the complex crosscorrelation Sathyaprakash and Dhurandhar (1991) between the detector output and the gravitational wave signal weighed inversely by the noise power spectral density of the detector Babak et al. (2013):
(2) 
where is the onesided noise power spectral density defined by , asterisk denotes the complex conjugation operator, the frequency range marks the effective bandwidth of the detector and is the time delay between these two signals. The signaltonoise ratio (SNR) after filtering is defined as:
(3) 
Without any loss of generality, we assume that all template waveforms are normalized such that . From Eq. (1) it is clear that this allows us to use the loglikelihood function (or equivalently, the SNR) maximized over the parameters, as the detection statistic.
The loglikelihood function can be maximized over all time lags () 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 loglikelihood 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 loglikelihood surface at each of these points and the maximum is suitably determined.
The template bank consisting of these discrete set of points on is constructed using a control parameter ; 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 intertemplate 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 . 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 , is given by
(4) 
From the normalization of waveforms, it is clear that . For high values of this parameter, the LHS can be Taylor series expanded upto leading order in the small parameter :
(5) 
where are the components of the vector and
(6) 
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 rearranged as , and identified to be the equation of an ellipsoid in 3D centered at a point . For higher dimensions it represents a hyperellipsoid. We shall refer to this as the minimalmatch ellipsoid elsewhere in the paper. In GW searches from CBC sources, template banks are usually constructed at Babak et al. (2013). This corresponds to a loss in detection rate of 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 initialLIGO 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 instead of component masses , 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 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 semianalytic 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 alignedspin 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 socalled 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 bottomup 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 topdown fashion. In this alternate implementation, one starts with an empty template list and a list of very large number of proposals distributed uniformly over the deemed parameter space. One picks a random point and appends it to , following which all points from that lie within the minimal match distance from are removed. The process continues in an iterative manner until all points from are exhausted. We call this the topdown approach as the templatebank 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 topdown approach is more useful in projecting a geometric structure over the stochastic template bank leading to the hybrid geometricrandom placement algorithm described in this paper. In Table 2 we demonstrate that the topdown approach is also faster by a factor over the traditional bottomup 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 bruteforce 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 alignedspin 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 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 geometricrandom template placement algorithm in three dimensions using a truncated octahedral lattice. Such lattices are the DirichletVoronoi polytope of bodycentred cubic lattice Schürmann and Vallentin (2006). The latter provide optimal coverage for conformally flat spaces where the metriccoefficients 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 14sided 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 using signal model. This signal model is constructed using postNewtonian (PN) template family of gravitational waveforms from inspiralling compact binaries with nonprecessing spins. Here, the spin effects are captured by the parameter defined as a massweighted, 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 threedimensional geometric random template bank, we require the metric on the parameter space. The metric for approximant varies rapidly in coordinate system where denotes the total mass and denotes the symmetric mass ratio. In order to enhance the efficiency of the algorithm, one places the templates in dimensionless chirp time coordinates , in which the metric components are slowly varying over the parameters. These new coordinates are defined as:
(7) 
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 nonflat space. In the following sections, we show that the truncatedoctahedral 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 lattice.
The method outlined in Algorithm (1) proceeds by initializing three lists as follows:

: a list of uniformly distributed random points sprayed over the deemed parameters space

: an empty list for template points

: an empty temporary list
At first is initialized or seeded with a random point chosen from which is immediately appended to the list of templates as its first element. The possible 14 TO neighbours of this initial point are then calculated and appended to followed by the removal of the first element. Only those TO neighbors which fall within the parameter space and are farther than from existing points in the and can be appended. This completes one iteration of the inner loop of the algorithm. We continue these steps until 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 . At this stage, all points from that are within the minimal match distance from the elements in are deleted; following which is reseeded with a new random point from . Termination of template placement algorithm occurs when is exhausted.
The advantage of casting the stochastic template placement algorithm (see Sec II.1) in a topdown fashion is clear when we consider the extreme case where no TO lattice neighbours can be found for any proposal point in . In this case, it is clear from Algorithm 1 that the geometricrandom algorithm naturally falls back to the vanillastochastic placement algorithm. This brings in the robustness of the latter to the proposed new algorithm. This contingency arises for a small fraction () of points in  due to boundary edge effects and small uncovered patches arising from curvature effects. By construction the geometricrandom placement method presented here is more efficient than vanillastochastic algorithm.
Also note that the TO, whose neighbors are added to in Lines 89 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 truncatedoctahedron. 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 lattice.
iii.1 Implementation
We now highlight some salient features for efficient implementation of the algorithm.
iii.1.1 Initializing 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 and for the range of parameters over which the bank is to be placed. A set of random points are generated in 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 which allows the calculation of the extremities in these coordinates. The extreme values of are evaluated using 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 upto a critical value, beyond which it does not change very much. Empirically, we find this number to be about 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 minimalmatch ellipsoid
We discuss an efficient way to delete points from 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:
(8) 
where ’s are eigenvalues of the metric evaluated at the ellipsoid centre. The orientation of the ellipsoid is defined by matrix composed of eigenvectors of the metric , where is component of the eigenvector. To visualize the deletion process, let us consider a cell centred around a template. The points in lying inside the cell can be determined by calculating the match between this template with the all other points in 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 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 .
This can be further refined by binning the template points in such that each bin contains about templates. For each such subset, the removal of points from proceeds by identifying those points that lie in the same bin (within some acceptable margin in ) 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 . Additional binning in other two coordinates may possibly improve the computational efficiency of removing the points from even more.
iii.1.3 Global coordinate transformation
In dimensionless chirp time coordinates, most of the ellipsoidal cells have semiaxes ratio around . 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 semimajor 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 semiaxes ratio to a large extent. We use the eigenvalues and eigenvectors of the metric at a fiducial point to transform the coordinates () as,
(9) 
where, is a diagonal scaling matrix and the rotational matrix is constructed from the component of the eigenvector of the metric calculated at the fiducial point. This is carried out by first binning the templates in . 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:
(10) 
It is trivial to check that the metric at the fiducial points transforms to a identity matrix.
This global coordinate transformation Eq. (9) causes other nearby cells to become almost spherical with semiaxes in the ratio . In Fig. 1 we show the effect of global coordinate transformation on nearby cells. This is doubly advantageous: not only does the ballpoint 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 , the coordinates of its 14 lattice neighbors are available in Table 4. Here is an index on the coordinates of the 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 can be calculated using the rotation and scaling matrices:
(11) 
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:

The point is inside the deemed parameter space and also satisfies . The latter corresponds to the condition that ’s can be inverted to yield physical masses.

The point is not within the minimal match distance of existing points in and .
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 ) 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 as the placement proceeds. At this stage, we need to reseed and continue the process iteratively until all points in 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 which is copied over as the first element of the template bank list .
One can start from any point in parameter space: in the present work we have started from “midpoint” corresponding to component masses with individual spins . 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 reseed .
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 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 templates with 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; and IMRPhenomB signal model Kalaghatgi et al. (2015). Our demonstration makes use of the model.
iv.1 Construction of the template banks for reducedspin binary system
We generate a template bank using metric in parameter space for aligned spin compact binary system. The range of parameters is chosen such that the mass of the first object lies between with dimensionless spin magnitude in the range . The mass of the second object is taken between with dimensionless spin magnitude in the range . The NSBH boundary mass condition is satisfied i.e any object with mass 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 neutronstarblackhole (NSBH) searches.
We construct several template banks by varying the number of uniformly sprayed random points over the parameter space, . We consder sizes of varying between points. In all the cases, the placement proceeds from a point corresponding to individual masses and individual spin magnitudes . The full specification of input arguments for template generation are given in SetI 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 and directions respectively. The 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 SetI 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 templates which is larger than the geometricrandom bank. The computational run time was also recorded and we observed that the geometricrandom method took minutes while the vanilla stochastic method took minutes to execute on a single, unloaded processor, which is times faster.
Bank Parameter  SetI  SetII 

Waveform model  
Noise model  
Lower cutoff frequency  
Higher cutoff frequency  
Mass of first object  
Mass of second object  
Spin of first object  
Spin of second object  
Size of  
Minimal Match  0.97  0.97 
iv.2 Validation of the template banks
We investigate the performance of both the geometricrandom and the vanilla stochastic template banks against a set of signal injections from the reduced spin 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” is defined as a measure of the maximum match over the template bank for a putative injected signal :
(12) 
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 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 ( 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 SetI of Table 1, we injected signals from reducedspin 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 edgeon 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 ), we found that signals were recovered with fitting factor for both geometricrandom 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 nonprecessing 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 geometricstochastic template bank) over various combinations of intrinsic parameters of the compact binary system where both the signals and template waveforms are generated from 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 highprecision test comparing the efficiency of the hybrid bank with that of the vanillastochastic 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 metric and parameters given in SetII of Table 1. The geometric random bank was constructed by initializing with uniform random points in dimensionless chirptime coordinates whereas the stochastic bank code was set to terminate when the rejection rate reached a value averaged over the last acceptances. The geometric random bank was found to contain templates whereas the stochastic bank contained ( more) template points. In this case, the geometricstochastic bank took times less time than the vanillastochastic 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: and TaylorF2. injections were made in both cases where the intrinsic parameters of the injected waveforms were chosen from SetII of Table 1 and other parameters chosen as earlier. For the case of geometricrandom template bank created using local TO lattice, we found signals were recovered below a fitting factor of 0.97 for injections, and 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 vanillastochastic bank are found to be ( injections) and (TaylorF2 injections) respectively. From these numbers it is evident that the geometricstochastic bank is equally efficient as the vanillastochastic bank. These results are summarized in Fig. 6.
As mentioned earlier, the 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 geometricrandom and vanillastochastic algorithms. At first, the deemed parameter space (corresponding to parameters in SetII of Table 1) was covered by a 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 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 templates. This is marginally () smaller than the vanillastochastic bank which has templates as reported above. The geometricrandom bank outperforms this seeded stochastic bank by a good margin of more than .
As remarked earlier, the vanillastochastic algorithm can be cast in two different ways. The traditional bottomup approach has been implementation in the LSC Algorithm Library and has been used in this work for comparison with the geometricrandom bank. In order to compare it headtohead with the topdown 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 topdown implementations gives nearly identical bank sizes ( difference in size) but takes less than half the time as the bottomup 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.
Bank Parameters  Placement Algorithm  Size of  Bank Size  Execution Time (min)  Comments  
SetI of Table 1 


Hybrid construction  
GeometricRandom  
VanillaStochastic (lalapps_cbc_sbank)  
SetII of Table 1  Hybrid construction GeometricRandom 


VanillaStochastic (lalapps_cbc_sbank)  
seeded Stochastic  
VanillaStochastic (top_down)  134, 426 

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 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 neutronstarblackhole compact binary systems where the waveform is described by two mass parameters and a massweighted spin magnitude parameter providing coverage for alignedspin systems.
The truncated octahedron which is a DirichletVoronoi polytope of the 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 geometricrandom algorithm proposed here is robust against such finetuning 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 . 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 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 geometricrandom algorithm will effectively fall back to the topdown version of the vanillastochastic 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 topdown 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 lattice.
Finally, the geometricrandom placement method presented here for 3D is generically scalable to signal manifolds in higher dimensions by using the appropriate lattices.
Appendix A The signal model for gravitational waves from inspiraling compact binary coalescence
The TaylorF2 âreduced spinâ waveform model in frequency domain is given byAjith (2011)
(13) 
where, the amplitude depends on the component masses, distance to the source, sky position, orientation of the binary’s plane, and is the instantaneous phase which can be explicitly written as:
(14) 
where is the time of arrival of the signal at the detector marking the epoch at which the instantaneous frequency takes a fiducial value, is the corresponding phase, is the instantaneous velocity, is the total mass and is the symmetric mass ratio of the binary and is the Euler gamma constant.
The spin effects are encoded through , and which appear at 1.5PN, 2PN and 2.5PN phase terms respectively and are given by:
(15) 
where the reduced spin parameter is defined as the weighted sum of individual spins and of the component masses as:
(16) 
The individual spins of the components are the projections of their spin vectors along the Newtonian orbital angular momentum vector and defined as:
(17) 
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 spacefilling polyhedra. Examples of such spacefilling polyhedron include cube, hexagonal prism etc. Solution to such spacefilling 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 spacefilling problem in curved space. The geometric properties of optimal spacefilling polyhedra can be understood from the volumetric quotient defined as
(18) 
where is the volume of polyhedron, is maximum distance from its center to any vertex and is the surface area of the polyhedron. Note that is the inverse of the thickness, a common measure for the efficiency of a given covering. The polyhedra with the highest is expected to provide the optimum coverage.
The TO is a 14faced 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 can be created by removing six pyramids of side from an octahedron of side . Fig. 7 shows the TO obtained from two pyramids. The geometric properties of a TO and pyramid are given in Table 3.
Name  Value 

Side length of Pyramid  
Side length of TO  
Height of Pyramid  
Height of TO along square face  
Height of TO along hexagonal face  
Height of TO along vertices  
Volume of TO  
Surface Area of TO  
Volumetric Quotient of TO 
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 Sneighbors and the others that share the hexagonal faces which we shall refer to as Hneighbors. Each TO has a maximum of six Sneighbors and eight Hneighbors. The distance between Sneighbors is twice the height of TO along square face and for Hneighbors, twice of height of TO along the hexagonal face. When a TO is inscribed in a sphere of radius such that the axis goes through the centre of one of the square faces, then, the sides of the squares and hexagons is given by , the distance from the centre to each of the Sneighbors is equal to , and the distance to the Hneighbors is . The coordinates of all the 14 neighbours are listed in Table 4.
Neighborhood Type  Position 

SNeighborhood  
square faces  , 
HNeighborhood  
hexagonal faces  
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.References
 Abbott et al. (2016a) B. P. Abbott et al. (LIGO Scientific Collaboration and Virgo Collaboration), Phys. Rev. Lett. 116, 061102 (2016a).
 Harry and the LIGO Scientific Collaboration (2010) G. M. Harry and the LIGO Scientific Collaboration, Classical and Quantum Gravity 27, 084006 (2010).
 Abbott et al. (2016b) B. P. Abbott et al. (LIGO Scientific Collaboration and Virgo Collaboration), Phys. Rev. Lett. 116, 241103 (2016b).
 Acernese et al. (2015) F. Acernese et al., Journal of Physics: Conference Series 610, 012014 (2015).
 Somiya (2012) K. Somiya, Classical and Quantum Gravity 29, 124007 (2012).
 Aso et al. (2013) Y. Aso, Y. Michimura, K. Somiya, M. Ando, O. Miyakawa, T. Sekiguchi, D. Tatsumi, and H. Yamamoto (The KAGRA Collaboration), Phys. Rev. D 88, 043007 (2013).
 Iyer et al. (2011) B. Iyer, T. Souradeep, C. Unnikrishnan, S. Dhurandhar, S. Raja, and A. Sengupta, LIGOIndia Technical Report (2011), https://dcc.ligo.org/cgibin/DocDB/ShowDocument?docid=75988.
 Unnikrishnan (2013) C. S. Unnikrishnan, International Journal of Modern Physics D 22, 1341010 (2013).
 Manca and Vallisneri (2010) G. M. Manca and M. Vallisneri, Phys. Rev. D 81, 024004 (2010).
 Messenger et al. (2009) C. Messenger, R. Prix, and M. A. Papa, Phys. Rev. D 79, 104017 (2009).
 Cokelaer (2007) T. Cokelaer, Phys. Rev. D 76, 102004 (2007).
 Babak et al. (2006) S. Babak, R. Balasubramanian, D. Churches, T. Cokelaer, and B. S. Sathyaprakash, Classical and Quantum Gravity 23, 5477 (2006).
 Harry et al. (2009) I. W. Harry, B. Allen, and B. S. Sathyaprakash, Phys. Rev. D 80, 104014 (2009).
 Brown et al. (2012) D. A. Brown, I. Harry, A. Lundgren, and A. H. Nitz, Phys. Rev. D 86, 084017 (2012).
 Buonanno et al. (2009) A. Buonanno, B. R. Iyer, E. Ochsner, Y. Pan, and B. S. Sathyaprakash, Phys. Rev. D 80, 084043 (2009).
 Poisson and Will (1995) E. Poisson and C. M. Will, Phys. Rev. D 52, 848 (1995).
 Ajith et al. (2014) P. Ajith, N. Fotopoulos, S. Privitera, A. Neunzert, N. Mazumder, and A. J. Weinstein, Phys. Rev. D 89, 084041 (2014).
 Capano et al. (2016) C. Capano, I. Harry, S. Privitera, and A. Buonanno, Phys. Rev. D 93, 124007 (2016).
 Ajith (2011) P. Ajith, Phys. Rev. D 84, 084037 (2011).
 Shoemaker (2010) D. Shoemaker (LIGO Scientific Collaboration), LIGO Document T0900288v3 (2010).
 (21) LIGO Data Analysis Software Working Group, “LALSuite: LSC Algorithm Library Suite,” .
 Finn (1992) L. S. Finn, Phys. Rev. D 46, 5236 (1992).
 Sathyaprakash and Dhurandhar (1991) B. S. Sathyaprakash and S. V. Dhurandhar, Phys. Rev. D 44, 3819 (1991).
 Babak et al. (2013) S. Babak et al., Phys. Rev. D 87, 024033 (2013).
 Conway and Sloane (1999) J. H. Conway and N. J. A. Sloane, Sphere packings, lattices and groups, Vol. 290 (Springer Science & Business Media, 1999).
 Prix (2007) R. Prix, Classical and Quantum Gravity 24, S481 (2007).
 Owen (1996) B. J. Owen, Phys. Rev. D 53, 6749 (1996).
 Owen and Sathyaprakash (1999) B. J. Owen and B. S. Sathyaprakash, Phys. Rev. D 60, 022002 (1999).
 Abbott et al. (2009a) B. P. Abbott et al. (LIGO Scientific Collaboration), Phys. Rev. D 79, 122001 (2009a).
 Abbott et al. (2009b) B. P. Abbott et al. (LIGO Scientific Collaboration), Phys. Rev. D 80, 047101 (2009b).
 Abadie et al. (2010) J. Abadie et al. (LIGO Scientific Collaboration and Virgo Collaboration), Phys. Rev. D 82, 102001 (2010).
 Keppel et al. (2013) D. Keppel, A. P. Lundgren, B. J. Owen, and H. Zhu, Phys. Rev. D 88, 063002 (2013).
 Harry et al. (2014) I. W. Harry, A. H. Nitz, D. A. Brown, A. P. Lundgren, E. Ochsner, and D. Keppel, Phys. Rev. D 89, 024010 (2014).
 Pai and Arun (2013) A. Pai and K. G. Arun, Classical and Quantum Gravity 30, 025011 (2013).
 Nitz et al. (2016) A. Nitz, I. W. Harry, C. M. Biwer, J. Willis, D. Brown, L. Pekowsky, T. D. Canton, T. Dent, A. R. Williamson, C. Capano, P. Kumar, Lenona, S. De, micamu, S. Fairhurst, tjma12, A. Nielsen, Shasvath, S. Babak, B. Machenschalk, L. Singer, D. Macleod, S. Reyes, C. Sugar, Couvares, B. Bockelman, A. Lundgren, V. Tewari, F. Ohme, and J. Veitch, “ligocbc/pycbc: Er10 production release 2,” (2016).
 Usman et al. (2016) S. A. Usman et al., Class. Quant. Grav. 33, 215004 (2016), arXiv:1508.02357 [grqc] .
 Dal Canton et al. (2014) T. Dal Canton et al., Phys. Rev. D90, 082004 (2014), arXiv:1405.6731 [grqc] .
 Allen et al. (2012) B. Allen, W. G. Anderson, P. R. Brady, D. A. Brown, and J. D. E. Creighton, Phys. Rev. D 85, 122006 (2012).
 Hibbard (1962) T. N. Hibbard, J. ACM 9, 13 (1962).
 Privitera et al. (2014) S. Privitera, S. R. P. Mohapatra, P. Ajith, K. Cannon, N. Fotopoulos, M. A. Frei, C. Hanna, A. J. Weinstein, and J. T. Whelan, Phys. Rev. D 89, 024003 (2014).
 Abbott et al. (2016c) B. P. Abbott et al. (LIGO Scientific Collaboration and Virgo Collaboration), Phys. Rev. D 93, 122003 (2016c).
 Schürmann and Vallentin (2006) A. Schürmann and F. Vallentin, Discrete & Computational Geometry 35, 73 (2006).
 Thomson (1887) S. W. Thomson, Philosophical Magazine Series 5 24, 503 (1887).
 Douglas (1959) A. S. Douglas, The Computer Journal 2, 1 (1959).
 Windley (1960) P. F. Windley, The Computer Journal 3, 84 (1960).
 Kalaghatgi et al. (2015) C. Kalaghatgi, P. Ajith, and K. G. Arun, Phys. Rev. D 91, 124042 (2015).
 Apostolatos (1995) T. A. Apostolatos, Phys. Rev. D 52, 605 (1995).
 Khan et al. (2016) S. Khan, S. Husa, M. Hannam, F. Ohme, M. Pürrer, X. J. Forteza, and A. Bohé, Phys. Rev. D 93, 044007 (2016).
 Husa et al. (2016) S. Husa, S. Khan, M. Hannam, M. Pürrer, F. Ohme, X. J. Forteza, and A. Bohé, Phys. Rev. D 93, 044006 (2016).
 Taracchini et al. (2014) A. Taracchini, A. Buonanno, Y. Pan, T. Hinderer, M. Boyle, et al., Phys. Rev. D 89, 061502 (2014).
 Pürrer (2014) M. Pürrer, Classical and Quantum Gravity 31, 195010 (2014).
 Alam and Haas (2006) S. M. N. Alam and Z. J. Haas, in Proceedings of the 12th Annual International Conference on Mobile Computing and Networking, MobiCom ’06 (ACM, New York, NY, USA, 2006) pp. 346–357.