Millimeter Flares and VLBI Visibilities from Relativistic Simulations of Magnetized Accretion onto the Galactic Center Black Hole
The recent VLBI observation of the Galactic center black hole candidate Sgr A* at 1.3mm shows source structure on event-horizon scales. This detection enables a direct comparison of the emission region with models of the accretion flow onto the black hole. We present the first results from time-dependent radiative transfer of general relativistic MHD simulation data, and compare simulated synchrotron images at black hole spin with the VLBI measurements. After tuning the accretion rate to match the millimeter flux, we find excellent agreement between predicted and observed visibilities, even when viewed face-on (). VLBI measurements on –km baselines should constrain the inclination. The data constrain the accretion rate to be (–) with confidence, consistent with but independent of prior estimates derived from spectroscopic and polarimetric measurements. Finally, we compute light curves, which show that magnetic turbulence can directly produce flaring events with hour rise times, – hour durations and – flux modulation, in agreement with observations of Sgr A* at millimeter wavelengths.
Subject headings:accretion, accretion disks — black hole physics — radiative transfer — relativity — galaxy: center
Due to its large angular size, the Galactic supermassive black hole candidate Sgr A* is a promising laboratory for precision black hole astrophysics using very long baseline interferometry (VLBI) at millimeter wavelengths. Previous measurements (Shen et al., 2005; Doeleman et al., 2001; Krichbaum et al., 1998) at 7 and 3.5mm found a small instrinsic source size for Sgr A* (100–200 as FWHM), but at those wavelengths interstellar scattering is the dominant contribution to the observed size. Recent measurements by Doeleman et al. (2008) at 1.3mm are the first at short enough wavelengths to avoid contamination due to interstellar scattering and at long enough baselines to achieve event-horizon scale resolution. These measurements have been used by Broderick et al. (2009) to constrain the black hole and accretion flow parameters in terms of a 1.5D, nonrelativistic, time-steady radiatively inefficient accretion flow model (RIAF; Yuan et al. 2003).
Sgr A* is also known to exhibit multiwavelength, polarized flares (e.g., Yusef-Zadeh et al. 2007; Marrone et al. 2008; Eckart et al. 2008a; Dodds-Eden et al. 2009). This activity has been modeled as hotspots orbiting in the inner radii (Broderick & Loeb, 2005), adiabatically expanding blobs (Yusef-Zadeh et al., 2006; Eckart et al., 2008b) and emission from a short, mildly relativistic jet (Falcke et al., 1993).
As a complementary approach to the relatively simple models outlined above, we perform relativistic radiative transfer on data from a recent 3D, general relativistic MHD (GRMHD) simulation of a black hole accretion disk (Fragile et al., 2007), and compare simulated synchrotron images to the VLBI observations. In addition, we produce light curves at 1.3mm as a first step towards explaining Sgr A*’s flaring activity in terms of magnetic turbulence.
In contrast to the results from RIAF modeling, we find that for the single black hole spin considered (), all inclinations provide excellent fits to the VLBI data. Visibility profiles from low (face-on) inclinations () resemble “ring” models (e.g., Doeleman et al. 2008), whereas those at higher (edge-on) inclinations are more Gaussian. Ongoing work with increased array sensitivity should be able to distinguish between these two types of visibility profiles.
Detection of a ring-like profile would provide an important constraint not only on the geometry of Sgr A*, but also on the necessary complexity of accurate accretion flow models. These results also demonstrate the power of mm VLBI to distinguish between accretion models, as well as the model-dependence of attempts to extract black hole parameters.
2.1. Simulation Data
Fragile et al. (2007) presented global GRMHD simulations of black hole accretion disks. They used a spherical grid with polar Kerr-Schild coordinates, and the main untilted simulation had an effective resolution of , except near the poles which are underresolved. The torus was initialized with the analytic, time-steady, axisymmetric solution of De Villiers & Hawley (2003) and threaded with a weak poloidal magnetic field with minimum initially. The magnetorotational instability (MRI) arose naturally from the initial conditions, and the disk quickly became fully turbulent. The simulation was evolved for M, or 40 orbits at M in units with . Only data from the final of the simulation are used in this letter, once the disk is fully turbulent and transient effects from the initial conditions have died down.
This simulation evolved an internal energy equation, and injected entropy at shocks. Such a formulation does not conserve energy, and produces a more slender, cooler torus than conservative simulations which capture the heat from numerical dissipation of magnetic fields (Fragile & Meier, 2009).
2.2. Ray Tracing
We performed relativistic radiative transfer on the simulation data via ray tracing. Starting from an observer’s camera, rays are traced backwards in time assuming they are null geodesics (geometric optics approximation), using the public code described in Dexter & Agol (2009). In the region where rays intersect the accretion flow, the radiative transfer equation is solved along the geodesic (Broderick, 2006) in the form given in Fuerst & Wu (2004), which then represents a pixel of the image. This procedure is repeated for many rays to produce an image, and at many time steps of the simulation to produce time-dependent images (movies). Light curves are computed by integrating over the individual images.
To calculate fluid properties at each point on a ray, the spacetime coordinates of the geodesic are transformed from Boyer-Lindquist to the Kerr-Schild coordinates used in the simulation. Since the accretion flow is dynamic, light travel time delays along the geodesic are taken into account. Data from the sixteen nearest zone centers (eight on the simulation grid over two time steps) were interpolated to each point on the geodesic. Between levels of resolution near the poles, data from the higher resolution layer were averaged to create synthetic lower resolution points, which were then interpolated. Little emission originates in the underresolved regions of the simulation.
In this work, we consider synchrotron emission from Sgr A* using the formulae from Mahadevan et al. (1996) and neglecting polarization. The simulation provided mass density, pressure, velocity and magnetic field in code units. These were converted into cgs units following the procedure described in Schnittman et al. (2006), except that we scaled the torus mass to match the observed flux of Sgr A* at 1.3mm. The length- and time-scales were set by the black hole mass, taken to be , at a distance of kpc.
The temperature was taken to be thermal, with the ion and electron temperatures equal. Many authors have used two temperature models of Sgr A*, with weak coupling between the dynamically important ions and the cooler electrons. However, our single temperature model fits the observed flux from Sgr A* at reasonable accretion rates and produces a spectrum near millimeter wavelengths that is consistent with observations. In this model, the gas is optically thin except at high inclinations during flares. The gas temperature is in the regions of peak synchrotron emissivity, which is consistent with peak electron temperatures in RIAF models. We also followed Goldston et al. (2005) and produced images with , keeping the total gas pressure fixed. These models became optically thick before reaching the observed 1.3mm flux, and were inconsistent with the location of the peak in the quiescent spectrum of Sgr A*.
3.1. Visibility Fitting
Our fitting procedure closely mirrors that described in Broderick et al. (2009). Simulated images were produced for a range of inclination angles, time-averaged accretion rates, and observer times (sets of simulation time steps). The images were then averaged over 10 minute observation periods (8 simulation time steps) to match the integration times used in Doeleman et al. (2008), padded sufficiently to resolve the shortest baselines, and Fourier transformed to 2D visibilities. The visibilities were rotated to a range of different orientations, , the projection of the black hole spin axis on the sky measured E of N. Interstellar scattering was then applied by multiplication by an elliptical Gaussian as described in Fish et al. (2009) using the fits from Bower et al. (2006). Sample images and scatter-broadened visibilities are shown in Fig. 1. At low inclinations, the emission region is essentially a ring around the shadow, which causes a minimum in the 2D visibilities. For higher inclinations , Doppler beaming causes the emission to be concentrated on one side of the black hole. The visibilities are then approximately elliptical, with the minor (major) axis corresponding to the major (minor) axis of the emission region.
Visibilities were interpolated to the detections in Doeleman et al. (2008), and an effective was computed as defined in Broderick et al. (2009). Grids of minimum effective, reduced are shown in Fig. 2 over inclination angle and sky orientation. At inclinations , our results are in at least qualitative agreement with Broderick et al. (2009). The orientation is well constrained at high inclinations, since the general image shape is nearly static and asymmetric, and the long baseline VLBI measurements pick out a specific orientation for the visibility ellipse. The two distinct bands of good fits are due to the (approximate) up-down symmetry of the image. However, RIAF models ruled out inclinations , whereas for our GRMHD models, good fits (reduced ) are possible at all inclinations. This is especially evident from the curve of () vs. overplotted in Fig. 2. Inclinations are less probable due to the prior and large variation of probability density with observer time. Although the best fit is roughly the same for low inclinations at all , at any given observer time good fits are restricted to less than half the range of sky orientations.
Visibility profiles from the best fits at many inclinations are shown in Fig. 3. The profiles represent the 2D visibility plotted along the line representing the average angle in the UV plane of the long baseline measurements from Doeleman et al. (2008). Also plotted are the data from that work. At low inclinations, the profile reaches a minimum as described above, whereas at , the profile decreases monotonically with baseline length. Our best fit face-on visibility profile is almost identical to that of the ring model from Doeleman et al. (2008).
Our images are strongly peaked at the circular photon orbit, as can be seen in Fig. 1. Low inclination () images produced from a single time step of the energy conserving simulation by McKinney & Blandford (2009) are similar to that in Fig. 1, suggesting that our results are not due to the nonconservative nature of the simulation. It is not surprising that the VLBI observations are well fit by face-on disks. In addition to the ring model, a delta function intensity profile at the circular photon orbit roughly matches the observations, as does the spherical accretion model from Dexter & Agol (2009) scaled to the mass of Sgr A*. In fact, any model sharply peaked at the circular photon orbit should approximately match the visibility data. This is because the ratio of the zeroth to first maximum in the visibility for such a model is roughly the same as that of the zeroth order Bessel function, and agrees with the ratio of the visibility amplitude at to the total flux.
If this GRMHD simulation approximates the correct model for Sgr A*, we can also constrain the accretion rate for the single spin value used. Upper and lower limits for and confidence from probability distributions marginalized over observer time and sky orientation are shown in Fig. 4. Doppler beaming leads to larger fluxes at high inclinations, and hence a lower accretion rate is necessary to match the observed flux. The limits are also more stringent at high inclinations, where the flux is less variable (see Fig. 5). Further marginalizing over inclination, we constrain the accretion rate to be (–) at confidence. This constraint on the accretion rate of Sgr A*, while model-dependent, only uses the VLBI visibility measurements.
A major advantage of using GRMHD simulations in place of RIAF models is the ability to probe variability from the same dynamical model used to produce spectra and visibilities. Previous work has used hotspots (Broderick & Loeb, 2005; Doeleman et al., 2009) to study variability, but this introduces extra free parameters, complicating attempts to extract black hole parameters. The light curves of our Sgr A* models at 1.3mm shown in Fig. 5 are the first from simulation data to consider synchrotron emission, and the first to include the time delay along geodesics.
The variability on short timescales is slightly more noticeable for the case. Interestingly, though, the longer time scale flaring behavior is more pronounced in the light curve. The power spectra are well described by red noise spectra with power law index , for , . These are both steeper than the observed power law index of in Mauerhan et al. (2005) at 3mm wavelength.
The twin flares occur simultaneously in both light curves, rise over half an hour and last – hours. The flux modulation is () at (), measured from the peak of the flare to the average of the light curve immediately preceding it. All of these features are consistent with mm flare observations of Sgr A* (Eckart et al., 2008b; Marrone et al., 2008; Zhao et al., 2003; Li et al., 2009). Since the flares are seen at low inclinations, they are not caused by Doppler shifts from hotspots. The peak intensity is attenuated at large inclination due to rising optical depth. The flares are caused by a rise in magnetic field strength near the midplane in the inner radii (–M). They are not due to heat from magnetic dissipation since this is not possible in the simulation used here. Thus heating from magnetic reconnection is not necessary to produce the mm flares of Sgr A*.
In this letter, simulated images of a single temperature synchrotron emission model from 3D, GRMHD simulations by Fragile et al. (2007) have been fit to recent 1.3mm VLBI observations of Sgr A*. For the single black hole spin considered, the high inclination results are in at least qualitative agreement with those presented for 1.5D, stationary, non-relativistic models used by Broderick et al. (2009). However, we find that in contrast to previous work, images at low inclination provide excellent fits to the data, which are relatively insensitive to sky orientation, and present a qualitatively different visibility profile due to approximate circular symmetry.
Radial intensity profiles from face-on images of recent RIAF models (Yuan et al., 2009) appear more extended than those from our GRMHD simulations. The minimum of the visibility profile for these images occurs at too short a baseline to match the VLBI data. The model in Broderick et al. (2009), on the other hand, has increasing number density up to the event horizon, in stark contrast to the results from the simulations, where the number density peaks at , well outside of the event horizon at and the marginally stable orbit at . Large number densities in the inner radii of the model from Broderick et al. (2009) cause the “shadow” to be much smaller than expected (Bardeen, 1973; Falcke et al., 2000), weakening the power on small scales (large baselines), and preventing a good fit to the VLBI data at low inclinations.
Increased 1.3mm array sensitivity should favor either the low or high inclination profiles in Fig. 3. Our predicted profiles would also be sensitive to measurements at 1.3mm along baselines of length –km. A future detection favoring a ring-like profile would constitute the first measurement of a black hole shadow.
The visibility analysis constrains the average accretion rate of the simulation to be (–) with confidence. While this constraint is model-dependent, it is independent of and consistent with constraints from spectral fitting (Yuan et al., 2003) and polarization measurements (Agol, 2000; Quataert & Gruzinov, 2000; Macquart et al., 2006; Marrone et al., 2006; Sharma et al., 2007), and in good agreement with the previous spectral fits of single temperature GRMHD synchrotron images by Noble et al. (2007).
We have also produced the first light curves of synchrotron emission from GRMHD modeling of Sgr A*. At 1.3mm, flares arise naturally from magnetic turbulence, and are consistent with observations. Unfortunately our analysis is poorly suited for studying variability in other wavebands such as the radio, NIR or X-ray. In the radio, most the emission originates outside of the initial torus pressure maximum at , while the NIR and X-ray luminosity is produced by some combination of synchrotron emission from non-thermal electrons and Compton scattering of synchrotron photons. Modeling of both is left for future work. In addition, polarization measurements from flares are an important observational clue. Adding this capability to our relativistic radiative transfer code is in progress.
There are many uncertainties in the analysis presented in this letter. While GRMHD allows for a self-consistent solution of the primitive variables, it is unknown whether the MHD approximation is valid for the low densities at the Galactic center. In addition, there are a range of possible solutions, depending on the energy prescription used. Current simulations also neglect the interplay between ions and electrons, which should be handled self-consistently.
If Sgr A* is radiatively inefficient, simulations which conserve total energy are probably more appropriate. We estimate the time-averaged energy dissipation rate from numerical losses in the simulation as ergs for an accretion rate of . This is larger than the bolometric luminosity of Sgr A*, and gives an effective radiative efficiency of . The accretion flow is advection dominated, in that an order of magnitude more energy is lost to the black hole than to numerical reconnection. Face on synchrotron images from a single timestep of the 3D, conservative code in McKinney & Blandford (2009) were similar to those found here, suggesting that our results may not be specific to nonconservative simulations.
Despite the uncertainties, relativistic radiative transfer is an important tool for connecting state of the art, 3D, global GRMHD simulations with observations. Our analysis shows that ongoing millimeter VLBI observations will not only constrain black hole parameters, but accretion models as well. Expanding this analysis to the full range of spins will help determine the model dependence of such future constraints. If Sgr A* is at a low inclination, mm VLBI may soon provide the first direct evidence of a black hole shadow.
In addition to the visibility fits and light curves considered here, the variability of Sgr A* can be probed through measurements of closure phases along triangles of baselines (e.g., Doeleman et al. 2009). Identifying salient features of these closure phases is a goal of future work. If the flares seen here are robust features of mm synchrotron light curves of GRMHD simulations, multiwavelength and polarization measurements can be used to constrain emission mechanisms and the non-thermal electron distribution.
- Agol (2000) Agol, E. 2000, ApJ, 538, L121
- Bardeen (1973) Bardeen, J. M. 1973, in Black holes (Les astres occlus), 215–239
- Bower et al. (2006) Bower, G. C., Goss, W. M., Falcke, H., Backer, D. C., & Lithwick, Y. 2006, ApJ, 648, L127
- Broderick (2006) Broderick, A. E. 2006, MNRAS, 366, L10
- Broderick et al. (2009) Broderick, A. E., Fish, V. L., Doeleman, S. S., & Loeb, A. 2009, ApJ, 697, 45
- Broderick & Loeb (2005) Broderick, A. E., & Loeb, A. 2005, MNRAS, 363, 353
- De Villiers & Hawley (2003) De Villiers, J.-P., & Hawley, J. F. 2003, ApJ, 589, 458
- Dexter & Agol (2009) Dexter, J., & Agol, E. 2009, ApJ, 696, 1616
- Dodds-Eden et al. (2009) Dodds-Eden, K., et al. 2009, ApJ, 698, 676
- Doeleman et al. (2009) Doeleman, S. S., Fish, V. L., Broderick, A. E., Loeb, A., & Rogers, A. E. E. 2009, ApJ, 695, 59
- Doeleman et al. (2001) Doeleman, S. S., et al. 2001, AJ, 121, 2610
- Doeleman et al. (2008) —. 2008, Nature, 455, 78
- Eckart et al. (2008a) Eckart, A., et al. 2008a, A&A, 479, 625
- Eckart et al. (2008b) —. 2008b, A&A, 492, 337
- Falcke et al. (1993) Falcke, H., Mannheim, K., & Biermann, P. L. 1993, A&A, 278, L1
- Falcke et al. (2000) Falcke, H., Melia, F., & Agol, E. 2000, ApJ, 528, L13
- Fish et al. (2009) Fish, V. L., Broderick, A. E., Doeleman, S. S., & Loeb, A. 2009, ApJ, 692, L14
- Fragile et al. (2007) Fragile, P. C., Blaes, O. M., Anninos, P., & Salmonson, J. D. 2007, ApJ, 668, 417
- Fragile & Meier (2009) Fragile, P. C., & Meier, D. L. 2009, ApJ, 693, 771
- Fuerst & Wu (2004) Fuerst, S. V., & Wu, K. 2004, A&A, 424, 733
- Goldston et al. (2005) Goldston, J. E., Quataert, E., & Igumenshchev, I. V. 2005, ApJ, 621, 785
- Krichbaum et al. (1998) Krichbaum, T. P., et al. 1998, A&A, 335, L106
- Li et al. (2009) Li, J., Shen, Z.-Q., Miyazaki, A., Huang, L., Sault, R. J., Miyoshi, M., Tsuboi, M., & Tsutsumi, T. 2009, ApJ, 700, 417
- Macquart et al. (2006) Macquart, J.-P., Bower, G. C., Wright, M. C. H., Backer, D. C., & Falcke, H. 2006, ApJ, 646, L111
- Mahadevan et al. (1996) Mahadevan, R., Narayan, R., & Yi, I. 1996, ApJ, 465, 327
- Marrone et al. (2006) Marrone, D. P., Moran, J. M., Zhao, J.-H., & Rao, R. 2006, J. Phys. Conference Series, 54, 354
- Marrone et al. (2008) Marrone, D. P., et al. 2008, ApJ, 682, 373
- Mauerhan et al. (2005) Mauerhan, J. C., Morris, M., Walter, F., & Baganoff, F. K. 2005, ApJ, 623, L25
- McKinney & Blandford (2009) McKinney, J. C., & Blandford, R. D. 2009, MNRAS, 394, L126
- Noble et al. (2007) Noble, S. C., Leung, P. K., Gammie, C. F., & Book, L. G. 2007, Class. and Quant. Gravity, 24, 259
- Quataert & Gruzinov (2000) Quataert, E., & Gruzinov, A. 2000, ApJ, 545, 842
- Schnittman et al. (2006) Schnittman, J. D., Krolik, J. H., & Hawley, J. F. 2006, ApJ, 651, 1031
- Sharma et al. (2007) Sharma, P., Quataert, E., & Stone, J. M. 2007, ApJ, 671, 1696
- Shen et al. (2005) Shen, Z.-Q., Lo, K. Y., Liang, M.-C., Ho, P. T. P., & Zhao, J.-H. 2005, Nature, 438, 62
- Yuan et al. (2003) Yuan, F., Quataert, E., & Narayan, R. 2003, ApJ, 598, 301
- Yuan et al. (2009) Yuan, Y.-F., Cao, X., Huang, L., & Shen, Z.-Q. 2009, ApJ, 699, 722
- Yusef-Zadeh et al. (2006) Yusef-Zadeh, F., Roberts, D., Wardle, M., Heinke, C. O., & Bower, G. C. 2006, ApJ, 650, 189
- Yusef-Zadeh et al. (2007) Yusef-Zadeh, F., Wardle, M., Cotton, W. D., Heinke, C. O., & Roberts, D. A. 2007, ApJ, 668, L47
- Zhao et al. (2003) Zhao, J.-H., Young, K. H., Herrnstein, R. M., Ho, P. T. P., Tsutsumi, T., Lo, K. Y., Goss, W. M., & Bower, G. C. 2003, ApJ, 586, L29