On the Formation of an Eccentric Nuclear Disk following the Gravitational Recoil Kick of a Supermassive Black Hole

The anisotropic emission of gravitational waves during the merger of two supermassive black holes can result in a recoil kick of the merged remnant. We show here that eccentric nuclear disks - stellar disks of eccentric, apse-aligned orbits - can directly form as a result. An initially circular disk of stars will align orthogonal to the black hole kick direction with a distinctive 'tick-mark' eccentricity distribution and a spiral pattern in mean anomaly.


INTRODUCTION
The double nucleus of the Andromeda galaxy is well-explained by a massive disk of apsealigned stellar orbits about a supermassive black hole (Tremaine 1995), a configuration we call an eccentric nuclear disk (END). Madigan et al. (2018) showed that ENDs are stabilized by inter-orbit, gravitational torques that limit differential apsidal precession and Gruzinov et al. (2020) demonstrated that they are a thermodynamic equilibrium configuration for rotating nuclear stellar clusters. These eccentric, lopsided stellar disks may turn out to be fairly common in galactic nuclei.
Tidal disruption events (TDEs) occur when stars are unbound by a black hole's tidal gravity (Hills 1975). TDEs are observed in poststarburst galaxies at higher rates than other galaxy types (Arcavi et al. 2014;Stone & van tatsuya.akiba@colorado.edu Velzen 2016; French et al. 2016) and these starbursts are likely linked to galaxy mergers (Barnes & Hernquist 1991;Cales & Brotherton 2015;Hammerstein et al. 2021). Stellar populations in ENDs can produce TDE rates 10 3 − 10 4 times that of a spherical, isotropic configuration (Madigan et al. 2018). We propose a formation mechanism for ENDs which links galaxy mergers to enhanced TDE rates.
Binary supermassive black holes can form as a result of galaxy mergers (Begelman et al. 1980). Gravitational waves are emitted during their inspiral and eventual coalescence (Einstein 1916). In an unequal-mass binary system, the less massive black hole has a higher orbital velocity and is thus more effective at forwardbeaming its gravitational wave radiation (Wiseman 1992). This results in anisotropic gravitational waves which carry away linear momentum and the remnant receives a recoil kick (Bekenstein 1973). The same effect is produced for spinning black hole binaries (Herrmann et al. 2007). In rare cases, the kick velocity may exceed the escape velocity of the galaxy leading to black hole ejection (Gualandris & Merritt 2008), but even for smaller kick velocities, there are important dynamical consequences (Blecha & Loeb 2008). Here we show that ENDs can form as a result of a recoil kick, explaining the enhanced TDE rates we observe in post-merger galaxies.
In Section 2, we introduce a toy model to analytically explain apsidal alignment that may result from a black hole recoil kick. In Section 3, we confirm this formation mechanism using Nbody simulations and discuss the orbital properties of the resultant ENDs. We summarize our findings and conclude in Section 4.

FORMATION MECHANISM
Here we set up a toy problem that allows us to predict the orbital structure of a disk of stars following a kick to the central black hole. We consider the simple case of a circular disk formed from a self-gravitating accretion disk (Paczynski 1978;Goodman 2003) and focus on an in-plane recoil kick motivated by Bogdanović et al. (2007) (but see also Lodato & Gerosa (2013)). To measure apsidal alignment we use the eccentricity vector, where r is the radius vector, v is the velocity vector, j is the angular momentum vector, and M BH is the mass of the central black hole. The eccentricity vector points from the apoapsis to the periapsis of an orbit and its magnitude is equal to the scalar eccentricity. We consider a circular orbit in the x-y plane with positive angular momentum (|j|= j z ) centered on the black hole located at the origin. Pre-kick, the star's initial position and velocity are given by r = r cos(Ωt)x + r sin(Ωt)ŷ , where Ω ≡ ( GM BH r 3 ) 1/2 is the angular velocity, and v circ ≡ Ωr is the circular velocity at this radius. Without loss of generality, we give the black hole a kick in the positive x-direction with magnitude v BH at t = t kick . We define the initial anomaly of the star at this instant as θ ≡ Ωt kick . Immediately following the kick, the position of the star remains unchanged but in the reference frame of the black hole, its velocity becomes Setting the star's pre-kick eccentricity vector to zero, and using trigonometric identities, the post-kick eccentricity vector becomes To approximate a population of stars in an initially circular disk, we assume that the initial anomalies are evenly distributed in [0, 2π) radians. The average eccentricity vector then becomes Figure 1(a) shows a schematic diagram representing two stars in this toy problem. We show the initial circular orbit in black with the kick direction indicated by the red arrow labeled v BH . At the time of the kick, the top (purple) star's instantaneous circular velocity is pointed in the negative x-direction. In the reference frame of the black hole the star takes  on a larger velocity as indicated by the top (purple) arrow, and moves onto a larger, eccentric orbit shown by the large (purple) ellipse. The bottom (green) star's instantaneous circular velocity is in the positive x-direction in the same direction as the kick. Thus, in the frame of the black hole the star slows down, takes on a smaller velocity as shown by the bottom (green) arrow, and moves onto a smaller, eccentric orbit indicated by the small (green) ellipse. Both ellipses have eccentricity vectors that point toward the positive y-direction, in agreement with Equation 6. This is analogous to the Hohmann transfer orbit, an orbital maneuver that is often used to transfer a spacecraft between two circular orbits of different radii (Palmore 1984). Black hole kicks induce apsidal alignment in nearby circular stellar orbits and, perhaps counter-intuitively, the eccentric disk forms in such a way that the net eccentricity vector direction is orthogonal to the kick direc-tion. We confirm our analytic approach with N -body simulations in the next section.

