Accretion disks around binary black holes of unequal mass: GRMHD simulations of postdecoupling and merger

Accretion disks around binary black holes of unequal mass:
GRMHD simulations of postdecoupling and merger

Roman Gold, Vasileios Paschalidis, Milton Ruiz, Stuart L. Shapiro Department of Physics, University of Illinois at Urbana-Champaign, Urbana, IL 61801
Department of Physics & Joint Space-Science Institute, University of Maryland, College Park, MD 20742
Department of Physics, Princeton University, Princeton, NJ 08544
Department of Astronomy & NCSA, University of Illinois at Urbana-Champaign, Urbana, IL 61801
   Zachariah B. Etienne and Harald P. Pfeiffer Department of Mathematics, West Virginia University, Morgantown, WV 26506
Canadian Institute for Theoretical Astrophysics, University of Toronto, Toronto, ON M5S 3H8, Canada
Canadian Institute for Advanced Research, Toronto, ON M5G 1Z8, Canada

We report results from simulations in general relativity of magnetized disks accreting onto merging black hole binaries, starting from relaxed disk initial data. The simulations feature an effective, rapid radiative cooling scheme as a limiting case of future treatments with radiative transfer. Here we evolve the systems after binary-disk decoupling through inspiral and merger, and analyze the dependence on the binary mass ratio with and . We find that the luminosity associated with local cooling is larger than the luminosity associated with matter kinetic outflows, while the electromagnetic (Poynting) luminosity associated with bulk transport of magnetic field energy is the smallest. The cooling luminosity around merger is only marginally smaller than that of a single, non-spinning black hole. Incipient jets are launched independently of the mass ratio, while the same initial disk accreting on a single non-spinning black hole does not lead to a jet, as expected. For all mass ratios we see a transient behavior in the collimated, magnetized outflows lasting after merger: the outflows become increasingly magnetically dominated and accelerated to higher velocities, boosting the Poynting luminosity. These sudden changes can alter the electromagnetic emission across the jet and potentially help distinguish mergers of black holes in AGNs from single accreting black holes based on jet morphology alone.

04.25.D-, 04.25.dg, 47.75.+f

I Introduction and Motivation

Mergers of supermassive black holes (SMBHs) and near-Eddington accretion of gas Soltan (1982) are both central ingredients in theoretical models of the assembly of the SMBH population (see e.g. Tanaka and Haiman (2009)). These models show a steadily growing consistency with data from quasar surveys Volonteri et al. (2003); Tanaka and Haiman (2009), indicating that SMBH mergers are not just a likely outcome of galaxy mergers, but necessary to explain the BH mass distribution in the Universe. When these mergers occur following the collision of galaxies, they are expected to be immersed in a magnetized plasma and surrounded by stars Begelman et al. (1980); Mayer (2013). Due to both their gravitational and electromagnetic (EM) radiation in their late evolutionary stages, such systems are unique probes of spacetime and relativistic plasmas, and are therefore interesting systems in General Relativity (GR), relativistic astrophysics and multi-messenger astronomy. Here we give only a brief introduction to this topic and refer to the more detailed discussion and additional references in Gold et al. (2014).

The first gravitational wave (GW) detectors LIG (), sensitive enough to observe GWs directly, will come online soon, but they will not be sensitive to SMBH binaries. However, pulsar timing arrays (PTAs) Hobbs et al. (2010) could detect GWs both from individual SMBH binaries Sesana et al. (2008a) and the stochastic background from unresolved SMBH binaries Sesana (2012); Sesana et al. (2008b) within this decade. Thanks to encouraging improvements in data analysis Taylor et al. (2014); Wang et al. (2014) and new discoveries in pulsar timing, PTAs may soon detect (sub-pc) SMBH binaries in the universe. Identifying EM signals from binary SMBHs will improve our understanding of the cosmic evolution of SMBHs, especially if a simultaneous GW signal from the same source is detected Sesana et al. (2011); Tanaka et al. (2012); Tanaka and Haiman (2013). Data from EM surveys like the Sloan Digital Sky Survey (SDSS) Volonteri et al. (2009); Ju et al. (2013) will probe larger redshifts with time. In addition, current and future EM detectors such as PanStarrs Kaiser et al. (2002), the LSST Abell et al. (2009) and WFIRST Green et al. (2012) will provide us with unprecedented data of transient phenomena.

In general, but especially until GWs from SMBH binaries are detected, it is crucial to have a thorough theoretical understanding of the full nonlinear dynamics and radiation properties of these systems to make the most out of EM observations. A key theoretical task is to predict observational features that distinguish accretion flows onto single versus binary SMBHs. From the point of view of transient signals, the most interesting regime is when the black hole binary merges. Modeling these systems through the late inspiral and merger phases requires a fully relativistic calculation, i.e. taking into account the dynamical black hole spacetime, as well as treating magnetized plasmas and radiative transfer relativistically.

