FAST VARIABILITY AND MILLIMETER/IR FLARES IN GRMHD MODELS OF Sgr A* FROM STRONG-FIELD GRAVITATIONAL LENSING

, , , , , , and

Published 2015 October 13 © 2015. The American Astronomical Society. All rights reserved.
, , Citation Chi-kwan Chan et al 2015 ApJ 812 103 DOI 10.1088/0004-637X/812/2/103

0004-637X/812/2/103

ABSTRACT

We explore the variability properties of long, high-cadence general relativistic magnetohydrodynamic (GRMHD) simulations across the electromagnetic spectrum using an efficient, GPU-based radiative transfer algorithm. We focus on both standard and normal evolution (SANE) and magnetically arrested disk (MAD) simulations with parameters that successfully reproduce the time-averaged spectral properties of Sgr A* and the size of its image at 1.3 mm. We find that the SANE models produce short-timescale variability with amplitudes and power spectra that closely resemble those inferred observationally. In contrast, MAD models generate only slow variability at lower flux levels. Neither set of models shows any X-ray flares, which most likely indicates that additional physics, such as particle acceleration mechanisms, need to be incorporated into the GRMHD simulations to account for them. The SANE models show strong, short-lived millimeter/infrared (IR) flares, with short (≲1 hr) time lags between the millimeter and IR wavelengths, that arise from the combination of short-lived magnetic flux tubes and strong-field gravitational lensing near the horizon. Such events provide a natural explanation for the observed IR flares with no X-ray counterparts.

Export citation and abstract BibTeX RIS

1. INTRODUCTION

The black hole in the center of the Milky Way, Sgr A*, is the optimal laboratory for studying accretion physics at very low mass accretion rates. The mass and distance to the black hole are well constrained from dynamical studies of the orbits of stars in its vicinity (see, e.g., Ghez et al. 2008; Gillessen et al. 2009) and its brightness and variability have been studied extensively across the electromagnetic spectrum during the past decades (see Baganoff et al. 2001 and Genzel et al. 2003 for early studies).

Long monitoring observations have revealed that Sgr A* is highly variable at timescales comparable to the dynamical time near its event horizon. In the infrared (IR; e.g., K band), Do et al. (2009) found that the power spectrum of variability can be modeled by a power law in the 1–15 minute range with an index that varies between −1.8 and −2.8. Similar power-law indices were obtained in other IR bands and at 230 GHz (Meyer et al. 2008; Yusef-Zadeh et al. 2009; Dexter et al. 2014). Moreover, Meyer et al. (2009) found evidence for a break in the submillimeter (submm) and IR power spectrum at ∼150 minute (see also Dexter et al. 2014).

The amplitude of variability depends strongly on the wavelength of observation and changes from ∼10% at 1.3 mm to ∼50% at 2 μm, and up to ≳100% in the X-rays in three-hour intervals (e.g., Marrone et al. 2008; Porquet et al. 2008; Do et al. 2009). In the IR, it is debated whether the variability comprises peaked flaring events on top of a more quiescent variability (see Dodds-Eden et al. 2011 and references therein) or if all the observed variability is due to a continuum of events with a range of amplitudes and durations (Witzel et al. 2012). Nevertheless, simultaneous IR/X-ray observations of Sgr A* have established that even though all strong X-ray flares have IR counterparts, the opposite is not true (Hornstein et al. 2007). This suggests that either there are two physical origins for the IR events or that the physical mechanism responsible for the flares has two different observational manifestations.

This wealth of observations has fueled substantial progress in modeling the accretion flow around Sgr A* (see Yuan & Narayan 2014 for a review). In particular, the combination of general relativistic magnetohydrodynamic (GRMHD) simulations and relativistic ray-tracing calculations has led to a relative convergence of numerical models that account for numerous observed characteristics of Sgr A*, such as its quiescent broadband spectrum and the size of its millimeter image (Mościbrodzka et al. 2009; Dexter et al. 2010; Shcherbakov & Baganoff 2010; Chan et al. 2015). However, a consensus has not been reached on the physical mechanisms responsible for the rich phenomenology of the variability observed from Sgr A* and, in particular, for the origin of the observed IR and X-ray flares.

MHD simulations of accretion flows are inherently turbulent and variable but, as we showed in Chan et al. (2009), the excursions in magnetic field and density generated by MHD turbulence are not large enough to generate the observed large amplitude of variability, unless additional physics are incorporated. In Chan et al. (2009), we investigated the effect of external perturbations on the accretion flow. Dodds-Eden et al. (2010) explored the possible connection between flares and the regions of potential magnetic reconnection in MHD simulations. Dexter & Fragile (2013) investigated the possibility of shock heating in a tilted accretion disk as an alternate avenue to generating large variability and flares.

In parallel to the numerical simulations, a number of phenomenological models of the variability have been explored, which involve, e.g., orbiting hotspots or sudden injection of non-thermal particles (e.g., Markoff et al. 2001; Yuan et al. 2004; Dodds-Eden et al. 2009; Eckart et al. 2009; Dibi et al. 2014). These models can reproduce some characteristics of the data, but their connection to the turbulence and thermodynamics of the accretion flow is unclear.

