A restless supermassive black hole in the galaxy J0437+2456
Dominic W. Pesce, Anil C. Seth, Jenny E. Greene, James A. Braatz, James J. Condon, Brian R. Kent, Davor Krajnovi?
DDraft version January 21, 2021
Typeset using L A TEX twocolumn style in AASTeX63
A restless supermassive black hole in the galaxy J0437+2456
Dominic W. Pesce,
1, 2
Anil C. Seth, Jenny E. Greene, James A. Braatz, James J. Condon, Brian R. Kent, and Davor Krajnovi´c Center for Astrophysics | Harvard & Smithsonian, 60 Garden Street, Cambridge, MA 02138, USA Black Hole Initiative at Harvard University, 20 Garden Street, Cambridge, MA 02138, USA Department of Physics and Astronomy, University of Utah, 115 South 1400 East, Salt Lake City, UT 84112, USA Department of Astrophysics, Princeton University, Princeton, NJ, USA National Radio Astronomy Observatory, 520 Edgemont Road, Charlottesville, VA 22903, USA Leibniz-Institut f¨ur Astrophysik Potsdam (AIP), An der Sternwarte 16, D-14482 Potsdam, Germany
ABSTRACTWe present the results from an observing campaign to confirm the peculiar motion of the super-massive black hole (SMBH) in J0437+2456 first reported in Pesce et al. (2018). Deep observationswith the Arecibo Observatory have yielded a detection of neutral hydrogen (H i ) emission, from whichwe measure a recession velocity of 4910 km s − for the galaxy as a whole. We have also obtainednear-infrared integral field spectroscopic observations of the galactic nucleus with the Gemini Northtelescope, yielding spatially resolved stellar and gas kinematics with a central velocity at the inner-most radii (0 . (cid:48)(cid:48) ≈
34 pc) of 4860 km s − . Both measurements differ significantly from the ∼ − H O megamaser velocity of the SMBH, supporting the prior indications of a velocity offset between theSMBH and its host galaxy. However, the two measurements also differ significantly from one another,and the galaxy as a whole exhibits a complex velocity structure that implies the system has recentlybeen dynamically disturbed. These results make it clear that the SMBH is not at rest with respectto the systemic velocity of the galaxy, though the specific nature of the mobile SMBH – i.e., whetherit traces an ongoing galaxy merger, a binary black hole system, or a gravitational wave recoil event –remains unclear. INTRODUCTIONGiven that nearly all galaxies are thought to harborcentral supermassive black holes (SMBHs; Magorrianet al. 1998), interactions between SMBHs have long beenrecognized as a natural and perhaps inevitable byprod-uct of galaxy mergers. The two primary dynamicalstates that result from such interactions are SMBH bi-naries (Begelman et al. 1980; Roos 1981) and gravita-tional recoil events (Fitchett 1983; Redmount & Rees1989), both of which predict substantial nonequilibrium(“peculiar”) motion of the SMBH with respect to itssurrounding environment. Yet despite much theoreticalattention and observational effort, direct dynamical ev-idence for SMBH peculiar motion has remained elusive(see, e.g., Eracleous et al. 2012; Popovi´c 2012; Komossa& Zensus 2016; Barack et al. 2019). In the absenceof recent interactions with comparable-mass objects, anSMBH is expected to be in kinetic equilibrium with its
Corresponding author: Dominic W. [email protected] surrounding environment (Merritt et al. 2007); for mostSMBHs, the equilibrium velocity is (cid:28) − with re-spect to the system barycenter.Pesce et al. (2018, hereafter P18) presented a tech-nique for using H O megamasers to measure SMBHpeculiar motions. The key idea is that masers resid-ing in the accretion disks around SMBHs (on scalesof ∼ (cid:46)
10 km s − ) measurements of the SMBH’sline-of-sight velocity (e.g., Miyoshi et al. 1995; Kuo et al.2011; Gao et al. 2017). P18 compared the maser-derivedSMBH velocity measurements for 10 systems with inde-pendent estimates of their host galaxy velocities to con-strain relative motions. One galaxy from the P18 sam-ple – SDSS J043703.67+245606.8, hereafter J0437+2456– showed a statistically significant ( > σ ) difference be-tween the SMBH and host galaxy line-of-sight velocities;P18 thus identified J0437+2456 as a promising candi-date for hosting either a recoiling or binary SMBH. a r X i v : . [ a s t r o - ph . GA ] J a n J0437+2456 is an approximately Sb-type spiral galaxylocated at a distance of ∼
70 Mpc (Greene et al. 2016;Pjanka et al. 2017). As measured by the Sloan DigitalSky Survey (SDSS) , J0437+2456 has an r -band abso-lute AB magnitude of M r = − .
37 and an estimatedstellar mass of 7 . × M (cid:12) . The megamaser systemin J0437+2456 was mapped by Gao et al. (2017, here-after G17), who also modeled the maser rotation curveand determined an SMBH velocity of 4818 ± . − .P18 used a SDSS spectrum to measure the recession ve-locity of J0437+2456 to be 4887 . ± . − . Theapparent 69 . ± . − blueshift of the SMBH withrespect to its host galaxy constitutes the putative pecu-liar motion. However, given the strong prior expectationfor zero peculiar motion and the possibility that system-atic effects such as SDSS fiber misalignment could plau-sibly account for a large fraction of the observed veloc-ity difference, P18 cautioned that the peculiar motionmeasurement should be regarded as tentative pendingcorroborating observations.In this paper we present the results from a followupobserving campaign to confirm the peculiar motion ofthe SMBH in J0437+2456. This paper is organized asfollows. In Section 2 we describe our observations andsubsequent data reduction procedures, and in Section 3we detail the velocity measurements made using thesedata. We discuss the results in Section 4, and we sum-marize and conclude in Section 5. Unless noted other-wise, all velocities quoted in this paper use the opticalconvention in the barycentric reference frame, and weassume a distance to J0437+2456 of 70 Mpc. OBSERVATIONS AND DATA REDUCTIONIn quiescent systems, H i provides an appealing reces-sion velocity tracer because it follows the global dynam-ics of the galaxy well outside of the SMBH sphere ofinfluence and does not suffer from reddening or extinc-tion. P18 targeted J0437+2456 with the Very LargeArray (VLA) to observe neutral hydrogen (H i ), but noemission was detected within the six-hour integrationtime. We have obtained followup H i observations ofJ0437+2456 using the Arecibo Observatory, which ismuch more sensitive than the VLA to low surface bright-ness emission but which lacks the ability to spatially re-solve the gas distribution (see Figure 1). Our Areciboobservations are presented in Section 2.1.Lacking H i data, P18 measured the recession veloc-ity for J0437+2456 using an SDSS spectrum. Like theArecibo spectrum, the kinematics contributing to the Here we quote the quantities compiled in the NASA-Sloan Atlas,http://nsatlas.org/. D e c li n a t i o n ( J ) Arecibo beamNIFS field of view
Figure 1.
False-color image of J0437+2456 made by com-bining the i -, r -, and g -band observations from the SDSSLegacy Survey (York et al. 2000), with the 2.4-arcminuteArecibo beam and the 3 (cid:48)(cid:48) × (cid:48)(cid:48) NIFS field of view overplottedin red and blue, respectively.
SDSS spectrum are spatially unresolved within the 3-arcsecond aperture of the optical fiber used to transportlight from the focal plane to the spectrograph (Gunnet al. 2006). However, because the aperture is smallerthan the region containing the emitting material, theSDSS measurement is subject to an unknown amountof systematic uncertainty associated with the relativeplacement of the fiber center and the galactic nucleus.We have thus obtained followup high-resolution integralfield spectra taken using the Gemini North telescope,which are able to spatially resolve the nuclear kinemat-ics. Additionally, dust absorption should be weaker inthe NIFS near-infrared waveband than at the SDSS opti-cal wavelengths, so any velocity errors caused by patchydust absorption will be smaller. Our Gemini observa-tions are presented in Section 2.2.2.1.
Arecibo data
We performed H i spectral-line observations ofJ0437+2456 over 6 nights using the Arecibo Observa-tory L-wide receiver. The observations were position-switched, with 5 minutes on and 5 minutes off source atmatched elevation. We used the Wideband Arecibo Pul-sar Processor (WAPP) spectrometer backend in single-polarization, 9-level autocorrelation mode with 4096channels spanning the bandwidth 1384.5—1409.5 MHz(i.e., ± − centered on the H i line). We usedtwo such boards, one per polarization. Calibration restless supermassive black hole in the galaxy J0437+2456 Table 1.
Arecibo observation details
Integration T sys GainDate (min.) (K) (K Jy − )2019 Jan 22 60 26.6 6.92019 Jan 23 65 26.8 6.92019 Jan 24 65 26.9 6.82019 Feb 11 60 27.1 7.02019 Feb 14 50 27.4 6.92019 Feb 17 50 27.0 6.9 Note —Observing dates, on-source integrationtimes, system temperatures, and gains for theArecibo observations. diodes were observed at the end of every scan to de-termine the flux density scale.Table 1 lists the on-source integration times for eachof the 6 nights. With a declination of +25 degrees,J0437+2456 passes through the Arecibo observing win-dow for ∼ . We firstconverted the flux scale from K to Jy using the regulargain curve monitoring scans performed by the obser-vatory (see Table 1). Each spectral scan was Hanningsmoothed to mitigate ringing, and a fourth-order poly-nomial fit to the emission-free regions of the spectrumwas subtracted off to remove low-frequency baseline rip-ples. We then combined each scan and both polariza-tions using an RMS-weighted average.Figure 2 shows the Arecibo spectrum, in which westrongly detect H i emission around the expected reces-sion velocity range. The spectrum peaks at ∼ . ∼ − , consis-tent with the non-detection reported in P18. Assumingthe H i is optically thin, the total H i mass is given by(Haynes et al. 2011; Condon & Ransom 2016) M H i = (cid:0) . × M (cid:12) (cid:1) (cid:18) D Mpc (cid:19) (cid:18) (cid:82) S ν ( v ) dv Jy km s − (cid:19) , (1)where D is the distance to the galaxy and S ν ( v ) is theflux density as a function of velocity v . For a distance of http://outreach.naic.edu/ao/scientist-user-portal/astronomy/IDL-Routines/Download-AO-IDL ∼ phil/sysperf/sysperfbymon.html )0.20.00.20.40.60.81.0 F l u x d e n s i t y ( m J y ) original spectrumsmoothed x4modelmodel component Figure 2. V velocity (see Section 3.1) is markedby a vertical green line, and the peak-to-peak velocity rangeis marked by a horizontal green line.
70 Mpc to J0437+2456, we estimate M H i ≈ . × M (cid:12) (see also Section 3).2.2. Gemini data
We obtained integral field spectra of a 3 (cid:48)(cid:48) × (cid:48)(cid:48) re-gion centered on the nucleus of J0437+2456 using theGemini North Near-Infrared Integral Field Spectrome-ter (NIFS) on 2018 November 21 in natural seeing mode.The spectrometer grating was set for K-band, with acentral wavelength of 2.18 µ m and spanning the range1.99–2.41 µ m. Nine 500-second exposures were taken,with five dithered exposures on-source and four offset toa blank sky location for subtraction. The observationswere performed at airmasses of 1.2–1.6 and seeing con-ditions corresponding to a zenith-corrected point spreadfunction of ∼ K -band.The NIFS data were reduced using the Gemini ver-sion 1.13 IRAF packages, with slight modifications toenable error array propagation and cube combinationas described in Ahn et al. (2018). The resulting finaldata cube has a central signal-to-noise of ∼
28, droppingto ∼
10 at 0 . (cid:48)(cid:48) emission lines are seen, and their veloc-ities are described in more detail in Section 3.2. D e c o ff s e t ( a r c s e c o n d s ) stars
100 pc 1.0 0.5 0.0 0.5 1.0RA offset (arcseconds) H V e l o c i t y ( k m s ) Figure 3.
Velocity maps derived from the Gemini NIFS data within the central 2 (cid:48)(cid:48) × (cid:48)(cid:48) region. Left : Velocity map for thestellar component, using Voronoi binning such that each bin has a signal-to-noise ratio of at least 25.
Right : H velocity map.In both panels, mean continuum contours are overplotted at 5%, 10%, 20%, 40%, and 80% of the peak value.3. ANALYSISIn this section we describe the analysis proceduresused to measure velocities from the Arecibo spectrum(Section 3.1) and the Gemini spectra (Section 3.2).3.1.
Neutral hydrogen spectral decomposition andvelocity measurements
Instead of the classic symmetric “double-horn” H i pro-file (Roberts 1978), J0437+2456 shows a more unusualtriple-peaked and asymmetric spectral structure. Be-cause our observations do not spatially resolve the H i kinematics, the association of individual spectral prop-erties with distinct dynamical components is ambiguous.While it is clear that a single double-horn componentcannot describe the observed spectral profile, there area variety of more complicated models that could poten-tially do so adequately. In Appendix A we explore threeplausible model extensions, from which we conclude thatthe observed spectral structure is most conservativelyand satisfactorily modeled using a sum of Gaussian com-ponents. Using the dynesty nested sampling routine(Speagle 2020) to explore the posterior distribution, wefind that N = 3 Gaussian components are sufficient tocapture the spectral structure and achieve a reduced- χ of ∼ i spectrum to make a measure-ment of V , defined to be to the midpoint between the two points on the profile that rise to 20% of the peak am-plitude (see, e.g., Fouque et al. 1990). V provides an es-timate of the galaxy recession velocity, and we find V =4909 . ± . − . For the associated width of the pro-file, W , we find W = 326 . ± . − , and for thetotal mass of H i we find M H i = (1 . ± . × M (cid:12) .3.2. Systemic velocities of the stellar and H components The NIFS K -band spectra show both the strong stel-lar absorption lines of CO at ∼ µ m and molecularhydrogen emission lines, including the strong H µ m.Stellar kinematics were derived by first Voronoi bin-ning the data cube to S/N ≥
25 (Cappellari & Copin2003), and then fitting the data with pPXF (Cappel-lari & Emsellem 2004) using high resolution stellar tem-plates from Wallace & Hinkle (1996). The resultingradial velocity map can be seen in Figure 3. Errorson individual bins are determined through Monte Carlosimulations and range from 5–10 km s − . The velocitymap was then analyzed using the Kinemetry code (Kra-jnovi´c et al. 2006) to determine the barycentric systemicvelocity as a function of radius from the observed photo-center; we note that the appearance of the galaxy in theNIFS data cubes is very symmetric. At the smallest ra-dius (0 . (cid:48)(cid:48) ± − .Kinemetry reveals that the velocity steadily declineswith radius – at 0 . (cid:48)(cid:48) ± − . These mea-surements are shown as orange dots in Figure 5. We notethat the quoted errors are the formal errors produced by restless supermassive black hole in the galaxy J0437+2456 Table 2.
Velocity measurements for different components in J0437+2456
Source of velocity Velocity (km s − ) Uncertainty (km s − ) Spatial scale (pc) Reference Maser rotation curve 4818.0 10.5 0.2 G17Maser rotation curve re-analysis 4809.3 10.0 0.2 this workNIFS stellar 4857–4844 ∼ ∼ ∼ . × Pesce et al. (2018)SDSS spectrum, stellar 4921.4 19.1 Pesce et al. (2018)SDSS spectrum, average 4887.6 7.1 Pesce et al. (2018)H i spectrum, first component 4774.6 3.2 4 . × this workH i spectrum, second component 4870.2 0.8 this workH i spectrum, third component 4989.4 1.9 this workH i spectrum, V Note —Velocity measurements considered in this paper and the spatial scales on which they are measured, assuming a distanceof 70 Mpc to J0437+2456. For the NIFS velocities, we quote the range of values corresponding to the systemic velocities at theinnermost and outermost annuli in which measurements were made; note that for the NIFS stellar measurements, the systemicvelocity measured from the outer annulus is smaller than that measured from the inner annulus. For the NIFS stellar and H velocity measurements, the uncertainties are dominated by an overall calibration systematic of ∼ − . the Kinemetry code and are smaller than the systematicerrors in our velocities discussed below.To check the veracity of the systemic velocity shiftwith radius, we also binned the spectra in circular an-nular bins, and we ran pPXF on the resulting spectra.The mean velocity of the innermost spectrum ( < . (cid:48)(cid:48)
1) is4865.4 ± − , while the annulus between 0 . (cid:48)(cid:48) . (cid:48)(cid:48) ± − . Thus it seemsquite clear that there is indeed a blueshift in the sys-temic velocity of ∼
15 km s − between the center of thegalaxy and the galaxy at radii of a few hundred parsecs.These “integrated light” measurements are shown as redpoints in Figure 5.We also determine the kinematics of H >
10 times the surrounding noise level.The result is shown in the right panel of Figure 3. Clearrotation with the same position angle as the stellar kine-matics is visible. However, Kinemetry reveals that whilethe systemic velocities of the stellar and H componentsare similar at the innermost radii, at larger radii the H systemic velocity is actually redshifted (not blueshiftedlike the stellar kinematics), as shown by the blue pointsin Figure 5. We note that the wavelength solution of the NIFS datawas verified through fitting of sky lines in the spectra; astandard deviation of 1.1 km/s from the mean velocitywas found from pixel to pixel, and an overall offset of-0.8 km/s was found. Thus the systematic errors onour velocity measurements are < The velocity of the SMBH
The velocity of the SMBH in J0437+2456 has previ-ously been measured to be 4818 ± . − by G17,who used VLBI measurements of H O megamasers inthe SMBH accretion disk to map out its rotation curvewell within the gravitational sphere of influence. Byfitting this rotation curve with a thin-disk Keplerianmodel, G17 were able to measure both the mass and ve-locity of the central SMBH. In this section, we re-analyzethe same VLBI dataset using an updated maser diskmodel, which relaxes several of the assumptions madeby G17 and thus permits an improved assessment of theassociated velocity uncertainty.The VLBI observations carried out by G17 resultedin position and velocity measurements for each of thedetected maser features, or “spots.” G17 fit their ro-tation curve using a two-step procedure, in which theentire VLBI map is first rotated and shifted such thatthe blueshifted and redshifted masers lie on the hori-zontal axis, and then the maser spot velocities are fitas a function of their measured one-dimensional posi-tions along this horizontal axis. The G17 rotation curvemodel contains three free parameters: the SMBH mass,the one-dimensional SMBH position along the horizon-tal axis, and the SMBH’s line-of-sight velocity.For the present analysis we employ a modified versionof the maser disk model described in Pesce et al. (2020)to fit the J0437+2456 VLBI data. The primary modifi-cation is the removal of acceleration measurements fromthe model likelihood, as the available VLBI dataset doesnot contain any such acceleration measurements. Ourfitting approach differs from that of G17 in several re-spects:1. We take the maser velocities, rather than their po-sitions, to be the “independent” quantities; i.e.,the model is essentially r ( v ) rather than v ( r ). Thisstrategy leverages the fact that the individual ve-locity measurements are uncertain at a level com-parable to a spectral channel width ( ∼ − )and therefore much smaller than the orbital ve-locities of several hundred km s − , while the posi-tion uncertainties are comparatively large ( ∼ x and y maser position measurements along-side the disk model parameters. These errorfloor parameters describe the additional uncer-tainty that it would be necessary to add into themeasurements to ensure that the data are consis-tent with the model; i.e., these parameters enforcea final reduced- χ value that is consistent withunity.Detailed descriptions of the model parameters, likeli-hood, and fitting procedure are provided in Pesce et al.(2020). Following G17, we fit only to the redshifted andblueshifted maser features because there are no avail-able acceleration measurements to constrain the sys-temic maser feature orbital radii. The final model con- Table 3.
Results from re-analysis of maser VLBI dataParameter Units Prior Best-fit value v km s − U (4500 , . ± . M M (cid:12) U (0 ,
30) 2 . ± . x mas U ( − . , .
5) 0 . ± . y mas U ( − . , .
5) 0 . ± . i deg. U (70 , deg. U (0 , . ± . deg. mas − U ( − , ± σ x µ as U (0 , ± . σ y µ as U (0 , < Note —Results from fitting a thin Keplerian disk model tothe J0437+2456 VLBI maser dataset from G17, as de-scribed in Section 3.3 and shown in Figure 4. The fittedmodel parameters are the SMBH velocity v , the SMBHmass M , the SMBH position ( x , y ), the disk inclination i ,the disk position angle Ω , a first-order warp in the disk po-sition angle with radius Ω , and two error floor parameters σ x and σ y for the maser x - and y -position measurements,respectively. The notation U ( a, b ) denotes a uniform distri-bution on the range ( a, b ). For most parameters we reportthe posterior mean and standard deviation, though we notethat the σ y parameter has a best-fit value that is consistentwith zero and so we report the 95% upper limit instead.For the i parameter, the posterior distribution matches theprior distribution, and so we do not report constraints onthis parameter. A more detailed description of the variousmodel parameters can be found in Pesce et al. (2020). tains 9 parameters, which are listed in Table 3 alongwith their priors and best-fit values.Our best-fit rotation curve and disk model are shownin Figure 4, from which we determine the SMBH velocityto be 4809 . ± . − . We find that the uncertaintyin the derived velocity matches well with the result fromG17, though our best-fit velocity itself is approximately9 km s − smaller. Because our disk model relies on fewerassumptions than that employed in G17, we hereafteradopt 4809 . ± . − as the velocity measurementof the SMBH in J0437+2456. DISCUSSIONThe recession velocity measurements considered inthis paper are listed in Table 2, and they are plottedagainst spatial scale in Figure 5. We find that all veloc-ity measurements fall within a ∼
100 km s − range span-ning ∼ − , and there is a general trend formeasurements made at larger spatial scales to recoverlarger recession velocities. In this section we discuss the restless supermassive black hole in the galaxy J0437+2456
200 250 300 350 400 450 500 550| v v | (km s )0.00.20.40.6 R a d i u s ( m a s ) v (km s )0.000.010.020.030.04 P r o b a b ili t y d e n s i t y ( k m s ) Gao et al. (2017)this work x (mas)0.40.20.00.20.4 y ( m a s ) Figure 4.
Results from fitting a thin Keplerian disk model to the J0437+2456 maser measurements from G17; the best-fitmodel parameters are listed in Table 3. The top left panel shows the on-sky projected radial separation from the SMBH versusorbital velocity for each of the maser spots, with the best-fit rotation curve plotted in black and 200 draws from the posteriordistribution plotted in gray. The points corresponding to individual maser features have been colored by velocity group, withblue points denoting blueshifted maser features and red points denoting redshifted maser features. The bottom left panelshows the posterior distribution for the SMBH velocity that we obtain from our fitting procedure (blue histogram) along witha Gaussian distribution with the mean and standard deviation reported in G17 (black line). The right panel shows the VLBImap of the maser system, with the best-fit warped disk midplane plotted as a dashed line and 200 draws from the posteriordistribution plotted as solid gray lines; the 1 σ and 2 σ contours for the SMBH location are shown as thick and thin black ellipses,respectively. various measurements and consider some possible inter-pretations.4.1. Velocity measurements in J0437+2456
The largest spatial scales are probed by the H i emis-sion, which traces gas throughout the galaxy and out tothe edge of the Arecibo beam (roughly ∼
50 kpc across).J0437+2456’s H i profile is atypical in that it shows threeprominent spectral peaks rather than the usual two thatare expected for a simply-rotating system. Similar pro-files have been classified as “anomalous” by previousauthors (e.g., UGC 2889 in Courtois et al. 2009), andthey are often attributed to spatial blending of galaxypairs in single-dish spectra, such as in the case of NGC876 and NGC 877 (Bottinelli et al. 1982; Lee-Waddellet al. 2014). However, for J0437+2456 we see neither ev- idence for a companion galaxy within the Arecibo beam(see Figure 1) nor obvious signs of morphological distur-bance in Hubble Space Telescope (HST) images (Pjankaet al. 2017). Nevertheless, the measured H i central ve-locity of V = 4910 km s − is in agreement with theSDSS stellar velocity measured by P18, supporting thenotion that both measurements trace the recession ve-locity of J0437+2456. Furthermore, the H i central ve-locity is in ∼ σ disagreement with the SMBH velocityas measured from the maser rotation curve (Section 3.3),indicating that the black hole is blueshifted by roughly100 km s − with respect to the galaxy’s recession veloc-ity.The NIFS measurements probe spatial scales of ∼ components agrees with that of the maser disk (G17), Spatial scale (pc)48004850490049505000 V e l o c i t y ( k m s ) maser diskNIFS H NIFS stellarNIFS integrated lightSDSS stellarSDSS emission linesHI spikeHI V HI velocity range
Figure 5.
The spatial scales on which the various velocitymeasurements considered in this paper are made. For the H i spike, we take the spatial scale to be ≥
600 pc as implied bythe brightness temperature limit given in Equation 2. Thefull peak-to-peak velocity range spanned by the H i emission(corresponding to the horizontal green line at the bottomof Figure 2) is shown as a vertical line that is horizontallyoffset from the V velocity for visual clarity. For the NIFSmeasurements, we plot the systemic velocities as a functionof annulus diameter and we include an overall 2 km s − cal-ibration systematic uncertainty on the error bars. All othervelocities are plotted with statistical error bars. though the maser disk has a position angle of ∼ ◦ whilethe outermost stellar and H components have positionangles of ∼ ◦ . We find that the systemic velocitiesof the stellar and H rotation curves agree with one an-other on the smallest scales ( ∼
30 pc), though they bothshow a ∼ σ redshift with the respect to the maservelocity. At larger radii the stellar and H systemic ve-locity measurements diverge, with the stellar systemicvelocity showing a ∼
30 km s − blueshift with respect tothe H on scales of ∼
200 pc. Such large variations in themeasured stellar systemic velocity as a function of radiusare rare; the typical dispersion of ATLAS galaxies be-tween the central and r = 500 pc velocities is only ∼ − (Appendix B; Krajnovi´c et al. 2011), and mostof the galaxies with substantially larger systemic veloc-ity gradients show evidence of interaction. We note, We note that this ∼ ◦ difference in position angle is consistentwith the offsets between maser disks and circumnuclear struc-tures seen in other galaxies (Greene et al. 2013) and comparableto the ∼ ◦ position angle difference between the J0437+2456maser disk and nuclear structure reported in Pjanka et al. (2017). however, that such a relative velocity offset could alsobe plausibly explained by a combination of geometricand obscuration effects (e.g., if the stellar and H emis-sion arose from two separate misaligned and mutuallyobscuring disks of material) while leaving the systemdynamically relaxed, and that the structure maps pro-duced by Pjanka et al. (2017) do show evidence of duston ∼ (cid:48)(cid:48) and larger scales.The outermost H emission ( ∼
200 pc) has a systemicvelocity of 4875 km s − that matches well with the emis-sion line velocity measured by P18 from the SDSS spec-trum, indicating that these two measurements may betracing similar material. These measurements are bothalso in agreement with the velocity of the “anomalous”central H i spike, which has a velocity of ∼ − and an amplitude of S ν ≈ . i spike represents a distinct dynamical subsystem(rather than, e.g., one “horn” of a double-horn profile),then we can set a lower limit on the area of the emissionregion by requiring that the H i brightness temperaturenot exceed its spin temperature of T s ≈
150 K (Condon& Ransom 2016), Ω ≥ S ν c kν T s . (2)Here, k is the Boltzmann constant, ν = 1 . i spike is ∼ . (cid:48)(cid:48) ≈
600 pc. This spatial scale is similar to thatprobed by the NIFS observations, and together with thecoincident velocities suggests that all three sources ofemission – i.e., the outermost H , the SDSS emissionlines, and the H i spike – may be originating from ma-terial with shared dynamics. The velocity of this mate-rial is significantly different from that of both the SDSSstellar and the central H i velocity (i.e., V ), perhapsindicating that there is a kinematically distinct sub-system located in the centermost few hundred parsecsof J0437+2456. However, we note that the observedFWHM of the H i spike of only ∼
55 km s − (see Table 4)is in tension with this interpretation, because at several-hundred parsec radii the material in this galaxy shoulddisplay a FWHM of ∼
200 km s − (Figure 3; see also No-ordermeer et al. 2007). It thus may not be viable tointerpret this H i spike as a distinct kinematic compo-nent. 4.2. Uncertainty in the black hole velocitymeasurement
Our measurement of the SMBH velocity (see Sec-tion 3.3) relies on accurate VLBI position measure-ments for each of the maser features, and if there are restless supermassive black hole in the galaxy J0437+2456 − usedas a reference feature) is known from phase-referencedVLBI measurements to a precision of better than 2 mas.G17 estimate that the expected additional positional un-certainties associated with this imperfectly-known refer-ence position, when propagated to the rest of the maserfeatures, should be (cid:46) µ as. This expectation is consis-tent with the magnitudes of the error floor parametersthat we recover from our model fitting (see Table 3).Additionally, we note that there are no obvious system-atic trends in the residual dispersion about the best fitsuch as would be expected if poor phase calibration werepresent at this level.4.3. An offset black hole
In our own Galactic Center, we have high-precision ev-idence that the SMBH is coincident with the dynamicalcenter of the Galaxy (Reid & Brunthaler 2020). Whilewe believe that a similar situation should generally holdfor other galaxies as well, a number of effects can atleast temporarily knock the SMBH out of this equilib-rium position. At very low galaxy mass, it is possiblethat SMBHs never settle at their galaxy center, giventhe very shallow galactic potential (e.g., Bellovary et al.2019; Reines et al. 2020). However, at higher galaxymasses, it is most likely that mergers are responsible forSMBH motions.Relative motions and spatial offsets between SMBHsand their host galaxies occur throughout the merger pro-cess. As galaxies merge, the SMBHs from each galaxywill be offset both spatially and in velocity from thecenter of the merger. This stage may be observableas velocity offset active galaxies (e.g., Comerford et al.2009; Comerford & Greene 2014) or as spatially resolvedpairs of active galactic nuclei (AGN; e.g., Komossa et al.2003; Gerke et al. 2007). Further along in the mergerprocess, when the two SMBHs become gravitationallybound, one may hope to observe the signatures of orbitalmotion for the bound pair (e.g., Eracleous et al. 2012;Shen et al. 2013; Ju et al. 2013, see also Appendix C).Finally, if an SMBH merger occurs, then any anisotropyin the radiated linear momentum will lead to a gravi-tational wave recoil (Fitchett 1983). These have beenmany observational recoil candidates proposed, but allhave their complications (see reviews in Komossa 2012and Blecha et al. 2016). The SMBH in the galaxy J0437+2456 is, to our knowl-edge, the most concrete case of an SMBH in motionwith respect to its galaxy. Because our initial searchfocused on megamaser disk galaxies (P18), the sourceswere all within 200 Mpc where detailed followup ob-servations are possible; luminous AGN that have beenidentified as recoil or binary SMBH candidates in thepast are often much more distant. Even in the case ofJ0437+2456, ambiguity remains about whether we areseeing an SMBH making its way to the galaxy centerfor the first time, SMBH binary orbital motion, or arecoil product. However, the fact that the galaxy onlarge scales is apparently out of equilibrium provides in-direct evidence that we are observing the aftermath ofa merger. SUMMARY AND CONCLUSIONFollowing the identification in P18 of the galaxyJ0437+2456 as a candidate for hosting a binary or re-coiling SMBH, we have obtained Arecibo and GeminiNIFS observations of the galaxy. Our new observationssupport the claim of a velocity offset between the SMBHand its host galaxy. Furthermore, the systemic velocityin J0437+2456 exhibits an apparent spatial scale de-pendence; the overall picture looks something like thefollowing:1. On the smallest spatial scales ( < O masers orbit witha central velocity of ∼ − . We associatethis velocity with the SMBH itself.2. At the photocenter of the galaxy, within the cen-tral ∼
30 pc and coincident with the location ofthe SMBH, both the stars and H gas emissionlines have a systemic velocity of ∼ − .However, on somewhat larger scales ( ∼ ∼
15 km s − in opposite directions.In all cases, these velocities are significantly offsetfrom the SMBH velocity as traced by the masers.3. On the largest spatial scales ( ∼ i emission is in agreement with theSDSS stellar velocity from P18. We find a cen-tral H i velocity of V ≈ − that we as-sociate with the recession velocity of the galaxyas a whole, though we note that the “anomalous”structure of the H i spectral profile complicates thisinterpretation.Multiple lines of evidence – including the differentinferred systemic velocities on different spatial scales,0the “anomalous” H i Facilities:
Arecibo Observatory, Gemini North
Software:
AOIDL, dynesty (Speagle 2020), Mon-tage , Gemini IRAF, Kinemetry (Krajnovi´c et al. 2006),pPXF (Cappellari & Emsellem 2004) http://montage.ipac.caltech.edu restless supermassive black hole in the galaxy J0437+2456 A. H i SPECTRAL MODELINGHere we describe three different models we use to fit the H i spectrum from Section 3.1. For each model, we use aGaussian likelihood given by ln ( L ) = − (cid:88) j (cid:32) S ν ( v j ) − ˆ S ν ( v j ) σ (cid:33) + ln (cid:0) πσ (cid:1) , (A1)where S ν ( v j ) is the model flux density for a spectral channel with velocity v j , ˆ S ν ( v j ) is the observed flux density in thatchannel, σ is the flux density uncertainty in a single channel, and the sum is taken over all channels. This likelihoodassumes that every spectral channel contains independent Gaussian-distributed noise with a standard deviation σ thatwe treat as a model parameter in each of our fits. We use the dynesty nested sampling code (Speagle 2020) forposterior exploration. The best-fit values and uncertainties for all model parameters are listed in Table 4.A.1. Modeling the profile using a sum of Gaussian components
The model we use in our primary analysis (Section 3.1) describes the H i spectral structure using a sum of Gaussiancomponents, S ν ( v ) = N (cid:88) i =1 A i exp (cid:34) − (cid:18) v − v i σ i (cid:19) (cid:35) , (A2)where the model parameters are the amplitude A i , central velocity v i , and width σ i for each component. The totalnumber of model parameters is 3 N + 1, where N is the number of Gaussian components; in this paper, we use N = 3.We impose uniform priors on all model parameters, in the range [0 ,
1] mJy for Gaussian component amplitudes,[0 , − for all Gaussian component standard deviations, [4500 , − for all Gaussian component centralvelocities, and [0 ,
1] mJy for σ . The posterior distribution is trivially multimodal upon pairwise swaps of Gaussiancomponents, but the modes are widely separated in parameter space and so we isolate a single mode when reportingparameter statistics. A plot of the resulting fit to the spectrum is shown in the left panel of Figure 6.A.2. Modeling the profile using a sum of double-horn components
Given that the galaxy J0437+2456 shows signs of dynamical disturbance (potentially indicating a recent merger)and that it exhibits an “anomalous” H i profile (Figure 2), it is natural to ask whether a combination of double-hornprofiles could give rise to the observed spectral structure. We have thus performed an alternative analysis using a sumof two double-horn components, each described using the parameterization developed by Stewart et al. (2014) for eachcomponent.The Stewart et al. (2014) model describes a double-horn profile using six parameters: the total flux, the centralvelocity, the velocity width, an asymmetry parameter, a parameter describing what fraction of the emission comesfrom solid-body rotation, and a velocity dispersion. Because we model the spectrum as a sum of N such double-horncomponents, and because we additionally model the channel uncertainty σ , the total number of model parameters is6 N + 1; in this paper, we use N = 2. We impose uniform priors on all model parameters, in the range [0 ,
1] Jy km s − for the total flux, [4500 , − for the central velocity, [0 , − for the velocity width, [ − ,
1] for theasymmetry parameter, [0 ,
1] for the solid-body fraction, [0 , − for the velocity dispersion, and [0 ,
1] mJy for σ .The results of fitting this alternative model to the H i data are shown in the central panel of Figure 6. We find thatthe best-fit model prefers only one of the two components to exhibit a standard double-horn profile, while the othercomponent is dominated by the solid-body contribution and so has only a single, wide spectral peak. This modelstruggles to fit the central H i spike, as evidenced by the large residual flux excess near ∼ − , so we disfavorit compared to the model composed of three Gaussian components.2 Table 4.
Results from H i spectral modelingModel description Parameter description Units Best-fit valuethree Gaussian components central velocity of first component km s − . ± . − . ± . − . ± . − . ± . − . ± . − . ± . . ± . . ± . . ± . . ± . − . ± . − . ± . − . ± . − . ± . − . ± . − . ± . . ± . . ± . . ± . . ± . − . ± . − . ± . . ± . − . ± . − . ± . − . ± . . ± . . ± . − . ± . − . ± . − . ± . . ± . . ± . Note —Results from fitting the three different models described in Appendix A to the Arecibo H i spectrum; these fitsare shown in Figure 6. For each parameter, we quote the posterior mean and standard deviation. The thermal noisehas been determined per 61 kHz ( ≈ − ) spectral channel. A.3.
Modeling the profile using a sum of double-horn and Gaussian components
Motivated by the appearance of the H i spectrum, we also attempt to model it using a sum of one double-horncomponent (parameterized as in Stewart et al. 2014 and Section A.2) and one Gaussian component. The resultingparameter values are listed in Table 4 and the best-fit spectrum is plotted in Figure 6. We again find that even thebest-fit model struggles to fit the observed spectral profile, with a substantial flux excess seen in the residuals around ∼ − . We thus disfavor this model compared to the model composed of three Gaussian components. restless supermassive black hole in the galaxy J0437+2456 F l u x d e n s i t y ( m J y ) three Gaussian components )0.250.000.25 R e s i d u a l s ( m J y ) two double-horn components ) one double-horn component andone Gaussian component ) Figure 6.
Similar to Figure 2, but showing the results of fits using the three different classes of model described in Appendix A;the best-fit parameter values for each model are listed in Table 4. In the left panel we show a more detailed breakdown of thefit from Figure 2 using three Gaussian components, with the individual best-fit Gaussian model components plotted in blue,green, and violet. Their corresponding best-fit velocities are marked by the vertical lines underneath each component. In thecenter panel we show a similar breakdown for the fit using two double-horn components, and in the right panel the fit using onedouble-horn component and one Gaussian component. In all panels, the spectrum is plotted at its native spectral resolution ingray, the spectrum after smoothing by a 4-channel boxcar is shown in black, and 1000 random posterior samples are overplottedin red. The bottom row of plots shows the residuals (i.e., the difference between the data and best-fitting model) for each fit.B.
ATLAS SYSTEMIC VELOCITY CURVESThe ATLAS project has collected integral field spectroscopic measurements for a sample of 260 early-type galaxiesin the local Universe (Cappellari et al. 2011). This sample provides a reference against which we can gauge the behaviorof the NIFS stellar systemic velocity measurements for J0437+2456, which show a systematic trend with radius (seeSection 3.2).Figure 7 shows the radial profile of the J0437+2456 stellar systemic velocity measurements plotted alongside thesame quantity measured for the “fast rotator” galaxies from ATLAS . The ATLAS sample is made up of early-typegalaxies, while J0437+2456 is a spiral, so for comparison we select only fast rotators from the ATLAS sample becausethey are galaxies with high angular momentum (Emsellem et al. 2011), stellar disks, and ordered (i.e., disk-like) stellarkinematics (Krajnovi´c et al. 2011, 2013). We note that unlike the ∼ observations were carried out under ∼ galaxies.At about 200 pc from the center the systemic velocity of a typical fast rotator deviates by only ∼ − from thesystemic velocity measured near the center. This trend does not change substantially with increasing radius.There are a few galaxies in the ATLAS sample that have systemic velocity deviations similar to or even largerthan those seen in J0437+2456, albeit at larger radii. The galaxies with the top four largest deviations are labeledin Figure 7: NGC 4753, UGC 09519, NGC 4342 and NGC 3665. NGC 4753, which has the largest difference in thesystemic velocity, also shows clear morphological evidence of a recent merger and contains complex dust filaments(Krajnovi´c et al. 2011; B´ılek et al. 2020), indicating that it is likely not in equilibrium. UGC 09519 might be dusty in4 | V s y s V s y s , | ( k m s ) NGC3665NGC4342UGC09519NGC4753J0437+2456
Figure 7.
The radial profile of the J0437+2456 systemic velocity as measured from the NIFS stellar emission (plotted in black;see Section 3.2) compared against similar profiles for “fast rotator” galaxies from ATLAS (plotted in blue). For each galaxy,we have subtracted off the systemic velocity measured at the smallest radii ( V sys , ) and then taken an absolute value of thedifference to aid comparison. the center, and it also has an unusual large-scale stellar disk characterised by blue colours and low surface brightness(Duc et al. 2015). NGC 3665 has a well defined nuclear dust and gas disk (Onishi et al. 2017), as well as asymmetricouter isophotes (B´ılek et al. 2020). NGC 4342 shows no evidence for disturbances in morphology or kinematics, exceptharbouring a central nuclear stellar disk (Scorza & van den Bosch 1998), though Cretton & van den Bosch (1999) notethat this galaxy has a remarkably large central velocity dispersion for its size and luminosity. The ATLAS sampleis perhaps not an ideal comparison sample, as it is made of early-type galaxies and the observations do not probe thesame spatial scales as the NIFS data of J0437+2456. Nevertheless, it is clear that the majority of ATLAS galaxiesdo not show strong variations in systemic velocity with radius, and there is some evidence that those with strongvariations tend to exhibit other indications of kinematic disturbance. C. OBSERVATIONAL CONSTRAINTS ON THE PROPERTIES OF A HYPOTHETICAL BINARY SMBHSYSTEM IN J0437+2456We observe the SMBH in J0437+2456 to have a velocity offset with respect to its host galaxy, as determined usingvarious different systemic velocity tracers (see Section 4.1). One possible explanation for this velocity offset is thatthe observed SMBH is part of a binary black hole system with a second, unseen SMBH. In this case, we have severalobservational constraints on the properties that such a binary system must have; these constraints are illustrated inFigure 8.Our first constraint comes from the fact that the observed SMBH in J0437+2456 is surrounded by an accretion disk,which is traced by H O maser emission to extend out to radii of ∼ , then its mass and separation from the observed SMBH mustbe such that it avoids tidally disrupting the accretion disk. This condition is roughly equivalent to requiring that the A second SMBH located within the innermost observed edge ofthe accretion disk would likely go undetected by the maser mea-surements (such a tight binary system would appear to the masersystem as a single SMBH with a mass equal to the combinedmasses of both SMBHs), but it would not by itself lead to an ob-served velocity offset between the maser measurements and thesystemic velocity of the host galaxy. restless supermassive black hole in the galaxy J0437+2456 Binary separation (pc)10 S e c o n d a r y m a ss ( M ) region excluded byHill sphere constraintregion excluded bySMBH velocity constraint region excludedby astrometry Figure 8.
Observational constraints on the space of secondary SMBH mass and binary separation for J0437+2456. Theblue shaded region is excluded by the requirement that the observed maser disk be tidally undisrupted, the red shaded regionis excluded by the requirement that the observed SMBH exhibit the measured velocity offset, and the gray shaded region isexcluded by the lack of an astrometric offset seen between the SMBH and the galactic center. The remaining unshaded regionindicates the permitted range of secondary SMBH mass and binary separation in the presence of these constraints. For thered and gray shaded regions, the solid and dashed lines represent 50% and 90% probability bounds, respectively, determined asdescribed in Appendix C. outer edge of the accretion disk lie within the Hill sphere of the observed SMBH. If we denote the mass of the observedSMBH as m , the mass of the second SMBH as m , their separation as r , and the Hill sphere radius as r H , then wecan cast this condition as an upper bound on m of m ≤ m r H (cid:34) r − r H ) − r (cid:35) − . (C3)The blue shaded region in Figure 8 shows the combinations of m and r that are excluded by this criterion. We usethe measured value of m = 2 . × M (cid:12) from G17 for the mass of the observed SMBH and the aforementioned valueof r H = 0 . v is related to the parameters of the binaryorbit via (see, e.g., Murray & Dermott 1999) v = m sin( i ) (cid:2) cos( ω + f ) + e cos( ω ) (cid:3)(cid:115) Gr (1 − e ) ( m + m ) . (C4)Here, i is the inclination of the orbital plane, ω is its argument of pericenter, f is the true anomaly of the observedSMBH, and e is the orbital eccentricity; m , m , and r are the same as in Equation C3. We do not currently have anyability to constrain the geometric parameters of the orbit, so we instead treat them probabilistically; we assume thatthe orbital plane is oriented randomly on the sphere (i.e., ω is distributed uniformly on [0 , π ] and cos( i ) is distributed6uniformly on [ − , f is oriented randomly on the circle, and that e is distributed uniformly in the range [0 , m and r given these assumptions about the distribution of possible orbit geometries. For the purposesof this constraint, we estimate the orbital velocity of the SMBH in J0437+2456 to be 48 km s − ≤ v ≤
101 km s − based on the measurements presented in this paper (see Table 2).Our third and final constraint comes from the apparent lack of an astrometric offset between the SMBH inJ0437+2456 and the center-of-light, determined in P18 to be (cid:46) Ahn, C. P., Seth, A. C., Cappellari, M., et al. 2018, ApJ,858, 102, doi: 10.3847/1538-4357/aabc57Barack, L., Cardoso, V., Nissanke, S., et al. 2019, Classicaland Quantum Gravity, 36, 143001,doi: 10.1088/1361-6382/ab0587Begelman, M. C., Blandford, R. D., & Rees, M. J. 1980,Nature, 287, 307, doi: 10.1038/287307a0Bellovary, J. M., Cleary, C. E., Munshi, F., et al. 2019,MNRAS, 482, 2913, doi: 10.1093/mnras/sty2842B´ılek, M., Duc, P.-A., Cuilland re, J.-C., et al. 2020,MNRAS, 498, 2138, doi: 10.1093/mnras/staa2248Blecha, L., Sijacki, D., Kelley, L. Z., et al. 2016, MNRAS,456, 961, doi: 10.1093/mnras/stv2646Bottinelli, L., Gouguenheim, L., & Paturel, G. 1982,A&AS, 50, 101Cappellari, M., & Copin, Y. 2003, MNRAS, 342, 345,doi: 10.1046/j.1365-8711.2003.06541.xCappellari, M., & Emsellem, E. 2004, PASP, 116, 138,doi: 10.1086/381875Cappellari, M., Emsellem, E., Krajnovi´c, D., et al. 2011,MNRAS, 413, 813, doi: 10.1111/j.1365-2966.2010.18174.xComerford, J. M., & Greene, J. E. 2014, ApJ, 789, 112,doi: 10.1088/0004-637X/789/2/112Comerford, J. M., Gerke, B. F., Newman, J. A., et al. 2009,ApJ, 698, 956, doi: 10.1088/0004-637X/698/1/956Condon, J. J., & Ransom, S. M. 2016, Essential RadioAstronomy (Princeton University Press)Courtois, H. M., Tully, R. B., Fisher, J. R., et al. 2009, AJ,138, 1938, doi: 10.1088/0004-6256/138/6/1938Cretton, N., & van den Bosch, F. C. 1999, ApJ, 514, 704,doi: 10.1086/306971Duc, P.-A., Cuillandre, J.-C., Karabal, E., et al. 2015,MNRAS, 446, 120, doi: 10.1093/mnras/stu2019Emsellem, E., Cappellari, M., Peletier, R. F., et al. 2004,MNRAS, 352, 721, doi: 10.1111/j.1365-2966.2004.07948.x Emsellem, E., Cappellari, M., Krajnovi´c, D., et al. 2011,MNRAS, 414, 888, doi: 10.1111/j.1365-2966.2011.18496.xEracleous, M., Boroson, T. A., Halpern, J. P., & Liu, J.2012, ApJS, 201, 23, doi: 10.1088/0067-0049/201/2/23Fitchett, M. J. 1983, MNRAS, 203, 1049,doi: 10.1093/mnras/203.4.1049Fouque, P., Durand, N., Bottinelli, L., Gouguenheim, L., &Paturel, G. 1990, A&AS, 86, 473Gao, F., Braatz, J. A., Reid, M. J., et al. 2017, ApJ, 834,52, doi: 10.3847/1538-4357/834/1/52Gerke, B. F., Newman, J. A., Lotz, J., et al. 2007, ApJL,660, L23, doi: 10.1086/517968Greene, J. E., Seth, A., den Brok, M., et al. 2013, ApJ, 771,121, doi: 10.1088/0004-637X/771/2/121Greene, J. E., Seth, A., Kim, M., et al. 2016, ApJL, 826,L32, doi: 10.3847/2041-8205/826/2/L32Gunn, J. E., Siegmund, W. A., Mannery, E. J., et al. 2006,AJ, 131, 2332, doi: 10.1086/500975Haynes, M. P., Giovanelli, R., Martin, A. M., et al. 2011,AJ, 142, 170, doi: 10.1088/0004-6256/142/5/170Ju, W., Greene, J. E., Rafikov, R. R., Bickerton, S. J., &Badenes, C. 2013, ApJ, 777, 44,doi: 10.1088/0004-637X/777/1/44Komossa, S. 2012, Advances in Astronomy, 2012, 364973,doi: 10.1155/2012/364973Komossa, S., Burwitz, V., Hasinger, G., et al. 2003, ApJL,582, L15, doi: 10.1086/346145Komossa, S., & Zensus, J. A. 2016, in IAU Symposium,Vol. 312, Star Clusters and Black Holes in Galaxiesacross Cosmic Time, ed. Y. Meiron, S. Li, F. K. Liu, &R. Spurzem, 13–25, doi: 10.1017/S1743921315007395Krajnovi´c, D., Cappellari, M., de Zeeuw, P. T., & Copin,Y. 2006, MNRAS, 366, 787,doi: 10.1111/j.1365-2966.2005.09902.x restless supermassive black hole in the galaxy J0437+245617