Theoretical modeling of circumbinary disks is still in its infancy and has shown that the evolution of a circumbinary disk is roughly composed of three phases: i) the early inspiral predecoupling phase, during which the disk viscous timescale () is shorter than the gravitational wave timescale (), and the disk relaxes to a quasi-equilibrium state as the BHBH slowly inspirals; ii) the postdecoupling phase, during which and the binary decouples from the disk before the disk can relax; iii) the post-merger or re-brightening/afterglow phase during which the disk begins to refill the hollow left behind by the BHBH and accretion ramps up onto the remnant BH. Simple, analytic considerations reveal that for geometrically thick disks, binary-disk decoupling occurs during the late stages of the binary inspiral. Analytic and semi-analytic models focus on the geometrically-thin, optically thick disk case (see e.g. Haiman et al. (2009); Tanaka and Menou (2010); Shapiro (2010); Liu and Shapiro (2010); Kocsis et al. (2012); Shapiro (2013) and references therein). These 1D (semi)analytic studies make simplifying assumptions such as the adoption of an azimuthally averaged formula for the binary tidal torques, which may overestimate the tidal-torque barrier Barausse et al. (2014), and also misses non-axisymmetric features such as the formation of accretion streams. Additional features, such as the presence of an inner cavity (hollow) of lowered density near the binary, were revealed in hydrodynamic studies in Newtonian gravity in 3D Artymowicz and Lubow (1994); Cuadra et al. (2008); Roedig et al. (2011, 2012) and 2D MacFadyen and Milosavljevic (2008); D’Orazio et al. (2012); Farris et al. (2014a). This picture has been refined by the first (Newtonian) ideal magnetohydrodynamic (MHD) treatments in 3D Shi et al. (2012) and in Post-Newtonian gravity Noble et al. (2012); Zilhão et al. (2014). Infalling clouds onto and the subsequent disk formation around BHBHs has recently been studied via Newtonian smoothed particle hydrodynamic simulations in Dunhill et al. (2014). The dynamics of EM fields in force-free electrodynamics in GR, but without modeling the disk itself, has been studied in Mosta et al. (2010); Neilsen et al. (2011); Palenzuela et al. (2010a, b); Alic et al. (2012). GR evolutions of geometrically thick disks have been achieved in Bode et al. (2010); Bogdanovic et al. (2011); Bode et al. (2012); Farris et al. (2011) and only quite recently magnetohydrodynamic simulations of these systems have been performed Farris et al. (2012); Gold et al. (2014) (see also Giacomazzo et al. (2012) for treatments of BHBHs inspiraling in magnetized gaseous clouds). The overall conclusions are that both the clearing of an inner cavity and the binary disk decoupling is at best only partial: Some gas remains near the BHs all the way through inspiral and merger. Note that these features in 3D MHD turbulent disks are drastically different from findings in hydrodynamical, geometrically thin disks MacFadyen and Milosavljevic (2008), but see Zilhão et al. (2014). The transient induced on the disk following the merger has been modeled approximately by boosting either the BH or the disk, with the BH mass suddenly reduced to mimic the energy loss through GWs, see e.g. Corrales et al. (2009); Rossi et al. (2010); Anderson et al. (2010); Megevand et al. (2009a); Zanotti et al. (2010); Ponce et al. (2012a). The work presented here constitutes a substantial advancement over the above treatments, because it (i) features a magnetized disk, which is self-consistently evolved through pre- and postdecoupling, and finally through merger and (ii) takes into account the dynamical spacetime in full GR, with no artificial boundary conditions imposed to mimic the role of a BH horizon.

In this paper we focus on the postdecoupling phase including the BHBH merger, as an extension of our results in the predecoupling regime Gold et al. (2014). We consider the BHBH binary mass ratios 1:1, 1:2 and 1:4. We use relaxed matter and magnetic field initial data, starting from the predecoupling epoch obtained in Gold et al. (2014). We consider geometrically thick disks resembling slim disks Abramowicz et al. (1988); Abramowicz and Fragile (2013). In our models no physical scale is set by microphysics, resulting in the scale freedom of our results both with binary ADM mass and the disk density. However, we have in mind disks that accrete near the Eddington limit and are dominated by radiation pressure. A key purpose of this paper is to develop and test computational machinery that will be required for a GRMHD treatment of the BHBH-disk problem with full radiative transport.

The paper is structured as follows: In Sec. II we summarize the adopted techniques before reporting our results in Sec. III. Section III.1 focuses on the dependence on mass ratio, while Sec. III.2 focuses on universal features independent of the binary mass ratio. Finally we discuss astrophysical implications of our results in Sec. IV and conclude in Sec. V. Geometrized units, where , are adopted throughout unless stated otherwise.

Ii Methods

Our BHBH-disk models adopt the following set of assumptions and simplifications: a) The non-spinning black hole binaries are initially in quasi-circular orbits, b) the disk self-gravity is neglected because we assume it is small compared to the gravity of the BHBH binary, c) ideal MHD describes well the plasma in the disk, d) the same effective emissivity employed in Gold et al. (2014) [Eq. (A2)] (see also Paschalidis et al. (2011)), with the same cooling time scale, is adopted to model rapid cooling as a limiting case of realistic cooling. See Gold et al. (2014) for a more detailed discussion and motivation for these simplifications.