The developments in the algorithms and the computational power in the last two years allow us to revisit the question of the broadband Sgr A* variability using long GRMHD simulations, with high-cadence ray tracing to compute observables, and more detailed descriptions of the thermodynamics of the accretion flow. We make use of GRMHD simulations of the accretion flow around Sgr A* performed with HARM (Narayan et al. 2012; Sa̧dowski et al. 2013) and radiative transfer calculations performed with the highly efficient GPU-based ray-tracing algorithm GRay (Chan et al. 2013). In earlier work, we used the combination of these two algorithms to carry out a comprehensive study of models with a wide range of black hole spins, magnetic field geometries, plasma physics models, and geometric parameters at high time, energy, and spatial resolutions (Chan et al. 2015). Among these, we identified five simulations with different characteristics that all reproduce the observed quiescent broadband spectrum and the 1.3 mm image size of Sgr A*. We use these as our five baseline numerical models.

In this paper, we calculate multi-wavelength light curves with a 10GMc−3 (≈3.5 minute for Sgr A*) time resolution that span approximately 10,000GMc−3 (60 hr) for each of the five simulations. We use these light curves to compute the amplitude of variability and the power spectrum over a wide range of wavelengths and to investigate the correlations and potential lags between different wavelengths. We find that the standard and normal evolution (SANE) models reproduce the short-timescale millimeter and IR variability with amplitudes and power spectra that are comparable to those inferred observationally, whereas magnetically arrested disk (MAD) models produce significantly less variability. We also find that none of the models produce X-ray flares, which most likely require additional physics, such as particle acceleration, that is not captured within the simulations. We trace the origins of the large-amplitude IR variability to short-lived, bright flux tubes and strong-field gravitational lensing near the black hole event horizon. We discuss the prospect of identifying these lensing events with future imaging experiments at millimeter and IR wavelengths.

2. THE GRMHD+RAY-TRACING SIMULATIONS

In Chan et al. (2015), we considered a large number of GRMHD models, with a wide range of black hole spins (0–0.9), accretion rates (which determine the scale of electron density in the accretion flow), initial configurations of the magnetic field (SANE and MAD in the terminology of Narayan et al. 2012), and with different prescriptions for the thermodynamic properties of the electrons in the accretion flow and in the magnetically dominated funnel. We then used GRay to compute average images and spectra from 100 snapshots from each simulation, with a time resolution of 10GMc−3, where G is the gravitational constant, M = 4.3 × 106 M is the mass of the Galactic black hole, and c is the speed of light. Hereafter, we will use gravitational units with G = c = 1. As described in Chan et al. (2015), we consider only thermal synchrotron and bremsstrahlung emission but neglect the possible contributions from non-thermal electrons and Compton scattering. We will address in the conclusions section whether these assumptions affect any of the results presented here.

We identified five GRMHD models that are able to account for the average spectrum of Sgr A* and the measured size of its emission region at 1.3 mm, as measured by the EHT. In particular, we required that each successful model reproduces (a) the correct average flux and spectral slope at 1011–1012 Hz, (b) a flux at ≃1014 Hz that falls within the observed range of highly variable IR fluxes, (c) an X-ray flux equal to 10% of the observed quiescent flux from Sgr A*, which has been attributed by Neilsen et al. (2013) to emission from the inner accretion flow, and (d) a scatter-broadened image at 1.3 mm with a size that is equal to the characteristic size measured by the EHT (Doeleman et al. 2008). We list the parameters of these five models in Table 1.

Table 1.  Summary of the Five Models from Chan et al. (2015)

Name a B0 Plasma Model tbegin tend
A 0.7 SANE Constant Te,funnel 89,760 99,990
B 0.9 SANE Constant Te,funnel 39,760 49,990
C 0.0 MAD Constant θfunnel 209,760 219,990
D 0.9 MAD Constant Te,funnel 38,520 48,750
E 0.9 MAD Constant θfunnel 38,520 48,750

Notes. A summary of the five best-fit GRMHD models of Sgr A* that we obtained in Chan et al. (2015). The first column lists the model names that we use throughout this paper. The parameters a and B0 are the black hole spin and the initial magnetic field configuration, respectively. We also specify the plasma models in the funnel, namely, constant electron temperature Te,funnel or constant electron–ion temperature ratio θfunnel, in the fourth column. The times listed in the last two columns correspond to the initial and final simulation times in gravitational units, GMc−3, that we use for the study of the variability presented here. Note that there are 1024 snapshots for each model separated by 10GMc−3, so that tendtbegin is always equal to 10,230GMc−3, which is approximately 60 hr for the mass of Sgr A*.

Download table as:  ASCIITypeset image

In this paper, in order to study the variability predicted by the GRMHD simulations, we use 1024 late-time snapshots of the best-fit models with the same time resolution (10M). For the mass of Sgr A*, this corresponds to a time resolution of ≈212 s with a span of ≈60 hr. For each snapshot, we compute four levels of images at each wavelength, with sizes equal to 1024M × 1024M, 256M × 256M, 64M × 64M, and 16M × 16M. Each level contains 512 × 512 pixels. Hence, the resolution improves by a factor of four in each direction going from one level to the next. We then calculate the total flux at each time and at each wavelength by composing all four levels of images together. In presenting the resulting flux, we often use the quantity $\nu {L}_{\nu }\equiv 4\pi {D}^{2}\nu {F}_{\nu },$ which would have been equal to the luminosity of the source had the emission been isotropic.

In the ray-tracing calculations, we use the fast light approximation, where the speed of light is assumed to be infinite, as was done in numerous previous studies (Dolence et al. 2009; Mościbrodzka et al. 2011; Shiokawa et al. 2012). This approximation is valid when the light crossing time tcrossing = R/c across a region of size R is much less than the dynamical timescale at that radius ${t}_{\mathrm{dynamical}}={({R}^{3}/{GM})}^{1/2},$ i.e., when RGM/c2. This is true for the emission from Sgr A* for wavelengths longer than approximately millimeter or shorter than ∼μm (IR). At millimeter to IR wavelengths, where the emitting region has a size that is comparable to the horizon of Sgr A*, the light crossing time is ≈20 s, which is a factor of 10 smaller than the time resolution of the snapshots. Because of this, the fastest timescales near the horizon of the black hole are not resolved and our power spectra will not be able to reveal any evidence for quasi-periodic variability (as seen in, e.g., Dolence et al. 2012), even if it were there.

