High accuracy binary black hole simulations with an extended wave zone

We present results from a new code for binary black hole evolutions using the moving-puncture approach, implementing finite differences in generalized coordinates, and allowing the spacetime to be covered with multiple communicating nonsingular coordinate patches. Here we consider a regular Cartesian near-zone, with adapted spherical grids covering the wave zone. The efficiencies resulting from the use of adapted coordinates allow us to maintain sufficient grid resolution to an artificial outer boundary location which is causally disconnected from the measurement. For the well-studied test case of the inspiral of an equal-mass nonspinning binary (evolved for more than 8 orbits before merger), we determine the phase and amplitude to numerical accuracies better than 0.010% and 0.090% during inspiral, respectively, and 0.003% and 0.153% during merger. The waveforms, including the resolved higher harmonics, are convergent and can be consistently extrapolated to $r\to \infty$ throughout the simulation, including the merger and ringdown. Ringdown frequencies for these modes (to $(\ell ,m)=(6,6)$) match perturbative calculations to within 0.01%, providing a strong confirmation that the remnant settles to a Kerr black hole with irreducible mass $M_{\mathrm{irr}}=0.884355\pm 20\times {10}^{-6}$ and spin $S_f/M_f^2=0.686923\pm 10\times {10}^{-6}$.

Figure 1

A schematic view of the $z=0$ slice of the grid setup that is used. The upper plot shows the central Cartesian grid surrounded by six “inflated-cube” patches (the four equatorial patches are shown, shaded). The boundaries of the nominal grids owned by each patch are indicated by thick lines. The lower plot shows an $r=\mathrm{\text{constant}}$ surface of the exterior patches, indicating their local coordinate lines.

Figure 2

The local radial coordinate, $R$ (solid line), can be stretched as a function of the global coordinate, $r$, in order to increase the effective size of the grid. The function specified by Eqs. (3) transitions between two almost constant radial resolutions over a distance $ϵ$ centered at $r_0$.

Figure 3

Schematic of interpolating patch boundaries in 1-dimension, assuming 4-point finite difference and interpolation stencils. Points in the nominal zones, ${\mathcal{N}}_{p,q}$, are indicated by filled circles, points in ghost zones, ${\mathcal{G}}_{p,q}$, by open squares, and points in overlap zones, ${\mathcal{O}}_{p,q}$, by closed squares. The vertical dotted line demarcates the boundary between nominal zones on each patch. Ghost points on patch $p$ are evaluated by centered interpolation operations from points in ${\mathcal{S}}_q$ on the overlapping patch (arrows) and vice versa.

Figure 4

Schematic of the causal propagation of information during the evolution. The gravitational wave source is located in the vicinity of $r=0$, with waves propagating outward at the speed of light $c=1$, and are measured at radius $r_d$ for a time of interest which would include the inspiral, merger and ringdown of a binary system. The unphysical outer boundary of the grid is located at $r_b$, which is chosen to be sufficiently far removed so that the future Cauchy horizon of the domain of dependence of the initial slice does not reach $r_b$ until the measurement is complete.

Figure 5

The dominant spherical harmonic modes of ${\psi }_4$ for $\ell =2$, 4, 6, 8, measured at $r=200M$ from the coordinate center. The plots on the right show amplitude and frequency evolution during the late-inspiral, merger and ringdown.

Figure 6

Convergence in amplitude (top) and phase (bottom) of the $(\ell ,m)=(2,2)$ mode of ${\psi }_4$ for detectors at $r=100M$ and $r=1000M$. The higher resolution difference, $h_{0.80}-h_{0.64}$, is scaled for 4th-order convergence.

Figure 7

Convergence in amplitude (top) and phase (bottom) of the $(\ell ,m)=(6,6)$ mode of ${\psi }_4$ for detector at $r=100M$ during the late through merger. The higher resolution difference, $h_{0.80}-h_{0.64}$, is scaled for 3rd-order convergence.

Figure 8

Absolute numerical error in the amplitude (top) and phase (bottom) accumulated over the course of the evolution for the highest resolution run, determined according to Eq. (35) for the pointwise differences in amplitude and phase between medium and high resolution runs. For the phase we assume the measured 8th-order convergence over the entire evolution, while for the amplitude we use 8th-order before $t\le -100$, and 4th-order thereafter (see text).

Figure 9

The ringdown frequencies for the dominant ${\psi }_4$ modes to $\ell =6$ of the merger remnant. From top to bottom, the plots show the frequencies of the $(\ell ,m)=(2,2)$, (4, 4) and (6, 6) modes, respectively, over a timescale from the (2, 2) waveform peak to $100M$ later, at which point the waveform amplitude is too small to measure an accurate frequency. The ${\psi }_4$ data measured at $r=1000M$ is plotted, in addition to the value extrapolated to $r\to \infty$, and the transformation to spheroidal harmonics. The expected quasinormal mode frequency is plotted as a dotted line, as well as a fit to the spheroidal harmonic data over the range $t\in [40M,80M]$, with error-bars determined by the standard deviation of the fit.

Figure 10

Trajectories of the two inspiralling punctures for a spinning configuration $a_1=-a_2=0.8$, with upwinded advection terms (solid lines) and without (dashed lines). In the case where no upwinding has been used, the black holes do not inspiral, due to the accumulation of numerical error.

Figure 11

Phase evolution of the $(\ell ,m)=(2,2)$ mode ${\psi }_4$ for the aligned-spin model with $a_1=-a_2=0.8$ $h=0.64M$. The 6th-order case at $h_{0.64}$ has a trajectory between the low resolution ($h_{0.80}$) and high resolution ($h_{0.64}$) 8th-order evolution.

Figure 12

Amplitude and phase evolution of the $(\ell ,m)=(2,2)$ mode of ${\psi }_4$ for the equal-mass aligned-spin model, comparing 8th-order spatial finite differencing with a scheme in which 8th-order is used only on the fine meshes surrounding the bodies, and 4th-order on the wave zone grids.

Figure 13

Differences in phase of a spinning configuration with resolution $h=0.80M$ and conformal variables $\varphi$ and $W$ against a simulation with $h=0.64M$ and conformal variable $W$. The dephasing is significant as we are on the coarse limit of resolution for this particular configuration.