ii.1 Initial data

ii.1.1 Metric initial data

For the initial black hole binary spacetime geometry we adopt conformal-thin-sandwich (CTS) solutions which correspond to quasi-equilibrium black hole binaries in quasicircular orbits  Pfeiffer and York (2003); Cook and Pfeiffer (2004); Caudill et al. (2006); Baumgarte and Shapiro (2010). These solutions possess a helical Killing vector. The CTS initial data have been generated using the spectral techniques described in Pfeiffer et al. (2003) as implemented in the Spectral Einstein Code (SpEC) SpE (); Szilagyi (2014) (see also Mroue et al. (2013)). We list the initial data parameters describing our spacetimes in Tab. 1.

0.98989 0.96865 0.50000 0.50000 0.42958 0.42958
0.99097 0.85933 0.66667 0.33333 0.60192 0.27312
0.99342 0.61603 0.80000 0.20000 0.75140 0.15832
Table 1: CTS initial data parameters for the BHBH vacuum spacetime. Columns show mass ratio (), ADM mass (), ADM angular momentum (), BH irreducible masses (), and apparent horizon radii () for the two black holes. Diagnostics generating these quantities, but computed from independent vacuum test simulations, agree with these values to within one part in .

We stress that the spectrally accurate, CTS initial data for the spacetime metric have been mapped directly onto our computational grids without requiring the lower-order interpolation from the spherical auxiliary grids used in Gold et al. (2014).

ii.1.2 Matter and B-field initial data

For the magnetized fluid we use as initial data the relaxed end state obtained in Gold et al. (2014), which started from equilibrium disk models around single BHs as in Chakrabarti (1985); De Villiers et al. (2003); Farris et al. (2011, 2012) with an adiabatic index , appropriate for thermal radiation pressure-dominated disks. These disks correspond to accretion flows driven by MHD turbulence, which is self-consistently triggered by the magnetorotational instability Balbus and Hawley (1991). These solutions are interpolated onto grids appropriate for a spacetime evolution, that have the same spatial extent but contain additional levels of refinement.

ii.2 Evolution equations and methods

We use the GRMHD AMR code developed by the Illinois Numerical Relativity Group Duez et al. (2005); Etienne et al. (2010, 2012a), which adopts the Cactus/Carpet infrastructure Goodale et al. (2003); Car (a, b), and includes an effective radiative cooling scheme. This code has been extensively tested and used in the past to study numerous systems involving compact objects and/or magnetic fields (see e.g. Paschalidis et al. (2011, 2011); Etienne et al. (2012b, c); Paschalidis et al. (2012); Etienne et al. (2013); Paschalidis et al. (2013)), including black hole binaries in gaseous media Farris et al. (2010, 2011, 2012) (see Gold et al. (2014) for additional references and details).

For the metric evolution we solve the BSSN equations Shibata and Nakamura (1995); Baumgarte and Shapiro (1999) coupled to the moving-puncture gauge conditions, see Eqs. (9)-(16) in Etienne et al. (2008). For the 1:4 case we use the spatially varying damping coefficient appearing in the shift condition, as was done in the case of the LEAN-code NRAR runs Hinder et al. (2014). See Lousto and Zlochower (2011); Müller et al. (2010) for a motivation of similar strategies.

We adopt a number of diagnostics to analyze accretion disks onto binary black holes. For brevity here we describe only those diagnostics that characterize the outgoing flow of energy which include: a) the Poynting luminosity , where is the EM stress-energy tensor. measures the outgoing flux of large scale EM energy. b) The cooling luminosity , where is the cooling emissivity, is the fluid four-velocity, and is the determinant of the 3-metric. measures the total thermal emission. c) The kinetic luminosity computed for unbound () material, where is the perfect fluid stress-energy tensor. measures the outgoing flux of kinetic energy carried by unbound matter. We also compute the gas Lorentz factors (measured by a normal observer) of the plasma in the funnel region before and after the outflow settles following merger. Here, is the lapse function. For definitions of all other diagnostics we adopt in this work see Farris et al. (2012); Gold et al. (2014).

Case name  levels(BH)  levels(bh) Grid hierarchy
1:1 1:1 9 9 ,
1:2 1:2 9 10 ,
1:4 1:4 9 11 ,
0 0 6 ,
Table 2: List of grid parameters for all models. Equatorial symmetry is imposed in all cases. The computational mesh consists of three sets of nested AMR grids, one centered on each BH and one in between (with levels of refinement for all cases), with the outer boundary at M in all cases. From left to right the columns indicate the case name, mass ratio , the coarsest grid spacing , number of AMR levels around the primary (BH) and the secondary (bh), and the half length of each AMR box centered on each BH. The grid spacing of all other levels is , where is the level number such that corresponds to the coarsest level. A dash “–” indicates “not applicable”.

The grids are similar to those in Gold et al. (2014), with additional finer AMR levels centered on each BH and increased resolution in between the BHs. The higher resolution is needed for a reliable metric evolution through inspiral and merger. The regridding procedure makes use of the Cactus/Carpet interpolation routines and is identical to the procedure used in Farris et al. (2012); Gold et al. (2014).