In Figure 1, we plot the average predicted spectrum for each model together with the observational constraints used in Chan et al. (2015). Even though we have expanded the time span of the snapshots used by a factor of 10, the average spectra do not violate any of the observational constraints.

Figure 1.

Figure 1. Mean spectra of the five best-fit GRMHD models to Sgr A* (see Table 1) and their observational constraints. The colorful curves are the mean spectra $\nu \langle {L}_{\nu }\rangle $ computed using all 1024 snapshots. The purple curve is plotted with a dashed line for clarity when it overlaps with the blue or green curves. The crosses, vertical arrows, and bowties are the observational data that we employed in Chan et al. (2015) to constrain the models. The gray band between ν ≈ 1011 Hz and ≈1012 Hz marks the frequency range over which we perform a least squares fit to the radio spectrum. The gray line at ν ≈ 1014 Hz marks the IR frequency at which we used the range of fluxes observed at different times to impose an upper and lower bound on the models. The gray line at ν ≈ 1018 Hz marks the X-ray frequency, at which we used 10% of the observed quiescent flux (the open circle below the bowtie) to fix the density normalization in the flow. The white dotted line inside the gray band marks $\lambda =1.3\;\mathrm{mm},$ at which we fit the models to the EHT measurement of the image size.

Standard image High-resolution image

3. THE VARIABILITY OF SIMULATED LIGHT CURVES

In this section, we study the light curves from each of the five models using a variety of time-series analysis techniques. We compute the power spectra and the fractional rms amplitudes as functions of photon frequency and investigate correlations between different frequencies.

3.1. Light Curves

We first compute the light curves Lν(t) from each of the five models at four representative frequencies. In each panel of Figure 2, we plot the radio ($\nu ={10}^{10}\;\mathrm{Hz}$), 1.3 mm (ν = 2.31 × 1011 Hz), IR ($\nu =1.38\times {10}^{14}\;\mathrm{Hz}$), and X-ray (ν = 1.04 × 1018 Hz) light curves in red, black, green, and blue, respectively. We will use this color scheme when we compare different wavelengths throughout the paper.

Figure 2.

Figure 2. Light curves of the best-fit models. In each panel, we plot the radio, 1.3 mm, IR, and X-ray light curves in red, black, green, and blue, respectively. We offset and normalize all light curves for easy comparison. In addition, we also detrend the X-ray light curves.

Standard image High-resolution image

The fluxes and the variability amplitudes span a wide dynamical range across the wavelengths of interest, as we will discuss below. In order to be able to compare the light curves at different wavelengths, we first subtract the mean flux at each wavelength from the corresponding light curve and normalize it to the standard deviation calculated over the entire interval of the simulation. As a result, the mean is zero and the standard deviation is equal to unity for the normalized flux at each wavelength. The X-ray light curves show a continuous overall dimming by about 15% from the beginning to the end of the simulations because of the depletion of the torus that is feeding the black hole. (This is pronounced in the X-rays because the emission originates from a large volume but is negligible at the other wavelengths.) To remove this artifact of the simulation setup, we also detrend the X-ray light curve by fitting a line to the X-ray flux as a function of time and subtracting this from the light curve.

We show in Figure 3 the fractional rms amplitudes over the entire length of the simulations as a function of photon frequency for the five different models, which are the standard deviations that we used to normalize the light curves shown in Figure 2. At low photon frequencies, the amplitudes increase with frequency, peaking near the optical at ≈1015 Hz. This is the transition between the synchrotron and the bremsstrahlung-dominated parts of the spectrum. At the wavelengths where the synchrotron emission dominates, a small change in the density, temperature, or the magnetic field can cause a large change in the flux and, therefore, give rise to a large amplitude of variability. In contrast, the bremsstrahlung emission is optically thin and originates from large volumes, averaging out fluctuations in density and temperature.

Figure 3.

Figure 3. Fractional rms amplitude computed using the full light curve. The color scheme of the amplitudes is identical to the one used in Figure 1, except that the purple curve is plotted as dashed line when it overlaps with the blue curves. The amplitudes peak near the visible frequencies, i.e., $\nu \sim {10}^{15}\;\mathrm{Hz},$ which is the transition region between the synchrotron- and bremsstrahlung-dominated regions.

Standard image High-resolution image

We can draw several conclusions from Figures 2 and 3. First, the MAD Models D and E overall have larger amplitudes of variability, but their radio-to-IR light curves are dominated by smoother, long-timescale flux changes. In contrast, the disk-dominated models (SANE Models A and B) have overall smaller amplitudes of variability, but their light curves are characterized by short-lived, high-amplitude excursions in their IR and millimeter flux. As we will discuss in more detail below, these excursions have a lot of the characteristics of the IR flares observed from Sgr A*. Model C, in which the long-wavelength emission originates primarily from the base of the jet, has characteristics of both SANE and MAD models and shows both long- and short-timescale variability. At 1.3 mm, the Model C light curve shows evidence for quasi-periodic dimming, the origin of which we will explore below. Finally, in the X-rays, all models show very weak and smooth variations in the flux, with no evidence for any flaring events such as those observed from Sgr A* by Chandra.

