{'fraction_non_alphanumeric': 0.030064790391559324, 'fraction_numerical': 0.011754462342185459, 'mean_word_length': 4.572998430141287, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Alternatives to standard puncture initial data for binary black hole evolution

Standard puncture initial data have been widely used for numerical binary black hole evolutions despite their shortcomings, most notably the inherent lack of gravitational radiation at the initial time that is later followed by a burst of spurious radiation. We study the evolution of three alternative initial data schemes. Two of the three alternatives are based on post-Newtonian expansions that contain realistic gravitational waves. The first scheme is based on a second-order post-Newtonian expansion in Arnowitt, Deser, and Misner transverse-traceless (ADMTT) gauge that has been resummed to approach standard puncture data at the black holes. The second scheme is based on asymptotic matching of the 4-metrics of two tidally perturbed Schwarzschild solutions to a first-order post-Newtonian expansion in ADMTT gauge away from the black holes. The final alternative is obtained through asymptotic matching of the 4-metrics of two tidally perturbed Schwarzschild solutions to a second-order post-Newtonian expansion in harmonic gauge away from the black holes. When evolved, the second scheme fails to produce quasicircular orbits (and instead leads to a nearly head-on collision). This failure can be traced back to inaccuracies in the extrinsic curvature due to low order matching. More encouraging is that the latter two alternatives lead to quasicircular orbits and show gravitational radiation from the onset of the evolution, as well as a reduction of spurious radiation. Current deficiencies compared to standard punctures data include more eccentric trajectories during the inspiral and larger constraint violations, since the alternative data sets are only approximate solutions of Einstein's equations. The eccentricity problem can be ameliorated by adjusting the initial momentum parameters.


I. INTRODUCTION
One of the best candidates for observation of gravitational radiation is the coalescence of two black holes. The successful detection of these gravitational waves by any ground or space-based interferometer (e.g. LIGO [1,2], VIRGO [3,4], GEO600 [5], eLISA/NGO [6]) will be confirmed by matched filtering of the observed signal against an extensive compilation of waveforms produced numerically from the binary black hole parameter space. Therefore, it is important that this collection of constructed waveforms represent the most physically viable scenarios of binary black hole evolution.
Full numerical evolution of the Einstein equations made possible since 2005 [7][8][9] is the method of choice for producing waveforms than span part of the inspiral, and the subsequent merger and ringdown phases. The waveform templates are dependent on the physical accuracy of the data being used. Currently, the most widely used approach for initial data construction is known as the "puncture" method. Developed by [10], this method is able to set up data for binary black holes with generic momenta and spins. These initial data are used almost exclusively by all groups that evolve black holes using the BSSNOK formulation [11] together with the moving punctures approach [8,9]. Standard puncture initial data have the disadvantage that they contain no realistic gravitational waves, because a conformally flat metric is used in this method. When evolved, these data lead to a gravitational wave signal that is zero for some time followed by a spurious burst. Only after this burst do we get an astrophysically realistic chirp signal, thereby making the prescription inherently unphysical. Other methods have begun to emerge whereby post-Newtonian (PN) approximations are used, either as the sole contributor [12][13][14] or in conjunction with analytical solutions near the black hole [15][16][17], to develop initial data. These methods allow for increased astrophysical accuracy by including, to leading PN order, gravitational waves on the spacetime from the onset of the simulation. In this paper, we investigate three methods of initial data construction and use them for full numerical evolution of binary black holes. These methods can be considered alternatives to the standard puncture technique for creating templates. We find that evolutions of two of these methods lead to improved gravitational waves. There is a chirp signal right from the start, and the (still present) burst of spurious radiation has a smaller amplitude than for evolutions starting from standard puncture data. The burst occurs when the built-in post-Newtonian waveform transitions to a fully numerically generated signal.
Throughout this paper, we will use units where G = c = 1. Latin indices such as i run from 1 to 3 and denote spatial indices, while Greek indices such as µ run from 0 to 3 and denote spacetime indices. The paper is organized as follows: In Section II, we briefly recall the 3+1 formulation for numerical relativity. In Section III, we summarize the standard puncture technique for creating initial data and discuss weaknesses in this approach. In Section IV, we describe our first alternative scheme; a post-Newtonian expansion of the 3-metric and extrinsic curvature in Arnowitt, Deser, and Misner transversetraceless (ADMTT) gauge [18] resummed to approach standard puncture data at the black holes. In Section V, we describe the technique of asymptotic matching the 4-metrics of tidally perturbed Schwarzschild black hole solutions to post-Newtonian expansions away from the black holes. Also within this section, we describe two more initial data schemes: asymptotic matching with a low-order post-Newtonian expansion of the 4-metric in ADMTT gauge away from the holes, and asymptotic matching of the 4-metrics with a second-order post-Newtonian expansion in harmonic gauge away from the black holes. In Section VI, we describe additional adjustments required for numerical evolution. In particular, we describe an algorithm we use to fill black hole interiors for initial data that contain singularities. We then show comparisons of evolutions starting from different initial data schemes. Section VII summarizes our results and includes discussion of future improvements to certain schemes.

II. 3+1 FORMULATION
For numerical evolutions of general relativity, it is useful to split the four-dimensional spacetime described in terms of the metric g µν into three-dimensional spatial hypersurfaces and time. The Arnowitt-Deser-Misner [19] decomposition of the Einstein equations defines the line element as where α is the lapse function, β i is the shift vector, and g ij is the 3-metric. The extrinsic curvature is defined as Under the Arnowitt-Deser-Misner decomposition, the Einstein equations yield evolution equations of the form and the Hamiltonian and momentum constraint equations Here, R ij and R are the Ricci tensor and scalar of the 3metric g ij , and D i is the derivative operator compatible with g ij . In binary black hole simulations, these equations simplify due to the vanishing of the source terms ρ, j i , S ij , and S = g ij S ij .

III. STANDARD PUNCTURE INITIAL DATA
The most common approach to producing initial data used currently is the standard puncture prescription of Brandt & Brügmann [10] using Bowen-York extrinsic curvature. This method uses a conformal transversetraceless decomposition of the extrinsic curvature [20,21] to simplify the constraint equations. A conformally flat three-metric, g ij , is defined as where ψ is the conformal factor. The extrinsic curvature is given by where is of Bowen-York form and where n a A is the radial normal vector, p a A is the linear momentum, and s a A is the spin of black hole A. This K ij is chosen so that it already satisfies the momentum constraint for arbitrary spin and momentum.
The exact solution for a black hole pair with zero momenta and spins is the Brill-Lindquist solution [22] where m A is the "bare" puncture mass of black hole A and r A = | x − x A | is the distance from the black hole at position x A . As one can see, this conformal factor diverges at the puncture location x A . Such divergences also occur for non-zero masses and spins. To avoid numerical problems with these divergences when solving the Hamiltonian constraint equation, one usually splits the conformal factor [10] ψ = ψ BL + u into a regular piece u and the divergent piece ψ BL . By doing this, the Hamiltonian constraint is reduced to a single elliptic equation (12) without any divergences, that can be easily evaluated by elliptic solvers such as in [23]. Notice that A ij BY = 0 for zero momentum and spin, so that u = 0. Furthermore, if one of the masses (e.g. m 2 ) is zero we obtain the Schwarzschild solution in isotropic coordinates.
The fact that data are built using the assumption of conformal flatness is an inherent problem. It is unrealistic to demand conformal flatness in the context of astrophysical constructions, because it implies that there are no gravitational waves on the initial time slice.
Often, gravitational waveforms are described in terms of the Newman-Penrose scalar, Ψ 4 , extracted at some distance from the sources. Usually this Ψ 4 is decomposed into individual ℓ and m modes by projecting onto spherical harmonics Y −2 ℓm of spin weight -2 [24]. Figure 1 shows the real part of the dominant ℓ = 2, m = 2 waveform mode of ψ 4 over time extracted at a distance of 90M . We observe a lack of any gravitational radiation until ∼ 80M ; this delay in the waveform signal at arbitrary extraction radii is common of standard punctures. Also, we see the beginning of an outward traveling wave profile only after about ∼ 150M and differing significantly from the initial burst; the waves up until this outward profile region have been commonly referred to as spurious radiation.
If such a simulation is used to make gravitational wave templates, the early part of the wave signal has to be discarded up to the time when the spurious radiation has become sufficiently small. How much exactly one has to discard depends on the desired accuracy, but in our experience it typically amounts to between 10% and 20% of the total simulation time. This kind of waste will only increase with future code upgrades that will allow us to extract the waves at larger distances, since then we have to wait longer for the spurious burst to even reach the extraction distance. These problems are the reason why we seek alternative constructions of binary black hole initial data. In [25] it was shown that spurious radiation can be reduced by using a 3-metric that is not conformally flat.

IV. INITIAL DATA BASED ON POST-NEWTONIAN EXPANSION IN ADMTT GAUGE
One alternative to standard puncture data has been introduced by Tichy et al. in [12] with further developments in [13,14]. This method uses post-Newtonian expansions in ADMTT gauge [18] for the 3-metric and extrinsic curvature on the initial time slice. While PN theory is strictly valid only away from regions of strong gravity such as black holes, it is possible to rewrite the 3metric and extrinsic curvature in ADMTT gauge in such a way that they approach the puncture form as in Eq. (7) and Eq. (9). In this way one can obtain globally valid data provided the initial black hole separation, r 12 , is within the PN regime. This is seen as a first step away from the standard puncture approach and a step towards constructing astrophysically realistic initial data. Using ADMTT gauge, the expressions for the 3-metric and the extrinsic curvature are taken up to O(v/c) 5 and O(v/c) 4 respectively, where v ∼ M/r 12 with M = m 1 + m 2 and c is the speed of light. A formal expansion parameter, ǫ ∼ (v/c), is used below to distinguish terms of different PN order. The 3-metric in ADMTT gauge is given as The full expression for h T T ij has been computed by Kelly et al. [13]. It denotes the gravitational wave contribution to the 3-metric. The conformal factor is very similar to that of standard punctures in Eqs. (10) and (11). The only difference is that it contains an additional term at order O(v/c) 4 and no u piece as in Eq. (11) since the constraints are satisfied only approximately. It must be noted that the PN expressions in Eqs. (13) and (15) are not pure Taylor expansions in (v/c) since terms such as ψ 4 in Eq. (13) are not expanded in powers of ǫ. This amounts to adding specific higher-order terms for the purpose of creating a metric that is similar to the standard puncture expression. This selective inclusion of terms does not improve the accuracy of the PN expansion, however it is performed to ensure the presence of black hole apparent horizons in the initial data [12]. The extrinsic curvature is where the leading order term contains the Bowen-York extrinsic curvature as in the case of standard punctures. The additional terms at order O(v/c) 5 can be found in [12]. The trace of the extrinsic curvature can be shown to vanish up to O(v/c) 6 .

V. ASYMPTOTIC MATCHING OF 4-METRICS: AN OVERVIEW
Another alternative to the standard puncture method involves a post-Newtonian expansion of the 4-metric away from the black holes. This expansion is then asymptotically matched to tidally perturbed black hole metrics used for, and in close proximity to, the black holes. This approach is believed to produce more astrophysically realistic initial data based on several reasons; the method is built on physical approximations and not on the potentially unphysical ansatz of assuming the punctures form close to each black hole, and the level of analytic control over the physical approximations allows for the method, in principle, to accommodate ever higher-order expansions if desired.
Matching requires dividing the spacetime into 4 zones; 2 inner zones, a near zone, and a far zone. The inner zones are considered as regions where perturbative Schwarzschild solutions are valid. The near zone is the region where PN theory will hold (provided (v/c) ≪ 1). The far zone is the region of spacetime where retardation effects will influence the system. Within each zone, 4-metric approximations are constructed in some coordinate system with corresponding parameters (mass, momentum, etc.). Matching these different approximations becomes possible because their regions of validity overlap in certain regions called buffer zones. Figure 2 shows a schematic sketch of the different regions. Asymptotic matching [26][27][28][29] involves comparing two asymptotic solutions inside of a 4-volume. This returns a map between the local coordinates and parameters of each of the different regions, and forces both solutions to be asymptotic to each other within the buffer zone. It is inside the buffer zone where these solutions can be merged into a smooth global metric by way of a transition function. We refer the reader to [15][16][17] for our choice of transition functions.
For a coordinate separation of r 12 the gravitational wavelength is given by A diagram of the regions is shown in Figure 2, and a summary of the regions' influences is seen in Table 1.  One of the first implementations that has used the above asymptotic matching approach is described in [16]. In this initial data scheme, a PN 4-metric in ADMTT coordinates [18,30,31] is used in the near zone. Since in this approach the Post-Newtonian metric is accurate only up to terms of order O(v/c) 2 , there are no retardation effects and the near zone results are also valid in the far zone. In the inner zones, a perturbed Schwarzschild metric in isotropic coordinates is used. A globally valid 4-metric is then obtained by asymptotic matching. The matching procedure is simplified by the fact that ADMTT and isotropic coordinates are similar. The matched result has errors of O(v/c) 4 away from the black holes and errors of order O(v/c) 2 O(r A /r 12 ) 3 near each black hole. The PN part of the matched results that is valid away from the black holes has the same 3-metric and extrinsic curvature as in Eqs. (13) and (15), but with all terms of order ǫ 4 and higher dropped. The perturbed black hole solutions that are valid near each hole can be found in [16]. Both metrics are matched in the buffer zones and smoothed through the transition function mentioned earlier.

B. Matching a Harmonic-Gauge Post-Newtonian
Metric to a Black Hole Perturbative Metric at Higher Order The second matching approach we consider here is from work by Johnson-McDaniel et al. [17]. It applies higher-order asymptotic matching of the 4-metric solutions in the buffer zones. For this approach the inner zones approximate the 4-metric in a quasi-Cartesian form of Cook-Scheel harmonic coordinates [32], the near zone now approximates the 4-metric in harmonic coordinates, and the far zone region is now considered through approximating the 4-metric in harmonic coordinates. The reader is referred to [17] for the complete construction of the 4-metric within the near zone regime. The inner zone solution is expressed to higher-order than what was used in the previous matching approach in ADMTT gauge and is found in its full form in [17]. The inclusion of a far zone contribution introduces an additional buffer zone, wherein another transition function will smoothly join the near and far zone approximations in the region where the errors from both approximations are comparable.
The main difference between this matching approach and the one prior is that it is now carried out to higher order; in both the PN expansion of the 4-metric away from the black holes and also the perturbation expansion of the inner zone black hole solutions. Also, the PN expansions are now in harmonic gauge. The matched result has errors of O(v/c) 5 away from the black holes and errors of order O(v/c) 5 O(r/r 12 ) 3 near each black hole. This result is then further augmented by additional (but formally unmatched) PN terms away from the black holes which include gravitational radiation, so that away from the black holes one has only an error of O(v/c) 6 . The choice of harmonic coordinates is one of convenience; matching has already been established between the near and far zones through [33] which simplifies the inclusion of the far zone region. This high-order matching yields more accurate results than the previous lower-order matching in ADMTT gauge. With the inclusion of initial gravitational radiation, we can expect improvements in the data set over the standard punctures approach.
In [34] the tidal deformations present in the inner zones of [17] (but not the near and far zones) are added to a superposition of two spinning black holes in Kerr-Schild coordinates. Evolutions display a reduction in the (2,0) mode of the outgoing gravitational waves. The expectation is that a full inclusion of outer zone phenomena, as in our work, will also reduce the (2,2) mode spurious signal.

VI. NUMERICAL EVOLUTION RESULTS
Evolutions from all initial data were done using the BAM code [24,[35][36][37]. The gravitational fields are evolved using the BSSNOK formalism [11,38,39] in the variation known as the "moving punctures" method [8,9]. Thus the 3-metric g ij is written as where the conformal metricγ ij has unit determinant. In addition, the extra variablẽ is introduced whereγ ij is the inverse of the conformal metric. Furthermore, the extrinsic curvature is split into its trace-free part,Ã ij , and its trace, K, which is given by The particulars of our BSSNOK implementation can be found in [24]. All our evolutions use equal-mass non-spinning black hole binaries at initial separation of 10M . Also, all sets are tested through evolutions using a resolution of ∼ 3M /224 near the black holes. Furthermore, all sets use the same expressions for α and β i at the initial time; To evolve the lapse and shift we use the "1+log" slicing condition [40] modified with the addition of an advection term as used in [9,41], and the "gamma freezing condition" [42] modified and used in [24,41], where it is labeled as the "shifting-shift case" in the former and the "000" shift choice in the latter. For the initial data obtained by matching or from PN approaches, the constraints are only approximately satisfied. In the case of standard punctures the constraints are solved numerically using the approach in [23].

A. Filling Black Hole Interiors
In the case of initial data obtained by asymptotic matching, each black hole is described by a perturbed black hole solution. This black hole solution can contain physical singularities inside the event horizons that have to be dealt with before a numerical evolution can be attempted. We have developed an algorithm (called BHfiller in our code) that fills a specified region inside the black hole apparent horizon with smooth data. The justification for this procedure is that the physics outside a black hole remains unaffected if we only change the inside of the horizon. Our algorithm modifies the BSS-NOK variables at the initial time as follows: For each black hole (centered at the point x A ) we pick a sphere of radius r f ill contained inside the black hole horizon. The radius r f ill is chosen such that this sphere contains all singular points. To set valid data at each point x inside this sphere we use a weighted average of linear extrapolation and a value that corresponds to standard puncture data. Let us define r = | x − x A |,n = ( x − x A )/r and denote a particular BSSNOK variable component at point x by u( x). If r < r f ill we set where otherwise. (24) and the weight is given by We apply this method of filling the region inside radius r f ill to all BSSNOK variables exceptΓ i , which we simply recompute from the filledγ ij using Eq. (18). We only need to use this algorithm for the case of asymptotically matched data sets, as all other data sets are free of any physical singularities.
Note that the constraints are generically not satisfied in the filled region. We have tested our filling algorithm by evolving standard puncture data with filling, applied at the initial time or some predetermined later time, and without filling and compared various quantities. Gauge invariant quantities, such as the gravitational wave amplitude as a function of gravitational wave frequency, are unchanged when we compare with and without black hole filling. We have also checked that in our evolutions no visible constraint violations are emitted by the black holes filled with this algorithm. The latter is expected since the BSSNOK system together with the gauge conditions used here has been shown [43] to lead to causal constraint propagation, so that any constraint violations introduced inside the black holes cannot affect the exterior spacetime. Notice, however, that the BSSNOK system has superluminal gauge modes [43]. Thus gauge dependent quantities such as the lapse are somewhat different, even outside the black holes. Hence the black hole trajectories (while qualitatively the same) are slightly different as well.

B. Shortcomings of Low-Order Asymptotic Matching
The real part of the dominant ℓ = 2, m = 2 waveform mode of ψ 4 illustrates whether a gravitational wave signal has been implemented starting from the t = 0 slice. Figure 3 shows the waveforms, extracted at radius r = 90M , for a simulation starting from standard puncture data and that produced from the first data set; asymptotic matching of a low-order PN expansion in ADMTT-Gauge to a perturbative black hole solution (hereafter called ADMTT-match). As with standard punctures, ADMTT-match fails to incorporate gravitational radiation at t = 0. This is not surprising as the 3-metric is still conformally flat.
Another undesired result is that the binary merges within a time of 250M , which is much too fast. This behavior indicates an accelerated evolution timescale that is not expected for a quasicircular inspiral. Figure 4 shows The (2,2) mode of ψ4 extracted at a radius of r = 90M for standard punctures and ADMTT-match initial data sets. Both waveforms show no initial wave signature at t = 0, due to their reliance on conformal flatness for the 3-metric expressions. It can also be seen that ADMTT-match merges much faster (within 250M ). This is evidence that the orbits were not quasicircular at all.
the orbit of the binary system in the ADMTT-match approach. We see that ADMTT-match leads to nearly head-on trajectories for the black holes. The individual black holes do not have enough momentum to complete even one orbit. One can track this problem back to the order of asymptotic matching used in the construction of these initial data. The terms in the extrinsic curvature that are associated with momentum are computed from the 4-metric perturbation expansions used near the black holes. However, as discussed in the Appendix, the matching procedure used in [16] is of an order that is too low to correctly obtain the the O(v/c) 3 piece of the extrinsic curvature that would lead to an orbiting black hole. Thus we see that this data set is merely a proof of concept. To obtain useful initial data, we need to use higher order matching.

C. Analysis of Higher-Order Methods
In this subsection, we compare evolutions starting from data obtained by the higher order matching approach discussed in Sec. V B and data from the pure PN approach of Sec. IV with standard puncture initial data. However, we first need to inform the reader of a modification to the transition function in [17] used for smoothing the far-zone and near-zone metric contributions in the outermost buffer zone. The parameters r 0 and w, which define where the transition begins and the width of the transition window respectively, have been adjusted. The analysis in [17] uses r 0 = λ GW /5 and w = λ GW [see Eqs. (8.4) and (8.6) in [17]], which then put the midpoint (where the transition function has a value of 1/2) at r ≈ 14λ GW . In the work here we choose values of r 0 = 0.044λ GW and w = 0.22λ GW , which lead to a midpoint at r ≈ λ GW . This adjustment is justified because one would expect the transition region (from near zone to far-zone) to be more suitable closer to λ GW than 14λ GW .

Waveforms
We first discuss the gravitational waveforms obtained from these data sets. In the figures that follow we denote our results as Harmonic-match or ADMTT-PN; depending on whether we evolve initial data constructed by asymptotic matching of a second-order post-Newtonian 4-metric expansion in Harmonic gauge to a perturbative black hole solution, or data from a post-Newtonian expansion in ADMTT gauge. We do not expect to see the same deficiencies as seen in ADMTT-match. Figure 5  shows the real part of the (ℓ = 2,m = 2) waveform mode of ψ 4 at an extraction radius of 90M for all three data sets. Harmonic-match produces gravitational waves starting from t = 0, more characteristic of an astrophysically realistic evolution. The diminished amplitude in the spurious radiation region for Harmonic-match is an additional positive outcome for this data set (although the ultimate goal is to completely eliminate this behavior from the waveform). Figure 5 also shows the waveform for ADMTT-PN and standard punctures. As with Harmonic-match, ADMTT-PN shows a non-zero waveform signature at t = 0 and comparable reduction of the spurious radiation amplitude. It is also noted that a comparative' reduction of the spurious radiation amplitude has been seen in [44] as well.

Eccentricity and Horizon Mass
We can see from Fig. 5 that the waveform behavior of both Harmonic-match and ADMTT-PN is nonmonotonic in amplitude. This differs from the standard chirp signal expected during the inspiral phase of the two black holes. The reason for this behavior is due to the black hole orbits being noticeably eccentric. The high eccentricity can be seen through plots of the coordinate separation over time (See Fig. 6). For a sense of scale, the coordinate separation for standard puncture initial data is included to indicate the small level of deviation for an evolution to be accepted as having low eccentricity. Note, however, that unlike in the case of standard punctures, the tangential momentum of the black holes has not been fine-tuned and is simply coming from PN expressions for circular orbits. For example, in the case of ADMTT-PN we use the momentum for circular PN orbits given by where M = m 1 + m 2 , µ = m 1 m 2 /M , η = µ/M and Ω is defined as (27) Increased eccentricity, with respect to standard punctures, is again apparent in Fig. 7 for both the Harmonicmatch and ADMTT-PN sets. The eccentricity calcula- The range between 100M and 600M is the range to analyze all sets simultaneously while guaranteeing no errors due to incomplete information (Before 100M our eccentricity estimator which averages over an orbital time is unreliable, and after 600M most sets have begun merger so that one cannot define an eccentricity.).
tion is performed as in [45], and is calculated through a time-average over the orbital period.
For the case of ADMTT-PN, we have added an additional parameter, C 3 , that can be used to adjust the initial tangential momenta of the black holes. The modified tangential momentum, p is given by Note that the term with C 3 comes into the expression at order O(v/c) 6 . Thus this additional term is of the order of one of the PN terms that were neglected in the expansion in Eq. (26). By adjusting C 3 we can then decrease the eccentricity. For C 3 = −396.4844 we change our tangential momentum to 98.3% of the original p P N . This choice effectively minimizes the eccentricity. In Figs. 6 and 7 we can compare the results with and without this modification to Harmonic-match and standard punctures.
A useful local definition for the mass of a black hole is the apparent horizon mass. An algorithm for finding apparent horizons without any symmetries is described in [46]. Over time, numerous variations on Gundlach's methodology have produced routines for many of the main computational groups involved in black hole simulation. The version used here is known as AHmod and was developed by Norbert Lages [47] as an improvement to the CACTUS thorn AHfinder developed by Miguel Alcubierre [48]. Figure 8 shows the apparent horizon mass of one of the (equal mass) black holes for all three data sets versus time. We see that ADMTT-PN and Harmonic-match show much more variation than standard puncture data. We see that the two alternative methods lead to a "mass loss" over the evolution time where they eventually approach the apparent horizon mass obtained from evolving standard puncture data. This "mass loss", albeit of low scale, is unphysical and we expect it is due to larger Hamiltonian constraint violations in the approximate initial data sets. For example, in [44] the Hamiltonian constraint (but not the momentum constraint) is solved numerically for ADMTT-PN, which leads to a horizon mass that is closer to the result for standard puncture data. Figure 9 shows the apparent horizons masses of the final black hole after merger for evolutions with each data set. The masses obtained with the alternatives tend to deviate from punctures, with Harmonic-match yielding an increase and ADMTT-PN yielding a decrease in the common horizon mass respectively over time. This decrease is certainly unphysical and is likely related to larger constraint violations.

Constraints
The violation of the Hamiltonian constraint equation at the beginning of the evolution, as seen in Fig. 10 is between 2 and 3 orders of magnitude larger for Harmonicmatch and ADMTT-PN than for standard punctures at t = 0. We expect that these violations are responsible for the observed drift in the apparent horizon masses. Over time, the evolution of ADMTT-PN relaxes the constraints to the regime of standard punctures where both are almost indistinguishable. This result brings the expectation that these initial data methods are stable over the life of the evolution, and a case for astrophysical relevance.  Figure 11 shows the L 2 -norm of the x-component of the Momentum constraint. It is seen that the initial discrepancies (compared to standard punctures) are of 2 orders of magnitude for ADMTT-PN and roughly 5 orders of magnitude for Harmonic-match. The y and z components are so similar that they would yield almost the same plots. Over time, all data sets relax to within an order of magnitude difference until they become almost indistinguishable during the merger and ringdown phases.

VII. DISCUSSION
We have investigated three alternative approaches of constructing initial data for a binary black hole system. Two of these alternatives involve asymptotic matching of post-Newtonian expansions to perturbative black hole solutions of the metric, and one uses a resummed post-Newtonian expansion of the metric everywhere. We numerically evolve all three alternatives using the BSS-NOK system [11] together with the moving punctures approach [8,9] and compare them with evolutions that start from standard puncture initial data. The two alternative data sets that were derived by asymptotic matching have not been evolved before. Since they contain singularities near the black hole centers, we have also developed a particular algorithm to fill a certain region inside each black hole with smooth data that approaches puncture data at the center. In this way we can avoid the use of black hole excision, so that we are able to use the same moving punctures approach to evolve all our data sets.
We find that the simplest data set, using low-order asymptotic matching, does not lead to physically accurate evolutions. Insufficient expansion of the extrinsic curvature leads to black hole momenta that are far too small for quasicircular orbits, so that we observe an almost head-on collision. Also, low-order expansion of the metric yields gravitational waveforms lacking an initial signature. Thus, this initial data set is not an improvement over standard punctures.
Higher-order expansion techniques resulted in the production of gravitational wave signals at the initial time as well as proper quasicircular inspiral orbits that subsequently lead to merger and ringdown. Eccentricity measures were higher for the Harmonic-match and ADMTT-PN sets than what is seen in standard punctures, when the puncture momentum parameters are fine-tuned to achieve low eccentricities. However, by tuning our 3PN order C 3 parameter (see Eq. (28)) we were able to bring the eccentricity of ADMTT-PN down as well. Both Harmonic-match and ADMTT-PN sets look promising in the sense that, unlike for standard puncture data, they have realistic gravitational wave signals built-in from the beginning of the evolution. In addition, they show a reduction of spurious radiation when we compare with standard puncture data.
There are, however, still problems with these alternative data sets that stem from the fact that they only approximately satisfy the constraint equations of general relativity. We have monitored the Hamiltonian and momentum constraint violations during the evolutions. We find that constraint violations of Harmonic-match and ADMTT-PN sets are above the usual level coming from numerical truncation errors observed when we evolve standard puncture data. In the case of ADMTT-PN these violations eventually decrease and reach a level similar to what is seen for standard puncture data evolutions. However, until this happens we see an unphysical drop in the apparent horizon masses. In the case of Harmonic-match, the constraint violations start at an even higher level and never quite fall to the level observed for standard puncture data evolutions, and again result in unphysical behavior of the apparent horizon masses.
Future improvements need to be focused on these constraint violations. An obvious remedy would be to project the alternative data sets onto the solution manifold of general relativity by solving the constraint equations, e.g. as in [12] by using the York-Lichnerowicz conformal decomposition [49]. If the constraints are satisfied we expect the apparent horizon masses to be strictly non-decreasing over time. Evolutions, horizon masses, and eccentricity could then be compared to the standard punctures approach and conclusions could be made on more stringent grounds.
Most astrophysical black holes are expected to be spinning. It is thus important to include spin in the initial data construction. As already mentioned in Sec. III, standard punctures contain spin parameters that enter the Bowen-York extrinsic curvature in Eq. (9), so that one can set up initial data with spins. Including spins in the alternative data sets discussed in this paper is possible as well. For any of the data sets that are produced from asymptotic matching one would have to redo the matching calculation as in [17] and match a post-Newtonian 4-metric with spin to a perturbed black hole metric with spin. This calculation, while certainly long and tedious, is in principle possible. The situation is much simpler for the post-Newtonian data set in ADMTT gauge discussed in Sec. IV. To leading order, the spin only enters the Bowen-York piece A ij BY in Eq. (15) in the same way as in the case of standard punctures. It would thus be easily possible to add these missing spin terms. This alone would already result in spinning black holes. However, in order to obtain the correct waveform at the initial time we have to ensure that h T T ij(4) in Eqs. (13) and (15) is computed for post-Newtonian particle trajectories where spin is included in the equation of motion. Since post-Newtonian trajectories for spinning particles are well known, and since h T T ij(4) is known for arbitrary trajectories [13], the inclusion of spin in this data set is not a complicated problem.
can only be fixed once matching has been performed two orders higher in v/c.
In the case of the lower order matching procedure used in [15,16], matching of the g ij and g 00 components of the two 4-metrics has only been performed up to O(v/c) 2 . This means that the piece (F βα ) 1 could not be computed. In fact, it was merely set to a value that simplified further calculations. This has important consequences. As shown in Sec. V.E. of [17] the only non-zero components of (F βα ) 1 are (F 02 ) 1  1/2 t in Y = X 2 (x β ), which corresponds to a boost. This boost is needed to transform the metric of the perturbed black hole solu-tion in inner zone 1 (from a frame where it is at rest to a frame where it is moving in the y-direction). In the case of lower order matching as in [15,16] this term is absent, so that the metric used near the black hole is missing this boost. The resulting extrinsic curvature (which is always smaller than the 3-metric by a factor of v/c) is thus correct only up to O(v/c) 2 and has errors of O(v/c) 3 , where the missing boost would enter. This leads to a black hole without sufficient momentum.
In general, this implies that if we want to obtain a 4-metric that is correct up to O(v/c) n with the matching procedure in [15][16][17] we need to find the matching coordinate transformation up to order O(v/c) n+2 .