In Tab. 2 we list the distinguishing characteristics and grid-hierarchy of the different cases we consider in this work. The labels are chosen to designate the mass ratio, e.g. the label 1:1 means mass ratio . We also evolve the same initial disk model with a single, non-spinning BH (case 0) to normalize some of our results and for comparisons to the binary cases.

We stress that the study of these systems over the duration of all epochs (from predecoupling to re-brightening) requires some of the longest GRMHD evolutions in full GR performed to date: The inner disk structure relaxes during the predecoupling epoch approximately on a viscous time scale at the inner disk edge, given by


Here is the disk (cylindrical) radius, is the disk scale height, and is the effective kinematic viscosity driven by MHD turbulence, which can be expressed as . We have assumed hydrostatic equilibrium in the vertical direction to derive an approximate relationship between and (see Shapiro and Teukolsky (1983)). The effective viscosity in our disks can be fit (approximately) to an “-disk” law for purposes of analytic estimates, and we use typical values found in our evolutions. In Gold et al. (2014) we empirically found that the relaxation at the inner edge of the same, geometrically thick disks takes , which is consistent with the order-of-magnitude estimate of Eq. (1) 111Due to the steep, inverse dependence of the viscous time scale on scale height, this time scale becomes prohibitively long for thinner disks. The subsequent inspiral occurs on a GW timescale Peters (1964)


where is the initial binary separation, and is the symmetric mass ratio. The normalization of is close to our initial binary coordinate separation and is chosen to be close to the decoupling radius as shown in Gold et al. (2014). Note that an equal-mass binary () has , while yields .

The inspiral epoch is the shortest epoch, but requires the highest resolution in order to track the inspiral reliably. The duration until the remnant disk viscously refills the inner cavity and accretes onto the merger remnant is largely determined by the radial matter distribution at decoupling and is expected to occur on a viscous time [Eq. (1)].

The disparity between the duration of the predecoupling (inspiral and merger) epoch and the dynamical (light-crossing) time scale across the horizon, where the latter determines the smallest time step in our explicit time integrations, makes these evolutions expensive and very time-consuming.

In summary, the minimum total simulation time for the computationally least expensive case (1:1) is . We have simulated the predecoupling (see Gold et al. (2014)), inspiral and merger epoch for a total of . The postdecoupling phase of our most expensive case (1:4) took approximately months (of wallclock time) to finish at a cost of CPU hours.

Iii Results

iii.1 Trend with mass ratio

In this section we discuss the dependence of our multiple diagnostics on the binary mass ratio . Results independent of are presented in Secs. III.2 and IV. For additional definitions of diagnostic quantities see Gold et al. (2014).

Case 222
0 333Near the funnel walls we find a wind-like outflow with substantially smaller velocities than the funnel regions in the binary cases.
Table 3: Columns show case name, the total accretion rate at merger normalized to the mean accretion rate for a single BH with the same cooling prescription , , , , and the 99th percentile of the gas Lorentz factors in the funnel region after the outflow settles following merger. and are computed through surface integrals over a spherical surface of coordinate radius . See also the description at the end of Sec. II. Values (except for ) are reported at merger. Based on the resolution study presented in Gold et al. (2014) we estimate the error of the quantities listed in the table to be .

iii.1.1 Evolution of the density

During the predecoupling phase, which furnishes the initial data for our postdecoupling evolutions, our disk models contain some matter in the inner cavity mainly in the form of dense accretion streams. The surface density depends on the mass ratio Gold et al. (2014). The distribution of material at decoupling determines the subsequent evolution through inspiral, merger and re-brightening. The evolution of is shown in Fig. 1, and the rest-mass density on the equatorial plane is plotted in Fig. 2. We observe accretion streams of dense gas attached to the BHs, as reported in Gold et al. (2014), during the inspiral and through merger. During inspiral there are always two diametrically opposite accretion streams as in the earlier predecoupling phase Gold et al. (2014). We therefore call this an mode by analogy to the terminology used for spiral density waves in other accretion disk studies, e.g. MacFadyen and Milosavljevic (2008). As already observed in Gold et al. (2014) there is an asymmetry between the two accretion streams (one stream is larger than the other) as departs from unity.

Spiral arms are observed throughout merger. However, the pre-merger mode ceases to dominate (see Fig. 2) for any case. Instead higher modes are excited which mostly decay over a few after merger. Following merger matter begins drifting inward for all binary cases. A “lump” feature reported in earlier work Shi et al. (2012); Noble et al. (2012); Gold et al. (2014) is also seen for the 1:1 case (see the densest regions in the upper panel in Fig. 2) but is hardly noticeable in the other cases.