Our multi-wavelength light curves allow us to look for correlations between the variability at different wavelengths. In the MAD models, all wavelengths (except X-rays) show quite similar trends in variability and are roughly synchronized. There is some weak evidence of the radio light curve lagging behind the 1.3 millimeter and the IR (see, e.g., the times between 6000M and 8000M in the lower panels of Figure 2). In the SANE models, it is only the 1.3 mm and the IR light curves that vary in sync, while the radio is decoupled. Moreover, the sharp IR flares have nearly coincident, albeit weaker, 1.3 mm counterparts with only very short (≲1 hr) time lags. In all models, the X-ray variability is decoupled from the long wavelengths, which is expected given their very different emission locations in the flow.

3.2. Flux Histograms

The range and character of variability observed from Sgr A* are often investigated via studies of flux histograms (see, e.g., Yusef-Zadeh et al. 2006, 2009; Do et al. 2009; Dodds-Eden et al. 2011; Witzel et al. 2012; Meyer et al. 2014). Figure 4 shows the flux histograms at 2.2 μm for our five models over the entire simulation time, with each flux corresponding to a time integration of 10M. (We chose this wavelength for this plot in order to compare it directly to the analysis of Dodds-Eden et al. 2011.) As in the case of the observational studies, we fit the low-flux ends of the histograms (up to 0.6 × 1036 erg s−1 for Models A and B and up to 0.03 × 1036 erg s−1 for Models C–E) with a lognormal distribution.

Figure 4.

Figure 4. IR flux histograms for the five models, evaluated at 2.2 μm (the green curves in Figure 2), together with the best-fit lognormal distributions. The SANE Models A and B are well described by the lognormal distribution, with marginal evidence for a small excess at large fluxes. Noting the different scales in the horizontal axes, the MAD Models C–E span a much narrower range in fluxes, but their distribution is significantly flatter than the lognormal curve.

Standard image High-resolution image

The flux histograms of the SANE Models A and B are well described by the lognormal distribution with marginal evidence for a small excess at large fluxes. In other words, the same physical mechanism appears to be responsible for most of the variability, with some evidence that the brightest of flares occur at larger than expected fluxes. This is consistent with our interpretation that the variability in these models is the result of the presence of short-lived hot, magnetically dominated flux tubes in the accretion flow and strong-field gravitational lensing.

The flux histograms of the MAD Models C–E are instead significantly flatter than the lognormal distribution at all but the smallest fluxes. This is another illustration of the fact that the IR light curves of these models show slow variability over a range of fluxes that cannot be described in terms of a range of short-lived flares.

3.3. Fractional rms Amplitudes

In Figure 3, we looked at the frequency dependence of the fractional rms amplitude using the entire ∼60 hr time span of the simulations. However, given the continuum of flux excursions and timescales visible in Figure 2, these amplitudes will depend on the time span used to calculate the rms and the particular interval for which it is computed. This is especially important when comparing theoretical amplitudes to observations, which usually span intervals that are much shorter.

In order to carry out the comparison to observations, we use the simultaneous multi-wavelength monitoring of Sgr A* of Marrone et al. (2008), which is one of the few data sets that spans the wavelengths from radio to X-rays. In particular, we use the 2006 July 17 light curve, during which an X-ray flare was detected. The longest time span of simultaneous observations is equal to 3 hr, and the range of rms amplitudes measured over a number of overlapping 3 hr intervals during this observation is shown as red data points with error bars in Figure 5. In the same figure, we plot as gray lines the predicted fractional rms amplitudes computed from many overlapping 510M (≈3 hr) intervals from the simulations.

Figure 5.

Figure 5. Gray curves show the fractional rms amplitudes computed using numerous overlapping 510M (≈3 hr) intervals from the simulations, starting at times ttbegin = 0, 100M, 200M,.... The blue curve in each panel, which shows the rms amplitude as a function of photon frequency computed using the entire span of the simulation, is identical to the corresponding curve in Figure 3. The red data points and error bars mark the range of rms amplitudes measured over a number of overlapping 3 hr intervals during the 2006 July 17 observations reported in Marrone et al. (2008).

Standard image High-resolution image

As expected, the simulations show a wide range of fractional amplitudes at each wavelength depending on the particular realization of the turbulent magnetic field and density in the accretion flow. However, all realizations have the same overall wavelength dependence and peak at optical frequencies. Furthermore, the data points in the millimeter and IR overlap with typical realizations of the simulated curves. On the other hand, the X-ray variability amplitude, which reflects the presence of the X-ray flare, is many orders of magnitudes larger than any of the predicted realizations.

This comparison between the observed and simulated variability amplitudes indicates that GRMHD simulations naturally produce IR flares and variability with no X-ray counterparts. Indeed, observations show two types of flaring events; IR flares are often not accompanied by X-ray counterparts, whereas X-ray flares always have IR counterparts (Hornstein et al. 2007). Interpreting the observations in the context of our simulations suggests that there are two distinct mechanisms that cause the two types of flaring activity. As we will discuss in detail in the next section, the turbulent GRMHD flows naturally produce IR flares without X-ray counterparts. However, additional physical processes must be incorporated into the simulations to generate the X-ray flares.