Numerical Set-Up
We use the open-source, N -body simulation package REBOUND with the IAS15 integrator, a high accuracy non-symplectic integrator with adaptive time-stepping (Rein & Liu 2012;Rein & Spiegel 2015). We use code units of G = 1, M BH = 1, and the inner edge of the disk a in = 1 such that the period of a circular orbit at the inner edge of the disk is P (a in ) = 2π.
At t kick = 10, we kick the black hole in the plane of the disk with varying magnitudes of v BH = (0.1 − 1.0) v circ (a out ) in increments of 0.1 v circ (a out ), where v circ (a out ) is the circular velocity of stars at the outer edge of the disk (a out = 2). Scaling our code units loosely to the M31 system (Lauer et al. 1993;Bender et al. 2005) with M BH = 10 8 M and a in = 1 pc means t = 1 translates to ∼ 1490 yr. The simulation in Figure 2 then runs until ∼ 75 Myr with a black hole kick of v BH = 0.4 v circ (a out ) ≈ 186 km/s.

Formation of an END
In Figure 1(b), we show the stellar orbits projected into the x-y plane before and after the recoil kick for a simulation with v BH = 0.3 v circ (a out ). The black hole is indicated by the red cross and the kick direction with a blue arrow. Immediately following the kick, the orbits become more eccentric and the eccentricity vectors point in the positive y-direction as predicted.
We quantify the net apsidal alignment using the average unit eccentricity vector, where N bound is the total number of bound stellar particles, andê i denotes the unit eccentricity vector of the i-th particle. In Figure 2, we show the time evolution of the magnitude of the average unit eccentricity vector of the disk stars for a simulation with v BH = 0.4 v circ (a out ). The black hole receives a kick at t kick = 10. The noise floor, calculated as 1/ √ N bound , is shown in the (orange) shaded region. The average unit eccentricity vector lies below the noise floor before the kick. After the kick, it jumps to a magnitude of close to one indicating near perfect apsidal alignment. It then undergoes a coherent, large-amplitude oscillation. This behavior is well-explained by the alignment/mis-alignment of stars separated in semi-major axis space, as explained below. In Figure 3, we show a polar plot of the eccentricity vector components (e x , e y ) after the kick, with the semi-major axis indicated with the color bar. This figure explains the origin of the oscillation in apsidal alignment seen in Figure 2. Immediately following the kick (a), stellar orbits become eccentric and aligned in the positive y-direction. Stars with small semimajor axes precess much faster than those with larger semi-major axes, and misalignment (b) and re-alignment (c) between these groups occur periodically.