Several studies Sesana et al. (2011); Farris et al. (2013); Yan et al. (2014); Farris et al. (2014b, a) investigated “mini disks” around each BH. While there is no universal definition of a mini disk, we consider a persistent mini disk to be a coherent density structure within the Hill sphere of each BH with an accretion time scale longer than a binary orbital period. For the systems under consideration here, we do not find persistent mini disks. Instead, occasionally matter piles up around the individual BHs before it is accreted. A necessary condition for persistent “mini disks” to form is that the Hill sphere (or Roche lobe) be significantly larger than the innermost stable circular orbit around the individual BH . A simple Newtonian estimate for the secondary BH yields and . These simple estimates demonstrate why we only see transient mini-disks. Given the mass ratio we expect “mini disks” around the (non-spinning) secondary BH to form at binary separations . The initial binary separations in our evolutions are and thus too close for persistent mini-disks to form. Note that for the geometrically thick disks we consider here, the decoupling separation () is given by Gold et al. (2014)


Thus, initial separations larger than the decoupling radius may be necessary for persistent mini disks to form in geometrically thick accretion flows.

Figure 1: Surface density at decoupling (yellow solid lines) and immediately after merger, time-averaged over (corresponding to the initial binary orbital period) (black solid line with circles). is normalized to – the maximum value of the surface density in the initial hydrostatic equilibrium solution used in Gold et al. (2014). In each panel the profile for the reference stationary single BH case is shown (gray solid line with dots).
Figure 2: Contours of rest-mass density normalized to the initial maximum (log scale) in the equatorial plane corresponding to the time at decoupling (left), an intermediate time prior to merger (middle), and the moment following merger (right). Upper panels: Mass ratio 1:1; . Middle panels: Mass ratio 1:2; . Lower panels: Mass ratio 1:4; . is the bolometric radiative luminosity at decoupling.
Figure 3: Total accretion rate, luminosities normalized to the single BH accretion rate, and GW strain ( polarization) as functions of time. For , and , is the retarded time. Red solid line: 1:1; yellow solid line: 1:2; black solid line: 1:4. The vertical lines indicate the corresponding merger times.
Figure 4: Contributions to the cooling luminosity from various spherical regions. Left panel: 1:1. Middle panel: 1:2. Right panel: 1:4.

iii.1.2 Accretion rates

We show accretion rates as a function of time, together with our luminosity estimates and the gravitational wave signal in Fig. 3. In the 1:1 case the total accretion rate drops as the inspiral proceeds, as expected. At merger the accretion rate is more than an order of magnitude below the value at decoupling. The 1:1 and 1:2 cases accrete at a similar rate near merger, see Tab. 3, which is not true over the entire evolution, see Fig. 3. Despite the stronger tidal torques in the 1:1 case, it can be seen from the snapshots in Fig. 2 that the density in the accretion streams, which give the dominant contribution to , reach higher values in the equal mass case. In Gold et al. (2014) we found that the 1:1 case accreted at a slightly larger rate than the 1:2 case during the predecoupling evolutions (see Tab. III in Gold et al. (2014)), but the difference in between the two cases decays close to merger.

Within the first after merger, spikes in the accretion rate from matter residing near the remnant appear in the 1:1 case. This is accompanied by a gradual rise over a longer time scale which will eventually lead to re-brightening for all mass ratio cases. In the 1:4 case the accretion rate at merger is . This value is the largest among the cases we study (see Tab. 3), but still substantially smaller than .

iii.1.3 Luminosities

We report the contribution to from within spheres of different radii and the total contribution from the disk in Fig. 4. While in the 1:1 case drops after about prior to merger, there is no such feature for the 1:2 case. As the tidal torques in the 1:1 case are stronger than the other mass ratios, this behavior is likely due to the decline in tidal heating in the 1:1 case as the inspiral proceeds. In the 1:2 case the inner cavity contribution to (see Fig. 4, middle panel) shows instead, a gradual rise, but without any prominent feature during merger. Also in the 1:4 case we observe no prominent feature in during merger in contrast to the thin disk case Shapiro (2013).

In contrast to , the binary-disk decoupling and the merger are reflected in (see also the merger aftermath feature in in Farris et al. (2012)). From the beginning of the inspiral slowly drops before rising after merger. In Fig. 3 we also plot the GW signal. One can compare the GW signal to different luminosity “light curves” and the accretion rates. For 1:1 we find a delay of between the peak in GWs and the rise in . For 1:2 this delay is significantly shorter and even shorter for 1:4, . The shortening of this delay may be explained by the fact that due to the decreasing tidal-torque barrier as decreases, there is more material near the BHs in the 1:2 and 1:4 cases, which is immediately available to be launched through the funnel. Despite the increase after merger, always remains lower than . We further report a “kinetic” luminosity associated with matter outflows, which includes only unbound material (), identical to used in Gold et al. (2014). We find in general . We give values at merger in Tab. 3. We normalize luminosities by the accretion rate of the single BH case, because it is not clear what a fair comparison to the time-dependent, instantaneous binary accretion rate would be. Note that actual efficiencies, i.e. luminosities normalized to the binary accretion rate would be much higher. The ratios range from to and are similar to the values found in Tab. 2 of McKinney and Gammie (2004), e.g. for the non-spinning case, where the same ratio is designated by . Even in the single BH case, differences with McKinney and Gammie (2004) are expected due to different disk models, our -based outflow diagnostic, their absence of radiative cooling, and possibly different locations where the ratio is evaluated (on the horizon in McKinney and Gammie (2004) vs far away from the black hole in our case).