Given the sensitivity of the variability amplitude on the time span used, we can use the simulations to explore this dependence in detail. In Figure 6, we show the fractional rms amplitude at two different millimeter and IR wavelengths each, as a function of the time interval over which it is evaluated. For each value of the time interval, the red points correspond to the average fractional rms amplitude for the same overlapping segments in the simulation discussed earlier, and the error bars represent the variance among them. As expected, the overall trend is that the amplitude increases with increasing time interval, with no evidence for a characteristic timescale in the range explored (i.e., 5–40 hr). This is typical of red noise processes, which we can explore in more detail using power spectral techniques. In addition, the fact that the fractional rms amplitude in the IR becomes larger than 100% points to a large number of short-lived, high-amplitude flux excursions (flares) being responsible for this variability. This observation also explains the fact that the range of variability amplitudes in different 3 hr segments is larger in the IR than at millimeter wavelengths, since the amplitude of any particular segment in the IR depends strongly on whether a strong flare has occurred within its duration.

Figure 6.

Figure 6. Fractional rms amplitude at 1.3 mm, 0.86 mm, 4.5 μm, and 2.1 μm, as a function of the time interval over which it is evaluated. For each value of the time interval, the red points correspond to the average fractional rms amplitude for all the overlapping segments in the simulation starting at times ttbegin = 0, 100M, 200M,... and the error bars represent the variance among them.

Standard image High-resolution image

3.4. Power Spectra

Red noise variability is best studied in Fourier space, using power spectral techniques. However, generating a power spectrum requires uniformly sampled light curves with fast sampling and minimal data gaps. This is very difficult to achieve with ground-based observations of Sgr A* (e.g., in the millimeter and in the IR), which would require data spanning many hours but sampled every few minutes. For this reason, very few attempts have been made to measure the power spectrum of Sgr A*, and the few efforts have primarily been in the IR via the measurement of the structure function. In all cases, the power spectrum was modeled by a power-law function (Meyer et al. 2008; Do et al. 2009; Zamaninasab et al. 2010, and Yusef-Zadeh et al. 2009 for near-IR, and Dexter et al. 2014 for submm). In a very detailed study of the systematic uncertainties that are inherent to the measurement of power spectra via structure functions, Do et al. (2009) reported indices for the K-band power spectra in the range −1.8 to −2.8, while Meyer et al. (2008) reported evidence for a break in the power spectrum at ${154}_{-87}^{+124}$ minutes, with the power-law index changing from −2.1 ± 0.5 above this break to $-{0.3}_{-0.4}^{+0.2}$ below it (see also Dexter et al. 2014 for a similar report of a break at ${8}_{-4}^{+3}$ hr for submm wavelengths).

To make a quantitative comparison of the variability predicted in the GRMHD simulations to the variability reported in these studies, we compute the power spectra for each of the five models at four different photon frequencies and plot them in Figure 7. To compute the power spectra, we divide the 1024 snapshots from each light curve into seven overlapping segments, each of which contains 256 data points. We multiply these segments by a symmetric "triangle" window function. We then compute the fast Fourier transform of the windowed segments and average them. Finally, we bin the spectra with a bin size equal to the Fibonacci sequence. Each panel of Figure 7 contains the power spectra for the five models A–E, color-coded as before, for the radio, 1.3 mm, IR, and the X-ray light curves.

Figure 7.

Figure 7. Power density spectra of the normalized light curves shown in Figure 2. In each panel, five power spectra are plotted, one each for Models A–E. The panels, from left to right, are the power spectra of the light curves in the radio, 1.3 mm, IR, and X-rays. The color scheme matches the one we used in Figure 1, with the exception that for the IR and X-rays (the third and fourth panels), we plot the purple curves as dashed lines as they overlap with the blue curves. In the third panel, the shape of the power spectrum inferred observationally by Meyer et al. (2008) is shown as a dashed black curve, with its uncertainty depicted by a shaded gray band. In the remaining three panels, the dashed lines show power-law dependencies that are similar to those exhibited by the simulation data.

Standard image High-resolution image

The power spectra in the X-rays (fourth panel in Figure 7) show an f−4 dependence, which is characteristic of long-wavelength variability over timescales that are longer than the span of the simulations. We can, therefore, only use them to infer that there is no evidence for fast X-ray variability in the simulations. In all three of the longer wavelengths, however, the models give rise to power spectra with a range of flatter slopes. In the radio, the slopes can be approximately characterized as f−3, with Models C and E having somewhat steeper power spectra. At 1.3 mm, all models show a flatter dependence, with slopes comparable to ≈−2. In addition, Model C has the flattest spectrum among all of the models, consistent with the fact that it has the strongest short-timescale variability at this wavelength (see Figure 2).

In the IR, the two SANE models (A and B) produce higher relative power at higher frequencies, as expected from the visual inspection of light curves that show high levels of short-timescale variability. These two models also exhibit some weak evidence for flattening toward smaller frequencies. To compare these power spectra to observations, we also plot with a dashed line a power spectrum that has a break at 154 minutes, with slopes −0.3 and −2.1 below and above the break and a gray band that corresponds to the formal uncertainties in these slopes (as inferred by Meyer et al. 2009). This figure demonstrates that not only the overall amplitude but also the detailed power spectrum of variability generated by our SANE simulations are consistent with Sgr A* observations. Therefore, the observed IR power spectra appear to favor the SANE models as opposed to the MAD ones.

3.5. Correlations

A higher-order Fourier statistic that is often used in the variability studies of Sgr A* is the cross-correlation function between different wavelengths that is used to search for time lags in the light curves. Such lags are often interpreted as signatures of short-lived, hot emitted regions that are advected with the accretion flow or are ejected in a jet or outflow (e.g., Marrone et al. 2008; Eckart et al. 2009; Brinkerink et al. 2015). We compute the cross-correlation between the 1.3 mm light curve and the light curves at other wavelengths to look for such signatures.