Eccentricity Gradient
Specifying v BH , r and its corresponding v circ , the orbital eccentricity of a star post-kick depends only on the anomaly of the star at the time of the kick, θ, which we see from taking the magnitude of the eccentricity vector in Equation 5. The semi-major axis of the post-kick orbit can be determined from the vis-viva equation, Substituting Equation 4 we obtain, In Figure 4(a), we show the post-kick eccentricities as a function of semi-major axes in an N -body simulation with v BH = 0.4 v circ (a out ). In solid (blue) and dash-dotted (red) lines, we show the analytic expectation for circular orbits that begin at the inner edge (r = 1) and the outer edge (r = 2) of the disk, respectively. The eccentricity profile has a distinctive 'tickmark' shape. There is a minimum at the radius corresponding to the initial semi-major axis of the stellar orbits and eccentricity values increase both inward or outward, such that the eccentricity gradient is negative for stars with smaller post-kick semi-major axes and positive for stars with larger post-kick semi-major axes. Whilst the orbits diffuse due to two-body relaxation, this pattern is maintained in our simulations for at least t = 50, 000.

Spiral in Mean Anomaly
We observe a clear pattern in the distribution of mean anomalies of the stars in the disk following the kick. We can see why a pattern might emerge from Figure 1(a). The top (purple) star, which takes on a large semi-major axis, is positioned such that it is located at the periapsis of its new orbit. On the other hand, the bottom (green) star, which takes on a small semi-major axis, is positioned at the apoapsis of its new orbit. This suggests that there is a pattern in anomalies following the kick that is correlated with the stars' semi-major axes. We calculate true anomaly as, We then convert to mean anomaly, M , by a standard series expansion (Smart 1953). Figure 4(b) shows results of an N -body simulation with v BH = 0.4 v circ (a out ) of mean anomalies as a function of semi-major axis immediately post-kick (orange circles) and shortly after at t = t kick + 10 (purple crosses). Once again, the expected distribution for initially circular orbits at r = 1 and r = 2 are shown by solid (blue) and dash-dotted (red) lines, respectively. A spiral structure emerges which continues to wind for several tens of orbital periods and dissipates in less than a precession period (≈ M BH /M disk P ).  In this letter we propose a new formation mechanism for eccentric nuclear disks (ENDs) via the gravitational wave recoil kick of the central black hole. We show that in-plane kicks cause surrounding circular stellar orbits to take on eccentricities orthogonal to the direction of the recoil kick. Post-kick, a transient spiral structure emerges in mean anomaly as a function of semi-major axis. The orbits display a distinctive and long-lived 'tick-mark' pattern in their eccentricity as a function of semi-major axis. On secular timescales we show that the average unit eccentricity vector undergoes a coherent oscillation due to the alignment/misalignment of stars separated in semi-major axis space (Figure 2).
In future work we will explore the stability of these aligned disks as a function of kick magnitude. Weak kicks will not result in significant apsidal alignment (Equation 6), while extreme kicks will leave few stars bound to the black hole. A recoil kick with a significant in-plane component is required to generate an END from a circular disk of stars. In future work we will consider recoil kicks with out-of-plane components, as expected for spinning Kerr black holes (Baker et al. 2008). We will explore other astrophysically-relevant stellar configurations including a scoured nuclear core (Begelman et al. 1980) and consider the effects of gas in the system. Stone & Loeb (2011 show that recoil kicks temporarily enhance TDE rates up to 0.1 yr −1 gal −1 by instantaneously filling the central black hole's loss cone (see also Komossa & Merritt (2008)). Our preliminary results suggest that if an END forms, elevated TDE rates are actually maintained for a long period of time (at least several precession periods). Such high rates could have a significant impact on the lifetime of ENDs formed via recoil kicks. However, much like bars in disk galaxies which are also composed of aligned stellar orbits, ENDs will evolve dynamically over time. They may dynamically attract new stellar orbits and grow in mass as they age. To this end, we will study the long-term evolution of ENDs formed in this manner.
Finally, this work applies to many other astrophysical scenarios, for example, when the recoil kick is due to anisotropic mass loss from a star or its non-spherical collapse to a stellar mass black hole or neutron star (which leads us to contemplate an END of planets). In the neutron star case, typical kicks are on the order ∼380 km/s (Faucher-Giguère & Kaspi 2006) which means that only the most tightly bound planets will be retained by the neutron star after an in-plane kick. White dwarf kicks on the other hand are of much smaller velocity. Kicks of order ∼1 km/s (Stone et al. 2015) may lead to END formation beyond the planetary regime at ∼ 1000 AU. In our Solar System, this relatively low mass region is populated by scattered disk objects on high eccentricity orbits. This low-kick-velocity, high orbital eccentricity configuration is one we will explore in future work.