iii.1.4 Outflows and jets

In Gold et al. (2014) we have identified collimated, magnetized outflows in the predecoupling epoch. As expected, no collimated outflows are observed for the non-spinning single BH case with the same initial disk and numerical parameters. In all binary cases we find that the incipient jets persist through merger and the immediate post-merger evolution; see Figs. 5 and 6 for all cases. The difference between single and binary cases as well as visualizations of B-field lines (Fig. 6) throughout the evolution lead us to attribute the outflows to magnetic winding and buildup of magnetic pressure above the poles of the orbiting black holes. Through accretion, B-field is accreted onto the black holes which can then tap the orbital kinetic energy as in a single spinning BH magnetic fields can tap the rotational kinetic energy of the BH, eventually giving rise to collimated, relativistic outflows, see Fig. 5.

After merger all cases reveal an increase in the Lorentz factor (measured by normal observers) of the flow in the funnel accompanied by an increase in (where is the magnetic pressure and the rest-mass density). Note that not only shows how dominant the B-field is over the inertia of the matter, but also equals the terminal Lorentz factor achieved by a steady-state, axisymmetric jet model Vlahakis and Königl (2003). Until merger we find mildly relativistic outflows with maximum at larger distance from the BHs. After merger and inside the funnel above the polar region increase to (), (), () and maximum .

Figure 5: Contours of , (log color scale) in a meridional slice, at merger (left panels), after merger (middle panels), after merger (right panels). Upper panels: 1:1. Middle panels: 1:2. Lower panels: 1:4.
Figure 6: Volume rendering of rest-mass density, normalized to its initial maximum value (see color coding), and magnetic field lines for the 1:1 case. Left panels: Halfway through the the inspiral . Right panels: Time after merger. Top panels: Global view out to . Bottom panels: Closeup view within . White field lines emanate from the BH apparent horizons. Green field lines emanate from the disk. The blue background indicates densities less than . Incipient jets are launched above each BH, merge at larger radii, and further collimate shortly after merger.

iii.2 Distinguishing pre-, postdecoupling and merger

In Gold et al. (2014) we have presented an analysis of the predecoupling phase and the relaxed disk properties. Here we report our results during the postdecoupling epoch, adopting the relaxed disk as initial data for a realistic calculation of the postdecoupling evolution.

iii.2.1 Postdecoupling

For most cases remains rather similar at merger and decoupling, at least in the bulk of the disk (see Fig. 1). This demonstrates that the response of the bulk of the disk material is slow compared to the merger time. Near the inner edge all cases show a small drift of matter inward. This behavior is expected, as the binary tidal torques decrease during inspiral while the outward angular momentum transport due to MHD turbulence persists. In all cases considered here, there is a significant reduction in density relative to the single BH case near the BHs.

For thin disks the postdecoupling evolution leads to a dimming of the source, as the binary inspirals in the nearly empty cavity while running away from the inner disk edge MacFadyen and Milosavljevic (2008); Shapiro (2013), but see Farris et al. (2014a). For thick disks the cavity contains a considerable amount of gas, which leads to a smearing of the classic decoupling picture for thin disks – in particular for mass ratios different from unity. The snapshots of the rest-mass density covering inspiral and merger, see Fig. 2, clearly demonstrate the persistence of two dense accretion streams threading the horizons and a gaseous environment, in which the BHs remain embedded in all the way through merger.

In the 1:1 case the onset of the postdecoupling evolution is signaled by a gradual decrease in the accretion rate (upper panel of Fig. 3; compare to in the predecoupling phase Gold et al. (2014)), and . In both the 1:2 and 1:4 cases there is no luminosity decrease during the inspiral of the binary. In all cases, the amplitude and frequency of the GW signal increase substantially, see Fig. 3.

iii.2.2 Early merger aftermath

A few hundred to after merger (depending on ), the Poynting luminosity (see Gold et al. (2014) for the definition) undergoes a sudden increase. For 1:1 increases by a factor of within following merger. By contrast, the cooling luminosity decreases by over a shorter time interval of before reaching a plateau at half the predecoupling value until the end of our simulation. Even so, the cooling luminosity still dominates: at the end of our simulation. We do not find a sudden, large increase in at merger as in Farris et al. (2012) (or as for thin disks Shapiro (2013)). This could be due to the different cooling emissivity or differences in the disk (adiabatic index) we employ here.

We confirm our findings in Farris et al. (2012) regarding a rapid change in the properties of the collimated, magnetized outflows in the polar regions during and shortly after merger in all cases. This can be seen, e.g., by comparing in the meridional plane at merger and shortly after merger, see Fig. 5. In all cases and for all epochs the magnetic pressure is sub-dominant relative to the rest-mass density in the bulk of the disk. The polar regions are dominated by magnetic pressure. Prior to merger we find relatively (compared to other single BH GRMHD accretion studies) small values at large distances from the BHs, but shortly after merger in a collimated cone, which quickly expands into the polar directions. This trend is shown in Fig. 5. The same collimation effect just after merger is obvious in 3D visualizations including the magnetic field lines, see Fig. 6. The simultaneous onset of the rapid change in the outflows with the increase in strongly suggests, that the increased magnetization and acceleration of the outflow is the main cause for the brightening. The collimation just after merger is observed in a similar way in all binary cases.