The first three panels in Figure 8, from left to right, show the cross correlations as functions of time lag between the 1.3 mm light curve and the radio, IR, and X-ray light curves, respectively. A positive lag implies that the light curve at the particular wavelength lags behind the 1.3 mm light curve, whereas a negative lag implies that it leads the 1.3 mm light curve. There are clear differences between the 1.3 mm–radio correlations and 1.3 mm–IR correlations with the 1.3 mm–X-ray correlations. This is not surprising because most of the X-ray emission comes from large radii, while the 1.3 mm emission originates near the event horizon. In order to focus on the 1.3 mm–X-ray correlations of the inner flow, we also compute X-ray light curves using the inner 16M × 16M images only and cross correlate the 1.3 mm light curves with these inner X-ray light curves. The result is plotted in the fourth panel of Figure 8.

Figure 8.

Figure 8. Cross correlations between the light curves at different wavelengths and the light curve at 1.3 mm for the five models shown in Figure 2. For the last panel (labeled as inner X-rays), we consider only the emission that originates from the inner 16M × 16M portion of each image.

Standard image High-resolution image

In the SANE models, there is clear evidence for correlated variability only between 1.3 mm and IR curves with small (≲1 hr) time lags. In Section 4.1.1, we use the simulated images to argue that variability due to strong-field gravitational lensing is the cause of the achromatic character of some flares seen in the simulations as well as of the lack of time lags. On the other hand, there is no correlated variability between radio, 1.3 mm, and X-ray wavelengths, as can also be seen in the light curves shown in Figure 2. The two peaks appearing at time lags of ∼10 hr for Models A and B between the radio and 1.3 mm as well as the two peaks at ∼−10 hr time lags between the inner X-rays and the 1.3 mm simply correspond to the time differences between physically unrelated peaks in the light curves. Such chance events emphasize the difficulty of using cross-correlation analysis in short observations to infer the time evolution of emission structures in the accretion flows.

In the MAD solutions, the radio, 1.3 mm, and IR light curves are highly correlated, with marginal evidence of radio light curves lagging and the IR leading the 1.3 mm light curves by ∼1 hr. There also seems to be a similar level of correlation between the inner X-ray and 1.3 mm light curves, with the former leading the latter by ∼1 hr. This ordering suggests that the physical mechanism driving the variability in these solutions causes first the X-rays (in the inner disk), then the IR, then the 1.3 mm, then finally the radio flux to increase. In Section 4.2.2, by studying the simulated images, we conjecture that magnetic reconnection is responsible for the variability in these models.

4. THE PHYSICAL ORIGINS OF THE VARIABILITY

In the previous section, we focused on calculating light curves of the simulated emission at different wavelengths and their statistical properties. In this section, we take a deeper look at the simulated time-dependent images in order to explore and understand the physical origins of the different types of variability that we see in the simulations.

4.1. SANE Models

The variability in the SANE models is characterized by a spectrum of intense, short-lived flares in the IR, with related (albeit weaker) flux excursions at 1.3 mm, and uncorrelated slow variability in the radio and in the X-rays.

We first discuss the origin of the IR/1.3 mm flares, which we attribute to hot flux tubes in the accretion flow and strong-field gravitational lensing. We will then explore the origin of the radio variability, which we attribute to the interplay between emission from a bright funnel and osculation by a cold accretion disk.

4.1.1. Strong-field Lensing of Magnetic Flux Tubes

In Figure 9, we plot the 1.3 mm and IR light curves of Model B as well as the corresponding images at a quiescent reference and during several representative flares with different amplitudes. It is clear from the achromatic character of the bright structures in the images, as well as from their characteristic shapes, that strong-field gravitational lensing plays an important role in these flares.

Figure 9.

Figure 9. 1.3 mm (black) and IR (green) light curves for Model B (neither rescaled nor normalized). The middle and lower rows show the 20M × 20M images at 1.3 mm and in the IR at different points in time that correspond to a quiescent reference and to several flares. All images use the same quasi-logarithmic color scale (see the adjacent color bars), where unity represents the minimum of the peak brightness value of the images in each row. The small-amplitude flare in the second column is generated by strong-field gravitational lensing of a small, magnetically dominated region that crossed a caustic behind the black hole and gave rise to a nearly perfect Einstein ring. The large-amplitude flares in the third and fourth columns are generated by sheared magnetically dominated regions, the light from which appears both above and below the equatorial plane, because of gravitational lensing.

Standard image High-resolution image

In some instances, such as the one shown in the second column from the left in Figure 9, a hot, magnetically dominated structure that has a small spatial extent crosses a caustic, giving rise to a nearly perfect Einstein ring and an accompanying flux increase. This is the phenomenon that was investigated by Rauch & Blandford (1994), who calculated the properties of caustics around rotating black holes, and is similar to the large increase of the brightness of stars during microlensing events. In other instances, a sheared flux tube of large azimuthal extent leads to very bright flares. Because of its large extent, the flux tube is easily lensed and appears both above and below the disk (see the right two columns of Figure 9). The only requirement to see these lensing effects is that the flow needs to be optically thin in order for the lensed images not to be obscured by intervening material.

The actual spectra of the flares depend strongly on our prescription for the plasma thermodynamics in the magnetically dominated regions. We recall that in both Models A and B, when the plasma-β parameter exceeds a threshold (here, β > 0.2), we use a thermal plasma model with a constant electron temperature (see the discussion in Mościbrodzka & Falcke 2013; Mościbrodzka et al. 2014; Chan et al. 2015). The plasmas in these short-lived, magnetically dominated filaments may not have time to thermalize completely, and this will affect the particular spectrum of the flares, predominantly in the IR. Nevertheless, within our simplifying assumptions, the IR luminosity at 2.1 μm becomes an order of magnitude larger than 1.3 mm luminosity during the strong flares, from which we can infer the spectral index between these two wavelengths as

Equation (1)

Even though this is not directly comparable to the $\alpha \simeq -0.6\pm 0.2$ value observationally inferred by Hornstein et al. (2007) for the index of the IR spectrum alone, it is significant that our simple thermal model for the electrons in the magnetic filaments generate a flare spectrum that is substantially bluer than that of the average emission.

4.1.2. Bright Radio Funnels Obscured by Optically Thick Disks

In Section 3.5, we asserted that in the SANE models, the large (≳10 hr) positive time lags obtained between the radio and the 1.3 mm light curves are artifacts of matching physically unrelated flux variations between these wavelengths. Figure 10 shows the radio (3 cm) images calculated during the global minimum of the light curve, the two flare events shown in Figure 9, and the global maximum of the light curve. A comparison of Figures 9 and 10 shows that the emission at the two wavelengths originates from two different, and disconnected locations in the flow. The radio emission comes from a hot funnel, which is partially blocked by a colder, optically and geometrically thick accretion disk. Note also that the scale of these images is much larger than the scales of the 1.3 mm and the IR images shown in Figure 9. In contrast, all magnetic filaments that are responsible for the 1.3 mm/IR flares and their lensed images lie very close to the horizon and are completely obscured by the cold accretion disk at radio wavelengths.

Figure 10.

Figure 10. Radio (3 cm) light curve for Model B (neither rescaled nor normalized) and the 256M × 256M images at different points in time corresponding to the global minimum of the light curve, the two flare events shown in Figure 9, and the global maximum. All images use the same quasi-logarithmic color scale (see the adjacent color bar), where unity represents the minimum of the peak brightness value of all images. The emission originates primarily in the hot, magnetically dominated funnel and is partially obscured by the colder accretion flow. Any gravitationally lensed structures, such as those shown in Figure 9 for the 1.3 mm and IR images, are completely obscured by the cold accretion flow.

Standard image High-resolution image

4.2. MAD Models

The light curves of the MAD models do not show any evidence for short-lived flaring events but rather slower, correlated flux variations across the radio-to-IR part of the spectrum. Moreover, there is marginal evidence for the IR light curve leading the 1.3 mm one by ≃1 hr and the 1.3 mm light curve leading the radio (3 cm) one by a comparable amount.

In this section, we attribute a set of quasi-regular dimming episodes in the 1.3 mm light curve of Model C to a time-dependent loading of the footpoint of the jet funnel. We also point out some evidence for variability generated by dissipation of magnetic energy in localized events.

4.2.1. Time-variable Loading of the Jet Funnel

The non-spinning MAD Model C shows a peculiar set of short-timescale, quasi-regular dimming episodes at 1.3 mm, marked by arrows in the light curves plotted in Figure 11. The same figure also shows representative images during the flux minima and recovering phases of two such dimming events. It is evident that during each dimming event, e.g., the images shown in the first and third panels from left, the bright funnel region becomes less active and more collimated. As the brightness recovers, the funnel becomes active again, with waves of hotter regions propagating along the funnel walls. Note that this effect is quite different from the radio/millimeter dimming of Sgr A* discussed in Yusef-Zadeh et al. (2010), which may be caused by obscuration of optically thick plasma.

Figure 11.

Figure 11. 1.3 mm light curve for Model C and the corresponding 20M × 20M images during representative short-timescale, quasi-regular dimming episodes. Because the variation in the flux during each event is small, we use a linear color scale for images (see the adjacent color bar, where unity represents the minimum of the peak brightness values of the images). Each dimming event occurs when the mass loading of the footpoint of the jet funnel is reduced and the funnel itself becomes narrower.

Standard image High-resolution image

Among the five best-fit models of Sgr A* that we are considering in this paper, the non-spinning MAD Model C is the only one for which the emission is dominated by the footpoint of the jet funnel (see the images in Chan et al. 2015). It is, therefore, not surprising that this is the only model that is characterized by this dimming behavior. The individual episodes are short-lived (≲1 hr) and quasi-regular (once every ∼1–3 hr). Even though the overall behavior of these events is well captured by our simulations, the dimming depth is affected by the fact that we excise the innermost grid cells near the rotational pole of the spacetime to avoid artificially hot pixels, as documented in Chan et al. (2015).

4.2.2. Magnetic Reconnection

In the rapidly spinning, MAD Models D and E, the relatively large extent of the emission region across the long-wavelength spectrum and the lack of well-defined, short-lived flaring events makes it difficult to identify the physical origin of the variability in the simulated light curves. Nevertheless, there is some evidence for variability generated by dissipation of magnetic energy in localized events, as shown in Figure 12.

Figure 12.

Figure 12. Radio (red), 1.3 mm (black), and X-ray (blue) light curves for Model E and the corresponding 64M × 64M images after a representative reconnection event. The fluxes in the total and inner X-ray light curves have been scaled by factors of 10 and 1000, respectively, for clarity. The color scales of the images are quasi-logarithmic (see the adjacent color bars, where unity represents the minimum of the peak brightness values of the images for each row). The first column shows a snapshot where the inner X-ray light curve reaches a local maximum. We attribute this to a reconnection event where magnetic energy is dissipated into thermal energy and the topology of the field changes. This generates a hot spot in the X-ray image at the left side of the black hole shadow. In the second column, a twisted feature develops in the left-hand side of the 1.3 mm image and moves upward in the third column, at which point the 1.3 mm light curve reaches a local maximum. Finally, the fourth column shows a snapshot where the inner X-ray light curve reaches a local minimum, but the first twisted feature appears in the funnel in the radio image above the cold, obscuring accretion disk. This set of events resembles coronal mass ejection episodes observed in the Sun.