All currently existing ideal MHD schemes (either relativistic or Newtonian) can accurately evolve regions only up to a certain critical value of the plasma parameter . Once the critical value is reached or exceeded (typically in the low-density atmosphere) certain inequalities must be imposed to continue the simulations, with their impact designed to be minimal. Based on previous results and tests with our code, we are convinced that the postmerger increase in the magnetization in the funnel is robust, but terminal values of in those regions may not be reliable.

The Blandford-Payne mechanism Blandford and Payne (1982) is probably not the cause for this transient behavior because of the special conditions under which the mechanism operates.

However, it is natural to attribute part of the increase in to the Blandford-Znajek Blandford and Znajek (1977) effect (see also Neilsen et al. (2011); Palenzuela et al. (2010b)), because all of our BH remnants are spinning with the funnel area above the remnant BH poles being nearly force-free. The BZ solution is known to describe the force-free regions in the funnel of magnetized, geometrically thick disks accreting onto single spinning BHs McKinney and Gammie (2004).

As in Farris et al. (2012), we can see the disk beginning to drift inwards towards the remnant BH by comparing at different times, see Fig. 1. By the time the bulk of the material will reach the remnant BH the system will likely undergo a re-brightening Milosavljevic and Phinney (2005).

Due to asymmetries in the gravitational radiation and momentum conservation the remnant BHs in the 1:2 and 1:4 cases, experience a recoil or “kick” (see Baker et al. (2006); Megevand et al. (2009b) and references therein). However, for initially non-spinning BHs the recoil velocity has a maximum value of km/s Gonzalez et al. (2007), and hence it is small compared to other characteristic velocities in the system (see also Bogdanovic et al. (2007); Dotti et al. (2010)) unlike Ponce et al. (2012b). Therefore, the remnant BH recoil in our simulations does not have any significant impact on the accretion flow.

Iv Astrophysical implications

Observational evidence for a SMBH binary near the postdecoupling regime remains elusive. There are two possible explanations: (I) There are too few sources. (II) Due to various reasons, identifying SMBH binaries is difficult. Point (II) includes the possible misinterpretation of a binary AGN as a single BH AGN. This confusion can arise because the second BH may not be an “active” AGN or may not be massive enough to alter the outgoing radiation at an observable level.

In our models several diagnostics reveal differences between the single and binary models. In the single BH more matter resides closer to the BH, which is visible in comparisons, see Fig. 1. In the binary case is reduced compared to the single BH case for the same disk. The decrease over time in and in the binary system after decoupling [ prior to merger], which we report here, is another signature which is absent in the single BH case. This reduction is expected to last until re-brightening ( after merger).

The magnetically-driven transient during merger and the resulting acceleration of matter along the polar regions is characteristic of the binary merger in all cases and is not observed (nor expected) in the single BH case. Instead, in the single non-spinning BH we observe little to no polar outflows. There are outflows from spinning BHs, but these are not likely to exhibit a one-time, dramatic transient behavior such as the one exhibited in the binary case. Based on our findings, the strongest evidence for the presence of a BH binary is the transient in the magnetized outflows/jets during and shortly after merger. The increased magnetic field strength and outflow velocity will likely lead to enhancements of the radio emission from the jet and perhaps the X-rays and additional brightening due to relativistic beaming.

While there are known effects that can cause a single BH accretion flow to flare (recurringly) – such as hotspots – the flare in the Poynting luminosity near a binary merger is a one-time event.

Our results motivate a search for binary SMBH candidates based on jet morphology. Even a merger event long in the past could be identified by a change in the collimation from the foot of the jet towards its head. In fact, observations in time of such jets could reveal that the stronger emission is propagating outwards. This is similar to the interpretation that X-shaped radio sources originate from a sudden spin-flip following a past BHBH merger. The difference here is that there is a transient feature in the jet even in the absence of a spin flip.

The effective temperature, magnetic energy density, and characteristic cyclotron frequencies during the inspiral phase remain similar to the values we reported in Gold et al. (2014) for the predecoupling phase:


However, a few hundred after merger we find an increase in the magnetic energy density in the funnel region by a factor of . This effect shifts the cutoff frequency of synchrotron emission, which may arise in these systems from the presence of relativistic electrons, towards higher frequencies; see e.g. Problem 4.2 in Rybicki and Lightman (1986). Therefore, a one-time frequency shift in the synchrotron emission could be detected in radio surveys Kaplan et al. (2011); O’Shaughnessy et al. (2011) and may reveal the presence of a BHBH merger.

Also blazar systems, similar to the SMBH binary candidate OJ-287 Sillanpaa et al. (1988); Valtonen and Ciprini (2011), constitute a promising class of systems where binaries might be identified through EM observations Decarli et al. (2013) based on variability studies Tanaka (2013).