Standard image High-resolution image

Visual inspection of the images during and after one of these events reveals that the increase in the emission starts as a localized hot spot in the X-rays, most probably due to magnetic reconnection and dissipation of the magnetic energy. This hot emitting region then spirals upward, and the change in its physical location is mirrored in the evolution of its temperature and, hence, in the peak wavelength of its emission. This is the origin of the ∼1 hr long time lags between the IR, 1.3 mm, and radio light curves discussed in Section 3.5. This is similar to the coronal mass ejection episodes observed in the Sun (see, e.g., Forbes 2000 for a review) and is the closest analogy within GRMHD simulations to the expanding plasmon model of van der Laan (1966) that has been used in several phenomenological description of time lags from Sgr A* (see, e.g., Yusef-Zadeh et al. 2008, 2009, 2010).

5. DISCUSSION AND CONCLUSIONS

In this paper, we performed an extensive exploration of the predicted variability properties of GRMHD numerical models of Sgr A*, with parameters that have been fixed in earlier work (Chan et al. 2015) by fitting their time-averaged properties to the observed spectrum of Sgr A* and to the size of its 1.3 mm image. In particular, we investigated the variability properties of five models, two of which were SANE models and three that were MAD models.

We found that the SANE models generated photon-frequency-dependent amplitudes of variability and IR power spectra that are in agreement with those inferred observationally. Their light curves also show a number of short-lived, large flaring events that have short (≲1 hr) time lags between the IR and 1.3 mm, but have no X-ray counterparts. The distribution of IR flare amplitudes is well described by a lognormal distribution, with potential excess at large amplitudes, in agreement with the observations. We showed that strong-field gravitational lensing of hot, magnetically dominated regions contributes to these large, achromatic flares. The reason that such flares were absent in our earlier work (Chan et al. 2009) is because we incorporated a more simplified model for plasma thermodynamics and did not include there the effects of gravitational lensing. However, we also confirmed our earlier result that MHD turbulence alone with thermal plasma processes does not produce large X-ray flares since the X-ray emission originates at large radii, where gravitational lensing is negligible.

We speculate that particle acceleration processes need to be incorporated in the plasma model in the GRMHD simulations to account for the observed X-ray flaring. We expect that the presence of a short-lived distribution of relativistic electrons in the inner accretion flow will generate a different class of significant flux excursions from the X-rays down to the IR. However, incorporating such non-thermal electrons in the model will not significantly affect the properties of the flares we report here. Indeed, because the large IR flares are the result of strong magnetic filaments and their lensed images, the total scattering optical depth (for Compton emission) or emission volume (for synchrotron emission) of non-thermal electrons in these filaments will be too small to generate significant variability in the X-ray flux (see also the estimates in Dodds-Eden et al. 2009, Section 6). This is true even if we were to include a non-thermal distribution of electrons, as demonstrated quantitatively by Dodds-Eden et al. (2010). In contrast to the SANE models, the variability of the MAD models shows only slow trends, with no evidence for short-lived flares at any of the wavelengths.

The non-spinning MAD model, where the 1.3 mm emission is localized at the footpoint of the funnel, shows quasi-regular dimming events that arise from the time-dependent loading of the funnel. The rapidly spinning jet models show weak evidence for time lags between wavelengths, with longer wavelengths lagging the shorter ones. These time lags seem to arise from the upward motion of hot regions heated by magnetic reconnection/dissipation events. Overall, the variability properties of the MAD simulations do not appear to agree with those inferred observationally for Sgr A*.

Our results indicate that in the near future, it will be possible to observe strong-field gravitational lensing effects around black holes with the Event Horizon Telescope at 1.3 mm and with GRAVITY in the IR. Even though the flare events are too short-lived to obtain full images during a single lensing event, their presence can be inferred from the evolution of the closure phases obtained along triangles with appropriate baselines (Doeleman et al. 2009; Vincent et al. 2014).

Finally, we remark that other physics mechanisms can also affect variability. For example, light-travel effects between the source and the observer will suppress, in principle, the fast variability at long wavelengths. However, the recent discovery of a magnetar in the direction of Sgr A* suggests a less efficient suppression than what was believed until now (Spitler et al. 2014). Diffractive (and refractive) effects in the scattering screen will introduce additional variability compared to what we calculate here (Johnson & Gwinn 2015). Incorporating these effects requires vast computational resources and a detailed model of the scattering screen to Sgr A*, which is beyond the scope of the current paper.

C.K.C., D.P., and F.O. were partially supported by NASA/NSF TCAN award NNX14AB48G and NSF grants AST 1108753 and AST 1312034. D.M. acknowledges support from NSF grant AST-1207752. R.N. acknowledges support from NASA/NSF TCAN awards NNX14AB47G. L.M. acknowledges support from NFS GRFP grant DGE 1144085. All ray-tracing calculations were performed with the El Gato GPU cluster at the University of Arizona that is funded by NSF award 1228509. We thank Jon McKinney, Heino Falcke, Jason Dexter, and Ann Mao for valuable comments and discussions.

Please wait… references are loading.
10.1088/0004-637X/812/2/103