The dimming of total luminosity observed in 1:1 signals the onset of the postdecoupling epoch and serves as a precursor for the upcoming merger with a lead time of . Such a one-time dimming is unique to the equal-mass SMBH binary and does not seem to occur either in a single BH accretion flow or for mass ratios significantly different from unity. Thus, a near-equal-mass binary AGN can potentially be distinguished from both a single BH AGN prior to merger and binary AGNs with mass ratios deviating from unity, even in the absence of a “sudden” EM feature during merger.

Force-free simulations in full GR have suggested that dual jets from BHBH systems may be detectable Palenzuela et al. (2010b). If such dual jets were detected they would be strong evidence for the presence of an accretion disk onto a BHBH. However, the force-free simulations of Moesta et al. (2012) argue that the power in these dual jets is only a small fraction of an otherwise more isotropic emission, suggesting that dual jets are likely not detectable. Our GRMHD simulations show, that the individual jets launched by each BH merge into one common jet structure (at least during the late inspiral). Therefore, we conclude that dual jets are unlikely to be detected from magnetized accretion disk-BHBH systems in which the BHs are slowly spinning and the disk orbital angular momentum is aligned with binary orbital angular momentum, see also O’Shaughnessy et al. (2011).

The merger of the two BHs poses a major change to the whole system. The results we report here are attributed to the effects immediately () following the merger event. By “immediate” we mean over time scales which involve the dynamics of material within a few M near the BH left over from the merger. This is in contrast to the independent effects of material from the bulk of the disk beginning to fall in on a (much longer) viscous time scale () after the merger has occurred.

V Conclusions

We have presented results from our follow-up study of Gold et al. (2014) by evolving relaxed GRMHD accretion flows through the binary inspiral and merger phases in full GR.

The key differences between signatures arising from the BHBH-disk at decoupling as compared to those near merger are:

  1. The mildly relativistic dual jets observed near decoupling and prior to merger, coalesce and form one common jet. The common jet further collimates near merger (see Figs. 5 and 6), while the outflow Lorentz factors are boosted following merger by .

  2. The Poynting luminosity increases shortly after merger by a factor of , and its value is sustained until the end of our simulations (see Fig. 3).

  3. The kinetic luminosity exhibits a large peak near merger (see Fig. 3) whose height is times larger than the values prior to merger.

  4. The cooling luminosity is largely insensitive to the dynamics during postdecoupling and merger (see Figs. 3 and 4). We find a decrease during postdecoupling only in the 1:1 case.

For decreasing mass ratio (1:1, vs. 1:2, vs. 1:4), the key trends of the BHBH-disk systems are:

  1. The increase in the Poynting and kinetic luminosities after merger begins earlier as decreases. For 1:4 the kinetic luminosity peaks almost simultaneously with the GW burst at merger.

  2. The boost in the Poynting luminosity after merger decreases with decreasing . This is most likely due to the fact that the spin of the remnant BH decreases with .

  3. There is significant variability in the accretion rate and the cooling luminosity in the 1:4 case, which is not observed for the other cases.

  4. The non-axisymmetric “lump” feature becomes weaker as decreases.

We find little decrease in nearly all luminosity diagnostics after decoupling, indicating that such sources may be bright. Aftermath EM signatures are more prominent than precursor EM signals. Generally, the dependence of EM signatures (the increase in the Poynting luminosity and its time lag after merger) on mass ratio is stronger after merger than before merger or in the predecoupling epoch. A robust acceleration and boost in magnetic energy density of the outflowing material is observed, which is an excellent candidate for a clear and pronounced, one-time EM signature for merging SMBH binaries. This transient and the reported features in the light curves are unlikely to occur in single BH disk systems.

In the future we intend to include realistic radiation processes and radiative transport, and refine our results by calculating an EM spectrum, in order to identify distinguishing features between single BH and binary BH AGN.

It is a pleasure to thank the Illinois Relativity group REU team (Brian R. Taylor, Lingyi Kong, Abid Khan, Sean E. Connelly, Albert Kim and Francis J. Walsh) for assistance in creating Fig. 6. We thank Charles Gammie, Brian Farris, Jonathan McKinney, and Alberto Sesana for useful discussions. This paper was supported in part by NSF Grants No. PHY-0963136 and No. PHY-1300903 as well as NASA Grants NNX11AE11G and NNX13AH44G at the University of Illinois at Urbana-Champaign. V.P. gratefully acknowledges support from a Fortner Fellowship at UIUC. H. P. P. acknowledges support by NSERC of Canada, the Canada Chairs Program, and the Canadian Institute for Advanced Research. The metric initial data was computed on the GPC supercomputer at the SciNet HPC Consortium Loken et al. (2010). SciNet is funded by: the Canada Foundation for Innovation under the auspices of Compute Canada; the Government of Ontario; Ontario Research Fund–Research Excellence; and the University of Toronto. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by NSF Grant number OCI-1053575. This research is part of the Blue Waters sustained-petascale computing project, which is supported by the National Science Foundation (Award No. OCI 07-25070) and the state of Illinois. Blue Waters is a joint effort of the University of Illinois at Urbana-Champaign and its National Center for Supercomputing Applications.


Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description